Speciation analysis the complexation of uranyl nitrate with tri-n-butyl phosphate in supercritical CO₂

Abstract

The complexation of solid uranyl nitrate with tri-n-butyl phosphate (TBP) in supercritical CO₂ is quite different from that of a liquid–liquid extraction system because fewer water molecules are involved. Here, the complexation mechanism was investigated by molecular dynamics simulation, emphasising on speciation distribution analysis. In the anhydrous uranyl nitrate system, poly-core uranyl-TBP species [UO₂(NO₃)₂]₂·3TBP and [UO₂(NO₃)₂]₃·3TBP were formed in addition to the predominant [UO₂(NO₃)₂]·1TBP and [UO₂(NO₃)₂]·2TBP species. The poly-core species was mainly constructed via the linkage of U=O⋯U contributed by pre-developed [UO₂(NO₃)₂]·1TBP species. However, in the hydrated uranyl nitrate system, TBP·[UO₂(NO₃)₂]·H₂O species form, preventing the formation of the poly-core species. The complexation developed differently depending on the TBP to the uranyl nitrate ratio, the solute densities and the participation of water. It suggested that the kinetically favoring species would gradually convert into the thermodynamically stable species [UO₂(NO₃)₂]·2TBP by ligand exchange.

---

1 Introduction

The spent nuclear fuel (SNF) discharged from nuclear power plants contains fissionable materials, such as uranium and plutonium, and other highly radioactive elements. Recycling uranium and plutonium from SNF is crucial for the sustainable development of nuclear power. This is achieved by the well-known hydrometallurgy process, Plutonium and Uranium Coextraction Process (PUREX), in which uranium and plutonium are extracted by tri-n-butyl phosphate (TBP) diluted in kerosene organic phase once SNF is dissolved in nitric acid. This leaves other radionuclides such as ²⁴¹Am, ^243/244Cm, ⁹⁰Sr, ¹³⁷Cs in the raffinate, thus generating a significant amount of secondary radioactive liquid waste. Supercritical fluid has unique properties, such as gas-like diffusion abilities and liquid-like solubilities. Using supercritical CO₂ for actinide extraction and separation has been extensively studied, and the main motivation is to reduce the secondary radioactive liquid waste compared to the hydrometallurgy process.^1,2 The extraction of uranyl nitrate with TBP into supercritical CO₂,³ or direct dissolution and extraction from uranium oxide using supercritical CO₂ containing the TBP–HNO₃ complex has been demonstrated.⁴ Conjunction with other newly investigated media such as ionic liquid or deep eluent solution was also explored.^5,6 Although utilizing supercritical fluid to extract actinides, particularly uranium, has been verified in diversified manners or application scenarios, the complexation mechanism of uranyl with TBP in supercritical CO₂ is still vague and may be quite different from a liquid–liquid extraction system. It has also been found that solid lanthanide nitrate can be easily extracted into supercritical CO₂ phase containing TBP,⁷ while lanthanide ions often remain in aqueous phase in a liquid–liquid extraction system.

Due to the high pressure of the supercritical fluid phase, it is difficult to explore the complexation reaction with analytical instrumentation. Molecular dynamic modelling is an alternative method, which has been actively developed for liquid–liquid extraction simulation.^8–11 In the liquid–liquid solvent extraction system, it is predicted that uranyl ion complexes with TBP form UO₂·NO₃·4TBP, UO₂·5TBP and UO₂·NO₃·3TBP·HNO₃ in the interface or in the organic phase.¹² Guilbaud et al.¹³ reported the aggregation of UO₂(NO₃)₂(TBP)₂ at a high uranium concentration. As for a supercritical CO₂ extraction system, only Wipff et al.^14,15 investigated the complexation of uranyl nitrate with TBP at the interface of the supercritical CO₂/aqueous phase. Since a large amount of H₂O molecules are involved in the simulation, the results are more comparable to the liquid–liquid extraction system. Thus, the complexation of uranyl nitrate with TBP in supercritical CO₂ with less water participation is still needed. Also, investigating the distribution of different species evolving with simulation time may be useful to track how the experimental conditions affect the results; however, no study has been reported yet because the number of uranyl nitrate molecules considered is usually limited. Therefore, we aim to investigate the complexation of anhydrous and hydrate uranyl nitrate with TBP in supercritical CO₂ and describe the species distribution evolved with the simulation procedure.

2 Computational methods

The force fields for uranyl nitrate, TBP, and CO₂ were carefully considered. The charge of O_P=O in the TBP molecule plays a very important role in determining the interaction strength of TBP with uranyl nitrate, but its value varies largely among the studies^8,9,11 from −0.54 e to −0.87 e. The original structure of TBP was adopted from ref. 16, and further optimized at the B3LYP/6-311G(d) level as implemented in the Gaussian 16 (ref. 17) code, with empirical dispersion corrections (GD3BJ).¹⁸ The IEFPCM solvation model was used, with a dielectric constant of 1.5, which is consistent with the experimental value of supercritical CO₂ at 323 K and 20 MPa.¹⁹ Also, by default, the solute cavity is defined as the UFF force field radii of the atoms with a scale factor of 1.1. After the optimization, the RESP charge of TBP was calculated by the Multiwfn software²⁰ version 3.7. The charge was 1.259 e and −0.670 e for P_P=O and O_P=O, respectively, close to 1.34 e for P_P=O and −0.707 e for O_P=O atoms reported by Mu et al.²¹ Other parameters for bonds, angles and non-bonded interactions of TBP were generated by Antechamber tools in Amber 18,²² by applying the build-in GAFF2 force field. The force field for CO₂ in the supercritical state was evaluated, as shown in Table S1.† The parameters developed by Zhu et al.²³ was used because the density was calculated to be 0.787 g cm⁻³, which is very close to that of the experimental value of 0.784 g cm⁻³ at 323 K and 20 MPa. The uranyl nitrate was treated as a bonded model to maintain the extracted species charge neutral.¹⁴ However, there are no exclusive force fields for uranyl nitrate, which is critical for the geometry of the uranyl complex. Here, the uranyl nitrate force field was generated through the Seminario²⁴ method by using MCPB.py.²⁵ The complete parameters for uranyl nitrate provide constrain to the coordination sites of uranyl ions, which ensures a complex a reasonable geometry. The established force field parameters of uranyl nitrate and TBP (Scheme 1) are listed in supplement information.

Scheme 1 The atom charge of UO₂(NO₃)₂ and TBP.

The simulation systems were constructed using the Packmol²⁶ software, as listed in Table S2.† The length of the simulation cubic box was either 98 or 150 Å. It was reported that the [UO₂(NO₃)₂]·2TBP complex had high solubility in supercritical CO₂,²⁷ which could be up to 0.4 mol L⁻¹. In this study, the concentration of uranyl nitrate was either 0.04 or 0.17 mol L⁻¹. In order to statistically calculate the number of species, up to 100 uranyl nitrates were included in the simulations. All the systems were first relaxed using steepest descent and conjugate gradient minimization algorithm to eliminate any bad contacts. Then, they were heated at gradually increased temperature up to 700 K in an NVT ensemble. At the production stage, an Langevin thermostat was employed for temperature control at 323 K and 20 MPa. Periodic boundary was imposed on the system during the calculation of non-bonded interactions. A cutoff radius was set at 10 Å for nonbonded interactions. The MD trajectories and velocities were calculated with a time step of 1 fs. The SHAKE constraints²⁸ was used for all the hydrogen atoms. Auto-image was employed when writing the coordinates to the restart and trajectory files. The production period of molecular dynamics lasted for 20 ns, and the radial distribution function of uranyl to TBP remained almost unchanged after 20 ns. As the system was relatively large, the GPU code^29,30 in Amber 18 was used. After the simulation, the trajectory was visualized by the VMD³¹ software.

The structures of complexes taken from the trajectories were optimized using B3LYP function. U atom relativistic effects were dealt using the quasi-relativistic effective core potentials (RECPs).³² The adopted small core RECPs had 60 electrons in the core of U (MWB60), and the affiliated segmental basis sets were applied for the valence electrons of U, while the def2-tzvp³³ basis set was used for N, H, O and P.

To describe the distribution of species with simulation, for example the uranyl-TBP species, the RDF of uranium with O_P=O was first calculated, and a cutoff of 3.5 Å was used to determine whether the ligand TBP complexed with uranyl nitrate. The coordinates in different simulation times were extracted to PDB files from the simulation trajectories. And then the residue IDs of the paired uranyl nitrates and the corresponding TBP molecules were generated by using “native contact” command in the cpptraj tool³⁴ from each PDB file. Lastly, discrete frequency counts of the obtained residue IDs of uranyl nitrate was performed twice using a shell script to obtain the number of species at different times. As for the species when uranyl nitrate complexed with TBP and water simultaneously, the contact analysis of uranyl with TBP and water were made separately, and then, the residue IDs of uranyl nitrate in the above output files were cross-indexed to obtain the number of species.

3 Results and discussion

3.1 Simulation without H₂O

Compared with liquid–liquid solvent extraction, typically there is less water in the supercritical CO₂ extraction system. The complexation of uranyl nitrate with TBP without water was first investigated. The RDF of U⋯O_P=O reached equilibrium after 20 ns, suggesting the simulation duration was adequate for the present system. Also, the cutoff of the complexation reaction of uranyl nitrate with TBP was 3.5 Å, as shown in Fig. 1.

Fig. 1 The RDF of U⋯O_P=O evolved with time in the system I.

Surprisingly, there were four kinds of uranyl-TBP species identified, [UO₂(NO₃)₂]·1TBP, [UO₂(NO₃)₂]·2TBP, [UO₂(NO₃)₂]₂·3TBP and [UO₂(NO₃)₂]₃·3TBP, denoted as 1 : 1, 1 : 2, 2 : 3 and 3 : 3, as illustrated in Fig. 2. The most predominant species were 1 : 1 and 1 : 2. The structure of UO₂(NO₃)₂·2TBP was similar to that reported in the literature,³⁵ in which two TBP molecules coordinated with uranyl moiety O=U=O at the equatorial plane. For the 1 : 1 species, TBP also coordinated to the uranium from the equatorial plane of the uranyl moiety, and the angle of N⋯U⋯N was approximately 153°, bending towards the side without TBP. The electrostatic potential (ESP) surfaces of the 1 : 1 and 1 : 2 species were also generated, as shown in Fig. 2. On the TBP-free side of the 1 : 1 species, more positive charges leaked out of uranium, and on the O_NO3 atom of the 1 : 2 species, there was more negative charge. The 1 : 1 and 1 : 2 species can interact with each other via electrostatic interactions, thus acting as the building blocks of the other species. The structure of 2 : 3 species was unique, with the two O=U=O moieties being linked by the cation–cation interactions.^36,37 Two types of 2 : 3 species were identified, of which type I had two N⋯U⋯N axes parallel to each other, while in type II, the N⋯U⋯N axis was arranged in a perpendicular manner. The structure of 3 : 3 species was constructed by the inner interaction of three O=U=O, with three N⋯U⋯N axes paralleled to each other. Also, the outer of three O=U=O clusters was covered by three TBP molecules. The parameters of the four identified species are listed in Table 1. In most cases, the length of the P=O bond in TBP was approximately 1.51 Å, but for type I of 2 : 3 species, there were two P=O bonds of 2.34 Å. The bond length of U⋯O_P=O was measured from 2.30 to 2.38 Å, except for the 1 : 1 species. It was noticed that some O=U=O bonds were enlarged because of the cation–cation interaction. Usually, the length was about 1.76 Å, while in the 2 : 3 species, the bond of U=O was enlarged to 1.82 Å. The average distances of U⋯O_NO3 were 2.47, 2.54, 2.51 and 2.48 Å for 1 : 1, 1 : 2, 2 : 3 and 3 : 3 species, respectively, among which, the 1 : 2 species had the longest U⋯O_NO3 distance. In the original uranyl nitrate, the distance of U⋯O_NO3 was 2.42 Å, and that of U=O was 1.75 Å; thus, the bond length would be enlarged upon complexing with TBP, and the longer distance of U⋯O_NO3 may suggest stability of the species.

Fig. 2 The optimized structures of uranyl species formed in supercritical CO₂ (the hydrogen atoms were hidden for clarity).

Table 1 The structural parameters of uranyl-TBP species

Species	Bond distance/Å					Angle/degree
	P=O	U⋯OP=O	O=U=O	U–ONO3	U⋯U	N⋯U⋯N
[UO2(NO3)2]·1TBP	1.52	2.29	1.77, 1.76	2.45, 2.47, 2.46, 2.48	—	153.3
[UO2(NO3)2]·2TBP	1.51, 1.51	2.37, 2.37	1.77, 1.77	2.53, 2.54, 2.53, 2.54	—	178.9
[UO2(NO3)2]2·3TBP type I	1.51, 2.34, 2.34	2.36, 2.34, 2.34	1.77, 1.77	2.50, 2.51, 2.51, 2.50	4.21	168.9, 176.3
			1.82, 1.77	2.50, 2.51, 2.52, 2.53
[UO2(NO3)2]2·3TBP type II	1.52, 1.52, 1.51	2.35, 2.34, 2.34	1.77, 1.81	2.51, 2.51, 2.51, 2.50	4.23	177.4, 174.1
			1.77, 1.76	2.51, 2.52, 2.52, 2.52
[UO2(NO3)2]3·3TBP	1.52, 1.52, 1.52	2.31, 2.31, 2.35	1.81, 1.77	2.48, 2.50, 2.47, 2.47	4.20, 4.21, 4.21	175.3, 173.2, 174.1
			1.76, 1.81	2.48, 2.51, 2.48, 2.46
			1.81, 1.76	2.48, 2.49, 2.49, 2.48

---

RDF_U⋯U was calculated from final 5 ns trajectories in different systems, and normalized by the number density of uranyl nitrate, as shown in Fig. 3. There were peaks at approximately 4.5 Å, coincided with the distance of U⋯U in the poly-core species in trajectories, thus indicating the probabilities of poly-core species. There were two factors affecting the formation of poly-core species, which were the ratio of TBP to uranyl nitrate (TBP : U), and the concentration of uranyl nitrate and TBP in the system. When the TBP : U ratio was low, the peak at 4.5 Å was high, indicating that a more poly-core complex formed. Even when the ratio of TBP : U was low, if the concentration of uranyl nitrate and TBP concentration was high, the formation of poly-core species also decreased dramatically. The diffusion constants of uranyl nitrate and TBP were larger for the system with lower solute concentration, as shown in Table 2; however, there was no direct clue to the formation of poly-core species from the diffusion constants.

Fig. 3 The RDF of uranium–uranium in systems I–IV.

Table 2 The diffusion constants of uranyl nitrate, TBP, CO₂

System	Number of molecules	Diffusion constants (×10^−5 cm^2 s^−1)
		Uranyl nitrate	TBP	CO2
I	100U + 200TBP + 10000CO2	0.4932	0.6024	7.1735
II	100U + 600TBP + 10000CO2	0.3193	0.5535	3.8410
III	100U + 200TBP + 40000CO2	0.7579	0.9831	8.7793
IV	100U + 600TBP + 40000CO2	0.5884	1.1091	7.7932

---

The speciation distribution is shown in Fig. 4, and the time 0 ns was the starting point of the production stage of simulation. It can be seen that, in different systems, the change of species with time was quite different. For the system with 600 TBP (system II and IV), the ultimate species was 1 : 2 species, and the condensed system IV formed 1 : 2 species more quickly. For the system without enough TBP, the species distribution was different. In the sparse system (system I), the 1 : 1 species increased first, and then decreased slowly, while for the condensed system (system III), the 1 : 1 species decreased continuously. Only system I had a large amount of 1 : 1 species after 20 ns, which was also easy to observe the poly-core uranyl species. This was also consistent with the RDF analysis. The value of RDF_U⋯U at 4.5 Å in system I was the highest among the systems. We retrieved the trajectories of the poly-core formation process. The 2 : 3 species was formed by the collision of 1 : 1 and 1 : 2 species. At first, the interaction was established by the attraction of O_NO3 in 1 : 2 species and U atom in 1 : 1 species, then a reconfiguration of the cluster and exchange of TBP molecules occurred, subsequently a stable 2 : 3 species formed. The 3 : 3 species was formed when three 1 : 1 species met almost in the same time. Therefore, the presence of 1 : 1 species was the key to the formation of poly-core species, and the side without TBP of 1 : 1 species had the ability to contact with other 1 : 1 species or 1 : 2 species. The system with a low TBP : U ratio and low TBP concentration produced more 1 : 1 species, which could further generate poly-core species.

Fig. 4 The speciation distribution in systems I to IV.

3.2 Simulation with H₂O

Anhydrous uranyl nitrate can be obtained by the reaction of U₃O₈ with N₂O₄ in a water-free container, whereas the readily available uranyl nitrate usually contains six water molecules.³⁸ Therefore, the system with six times of water molecules was examined with the TIP3P water model. The main species in systems V–VIII were [UO₂(NO₃)₂]·2H₂O, TBP·[UO₂(NO₃)₂]·H₂O, [UO₂(NO₃)₂]·2TBP and TBP·H₂O. Other species such as [UO₂(NO₃)₂]·H₂O, [UO₂(NO₃)₂]·1TBP were seldom observed, which only appeared at the beginning of production stage. Interestingly, the TBP-H₂O species could coordinate with TBP-uranyl species to obtain a [UO₂(NO₃)·H₂O]·2TBP species, as shown in Fig. 5.

Fig. 5 The species of TBP·[UO₂(NO₃)₂]·H₂O (left) and [UO₂(NO₃)·H₂O]·2TBP (right).

The peak at 4.5 Å of the RDF_U⋯U curve disappeared, as shown in Fig. 6. It was difficult to form the poly-core complex in the presence of water because the coordination of H₂O with uranium atoms would prevent the approaching of another uranium atom. The RDF of U⋯O_H2O was about 2.7 Å, and when the amount of TBP molecules was high, the peak intensity at 2.7 Å decreased, suggesting that more number of water molecules were replaced by TBP. This coincides with solvent extraction situation, that the complexation of uranyl with TBP is achieved by the replacement of the water molecules of the hydrated uranyl ions by TBP.

Fig. 6 The RDF of uranium–uranium when water present at last 5 ns.

The speciation distribution is shown in Fig. 7. The species TBP·H₂O fluctuated largely in all systems, which suggested that this species was less stable. Also, TBP·H₂O decreased gradually over time when TBP was insufficient, while it remained constant in the systems with 600 TBP. [UO₂(NO₃)₂]·2TBP species increased gradually, while [UO₂(NO₃)₂]·2H₂O and TBP·[UO₂(NO₃)₂]·H₂O increased first, and then decreased gradually, except for system V. Therefore, it is assumed that, at the beginning of the simulation, uranyl nitrate coordinated with water and TBP quickly at its two coordination sites, resulting in the formation of [UO₂(NO₃)₂]·2H₂O or TBP·[UO₂(NO₃)₂]·H₂O, rather than one site coordination species, [UO₂(NO₃)₂]·H₂O or TBP·[UO₂(NO₃)₂]. Moreover, excess water molecules complexed with TBP, forming TBP·H₂O species. Then, the thermodynamically stable species [UO₂(NO₃)₂]·2TBP developed gradually, and kinetically favoured species decreased accordingly. We also compared system V and VII. Because the TBP was not sufficient in system V, TBP·H₂O species decreased, suggesting that TBP·H₂O species would dissociate, and that the released TBP molecules were captured by uranyl nitrate or uranyl complex. There was a clear tendency that TBP·[UO₂(NO₃)₂]·H₂O species converted into [UO₂(NO₃)₂]·2TBP species. As shown by systems VI, VII and VIII, when the amount of [UO₂(NO₃)₂]·2H₂O was already very low, the number of TBP·[UO₂(NO₃)₂]·H₂O molecules decreased, while that of [UO₂(NO₃)₂]·2TBP increased. An exception was noticed in system V, in which the amount of TBP·[UO₂(NO₃)₂]·H₂O remained constant in the closure of simulation because of the low density of free TBP; this suggested that free TBP is needed for the conversion reaction to occur. As for the system with 600 TBP, the number of [UO₂(NO₃)₂]·2TBP at 20 ns was much higher, and even when there were free TBP, the TBP·[UO₂(NO₃)₂]·H₂O species still coexisted.

Fig. 7 The speciation distribution in systems V to VIII (U represents UO₂(NO₃)₂ in the figure).

3.3 Simulation at the interface

In experimental situation, the complexation would first take place in the interface of uranyl nitrate and supercritical CO₂ phases. Although it is difficult to simulate the system starting from the phase interface, we preliminary attempted to explore the situation. A system containing uranyl nitrate layer of 25 Å thickness was constructed using the Packmol software, in which 3000 uranyl nitrates were packed in the box of 150 × 150 × 25 Å, while TBP and CO₂ were packed in the rest of the cube of 150 Å length, as shown in Fig. 8. The solid phase present for a long time after a mixing period at 700 K was applied. A large proportion of TBP complexed with uranyl nitrate at the surface of uranyl nitrate layer, leading to a further decline in the concentration of TBP in supercritical CO₂. There were about 170 TBP and 60 uranyl nitrates in the bulk of supercritical CO₂ phase at 1 ns, whereas the number decreased to 70 and 50, respectively, at 20 ns. In the supercritical CO₂ phase, 1 : 1, 1 : 2 and 2 : 3 species were also identified, which further verified our results that poly-core species would form in supercritical CO₂.

Fig. 8 Complexation at the interface of anhydrate uranyl nitrate/supercritical CO₂ (CO₂ was hidden for clarity).

As for the system of hydrated uranyl nitrate, uranyl nitrate and water molecules were packed in a layer of 25 Å thickness. After heated at 700 K, some uranyl nitrate and water entered the supercritical CO₂ phase. The solid layer remained aggregated at the stage of production simulation. Also, more TBP aggregated at the surface of the solid phase, further reducing the concentration of TBP in the supercritical phase. In the bulk of the supercritical phase, the species were similar to that of the perfect mixing system containing water. No poly-core species was formed. Presumably, when the solid phase was “disturbed” by heating, the uranyl nitrate and water entered the supercritical CO₂ phase simultaneously, holding a proportion consistent with the solid phase.

4 Conclusions

The complexation mechanism of uranyl with TBP in a water-deficient or non-water system, such as supercritical CO₂ extraction or ionic liquid extraction, is quite different from that of a liquid–liquid extraction system. Here, the complexation of anhydrous and hydrate uranyl nitrate with TBP in supercritical CO₂ was investigated. Poly-core uranyl species were present in the system without water once 1 : 1 species was formed. However, in the system with water, the 1 : 1 species quickly coordinated by one water molecule, resulting in the formation of TBP·[UO₂(NO₃)₂]·H₂O, preventing the formation of poly-core species. It also revealed that uranyl nitrate would complex with TBP and H₂O quickly, and then gradually changed to a thermodynamically stable 1 : 2 complex. The attempt of describing the species distribution with time was successful, which helped to understand the path of complexation, and how the experimental conditions take effect. The results facilitate the understanding of complexation of uranyl with TBP in supercritical CO₂, and the developed speciation distribution analysis would be beneficial to the simulation research on the liquid–liquid extraction system.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was supported by the National Science Foundation of China, No. 21790371. The authors thank Prof. Guoxin Tian (CIAE) for the deep discussion of actinide coordination. Liyang Zhu thanks Prof. WA De Jong (LBNL) for the fruitful discussion of molecular dynamics and DFT calculation.

References

1. Y. Lin, R. D. Brauer, K. E. Laintz and C. M. Wai, Anal. Chem., 1993, 65, 2549–2551.
2. A. Leybros, L. Hung, A. Hertz, D. Hartmann, A. Grandjean and O. Boutin, Chem. Eng. J., 2017, 316, 196–203.
3. L. Zhu, W. Duan, J. Xu and Y. Zhu, J. Hazard. Mater., 2012, 241–242, 456–462.
4. M. D. Samsonov, C. M. Wai, S. C. Lee, Y. Kulyako and N. G. Smart, Chem. Commun., 2001, 1868–1869.
5. J. S. Wang, C. N. Sheaff, B. Yoon, R. S. Addleman and C. M. Wai, Chem.–Eur J., 2009, 15, 4458–4463.
6. A. Rao and A. Srivastava, Sep. Purif. Technol., 2021, 257, 117950.
7. L. Zhu, Z. Wang, H. He and G. Tian, RSC Adv., 2016, 6, 96531–96537.
8. P. Sahu, A. K. Deb, S. M. Ali and K. T. Shenoy, J. Mol. Liq., 2021, 330, 115621.
9. G. Benay and G. Wipff, J. Phys. Chem. B, 2014, 118, 3133–3149.
10. M. J. Servis, D. T. Wu, J. C. Shafer and A. E. Clark, Chem. Commun., 2018, 54, 10064–10067.
11. S. Cui, V. F. de Almeida and B. Khomami, J. Phys. Chem. B, 2014, 118, 10750–10760.
12. P. Sahu, S. M. Ali and K. T. Shenoy, Phys. Chem. Chem. Phys., 2016, 18, 23769–23784.
13. P. Guilbaud, L. Berthon, W. Louisfrema, O. Diat and N. Zorz, Chem.–Eur J., 2017, 23, 16660–16670.
14. M. Baaden, R. Schurhammer and G. Wipff, J. Phys. Chem. B, 2002, 106, 434–441.
15. R. Schurhammer and G. Wipff, J. Phys. Chem. A, 2005, 109, 5208–5216.
16. G. Gopakumar, B. Sreenivasulu, A. Suresh, C. V. S. Brahmmananda Rao, N. Sivaraman, M. Joseph and A. Anoop, J. Phys. Chem. A, 2016, 120, 4201–4210.
17. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. V. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery Jr, J. E. Peralta, F. Ogliaro, M. J. Bearpark, J. J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, T. A. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman and D. J. Fox, Gaussian, Inc., Wallingford CT, 2016.
18. S. Grimme, S. Ehrlich and L. Goerigk, J. Comput. Chem., 2011, 32, 1456–1465.
19. T. Moriyoshi, T. Kita and Y. Uosaki, Ber. Bunsenges. Physik. Chem., 1993, 97, 589–596.
20. T. Lu and F. Chen, J. Comput. Chem., 2012, 33, 580–592.
21. J. Mu, R. Motokawa, C. D. Williams, K. Akutsu, S. Nishitsuji and A. J. Masters, J. Phys. Chem. B, 2016, 120, 5183–5193.
22. D. A. Case, I. Y. Ben-Shalom, S. R. Brozell, D. S. Cerutti, T. E. Cheatham III, V. W. D. Cruzeiro, T. A. Darden, R. E. Duke, D. Ghoreishi, M. K. Gilson, H. Gohlke, A. W. Goetz, D. Greene, R. Harris, N. Homeyer, Y. Huang, S. Izadi, A. Kovalenko, T. Kurtzman, T. S. Lee, S. LeGrand, P. Li, C. Lin, J. Liu, T. Luchko, R. Luo, D. J. Mermelstein, K. M. Merz, Y. Miao, G. Monard, C. Nguyen, H. Nguyen, I. Omelyan, A. Onufriev, F. Pan, R. Qi, D. R. Roe, A. Roitberg, C. Sagui, S. Schott-Verdugo, J. Shen, C. L. Simmerling, J. Smith, R. Salomon-Ferrer, J. Swails, R. C. Walker, J. Wang, H. Wei, R. M. Wolf, X. Wu, L. Xiao, D. M. York and P. A. Kollman, AMBER 2018, University of California, San Francisco, 2018.
23. A. Zhu, X. Zhang, Q. Liu and Q. Zhang, Chin. J. Chem. Eng., 2009, 17, 268–272.
24. J. M. Seminario, Int. J. Quantum Chem., 1996, 60, 1271–1277.
25. P. Li and K. M. Merz, J. Chem. Inf. Model., 2016, 56, 599–604.
26. L. Martínez, R. Andrade, E. G. Birgin and J. M. Martínez, J. Comput. Chem., 2009, 30, 2157–2164.
27. M. J. Carrott, B. E. Waller, C. M. Wai, M. J. Carrott and N. G. Smart, Chem. Commun., 1998, 373–374.
28. J.-P. Ryckaert, G. Ciccotti and H. J. C. Berendsen, J. Comput. Phys., 1977, 23, 327–341.
29. A. W. Götz, M. J. Williamson, D. Xu, D. Poole, S. Le Grand and R. C. Walker, J. Chem. Theory Comput., 2012, 8, 1542–1555.
30. R. Salomon-Ferrer, A. W. Götz, D. Poole, S. Le Grand and R. C. Walker, J. Chem. Theory Comput., 2013, 9, 3878–3888.
31. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38.
32. X. Cao and M. Dolg, J. Mol. Struct.: THEOCHEM, 2004, 673, 203–209.
33. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
34. D. R. Roe and T. E. Cheatham, J. Chem. Theory Comput., 2013, 9, 3084–3095.
35. P. K. Verma, N. Kumari, P. N. Pathak, B. Sadhu, M. Sundararajan, V. K. Aswal and P. K. Mohapatra, J. Phys. Chem. A, 2014, 118, 3996–4004.
36. V. Serezhkin, G. Sidorenko, D. Pushkin and L. Serezhkina, Radiochemistry, 2014, 56, 115–133.
37. P. Tecmer, F. Schindler, A. Leszczyk and K. Boguslawski, Phys. Chem. Chem. Phys., 2020, 22, 10845–10852.
38. L. Zhu, W. Duan, J. Xu and Y. Zhu, Ind. Eng. Chem. Res., 2010, 49, 11195–11199.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/d1ra06512b

---