Liquid–liquid equilibrium measurements and computational study of salt–polymer aqueous two phase system for extraction of analgesic drugs

In recent decades, aqueous two phase systems have gained a lot of attention for extraction of different materials. In this work, an aqueous two phase system was made by polyethylene glycol 600 and potassium hydroxide and phase diagram were determined for this system. The experimental binodal data were described using two empirical nonlinear three parameter expressions developed by Merchuk and Zafarani-Moattar. The consistency of the experimental tie-line data was determined by utilizing the Othmer-Tobias, Bancraft, and Setschenow correlations. Also, the extraction of two analgesic drugs, namely ibuprofen and acetaminophen were investigated by the mentioned ATPS. For this purpose, partition coefficients and extraction efficiencies of each drug were calculated. The trend of extraction efficiencies indicated that the responsibility of extraction of the mentioned drugs into the polymer-rich top phase is related to their hydrophobicity. The Diamond-Hsu equation and its modified version were used to correlate the drugs experimental partition coefficients. Furthermore, the interaction of mentioned drugs with polyethylene glycol was investigated employing quantum computing techniques based upon density functional theory (DFT). These results were in good agreement with the trend of extraction efficiencies of studied drugs.


Scientific Reports
| (2022) 12:13848 | https://doi.org/10.1038/s41598-022-18122-x www.nature.com/scientificreports/ Also, the application of this ATPS was studied for partitioning of ibuprofen and acetaminophen as widely used analgesic drugs. The partitioning coefficients, K, and the corresponding extraction efficiency, EE %, at different tie-lines were calculated to the investigation of effect of drug nature on the separation of mentioned drugs. The Diamond-Hsu 17 equation and its modified form were used to fit the partitioning coefficient values.
Finally, this paper has investigated the aspects of binary systems, including PEG + Ibuprofen, PEG + Acetaminophen, KOH + Ibuprofen, and KOH + Acetaminophen. In this case, the geometry optimization of all possible configurations of systems incorporating these systems was computed utilizing the density functional theory (DFT) approach. From an experimental and comprehensive theoretical standpoint, including natural bond orbital (NBO), quantum theory of atoms in molecules (QTAIM), and non-covalent interaction (NCI) studies, have been applied to elucidate the posture of extraction and phase separation of these pharmaceutical compounds.

Experimental section
Materials. All data corresponding to the chemicals used in this study have been shown in Table 1. The poly(ethylene glycol) with a molar mass of 600 g mol −1 and potassium hydroxide with a minimum purity of 99% were purchased from Merck. Ibuprofen and acetaminophen (> 99.5 wt%) were obtained from Danna pharmaceutical company (Iran). All materials were used without further purification and double distilled and deionized water was used in preparation of solutions.

Methods. Determination of phase diagrams and tie-lines.
The cloud point titration method 18,19 was used to create the phase diagrams at 298 K and under atmospheric pressure (85 kPa). In this method dropwise additions of potassium hydroxide aqueous solution to PEG 600 solution were repeated until a cloudy solution (biphasic region) was detected, then water was added dropwise until a clear and monophasic region was observed. To determine the binodal curves, aqueous solutions of PEG 600 at 60 wt% and aqueous solutions of potassium hydroxide at 50 wt% were used. Throughout the experiment, the solution was constantly stirred. During the preparation of the solutions water contents of chemicals, which were measured by Karl Fischer coulometer (Metrohm 751 GPD) were taken into account. An analytical balance (Shimadzu, 321-34553, Shimadzu Co., Japan) with a precision of ± 1•10 -7 kg was used to determine the mass fraction of the components. The maximum uncertainty in determining the mass fraction of polymer and salt was found to be 0.002.
Five different solutions of studied ATPS consist of {PEG 600 + KOH + water} were prepared to determination of tie-lines. For each tie-line, a mixture of salt, polymer and water at the biphasic region was gravimetrically prepared within ± 10 -7 kg. These solutions were vigorously stirred for 30 min, centrifuged and placed in water bath at 298 K to reach equilibrium. The compositions of coexisting phases were analytically determined. After separation of the two phases, the concentrations of KOH in the top and bottom phases were determined by flame photometer (JENVEY model PFP7, England). The concentration of PEG in both phases was determined by refractive index measurements 20 performed at T = 298.15 K using a refractometer (ATAGO DR-A1, Japan) with a precision of ± 0.0001. The uncertainty in refractive index measurement is ± 0.0002. For dilute aqueous solutions containing a polymer and a salt, the relation between the refractive index, n D , and the mass fractions of polymer, w p , and salt, w s is given by: here n 0 is the refractive index of pure water for which the value 1.3325 at T = 298.15 K was obtained. The values of constants a p , and a s corresponding to polymer, and salt respectively were obtained from the linear calibration plots of the corresponding refractive index of the solutions diluted in mass fraction range (C Range(w/w)) as showed in Table S1. These constants together with the coefficient of determination values, R 2 , are reported in Table S1.
Partition coefficient and the corresponding extraction efficiency. To assess the performance of the investigated ATPS in terms of drug partitioning and extraction efficiency, the five mixture compositions (reported in Table 2) were chosen based on the phase diagrams determined before for {PEG 600 + KOH + H 2 O} system. These mixtures had similar overall compositions used to obtain tie-lines. The prepared mixtures were vigorously stirred for 30 min before being placed in a water bath at the desired temperature. After observation of phase separation 1 mL of each phase was carefully separated, mixed again and then 0.002 mass fractions of drugs (ibuprofen and acetaminophen) was added to each sample. The samples were centrifuged at 2,000 rpm for 10 min and kept in a water bath for 24 h to ensure complete phase separation and equilibrium. Finally, after careful separation of phases in equilibrium, the concentrations of acetaminophen in both phases were determined by using a UV (1) n D = n 0 + a p w p + a s w s www.nature.com/scientificreports/ spectroscopy using a spectrophotometer mentioned above at the wavelengths of (284, 261 and 303) nm. Also, for analysis of ibuprofen a fluorescence spectrophotometer (F-4500, Hitachi, Tokyo, Japan) was used. To avoid interference from the phase components, the samples were diluted and analyzed against the blanks containing the same phase components but without drugs 21 . The partitioning coefficient, K, and extraction efficiency, EE %, were calculated respectively by Eqs.
(2) and (3) as below: here wtop drug and wbot drug are the mass fractions of drug in the top and bottom phases respectively.
Determination of the tie-lines length. The tie-line length (TLL) has unit of %w/w, same as the component concentrations. The TLL can be related to the equilibrium phase composition as follows: where w p and w s are equilibrium compositions of PEG 600 , and salt, respectively. Superscripts "bot" and "top" designate the PPG 600 -rich phase (top phase) and salt-rich phase (bottom phase), respectively.

Computational
The interaction of Ibuprofen and acetaminophen with PEG was investigated employing quantum computing techniques based upon density functional theory (DFT). Noteworthy to mention that the investigation of intramolecular interactions in systems including pharmaceutical compounds with KOH and PEG indicates that drugs with PEG have stronger intermolecular interactions, including hydrogen interactions toward drugs with KOH, which would be the results of extraction and separation of these compounds confirm this pattern. In this regard, a computational investigation of binary systems comprising drugs and PEG is being investigated in this research effort. All electrical calculations have been performed using the Gaussian 09 set of tools. The natural bond orbitals (NBO), quantum theory of atoms in molecules (QTAIM), non-covalent interaction (NCI), and energetic parameters were estimated following the optimization of the ibuprofen, acetaminophen and PEG structures. The NBO function in the Gaussian 09 software package was used for all NBO computations. QTAIM investigations were also carried out using Multiwfn software.
The B3LYP functional, a hybrid functional with reliable and consistent performance in investigating diverse compounds' electronic structures and characteristics, was utilized to improve the systems, including Ibuprofen, acetaminophen and PEG. In addition, Grimm's third version of the empirical correction with the Bijerlin effect was included to account for the effect dispersion effect, and the B3LYP-D3 (BJ) functional was employed to optimize the structures. The Multiwfn software has been employed to optimize the most stable structures, further analyzed using the QTAIM approach 22 .
Based on the sign of the eigenvalues of the diagonalited hessian matrix of the Laplacian of electron density ( ∇ 2 ρ(r) ), the critical points (i.e. where the divergence of ρ is equal to zero) were calculated and used to characterize the critical points as the atom, bond, ring, and cage critical points in QTAIM analysis 23,24 . Bond critical points (BCP) are nearly (3,−1) critical points that are used to detect and characterize potential bond interactions. In this analysis, different factors can be used to decide about the magnitude of ρ the magnitude and the sign of ∇ 2 ρ(r) , the ratio of kinetic (G(r)) to potential energy (|V (r)|) and the value of the total energy (H(r)) and also the ellipticity of the electron density ( ε ) are some of the applicable parameters.
Because the B3LYP functional is somewhat unsuccessful in determining weak van der Walls interactions due to its asymptotic behavior, single-point energy calculations on optimized structures were performed using other functional such as M06-2X and WB97X-D3 to explore the interaction energies more consistently. which E int represents the total interaction of the (Drugs-PEG) systems, E(A) denotes the Drugs energy, and E(B) includes the energy of systems with PEG. The compounds PEG + Ibo and PEG + ACE with the most energetically stable structure were selected. Furthermore, at the B3LYP/6-311 G (d) level of theory, some quantitative conceptual parameters of natural bond orbitals (NBO), quantum theory of atoms in molecules (QTAIM), and Reduced Density Gradient (RDG) analysis were used to confirm the presence of intramolecular interaction, including hydrogen bonds. The Wiberg bond index (WBI), which is defined using the following equation, is another idea that may be used to analyze bond orders 25 : The density matrix elements are denoted by pjk, whereas the atomic orbital charge density is shown by pjj. Furthermore, the energy of donor-acceptor interactions in the NBO basis might be estimated using the secondorder perturbation technique as follows 26 : in which b ∧ F is the effective orbital Hamiltonian, and ε i and εj* are defined as

Results and discussion
Phase diagram. The phase diagrams of studied ATPS composed of PEG 600 , KOH and water at 298.15 K are plotted in Fig. 1, and the experimental mass fraction data are also listed in Table S2.
In above equations w p and w s are mass fractions of polymer and salt, respectively. The adjustable parameters of Eqs. (5) and (6) {(a, b, and c) and ( α , β and γ )} together with the corresponding standard deviation, sd, are presented in Table S3. The obtained values of sd indicate a good performance of both equations in regeneration of binodal data.
Analysis of tie-lines. The tie-lines together with binodal curves are presented in Fig. 1 for ATPS composed of PEG 600 , KOH, and water. Five experimental tie-lines and the corresponding length (TLL) were determined  Tie-line correlation. In this work for the correlation of LLE data of (PEG 600 + KOH + water) system several models with different thermodynamics basis were used and compared their performances in the correlation of the tie-lines.
Othmer-Tobias and Bancraft equations. The correlation equations given by Othmer-Tobias 14 , Eq. (10), and Bancraft 15 , Eq. (11), have been used to correlate the tie-lines composition: where k,n , k 1 , and r represented fit parameters. These equations have also been used to assess the reliability of LLE data. Using the tie-line data reported in Table 2 Table S4. On the basis of the obtained standard deviations (sd), we conclude that Eqs. (10) and (11) can be satisfactorily used to correlate the tie-line data of the investigated system. Setschenow equation. The following Setschenow equation proposed by Hey and coworkers 16 has also been used for the correlation of tie-line data of studied systems: in which the k s is the salting-out coefficient, k p is a constant, and C p and C s are the molality of polymer and salt respectively. The parameters of the Eq. (12) which were obtained from the correlation of the experimental LLE data are also given in Table S4 along with the corresponding deviations. According to the sd values in Table S4, it could be estimated that correlation with Setschenow equation is satisfactory.
The performance of salt -PEG ATPS in the partition behavior of drugs. The experimental partition coefficients, K, and the percentage extraction efficiencies, EE %, of investigated drugs (ibuprofen and acetaminophen) are presented in Table 3 and their variation with TLL can be seen in Figs. 2 and 3, respectively. The  www.nature.com/scientificreports/ observed trend is as follows: K (ibuprofen) > K (acetaminophen). It is observed from Table 3 that drug effectively partitioned to the PPG-rich phase (more hydrophobic phase).
The effect of drug structure on the drug partitioning behavior. The chemical structures and polarities of drug are responsible for the observed trend of the K and EE % values ( Table 3). The studied drugs because of their non-polar properties have favorable interactions with hydrophobic phase (PEG-rich phase), and therefore, partitioning of these drugs to the salt-rich phase occurs significantly less frequently. In fact, the presence of nonpolar groups in drug structure increases its tendency to PEG-rich top phase. The obtained high values of EE % collected in Table 3 imply that both of the investigated drugs partitioned mostly in the PEG-rich top phase. The observed trend for the drugs K values is in good agreement with their hydrophobicity values defined by the logarithm of octanol/water partition coefficient (log K ow ); so that higher log K ow values show that the solute have more hydrophobicity and lower tendency to water molecules 18 . The log K ow values of drugs are 3.97 27 and 2.34 28 for ibuprofen and acetaminophen, respectively. Therefore the drug extraction efficiency is in the order: ibuprofen > acetaminophen. This study and previous one 29 indicate that the hydrophobicity of drugs can be an important criteria to predict of drug partitioning between two-phases of an ATPS.  13) and (14). The symbol �w(PEG 600 ) is used for the mass fraction difference of PEG 600 in the top and bottom phase.
For the studied drugs the results of fitting experimental partition coefficients values to Eqs. (13) and (14) together with the corresponding standard deviations, sd are reported in Table S5. According to the sd values in Table S5, it could be estimated that correlation with both of the equations is satisfactory.

Computational section. The optimized structures and energetic parameters.
The initial structures of the drugs and PEG have been optimized in various configurations to determine the most stable compounds with the minimum relative energy in binary systems based on PEG + acetaminophen and PEG + ibuprofen, which can be seen in Table 4. Furthermore, the energy of various configurations of these binary systems has been provided in the supplementary information from Tables S6 and S7.
The various configurations for each binary system (6 and 8) were considered. In this regard, a total of 6 and 8 distinct binary complexes, comprising drugs and PEG, were evaluated and geometrically optimized.
The obtained results show that the PEG + ibuprofen systems have the lowest structural energy (the most negative value) in comparison to the other compounds in all three functional studies, indicating that the investigated drugs + PEG are the most stable in terms of energy and structure.
Besides addition, based on the QTAIM, NBO, and NCI analyses, the most stable acetaminophen + PEG and ibuprofen + PEG complexes were utilized to evaluate.
The (QTAIM) analyses. The quantum theory of atoms in molecules (QTAIM) would be another critical and prevalent technique. By employing Multiwfn software, the most stable optimized structures were further evaluated using the QTAIM approach. In investigating the nature of possible intermolecular interactions, the QTAIM was also researched on the studied binary system.
Bader's QTAIM method 30,31 , which uses specific definitions including electron density ( ρ(r) ) and Laplacian of electron density ( ∇ 2 ρ(r) ), is a sophisticated analysis method for determining the nature of various sorts of interactions 32 . Critical points (i.e., a minimum, maximum, or saddle point in the spatial curvature of ( ρ(r) ) are discovered in QTAIM analysis when the gradient of electron density ( ∇ 2 ρ(r) ) is identically zero, and the Hessian matrix of ( ρ(r) ) is calculated and utilized to characterize critical points. The eigenvalues ( 1 ), ( 2 ), and ( 3 ) of the diagonalized Hessian matrix would be used to categorize the critical locations (CPs).
To evaluate whether a weak or strong attractive intermolecular or interatomic bonding interaction exists, the presence of a bond critical point (BCP) is crucial 23,24 .
Chemical interactions could be classified into various classes based on the ρ(r) and ∇ 2 ρ(r) amounts obtained by the QTAIM analysis at the BCP. As a outcomes, a negative ∇ 2 ρ(r) value and a ρ(r) value more than 0.1 a.u. at a BCP identify the BCP as a shared (covalent) bonding 33 . However, the bonding interactions can be classified as non-covalent, strong ionic, or weak Van der Waals for positive ∇ 2 ρ(r) values.
The quantum virial theorem [33][34][35] could be a significant relationship between energetic QTAIM parameters and P1 at CPs, as evidenced by the following relationship: where V (r) are G(r) the potential and kinetic energy densities of electrons, respectively.
Applying Multiwfn software, the most stable optimized structures were further analyzed using the QTAIM method. Table 5 exhibit binary systems including acetaminophen + PEG and ibuprofen + PEG complexes,  Table 5. The AIM topological parameters, including electron density (ρ), Laplacian of electron density (∇ 2 ρ(r)), the kinetic electron density G(r), potential electron density V(r), eigenvalues of Hessian matrix (λ) and bond ellipticity index (ε) of binary systems containing drugs and PEG.

System
Interaction ρ(r) The binary interaction between ibuprofen-PEG and Ace PEG has 0.03 ≤ ρ(r) ≤ 0.05, positive ∇ 2 ρ(r) and 0.5 ≤ G(r) |V (r)| ≤ 1 and negative values, indicating that all interactions are mixed ionic-covalent. As an outcome, the results for ibuprofen-PEG and acetaminophen-PEG were reported as 0.0.0395 and 0.0342 kg/m 3 , respectively, as shown in Table 5. Additionally, the values are also used to determine the stability of the intramolecular interaction, with ε 1 indicating the most stable structure. In binary systems, the ε values for ibuprofen-PEG and Acetaminophen-PEG, for instance, have been reported as 0.0373 and 0.0543, respectively. The obtained G(r) |V (r)| values are shown in Table 5 and Fig. 4, reported as 0.9300 and 1.0101 for the ibuprofen-PEG and acetaminophen-PEG, respectively.

The (NBO), (WBI) values and (N.C.I) analyses.
Weinhold et al. 36 outline the scope of the natural bond orbital (NBO) analysis applied in this work. NBO methods have a mathematical and historical background that could be discovered elsewhere. The second-order perturbation theory defined H-bond strength as delocalization energy E2 (e.g., electron transfer from donor to acceptor orbital) because of the NBO conception. The most stable structures of binary systems were also applied to NBO analysis, with the outcomes presented in Table 6 for some donner-acceptor interactions. It is noteworthy that the interactions with the E 2 values are included in Table 6. In all circumstances, the most significant interaction is due to the transfer of lone pair electrons (LP) from the π orbitals of O of drugs to the antibonding H-O of PEG substances, as seen in this table. According to Table 6, the ibuprofen-PEG with the highest E 2 values is the most stable structure compared to the other binary systems. Second-order perturbation energies (E 2 ) for ibuprofen-PEG. For instance, the data has been reported as being of the order of 16.12 kcal/mol for the ibuprofen-PEG binary system. The natural charge (N.C.I) and Wiberg bond indices (WBI) values were also employed to investigate the interactions further. Furthermore, the natural charge (N.C.I) and Wiberg bond indices (WBI) values were utilized to evaluate the interactions further. The Table 7 indicates the values of the WBI (bond length (re) and Wiberg bond index of the most significant    34,37,38 . The values of WBI are in the region of strong non-covalent interactions, according to the reported data in Table 7. Furthermore, the findings of the NCI analysis are as follows: This investigation utilized NCI analysis to investigate the weak non-covalent interactions outlined in the QTAIM and NBO sections. The reduced density gradient (RDG) analysis is a method for visualizing and evaluating noncovalent interactions in real space, including as van der Waals contacts, hydrogen bonds, and steric effects, utilizing electron density and its derivatives 39 . One of the other descriptor for evaluating intermolecular interactions appears to be the reduced density gradient (RDG), which is defined as follows 39,40 : Realizing that the strength of the interaction illustrated the good correlation with sign 2 (r)ρ(r) , a scatter plot of the RDG as a function of sign 2 (r)ρ(r) could be applied to recognize the nature of the polymer and drugs interaction.
NCI analysis results in this respect the scatter plot of RDG values were calculated using Multiwfn software and the RDGs versus sign( 2 )ρ(r) were plotted. The results were illustrated for binary complexes in Fig. 5. The existence of narrow spikes in the region of sin (r) 0 , for is surfaces of RDG ˂0.5, the existence of reveals strong non-covalent H-O and hydrogen bonding in all cases. The results also are in full agreement with QTAIM results, the Van der Waals interaction (i.e. where ρ ≈ 0,∇ 2 ρ ≈ 0 ) and the repulsive steric effects (i.e. where 2 (r) 0 ). This conclusion can also be made from QTAIM and NBO analysis, where in both cases the antisymmetry of the similar interactions is obvious.

Conclusions
Phase diagrams of an ATPS composed of polyethylene glycol with a molar mass of 600 g mol −1 , and potassium hydroxide at 298.15 K were studied. The experimental binodal values were satisfactorily adjusted by Merchuk and Zafarani-Moattar et al. two semi-empirical equations. To correlation of tie-lines in the studied system various equations (Othmer-Tobias, Bancraft and Setschenow) were used. Based on the obtained outcomes, it was determined that the performances of all the analyzed equations in the correlation are satisfactory. Additionally, partition coefficient and extraction efficiency for two drugs namely ibuprofen and acetaminophen on the studied ATP were obtained. The trend of the extraction efficiency showed that these drugs due to their hydrophobic properties preferentially partition into the top phase (polymer-rich).
Finally, density functional theory (DFT), QTAIM and NBO calculations were employed to study binary systems' intramolecular and structural properties. Different kinds of DFT functional, including B3LYP-D3 (BJ), M062X and WB97X-D3 combined with two 6-311G (d,p) basis sets were used. The formation energies confirmed that B3LYP-D3 (BJ) is the most suitable functional group among the other functional groups. Also, it could be seen that the ibuprofen-PEG is more stable than the acetaminophen-PEG. The interactions with the highest values of E 2 obtained from NBO analysis showed the most significant interaction in all cases is due to the donation of the lone pair electrons (LP) from π orbitals of O to the antibonding H of PEG. QTAIM topology analysis showed the positive values for the Laplacian of electron density, which confirms the strong interactions are due to the interaction of the O atom of drugs with H (-O) of PEG. These changes are more intensive in the ibuprofen-PEG binary system than in acetaminophen-PEG. Since the interaction ibuprofen-PEG could change www.nature.com/scientificreports/ the intermolecular interaction trends more than acetaminophen-PEG, it could be more beneficial for drug separation in the pharmaceutical and chemical industries.