Comparison of prominent methods for computational studies of lanthanoid cation complexation

We compared different computational methods (quantum chemical and DFT) for calculations of binding energies of 8and 9-coordinated lanthanoid–aqua complexes. We used nine computational methods and compared the results with those obtained by the CCSD(T) method. All the nine methods provided relatively similar results and calculated energies correlated very well with the CCSD(T) obtained energies for complexes of this type. The comparison of basis sets revealed that combination of Dolg’s (5s5p4d)/[4s4p3d] + 2s1p1d basis set for lanthanoids and the cc-pvdz basis set for non-lanthanoids can be suggested as optimal for further studies of lanthanoids cation complexation.


INTRODUCTION
Lanthanoid compounds have been extensively used as luminescent chemosensors for medical diagnostics and optical cell imaging, contrast reagents for magnetic resonance imaging, shift reagents for NMR spectroscopy, as well as for applications in fundamental and applied sciences such as organic synthesis, bioorganic chemistry, and catalysis. Thus, the studies of lanthanoid complexes have shown a rising trend during the past few decades. Most of these studies mainly concern compounds where the metal ions have the common +3 oxidation state (Ln 3+ ) and behave as hard acids with a strong affinity for hard bases such as O-donors (neutral or negatively charged) or N-donors with which they form essentially non-directional bonds of a predominant ionic nature [1].
Research has also been focused on lanthanoids complexes and their applications in distinct devices such as chemical sensors, diagnostic systems, luminescent materials, and liquid crystals. The main advantages of luminescent lanthanoids complexes with chelating ligands in comparison to traditional organic fluorophores are the long luminescence lifetimes in addition to narrow emission bands and low concentration quenching [2].
Several reviews have been published on lanthanoids chemistry. In 1985 Rard [3] published a comprehensive overview on the chemistry and thermodynamics of europium. More recent reviews concentrate on organometallic complexes of lanthanoids [1,4] as well as on their luminescence properties [5]. In 2012 D'Angelo and Spezia [6] published a review of the hydration of lanthanoids and actinides from both theoretical and experimental aspects.
There are also theoretical studies on Eu 3+ -ligand complex formation [7,8,9]. Ricca et al. [7] performed ab initio calculations on Eu 3+ -L complexes, which can be used as benchmarks to assess the accuracy of calculations done at a lower level of theory on larger complexes. They focused on the interaction of Eu 3+ with CH 3 NH 2 as a model of lysine, S(CH 3 ) 2 as a model of methionine, and imidazole as a model of histidine. Respective results were compared with those obtained for Eu 3+ -H 2 O and Eu 3+ -H 2 S. It was shown that the binding energies follow the trends in dipole moment and polarizabilities as expected for electrostatic bonding.
Freire et al. [8] compared the efficacy of the semiempirical sparkle model (SMLC II, MOPAC 93r2) to ECP ab initio calculations for the prediction of ligand field parameters of europium(III) complexes using the quasi-relativistic effective core potential of Dolg et al. [10] with related [5s4p3d]-GTO basis set for Eu and HF/STO-3G and HF/3-21G basis sets for all other atoms. They carried out a comparison between the predicted and crystallographic ground state geometries and concluded that SMLC II is applicable in molecular biology where biotechnological applications of very large systems are addressed.
Freire et al. [9] also studied the structure and emission spectrum of the europium cryptate containing 3,3-biisoquinoline-2,2-dioxide with the framework of the Sparkle/AM1 model and Judd-Ofelt theory. They found a satisfactory similarity between the theoretical study and X-ray measurements found in the literature but also stated that some differences can be attributed to the fact that the optimization was performed in vacuum. However, they once again demonstrated that the Sparkle/ AM1 model can be used for simulating the larger Eu 3+ molecular systems.
Albuquerque et al. [11] combined an empirical relation involving the energy of the lowest triplet state of Tb 3+ complexes with semi-empirical quantum mechanical calculations and studied the relation between the structure and luminescent properties of Tb 3+ complexes containing β-diketonate.
Studies on the complexation of actinides by mono-, di-, and hydroxycarboxylic acids also show a rising trend, mainly for two reasons. First, carboxylic acids present in the natural waters may play an important role in influencing the migration of actinides in the aquatic environment around a deep geological repository. They are also used as model compounds to investigate the mechanism of the binding of actinides by humic substances in the natural waters. Second, hybrids containing rare earths and organic moieties have been found to enhance the optical properties of rare earth elements. For instance, the luminescence of Eu(III) was found to be sensitized by complexation with carboxylic acids [12].
A number of studies have been carried out for computational simulation of the hydration of lanthanoids. Terrier et al. [13] studied structural and electronic properties of La 3+ immersed in bulk water by means of DFT-based CPMD simulations. They reported firstprinciple molecular dynamics simulations of La 3+ in bulk liquid water that are able to reproduce structural properties of this system, which includes not only the La-water distance in the first hydration shell but also a stable coordination number equal to nine for the time spanned by the simulation. Zhang et al. [14] studied the hydration of all trivalent lanthanum ions theoretically from two aspects: energy and wave function. In comparing two computational methods, they stated that SCS-MP2 [15] is nearly as accurate as the computationally more demanding CCSD [16] and enables giving the most accurate first principle hydration Gibbs free energies and reliable preferred coordination numbers of lanthanoid(III)aqua complexes: 9, 8, and both, for light, heavy, and intermediate lanthanoids, respectively. In the current study, we performed calculations for 15 eight-and nine-coordinated lanthanoids by using nine computational chemistry methods.
Full geometry optimizations and vibrational analyses were performed for all reactants and products of the studied hydration reactions at the BP level with the def-TZVP basis set. All stationary points were found to be true minima (NImag = 0). For the rest of the calculations only geometry optimizations were performed and the correction terms calculated at the BP/def-TZVP level were applied.

RESULTS AND DISCUSSION
Hydration binding energies were calculated for all lanthanoids (Ln) by Eq. (1): where n is 8 or 9.
The hydration energies for eight-and nine-coordinated aqua complexes calculated with the def-TZVP basis as described in Methods are presented in Tables 1  and 2.
With this basis the best correlation between the observed and the calculated hydration energies was obtained with the B3LYP and B3LYP-D3 functionals (R 2 = 0.97 for both functionals, Fig. 1). These functionals are closely followed by PBE0 with the R 2 value of 0.97. The order of functionals based on R 2 values is as follows: B3LYP, B3LYP-D3 > PBE0 > TPSSh > TPSS, PBE > BP, BP-D3 > MP2. However, the slope, intercept, and standard error values should also be considered in the comparison of functionals presented in Tables 1 and 2 Table 1 reveal that this is only in the case of lanthanum. Similar effects can be seen for all other tested computational methods. The reason is illustrated in Fig. 1, as the hydration energies for La are smaller than calculated by Zhang et al. [14], which has a great effect on the    [14].  Table 1. The deviation of calculated results for lanthanum is most likely caused by its somewhat different ECP, which replaces 46 inner electrons instead of 28 as for the rest of lanthanoids.
Removal of La from the data set greatly improved the correlation between previously calculated hydration energies and the results presented here. The best results were obtained with the MP2 method as it has the best R 2 value, the intercept closest to one, and the smallest standard error values. However, B3LYP and B3LYP-D3 still offer the best results with the highest R 2 , the intercept closer to zero, and the lowest standard error values among the DFT methods. The studied DFT functionals can be considered as equally good for the differences in correlation results are very small. It must be noted that the values of hydration energies are relatively large and most of the standard error comes from the intercept. Therefore, a small difference in intercept values between functionals is not the best measure to compare different functionals.
In the case of nine-coordinated lanthanoid complex formation energies, the whole picture is pretty much the same. Calculated energies are presented in Table 2 and illustrated in Fig. 2. Updated statistical parameters after the removal of La are also presented in Table 2. Similarly to eight-coordinated complexes, the MP2 was by far the best method when considering the R 2 , intercept, and standard error values. However, B3LYP and B3LYP-D3 were the best DFT methods because of small intercept and standard error values.

Comparison of basis sets
As demonstrated, all functionals provide rather similar results for water binding energy calculations. Therefore, for comparison of basis sets, we chose the BP functional as the fastest. We tested the binding energies for eightand nine-coordinated lanthanoid complexes using the BP functional and the second set of basis sets described in the methodology section (pseudopotentials where the core-shell electrons were substituted by the scalarrelativistic 4f-in-core pseudopotential [10,33] and the corresponding (7s6p5d)/[5s4p3d] [10], (5s5p4d)/[4s4p3d] + 2s1p1d [10,34], and (6s6p5d)/[4s4p4d] + 2s1p1d [10,34,35] basis sets for lanthanoids, and cc-pvdz, aug-cc-pvdz, or aug-cc-pvtz [36,37] basis sets (altogether 9 combinations) for the other atoms). The best results obtained with the (5s5p4d)/[4s4p3d] + 2s1p1d [10,34] basis set for lanthanoids are given in Table 3 and illustrated in Fig. 3. The results with two other basis sets for lanthanoids are given in Supporting Information.
Comparison of the results in Tables 1, 2, and 3 clearly indicates that the (5s5p4d)/[4s4p3d] + 2s1p1d [10,34] basis, especially with cc-pvdz or aug-cc-pvdz for nonlanthanoid atoms, performs much better than the def-TZVP basis. All lanthanoids are treated equally well and there are no clear outliers as was La in the case of def-TZVP. We are tempted to suggest the (5s5p4d)/[4s4p3d] + 2s1p1d basis set for lanthanoids and cc-pvdz basis set for non-lanthanoids for further studies of lanthanoid complexation as a good compromise between cost and accuracy. However, one should also keep in mind that this suggestion is based only on hydration energies and in the case of softer bases (like N-or S-bases) the inclusion of f-electrons in ECP might lead to some errors.

CONCLUSIONS
In the current study we demonstrated that nine quantum chemical methods could equally well be used alternatively to the computationally costly CCSD(T) method in the calculation of binding energies for eight-and ninecoordinated water complexes of lanthanoid cations with the def-TZVP basis. All the used computational methods -BP, BP-D3, B3LYP, B3LYP-D3, PBE, PBE0, TPSS, TPSSh, and MP2 -provided relatively similar results and correlated well with the respective results from the CCSD(T) method.
Comparison of Dolg's f-electrons in-core pseudopotentials and the corresponding basis sets for lanthanoids together with Dunning basis sets on non-lanthanoid atoms revealed that these basis set combinations may be even more accurate than the Stuttgart-Dresden pseudopotentials, which replaced 28 core electrons (46 in the CCSD(T) B3LYP case of La) with the corresponding def-TZVP basis sets on lanthanoids and the TZVP basis on other atoms.
The combination of the (5s5p4d)/[4s4p3d] + 2s1p1d basis set for lanthanoids and the cc-pvdz basis set for non-lanthanoids can be suggested as optimal for further studies of lanthanoid cation complexation as it is a good compromise between cost and accuracy. However, as the current evaluation is based only on hydration energies, one must be careful when extending these conclusions to softer bases. Table 3. Gas phase binding energies for eight-and nine-coordinated lanthanoids (kJ/mol) using the scalar-relativistic 4f-in-core pseudopotential [10,33] and the (5s5p4d)/[4s4p3d] + 2s1p1d [10,34]