Numerical modeling of mass transfer in the supercritical antisolvent process

Numerical modeling of mass transfer in the supercritical antisolvent process

Journal of Supercritical Fluids 16 (1999) 167–181 Numerical modeling of mass transfer in the supercritical antisolvent...

2MB Sizes 0 Downloads 17 Views

Journal of Supercritical Fluids 16 (1999) 167–181

Numerical modeling of mass transfer in the supercritical antisolvent process Jane O. Werling, Pablo G. Debenedetti * Department of Chemical Engineering, Princeton University, Princeton, NJ 08544, USA Received 19 February 1999; received in revised form 18 June 1999; accepted 18 June 1999

Abstract A mathematical model for mass transfer between a droplet of organic solvent and a compressed antisolvent is presented. The model, which allows for two-way mass transfer both into the droplet and into the bulk antisolvent, is applicable to the supercritical antisolvent (SAS) method of particle formation. In this process, solute particles are formed by precipitation from an organic solution upon contacting with a supercritical fluid antisolvent which is miscible with the organic solvent, but immiscible with the solute. The mass transfer behavior of the droplets is thought to be a key factor affecting particle morphology. In this work, conditions that are supercritical with respect to the pure antisolvent, but subcritical with respect to the solvent–antisolvent mixture, are considered, meaning that two phases are always present. Radial profiles of composition and density as a function of time, both inside and outside the droplet, are generated by numerical solution of the governing equations. The radius of the droplet as a function of time is also determined. Calculations with toluene as the organic solvent and carbon dioxide as the compressed antisolvent show that the initial interfacial flux is always into the droplet, causing droplet swelling. Trends in the extent of droplet swelling as a function of operating conditions correlate with trends in the initial interfacial flux, which shows a non-monotonic pressure and temperature dependence and is especially sensitive near the critical locus of the mixture. At the mixture critical point, droplet lifetime diverges. Away from the critical point, the lifetime decreases as the pressure is increased, but shows a non-monotonic temperature dependence. However, the time required for a droplet to reach saturation at its center does not correlate with droplet lifetime, so that droplets with short lifetimes are not necessarily the fastest to reach saturation. © 1999 Elsevier Science B.V. All rights reserved. Keywords: Counterdiffusion model; Droplet swelling; Particle formation; Supercritical antisolvent; Two-way mass transfer

1. Introduction There is currently a growing interest in using supercritical fluids in materials processing applications [1]. By taking advantage of the distinguishing * Corresponding author. Tel.: +1-609-258-5480; fax: +1-609-258-0211. E-mail address: [email protected] (P.G. Debenedetti)

thermophysical properties of supercritical fluids [2], such alternative routes can offer important advantages with respect to conventional processing methods. In addition to the environmentally benign character of carbon dioxide, the most common supercritical solvent, materials processing applications can benefit from the sensitivity of solvent power to small changes in pressure, and the favorable mass transfer characteristics (comparatively low viscosity and high diffusivity relative

0896-8446/99/$ – see front matter © 1999 Elsevier Science B.V. All rights reserved. PII: S0 8 9 6 -8 4 4 6 ( 9 9 ) 0 0 02 7 - 3


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

Nomenclature a parameter in Peng–Robinson equation of state [mass length5/(time2 mol2)] b parameter in Peng–Robinson equation of state ( length3/mol ) D diffusion coefficient from Fick’s Law ( length2/time) D modified diffusion coefficient containing thermodynamic correction ( length2/time) f fugacity [mass/( length time2)] k, l binary interaction parameters n number of moles N flux due to composition gradients [mol/( length2 time)] ˜ N flux at interface [mol/( length2 time)] P pressure [mass/( length time2)] P vapor pressure [mass/( length time2)] vap r radial coordinate ( length) R gas constant [mass length2/(time2 mol K )] R instantaneous droplet radius ( length) R initial droplet radius ( length) 0 t time T temperature v molar volume ( length3/mol ) x mole fraction in the liquid phase y mole fraction in the gas phase Z compressibility factor Greek letters w fugacity coefficient r density (mol/length3) j liquid phase radial coordinate in frame of reference moving with the interface j∞ gas phase radial coordinate in frame of reference moving with the interface Subscripts/superscripts A antisolvent (CO ) 2 B solvent (toluene) c critical value g gas phase i, j component indices l liquid phase m mixture r reduced value s solid phase 1 value at saturation 0 value at infinite dilution to liquids) of fluids in the relative vicinity of the critical point [1,2]. One such application involves precipitation of a solute from an organic solution using a compressed or supercritical fluid as an antisolvent. This method is known as the supercritical anti-

solvent (SAS ) or precipitation with a compressed antisolvent (PCA) process [3–6 ]. In this approach, the solute of interest is dissolved in the organic solvent, and the solution is cocurrently sprayed into an excess of the antisolvent. The organic solvent is miscible with the antisolvent, but the

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

solubility of the solute in the antisolvent is limited. Under appropriate contacting conditions, high mass transfer rates occur, leading to the formation of microparticulate solids. SAS offers several advantages as a particle formation process. Using supercritical carbon dioxide as an antisolvent (T =31°C ) allows pure c products to be generated at mild operating conditions. Carbon dioxide is a common choice for supercritical applications due to the fact that it is non-flammable, non-toxic, inexpensive, and environmentally acceptable. In addition, the ability to tune the density of a supercritical fluid with pressure allows for easy separation and recovery of solvent and antisolvent by a simple downstream depressurization step. These advantages make particle formation with supercritical fluids especially attractive for pharmaceutical applications, where the production of high-value-added products can make operating at supercritical conditions economically viable [1]. Two noteworthy features of SAS are the wide range of materials that can be so processed and the variety of particle morphologies that have been observed. Applications of SAS include the precipitation of explosive materials that are difficult to recrystallize by conventional methods [7], the formation of microparticulate protein powders that are biologically active upon reconstitution [4,8], precipitation of pharmaceutical compounds [9,10], pigments [11], and superconductor precursors [12,13], and the formation of polyamide fibers [5] and biocompatible polymers [14–16 ]. Coprecipitation of drug with a biopolymer has also been investigated [17,18]. These and other studies have yielded interesting and varied precipitate morphologies, including microspheres [3,4,6,14,19], fibers [3,5,6,19], microballoons (hollow microspheres) [20], and interconnected networks [6 ], all with varying degrees of porosity. Although some trends in particle size and morphology as a function of process conditions have been reported [3,4,6,9,11,14,19,21], conflicting observations abound in the literature. For example, increasing the precipitation pressure has been reported to decrease particle size [11,21], increase particle size [14], and have no effect on particle size [9]. In all cases, the theoretical basis for


experimental observations has been largely unexplored. Although currently deficient, the rational understanding of the relationship between particle size and morphology and process conditions should be valuable for the design and interpretation of experiments, and will be needed for the engineering design of industrial-scale SAS units. As a first important step towards a theoretical understanding of SAS, Dixon et al. [3] proposed a framework for rationalizing the relationship between phase behavior, mass transfer, and particle morphology. Their diagram for an antisolvent–solvent–polymer system is shown in Fig. 1. The diagram shows the binodal and spinodal curves for phase separation, while the arrows represent the mass transfer pathways available to the system. If we consider a droplet of organic solvent with dissolved polymer, each arrow can be thought of as showing the change in composition with time at a fixed point in the droplet or, alternatively, the change in composition from the center of the droplet to the droplet–fluid interface at a fixed time. In path A, the net mass transfer is primarily out of the droplet. The organic solvent evaporates with little antisolvent influx, thus avoiding phase separation. As a result, the polymer vitrifies as a dense particle. In path B, enough antisolvent penetrates the droplet to cause the system to cross the binodal curve on the polymer-rich side. As the amount of antisolvent in the drop increases, antisolvent-rich bubbles nucleate and grow in the

Fig. 1. Schematic diagram showing phase behavior and potential mass transfer pathways in SAS. Solid line=binodal curve; dashed line=spinodal curve. (Adapted from Dixon et al. [3].)


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

polymer-rich phase until the polymer vitrifies as a porous solid. The greater the extent of antisolvent penetration into the droplet, the greater the porosity of the resulting particles. If the nucleation kinetics of the system are sufficiently slow, as may be the case for systems with polymeric solutes, the system may avoid nucleation and cross the spinodal curve. An example is shown in path C, where the mass transfer path causes the system to pass close to the critical point. If the spinodal curve is crossed, phase separation occurs by spinodal decomposition [22]. This path is consistent with the observation of a bicontinuous morphology with fine structure, instead of discrete particles [6 ]. Finally, if the extent of antisolvent diffusion into the droplet is large compared to organic evaporation (path D), the droplets cross the binodal curve on the antisolvent-rich side. In this regime, the antisolvent swells the droplet until the polymer nucleates and grows in the antisolvent-rich phase. These pathways show that knowledge of phase behavior and mass transfer is essential for complete understanding and eventual control of system behavior and, ultimately, particle morphology. To date, theoretical modeling of SAS has consisted primarily of phase equilibria calculations. Dixon and Johnston [23] used the Peng–Robinson equation of state to predict the solubilities of phenanthrene and naphthalene in toluene expanded with CO over a range of pressures at 2 25°C. Their model predicts solubilities at low and high pressures, and accurately models the large drop in solubility at high pressures observed for solutes in toluene. Kikic et al. [24] formulated the equations for four-phase (solid–liquid–liquid– vapor) equilibria in the CO –toluene–naphthalene 2 system. The Peng–Robinson equation of state with van der Waals mixing rules was used to model the system’s behavior. Results of their model were compared to experimental data, and schematic phase diagrams illustrating the range of possible phase behavior were presented for temperatures between 315 and 325 K and pressures between 82 and 98 bar. Although the thermodynamics of SAS have received attention, the overall process is considerably more complex, involving jet hydrodynamics, mass transfer, and nucleation kinetics in addition

to phase equilibria. Kikic et al. [25] have recently developed a model that incorporates several of these processes. These workers constructed a pseudo-steady-state model that accounts for the thermodynamics, hydrodynamics, and mass transfer of the semi-continuous process. In addition to fugacity relationships, overall material balances and a force balance over the droplet were solved. In the mass transfer aspect of the model, an empirical correlation originally developed for spray towers was used to calculate mass transfer coefficients. Their calculations predict that droplets experience an initial increase in diameter followed by evaporation, and that larger drops require more time to evaporate. By defining precipitation to be the time at which the mole fraction of solute in the droplet reaches the equilibrium value, they found that higher gas-to-liquid flow rate ratios cause earlier precipitation and that the sensitivity of the process to the flow rate ratio is dependent on the solute. The important work of Kikic and coworkers does not, however, address the details of mass transfer within the droplet. This information is required in order to understand and predict the morphology of solute particles formed during SAS precipitation. The complexity of the SAS process suggests that a fundamental understanding may be gained by breaking down the whole into its constituent parts. A detailed analysis of the hydrodynamics of jet breakup in a supercritical fluid, the mass transfer between an individual droplet and its environment, the range of phase behavior that a solvent–solute–antisolvent system can attain, and the nucleation kinetics within the droplet will all be valuable steps in the direction of developing a comprehensive model and improved understanding of SAS. Among these constituent parts, we choose here to model the problem of mass transfer between an isolated droplet and an antisolvent continuum. The objective of this work is to gain a detailed understanding of this aspect of SAS through numerical solution of the governing equations for isothermal mass transfer between a drop of solvent and an excess of compressed antisolvent. We take the existence of the droplet as a starting point. Mass transfer occurs in both directions, from the droplet to the continuous phase and vice

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

versa. The change in droplet composition and density as a function of time and radial position is calculated, and the droplet radius as a function of time is also determined. In addition, the effect of temperature, pressure, and initial droplet radius on mass transfer is investigated. An area of research closely related to mass transfer in SAS is the evaporation of fuel droplets in spray combustion processes. In spray combustion, a cold fuel droplet is injected into a gaseous or supercritical medium at high temperatures and pressures. Mass, species, and energy balances are solved to track the change in droplet radius and the extent of heat and mass transfer between the droplet and its environment as a function of time. Hsieh et al. [26 ], Delplanque and Sirignano [27], Lage et al. [28], and Jia and Gogos [29] have presented detailed models for droplet vaporization near the mixture critical point. The present mass transfer model is related to droplet combustion studies, but with a very significant difference. In the systems which are typically of interest in combustion, the solubility of the gas in the liquid phase is relatively low. For example, Jia and Gogos [29] and Hsieh et al. [26 ] considered the evaporation of n-pentane droplets in a nitrogen environment. The mole fraction of nitrogen in n-pentane is limited at low pressures (P <1) and ranges from 0.175 to 0.375 at P =3, r r depending on the temperature (where P =reduced r pressure of pure n-pentane) [26 ]. Due to the low solubility of the gas in the liquid phase, it is often assumed that the density of the liquid phase is spatially uniform within the droplet and that convection in the liquid phase can be neglected [27,29]. Lage et al. [28] investigated the vaporization of n-heptane, methanol, and n-heptane/ methanol bicomponent droplets in air. They included the continuity equation in the liquid phase to account for the change in liquid phase density with time due to temperature changes, but assumed that the gas phase was insoluble in the liquid. In sharp contrast, the solubility of CO in organic 2 liquids is very high at conditions typical for SAS experiments. For example, at 310.1 K and 74.4 bar, the equilibrium mole fraction of CO in toluene is 2 0.913 [30]. Therefore, the droplets experience a very appreciable change in composition as they


undergo mass transfer with the CO environment. 2 Furthermore, the liquid phase density is extremely sensitive to the antisolvent content. Hence the density changes dramatically in the droplet’s radial direction, and accurate accounting of this behavior is crucial to understanding the evolution of the environment to which a solute is exposed on its path to precipitation.

2. Model formulation A fundamental analysis of mass transfer between the solvent droplet and the antisolvent continuum is made using time-dependent conservation equations and an equation of state. We assume spherical symmetry, equilibrium at the vapor– liquid interface, and a stagnant droplet. The only convective motion, therefore, is that induced by mass transfer. In this paper, we focus our attention on the subcritical regime. We use the term ‘subcritical’ to refer to temperatures and pressures at which the pure antisolvent is supercritical, but the solvent–antisolvent mixture is in the two-phase region. At these conditions, an interface always exists between the droplet and antisolvent phases. For the toluene (solvent)–CO (antisolvent) system 2 considered here, this means that the mixture is below its critical pressure at the given temperature, and above its critical temperature at the given pressure. This particular solvent–antisolvent combination was selected because it has been extensively used in particle formation experiments [3,6,20,31]. In addition, experimental evidence suggests that thermal effects are not important under practical SAS conditions [3,6,8,9]. Therefore, as a first step in our ongoing study of mass transfer in SAS, we consider the isothermal case here; work is in progress on extending the model by coupling it to an energy balance. The present work is therefore a base case against which more detailed, but necessarily more complex, calculations will be compared. Special care must be taken in describing the diffusion process. Mass transfer can occur in two directions: toluene can evaporate from the droplet, and carbon dioxide can diffuse into the droplet. Under typical SAS temperatures and pressures,


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

carbon dioxide is a dense fluid which is very soluble in toluene and will alter the density of the liquid droplet in a highly non-ideal way. Therefore, thermodynamics and mass transfer are inherently coupled. In addition, the rapid mass transfer and the large changes in density will induce corresponding changes in the droplet radius, which introduces an additional convective flux into the system. Furthermore, the fact that typical SAS operating conditions involve temperatures and pressures in the relative proximity of mixture critical points implies that care must also be taken in describing molecular diffusivity in these mixtures. At conditions far removed from the critical point, the diffusive flux can be assumed to be proportional to the concentration gradient as described by Fick’s Law: N =−rDVx +x N. A A A


In reality, however, the flux is proportional to the chemical potential gradient, and we write [32]: N =−rDVx +x N A A A


∂ ln wˆ l A D =D 1+ l l ∂ ln x T,P A


∂ ln wˆ g A D =D 1+ . g g ∂ ln y T,P A





This thermodynamic correction allows the modified Fick’s Law to reproduce the experimentally observed behavior at the critical point, where the diffusion coefficient goes to zero [33]. The equations describing the mole fraction profiles in each phase are obtained from the continuity equation and the modified Fick’s Law. All equations are in terms of molar quantities. The subscript A denotes carbon dioxide and B denotes toluene, and l and g represent the liquid and gas phases, respectively. Then, with R(t) denoting the droplet radius, we write, for the liquid phase [r
(r x )+V · (−r D Vx +x N )=0 l A l l A A l


and for the gas phase [r>R(t)]: ∂ ∂t

(r y )+V · (−r D Vy +y N )=0. g B g g B B g


Initially, the bulk droplet contains no CO and 2 the gas phase contains no toluene, but the interface is saturated; therefore, initial mole fraction profiles must be assumed in order to avoid convergence problems associated with the presence of discontinuities. In order to ensure that the simulation results were insensitive to the choice of initial profile, three types of initial profiles were tested: linear (ramp), hyperbolic tangent (sigmoidal ), and exponential. For ease of comparison, the profiles were chosen such that all three contained the same amount of minority component in the respective majority phase. Provided that the initial loading was kept below a critical value (5.4 mol% CO in 2 the droplet and 6.7×10−5 mol% toluene in the continuous phase), the results were insensitive (to within 1%) to the choice of initial profile. Note that this is a very conservative criterion, since in principle the initial loading should be exactly zero. In order to solve the species equations, it is necessary to solve for the convective flux due to density variations. Therefore, the continuity equation is used in each phase: ∂r l +V · N =0 l ∂t


∂r g +V · N =0. g ∂t


A molar balance over the droplet gives the expression for the change in droplet radius with time: ∂




˜ 4pr2r dr =−4pR2N l

(8) 0 ˜ is the interfacial flux. The expression for where N ˜ is derived from the physical requirement that N the total flux in each phase be equal at the interface. This is achieved by equating the modified Fick’s Law expressions in each phase at the interface: ∂t

˜ ] =[−r D Vy +y N ˜] . [−r D Vx +x N l l A A R g g A A R



J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

˜ is comEq. (8) can be rewritten to show that N posed of a contribution from the convective flux at the interface and a contribution from the rate of change of the droplet radius: ˜ =N N




. −r1 l dt


˜ can be defined using the gas phase Alternatively, N characteristics: ˜ =N N



dR . −r1 g dt


In order to describe the composition dependence of the gas and liquid densities, an equation of state is needed. In this work, the Peng–Robinson equation of state [34] was chosen because of its ability to predict phase equilibria, its relatively accurate density predictions at high pressures, and its ease of application [2]. The Peng–Robinson equation can be written as: RT


m (11) v(v+b )+b (v−b ) m m m m where the component a and b parameters are functions of the critical properties of the individual species and the mixture parameters a and b are m m obtained from the following mixing rules [32]: P=


a =∑ ∑ x x a m i j ij i j a =(a a )0.5(1−k ) ij ii jj ij b =∑ ∑ x x b m i j ij i j (b +b ) jj (1−l ). b = ii ij ij 2

of the mixture. Cussler [38] and Reid et al. [32] recommend an empirical relation proposed by Vignes for concentrated liquid mixtures [39]: D =(D0 )xA (D0 )1−xA (14) l g l Eq. (14) is thus an empirical mixing rule, whereas Eqs. (3a,b) ensure thermodynamic consistency [32]. Since the amount of toluene entering the gas phase is small, the gas phase may be considered dilute. In order to determine the boundary conditions for the species equations, it is necessary to solve for the composition of the saturated system at the temperature and pressure of interest. This is achieved by equating the fugacity of each component in the liquid and gas phases and solving for the liquid and gas phase mole fractions. fˆ l =fˆ g (15a) A A ˆf l =fˆ g (15b) B B x wˆ l =y wˆ g (16a) A A A A (1−x )wˆ l =(1−y )wˆ g . (16b) A B A B Eqs. (16a,b) are the ones used in actual calculations. The fugacity coefficient wˆ is a function of i temperature, pressure and composition. It can be related to the volumetric properties of the mixture by the following expression [40]:

(12a) (12b)

RT ln wˆ =− i

(13a) (13b)

The critical properties and the binary interaction parameters for CO and toluene used in this work 2 are those reported by Kikic et al. [24]. The liquid phase diffusion coefficient at infinite dilution was calculated using the Tyn–Calus method [35]. The gas phase diffusion coefficient at infinite dilution was calculated using the method of Fuller et al. [36 ] with the Takahashi correction for high pressures [37]. Since the droplet contains a high concentration of CO at saturation, an 2 equation is needed to approximate the diffusivity





i T,V,nj



dV−RT ln Z (17)

where Z is the compressibility factor of the mixture. The Peng–Robinson equation is used to evaluate this expression.

3. Numerics Eqs. (4)–(9) and (11) make up the mass transfer model. These equations were solved numerically to obtain the mole fraction and density profiles both inside and outside the droplet and the change in droplet radius with time. The calculations were performed at temperatures and pressures typical for SAS operation: 310–334 K and 75–110 bar.


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

The partial differential equations were reduced to ordinary differential equations by replacing the radial dependence with second-order finite difference approximations. The FORTRAN subroutine DDASSL (Double-precision DifferentialAlgebraic System SoLver) was then used to solve the system of equations. DDASSL uses a finite difference method in which the step size and the order of the finite difference approximation are chosen by the code based on the behavior of the solution. A detailed discussion of DDASSL can be found in the book by Brenan et al. [41]. For ease of solution, the equations were nondimensionalized and the mole fractions in each phase were scaled with their saturation values. Because boundary conditions are applied at the droplet interface, which itself is a function of time, a coordinate transformation to a frame of reference moving with the interface was applied using the following relations. In the liquid phase [r
r R(t)



In the gas phase [r>R(t)]: j∞=ln

r R(t)



With this transformation, the interface is always at j=1 and j∞=0. The finite difference meshes in both the liquid and gas phases consisted of 100 radial nodes. A logarithmic mesh was used in the gas phase in order to cluster the grid points near the droplet interface where gradients are steep, while conserving grid points far away from the droplet for computing efficiency. A mass balance was performed over the entire system after each time step to ensure that mass was conserved. The boundary conditions at infinite distance from the drop were applied at r=20R , where R is the initial radius 0 0 of the droplet. Test runs with different values of mesh size and infinity setting were performed in order to choose the optimal values for efficient and accurate computing.

4. Results and discussion In order to predict the conditions under which evaporation or droplet swelling dominate, it is ˜ , the most convenient to examine the variable N net molar flux at the interface. Solving Eq. (9) for ˜ gives: N ˜= N



r1 D1 ∂y g g B x1 +y1 −1 ∂r r=R A B r1 D1 ∂x l l A . + R x1 +y1 −1 ∂r r= A B



The two coefficients of the derivative terms in Eq. (19) depend on the properties of the liquid and gas phases at saturation. Furthermore, the mole fraction gradient is always negative in the gas phase and positive in the liquid phase. For most organic solvent–CO mixtures, the solubility 2 of the CO in the organic liquid (x1 ) is much 2 A greater than the solubility of the organic in CO 2 (y1 ). For example, for toluene–CO at 318 K and B 2 82.5 bar, x1 =0.89 and y1 =0.01. Therefore, A B y1 %x1 in the denominator of the coefficients and B A we conclude that (x1 +y1 −1) is always negative. A B This causes the first term in Eq. (19) to be positive, favoring evaporation, whereas the second term will be negative, favoring swelling. ˜ (0) as Fig. 2 shows the initial interfacial flux N a function of pressure at several different temperatures. The flux has been scaled with [2500r1 D1 /R ], where r1 is the saturated molar l l 0 l density, D1 is the saturated liquid phase diffusivity, l and R is the initial droplet radius at 318 K and 0 82.5 bar. For the toluene–CO system over the 2 ˜ range of conditions considered, N (0) is always negative. This is primarily due to the difference in solubility, which causes the mole fraction gradient in the liquid phase to be much greater than the mole fraction gradient in the gas phase. Therefore, the net mass transfer is into the droplet, and ˜ (0) swelling results. The pressure dependence of N is non-monotonic: each isotherm passes through a sharp minimum, indicating a maximum initial flux into the droplet. The sensitivity to pressure is greatest near the mixture’s critical locus, and the model correctly predicts that the flux becomes zero

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181


˜ is a useful The net molar flux at the interface N measure of system behavior. It is also important to know how the size of the droplet is changing with time. Rewriting Eq. (8) and solving for dR/dt gives: dR dt

Fig. 2. Effect of temperature and pressure on the initial interfacial molar flux for toluene droplets exposed to an excess of supercritical carbon dioxide. Negative values indicate net flux into droplet. The molar flux has been scaled with [2500r1 D1 /R ], where r1 is the saturated molar density, D1 is l l 0 l l the liquid phase diffusivity, and R is the initial droplet radius. 0

(and hence the droplet lifetime diverges) when the two phases become indistinguishable. The minimum in each isotherm is a result of competing effects: the increased density of CO at higher 2 pressures results in higher solubilities of CO in 2 the liquid phase, causing the flux to become increasingly negative. Near the critical pressure, however, the requirement that the diffusion coefficient approaches zero causes the flux also to approach zero. It is important to note that a calculated value ˜ at time=0 depends on the initial mole fraction of N gradients in each phase. As discussed in the model formulation, initial mole fraction profiles must be assumed in order to avoid convergence problems associated with discontinuities. Although the ˜ (0) depend on the initial profiles, actual values of N ˜ (0) do not, for all the initial loadings the trends in N ˜ (0) are tested. Therefore, the actual values of N ˜ (0) is not important; instead, the usefulness of N in comparing droplets under different conditions and elucidating mass transfer trends. In generating ˜ (0) profiles in this work, the same initial the N mole fraction profiles were used in all cases.


1 R2r1 l


R 0


˜ ∂r N l dr− . ∂t r1 l


It can be seen from Eq. (20) that dR/dt depends on the rate of change in the droplet density, as ˜ . Therefore, if the net flux is into the well as on N ˜ <0) and the molar density of the solvent droplet (N decreases upon addition of antisolvent, then dR/dt will be positive and the droplet will grow. However, it is possible for the molar density of the liquid phase to increase with the addition of antisolvent, which would tend to decrease the droplet radius. Kordikowski et al. [42] investigated volume expansions of several polar solvents with ethane, ethene, and carbon dioxide in order to identify potential solvent–antisolvent systems. They found that addition of the antisolvents ethene and ethane causes a steep decrease in liquid mass density with increasing pressure; however, addition of CO results in a maximum in liquid mass density 2 with respect to pressure. According to Eq. (20), an increase in molar density tends to decrease the droplet radius. In the case of toluene, however, which shows a mass density [30] and a molar density maximum with addition of CO , the value 2 ˜ (0) outweighs the other effects. That of N is to say, the amount of mass entering the droplet is large enough that the increase in density is overwhelmed and the droplet radius increases nevertheless. Thus, we conclude that the swelling/ evaporation behavior can change dramatically from one system to another, being as it is the result of several competing factors. Figs. 3 and 4 show the time evolution of several toluene droplets as a function of temperature and pressure. Droplet lifetimes are characterized by the initial rapid diffusion of CO into the droplet, due 2 to the large mole fraction gradient in the liquid phase. As diffusion occurs, the gradient decreases ˜ balance. until the two terms in the expression for N After this point, the net diffusion is out of the droplet and the droplet shrinks. Eventually, the


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

Fig. 3. Evolution of toluene droplet radius at 82.5 bar and different temperatures. Initial radius of all droplets is 50 mm. The mixture’s critical temperature at this pressure is 313.49 K.

Fig. 4. Evolution of toluene droplet radius at 318 K and different pressures. Initial radius of all droplets is 50 mm. The mixture’s critical pressure at this temperature is 87.405 bar.

droplet becomes fully saturated and mass transfer is completely out of the droplet until the droplet disappears. Figs. 3 and 4 show that the droplet swelling can be appreciable. At 82.5 bar and 314 K the droplet swells from 50 to 192 mm, almost a four-fold increase in radius. In addition, the extent of swelling for each droplet correlates with the trends in

˜ (0) shown in Fig. 2: at a constant pressure of N 82.5 bar, the extent of droplet swelling decreases with increasing temperature, except in the immediate vicinity of the critical point, T=313.49 K (not shown), where the swelling decreases sharply and the lifetime diverges. At a constant temperature of 318 K, swelling increases with pressure, except in the immediate vicinity of the critical point, P= 87.405 bar, where, as shown, the swelling decreases sharply and the lifetime diverges. The effect of pressure is again shown in Fig. 5, which compares the evolution of droplets at 75 and 87 bar, at 318 K. Initially, both droplets experience a rapid change in composition, but the swelling is significantly greater at 87 bar than at 75 bar, due to the greater solubility of CO in 2 toluene at the higher pressure. However, after 1.10 s, the low pressure droplet is still larger than its original size, but the high pressure droplet has shrunk considerably. This example illustrates the complexity of the mass transfer process. From Eq. (19) it can be seen that several factors contribute to the evolution of a droplet during SAS, and interact in a complicated way. The compositions, densities, and diffusivities of the two phases, and how these variables change with time, all combine to affect mass transfer and the change in droplet radius as a function of time. Therefore, it is not in general straightforward to isolate a single cause underlying a given trend in droplet behavior. The temperatures and pressures considered in this work range from just above the critical point of pure CO to the critical points of the 2 CO –toluene mixture. The critical point of the 2 mixture was selected as a stopping point for the simulation because above the critical point, the vapor–liquid interface ceases to exist, and the equations derived in this work no longer apply. Subcritical conditions are not uncommon for SAS operation: for example, experiments using dimethyl sulfoxide as the organic solvent have been conducted over a wide range of temperatures and pressures, both in the subcritical and supercritical regimes [4,5,8,16 ]. The estimation of time scales is another important use of the mass transfer model. Process scale-up, for example, requires a knowledge of the

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181


Fig. 5. Evolution of droplet composition and radius at 75 and 87 bar at 318 K. Initial droplet radius is 50 mm. On the color continuum, 0=pure toluene and 100=toluene saturated with CO at the respective pressure. 2

characteristic times for mass transfer and an understanding of how this time changes with operating conditions. The times required for toluene droplets to reach their maximum size, for saturation at the droplet center to occur, and for complete evaporation (total lifetime) were examined as a function of temperature at 82.5 bar, and as a function of pressure at 318 K. The results are shown in Figs. 6 and 7. Droplet lifetimes show two regimes. In the nearcritical region, all measures of lifetime diverge due to the vanishing interfacial flux. Away from this region, the total lifetime decreases with increasing pressure and is non-monotonic in temperature. The time to saturation is non-monotonic in both temperature and pressure. The time required for maximum swelling is relatively insensitive to both temperature and pressure. All calculations presented in this work assume

an initial droplet radius of 50 mm. The equations in the mass transfer model are scaled by the initial droplet radius, and since we assume a spherical droplet, the effect of droplet radius is determined exclusively by the characteristic time of the system, R2 /D1. Therefore, doubling the initial droplet 0 l radius causes a four-fold increase in the amount of time required for these mass transfer events. The dependence of the characteristic time on the liquid phase diffusion coefficient implies that the time to the mass transfer events will also depend on the solvent used. To see how the mass transfer model is useful for understanding basic trends in SAS, the results can be analyzed in conjunction with a ternary phase diagram. The phase diagram for the model ternary system CO , toluene, naphthalene at 325 K 2 and 85 bar is shown in Fig. 8. The coexistence curves delineating the liquid, liquid–solid, solid–


J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

Fig. 6. Time required for a toluene droplet to reach maximum radius, time to saturation at the center of the droplet, and total lifetime as a function of temperature at 82.5 bar. Curves shown are for droplets at 314, 316, 318, 320, 322, 326, and 332 K. Initial radius of droplets is 50 mm.

Fig. 8. Calculated phase behavior of carbon dioxide, toluene, and naphthalene at 325 K and 85 bar, including tie lines. L= liquid phase, S=solid phase, V=vapor phase. The Peng– Robinson equation was used in the calculations. Details are given in the text.

solid phase is assumed to be pure naphthalene, for which the fugacity can be written as: f s =P exp naph vap

Fig. 7. Time required for a toluene droplet to reach maximum radius, saturation at the center, and total lifetime as a function of pressure at 318 K. Curves are for droplets at 76, 78, 80, 82, 84, and 85.6 bar. Initial radius of droplets is 50 mm.

liquid–vapor, solid–vapor and liquid–vapor regions of the phase diagram were constructed as described in the model formulation section. The



v (P−P ) s vap . RT


The phase equilibrium conditions were solved using the Peng–Robinson equation of state. The component critical constants and naphthalene vapor pressure and molar volume are those reported by Dixon and Johnston [23], and the interaction parameters used were those reported by Kikic et al. [24]. In cases where the coexistence curves created overlapping regions such that one phase envelope was completely contained within another one, a point was selected in the overlapping region and the Gibbs free energy was calculated assuming separation along each set of tie lines. The coexisting system with the lower Gibbs free energy is stable, while the other phases are metastable. Since the mass transfer model presented in this work considers the diffusion of solvent and antisolvent only, discussion of trends with respect to the phase diagram in Fig. 8 is limited to the CO –toluene axis. When the droplet is first intro2 duced to the CO environment, the center of the 2 droplet is pure toluene. As mass transfer occurs,

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181

the composition of the center of the droplet moves along the CO –toluene axis towards the toluene2 rich side of the liquid-vapor binodal curve. Simultaneously, the vapor phase composition approaches the CO -rich side of the binodal from 2 the CO vertex. When the center of the droplet 2 reaches the toluene-rich side of the binodal, the droplet is fully saturated.

5. Conclusion A mathematical model for the mass transfer between an isolated solvent droplet and an antisolvent continuum is presented. The model accounts for two-way mass transfer of solvent and antisolvent, and includes the composition dependence of the diffusion coefficients, non-idealities in the liquid phase density due to addition of antisolvent, and corresponding changes in the droplet radius. Due primarily to the high solubility of CO in toluene, the initial interfacial flux is always 2 into the droplet for this system, resulting in an increase in droplet radius. Furthermore, the initial interfacial flux can be used to relate the operating temperature and pressure to the extent of droplet swelling. In addition, investigation of the effect of process conditions on the time scales for mass transfer shows that lower pressures and larger initial radii result in longer droplet lifetimes, and that lifetimes increase sharply at near-critical conditions. However, time for a droplet center to reach saturation does not correlate with droplet lifetime, meaning that droplets with short lifetimes are not necessarily the fastest to reach saturation. The mass transfer model developed in this work is also useful for describing the condition of the droplet at saturation. Whether the center of the droplet reaches saturation as a result of solvent evaporation or swelling will significantly affect the concentration of the solute in the droplet. In the case of a low molecular weight solute such as naphthalene, it is unlikely that concentration effects play a role in particle morphology, since the crystallization kinetics of low molecular weight solutes are orders of magnitude faster than the mass transfer process. However, for more complex, high molecular weight solutes such as polymers


and proteins, the nucleation kinetics may be appreciably slower. In such cases the solute concentration dictates the region of the phase diagram in which precipitation occurs, therefore affecting the morphology of the resulting particles. This work represents the first attempt at a fundamental understanding of the detailed mass transfer between a solvent droplet and an antisolvent continuum in the subcritical regime. This work is complementary to the global model developed by Kikic et al. [25] and to previous thermodynamic modeling of SAS-related systems [23,24]. In future work, we will extend the mass transfer model to the supercritical regime, and investigate thermal effects by including an energy balance in the model equations. In addition, calculating ternary phase diagrams for a complex solute and incorporating the nucleation kinetics of the solute into the mass transfer model will allow for more accurate prediction of mass transfer pathways and precipitation conditions, leading to a more comprehensive understanding of the SAS process. Such a model should provide valuable insight into the environment to which a solute is exposed upon precipitation, which can have significant consequences for particle morphology.

Acknowledgements The authors would like to thank Y.G. Kevrekidis of Princeton University and G. Gogos of the University of Nebraska at Lincoln for helpful discussions on the numerics of the mass transfer model. P.G.D. gratefully acknowledges the support of the Air Force Office of Scientific Research through grant F49620-96-1-0169. J.O.W. gratefully acknowledges the National Science Foundation for a Graduate Fellowship.

References [1] C.A. Eckert, B.L. Knutson, P.G. Debenedetti, Supercritical fluids as solvents for chemical and materials processing, Nature 383 (1996) 313. [2] M.A. McHugh, V.J. Krukonis, Supercritical Fluid Extrac-





[6 ]










[16 ]



J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181 tion: Principles and Practice, Butterworth-Heinemann, Stoneham, MA, 1994. D.J. Dixon, K.P. Johnston, R.A. Bodmeier, Polymeric materials formed by precipitation with a compressed fluid antisolvent, AIChE J. 39 (1993) 127. S.-D. Yeo, G.-B. Lim, P.G. Debenedetti, H. Bernstein, Formation of microparticulate protein powders using a supercritical fluid antisolvent, Biotech. Bioeng. 41 (1993) 341. S.-D. Yeo, P.G. Debenedetti, M. Radosz, H.-W. Schmidt, Supercritical antisolvent process for substituted paralinked aromatic polyamides: phase equilibrium and morphology study, Macromolecules 26 (1993) 6207. D.J. Dixon, K.P. Johnston, Formation of microporous polymer fibers and oriented fibrils by precipitation with a compressed fluid antisolvent, J. Appl. Polym. Sci. 50 (1993) 1929. P.M. Gallagher, M.P. Coffey, V.J. Krukonis, W.W. Hillstrom, Gas anti-solvent recrystallization of RDX: formation of ultra-fine particles of a difficult-to-comminute explosive, J. Supercrit. Fluids 5 (1992) 130. M.A. Winters, B.L. Knutson, P.G. Debenedetti, H.G. Sparks, T.M. Przybycien, C.L. Stevenson, S.J. Prestrelski, Precipitation of proteins in supercritical carbon dioxide, J. Pharm. Sci. 85 (1996) 586. W.J. Schmitt, M.C. Salada, G.G. Shook, S.M. Speaker, Finely-divided powders by carrier solution injection into a near or supercritical fluid, AIChE J. 41 (1995) 2476. M. Hanna, P. York, B.Yu. Shekunov, Control of the polymorphic forms of a drug substance by solution enhanced dispersion by supercritical fluids (SEDS ), in: M. Perrut, P. Subra ( Eds.), Proc. Fifth Meeting on Supercritical Fluids Vol. 1 (1998) 325. Y. Gao, T.K. Mulenda, Y.-F. Shi, W.-K. Yuan, Fine particles preparation of Red Lake C pigment by supercritical fluid, J. Supercrit. Fluids 13 (1998) 369. E. Reverchon, G. Della Porta, C. Celano, S. Pace, A. Di Trolio, Supercritical antisolvent precipitation: a new technique for preparing submicronic yttrium powders to improve YBCO superconductors, J. Mater. Res. 13 (1998) 284. E. Reverchon, G. Della Porta, A. Di Trolio, S. Pace, Supercritical antisolvent precipitation of nanoparticles of superconductor precursors, Ind. Eng. Chem. Res. 37 (1998) 952. T.W. Randolph, A.D. Randolph, M. Mebes, S. Yeung, Sub-micrometer sized biodegradable particles of poly(Llactic acid) via the gas antisolvent spray precipitation process, Biotechnol. Progr. 9 (1993) 429. J. Bleich, B.W. Mu¨ller, W. Waßmus, Aerosol solvent extraction system — a new microparticle production technique, Int. J. Pharm. 97 (1993) 111. L. Benedetti, A. Bertucco, P. Pallado, Production of micronic particles of biocompatible polymer using supercritical carbon dioxide, Biotech. Bioeng. 53 (1997) 232. Y.-H. Chou, D.L. Tomasko, Gas crystallization of polymer–pharmaceutical composite particles, in: K. Arai (Ed.), Proc. Fourth Int. Symp. on Supercritical Fluids (1997) 55. L. Sze Tu, F. Dehghani, A.K. Dillow, N.R. Foster, Appli-




[22] [23]



[26 ]








cations of dense gases in pharmaceutical processing, in: M. Perrut, P. Subra ( Eds.), Proc. Fifth Meeting on Supercritical Fluids Vol. 1 (1998) 263. R. Bodmeier, H. Wang, D.J. Dixon, S. Mawson, K.P. Johnston, Polymeric microspheres prepared by spraying into compressed carbon dioxide, Pharm. Res. 12 (1995) 1211. D.J. Dixon, G. Luna-Barcenas, K.P. Johnston, Microcellular microspheres and microballoons by precipitation with a vapour–liquid compressed fluid antisolvent, Polymer 35 (1994) 3998. J. Robertson, M.B. King, J.P.K. Seville, D.R. Merrifield, P.C. Buxton, Recrystallisation of organic compounds using near critical carbon dioxide, in: K. Arai (Ed.), Proc. Fourth Int. Symp. on Supercritical Fluids (1997) 47. P.G. Debenedetti, Metastable Liquids: Concepts and Principles, Princeton University Press, Princeton, NJ, 1996. D.J. Dixon, K.P. Johnston, Molecular thermodynamics of solubilities in gas antisolvent crystallization, AIChE J. 37 (1991) 1441. I. Kikic, M. Lora, A. Bertucco, A thermodynamic analysis of three-phase equilibria in binary and ternary systems for applications in rapid expansion of a supercritical solution (RESS), particles from gas-saturated solutions (PGSS) and supercritical antisolvent (SAS ), Ind. Eng. Chem. Res. 36 (1997) 5507. I. Kikic, A. Bertucco, M. Lora, Thermodynamics and mass transfer for the simulation of recrystallization processes with a supercritical antisolvent, in: E. Reverchon ( Ed.), Proc. Fourth Italian Conf. on Supercritical Fluids and their Applications (1997) 299. K.C. Hsieh, J.S. Shuen, V. Yang, Droplet vaporization in high pressure environments I: near-critical conditions, Combust. Sci. Technol. 76 (1991) 111. J.-P. Delplanque, W.A. Sirignano, Numerical study of the transient vaporization of an oxygen droplet at sub- and supercritical conditions, Int. J. Heat Mass Transfer 36 (1993) 303. P.L.C. Lage, C.M. Hackenberg, R.H. Rangel, Non-ideal vaporization of dilating binary droplets with variable properties, Int. J. Heat Mass Transfer 36 (1993) 3731. H. Jia, G. Gogos, High pressure droplet vaporization; effects of liquid-phase gas solubility, Int. J. Heat Mass Transfer 36 (1993) 4419. C.J. Chang, C.-Y. Chen, H.-C. Lin, Solubilities of carbon dioxide and nitrous oxide in cyclohexanone, toluene, and N,N-dimethylformamide at elevated pressures, J. Chem. Eng. Data 40 (1995) 850. E.M. Berends, O.S.L. Bruinsma, J. de Graauw, G.M. van Rosmalen, Crystallization of phenanthrene from toluene with carbon dioxide by the GAS process, AIChE J. 42 (1996) 431. R.C. Reid, J.M. Prausnitz, B.E. Poling, The Properties of Gases and Liquids, 4th edn., McGraw-Hill, New York, 1987. V.G. Levich, Physicochemical Hydrodynamics, 2nd edn., Prentice-Hall, Englewood Cliffs, NJ, 1962.

J.O. Werling, P.G. Debenedetti / Journal of Supercritical Fluids 16 (1999) 167–181 [34] D.Y. Peng, D.B. Robinson, A new two-constant equation of state, Ind. Eng. Chem. Fundam. 15 (1976) 59. [35] M.T. Tyn, W.F. Calus, Diffusion coefficients in dilute binary liquid mixtures, J. Chem. Eng. Data 20 (1975) 106. [36 ] E.N. Fuller, K. Ensley, J.C. Giddings, Diffusion of halogenated hydrocarbons in helium: the effect of structure on collision cross sections, J. Phys. Chem. 75 (1969) 3679. [37] S. Takahashi, Preparation of a generalized chart for the diffusion coefficients of gases at high pressures, J. Chem. Eng. Jpn. 7 (1974) 417. [38] E.L. Cussler, Diffusion: Mass Transfer in Fluid Systems, Cambridge University Press, Cambridge, 1984. [39] A. Vignes, Diffusion in binary solutions: variation of


diffusion coefficient with composition, Ind. Eng. Chem. Fundam. 5 (1966) 189. [40] M. Modell, R.C. Reid, Thermodynamics and Its Applications, 2nd edn., Prentice-Hall, Englewood Cliffs, NJ, 1983. [41] K.E. Brenan, S.L. Campbell, L.R. Petzold, Numerical Solution of Initial-Value Problems in Differential-Algebraic Equations, Society for Industrial and Applied Mathematics, Philadelphia, PA, 1996. [42] A. Kordikowski, A.P. Schenk, R.M. Van Nielen, C.J. Peters, Volume expansions and vapor–liquid equilibria of binary mixtures for a variety of polar solvents and certain near-critical solvents, J. Supercrit. Fluids 8 (1995) 205.