Hydration of trivalent lanthanum revisited - An ab initio QMCF-MD approach.

The previously investigated La3+-hydrate has been re-evaluated by means of the quantum mechanical charge field (QMCF) molecular dynamics (MD) approach. Improved description of the hydration characteristics has been realised by including the full second hydration shell into the quantum mechanically treated region and by introducing the influence of the surrounding bulk via an electrostatic embedding technique. Analytical tools such as the ligand angular radial distribution analysis have been employed to gain deeper insight into the structural features of the hydrate. La3+ simultaneously forms nona- and decahydrates with capped trigonal and quadratic prismatic structure, besides small amounts of an octahydrate.

Introduction

Ionic lanthanoids with their highly coordinating character are currently of particular interest as understanding of detailed structural and dynamical properties is the key for a wide variety of applications with increasing economical interest. Their behaviour in aqueous environment is notably important in medical diagnostics as luminescent probe or as contrast agent in magnetic resonance imaging (MRI) [1]. Other main fields of interest are the development of novel catalysts in organic synthesis [2] and the application in nuclear waste repositories. Lanthanum oxides are of high importance in the optical industry [3].

Characterisation of hydrated trivalent lanthanum has been carried out by extended X-ray absorption fine structure (EXAFS) experiments [4-6] and via mid-infrared/Raman spectroscopy [7]. However, as the capturing of ultrafast dynamical inter- and intramolecular processes and an in-detail structural description are still not feasible with modern experimental techniques, the application of theoretical simulations for this task is promising. Available computational investigations of the hydrate are based on QM/MM-MD [8], MM-MD [9-11] and Car-Parinello-MD [12] approaches. As previous applications of the QMCF-MD framework [13-16] attest strongly improved descriptions of the studied systems over conventional QM/MM-MD based studies of hydrated ions [17-19], the hydration of La3+ was investigated with this new ab initio approach as well.

Methods

Simulation method and setup

Similar as the QM/MM-MD approach [20,21], the QMCF-MD ansatz [13-16] is also based on a partitioning of the simulation box into quantum mechanically and molecular mechanically treated regions. The QM region is further split into two sub-regions, in this case resembling the chemically relevant first and second hydration spheres. The difficult and time consuming construction of solute–solvent potentials is no longer a necessity with this approach, as the radius of the quantum mechanically treated region usually is 50–100% larger than in the case of a simple QM/MM-MD study. Detailed information on the employed framework is given in the corresponding works [13-16], describing the setup of interaction forces, the employment of calculation formalisms, the realisation of smooth particle migration between QM and MM region and also the versatile applicability of this approach.

Prior to the the simulation, La3+–(H2O) (n = 1–9) geometry optimisations were conducted with the Gaussian09 software package [22] in order to confirm the quality of the chosen LANL2DZ effective core potential (ECP) basis-set for lanthanum provided by Hay and Wadt [23]. The quality of the employed calculation method (Hartree–Fock) was also proven along the mentioned optimisations, yielding only minimal electron correlation contributions when compared to MP/2, CC-SD and B3LYP based calculations. The resulting equilibrium geometries of the ion–water clusters (see Supporting information) are in agreement with the reported experimental [5,6] and theoretical results [8-12], thus underlining the quality of the chosen method and the basis set further. For water, the well established and widely used Dunning double-ζ basis set including polarising functions for oxygen and hydrogen was employed [24-26]. The BJH-CF2 water model [27,28], enabling explicit hydrogen movement, was utilised for the MM-treated water molecules.

Simulation protocol

Compared to the previously conducted QM/MM-MD based study, the simulation box was enlarged substantially (side length of 31.15 Å), now containing 1000 water molecules, however resembling the same density of 0.997 g/cm3. The radii of the QM core and layer zones were set to 3.5 and 6.2 Å, respectively, already including a 0.2 Å thick smoothing region. The canonical ensemble (NVT) was chosen for the simulation and the constant temperature of 298 K was maintained by the Berendsen weak-coupling algorithm with a 0.1 ps relaxation time [29]. The long-range Coulombic cut-off at 15 Å was supplemented by the reaction field approach [30] with a permittivity of , and the integration of the equations of motion was realised with an Adams–Bashforth predictor–corrector algorithm. In order to create a bulk-like environment and thus avoiding surface conditions, periodic boundary condition and minimum image convention were employed. After an initial equilibration period of 2 ps, the hydrate was heated to 600 K and subsequently cooled and re-equilibrated for 5 ps. The actual sampling was conducted for 40 ps and the simulation time step was chosen as 0.2 fs. The utilised 16 CPU core platform with 32 GB of installed RAM allowed for a loop time of roughly 120 s per simulation step, resulting in approximately 1 year of calculation time for the simulation. For calculation of the quantum mechanical part of the simulation, the Turbomole 5.10 software package [31] was used.

Results and discussion

La3+ in aqueous solution proved to be much more dynamic than suggested by the previous QM/MM-MD study [8] – numerous ligand exchange reactions were observed between the first and the second hydration sphere; however, during the heating period, where the solution was heated to 600 K, no hydrolysis reactions occurred.

Figure 1 shows the respective lanthanum oxygen and lanthanum hydrogen pair distribution functions and their corresponding integrations. The hydration spheres are similarly well defined as previously reported [8], but the non-zero minimum between the first and the second hydration sphere is higher than in the QM/MM-MD based simulation, indicating more ligand exchange reactions. The region between 6.2 and 8.5 Å indicates a weak third hydration sphere and the plots [32-34] in Figure 2 confirm the existence of this ill-defined sphere, being in agreement with the previous simulation. The two innermost hydration spheres range from 2.3 to 3.3 Å and 3.9 to 6.1 Å with their maxima being located at 2.61 and 4.70 Å, respectively. Table 1 shows the structural parameters obtained from other approaches and compares them to the QMCF-MD related results, attesting a correct trend towards experimental values. The deviation from the experiments can be explained by a variation of the hydration due to the presence of eight-, nine- and tenfold coordination, which leads to an average first shell cordination number of 9.5. This trend also agrees with a recent EXAFS experiment [6], reporting a coordination number of 9.2, indicating that nine- and tenfold coordination occurs in solution. Although the first shell ion–oxygen distance reported from a MM-MD simulation [10] is in better agreement with experimental data, the resulting coordination number of 8.9 shows a different trend. The QMCF-MD obtained shell boundaries are very similar to the ones retrieved from the QM/MM-MD simulation [8] and the slightly shifted second shell maximum can be attributed to the larger QM-treated region (6.2 vs. previously 4.2 Å) in the QMCF-MD approach, which implies that all hydrogen bonds between these shells were treated quantum mechanically.

Figure 1

La3+–O and La3+–H RDFs and their running integrations.

Figure 2

Local density corrected three-body distribution functions for the first, second and third shell of hydration. Overlay of the O–O pair distribution function for pure water (dashed line) is given for comparison [44].

Table 1

Hydration sphere radii (rmin,, rmax,) and maximum shell peaks of the La3+-hydrate in Å – all values refer to oxygen.

Study	rmin,1	rmax,1	rmin,2	rmax,2	r¯_1	r¯_2
QMCF-MD	2.3	3.3	3.9	6.1	2.61	4.70
QM/MM-MD [8]	2.4	3.2	4.1	5.7	2.65	5.0
MM-MD [9]	2.0	3.0	3.8	5.4	2.66	5.0
MM-MD [10]	2.2	3.1	3.6	5.8	2.56	4.68
MM-MD [11]	n.a.	n.a.	n.a.	n.a.	2.52	4.65
CP-MD [12]	2.3	3.0	n.a.	n.a.	2.53	n.a.
EXAFS [4]	n.a.	n.a.	n.a.	n.a.	2.58	n.a.
EXAFS [5]	n.a.	n.a.	n.a.	n.a.	2.54	n.a.
EXAFS [6]	n.a.	n.a.	n.a.	n.a.	2.54	n.a.

The dynamic yet well defined first hydration sphere is characterised by ninefold (43% relative occurence) and tenfold (56% relative occurrence) coordination, an octahydrate is rarely observed. The coordination number distribution of the second hydration shell is not as sharply defined as in the QM/MM-MD simulation [8] whose mean coordination number of 23.4 water molecules is slightly smaller than the coordination number observed in the QMCF-MD study (25.6). The mean first shell coordination number of 9.5 is almost identical with the QM/MM-MD simulation, which however lacks the rare events of an eightfold coordination (1% relative occurrence). Table 2 compares the obtained data with other theoretical studies [8-12] and the reported X-ray absorption measurements [4-6]. The out-dated twelvefold coordination from the 1995 EXAFS experiment reported by Solera et al. [4] has been corrected by the works of Helm and Merbach [35-37], Persson et al. [5] and Allen et al. [6]. Today, it is generally accepted that lighter lanthanoid ions are predominantly ninefold coordinated while the higher species from gadolinium onwards are mostly eightfold coordinated.

Table 2

Minimum (CNmin,) and maximum (CNmax,) coordination numbers and mean coordination numbers (CN) of the La3+-hydrate.

Study	CNmin,1	CNmax,1	CNmin,2	CNmax,2	CN1	CN2
QMCF-MD	8	10	19	32	9.5	25.6
QM/MM-MD [8]	9	10	18	28	9.6	23.4
MM-MD [9]	8	10	n.a.	n.a.	9.1	19.1
MM-MD [10]	8	10	n.a.	n.a.	8.9	15.9
MM-MD [11]	9	10	n.a.	n.a.	9.02	18.8
CP-MD [12]	8	9	n.a.	n.a.	n.a.	n.a.
EXAFS [4]	n.a.	n.a.	n.a.	n.a.	12	n.a.
EXAFS [5]	n.a.	n.a.	n.a.	n.a.	9	n.a.
EXAFS [6]	n.a.	n.a.	n.a.	n.a.	9.2	n.a.

The O–La–O angular distribution function is depicted in Figure 3, yielding maxima in the regions of 70.5° and 136.5° – the valley between the two peaks ranges from 84° to 121°. The narrowness of the peaks reflects a strong solute–solvent interaction and the valley indicates flexibility of the hydrate resulting from structural re-organisation. While the ninefold coordinated species can form either a tri-capped trigonal prism (Figure 4, top) or a capped square antiprism, the tenfold coordinated complex mostly resembles a decagon with a capped pentagonal planar substructure and a square arrangement of the remaining four ligands below the pentagonal plane (Figure 4, bottom). The eightfold coordinated complex, only being stable for 100–200 fs, is characterised by a more or less well defined square antiprism surrounding the solute atom. Analysis of the angular radial distribution (ARD) of ligands [38,39] is shown in Figure 4, underlining the previous statements. For the ninefold coordinate complex, the diffuse first shell ligand localisation observed from the ARD [38,39] in Figure 4 confirms the previously reported [11,12] easily accessible interconversion between capped square antiprismatic and tri-capped trigonal prismatic structure. The single cap water molecule in the decagonal structure residing in the axial position is bound in a more rigid manner than it is the case for the remaining ligands (Figure 4, bottom). The concluded lack of axial ligand density until the third hydration sphere in the case of the Ce4+-hydrate [40] does not hold true for the isoelectronic La3+-hydrate. Significant preferential location of second shell water molecules in the axial region becomes obvious in both the nine- and the tenfold coordinated lanthanum hydrate, in contrast to the observations made in the hydration structure of Ce4+ [40].

Figure 3

O–La3+–O ADF evaluated for the first hydration sphere.

Figure 4

La3+–O ARD plots drawn over 3 ps each.

The mean ligand residence times have been computed and compared to the previously obtained values and to those retrieved for other highly charged metal ions in Table 3. The results suggest that the QMCF-MD approach seems to describe hydrated La3+ in a more dynamic system than the QM/MM-MD framework. For the previous simulation of La3+, first shell water exchange reactions were rarely observed (0.4 vs. now 5.8 exchanges per 10 ps), indicating a less rigid first hydration sphere than previously assumed [8]. Four-hundred and forty one ligand displacement reactions with at least τ ⩾ 0.5 ps and 1016 attempted exchange reactions from the second hydration sphere were observed during the 40 ps of sampling, further underlining the more dynamic character of the hydrate. The REX,x values denote the required number of ligand exchange attempts until a successful event is registered, allowing direct comparison of various highly charged ions.

Table 3

First shell and second shell MRTs and corresponding REX, values obtained for various polarising ions in aqueous solution.

Ion	τ10.5	τ20.5	REX,1	REX,2
	(ps)	(ps)
AlQMCF3+[17]	n.a.	17.7	n.a.	21. 1
ZrQMCF4+[41]	n.a.	5.5	n.a.	6.8
LaQM/MM3+[8]	⩾240	8.4	6.0	13.2
Lathis work3+	16.6	2.3	2.4	2.3
CeQMCF4+[40]	n.a.	6.0	n.a.	8.5
HfQMCF4+[42]	n.a.	15.4	n.a.	10.6
UQMCF4+[43]	n.a.	8.1	n.a.	5.9
H_2O_QM/MM[44]	1.7	n.a.	11.2	n.a.

The well equilibrated, 40 ps long, simulation trajectory represents an ideal basis for the calculation of the lanthanum–oxygen stretching motion and its corresponding force constant, enabling a direct comparison with the QM/MM-MD simulation and other previously simulated ions (Table 4). As vibrational motions can be excited with infrared light, comparison with experimentally derived results is possible. Trivalent lanthanum, despite its much more dynamic character than previously suggested [8], belongs to a series of highly charged ionic species with notable structure forming capabilities. The reported wave number of 360 cm−1 for the ion–oxygen stretching frequency is in excellent agreement with the value reported by Kanno and Hiraishi (354 cm−1) [7], once more showing the improved dynamical description of the QMCF-MD approach. The force constant of the ion–oxygen bond is a good criterion for the relative binding strength of water ligands to the metal ion. Due to the dependence of the ion–oxygen stretching motion’s frequency on the second derivative of the energy with respect to the nuclear coordinates, it serves as a sensitive probe of the accuracy of the simulation. The value obtained for the vibrational frequency highlights the improvement achieved when the QMCF-MD methodology is applied, especially when comparing the obtained value (360 cm−1) to the previously conducted QM/MM-MD simulation (253 cm−1) [8].

Table 4

Comparison of theoretically derived stretching frequencies (top) with experimental results (bottom).

Ion	QIon-O	kIon-O
	(cm−1)	(Nm−1)
Al3+a[17]	560	194
Fe2+[45]	357	93
Fe3+[45]	513	193
Zr4+[41]	484	188
Hf4+[42]	n.a.	212
Ir3+[40]	n.a.	260
Lathiswork3+	360	110
Ce3+[40]	354	106
Ce4+[40]	420	149
La3+[46,47]	354	106
Ce3+[46]	359	109
Ce4+[46,48]	408	141

Values have been unscaled prior to comparison [49,50].

Conclusion

This Letter presents a revised ab initio treatment of trivalent lanthanum in aqueous environment. Due to the significant sampling period and extended QM region, valuable insight into structural and dynamical features was gained. With the ab initio investigation of highly charged ions being a rather delicate task, the quality of the employed QMCF-MD ansatz [13-16] was proven also by comparing the obtained results with the ones from a previously conducted QM/MM-MD simulation of hydrated La3+ [8]. The charge flexibility in the QM region and the influence of QM particles’ partial charges on MM particles contribute their share to the improved quality of the QMCF-MD ansatz. Minor structural improvements and major improvements with regard to the dynamics of the system are ascribed to the quantum mechanically treated second hydration shell, facilitating an accurate description of the inter-shell hydrogen bonds. Comparison of the theoretically derived results with experimental data delivered very satisfactory results. The employment of analytical tools like the evaluation of ligand angular radial distributions proved as valuable for the identification of distinct structural motifs. The local density corrected three-body functions have identified a weak third hydration sphere, indicating similarities in hydration structure of highly charged lanthanoid species. Comparison of the lanthanum–oxygen stretching motion’s wave number yielded excellent agreement with experimental results, which was also observed for structural properties such as ion–oxygen bond lengths and mean coordination numbers.