Polarizable and Non-Polarizable Force Field Representations of Ferric

Article pubs.acs.org/JPCB

Polarizable and Non-Polarizable Force Field Representations of Ferric Cation and Validations Miaoren Xia,† Zhifang Chai,†,‡ and Dongqi Wang*,† †

Multidisciplinary Initiative Center, Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China School of Radiation Medicine and Interdisciplinary Sciences (RAD-X), Soochow University, Suzhou 215123, China

‡

S Supporting Information *

ABSTRACT: The AMOEBA polarizable force field of ferric ion was optimized and applied to study the hydration of ferric ion and its complexation with porphine in the aqueous phase. The nonpolarizable force field was also optimized for comparison. The AMOEBA force field was found to give a more accurate hydration free energy than the nonpolarizable force field with respect to experimental data, and correctly predict the most stable electronic state of hydrated Fe3+, which is the sextet state, and of the Fe(III)− Por complex, which is the quartet state, consistent with the literature that was carried out using the DFT method. The explicit inclusion of charge transfer between Fe3+ and ligand was found to be important in order to obtain a precise picture of polarization energy and van der Waals energy, which otherwise deviate from the corresponding energy components derived from ab initio calculations. The successful application of the AMOEBA force field in the characterization of aquo Fe(III)−Por complexes suggests that its use may be extended to the study of the dynamics of metalloenzyme containing highly charged metal ions in the condensed phase with reliable treatment of the interactions between metal atom and protein.

■

INTRODUCTION Classical molecular dynamics simulation (MD) has evolved into a ubiquitous tool in studying the dynamic behavior of the condensed phase after decades of its birth. It applies force field parameters optimized in advance to describe the interactions between particles, and has contributed to the enrichment of our knowledge on the molecular mechanisms in the condensed phase. The widely used force fields, such as AMBER,1−3 CHARMM,4 GROMOS,5 and OPLS,6 handle individual atoms as a point charge with fixed atomic charge, and do not explicitly consider the polarization of the atoms as a response to external electric fields, exerted either by the other atoms in the surrounding or for special purposes. Such polarization becomes even more significant in the systems including transition metal ions,7 which in general carry highly positive charge concentrated on the metallic atoms, require sufficiently accurate treatment of electrostatics, and cast doubt on the application of the nonpolarizable force fields on these model systems. A remediation to classical force fields is the explicit treatment of the polarization effect to improve the calculation of electrostatics, which is among the most important issues in the simulation of realistic modeling.8 Various strategies have been proposed, and the development of a polarizable force field has been reviewed.8−13 Depending on the treatment of the permanent component of electrostatics, these approaches may be divided into two families:11 the first group, such as fluctuating charge models,14 methods based on Drude oscillators,15 and the induced-point polarization model,13 © XXXX American Chemical Society

calculates electrostatics on the basis of atomic point charge models (monopoles), which is isotropic, while the second one introduces multipoles, and represents a model system intrinsically with the ability to handle anisotropic properties, such as lone pairs and π clouds. The AMOEBA (atomic multipole optimized energetics for biomolecular applications)13 polarizable force field belongs to the second group, and it includes multipole (up to quadrupole) electrostatic interaction terms to take into account the many-body effect in polarization and aims at chemical accuracy in the prediction of conformational and interaction energies.16 Concerning the importance of ferric cation in biology,17−24 environmental science,25−28 and catalysis,29−35 in this work, we extended the AMOEBA polarizable force field to Fe3+ cation. Iron is known as the most abundant metal element next to aluminum in the Earth’s crust.36 Its significance in biology, chemistry, physics, environment, etc., has driven numerous studies to unravel its role in the complex processes for better understanding and manipulation of its behavior toward specific purposes. It is known that iron appears as a key component in metalloproteins, where its presence and its valence state determine the functionality of the target proteins, such as hemoproteins,37 transferrins,38,39 and ferritin,40,41 and thus plays critical roles in bioorganisms.42 Received: March 1, 2017 Revised: May 5, 2017 Published: May 16, 2017 A

DOI: 10.1021/acs.jpcb.7b02010 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

where Kb (kcal·mol−1·Å−2), Kθ (kcal·mol−1·rad−2), Kχ (kcal· mol−1·rad−2), and Knϕ (kcal/mol) are the force constants of bond length, angle, out-of-plane bending, and torsion; b0 (Å), θ0 (rad), and δ (rad) represent the equilibrium bond length, angle, and torsional phase angle. The nonbonded part contains the contributions of van der Waals (vdW) and electrostatics. The vdW term takes the buffered 14-7 potential form63

In the condensed phase, Fe often exists in its ferric and ferrous states. This governs its chemical behavior and biological effect. We note that, in earlier work, Gagliardi et al. reported the AMOEBA force field of the ferrous ion.43 Concerning the importance of ferric ion, we have optimized its AMOEBA force field in its doublet, quartet, and sextet spin states. The parameters were then validated and applied to study the aquo structure of ferric ion and its complex with porphine and their dynamics in aqueous solution. Apart from polarization, charge transfer (CT) is also considered important for accurate modeling of solvation of ions in the aqueous phase,44,45 while it is often ignored for simplicity in modeling. In order to evaluate the significance of the CT contribution at the AMOEBA level, we adopted the discrete charge-transfer (DCT) method according to Lee and Rick46 in the energy decomposition analysis and hydration free energy calculations of ferric cation.

Uvdw

where rij is the distance between atom i and atom j. The parameters εij (kcal/mol) and Rij0 (Å) are calculated by combination rules below:

■

εij =

COMPUTATIONAL DETAILS Nonpolarizable AMBER Force Field of Ferric Ion and Porphine. The Amber force field parameters of porphine were generated by Antechamber47,48 included in AmberTools,49 and Acpype50 was used to convert the topology file to Tinker format. Water was considered as solvent to optimize the parameters. The parameters of ferric ion were optimized by fitting the [Fe(H2O)5]3+−H2O dissociation energy profile to the energies obtained from ab initio (MP2) calculations. Polarizable AMOEBA Force Field of Ferric Ion and Porphine. Reference data of porphine were obtained on the basis of QM calculations employing the well-tested B3LYP functional,51−53 and Poltype16 was used to generate the AMOEBA force field parameters of porphine. In the QM calculations, the basis sets used for structure optimization, distributed multipole analysis (DMA),54 and electrostatic potential calculations are 6-31G(d), 6-311G(d,p), and 6-311+ +G(d,p),55−61 respectively. The parameters of Fe3+ were fitted using an in-house program to the BSSE-corrected QM reference energies which will be detailed below. The philosophy of the AMOBEA force field has been described previously by Ren and Ponder.13,16,62 Here we only give a brief description relevant to this work. Interactions between two particles are represented by bonded and nonbonded terms in the AMOEBA potential. The former contains bond stretching, angle bending, out-of-plane bending, torsion, and cross-terms between particles, which reads

4εiiεjj ( εii +

εjj )2

R ij0

=

(R ii0)3 + (R jj0)3 (R ii0)2 + (R jj0)2

The electrostatics contain the contributions of both permanent and induced components. The former is represented by multipoles up to the quadrupole moment, which are obtained from distributed multipole analysis according to Stone64 and fitted to the MP2-derived potential. The orientations of multipoles are defined by local frames through neighboring atoms. The latter, i.e., the polarization, is calculated via a self-consistent procedure of induced dipole. In the calculation of polarization energy, the damping factor of the Tholé model65 is introduced to avoid the divergence at short range, the so-called polarization catastrophe. Implementation of the Charge Transfer Term in the AMOEBA Force Field. In the original formula of the AMOEBA force field, atomic charge is a fixed quantity; thus, its fluctuation due to the charge transfer (CT) is not explicitly handled. In transition metal coordination chemistry, the charge transfer becomes important and may contribute to the interaction between the center metal atom and the ligand significantly, especially when the metal atom is highly positively charged. In this work, to take into account the charge transfer explicitly, we have also implemented in Tinker the discrete charge-transfer (DCT) model according to Lee and Rick,46 and the CT potential is defined as below N+M−1

UCT =

⎛

∑ ∑ ⎝−μijCT |qijCT| + ⎜

i=1

Ubonded = Ubond + Uangle + Uoop + Utorsion + Ucross ‐ terms

j R 2CT ⎩

where QCT is the maximal amount of charge transferred ij CT between each pair. The distances RCT 1 and R2 define the range of the switch function. As the charge transfer strongly depends on the two atoms that are bound to each other, the

Uoop = K χ χ 2 Utorsion =

⎞ ⎛ ⎞7 ⎛ ⎟ ⎜ ⎟⎜ 1.07 1.12 ⎟ ⎟⎜ 2 = εij⎜ − ⎟ ⎜ ⎛⎜ rij ⎞⎟ ⎟ ⎜ ⎛ rij ⎞7 ⎜ 0 + 0.07 ⎟ ⎜ ⎜ 0 ⎟ + 0.12 ⎟ ⎝ ⎝ R ij ⎠ ⎠ ⎝ ⎝ R ij ⎠ ⎠

∑ K nϕ[1 + cos(nϕ ± δ)] n

B


Article


to explore the nature of the interaction between [Fe(H2O)5]3+ and H2O, and to evaluate the quality of the optimized force field parameters. The calculations were done at the MP2/ccpVTZ level, and all energies were corrected for BSSE. In LMOEDA, the total interaction energy Eint consists of five terms:

parametrization was done on the basis of atom pairs. In order to be compatible with the original AMOEBA water parameters, we only considered energy with force for the CT term and do not change the charge on an atom. Atomic Polarizabilities and Reference Energies from QM Calculations. The polarizability of Fe3+ used in the AMOEBA force field was determined using the second-order Møller−Plesset perturbation (MP2)66−68 method with all atoms treated by the aug-cc-pVQZ basis set,69 as implemented in the Gaussian 09 software package.70 All electrons were taken into account in the dynamic electronic correlation energy calculation to get accurate polarizabilities.71 According to our calculations, the polarizabilities of Fe3+ for doublet, quartet, and sextet states are 0.297, 0.273, and 0.258 Å3, respectively, and a dependence on spin states is observed. These values are smaller than those of Fe2+, which were reported to be 0.547, 0.549, and 0.550 Å3 in its singlet, triplet, and quintet states, respectively.43 This is consistent with the feature of Fe3+ being a “harder” Lewis acid than Fe2+ cation. The interaction energy profiles between the [Fe(H2O)5]3+ cluster and one water molecule in each spin state were calculated using the Gaussian 09 program. The structures of [Fe(H2O)6]3+ clusters were first optimized in the Th symmetry using the M06L method72 with Fe treated by the Stuttgart relativistic pseudopotential ECP10MDF73 and the corresponding valence basis set contracted as [8s7p6d1f|6s5p3d1f], and 631G(d)55−58,74 for H and O. The targeting H2O was then moved along the Fe−O vector, as shown in Figure 1 with other

E int = Eele + Epol + Eex + Erep + Edisp

where Eele is the electrostatic energy contributed by the permanent component, Epol the polarization, Eex the exchange, Erep the repulsion, and Edisp the dispersion term. The vdW energy is the sum of exchange, repulsion, and dispersion terms, i.e., EvdW = Eex + Erep + Edisp

Note that these energy components are not observables, and the decomposition of the total interaction energy into the individual components is to understand the nature of the interactions between two fragments, and to evaluate how much the optimized force field parameters can reproduce the data obtained at the QM level. Simulation Methodologies. All MD simulations were performed by using the Tinker 7.1 package.78 For the DCT enhanced AMOEBA force field, a modified version was used. The AMOEBA03 water model79 was used in parametrization and simulations. In the simulation of hydrated ferric cation, the model system contains one ferric cation and 500 water molecules, with a box size of 24 Å. The Beeman algorithm implemented in Tinker was used for integration with a time step of 1 fs. The simulations were conducted with an isothermal isobaric ensemble (NPT), where the temperature was kept at 298 K by Berendsen thermostat with a coupling time (τT) of 0.1 ps, and the pressure was maintained at 1 bar using a Berendsen barostat with a coupling time (τP) of 2 ps. Long-range electrostatic interactions were calculated using particle-mesh Ewald summation of atomic multipoles with a cutoff of 7 Å in real space. Periodic boundary condition (PBC) was applied in all simulations. The van der Waals (vdW) interactions were calculated with a cutoff of 10 Å. The samplings were run for 2 ns after being equilibrated for 50 ps, and the trajectories were saved every 0.1 ps. Radial distribution functions (RDFs) were computed from 2 ns production run trajectories. For the Fe(III)−porphine model system, plain simulation was done to obtain the solvation structure of the Fe(III)−porphine complex in the aqueous phase. It contains one chloride, one dianionic porphine, one ferric cation, and 1500 water molecules. The model system was first equilibrated for 200 ps followed by a 2 ns production run under the conditions similar to the simulations of ferric cation for analysis. The free energy perturbation method (FEP)80 was used to evaluate the hydration free energies of Fe3+ cation in its free state and in its bound state with porphine. All of its three spin states, i.e., doublet, quartet, and sextet states, were considered. The simulations were done with the NVT ensemble at both nonpolarizable and polarizable force field levels, with a cubic box size of ca. 24.6 Å, containing one Fe3+ cation and 500 water molecules. Soft-core approximation was used during the simulation to avoid the presence of singularity.81−83 The Bennett acceptance ratio84 was used to calculate the Helmholtz free energy between two perturbation states λi and λi+1

Figure 1. Model used to derive the Fe−OH2 dissociation energy profile.

geometric parameters of the fragments fixed, to generate a series of structures along the Fe−OH2 dissociation path for all spin multiplicities. For each structure, single point energies was calculated at the MP2 level of theory66−68 with core electrons frozen, and Fe was treated by the Stuttgart ECP10MDF pseudopotential basis set and other elements by def2-TZVP.75 For each single point energy calculation, basis set superposition error (BSSE) was obtained by using the counterpoise method of Boys and Bernardi.76 The corrected QM energy profiles were used to optimize the AMOEBA (R0, ε0, and damping factor) and AMBER force field parameters of ferric cation. Energy Decomposition Analysis. The localized molecular orbital energy decomposition analysis (LMOEDA) method implemented in GAMESS_US (2014-R1 version)77 was used C


Article

The Journal of Physical Chemistry B Table 1. Optimal AMOEBA and AMBER Force Field Parameters of Fe3+ AMOEBA

AMBER

spin state

R0Fe3+ (Å)

ε0Fe3+ (kcal/mol)

α (Å3)

damping factor

RFe3+ (Å)

εFe3+ (kcal/mol)

doublet quartet sextet

2.050 3.354 2.066

0.271 0.354 0.833

0.297 0.273 0.258

0.109 0.390 0.052

2.787 3.182 3.329

0.00126 0.00076 0.00122

λ

ΔA λi , λi+1 = RT ln

∑j i λ

∑ j ′i+1

1 f (U (j , λi + 1) − U (j , λi) − C) 1 f (U (j ′ , λi) − U (j ′ , λi + 1) + C)

+C

where the function f has an exponential form of gas constant R and temperature T: f (x ) =

1 1 + e x / RT

The simulations first considered the perturbation of vdW interaction, and the whole procedure was divided into 14 windows with λ values varied from 0 to 1, i.e., 0.0, 0.0001, 0.001, 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, and 1.0. Then, 41 additional FEP steps were taken for electrostatic (both permanent and polarization terms) with the λ′ values changed from 0 to 1 with a spacing of 0.025 by linearly scaling the atomic charge qi. For each perturbed state, the simulation ran for 500 ps with the first 50 ps for equilibration and the remaining 450 ps for sampling used in free energy calculations. In the calculation of the hydration free energy of Fe3+ with CT explicitly treated, the DCT term was scaled by the λi values as electrostatics. We note that the finite size effect may bring a significant influence on the hydration free energy of cations, and potential correction schemes have been proposed by several groups,85−88 which should be taken into account when making a direct comparison between the calculated values and the experimental ones. In the present work, the simulations were done with a canonical ensemble; thus, the calculated values are for the end states that are not fully equivalent with the ones measured in experiments under constant pressure. Therefore, in the present work, the raw data of the hydration free energy of ferric cation was reported for a qualitative comparison with experimental results and evaluation of force field performance.

■

RESULTS AND DISCUSSION A. Parameterization of AMOEBA and Nonpolarizable Force Fields. The optimized parameters derived from a nonlinear least-squares fitting procedure of the intercluster interaction energies at the force field level fitted to the MP2 reference energies are shown in Table 1, and the BSSEcorrected [Fe(H2O)5]3+−H2O dissociation energy profiles from calculations at MP2, AMOEBA, and AMBER force field levels are plotted for doublet, quartet, and sextet spin states in Figure 2. As seen in Figure 2, the AMOEBA force field reproduces the MP2 energy profiles well, with optimal Fe−OH2 distances of 1.89, 1.90, and 2.05 Å in the three spin states of doublet, quartet, and sextet, respectively, corresponding to well depths of −71.88, −57.05, and −63.72 kcal/mol. In contrast, the nonpolarizable force field overestimates the Fe−OH2 interaction strength in the medium- and long-range regions. B. Nature of Interaction. In order to understand the nature of the interaction between the [Fe(H2O)5]3+ fragment

Figure 2. [Fe(H2O)5]3+−H2O dissociation profiles: (a) doublet; (b) quartet; (c) sextet. The AMOEBA03 and TIP3P water models were used in the calculations of energy profiles at AMOEBA and AMBER force field levels, respectively. In the nonpolarizable force field calculations, the SPC/E model was also used and similar results were obtained as those using the TIP3P model.

and H2O and evaluate the accuracy of the optimized AMOEBA and nonpolarizable parameters, LMOEDA was performed for D


Article


Figure 3. Energy decomposition for the [Fe(H2O)5]3+−H2O dissociation profile for sextet spin states from ab initio based LMOEDA (black squares), AMBER force field (red triangles), AMOEBA force field (blue diamonds), and AMOEBA with CT term included (green stars): (a) van der Waals energy; (b) polarization energy; (c) electrostatic energy (permanent component); (d) total binding energy.

interactions. In order to find out the origin of this deviation, we implemented the charge transfer term (CT) in the AMOEBA Hamiltonian and reoptimized the force field parameters with CT taken into account for the sextet state, which are tabulated in Table 2. With the explicit treatment of the CT term, we adopted the standard polarization damping factor in AMOEBA, a = 0.39. The decomposed energy components are plotted in Figure 3. As seen in Figure 3, for the [Fe(H2O)5]3+−H2O intercluster dissociation energy profile, taking into account CT improves the agreement of the vdW and polarization energies with the LMOEDA values, suggesting the importance of an appropriate treatment of charge transfer in order to obtain a correct

each spin state. The data for the sextet state are shown in Figure 3, and those of the other states are provided in the Supporting Information. Figure 3a shows van der Waals (vdW) energy along the dissociation pathway. It is clear that the vdW energy quickly decays in the short range and fades beyond 3.0 Å. With respect to the LMOEDA data, the AMOEBA force field was systematically overestimated for the vdW energy especially in short distance. In the AMOEBA force field, the overestimation was corrected by the inclusion of the polarization energy, and a good agreement with the LMOEDA result in total energy was achieved. For comparison, the nonpolarizable force field (ca. 80 kcal/mol) overestimated the interaction between Fe3+ and water by more than 10 kcal/mol. We note that there is large deviation between the polarization energy values calculated at the AMOEBA force field level and the ab initio reference data. The LMOEDA values decrease monotonously as the interatomic distance increases, suggesting a diminishing polarization effect when two atoms moved away from each other. At the AMOEBA force field level, a shallow well of −4.98 kcal/mol at an Fe−Ow distance of 2.25 Å is observed on the energy profile. C. Charge Transfer. As seen in Figure 3, the polarization energy calculated at the AMOEBA force field level is underestimated with respect to the ab initio values in the range of short Fe−Ow distance. At the AMOEBA force field level, this deviation is compensated by the overestimated vdW

Table 2. Optimal AMOEBA Parameters with the CT Term Included for the Fe3+−Ow Pair in Its Sextet State van der Waals polarization term charge transfer term

E

R (Å) ε (kcal/mol) α (Å3) a μCT (kcal/mol/e) ηCT (kcal/mol/e2) QCT (e) R1 (Å) R2 (Å)

4.983 0.021 0.258 0.390 230.0 100.0 0.100 1.300 4.000


Article


AMOEBA with CT term included, and nonpolarizable force field (AMBER) were calculated and plotted in Figure 4. When going from doublet to sextet states, the Fe−Ow distance increases for the first coordination shell, while it does not change much for the second shell. In the sextet state, the first maximum of RDF is located at a Fe−Ow distance of about 2.0 Å, which agrees well with experimental results ranging from 1.99 to 2.10 Å.90,91,93 The integration of the radial distribution functions gives a coordinate number of the first solvation shell of 6.0 for all the three spin states, which agrees with experimental results.89 These results show that all force field treatment of Fe3+ gives reasonable solvation structures of ferric ion in the aqueous phase. The diffusion coefficients of Fe3+(aq) were calculated with the optimized AMOEBA force field parameters in its three spin states, which are 0.72 × 10−10, 2.81 × 10−10, and 6.16 × 10−10 m2/s (Table 4). The calculated D value in the sextet state shows a good agreement with the experimentally measured data98 of 6.04 × 10−10 m2/s. The smaller D values in the lower spin states could be due to the closer Fe−Ow distance, leading to a more structured solvent environment in the vicinity. D.3. Free Energy of Hydration. In dilute solution, ferric cation exists in its hydrated form. By using the free energy perturbation method, the hydration free energy of ferric ion was calculated at both the AMOEBA and nonpolarizable force field levels in the three spin states and compared to experimentally determined data. In its sextet state, the value is −916.0 kcal/mol using the AMOEBA force field and −873.7 kcal/mol with the nonpolarizable force field. These values are less exothermic than the experimental data of −1036.499 or −1019.8 kcal/mol100 with more underestimation of the hydration of ferric ion at the nonpolarizable FF level than at the PFF level (11.6%). Taking into account the CT explicitly improved the calculated hydration energy to −940.4 kcal/mol. This suggests that charge transfer contributes to the Fe−water interaction to form a stable aqua complex of ferric ion when solvated in aqueous solution, as in the case of Fe2+(aq).101 In order to evaluate the relative thermodynamic stabilities of the three spin states, we calculated their relative free energy of hydration following the thermodynamic cycle in Scheme 1,

physical picture for the interaction in a transition metal−ligand coordination system. D. Validation and Application of Force Field Parameters. D.1. Structures. To evaluate the quality of the force field parameters obtained above, we first used them to optimize the geometry of [Fe(H2O)6]3+, and compared with the data from M06L calculations. The hydrated structure of ferric cation was optimized to Th symmetry at the M06L/[Fe: ECP10MDF; C,H:6-31G(d)] level with Fe3+ in its doublet, quartet, and sextet spin states, and the Fe−Ow distances in the first coordination shell are given in Table 3. Table 3. Fe−Ow Distances (Å) in the Optimized Structures of [Fe(H2O)6]3+ Using QM(M06-L) and Force Field (AMOEBA and AMBER) QM (M06-L) AMOEBA AMBER AMOEBA with CT

doublet

quartet

sextet

1.94 1.92 1.86

2.00 1.95 1.91

2.05 2.06 1.96 2.05

The [Fe(H2O)6]3+ complex in the higher spin states has longer Fe−Ow distances than in the low spin state. The calculations using both AMOEBA and nonpolarizable force fields reproduced the tendency. The optimized Fe−Ow distance in the sextet spin state has the best agreement between M06-L and AMOEBA values, and that in the quartet state has the largest deviation of 0.05 Å from the MP2 value of 2.00 Å. The deviation in the doublet spin state cluster is 0.02 Å. When using the nonpolarizable force field, the Fe−Ow distance is predicted to always be shorter than the corresponding QM and AMOEBA values in all three spin states, suggesting a stronger ferric−water interaction. These results show that the AMOEBA force field gives a more precise description for the interaction of ferric cation and water than the nonpolarizable force field with respect to QM reference. D.2. Hydration of Ferric Cation (Fe3+). To understand the hydrated structure of ferric ion, the radial distribution function (RDF) and its integral of the hydrated ferric cation from simulations with the model system described by AMOEBA, and

Figure 4. Radial distribution function (solid line) and its integral (dashed line) of water molecules around Fe3+ in the sextet state during simulations using (a) AMOEBA, (b) AMOEBA with CT included, and (c) AMBER force fields. F


Article


Table 4. Distances between Fe3+ and the Ow Atoms of Water in the First (dFe−Ow(1), in Å) and Second (dFe−Ow(2), in Å) Coordination Shells and the Corresponding Coordination Numbers (C.N.(1) and C.N.(2)) and the Diffusion Coefficient (D, in 10−10 m2/s) of Fe3+ Cation in Water AMOEBA dFe−Ow(1) C.N.(1) dFe−Ow(2) C.N.(2) D a

AMBER a

doublet

quartet

sextet

1.86 6.0 3.99 12.3 0.72

1.90 6.0 4.03 11.0 2.81

1.99 (2.01) 6.0 (6.0) 4.16 (4.25) 12.2 (12.2) 6.16 (5.2)

doublet

quartet

sextet

Exp.

1.84 6.0 4.04 13.1 7.05

1.89 6.0 3.96 11.4 3.77

1.94 6.0 4.02 11.3 21.05

1.99−2.1089−96 6.092−96 4.1596 12.096 6.0497

Data in parentheses were obtained by using the AMOEBA force field with CT included.

Scheme 1. Thermodynamic Cycle to Calculate the Solvation Free Energy Difference of Ferric Ion in Its Three Spin States

D.4. Fe(III)−Porphine Complex. The iron−porphine complex, which represents the core unit of heme, is vital to bioorganisms by functioning as the active center in monooxygenase and peroxooxygenase proteins.102 In these proteins, the distal side of heme is in general bound to cysteine or histidine residue, and its proximal side offers a catalytic site for reactions to happen. Previous studies on heme, if carried out at force field level, always handled heme with Fe and porphine (Por) moiety as an integrated residue. Such a treatment is to circumvent the difficulty in the description of the interaction between the highly charged Fe cation and the macrocyclic ligand at the cost of losing the message on the Fe−N coordination interaction. In order to delineate the nature of the Fe−N dative interaction and validate the use of the ferric force field parameters in organometallic and biomolecular systems, here we chose the Fe(III)−Por complex as the test model, which is shown in Scheme 2.

where the energy differences in the gas phase between the sextet state and the doublet (ΔEsextet→doublet = 174.3 kcal/mol) or the quartet (ΔEsextet→quartet = 107.9 kcal/mol) state were determined at the MP2(full)/aug-cc-pVQZ level, and the results are collected in Table 5. Table 5. Binding Free Energies (ΔF, kcal/mol) of Ferric Ion with Water and Porphrin in Aqueous Solution Fe3+ (aq) ΔFdoublet ΔFquartet ΔFsextet ΔFsextet→doublet ΔFsextet→quartet

FePor+ (aq)

AMOEBA

AMBER

AMOEBA

AMBER

−1007.8 −975.6 −916.0 82.5 48.3

−935.1 −907.4 −873.7 112.9 74.2

−1018.6 −1014.1 −874.5 30.0 −31.4

−807.3 −796.9 −781.4 148.4 92.4

As seen in Table 5, in aqueous solution, the hydrated ferric ion in its sextet state is more stable than that in its doublet and quartet states with a relative free energy difference of 82.5 and 48.3 kcal/mol, respectively, at the AMOEBA level. Calculations at the AMBER FF level gave similar predictions but with a larger free energy difference. These results qualitatively agree with the recent multiconfiguration electronic structure study101 that identified an 6Ag ground state for the [Fe(H2O)6]3+ cluster with Th symmetry, and are consistent with the coordination chemistry of ferric ion that high spin state is favored when coordinated to weak field ligands, here water, according to coordination field theory. We have also tested hydrolyzed ferric ionic species, [Fe(H2O)5(OH)]2+ and [Fe(H2O)4(OH)2]+, in their doublet and sextet states. At the nonpolarizable force field level, the free energy differences (ΔFsextet→doublet) are 103.7 and 80.5 kcal/mol, respectively, and at the AMOEBA force field level, the values were calculated to be 52.7 and 66.3 kcal/mol, respectively. This shows that the simulations at the force field levels predict that the hydrated ferric hydroxide species are more stable in its sextet state than in its doublet state.

Scheme 2. Fe(III)−Porphine Complex

In this complex, in its equatorial plane, Fe(III) is coordinated to four N atoms from four pyrrolic rings linked to each other via a methylene group. When solvated in aqueous solution, its axial coordination sites were accessible for water molecules. On the basis of MD simulations, the binding of ferric ion with Por G


Article


Figure 5. Radial distribution function of water molecules around Fe in the Fe(III)−Por complex in its quartet state and its integral calculated from the geometric trajectories obtained at the AMOEBA and AMBER FF levels.

Table 6. Fe−Ow Distances (Å) in the Optimized Structures of [Fe(Por)]+ Using QM(M06-L) and Force Field (AMOEBA and AMBER) doublet AMOEBA Non-Pol QM

quartet

sextet

Fe−N

Fe−Ow

Fe−N

Fe−Ow

Fe−N

Fe−Ow

2.04 ± 0.05 2.18 ± 0.05 2.01

1.84 1.80 1.94

2.01 ± 0.09 2.19 ± 0.04 2.00

1.93 1.86 2.00

2.23 ± 0.07 2.20 ± 0.04 2.05

2.00 1.93 2.05

Scheme 3. Thermodynamic Cycle to Evaluate the Relative Binding Strength of Fe3+ with Por in Its Three Spin States

which is about 0.20 Å more than the latter ones, suggesting Fe3+ ion is more weakly bound to Por in its sextet state. We note that there are five water molecules in the second coordination shell of FeIII in the FeIII−Por complex in aqueous solution with the AMBER force field while four with the AMOEBA force field. This shows that the explicit treatment of polarization favors a less compact solvation structure for the Fe−Por complex. This coordination complex is stable, and during the 2.0 ns of sampling, no exchange was observed between the bound water molecules and the bulk ones in the doublet and quartet states. The higher stability of the complex in the doublet and quartet states than in the sextet state is also reflected by the observation that, during the simulation, ferric ion remained in the center of porphine with the Por ligand in its planar conformation in the former case, while, in the sextet state, the ferric ion migrates away the ring center and favors staying above the Por ligand. This is consistent with the fact that in its sextet state ferric ion favors being maximally hydrated. As discussed above, in the aqueous phase, the ferric ion prefers to exist in its sextet state. When coordinating with Por, the preference may change, since the N atom of pyrrol ring of Por is a “softer” electron donor than a water molecule, and thus

was analyzed in view of both geometric and energetic properties, and compared to available data from QM studies. The analysis of the geometric trajectory shows that, during the simulation, Fe3+ remained bound at the center position of the porphrin ligand, with its proximal and distal sites occupied by two water molecules. The radial distribution functions of water in the surrounding region of Fe3+ in its three spin states were calculated; the one in the quartet state is plotted in Figure 5, and the averaged Fe−N and Fe−Ow distances are collected in Table 6. The calculations give Fe−Ow distances of 1.84 (AMOEBA) or 1.80 (AMBER) for the first shell in the doublet state. These values are smaller than those in the [Fe(H2O)6]3+ complex. Going to higher spin states leads to longer Fe−Ow distances, 1.93 Å in the quartet and 2.00 Å in the sextet states. This trend is qualitatively consistent with the QM results. The longer Fe−Ow distance in the sextet state suggests weaker interaction. This agrees with the analysis of two sequentially collected 2 ns long trajectories at the AMOEBA level that the exchange between the bound water and the bulk water happened only for the sextet states (residence time of 14−20 ps) while not for the lower spin states. For the Fe−N distances, the AMOEBA force filed distinguishes the sextet state from the doublet and quartet states by a larger value for the former H



■

may exert a relatively stronger field on Fe3+ than water and modulate its electronic configuration. To evaluate the relative binding strength of Fe3+ with Por in its three spin states, the relative free energies were calculated following the thermodynamic cycle in Scheme 3 and the data were tabulated in Table 5. As shown in Table 5, with the AMOEBA force field, the binding of Fe3+ with Por ligand in water is exothermic by 1018.6, 1014.1, and 874.5 kcal/mol in doublet, quartet, and sextet states. Compared with the hydration free energies of ferric ion in water, higher exothermicity is observed for the doublet (−10.8 kcal/mol) and quartet (−38.5 kcal/mol) states while lower for the sextet state (41.5 kcal/mol) when binding with Por ligand, suggesting that Fe3+ favors coordinating to Por ligand in its low spin state while it prefers staying in the aqueous phase in its high spin state. With a nonpolarizable force field, the ligand change process of water ligands by Por is estimated to be an endothermal process in all of the three spin states. This shows that a nonpolarizable force field tends to predict stronger interaction between ferric ion with water than with Por ligand. The calculations also show that, at the AMOEBA force field level, the FeIII−Por complex was predicted to exist in its quartet state, with a free energy preference of 30.0 kcal/mol over doublet and 31.4 kcal/mol over sextet states. This is qualitatively consistent with the recently reported density functional theory study.103 The nonpolarizable force field fails to capture this spin inversion process of ferric ion upon the substitution of Por for the water ligands, and favors a sextet state always over doublet and quartet states. Replacing a CH3S− ligand at the distal position of Fe(III) to mimic the situation in P450 predicts stable complexes in the doublet and quartet states at the AMOEBA level, while, in sextet, the CH3S− cannot build a strong enough interaction to resist the competition of water.

Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.7b02010. ELF basins for the sextet [Fe(H2O)6]3+ complex (Th point group); energy components of the [Fe(H2O)5]3+− H2O interfragment energy at the levels of AMOEBA, AMBER, AMOEBA with CT, and LMOEDA (MP2); energy components of the Fe3+−(H2O)6 interfragment energy at the levels of AMOEBA, AMBER, AMOEBA with CT, and LMOEBA (MP2); RDF of water around Fe3+ in water in the doublet and quartet spin states at AMOEBA and AMBER force field levels; RDF of water oxygen around Fe3+ of FeIII−Por at the AMOEBA and AMBER force field levels in the doublet and sextet spin states; the distribution of the angle between the water dipole moment and the centroid of the pyrrol ring nearby; AMOEBA and AMBER force field parameters of porphine in Tinker format; coordinates of optimized [Fe(H2O)6]3+ by Gaussian (M06L, Th) in the doublet, quartet, and sextet states; and Fe3+−OW dissociation profile from MP2 calculations (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Miaoren Xia: 0000-0002-9970-0802 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was financially supported by the National Natural Science Foundation of China to Z.C. (91026000) and to D.W. (21473206, 91226105) and by the CAS Hundred Talents Program (Y2291810S3, Y1515540U1), which are gratefully acknowledged. The Special Program for Applied Research on Super Computation of the NSFC-Guangdong Joint Fund (the second phase) and the National Supercomputing Center in Tianjin (NSCC-TJ) are acknowledged for providing computational grids.

■

CONCLUSIONS We have optimized the AMOEBA force fields for iron(III) cation in its three spin states based on ab initio calculations. The localized molecular orbital energy decomposition analysis (LMOEDA) approach implemented in GAMESS-US was used to understand the nature of intercluster interaction of [Fe(H 2O) 5]3+−H 2O and calibrate the accuracy of the optimized AMOEBA and AMBER force field parameters. The discrete charge-transfer (DCT) term proposed by W. Rick was implemented in Tinker to evaluate the contribution of charge transfer between the ferric ion and the ligands, which showed that CT was important for a correct picture of the van der Waals and polarization energies of the model system. The optimized AMOEBA and AMBER force fields of ferric ion were validated and applied to study the aquo structure of Fe3+ and its complex with porphine. Radial distribution functions (RDFs) show good agreement with experimental data for coordinate structure in the bulk phase. The hydration free energy of ferric ion and the binding free energy of Fe(III) with porphine ligand were calculated by using the free energy perturbation (FEP) method, and the most stable spin states of [Fe(H2O)6]3+ and [Fe(Por)(H2O)2]+, i.e., the sextet and quartet states, respectively, were identified correctly at the AMOEBA level, which is consistent with the available QM results. This suggests that the optimized AMOEBA parameters of ferric ion may be used in the study of its dynamics in the condensed phase and in biomolecular systems.

■

REFERENCES

(1) Corongiu, G. Molecular Dynamics Simulation for Liquid Water Using a Polarizable and Flexible Potential. Int. J. Quantum Chem. 1992, 42, 1209−1235. (2) Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T.; et al. A Point-Charge Force Field for Molecular Mechanics Simulations of Proteins Based on Condensed-Phase Quantum Mechanical Calculations. J. Comput. Chem. 2003, 24, 1999−2012. (3) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Comparison of Multiple Amber Force Fields and Development of Improved Protein Backbone Parameters. Proteins: Struct., Funct., Genet. 2006, 65, 712−725. (4) MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586−3616. (5) Oostenbrink, C.; Villa, A.; Mark, A. E.; Gunsteren, W. F. V. A Biomolecular Force Field Based on the Free Enthalpy of Hydration and Solvation: The GROMOS Force-Field Parameter Sets 53A5 and 53A6. J. Comput. Chem. 2004, 25, 1656−1676.

I


Article

The Journal of Physical Chemistry B (6) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (7) Li, P.; Merz, K. M., Jr. Metal Ion Modeling Using Classical Mechanics. Chem. Rev. 2017, 117, 1564−1686. (8) Warshel, A.; Kato, M.; Pisliakov, A. V. Polarizable Force Fields: History, Test Cases, and Prospects. J. Chem. Theory Comput. 2007, 3, 2034−2045. (9) Halgren, T. A.; Damm, W. Polarizable Force Fields. Curr. Opin. Struct. Biol. 2001, 11, 236−242. (10) Ji, C.; Mei, Y. Some Practical Approaches to Treating Electrostatic Polarization of Proteins. Acc. Chem. Res. 2014, 47, 2795−2803. (11) Baker, C. M. Polarizable Force Fields for Molecular Dynamics Simulations of Biomolecules. WIREs Comput. Mol. Sci. 2015, 5, 241− 254. (12) Xu, P.; Wang, J.; Xu, Y.; Chu, H.; Liu, J.; Zhao, M.; Zhang, D.; Mao, Y.; Li, B.; Ding, Y.; et al. Advancement of Polarizable Force Field and Its Use for Molecular Modeling and Design. Adv. Exp. Med. Biol. 2015, 827, 19−32. (13) Ponder, J. W.; Wu, C. Current Status of the AMOEBA Polarizable Force Field. J. Phys. Chem. B 2010, 114, 2549−2564. (14) Chen, J. Theory and Applications of Fluctuating−Charge Models; Proquest, Umi Dissertation Publishing: Urbana, IL, 2011. (15) Lemkul, J. A.; Huang, J.; Roux, B.; MacKerell, A. D. An Empirical Polarizable Force Field Based on the Classical Drude Oscillator Model: Development History and Recent Applications. Chem. Rev. 2016, 116, 4983−5013. (16) Wu, J. C.; Chattree, G.; Ren, P. Automation of AMOEBA Polarizable Force Field Parameterization for Small Molecules. Theor. Chem. Acc. 2012, 131, 1138. (17) Dixon, S. J.; Stockwell, B. R. The Role of Iron and Reactive Oxygen Species in Cell Death. Nat. Chem. Biol. 2014, 10, 9−17. (18) Pain, D.; Dancis, A. Roles of Fe-S Proteins: From Cofactor Synthesis to Iron Homeostasis to Protein Synthesis. Curr. Opin. Genet. Dev. 2016, 38, 45−51. (19) Milto, I. V.; Suhodolo, I. V.; Prokopieva, V. D.; Klimenteva, T. K. Molecular and Cellular Bases of Iron Metabolism in Humans. Biochemistry (Moscow) 2016, 81, 549−564. (20) Stiban, J.; So, M.; Kaguni, L. S. Iron-Sulfur Clusters in Mitochondrial Metabolism: Multifaceted Roles of a Simple Cofactor. Biochemistry (Moscow) 2016, 81, 1066−1080. (21) Pietrangelo, A. Iron and the Liver. Liver Int. 2016, 36 (Suppl. 1), 116−123. (22) Brear, E. M.; Day, D. A.; Smith, P. M. Iron: An Essential Micronutrient for the Legume-Rhizobium Symbiosis. Front. Plant Sci. 2013, 4, 359. (23) Pra, D.; Franke, S. I.; Henriques, J. A.; Fenech, M. Iron and Genome Stability: An Update. Mutat. Res., Fundam. Mol. Mech. Mutagen. 2012, 733, 92−99. (24) Ward, R. J.; Crichton, R. R.; Taylor, D. L.; Della Corte, L.; Srai, S. K.; Dexter, D. T. Iron and the Immune System. J. Neural Transm. 2011, 118, 315−328. (25) Chen, A.; Shang, C.; Shao, J.; Zhang, J.; Huang, H. The Application of Iron-Based Technologies in Uranium Remediation: A Review. Sci. Total Environ. 2017, 575, 1291−1306. (26) Gong, Y.; Tang, J.; Zhao, D. Application of Iron Sulfide Particles for Groundwater and Soil Remediation: A Review. Water Res. 2016, 89, 309−320. (27) Jack, R. S.; Ayoko, G. A.; Adebajo, M. O.; Frost, R. L. A Review of Iron Species for Visible-Light Photocatalytic Water Purification. Environ. Sci. Pollut. Res. 2015, 22, 7439−7449. (28) Feng, W.; Nansheng, D. Photochemistry of Hydrolytic Iron (III) Species and Photoinduced Degradation of Organic Compounds. A Minireview. Chemosphere 2000, 41, 1137−1147. (29) Sun, C. L.; Li, B. J.; Shi, Z. J. Direct C-H Transformation via Iron Catalysis. Chem. Rev. 2011, 111, 1293−1314.

(30) Sarhan, A. A.; Bolm, C. Iron(III) Chloride in Oxidative C-C Coupling Reactions. Chem. Soc. Rev. 2009, 38, 2730−2744. (31) Correa, A.; Garcia Mancheno, O.; Bolm, C. Iron-Catalysed Carbon-Heteroatom and Heteroatom-Heteroatom Bond Forming Processes. Chem. Soc. Rev. 2008, 37, 1108−1117. (32) Cera, G.; Ackermann, L.; Iron-Catalyzed, C.-H. Functionalization Processes. Top. Curr. Chem. 2016, 374, 57. (33) Lv, L.; Li, Z. Fe-Catalyzed Cross-Dehydrogenative Coupling Reactions. Top. Curr. Chem. 2016, 374, 38. (34) Mérel, D. S.; Do, M. L. T.; Gaillard, S.; Dupau, P.; Renaud, J.-L. Iron-Catalyzed Reduction of Carboxylic and Carbonic Acid Derivatives. Coord. Chem. Rev. 2015, 288, 50−68. (35) Jia, F.; Li, Z. Iron-Catalyzed/Mediated Oxidative Transformation of C−H Bonds. Org. Chem. Front. 2014, 1, 194. (36) Yaroshevsky, A. A. Abundances of Chemical Elements in the Earth’s Crust. Geochem. Int. 2006, 44, 48−55. (37) Nastri, F.; Lombardi, A.; D’Andrea, L. D.; Sanseverino, M.; Maglio, O.; Pavone, V. Miniaturized Hemoproteins. Biopolymers 1998, 47, 5−22. (38) Gkouvatsos, K.; Papanikolaou, G.; Pantopoulos, K. Regulation of Iron Transport and the Role of Transferrin. Biochim. Biophys. Acta, Gen. Subj. 2012, 1820, 188−202. (39) Irie, S.; Tavassoli, M. Review: Transferrin-Mediated Cellular Iron Uptake. Am. J. Med. Sci. 1987, 293, 103−111. (40) Theil, E. C. Iron, Ferritin, and Nutrition. Annu. Rev. Nutr. 2004, 24, 327−343. (41) Wang, W.; Knovich, M. A.; Coffman, L. G.; Torti, F. M.; Torti, S. V. Serum Ferritin: Past, Present and Future. Biochim. Biophys. Acta, Gen. Subj. 2010, 1800, 760−769. (42) Abbaspour, N.; Hurrell, R.; Kelishadi, R. Review on Iron and Its Importance for Human Health. J. Res. Med. Sci. 2014, 19, 164−174. (43) Semrouni, D.; Isley, W. C., III; Clavaguera, C.; Dognon, J.-P.; Cramer, C. J.; Gagliardi, L. Ab Initio Extension of the AMOEBA Polarizable Force Field to Fe2+. J. Chem. Theory Comput. 2013, 9, 3062−3071. (44) Soniat, M.; Rick, S. W. The Effects of Charge Transfer on the Aqueous Solvation of Ions. J. Chem. Phys. 2012, 137, 044511. (45) Soniat, M.; Hartman, L.; Rick, S. W. Charge Transfer Models of Zinc and Magnesium in Water. J. Chem. Theory Comput. 2015, 11, 1658−1667. (46) Lee, A. J.; Rick, S. W. The Effects of Charge Transfer on the Properties of Liquid Water. J. Chem. Phys. 2011, 134, 184507. (47) Wang, J.; Wang, W.; Kollman, P. A.; Case, D. A. Automatic Atom Type and Bond Type Perception in Molecular Mechanical Calculations. J. Mol. Graphics Modell. 2006, 25, 247−260. (48) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General Amber Force Field. J. Comput. Chem. 2004, 25, 1157−1174. (49) Case, D. A.; Betz, R. M.; Botello-Smith, W.; Cerutti, D. S.; Cheatham, T. E., III; Darden, T. A.; Duke, R. E.; Giese, T. J.; Gohlke, H.; Goetz, A. W.; et al. Amber, version 2016; University of California: San Francisco, CA, 2016. (50) Sousa da Silva, A. W.; Vranken, W. F. ACPYPE - Antechamber Python Parser Interface. BMC Res. Notes 2012, 5, 367. (51) Parr, R. G.; Yang, W. Density Functional Theory of Atoms and Molecules; Oxford University Press: New York, 1989. (52) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (53) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98, 11623−11627. (54) Stone, A. J. Distributed Multipole Analysis: Stability for Large Basis Sets. J. Chem. Theory Comput. 2005, 1, 1128−1132. (55) Ditchfield, R.; Hehre, W. J.; Pople, J. A. Self-Consistent Molecular-Orbital Methods. IX. An Extended Gaussian-Type Basis for Molecular-Orbital Studies of Organic Molecules. J. Chem. Phys. 1971, 54, 724−728. J


Article


(77) Dykstra, C. E.; Frenking, G.; Kim, K. S.; Scuseria, G. E. Theory and Applications of Computational Chemistry: The First Forty Years; Elsevier: Amsterdam, The Netherlands, 2005. (78) Ponder, J. W. TINKER: Software Tools for Molecular Desgin, version 7.1; Washington University, School of Medicine: Saint Louis, MO, 2015. (79) Ren, P. Y.; Ponder, J. W. Polarizable Atomic Multipole Water Model for Molecular Mechanics Simulation. J. Phys. Chem. B 2003, 107, 5933−5947. (80) Zwanzig, R. W. High-Temperature Equation of State by a Perturbation Method 0.1. Nonpolar Gases. J. Chem. Phys. 1954, 22, 1420−1426. (81) Beutler, T. C.; Mark, A. E.; van Schaik, R. C.; Gerber, P. R.; van Gunsteren, W. F. Avoiding Singularities and Numerical Instabilities in Free Energy Calculations Based on Molecular Simulations. Chem. Phys. Lett. 1994, 222, 529−539. (82) Jiao, D.; Golubkov, P. A.; Darden, T. A.; Ren, P. Calculation of Protein−Ligand Binding Free Energy by Using a Polarizable Potential. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 6290−6295. (83) Schnieders, M. J.; Baltrusaitis, J.; Shi, Y.; Chattree, G.; Zheng, L.; Yang, W.; Ren, P. The Structure, Thermodynamics, and Solubility of Organic Crystals from Simulation with a Polarizable Force Field. J. Chem. Theory Comput. 2012, 8, 1721−1736. (84) Bennett, C. H. Efficient Estimation of Free-Energy Differences from Monte-Carlo Data. J. Comput. Phys. 1976, 22, 245−268. (85) Hummer, G.; Pratt, L. R.; Garcia, A. E. Free Energy of Ionic Hydration. J. Phys. Chem. 1996, 100, 1206−1215. (86) Figueirido, F.; DelBueno, G. S.; Levy, R. M. On Finite-Size Corrections to the Free Energy of Ionic Hydration. J. Phys. Chem. B 1997, 101, 5622−5623. (87) Reif, M. M.; Huenenbergera, P. H. Computation of Methodology-Independent Single-Ion Solvation Properties from Molecular Simulations. III. Correction Terms for the Solvation Free Energies, Enthalpies, Entropies, Heat Capacities, Volumes, Compressibilities, and Expansivities of Solvated Ions. J. Chem. Phys. 2011, 134, 144103. (88) Rocklin, G. J.; Mobley, D. L.; Dill, K. A.; Huenenberger, P. H. Calculating the Binding Free Energies of Charged Species Based on Explicit-Solvent Simulations Employing Lattice-Sum Methods: An Accurate Correction Scheme for Electrostatic Finite-Size Effects. J. Chem. Phys. 2013, 139, 184103. (89) Sham, T. K.; Hastings, J. B.; Perlman, M. L. Structure and Dynamic Behavior of Transition-Metal Ions in Aqueous Solution: An EXAFS Study of Electron-Exchange Reactions. J. Am. Chem. Soc. 1980, 102, 5904−5906. (90) Asakura, K.; Nomura, M.; Kuroda, H. Fe K-Edge XANES and EXAFS of the X-Ray Absorption Spectra of FeCl3 Aqueous Solutions A Structural Study of the Solute, Iron(III) Chloro Complexes. Bull. Chem. Soc. Jpn. 1985, 58, 1543−1550. (91) Magini, M.; Caminiti, R. Structure of Highly Concentrated Iron(III) Salt-Solutions. J. Inorg. Nucl. Chem. 1977, 39, 91−94. (92) Apted, M. J.; Waychunas, G. A.; Brown, G. E. Structure and Specification of Iron Complexes in Aqueous Solutions Determined by X-Ray Absorption Spectroscopy. Geochim. Cosmochim. Acta 1985, 49, 2081−2089. (93) Herdman, G. J.; Neilson, G. W. Ferrous Fe(II) Hydration in a 1 Molal Heavy Water Solution of Iron Chloride. J. Phys.: Condens. Matter 1992, 4, 649−653. (94) Ohtaki, H.; Radnai, T. Structure and Dynamics of Hydrated Ions. Chem. Rev. 1993, 93, 1157−1204. (95) Herdman, G. J.; Neilson, G. W. Ferric Ion (Fe(III)) Coordination in Concentrated Aqueous Electrolyte Solutions. J. Phys.: Condens. Matter 1992, 4, 627−638. (96) Caminiti, R.; Magini, M. An X-Ray Diffraction Study on the First and the Second Hydration Shell of the Fe(III) Ion in Nitrate Solutions. Chem. Phys. Lett. 1979, 61, 40−44. (97) Buffle, J.; Zhang, Z.; Startchev, K. Metal Flux and Dynamic Speciation at (Bio)interfaces. Part I: Critical Evaluation and Compilation of Physicochemical Parameters for Complexes with

(56) Hehre, W. J.; Ditchfield, R.; Pople, J. A. Self-Consistent Molecular Orbital Methods. XII. Further Extensions of Gaussian-Type Basis Sets for Use in Molecular Orbital Studies of Organic Molecules. J. Chem. Phys. 1972, 56, 2257−2261. (57) Hariharan, P. C.; Pople, J. A. Accuracy of AHn Equilibrium Geometries by Single Determinant Molecular Orbital Theory. Mol. Phys. 1974, 27, 209−214. (58) Hariharan, P. C.; Pople, J. A. The Influence of Polarization Functions on Molecular Orbital Hydrogenation Energies. Theor. Chem. Acc. 1973, 28, 213−222. (59) Spitznagel, G. W.; Clark, T.; Schleyer, P. V.; Hehre, W. J. Efficient Diffuse Function-Augmented Basis-Sets for Anion Calculation 0.4. An Evaluation of the Performance of Diffuse FunctionAugmented Basis-Sets for 2nd-Row Elements, Na-Cl. J. Comput. Chem. 1987, 8, 1109−1116. (60) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. V. Efficient Diffuse Function-Augmented Basis Sets for Anion Calculations. III. The 3-21+G Basis Set for First-Row Elements, Li-F. J. Comput. Chem. 1983, 4, 294−301. (61) Frisch, M. J.; Pople, J. A.; Binkley, J. S. Self-Consistent Molecular Orbital Methods 25. Supplementary Functions for Gaussian Basis Sets. J. Chem. Phys. 1984, 80, 3265−3269. (62) Zhou, R. Molecular Modeling at the Atomic Scale: Methods and Applications in Quantitative Biology; CRC Press: New York, 2015. (63) Halgren, T. A. Representation of van der Waals (vdW) Interactions in Molecular Mechanics Force Fields: Potential Form, Combination Rules, and vdW Parameters. J. Am. Chem. Soc. 1992, 114, 7827−7843. (64) Stone, A. The Theory of Intermolecular Forces; Oxford University Press: Oxford, U.K., 1997. (65) Thole, B. T. Molecular Polarizabilities Calculated with a Modified Dipole Interaction. Chem. Phys. 1981, 59, 341−350. (66) Friesner, R. A. Ab Initio Quantum Chemistry: Methodology and Applications. Proc. Natl. Acad. Sci. U. S. A. 2005, 102, 6648−6653. (67) Headgordon, M.; Pople, J. A.; Frisch, M. J. MP2 Energy Evaluation by Direct Methods. Chem. Phys. Lett. 1988, 153, 503−506. (68) Frisch, M. J.; Head-Gordon, M.; Pople, J. A. Semi-Direct Algorithms for the MP2 Energy and Gradient. Chem. Phys. Lett. 1990, 166, 281−289. (69) Balabanov, N. B.; Peterson, K. A. Systematically Convergent Basis Sets for Transition Metals. I. All-Electron Correlation Consistent Basis Sets for the 3d Elements Sc−Zn. J. Chem. Phys. 2005, 123, 064107. (70) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision C.01; Gaussian, Inc.: Wallingford, CT, 2009. (71) Maroulis, G. Atoms, Molecules and Clusters in Electric Fields: Theoretical Approaches to the Calculation of Electric Polarizability; Imperial College Press: London, United Kingdom, 2006. (72) Zhao, Y.; Truhlar, D. G. A New Local Density Functional for Main-Group Thermochemistry, Transition Metal Bonding, Thermochemical Kinetics, and Noncovalent Interactions. J. Chem. Phys. 2006, 125, 194101. (73) Dolg, M.; Wedig, U.; Stoll, H.; Preuss, H. Energy-Adjusted Ab Initio Pseudopotentials for the First Row Transition Elements. J. Chem. Phys. 1987, 86, 866−872. (74) Francl, M. M.; Pietro, W. J.; Hehre, W. J.; Binkley, J. S.; Gordon, M. S.; DeFrees, D. J.; Pople, J. A. Self-Consistent Molecular Orbital Methods. XXIII. A Polarization-Type Basis Set for Second-Row Elements. J. Chem. Phys. 1982, 77, 3654−3665. (75) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297−3305. (76) Boys, S. F.; Bernardi, F. Calculation of Small Molecular Interactions by Differences of Separate Total Energies - Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553−566. K


Article

The Journal of Physical Chemistry B Simple Ligands and Fulvic/Humic Substances. Environ. Sci. Technol. 2007, 41, 7609−7620. (98) Lide, D. R.; Kehiaian, H. V. CRC Handbook of Thermophysical and Thermochemical Data; CRC Press Inc: Boca Raton, FL, 1994. (99) Friedman, H. L.; Krishnan, V. V. In Water: A Comprehensive Treatise; Franks, F., Ed.; Plenum Press: New York, 1973. (100) Marcus, Y. Thermodynamics of Solvation of Ions. Part 5.Gibbs Free Energy of Hydration at 298.15 K. J. Chem. Soc., Faraday Trans. 1991, 87, 2995−2999. (101) Miliordos, E.; Xantheas, S. S. Ground and Excited States of the [Fe(H(2)O)(6)](2)(+) and [Fe(H(2)O)(6)](3)(+) Clusters: Insight into the Electronic Structure of the [Fe(H(2)O)(6)](2)(+)-[Fe(H(2)O)(6)](3)(+) Complex. J. Chem. Theory Comput. 2015, 11, 1549−1563. (102) Sligar, S. G.; Makris, T. M.; Denisov, I. G. Thirty Years of Microbial P450 Monooxygenase Research: Peroxo-Heme Intermediates–The Central Bus Station in Heme Oxygenase Catalysis. Biochem. Biophys. Res. Commun. 2005, 338, 346−354. (103) Wondimagegn, T.; Rauk, A. The Structures and Stabilities of the Complexes of Biologically Available Ligands with Fe(III)− Porphine: An Ab Initio Study. J. Phys. Chem. B 2011, 115, 569−579.

L


Polarizable and Non-Polarizable Force Field Representations of Ferric

Recommend Documents