Quantum theory of atoms in molecules (QT-AIM) allows detailed insight into the electronic structure of molecules by analysis of the gradient vector field of the electronic density distribution function. First results of a QT-AIM analysis of neutral pterin as well as its anionic and cationic forms in aqueous solution are reported based on density functional theory using B3LYP/6–311+G(2d,p)//6–31G(d) level of theory. Besides reporting QT-AIM results of the atomic partial charges and bond orders of the molecules, their electron density functions, Laplacians and electrostatic potential functions are also shown. The results demonstrate the rather extensive delocalization of both negative and positive extra charges predominantly over the pyrimidine moiety of the pterin ring system and allow a precise quantitation of the effects of addition or elimination of a proton on the bonding structure. Quantum chemical techniques in combination with QT-AIM procedures are thus well-suited to describing the bonding structure as well as the major atomic and bonding features of neutral pterin and its cationic and anionic forms. Moreover, the data demonstrate that the chosen level of theory yields chemically reasonable results while not being computationally too expensive, thus providing a sound basis for further theoretical chemical work on pteridines.
Pteridine derivatives are widespread in biological systems, exhibiting a broad spectrum of interesting functions. From a chemical point of view, their character as pyrazino[2,3-d]pyrimidine, containing four ring nitrogen atoms, leads to remarkably complex chemical behavior. Thus, they have early attracted the interest of theoretically minded chemists, and during the last two decades of the 20th century, quite a number of papers described results of quantum chemical work on various biologically and/or chemically interesting pteridine derivatives [1–10]. Given the state of theoretical chemistry and particularly the availability of computing power during those times, it is not surprising that most of these earlier papers employed ab initio quantum chemical procedures confined to the Hartee-Fock method using rather small basis function sets or even resorting to using semi-empirical methods.
More recently, studies employing Density Functional Theory with moderately large basis function sets reported on spectroscopic properties of pterin in aqueous environments  and investigated the energetics of the plethora of possible tautomeric structures of neutral 6-methyl-pterin and its anionic and cationic forms in aqueous solutions [12–14].
Using pterin structures concordant with the most stable tautomeric forms of 6-methyl-pterin according to [12–14], I report here on detailed analyses of the neutral as well as the anionic and the cationic forms of pterin (2-aminopteridine-4(3H)-one, ), obtained at the B3LYP/6–311+G(2d,p)//6-31G(d) level of theory for the gas phase as well as for aqueous solutions (B3LYP is Becke’s 3-parametric density exchange functional using the correlation function by Lee, Yang and Parr) [16, 17]. Specifically, the Quantum Theory of Atoms in Molecules (QT-AIM)  is applied to the resulting wave functions yielding, among others, realistic atomic partial charges and bond orders. Additionally, for each of the three forms of pterin, I present visualizations of the electron density function, its gradient vector field topology and its Laplacian along with the electrostatic potential (ESP) function, in order to provide enhanced insights into the intimate bonding structures of the molecules.
This aim of this study is to evaluate, for the first time, the application of QT-AIM analysis on good quality density functional theory wave functions of pteridine systems with particular focus on binding characteristics of such molecules and to lay a firm foundation for further theoretical work on biologically and chemically interesting pteridine derivatives.
Traditional quantum chemistry faced some difficulty with the concept of “an atom within the context of a molecule”, i.e., neither the molecular wave function nor the electron density function computed from the latter yield a direct perception of “atoms” or of such age-honoured and incredibly important and fruitful chemical ideas as “electron pairs” (bonding as well as lone pairs). It was the outstanding work of Richard WF Bader and his co-workers over decades, excellently summarized in , that paved the way for a physically correct and unique quantum mechanical description of “atoms in molecules” and their properties. Briefly, an atom within the context of a molecular system is defined in QT-AIM on the basis of a topological analysis of the gradient vector field associated with the electron density function ρ(x, y, z), a scalar function in physical space which is an “observable” in the strict sense of quantum mechanics. The major topological feature of this function is that its local maximum values occur at the nuclear positions. All gradient paths of the associated gradient vector field ∇ρ(x, y, z) (the “nabla” operator ∇ is defined as:
with unit vectors and ), with very few and special exceptions, are “attracted” by these maximum positions where, of course, the gradient vectors become zero. Moreover, besides the nuclear positions (atomic critical points), there are also other “critical points” at which the gradient vector vanishes ( ). The true nature of such a critical point is determined by the behavior of the tensor ∇∇ρ leading to the so-called “Hessian matrix” of the nine possible second derivatives of ρ. The latter can be diagonalized and their “eigenvalues” λ1, λ2 and λ3 can be determined. In 3D space, four types of critical points can be distinguished and characterized as pair (rank, signature), according to their rank (the number of the non-zero “curvatures”, i.e., the non-zero second derivatives) and signature (the sum of the signs of the eigenvalues):
(3, -3): atomic critical points (ACP): local maximum points of ρ (three negative curvatures)
(3, -1): bond critical points (BCP): saddle points between ACP (two negative curvatures, one positive)
(3, +1): ring critical points (RCP): saddle points characteristic for ring system (one negative, two positive curvatures)
(3, +3): cage critical points (CCP): local minimum points, characteristic for caged molecules (three positive curvatures).
Each gradient path originates or terminates either at a critical point, or at infinity. There are “special” gradient paths; for example, those originating at a BCP and terminating at the two neighboring ACPs. Their entirety defines the framework of the chemical bonds of a molecule. Other special gradient paths originate at infinity and terminate at one BCP, or they connect an RCP with a BCP; in their entirety, these paths define the interatomic surfaces.
In QT-AIM, the topological analysis of ρ(x, y, z) provides a means of partitioning the physical space of a molecule into non-overlapping domains (“basins of attraction”) which consist of the entirety of all gradient paths terminating at one specific ACP, i.e., at one specific nucleus. The surfaces between adjacent basins of attraction are defined by a “zero-flux” condition. This means that no gradient path can cross such a surface, and this condition guarantees the identification of the unequivocally defined basins of attractions as “atoms in a molecule” in a strictly physical sense, allowing the computation of many atomic properties. For example, the electronic charge of an atom can be obtained by integrating the electronic density function over the whole basin of attraction of that atom. Similarly, integration of the electron density function over the interatomic surface yields a physically sound measure of the bond order between the involved atoms. Importantly, QT-AIM properties of an atom in a molecule are derived from the electron density function (which is obtainable either by quantum chemistry or by X-ray diffraction experiments) and are not, in contrast to other quantum chemical techniques, estimated by more or less arbitrary partitioning schemes between atoms applied to molecular orbitals (which are not “observables” in a strict sense).
QT-AIM also emphasizes the importance of the so-called Laplacian Δ(ρ), which is obtained by computing, at each point of space, the sum of the diagonal elements of the diagonalized Hessian matrix:
a scalar field which is useful for visualizing the local concentration (Δ(ρ)<0) or depletion of electron density (Δ(ρ)>0) thus providing chemically meaningful information regarding chemical bonds and/or nonbonding electron pairs.
Of the many possible tautomeric structures of neutral, cationic and anionic pterin, I choose analogous structures according to the minimum energy tautomers of 6-methyl-pterin [12–14]. (Note: In [12, 13] 6-methyl-pterin erroneously is designated as “pterin”.) For each of the three studied structures, a reasonable starting geometry is obtained by Avogadro software . These structures are then employed as input geometries for a geometry optimization at the B3LYP/6–31G(d) level of theory, using the GAUSSIAN suite of quantum chemistry programs (Gaussian G09W, version 9.5, Gaussian Inc., Pittsburgh, PA, USA). In order to ensure that the geometries found are indeed minimum energy structures, a frequency calculation is then performed at the same level of theory. Notably, when strictly planar starting geometries are used, in the case of neutral and anionic pterin, the first of these two steps indeed yields planar optimized structures which, however, by the second step are recognized as transition states. The frequency analysis reveals that the amino group at C2 is not strictly planar in a true minimum geometry (the two hydrogen atoms are slightly bent off the molecular plane). In contrast, in the cationic form, even when starting with a non-planar geometry, the optimum geometry obtained is strictly planar; obviously a result of the positive charge, the delocalization of the lone pair at the nitrogen atom is substantially stronger, thus favouring the planar structure. All results shown in this paper refer to the true minimum energy geometries according to frequency analysis.
The remaining computations of the final electronic energy, the molecular wave function, the electron density distribution, the electrostatic potential (ESP), and the atomic partial charges and bond orders according to QT-AIM, for each structure were done using G09W software at the B3LYP/6–311+G(2d,p) level of theory. As shown in , the B3LYP/6–311+G(2d,p)//6–31G(d) level of theory yields quite accurate results while at the same time being computationally feasible for medium-sized molecules as well. The choice of a relatively extended basis set including polarization and diffuse functions for the final computations is motivated by the aim to also get reliable results for the anionic form of pterin. Notably, all quantum chemical computations were performed in the gas phase as well as in solution phase. In this work, the SMD model  is employed. The SMD model is a continuum solvation model based on quantum mechanical charge density of a solute molecule interacting with a continuum description of the solvent. Full electron density is employed without defining partial atomic charges. The model is universally applicable to charged or uncharged solute molecules in any solvent for which a few key descriptors such as the dielectric constant, the refractive index, the bulk surface tension and the acidity and basicity parameters are known.
In-depth QT-AIM analysis (computation of the critical points of the electron density function as well as the molecular graph and the gradient paths) and visualization of these results as well as of the atomic basins of attraction is done using the programs AIM2000 version 2.0 (Innovative Software, F. Biegler-König, J. Schönbohm, Bielefeld, Germany) and AIMStudio Version 13.11.04 Professional (TK Gristmill Software, Todd A. Keith, Overland Park, KS, USA). These programs are also used to verify the QT-AIM charges obtained with G09W. Finally, visualizations of the electron density function, the Laplacian and the ESP of the molecules are done using AVS Express 5.3 software (Advanced Visual Systems Inc., Waltham, MA, USA).
Figure 1 shows the chemical formulae of the three molecular systems studied (left column) and first results of the QT-AIM analyses (right column). The electron density functions, computed in the molecular planes, are depicted as isodensity lines. Perpendicular to these, some gradient paths are shown. There are three types of critical points (ACP, BCP and RCP); the “special” gradient paths connecting these are emphasized in red color. The set of the critical points together with these special gradient paths define the “molecular graph”. Notably, the molecular graphs are calculated exclusively from the properties of the electron density functions; however, as the pictures demonstrate, they show perfect analogy with the chemical formulae.
In order to visualize the 3-dimensional “basins of attraction” which one can identify as “atoms” in the sense of QT-AIM, Figure 2 shows selected atoms for the case of pterin. Obviously the shapes of these basins of attraction are quite different for different elements. For example, the strongly electronegative oxygen and nitrogen atoms appear more bulky than the positively charged carbon atoms, and the hydrogen atoms appear as “caps” sitting on their respective bond partners. Importantly, all atomic properties according to the QT-AIM framework are obtained by integration over the respective atomic basin of attraction.
Thus, by integrating the electron density function of the atomic basins of attraction, one obtains the atomic partial charges (Figure 3, left column). Notably, charges obtained by QT-AIM are generally considerably larger (in terms of absolute values) than the charges obtained by other quantum chemical techniques, including Mulliken as well as natural bond orbital (NBO) population analyses (Table 1 shows the comparison of the atomic partial charges obtained by the three methods for the atoms of neutral pterin). It is obvious from a comparison of the three forms of pterin that the largest changes of the charge distribution among the ring atoms occur in the pyrimidine moiety of the pterin ring system. In contrast, the pyrazine part remains nearly unaffected from addition or elimination of a proton. The strongest changes occur – in this order and with the expected direction – at carbon C2, oxygen and carbon C4 (observe Figure 1, top, for the numbering scheme).
|Atom||Mulliken charge||Natural charge||QT-AIM charge|
|H at N3||0.352||0.449||0.463|
|H at C6||0.151||0.216||0.096|
|H at C7||0.155||0.216||0.100|
|H at N2′ (N1)a||0.305||0.432||0.446|
|H at N2′ (N3)a||0.305||0.428||0.445|
aThe two hydrogens of the NH2-group are oriented towards N1 and N3, respectively.
Figure 3 (right column) shows isosurfaces of the electron density functions of pterin and its cation and anion mapping, by colors, of the local ESP. Additionally, the ESP functions obtained in the molecular planes are visualized as contour diagrams in order to demonstrate how the molecules would interact with a positive probe charge. Red contour lines (positive ESP) define parts of the molecules’ neighborhoods where a positive probe charge would be repelled, and blue contour lines (negative ESP) denote regions where a positive probe charge would be attracted. While neutral pterin shows negative ESP along the axis of the oxygen atom and nitrogen atom N5, and along nitrogen atoms N3 and N8, a pterin cation is surrounded nearly exclusively by positive ESP. In contrast, a pterin anion is surrounded by negative ESP.
Similarly, Figure 4 (left column) shows the bond orders obtained by integrating the electron density along the interatomic surfaces and visualizations of the Laplacian functions of the molecules, overlaid on an electron density isosurface (right column). When comparing the cation as well as the anion with the neutral form of pterin, it is obvious that the largest differences are again seen in the bonding situation of the pyrimidine moiety, i.e., the pyrazine part remains essentially unchanged. The strongest variation occurs along the N1–C2 bond, and, interestingly, bond order (and, thus, the partial double bond character) decreases for the cation as well as for the anion, when compared to neutral pterin. The second-largest changes occur along the C2–N3 bond, where bond order increases for both the cation and the anion. For the bonds between N3–C4 and N1–C8a, bond order increases from the cation via neutral pterin to the anion, and for the C4–O and the C2–NH2 bond, bond order decreases in that direction. These changes are compatible with the formulation of mesomeric resonance structures that are shown in Figure 5. In the cationic form, strong delocalization of the positive charge between the nitrogen atoms bonded to C2 is found, thus increasing bond strengths along the C2–N3 and the C2–NH2 bonds and decreasing bond strength along the N1–C2 bond. In the anion, delocalization of the extra negative charge is observed mainly from N3 to the oxygen atom, thus increasing the aromatic character of the pyrimidine part of the system.
These changes of the bond orders are nicely reflected by the Laplacian functions also shown in Figure 4. As an example, at the chosen isovalue of –2.6 a.u. the negative charge concentration (yellow isosurfaces) in neutral pterin is interrupted along the C2-N3 bond, and a contiguous isosurface is found along this bond in the cation as well as in the anion. On the other hand, in the cation, an interruption of this isosurface is seen along the N1-C8a bond, while in the anion a nearly contiguous charge concentration in the whole pyrimidine moiety is found. Only between C4 and C4a is there a gap in the respective isosurface. In contrast, the negative charge density along the C2-NH2 bond is markedly “thinner” in the anion, pointing to a decrease of the delocalization of the lone pair at the NH2 nitrogen atom. In accordance with the atomic partial charges as well as the bond orders, in the pyrazine ring, in all cases a contiguous isosurface is observed, demonstrating its rather undisturbed aromatic character in all three forms of the molecule.
Particularly for the more electropositive carbon atoms, at the chosen isovalue of –2.6 “holes” in the yellow surfaces at the atomic positions allow a glance at the “core” regions of electronic charge depletion (in red) which mark the “gap” in electronic density between the core and the valence regions of the respective atoms. Inside these core depletion regions one finds even smaller yellow regions of the core regions of electronic charge accumulation, corresponding in principle to the 1s2 electrons.
All results shown above are obtained using the SMD model simulating solution of the molecules in water. As stated in the Methods section, all molecules were also treated by the same quantum chemical procedures in the gas phase. There are no large differences in the QT-AIM atomic charges, which may be best summarized by the molecular dipole moments. While for neutral pterin, solvation caused an increase of the dipole moment from 1.555 to 1.678 D (debye), for cationic pterin (1.656 in gas phase vs. 1.633 D in solution) and for anionic pterin (1.240 vs. 1.226 D, respectively) solvation in water does not alter the polarisation markedly. In contrast, the solvation enthalpy change (ΔHsolv) is markedly larger in absolute size for the two charged species: ΔHsolv for neutral pterin is –25.2 kcal/mol; for cationic pterin, – 76.7 kcal/mol; and for anionic pterin, –69.1 kcal/mol.
In this paper I analyse the electronic structure of neutral as well as cationic and anionic pterin using quantum chemistry at the Density Functional level of theory, combined with – to the best of my knowledge – the first application of QT-AIM analysis in the field of pteridine research. In order to get a reliable description of the electronic structure, I employ geometry optimization with the 6–31G(d) basis set and computation of the final wave function with the extended 6–311+G(2d,p) basis set. The use of polarization functions as well as diffuse functions is mainly dictated by the desire to also get a good representation of the electron-rich anion of pterin.
One might criticize that I have not considered all the possible tautomeric forms of neutral, cationic and anionic pterin. However, here I aim to demonstrate the basic usefulness of QT-AIM for the analysis of the electronic structures of pteridines in general, and I have chosen starting molecular geometries which are in accordance with “chemical” knowledge as well as with recent results of quantum chemical calculations on 6-methyl pterin using density functional level of theory, albeit at smaller basis set sizes [12–14].
QT-AIM in combination with high-level ab initio quantum chemistry enables detailed insight into the electronic structure and the bonding situation of neutral as well as charged forms of pterin. In particular, the changes occurring when adding or eliminating a proton are well described and are in excellent accordance with conventional chemical wisdom and concepts. The particular charm of QT-AIM is the ready transformation of more abstract quantum chemical constructs, such as molecular orbitals, into a picture of the bond structures which is much more similar to traditional chemical concepts such as bonding or non-bonding localized electron pairs according to Lewis’ classical model . In the particular investigated forms of pterin, these classical ideas can be nicely identified, employing – in a physically correct way – only the topological properties of the electron density function which by itself is an observable in the strict sense of quantum mechanics.
The present results are encouraging with regards to further work on interesting aspects of pteridine chemistry and biochemistry. Theoretical studies of this type may provide deeper insight into, for instance, the field of spectroscopic properties of pteridines or even the quite complicated and by no means fully understood interactions of various pteridine derivatives with free radicals – a topic most relevant to the plethora of biological roles and functions of these fascinating molecules.
1. Bieri J, Geiger RE. Über Pterinchemie. 51. Mitteilung. CNDO-Rechnungen an pterin, 6,7-dimethyl-7,8-dihydropterin and 5-formyl-6,7-dimethyl-5,6,7,8-tetrahydropterin. Helv Chim Acta 1975;58:1201–5. Search in Google Scholar
2. Gready JE. The value of the π-bond order – bond length relationship in geometry prediction and chemical bonding interpretation. J Comput Chem 1984;5:411–26. Search in Google Scholar
3. Gready JE. Theoretical studies on pteridines. 2. Geometries, tautomers, ionization and reduction energies of substrates and inhibitors of dihydrofolate reductase. J Comput Chem 1985;6:377–400. Search in Google Scholar
4. Gready JE. Analysis of the inter-ring coupling effects in N-heterobicyclic π-systems using a structural definition of aromaticity. Int J Quant Chem 1987;31:369–82. Search in Google Scholar
5. Uchimaru T, Tsuzuki S, Tanabe K, Benkovic SJ, Furukawa K, Taira K. Computational studies on pterins and speculations on the mechanism of action of dihydrofolate reductase. Biochem Biophys Res Commun 1989;161:64–8. Search in Google Scholar
6. Katoh S, Sueoka T, Kurihara T. Computer studies on the stereostructure and quantum chemical properties of 6-pyruvoyl tetrahydropterin, the key intermediate of tetrahydrobiopterin biosynthesis. Biochem Biophys Res Commun 1991;176:52–8. Search in Google Scholar
7. Reibnegger G, Denny BJ, Wachter H. Ab initio quantum chemical calculations on the stability of different tautomers of 6- and 7-phenacetyl pterins. Pteridines 1993;4:23–6. Search in Google Scholar
8. Estelberger W, Fuchs D, Murr C, Wachter H, Reibnegger G. Conformational investigation of the cofactor (6R,1′R,2′S-)-5,6,7,8-tetrahydrobiopterin. BBA 1995;1249:23–8. Search in Google Scholar
9. Reibnegger G, Horejsi R, Oettl K, Mlekusch W. Electronic charge density and electrostatic potential of pterin, 7,8-dihydropterin and 5,6,7,8-tetrahydropterin – an ab initio quantum chemical study. Pteridines 1998;9:85–90. Search in Google Scholar
10. Reibnegger G, Pauschenwein J, Werner ER. Electronic structure of tetrahydropteridine derivatives. Pteridines 2000;10:91–4. Search in Google Scholar
11. Jaramillo P, Coutinho K, Canuto S. Continuum, discrete and explicit solvation models for describing the low-lying absorption spectrum of the pterin acid in aqueous environment. Int J Quantum Chem 2010;110:2371–7. Search in Google Scholar
12. Soniat M, Martin CB. Theoretical study on the relative energies of neutral pterin tautomers. Pteridines 2008;19:120–4. Search in Google Scholar
13. Soniat M, Martin CB. Theoretical study on the relative energies of anionic pterin tautomers. Pteridines 2009;20:124–9. Search in Google Scholar
14. Soniat M, Martin CB. Theoretical study on the relative energies of cationic pterin tautomers. Pteridines, in press. Search in Google Scholar
15. Ferré J, Jacobson KB, Pfleiderer W. Proposal towards a normalization of pteridine nomenclature. Pteridines 1990;2:129–32. Search in Google Scholar
16. Becke AD. A new mixing of Hartree-Fock and local density-functional theories. J Chem Phys 1993;98:1372–7. Search in Google Scholar
17. Lee C, Yang W, Parr RG. Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density. Phys Rev B 1988;37:785–9. Search in Google Scholar
18. Bader RW. Atoms in molecules: a quantum theory. Oxford, UK: Oxford University Press, 1990. Search in Google Scholar
19. Hanwell MD, Curtis DE. Lonie DC, Vandermeersch T, Zurek E, Hutchison GR. Avogadro: an advanced semantic chemical editor, visualization and analysis platform. J Cheminformatics 2012;4:17. http://www.jcheminf.com/content/4/1/17. Search in Google Scholar
20. Foresman JB, Frisch AE. Exploring chemistry with electronic structure methods, 2nd Ed. Pittsburgh, PA: Gaussian Inc., 1996. Search in Google Scholar
21. Marenich AV, Cramer CJ, Truhlar DG. Universal solvation model based on solute electron density and on a continuum model of the solvent defined by the bulk dielectric constant and atomic surface tensions. J Phys Chem B 2009;113:6378–96. Search in Google Scholar
22. Lewis GN. The atom and the molecule. J Am Chem Soc 1916;38:762–85. Search in Google Scholar
©2014 by Walter de Gruyter Berlin/Boston
This article is distributed under the terms of the Creative Commons Attribution Non-Commercial License, which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.