Molecular dynamics simulations of Y(III) coordination and hydration properties

Abstract

Y mainly exists in ionic rare-earth resources. During rare-earth carbonate precipitation, rare-earth ion loss in the precipitated rare-earth mother liquor often occurs due to CO₃²⁻ coordination and Y(III) hydration. Microscopic information on the coordination and hydration of CO₃²⁻ and H₂O to Y(III) has not yet been elucidated. Therefore, in this study, the macroscopic dissolution of Y(III) in different aqueous solutions of Na₂CO₃ was studied. The radial distribution function and coordination number of Y(III) by CO₃²⁻ and H₂O were systematically analyzed using molecular dynamics (MD) simulations to obtain the complex ion form of Y(III) in carbonate solutions. Density functional theory (DFT) was used to geometrically optimize and calculate the UV spectrum of Y(III) complex ions. This spectrum was then analyzed and compared with experimentally determined ultraviolet-visible spectra to verify the reliability of the MD simulation results. Results showed that Y(III) in aqueous solution exists in the form of [Y·3H₂O]³⁺ and that CO₃²⁻ is present in the bidentate coordination form. In 0–0.8 mol L⁻¹ CO₃²⁻ solutions, Y(III) was mainly present as the 5-coordinated complex [YCO₃·3H₂O]⁺. When the concentration of CO₃²⁻ was increased to 1.2 mol L⁻¹, [YCO₃·3H₂O]⁺ was converted into a 6-coordinated complex [Y(CO₃)₂·2H₂O]⁻. Further increases in CO₃²⁻ concentration promoted Y(III) dissolution in solution in the form of complex ions. These findings can be used to explain the problem of incomplete precipitation of rare earths in carbonate solutions.

---

1. Introduction

The demand for elemental Y has increased yearly due to its excellent optical applications. Rare-earth carbonates are important inorganic salts used in the extraction and application of rare-earth elements.¹ Y mainly exists in ionic rare-earth resources and xenotime, and the solubility of Y and Ln carbonates in water is merely 10⁻⁵ to 10⁻⁷ mol L⁻¹; therefore, carbonate is used in the industry to precipitate rare-earth ions.^2–4 However, rare-earth ions in the mother liquor are considerably lost due to the coordination and hydration in the carbonate solution. The microscopic information of the coordination and hydration by CO₃²⁻ and H₂O on rare-earth ions has not been elucidated.^5,6

From the 1970s to the 1990s, many scholars have researched the solubility and complex ion morphology of Ln elements in carbonate solution. The amount of Ln dissolved in Na₂CO₃ solution and the ion morphology are mainly related to the concentration of CO₃²⁻. Taketatsu^7,8 studied the solubility of rare-earth elements in different concentrations of K₂CO₃ solution and showed that the dissolved amount of rare-earth elements in K₂CO₃ solution increases with the carbonate concentration and the atomic number (except for Ce, Y, and Sc). Tran⁹ investigated the geological migration of Ce in fractured carbonate rock formations and found that in a groundwater environment with pH 7–9, the main rare-earth carbonate complex ions are RECO₃⁺ and RE(CO₃)₂⁻; complex ions, such as REHCO₃²⁺ and REOH²⁺, are nearly negligible. In a study on low concentration CO₃²⁻ natural groundwater, Noack et al.¹⁰ found that heavy rare-earth elements (HREEs) are more easily coordinated with CO₃²⁻ than light rare-earth elements (LREEs) and produce RE(CO₃)₂⁻. Johannesson^11,12 researched the geochemical behavior of CO₃²⁻ and RE³⁺ in Mono Lake water; the forms of rare-earth carbonate complex ions are RE(CO₃)₂⁻ at 0.27 M CO₃²⁻ in lake water and RECO₃⁺ in lower-concentration CO₃²⁻ groundwater. Philippini^13,14 found that RE(CO₃)₂⁻, RE(CO₃)₃³⁻, and other complex ions exist in aqueous solutions with different CO²⁻₃ concentration gradients (0.1–1.5 mol L⁻¹). In summary, the complexation mode and morphology of rare-earth ions change with the concentration of CO₃²⁻.

Reiller et al.^15,16 performed inductive analysis of the existing migration law of Am and associated it with rare-earth elements; they concluded that when the pH of the solution is less than 7, between 7 and 8.5, and above 8.5, the main complex ions are RE³⁺, RECO₃⁺ and RE(CO₃)₂⁻, and possibly RE(CO₃)₃³⁻ or higher-valence-state complex ions, respectively. Given similarities in the physicochemical properties of Am and rare-earth elements, Maloubier et al.¹⁷ investigated the results of Am migration studies and concluded that the main forms of rare-earth carbonate ions in a water environment with pH 7–9 are RECO₃⁺ and RE(CO₃)₂⁻. This work thus indicated that the morphology of rare-earth complex ions is also related to the pH of the solution.

Using available crystallographic and spectral data, including ultraviolet-visible [UV-vis], near-infrared, and infrared spectra, Jeanvoine¹⁸ stated that all Ln(III) ions form tetra-carbonate ions in aqueous solution when carbonate ions are not restricted. Philippini et al.^13,14 used electrophoretic mobility shift and time-resolved laser-induced fluorescence spectroscopy (TRLFS) to study several Ln(III)–carbonate complexes and concluded that light rare-earth Ln(III) ions can coordinate 4 carbonate ligands, whereas heavy rare-earth Ln(III) can only coordinate with 3. Modern research methods have enabled significant advances in the analysis of the forms of ions in rare-earth carbonate solutions,¹⁹ but nearly all available studies involved are extrapolated analyses of different rare-earth carbonates.²⁰ The structures of rare-earth complex ions in carbonate aqueous solutions are complex and diverse, and their specific forms remain unclear; thus, difficulties have been encountered when directly using UV-vis, near-infrared, infrared, X-ray, and TRLFS to clarify the structure of rare earth-complex ions.^21,22

Martelli et al.²³ achieved remarkable results through polarization MD simulation studies and examined the structure and kinetics of tri-and tetra-carbonate complexes. Jeanvoine et al.¹⁸ used the density functional theory (DFT) method to study the bonding theory of Ln(III) tri- and tetra-carbonate ions and provided the following conclusions: (1) in the gas phase, the most stable structures are full bidentate and monodentate for tri-carbonate and tetra-carbonate complexes, respectively; (2) in the aqueous phase, when the interaction of water molecules is not explicitly considered in the continuous model, the full bidentate complexes are the most stable structures; and (3) Ln(III)–carbonate interactions mainly occur through ionic bonds.

After decades of development in the field of rare-earth metals, many scholars have studied the macroscopic dissolution of aqueous solutions in rare-earth elements in carbonates. However, studies on the microstructure and coordination behavior with CO²⁻₃ and H₂O in carbonate solution are scarce. Therefore, the coordination and hydration microscopic information of Ln(III) ions must be further examined.

Y and Ln do not belong to the same period in the periodic table of elements, and their atomic and ionic electronic layer structures considerably differ. However, the properties of Y are substantially similar to those of heavy rare-earth elements, and this similarity requires deep understanding. Therefore, the dissolution law of Y and the microscopic information of its coordination and hydration should be investigated. Motivated by previous studies on the dissolution of Ln elements and complex ion morphology, in the present work, we studied the microscopic interaction of Y with CO₃²⁻ and H₂O in aqueous carbonate solutions. We performed a dissolution experiment and considered the effect of CO₃²⁻ coordination and hydration in a real solution using MD simulation. Then, we performed DFT quantum chemistry calculations, UV-vis spectral analysis, inversion, and in-depth analysis [Y(CO₃)_m·nH₂O]^3−2m of the morphology and ionic structure of Y to clarify its evolution mechanism in carbonate solution. The results of this work could explain the problem of rare-earth yield in ionic rare-earth precipitation.

2. Experimental materials and research methods

2.1 Experimental materials and equipment

The raw materials used in this experiment were obtained from a high-purity YCl₃ solution produced by a rare-earth separation plant in Longnan county, Ganzhou. The main component of the YCl₃ solution was prepared by oxalic acid precipitation, burning, HCl dissolution, and determined by ICP-MS and are shown in Table 1. Other chemical reagents used, such as Na₂CO₃ and H₂O, were of analytical purity. The experimental equipment or instruments used in this work are shown in Table 2.

Table 1 Contents of each element and REE distribution in high purity YCl₃

Y^3+, mol L^−1	H^+, mol L^−1	Density, g mL^−1	Rare earth impurities/REO (μg mL^−1)
0.3419	<0.10	1.0568	La2O3	CeO2	Pr6O11	Nd2O3	Sm2O3	Eu2O3	Gd2O3
			<10	<10	<10	<10	<10	<10	<10
			Tb2O3	Dy2O3	Ho2O3	Er2O3	Tm2O3	Yb2O3	Lu2O3
			<10	<10	<10	<10	<10	<10	<10

---

Table 2 Experimental equipment specific information table

Equipment	Type specification	Manufactures
High-speed centrifuge	TGL16MS	Yancheng Anxin Experimental Instrument Co., Ltd.
Computer server	IBM System X3850 X5	International Business Machines Corporation
UV-visible spectrophotometer (UV-vis)	UV-5500PC	Shanghai yoke instrument Co., Ltd.
Inductively coupled plasma-optical emission spectroscopy (ICP-OES)	ULTIMA2	HORIBA Jobin Yvon
pH meter	KL-009	Xuzhou Yaming Instrument Co., Ltd.

---

2.2 Research method

2.2.1 Instantaneous saturated dissolution of Y(III). The calibrated YCl₃ solution was added dropwise to a series of 50 mL solutions containing 0.4, 0.8, 1.2, 1.6, or 2.0 mol L⁻¹ Na₂CO₃ keeping in a conical flask at 35 °C under water-bath conditions. The conical flask was then removed from the water bath and shaken. After drop-wise addition of the YCl₃ solution to the carbonate solutions, a white precipitate was produced. As the conical flask oscillated, the white precipitate dissolved. Addition of the rare-earth liquid was continued until dissolution of the white precipitate could no longer be completed. Next, the solution was centrifuged for 10 min at 35 ± 1 °C and 9000 rpm. The supernatant was titrated with EDTA to determine the molar concentration of Y(III), which was defined as the instantaneous saturated dissolution amount of Y(III). The instantaneous saturated molar concentrations of Y(III) in different concentrations of Na₂CO₃ solution was calculated at 35 ± 1 °C.

2.2.2 Molecular dynamics calculations. Molecular Dynamics (MD) is a method of analyzing the movement of ions or molecules in a system from the microscopic level; here, theoretical simulation of various characteristics of the system is achieved through correlation calculation of statistical mechanics.

For the microscopic behavior of Y element in carbonate solution system, according to the instantaneous saturated solution of Y(III) in different concentrations of Na₂CO₃ solution, the SPC/E water model²⁴ was used to calculate the MD of its real solution system. The radial distribution function (RDF) and coordination number of each ion pair in the carbonate solution system were calculated and analyzed, and the Y(III) coordination form considering CO₃²⁻ was determined.

Geometric and energy optimization of Y³⁺, Cl⁻, H₂O, and CO₃²⁻ contained in the solution system were achieved by using the Forcite module in the calculation platform Materials Studio 8.0.^25,26 The initial simulation model was a cube system established by the Amorphous cell module in Materials Studio 8.0. The water molecules and ions were randomly distributed in the system. The universal force field was used,²⁷ and the system was geometrically optimized by the Smart method. The parameters were preset in the operation, the calculation step was 1 fs, the total calculation time was 500 ps, the temperature control setting mode was set to Andersen,²⁸ the system temperature was based on 308 K, and the NVT ensemble was selected. The accuracy of the simulation calculations was set to fine. Calculation results were deemed credible when the temperature was stable within 308 K ± 10% and no large temperature disturbance could be observed.^29,30 The YCl₃–H₂O solution system was simulated, and the parameters of each particle are listed under model a in Table 3. YCl₃–Na₂CO₃ solution systems with Na₂CO₃ concentrations of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ were simulated. The parameters of the corresponding solution particle components are shown in models b, c, d, e, and f in Table 3.

Table 3 Modeling parameters of solution components

Model	Solution system	Projects	Composition of each particle of the model
a	YCl3–H2O solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
		Number	5	15	—	—	4444
	ρ: 1.009 g mL^−1	Weight%	0.5	0.7	—	—	98.8
b	YCl3–Na2CO3 solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
	0.4 mol L^−1 Na2CO3	Number	1	3	24	48	4388
	ρ: 1.075 g mL^−1	Weight%	0.1	0.1	1.8	1.3	96.7
c	YCl3–Na2CO3 solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
	0.8 mol L^−1 Na2CO3	Number	2	6	48	96	4312
	ρ: 1.095 g mL^−1	Weight%	0.2	0.3	3.5	2.7	93.3
d	YCl3–Na2CO3 solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
	1.2 mol L^−1 Na2CO3	Number	3	9	72	144	4236
	ρ: 1.115 g mL^−1	Weight%	0.3	0.4	5.2	4.0	90.1
e	YCl3–Na2CO3 solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
	1.6 mol L^−1 Na2CO3	Number	4	12	96	192	4160
	ρ: 1.136 g mL^−1	Weight%	0.4	0.5	6.8	5.2	87.1
f	YCl3–Na2CO3 solution	Component	Y^3+	Cl^−	CO3^2−	Na^+	H2O
	2.0 mol L^−1 Na2CO3	Number	5	15	120	240	4084
	ρ: 1.156 g mL^−1	Weight%	0.5	0.6	8.4	6.4	84.1

---

2.2.3 UV-vis full wavelength scanning analysis. The coordination form of Y(III) in carbonate solution was quickly and economically analyzed³¹ by performing UV-vis full-wavelength scanning on different carbonate concentration solution systems. Differences in the absorption spectrum of each solution system after addition of Y(III) to the different carbonate concentration solutions were analyzed. Using a burette, 1 mL of YCl₃ was added to 25 mL of 0.4, 0.8, 1.2, 1.6, or 2.0 mol L⁻¹ Na₂CO₃ solution by using a burette. The mixtures were shaken well at room temperature (30–35 °C) and then centrifuged for 10 min at 9000 rpm. Then, the supernatant was subjected to UV-vis full-wavelength scanning. The wavelength range was set to 190–1100 nm, deionized water was used as the reference sample, and the optical path was a 1 cm quartz glass cuvette. Given the free CO₃²⁻ in solution, 1 mL of pure water was added dropwise to 25 mL of 0.4, 0.8, 1.2, 1.6, or 2.0 mol L⁻¹ Na₂CO₃ solution for experimental comparison. The solution with pure water was utilized as the blank control, and all UV-vis absorption spectra were measured by the same test method.

2.2.4 DFT geometry optimization. Based on the chemical composition of the Y(III) complex ions obtained in the MD simulations, geometric optimization of the Y(III) complex ion structure and calculation of the characteristic UV absorption spectrum were performed by DFT.³² The DMol3 module in Materials Studio 8.0 was used for DFT computation.³³ The GGA/PBE method was selected,³⁴ and the base was set to DND 4.4. TD-DFT calculation was conducted using the COSMO solvation model. The integration precision, SCF tolerance, and track cutoff were set to fine. During MD simulations, the possible ion structures were obtained, geometric optimization and total energy calculation of the system were conducted, and excited-state UV spectra were calculated using the same basis set and method.

3. Results and discussion

3.1 Instantaneous saturated dissolution of Y(III)

The instantaneous saturated molar concentrations of Y(III) in 50 mL of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃ solution at 35 °C are shown in Fig. 1.

Fig. 1 Concentration relationship between instantaneous saturation solubility of Y element with CO₃²⁻.

The experiments confirm that the dissolution process of Y occurs in the Na₂CO₃ solution. A precipitate form when YCl₃ is added to the Na₂CO₃ solution, but this solid is dissolved as the reaction proceeds, thus leaving a clear solution. At this time, Y(III) in the solution does not precipitate as a rare-earth carbonate but is dissolved in the Na₂CO₃ solution. However, when the added amount of Y(III) exceeds the instantaneous saturated dissolved amount, the solution no longer clears; instead, it assumes the form of a suspension. In this case, part of Y(III) in the solution is precipitated; the amount of dissolved rare-earth ions when the flocculent turbidity is produced is considered the instantaneous saturated dissolved amount. Fig. 1 shows that the instantaneous saturated dissolved amounts of Y in 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃ solutions are 0.0041, 0.0163, 0.0377, 0.0693, and 0.0909 mol L⁻¹, respectively. The relationship between the instantaneous saturated dissolved amount of Y and the concentration of CO₃²⁻ in the base solution was analyzed, and Fig. 1 shows that the Y(III) instantaneous saturated dissolved amount has a linear relation with the concentration of CO₃²⁻. The equation describing this relationship is y = 0.0416x + 0.0074, and the R² is 0.985. This result provides a basis for designing a model of a real solution with different concentrations of CO₃²⁻.

3.2 MD

3.2.1 Radial distribution and coordination number. Radial distribution function (RDF) analysis is conducted on the MD calculation results, and the RDF between the components in the system is associated with the coordination number. Eqn (2)³⁵ is used to calculate the average coordination number between components. The definitions of g(r) and Cn are shown in Fig. 2.

Fig. 2 RDF schematic.

The calculation formulas are as follows:

where V is the total volume of the system; N is the total number of particles in the system; n(r, Δr) is the number of particles from r to r + Δr; ρ is the bulk density (ratio of the number of particles to the bulk volume); Cn(r) is the number of coordination atoms (molecules) in the range of 0–r spherical shell around the target atom; ρ is the number density of the coordination atom (or molecular), that is, the ratio of the number of atoms (molecule) to the volume of space; g(r) is the value of the radial distribution function, indicating the probability of occurrence of a coordinating atom (or molecular) within a certain distance; r is the intercept radius.

Before the results are described and discussed, a representation the RDF for the studied simulation system should be defined. O_C(n) is an oxygen atom of carbonates in the first coordination shell of Y(III), and O_C(f) is an oxygen atom farther away from Y(III) in the carbonates. In order to directly investigate the distribution of H₂O around Y(III), H₂O is considered a whole particle to analyze the RDF of H₂O for Y(III) and calculate the coordination number. After each solution simulation is performed and the MD simulations are converged, the RDF and coordination number of Y³⁺–O, Y³⁺–H₂O, and Y³⁺–Cl⁻ are calculated separately. The simulation results are shown in Fig. 3–5.

Fig. 3 RDF and coordination number of each ligand for Y(III) in YCl₃–H₂O solution system. (Model a.)

Fig. 4 RDF and coordination number of each ligand for Y(III) in YCl₃–H₂O solution system. (Model b: b1 and b2; Model c: c1 and c2.)

Fig. 5 RDF and coordination number of each ligand for Y(III) in YCl₃–H₂O solution system. (Model d: d1 and d2; Model e: e1 and e2; Model f: f1 and f2.)

Fig. 3–5 show the RDFs of Y³⁺ and O_C(n), O_C(f), H₂O, and Cl⁻; the RDF directly reflects the microscopic information of the ions in the rare-earth carbonate solution system. According to the figures (a, b2, c2, d2, e2, and f2), the average coordination of Y³⁺ and Cl⁻ has zero RDF intensity, and Y³⁺ and Cl⁻ coordination does not occur. These results indicate that the interaction of Y³⁺ with Cl⁻ in the solution may be negligible. The Y(III) short-range has two hydrated shells in the YCl₃–H₂O solution and the YCl₃–Na₂CO₃ solution system; in the central ion close-range ordered layer, RDF counts all ordered short-range particles.³⁶ In existing research and literature, we can conclude that the hydration number of Ln elements decreases with the increase of atomic number, their effective hydration number decreases with the ionic radius reduced. Especially in the diluted LaCl₃ and CeCl₃ aqueous solution, the hydration number of La and Ce can reach 8–9. While Y element is one cycle shorter than the Ln elements, and its ionic radius is almost smaller than that of the Ln elements. In this paper, in the pure YCl₃ aqueous solution in simulation a, the number of water molecules in the first shell of H₂O to Y(III) is 3.0, when the g(r) was integrated to the second RDF peak valley, and the short-range H₂O molecule number is 9.1, but we consider that only 3.0 H₂O molecules in the first shell coordinate with Y(III), and the H₂O molecules in the second shell are considered to be under the action of electrostatic attraction. With increasing CO₃²⁻ concentration in the YCl₃–Na₂CO₃ solution system, eqn (2) is integrated to the second peak valley of the Y(III) short-range order,^37,38 that is, in the short-range order, the total H₂O molecular number decreases from 9.1 to 7.0 while the coordination number of H₂O for Y(III) decreases from 3.0 to 2.2 within r < 2.6 Å. Interestingly, research results reported by Migliorati et al.³⁹ shows, analysis of the simulation results has shown that both Lennard-Jones and Buckingham potentials are able to properly describe the radial distribution of water molecules around the Ln³⁺ ions, the smooth decrease of the hydration number along the lanthanoid series. In model a, within the first shell layer of Y(III), there are only 3 hydration numbers. In fact, the ionic radius of Y(III) is almost smaller than the ionic radius of the lanthanide, resulting in less coordination number of H₂O in first shell layer than that of lanthanides. This does not violate Migliorati's conclusion.

We consider the coordination of CO₃²⁻ with Y(III). Fig. 3–5 show that in the first coordination shell, formula (2) is integrated to r < 2.6 Å; that is, models b, c, d, f, and e, are integrated to the first peak valleys of the RDF of CO₃²⁻ for Y(III) coordination. While eqn (2) is integrated to r < 4.0 Å, the five models are integrated to the second peak valleys of the RDF; thus, all the atoms in the short-range order are statistically analyzed. The micro information of the MD simulation results of models a, b, c, d, e, and f is summarized in Table 4 for analysis and comparison.

Table 4 MD simulation calculation results

Model	r Y–OC(n) (Å)	CnY–O (r < 2.6 Å)	r Y–OC(f) (Å)	CnY–O (r < 4.0 Å)	N OC(n)	N OC(f)	N OC(n)/NOC(f)	r Y–H2O (Å)	CnY–H2O (r < 2.6 Å)
a	—	—	—	—	—	—	—	2.27	3.0
b	2.15	2.0	3.75	3.0	2.0	1.0	2.0	2.35	3.0
c	2.15	2.1	3.75	3.1	2.1	1.0	2.1	2.35	2.8
d	2.15	4.0	3.75	6.0	4.0	2.0	2.0	2.37	2.6
e	2.15	4.0	3.75	6.0	4.0	2.0	2.0	2.40	2.4
f	2.25	4.1	3.75	6.1	4.1	2.0	2.0	2.40	2.2

---

Fig. 3–5 show that the O atom of the carbonate is short-range ordered. Moreover, O_C(n) is in the first shell, and O_C(f) is in the second shell. According to Table 4, when CO₃²⁻ is 0.4–1.6 mol L⁻¹, r_Y–OC(n) is 2.15 Å, and only when CO₃²⁻ increases to 2.0 mol L⁻¹ does r_Y–OC(n) increase to 2.25 Å; in addition, r_Y–OC(f) of all MD results is 3.75 Å. Therefore, in the YCl₃–Na₂CO₃ solution system, the change in CO₃²⁻ concentration has nearly no effect on the peak position of the RDF; it affects the peak of g(r) only.

According to the specific simulation values in Table 4, we initially predict, the ionic composition of Y(III) in the YCl₃ solution is [Y·3H₂O]³⁺. Y(III) mainly exists in the form of a complex [YCO₃·3H₂O]⁺ in 0–0.8 mol L⁻¹ CO₃²⁻ solution and in the form of a complex [Y(CO₃)₂·2H₂O]⁻ in 1.2–2.0 mol L⁻¹ CO₃²⁻ solution. In the solution with high CO₃²⁻ concentration, H₂O in Y(III) short-range is crowded by the CO₃²⁻ space. In the first water molecule shell, the coordination number of H₂O for Y(III) is reduced to 2.2, consistent with the conclusion reported by Martelli.²³ In the high-CO₃²⁻-concentration solution, the hydration number of heavy rare-earth carbonate Ln(III) is 4–2. When CO₃²⁻ ≤ 0.8 mol L⁻¹, one CO₃²⁻ passes through the second shell of the H₂O molecules from the outer layer and then enters the first shell to coordinate with Y(III). At this point, the effect of H₂O molecule crowding is not noticeable. However, when CO₃²⁻ is further increased, the second CO₃²⁻ passes through the second H₂O ordered layer from the outer layer then and enters the first shell to coordinate with Y(III). Under this condition, short-range CO₃²⁻ space crowding results in a decrease of 1–2 H₂O molecules around the Y(III) short-range. This result indicates that the interaction between CO₃²⁻ and Y(III), is that CO₃²⁻ is originally due to electrostatic attraction at a distant location, passing through the second layer of H₂O molecules, and entering the inner shell. Then, a coordination bond with Y(III) occurs.

It is worth mentioning, carbonate not only complexes rare-earth ions but also precipitates with rare-earth ions in carbonate solutions. In the high concentration CO₃²⁻ solution system of this paper, due to the higher CO₃²⁻ ion strength, CO₃²⁻ will be squeezed into the short-range. It has been reported that with the increase of CO₃²⁻ concentration, in short-range of H₂O will replace 1–2, a CO₃²⁻ replaces a H₂O molecule, resulting in increase the coordination number 1 because there is two O elements of CO₃²⁻ coordinated with Y(III). In a certain sense, high concentration of CO₃²⁻ seem to have interactions like counterions on CO₃²⁻ coordination, related theory of counterions on the hydration structure of lanthanide ions can be used to study this problem.⁴⁰

3.2.2 Coordination mode of CO₃²⁻ for Y(III). The MD simulation results show the change trend of the average coordination number of H₂O and CO₃²⁻ for Y(III), and the approximate evolution process of Y(III) complex ions is obtained. According to the simulation results of b, c, d, e, and f, two O atoms in the carbonate are arranged in the first shell, and the one remaining O atom is ordered in the second shell far away from Y(III), and their relative distance is stable. Considering the interaction distance of the coordination bond in the range of the overall chemical bond (r < 2.6 Å),⁴¹ we believe that two O_C(n) atoms in the carbonate solution are combined with Y(III) in a coordinated manner. However, far away from Y(III), the O_C(f) atom is combined with the C atom, the relative distance between r_Y–OC(n) and r_Y–OC(f) is also stable, and the characteristic peak positions of O_C(n) and O_C(f) are at a consistent ratio. According to the actual geometry of CO₃²⁻, in the solution, CO₃²⁻ takes the O_C(n) atoms to the Y atom and combines with Y(III) via bidentate coordination.

In order to demonstrate CO₃²⁻ for the Y(III) coordination mode, we assume that the coordination mode of each CO₃²⁻ to Y atom is bidentate form, r_Y–OC(n) is the distance between the Y(III) short-range Y atom and the CO₃²⁻ first-coordination-layer O atom, and r_Y–OC(f) is the distance between Y(III) and O_C(n). The distance between the other atoms is defined as shown in Fig. 6. In these MD simulation studies, CO₃²⁻ is a rigid particle, where r_C–O and r_C=O are 1.269 Å, ∠OCO = 120° according to the principle of minimum energy, and ∠YCO = 60°.

Fig. 6 Schematic diagram of coordination of CO₃²⁻ to Y(III).

On the basis of these known conditions, we infer the relationship between r_Y–OC(n) and r_Y–OC(f)via the following reasoning:

x₁ = sin 60°r_C–O (3)

x₂ = cos 60°r_C–O (4)

r_Y–C = x₂ + x₃ (6)

r_Y–OC(f) = x₂ + x₃ + r_C=O (7)

By substituting the values of r_C–O and r_C=O to eqn (8), we calculate r_Y–OC(f).

According to the definition of Fig. 6, the following relationship can also be determined: N_OC(n) : N_OC(f) = 2 : 1.

N_OC(f) = Cn_Y–O(r < 4.0 Å) − Cn_Y–O (r < 4.0 Å) (10)

According to the calculation results in Table 4, N_OC(n) : N_OC(f) are all 2 : 1 among models b, c, d, e, and f. The relationship between their r_Y–OC(n) and r_Y–OC(f) agrees well with eqn (9). Therefore, our assumptions are true. Regardless of CO₃²⁻ concentration, CO₃²⁻ and Y(III) are combined in the form of bidentate coordination. In Na₂CO₃ solution systems of different concentrations, only the Y(III) ion CO₃²⁻ average coordination number and hydration number change. With increasing CO₃²⁻, the coordination number of the O atom to Y(III) changes. The O_C(n) atoms are directly coordinated with Y(III), that is, Cn_Y–O (r < 2.6 Å) is the average coordination number of O.³³ When CO₃²⁻ concentration is 0.4–0.8 mol L⁻¹, the coordination number of O atom to Y(III) is 2, the coordination number of H₂O molecule to Y(III) is 3; while the concentration of CO₃²⁻ is 1.2–2.0 mol L⁻¹, the coordination number of O atom to Y(III) is 4, and the coordination number of H₂O molecule to Y(III) is 2. Therefore, the total coordination number of Y(III) in YCl₃ aqueous solution is 3, and the total coordination number of Y(III) is 5 when CO₃²⁻ concentration is 0.4–0.8 mol L⁻¹, and the total coordination number of Y(III) is 6 when the concentration of CO₃²⁻ increases to 1.2–2.0 mol L⁻¹. In a high-concentration (CO₃²⁻ > 1.2 mol L⁻¹) carbonate solution, Y(III) mainly exists in the bidentate coordination form of [Y(CO₃)₂·2H₂O]⁻. This finding is consistent with the conclusion by Jeanvoine¹⁸ that the carbonate of Ln(III) in aqueous solution shows bidentate coordination.

3.3 UV-vis full wavelength scanning absorption spectroscopy

Before the results are analyzed and discussed further, the characteristics of some experimental results is explained. A0, which is the UV-vis absorbance curve, is determined by adding 1 mL of YCl₃ to 25 mL of pure water. A0.4, A0.8, A1.2, A1.6, and A2.0 are the UV-vis absorbance curves determined by adding 1 mL of YCl₃ to 25 mL solutions of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃, respectively. After mixing, the mixture is centrifuged, and the supernatant is subjected to UV-vis full-wavelength scanning to obtain the corresponding UV-vis absorbance curves. B0.4, B0.8, B1.2, B1.6, and B2.0 are the UV-vis absorbance curves determined by adding 1 mL of pure water to 25 mL solutions of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃, respectively. After mixing, the UV-vis full-wavelength scanning absorption curve is measured by the same method as a blank control spectrum. The results of UV-vis full-wavelength scanning analysis are shown in Fig. 7.

Fig. 7 UV-vis spectrum of YCl₃–H₂O solution system (a0) and YCl₃–Na₂CO₃ solution system (b0, c0, d0, e0, f0). (A0.4, 0.8, 1.2, 1.6, and A2.0 are determined by adding 1 mL of YCl₃ to 25 mL solutions of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃, respectively. After mixing, the mixture is centrifuged, and the supernatant is subjected to UV-vis full-wavelength scanning to obtain the corresponding UV-vis absorbance curves. B0.4, 0.8, 1.2, 1.6, and 2.0 are determined by adding 1 mL of pure water to 25 mL solutions of 0.4, 0.8, 1.2, 1.6, and 2.0 mol L⁻¹ Na₂CO₃.)

Fig. 7(a0) illustrates that the UV-vis characteristic absorption peak of the YCl₃–H₂O solution system is λ = 203 nm. Fig. 7(b0 and c0) show the UV absorption peaks of the YCl₃–Na₂CO₃ solutions (0.4 and 0.8 mol L⁻¹ CO₃²⁻) and the corresponding blank control spectrum. The characteristic absorption wavelengths of A0.4 and B0.4 is close, A0.8 and B0.8 is also close. But the absorbance of A0.4 at λ = 224 nm is significantly stronger than that of B0.4, and the absorbance of A0.8 at λ = 225 nm is significantly stronger than B0.8 at peak position. Fig. 7(d0, e0 and f0) show the UV absorption peaks of the YCl₃–Na₂CO₃ solutions (1.2, 1.6, and 2.0 mol L⁻¹ CO₃²⁻) and the corresponding blank control spectrum. The UV-vis characteristic absorption peaks at λ = 217–218 nm for A1.2, A1.6, and A2.0 and their peak absorbances are remarkably stronger than those of the corresponding B1.2, B1.6, and B2.0 systems. The results of MD calculation should be verified via quantum chemical calculation to further analyze and demonstrate the Y(III) coordination morphology in the YCl₃–Na₂CO₃ solution systems of different CO₃²⁻.

3.4 DFT calculation results

3.4.1 Ion structure and theoretical calculation UV spectroscopy. Based on the complex ion composition obtained by the MD simulations and considering the first shell of Y(III) in the Y–carbonate complex solution system, we optimize the Y(III) complex structure and calculate the UV spectrum via the DFT method. Fig. 8 shows the quantum chemical calculation UV spectrum and optimized structure of [Y·3H₂O]³⁺, [YCO₃·3H₂O]⁺, and [Y(CO₃)₂·2H₂O]⁻.

Fig. 8 [Y·3H₂O]³⁺, [YCO₃·3H₂O]⁺ and [Y(CO₃)₂·2H₂O]⁻ calculated UV spectrum and structure (S1, S2, S3).

In this section, the absorption curve measured by UV-vis full-wavelength scanning is obtained by adding 1 mL of YCl₃ to several Na₂CO₃ solutions; in fact, it is the overlap of the absorption curve generated by Y(III) complex ions (including [Y·3H₂O]³⁺, [YCO₃·3H₂O]⁺, and [Y(CO₃)₂·2H₂O]⁻) and the absorption bands of the solvent (Na₂CO₃, H₂O). Although the solvent itself has its own UV absorption band, the characteristic peak wavelength of the UV-vis full-wavelength scan of Fig. 7(a0) is close to the wavelength of the calculated UV spectral characteristic absorption of Fig. 8(S1). Moreover, the UV-vis characteristic absorption peaks at λ = 217–218 nm for A1.2, A1.6, and A2.0 are significantly stronger than those for the corresponding B1.2, B1.6, and B2.0 curves. This finding verifies that Y(III) exists in the form of [Y·3H₂O]³⁺ in the YCl₃ aqueous solution. The peak positions of the UV-vis full-wavelength scans of Fig. 7(d0, e0 and f0) are consistent with the absorption peak position of the UV spectrum characteristic of Fig. 8(S3), that is, the theoretically calculated UV spectral absorption characteristic peak position and the experiment measured absorption peaks of the UV-vis full-wavelength scan are consistent. Thus, Y(III) is a 6-coordinated complex [Y(CO₃)₂·2H₂O]⁻ in the 1.2–2.0 mol L⁻¹ Na₂CO₃ solution. In Fig. 7, the results of UV-vis full-wavelength scanning analysis of B0.4, B0.8, B1.2, B1.6, and B2.0 show a certain absorption band in the blank control sample, which means the Na₂CO₃ solution itself presents some UV absorption, and the absorption curves of different concentrations of Na₂CO₃ solution are basically the same.

However, comparison between the UV test results of Fig. 7(b0 and c0) and the calculated spectra of Fig. 8(S2) shows that the absorption band of Na₂CO₃ solution considerably influences the absorption intensity and absorption peak wavelength of A0.4 and A0.8. That is, in the Na₂CO₃ solution with low CO₃²⁻ concentration, the UV absorption band of the solvent itself greatly affects the peak position and absorption intensity of the UV absorption curve of [YCO₃·3H₂O]⁺. It can be seen from Fig. 8(S2) and Fig. 7(b0, c0), the characteristic absorption peak wavelength of [YCO₃·3H₂O]⁺ is λ = 203.7 nm, which is not seen in the experimental analysis of the low-CO₃²⁻-concentration Na₂CO₃ solutions (in absorbance curve of A0.4 and A0.8), because the characteristic peak of the Y–O_w bond at λ = 203.7 nm is not revealed in such solutions. The absorption bands B0.4 and B0.8 of the Na₂CO₃ solution itself are at a low absorption level of λ = 203.7 nm. Therefore, the UV spectrum obtained by the experimental test (A0.4 and A0.8) is the result of superimposing the UV absorption band of [YCO₃·3H₂O]⁺ as the solute in the solution and the absorption band of the Na₂CO₃ solution.⁴² This process results in the theoretical UV absorption peak of λ = 203.7 nm, which is not reflected in the experimental analysis of the UV spectrum. In Fig. 8(S2), the characteristic absorption band of [YCO₃·3H₂O]⁺ at λ = 220–240 nm, which resembles a shoulder peak, is enhanced in the experimental UV spectrum.⁴³ Thereby, the absorption band of [YCO₃·3H₂O]⁺ calculated UV spectrum (Fig. 8(S2)) at λ = 220–240 nm, is exhibited in Fig. 7(b0, c0) at λ = 224–225 nm.

We analyze the bond length information of Y(III) complex ions by DFT geometry optimization. The structural results of carbonate complex ion S1, S2, S3 complexes obtained from DFT calculation are shown in Table 5. Here, Y–O_C(n) represents the distance between the Y atom and the nearer O atom in the carbonate, and Y–O_W represents the distance between the Y atom and the O atom in the water molecule.

Table 5 Structural results for S1, S2, S3 complexes in bulk water obtained from DFT calculations

Structure	Complex ion composition	Y–OC(n) (Å)	Y–OW (Å)
S1	[Y·3H2O]^3+	—	2.24
S2	[YCO3·3H2O]^+	2.17	2.32
S3	[Y(CO3)2·2H2O]^−	2.25	2.43

---

In the RDF analysis of the MD simulation, the first shell of the O atoms interacting with Y(III) and the first shell of H₂O are analyzed, that is, Y(III) and each ligand the distance between them. Comparison of the values in Tables 4 and 5 shows that in model a, r_Y–H2O is close to the value of Y–O_W; in models b and c, r_Y–Oc(n) and the values of Y–O_C(n) are consistent, and r_Y–H2O is close to the value of Y–O_C(n); in models d, e, and f, r_Y–H2O is close to the value of Y–O_C(n), and the values of r_Y–H2O and Y–O_W are close. The comparison also shows that the coordination around Y(III) is CO₃²⁻. When the coordination number of CO₃²⁻ is increased from 1 to 2, the length of Y–O_W in the first shell of the H₂O molecule increases. Therefore, the distance between Y(III) and O_C(n) obtained by MD simulation is consistent with the bond length Y–O_C(n) calculated by DFT. The distance r_Y–H2O of the first shell of the H₂O molecule is also consistent with the bond length of Y–O_W calculated by DFT. In summary, the DFT calculation results and the experimental UV-vis analysis are well-validated by the MD simulation results.

According to the results of the MD simulation, dissolve experiment and UV-vis full-wavelength scanning experiment, three complex ions of Y(III) can be gradually transformed into each other increases with the concentration of CO₃²⁻ increasing, the process of which is shown in Fig. 9.

Fig. 9 Schematic diagram of morphological transformation of S1, S2 and S3 complex ions in YCl₃–Na₂CO₃ solution system.

With increasing of the CO₃²⁻ concentration, the following dynamic balance exists in the YCl₃–Na₂CO₃ solution system:

YCl₃ ⇌ Y³⁺ + 3Cl⁻ (11)

Y³⁺ + 3H₂O ⇌ [Y·3H₂O]³⁺ (12)

[Y·3H₂O]³⁺ + CO₃²⁻ ⇌ [YCO₃·3H₂O]⁺ (13)

[YCO₃·3H₂O]⁺ + CO₃²⁻ ⇌ [Y(CO₃)₂·2H₂O]⁻ + H₂O (14)

Fig. 9 illustrates that for [Y·3H₂O]³⁺ has one coordination point, with the concentration of CO₃²⁻ increasing, preferentially combines one CO₃²⁻ at the coordination point of [Y·3H₂O]³⁺, and converts to [YCO₃·3H₂O]⁺ complex ions, while the concentration of CO₃²⁻ further increase, there will be a second CO₃²⁻ to replace one of the H₂O molecules in the first shell, as in the reaction of (14) process, [YCO₃·3H₂O]⁺ transforms to [Y(CO₃)₂·2H₂O]⁻. Under different CO₃²⁻ concentration of the YCl₃–Na₂CO₃ solution system, the complex ion morphology of Y(III) is consistent with the conclusion reported by Jeanvoine.¹⁸ This result is also in line with Philippini's¹³ conclusion that Ln(III) coordinates with additional CO₃²⁻ ions with increasing concentration of CO₃²⁻.

In our analysis of the instantaneous saturated dissolution of Y(III) in Na₂CO₃ solution, we believe that, given the coordination effect of a high concentration of CO₃²⁻ ions on Y(III) along with the CO₃²⁻ concentration increase, reaction eqn (13) and (14) are promoted in the positive direction, thereby increasing the instantaneous saturated solubility of Y(III). Therefore, according to the Y(III) dissolution experiment results, MD simulation studies of different concentrations of the YCl₃–Na₂CO₃ real solution system are performed. By combining the DFT calculation results and experimental UV-vis analysis, we can accurately resolve the microstructure and behavior of CO₃²⁻ coordination and hydration for Y(III). We reveal that the large mother liquor ion loss in rare-earth carbonate precipitation is caused by CO₃²⁻ coordination and hydration; consequently, REE is dissolved in the carbonate solution in the form of complex ions, resulting in incomplete precipitation of rare-earth.⁴⁴

4. Conclusions

In this study, a Y(III) dissolution experiment in Na₂CO₃ solution was used to determine the instantaneous saturated solubility of Y(III) in Na₂CO₃ solution, and the macroscopic dissolution law of Y(III) with CO₃²⁻ concentration was obtained. The RDF and coordination number of Y(III) by CO₃²⁻ and H₂O were systematically analyzed using MD simulations to obtain the complex ion form of Y(III) in carbonate solution. DFT was adopted to geometrically optimize and calculate the UV spectrum of the Y(III) complex ions. Findings are then compared with the experimentally determined UV-vis spectra and analyzed to further verify the reliability of the MD simulation results. The main conclusions are as follows.

(1) YCl₃ is added to Na₂CO₃ solution with instantaneous saturated dissolution, and a coordination effect on Y(III) is observed. The higher the CO₃²⁻ concentration, the easier the conversion to [Y(CO₃)₂·2H₂O]⁻ and the higher the instantaneous saturation dissolution of Y(III). When the initial concentration of Na₂CO₃ is 2 M, the instantaneous saturated dissolution of Y(III) in the carbonate reaches 0.0909 mol L⁻¹.

(2) In the aqueous solution of YCl₃, Y(III) exists in the form of [Y·3H₂O]³⁺. Each CO₃²⁻ with O atoms is close to the Y atom and combines via bidentate coordination. In 0–0.8 mol L⁻¹ CO₃²⁻ aqueous solutions, Y(III) is mainly 5-coordinated [YCO₃·3H₂O]⁺; with increasing CO₃²⁻ concentration, [YCO₃·3H₂O]⁺ is converted into 6-coordinated [Y(CO₃)₂·2H₂O]⁻.

(3) When Y(III) is added to aqueous solutions, hydrated [Y·3H₂O]³⁺ ions are preferentially formed, and [Y·3H₂O]³⁺ exists as a coordination point. When CO₃²⁻ ≤ 0.8 mol L⁻¹, one CO₃²⁻ passes through the second shell of the H₂O molecules from the outer layer and enters the first shell to coordinate with Y(III). In this case, the effect of H₂O molecule crowding is not noticeable. As CO₃²⁻ is further increased, the second CO₃²⁻ passes through the second shell and enters the first shell of Y(III), thereby replacing one of three hydrated H₂O molecules in the form of a substitution reaction and initiating a bidentate coordination bond with Y(III). In this case, short-range CO₃²⁻ space crowding causes the H₂O molecule number around the short range of Y(III) to be reduced by 1–2.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work is supported by the National Natural Science Foundation of China (NSFC) (No. 51774155 and 51564017), and thanks to the outstanding doctoral training program of Jiangxi University of Science and Technology (No. 3105500030).

References

1. M. K. Jha, A. Kumari, R. Panda, J. R. Kumar, K. Yoo and J. Y. Lee, Hydrometallurgy, 2016, 165, 2–26.
2. Y. Xiong, W. Kuang, J. Zhao and H. Z. Liu, Sep. Purif. Technol., 2017, 179, 349.
3. M. Vogit, J. D. Rodriguez, B. Vallina, L. G. Benning and E. H. Oelkers, Chem. Geol., 2016, 43, 70.
4. R. Safarbali, M. R. Yaftian and A. Zamani, J. Rare Earths, 2016, 34, 91–98.
5. X. L. Shen, M. M. Xing, Y. Tian, Y. Fu, Y. Peng and X. X. Luo, J. Rare Earths, 2016, 34(5), 458.
6. Y. Xiao, Z. Long, X. Huang, Z. Feng, D. Cui and L. Wang, J. Rare Earths, 2013, 31(5), 512.
7. T. Taketatsu, Anal. Chim. Acta, 1965, 32, 40.
8. T. Taketatsu, Bull. Chem. Soc. Jpn., 1963, 36(5), 549.
9. E. L. Tran, O. Klein-Bendavid, N. Teutsch and N. Weisbrod, Environ. Sci. Technol., 2015, 49(22), 13275–13282.
10. C. W. Noack, D. A. Dzombak and A. K. Karamalidis, Environ. Sci. Technol., 2014, 48(8), 4317.
11. K. H. Johannesson and W. B. Lyons, Limnol. Oceanogr., 1994, 39(5), 1141–1154.
12. M. E. Vasconcellos, S. M. R. Rocha, W. R. Pedreira, C. A. Queiroz and A. Abrao, J. Alloys Compd., 2008, 451(1), 426–428.
13. V. Philippini, T. Vercouter, J. Aupiais, S. Topin, C. Ambard, A. Chausse and P. Vitorge, Electrophoresis, 2008, 29(10), 2041–2050.
14. V. Philippini, T. Vercouter and P. Vitorge, J. Solution Chem., 2010, 39(6), 747–769.
15. P. E. Reiller, T. Vercouter, L. Duro and C. Ekberg, Appl. Geochem., 2012, 27(2), 414–426.
16. T. Vercouter, P. Vitorge, N. Trigoulet, E. Giffaut and C. Moulin, New J. Chem., 2005, 29(4), 544.
17. M. Maloubier, H. Michel, P. L. Solari, P. Moisy, M. A. Tribalat, F. R. Oberhaensli, M. Y. D. Bottein, O. P. Thomas, M. Monfort, C. Moulin and C. D. Auwer, Dalton Trans., 2015, 44(47), 20584–20596.
18. Y. Jeanvoine, P. Miro, F. Martelli, C. J. Cramer and R. Spezia, Phys. Chem. Chem. Phys., 2012, 14, 14822–14831.
19. L. Rao, D. Rai, A. R. Felmy and C. F. Novak, Actinide Speciation in High Ionic Strength Media, 1999, pp. 153–169.
20. P. Thakur, Y. Xiong and M. Borkowski, Chem. Geol., 2015, 413, 7.
21. I. Perssin, E. D. Risberg, P. D'Angelo, S. De Panfilis, M. Sandstrm and A. Abbasi, Inorg. Chem., 2007, 46, 7742–7748.
22. O. Seo, J. Kim, A. Tayal, C. Song and L. S. R. Kumara, RSC Adv., 2019, 9, 21311–21317.
23. F. Martelli, Y. Jeanvoine, T. Vercouter, C. Beuchat, R. Vuilleumier and R. Spezia, Phys. Chem. Chem. Phys., 2014, 16, 3693–3705.
24. Z. F. Tong, Y. B. Xie and Y. H. Zhang, J. Mol. Liq., 2018, 259, 65–75.
25. J. Tang and K. H. Johannesson, Chem. Geol., 2010, 279(3), 120.
26. Y. Cheng, Z. Yang, J. Liao, J. Qiu, Z. Song and Y. Yang, J. Rare Earths, 2015, 33, 599–603.
27. Y. M. Yang, X. L. Zhang, L. Li, T. M. Wei and K. Z. Li, ACS Omega, 2019, 4, 9160–9168.
28. M. T. Alhaffar, S. A. Umoren, I. B. Obot and S. A. Ali, RSC Adv., 2018, 8, 1764.
29. P. Liang, C. Du, G. Li, X. Chen and Z. Liang, Int. J. Miner., Metall. Mater., 2009, 16(4), 407–413.
30. A. Saxena and D. Sept, J. Chem. Theory Comput., 2013, 9, 3538.
31. G. Hu, Y. Sun, Y. Xie, S. Wu, X. Zhang, J. Zhuang, C. Hu, B. Lei and Y. Liu, ACS Appl. Mater. Interfaces, 2019, 11, 6634–6643.
32. E. E. El-katori and A. S. Abousdem, RSC Adv., 2019, 9, 20760–20777.
33. M. Zhou, T. Jiang, S. T. Yang and X. X. Xue, Int. J. Miner., Metall. Mater., 2015, 22(9), 917.
34. H. Tomono, H. Nada, F. Zhu, T. Sakamoto, T. Nishimura and T. Kato, J. Phys. Chem. B, 2013, 117, 14849.
35. D. Liu, J. Zhang, X. Xue, G. Wang, K. Jiang and Z. Liu, Int. J. Miner., Metall. Mater., 2016, 23(6), 618.
36. R. Hayes, G. G. Warr and R. Atlin, Chem. Rev., 2015, 115(13), 6657.
37. H. Zhang, Y. Jiang, H. Yan, C. Yin, T. Tan and D. V. Spoel, J. Phys. Chem. Lett., 2017, 8, 2705–2712.
38. A. Bankura, V. Carnevale and M. L. Klein, J. Chem. Phys., 2013, 138, 014501.
39. V. Migliorati, A. Serva, F. M. Terenzio and P. D'Angelo, Inorg. Chem., 2017, 56, 6214–6224.
40. V. Migliorati, A. Serva, F. Sessa, A. Lapi and P. D'Angelo, J. Phys. Chem. B, 2018, 122, 2779–2791.
41. I. C. Yeh and G. Hummer, J. Phys. Chem. B, 2004, 108(40), 15873.
42. S. C. Hess, F. A. Permatasari, H. Fukazawa, E. M. Schneider, R. Balgis, T. Ogi, K. Okuyama and W. J. Stark, J. Mater. Chem. A, 2017, 5, 5187–5194.
43. S. Soumya, V. N. Sheemol, P. Amba, A. P. Mohamed and A. Solaiappan, Sol. Energy Mater. Sol. Cells, 2018, 174, 554–564.
44. F. Lai, G. Gao, L. Huang, Y. Xiao, R. Yang and K. Li, Hydrometallurgy, 2018, 179, 25–35.

---