- Email: [email protected]
Thulium(III) and ytterbium(III) in aqueous solution ab initio quantum mechanical charge field molecular dynamics studies
Chemical Physics Letters 638 (2015) 128–132
Contents lists available at ScienceDirect
Chemical Physics Letters journal homepage: www.elsevier.com/locate/cplett
Thulium(III) and ytterbium(III) in aqueous solution ab initio quantum mechanical charge field molecular dynamics studies Peter P. Passler, Bernd M. Rode ∗ Theoretical Chemistry Division, Institute of General, Inorganic and Theoretical Chemistry, University of Innsbruck, Innrain 80-82, A-6020 Innsbruck, Austria
a r t i c l e
i n f o
Article history: Received 10 July 2015 In final form 13 August 2015 Available online 20 August 2015
a b s t r a c t Hydration properties of trivalent thulium and ytterbium ions in aqueous solution are investigated via quantum mechanical charge field molecular dynamics (QMCF-MD) simulations. The QMCF-MD formalism is a special type of QM/MM simulation, where the chemically most relevant part of the system – in this case the ion with its first and second hydration shells – is treated by quantum mechanics. The mean ion O distances and the average coordination numbers of the first hydration shells are compared with experimental EXAFS data. Mean ligand residence times, vibrational frequencies and force constants were evaluated to characterise the dynamics of the systems. © 2015 Elsevier B.V. All rights reserved.
1. Introduction The importance of lanthanoids and their ions is steadily increasing because of their manifold use in materials and industrial processes. Thulium and ytterbium are used in lasers, optical spectroscopy [1–7] and in various medical diagnostic methods such as X-ray, magnetic resonance imaging (MRI) and positron emission tomography (PET) technology [8–15]. Production and disposal of such apparatus often involve the presence of the trivalent ions of these elements in aqueous medium and hence detailed knowledge of the chemical and physical behavior of these ions has gained much attention stimulating experimental studies by extended Xray absorption fine structure (EXAFS) [16] as well as theoretical investigations [17]. However, further investigations have revealed that a higher level of theory is required to deal with these ions in water [18], and that the experiments might need some further analysis and interpretation. A particular property of the trivalent lanthanoid ions in water is their extremely rapid dynamics compared to trivalent transition metal and main group ions [18], leading to mean ligand residence times in the order of 15 (Ce) to 150 (Lu) ps and to very fast structural reorientation processes even without change of coordination number. Previous investigations of some lanthanoid ions, namely the trivalent ions of La [19], Ce [20], Gd [21], Tb [21], Dy [22], Ho [22], Er [23] and Lu [24,25], have assured that the QMCFMD methodology which includes two full hydration shells in the
∗ Corresponding author. E-mail address: [email protected] (B.M. Rode). http://dx.doi.org/10.1016/j.cplett.2015.08.041 0009-2614/© 2015 Elsevier B.V. All rights reserved.
quantum mechanical description of the system and a number of other improvements (vide infra) is a suitable tool for the study of such ions. This methodology has been employed, therefore, also for the present investigation of Tm(III) and Yb(III). The scope of this study was to reveal differences and similarities between these neighboring elements and to compare them on the other hand with their other neighbors Er(III) [23] and Lu(III) [24,25], respectively. 2. Methods In the QMCF molecular dynamics the quantum mechanically (QM) calculated part of the system, the QM region, consists of the solute and two complete hydration shells, with the advantages of describing all hydrogen bonds between 1st and 2nd shell by ab initio quantum mechanics and by shifting the smoothing zone between QM and molecular mechanics (MM) region quite far from the central solute, thus reducing the possible errors in the smoothing procedure [26]. The main advantage, however, is that due to the large distance between solute and MM region no potential functions for solute–solvent interactions are needed, as at that distance all non-Coulombic interactions have converged to zero, and a potential function is only required for solvent–solvent interactions. This is of particular importance in the case of highly charged ions, for which the construction of ion–solvent potentials on the basis of quantum mechanical calculations is difficult due to artificial charge transfers in pair potential evaluations. In practice, the QM region consists of two sub-regions, the core zone with ion and first hydration shell included, and the layer region, containing only solvent molecules. Further advantages of the QMCF-MD formalism are the inclusion of the MM charges as perturbation into the core
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
129
Table 1 ˚ Average Ln3+ –O1stshell distance in A. Ion
Ln3+ –O
Method
Reference
Tm(III)
2.44 2.33 2.42 2.32
QMCF-MD EXAFS QMCF-MD EXAFS
This work Persson et al. [16] This work Persson et al. [16]
Yb(III)
Figure 1. Definitions of the QM regions in QMCF-MD simulations.
Hamiltonian of the QM calculation of energies and forces, and the use of continuously newly determined QM charges in the evaluation of Coulombic interactions between QM and MM region. It is evident that the size of the two QM zones leads to a much higher computational effort than a conventional QM/MM approach with only one hydration shell. The increase in accuracy, especially for the dynamics of the system, fully justifies the cost of it. Further details and a mathematical description of the formalism are available in the literature [27–29]. The radius of the QM core region was ˚ that of the total QM region (core plus layer) to 5.7 A, ˚ set to 3.3 A, ˚ 1). smoothing being applied within 5.5 and 5.7 A(Figure For the quantum mechanical calculation of energies and forces ab initio calculations at Hartree–Fock (HF) level with dunning double zeta-polarisation (DZP) basis sets [30] for water and Stuttgart effective core potential (ECP) valence basis sets [31] for the lanthanoid ions have been employed. All quantum mechanical calculations were performed by Turbomole 6.4 [32]. Comparison of the results for ion–water clusters of other ions [19–24] at HF level with those of correlated methods have shown that the influence of electron correlation in such systems is very minor. Some previous results of similar simulations [19–24] have indicated that the relativistic correction via the ECP might be insufficient to take into account the relativistic contraction of the inner shells thus leading to an increase of the ion O distances compared to experiment [16]. On the other hand, the employment of correlated methods and even more, a full relativistic treatment of the systems in the simulations would lead to unrealistic computation times, as even the HF level
requires computation times of over 8000 h per system. Water in the MM-region was treated by the flexible Bopp–Jancsó–Heinzinger Center Force 2 model (BJH-CF2) [33,34], as this model proved most compatible with the Mullikan charges of the QM calculations. The simulations were carried out for a NVT ensemble consisting of one ion and 1000 water molecules at the experimental density of pure water and T = 298.15 K, thermostatised by the Berendsen algorithm [35], with a time step of 0.2 fs starting from a pre-equilibrated configuration of a simulation of Lu(III) in water. Equilibration required 20 000 steps, i.e. 4 ps, and sampling was performed for further 180 000 steps (36 ps). Hence the total computation time for the simulations presented here amounted to more than 200 000 CPUhours on a high-performance computer cluster with 12 dedicated 3.0 GHz processors, making use of 5 GB RAM disk. The trajectories have been analyzed with regard to radial distribution functions (RDF), angular distribution functions (ADF) and coordination number distribution (CND) for structural data. For the dynamics ligand mean residence times (MRT) and vibrational spectra, in particular frequencies and force constants of the Ln(III)–O stretching vibration were computed. The determination of the MRT was made by the ‘direct’ method [36] with a minimum displacement time from a hydration shell of t* = 0.5 ps. For the evaluation of the vibrational spectra the velocity autocorrelation functions (VACF) were Fourier-transformed in analogy to previous works [19–24]. 3. Results and discussion 3.1. Structural data Table 1 displays the characteristic data of the Ln(III)–O and Ln(III)–H RDFs of both ions, which are very similar. The absolute ˚ values for the first maxima differ by ∼0.1 Afrom the EXAFS data, ˚ but the small shortening of ∼0.01 Agoing from Tm(III) to Yb(III) is
Figure 2. (a) Ion O (black) and ion H (red) RDFs and their integrations of Tm3+ (upper) and Yb3+ (lower). The zoomed first peak is shown in the inserts, illustrating the presence of different hydrate structures. (b) Three-body correlation functions for the hydration shells of Tm3+ (left) and Yb3+ (right).
130
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
Figure 3. CNDs for the first and second hydration spheres of the hydrated Tm3+ (upper) and Yb3+ ion (lower).
well reproduced. The EXAFS data have been fitted to a model with six shorter (prisma) and three longer (caps) distances, and the values given in Table 1 are the average of these lengths. However, the RDFs of both ions (Figure 2) show upon zooming in tiny shoulder peaks near the first maximum of the ion O RDF, which could be an indication of nonequivalent bond lengths and will be discussed later in connection with the vibrational spectra. The too long absolute distances are another indication of an insufficient relativistic correction as observed in previous studies [19–24]. A closer inspection of the trajectories reveals that even without ligand exchange, reorientation processes of the first hydration shell occur on the picosecond scale, preventing the assignment of a unique structure to the hydrate. This seems to be a common characteristic of Ln(III) ions [19–24]. The ADFs of the first coordination shell are very similar again and are most compatible with a squared antiprismatic, a bi- and a tricapped trigonal prismatic structure (see Figure 4). The average coordination numbers of this shell resulting from the RDFs’ integration are 8.3 for Tm(III) and 8.4 for Yb(III), and the CND diagrams in Figure 3 display the slightly higher amount of 9-coordinated Yb(III) ion. The inserted snapshots in this figure are typical examples of the 8- and the 9-coordinated species. The dominance of coordination number 8 is not in full agreement with the EXAFS results [16] of 8.8 and 8.7 for Tm(III) and Yb(III), respectively. At this point is has to be mentioned, however, that fitting of a model for systems with very fast exchange and reorientation mechanisms seems an extremely difficult task in scattering experiments, and other literature data of coordination numbers ranging from 8 to 12 [37] seem to confirm this statement. The aforementioned EXAFS study concludes that almost all lanthanoid ions are up to Er(III) exclusively and even beyond dominantly 9-coordinated with only exception of Lu(III) (CN = 8.4), while other studies have discussed the so-called ‘Gadolinium break’, meaning a change of different physicochemical properties, possibly associated with a change of the coordination number from 9 to 8 at this element [38–40]. QMCF-MD simulations of the other immediate neighbors of Tm and Yb, Er and Lu have determined the average coordination numbers of 8.4 [23] and 8.0 [24,25] for these ions. The second hydration shells of Tm(III) and Yb(III) contain in average 20.1 and 19.4 water molecules, both extending over a wide range (cf. Figure 3), with a maximum at 19 and 20 ligands, respectively. The MRTs of these ligands are 3.1 and 2.7 ps, which significantly differs from the corresponding bulk value of water, being 1.6 ps [41]. While the shape of the RDFs (cf. Figure 2) shows a
Figure 4. ADF for the first hydration spheres of the hydrated Tm3+ (black) and Yb3+ ion (red), compared with the ADFs (blue lines) of an ideal (a) squared antiprism (CN = 8), (b) bicapped trigonal prism (CN = 8) and (c) tricapped trigonal prism (CN = 9).
distinct second shell connected to the first one mainly by H bonds according to the orientation of ligands, a third shell is hardly recognized. In order to confirm this behavior, the three-body analysis [42,43] was performed. Figure 2b proves that the first and second shell show clear deviations from the bulk structure, whereas the third shell shows identical features as the bulk. Evaluation of the second shell ligands’ MRT with t* = 0.0 ps provides an estimate of the average lifetime of the H bonds between first and second shell (each removal of a ligand from the second shell means breaking a hydrogen bond of this ligand to one in the first shell). The resulting values are 0.76 ps for Yb(III) and 0.68 ps for Tm(III), respectively. Consequently, the stability of the H bonds between first and second shell of both ions is considerably higher than between bulk water molecules, for which QMCF-MD simulations produced a value of 0.36 ps [41]. 3.2. Dynamics The mean residence times of water molecules in the first shell are 37 and 25 ps for Tm(III) and Yb(III), both lower than the corresponding value of Er(III) (41 ps), but also much lower than the value for the f14 ion Lu(III) (>160 ps). Apparently the MRTs are not or only loosely correlated to ion–ligand distance and coordination number. On the other hand, it seems that MRTs of trivalent lanthanoid ions are by many orders of magnitude shorter than those of trivalent d-transition elements and main group elements [44], almost comparable to monovalent and some divalent ions. This will have a considerable influence on their complex chemistry, e.g. in the binding of proteins utilized in specific separation processes recently [45]. The distance plots in Figure 5 provide some insight into the ligand exchange reaction mechanism. Although only a few
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
131
Table 2 First and second shell average coordination numbers, MRTs in ps, ion O stretching frequencies and corresponding force constants of Tm3+ , Yb3+ and other trivalent ions, all obtained by QMCF-MD studies. Ion 3+
Al Fe3+ Ir3+ La3+ Ce3+ Gd3+ Tb3+ Dy3+ Ho3+ Er3+ Tm3+ Yb3+ Lu3+
CN1
MRT1
CN2
MRT2
Qion
6 6 6 9.5 9.1 8.5 8.4 8.1 8.1 8.4 8.3 8.4 8
n.o. n.o. n.o. 16.6 >118 32.6 44.5 36.0 22.0 40.5 27.1 25.1 >160
11.8 13.6 13.5 25.6 21.4 19.8 20.1 19.2 19.5 18.8 20.1 19.4 18.7
17.7 3.1 3.6 2.3 2.6 3.0 2.7 3.0 2.1 3.0 3.1 2.7 2.6
560 513 n.a. 360 354 340 336 334 334 360 356 352 360
O
(cm−1 )
kion
O
(Nm−1 )
194 193 260 110 106 99 97 96 96 110 109 107 112
Reference Hofer et al. [46] Moin et al. [47] Pedevilla et al. [48] Lutz et al. [19] Lutz et al. [20] Canaval et al. [21] Canaval et al. [21] Tirler et al. [22] Tirler et al. [22] Canaval et al. [23] This work This work Hitzenberger et al. [24,25]
n.o.: no exchange processes observed during simulation time; n.a.: not available.
Figure 5. Distance plot of all Tm3+ –O (upper) and Yb3+ –O (lower) distances.
exchange events have been observed, they all seem to follow dominantly the associative (A) and less the dissociative (D) type, the latter exclusively starting from a 9-coordinated structure. The distance plots in Figure 5 also reflect the vigorous dynamics outside the first hydration shell, but the reduced flexibility of the second
shell is still recognized by a slightly less populated region between ˚ second shell and bulk around 6 A. MRTs are related to the strength of the ion–ligand binding, which is reflected in the ion O stretching vibration frequency and force constant. These values have been evaluated, therefore, via the VACFs and subsequent Fourier transformation, and they are collected in Table 2, together with some characteristic values of transition metal and main group elements. Figure 6 displays the vibrational spectra obtained from the simulations, separated for both ions for 8- and 9- coordinated species. Similar to the RDFs the main peak is accompanied by a smaller one or shoulder at lower wave numbers indicating a fraction of weakly bound water ligands, in agreement with the EXAFS result of two different types of water ligands coordinated to the ion [16]. Table 2 also lists the frequencies and force constants for peaks and shoulders resulting for Tm(III) and Yb(III) and compares them to characteristic values for other trivalent ions. The observed weakness of the ion O binding is in full agreement with the short MRT values and the generally much higher lability of water ligands bound to lanthanoid ions. However, the differences between Tm(III) and Yb(III) are very small and hence within the methodical accuracy limits for the main peak. The MRT differences can hence be explained only by the differences in the secondary peaks. As they result from the weaker bound ligands, which should dominate the exchange processes, these findings appear consistent not only with the interpretation of the EXAFS studies predicting two types of ligands, but also with our RDFs
Figure 6. Vibrational spectra of Tm3+ –O (a) and Yb3+ –O (b). Black: total, red: CN = 8 only, blue: CN = 9 only.
132
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
and the MRT data: Yb(III) with a considerably shorter MRT shows the frequency of the secondary peak at a much more red-shifted position (285 cm−1 for Tm(III) and 228 cm−1 for Yb(III) cf. Figure 6). The dynamics of the second shell have been discussed already in the previous chapter. 4. Conclusions The present work has once more demonstrated the suitability of the ab initio QMCF-MD methodology for the simulation of ions in aqueous environment, even for the fastidious treatment of the lanthanoid ions, where subtle differences of very similar but chemically often different elements have to be worked out. The use of HF level and relativistically corrected ECP appear a reasonable although not fully satisfactory compromise with respect to computational demand. The ions treated here assume hydration numbers of 8 and 9 in the first shell, and vigorous exchange processes on the picosecond scale lead to rapid changes between structure of the aquo-complexes. The mean residence times of the water ligands are considerably shorter than for their neighbor Er(III), despite a very similar force constant of the ion O bond. The striking differences between Ln(III) ions and transition element and main group ions of the same charge are also observed for the examples studied here. Acknowledgement Generous supply of computer time and a scholarship for Peter P. Passler by the University of Innsbruck is gratefully acknowledged. References [1] W. Koechner, Solid-State Laser Engineering; Springer Series in Optical Sciences, Springer, 2006, pp. p.49. [2] F.J. Duarte, Tunable Laser Applications, CRC press, 2010. [3] P. Lacovara, H.K. Choi, C.A. Wang, R.L. Aggarwal, T.Y. Fan, Opt. Lett. 16 (1991) 1089. [4] J. Koponen, M. Söderlund, H. Hoffman, S. Tammela, Opt. Express 14 (2006) 11539. [5] J.-F. Bisson, D. Kouznetsov, K.-I. Ueda, S. Fredrich-Thornton, K. Petermann, G. Huber, Appl. Phys. Lett. 90 (2007) 201901. [6] N.V. Sochinskii, M. Abellán, J. Rodríguez-Fernández, E. Saucedo, C.M. Ruiz, V. Bermúdez, Appl. Phys. Lett. 91 (2007) 202112. [7] D.A. Grukh, V.A. Bogatyrev, A.A. Sysolyatin, V.M. Paramonov, A.S. Kurkov, E.M. Dianov, IEEE J. Quantum Electron. 34 (2004) 247. [8] N. Krishnamurthy, C. Gupta, Extractive Metallurgy of Rare Earths, Taylor & Francis, 2004, pp. 32. [9] D. Krishnamurthy, V. Weinberg, J.A.M. Cunha, I.-C. Hsu, J. Pouliot, Brachytherapy 10 (2011) 461.
[10] A. Ayoub, G. Shani, World Congress on Medical Physics and Biomedical Engineering, September 7–12, Munich, Germany, 2009, p. 1. [11] B. Raj, B. Venkataraman, Practical Radiography, Alpha Science International Limited, 2004. [12] P. Enghag, Encyclopedia of the Elements: Technical Data - History - Processing - Applications, Wiley, 2008, pp. 448. [13] R. Halmshaw, Industrial Radiology: Theory and Practice, Non-Destructive Evaluation Series Springer, 1995, pp. 168. [14] Y. Sun, M. Yu, S. Liang, Y. Zhang, C. Li, T. Mou, W. Yang, X. Zhang, B. Li, C. Huang, F. Li, Biomaterials 32 (2011) 2999. [15] S. Aime, A. Barge, D. Delli Castelli, F. Fedeli, A. Mortillaro, F.U. Nielsen, E. Terreno, Magn. Reson. Med. 47 (2002) 639. [16] I. Persson, P. D’Angelo, S. De Panfilis, M. Sandström, L. Eriksson, Chem.-Eur. J 14 (2008) 3056. [17] M. Duvail, R. Spezia, P. Vitorge, ChemPhysChem 9 (2008) 693. [18] T.S. Hofer, A.K. Weiss, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 512 (2011) 139. [19] O. Lutz, T.S. Hofer, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 536 (2012) 50. [20] O. Lutz, T.S. Hofer, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 539 (2012) 50. [21] L.R. Canaval, P.P. Passler, B.M. Rode, Chem. Phys. Lett. 625 (2015) 116. [22] A.O. Tirler, P.P. Passler, B.M. Rode, Chem. Phys. Lett. 635 (2015) 120. [23] L.R. Canaval, T. Sakwarathorn, B.M. Rode, C.B. Messner, O.M. Lutz, G.K. Bonn, J. Phys. Chem. B 117 (2013) 15151. [24] M. Hitzenberger, T.S. Hofer, A.K. Weiss, J. Chem. Phys. 139 (2013) 114306. [25] M. Hitzenberger, Solvation Properties and Behaviour of Lutetium (III) in Aqueous Solution–A Quantum Mechanical Charge Field (QMCF) Study, University of Innsbruck, 2012 (M.Sc. thesis). [26] H. Lin, D.G. Truhlar, Theor. Chem. Acc. 117 (2007) 185. [27] T.S. Hofer, A.B. Pribil, B.R. Randolf, B.M. Rode, in: J.R. Sabin, E. Brändas (Eds.), Combining Quantum Mechanics and Molecular Mechanics. Some Recent Progresses in QM/MM Methods, vol. 59, Academic Press, 2010, p. 213 (Chapter 7). [28] T.S. Hofer, B.M. Rode, A.B. Pribil, B.R. Randolf, in: R. van Eldik, J. Harvey (Eds.), Theoretical and Computational Inorganic Chemistry, vol. 62, Academic Press, 2010, p. 143. [29] B.M. Rode, T.S. Hofer, B.R. Randolf, C.F. Schwenk, D. Xenides, V. Vchirawongkwin, Theor. Chem. Acc. 115 (2006) 77. [30] T.H. Dunning Jr., J. Chem. Phys. 53 (1970) 2823. [31] X. Cao, M. Dolg, J. Chem. Phys. 115 (2001) 7348. [32] TURBOMOLE V6.4 2012, A Development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989–2007, TURBOMOLE GmbH, since 2007. Available from http://www.turbomole.com [33] P. Bopp, G. Jancsó, K. Heinzinger, Chem. Phys. Lett. 98 (1983) 129. [34] F. Stillinger, A. Rahman, J. Chem. Phys. 68 (1978) 666. [35] H.J. Berendsen, J.P.M. Postma, W.F. van Gunsteren, A. DiNola, J. Haak, J. Chem. Phys. 81 (1984) 3684. [36] T.S. Hofer, H.T. Tran, C.F. Schwenk, B.M. Rode, J. Comput. Chem. 25 (2004) 211. [37] C.-H. Huang, Rare Earth Coordination Chemistry: Fundamentals and Applications;, John Wiley & Sons, 2010. [38] F. Spedding, M. Pikal, B. Ayers, J. Phys. Chem. 70 (1966) 2440. [39] F. Spedding, K. Jones, J. Phys. Chem. 70 (1966) 2450. [40] F. Spedding, M. Pikal, J. Phys. Chem. 70 (1966) 2430. [41] D. Xenides, B.R. Randolf, B.M. Rode, J. Chem. Phys. 122 (2005) 174506. [42] A. Bhattacharjee, T.S. Hofer, B.M. Rode, Phys. Chem. Chem. Phys. 10 (2008) 6653. [43] A.K. Weiss, T.S. Hofer, RSC Adv. 3 (2013) 1606. [44] F.A. Dunand, L. Helm, A.E. Merbach, Adv. Inorg. Chem. 54 (2003) 1. [45] M.R. Mirza, M. Rainer, C.B. Messner, Y. Güzel, D. Schemeth, T. Stasyk, M.I. Choudhary, L.A. Huber, B.M. Rode, G.K. Bonn, Analyst 138 (2013) 2995. [46] T.S. Hofer, B.R. Randolf, B.M. Rode, J. Phys. Chem. B 112 (2008) 11726. [47] S.T. Moin, T.S. Hofer, A.B. Pribil, B.R. Randolf, B.M. Rode, Inorg. Chem. 49 (2010) 5101. [48] P.A. Pedevilla, T.S. Hofer, B.R. Randolf, B.M. Rode, Aust. J. Chem. 65 (2013) 1582.
Contents lists available at ScienceDirect
Chemical Physics Letters journal homepage: www.elsevier.com/locate/cplett
Thulium(III) and ytterbium(III) in aqueous solution ab initio quantum mechanical charge field molecular dynamics studies Peter P. Passler, Bernd M. Rode ∗ Theoretical Chemistry Division, Institute of General, Inorganic and Theoretical Chemistry, University of Innsbruck, Innrain 80-82, A-6020 Innsbruck, Austria
a r t i c l e
i n f o
Article history: Received 10 July 2015 In final form 13 August 2015 Available online 20 August 2015
a b s t r a c t Hydration properties of trivalent thulium and ytterbium ions in aqueous solution are investigated via quantum mechanical charge field molecular dynamics (QMCF-MD) simulations. The QMCF-MD formalism is a special type of QM/MM simulation, where the chemically most relevant part of the system – in this case the ion with its first and second hydration shells – is treated by quantum mechanics. The mean ion O distances and the average coordination numbers of the first hydration shells are compared with experimental EXAFS data. Mean ligand residence times, vibrational frequencies and force constants were evaluated to characterise the dynamics of the systems. © 2015 Elsevier B.V. All rights reserved.
1. Introduction The importance of lanthanoids and their ions is steadily increasing because of their manifold use in materials and industrial processes. Thulium and ytterbium are used in lasers, optical spectroscopy [1–7] and in various medical diagnostic methods such as X-ray, magnetic resonance imaging (MRI) and positron emission tomography (PET) technology [8–15]. Production and disposal of such apparatus often involve the presence of the trivalent ions of these elements in aqueous medium and hence detailed knowledge of the chemical and physical behavior of these ions has gained much attention stimulating experimental studies by extended Xray absorption fine structure (EXAFS) [16] as well as theoretical investigations [17]. However, further investigations have revealed that a higher level of theory is required to deal with these ions in water [18], and that the experiments might need some further analysis and interpretation. A particular property of the trivalent lanthanoid ions in water is their extremely rapid dynamics compared to trivalent transition metal and main group ions [18], leading to mean ligand residence times in the order of 15 (Ce) to 150 (Lu) ps and to very fast structural reorientation processes even without change of coordination number. Previous investigations of some lanthanoid ions, namely the trivalent ions of La [19], Ce [20], Gd [21], Tb [21], Dy [22], Ho [22], Er [23] and Lu [24,25], have assured that the QMCFMD methodology which includes two full hydration shells in the
∗ Corresponding author. E-mail address: [email protected] (B.M. Rode). http://dx.doi.org/10.1016/j.cplett.2015.08.041 0009-2614/© 2015 Elsevier B.V. All rights reserved.
quantum mechanical description of the system and a number of other improvements (vide infra) is a suitable tool for the study of such ions. This methodology has been employed, therefore, also for the present investigation of Tm(III) and Yb(III). The scope of this study was to reveal differences and similarities between these neighboring elements and to compare them on the other hand with their other neighbors Er(III) [23] and Lu(III) [24,25], respectively. 2. Methods In the QMCF molecular dynamics the quantum mechanically (QM) calculated part of the system, the QM region, consists of the solute and two complete hydration shells, with the advantages of describing all hydrogen bonds between 1st and 2nd shell by ab initio quantum mechanics and by shifting the smoothing zone between QM and molecular mechanics (MM) region quite far from the central solute, thus reducing the possible errors in the smoothing procedure [26]. The main advantage, however, is that due to the large distance between solute and MM region no potential functions for solute–solvent interactions are needed, as at that distance all non-Coulombic interactions have converged to zero, and a potential function is only required for solvent–solvent interactions. This is of particular importance in the case of highly charged ions, for which the construction of ion–solvent potentials on the basis of quantum mechanical calculations is difficult due to artificial charge transfers in pair potential evaluations. In practice, the QM region consists of two sub-regions, the core zone with ion and first hydration shell included, and the layer region, containing only solvent molecules. Further advantages of the QMCF-MD formalism are the inclusion of the MM charges as perturbation into the core
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
129
Table 1 ˚ Average Ln3+ –O1stshell distance in A. Ion
Ln3+ –O
Method
Reference
Tm(III)
2.44 2.33 2.42 2.32
QMCF-MD EXAFS QMCF-MD EXAFS
This work Persson et al. [16] This work Persson et al. [16]
Yb(III)
Figure 1. Definitions of the QM regions in QMCF-MD simulations.
Hamiltonian of the QM calculation of energies and forces, and the use of continuously newly determined QM charges in the evaluation of Coulombic interactions between QM and MM region. It is evident that the size of the two QM zones leads to a much higher computational effort than a conventional QM/MM approach with only one hydration shell. The increase in accuracy, especially for the dynamics of the system, fully justifies the cost of it. Further details and a mathematical description of the formalism are available in the literature [27–29]. The radius of the QM core region was ˚ that of the total QM region (core plus layer) to 5.7 A, ˚ set to 3.3 A, ˚ 1). smoothing being applied within 5.5 and 5.7 A(Figure For the quantum mechanical calculation of energies and forces ab initio calculations at Hartree–Fock (HF) level with dunning double zeta-polarisation (DZP) basis sets [30] for water and Stuttgart effective core potential (ECP) valence basis sets [31] for the lanthanoid ions have been employed. All quantum mechanical calculations were performed by Turbomole 6.4 [32]. Comparison of the results for ion–water clusters of other ions [19–24] at HF level with those of correlated methods have shown that the influence of electron correlation in such systems is very minor. Some previous results of similar simulations [19–24] have indicated that the relativistic correction via the ECP might be insufficient to take into account the relativistic contraction of the inner shells thus leading to an increase of the ion O distances compared to experiment [16]. On the other hand, the employment of correlated methods and even more, a full relativistic treatment of the systems in the simulations would lead to unrealistic computation times, as even the HF level
requires computation times of over 8000 h per system. Water in the MM-region was treated by the flexible Bopp–Jancsó–Heinzinger Center Force 2 model (BJH-CF2) [33,34], as this model proved most compatible with the Mullikan charges of the QM calculations. The simulations were carried out for a NVT ensemble consisting of one ion and 1000 water molecules at the experimental density of pure water and T = 298.15 K, thermostatised by the Berendsen algorithm [35], with a time step of 0.2 fs starting from a pre-equilibrated configuration of a simulation of Lu(III) in water. Equilibration required 20 000 steps, i.e. 4 ps, and sampling was performed for further 180 000 steps (36 ps). Hence the total computation time for the simulations presented here amounted to more than 200 000 CPUhours on a high-performance computer cluster with 12 dedicated 3.0 GHz processors, making use of 5 GB RAM disk. The trajectories have been analyzed with regard to radial distribution functions (RDF), angular distribution functions (ADF) and coordination number distribution (CND) for structural data. For the dynamics ligand mean residence times (MRT) and vibrational spectra, in particular frequencies and force constants of the Ln(III)–O stretching vibration were computed. The determination of the MRT was made by the ‘direct’ method [36] with a minimum displacement time from a hydration shell of t* = 0.5 ps. For the evaluation of the vibrational spectra the velocity autocorrelation functions (VACF) were Fourier-transformed in analogy to previous works [19–24]. 3. Results and discussion 3.1. Structural data Table 1 displays the characteristic data of the Ln(III)–O and Ln(III)–H RDFs of both ions, which are very similar. The absolute ˚ values for the first maxima differ by ∼0.1 Afrom the EXAFS data, ˚ but the small shortening of ∼0.01 Agoing from Tm(III) to Yb(III) is
Figure 2. (a) Ion O (black) and ion H (red) RDFs and their integrations of Tm3+ (upper) and Yb3+ (lower). The zoomed first peak is shown in the inserts, illustrating the presence of different hydrate structures. (b) Three-body correlation functions for the hydration shells of Tm3+ (left) and Yb3+ (right).
130
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
Figure 3. CNDs for the first and second hydration spheres of the hydrated Tm3+ (upper) and Yb3+ ion (lower).
well reproduced. The EXAFS data have been fitted to a model with six shorter (prisma) and three longer (caps) distances, and the values given in Table 1 are the average of these lengths. However, the RDFs of both ions (Figure 2) show upon zooming in tiny shoulder peaks near the first maximum of the ion O RDF, which could be an indication of nonequivalent bond lengths and will be discussed later in connection with the vibrational spectra. The too long absolute distances are another indication of an insufficient relativistic correction as observed in previous studies [19–24]. A closer inspection of the trajectories reveals that even without ligand exchange, reorientation processes of the first hydration shell occur on the picosecond scale, preventing the assignment of a unique structure to the hydrate. This seems to be a common characteristic of Ln(III) ions [19–24]. The ADFs of the first coordination shell are very similar again and are most compatible with a squared antiprismatic, a bi- and a tricapped trigonal prismatic structure (see Figure 4). The average coordination numbers of this shell resulting from the RDFs’ integration are 8.3 for Tm(III) and 8.4 for Yb(III), and the CND diagrams in Figure 3 display the slightly higher amount of 9-coordinated Yb(III) ion. The inserted snapshots in this figure are typical examples of the 8- and the 9-coordinated species. The dominance of coordination number 8 is not in full agreement with the EXAFS results [16] of 8.8 and 8.7 for Tm(III) and Yb(III), respectively. At this point is has to be mentioned, however, that fitting of a model for systems with very fast exchange and reorientation mechanisms seems an extremely difficult task in scattering experiments, and other literature data of coordination numbers ranging from 8 to 12 [37] seem to confirm this statement. The aforementioned EXAFS study concludes that almost all lanthanoid ions are up to Er(III) exclusively and even beyond dominantly 9-coordinated with only exception of Lu(III) (CN = 8.4), while other studies have discussed the so-called ‘Gadolinium break’, meaning a change of different physicochemical properties, possibly associated with a change of the coordination number from 9 to 8 at this element [38–40]. QMCF-MD simulations of the other immediate neighbors of Tm and Yb, Er and Lu have determined the average coordination numbers of 8.4 [23] and 8.0 [24,25] for these ions. The second hydration shells of Tm(III) and Yb(III) contain in average 20.1 and 19.4 water molecules, both extending over a wide range (cf. Figure 3), with a maximum at 19 and 20 ligands, respectively. The MRTs of these ligands are 3.1 and 2.7 ps, which significantly differs from the corresponding bulk value of water, being 1.6 ps [41]. While the shape of the RDFs (cf. Figure 2) shows a
Figure 4. ADF for the first hydration spheres of the hydrated Tm3+ (black) and Yb3+ ion (red), compared with the ADFs (blue lines) of an ideal (a) squared antiprism (CN = 8), (b) bicapped trigonal prism (CN = 8) and (c) tricapped trigonal prism (CN = 9).
distinct second shell connected to the first one mainly by H bonds according to the orientation of ligands, a third shell is hardly recognized. In order to confirm this behavior, the three-body analysis [42,43] was performed. Figure 2b proves that the first and second shell show clear deviations from the bulk structure, whereas the third shell shows identical features as the bulk. Evaluation of the second shell ligands’ MRT with t* = 0.0 ps provides an estimate of the average lifetime of the H bonds between first and second shell (each removal of a ligand from the second shell means breaking a hydrogen bond of this ligand to one in the first shell). The resulting values are 0.76 ps for Yb(III) and 0.68 ps for Tm(III), respectively. Consequently, the stability of the H bonds between first and second shell of both ions is considerably higher than between bulk water molecules, for which QMCF-MD simulations produced a value of 0.36 ps [41]. 3.2. Dynamics The mean residence times of water molecules in the first shell are 37 and 25 ps for Tm(III) and Yb(III), both lower than the corresponding value of Er(III) (41 ps), but also much lower than the value for the f14 ion Lu(III) (>160 ps). Apparently the MRTs are not or only loosely correlated to ion–ligand distance and coordination number. On the other hand, it seems that MRTs of trivalent lanthanoid ions are by many orders of magnitude shorter than those of trivalent d-transition elements and main group elements [44], almost comparable to monovalent and some divalent ions. This will have a considerable influence on their complex chemistry, e.g. in the binding of proteins utilized in specific separation processes recently [45]. The distance plots in Figure 5 provide some insight into the ligand exchange reaction mechanism. Although only a few
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
131
Table 2 First and second shell average coordination numbers, MRTs in ps, ion O stretching frequencies and corresponding force constants of Tm3+ , Yb3+ and other trivalent ions, all obtained by QMCF-MD studies. Ion 3+
Al Fe3+ Ir3+ La3+ Ce3+ Gd3+ Tb3+ Dy3+ Ho3+ Er3+ Tm3+ Yb3+ Lu3+
CN1
MRT1
CN2
MRT2
Qion
6 6 6 9.5 9.1 8.5 8.4 8.1 8.1 8.4 8.3 8.4 8
n.o. n.o. n.o. 16.6 >118 32.6 44.5 36.0 22.0 40.5 27.1 25.1 >160
11.8 13.6 13.5 25.6 21.4 19.8 20.1 19.2 19.5 18.8 20.1 19.4 18.7
17.7 3.1 3.6 2.3 2.6 3.0 2.7 3.0 2.1 3.0 3.1 2.7 2.6
560 513 n.a. 360 354 340 336 334 334 360 356 352 360
O
(cm−1 )
kion
O
(Nm−1 )
194 193 260 110 106 99 97 96 96 110 109 107 112
Reference Hofer et al. [46] Moin et al. [47] Pedevilla et al. [48] Lutz et al. [19] Lutz et al. [20] Canaval et al. [21] Canaval et al. [21] Tirler et al. [22] Tirler et al. [22] Canaval et al. [23] This work This work Hitzenberger et al. [24,25]
n.o.: no exchange processes observed during simulation time; n.a.: not available.
Figure 5. Distance plot of all Tm3+ –O (upper) and Yb3+ –O (lower) distances.
exchange events have been observed, they all seem to follow dominantly the associative (A) and less the dissociative (D) type, the latter exclusively starting from a 9-coordinated structure. The distance plots in Figure 5 also reflect the vigorous dynamics outside the first hydration shell, but the reduced flexibility of the second
shell is still recognized by a slightly less populated region between ˚ second shell and bulk around 6 A. MRTs are related to the strength of the ion–ligand binding, which is reflected in the ion O stretching vibration frequency and force constant. These values have been evaluated, therefore, via the VACFs and subsequent Fourier transformation, and they are collected in Table 2, together with some characteristic values of transition metal and main group elements. Figure 6 displays the vibrational spectra obtained from the simulations, separated for both ions for 8- and 9- coordinated species. Similar to the RDFs the main peak is accompanied by a smaller one or shoulder at lower wave numbers indicating a fraction of weakly bound water ligands, in agreement with the EXAFS result of two different types of water ligands coordinated to the ion [16]. Table 2 also lists the frequencies and force constants for peaks and shoulders resulting for Tm(III) and Yb(III) and compares them to characteristic values for other trivalent ions. The observed weakness of the ion O binding is in full agreement with the short MRT values and the generally much higher lability of water ligands bound to lanthanoid ions. However, the differences between Tm(III) and Yb(III) are very small and hence within the methodical accuracy limits for the main peak. The MRT differences can hence be explained only by the differences in the secondary peaks. As they result from the weaker bound ligands, which should dominate the exchange processes, these findings appear consistent not only with the interpretation of the EXAFS studies predicting two types of ligands, but also with our RDFs
Figure 6. Vibrational spectra of Tm3+ –O (a) and Yb3+ –O (b). Black: total, red: CN = 8 only, blue: CN = 9 only.
132
P.P. Passler, B.M. Rode / Chemical Physics Letters 638 (2015) 128–132
and the MRT data: Yb(III) with a considerably shorter MRT shows the frequency of the secondary peak at a much more red-shifted position (285 cm−1 for Tm(III) and 228 cm−1 for Yb(III) cf. Figure 6). The dynamics of the second shell have been discussed already in the previous chapter. 4. Conclusions The present work has once more demonstrated the suitability of the ab initio QMCF-MD methodology for the simulation of ions in aqueous environment, even for the fastidious treatment of the lanthanoid ions, where subtle differences of very similar but chemically often different elements have to be worked out. The use of HF level and relativistically corrected ECP appear a reasonable although not fully satisfactory compromise with respect to computational demand. The ions treated here assume hydration numbers of 8 and 9 in the first shell, and vigorous exchange processes on the picosecond scale lead to rapid changes between structure of the aquo-complexes. The mean residence times of the water ligands are considerably shorter than for their neighbor Er(III), despite a very similar force constant of the ion O bond. The striking differences between Ln(III) ions and transition element and main group ions of the same charge are also observed for the examples studied here. Acknowledgement Generous supply of computer time and a scholarship for Peter P. Passler by the University of Innsbruck is gratefully acknowledged. References [1] W. Koechner, Solid-State Laser Engineering; Springer Series in Optical Sciences, Springer, 2006, pp. p.49. [2] F.J. Duarte, Tunable Laser Applications, CRC press, 2010. [3] P. Lacovara, H.K. Choi, C.A. Wang, R.L. Aggarwal, T.Y. Fan, Opt. Lett. 16 (1991) 1089. [4] J. Koponen, M. Söderlund, H. Hoffman, S. Tammela, Opt. Express 14 (2006) 11539. [5] J.-F. Bisson, D. Kouznetsov, K.-I. Ueda, S. Fredrich-Thornton, K. Petermann, G. Huber, Appl. Phys. Lett. 90 (2007) 201901. [6] N.V. Sochinskii, M. Abellán, J. Rodríguez-Fernández, E. Saucedo, C.M. Ruiz, V. Bermúdez, Appl. Phys. Lett. 91 (2007) 202112. [7] D.A. Grukh, V.A. Bogatyrev, A.A. Sysolyatin, V.M. Paramonov, A.S. Kurkov, E.M. Dianov, IEEE J. Quantum Electron. 34 (2004) 247. [8] N. Krishnamurthy, C. Gupta, Extractive Metallurgy of Rare Earths, Taylor & Francis, 2004, pp. 32. [9] D. Krishnamurthy, V. Weinberg, J.A.M. Cunha, I.-C. Hsu, J. Pouliot, Brachytherapy 10 (2011) 461.
[10] A. Ayoub, G. Shani, World Congress on Medical Physics and Biomedical Engineering, September 7–12, Munich, Germany, 2009, p. 1. [11] B. Raj, B. Venkataraman, Practical Radiography, Alpha Science International Limited, 2004. [12] P. Enghag, Encyclopedia of the Elements: Technical Data - History - Processing - Applications, Wiley, 2008, pp. 448. [13] R. Halmshaw, Industrial Radiology: Theory and Practice, Non-Destructive Evaluation Series Springer, 1995, pp. 168. [14] Y. Sun, M. Yu, S. Liang, Y. Zhang, C. Li, T. Mou, W. Yang, X. Zhang, B. Li, C. Huang, F. Li, Biomaterials 32 (2011) 2999. [15] S. Aime, A. Barge, D. Delli Castelli, F. Fedeli, A. Mortillaro, F.U. Nielsen, E. Terreno, Magn. Reson. Med. 47 (2002) 639. [16] I. Persson, P. D’Angelo, S. De Panfilis, M. Sandström, L. Eriksson, Chem.-Eur. J 14 (2008) 3056. [17] M. Duvail, R. Spezia, P. Vitorge, ChemPhysChem 9 (2008) 693. [18] T.S. Hofer, A.K. Weiss, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 512 (2011) 139. [19] O. Lutz, T.S. Hofer, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 536 (2012) 50. [20] O. Lutz, T.S. Hofer, B.R. Randolf, B.M. Rode, Chem. Phys. Lett. 539 (2012) 50. [21] L.R. Canaval, P.P. Passler, B.M. Rode, Chem. Phys. Lett. 625 (2015) 116. [22] A.O. Tirler, P.P. Passler, B.M. Rode, Chem. Phys. Lett. 635 (2015) 120. [23] L.R. Canaval, T. Sakwarathorn, B.M. Rode, C.B. Messner, O.M. Lutz, G.K. Bonn, J. Phys. Chem. B 117 (2013) 15151. [24] M. Hitzenberger, T.S. Hofer, A.K. Weiss, J. Chem. Phys. 139 (2013) 114306. [25] M. Hitzenberger, Solvation Properties and Behaviour of Lutetium (III) in Aqueous Solution–A Quantum Mechanical Charge Field (QMCF) Study, University of Innsbruck, 2012 (M.Sc. thesis). [26] H. Lin, D.G. Truhlar, Theor. Chem. Acc. 117 (2007) 185. [27] T.S. Hofer, A.B. Pribil, B.R. Randolf, B.M. Rode, in: J.R. Sabin, E. Brändas (Eds.), Combining Quantum Mechanics and Molecular Mechanics. Some Recent Progresses in QM/MM Methods, vol. 59, Academic Press, 2010, p. 213 (Chapter 7). [28] T.S. Hofer, B.M. Rode, A.B. Pribil, B.R. Randolf, in: R. van Eldik, J. Harvey (Eds.), Theoretical and Computational Inorganic Chemistry, vol. 62, Academic Press, 2010, p. 143. [29] B.M. Rode, T.S. Hofer, B.R. Randolf, C.F. Schwenk, D. Xenides, V. Vchirawongkwin, Theor. Chem. Acc. 115 (2006) 77. [30] T.H. Dunning Jr., J. Chem. Phys. 53 (1970) 2823. [31] X. Cao, M. Dolg, J. Chem. Phys. 115 (2001) 7348. [32] TURBOMOLE V6.4 2012, A Development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989–2007, TURBOMOLE GmbH, since 2007. Available from http://www.turbomole.com [33] P. Bopp, G. Jancsó, K. Heinzinger, Chem. Phys. Lett. 98 (1983) 129. [34] F. Stillinger, A. Rahman, J. Chem. Phys. 68 (1978) 666. [35] H.J. Berendsen, J.P.M. Postma, W.F. van Gunsteren, A. DiNola, J. Haak, J. Chem. Phys. 81 (1984) 3684. [36] T.S. Hofer, H.T. Tran, C.F. Schwenk, B.M. Rode, J. Comput. Chem. 25 (2004) 211. [37] C.-H. Huang, Rare Earth Coordination Chemistry: Fundamentals and Applications;, John Wiley & Sons, 2010. [38] F. Spedding, M. Pikal, B. Ayers, J. Phys. Chem. 70 (1966) 2440. [39] F. Spedding, K. Jones, J. Phys. Chem. 70 (1966) 2450. [40] F. Spedding, M. Pikal, J. Phys. Chem. 70 (1966) 2430. [41] D. Xenides, B.R. Randolf, B.M. Rode, J. Chem. Phys. 122 (2005) 174506. [42] A. Bhattacharjee, T.S. Hofer, B.M. Rode, Phys. Chem. Chem. Phys. 10 (2008) 6653. [43] A.K. Weiss, T.S. Hofer, RSC Adv. 3 (2013) 1606. [44] F.A. Dunand, L. Helm, A.E. Merbach, Adv. Inorg. Chem. 54 (2003) 1. [45] M.R. Mirza, M. Rainer, C.B. Messner, Y. Güzel, D. Schemeth, T. Stasyk, M.I. Choudhary, L.A. Huber, B.M. Rode, G.K. Bonn, Analyst 138 (2013) 2995. [46] T.S. Hofer, B.R. Randolf, B.M. Rode, J. Phys. Chem. B 112 (2008) 11726. [47] S.T. Moin, T.S. Hofer, A.B. Pribil, B.R. Randolf, B.M. Rode, Inorg. Chem. 49 (2010) 5101. [48] P.A. Pedevilla, T.S. Hofer, B.R. Randolf, B.M. Rode, Aust. J. Chem. 65 (2013) 1582.