A numerical approach to determine the optimal condition of the gas anti-solvent supercritical process for nanoparticles production

Supercritical gas antisolvent (GAS) process is an efficient method for nanoparticles production, in which accurate selection of operational condition is essential. Thermodynamic models can be applied for evaluation the phase equilibrium behavior and determination the required precipitation pressure of GAS process. In this research, thermodynamic behavior of (CO2—dimethyl sulfoxide (DMSO)) binary system and both of (CO2–DMSO-anthraquinone Violet 3RN (AV3RN)) and (CO2–DMSO-solvent Yellow 33 (SY33)) ternary systems in the GAS process were studied at different temperatures (308, 318, 328 and 338) K and pressures (1.0–14.0) MPa, using Peng–Robinson equation of state (PR-EoS). The minimum precipitation pressure of AV3RN and SY33 at 308, 318, 328 and 338 K were 7.80, 8.57, 9.78 and 11 MPa and 8, 8.63, 9.5 and 10.77 MPa, respectively. Also, the mole fraction of substances in liquid phase of ternary systems were determined by PR-EoS, at 328 K versus pressure. The accuracy of the obtained results were investigated using the experimental data reported in the literatures.


List of symbols
Heat of fusion P min Minimum required pressure for precipitation in GAS process AV3RN Anthraquinone Violet 3RN SY33 Solvent Yellow 33 Greek symbols α(T r ω), Temperature-dependent function for the considered parameter of Peng Robinson EoS φ Fugacity coefficient ω Acentric factor Δ Property change www.nature.com/scientificreports/ SY33 (Quinoline Yellow) is an anionic dye applied as a colorant in textile, plastic, medicine, cosmetic, rubber and many other industries. Generally, textile dyeing is known as a process with remarkable water consumption and also it is the second largest polluter of water globally 21 . Alternatively, utilization of sc-CO 2 as dyeing medium, instead of water, provide great opportunities for significant water savings along with reducing the amount of discharged wastewater to the environment leading to decrement the environmental pollution. Due to necessity of dissolving dyes in sc-CO 2 phase for penetration to textile fibers, supercritical dyeing (SCD) process efficiency is significantly depends on dye solubility in sc-CO 2 . It has been confirmed in several researches that production of nanoparticles with high surface area is the most efficient approach for solubility enhancement of different components [7][8][9]22,23 . Solubility of AV3RN in sc-CO 2 with/without methanol as co-solvent has been measured in temperatures range of 308-348 K and pressures range of 100-350 bar in our previous work 24 . It was shown that AV3RN solubility in the binary (Sc-CO 2 + AV3RN) and ternary (Sc-CO 2 + AV3RN + methanol) systems were in the range of 4.7 × 10 −7 to 5.5 × 10 −6 and 4.4 × 10 −6 to 5.8 × 10 −5 mol fraction, respectively. Solubility of SY33 in sc-CO 2 has not yet been reported, but according to available data for solubility of some other similar structure dyes 25 , it can be concluded that solubility of SY33 in sc-CO 2 is also poor. Based on poor solubility of both substances in sc-CO 2 , GAS anti solvent process can be one of the best techniques for production of AV3RN and SY33 nanoparticles.
In this semi batch process, first the solute is completely dissolved in an organic solvent, then sc-CO 2 is injected into this stationary liquid phase with special flow rate. Sc-CO 2 dissolution into organic solvent leads to a large volume expansion and reduction of liquid density. This is associated with a decrement in solvation power of organic solvent and a sharp increment of supersaturation degree within the liquid phase. Degree of supersaturation which defines as the ratio of solute concentration to its equilibrium solubility, is the driving force of precipitation in GAS process, impressing the precipitates nucleation and growth rates 26 . Therefore, severe reduction of solid solute solubility causes to precipitation of pure solute nanoparticles which can be separated from the expanded liquid phase by high pressure filtration. More purity, narrow particle size distribution and easier quality control in terms of size, morphology and crystallinity of deposited nanoparticles along with lower operational temperature and negligible levels of particles contamination to organic solvent are the most important preferences of GAS process 12,27 . Due to direct effect of operational condition on equilibrium solubility and supersaturation of solution and consequently on quality of deposited particles, it can be concluded that nanoparticles with desired characteristics are not produced at arbitrary conditions and selection the appropriate operational conditions is so important.
Generally, thermodynamic modelling can be helpful in reduction of required time and energy of experimental studies. Simultaneous presence of solid solute, organic liquid solvent and sc-CO 2 anti solvent with a certain solubility of solute and sc-CO 2 in the solvent is one of the thermodynamic aspects of the GAS process. Additionally, phase equilibrium of antisolvent-solvent mixture has determinative role in controlling the particles deposition in this process. So, the optimal conditions of nanoparticles formation in the GAS process should be determined by volume expansion and thermodynamic modeling. Also, thermodynamic modelling can be helpful for detailed understanding the phase behavior of the ternary system (gas/liquid/ solid), analyzing the effect of operational conditions and finally evaluation the feasibility and optimization of GAS process. Up to now, different thermodynamic models using an equation of state (EoS) for calculation the fugacity of three phases, have been applied for determination the solute and sc-CO 2 solubilities in liquid solvent and calculation the volume expansion of liquid phase. Peng-Robinson (PR) with conventional quadratic mixing rules (vdW2) is one of the most popular EoS used for this purpose 28,29 .
In current work, for the first time, the optimal operational conditions of AV3RN and SY33 nanoparticles precipitation via GAS process specified through modeling the volume expansion and phase equilibrium via PR-EoS with vdW2 mixing role. In this system, AV3RN and SY33, Dimethyl sulfoxide (DMSO) and sc-CO 2 were considered as solutes, solvent and antisolvent, respectively. DMSO is one of the commonly used solvent in supercritical antisolvent processes, due to its high miscibility with sc-CO 2 [30][31][32] . Also, remarkable high volumetric expansion of DMSO with sc-CO 2 is reported in several researches 29,33 . Theoretical studies Thermodynamic framework. For feasiblity analysis of AV3RN and SY33 nanoparticles production via GAS process, determination the optimum process condition and understanding the solute precipitation mechanism are performed by thermodynamic modeling. In this research, the Peng-Robinson equation of state (PR-EoS) was utilized for investigation the phase equilibrium and volume expansion of the liquid phase.
Phase equilibrium and volume expansion analysis. In this study, sc-CO 2 as antisolvent (1) in gas phase (g), DMSO as solvent (2) in liquid phase (l) and AV3RN and SY33 as solutes (3) in solid phase (s) are in equilibrium with each other. The equilibrium concentration of each of these components is an effective thermodynamic parameter which determines the maximum accessible supersaturation in the system. Equality of pressure, temperature and fugacity of these components in all the phases are the equilibrium criteria of these phases. Also, due to mixing of the liquid and gas phases, the mass transfer resistance between these two phases can be ignored 30 . Thus, the equilibrium condition can be explained by the following equilibrium equations: www.nature.com/scientificreports/ where φ i , y i and x i are the fugacity coefficient and equilibrium mole fraction of component "i" in gas and liquid phases, respectively. In Eq. (4), it was assumed that the solid phase is pure and the solubility of the solvent and the antisolvent in this phase is insignificant. Obviously, the sum of components mole fractions in both of the liquid (x i ) and gas (y i ) phases are equal to one 34 . Accordingly, a system of six equations with six unknowns of molar compositions in liquid and gas phases is obtained.
The following relationship along with an appropriate equation of state can be used for calculation the fugacity coefficient of each component in the liquid and gas phases: Peng-Robinson equation of state (PR-EoS) have been used as follows: where P, T and R are the absolute pressure (MPa), temperature (K) and the universal gas constant (8.314 J mol −1 K −1 ), respectively. Also, ν shows the solute molar volume (m 3 mol -1 ) which was computed by Fedors method 35 . Moreover, a(T) and b terms account for interactions between the species in the mixture and excluded volume of the mixture species, respectively. Generally, these two parameters are obtained via quadratic mixing rules in terms of the equilibrium compositions. The van der Waals mixing rules for a(T) and b is given as: where x is the mole fraction. Also, a ij (T), the cross energetic parameter, and b ij can be calculated as follows: where k ij and l ij are the binary interaction parameters; a(T) and b parameters are dependent on the critical and physical properties of the components: where where ω is acentric factor and T r is the reduced temperature (T r = T/T c ). However, the PR-EoS is inappropriate for modeling the phase behavior of solid phase. So, modified form of the proposed equation by de la Fuente Badilla et al. 36 is applied for calculation of solid fugacity coefficient ( ϕ s 3 ): This equation relates the fugacity coefficient of the solid solute ( ϕ s 3 ) to the fugacity of the sub-cooled liquid ( ϕ l 3 ) at temperature T and pressure P, which can be calculated by the PR-EOS. For calculation of ϕ s 3 , the heat of fusion at the triple point ( H tp ), the triple point temperature (T tp ), the triple point pressure (P tp ) and the molar volume of the solute at the triple point (v tp ) should be specified. These characteristics for AV3RN and SY33 components are listed in Table 2.
The volume expansion of liquid phase has a decisive role in GAS process and computation of this parameter is essential for determination of optimum process condition. The definition of the relative molar volume change is given as follows 30,36 : www.nature.com/scientificreports/ In this relation, v(T, P) is the molar volume of the liquid solution at the system temperature and pressure and v 0 (T, P 0 ) stands for the molar volume of the pure solvent at the system temperature and reference pressure (P 0 ) (usually atmospheric pressure). If CO 2 dissolution accompanied with reducing the molar volume, the sign of this parameter became negative and vice versa.
Optimization algorithm: particle swarm optimization (PSO). Process optimization is an important part of a process design where input process parameters are optimized to reach satisfactory output parameters. Numerous optimization techniques have been developed for different optimization problems. Particle Swarm Optimization, commonly referred as PSO, is one of the most commonly utilized optimization algorithms for solving continuous nonlinear optimization problems, firstly introduced by Kennedy and Eberhart 37 . Generally, it is a computational method that optimizes an objective function of mathematical problem by iteratively trying to improve a candidate solution with regard to a given measure of quality 38 . PSO is a popular technique with very easy execution only needs a few programming code lines. Also, simplicity of required mathematical operators makes it computationally economical in terms of both memory requirements and speed. PSO algorithm is created by natural swarms that are formed of volume-less particles with random velocities, each of which shows a feasible solution. This algorithm finds the optimum solution by moving the particles in the solution space 39 . PSO performs searching via a swarm of particles that updates from iteration to iteration. To find the optimum solution, each particle moves in the direction to its previously best position and the global best position in the swarm 40 .
Determination the physicochemical properties of AV3RN and SY33. Physicochemical properties of CO 2 , organic solvent (DMSO) and solid solutes (AV3RN and SY33) are needed parameters in modelling the supercritical processes by cubic equation of states. Required characteristics of CO 2 and DMSO including critical properties, boiling and melting point, acentric factor and molar volume are reported in literatures. However, these properties for complex molecules such as pigments and dyes are not available and their experimental determination is so difficult. In these situations, suitable group contribution methods can be applied for computation of required solids properties.
In this work, critical pressure (P c ), critical temperature (T c ) and normal boiling point (T b ) are estimated by Marrero and Gani method 41 . The acentric factor (ω) and solid molar volume (ν) are calculated by Constantinou-Gani 42 and Fedors 35 methods, respectively. Also, the heat of fusion at the triple point ( H tp ) and the triple point temperature (T tp ), which can be well estimated to the melting point 26 , of AV3RN and SY33 molecules are obtained from differential scanning calorimetry (DSC) analysis. The critical and physicochemical properties of AV3RN and SY33, DMSO, and CO 2 are presented in Table 1. Also, required physical characteristics of Eq. (14) for AV3RN and SY33 molecules are presented in Table 2.

Result and discussion
Analyzing the phase behavior through thermodynamic modeling is a suitable approach for assessing the optimal conditions of nanoparticles production via GAS process. In current work, the fluid phase behavior in binary (CO 2 -DMSO) and ternary (CO 2 -DMSO-AV3RN) and (CO 2 -DMSO-SY33) systems was investigated using PR-EoS.
Binary system (CO 2 -DMSO). Degree of CO 2 solubility in organic solvent, as a function of temperature and pressure, is a key parameter in fine particles precipitation by the GAS process. The vapor-liquid equilibrium (VLE) behavior of CO 2 -DMSO binary system as a function of temperature and pressure is shown in Fig. 1. For validation the accuracy of the considered thermodynamic model, calculated equilibrium mole fraction of dissolved sc-CO 2 in DMSO (x 1 ) is compared with the experimental data reported by Gonzalez et al. 43 and Lee et al. 17 , in Fig. 1a. As is evident, model results are well consistent with the reported experimental data.
As expected, CO 2 solubility in organic solvent (x 1 ) increases by pressure increment and temperature decrement (Fig. 1). Bakhshi et al. 44 have presented similar results with calculation the solubility of CO 2 in toluene and ethanol with SRK-vdW2 and PR-HV models, respectively. Also, the same trend was reported by Yao et al. 45 for CO 2 -acetone system and Ghoreishi et al. 30 for CO 2 -DMSO system, in which the VLE behavior of both systems is calculated by PR-EoS. www.nature.com/scientificreports/ For determination the optimum operating pressure, the relative molar volume change is plotted as a function of pressure at the different temperatures. The optimum operating pressure of the GAS process has to be selected in such a way that the relative molar volume change of the CO 2 -solvent mixture shows a minimum value. So, a minimum point of pressure (P min ) along with a sharp increment of volume expansion in this plot is an accepted criterion for this purpose. Based on reported results in de la Fuente et al. research 36 , the P min value indicates the minimum required pressure in GAS process for precipitation of about 95% of the solute in the form of fine particles. Thus, the operating pressure should be above the pressure of this minimum point (P min ). Variation of relative molar volume for CO 2 -DMSO binary system as a function of pressure at constant temperature (328 K) is presented in Fig. 2. P min value of 7.95 MPa is evident in this graph. As shown in Fig. 3, trend of relative molar volume change of this binary system versus the dissolved CO 2 mole fraction in liquid phase at 328 K is completely similar to Fig. 2. Molar volume contraction was observed by increasing dissolved CO 2 concentration in DMSO (mole fraction) up to a specified value (x 1min = 0.52). Afterwards, further CO 2 dissolution in DMSO change the molar volume status to molar volume expansion.
Observed trend can be explained based on definition of liquid molar volume of solution (ν) 46 : where ν 1 and ν 2 are partial molar volume of component 1 (CO 2 ) and 2 (DMSO) in a binary mixture of them. Regardless of whether ν 2 is bigger or smaller than ν 1 , generally, the molar volume of the solution (ν) is frequently Given that liquid molar volume of solution is a function of pressure, temperature and CO 2 mole fraction, dν dx 1 can be defined as follows:  dP dx 1 )) is a negative and small value, whereas ( ∂ν ∂x 1 ) P,T is negative at the first and becomes positive with increasing CO 2 dissolution in organic solvent. Accordingly, with initiation of CO 2 dissolution in DMSO, the liquid molar volume of solution reduces at the first (as ν 1 ≺ ν 2 ). After that, by increasing dissolved CO 2 concentration in DMSO, ν 1 becomes equal with ν 2 and the liquid molar volume of solution reaches to its lowest value (the minimum point of corresponding graph). Subsequently, further CO 2 dissolution leading ν 1 becomes larger than ν 2 and volume expansion be revealed. These trends are confirmed by other researchers, too. Chen et al. 31 applied Volume-Translated Peng-Robinson (VTPR) EoS for calculation the relative molar volume change of liquid phase in binary systems of CO 2 and different solvents including acetone, DMSO, ethanol, and ethyl acetate. Different P min values obtained at 308 K as 5.40, 7.65, 7.10, and 6.15 MPa K for acetone, DMSO, ethanol, and ethyl acetate, respectively. These values clearly show the effect of solvent type on solute precipitation in GAS process. Ghoreishi et al. 30 have also reported the P min values of 7, 7.74 and 8.5 MPa for CO 2 -DMSO binary system at different temperatures (308, 313 and 319 K). They calculated the relative molar volume variation vs. pressure using PR-EoS with linear combination of Vidal and Michelsen mixing rules (PR-LCVM). Reported results in Mamata Mukhopadhyay research 28 show that partial molar volume of DMSO and CO 2 are descending and ascending functions of dissolved CO 2 mole fraction, respectively.
Variations of relative molar volume against pressure for CO 2 -DMSO binary system at various temperatures (308, 318, 328 and 338 K) were shown in Fig. 4. Accordingly, the computed P min values were 7.27, 7.61, 7.95 and 8.29 MPa at 308, 318, 328 and 338 K, respectively. As is evident, required P min for solute precipitation is increased with temperature increment which is related to reducing the solubility of sc-CO 2 gas in DMSO due to temperature increment.
The binary interaction parameters of Eqs. 9 and 10 (k ij and l ij ) have been optimized by PSO algorithm for each temperature and optimally fitted parameters are shown in Table 3.

Ternary systems of (CO 2 -DMSO-AV3RN) and (CO 2 -DMSO-SY33).
For determination the optimum operational conditions of AV3RN (pigment) and SY33 (dye) nanoparticles precipitation via GAS process, the phase behavior of (CO 2 -DMSO-AV3RN) and (CO 2 -DMSO-SY33) ternary systems was investigated by PR-EoS with vdW2 mixing rule. For validation the model results in a ternary system, relative molar volume change  Table 3. Binary interaction parameters from PR-EoS with the vdW2 mixing rules for the binary system (CO 2 (1)-DMSO (2)). www.nature.com/scientificreports/ and variation of solute (naphthalene) mole fraction in CO 2 -toluene-naphthalene ternary system were determined by this model and compared with experimental data reported by Fuente Badilla 36 . As shown in Fig. 5a,b, good consistence between model results and experimental data confirms the model accuracy.

System CO 2 -DMSO
The binary interaction parameters of the (CO 2 -DMSO-AV3RN) and (CO 2 -DMSO-SY33) ternary systems, based on PR-EoS, are reported in Tables 4 and 5, respectively. The relative molar volume change of both of these ternary systems at 328 K was shown in Fig. 6. As can be seen, the lowest relative molar volume change of the liquid phase is appeared at 9.8 and 9.5 MPa for AV3RN and SY33, respectively. So, AV3RN and SY33 nanoparticles can optimally precipitated above the relevant P min values and at the same temperature, SY33 precipitates a bit sooner. The behavior of the relative molar volume change versus pressure at different temperatures (308, 318, 328 and 338 K) is shown in Fig. 7a and b for CO 2 -DMSO-AV3RN and CO 2 -DMSO-SY33 ternary systems, respectively. Also, the obtained results were validated with the experimental data of CO 2 -DMSO-phthalocyanine green (Pc-G) ternary system, reported by Sodeifian et al. 47 . As is evident, model results and the experimental data are in good agreement with each other, which approving the model precision. As previously described, increasing temperature leads to higher minimum pressure (P min ) value. The obtained P min values at 308, 318, 328 and 338 K were 7.80, 8.57,9.78 and 11 MPa and 8.0,8.63,9.50 and 10.77 MPa for AV3RN and SY33, respectively. Comparision the Figs. 7 and 3 indicates that the calculated P min values for both of ternary systems (Fig. 7) were higher than that of the CO 2 -DMSO binary system (Fig. 3) at corresponding temperatures. Difference between www.nature.com/scientificreports/ P min values of binary and ternary systems highlights the significance of selection the optimum operating conditions for each individual ternary system. At pressures higher than P min , a high supersaturation status of solute is expected for successful precipitation of fine solid particles. So, the operation pressure should be selected above the P min value, which was confirmed by others, too 31,48,49 .
Variation the concentration (mole fraction) of antisolvent (CO 2 ), solvent (DMSO) and both of solutes in the GAS process at 328 K are calculated by PR-EoS and the results are shown in Figs. 8 and 9 for AV3RN and SY33, respectively. As seen in these figures, trend of mole fractions variations in both of ternary systems are the same, where in increment CO 2 mole fraction and decrement DMSO mole fraction with increasing pressure were observed, simultaneously. Also for detailed interpretation, variation the mole fraction of dissolved AV3RN and SY33 in DMSO solution along with the change of relative molar volume of this liquid phase at 328 K is shown in Figs. 8b and 9b, respectively. Solutes solubility (mole fraction) in liquid phase gradually reduces by initiation CO 2 injection into DMSO solution. By further pressure enhancement, severe drop in solutes solubility is observed Table 4. Binary interaction parameters from PR-EoS with the vdW2 mixing rule for the ternary system (CO 2 (1)-DMSO (2)-AV3RN (3)).  Table 5. Binary interaction parameters from PR-EoS with the vdW2 mixing rule for the ternary system (CO 2 (1)-DMSO (2)-SY33 (3)). www.nature.com/scientificreports/ in the P min of each system (9.8 MPa for AV3RN and 9.5 MPa for SY33). Consequently, solutes mole fraction decreases to about zero, which confirms the precipitation of almost all of the solutes. This clearly approves that operational pressure of GAS process for each ternary system should be selected higher than its related P min value. It is worth noting that trend of variations the mole fractions of CO 2, DMSO and both of pigments in the liquid phase of both of ternary systems in other temperatures is the same of 328 K, which are not shown here. Experimental studies of Dixon and Johnston 50 on determination the mole fraction changes of CO 2 as antisolvent, toluene as solvent and naphthalene and phenanthrene as solutes have also the same results of our research. Thermodynamic modelling with regular solution theory and expanded liquid equation of state models confirmed their experimental results, too. Peters et al. 36 calculated the change of naphthalene and phenanthrene solutes of (CO 2 -toluene-naphthalene) and (CO 2 -toluene-phenanthrene) ternary systems with PR-EoS. They simultaneously presented the variation of relative molar volume change of the solution and solutes mole fraction vs. pressure in one chart. It is completely obvious that in pressures higher than reported P min , the solutes concentration approaches to zero. Variation the solubility of ampicillin as a function of pressure in a ternary system of CO 2 -DMSO-ampicillin was calculated via PR EoS, by Ghoreishi et al. 30 . A slight increase in ampicillin solubility (mole fraction) in the liquid phase in pressure increment to about 4.5 MPa, and after that, a sharp decrease of its solubility have observed in the pressure range of 4.5-7.3 MPa. The minimum ampicillin solubility was indicated www.nature.com/scientificreports/ at the pressure above 7.3 MPa, which is equal to calculated P min value of this ternary system. Change of CO 2 , toluene and naphthalene in liquid phase of corresponding ternary system was experimentally measured and thermodynamically predicted (PR-EoS and SRK-EoS) by Pahlavanzadeh et al. 26 . Also, they did the same work for CO 2 , ethanol and acridine system. In both of systems, sharp decrement of solutes (naphthalene and acridine) and solvents (toluene and ethanol) solubility were observed at calculated P min value of each system. Reduction of salicylic acid mole fraction in CO 2 -1-propanol-salicylic acid ternary system at 288 K was experimentally determined and thermodynamically predicted by Peters et al. 51 . They also reported the same trend for salicylic acid and benzoic acid mole fraction variation at 313 K in CO 2 -1-propanol-salicylic acid and CO 2 -acetone-benzoic acid ternary systems, respectively 52 .

Conclusion
In current study, the Peng-Robinson equation of state (PR EoS) was utilized for determination the phase equilibrium of the binary (CO 2 -DMSO) and both of ternary (CO 2 -DMSO-anthraquinone Violet 3RN (AV3RN)) and (CO 2 -DMSO-solvent yellow 33 (SY33)) systems. The critical properties of AV3RN and SY33 were calcuated using group contribution methods. The optimal operational condition (pressure and temperature) for precipitaion of AV3RN and SY33 nanoparticles via the GAS process were specified by the plot of relative liquid phase molar volume variations (computed with the definition presented by de la Fuente Badilla et al.) vs. pressure at the operational temperature. The minimum pressure of the (CO 2 -DMSO-AV3RN) and (CO 2 -DMSO-SY33) ternary systems were determined 7.80, 8.57,9.78,10.46 MPa and 8,8.63,9.5 and 10.77 MPa at 308, 318, 328 and 338 K, respectively. As is evident, temperature increment causes to increasing the minimum pressure value. Finally, variations the mole fractions of CO 2 , DMSO, AV3RN and SY33 of each ternary system were determined by the PR-EoS. By increasing dissolved CO 2 concentration in the DMSO solution, solutes solubility reduces to about zero, which confirms the precipitation of almost all of the solutes. Sharp decrement of solutes solubility at pressures higher than P min value of each ternary system approves that operational pressure of GAS process should be selected above this pressure.