Baixe Study Material on Molecular Structures and Physics from Various Sources e outras Notas de estudo em PDF para Engenharia Elétrica, somente na Docsity!
Density Functional Theory and Isodesmic Reaction Based Prediction of Four Stepwise
Protonation Constants, as log K H( n ), for Nitrilotriacetic Acid. The Importance of a Kind and
Protonated Form of a Reference Molecule Used
Krishna K. Govender and Ignacy Cukrowski*
Department of Chemistry, Faculty of Natural Sciences, Uni V ersity of Pretoria, Lynnwood Road, Hillcrest, Pretoria 0002, South Africa
Recei V ed: July 29, 2009; Re V ised Manuscript Recei V ed: No V ember 25, 2009
An explicit application of isodesmic reaction (a proton exchange between the studied and structurally similar reference molecule), where the free energy change of the protonation reaction in water was obtained using the free energies in solution from a single continuum model, was used to predict stepwise protonation constants of nitrilotriacetic acid. Calculations were performed at the RB3LYP/6-311+G(d,p) level of theory in conjunction with the PCM-UA0 solvation model. Five reference molecules were investigated. It has been established that one must pay special attention to structural similarities between the studied and reference molecules and selection of a protonated form of the reference molecule. The protonation reactions in which the studied and reference molecule are involved in must be (if possible) of the same order; e.g., the first (or generally n th) protonation reaction of the reference molecule must be used to compute the first (or n th) protonation constant of the studied molecule. The lowest energy conformer must always be used. The first, second, third, and fourth computed protonation constants differed, on average, from experimental values by 3.3, 0.8, 0.2, and 0.2 log units, respectively. It appears that the charge on the reference molecule has more decisive influence on the accuracy of computed protonation constants than its structural differences when compared with the studied molecule. Results reported can be used as a guide in constructing isodesmic reactions useful for the theoretical prediction of protonation constants by use of methodology described in this work.
1. Introduction
Knowledge of protonation, K H, and dissociation, K a, constants is of special interest to many chemists and life scientists 1 as they constitute important thermodynamic property of a com- pound that might be of biological, medicinal, or industrial (just to mention few) importance. Although a number of experimental techniques has been developed to measure protonation/dissocia- tion constants under various experimental conditions, many of the chemical species are not easily amenable to a full experi- mental characterization. 2 A number of papers has been reported 1 -^54 on theoretical prediction of dissociation constants. Most of them employed thermodynamic cycles (TC) to compute the free energies of dissociation reaction. Often, high-level theories were used in the gas-phase calculations where they are known to be accurate. The solution-phase calculations were used to provide the solvation energies ( G sol ); usually low-level continuum models were employed for the purpose. When the above protocol is used, the absolute p K a value is obtained. To avoid uncertainties related to the solvation energies of either H+^ or H 3 O+^ ions, isodesmic reaction (IRn) was incorporated within TC;1,3-^7 this protocol of calculation results in relative^3 p K a values. Results reported to date predominantly describe the calculations of singly charged molecules, either anions^1 -^29 (a study of doubly charged anions is very rare) or cations.^30 -^35 This is most likely due to the fact that (i) it is very difficult for DFT methods to properly describe anions (with multiple negative charges) in the gas phase because in the absence of an external stabilization of the charge (e.g., solvent) DFT methods have a bias toward “over-delocalization” of the charge (one
might observe bonds that are longer than expected and significant reduction on the HOMO-LUMO gap)^55 and (ii) computational evaluations of ionic solvation free energies for highly charged anions are inaccurate (these energies are highly dependent on the solvation model used due to different models chosen to generate the “best” electrostatic cavity).^56 Accuracies achieved thus far for computed dissociation constants (for a singly dissociable organic acids) are often within (1.0 log unit, on average, when compared with experimentally available values, but differences of several log units are not uncommon.^3 - 5, Recently we reported the DFT-predicted four stepwise protonation constants, expressed as log K H( n )^ , for a highly charged molecule nitrilotripropanoic acid (NTPA).^54 An explicit application of an isodesmic reaction involving two structurally similar ligands, where the free energy change of the protonation reaction in solution was obtained using the free energies in solution from a single continuum model, resulted in the average difference between predicted and experimental stepwise pro- tonation constants being (0.5 log unit. This suggested that, in principle, even though serious concerns were expressed,55, accurate determination of stepwise protonation constants for highly negatively charged molecules is possible. In this paper our focus is on parameters that influence the accuracy in predicting four consecutive protonation constants when the IRn-based procedure, as reported by us recently, 54 is employed. It is important to investigate a wide range of polycharged compounds (with negative and positive charges) to establish (i) whether implementation of the protocol can indeed produce consistently good predictions, (ii) how significant the selection of a reference molecule is from the point of view
- Corresponding author. E-mail: [email protected]. of its structural similarity to the studied compound, (iii) to what
J. Phys. Chem. A XXXX, xxx, 000 A
10.1021/jp9092964 XXXX American Chemical Society
degree the selection of a conformer impacts the accuracy in computed protonation constants, (iv) if the selection of a different protonated form of the reference molecule has an influence on the accuracy of computed protonation constants, (v) which one, structural similarity or the charge on a protonated form of a reference molecule, plays a more important role when the IRn-based methodology described in this work is employed, (vi) if it is possible to predict protonation constants in correct order, as determined from an experiment, (vii) to what extent the selection of level of theory and a basis set has an influence on the accuracy in computed protonation constants, and (viii) how significant is the kind of solvation model employed. These are just a few of the important questions, and this work will address only some of them. We have chosen the ligand nitrilotriacetic acid (NTA) because it is an important derivative of glycine that is widely studied due to its excellent chelating abilities.^57 This is a ligand that has enjoyed numerous applications in medicine,58,59^ biochemistry,58,59^ and industry.^60 -^62 In medicinal and biological studies it was shown that aliphatic amine salts of NTA inhibit the growth of bacteria and fungi and have herbicidal activity.^63 NTA has also been used as a transient phytoextraction agent that combines high biodegradability and low phytotoxicity with chelating strength.^64
The application of IRn requires a reference molecule and its selection appears to be crucial. Molecules shown in Figure 1 (iminodiacetic acid (IDA), methyliminodiacetic acid (MIDA), ethyliminodiacetic acid (EIDA), propyliminodiacetic acid (PIDA), and hydroxyethyliminodiacetic acid (HIDA)) were used here as they have many structural similarities with the compound of interest (NTA) and their experimental stepwise protonation constants are well-known. 57
2. Computational Details
All calculations were performed using GAUSSIAN 03, revision D.01,^65 on a 64-bit Linux workstation in parallel environment (Opensuse 10.3). Molecular visualizations were accomplished with the aid of GaussView 4.^66 Since it is of paramount importance to include diffuse functions for anions,^19 both gas-phase and solvent (water, ε ) 78.39) optimizations were performed at the RB3LYP level of theory^67 in conjunction with a 6-311+G(d,p) basis set. Full solvent-optimization was performed with the default solvation model provided by Gaussian, i.e., Tomasi’s polarized continuum model (PCM),^68 -^70 and UA0 radii (united atom topological model). We have chosen this model because it generated acceptable results in the
Figure 1. Top view of the ligands (in fully protonated forms) discussed in this work.
B J. Phys. Chem. A, Vol. xxx, No. xx, XXXX Govender and Cukrowski
Preliminary Investigations. Initially, we have used TCs to compute absolute protonation constants. A detailed description of TCs used and results obtained are provided in Supporting Information (see Tables S5-S7). Unfortunately, the methodol- ogy based on TC principles could not be fully applied in our studies due to significant structural differences of molecules in the gas phase and solvent. Migration of a proton from the central N-atom to the - COO-^ group took place when optimization was performed in the gas phase; a similar observation was also reported elsewhere.13, We have also tested whether a higher level of theory could result in preserving the solvent structure. The H 3 L form of NTA was subjected to the full gas-phase optimization at the RMP2/ 6-311+G(d,p) level of theory. Unfortunately, the proton mi- grated again from the N-atom to the carboxylic group. Clearly, even with a higher level of theory that generates accurate gas- phase free energies, we were unable to apply commonly used TCs in the study of molecules considered in this work. Isodesmic Reaction. To date, isodesmic reaction principles (where two structurally similar compounds were used, inves- tigated and reference molecule) have been extensively utilized in the prediction of enthalpies of formation. 76 -^84 In some cases a reference molecule has been incorporated within TCs1,3-^7 to eliminate uncertainties related to either H+^ or H 3 O+^ ions. Dissociation constants of a number of compounds with a single
dissociable proton have also been computed directly in solvent from ∆ G diss for the dissociation reaction HL + H 2 O ) L-^ + H 3 O+.^23 The implementation of IRn has an advantage, when the total free energies in solution from a single continuum model calculation are used, because it should minimize (or systemati- cally eliminate) errors related to the solvation models used, provided that the same level of theory, basis set, and solvation model are used for each component involved in the reaction of interest.^56 The main challenge associated with the use of IRn, however, is the selection of appropriate reference molecule.56, The ligand NTA can be seen as a set of five molecules that differ in (i) a number of protons (from 0 to 4) and (ii) charges on the molecule (from +1 to - 3). Clearly, a careful selection must be made to find the most appropriate protonated form of the reference molecule that has to be included in each of IRn needed to compute four protonation constants of NTA. By considering the structural properties of NTA (called further L(1)), we opted for IDA, MIDA, EIDA, PIDA, and HIDA as reference compounds (called further L(2) ) because each of them has two acetate groups (there are three in NTA) and the same kind of electron donor atoms (-COO-^ and R 3 N:) that can be protonated in a solution. In addition, protonation constants for all of the chosen reference molecules are well-known as these ligands are widely studied.^57
Figure 4. Self-constructed protonated forms of NTA and a crystal structure H 3 L* fully optimized at the RB3LYP/6-311+G(d,p) level of theory in solvent (PCM/UA0).
D J. Phys. Chem. A, Vol. xxx, No. xx, XXXX Govender and Cukrowski
IRn employed here can be seen as a competition reaction between two ligands for a proton (proton transfer reaction), and for the first protonation constant of NTA it can be written as
To investigate the impact the kind and protonated form of the reference molecule has on the theoretically generated protonation constants of NTA, a large number of isodesmic reactions was tested. Here, each reference molecule has three protonation constants (NTA has four), hence for each pair of ligands (L(1) and L (^) (2) ) 12 isodesmic reactions (such as eq 1, but involving different protonated forms of the ligands) had to be considered. For simplicity, only the first protonation reaction (PRn) in which each of the two ligands (NTA and a reference molecule) is involved is shown as eqs 2 and 3
A protonation reaction is the reverse of a weak acid dissociation reaction (DRn) and in the case of stepwise reactions the following holds
where k ) 1 + m - n , m and n represent the highest dissociation constant (here m ) 4) and an n th consecutive dissociation constant (1 e n e m ), respectively, and k applies to a k th consecutive protonation constant, 1 e k e m. Note that the ligand NTA has three acidic groups and only three dissociation constants would be reported when, for example, the TC-based methodology was employed. However, due to the protonation/ deprotonation of the N-atom in NTA, it is of paramount importance to consider also the first protonation constant, log K H(1)^. From this it follows that the fourth dissociation constant of NTA is linked through eq 4 with the first protonation constant of this ligand. In the case of the reference molecule, the third dissociation constant is linked through eq 4 with the first protonation constant of this reference molecule. (The above is provided here for convenience and to ensure clarity in nomen- clature used in this work.) The relationships between the change in Gibbs energies for protonation and dissociation reactions applicable to eqs 2 and 3 can be written as
where ∆ G PRn (aq)L(1)( k )^ and ∆ G DRn (aq)L(1)( n )^ refer to the k th (here k ) 1) stepwise protonation reaction and relevant n th (here n ) 4 for NTA) stepwise dissociation reaction, respectively (the same applies to the reference molecule L (^) (2) , but n ) 3). The change in Gibbs energies for each PRn, eqs 2 and 3, can be written, respectively, as
The isodesmic reaction of interest (eq 1) can be obtained by subtracting eq 3 from eq 2, and hence by subtracting eq 8 from eq 7 one obtains expressions for the change in Gibbs energy applicable to this isodesmic reaction
where the uncertainty related to G aq(H+) is no longer applicable as this term cancels out (this eliminates any error that might have been introduced by the use of an experimental value for this quantity). Equation 10 was used to calculate ∆ G IRn (aq) of IRn (eq 1) from appropriate Gibbs energies obtained for relevant and fully solvent-optimized structures of the ligand NTA and reference molecule L (^) (2). Table S8 in Supporting Information provides the ZPVE-corrected minimum energies, E min , as well as the Gibbs free energies of NTA and all of the reference molecules studied here. The value for ∆ G PRn (aq)L(2) was obtained from the well- known relationship
using the reported protonation constants^57 (at 20 and 25 °C, μ ) 0.0 and 0.1 M) of the reference molecule L(2). Once ∆ G IRn(aq) and ∆ G PRn(aq)L(2) have been calculated, the value of ∆ G PRn(aq)L(1), which is needed to calculate the protonation constants of NTA from eq 11, was obtained from eq 9. Table 1 provides values for the functions required to calculate protonation constants, calculated and experimental protonation constants of NTA, along with differences between calculated and experimental protonation constants (δ). Values for ∆ G PRn (aq)L(2)( k )^ were calculated from the experimentally available stepwise protonation constants, 57 which have been reproduced in Table S9, Supporting Information. Only isodesmic reactions that produced the best results are shown in Table 1; the remaining results are provided in Table S10, Supporting Information. There are several interesting observations one can make from the analysis of data seen in Table 1. (i) All isodesmic reactions seen in Table 1 predicted four protonation constants in correct order, log K H(1)^ > log K H(2)^ > log K H(3)^ > log K H(4), as observed from an experiment. One must realize that the experimental values of the second, third, and fourth protonation constants of NTA differ between each other only by one log unit (or less) and this is a typical error reported in theoretically predicted values of dissociation constants reported to date for organic acids containing only one carboxylic group. Clearly, results obtained here should be seen as satisfac- tory and encouraging for further studies aimed at improvements in accuracy of theoretical predictions. (ii) As one goes from IDA to HIDA, it is seen that the prediction of protonation constants becomes more accurate, with HIDA yielding the most accurate estimates. The observed trend
L(1)(aq) + HL(2) (aq) ) HL(1)(aq) + L(2)(aq) ∆ G IRn(aq) (1)
H+^ + L(1) T HL(1) ∆ G PRn (aq)L(1) (2)
H+^ + L(2) T HL(2) ∆ G PRn (aq)L(2) (3)
log K H( k )^ ) log
K a( n )^
) p K a( n )^ (4)
∆ G PRn(aq)L(1)(1)^ ) -∆ G DRn(aq)L(1)(4)^ (5)
∆ G PRn(aq)L(2)(1)^ ) -∆ G DRn(aq)L(2)(3)^ (6)
∆ G PRn(aq)L(1) ) G aq (HL(1)) - G aq (H+) - G aq (L(1) ) (7)
∆ G PRn(aq)L(2) ) G aq (HL(2)) - G aq (H+) - G aq (L(2) ) (8)
∆ G IRn(aq) ) ∆ G PRn(aq)L(1) - ∆ G PRn(aq)L(2) (9)
∆ G IRn (aq) ) G aq (HL(1)) - G aq (L(1) ) - G aq(HL(2)) + G aq (L(2)) (10)
∆ G (aq) ) - RT ln K (11)
Protonation Constants for Nitrilotriacetic Acid J. Phys. Chem. A, Vol. xxx, No. xx, XXXX E
molecules. Also, careful attention needs to be paid to the positioning and presence of atoms, especially the heteroatom, as the additional - OH group present in HIDA, which is not present in IDA, MIDA, EIDA, or PIDA, seems to make a huge difference, as far as prediction of protonation constants is concerned. Additional Test. All four theoretically predicted protonation constants of NTA seen in Table 1 differ from experimental values by less than a single log unit when HIDA was used; results of this accuracy are often referred to as excellent in the literature. On the other hand, the absolute protonation constants obtained from TCs carried significantly larger errors (see Table S7, Supporting Information). One might argue that accuracy obtained from straight continuum model free energies in solution calculation is not directly comparable with that obtained from TC as totally different chemical reactions are involved. Since some solvent structures were not preserved in the gas phase, we performed a single point frequency calculation in the gas phase on the solvent-optimized structures. We selected IDA for the test, as inclusion of this reference molecule in the proton competition reaction (IRn) resulted in the worst results; hence one might expect significant improvement in protonation constants calculations. The test was performed for the reaction
where L (^) (1) and L (^) (2) are NTA and the reference molecule IDA, respectively. Our aim here was to find out whether improvement of the G gas component to a better (higher) level of theory improves the computed protonation constant when the TC-based protocol is implemented. The following protocol was imple- mented. Let G sol stand for the Gibbs free energy in solution obtained from a single continuum model calculation at B3LYP/ 6-311+G(d,p) in conjunction with the PCM/UA0 solvation model. For each component involved in the above IRn (eq 12) one can compute the free energy of solvation from ∆ G sol ) G sol - G gas(1)^ , where G gas(1)^ is computed as a single point on the solution-phase geometries at the same B3LYP/6-311+G(d,p) level of theory. Let G gas(2)^ stand for the improved estimate obtained at the higher level of theory, in this case also calculated as a single point on the solution-phase geometry, then G sol(TC), the improved Gibbs free energy obtained from thermodynamic cycle calculations, can be expressed as G sol (TC) ) G gas(2)^ + ∆ G sol. The values of G gas(2)^ obtained at the RMP2/6-311+G(d,p) level of theory along with other necessary components are presented in Table 2. From G sol (TC) values we calculated the change in the Gibbs free energy for reaction 12, ∆ G sol (TC) ) 6.137 kcal mol-^1 , followed by the fourth protonation constant of NTA, log K H(4)^ ) -4.50. This is a much worse estimate when compared with the value of - 1.19 obtained directly from implementation of IRn (see Table 1). An attempt was made to obtain G gas(2)^ values using the G3/ 6-311+G(d,p) level of theory for the single point frequency
calculations in the gas phase on the structures optimized at the B3LYP/6-311+G(d,p) level of theory in conjunction with the PCM/UA0 solvation model. Unfortunately, the solvent structures of H3 L(1) and H 2 L(2) have not been preserved. The H-atom has moved away from nitrogen and protonated the - COO-^ group; see Figure S8, Supporting Information. Conformational Considerations. All protonated forms of NTA and the reference molecules were subjected to the Schro¨dinger’s Maestro^85 conformational analysis in solvent. Generated, on the basis of molecular mechanics/molecular dynamics (MM/MD) principles, lowest energy conformers C- are shown in Figures S9-S14, Supporting Information. Table 3a provides energies E C-1 to E C-5 (in kJ/mol) for the five lowest energy conformers of all the protonated forms of NTA. MM/ MD-based SPC in solvent was also performed on the NTA DFT- structures shown in Figure 4; results are shown as E SPC in Table 3a. It is seen that all E SPC are larger than energies of the MM/ MD-generated conformers. Initially, all of the C-1 conformers of NTA were fully DFT-optimized in solvent; results obtained are shown in Table 3b. It was gratifying to note that the differences δ G in Table 3b between the relevant structures became almost negligibly small. For all protonated forms of NTA, except H 4 L+, the value of δ G is about (0.1 kcal mol-^1. This is equivalent to about (0.07 log unit of the computed protonation constant, a typical experimental uncertainty. A similar procedure was applied to all the protonated forms of IDA, MIDA, EIDA, PIDA, and HIDA (Table S11, Supporting Information). It is seen in Table 3b that only two out of six of the DFT-optimized C-1 conformers have lower energies (HL^2 - and H 2 L-) when compared with the energies of the original structures seen in Figure 4. A similar observation applies to energies of the DFT-optimized C-1 conformers of the reference molecules. The lowest Gibbs free energies G aq were used to calculate protonation constants; results are provided in Table S12, Supporting Information. Because the use of C-1 conformers has not resulted in any significant change in the computed protonation constants and since there was no direct correlation between energies of the DFT-optimized C-1 conformers and original structures, we have decided to use six lowest energy MM/MD-generated conformers of NTA and HIDA for further studies. HIDA was selected because this reference molecule, when used in IRn’s seen in Table 1, generated best estimates in protonation constants. Energy-minimized in solvent conformers C-1 to C-6 of all the protonated forms of NTA and HIDA are shown in Figures S15-S23, Supporting Information. Their ZPVE-corrected ener- gies ( E min ) and Gibbs free energies ( G aq ) along with the values of the self-constructed (S-c) structures are shown in Table 4; for each protonated form of L(1) and L(2) , the lowest E min and G aq are printed in bold. Implementation of an extended con- formational analysis resulted in three new lower in energy conformers of NTA (H 2 L-, H 3 L, and H 4 L+) as well as HIDA (L^2 - , H 2 L, and H 3 L+).
TABLE 2: Thermochemical Data Used To Calculate the Fourth Protonation Constant of NTA Involving IRn in Eq 12 a
molecule G sol , au G gas(1)^ , au ∆ G sol , au G gas(2)^ , au G sol (TC), au H 3 L (^) (1) - 740.332226 - 740.252149 - 0.080077 - 738.373664 - 738. H 3 L (^) (2)+^ - 512.862287 - 512.727315 - 0.134972 - 511.414007 - 511. H 4 L (^) (1)+^ - 740.764545 - 740.625849 - 0.138696 - 738.746987 - 738. H 2 L (^) (2) - 512.423365 - 512.349858 - 0.073507 - 511.03375 - 511. a (^) G gas(1) (^) values were obtained at the B3LYP/6-311+G(d,p) level of theory. G gas(2) (^) values were obtained at the RMP2/6-311+G(d,p) level of
theory. G sol values were obtained at the B3LYP/6-311+G(d,p) level of theory in solvent (PCM/UA0). L (^) (1) ) NTA, L (^) (2) ) IDA. For details see the text.
H 3 L(1) + H 3 L(2)^ +^ + H 4 L(1)^ +^ + H 2 L(2) (12)
Protonation Constants for Nitrilotriacetic Acid J. Phys. Chem. A, Vol. xxx, No. xx, XXXX G
All lowest DFT-computed G aq values for NTA and HIDA were used to calculate protonation constants; see the first row in Table 5. One can note from Table 4 that in two cases involving HIDA (L^2 -^ and H 2 L) the conformers with the lowest energy, E min , do not have the lowest G aq value. It has been decided to include these G aq values (printed in Italic in Table
- in computing protonation constants to estimate to what degree they influence the final result; see second row in Table 5. When one assumes a study of a new ligand, conformational search most likely would be the first analysis performed to find the lowest energy conformer. To simulate this approach, the self- constructed molecules were excluded in search of lowest energy
TABLE 3: (a) Minimum Energies of MM/MD-Generated Conformers, C-1 to C-5, in Solvent ( E C-1 to E C-5) and Energies Obtained from MM-Based SPC ( E SPC) Performed on the NTA Structures Seen in Figure 4 and (b) DFT-Calculated Solvent-Optimized Energies ( E min ) ZPVE-Corrected Energy) of All Protonated Forms of the Ligand NTA Seen in Figure 4 and Lowest Energy C-1 Conformers
(A) Minimum Energies (kJ mol-^1 ) L ) NTA E SPC , kJ/mol E C-1 , kJ/mol δ E , a^ kJ/mol δ E , a^ kcal/mol E C-2 , kJ/mol E C-3 , kJ/mol E C-4 , kJ/mol E C-5 , kJ/mol L 3 -^ - 770.81 - 799.18 28.37 6.78 - 799.17 - 798.69 - 798.69 - 796. HL 2 -^ - 1190.19 - 1216.05 25.86 6.18 - 1216.05 - 1201.93 - 1201.93 - 1201. H 2 L-^ - 991.30 - 1027.95 36.65 8.76 - 1027.95 - 1027.67 - 1027.67 - 1025. H 3 L - 758.87 - 803.94 45.07 10.77 - 803.94 - 802.89 - 802.89 - 802. H 3 L* - 750.22 - 803.95 54.38 13.00 - 803.94 - 802.89 - 802.02 - 802. H 4 L+^ - 489.11 - 543.49 53.73 12.84 - 543.49 - 539.35 - 538.46 - 538. (b) DFT-Calculated Solvent-Optimized Energies structures seen in Figure 4 C-1 structures L ) NTA E min , hartree G aq , hartree E min , hartree G aq , hartree δ G (aq), b^ hartree δ G (aq), b^ kcal/mol L 3 -^ - 738.936972 - 738.978069 - 738.936883 - 738.977928 - 0.000141 - 0. HL 2 -^ - 739.410772 - 739.451765 - 739.410827 - 739.451947 0.000182 0. H 2 L-^ - 739.852748 - 739.893949 - 739.852711 - 739.894128 0.000179 0. H 3 L - 740.290659 - 740.332226 - 740.290552 - 740.332033 - 0.000193 - 0. H 3 L* - 740.293694 - 740.335141 - 740.290657 - 740.332261 - 0.000182 - 0. H 4 L+^ - 740.722681 - 740.764545 - 740.722433 - 740.764363 - 0.002880 - 1. a (^) δ E ) E SPC - E C-1. b (^) δ G aq ) G aq (self-constructed structure) - G aq (C-1).
TABLE 4: Energies E min ( E min ) ZPVE-Corrected Energy) and G aq Obtained for All Protonated Forms of the Self-Constructed Structures and Conformers C-1 to C-6 That Were Fully Energy-Optimized in Gaussian at the RB3LYP/6-311 + G(d,p) Level of Theory in Solvent (PCM-UA0)
(a) L ) NTA L 3 -^ HL 2 -^ H 2 L-^ H 3 L H 4 L+ E min , au G aq , au E min , au G aq , au E min , au G aq , au E min , au G aq , au E min , au G aq , au S-c - 738.936972 - 738.978069 - 739.410772 - 739.451765 - 739.852748 - 739.893949 - 740.290659 - 740.332226 - 740.722681 - 740. C1 - 738.936883 - 738.977928 - 739.410827 - 739.451947 - 739.852711 - 739.894128 - 740.290552 - 740.332033 - 740.722433 - 740. C2 - 738.936941 - 738.977965 - 739.410753 - 739.451937 - 739.852772 - 739.894052 - 740.290599 - 740.332198 - 740.722408 - 740. C3 - 738.934628 - 738.976228 - 739.407682 - 739.448582 - 739.854620 - 739.895947 - 740.292742 - 740.334330 - 740.722671 - 740. C4 - 738.934684 - 738.976365 - 739.410791 - 739.451885 - 739.852733 - 739.893980 - 740.290717 - 740.332398 - 740.727700 - 740. C5 - 738.931462 - 738.973284 - 739.402567 - 739.447349 - 739.846053 - 739.888621 - 740.294578 - 740.336185 - 740.729966 - 740. C6 - 738.934569 - 738.976095 - 739.852702 - 739.894193 - 740.292501 - 740.334200 - 740.723142 - 740. (b) L ) HIDA L 2 -^ HL-^ H 2 L H 3 L+ E min , hartree G aq , hartree E min , hartree G aq , hartree E min , hartree G aq , hartree E min , hartree G aq , hartree S-c - 665.299655 - 665.340516 - 665.767828 - 665.808578 - 666.207584 - 666.248534 - 666.642422 - 666. C1 - 665.305277 - 665.344936 - 665.766303 - 665.807676 - 666.205011 - 666.246706 - 666.642406 - 666. C2 - 665.305221 - 665.344973 - 665.765839 - 665.807813 - 666.204984 - 666.246766 - 666.642341 - 666. C3 - 665.300319 - 665.339783 - 665.765002 - 665.805790 - 666.208044 - 666.248784 - 666.642415 - 666. C4 - 665.298207 - 665.340665 - 665.765071 - 665.805962 - 666.206880 - 666.248997 - 666.644447 - 666. C5 - 665.301309 - 665.341919 - 665.765305 - 665.806013 - 666.205287 - 666.246723 - 666.641862 - 666. C6 - 665.298902 - 665.340283 - 665.765156 - 665.805658 - 666.204193 - 666.245406 - 666.643550 - 666.
TABLE 5: Computed Stepwise Protonation Constants (as log K H( n )) of NTA Obtained from Isodesmic Reactions Involving HIDA (Table 1) a
structures considered K H(1)^ δ K H(2)^ δ K H(3)^ δ K H(4)^ δ C-1 to C-6 and S-c b^ 13.40 3.07 3.85 0.91 2.12 0.12 1.26 0. C-1 to C-6 and S-c c^ 13.39 3.06 3.94 1.00 2.21 0.21 1.16 0. C-1 to C-6 b^ 13.80 3.47 3.50 0.56 1.76 - 0.24 1.26 0. C-1 to C-6 c^ 13.85 3.52 3.53 0.59 1.80 - 0.20 1.16 0. S-c b^ 11.27 0.94 3.22 0.28 1.43 - 0.57 0.43 - 0. a (^) δ is the difference between the computed and experimental value. S-c stands for self-constructed structures. C-1 to C-6 represent MM/ MD-generated conformers that were fully energy optimized in Gaussian at the RB3LYP/6-311+G(d,p) level of theory in solvent (PCM-UA0). b (^) Minimum G aq were used. c (^) G aq for E min were used.
H J. Phys. Chem. A, Vol. xxx, No. xx, XXXX Govender and Cukrowski
involves highly charged anions (-3 and - 2 for the studied and the reference molecule, respectively), a significantly larger than experimental value (by 3 log units) was obtained. Clearly, more work has to be done to improve this value. Nevertheless, we are convinced that results reported here can be used as a guide in constructing isodesmic reactions that are useful for the theoretical prediction of protonation/dissociation constants, particularly for compounds for which solvent and gas structures differ significantly.
Acknowledgment. Financial support of the National Research Foundation of South Africa and the University of Pretoria is highly appreciated.
Supporting Information Available: Reported crystal struc- tures of IDA, MIDA, and HIDA; structural matrix of the solvent- optimized H 3 L* and H 3 L forms of NTA, H 2 L form of IDA and MIDA, and H 3 L form of HIDA; energy-minimized structures in solvent (PCM/UA0) at the RB3LYP/6-311+G(d,p) level of theory for all self-constructed protonated forms of IDA, MIDA, EIDA, PIDA, and HIDA; results obtained from computations involving thermodynamic cycles; ZPVE-corrected minimum and Gibbs free energies of all protonated forms of NTA and reference molecules (IDA, MIDA, EIDA, PIDA, and HIDA) obtained from full energy optimization at the RB3LYP/6- 311 +G(d,p) level of theory in conjunction with the PCM/UA solvation model; experimental stepwise protonation constants of IDA, MIDA, EIDA, PIDA, and HIDA at μ ) 0.0 and 0.1 M and 25 °C; comparison of experimental and calculated stepwise protonation constants of NTA using protonation constants of the reference molecules IDA, MIDA, EIDA, PIDA, and HIDA; output structures (H 3 L form of NTA and H 2 L form of IDA) obtained from the single point frequency calculations at the G3/ 6-311+G(d,p) level of theory on the solvent optimized mol- ecules optimized at the RB3LYP/6-311+G(d,p) level of theory in conjunction with the PCM/UA0 solvation model; energy- minimized structures, at the RB3LYP/6-311+G(d,p) level of theory in conjunction with the PCM/UA0 solvation model, for all protonated forms of C-1 conformers of NTA, IDA, MIDA, EIDA, PIDA, and HIDA; minimum energies of MM/MD- generated conformers in solvent and energies obtained from MM-based SPC performed on the IDA, MIDA, EIDA, PIDA, and HIDA structures; DFT-calculated solvent-optimized energies of all self-constructed protonated forms of IDA, MIDA, EIDA, PIDA, and HIDA and lowest energy MM/MD-generated C- conformers; comparison of experimental and calculated stepwise protonation constants of NTA obtained from the lowest energy MM/MD conformers of IDA, MIDA, EIDA, PIDA, and HIDA; energy-minimized structures in solvent (PCM/UA0) at the RB3LYP/6-311+G(d,p) level of theory for MM/MD-generated conformers of the protonated forms of NTA and HIDA. This material is available free of charge via the Internet at http:// pubs.acs.org.
References and Notes
(1) Namazian, M.; Kalantary-Fotooh, F.; Noorbala, M. R.; Searles, D. J.; Coote, M. C. J. Mol. Struct. (THEOCHEM) 2006 , 758 , 275. (2) Jacquemin, D.; Perpete, E. A.; Ciofini, I.; Adamo, C. J. Phys. Chem. A 2008 , 112 , 794. (3) Lim, C.; Bashford, D.; Karplus, M. J. Phys. Chem. 1991 , 95 , 5610. (4) Namazian, M.; Heidary, H. J. Mol. Struct. (THEOCHEM) 2003 , 620 , 257. (5) Namazian, M.; Halvani, S.; Noorbala, M. R. J. Mol. Struct. (THEOCHEM) 2004 , 711 , 13. (6) Namazian, M.; Halvani, S. J. Chem. Thermodyn. 2006 , 38 , 1495. (7) Namazian, M.; Zakery, M.; Noorbala, M. R.; Coote, M. L. Chem. Phys. Lett. 2008 , 451 , 163.
(8) Charif, I. E.; Mekelleche, S. M.; Villemin, D.; Mora-Diez, N. J. Mol. Struct. (THEOCHEM) 2007 , 818 , 1. (9) Liptak, M. D.; Shields, G. C. Int. J. Quantum Chem. 2001 , 85 ,
(10) Liptak, M. D.; Shields, G. C. J. Am. Chem. Soc. 2001 , 123 , 7314. (11) Liptak, M. D.; Gross, K. C.; Seybold, P. G.; Feldgus, S.; Shields, G. C. J. Am. Chem. Soc. 2002 , 124 , 6421. (12) Bensen, M. T.; Moser, M. L.; Peterman, D. R.; Dinescu, A. J. Mol. Struct. (THEOCHEM) 2008 , 867 , 71. (13) Sang-Aroon, W.; Ruangpornvisuti, V. Int. J. Quantum Chem. 2008 , 108 , 1181. (14) Jang, Y. H.; Goddard, W. A., III; Noyes, K. T.; Sowers, L. C.; Hwang, S.; Chung, D. S. J. Phys. Chem. B 2003 , 107 , 344. (15) Jitariu, L. C.; Masters, A. J.; Hiller, I. H. J. Chem. Phys. 2004 , 121 , 7795. (16) Bryantsev, V. S.; Diallo, M. S.; Goddard, W. A., III. J. Phys. Chem. A 2007 , 111 , 4422. (17) Gao, D.; Svoronos, P.; Wong, P. K.; Maddalena, D.; Hwang, J.; Walker, H. J. Phys. Chem. A 2005 , 109 , 10776. (18) Murlowska, K.; Sadlej-Sosnowska, N. J. Phys. Chem. A 2005 , 109 ,
(19) Saracino, G. A. A.; Improta, R.; Barone, V. Chem. Phys. Lett. 2003 , 373 , 411. (20) da Silva, C. O.; da Silva, E. C.; Nascimento, M. A. C. J. Phys. Chem. A 1999 , 103 , 11194. (21) Li, G.; Cui, Q. J. Phys. Chem. B 2003 , 107 , 14521. (22) Dahlke, E. E.; Cramer, C. J. J. Phys. Org. Chem. 2003 , 16 , 336. (23) Klamt, A.; Eckert, F.; Diedenhofen, M.; Beck, M. E. J. Phys. Chem. A 2003 , 107 , 9380. (24) Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. J. Phys. Chem. A 2006 , 110 , 2493. (25) Pliego, J. R., Jr.; Riveros, J. M. J. Phys. Chem. A 2002 , 106 , 7434. (26) Klicic, J. J.; Friesner, R. A.; Liu, S.; Guida, W. C. J. Phys. Chem. A 2002 , 106 , 1327. (27) Jang, Y. H.; Sowers, L. C.; Cagin, T.; Goddard, W. A., III. J. Phys. Chem. A 2001 , 105 , 274. (28) Himo, F.; Eriksson, L. A.; Blomberg, M. R. A.; Siegbahn, P. E. M. Int. J. Quantum Chem. 2002 , 76 , 714. (29) Himo, F.; Noodleman, L.; Blomberg, M. R. A.; Siegbahn, P. E. M. J. Phys. Chem. A 2002 , 106 , 8757. (30) Chipman, D. M. J. Phys. Chem. A 2002 , 106 , 7413. (31) Mujika, J. I.; Mercero, J. M.; Lopez, X. J. Phys. Chem. A 2003 , 107 , 6099. (32) Scharnagl, C.; Raupp-Kossmann., R. A. J. Phys. Chem. B 2004 , 108 , 477. (33) Lee, I.; Kim, C. K.; Lee, I. Y.; Kim, C. K. J. Phys. Chem. A 2000 , 104 , 6332. (34) Maksic, Z. B.; Vianello, R. J. Phys. Chem. A 2002 , 106 , 419. (35) Magill, A. M.; Cavell, K. J.; Yates, B. F. J. Am. Chem. Soc. 2004 , 126 , 8717. (36) Li, H.; Hains, A. W.; Everts, J. E.; Robertson, A. D.; Jensen, J. H. J. Phys. Chem. B 2002 , 106 , 3486. (37) da Silva, C. O.; da Silva, E. C.; Nascimento, M. A. C. J. Phys. Chem. A 2000 , 104 , 2402. (38) da Silva, C. O.; da Silva, E. C.; Nascimento, M. A. C. Chem. Phys. Lett. 2003 , 381 , 244. (39) Alvarado-Gonzalez, M.; Orrantia-Borunda, E.; Glossman-Mitnik, D. J. Mol. Struct. (THEOCHEM) 2008 , 869 , 105. (40) Sulpizi, M.; Sprik, M. Phys. Chem. Chem. Phys. 2008 , 10 , 5238. (41) Tao, L.; Han, J.; Tao, F. M. J. Phys. Chem. A 2008 , 112 , 775. (42) Kieseritzky, G.; Knapp, E. W. J. Comput. Chem. 2008 , 29 , 2575. (43) da Silva, G.; Kennedy, E. M.; Dlugogorski, B. Z. J. Phys. Chem. A 2006 , 110 , 11371. (44) Ulander, J.; Broo, A. Int. J. Quantum Chem. 2005 , 105 , 866. (45) Gross, K. C.; Seybold, P. G. Int. J. Quantum Chem. 2000 , 80 , 1107. (46) Gross, K. C.; Seybold, P. G. Int. J. Quantum Chem. 2001 , 85 , 569. (47) Pliego, J. R., Jr.; Riveros, J. M. Chem. Phys. Lett. 2000 , 332 , 597. (48) Pliego, J. R., Jr. Chem. Phys. Lett. 2003 , 367 , 145. (49) Adam, K. R. J. Phys. Chem. A 2002 , 106 , 11963. (50) Kim, Y. J.; Streitwieser, A. J. Am. Chem. Soc. 2002 , 124 , 5757. (51) Schuurmann, G.; Cossi, M.; Barone, V.; Tomasi, J. J. Phys. Chem. A 1998 , 102 , 6706. (52) Soriano, E.; Cerdan, S.; Ballesteros, P. J. Mol. Struct. (THEOCHEM) 2004 , 684 , 121. (53) Sadlej-Sosnowska, N. Theor. Chem. Acc. 2007 , 118 , 281. (54) Govender, K. K.; Cukrowski, I. J. Phys. Chem. A 2009 , 113 , 3639. (55) Private communication. Gaussian Technical Support. (56) Cramer, C. J. In Essentials of Computational Chemistry ; John Wiley & Sons Ltd.: New York, 2002; pp 371-372. (57) NIST Standard Reference Database 46. NIST Critically Selected Stability Constants of Metal Complexes Database; Version 8.0; Smith, R. M., Martell, A. E., Eds.; US Department of Commerce, National Institute of Standards and Technology: Gaithersburg, MD, 2004.
J J. Phys. Chem. A, Vol. xxx, No. xx, XXXX Govender and Cukrowski
(58) Grigoriev, M. S.; Auwer, C. D.; Meyer, D.; Moisy, P. Acta Crystallogr. 2006 , C62 , m163. (59) Yu, L. C.; Liu, S. L.; Liang, E. X.; Wen, C. L. J. Coord. Chem. 2007 , 60 , 2097. (60) Malevich, D.; Wang, Z.; Tremaine, P. R. J. Sol. Chem. 2006 , 35 ,
(61) Vohra, M. S.; Davis, A. P. J. Colloid Interface Sci. 1997 , 194 , 59. (62) Karhu, J.; Harju, L.; Ivaska, A. Anal. Chim. Acta 1999 , 380 , 105. (63) Souaya, E. R.; Hanna, W. G.; Ismail, E. H.; Milad, N. E. J. Coord. Chem. 2004 , 57 , 825. (64) Quartacci, M. F.; Irtelli, B.; Baker, A. J. M.; Navari-Izzo, F. Chemosphere 2007 , 68 , 1920. (65) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Pettersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03 , revision D.01; Gaussian, Inc.: Wallingford, CT, 2004. (66) GaussView 4.1.2 ; Gaussian Inc.; Wallingford, CT, 2004. (67) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994 , 98 , 11623.
(68) Cances, E.; Mennucci, B.; Tomasi, J. J. Chem. Phys. 1997 , 107 ,
(69) Cossi, M.; Barone, V.; Cammi, R.; Tomasi, J. Chem. Phys. Lett. 1996 , 255 , 327. (70) Miertus, S.; Scrocco, E.; Tomasi, J. Chem. Phys. 1981 , 55 , 117. (71) Barone, V.; Cossi, M. J. Phys. Chem. A 1998 , 102 , 1995. (72) Cossi, M.; Rega, N.; Scalmani, G.; Barone, V. J. Comput. Chem. 2003 , 24 , 669. (73) Allen, F. H. Acta Crystallogr. 2002 , B58 , 380. (74) Skrzypczak-Jankun, E.; Smith, D. A.; Maluszynska, H. Acta Crystallogr. 1994 , C50 , 1097. (75) Kaneyoshi, M.; Bond, A.; Jones, W. Acta Crystallogr. 1999 , C55 ,
(76) da Silva, G.; Bozzelli, J. W. J. Phys. Chem. A 2006 , 110 , 13058. (77) da Silva, G.; Bozzelli, J. W. J. Phys. Chem. C 2007 , 111 , 5760. (78) Guo, Y.; Gao, H.; Twamley, B.; Shreeve, J. M. Ad V_. Mater._ 2007 , 19 , 2884. (79) Espinosa-Garcia, J.; Garcia-Bernaldez, J. C. Phys. Chem. Chem. Phys. 2002 , 4 , 4096. (80) Krossing, I.; Raabe, I. Chem. s Eur. J. 2004 , 10 , 5017. (81) Kutt, A.; Movchun, V.; Rodima, T.; Dansauer, T.; Rusanov, E. B.; Leito, I.; Kaljurand, I.; Koppel, J.; Pihl, V.; Koppel, I.; Ovsjannikov, G.; Toom, L.; Mishima, M.; Medebielle, M.; Lork, E.; Roschenthaler, G. V.; Koppel, I. A.; Kolomeitsev, A. A. J. Org. Chem. 2008 , 73 , 2607. (82) Nolan, E. M.; Linck, R. G. J. Am. Chem. Soc. 2000 , 122 , 11497. (83) Peeters, D.; Leroy, G.; Wilante, C. J. Mol. Struct. 1997 , 416 , 21. (84) Wang, L.; Heard, D. E.; Pilling, M. J.; Seakins, P. J. Phys. Chem. A 2008 , 112 , 1832. (85) Schrodinger, L. L. C. Schrodinger Maestro ; 32nd Floor, Tower 45, 120 West Forty-Fifth Street, New York, 10036, 2003. JP
Protonation Constants for Nitrilotriacetic Acid J. Phys. Chem. A, Vol. xxx, No. xx, XXXX K
