Molecular Insights into Water Clusters Formed in Tributylphosphate

Publication Date (Web): February 7, 2019. Copyright © 2019 American Chemical Society. (V.H.D.) E-mail: [email protected]; ...
0 downloads 0 Views 4MB Size
Article Cite This: J. Phys. Chem. B XXXX, XXX, XXX−XXX

pubs.acs.org/JPCB

Molecular Insights into Water Clusters Formed in Tributylphosphate−Di-(2-ethylhexyl)phosphoric Acid Extractant Systems from Experiments and Molecular Dynamics Simulations Deepak U. Bapat and Vishwanath H. Dalvi* Department of Chemical Engineering, Institute of Chemical Technology, Mumbai 400019, India

Downloaded via WEBSTER UNIV on February 8, 2019 at 05:18:46 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Di-(2-ethylhexyl)phosphoric acid (D2EHPA) and tributylphosphate (TBP) are two of the most studied and researched organophosphorous extractants. D2EHPA is an acidic extractant, offering both hydrogen bond donor and acceptor sites while TBP, a neutral extractant, only offers a single acceptor site per molecule. In spite of this, it is observed that 1 M D2EHPA in dodecane is a poorer extractant for water than 1 M TBP in dodecane. The objective of present work is to look into the molecular interactions that cause such behavior. Experiments were carried out with varying molar ratios of TBP and D2EHPA in the organic dodecane phase. Total extractant concentration was kept constant at 1 M with dodecane as diluent. Water extraction was quantified by measuring the moisture content of the organic phase after equilibration. 1H and 31P NMR spectra of the organic phase samples were recorded to study the change in the chemical environment upon extraction. Small angle X-ray scattering data of water saturated extractant phases were analyzed for the possibility of a reverse micellar aggregate formation. Molecular dynamics simulations could calculate free energies in quantitative agreement with experiments. Experimental and simulation studies showed that aggregation in the organic phase was promoted by the presence of water. This combined approach, of experiments and simulation, has shown that water is indispensable for the formation of ordered aggregates of extractants in nonpolar organic solvents. It is seen that, in the organic phase, around 80% of water’s hydrogen bonds are with extractant molecules rather than with itself. The analysis clearly indicates that, rather than forming an aqueous core surrounded by extractant, water acts as a bridge between extractant molecules.

1. INTRODUCTION Solvent extraction is one of most important unit operations used in the chemical industry. Use of a ligand or an extractant for preferentially isolating a solute is widely practiced. To this end, various extractants such as ketones, organophosphorous compounds, long chain alkyl amines, diglycolyamides, etc. have been used. Partitioning of the solute is driven by the preferential affinity or binding of the solute with the extractant or ligand. Very often a mixture of extractants is deployed to get enhanced extraction; a phenomenon usually called synergistic extraction. Organophosphorous compounds represent a very important class of compounds which are widely used as extractants during liquid liquid extraction in the nuclear and metallurgical industries. Organophosphorous extractants can be broadly classified into acidic and neutral extractants. Di-(2-ethylhexyl)phosphoric acid (D2EHPA), one of the most important extractants belonging to an acidic category, has long been used in the extraction of uranium. Use of neutral organophosphorous extractants along with D2EHPA was found to give a better extraction coefficient for uranium extraction than would be suggested by a linear combination of contributions from the separate extractants.1 D2EHPA with tri-n-octylphosphine oxide (TOPO) is known to give synergistic extraction of uranium.2 A mixture of D2EHPA with cyanex 302 was utilized for the separation of cobalt and nickel.3 D2EHPA-cynex-302 © XXXX American Chemical Society

combination gave better performance over D2EHPA-cynex272 mixture in terms of selective extraction of cobalt over nickel. The mixture of N,N′-dimethyl-N,N′-dioctylhexylethoxymalonamide (DMDOHEMA) and D2EHPA shows synergy for lanthanide extraction.4 The mixture showed synergy in the presence of high metal concentration or at 1 M HNO3. But the same mixture showed antagonistic behavior at pH > 2 or when the metal concentration was in traces. Swami et al.5 used a mixture of D2EHPA and N,N-dioctyl-2-hydroxyacetamide for isolating Eu (III) and Am(III) from nitric acid and achieved a separation factor of about 10 for Eu(III) over Am(III). A combination of D2EHPA with tributyl phosphate (TBP) was found to demonstrate high selectivity for extracting lithium from sea bittern over sodium and potassium ions.6 Separation factor for Li over Na or K ions was found to decrease with increase in the concentration of Na and K ions. TBP is an important neutral organophosphorous extractant which has been investigated by many researchers due to its ability to extract a variety of solutes such as metals, organic carboxylic acids, phenols, etc. TBP is also known to form synergistic mixtures along with various extractants. TBP along with its degradation product, di-n-butylphosphoric acid Received: November 7, 2018 Revised: January 10, 2019

A

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Various analytical techniques such as SANS and SAXS studies have been used to study the synergistic behavior of extractant systems. It has been observed that the extractants assemble to form aggregates in the organic phase. Synergy has been studied with the help of various experimental techniques.18−21 In addition to these analytical techniques, spectroscopic methods such as FTIR, NMR, etc., are also useful tools to probe the chemical environment of the system. For instance, NMR has been effectively used to study the nature of molecular aggregates by surfactant and hydrotropes in organic media.22−25 Marie et al.26 used 31P NMR and electrospray ionization mass spectrometry to study the complexation behavior of D2EHPA with lanthanides (La and Sm) in 1,3-diisopropylbenzene diluent. The solvent extraction system was analyzed in the presence and absence of lactic acid. On the basis of the NMR spectroscopy authors identified the presence of different phosphorus environments existing in the presence of lactic acid. In another study, the authors were able to deduce 1:1 adduct between D2EHPA and a neutral extractant octyl(phenyl)-N,N-diisobutylcarbamoylmethylphosphine oxide (CMPO) from 31P NMR studies.27 Molecular dynamics (MD) is a versatile tool for studying molecular interactions. Reliable molecular level description of TBP diluent solute system is essential considering widespread applications of TBP. MD simulations of tributyl phosphate have been tried by several authors. Some of the earliest works on TBP simulations were undertaken by Beudaert et al.28 to study TBP conformation in vacuum and at the water chloroform interface. Cui et al.29 have compared the performance of different force field parameters on physical and electrical properties such as mass density, self-diffusion coefficient and electrical dipole moment. Efforts have also been directed toward quantification of water uptake by 30 vol % TBP organic phase. It was found that the dipole moment had a significant impact on the prediction of water solubility. After scaling down the atomic charges by a factor of 0.9, the authors managed to obtain simulated solubility very close to experimental value. Vo et al.30 have carried out MD simulations of TBP in dodecane. Enthalpy of vaporization, self-diffusion coefficient, and dipole moment were used to validate models. TBP is known to form aggregates in nonpolar solvents such as dodecane. The authors also studied selfassociation in TBP molecules using umbrella sampling. Potential of mean force was studied as a function of the separation between phosphorus of one TBP molecule and double bonded oxygen of the other molecule. The energetically favorable structure was found with the separation of around 0.525 nm and anti parallel orientation of molecules. Similar studies have been carried out to study TBP aggregates in nonpolar solvents.31,32 At higher concentrations, trimers of TBP were also found to exist along with dimers and monomers in MD simulations. All studied models considerably underpredicted the self-diffusion coefficient of TBP. Water extracted along with TBP is an important factor during LLE processes involving TBP. Mu et al.33 have used the OPLS force field to simulate TBP. Experimental scattering data of TBP solutions in dodecane obtained by Motokawa et al.15 was used to validate the model. The mixtures of TBP water, as well as TBP alkane, were studied. TBP water mixtures were studied for mixture density and phase separations. A refined model developed by the authors did not show phase splitting in TBP water system when water was within solubility limits. The model developed was able to predict TBP alkane mixture

(HDBP) was shown to have a synergistic effect toward the extraction of lanthanide from nitric acid solution. This synergistic behavior was observed at a molar ratio of a 0.75 M TBP to 0.25 M HDBP.7 TBP along with 1,1,1-trifluoro-3,2′thenoylacetone (TTA) shows synergistic behavior toward uranium extraction from nitric acid solutions.8 The mixture of TBP and D2EHPA has been advantageously used to separate Fe(III) and Zn(II).9 4-(4-fluorobenzoyl)-3-methyl-1-phenylpyrazol-5-one along with neutral organophoshorous extractants such as TBP, triphenylphosphine oxide (TPPO), TOPO, tributylphosphine oxide (TBPO) were investigated for the extraction of lanthanoids.10 The synergistic effect as well separation factors between various lanthanoids were found to be dependent upon the basicity of the neutral organophosphorous extractant. TBP along with a basic extractant tri-n-octylamine was shown to form a synergistic mixture for carboxylic acid extraction.11 But it is worth noting that the synergistic behavior was not equally effective toward all acids. In the case of glycolic, fumaric, lactic and succinic acids this effect was much more pronounced than other carboxylic acids such as acetic, propionic, L-maleic, and itaconic acids. This list of references, not exhaustive by any means, demonstrates that synergy is a prevalent phenomenon in the solvent extraction field. At the same time, we find that this phenomenon is very complex and often dependent on various parameters such as the chemical nature of compounds under consideration, their concentration, pH, etc. Hence it is important to understand molecular origins of synergism considering its widespread occurrence in the solvent extraction processes. Osseo-Asare 12 reviewed the behavior of TBP-based extraction systems. In this seminal work, the formation of TBP reverse micelle or microemulsion in the organic phase was proposed. In a significant deviation from the conventional approach of solute solvent interaction, a colloidal approach was presented. Since this work, various experimental results have been obtained corroborating this claim. Chiarizia et al.13 have used small angle neutron scattering (SANS) to study the extraction system consisting of TBP, dodecane, and uranyl nitrate in nitric acid. It was deduced that TBP formed reverse micelles in the organic phase in the presence of water. The same system was further investigated14 using SANS, smallangle X-ray scattering (SAXS), vapor phase osmometry, tensiometry, and conductivity measurements. Experimental evidence and modeling of aggregates using Baxter’s hard sphere model indicated TBP aggregates consisting of two to five TBP molecules around an aqueous core. The size of the aggregates was found to be dependent on extractant concentration, acid nature, and concentration. However, the shape of the aggregates was found to remain spherical. TBP in dry octane phase was investigated using scattering data, and dimeric TBP aggregates were found to exist.15 It has been generally accepted that D2EHPA forms dimeric structures in nonpolar diluents such as kerosene owing to the presence of an acidic proton. Some of the earliest studies by Sato16 have tried to explain the synergism mechanism of uranium extraction by IR studies. The synergy of D2EHPA with another neutral extractant, trioctyl phosphine oxide (TOPO) was studied using SAXS data.17,18 In this study, the role of preorganized structures was investigated as a prerequisite for synergistic behavior toward uranium extraction. It was found that extractants did not form aggregates in absence of an aqueous phase containing the uranium cation. B

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B densities satisfactorily. Das et al.34 have used OPLS force field with their set of partial atomic charges to study neutral trialkyl phosphates. Looking at various models available we can see that a complete description of “all” thermophysical properties of TBP using any of the above-said models is very difficult. Accurate predictions of interfacial properties of organic and aqueous phases require good predictions of both TBP−organic as well as TBP−aqueous mixtures. Hence the model developed by Mu et al.33 was used in subsequent studies considering the model’s ability to predict TBP−water and TBP−alkane interactions adequately. Pecheur et al.18 have carried D2EHPA simulations using AMBER force field to investigate synergism between D2EHPA and trioctyl phosphine oxide (TOPO). The organic phase containing D2EHPA and TOPO was simulated to study the aggregation in absence of water. The average aggregation number was found to be 2−3, which did not show deviation even at the synergistic ratio of 1 TOPO: 4 D2EHPA. Aggregation numbers obtained are too small to classify it as a micelle. It would be more appropriate to term these aggregates as dimers or trimers. Murzyn et al.35 have presented refined OPLS parameters for simulating dimethyl phosphate in n-pentadecane. Model parameters were validated against liquid phase density, enthalpy of vaporization, and hydration free energy. Although sodium salt of D2EHPA has been studied using molecular dynamics, parametrization of D2EHPA has not been undertaken systematically. Inert nonpolar solvents such as dodecane, heavy normal paraffin, etc. have been used as diluents for organo-phosphorus extractants. It is important to select a proper model describing the diluent system since it is known to affect the extraction properties significantly. Use of default OPLS parameters which are adequate to describe short chain hydrocarbons, when applied to long chain hydrocarbons, showed undesirable liquid to gel type transition. Siu et al.36 have proposed a new set of OPLS parameters for alkanes and alkenes. The developed parameters were found to give improved predictions of liquid phase properties over previous studies. Although synergism has been found in a variety of extraction systems, molecular understanding of these systems remains far from complete. It has been thought that the extractant molecules form reverse micellar structures and the synergism coincides with the tendency of extractants to form aggregates in the organic phase. Baldwin et al.37 have suggested that TBP solutions are molecular solutions questioning the use of a colloidal approach. In the present work, molecular level description of TBP and D2EHPA and the role of water in this extraction system is investigated. Role of water in the synergism is studied by using 1 H NMR, 31P NMR spectroscopy, and SAXS. Molecular dynamics simulations have been performed to gain molecular level insights into the system. Free energy calculations using MD simulations have been undertaken to predict the partitioning of water in the extractant phase. It was seen that water formed majority of its hydrogen bonds with extractant molecules rather than with itself: effectively acting as a bridge between extractant molecules rather than forming an aqueous core surrounded by extractant molecules. Results show high fidelity between simulated and experimental values of free energies. Hydrogen bonding analysis shows that water forms a greater number of hydrogen bonds per molecule in D2EHPA solutions than in TBP solutions but still is extracted to a lesser degree in the former compared to the latter. This may indicate

a phenomenon analogous to the hydrophobic effect at work where strong enthalpic interactions carry a substantial entropic penalty.38

2. MATERIALS AND METHODS 2.1. Experimental Section. Tributyl phosphate (97%), dodecane were purchased from S. D. Fine Chemicals India. D2EHPA (98%) was purchased from Spectrochem India. DI water from Millipore system was used throughout the experiments. Materials were used as received without further purification. Organic extractant phases were prepared by dissolving appropriate quantities of extractants in dodecane. The total concentration of extractant was kept constant at 1 M (TBP + D2EHPA). Solutions with varying molar ratios (0:1, 0.25:0.75, 0.5:0.5, 0.75:0.25 and 1:0 M; TBP:D2EPHA) were prepared. In order to study the amount of water extracted in the organic phase, equal volumes of organic and aqueous phases were thoroughly mixed at 300 K. After ensuring equilibrium, both phases were separated and the organic phase was analyzed for moisture content. The water content of the organic phase was measured using Karl Fischer titration. Sample volumes of organic phases were adjusted such that analyzed water content was greater than 15 mg. Karl Fischer readings were done in triplicate, and an average value of three readings is reported. NMR analysis of samples was carried out at Sophisticated Analytical Instrument Facility (SAIF), IIT, Bombay, on a JEOL, 600 MHz instrument. The coaxial insert technique was used to record spectra. CDCl3 was used as a lock solvent in the inner capillary tube. 31P NMR shifts are reported against external 85% H3PO4. All measurements were done at 298 K. To study the effect of water content in organic phase one set of extractant mixture was analyzed with varying degree of water content in it. An extractant mixture with 0.5 M TBP and 0.5 M D2EHPA was selected for this purpose. Organic mixtures were prepared with varying water content (0, 0.14, 0.28, 0.44, 0.61, 0.71, and 1.03 wt %). The shift of D2EHPA acidic protons was monitored in these mixtures. 1 M organic solutions equilibrated with water were used for SAXS studies. Small angle X-ray scattering experiments were carried out at CSIR-CCMB SAXS facility in Hyderabad by using BioSAXS-2000 (Rigaku Systems) X-ray generator MicroMax-007 HF with a microfocus rotating anode using X-ray wavelength of 1.5418 Å. The scattering data were recorded on a HyPix3000 detector with a pixel size of 100 × 100 mm. Scattering was measured with a q-range of 0.008−0.3 Å−1 at an exposure time of 30 min for each sample. The scattering vector (q), was defined as q = 4π(sin θ)λ−1, where 2θ is the scattering angle and λ is the X-ray wavelength. Intensities were converted to absolute scale by measuring scattering of water as the primary standard. The scattering intensity I(q) (cm−1) is expressed as I(q) = npVp 2Δρ2 P(q)S(q)

(1)

where np is the concentration of particles in the solution (number.cm−3), Vp is the volume of scattering particle, and P(q) and S(q) represent the form factor and structure factor of the scattering particle, respectively. Δρ (cm−2) is the difference between scattering length densities of solvent and scattering particle. In the present work, factor npVp2Δρ2 from eq 1 was modeled by using a single constant, Ciq, having units of cm−1. The rationale for not calculating np explicitly is provided in the C

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

dihedral angle between the atoms i, j, k, and l. V1, V2, V3 are constants. Nonbonded interactions are represented by a fourth term, which takes into account electrostatic and Lennard-Jones interactions between atoms i and j. TBP partial atomic charges were adopted from Mu et al.33 Chemical structures and atomic names used in simulation are shown in Figure 1. In the present work, equilibrium bond

subsequent Results and Discussion section describing the SAXS data. The form factor (P(q)) of a spherical particle with radius R is given by ÄÅ É ÅÅ 3(sin(qR ) − qR cos(qR )) ÑÑÑ2 Å ÑÑ Å P(q) = ÅÅ ÑÑ 3 ÅÅÅ ÑÑÑ ( qR ) (2) Ç Ö It is known that the monodispersed spherical form factor of sphere shows oscillations at high q range. These oscillations smear out as the polydispersity increases. The effect of polydispersity on form factor is represented by ̅

P(q) =

∫ n(R)P(q) dR

(3)

Here, n(R) represents the particle size distribution. In line with most commonly used approach for modeling particle size distribution, the Schultz−Zimm function is used to represent n(R).39 The functional form of the distribution with numberaveraged radius R0 and root mean squared deviation σ is shown in the Supporting Information. Interparticle attractions were modeled using the Baxter’s surface adhesion model.40 The potential (U(r)) in kbT (kb = Boltzmann constant) acting between particles separated by distance r, obeying Baxter’s surface adhesion model is given by ÅÄÅ i ÅÅ j δ − Å j U (r ) = lim lnÅÅÅÅ12τ jjjj δ→ dhs /2 Å ÅÅ jj δ ÅÇ k U (r ) = ∞

U (r ) = 0

Figure 1. Chemical structure of extractants is to the left and the corresponding naming convention used in the MD simulations is shown on the right.

dhs 2

É dhs yÑ Ñ zzÑÑÑ d 2 z zÑÑ hs < r < δ zzÑÑ zzÑÑ 2 {ÑÑÖ r
δ

length of P-OS was modified from the paper. In the paper, the equilibrium length of the P-OS bond was quoted as 0.1697 nm, which is different from 0.161 nm used by all models. Hence, equilibrium bond length of 0.161 nm, consistent with all other models was used. D2EHPA parameters were adopted from OPLS simulations of dimethyl phosphoric acid (DMP).35 In this work, simulation, charges on C1 protons of D2EHPA were adjusted. Charges of methyl group of DMP were redistributed on methylene group to make it neutral. All other parameters for the alkyl chain of D2EHPA were borrowed from the standard OPLS library. Partial atomic charges and LennardJones parameters used in the simulations are shown in Table 1 and Table 2 respectively.

(4)

where dhs is the hard sphere diameter and δ − dhs is the range of the small attractive well. For particles interacting via a potential defined by Baxter’s model, an analytical expression for the Ornstein−Zernike equation in the Percus−Yevic closure has been used. The expression for structure factor (S(q)) for such systems with η volume fraction of particles has been given previously but is reproduced in the Supporting Information for the sake of clarity.41 2.2. Simulation Details. All simulations were performed by using GROMACS package (v.5.1.2).42 The OPLS force field was used for calculating the interaction energy.43 The potential energy of interaction E is given by E=

Table 1. Partial Atomic Charges and Lennard-Jones Parameters of TBP

Kr K (rij − req)2 + ∑ θ (θijk − θeq)2 2 bonds angles 2 ÄÅ ÅV V + ∑ ÅÅÅÅ 1 (1 + cos Φijkl ) + 2 (1 + 2 cos Φijkl ) ÅÅÇ 2 2 dihedrals Ä ÉÑ ÉÑ ÅÅÅÅ yÑÑÑ 6z qiqje 2 jij σ 12 ÑÑ ÅÅ V3 σ Ñ z + (1 + cos 3Φijkl )ÑÑÑ + ÅÅ∑ ∑ + 4ϵijjjj 12 − 6 zzzÑÑÑ zzÑÑ jj r ÑÑÖ ÅÅÅ 2 4π ϵ0rij r ij ij i j {ÑÑÖ k ÇÅ (5)



where, first term represents bond stretching potential between atoms i and j. rij and req represent the distance between the atoms i and j and equilibrium distance, respectively. Kr is the force constant. Second term is angle bending potential between atoms i, j, and k with θijk being the actual and θeq the equilibrium bond angle. Kθ is the corresponding force constant for bending. Dihedral torsion potential between atoms i, j, j, and k is represented by the third term. ϕijkl is the

atom type

charge

ϵ (kJ mol−1)

σ (nm)

O2 P OS C1 C2 C3 C4 H1 H2 H3 H4

−0.7070 1.34 −0.5047 0.1222 0.0297 0.1349 −0.2534 0.0361 0.0232 −0.0175 0.0589

0.8368 0.8368 0.5858 0.2716 0.2716 0.2716 0.2716 0.1255 0.1255 0.1255 0.1255

0.2980 0.374 0.2850 0.3500 0.3500 0.3500 0.3500 0.2500 0.2500 0.2500 0.2500

The cut off radius for short-range electrostatic and LennardJones interactions was 1.0 nm. Long range electrostatic interactions were computed by using particle-mesh Ewald summation.44 The starting point for each simulation was a random insertion of molecules in a cubic box of large volume to avoid molecular overlap. Energy of the system was minimized using the steepest descent method. The energy D

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

simulation of 25 ns, water molecules were added to the box which was extended in z direction so as to ensure equal volumes of aqueous and organic phases. Data was collected from the procedure similar to pure component simulations. This biphasic data was used to calculate interfacial tension. Interfacial tension (IFT) was calculated by

Table 2. Partial Atomic Charges and Lennard-Jones Parameters of D2EHPA atom type

charge

ϵ (kJ mol−1)

σ (nm)

O2 P OS OH HO C1 H1 C2 C3 C4 HC

−0.92 1.62 −0.60 −0.53 0.45 0.38 −0.045 −0.18 −0.12 −0.06 0.06

0.8368 0.8368 0.585 76 0.711 28 0.0 0.2716 0.1255 0.2716 0.2716 0.2716 0.125 52

0.315 0.374 0.29 0.300 0.0 0.3500 0.2500 0.3500 0.3500 0.3500 0.2500

l Pxx(z , t ) + Pyy(z , t ) | o o o o dz ( , ) − P z t m } zz o o o o 2 n ~ | Pxx(t ) + Pyy(t ) o o L l o o = zm Pzz(t ) − } o o o 2 no (6) n ~

γ (t ) =

1 n

∫0

Lz

where γ(t) is the average surface tension, n is number of surfaces (2 in all cases run here), Lz is the height of the box in the z direction. Pxx, Pyy, and Pzz denote normal components of pressure in x, y, and z directions, respectively. The term in the braces is the difference between normal and lateral pressures. Organic phases saturated with water were also simulated. Amount of water was determined from the Karl Fischer experiments. Number of water molecules were inserted to ensure agreement with experimental results. In Table 3, number of molecules used for each organic phase simulation is tabulated. The cubic cell dimension for this simulation setup was around 10 nm. Hydrogen bond (H-bond) analysis on the generated trajectory was performed. In this work, double bonded oxygen of TBP, double bonded oxygen of D2EHPA, hydroxy oxygen of D2EHPA and oxygen of water were identified as hydrogen bonding centers. A H-bond between two such oxygens is defined when O−O distance is less than 0.35 nm and HOO bond angle is less than 30°. In our script, nearest neighbors are obtained by using a Voronoi tessellation. The in-house code for hydrogen bonding analysis was checked against the in built program utility from GROMACS. Both the procedures gave identical results, often total number of hydrogen bonds identified not differing by more than 10 counts. Hence depending on the ease of operation either of the codes was used interchangeably. The cluster analysis was performed on the extractant mixture simulations. The analysis was based on the hydrogen bonded pairs formed in the simulation. Well

minimization step was followed by a canonical ensemble (NVT) run for 100 to 300 ps. The velocity-rescale modified Berendsen thermostat was used for keeping the temperature of the system at 300 K during NVT simulations. After completing the NVT simulation, system was equilibrated for 50 to 100 ps in isobaric isothermal NPT ensemble to ensure close agreement with the experimental density value. Once system showed no significant deviations in the density values, data was collected from long production runs. A production run was carried out within isobaric isothermal NPT ensemble. The Parrinello−Rahman pressure coupling was used for keeping the pressure of the system to 1 bar. The time constant for barostat was 2 ps. Temperature of the system was maintained at 300 K using the Nose-Hoover thermostat. The equations of motion were integrated with the leapfrog algorithm with 1 fs time step. All hydrogen bonds were constrained by the LINCS algorithm.45 All long-run simulations were performed for 20 to 30 ns. Number of molecules used for simulations along with the corresponding length of the cubic box is shown in Table 3. TIP4P water model was used throughout the work.46 The force field parameters for TIP4P water model can be found in the Supporting Information, Table S1. The starting point for biphasic simulations was the pre-equilibrated box of the compound which was obtained after NPT simulation. After getting desired density of the compound and ensuring long

Table 3. Number of Molecules Taken for Organic Phase Simulation organic phasea

TBP

D2EHPA

dodecane

water

L (nm)

1M-TBP-N 1M-TBP-W 0.75-TBP-0.25-D2EHPA-N 0.75-TBP-0.25-D2EHPA-W 0.5-TBP-0.5-D2EHPA-N 0.5-TBP-0.5-D2EHPA-W 0.25-TBP-0.75-D2EHPA-N 0.25-TBP-0.75-D2EHPA-W 1M-D2EHPA-N 1M-D2EHPA-W TBP TBP-W D2EHPA dodecane

450 450 354 354 234 234 114 114 − − 500 500 − −

− − 114 114 234 234 354 354 450 450 − − 500 −

1470 1470 1398 1398 1398 1398 1398 1398 1350 1350 − − − 500

− 204 − 204 − 204 − 150 − 96 − 500 − −

9.13 9.14 9.09 9.11 9.13 9.14 9.17 9.18 9.11 9.11 6.06 6.18 6.46 5.74

Organic phases containing water have been represented by names ending with -W whereas “dry” organic phases are presented by names ending with -N. Pure component is represented as it is. Free energy calculations for mixtures were performed with one thirds of the molecules used for bulk simulations. The length of cubic box for the free energy calculations was in the range of 6.1 to 6.3 nm.

a

E

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

configuration starting with a neutral water molecule was alchemically charged to TIP4P water charges. For steps two and three, simulations were carried out at discrete points and the results at these points were numerically integrated to give the free energy (ΔG) change in every step.

spaced frames (around 35 to 40 frames per simulation) from simulation were collected and a procedure similar to the one used by Servis et al.47 was used to identify the clusters formed. On the basis of the cluster analysis the cluster size distribution histograms were obtained. Radial distribution function (RDF) was used to analyze the spatial distribution of the extractants in the system. The built in GROMACS utility was used calculate the RDF. Coordination number (CN) analysis was done to quantify the trajectory averaged neighbor distribution in the system. n(r′) = 4πρ

∫0

r′

g (r )r 2 d r

ΔG =

∫0

1

d H (λ ) λ



(8)

where, H is the Hamiltonian of the system which is dependent on the “state” of the system defined by λ. Steps 2 and 3 were carried out at 15 Kronrod Gauss nodes. Data collection for free energy was done for 3 to 5 ns after doing energy minimization and system equilibration. Data of first 1 ns was discarded and the remaining data was used for calculations. Error associated with discrete intermediate point was calculated by using block averaging technique. This procedure can be carried out both in the aqueous and the organic phase to calculate the excess chemical potential of water. On the basis of free energy estimation, the distribution coefficient of water (KD) in terms of excess chemical potential sol of water in aqueous (μhyd excess) and solvent phase (μexcess) can be calculated as49,50

(7)

Here, ρ is the particle density, g(r) is the RDF between the particles under consideration. n(r′) gives the number of neighbors around the particle in a spherical shell of radius r′. Free energy calculations were used for predicting distribution of water between organic extractant phase and pure water phase. In order to calculate partitioning, free energy of inserting a water molecule from vacuum into the concerned phase was calculated. The process of solvation free energy calculation was broken down in three alchemical intermediate steps. This procedure is represented schematically in Figure 2.

ln KD =

hyd sol μexcess − μexcess

RT

=−

excess ΔGtransfer

RT

(9)

ΔGexcess transfer

where is also called as the excess transfer free energy between two phases, R is the universal gas constant and T is the system temperature in K. The experimental excess free excess energy of this transfer ΔGexpt is evaluated from the experimental distribution coefficient (KD) as follows. excess ΔGexp t = RT ln KD

(10)

3. RESULTS AND DISCUSSION 3.1. Water Content. When only TBP is equilibrated with water, around 68 g L−1 water is extracted by TBP phase.51 This water content corresponds to a molar ratio of TBP to water of around 1. However, TBP is rarely used without the diluent for extraction purpose. Hence it is important to study the water extracted by the solvent−TBP organic phase. The water content of the organic phase as a function of TBP volume fraction is shown in Figure S1. With increase in TBP content of organic phase, the amount of water extracted was found to increase in good agreement with the water extraction data reported earlier.51,52 In Figure 3, the amount of water extracted by the mixture of extractants as a function of extractant composition is shown. It is evident that the extent of water extraction by 1 M D2EHPA is lesser than that of 1 M TBP. When TBP is present along with D2EHPA in the extractant mixture, the amount of extracted water increases. Water extracted by the equimolar mixture of 0.5 M TBP and 0.5 M D2EHPA approaches the values shown by 1 M TBP. A line connecting the two extremes is shown in the figure. This line represents a hypothetical linear variation of extractant water content as a function of TBP to D2EHPA ratio. The experimental water content points for all mixtures always lie above the straight line representing a linear variation and positive deviation is seen at the equimolar ratio. Similar observations were made in the extraction system of TBP and di-n-butyl phosphoric acid (HDBP).19 Non linear variation of water extraction as well as dysprosium(III) (Dy(III)) from nitric acid mixture as a function of TBP to

Figure 2. Schematic representation of alchemical transformation to calculate free energy of solvation. State A to state B represents end states i.e. inserting solute from vacuum state (A) into a solvent (B). This procedure is subdivided in three steps shown in the figure. First step involves Widom insertion of small neutral test particle in the solvent (A−I). Step two includes increasing the Lennard-Jones parameters on the test particle to that of solute particle (B−I). The last step involves turning on the electrostatic interactions.

The first step involved calculation of excess chemical potential of a single neutral small-diameter Lennard-Jones (LJ) sphere in the fluid. The diameter was small enough to fit into interstices of extractant molecules. Hence, Widom’s test particle insertion method was used for the same.48 The LJ parameters were chosen to ensure high probability of successful insertion. In the present work, σ and ϵ parameters of the test particle were 0.09 nm and 0.065 kJ mol−1 respectively. In the second step, test particle LJ parameters were transformed to LJ parameters of TIP4P oxygen molecule. Since water modeled using TIP4P model has all LJ parameters situated on the oxygen atom, steps 1 and 2 give free energy of inserting a neutral water like particle in the system. For the organic phase simulations, widom insertion of a particle with LJ parameters of σ = 0.315365 nm and ϵ = 0.64852 kJ mol−1 was done. In the last step, F

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

observed at 3.92 ppm. This peak occurs considerably upfield compared to the bulk water peak which is found around 4.8 ppm. It is evident from this upfield shift that the water is closely associated with the TBP molecules. The peak proton peak integration yields the ratio of water protons to TBP proton as 1:15. This ratio is in excellent agreement with experimental findings of around 1 mol water extraction per mole of TBP. The upfield shift of water protons indicates that TBP−water bond is weaker than the water−water hydrogen bond. This upfield shift also suggests that the water is not able to form its hydrogen-bonded network as it would in the bulk water.53 In Figure 4, NMR spectrum of D2EHPA containing extractant mixtures is shown. The alkoxy proton for D2EHPA are found as multiplets around 4.18 to 4.24 ppm. Alkyl protons are found at 1.48 to 1.62 and 0.98 to 1.10 ppm. The NMR spectrum of pure D2EHPA is characterized by the lone acidic proton which occurs downfield due to acidic nature and hydrogen bonding. The acidic proton of D2EHPA appears at 12.72 ppm. In Figure 4, we find that, in going from pure D2EHPA to a 0.75 M TBP + 0.25 M D2EHPA mixture equilibrated with water, the acidic proton moves upfield from 12.72 to 6.91 ppm. The resonance frequency of acidic protons is dependent on the nature of solvent as well as dilution. In the present work, the presence of water and TBP contributes to the change in chemical environment experienced by D2EHPA molecules. The presence of water in the organic phase results in upfield shift of acidic proton. The effect of water content on the NMR pattern for one of the mixtures was studied. The water content of the equimolar mixture of TBP and D2EHPA was increased gradually. Acidic

Figure 3. Water content of organic phase as a function of TBP to total extractant concentration ratio. The dotted line connects two extreme points corresponding to 1 M TBP and 1 M D2EHPA mixture water contents.

HDBP molar ratio was seen. For Dy(III), maximum of distribution coefficient peaked at 25:75 TBP to HDBP whereas for water the maximum was observed at equimolar mixture. In that study this increase in water content was seen as indicative of a reverse micellear system. HDBP is structurally similar to D2EHPA, barring slight differences in the alkyl chain. 3.2. NMR of Mixtures. The NMR spectrum of TBP solution saturated with water is shown in supporting figure S2. When the NMR of pure TBP saturated with water is recorded using the coaxial insert technique, a sharp water peak is

Figure 4. 1H NMR spectrum of D2EHPA and extractant mixtures: (1) neat D2EHPA; (2) 1 M D2EHPA in dodecane; (3) mixture of 0.25 M TBP and 0.75 M D2EHPA in dodecane; (4) mixture of 0.5 M TBP and 0.5 M D2EHPA in dodecane; (5) 1 M D2EHPA equilibrated with water; (6) mixture of 0.75 M TBP and 0.25 M D2EHPA in dodecane; (7) D2EHPA equilibrated with water; (8) mixture of 0.25 M TBP and 0.75 M D2EHPA equilibrated with water; (9) mixture of 0.5 M TBP and 0.5 M D2EHPA equilibrated with water; (10) mixture of 0.75 M TBP and 0.25 M D2EHPA equilibrated with water. In the inset, a spectrum showing the shift of acidic proton is highlighted. Numbering is same as in the main figure. Stacking is done in the order of increasing upfield shift of acidic proton. G

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

TBP is susceptible to hydrogen bonding of TBP with other molecules. In case of TBP complexed with uranyl nitrate, alkoxy protons of TBP were found to resonate around 0.5 ppm downfield from alkoxy protons of uncomplexed TBP protons.54 The behavior shown by D2EHPA is completely opposite to the one shown by TBP. The alkoxy protons move gradually upfield with increase in water content in the mixture. This behavior was evident even at other extractant ratios which were investigated in the present work. In all extractant mixtures, TBP alkoxy protons are found at downfield shifts compared to the corresponding neat extractant mixtures. D2EHPA protons show exactly the reverse trend. TBP and D2EHPA contain two nuclei useful for NMR. 31P peak for pure TBP and pure D2EHPA are found at −2.061 and −2.116 ppm, respectively. The corresponding mixtures saturated with water show peaks at −2.047 and −2.306 ppm, respectively. Variation of peak shift in extractant mixtures in dodecane is depicted in Figure 7. 31P shifts for D2EHPA and

D2EHPA proton NMR shift in equimolar mixture of TBP and D2EHPA as a function of water content in organic phase is plotted in Figure 5. The progress in upshift of proton peak

Figure 5. Effect of water content in organic phase on NMR shift of acidic proton in D2EHPA. Extractant mixture: equimolar mixture of TBP and D2EHPA.

with increase in water content appears gradual. Increased hydrogen bonding is usually associated with downfield shift. From these data, it can be deduced that dimeric structure of D2EHPA is affected by the presence of water along with TBP. In this mixture we do not find a separate peak for water protons. This indicates a fast proton exchange between water and D2EHPA exchangeable protons. Hence the peak could be attributed to the average peak of D2EHPA and water hydrogen exchange. In Figure 6, effect of gradual increase in water content on the NMR shift of alkoxy protons in the equimolar mixture of extractants is depicted. In this figure, the chemical shift variation of the largest peak of alkoxy protons is plotted. With increase in water content, TBP alkoxy protons show downfield shift as opposed to the upfield shift shown by D2EHPA alkoxy protons. It is known that the NMR shift of alkoxy protons of

Figure 7. Phosphorus chemical shifts in various TBP-D2EHPA extractant mixtures. Total extractant concentration in all mixtures is 1M. First name in the legend entry is used to tag the corresponding P atom in the moiety. Legends ending with -N indicate neat extractant mixtures. Legends ending with -W indicate extractant mixtures saturated with water. Continuous lines are shown only to guide the eye.

TBP molecules are plotted on different axes. The trend in chemical shifts of mixtures with change in TBP to D2EHPA ratio is clearly noticeable in the shown figure. Going from neat mixture containing 0.25 M TBP to 0.75 M TBP, phosphorus peak of TBP as well as D2EHPA move sequentially downfield. This clearly indicates that in the neat extractant mixtures, with increase in TBP content, the electron density at P atom in both the extractant molecules decreases. The peak for the water-saturated mixture always appears at upfield shifts (becomes more negative) compared to the corresponding neat mixture peak. The P shift of TBP molecules in water saturated mixture moves downfield with increase in concentration of TBP; a trend which is similar to that shown by phosphorus peak of neat extractant mixtures. In contrast to this, D2EHPA phosphorus peak in water saturated mixture gradually moves upfield with increase in TBP concentration. We can quantify the difference between the neat and water saturated peak position as Δ = δN − δW, where δN and δW represent the shift in the neat and water saturated mixture, respectively. The Δi value increases with decrease in the concentration of the species i (i = TBP or D2EHPA) the

Figure 6. Effect of water content in organic phase on NMR shift of alkoxy protons in TBP and D2EHPA. Extractant mixture: equimolar mixture of TBP and D2EHPA. H

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 8. SAXS data of extractant mixtures. (a) Discrete data points represent the experimental points. Red line indicates SAXS pattern obtained by using fitted parameters. Curves are shifted for better visualization. (b-i) Enlarged view of the SAXS curve representing characteristic peak region is shown. In this figure intensities are reported without scaling factors. (b-ii) Structure factor variation in extractant mixtures according to model; the mixture can be identified with the color of discrete data points from part a.

Table 4. Fitted Parameters Obtained from Least Square Regressiona sample name

R0 (Å)

σ (Å)

τ

Ciq (cm−1)

d (Å)

n

1M-TBP-W 0.75-TBP-0.25-D2EHPA-W 0.5-TBP-0.5-D2EHPA-W 0.25-TBP-0.75-D2EHPA-W 1M-D2EHPA-W

6.21 6.06 5.97 6.01 7.12

3.90 3.86 3.21 2.87 4.87

0.23 0.31 0.47 0.81 0.91

0.24 0.24 0.27 0.33 0.33

10.57 10.51 10.78 11.42 13.63

0.14 0.16 0.16 0.17 0.17

a

The nomenclature is similar to the one used earlier in the text.

mixture. The magnitude of ΔD2EHPA is always higher than that of ΔTBP at all concentration ranges. From this observation we can conclude that the magnitude of change in chemical environment experienced by D2EHPA molecules is higher than that experienced by TBP in the presence of water. TBP as well as D2EHPA possess one hydrogen bond accepting center in the form of double bonded oxygen with phosphorus. This could imply that the electron density at the phosphorus center is affected to a greater extent by the hydrogen bonding due to −OH than it is affected by the −PO linkage. 3.3. SAXS. Extractant mixtures equilibrated with water were analyzed with SAXS. The SAXS profile of mixtures is shown in Figure 8. It is evident from figure that going from mixture containing 1 M TBP to 1 M D2EHPA, scattering pattern shows characteristic changes in peak positions. We see that 1 M D2EHPA mixture is characterized by a peak around 0.37 Å−1. The intensity of this peak reduces with decrease in D2EHPA proportion in the mixture and finally vanishing at 1 M TBP. It is also interesting to note an upturn observed in case of 1 M D2EHPA solution below q range of 0.02 Å−1. Similar trends are also observed in case of other extractant mixtures at lower q values except equimolar mixture. This upturn at low q regions was not noticed in some of the similar studies reported earlier. Earlier studies have not recorded data below q range of 0.01 Å−1 and hence probably behavior obtained with 1 M TBP mixture as well has not been reported earlier. An earlier study with TBP and HDBP which reported mixture SAXS data also did not have any indication of an upturn at low q data.19 The study on diblock polymeric materials has shown such increase in behavior at low q values.55 Such behavior was attributed to the presence of large micelles either interconnected through polymeric chains or micelle packings with strong correlations. In case of D2EHPA,

long chains of D2EHPA, interconnected via hydrogen bonds could give rise to such structures. The SAXS behavior was not modeled using this data with the upturn. The data range selected for modeling was with q range grater than 0.05 Å−1. As mentioned in the modeling section, np, i.e., the number of particles responsible for scattering, was not calculated explicitly. One of problems which made this calculation difficult was lack of knowledge about exact stoichiometry of extractant molecules participating to form one “micelle” or cluster. It has been generally accepted that water molecules form the inner core of the aggregate, but ratio in which TBP and D2EHPA contribute in micelle corona in mixtures was not clear. Also, the scattering length density difference of TBP and D2EHPA from the solvent is not dramatically different. Hence considering these constraints the entire term of npVP2Δρ was modeled as a single constant Ciq with the units of cm−1. The parameter optimization was done using least-squares fitting. The fitted parameters obtained are shown in Table 4. The Figure 8b shows the variation of structure factor as a function extractant concentration in the mixture. Structure factor gives insight into the inter particle separation obtained from the first peak in S(q). This separation (d) is obtained with d = 2π/q*, with q* being the position of first maximum in S(q). From Table 4, we see that going from 1 M TBP to 1 M D2EHPA, the value of d increases. This may be attributed to increased alkyl chain length of D2EHPA. Also, we see that fitted parameters give R0 values which increase with increase in D2EHPA concentration in the solution. We find that fitted particle diameter is always grater than d values by 2 to 3 Å. This indicates that length scales probed by structure factor are falling within particle hard sphere diameter. This might be attributed to alkyl chains of extractant mixtures inter penetrating. Similar observations were made by Chiarizia et I

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B al.13 for the aggregate system consisting of U(VI)-HNO3-TBP aggregate system in n-dodecane. τ−1, often called as stickiness parameter is a measure of attractive interaction between particles. We see this value is smaller for extractants with higher D2EHPA concentration. The diminishing attractions between micelles might me attributed to increased chain length of D2EHPA. We also see that polydispersity obtained from fitted parameters as evident from values of σ, is on the higher side. Polydispersity smears out spherical form factor oscillations at high q range. The exercise of fitting mixture data with a single constant term and assuming spherical sticky hard sphere model with continuous distribution might appear to be an oversimplification of the problem. However, despite this simplification and assumptions stated, the overall trend is captured adequately as is evident from Figure 8b. Further improvements in the model could be achieved with the help of more complex terms, defining the form and structure factor of the particles.56 In SAXS studies, the contrast between scattering particles and the solvent is based on the electron density difference, which offers minimal freedom in manipulating the contrast term. In mixed extractant system such as the one studied in the present work, SANS studies might prove more useful. SANS offers the possibility of changing contrast term based on the hydrogen deuterium ratio of solvent (and) or scattering particle. The scattering length density matching could be done to “hide” one of the extractants to gain more insights into the exact stoichiometry of the aggregate. 3.4. MD Simulations. Dodecane is an inert solvent which has been commonly used as diluent for TBP and D2EHPA. Use of modified parameters avoided the liquid to gel transition reported earlier. Interfacial tension of dodecane with water was evaluated using eq 6 from biphasic simulations. Similarly, pure component physical properties were evaluated for TBP and D2EHPA. In Table 5, values obtained from simulations are

Molecular dynamic simulations were carried out on organic extractant mixtures. Radial distribution function (RDF) between O2 atoms of TBP is plotted in Figure 9a. RDF clearly shows difference between organic phases devoid of water and organic phases saturated with water. In case of dry organic mixtures, the first peak is at 0.65 nm. The first maximum of this peak position is in good agreement with earlier results.61 Vo et al.30 reported the first peak at 0.49 nm. However, this difference may be attributed to the use of different force field used by the authors. In a force field comparison study made earlier, authors found that RDF peak positions depend on the corresponding model used.29 All organic phases saturated with water show one additional peak at 0.43 nm. Presence of the additional peak indicates higher degree of ordering in TBP. Earlier MD studies of TBP in nonpolar solvents have identified TBP dimers.29,30 The fact that the new peak is apparent only in the presence of water indicates that water is helping TBP to form more ordered structure. The electron rich center of TBP, i.e., O2 can form a hydrogen bond with water. The RDF between O2 atoms of TBP and oxygen atoms of TIP4P water is shown in Figure S3. This RDF is characterized by a sharp peak at 0.27 nm followed by a second peak at 0.45 nm. D2EHPA is an acidic extractant. RDF between OH atoms of D2EHPA and O2 atoms of D2EHPA shows a sharp peak situated around 0.26 nm highlighting intermolecular hydrogen bonding between D2EHPA molecules (Supporting Information, Figure S4). D2EHPA when present alone in the system, forms intermolecular hydrogen bonds owing to acidic proton present in the moiety. Hence it predominantly exists in dimeric form in organic nonpolar solvent. RDFs between O2 atoms of D2EHPA is shown in Figure 9b. We see that RDF pattern changes with incorporation of water in the organic phase. D2EHPA can form hydrogen bonds with water and TBP as well. The hydrogen bonding between D2EHPA and TBP is also confirmed with the RDF between OH atoms of D2EHPA and O2 atoms of TBP which has a sharp peak at 0.26 nm (Supporting Information, Figure S5). The peak intensity decreases with incorporation of water in the system. This clearly indicates that some of the TBP D2EHPA hydrogen bonds are broken by water. The RDF between O2 atoms of D2EHPA and OW atoms of TIP4P clearly indicates the intermolecular hydrogen bonding (Supporting Information, Figure S6). The RDF structure between oxygen atoms of TIP4P molecules can shed some light on the hydrogen bonding pattern in the extractant mixture. Figure 10 compares RDF of bulk water and water in extractant mixtures. The bulk water RDF is characterized by first hydration shell followed by two distinct peaks. The water hydrogen bonding network which characterizes the bulk water is disturbed in the presence of extractant mixtures. Although the position of the first peak is never altered, subsequent peaks are broadened. In bulk TIP4P water simulations, we find around 4.7 water molecules in the first coordination shell. In the simulations of organic phases saturated with water in the present work, we find that this number is lowered considerably, and the first coordination shell contains around 0.8 to 1.2 water molecules. At the same time, owing to hydrogen bonding, we find around 1.2 O2 atoms of TBP in the first coordination sphere of water in the case of 1 M TBP. This number gradually decreases to around 0.3 with increase in concentration of D2EHPA in the extractant mixture. Number of D2EHPA O2 atoms in the

Table 5. Physical Properties of Pure Components Obtained with MD Simulations and Comparison with the Corresponding Experimental Values density (kg m−3)

IFT (mN m−1)

compound

simulated

experimental

simulated

experimental

dodecane TBP D2EHPA

752 ± 3 994 ± 3 987 ± 4

745.3157 972.4959 97060

50.1 ± 0.7 7.3 ± 0.3 26 ± 0.5

5358 751 14a

a

IFT calculated by using drop volume method. Value obtained within ±2 mN m−1.

compared to the corresponding experimental values. We find excellent agreement between experimental and simulated values for TBP and dodecane. The deviation of interfacial tension (IFT) values in case of D2EHPA is slightly higher compared to the other two components. Considering the effect of water in extraction process, it was essential to use TBP model representing TBP−water interactions as closely as possible. When biphasic simulations were carried out with 1 M TBP, TBP molecules were found to preferentially assemble at the aqueous organic interface in accordance with the known surface active nature of TBP molecules. D2EHPA simulations predicted a liquid phase density of 987 kg m−3, in very good agreement with the experimental value of 970 kg m−3. IFT between D2EHPA and water was found to be 26 mN m−1. J

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 9. (a) Radial distribution function of O2 atoms in TBP. (b) Radial distribution function of O2 atoms in D2EHPA. Dry organic solutions are displayed with continuous lines and the corresponding water containing solution RDFs are shown with dotted lines. Legend entries follow nomenclature similar to that used in Table 3.

tetrahedral hydrogen bonded structure of water is broken in the organic phases. Polar head groups of extractant molecules provide hydrogen bonding sites for water. Because of the bulky nature of extractant molecules owing to dangling alkyl chains, water has on an average around two to three sites for hydrogen bonding in the first shell. It could be said that the polar head groups mimic the hydrogen bonding of water molecules, which keeps the water molecules well dispersed in the organic phase. In all extractant mixtures an additional small peak is observed in the range from 1 to 1.8 nm. The presence of this additional peak in conjunction with interparticle separation (d), predicted by SAXS fitting exercise could be a useful tool in describing water behavior in the mixture. It could be deduced from these studies that discrete water pockets are present inside the extractant system, screened by alkyl chains of the extractant molecules. The representative images from MD simulation setup are shown in Figure 11. It is clearly evident that water remains well dispersed in the system without agglomerating into a single large water droplet in the organic mixture. Water molecules form very small aggregates of two to three water molecules incorporating polar head groups of extractant molecules. The RDF was used to evaluate the composition of the clusters or aggregates formed. The coordination number gives the trajectory averaged nearest neighbor count around the

Figure 10. (a) Radial distribution function between oxygen atoms of TIP4P molecules. (b) Coordination number (CN) of the corresponding RDFs shown in part a.

first coordination shell of water oxygen atoms decreases from around 1.5 for 1 M D2EHPA to around 0.5 for a mixture containing 0.25 M D2EHPA. This clearly indicates that

Figure 11. Snapshot from MD simulations. The dodecane molecules are not shown for clarity. Also, hydrogen atoms on the alkyl chain are not shown. TBP (green), D2EHPA (pink), and water (blue) are depicted in line representations. Key: (a) 0.75 M TBP−0.25 M D2EHPA TIP4P water mixture; (b) 0.50 M TBP−0.50 M D2EHPA TIP4P water mixture; (c) 0.25 M TBP−0.75 M D2EHPA TIP4P water mixture. Structural visualization was done with the help of VMD. K

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B atom/molecule. The coordination number analysis of parts a and b of Figure 9 is shown in the Supporting Information, Figures S7 and S8, respectively. D2EHPA and TBP molecules were identified with the phosphorus atom center. Oxygen atoms were used to tag TIP4P water molecules. The RDF and coordination number analysis of P atoms in TBP molecules is shown in Figure 12. As expected, the number of TBP P atoms

Figure 13. Variation of number of hydrogen bonds per molecule as a function of TBP to total extractant concentration ratio. First name in the legend entry is used to identify the hydrogen bond per molecule of the legend. Legend entries beginning with “Total” represent total number of hydrogen bonds formed averaged over all hydrogen bonding centers in the system. (Dodecane molecules are not considered in the denominator while calculating this value.) Continuous lines have been used for visual aid only. Dashed lines are used to represent the variation in the neat extractant mixture. Dotted lines are used for mixtures containing water. Figure 12. (a) Radial distribution function between P atoms of TBP molecules. (b) Coordination number (CN) of the corresponding RDFs shown in part a. Dry organic extractant mixtures are depicted with continuous lines. Dashed lines are used to show the corresponding organic phases saturated with water.

The H-bonds formed were categorized into D2EHPA− D2EHPA, D2EHPA−TBP, D2EHPA−water, TBP−water and water−water fractions. Figure 14a, shows the distribution of hydrogen bonds involving D2EHPA molecules in the extractant mixtures. Detailed breakup of the bonds can be found in the Supporting Information, Table S2. In case of the neat extractant mixtures, we see that the fraction of D2EHPA molecules involved in dimeric hydrogen bonding with itself decreases with the increase in TBP concentration in the mixture. This clearly indicates that the dimeric form of D2EHPA is affected in the presence of TBP molecules. We see that in 1 M D2EHPA mixture containing water, water almost exclusively forms H-bonds with D2EHPA and the number of H-bonds per water molecule is comparable to that in the bulk water (see Figure 14b). However, in the extractant mixtures containing TBP we find that the water has a freer choice of Hbonding partners and is therefore much less constrained. Parts b and c of Figure 14 depict the distribution of hydrogen bonds formed by water and TBP molecules, respectively. The knowledge of H-bonds formed was used to identify the clustering in the extractant mixtures. The cluster definition was based on the hydrogen bonding pattern. The cluster analysis of water containing extractant mixtures demonstrated a wide range of cluster size distribution. Figure 15 shows the cluster size distribution in 0.75 M TBP−0.25 M D2EHPA solution in the presence and absence of water. It is clearly evident from this figure that this mixture, rich in TBP content, shows a high portion of extractants remaining in the free form in the absence of water. However, with the introduction of water in the system, we find that other cluster structures with differing cluster size are found in the simulation box. This again highlights the crucial role that water plays in forming these extractant aggregates. The cluster size distribution in other mixed extractant mixtures can be found in the Supporting Information, Figures S11, S12, and S13. Figure 16 shows the cluster size distribution in 1 M D2EHPA solution in the

in the first coordination shell is dependent on the TBP concentration in the mixture. With incorporation of water in the mixture, the coordination number increases. This is driven by the hydrogen bonding between water and TBP molecules. If we integrate the RDF between P atoms of TBP and OW atoms of TIP4P until 0.85 nm we get the number of water molecules encapsulated by polar headgroup of TBP. We can make out that around 1.5 water molecules reside within this shell. This number decreases with increase in concentration of D2EHPA in the mixture. Figures S9 and S10 of the Supporting Information give RDF analysis of D2EHPA P atoms with TBP P atoms and water oxygen atoms, respectively. This combined with the coordination number analysis shown in Figure 12 suggests that water molecules are surrounded by around two to four extractant molecules, depending on the extractant concentration. H-bond analysis was also performed on the simulation data obtained. Well spaced trajectory coordinates were extracted and analyzed for hydrogen bonding in the system. Figure 13 depicts the variation of average number of hydrogen bonds formed per molecule in the extractant mixture. It is evident that hydrogen bonding increases with increase in D2EHPA fraction in the mixture. TBP is incapable of forming H bonds with itself. But it can form H-bonds with D2EHPA and hence we see that in the neat extractant mixture number of H-bonds per TBP molecule increase with increase in D2EHPA concentration. We find that in the neat extractant mixtures average number of H-bonds per D2EHPA molecule hovers around 1. But with incorporation of water in the system it increases and shows around two H-bonds per D2EHPA molecule in 0.75M-TBP-0.25M-D2EHPA mixture. L

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 14. (a) Hydrogen bond distribution of D2EHPA molecules in the extractant mixtures. (b) Hydrogen bond of distribution water molecules in the extractant mixtures. (c) Hydrogen bond of distribution TBP molecules in the extractant mixtures. Note in part c, TBP H-bond distribution in neat mixtures is not depicted because in the absence of water TBP forms all hydrogen bonds with D2EHPA only.

Figure 15. (a) Cluster size distribution in the extractant mixture of 0.75 M TBP−0.25 M D2EHPA without water. (b) Cluster size distribution in the extractant mixture of 0.75 M TBP−0.25 M D2EHPA saturated with water.

presence and absence of water. Similar to the TBP-rich solutions, here we find the cluster size distribution of water rich mixtures also shows a marked change from the neat 1 M D2EHPA mixture. It is clear that neat D2EHPA primarily exists as dimers. However, in water saturated D2EHPA

solutions the dimer fraction decreases considerably. Further when TBP and water are “combined” together the D2EHPA− D2EHPA dimeric fraction is reduced still further (see Figures 15 and 16). Owing to this D2EHPA can be involved in more hydrogen bonds per D2EHPA molecule; a trend which is M

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 16. (a) Cluster size distribution of 1 M D2EHPA without water. (b) Cluster size distribution of 1 M D2EHPA saturated with water.

consistent with the one observed in Figure 13. These findings corroborate the observations of NMR analysis where it was observed that the D2EHPA molecules experience a change in chemical environment around the phosphorus center in response to the combined effect of TBP and water in the system. It is also interesting to note that the theoretical DFT calculations done by Gullekson et al.62 show that water interaction with D2EHPA results in significant change in the chemical nature and dipole around the phosphorus center: in line with our NMR observations(see Figure 7). Figure 17 shows representative aggregates formed in the organic phase. The structure visualization was done with the help of VMD.63 Other examples of typical clusters observed during the course of simulations can be found in the Supporting Information Figure S14. We can clearly see that water acts a bridge connecting the extractant molecules. Counterintuitively, in spite of greater hydrogen-bonding sites in D2EHPA compared to TBP and, indeed, greater hydrogen bonding (see Figure 16), it is observed that 1 M TBP in dodecane is a better extractant for water than 1 M D2EHPA in dodecane. It may indicate that increased Hbonding after a threshold imposes an entropic penalty on the system (analogous to the hydrophobic effect) that offsets the gain due to enthalpically favorable interactions. Further simulations might be required to quantitatively establish this. In Figure 18, uncorrected self-diffusion coefficient values of TBP obtained by MD simulations are presented. All TBP models which have been used until date have always under predicted self-diffusion coefficient of TBP. The self-diffusion constant of pure TBP found in present work, was around 10 times smaller than experimental value of (2.29 ± 0.01) × 10−10 m2 s−1.51 This deviation from experimental findings is somewhat smaller when TBP in dodecane is considered. In 1 M TBP simulated value is around 2.5 times smaller than experimental value. Although values shown by simulations are not very accurate, these figures could still be used to draw qualitative conclusions. One trend which is very clearly noticeable is the reduction in the diffusivity of TBP in neat mixtures with increase in the concentration of D2EHPA. This could be explained with the help of discussions about hydrogen bonding between TBP and D2EHPA as demonstrated from RDF and H-bonding analysis. Also, all extractants saturated

Figure 17. Representative figures of clusters formed are shown in parts a and b. Extractants are depicted as VDW spheres. Water, D2EHPA, and TBP molecules are colored in red, blue, and green, respectively. In structure c, the cluster depicted in part b is shown without alkyl chain carbons and hydrogens. The cluster is rotated favorably to show the hydrogen bonding centers in the cluster. Red spheres represent the oxygens on the extractants, dark green sphere represents the phosphorus atom, and white spheres represent the hydroxyl hydrogens. Pink spheres represent water molecule oxygen (virtual site and oxygen atom).

with water show reduced TBP diffusivity than the corresponding dry mixture. This dip in diffusivity may be attributed to increased ordering or clustering of TBP molecules in organic phase. This trend is consistent with RDF and H-bonding findings of organic phases containing water. 3.4.1. Free Energy Calculations for Distribution Coefficient Prediction. Free energy of inserting a water molecule from vacuum phase to pure TIP4P water phase at 1 atm and 300 K was calculated by alchemical transformation. The free energy value obtained was −27.14 kJ mol−1, which can be compared with −26.56 kJ mol −1 reported earlier.64 The values obtained in the present work are comparable with some of the earlier studies reported on different water models.65,66 Similarly, free energy of inserting water in pure dodecane phase was calculated. One molecule was inserted in box containing 500 dodecane molecules. For dodecane a value of −0.27 kJ mol −1 was obtained. Owing to negligible solubility of water in long chain hydrocarbon solvents, one would expect very small free energy change associated with this transfer. This N

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

inappropriate. In the present work, water distribution coefficients were measured by saturating organic phases with water. So it was more pertinent to compare free energy results with organic mixtures containing water. Composition of the organic phases was stated earlier in Table 3. The free energy simulations were carried out in which organic phase water content was adapted from experimentally obtained moisture content of the organic phase. The excess free energy of transferring one water molecule from aqueous phase to this water saturated system was evaluated. Hence we are effectively reporting lnKD in the manner of earlier workers.67−71 In Figure 19, the excess transfer free energy

Figure 18. Variation of diffusivity of TBP molecules in extractant mixture. Unfilled bars indicate values from simulations of dry organic phase. Filled bars depict the diffusivity after incorporating water in the mixture.

value is close to the experimental value of −1.84 kJ mol −1. Using these values one obtains a solubility of water in pure dodecane phase with the help of eq 9. We see that the method employed here is able to predict extremely low solubility of water in dodecane. Similar exercise was undertaken in pure TBP. One water molecule was inserted in 500 TBP molecules. The free energy associated with this step was found out to be −24.31 kJ mol −1. Mu et al.33 have reported a value −17.42 kJ mol −1 for TIP3P water model with atomic charges used in current work. The variation between free energy values reported earlier and obtained in the present work is significant. Mu et al.33 have used TIP3P and SPC/E water models, whereas in the present work TIP4P water model was used. Another reason for the variation might be the different set of cut off parameters used in the present simulations as opposed to 1.2 nm used by Mu et al.33 To check the sensitivity of used partial atomic charges on the free energy estimate, one set of simulations was carried out with the default OPLS-2005 partial atomic charges for inserting water in pure TBP. These simulations yielded a free energy change of −11.31 kJ mol −1. Mu et al.33 obtained a value of −6.7 kJ mol−1 for TIP3P water with their methodology. Both sets of calculations show under prediction of free energy change, but nevertheless, the difference between the present work and the earlier work is of similar magnitude. These results highlight the crucial role of electrostatic contributions to water TBP interactions. Since it was not the objective of the current work to parametrize TBP molecule, partial atomic charges were used without further modifications. The free energy of inserting a solute from vacuum to pure organic or aqueous phase represents a case where partitioning is compared at infinite dilution. It should be noted that calculating partition coefficient by considering free energy change in pure solvent phase assumes complete immiscibility of two phases. The comparison of these values with experimental measurements carried out at concentrations deviating significantly from infinite dilution would be

Figure 19. Comparison of excess transfer free energy of different extractant solutions computed from MD simulations and experiments.

obtained from MD simulations is compared with the experimental values for different extractant mixtures. We see very good agreement between experiments and simulations for 1 M TBP and 1 M D2EHPA: indicative of the fidelity of simulations to experiments. However, the results show deviation from experimental values for the mixture of extractants. The reasons for deviations could be attributed to the dominant effect of electrostatic interactions in the free energy calculations. An earlier example of TBP simulations with different electrostatic charges underscore the contribution of electrostatics in free energy calculations (a difference of around 12 kJ mol−1 between two sets of partial atomic charges on TBP). Another source of deviation to the free energy calculations might arise out of failure to account for change in the chemical environment experienced by D2EHPA molecules around the phosphorus center in response to the presence of TBP and water. 1H NMR and 31 P NMR studies clearly hint toward these variations. So any variations arising due to change in the dipole may not be adequately captured in the present simulation setup. Also, the partial atomic charges are kept constant so the variation arising out of change of charge distribution will not be captured in the current simulation methodology. But nevertheless qualitative trends predicted by MD simulations are promising. Improved partial atomic charges might result in better prediction of system properties. Bannan et al.71 have used MD simulations for predicting distribution coefficient of various solutes between water and octanol or cyclohexane. In case of very large solutes, it was found that results were dependent on the starting configuration. In such situations, results were O

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B susceptible to poor sampling and slow convergence. The organic phase was considered completely immiscible with water. Concerns were raised in usage of water in the organic phase, citing the possibility of aggregate formation and hence slower convergence of results. In the present work, the solute molecule which has been tested is not very bulky, but the organic system used in the present system is very complex, consisting of extractant molecules for which complete MD parametrization has not been done. We see that the effect of electrostatics is significant, as observed from the large difference in values obtained with two sets of partial atomic charges. The simulations in the relatively simple system consisting of pure dodecane gave good results. Hence although results obtained capture qualitative behavior, for accurate quantitative description better models will have to be tried.



AUTHOR INFORMATION

Corresponding Author

(V.H.D.) E-mail: [email protected]; [email protected]. ORCID

Deepak U. Bapat: 0000-0002-7892-5340 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS D.U.B. thanks Pidilite Industries, India, for providing research fellowship under the aegis of Professor M. M. SharmaPidilite Industries Fellowship scheme. The authors also thank DST for providing the high computing facility at ICT. D.U.B. appreciates the kind help from K. Mallesham and R. Rukmini at CCMB Hyderbad for SAXS data collection. The authors thank Prof. Ashwin Patwardhan for fruitful discussions and helpful comments.

4. CONCLUSIONS Water extraction was quantified by measuring moisture content of the organic phase after equilibration. Pure TBP extracted around 1 mol water per mole of TBP. Acidic proton of D2EHPA showed wide range of peak shift, ranging from 12.72 ppm for pure D2EHPA to 6.91 ppm for extractant mixtures. Presence of TBP and water resulted in significant chemical environment change around the phosphorus center of D2EHPA. Molecular dynamics simulations were carried out to study interactions between extractants. MD simulations showed that aggregation of both extractants in organic phase was promoted by the presence of water. The presence of water in a mixed extractant system of TBP and D2EHPA promotes ordered aggregate formation. SAXS studies in conjunction with MD simulations and NMR studies indicate the presence of small aggregates being formed in the organic phases. Molecular simulations show that, rather than a structure where the inner core of these aggregates consists of water molecules surrounded by polar head groups of organophosphorous extractants, the aggregates are actually the extractants bridged by water molecules. This assembly of molecules is driven by the hydrogen-bonding pattern demonstrated by the mixed extractant system with water. Our simulations are able to reproduce, with high fidelity, the excess free energy change associated with transferring a water molecule from aqueous to saturated organic phase. A striking observation is that D2EHPA solutions, in spite of presenting both H-bonding acceptor and donor sites, extract less water than TBP solutions with only acceptor sites. This may be due to an effect analogous to the hydrophobic effect where favorable enthalpic interactions are substantially offset by unfavorable entropic contributions. In other words, water finds enthalpically very favorable clusters in D2EHPA solutions, but these have very few accessible configurations. Finally, we believe that the methodology used in the present work, can in principle be extended to gain insights into any LLE system.



extractant mixtures, and representative cluster snapshots from simulations (PDF)



REFERENCES

(1) Blake, C.; Baes, C.; Brown, K. Solvent Extraction with Alkyl Phosphoric Compounds. Ind. Eng. Chem. 1958, 50, 1763−1767. (2) Bunus, F. T.; Domocoş, V. C.; Dumitrescu, P. Synergic Extraction of Uranium from Phosphate Solutions with Di-(2Ethylhexyl) Phosphotic Acid and Tri-N-Octylphsphine Oxide. J. Inorg. Nucl. Chem. 1978, 40, 117−121. (3) Darvishi, D.; Haghshenas, D.; Alamdari, E. K.; Sadrnezhaad, S.; Halali, M. Synergistic Effect of Cyanex 272 and Cyanex 302 on Separation of Cobalt and Nickel by D2EHPA. Hydrometallurgy 2005, 77, 227−238. (4) Muller, J. M.; Berthon, C.; Couston, L.; Zorz, N.; Simonin, J.-P.; Berthon, L. Extraction of Lanthanides (III) by a Mixture of a Malonamide and a Dialkyl Phosphoric Acid. Solvent Extr. Ion Exch. 2016, 34, 141−160. (5) Rama Swami, K.; Kumaresan, R.; Venkatesan, K.; Antony, M. Synergic Extraction of Am (III) and Eu (III) in N, N-Dioctyl-2hydroxyacetamide-bis (2-Ethylhexyl) Phosphoric acid Solvent System. J. Mol. Liq. 2017, 232, 507−515. (6) Sharma, A. D.; Patil, N. D.; Patwardhan, A. W.; Moorthy, R. K.; Ghosh, P. K. Synergistic Interplay between D2EHPA and TBP towards the Extraction of Lithium Using Hollow Fiber Supported Liquid Membrane. Sep. Sci. Technol. 2016, 51, 2242−2254. (7) Braatz, A. D.; Antonio, M. R.; Nilsson, M. Structural Study of Complexes Formed by Acidic and Neutral Organophosphorus Reagents. Dalton Trans. 2017, 46, 1194−1206. (8) Irving, H.; Edgington, D. Synergic Effects in the Solvent Extraction of the Actinides-I Uranium (VI). J. Inorg. Nucl. Chem. 1960, 15, 158−170. (9) Azizitorghabeh, A.; Rashchi, F.; Babakhani, A. Stoichiometry and Structural Studies of Fe (III) and Zn (II) Solvent Extraction Using D2EHPA/TBP. Sep. Purif. Technol. 2016, 171, 197−205. (10) Petrova, M. A.; Kurteva, V. B.; Lubenov, L. A. Synergistic Effect in the Solvent Extraction and Separation of Lanthanoids by 4-(4Fluorobenzoyl)-3-methyl-1-phenyl-pyrazol5-one in the Presence of Monofunctional Neutral Organophosphorus Extractants. Ind. Eng. Chem. Res. 2011, 50, 12170−12176. (11) Matsumoto, M.; Otono, T.; Kondo, K. Synergistic Extraction of Organic Acids with Tri-N-Octylamine and Tri-N-Butylphosphate. Sep. Purif. Technol. 2001, 24, 337−342. (12) Osseo-Asare, K. Aggregation, Reversed Micelles, and Microemulsions in Liquid-Liquid Extraction: The Tri-N-butyl PhosphateDiluent-Water-Electrolyte System. Adv. Colloid Interface Sci. 1991, 37, 123−173.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.8b10831. Calculation of structure factor and Zimm distribution function, NMR spectrum of TBP saturated with water, RDFs obtained from simulations, TIP4P water model parameters, cluster size distribution of extractant mixtures, number of hydrogen bonded pairs in different P

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Association of Tri-N-Butyl Phosphates in N-Dodecane. J. Phys. Chem. B 2015, 119, 1588−1597. (31) Vo, Q. N.; Dang, L. X.; Nilsson, M.; Nguyen, H. D. Quantifying Dimer and Trimer Formation by Tri-N-Butyl Phosphates in NDodecane: Molecular Dynamics Simulations. J. Phys. Chem. B 2016, 120, 6985−6994. (32) Servis, M. J.; Tormey, C. A.; Wu, D. T.; Braley, J. C. A Molecular Dynamics Study of Tributyl Phosphate and Diamyl Amyl Phosphonate Self-Aggregation in Dodecane and Octane. J. Phys. Chem. B 2016, 120, 2796−2806. (33) Mu, J.; Motokawa, R.; Williams, C. D.; Akutsu, K.; Nishitsuji, S.; Masters, A. J. Comparative Molecular Dynamics Study on Tri-NButyl Phosphate in Organic and Aqueous Environments and Its Relevance to Nuclear Extraction Processes. J. Phys. Chem. B 2016, 120, 5183−5193. (34) Das, A.; Sahu, P.; Ali, S. M. Molecular Dynamics Simulation for the Calibration of the OPLS Force Field Using DFT Derived Partial Charges for the Screening of Alkyl Phosphate Ligands by Studying Structure, Dynamics, and Thermodynamics. J. Chem. Eng. Data 2017, 62, 2280−2295. (35) Murzyn, K.; Bratek, M.; Pasenkiewicz-Gierula, M. Refined OPLS All-Atom Force Field Parameters for N-Pentadecane, Methyl Acetate, and Dimethyl Phosphate. J. Phys. Chem. B 2013, 117, 16388−16396. (36) Siu, S. W. I.; Pluhackova, K.; Bockmann, R. Optimization of the OPLS-AA Force Field for Long Hydrocarbons. J. Chem. Theory Comput. 2012, 8, 1459−1470. (37) Baldwin, A. G.; Servis, M. J.; Yang, Y.; Bridges, N. J.; Wu, D. T.; Shafer, J. C. The Structure of Tributyl Phosphate Solutions: Nitric acid, Uranium (VI), and Zirconium (IV). J. Mol. Liq. 2017, 246, 225− 235. (38) Chandler, D. Interfaces and the Driving Force of Hydrophobic Assembly. Nature 2005, 437, 640. (39) Li, T.; Senesi, A. J.; Lee, B. Small Angle X-ray Scattering for Nanoparticle Research. Chem. Rev. 2016, 116, 11128−11180. (40) Baxter, R. Percus−Yevick Equation for Hard Spheres with Surface Adhesion. J. Chem. Phys. 1968, 49, 2770−2774. (41) Menon, S.; Kelkar, V.; Manohar, C. Application of Baxter’s Model to the Theory of Cloud Points of Nonionic Surfactant Solutions. Phys. Rev. A: At., Mol., Opt. Phys. 1991, 43, 1130. (42) Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.; Hess, B.; Lindahl, E. GROMACS: High Performance Molecular Simulations Through Multi-Level Parallelism from Laptops to Supercomputers. SoftwareX 2015, 1, 19−25. (43) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (44) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N log (N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (45) Hess, B. P-LINCS: A Parallel Linear Constraint Solver for Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 116−122. (46) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935. (47) Servis, M. J.; Wu, D. T.; Braley, J. C. Network Analysis and Percolation Transition in Hydrogen Bonded Clusters: Nitric Acid and Water Extracted by Tributyl Phosphate. Phys. Chem. Chem. Phys. 2017, 19, 11326−11339. (48) Widom, B. Some Topics in the Theory of Fluids. J. Chem. Phys. 1963, 39, 2808−2812. (49) Wand, C. R.; Totton, T. S.; Frenkel, D. Addressing Hysteresis and Slow Equilibration Issues in Cavity-Based Calculation of Chemical Potentials. J. Chem. Phys. 2018, 149, No. 014105. (50) Espinosa, J. R.; Wand, C. R.; Vega, C.; Sanz, E.; Frenkel, D. Calculation of the Water-Octanol Partition Coefficient of Cholesterol for SPC, TIP3P, and TIP4P Water. J. Chem. Phys. 2018, 149, 224501.

(13) Chiarizia, R.; Nash, K. L.; Jensen, M. P.; Thiyagarajan, P.; Littrell, K. C. Application of the Baxter Model for Hard Spheres with Surface Adhesion to SANS Data for the U (VI)- HNO3, TBP-NDodecane System. Langmuir 2003, 19, 9592−9599. (14) Nave, S.; Mandin, C.; Martinet, L.; Berthon, L.; Testard, F.; Madic, C.; Zemb, T. Supramolecular Organisation of Tri-N-butyl Phosphate in Organic Diluent on Approaching Third Phase Transition. Phys. Chem. Chem. Phys. 2004, 6, 799−808. (15) Motokawa, R.; Suzuki, S.; Ogawa, H.; Antonio, M. R.; Yaita, T. Microscopic Structures of Tri-N-Butyl Phosphate/N-Octane Mixtures by X-ray and Neutron Scattering in a Wide q Range. J. Phys. Chem. B 2012, 116, 1319−1327. (16) Sato, T. The Synergic Effect of Tri-N-Butyl Phosphate in the Extraction of Uranium (VI) from Sulphuric Acid Solutions by Di-(2ethylhexyl)-Phosphoric Acid. J. Inorg. Nucl. Chem. 1964, 26, 311−319. (17) Dourdain, S.; Hofmeister, I.; Pecheur, O.; Dufreche, J.; Turgis, R.; Leydier, A.; Jestin, J.; Testard, F.; Pellet-Rostaing, S.; Zemb, T. Synergism by Coassembly at the Origin of Ion Selectivity in Liquid− Liquid Extraction. Langmuir 2012, 28, 11319−11328. (18) Pecheur, O.; Dourdain, S.; Guillaumont, D.; Rey, J.; Guilbaud, P.; Berthon, L.; Charbonnel, M.; Pellet-Rostaing, S.; Testard, F. Synergism in a HDEHP/TOPO Liquid− Liquid Extraction System: An Intrinsic Ligands Property? J. Phys. Chem. B 2016, 120, 2814− 2823. (19) Ellis, R. J.; Anderson, T. L.; Antonio, M. R.; Braatz, A.; Nilsson, M. A SAXS Study of Aggregation in the Synergistic TBP−HDBP Solvent Extraction System. J. Phys. Chem. B 2013, 117, 5916−5924. (20) Anderson, T. L.; Braatz, A.; Ellis, R. J.; Antonio, M. R.; Nilsson, M. Synergistic Extraction of Dysprosium and Aggregate Formation in Solvent Extraction Systems Combining TBP and HDBP. Solvent Extr. Ion Exch. 2013, 31, 617−633. (21) Rey, J.; Dourdain, S.; Dufrêche, J.-F.; Berthon, L.; Muller, J.; Pellet-Rostaing, S.; Zemb, T. Thermodynamic Description of Synergy in Solvent Extraction: I. Enthalpy of Mixing at the Origin of Synergistic Aggregation. Langmuir 2016, 32, 13095−13105. (22) Flores, M. E.; Martínez, F.; Olea, A. F.; Shibue, T.; Sugimura, N.; Nishide, H.; MorenoVilloslada, I. Water-Induced Phase Transition in Cyclohexane/N-Hexanol/Triton X-100 Mixtures at a Molar Composition of 1/16/74 Studied by NMR. J. Phys. Chem. B 2017, 121, 876−882. (23) Lemyre, J.-L.; Ritcey, A. M. Characterization of a Reverse Micellar System by 1H NMR. Langmuir 2010, 26, 6250−6255. (24) Calandra, P.; Nicotera, I.; Rossi, C. O.; Liveri, V. T. Dynamical Properties of Self-Assembled Surfactant-Based Mixtures: Triggering of One-Dimensional Anomalous Diffusion in Bis-(2-Ethylhexyl) Phosphoric Acid/N-Octylamine Systems. Langmuir 2013, 29, 14848− 14854. (25) Knight, A. W.; Qiao, B.; Chiarizia, R.; Ferru, G.; Forbes, T.; Ellis, R. J.; Soderholm, L. Subtle Effects of Aliphatic Alcohol Structure on Water Extraction and Solute Aggregation in Biphasic Water/NDodecane. Langmuir 2017, 33, 3776−3786. (26) Marie, C.; Hiscox, B.; Nash, K. L. Characterization of HDEHPLanthanide Complexes Formed in a Non-Polar Organic Phase Using 31P NMR and ESI-MS. Dalton Trans. 2012, 41, 1054−1064. (27) Lumetta, G. J.; Neiner, D.; Sinkov, S. I.; Carter, J. C.; Braley, J. C.; Latesky, S. L.; Gelis, A. V. Combining Neutral and Acidic Extractants for Recovering Transuranic Elements from Nuclear Fuel; International Solvent Extraction Conference, Santiago, Chile, Oct 3− 7, 2011. (28) Beudaert, P.; Lamare, V.; Dozol, J.-F.; Troxler, L.; Wipff, G. Theoretical Studies on Tri-N-Butyl Phosphate: MD simulations in Vacuo, in Water, in Chloroform, and at a Water/Chloroform Interface. Solvent Extr. Ion Exch. 1998, 16, 597−618. (29) Cui, S.; de Almeida, V. F.; Hay, B. P.; Ye, X.; Khomami, B. Molecular Dynamics Simulation of Tri-N-Butyl-phosphate Liquid: A Force Field Comparative Study. J. Phys. Chem. B 2012, 116, 305−313. (30) Vo, Q. N.; Hawkins, C. A.; Dang, L. X.; Nilsson, M.; Nguyen, H. D. Computational Study of Molecular Structure and SelfQ

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (51) Johnson, W. F.; Dillon, R. Physical properties of Tributylphosphate-Diluent Solutions; U. S. Atomic Energy Commission Report HW-29086, 1953. (52) Naganawa, H.; Tachimori, S. Highly Hydrated Associate of Tributyl Phosphate in Dodecane. Anal. Sci. 1994, 10, 607−613. (53) Murray, B.; Axtmann, R. Water Content of Tributyl Phosphate by High Resolution Nuclear Magnetic Resonance. Anal. Chem. 1959, 31, 450−451. (54) Siddall, T.; Stewart, W. Metal Nitrate-Tributyl Phosphate and Related Exchange Studies by Means of Nuclear Magnetic Resonance. Inorg. Nucl. Chem. Lett. 1967, 3, 279−284. (55) Srivastava, S.; Andreev, M.; Levi, A. E.; Goldfeld, D. J.; Mao, J.; Heller, W. T.; Prabhu, V. M.; De Pablo, J. J.; Tirrell, M. V. Gel Phase Formation in Dilute Triblock Copolyelectrolyte Complexes. Nat. Commun. 2017, 8, 14131. (56) Pedersen, J. S.; Svaneborg, C.; Almdal, K.; Hamley, I. W.; Young, R. N. A Small-Angle Neutron and X-ray Contrast Variation Scattering Study of the Structure of Block Copolymer Micelles: Corona Shape and Excluded Volume Interactions. Macromolecules 2003, 36, 416−433. (57) Pardo, J. M.; Tovar, C. A.; González, D.; Carballo, E.; Romaní, L. Thermophysical Properties of the Binary Mixtures Diethyl Carbonate+(N-dodecane or N-tetradecane) at Several Temperatures. J. Chem. Eng. Data 2001, 46, 212−216. (58) Zeppieri, S.; Rodríguez, J.; López de Ramos, A. Interfacial Tension of Alkane+ Water Systems. J. Chem. Eng. Data 2001, 46, 1086−1088. (59) Tian, Q.; Liu, H. Densities and Viscosities of Binary Mixtures of Tributyl Phosphate with Hexane and Dodecane from (298.15 to 328.15) K. J. Chem. Eng. Data 2007, 52, 892−897. (60) Biswas, R.; Banu, R.; Islam, M. Some Physico-Chemical Properties of D2EHPA: Part 2. Distribution, Dimerization and Acid Dissociation Constants in N-hexane/1 M (Na+, H+) SO42- System, Interfacial Adsorption and Excess Properties. Hydrometallurgy 2003, 69, 157−168. (61) Cui, S.; de Almeida, V. F.; Khomami, B. Molecular Dynamics Simulations of Tri-N-Butyl Phosphate/N-Dodecane Mixture: Thermophysical Properties and Molecular Structure. J. Phys. Chem. B 2014, 118, 10750−10760. (62) Gullekson, B. J.; Breshears, A. T.; Brown, M. A.; Essner, J. B.; Baker, G. A.; Walensky, J. R.; Paulenova, A.; Gelis, A. V. Extraction of Water and Speciation of Trivalent Lanthanides and Americium in Organophosphorus Extractants. Inorg. Chem. 2016, 55, 12675−12685. (63) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (64) Jorgensen, W. L.; Blake, J. F.; Buckner, J. K. Free Energy of TIP4P Water and the Free Energies of Hydration of CH4 and ClFrom Statistical Perturbation Theory. Chem. Phys. 1989, 129, 193− 200. (65) Hermans, J.; Pathiaseril, A.; Anderson, A. Excess Free Energy of Liquids from Molecular Dynamics Simulations. Application to Water Models. J. Am. Chem. Soc. 1988, 110, 5982−5986. (66) Deschamps, J.; Costa Gomes, M.; Pádua, A. A. Solubility of Oxygen, Carbon Dioxide and Water in Semifluorinated Alkanes and in Perfluorooctylbromide by Molecular Simulation. J. Fluorine Chem. 2004, 125, 409−413. (67) Yang, L.; Ahmed, A.; Sandler, S. I. Comparison of Two Simulation Methods to Compute Solvation Free Energies and Partition Coefficients. J. Comput. Chem. 2013, 34, 284−293. (68) Bhatnagar, N.; Kamath, G.; Chelst, I.; Potoff, J. J. Direct Calculation of 1-Octanol-Water Partition Coefficients from Adaptive Biasing Force Molecular Dynamics Simulations. J. Chem. Phys. 2012, 137, No. 014502. (69) Garrido, N. M.; Economou, I. G.; Queimada, A. J.; Jorge, M.; Macedo, E. A. Prediction of the N-Hexane/Water and 1-Octanol/ Water Partition Coefficients for Environmentally Relevant Compounds Using Molecular Simulation. AIChE J. 2012, 58, 1929−1938.

(70) Ogata, K.; Hatakeyama, M.; Nakamura, S. Effect of Atomic Charges on Octanol−Water Partition Coefficient Using Alchemical Free Energy Calculation. Molecules 2018, 23, 425. (71) Bannan, C. C.; Calabró, G.; Kyu, D. Y.; Mobley, D. L. Calculating Partition Coefficients of Small Molecules in Octanol/ Water and Cyclohexane/Water. J. Chem. Theory Comput. 2016, 12, 4015−4024.

R

DOI: 10.1021/acs.jpcb.8b10831 J. Phys. Chem. B XXXX, XXX, XXX−XXX