PCCP View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PAPER

Cite this: Phys. Chem. Chem. Phys., 2014, 16, 3693

View Journal | View Issue

Hydration properties of lanthanoid(III) carbonate complexes in liquid water determined by polarizable molecular dynamics simulations† e Fausto Martelli,ab Yannick Jeanvoine,ac Thomas Vercouter,d Ce ´sar Beuchat, fgh ac Rodolphe Vuilleumier and Riccardo Spezia*

In this work we have studied the structure and dynamics of complexes formed by three and four carbonates and a central lanthanoid(III) ion in liquid water by means of polarizable molecular dynamics simulations. With this aim we have developed a force field employing an extrapolation procedure that was previously developed for lanthanoid(III) aqua ions and then we have validated it against DFT-based data. In this way we were able to shed light on properties of the whole series, finding some similarities and differences across the series, and to help in interpreting experiments on those systems. We found that the bi-dentate tri-carbonate complexes are the most stable for all the atoms, but a variation of the number of water molecules in the first ion shell, and the associated exchange dynamics, is observed from lighter to heavier elements. On the other hand, for four-carbonate systems only one water molecule is observed in the first shell, with 10–20% probability, for La(III) and Ce(III), while for the rest of Received 21st September 2013, Accepted 10th December 2013 DOI: 10.1039/c3cp54001d

the series it seems impossible for a water molecule to enter the first ion shell in the presence of such an excess of carbonate ligands. Finally, the good performance of our extrapolation procedure, based on ionic radii, makes us confident in extending such approaches to study the structure and dynamics of other systems in solution containing Ln(III) and An(III) ions. This parametrization method results particularly useful since it does not need expensive quantum chemistry calculations for all the atoms in

www.rsc.org/pccp

the series.

1. Introduction The structure and dynamics of water around lanthanoid (Ln) ions, mainly at oxidation state III, have been widely studied both experimentally and theoretically.1–4 Behavior of Ln(III) in water is in fact related to various fields such as rare earth elements chemistry, environmental impact of pollutants, a

CNRS, Laboratoire Analyse et Mode´lisation pour la Biologie et l’Environnement, UMR 8587, France b Frick Chemistry Laboratory, Department of Chemistry, Princeton University, Princeton, 08544 USA c Universite´ d’Evry Val d’Essonne, UMR 8587 LAMBE, Boulevard F. Mitterrand, 91025 Evry Cedex, France. E-mail: [email protected]; Tel: +33 (0)169 47 01 41 d CEA, DEN, DANS, DPC, SEARS, LANIE, F-91191 Gif sur Yvette cedex, France e Department of Physical Chemistry, University of Geneva, 30 Quai Ernest Ansermet, CH-1211 Geneva, Switzerland f Ecole Normale Supe´rieure, De´partement de Chimie, 24, rue Lhomond, 75005 Paris, France g UPMC Univ Paris 06, 4, Place Jussieu, 75005 Paris, France h UMR 8640 CNRS-ENS-UPMC, France † Electronic supplementary information (ESI) available. See DOI: 10.1039/ c3cp54001d

This journal is © the Owner Societies 2014

radioactive waste remediation and medical imaging.5–7 In the context of radioactive waste management, these elements are used as chemical analogues of actinide ions (An).8 Their interactions with inorganic anions may control their speciation in solution with the formation of aqueous complexes. The Ln(III) ions can accept several ligands in their first coordination layer by replacing all the hydration molecules. Carbonate ligands have gained major attention due to environmental and process applications. Interestingly, reliance on differential lanthanide– carbonate interactions has been proposed as a possible separation procedure for Ln(III) and An(III) ions in solution.9 Up to three or four carbonate anions can bind to Ln(III) trications in concentrated carbonate solutions. Several experimental studies have been carried out in order to determine the exact stoichiometry though a clear consensus has not yet been reached. Crystallographic data for Ln(III) carbonate hydrates are available for tri-carbonate ligands,10 and for Nd(III). Runde et al.11 have suggested the formation of [Nd(CO3)4H2O]5 structure at high carbonate concentrations. Recently several Ln(III)–carbonate complexes have been studied in solution using electrophoretic mobility measurements and time-resolved laser-induced fluorescence spectroscopy (TRLFS).12–14 They concluded that light Ln(III) ions

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3693

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

Paper

can coordinate up to four carbonate ligands while heavier ones coordinate only up to three ligands, in line with the effect of the lanthanoid ion contraction. In contrast, considering available crystallographic and spectroscopic data (including UV-Vis, near infrared, and infrared), Janicki et al. concluded that in aqueous solution all Ln(III) ions form tetra-carbonates when carbonate ions are not limited.15 Another recent theoretical contribution in this area was a report by Sinha et al. on [Nd(CO3)4]5 using the Parameterized Model 3 (PM3) semi-empirical method.16 Notwithstanding these two studies, the first quantitative and systematic theoretical study undertaken to characterize the structures and bonding of lanthanoid(III) tri- and tetra-carbonates was reported recently by Jeanvoine et al.17 They studied different clusters by means of Density Functional Theory (DFT) pointing out that: (1) in the gas phase the most stable structure is the full bi-dentate for tri-carbonate complexes, while for the tetra-carbonate complexes it is the full mono-dentate; (2) in water, modeled there as a continuum – i.e. no explicit water molecules were considered – in both cases the full bi-dentate complexes are the most stable structures; (3) the tri-carbonate structure is more stable than the tetra-carbonate one; (4) the Ln(III)–carbonate interaction is mainly ionic, thus allowing use of a simple, physics based, interaction potential to study such complexes in liquid water by means of classical molecular dynamics simulations. Molecular dynamics (MD) simulations were performed recently with success to study hydration of lanthanoids(III) in water or in salts, in particular when polarizable force fields are adopted.18–23 In aqueous solution with non-coordinating counterions, the difference in ionic radius for two elements gives rise to a difference in hydration number across the series (9-fold vs. 8-fold for La and Lu, respectively).24 Ln(III)–water interactions have been determined to be mainly electrostatic in nature, as one might expect given the ‘‘hard’’ characters of both Ln(III) ions and water. As such, the variation in ionic radius is the main physical quantity that affects hydration properties.25 The fact that ionic radii can dictate the complexation properties has also been pointed out for the case of ligands that are potentially less hard than water, like hexacyanoferrate.26 Albeit carbonates are ligands softer than water, and it is also possible that the metal–ligand interaction may change across the spectrum of the lanthanoid series, we have recently shown that also in this case the interaction is mainly electrostatic and the contribution of 4f orbitals to Ln/carbonate interaction is negligible.17 Since 4f orbitals are compact around lanthanoids, they rarely contribute to valence bonding; indeed this behavior rationalizes the importance of the role that ionic radius plays in dictating interactions with water as a ligand.27 As we will show in the present study, the same physical picture can be evoked and used to successfully develop a polarizable potential in order to study the complexes between all the Ln(III) ions in the series and carbonates in water. MD simulations allow for an explicit description of the solvent and inclusion of temperature and entropic effects that are so important in the liquid state. Recently, the presence of Cl, ClO4 and NO3 anions in the first shell of Ln(III) was studied using polarizable force fields,19,28,29 but, with the

3694 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

PCCP

aforementioned exceptions, few studies were devoted to study complexes in liquid water using molecular dynamics simulations. In the present work we have thus developed a force field for classical molecular dynamics of such complexes in liquid water. In particular, we have developed a polarizable force field for La(III) (the first ion in the series) interacting with carbonate ions and then extrapolated it to Lu(III) following the same procedure developed for Ln(III)– and An(III)–water interaction.30,31 Results obtained from classical molecular dynamics simulations were validated against DFT results on clusters17 and discussed with respect to DFT-based molecular dynamics. The good performance of the extrapolation procedure to Lu(III), which is based on the picture obtained from a previous study of electron density, pointing out the ionic nature of Ln–carbonate interaction,17 gave us confidence to use the same procedure to study the whole series. We have thus extended the polarizable force field to the whole series and performed molecular dynamics simulations of Ln(III)–carbonate clusters in liquid water using the water ionic radii recently reported by D’Angelo et al.18 We were then able to study structural and energetic properties of Ln(III)–carbonate complexes in liquid water. In particular, we found that the tri-carbonate structures are the most stable structures in both the gas phase and solution and that the bi-dentate complexes are spontaneously formed in liquid water when starting from mono-dentate ones, thus confirming results obtained from calculations done on clusters. Furthermore, we point out that for tri-carbonate structures, a variable number of water molecules can enter the first Ln(III) hydration shell: from 4 to 2 moving from lighter to heavier elements. In the case of tetra-carbonate structures, only one water molecule can come into the first hydration shell, as well as for the first two elements of the series, La(III) and Ce(III).

2. Methods 2.1

Classical force field

Our main goal is to provide a classical force field in order to study liquid systems composed of Ln(III) and carbonate ions in water by means of classical molecular dynamics simulations (CLMD). With this aim, our classical potential energy function consists of the sum of different terms: LJ LJ Vtot = Velec + V LJ water–water + V water–carb + V carb–carb + VLn–water + VLn–carb (1)

where Velec is the electrostatic energy term composed of a Coulomb and a polarization term (when activated) following LJ Thole’s induced dipole model,32 V LJ water–water, V water–carb and LJ V carb–carb are the 12-6 Lennard-Jones potential terms describing the interactions between water and carbonates. In Thole’s model, each site carries one permanent charge and an induced dipole associated with an isotropic polarizability tensor. To avoid the polarization catastrophe in Thole’s model, a screening function is employed for dipole–dipole interaction at short distances, using the original exponential function, with the

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

Paper

same parameters used in our previous studies.33 Carbonate Lennard-Jones parameters were obtained from the literature34 and for water we used those of the TIP3P water model.35 VLn–water and VLn–carb account for the non-electrostatic Ln–water and Ln–carbonate interaction potential for which we have used a potential composed of a long range attractive part with a 1/r6 behavior and a short range repulsive part modeled via an exponential function, dealing with the well-known Buckingham potential V Buck ðrÞ ¼ A expðBrÞ 

C r6

(2)

where r is the distance between Ln and O atoms of water and carbonate. This expression applied to Ln(III) and An(III) in water was shown to be able to correctly reproduce various properties.30,31 This interaction potential, in conjunction with the TIP3P/P model for water, was able, for example, to accurately reproduce EXAFS spectra both for Ln(III) and An(III)18,36 and thus for Ln(III)–water interactions we used the same parameters developed previously.18 Parameters of the Buckingham potential can be related to ionic radii, as shown recently by our group31 and in previous work of Madden and co-workers in the case of molten salts.37–39 In particular, we have shown that starting from B and C reference values of an atom of the Ln (or An) series, we can extrapolate values for the whole series. To this end, we have first obtained B and C parameters for La(III) in interaction with carbonates. Here, as previously, we constrained the fit by keeping A constant since it modulates the height of the barrier that is never reached for liquid systems at room temperature, and we used the same value,30,33 i.e. 240 000 kcal mol1. The B and C parameters were obtained in order to reproduce DFT results in gas and liquid phases. We have explored the performances of both polarizable and non-polarizable potentials, such that we consider two force fields labeled POL-ff and UPOL-ff, respectively, and corresponding molecular dynamics (POL-CLMD and UPOL-CLMD). Then, assuming the same extrapolation rules for water–Ln(III) and water–An(III) interaction,18,30,31,36 we have obtained values for Lu(III), and then for the whole series. We should point out that this is not a strict derivation, but an extrapolation based on the assumption that the ionic radii decrease across the series similarly to what happens in water. The validity of the extrapolation is checked by comparing results of the extrapolated force fields with DFT-based calculations. Extrapolation was done for both POL-ff and UPOL-ff. The electrostatic term is composed of partial charges and atomic polarizabilities for POL-ff; for UPOL-ff only partial charges were used. For water we used TIP3P35 charges and TIP3P/P40 charges and atomic polarizabilities for UPOL-ff and POL-ff, respectively. For Ln(III) we used: charge q = +3; for ionic polarizability we employed the values reported and used in our force field developed for hydration.30,41 For carbonates, we used partial charges reported in the literature34 for the UPOL-ff, while for POL-ff we obtained them via the Lo-PROP procedure42 based on MP2 calculations using an ANO-RCC-VQZP basis

This journal is © the Owner Societies 2014

Table 1 Charges and atomic polarizabilities used in both non-polarizable (UPOL-ff) and polarizable (POL-ff) force fields. CC and OC are C and O of carbonates, while OW and HW are O and H of water. We also report intramolecular parameters for carbonates used

La Lu OW HW CC OC

C–OC

OC–C–OC

OC–C–OC–OC

UPOL-ff

POL-ff

q

Q

a (Å3)

+3 +3 0.82 0.41 0.68 0.89

+3 +3 0.66 0.33 0.78 0.93

1.41 0.77 0.85 0.41 0.52 0.71

req (Å)

Kr (kcal mol1 Å2)

1.31

383.3

yeq (1)

Ky (kcal mol1 rad2)

120.0

200.0

xeq

Kx (kcal mol1 rad2)

180.0

300.0

set.43 These calculations were done using the MOLCAS package.44 The same Lo-PROP procedure was used to obtain atomic site polarizabilities of carbonates. Partial charges and atomic site polarizabilities used in both UPOL and POL force fields are summarized Table 1. For intermolecular potential of carbonates we used values reported by Bruneval et al.,45 where the Morse potential was replaced by harmonic potential and for the equilibrium distance we used the value suggested for liquid water.46 These parameters are reported in Table 1. To check the effect of modifying these parameters on results we have performed some simulations with higher values of Kr and Ky, 570 kcal mol1 Å2 and 400 kcal mol1 rad2 respectively. 2.2

Classical molecular dynamics simulations set-up

CLMD simulations in liquid water were done as follows. A cluster composed of one Ln(III) and 3 or 4 carbonates (CO3)2 was immersed in a box composed of 216 or 512 water molecules. Box sizes have been obtained with 5 ps of equilibration in the NpT ensemble with the DL_POLY classic package47 using the UPOL force field. The Berendsen thermostat and barostat were employed48 resulting in box edges of 19.5 and 25.05 Å (for 216 and 512 water molecule boxes respectively) for T = 300 K and p = 1 atm. Periodic boundary conditions (PBC) were applied to each simulation box in order to mimic bulk conditions. Long-range interactions were calculated by using the smooth particle mesh Ewald method.49 Once the system equilibrated we switch to the NVE ensemble and use a velocityVerlet algorithm with a time step of 1 fs. The SHAKE algorithm was used to constrain bonds and angles of water molecules.50 Trajectories were thus 3 ns long for each system. POL-CLMD simulations were started after UPOL NpT equilibration as described previously. In these simulations we employed the extended Lagrangian method to propagate induced dipoles in time51 that was shown to work properly with Thole’s induced dipole model for similar systems. Note that the dynamics of the

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3695

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

Paper

PCCP

induced dipole degree of freedom is fictitious, it only serves the purpose of keeping the induced dipoles close to their values at minimum energy, which would be obtained from the exact solution of the self-consistent equations (SCF) at each step. Dipoles and temperature stability were verified along the present simulations (induced dipoles obtained with SCF are within the oscillations observed using the dipole dynamics) while a more detailed report of the performances of the extended Lagrangian implementation is reported elsewhere.52 Equations of motion were numerically integrated using a 1 fs time step, with the SHAKE algorithm50 to constrain bonds and angles of water molecules. The system was equilibrated at 298 K for 5 ps. Production runs were subsequently collected for 3 ns. The 216 water molecule box was used to develop and validate the force fields and then final results were checked with the 512 water molecules box. No big differences were noticed in results obtained with two boxes, as previously observed.33,53,54 We have studied the complexation between Ln(III) ions and 3 or 4 carbonate ions. The starting points of the dynamical simulations are all the possible arrangements of 3 and 4 carbonate ions, i.e. mono-dentate (Z1) or bi-dentate (Z2). The general formula we adopted to identify those clusters is: [Ln(Z1-CO3)l(Z2-CO3)m]k, where l + m = 3 or 4 and k = 3- or 5- respectively. Representative illustrations of the clusters considered are reported in Fig. 1. As initial structures in MD simulations we have used clusters previously optimized at the DFT level in our group17 and then immersed in explicit water boxes. The code MDVRY52 was used for production dynamics with both UPOL and POL force fields. 2.3

Car–Parrinello molecular dynamics

In order to have a DFT comparison where solvent water molecules are explicitly taken into account, we performed Car–Parrinello molecular dynamics55 (CPMD) simulations of La(III) and Lu(III) clusters with 3 and 4 carbonates in full bi-dentate binding conformation. We immersed the clusters previously optimized with DFT17 in a box of 120 water molecules and after an NpT equilibration using classical dynamics (as reported in Section 2.1) in order to set the correct box dimensions for T = 300 K and p = 1 atm; then we performed NVE simulations. Periodic boundary conditions were employed with a box length of 15.5 Å as obtained by NpT equilibration. The BLYP functional56,57 was employed with plane waves basis set and a cut-off of 110 Ry. For C, O and H we used standard Troullier–Martins semi-core pseudo-potentials58 used previously in similar CPMD simulations.59–61 For La(III), we used a semi-core Troullier–Martins pseudo potential developed recently by some of us.27 For Lu(III) we developed a Troullier– Martins pseudo-potential as follows. The reference configuration used to generate the PP for Lu was Lu(III). The orbitals 5s, 5p, and 5d were included in the PP with cutoffs of 1.01, 1.15 and 1.163 au, respectively. 4f orbitals were considered as core orbitals, and a non-linear core correction was thus employed.

3696 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

Fig. 1

Chemical structures of different Ln(III)–carbonate clusters.

When using pseudo-potentials with plane waves, the semilocal Kleinman–Bylander form was used with the p channel as the local channel.62 Molecular dynamics simulations were done

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

using a fictious mass of 150 au and a time-step for numerical integration of 2 au (= 0.0484 fs). For La(III) with three carbonates we have performed two simulations with different initial conditions: one obtained from a classical equilibration keeping the complex fixed and relaxing water molecules around; the other was obtained after a classical equilibration where the cluster was allowed to move freely (using the Pol-ff described previously). In the first case we do not have any water molecule at t = 0 in DFT-based dynamics while in the second case we have four molecules – in agreement with Pol-ff results shown in the following section. Thus we label these two simulations CPMD(0) and CPMD(4). In the case of Lu(III) with three carbonates we have performed two simulations with one water molecule in the first shell at t = 0, called CPMD(1)S1 and CPMD(1)S2. We have also done simulations with two water molecules at t = 0 in the first shell taken from the Pol-ff equilibration. In this case, when using the CP Lagrangian, the kinetic energy of electrons quickly diverged, and thus for this specific simulation we performed a Born–Oppenheimer dynamics simulation, with an optimization convergence of 106 and an extrapolation scheme of order 3 and a time step of 20 au (= 0.48378 fs). This simulation is labeled BOMD(2). Trajectories, after an equilibration of 0.5 ps, were then accumulated for 2–8 ps depending on the system. CPMD package was employed for these simulations (CPMD).63 2.4

DFT optimizations

Finally, based on the structures obtained from simulations in liquid water of La(III) and Lu(III) complexes we have optimized Ln–carbonate complexes where the first shell water molecules were treated explicitly and the bulk as a continuum. With this aim we have used the PBE functional64,65 with an all-electron triple-z plus two polarization function basis set on all atoms. Relativistic corrections were introduced by the scalar relativistic zero-order regular approximation (ZORA).66,67 Implicit solvation was considered via the COSMO68 solvent model with standard radii except for La (R = 2.42 Å) and Lu (R = 2.24 Å) centers.69 Dispersion corrections were added using the Grimme method with D3 parameters.70 We have optimized structures of La(III) with three carbonates and four or three water molecules in the first shell, labeled PBE/COSMO(4) and PBE/COSMO(3), Lu(III) with three carbonates and two or one water molecules in the first shell, labeled PBE/COSMO(2) and PBE/COSMO(1) and La(III) with four carbonates and one water molecule in first shell, labeled PBE/COSMO(1). These structures were selected on the basis of results obtained in bulk water. These calculations were performed using the Amsterdam Density Functional program (ADF 2012.01) developed by Baerends and co-workers.71

3. Results and discussion Before describing and discussing our results it is useful to define a nomenclature used here and hereafter for atoms of the

This journal is © the Owner Societies 2014

Paper

systems studied. Ln is for the generic lanthanoid(III) ion (if not specified), OW and HW for oxygen and hydrogen atoms of the water molecules, OC(n) for oxygen atoms of carbonates that are in the first shell of Ln(III), OC(f) for the oxygen atoms of carbonates that are exposed to the solvent and C for the carbon atoms of carbonates.

3.1

La(III) and Lu(III) force field development and validation

First we have obtained non-electrostatic parameters for the La(III)–carbonate interaction. With this aim we obtained La(III)–OC B and C parameters (see eqn (2)) in order to better reproduce DFT results obtained with the B3LYP functional.17 As done in our previous studies,30,31 we kept fixed the A parameter. Starting from values developed for La(III)–water, we adjusted the B and C parameters in a systematic way in order to provide structural properties for the [La(Z2-CO3)3]3 complex in liquid water as close as possible to Ln–OC(n) distances resulting from DFT-optimization and DFT-based dynamics. In this way we obtained UPOL and POL parameters. In Table 2 we report the final parameters that best reproduce known Ln–OC(n) distances of the [La(Z2-CO3)3]3 complex. We should note that once we obtained B and C values for UPOL-ff, we only had to adjust the C parameter to obtain a POL-ff that is in agreement with the DFT calculations. We have then extrapolated the Lu(III) parameters following the procedure developed (and justified) in our previous work:30,31 A is kept fixed across the series, B is increased as a function of the ionic radius differences and C is decreased linearly as a function of the ionic radius behavior across the series. With this aim, we have used effective ionic radii obtained recently18 and reported in Table 2. Results provided by UPOL-ff and POL-ff for [La(Z2-CO3)3]3, [Lu(Z2-CO3)3]3, [La(Z2-CO3)4]5 and [Lu(Z2-CO3)4]5 complexes are reported in Table 3, where we compare them with structural results from DFT optimization,17 DFT-based dynamics with different initial conditions (i.e. different number of water molecules in the first shell of the cation) and geometry optimizations of clusters with explicit water molecules in the first shell and the implicit solvation model (simulations labeled PBE/COSMO). We can see that our force field is able to correctly reproduce in particular Ln–OC(n) distances of the different complexes; that was the principal aim of our parametrization. Ln–C and Ln–OC(f) distances are underestimated in classical simulations. This is due mainly to different OC(n)COC(n) angles and COC(n) distances of the carbonates, while COC(f) distances are almost the same in both classical and DFT results, as reported in detail in ESI†

Table 2 Interaction parameters of UPOL and POL force fields for the Buckingham potential (eqn (2)) between La3+, Lu3+ and carbonate oxygen atoms. A is fixed at 240 000 kcal mol1 in all cases. We report also effective ionic radii as reported by D’Angelo et al.18

3+

La (UPOL) La3+ (POL) Lu3+ (UPOL) Lu3+ (POL)

B (Å1)

C (Å6 kcal mol1)

Ri (Å)

3.47 3.47 3.725 3.725

9257.0 9157.0 6135.5 6035.5

1.250 1.250 0.995 0.995

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3697

View Article Online

Paper

PCCP

Table 3 Structural results for [La(Z2-CO3)3]3, [La(Z2-CO3)4]5, [Lu(Z2-CO3)3]3, [Lu(Z2-CO3)4]5 complexes in bulk water obtained from UPOL-ff, POL-ff, PBE/COSMO and CPMD calculations. We report the maximum of radial distribution functions. Ln–OC(n) is the distance from the oxygen atoms of carbonates in the first shell of Ln, while Ln–OC(f) is the remaining oxygen atom to Ln distance; Ln–OW is the distance of Ln from water molecules and in parentheses we report the number of water molecules in the first shell (if any)

Method

Ln–OC(n)

Ln–OC(f)

Ln–C

Ln–OW

[La(Z -CO3)3] [La(Z2-CO3)3]3 [La(Z2-CO3)3]3 [La(Z2-CO3)3]3 [La(Z2-CO3)3]3 [La(Z2-CO3)3]3 [La(Z2-CO3)3]3 [La(Z2-CO3)3]3

UPOL-ff POL-ff POL-ff (Ky = 400) POL-ff (Kr = 570) DFTa CPMD(0) CPMD(4) PBE/COSMO(4)

2.50 2.49 2.49 2.49 2.47 2.485 2.54 2.52

4.09 4.09 4.07 4.06 4.19 4.19 4.27 4.23

2.82 2.77 2.80 2.79 2.93 2.93 3.01 2.97

2.75 2.75 2.76 2.76 — 2.63 2.68 2.73

[La(Z2-CO3)3]3

PBE/COSMO(3)

2.52

4.22

2.95

2.74 (3)

[La(Z2-CO3)4]5 [La(Z2-CO3)4]5 [La(Z2-CO3)4]5 [La(Z2-CO3)4]5 [La(Z2-CO3)4]5

UPOL-ff POL-ff DFTa CPMD PBE/COSMO

2.53 2.48 2.64 2.49 2.55

4.09 4.01 4.42 4.20 4.24

2.83 2.75 3.12 2.97 2.98

2.62 (1.98) 2.72 (0.17) — — 2.71 (1.00)

[Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)3]3

UPOL-ff POL-ff DFTa CPMD BOMD(2) PBE/COSMO(2) PBE/COSMO(1)

2.32 2.29 2.27 2.29 2.29 2.30 2.28

3.89 3.82 3.97 4.00 4.03 3.99 3.96

2.62 2.53 2.71 2.67 2.73 2.73 2.70

2.42 2.50 — 2.30 2.38 2.50 2.38

[Lu(Z2-CO3)4]5 [Lu(Z2-CO3)4]5 [Lu(Z2-CO3)4]5 [Lu(Z2-CO3)4]5

UPOL-ff POL-ff DFTa CPMD

2.33 2.33 2.48 2.30

3.89 3.86 4.20 4.04

2.63 2.60 2.90 2.75

2.45 (0.86) — — —

System

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

2

a

3

(3.91) (3.76) (3.78) (3.79) (1.20) (3.0) (3) 3.79 (1)

(2.99) (2.00) (1.00) (2.00) (2.00) (1.00)

DFT B3LYP optimization results from Jeanvoine et al.17 All distances are in Å.

(Tables S1 and S2). Furthermore the LnCOC(f) angles seem to deviate more from linearity in DFT-based simulations while in classical and DFT optimized geometries they are linear. In particular, what is promising for our parametrization procedure is that the Lu–OC(n) distances obtained from classical MD simulations are in good agreement with DFT ones. We should remark that the Lu(III)–carbonate force field is extrapolated and not adjusted to reproduce DFT results: parameters were adjusted only to reproduce the structures of tri-carbonate complexes with La(III), and they give reasonable results (i.e. without any further adjustment), as well as for tetra-carbonate systems and for tri- and tetra-carbonate Lu(III) complexes. We now move to analyze in detail results obtained from MD simulations. In Fig. 2 and 3 we report radial distribution functions (RDFs) between the ions and the oxygen atoms of carbonates and water molecules obtained from POL-ff and CPMD simulations for simulations containing three carbonates. The La–OC first peak is very well reproduced while the second La–OC peak, corresponding to the distance with the ‘‘far’’ oxygen atom of the carbonate, is shorter in the classical simulations than in the CPMD ones. Note that in classical dynamics the C–OC bond is modeled using a harmonic potential where the equilibrium distance is 1.31 Å and all the bonds are equal, while in DFT simulations this is not the case. This is an intrinsic limit of classical simulations, and we used an equilibrium value in the classical force field of 1.31 Å,

3698 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

Fig. 2 La–OC and La–OW radial distribution functions obtained from POL-ff, CPMD(0) and CPMD(4) simulations for [La(Z2-CO3)3]3 simulations.

which is reported for C–O distances of carbonates in explicit water at the Hartree–Fock level;46 this geometry seems the most adapted for simulating free carbonates in water and, in perspective, more complex mixtures of lanthanoid/carbonate solutions or even crystal growth can be studied with the same (or similar) force field. When the force constants of CQO

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

Fig. 3 Lu–OC and Lu–OW radial distribution functions obtained from POL-ff, CPMD(0) and BOMD(2) simulations for [Lu(Z2-CO3)3]3 simulations.

stretching and OCO bending are changed, the results are not largely modified (as reported in Table 3, Tables S1 and S2, ESI†). We should also note that the force field gives COC(f) distances that are in good agreement with DFT ones. The difference occurs when oxygen atoms (partially negatively charged) are in the vicinity of a 3+ ion. As shown by results obtained by increasing Kr and Ky that are very similar to what was obtained with other parameters, the discrepancy in the COC(n) distance and the OC(n)COC(n) angle should come from intrinsic limitations of harmonic potentials or more in general classical models that do not take into account resonance hybrids. Explicit polarization through atomic site induced dipoles provide more physical results than a simple fixed charge model, but it is probably not able to catch all aspects of the interactions involved. This is reflected in differences in hydration of OC(f), as reported in ESI† (Fig. S11). Similarly to what was obtained by Leenders et al.,72 the classical force fields tend to form more structuration between the water molecules and the carbonates. Furthermore, even if our polarizable potential gives an improved description of electrostatics with respect to non-polarizable models, surely it misses the full charge transfer effect that is only partially taken into account by Buckingham parameters (implicitly) and by means of induced dipoles. NBO analysis found that, even if small, there is a partial charge transfer contribution in Ln–carbonate binding.17 Here our aim was to provide a first polarizable classical potential able to reliably describe Ln(III)–carbonate systems that can be used for further studies with typical advantages over the limitations of classical MD simulations. We now move to analyze the structure of (possible) water molecules in the first shell of Ln(III)–carbonate complexes. In Fig. 2 and 3, we report the RDFs between La(III) and Lu(III) and oxygen atoms of water molecules. In the case of [La(Z2-CO3)3]3 complexes (Fig. 2), from POL-ff simulations we obtained a dynamical equilibrium between three and four water molecules at a distance of 2.75 Å (maximum of RDF). When optimizing clusters with explicit water molecules, we have results very similar in distances (average distances are shown in Table 3).

This journal is © the Owner Societies 2014

Paper

Interestingly we saw that when four waters are in the first shell, three are at 2.73 Å and one is far away – nearly in a second shell. We then performed CPMD simulations with different initial conditions: CPMD(0) where there are no water molecules in the first shell at t = 0 and CPMD(4) with 4 water molecules at t = 0. In the CPMD(0) simulations the cluster needs some time to open and in 3 ps we have two water molecules that are able to come into the first shell. Thus, from CPMD(0) simulation we obtained a RDF that is different than what was obtained from Pol-ff simulations, but also distances that are much shorter than what resulted from DFT optimizations. To let water molecules enter, the cluster must open (as shown by CLaC angles in ESI,† Fig. S1). We thus run a CPMD simulation with four water molecules at the beginning. In this case one water molecule leaves quickly and then three water molecules are stable in the first shell for 8 ps. Now the RDF is more in agreement with the Pol-ff in terms of height, integral and also the position is now closer (2.68 Å) to Pol-ff and PBE/COSMO results. A similar situation is obtained in the case of [Lu(Z2-CO3)3]3 complexes for which RDFs are shown in Fig. 3. In the case of Pol-ff simulations we have two water molecules at 2.50 Å, in agreement with DFT optimization with two water molecules – PBE/COSMO(2) results in Table 3. In DFT-based dynamics by increasing the number of water molecules in the first shell the Lu–OW distance increases, in agreement with PBE/COSMO results. When comparing Pol-ff and DFT-based dynamics with two water molecules in first shell, we obtain that DFT tends to provide shorter distances. This discrepancy should result from differences in Lu(III)–water interaction. We should remark that the POL-ff Lu(III) water interaction potential was shown to very accurately reproduce structural and thermodynamical experiments,53 and PBE/COSMO optimizations provide the same Lu(III)–water distance. Furthermore, CLnC angles obtained from POL-ff dynamics and PBE/COSMO optimizations (shown in ESI,† Tables S1 and S2 and Fig. S2 and S4) are similar. The optimized geometries are reported in ESI† (Table S3) and can be used to have a simple static picture of water–carbonate arrangement – more details will be given in Section 3.3. It could be that BLYP-based dynamics, while describing well structures of complexes, are not highly accurate in describing ion–water interaction (that can be better described by a well parameterized force field), as we have pointed out in the case of cobalt–cysteine complexes.60,73,74 Furthermore, it is well known that BLYP-based dynamics at T = 300 K give an overstructurated water.75 Here we decided to perform dynamics at this temperature since increasing it (as it is sometimes done) can result in a problematic description of Ln–carbonate complexes, which is the main subject of this work. Finally, we inspect differences between POL-ff and UPOL-ff results. The Ln–OC(n), Ln–OC(f) and Ln–C distances are very similar, in particular for [La(Z2-CO3)3]3, while they are slightly longer in UPOL-ff results for the other complexes. We should also note that the extrapolation procedure used to obtain the Lu(III)–carbonate force field was developed using a polarizable potential and thus it is not surprising that it is less accurate in

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3699

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

Paper

the case of a non-polarizable potential. The carbonate structure is also similar (see Tables S1 and S2, ESI†), while the most noticeable difference is in tri-dimensional arrangement of carbonates around the ion, as shown by the CLnC angle dynamics (reported in ESI,† Fig. S2 and S4). As shown previously, these angles are related to the water molecules in the first shell of the Ln(III) ion. The main difference between POL-ff and UPOL-ff results was found mainly in the water structure and dynamics. In particular the UPOL-ff simulations provide higher Ln(III) first shell coordination numbers (CN) than POL-ff, as shown in Table 3: 3.91 vs. 3.76 for [La(Z2-CO3)3]3, 2.99 vs. 2.0 for [Lu(Z2-CO3)3]3, 1.98 vs. 0.17 for [La(Z2-CO3)4]5 and 0.86 vs. 0.0 for [Lu(Z2-CO3)4]5. This is in agreement to what was observed in La(III) hydration, where UPOL-ff overestimates the number of water molecules in the first shell.33 Furthermore, the dynamics of the first shell water molecules is a property very sensitive to the potential, as largely shown in the case of bare ions in liquid water.22,23 In Fig. 4 we show the number of water molecules in the first shell of La(III) and Lu(III) as a function of time as obtained from UPOL-ff and POL-ff. We note that UPOL-ff provides much less exchange than POL-ff and higher coordination numbers (as obtained for La(III) in bulk water). Thus, while UPOL-ff provides results for Ln(III)–carbonates similar to POL-ff results, it gives different results for structure

Fig. 4 Number of water molecules in the first shell of Ln(III) ion obtained from UPOL-ff (in red) and POL-ff (in black) for [La(Z2-CO3)3]3 (panel a) and [Lu(Z2-CO3)3]3 (panel b) simulations.

3700 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

PCCP

and dynamics of water around the central ion. This difference in the behavior of polarizable and non-polarizable force fields has been already noticed in other situations.22,23,33 This should be taken into account if one wanted to use the UPOL-ff model for simulating bigger systems in the future. 3.2

Energetic properties

Before studying the energetic properties in the liquid phase, we determined using both UPOL-ff and POL-ff the formation energies of different clusters (i.e. for the different coordination motifs shown in Fig. 1) in the gas phase, and compared them to previous DFT results.17 In Fig. 5 we show results for La(III) and Lu(III) as a function of the cluster structures for both tri- and tetra-carbonate complexes. Both POL-ff and UPOL-ff reproduce correctly the trend obtained from DFT: in the gas phase for three carbonate complexes the full bi-dentate is the most stable structure, while for four carbonate complexes the full monodentate is the most stable one. The same behavior is found for La(III) and Lu(III). POL-ff better reproduces the formation energy than UPOL-ff and in particular in the extrapolation to Lu(III). Note that the parameters for Lu(III) were not obtained in order to reproduce structural properties of Lu(III)–carbonate complexes. Parameters were only adjusted for La(III) in order to match structural properties of the [La(Z2-CO3)3]3 complex with a ‘minimal’ modification from the La(III)–water parameters. Concerning the stability of clusters in water, first, we should remark that when performing simulations where the initial structure is a mono-dentate species (full or partial, for both three of four carbonate complexes) we always obtain after a few picoseconds a full bidentate complex (see a prototypical movie in ESI†). We then calculate hydration and formation enthalpies of the stable full bi-dentate species, i.e. [La(Z2-CO3)3]3, [Lu(Z2-CO3)3]3, [La(Z2-CO3)4]5 and [Lu(Z2-CO3)4]5 complexes. To obtain formation enthalpies in water, we used a thermodynamic cycle shown in Fig. 6. Hydration energies are obtained from simulations in the liquid phase, using the

Fig. 5 Gas phase formation energies of [Ln(CO3)3]3 (left) and [Ln(CO3)4]5 (right) clusters (Ln = La, Lu) as a function of the number of mono-dentate ions (from left to right we move from (Z2-CO3)3 and (Z2CO3)4 to (Z1-CO3)3 and (Z2-CO3)4, respectively). Energies calculated with UPOL-ff and POL-ff are compared with values obtained from DFT.17

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

Paper

Fig. 6 The thermodynamic cycle used to evaluate formation enthalpies in aqueous solution from hydration enthalpies and gas phase formation enthalpies.

correction terms described by Garcia and co-workers,76–79 as recently done by us to obtain hydration enthalpies of Ln(III) and An(III).31,53,54 Results for the reaction Ln3+(aq) + nCO32(aq) - [Ln(Z2-CO3)n]32n

(3)

are reported in Table 4. In the same table we report values obtained for the reactions in the gas phase. Note that the gas phase DH overestimates the interaction, but by less than 4% of the interaction energy; that is quite surprising since the parameters were obtained by modifying Ln–water parameters in order to reproduce La–OC distances in tri-carbonate complexes. Solvation has a small effect on formation enthalpies and the DH POL-ff –DH POL-ff is in the 7/+51 kcal mol1 range – small aq gas compared to formation enthalpy values of the order of thousands of kcal mol1. Following the simple Born model one should expect a strong contribution in solvation for small ions, but, as we have recently shown,54 effective Born radii of highly charged cations like lanthanoids(III) are much bigger (in the 7.2–9.5 Å range) than corresponding ionic radii. This corresponds to the fact that a Ln(III) ion is, from a Born point of view, of a size similar to the complex because of the saturation of the dipoles of structured water molecules. On the other hand, the formation of such complexes in liquid water corresponds to a decreasing number of charged atoms exposed to the solvent with respect to free ions. Furthermore, for both La(III) and Lu(III) systems we found that the tri-carbonate complexes are more stable than the tetracarbonate ones in solution, as was similarly found in the gas phase. This strengthens the picture that the tri-carbonate

Table 4 Formation enthalpies (see eqn (3)) obtained from POL-ff simulations. In water values (aq) are results from the thermodynamic cycle of Fig. 6

[La(Z2-CO3)3]3 [La(Z2-CO3)4]5 [Lu(Z2-CO3)3]3 [Lu(Z2-CO3)4]5 DDH(La) DDH(Lu)

a DHDFT (gas)

DH POL-ff (gas)

DH POL-ff (aq)

1228 957 1346 1036 271 310

1272 962 1355 1025 310 330

1258 969 1304 1011 289 293

a DFT B3LYP results from Jeanvoine et al.17 DDH(Ln) are the differences between tri- and tetra-carbonate complex formation enthalpies. All values are in kcal mol1.

This journal is © the Owner Societies 2014

species is the most stable complex in solution. Of course other effects like ionic strength can change the relative stability; here we just provide enthalpy differences under ideal conditions that can be used as a starting point for further evaluation of thermodynamics quantities under different conditions. Another interesting effect is the difference between tri- and tetra-carbonate complexes for La and Lu, represented by the difference between the two formation enthalpies – DDH(Ln) in the same Table 4. Moving from gas to solution, the difference between the two lanthanoids(III) decreases from 20 kcal mol1 (40 kcal mol1 using DFT energies) to only 4 kcal mol1, that is largely within statistical fluctuations. We should note that in the gas phase we have a difference of 20 kcal mol1 between DFT and POL-ff energies; the force field could underestimate the relative energy in the liquid phase, too. 3.3 Extension to the whole Ln(III) series and water exchange dynamics On the ground of the good performances in the extrapolation of the POL-ff from La(III) to Lu(III), we extended it to the whole series, and thus performed simulations of [Ln(Z2-CO3)3]3 and [Ln(Z2-CO3)4]5 complexes in liquid water for the whole Ln(III) series (except Pm that is not naturally available and for which accurate experimental ionic radii are not present). Extrapolated parameters and results for [Ln(Z2-CO3)3]3 complexes are shown in Table 5. The Ln–OC distance (i.e. the maximum of Ln–OC RDF) decreases as a function of the ionic radius as well as the Ln–OW one. In Fig. 7 we show how the Ln–OC distance decreases across the series – here we represent the series by the inverse of the ionic radius – for both complexes. We note that, while for light atoms the Ln–OC distances of tri- and tetracarbonate complexes are very similar, decreasing the ionic radius, in the case of tri-carbonate complexes, the ‘contraction’ is bigger than for tetra-carbonate ones. This is likely due to the steric hindrance that is higher for tetra-carbonate complexes so that a plateau in the distance is obtained for smaller radii. It is interesting to note that while in both series we employed B and C parameters corresponding to the same decrease of ionic radii

Table 5 Buckingham parameters (see eqn (2)) and structural results for the [Ln(Z2-CO3)3]3 series in liquid water as obtained from POL-ff simulations. Ri is the effective ionic radius as reported by D’Angelo et al.18 and CN is the number of water molecules in the Ln(III) first shell

La Ce Pr Nd Sm Eu Gd Tb Dy Ho Er Tm Yb Lu

Ri (Å) B (Å1) C (kcal mol1 Å6) rLn–OC (Å) rLn–OW (Å) CN

t (ps)

1.250 1.220 1.200 1.175 1.140 1.120 1.105 1.090 1.075 1.055 1.040 1.025 1.010 0.995

147.7 203.6 92.1 165.9 191.1 110.2 97.7 104.7 69.6 143.1 172.9 177.0 361.9 591.0

3.470 3.475 3.520 3.545 3.580 3.600 3.615 3.630 3.645 3.665 3.680 3.695 3.710 3.725

9157.0 8789.5 8544.5 8238.2 7809.4 7564.4 7380.6 7196.8 7013.1 6768.0 6584.3 6400.5 6216.8 6033.0

2.49 2.50 2.43 2.42 2.40 2.37 2.37 2.35 2.34 2.32 2.32 2.30 2.29 2.29

2.75 2.75 2.70 2.68 2.64 2.64 2.62 2.61 2.59 2.56 2.57 2.56 2.54 2.53

3.76 3.83 3.25 3.15 3.01 2.93 2.89 2.73 2.60 2.28 2.21 2.14 2.04 2.04

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3701

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

Paper

PCCP

Fig. 7 The Ln–Oc(n) distance across the series (as a function of the inverse of the ionic radius), as obtained from POL-ff simulations in liquid water for both [Ln(Z2-CO3)3]3 and [Ln(Z2-CO3)4]5 complexes.

Fig. 8 Structures of [Ln(Z2-CO3)3]3 and [Ln(Z2-CO3)4]5 complexes obtained from PBE/COSMO geometry optimizations.

Table 6 Structural results for the [Ln(Z2-CO3)4]5 complexes across the Ln(III) series in liquid water as obtained from POL-ff simulations

La Ce Pr Nd Sm Eu Gd Tb Dy Ho Er Tm Yb Lu

rLn–OC (Å)

rLn–OW (Å)

CN

2.48 2.49 2.44 2.43 2.40 2.39 2.38 2.37 2.36 2.35 2.35 2.34 2.34 2.33

2.75 2.75 — — — — — — — — — — — —

0.17 0.27 0.01 — — — — — — — — — — —

along the series, the resulting Ln–OC(n) distances do not decrease with the same slope for tri- and tetra-carbonate complexes. We should remark that our Nd–O distance in the [Nd(Z2-CO3)4]5 complex, 2.42 Å, is not far from the reported X-ray structure of the hydrated crystal, where a value between 2.55 and 2.45 Å in the presence of one water molecule is reported.11 Intriguingly, this Ln–O distance is even closer to La–O or Ce–O distances corresponding to tetra-carbonate complexes where one water molecule can come into the first ion hydration shell (Table 6). The structures of carbonates, in terms of bond distances and angles, are similar for all the complexes (see Table S3 in ESI†). Concerning the tridimensional arrangement around the ion, the CLnC angles show a slightly different dynamics, reflecting different water exchange dynamics as we will discuss in the following paragraph. [Ln(Z2-CO3)3]3 complexes fluctuate more from Pr(III) to Ho(III), while they do less at the end of the series. Similarly for [Ln(Z2-CO3)4]5, the CLnC angles fluctuate for La(III) and Ce(III) while they do not (or slightly) for other ions. As already remarked, the values are similar to what was obtained from PBE/COSMO optimization with explicit water molecules, and thus these geometries can provide a simple picture of how carbonates and water are

3702 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

Fig. 9 Number of water molecules in the first shell of Ln(III) for [Ln(Z2-CO3)3]3 complexes as a function of time from POL-ff simulations.

disposed around the ion. In Fig. 8 we show the optimized structures with CN = 4, 3, 2, 1 for prototypical La(III) and Lu(III) complexes – geometries are given in ESI,† Table S3. We now move to study the number of water molecules (CN) in the first shell of Ln(III) and their exchange dynamics across the series. Results for [Ln(Z2-CO3)3]3 complexes simulations are reported in Table 5 in terms of distances and average CN. They both decrease across the series, as expected from the lanthanoid contraction. Details of the first shell water molecules dynamics for [Ln(Z2-CO3)3]3 complexes across the series are reported in Fig. 9, where we show CN as a function of time for Ln(III) in the beginning, middle and end of the series. We should note, with regard to the first shell, that moving from light to heavy Ln(III), we first have an equilibrium between 4 and 3 water molecules, then in the middle of the series we have 3 water molecules as an average and the possibility of both adding or removing one, and for heavy elements we have 2 water molecules with possibly a third molecule coming in. These results may be compared with experimental estimations of the number of water molecules in Eu(III) and Dy(III) first coordination shells by time-resolved laser-induced fluorescence spectroscopy (TRLFS). An average CN of 1.8 was obtained by fluorescence data fitting for Eu(III) in 1 M Na2CO3,12 while values of 2.1 and 1.9 were obtained for Eu(III) and Dy(III) respectively, in 3 M K2CO3.14

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

Further information on such first shell water molecules dynamics is obtained by calculating the mean residence time (MRT) of water molecules in the first shell of the Ln(III). Here we adopted the so-called Impey method;80 we used here a t* of 2 ps as in Impey’s original work. Results are given in Table 5. We note that from La(III) to Dy(III) we do not have a clear behavior of MRT across the series as was obtained in the case of bare ions in water.24 This is probably due to the fact that while for the light elements we have an equilibrium between CN = 4 and CN = 3, then in the middle we have equilibrium between CN = 4, CN = 3 and CN = 2 and, at the end, we have the equilibrium only between CN = 3 and CN = 2. Upon increasing the atomic number (i.e. moving from light to heavy ions) this low coordination number becomes much more probable and thus the time spent by water molecules in the first shell of Ln(III) is longer. In fact, from Dy(III) to Lu(III) we have a steady increase of MRT (Table 5). In the case of [Ln(Z2-CO3)4]5 complexes the situation is much simpler: only one water molecule is able to enter the first shell (Table 6). This was only observed for La(III), Ce(III) and, with very low probability, Pr(III), while for the rest of the series the carbonates fulfill the coordination sphere of Ln(III) and there is no room for water molecules to enter. In Fig. 10 we report some prototypical structures (snapshots from the dynamics) of [Ln(Z2-CO3)3]3 and [Ln(Z2-CO3)4]5 clusters to show how such water molecules are spatially arranged in the first shell in the presence of three or four carbonates. We should note that they are in-between carbonate molecules, similarly to what was obtained from PBE/COSMO optimizations, but the hydrogen atoms now point towards bulk while it was not always the case in cluster optimization. This clearly comes from the explicit inclusion of solvent water molecules outside the cluster.

Fig. 10 Snapshots of Ln–carbonate structures obtained from POL-ff simulations in liquid water. From left to right and top to bottom: [Ln(Z2-CO3)3]3 with 4, 3 and 2 water molecules in the first shell and [Ln(Z2-CO3)4]5 with one water molecule in the first shell.

This journal is © the Owner Societies 2014

Paper

4. Conclusions In this paper we have studied the complexes formed by carbonates with ions of the lanthanoids(III) series by means of polarizable molecular dynamics simulations, DFT-based molecular dynamics simulations and DFT optimizations with explicit water molecules in the first shell and continuum method (COSMO) to mimic the bulk effect. With this aim we have first built a polarizable potential able to correctly reproduce DFTbased structure and energetics. Starting from parameters adjusted on La(III), the first atom of the series, we were then able to extrapolate parameters to the whole series without the need of electronic structure calculations that are slow and problematic for open shell Ln(III) of the series. No classical force field, with or without polarization, between lanthanoid(III) ions and carbonates is reported in the literature, and thus the present work can serve as a reference for future work on the application of classical molecular dynamics on systems containing such species under different conditions of concentration and ionic strength. Here, we studied complexes formed with three and four carbonates. When considering these complexes, different coordination motifs are possible: we found that bi-dentate structures are spontaneously formed when starting from monodentate ones, thus confirming the picture obtained from previous PCM calculations.17 Further, the tri-carbonate complexes seem to be preferred from an enthalpy point of view in the limit of infinite dilution conditions. Of course this result holds under ideal conditions, and, for example, the high ionic strength of solutions needed to form these complexes can modify the equilibrium. Thanks to the possibility of performing long simulations with explicit water molecules provided by the developed force field and subsequent classical MD simulations, we were able to study whether and how water molecules enter into the first shell of the Ln(III) ion in the presence of three or four carbonates surrounding the ion. We obtained that for tri-carbonate complexes there is room for up to four water molecules for light Ln(III) ions, while this decreases to two for heavy ions. In particular, MD simulations provided details of water exchange dynamics and how the first shell water molecules residence time evolves across the series, which could help in interpreting experimental data like TRLFS data. On the other hand, for tetracarbonate complexes it is very difficult for a water molecule to enter into the ion’s first shell, and this possibility is only allowed for light ions, even with relatively low probability (less than 20%). Finally, the ability of the extrapolation procedure of parameters for the whole series starting from those for La(III) to reproduce DFT-based structure and energetics for Lu(III) and TRLFS of Eu(III) and Dy(III) paves the way for using such approaches for other ligands of interest like silicates or nitrates.

Acknowledgements ` for useful We thank Prof. Laura Gagliardi and Dr Pere Miro discussions and Mrs Elizabeth A. Kish for careful reading of the

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3703

View Article Online

Paper

manuscript. This work was supported by the French National Research Agency (ANR) on project ACLASOLV (ANR-10-JCJC0807-01) (F. M., Y. J, T. V and R. S.). We acknowledge GENCI (grants x2012071870 and x2013071870) for computing time.

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

References 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

19 20 21 22 23 24 25

26

P. D’Angelo and R. Spezia, Chem.–Eur. J., 2012, 18, 11162. L. Helm and A. E. Merbach, Chem. Rev., 2005, 105, 1923. H. Ohtaki and T. Radnai, Chem. Rev., 1993, 93, 1157. Y. Marcus, Chem. Rev., 1988, 88, 1475. R. Pollet and D. Marx, J. Chem. Phys., 2007, 126, 181102. R. Pollet, N. N. Nair and D. Marx, Inorg. Chem., 2011, 50, 4791. O. V. Yazyev and L. Helm, Eur. J. Inorg. Chem., 2008, 201. N. Kaltsoyannis and P. Scott, The f Elements, 2007, Oxford Chemistry Primers Oxford University Press, Oxford. H. S. Sherry and J. A. Marinsky, Inorg. Chem., 1963, 3, 330. M. S. Wickleder, Chem. Rev., 2002, 102, 2011. W. Runde, M. P. Neu, C. Van Pelt and B. L. Scott, Inorg. Chem., 2000, 39, 1050–1051. T. Vercouter, P. Vitorge, N. Trigoulet, E. Giffaut and C. Moulin, New J. Chem., 2005, 29, 544. V. Philippini, T. Vercouter, J. Aupiais, S. Topin, C. Ambard, ´ and P. Vitorge, Electrophoresis, 2008, 29, 2041–2050. A. Chausse V. Philippini, T. Vercouter and P. Vitorge, J. Solution Chem., 2010, 39, 747. R. Janicki, P. Starynowicz and A. Mondry, Eur. J. Inorg. Chem., 2011, 3601. S. P. Sinha, A. M. Simas and G. L. C. Moura, J. Rare Earths, 2010, 28, 847. Y. Jeanvoine, P. Miro, F. Martelli, C. J. Cramer and R. Spezia, Phys. Chem. Chem. Phys., 2012, 14, 14822. P. D’Angelo, A. Zitolo, V. Migliorati, G. Chillemi, M. Duvail, P. Vitorge, S. Abadie and R. Spezia, Inorg. Chem., 2011, 50, 4572. C. Beuchat, D. Hagberg, R. Spezia and L. Gagliardi, J. Phys. Chem. B, 2010, 114, 15590. ´ra, F. Calvo and J.-P. Dognon, J. Chem. Phys., C. Clavague 2006, 124, 074505. ´ra, R. Pollet, J. M. Soudan, V. Brenner and J.-P. C. Clavague Dognon, J. Phys. Chem. B, 2005, 109, 7614. T. Kowall, F. Foglia, L. Helm and A. E. Merbach, J. Am. Chem. Soc., 1995, 117, 3790. T. Kowall, F. Foglia, L. Helm and A. E. Merbach, J. Phys. Chem., 1995, 99, 13078. M. Duvail, R. Spezia and P. Vitorge, ChemPhysChem, 2008, 9, 693. C. Apostolidis, B. Schimmelpfennig, N. Magnani, P. LidqvistReis, O. Walter, R. Sykora, A. Morgenstern, R. C. Colineau, R. Klenze, R. G. Haire, J. Rebizant, F. Bruchertseifer and T. Fanghanel, Angew. Chem., Int. Ed., 2010, 49, 6343. G. Dupouy, I. Bonhoure, S. D. Conradson, T. Dumas, C. Hennig, C. Le Naour, P. Moisy, S. Petit, A. Scheinost, E. Simoni and C. Den Auwer, Eur. J. Inorg. Chem., 2011, 1560.

3704 | Phys. Chem. Chem. Phys., 2014, 16, 3693--3705

PCCP

27 C. Terrier, P. Vitorge, M.-P. Gaigeot, R. Spezia and R. Vuilleumier, J. Chem. Phys., 2010, 133, 044509. 28 M. Duvail and P. Guilbaud, Phys. Chem. Chem. Phys., 2011, 13, 5840. 29 M. Duvail, A. Ruas, L. Venault, P. Moisy and P. Guilbaud, Inorg. Chem., 2010, 49, 519. 30 M. Duvail, P. Vitorge and R. Spezia, J. Chem. Phys., 2009, 130, 104501. 31 M. Duvail, F. Martelli, P. Vitorge and R. Spezia, J. Chem. Phys., 2011, 135, 044503. 32 B. T. Thole, Chem. Phys., 1981, 59, 341. 33 M. Duvail, M. Souaille, R. Spezia, T. Cartailler and P. Vitorge, J. Chem. Phys., 2007, 127, 034504. 34 S. J. Weiner, P. A. Kollman, D. A. Case, C. Singh, C. Ghio, G. Alagona, S. Profeta Jr and P. A. Weine, J. Am. Chem. Soc., 1984, 106, 765. 35 W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, J. Chem. Phys., 1983, 79, 926. 36 P. D’Angelo, F. Martelli, R. Spezia, A. Filipponi and M. A. Denecke, Inorg. Chem., 2013, 52, 10318. 37 F. Hutchinson, M. Wilson and P. A. Madden, Mol. Phys., 2001, 99, 811. 38 M. Salanne, C. Simon, P. Turq and P. A. Madden, J. Phys. Chem. B, 2008, 112, 1177. 39 Y. Okamoto, S. Suzuki, H. Shiwaku, A. Ikeda-Ohno, T. Yaita and P. A. Madden, J. Phys. Chem. A, 2010, 114, 4664. 40 P. van Duijnen and M. Swart, J. Phys. Chem. A, 1998, 102, 2399. 41 Handbook of Physics and Chemistry, 1996, CRC, Boca Raton, FL. 42 L. Gagliardi, R. Lindh and G. Karlstrom, J. Chem. Phys., 2004, 121, 4494. 43 B. O. Roos, R. Lindh, P.-Å. Malmqvist, V. Veryazov and P.-O. Widmark, J. Phys. Chem. A, 2004, 108, 2851. ¨m, R. Lindh, P.-A. Malmqvist, B.-O. Roos, U. Ryde, 44 G. Karlstro V. Veryazov, P.-O. Widmark, M. Cossi, B. Schimmelpfennig, P. Neogrady and L. Seijo, Comput. Mater. Sci., 2003, 287, 222. 45 F. Bruneval, D. Donadio and M. Parrinello, J. Phys. Chem. B, 2007, 111, 12219. 46 V. Vchirawongkwin, C. Kritayakornupong, A. Tongraar and B. M. Rode, J. Phys. Chem. B, 2011, 115, 12527. 47 W. Smith, T. R. Forester and I. T. Todorov, DL_POLY Classic, Version 1.0, STFC Daresbury Laboratory Daresbury, Warrington WA4 4AD Cheshire, UK, 2010. 48 H. J. C. Berendsen, J. P. M. Postma, W. van Gunsteren, A. Di Nola and J. R. Haak, J. Chem. Phys., 1984, 81, 3684. 49 U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577. 50 J.-P. Ryckaert, G. Ciccotti and H. J. C. Berendsen, J. Comput. Phys., 1977, 23, 327. 51 M. Sprik, J. Phys. Chem., 1991, 95, 2283. 52 M. Souaille, H. Loirat, D. Borgis and M.-P. Gaigeot, Comput. Phys. Commun., 2009, 180, 276. 53 F. Martelli, S. Abadie, J.-P. Simonin, R. Vuilleumier and R. Spezia, Pure Appl. Chem., 2013, 85, 237. 54 F. Martelli, R. Vuilleumier, J.-P. Simonin and R. Spezia, J. Chem. Phys., 2012, 137, 164501. 55 R. Car and M. Parrinello, Phys. Rev. Lett., 1985, 55, 2471.

This journal is © the Owner Societies 2014

View Article Online

Published on 13 December 2013. Downloaded by University of Windsor on 28/10/2014 15:54:36.

PCCP

56 A. D. Becke, Phys. Rev. A, 1998, 38, 3098. 57 C. Lee, W. Yang and R. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785. 58 N. Troullier and J. L. Martins, Phys. Rev. B: Condens. Matter Mater. Phys., 1991, 43, 1993. 59 R. Spezia, G. Tournois, J. Tortajada, T. Cartailler and M.-P. Gaigeot, Phys. Chem. Chem. Phys., 2006, 8, 2040. 60 R. Spezia, C. Bresson, C. Den Auwer and M.-P. Gaigeot, J. Phys. Chem. B, 2008, 112, 6490. 61 R. Ayala, R. Spezia, R. Vuilleumier, J. M. Martinez, R. R. Pappalardo and E. Sanchez Marcos, J. Phys. Chem. B, 2010, 114, 12866. 62 L. Kleinman and D. M. Bylander, Phys. Rev. Lett., 1982, 48, 1425. 63 CPMD, http://www.cpmd.org/, Copyright IBM Corp. 1990-2008, ¨rperforschung Stuttgard 1997-2001. ¨r Festko Copyright MPI fu 64 J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1995, 77, 3865. 65 J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1997, 78, 1396. 66 E. van Lenthe, E. J. Baerends and J. G. Snijders, J. Chem. Phys., 1994, 101, 9783. 67 E. van Lenthe, E. J. Baerends and J. G. Snijders, J. Chem. Phys., 1993, 99, 4597–4610. 68 A. Klamt and G. Schuurmann, J. Chem. Soc., Perkin Trans. 2, 1993, 799.

This journal is © the Owner Societies 2014

Paper

69 A. V. Marenich, C. J. Cramer and D. G. Truhlar, J. Phys. Chem. B, 2009, 113, 6378. 70 S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104. 71 G. T. Velde, F. M. Bickelhaupt, E. J. Baerends, C. F. Guerra, S. J. A. Van Gisbergen, J. G. Snijders and T. Ziegler, J. Comput. Chem., 2001, 22, 931. 72 E. J. M. Leenders, P. G. Bolhuis and E. J. Meijer, J. Chem. Theory Comput., 2008, 4, 898. 73 R. Spezia, M. Duvail, P. Vitorge, T. Cartailler, J. Tortajada, G. Chillemi, P. D’Angelo and M.-P. Gaigeot, J. Phys. Chem. A, 2006, 110, 13081. 74 C. Bresson, R. Spezia, S. Esnouf, S. Coantic, P. L. Solari and C. Den Auwer, New J. Chem., 2007, 31, 1789. 75 S. Izvekov and G. A. Voth, J. Chem. Phys., 2005, 123, 044505. 76 G. Hummer, L. R. Pratt and A. E. Garcia, J. Phys. Chem., 1996, 100, 1206. 77 P. H. Hunenberger and J. A. McCammon, J. Chem. Phys., 1999, 110, 1856. 78 M. A. Kastenholz and P. H. Hunenberger, J. Chem. Phys., 2006, 124, 124106. 79 M. A. Kastenholz and P. H. Hunenberger, J. Chem. Phys., 2006, 124, 224501. 80 R. W. Impey, P. A. Madden and I. R. McDonald, J. Phys. Chem., 1983, 87, 507.

Phys. Chem. Chem. Phys., 2014, 16, 3693--3705 | 3705

Hydration properties of lanthanoid(III) carbonate complexes in liquid water determined by polarizable molecular dynamics simulations.

In this work we have studied the structure and dynamics of complexes formed by three and four carbonates and a central lanthanoid(III) ion in liquid w...
2MB Sizes 0 Downloads 0 Views