Surface Energy Evolution in Pharmaceutical Powder Micronization

Article pubs.acs.org/IECR

Surface Energy Evolution in Pharmaceutical Powder Micronization Using Compressed Gas Antisolvent (Re-)Precipitation Daniel E. Rosner*,† and Manuel Arias-Zugasti‡ Department of Chemical & Environmental Engineering, Sol Reaction Engineering Group, Yale University, New Haven, Connecticut 06520-8286, United States ABSTRACT: We illustrate the importance of environment-dependent surface energy changes in predicting the micronization of active pharmaceutical ingredients (APIs) in gas antisolvent precipitation (GASP) processes. This size-reduction scheme exploits compressed CO2(g) as antisolvent (AS) at near-ambient temperatures. Ordinary API-loaded solvents (often sprays) are contacted with dense CO2, and during CO2 uptake in an evolving expanding liquid API + solvent + CO2 solution droplet, particle nucleation (N) sets in, continuing along with growth (G) and, ultimately, coagulation. A rational method [due to Nielsen and Sohnel (J. Cryst. Growth 1971, 11, 233) and Mersmann (J. Cryst. Growth 1990, 102, 841)) is used to estimate the changing embryonic solid/ternary solution interfacial energy, γ. We demonstrate the dramatic yield and crystal size distribution (CSD) consequences of surface energy evolution (SEE) by carrying out N/G calculations for the surrogate organic API: phenanthrene dissolved in representative well-mixed micrometer-sized toluene droplets (sprayed into 298 K CO2 for p < 60 bar). To solve the population balance partial differential equation, we exploit the method of characteristics. Our results demonstrate that assuming constant surface energy, sometimes reasonable for API precipitation via the rapid expansion of supercritical-CO2 solvent (i.e.: relatively dilute rapid expansion of a supercritical solution conditions), fails for GASP-process modeling. When the crystal growth kinetics are sufficiently rapid, SEE also modifies performance via the Gibbs−Kelvin reduction of small particle growth rates. Rational yet tractable methods to incorporate both systematic effects in future design/optimization/parameter estimation calculations are suggested. called RESS applications,1,15 in which an API “surrogate”, phenanthrene, particles are precipitated from dilute solutions in dense CO2 which is rapidly expanded, this effective interfacial energy has been estimated to be only about 19 mJ/m2,13,14 and, of particular concern here, treated as a constant. We demonstrate below that while this latter assumption may be acceptable in some (RESS-type) applications, it is inappropriate to realistically describe API-particle inception and growth in the more flexible GASP processwhich starts with nondilute API solutions in conventional (effective) liquid solvents. The need to include surface energy evolution (SEE) for nucleation was apparently first recognized by Dodds et al.16 in their modeling of a semibatch GASP process for the anti-inflammatory corticosteroid BDP dissolved in liquid acetone. However, these authors did not focus on this particular effect and also treated the small particle growth process as unmodified by SEE. By exploiting our present well-mixed-droplet (WMD) GASP model, along with the N/G laws subsequently described, we are able to unambiguously demonstrate here the expected role of SEE in dramatically modifying product yield and powder size distributions. However, more detailed experiments and analysis will probably be necessary to establish the adequacy of this present approach and/or point the way to further systematic improvements. In this regard, one possible alternative route to

1. INTRODUCTION AND OBJECTIVES An emerging and versatile process for the size reduction of pharmaceutically active powders under mild processing conditions exploits compressed CO2 gas as an antisolvent (AS), into which a convenient active pharmaceutical ingredient (API)-containing liquid solvent is injected .1−7 This results in API reprecipitation as an ultrafine powder. Often abbreviated with the acronym GASP (for gas(-induced) antisolvent precipitation) or PCA (for precipitation via compressed antisolvent), this process presents many thermodynamic, kinetic (nucleation, growth), and convective-diffusion transport challenges to the theorist seeking a realistic yet tractable process model.8−11 The present contribution12 focuses on relaxing some of the troubling assumptions of relevant earlier work and is intended to ultimately improve the ability of process engineers to predict performance and select optimum GASP operating conditionsi.e., to produce, at higher yields, active pharmaceutical ingredients of desired mean particle size with narrow population spread. Previous work on modeling near-critical CO2 antisolventinduced precipitation of “micronized” pharmaceuticals9,10 has made extensive use of “classical nucleation theory” (CNT) to estimate particle birth rates.13,14 This convenient (yet much maligned!) theory, in which the creation of interfacial energy appears as an activation barrier that must be overcome, eliminates the need to perform a separate population balance on small embryonic precursor particles. However, as is now well-known, it introduces the effective surface free energy, γ, as a decisive “material property” appearing in the nucleation rate expression cubed in the dimensionless argument of the exponential function (see section 2.2). For example, in so© 2014 American Chemical Society

Received: Revised: Accepted: Published: 4489

August 7, 2013 December 26, 2013 February 7, 2014 February 17, 2014 dx.doi.org/10.1021/ie4025853 | Ind. Eng. Chem. Res. 2014, 53, 4489−4498

Industrial & Engineering Chemistry Research

Article

anticipate/incorporate SEE effects17 is briefly considered in section 4.3. Historically, especially in the field of industrial crystallization, effective constant surface energies have been “back-calculated” based on direct or indirect measurements of sparingly soluble particle homogeneous nucleation rates. But, in the ASP applications of interest here, the API solute is initially often present in an effective organic solvent, implying the likelihood (if not desirability) of initial nondiluteness prior to extensive CO2 addition. For example, near 300 K the presently considered API surrogate phenanthrene (C14H10) exhibits a saturation mole fraction in toluene (or benzene) of about 20%, corresponding to a saturation mass fraction near 1/3. With the introduction of the antisolvent (here CO2) we demonstrate in later text that the effective surface energy of the solid/fluid interface will necessarily “evolve”, increasing more than enough to dramatically reduce the predicted particle birth rate: Ḃ ‴, in the time-dependent supersaturated “expanding liquid” solution. This reduction will diminish with the additional uptake of antisolvent (CO2), as shown in section 2.2. These expectations, along with significant, if less dramatic, effects on recently nucleated particle growth rates, motivated the following analysis to provide a rational yet tractable way to anticipate the N + G consequences of such effective surface energy evolution.

manner). Then, exploiting our convenient idealized mathematical model for an isobaric GASP process,20 we can compare crystal size distribution (CSD) predictions using either constant surface energy or eq 1. In this way it is possible to unambiguously and economically demonstrate the role of SEE in modifying product powder distributions and precipitate yield. While instructive, more detailed experiments and analysis will probably be necessary to either establish the adequacy of our present approach, and/or point the way to further systematic improvements (see, also, section 4.3). 2.2. Surface Energy Evolution for CO2-Saturated Solvent + API Solutions. As a useful preliminary exercise, consider a saturated CO2 + T + Ph solution at 1 bar, 298 K. We now ask, how would Γ (computed using the method of section 2.1) change if the VSE level of Ph decreased as this liquid became saturated with the antisolvent CO2 at higher pressures (here up to 64 bar) while disallowing solvent evaporation or Ph precipitation? Using our thermodynamic model (section 3), the results for this reference state Γ are shown in Figure 1also

2. METHODS: QUASI-STEADY ESTIMATION OF RELEVANT TIME-DEPENDENT SURFACE ENERGY IN GASP PROCESSES Our premise is that while there is some intrinsic API (surrogate phenanthrene) crystal surface energy (relative to a vacuum or its own low sublimation pressure at its triple point), the effective value, γeff (when “wetted” by a liquidlike solution containing the solvent (T), the solute (Ph), and the antisolvent (CO2)), will be reduced significantly, in accord with the equations developed below. For simplicity, in what follows we neglect the presence of contaminants (e.g., solvent) in the solid phase. 2.1. Estimation of Evolving Effective Surface Energy. A useful precedent and starting point is provided by the laboratory-scale semibatch GASP-process study of Dodds et al.,16 who recognized the existence of SEE in modeling GASP nucleation rates. They incorporated a Nielsen−Sohnel− Mersmann (NSM) correlation of the following form:18,19 ⎛ ns ⎞1/3 ⎟ Γ ∝ ln⎜ API L ⎝ nAPI ⎠

Figure 1. Dimensionless surface energy vs p for an initially vapor− liquid-equilibrium (VLE) CO2 + T + Ph mixture at p = 1 bar and T = 298.15 K, with xPh/xsat Ph = 1 (lower dotted curve) and 0.1 (upper dotted curve). Solid curves show the corresponding Γ found for this APIdilute mixture after CO2 addition up to the saturation level, but without Ph/T precipitation or vaporization. The gray curve corresponds to vapor−liquid−solid-equilibrium (VLSE) conditions.

including the effects of initial Ph undersaturation. While these Γ-values are found to be less than 0.6 below ca. 40 bar, a rapid rise is observed above ca. 50 bar, with values well above unity when the pressures are such that the expanded liquid has a CO2 mole fraction above ca. 0.8, with a corresponding reduction in the saturation number density nLAPI (see eq 1). 2.3. Effects of Surface Energy Evolution on both Nucleation and Growth Kinetics. Because classical nucleation theory (CNT) predicts an API-particle “birth” rate, Ḃ ‴ (particles per unit time and volume), proportional to

(1)

(where Γ ≡ γvm /(kBT) is the dimensionless surface energy, nAPI is the number density of API in the droplet [superscript L] or crystal [superscript s], vm is the molecular volume of API in the solid state, and kBT is the Boltzmann constant times absolute temperature) to evaluate local nucleation rates via classical nucleation theory (CNT) while leaving their assumed crystal growth rate law unaltered in their numerical analysis of the relevant univariate population balance. Because of its simplicity and our present goals, we also adopt this NSM correlation (eq 1). We also suggest (section 2.4) a method to incorporate SEE in a more appropriate law for crystal growth kinetics in transcritical environments.20 For the present we select the “constant” (proportional to kBT/vm2/3) in eq 1 to agree with previous estimates of γ for the Ph/CO2 system (e.g., 19 mJ/m2 at T = 343 K and p = 260 bar,21 in the absence of the solvent tolueneobtaining the value 1.06 in this 2/3

3⎞ ⎛ ̇ ∝ Γ1/2 exp⎜ −16π Γ ⎟ B‴ ⎝ 3 ln 2 S ⎠

(2)

where S is the prevailing (activity-based) supersaturation,13−15,21,22 the preliminary results of section 2.2 immediately lead us to expect dramatic modifications in the expected output particle number densities and mean sizes. This is confirmed when allowance is made (section 2.4) for the simultaneous evolution of supersaturation using appropriate solute mass balance equationsand we provide comparable 4490

dx.doi.org/10.1021/ie4025853 | Ind. Eng. Chem. Res. 2014, 53, 4489−4498


Article

results (shown dashed) when Γ is (unrealistically) forced to remain constant (at its initial value at the outset of CO2 influx). The reader may have noticed that many published descriptions of classical nucleation theory refer to γ as the “surface energy of the new (precipitated) phase”as if this were a property of the new phase aloneindependent of the host environment. This view is quite misleading, as convincingly illustrated in this study. 2.3.1. Remarks on “Secondary” Nucleation and Heterogeneous Nucleation. Several earlier treatments of GASP/PCA processes have explicitly included what is called “secondary” nucleation and its parametrization. This class of mechanisms, associated with nuclei formation from existing particles, is deliberately omitted here for simplicity. Dependent upon the solvent/antisolvent contacting scheme, their relative importance is indeed often “secondary”. However, when such mechanisms cannot be neglected (see, e.g., Jarmer et al.23), SEE remains likely to play an important role. [Note: For example, in the “mixed suspension mixed product removal” (MSMPR) study of poly(lactic acid) (PLLA) precipitation from methylene chloride in liquid CO2 (see Jarmer et al.23), the rate of secondary nucleation was empirically inferred to be nearly linear in the product of the single particle growth rate, Ġ , and the suspension mass density (their MT). But in section 2.4 we show how SEE can influence monomer-based growth rates via a particle-size-dependent Seff. Keeping in mind that MT is proportional to our particle volume fraction, ϕp (section 3.1), then our Figure 3 suggests how SEE would modify this contribution to secondary nucleation. On this basis alone we are led to expect that when AS/S contacting conditions are such that secondary nucleation is important, SEE (section 2.1) is also likely to be consequential.] This should also be true if foreign nanoparticle nuclei were present in the initial solution sprayed into the CO2-containing chamber. Here we simply assume that any such preexisting nuclei cannot compete with the embryos formed by homogeneous nucleation of the nondilute solute. 2.4. Effects of Surface Energy Evolution on Particle Growth Kinetics. Coupled component balance equations that enable the actual AS-initiated supersaturation “history” to be calculated must include the growth of previously nucleated particles by API (supersaturated solute) “monomer scavenging” as well as current homogeneous nucleation. Thus, another systematic effect of surface energy evolution, surprisingly omitted in all earlier GASP modeling, is that associated with the reduced growth rate of recently born particles. In our present GASP-process model (described more completely in ref 20) we incorporate this additional effect by introducing a Gibbs−Kelvin−Ostwald factor, FGKO, presuming that the volumetric growth rate of a particle of volume v will be approximately proportional to v2/3(1 − (v*/v)1/3), where v* is the CNT critical volume at the prevailing instantaneous supersaturation, S; i.e., v ⎛ ⎞3 * = 32π ⎜ Γ ⎟ vm 3 ⎝ ln S ⎠

where the grouping k G (T,...) ln n−1 S, containing the phenomenological coefficients kG(T,...) and n, plays the role of an overall incorporation probability (ε), Ż API ″ is the prevailing 1/3 API molecular impingement flux, and FGKO ≃ S(v*/v) . [Note: Here, and in eq 6 in subsequent text, we make use of a rational interpolation formula between the ideal gas limit (for ϕ/ϕL ≪ 1) and liquid environments using a normalized molecular volume fraction: ϕ/ϕL, where ϕL ≈ 0.62. For details, see ref 20.] Note that, for v ≫ v*, this law reduces to the frequently observed behavior: Ġ ∼ lnn S in liquid solvents.24 It is also noteworthy that when n = 1, the dimensionless rate constant kG(T,...) and the corresponding overall incorporation probability, ε, become identical. Further details of these rate laws and our mathematical/numerical methods, beyond the scope of this present work, are contained in ref 20. Because of the CNT self-consistency requirement v* > vm, an interesting corollary of eq 2 is that CNT will no longer be selfconsistent if the local supersaturation S exceeds exp(3.224Γ). This necessary condition was satisfied (see Figure 2) for the simulations discussed in sections 4.1 and 4.2.

Figure 2. Time history of dimensionless critical particle volume, v*, for an initially VLE CO2 + T + Ph mixture at p = 1 bar and T = 298.15 K, after instantaneous compression to p = 56 bar and subsequent CO2 addition. The solid curve shows results based on prevailing dimensionless surface energy Γ given by eq 1. The dotted curve shows corresponding results when Γ is assumed to be constant. Present results are based on a constant AS-uptake rate assumption, which becomes inaccurate toward the end of our numerical calculations, when droplet AS saturation level is reached. As a consequence the plateau followed by a sharp variation observed toward the end of the present calculations (constant Γ case), and incipient leveling-off of v* (variable-Γ case), are not expected to be seen in a real experiment.

3. GASP-PROCESS MODEL TO INVESTIGATE “SEE” CONSEQUENCES 3.1. Well-Mixed Droplet Model. As a tractable “platform” for examining the GASP-process consequences of surface energy evolution, we have adopted a well-mixed (single) droplet (WMD) model in which we focus on crystal size distributions resulting from homogeneous nucleation and interface-controlled particle growth during the stage in which isobaric CO2 uptake at constant (impingement-controlled) rate occurs with negligible loss of solvent from the swelling droplet. As indicated previously, we make use of classical nucleation theory for the birth rate, and, for the growth kinetics in these transcritical environments, we introduce a collision theory

(3)

More explicitly, we have implemented a droplet growth rate law (2Ġ = ddp/dt, where dp is the instantaneous particle diameter and t is time) of the factorable form ″̇ k G(T , ...)(ln n − 1 S) Ġ = vm,APIZAPI

FGKO ⎛ S ⎞ ln⎜ ⎟ S ⎝ FGKO ⎠

(4) 4491



Article

Table 1. Pure Component Properties (AS, S, and API Surrogate)a component

M (kg/kmol)

Tc (K)

pc (bar)

Vc (m3/kmol)

ω

Zc

CO2 (AS) toluene (S) phenanthrene (API surrogate)

44.01 92.14 178.23

304.2 591.8 869.3

73.8 41.1 29.0

0.0940 0.3155 0.554

0.225 0.257 0.495

0.274 0.263 0.222

a The value used here for the sublimation enthalpy of phenanthrene at the triple point (TAPI,tp = 372.4 K) was ΔHsublim = 90.34 MJ/kmol, and the corresponding saturation pressure was ptp = 24.04 Pa.

perspective (see ref 25 and preceding section 2.4). For the relevant solution thermodynamics already implicit in the results shown in Figure 1, we adopted a Peng−Robinson cubic EOS incorporating binary interaction coefficients (for the a and b parameters) specific to the ternary Ph + T + CO2 system26−28 (see Table 1). Thus, the a and b parameters were computed according to the usual mixing rules, a = ∑i∑j xixjaij and b = ∑i∑j xixjbij, with aij = (1 − kij)(aiaj)1/2 ,

bij = (1 − lij)

We consider here the individual component balance equations appropriate to our present GASP mathematical model. Recall that we focus on the stage of isobaric, isothermal antisolvent AS uptake with negligible loss of solvent S or solute APItreating the hypothetical expanding liquid droplet as well-mixed at each instant. For the solvent S (toluene) we impose the condition of negligible loss during CO2 uptake; i.e., ⎡ ⎤ d ⎢ xS(1 − ϕp) ⎥ = =0 dt ⎢⎣ V (x ; p , T ) ⎥⎦

bi + bj 2

(5)

where the values of the binary interaction parameters kij and lij are as follows: for CO2−T, kij = 0.090 and lij = 0.0; for CO2− Ph, kij = 0.078 and lij = −0.030; for T−Ph, kij = −0.004 and lij = −0.059. Molar balance equations25 are written for each of the solution components, and a population balance equation (a linear PDE containing the above-mentioned rate laws and host fluid properties29) is written to track the univariate particle number density distribution function n(v,t)assuming negligible “growth” by Brownian coagulation. 20 This latter approximation, defended a posteriori, enables our use of the method of characteristics (MOC) to numerically generate complete particle size distributions (PSD), free of any imposed PSD-“shape” constraints (see section 3.3 and refs 29 and 30). While several volume-based moments of n(v,t) are of practical and theoretical interest,29 we focus subsequent discussion on a dimensionless Sauter mean diameter, SMD, and a dimensionless relative spread parameter, σ/SMD (where σ is the particle diameter standard deviation of the precipitated phase), when evaluated at times which are specified multiples of the characteristic uptake time: tAS ≡

sat d0xAS = sat(p ,T ) ″̇ 6V LsatZAS,w =0

(7)

where V(x;p,T) is the molar volume of the prevailing ternary fluid mixture of composition x ≡ {xS, xAPI, xAS}T, = is the present value of the (expanded) droplet total volume (see subsequent text), and ϕp is the present value of the API particle volume f raction (initially zero). For the antisolvent (AS) we consider the molar inflow rate (Ṅ AS) to be constant at the initial value: ̇ = NAS,0 ̇ = πd0 2ZAS ″̇ NAS

(8)

In that case the AS-component balance ODE can be written ⎡ ⎤ d ⎢ xAS(1 − ϕp) ⎥ ̇ = = NAS dt ⎢⎣ V (x ; p , T ) ⎥⎦

(9)

Perhaps most interesting is the balance equation for the component API which is ultimately distributed between the initial solute (monomer in solution) and the precipitated particles. Thus, we can impose the condition of constancy of NAPI, where ⎡ xAPI(1 − ϕ ) ϕp ⎤ p ⎥= = NAPI,0 NAPI = ⎢ + NAvm ⎥⎦ ⎢⎣ V (x ; p , T )

(6)

where =0 = (π/6)d03 is the initial droplet volume, superscript sat refers to saturation conditions, subscript 0 refers to initial conditions, xAS is the AS mole fraction, and VL is the molar volume of the liquid phase. In text that follows (sections 4.1 and 4.2), using this WMD physical model, we illustrate the predicted dramatic consequences of surface energy evolution on the crystal size distribution function, as well as the corresponding dimensionless API powder mean diameteri.e., SMD/d ref and dimensionless powder population spread (i.e., σ/SMD) for initially Ph-saturated toluene droplets of 1 μm diameter sprayed into CO2 at 56 bar. While the effects of SEE on particle nucleation rates via CNT are decisive (eq 2 and section 4.3), we also incorporate (sections 2.3 and 4.2) the significant effects of SEE on growth rate reductions for recently “born” particlesbut these growth rate effects are expected to be secondary under these particular conditions for a substance such as phenanthrene near 298 K.20

(10)

noting that ϕp,0 = 0, where NA is Avogadro’s number (0.6023 × 1027 molecules/(kg-mol)). By assumption, in the expanded solvent all mole fractions will sum to unity at any instant, and the total “droplet volume” = must satisfy the overall molar balance: ⎡ ϕp ⎞ ⎤ d ⎢⎛ 1 − ϕp ̇ ⎟⎟=⎥ = NAS ⎜⎜ + dt ⎢⎣⎝ V (x ; p , T ) NAvm ⎠ ⎥⎦

(11)

where V(x;p,T) is fixed by the prevailing composition and the Peng−Robinson EOS. Of course, the initial value of = is simply the (specified) value of (π/6)d03. 3.2. Particle Population Balance Equation. Before the onset of appreciable Brownian coagulation, the univariate29 population balance equation (PBE) for the condensate API population (particles per unit droplet volume n(v,t)) is31 4492



Article

∂n ∂ ̇ δ(v − v ) − n d ln[(1 − ϕ )=] + (vṅ ) = B‴ p * ∂t ∂v dt

n+ = n− + (12)

(13)

and δ(v − v*) is the Dirac δ function centered at the critical size v*. The first term on the right-hand side (RHS) of eq 12 represents the time variation of n(v,t) as a consequence of homogeneous nucleation (see eq 4), and the second term represents the time variation due to changes in the available droplet volume (1 − ϕp)= . We note that (when (1 − ϕp)= is increasing) this term has the character of a first-order homogeneous sink of particles, with a size-independent effective rate constant (d/dt) ln[(1 − ϕp)= ]. While the former PBE (eq 12) is clearly coupled to the (time-dependent) environment conditions through the relevant quantities S, Γ, v*, Ḃ ‴, and Ġ , the environment evolution equations (given by eqs 7−11), together with the Peng− Robinson EOS, are coupled to the particle number density distribution function n(v,t) only through the particle volume fraction ϕp, given in terms of n by ϕp =

∫0

∞

vn(v ,t ) dv

(14)

The evolution equation for ϕp is found by integrating (v × eq 12). Thus we find dϕp dt

=

∫0

∞

̇ − ϕ d ln[(1 − ϕ )=] vṅ (v ,t ) dv + v B‴ p p * dt (15)

which, together with eq 12, eqs7−11, and the Peng−Robinson EOS, closes the system of evolution equations that define our mathematical model of this physical system. 3.3. Numerical Methods Employed: Method of the Characteristics Formulation. The first-order linear PDE for n(v,t) allows for a formal solution via the MOC, using the particle volume as the coordinate along the characteristic paths. The characteristic path v = v(t) is considered, given by the solution of the first-order ODE dv/dt = v̇, v(t0) = v0 (where the starting time t0 and the initial value v0 have arbitrary values). As long as the characteristic path v(t) does not cross the timedependent value of v*, the particle nucleation term is identically zero, allowing for the following formal quadrature of eq 12: n[v(t ),t ] = n[v0 ,t0]

∫t

t 0

dv ̇ dv

⎡ exp⎢ − (1 − ϕp(t ))=(t ) ⎢⎣

(1 − ϕp(t0))=(t0)

v(t )

⎤ dt ⎥ ⎥⎦

(17)

Hence, whenever a characteristic line crosses the nucleation size, the corresponding value of n[v(t),t] undergoes a discontinuity, with intensity given by Ḃ ‴(t*)/|v̇ − v̇*|t*. Although the former quadrature (eq 16) and boundary condition (eq 17) provide a formal solution of the PBE, it is not straightforward to extract numerical information from this solution without further simplifications. This is because the numerical evaluation of eqs 16 and 17 needs to be carried out along with the numerical integration of the system of ODEs that determine the evolution of the environment and the integro-ODE that determines the value of ϕp. The solution of this last equation (eq 15), or alternatively the immediate result eq 14, is difficult to evaluate in a numerically accurate way using the information provided by eqs 16 and 17, partly because of the discontinuous nature of n[v(t),t] at each crossing time t*, which forces one to consider relatively large numbers of characteristic paths to get acceptable results. In the present work this difficulty has been solved by means of a perturbation scheme, in which particle growth acts as a small correction to the time evolution dictated by droplet swelling and particle nucleation. This perturbation scheme is based on the large disparity observed, under present conditions, between the characteristic growth and nucleation times, both of them determined by an order of magnitude analysis of the former system of evolution equations (eqs 7−11, eq 12). In the results shown as follows, we have used the following dimensionless variables. The dimensionless particle volume is given by w ≡ v/vref, where vref = 6.44vm is defined as the critical nucleation size corresponding to reference values of surface tension and supersaturation. In the present work these reference values have been defined as the corresponding Γ and S found at current pressure (56 bar) and temperature (298 K), and reference droplet composition. The latter has been defined as that corresponding to instantaneous compression up to 56 bar of a 1 bar-saturated VLE mixture, after addition of CO2 up to the saturation level xsat CO2 (56 bar, 298 K), assuming the (hypothetical) absence of precipitation despite saturation with the antisolvent CO2, and negligible solvent loss. This convenient reference state has also been used to define the corresponding reference values of the particle nucleation and growth rates. Finally, a dimensionless particle number distribution function f has been defined by f ≡ nvref/nref, where the reference particle number density is nref ≡ ϕp,∞/vref, with ϕp,∞ being the limiting API precipitate volume fraction, found by assuming complete precipitation of API above xsat Ph after the AS-uptake stage

where v̇ is the single particle growth rate law in terms of particle volume, v ̇ = a(v)Ġ = (36π )1/3 v 2/3Ġ

̇ (t ) B‴ * |v ̇ − v ̇ |t **

ϕp, ∞ =

sat sat x Ph,0(1 − xCO ) − x Ph 2

x Ph,0 + (V Lsat /VPh,s)

(18)

To facilitate comparison between the constant- and variable-Γ results, the same dimensionless units (based on variable Γ) have been used in both cases.

(16)

On the other hand, if at a certain value of time (say t*) this characteristic line crosses the time-dependent critical volume (v(t*) = v*(t*)), then the nucleation term sets in. In this case integration of eq 12 over the infinitesimal time interval (t* − dt,t* + dt) leads to the following boundary (i.e., “jump”) condition25 that connects the values of n[v(t), t] before (n[v*, t*−] = n−) and after (n[v*, t*+] = n+) the crossing

4. RESULTS: EFFECTS OF SURFACE ENERGY EVOLUTION USING THE WELL-MIXED-DROPLET MODEL Illustrative results will be displayed as a function of t/tAS, for the typical cases: n = 2, kG(298K) = 10−4, tAS = 0.0056 ms, and d0 = 4493



Article

Figure 3. Precipitated API volume fraction, ϕp (left), and fraction of API precipitated, ξ (right), vs dimensionless time. Dotted curves correspond to constant surface energy; solid curves correspond to Γ given by eq 1. (T = 298 K, p = 56 bar, and d0 = 1 μm.)

Figure 4. Predicted dimensionless particle Sauter mean diameter (SMD, left) and particle-diameter-based standard deviation (σ) divided by the population SMD (right) vs t. Dotted curves correspond to constant surface energy; solid curves correspond to Γ given by eq 1. (T, p, and d0, as in Figure 3.)

Figure 5. Dimensionless precipitated particle number density distribution function (NDDF) wf vs w ≡ v/vref (log scale) for 10 equispaced values of time between 0 and 0.2tAS. The left plot corresponds to constant surface energy, and the right plot corresponds to Γ given by eq 1. (T, p, and d0, as in Figure 3.) The peaks that can be seen in the NDDF for values of time close to the final time are motivated by the behavior of v* in that limit (see Figure 2) and are beyond the domain of applicability of our present simplified constant AS-uptake assumption, and as a consequence are not expected to be accurate.

regard to our choice of (only) a 1 μm diameter toluene solvent droplet, used for illustration purposes, there are two reasons for this. The first one is the self-consistency of our WMD mathematical model. Second, there is an incentive to move

1 μm: however, we remark that our present WMD theoretical model also leads to the identification of a potentially useful rate-based dimensionless scaling (similitude or correlation) parameter for this class of GASP processes.20 [Note: With 4494



Article

ppm (Figure 5), we have also formally run this (upper limit) case to illustrate the consequences of SEE via the Gibbs− Kelvin−Ostwald effect on embryonic particle growth rates. [Note: In practice, size-dependent fluid phase diffusion would insert itself to limit the Ph(s) growth rate in such cases.] 4.3. Improved Quasi-Steady Methods To Anticipate/ Incorporate Surface Energy Evolution. Hopefully, to successfully model future GASP processes it will not be necessary to introduce a much more detailed (e.g., selective, nonquasi-steady adsorption-based) computational method. However, it may very well prove necessary to go beyond the presently employed Nielsen−Sohnel−Mersmann (NSM) correlation (eq 1), which has the remarkable property that the only relevant property of the multicomponent “solvent” is its ability to alter the saturation number density of the solute. An interesting possible alternative will be the Girifalco− Good relationship,17 which relates the free energy γαβ of a unit area of the αβ interface to the individual liquid adjacent phase energies γα and γβ, allowing for phases α and β to have unequal molar volumes. In the present case γα would presumably be the surface free energy of phenanthrene solid and γβ the surface tension of the prevailing ternary host fluid, i.e., here T + Ph + CO2perhaps estimated via the constituent pure saturated fluid surface tensions (at the prevailing temperature T and molar composition xi) via a semiempirical “parachor-based” mixing rule. While under investigation, this route to estimating SEE is understandably quite sensitive to the choice of individual surface energies−which must also be fully self-consistent. With these sensitivities in mind, the simpler, yet remarkably successful, NSM correlation was considered to be adequate for our present purposes. Of course, in retrospect, “solvent effects” on crystal morphology (via relative [crystal face] growth rates in supersaturated liquid solvents) are quite familiar,34 so the importance of evolving surface energy we find for GASPprocess modeling is not without precedents. Ultimately, surface energy anisotropy could be introduced, opening the door to useful predictions of precipitated crystal shape.35

toward smaller injector droplet sizes to minimize crystal growth effects, and, perhaps, also delay the importance of Brownian coagulation.32] 4.1. Crystal Size Distribution Consequences of Surface Energy Evolution. The main consequence of SEE, due to AS-induced reduction in the local API concentration (via eq 1) is a significant reduction (over 50%) in the nucleationkinetics-limited amount of precipitated material. This result is illustrated in Figure 3, where the precipitated API volume fraction is shown (left), together with the relative amount of precipitated API (right), as compared to the initial inventory in the droplet. Figure 3 clearly shows that the surface energy increase resulting from the time evolution of the droplet composition produces a remarkable reduction in the total amount of precipitated (as a fraction of the maximum possible extent of precipitation) material (solid curves), as compared to the results found when Γ is assumed to remain equal to its initial value (dotted lines). Because of the dependence on Γ of the critical nucleation size (see Figure 2), surface energy evolution also affects the shape of the CSD. This is shown in Figure 4 (left), where it can be seen that SMD corresponding to constant Γ is about 20% higher than the SMD computed retaining SEE effects. Figure 4 (right) shows that the constant Γ assumption leads to narrower CSDs, which could have also been anticipated by inspection of the time dependence of v* (Figure 2). The CSDs of the precipitated API are shown in Figure 5 (see Figure 5 (left) for constant Γ and Figure 5 (right) for SEE). Because under present conditions particle growth has been found to play a minor role, the shape of n(v,t) responds primarily to the time history of v*, which depends strongly on Γ. As a consequence the results for n(v,t) based on the constant Γ assumption show a “plateau” at small values of v (v/vref between 6 and 7), which corresponds to the corresponding (near-)plateau of v* during a considerable interval of time. This plateau is absent in the results found when taking into account environment-dependent SEE. On the other hand, also because under present conditions particle growth is weak, n(v,t) shows a strong peak at the values of v where v* reaches a relative extremum (either minimum or maximum). As a consequence we may find in both cases a strong peak reached at values close to the minimum v* as a function of time (see Figure 5 and Figure 2). For times close to the final time we may also find a strong peak in n at intermediate values of v, which corresponds to the maximum reached by v* as a function of time when the droplet mixture becomes transcriticalat which point the numerical integration is stopped. In this regard, the very last portion of our numerically computed time-history results (for times close to the final time) serves the purpose of illustrating the correspondence between the v* time history and precipitate NDDF. However, it is important to remark that these results are based on our constant AS-uptake rate assumption, which is expected to become inaccurate as the droplet AS level approaches the saturation level. According to our numerical results, this condition is met toward the end of our present numerical integrations, and as a consequence the last portion of our present numerical results will not be observed in a real experiment. 4.2. Effects of Surface Energy Evolution on Particle Growth. Because, in principle, the incorporation probability, ε, for phenanthrene crystal growth could be as high as 1/2 (see also ref 33), rather than our present 298 K expectation of only ca. 1 ppm, and our previous generous assignment of kG = 100

5. CONCLUSION AND RECOMMENDATIONS 5.1. “Particle Design” Implications of Surface Energy Evolution for GASP Modeling. While each of our present property estimation methods (sections 2 and 3) may require refinement in the light of new GASP-performance data, it seems clear from the well-mixed-droplet idealized examples selected above that useful predictions of API yield, mean particle size, and API powder population spread will have to be made by allowing for an environment-dependent effective surface energy. This API solid/local fluid “surface energy evolution” (SEE) dramatically alters the instantaneous particle nucleation rates and also reduces the growth rates of newly born API particles in the prevailing supersaturated solutions. Apparently, with the notable exception (for nucleation) of ref 16, neither of these systematic effects has been accounted for in earlier GASP modeling studies, calling into question the use of “constant-γ” models to infer meaningful physicochemical parameters (such as γ, kG, and n) from laboratory scale GASP experiments, and probably also compromising their ability to predict the performance of scaled-up GASP systems.36 Experiments needed to more directly test our present isobaric, small solvent droplet predictions are not readily available, but the lessons learned here should facilitate presently needed extensions (section 5.2). 4495



Article

5.2. Necessary Extensions and Work in Progress. Unless internal circulation and/or liquid-phase turbulence “come to the rescue”, the well-mixed-droplet-based GASPprocess model exploited here is probably not realistic for the (≫10 μm) solvent droplet sizes that were actually used in the earliest laboratory-scale experimentsi.e., the characteristic intradroplet antisolvent-diffusion time: d02/(4DAS−mix) is not much smaller than the AS-uptake time, tAS, above (see, also, Chavez et al.37). However, this convenient WMD theoretical platform has been shown here to provide access to rate-limited precipitate yields, unusual precipitate particle size distributions (including SMD and spread information), and sensitivity to important modeling approximations (e.g., surface energy evolution, Ġ law, ...). In our more recent studies (ref 20), we have exercised this mathematical model to include predicted sensitivity to initial solvent droplet size, d0 (and DSD), and also propose scaling (similitude) parameter(s). We view these as useful steps toward a much-needed, more comprehensive ratebased GASP-process modela “holy grail” which has now proven to be elusive for well over 1 decade (see refs 2, 7, 36, and 37.). Such a model should not only account for intradroplet concentration nonuniformities and surface energy evolution but also ultimately include surface energy anisotropyopening the possibility of also predicting precipitate particle morphology.

■

Mi n n n(v,t) NA ni Ni Ṅ i N‴ p p R S t tAS T Tc v V Vc Vi vm V x xi yi Z Zc,i Ż ″i

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] Notes

The authors declare no competing financial interest. † D.E.R. is the L. W. Jones Professor of Chemical Engineering, Director of the Yale Sol Reaction Engineering Group. ‡ M.A.-Z. is a visiting Associate Professor. Permanent address: ́ Departamento Fisica Matemática y de Fluidos, UNED, Apdo 60141, 28080 Madrid, Spain.

molar mass of component i (e.g., kg/kmol) kinetic “order” parameter in Ġ (S) expression total number (or molar) density of mixture particle number density distribution function Avogadro’s number (0.6023 × 1027 molecules/(kgmol)) number density of component i (molar or molecular) number of moles of component i molar flow rate of component i particle number density pressure (bar) universal gas constant supersaturation time (reckoned from the onset of AS addition) characteristic time associated with antisolvent (CO2) addition rate absolute temperature (K) critical temperature volume of an individual precipitate particle molar volume of mixture critical molar volume molar volume of component i molecular volume (for API in solid state) total volume of droplet (containing S + API + AS) composition vector {x1,x2,x3}T mole fraction of component i (liquid phase) mole fraction of component i (vapor phase) compression factor (Z ≡ pV/(RT)) critical value of compression factor for component i molar (or molecular) impingement flux of component i

Greek Letters

α αβ β γ Γ ε ξ ϕp

pertaining to phase α pertaining to interfacial phase αβ pertaining to phase β prevailing surface free energy (environment-dependent) dimensionless surface energy (γv2/3 m /(kBT)) incorporation probability for growth species (API) fraction of API precipitated (compared to VSE maximum) API precipitate volume fraction (first volume moment of n(v, t)) ω component acentric factor

■

ACKNOWLEDGMENTS This work was supported by Industrial Affiliates of the Sol Reaction Engineering Research Group at Yale, Yale ChE graduate alumni, and the Yale School of Engineering and Applied Science. Much of the text was prepared while the senior author (D.E.R.) enjoyed the hospitality of Stanford University/HTGL during his Fall 2012 academic leave. M.A.-Z. also acknowledges the support of Ministerio de Ciencia e Innovación (Grant Nos. CSD2010-00011 and ENE201126868) and Comunidad de Madrid (Grant No. S2009/ENE1597) at UNED (Madrid).

Subscripts and Superscripts

* 0 ∞ AS API G i, j L N Ph ref s S T

■

NOMENCLATURE Ḃ ‴ birth rate of particles per unit volume and time d diameter (of particle or droplet) D Fick molecular diffusion coefficient F(GK) correction factor due to interface curvature Ġ linear growth rate of particle or crystal kB Boltzmann constant kG(T,...) dimensionless growth rate constant kij binary interaction parameter for mixture a in Peng− Robinsion EOS lij binary interaction parameter for mixture b in Peng− Robinson EOS mi molecular mass of component i

pertaining to CNT critical size initial value (subscript) final value antisolvent (here CO2) active pharmaceutical ingredient (here Ph (surrogate)) growth components i, j liquid phase nucleation phenanthrene (API surrogate) reference value solid phase solvent (here toluene) toluene (C7H8)

Acronyms/Abbreviations

API AS ASP 4496

active pharmaceutical ingredient (or Ph surrogate) antisolvent (here CO2) antisolvent precipitation dx.doi.org/10.1021/ie4025853 | Ind. Eng. Chem. Res. 2014, 53, 4489−4498

Industrial & Engineering Chemistry Research BDP CNT CSD EOS G GASP GKO IGKT MOC N NDDF N/G NSM ODE PDE Ph PR PSD RESS S SASP SEE SMD T VLSE WMD

■

Article

following: Fundamentals Session, 13th European Aerosol Conference, Prague, Czech Republic, Sep. 1, 2013; AIChE Annual Meeting, Paper/ Poster No. 394c, San Francisco, CA, USA, Nov. 5, 2013. (13) Debenedetti, P. G. Homogeneous nucleation in supercritical fluids. AIChE J. 1990, 36, 1289−1298. (14) Kwauk, X.; Debenedetti, P. G. Mathematical modeling of aerosol formation by rapid expansion of supercritical solutions in a converging nozzle. J. Aerosol Sci. 1993, 24, 445−469. (15) Turk, M. Influence of thermodynamic behaviour and solute properties on homogenerous nucleation in supercritical solutions. J. Supercrit. Fluids 2000, 18, 169−184. (16) Dodds, S.; Wood, J. A.; Charpentier, P. A. Modeling of the GasAntisolvent (GAS) Process for Crystallization of Beclomethasone Dipropionate using Carbon Dioxide. Ind. Eng. Chem. Res. 2007, 46, 8009−8017. (17) Girifalco, L. A.; Good, R. J. A Theory for the Estimation of Surface and Interfacial Energies. 1. Derivation and Application to Interfacial Tension. J. Phys. Chem. 1957, 61, 904−909. (18) Nielsen, A. E.; Sohnel, O. Interfacial Tensions Electrolyte Crystal-Aqueous Solution, from Nucleation Data. J. Cryst. Growth 1971, 11, 233−242. (19) Mersmann, A. Calculation of Interfacial-Tensions. J. Cryst. Growth 1990, 102, 841−847. (20) A more comprehensive account of this mathematical model of the isobaric GASP process and its consequences is: Rosner, D. E.; Arias-Zugasti, M. Theory of Pharmaceutical Powder “Micronization” Using Compressed Gas Anti-solvent (Re-)precipitation. J. Supercrit. Fluids 2014, submitted for publication. In this connection, an overview of our reformulation of crystal growth laws for application to transcritical environments was presented: AIChE Annual Meeting, Paper No. 300f, San Francisco, CA, USA, Nov. 5, 2013, submitted for publication in Cryst. Growth Des. (21) Weber, M.; Russell, L. M.; Debenedetti, P. G. Mathematical modeling of nucleation and growth of particles formed by the rapid expansion of a supercritical solution under subsonic conditions. J. Supercrit. Fluids 2002, 23, 65−80. (22) Kashchiev, D. Nucleation: Basic Theory with Applications; Butterworth-Heinemann: Boston MA, USA, 2000. (23) Jarmer, D. J.; Lengsfeld, C. S.; Randolph, T. W. Nucleation and growth rates of poly(L-lactic acid) microparticles during precipitation with a compressed-fluid antisolvent. Langmuir 2004, 20, 7254−7264. (24) Mahajan, A. J.; Kirwan, D. J. Nucleation and growth-kinetics of biochemicals measured at high supersaturations. J. Cryst. Growth 1994, 144, 281−290. (25) Rosner, D. E. Transport Processes in Chemically Reacting Flow Systems; DOVER: Mineola, NY, USA, 2000. (26) de la Fuente Badilla, J.; Peters, C. J.; de Swaan Arons, J. Volume Expansion in Relation to the Gas-Antisolvent Process. J. Supercrit. Fluids 2000, 17, 13−23. (27) Dixon, D. J.; Johnston, K. P. Molecular Thermodynamics of Solubilities in Gas Antisolvent Crystallization. AIChE J. 1991, 37, 1441−1449. (28) Reverchon, E.; Caputo, G.; De Marco, I. Role of phase behavior and atomization in the supercritical anti-solvent precipitation. Ind. Eng. Chem. Res. 2003, 42, 6406−6414. (29) Rosner, D. E.; McGraw, R.; Tandon, P. Multivariate population balances via moment and Monte Carlo simulation methods: An important sol reaction engineering bivariate example and “mixed” moments for the estimation of deposition, scavenging, and optical properties for populations of nonspherical suspended particles. Ind. Eng. Chem. Res. 2003, 42, 2699−2711. (30) Rosner, D. E.; Arias-Zugasti, M. Bi-variate Population Balance Model of Ethanol-Fueled Spray Combustors. AIChE J. 2011, 57, 3534−3554. (31) Hulburt, H. M.; Katz, S. Some Problems in Particle Technology. A Statistical Mechanical Formulation. Chem. Eng. Sci. 1964, 19, 555− 574.

beclonemethasone-17,21-diproportionate (anti-inflammatory) classical nucleation theory crystal size distribution (function) equation of state growth gas(-induced) antisolvent precipitation Gibbs-Kelvin-Ostwald ideal gas kinetic theory method of characteristics nucleation number density distribution function nucleation/growth Nielsen−Sohnel−Mersmann (eq 1) ordinary differential equation partial differential equation phenanthrene (C14H10; surrogate API) Peng−Robinson (EOS) particle size distribution (function) rapid expansion of a supercritical solution solvent (for API or API surrogate) supercritical (vapor) antisolvent precipitation surface energy evolution (see eq 1) Sauter mean diameter toluene (C7H8) solvent for API or API surrogate vapor−liquid−solid equilibrium well-mixed droplet20,37

REFERENCES

(1) Chang, C. J.; Randolph, A. D. Precipitation of microsize organic particles from supercritical fluids. AIChE J. 1989, 35, 1876−1882. (2) Fages, J.; Lochard, H.; Letourneau, J.-J.; Sauceau, M.; Rodier, E. Particle generation for pharmaceutical application using supercritical fluid technology. Powder Technol. 2004, 141, 219−226. (3) Gallagher, P. M.; Coffey, M. P.; Krukonis, V. J.; Klasutis, N. Gas anti-solvent recrystallization: New process to recrystallize compounds insoluble in supercritical fluids. In Supercritical Fluid Science and Technology, American Chemical Society, 1989; Chapter 22, pp. 334− 356. (4) Chang, S.-C.; Lee, M.-J.; Lin, H.-M. Role of phase behavior in micronization of lysozyme via a supercritical anti-solvent process. Chem. Eng. J. 2008, 139, 416−425. (5) Kikic, I.; Sist, P. Applications of supercritical fluids to pharmaceuticals: Controlled drug release systems. In Supercritical Fluids: Fundamentals and Applications; Kiran, E., Debenedetti, P. G., Peters, C. J., Eds.; NATO Science Series; Kluwer: Dordrecht, The Netherlands, 2000; Vol. E366, pp 291−306. (6) Jung, J.; Perrut, M. Particle design using supercritical fluids: Literature and patent survey. J. Supercrit. Fluids 2001, 20, 179−219. (7) Martin, A.; Cocero, M. J. Micronization Processes with Supercritical Fluids: Fundamentals and Mechanisms. Adv. Drug Delivery Rev. 2008, 60, 339−350. (8) Werling, J.; Debenedetti, P. G. Numerical modeling of mass transfer in the supercritical anti-solvent process. J. Supercrit. Fluids 1999, 16, 167−181. (9) Elvassore, N.; Cozzi, F.; Bertucco, A. Mass Transport Modeling in a Gas Antisolvent Process. Ind. Eng. Chem. Res. 2004, 43, 4935− 4943. (10) Muhrer, G.; Lin, C.; Mazzotti, M. Modeling the Gas Antisolvent Recrystallization Process. Ind. Eng. Chem. Res. 2002, 41, 3566−3579. (11) Bakhbakhi, Y.; Rohani, S.; Charpentier, P. A. Micronization of Phenanthrene Using the Gas Antisolvent Process: Part 2. Theoretical Study. Ind. Eng. Chem. Res. 2005, 44, 7345−7351. (12) This paper is based, in part, on the following: 12th European Aerosol Conference, Paper No. WG10S1O, Granada, Spain, Sep. 6, 2012. A more recent account of this work was presented at the 4497



Article

(32) Rosner, D. E.; Arias-Zugasti, M. Coupling between homogeneous rate processes and fluid deformation rate: Brownian particle coagulation in a rapidly dilating solvent. AIChE J. 2011, 57, 307−318. (33) Cammenga, H. K.; Petrick, H. J.; Schulze, F. W. Sublimation Kinetics of Organic Molecular-Crystals. J. Cryst. Growth 1981, 55, 351−362. (34) Mullin, J. W. Crystal Growth. Crystallization, 4th ed.; Elsevier/ Butterworth-Heinemann: Oxford, U.K., 2001; Chapter 6, pp 224−296. (35) Snyder, R. C.; Doherty, M. F. Predicting crystal growth by spiral motion. Proc. R. Soc., A 2009, 465, 1145−1171. (36) Thiering, R.; Dehghani, F.; Foster, N. R. Current issues relating to anti-solvent micronisation techniques and their extension to industrial scales. J. Supercrit. Fluids 2001, 21, 159−177. (37) Chavez, F.; Dibenedetti, P. G.; Luo, J. J.; Dave, R. N.; Pfeffer, R. Estimation of the characteristic time scales in the supercritical antisolvent process. Ind. Eng. Chem. Res. 2003, 42, 3156−3162.

4498


Surface Energy Evolution in Pharmaceutical Powder Micronization

Recommend Documents