Znd. Eng. Chem. Res. 1993,32, 2118-2128
2118
Integral Equation Study of Microstructure and Solvation in Model Attractive and Repulsive Supercritical Mixtures Jean W.Tom and Pablo G . Debenedetti’ Department of Chemical Engineering, Princeton University, Princeton, New Jersey 08544
The Ornstein-Zernike equation with Percus-Yevick closure was used to investigate the local environment around solute molecules in dilute supercritical mixtures and its relationship t o solubility. Two binary Lennard-Jones mixtures were studied: one, attractive; the other, repulsive. Both systems were studied a t high solute dilution, at supercritical temperature, and over a broad range of reduced densities (0.33 C p / p c < 1.6 for the attractive system; 0.33 C p / p c 2.3 for the repulsive one). The attractive system exhibited significant short-ranged solvent enrichment around the solute. The difference between the average solvent density within the first, second, and third solvation shells and the bulk solvent density was found t o be more pronounced over the approximate range of reduced densities 0.5 C p/pc < 0.8. In contrast, the repulsive system exhibited local solvent depletion, but this effect was particularly pronounced at near-critical density. A simple thermodynamic analysis shows that in dilute attractive supercritical systems the solute’s chemical potential is necessarily a weak function of bulk density a t constant temperature. This was indeed found to be the case. In contrast, the solute’s fugacity coefficient (whose reciprocal is a direct measure of solubility enhancement for attractive systems) was found to be a very sensitive function of density. The solute’s fugacity coefficient was calculated as a function of a variable distance cutoff, beyond which the mixture was assumed to be uniform. The value of the cutoff beyond which no further change in the fugacity coefficient resulted is thus an unambiguous measure of the size of the local region around solute molecules that is important for solubility. Over the range of densities studied here, this quantity ranged from 3 to 5 solvent diameters for the attractive mixture.
I. Introduction More than a century after Hannay and Hogarth (1879) discovered that nonvolatile solids can dissolve in supercritical fluids, the high level of current interest in such solvents is reflected in the existence of a journal (the Journal of Supercritical F l u i d s ) , a biennial symposium at Fall AIChE meetings, and a triennial international symposium, all devoted exclusively to the topic. During this century, the first major study to address solubility in supercritical fluids was done at Delft University by Scheffer and co-workers [e.g., Diepen and Scheffer, 1948a,b; see McHugh and Krukonis (1986) for additional references]. A decade later, solubility in supercritical fluids was studied by Elgin and co-workers at Princeton (e.g., Todd and Elgin, 1955; Chappelear and Elgin, 1961). Elgin’s work was followed by investigations in the former Sovient Union (see, e.g., Tsekhanskaya et al., 1962). Widespread interest in the West, starting in the early-to-mid 19709,has grown continuously ever since. Most of the work on supercritical fluids has traditionally been geared toward extraction and separation applications (McHugh and Krukonis, 1986). Since the mid 19809, however, attention has steadily shifted toward ever more complex and challenging systems (Johnston, 1989). This increased complexity encompasses the chemical structure of the solute molecules (e.g., Wong and Johnston, 19861, the self-assembled nature of some supercritical phases of interest [e.g., reverse micelles (Gale et al., 1987; Smith et al., 1988)], the sophistication of the instrumental and analytical techniques employed in experiments (see, e.g., Bright and McNally, 19921, the nature of the physical transformations now being investigated [e.g., particle formation (Lele and Shine, 1992; Tom and Debenedetti, 1991)1,and the specificityof the chemical transformations
* Author to whom correspondence should be addressed.
that are carried out in a supercritical medium (e.g., Randolph et al., 1988a,b). As is often the case, experiments and applications in this field have far outpaced theoretical studies. The supercritical systems that are currently being investigated by simulations (e.g., Chialvo and Debenedetti, 1992a) or by integral equations (e.g., Wu et al., 1990) are at best reductionist models of reality (typically Lennard-Jones mixtures). In discussing our own work, we have described these models as caricatures (Knutaon et al., 1992). Yet such calculations represent in many cases (simulations,in particular) the current state of the art. This is because the study of solute-solute interactions in dilute supercritical mixtures via molecular dynamics, for example (Chialvoand Debenedetti, 1992a),requires the use of very large systems and long simulations in order to attain statistical significance. Even using codes designed specifically for large systems (Chialvo and Debenedetti, 1992b), such studies would be virtually impossible if more realistic (e.g., multisite) representations of intermolecular potentials were used instead of Lennard-Jones-typeinteractions. Some computational studies of dilute supercritical mixtures involving more realistic pair potentials exist (e.g., Cummings et al., 1991); however, the difference in complexity between real and simulated supercritical systems remains generally large. As for integral equations, calculations for molecular fluids have been done, mostly (Lee, 1988) using the reference interaction site model [RISM (Chandler and Andersen, 197211. Some examples include studies of small-molecule liquids (e.g., Chandler and Andersen, 1972;Chandler, 1988);polymer melts and blends (e.g., Curro and Schweitzer, 1987; 1990), and model supercritical mixtures (Wu et al., 1990). In spite of these limitations, calculations on idealized models play a very important role in supercritical fluids research. Both simulationsand integral equations provide a level of microscopic detail that is inaccessible to even
0888-5885/93/2632-2118$04.00/00 1993 American Chemical Society
Ind. Eng. Chem. Res., Vol. 32, No. 9,1993 2119 the most sophisticated of experimental techniques. Thus, what is lost in predictive power by virtue of simplifications in the representation of intermolecular forces is gained in fundamental understanding of the underlying phenomena through judicious choice of the potential parameters [for example, invoking the attractiverepulsive criteria (Petsche and Debenedetti, 1991)l. This improved microscopic understanding of solvation mechanisms is important practically as well as scientifically. A large body of experimental (Kim and Johnston, 1987a,b; Johnston et al., 1989;Brennecke and Eckert, 1988,1989; Brennecke et al., 1990a,b; Betts et al., 1992a,b; Sun et al., 19921, theoretical (Wu et al., 1990;Munoz and Chimowitz, 1992), and computational (Petache and Debenedetti, 1989; Knutaon et al., 1992;O’Brien et al., 1993)evidence appears to suggest that in dilute supercritical mixtures of practical interest the local environment surrounding solute molecules (microstructure) differs appreciably from the bulk. For example, the remarkable solubility enhancements brought about by the addition of small amounts of cosolventa (Wong and Johnston, 1986) appear to be correlated with an enhancement in the local concentration of cosolvent around the solute (Kim and Johnston, 1987a). Rate and selectivity enhancements in a bimolecular photochemical reaction in a supercritical solvent (Combes et al., 1992) have been interpreted in terms of enhanced solute-solute interactions at high dilutions. Likewise, cholesterol aggregation plays an important role during the substance’s enzyme-catalyzed oxidation in supercritical carbon dioxide under conditions where the mole fraction in the supercritical fluid does not exceed 10-4 (Randolph et al., 1988a,b). Furthermore, appreciable augmentation in the solvent density around solute molecules has been reported in a recent comparative study of the system pyrene-carbon dioxide by fluorescence spectroscopy and molecular dynamic simulations (Knutson et al., 1992). These examples suggest that an improved understanding of the local environment around the solute species (preciselythe type of information provided by simulations and integral equation calculations) could well lead to the tailoring of specific solvent environments for targeted separations or reactions (Brennecke et al., 1990~). While the microstructure has received much attention in recent theoretical and computational work, the emphasis has been almost exclusively descriptive [e.g., radial distribution functions (Chialvo and Debenedetti, 1992); structure and dynamics in the vicinity of solute molecules (Petache and Debenedetti, 1989)l. Comparatively little effort has been devoted to relating the microstructure to bulk thermodynamic behavior. The main purpose of this paper is precisely to discuss quantitatively the relationship between microstructure and thermodynamics. This is done by using integral equations to calculate the extent of the region around the solute beyond which the mixture can be assumed to be homogeneous without this leading to changes in the solute’s fugacity coefficient. This is an unambiguous measure of what is meant by local environment, and we calculate the dependence of its size upon bulk density. Not surprisingly, we find that solvation is entirely a short-ranged phenomenon, and that it is unrelated to near-critical effects, such as long tails in the correlation functions. In this, supercritical solvents are of course no different from ordinary liquids. However, the coexistence of widely different relevant length scales, with microstructure controlling solubility and long-ranged correlations determiningpartial molar property anomalies, is unique to supercritical solutions.
Our work is closely related to the recent fundamental study of Munoz and Chimowitz (1992). These authors investigated microstructure and chemical potentials in an attractive Lennard-Jones system similar to the one studied here, although we study a different repulsive system. However, we differ from these authors in our analysis of the relation between microstructure and bulk behavior. Where Munoz and Chimowitz (1992) focus on the insensitivity of the solute’s chemical potential to changes in microstructure, we show that the fugacity coefficient, in contrast, is extremely sensitive to the nature of the local environment around the solute. This is an important distinction because the relative constancy of the chemical potential is unrelated to solubility enhancement and can be predicted a priori from simple thermodynamic argumenta. This paper is organized as follows. Section I1 discusses integral equations and the numerical scheme used here to solvethem, it presents equations for the different methods used to calculate the solute’s chemical potential, and it gives expressions for several thermodynamic quantities in terms of pair correlation functions. In section 111we contrast the molecular distribution functions in attractive and repulsive mixtures, we relate the microstructure around the solute molecule to its fugacity coefficient (and hence to ita solubility), and we calculate the density dependence of the size of the solute’s local environment beyond which the mixture can be considered homogeneous. The main conclusions are summarized in section IV.
11. Integral Equation Theory and Calculations Theory. An important goal of this work is to calculate the average distribution of pairs of molecules (the pair correlation functions). The starting point is the OrnsteinZernike (OZ) equation for mixtures
where hij(r) and cij(r) are the total and direct correlation functions for a pair of molecules of species i and j , r is the distance between the centers of two molecules, P k is the number density of species k, and n is the number of components in the mixture (n= 2 throughout this study). Equation 1is simply the definition of the direct correlation function. To solve the OZ equations, an independent closure relating the total and direct correlation functions is needed. We use the Percus-Yevick (PY)approximation (Percus and Yevick, 1958) which yields long-ranged correlations in the critical region (McGuigan and Monson, 1990). Incontrast,aclosuresuchaathehypemettedchain (HNC) approximation does not exhibit critical behavior since it has no solutions in the near-critical region (Brey and Santos, 1986). The PY closure reads
+
{ [uf;)])
cij(r)= [l hij(r)l 1-exp -
(2)
where uij(r) is the interaction potential for an ij’pair. The total correlation function, hij(r), is related to the pair correlation function, gij(r),by
gij(r)= h,(r) + 1
(3)
Thermodynamic quantities can be calculated from knowledge of the pair correlation functions. In this work, we calculate the pressure from the pressure equation (Lee, 1988)
2120 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993
where k is Boltzmann’s constant, and uij’(r) is the derivative of uij with respect to r. Alternatively, one can use the compressibility equation, whose validity does not depend on the assumption of a pairwise additive potential (McQuarrie, 1976). The statistical mechanical theory of fluctuations in mixtures (Kirkwood and Buff, 1951)yields an exact relationship between thermodynamic derivatives (isothermal compressibility, partial molar volumes, chemical potential compositional derivatives) and so-called fluctuation integrals Hij = Jhij(r) d3r. The isothermal compressibility KT [=(a In p / d P ) T ~ l ~ 2of1a binary mixture is given by (0 = l / k n
An explicit expression for component i’s chemical potential in a mixture was first developed by Lee (1974) for the PY closure at low densities; this expression was recently extended to all densities by Kjellander and Sarman (19891, and simultaneously by Kiselyov and Martynov (1990) for the PY closure. The latter expression reads Bp; = Bpi - ln[piA~l=
Labik et al., 1985). Here we simply outline the essential features of the method. A detailed discussion can be found in Labik et al. (1985). In the Gillan-Labik approach, the OZ equations are transformed into a set of nonlinear algebraic equations in terms of Fourier transforms of the total and direct correlation functions. The technique combines Newton-Raphson (NR) and direct iteration (DI) methods to solve the system of nonlinear algebraic equations. The variable Yij(r) hij(r) -cij(r) is introduced, and Yij(r), hij(r), and cij(r) are then discretized to yij(ri), hij(ri), and cij(ri) where ri = iAr (i 1, 2, ..., N), with Ar and N representing the grid size and the number of grid points. The total correlation function, hij(r), is truncated to 0 for r L NAr. This algorithm diverges beyond a dimensionless isothermal compressibility, pk TKT, of -15 (Wu et al., 1990)because of the truncation of hij(r), which becomes long-ranged in the highly compressible region. However, all state points used in this study had PkTKT values below 15. Our solution utilized agrid of Ar = 0.02~2 with a potential cutoff at 20.48~2or N = 1024 (here u2 is the solvent’s Lennard-Jones size parameter; see section I11 for details on the intermolecular potentials). The numerical solution of gij(r), hij(r), and cij(r) is used to calculate the fluctuation integrals and thermodynamic properties. Since hij(r) is long-ranged near the critical point, direct iteration of hij(r) to obtain fluctuation integrals is not possible (Wu et al., 1990). An indirect method is then used: it requires calculating the Fourier transform of hij(r) and cij(r)to obtain fiij(k) and ?ij(k),and taking the k 0 limit to obtain
-
1
fiij(0) = 4nKr2hij(r)dr eij(0) = 4r1r2cij(r)dr
(6) whererij = hij(r)-Cij(r), hi = {h2/(2rmikT))’/2, Pik = (pi~k)O.~, h is Planck’s constant, mi is the mass of a molecule of type i, and pi*is the residual chemical potential, defined as the differencebetween component i’s actual chemical potential and its chemical potential in an ideal gas mixture at the same density, temperature, and composition. In the Appendix, chemical potentials calculated with this expression are compared to results based on Kirkwood’s coupling parameter equation (Kirkwood, 1935, 1936) which, for a binary mixture reads
(1- ~ ~ ) u g12(r,t)l ~ ~ ( r dr ) (7) where y1 is the solute’s mole fraction. At infinite dilution (y1- 01, eq 7 can be approximated as
In the Kirkwood equation, is the coupling parameter (0 5 5 1); it has the effect of replacing the interaction of a particular solute molecule with another molecule j by tuij (where j can represent either a solute or a solvent molecule). As shown in the Appendix, solute chemical potentials calculated via eqs 6 and 7 differ by roughly 0.2%. Numerical Solution of the Ornstein-Zernike Equations. To solve the OZ equations with PY closure, we used the Gillan-Labik numerical technique (Gillan, 1979;
(10)
Since cij(r) is short-ranged, eq 10 is integrated using Simpson’s method directly. fiij(0) can then be obtained from ?ij(O) by solving the OZ equations in Fourier space (Wu et al., 1990). Equation 1is reformulated in Fourier space to
A = (I - e)-%
(11) where A is the matrix [fiijl, I is the identity matrix, and is the matrix [eij].
c
111. Results and Discussion Two highly dilute binary Lennard-Jones mixtures, one attractive and one repulsive, were studied. Dilute supercritical solutions have been classified as attractive, weakly attractive, and repulsive by Debenedetti and co-workers (Debenedetti and Mohamed, 1989; Petache and Debenedetti, 1989,1991). The first are of technological interest, as nearly all applications of supercritical fluids involve attractive solutions. In such a class of binary mixture, there is long-ranged solvent enrichment around solute molecules in the neighborhood of the solvent’s critical point, solubility enhancement, and the solute’s partial molar volume and enthalpy are large and negative, diverging to --m at the solvent’s critical point in the limit of infinite dilution. Attractive behavior occurs generally when the solute is larger than the solvent and has a correspondingly larger characteristic interaction energy. Repulsive solutions are mainly of theoretical interest, although experiments with repulsivesupercritical solutions have been reported (Biggerstaff and Wood, 1988a,b). In a repulsive mixture, there is long-ranged solvent depletion
Ind. Eng. Chem. Res.. Vol. 32, No. 9,1993 2121 Table I. Lannard-Jones Parameters for Naphthalene in
a)
co.
naphthalene in COz solutesolute
solventaolvent salutesolvent
P* ij 11 22 12
(A)
qlm
6.199 3.794 4.996
1.634 1.ooO 1.317
ajj
q l k (K) ejjle22 554.4 2.458 225.5 1.OOO 353.4 1.567
around the solute molecules in the neighborhood of the solvent's critical point, and the solute's partial molar volume and enthalpy are large and positive, diverging to at the solvent's critical point in the limit of infinite dilution. Repulsive behavior occurs when the solvent is larger than the solute and has a correspondingly higher characteristic interaction energy. The microstructures around solute molecules tend to be solvent-rich in the attractive case and solvent-lean for repulsive binaries. The mixtures used here have been called naphthalene in carbon dioxide [attractive case; see Cochran and Lee (1989)l and neon in xenon [repulsive case; see Petsche and Debenedetti (1989)l. We will use this terminology, although it should be understood that, except for yielding reasonable values for the pure component critical temperatures and densities with the choice of potential parameters listed below, the Lennard-Jones potential is a poor representation of the actual interaction potentials, especially in the attractive case. The basic features of attractive and repulsive behavior to he discussed below, however, are quite general and hence independent of the details of the intermolecular potentials. Attractive Mixture. The Lennard-Jones parameters of a model attractive mixture, naphthalene in carbon dioxide,are shown in Table I. In this paper, dimensionless temperatures, densities, and distances are denoted by a 2 superscript *, and are given by T* = kTlrz, p* = ~ 0 and r* = riaz9respectively (1= solute; 2 = solvent). We begin by studying the pair correlation functions. This information is needed both for the calculation of thermodynamic properties and to gain a qualitative understanding of the microstructure. The solute-solute and solute solvent pair correlation functions at a solute mole fraction and T* = 1.4, for p* ranging from 0.1 to 0.5 are of 1X shown in Figures 1and 2. The conditions studied here are therefore TiT, = 1.07 and p / p c from 0.32 to 1.61, using the best estimate of the critical point for a Lennard-Jones fluid via direct simulation of phase coexistence,pe* = 0.31 and T,* = 1.31(Smit et al., 1989). Note, however, that the critical constants of the Lennard-Jones fluid in the PY approximation have been estimated to range from T,* = 1.291, pc* = 0.27 (Wu et al., 1990)to T,* = 1.31, pe* = 0.28 (Rowlinson and Swinton, 1982). The first peak in the solute-solute pair correlation function (Figure 1)occurs at r* = 1.8, and represents the direct contact between naphthalene molecules. The distance r* = 1.8 corresponds to the well of the solutesolute potential, which occurs at a separation of 1.83. This peakismaximum(gl1=7,455)atp* =0.18(pr=0.58).The second interesting feature of the solute-solute pair correlation is a broad shoulder occurring approximately between r* = 2.3 and 3.0 at subcritical densities. For densities above p* = 0.25, the shoulder becomes a (second) peak. The location of this second peak, r* = 2.76, suggests a contribution from naphthalene-C0z-naphthalenetriplets. The minimum between the first two peaks decreases with density; it is as high as 2.5 at p* = 0.25 and falls below unity at p* = 0.4. These general features of the solutesolute pair correlation function agree with those found via integral equation calculations (Wu et al., 1990) and molecular dynamics simulations(Chialvoand Debenedetti,
0.3
0.5 8
&
" N
+-
b)
P* 0.5 0 8
-
6
Q
I
4
.G
4 2 0
3
4
r/o2 Figure 1. Density dependenee of the solutesolute pair correlation function for the naphthalene402 system at F = 1.4 and y1 = 1.0 x 10-4
P* ~
~
0
-.
0.3 O.'
/I
r/o2
Figure 2. Density dependenceofthe solutesolvent pair correlation function for the naphthalene402 system at P ' = 1.4 and y1 = 1.0 x 1w.
1992a) for another attractive binary mixture at dilute supercritical conditions. In particular, the solute-solute pair correlations described here provide further evidence of high solute-solute peaks (especially at subcritical densities). The solutesolvent pair correlation (Figure 2) shows a first peak at r* = 1.44, representing the direct contact between naphthalene and COz molecules. The well of the solutesolvent potential occurs at a separation of 1.48. The first peak attains its highest value (g&) = 3.07) at the lowest density used in this study (p* = 0.1). This is in contrast to the solute-solute peak, whose maximum occurs a t p* = 0.18. As discussed in the Introduction, the local distribution of solvent molecules around the solute species is of considerable interest in light of numerous experimental, theoretical, and computational studies (e.g., Kim and Johnston, 1987a,b; Brennecke et al., 1990a,b; Wu et al., 1990; Munoz and Chimowitz, 1992; Betts et al., 1992a,b; Randolph and Carlier, 1992;Carlier and Randolph, 1993; OBrien et al., 1993; Knutson et al., 1992)that suggest the existence of asignificant enhancement in the local solvent density with respect to bulk conditions. The solvation
2122 Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 2 x
251 o,*/,,[
:/-:-- - _ - - - --...-------* , , ,
* , , ,
3
w
1
0.2 ........... r*=3.723
L
0 0.1
0.2
0.3
Bulk Density
0.4
I
, , , ,
1.5
2
2.5
3
3.5
4
r/02
0.5
0.1
,,,,
0
Figure 3. Location of solvationshells used to compute local solvent densities for the naphthalene402 system. Molecular sizes are characterized by the location of the well depth.
s 4
0.4
I4
3.723
7
I
--.....
p*=0.1--p*=0.4 p*=o.2- - - p*=0.5 ........... p 4 . 3
c8
0"
1 I I
t
0.5
Figure 4. Relationshipbetween local and bulk solvent density for the naphthalene402 system at Tu = 1.4 and y 1 = 1.0 X 1O-g. Local and bulk densities are equal along the full line.
environment around a solute molecule determines its solubility in a particular solvent, and can also influence the rates of reactions in which the solute participates. It is therefore a quantity of obvious significance in any molecular-based interpretation of supercritical behavior. In this work, we use the pair correlation functions calculated from integral equation solutions to quantify the asymmetry between local and bulk solvent densities. The solute-solvent pair correlation can be integrated to obtain (Nzl(r)) ,the average number of solvent molecules within a sphere of radius r centered around a solute molecule
The local solvent density withina sphere of radius r around a central solute molecule can then be defined as (Knutson et al., 1992)
Here, the local solvent density is calculated over a range of densities (p* = 0.1-0.5), and at radial locations r* = 2.039,2.600,and 3.723. These represent, respectively, the boundary of the first solvation shell and the midpoints of the second and third shells (using the potential well, rather than the point of zero potential, to characterize molecular size). The geometric representation of these radial distances is depicted in Figure 3. Figure 4 compares local and bulk solvent densities. Deviations from bulk conditions are most pronounced over a range of bulk densities extending approximately between 0.15 and 0.25. This is in agreement with recent electron paramagnetic resonance (EPR)spectroscopy of dilute di-tert-butyl nitroxide radicals in near-critical and supercritical ethane (Carlier
Figure 5. Radial dependence of the ratio of local to bulk density for the naphthalend02 system at !P = 1.4 and y1 = 1.0 X 1o-B.
and Randolph, 1993),where the maximum deviation from bulk density was found at 50 5% of the critical density. The gradual deepening of the first minimum in g12 around r* = 2 (Figure 2) a t higher densities indicates the beginning of liquidlike packing and manifests itself in the crossing of the first and second solvation shell curves at p* = 0.3, and in the crossing of the first solvation shell and bulk lines at p* = 0.35 (Figure 4). Figure 4 is in excellent qualitative agreement with analogous calculations by molecular dynamics for another model attractive LennardJones binary (Knutson et al., 1992). The ratio of local-to-bulk density provides a measure of solvent density enhancement. Figure 5 shows how the average local density (normalized by the bulk density) depends on the distance over which the averaging is performed. The radius of the partially hollow sphere within which density augmentation is largest decreases monotonically with density, from r* = 1.84 at p* = 0.1 to r* = 1.72 at p* = 0.5. The maximum density enhancement always occurs within the first solvation shell. However, significant density augmentation in the attractive mixture occurs even up to three solvation shells. At low and nearcritical densities, the density enhancement persists over large distances away from the solute molecule (for example, = 1.045 at r* = 8.0). However, at a t p* = 0.2, higher densities, the enhancement decays much more rapidly (e.g., at p* = 0.5, p l d p b u l k = 1.01 a t r* = 5.0). This is due to a cancellation between peaks and valleys in the correlation function is liquidlike densities are approached. Thus, we see an augmentation in the solvent's local density with respect to bulk conditions that is particularly pronounced at subcritical densities and decreases rapidly at liquidlike densities. The present integral equation study provides a detailed picture of the microstructure of dilute attractive supercritical solutions that is consistent with recent spectroscopicinvestigations (Brennecke et al., 1990a,b),molecular dynamic calculations (Knutson et al., 19921, as well as other integral equation studies (Wu et al., 1990; Munoz and Chimowitz, 1992). In addition, the use of integral equations provides the means of relating microstructure to bulk thermodynamic behavior through calculation of chemical potentials. We are specifically interested in understanding solubility at supercritical conditions, and how it depends on the microstructural features discussed so far. We begin by recalling that the chemical potential, pi, fugacity, fi, and fugacity coefficient, @i,of component i in a mixture are defined as follows pi
= kT In f i fi
+ kT ln(@A:)
= YiP@i
(14) (15)
where yi is the mole fraction of component i, and Ai has
Ind. Eng. Chem. Res., Vol. 32, No. 9,1993 2123 already been defined. For a pure, incompressible, nonvolatile solid a t a temperature T and pressure P, the fugacity can be written as (Model1 and Reid, 1983)
where Pvap is the vapor pressure at a temperature T,Usolid is the solid‘s molecular volume, and we assume that the vapor phase in equilibrium with the solid at Pvap and T behaves ideally. The term exp [umlid(P-Pvap)/km is called the Poynting correction factor. When a pure nonvolatile solid 1 (e.g., naphthalene) is in equilibrium with a binary supercritical solution [e.g., carbon dioxide (2) plus naphthalene (111, we have the condition f1 = fsolid. Invoking eqs 15 and 16, and noting furthermore that in such situations P >> Pvap, we have
0.1
0.2
0.3
0.4
0.5
P* Figure 6. Density dependence of the solute’s chemical potential for the naphthalene402 system at r* = 1.4 and y1 = 1.0 X 1O-e. 10’
lo4
2
op
Because of the small value of Usolid, the chemical potential of the solid changes very little with applied pressure. For example, increasing the pressure from 1to 200 bar at 300 K causes a relative change of only 11%in the chemical m3 potential of a substance with Usolid = 1.66 X (corresponding to a density and molecular weight of 1g/cm3 and 100 g/mol, respectively). Because at equilibrium the chemical potential of component 1 in the solid and supercritical phase are equal, one would expect any realistic model of dilute attractive supercritical solutions to yield solute chemical potentials that change very little with bulk density, regardless of changes in microstructure. Note that this must be so even though the solid phase is not included in the model. This implies that the solute’s chemical potential, being virtually fixed, is not a very informative quantity in analyzing solvation and solubility at supercritical conditions. Instead, we rearrange eq 17 to read
The right-hand-side numerator is the ratio of actual to ideal-gas-predicted solubility at fixed temperature and pressure: this quantity is usually referred to as solubility enhancement or enhancement factor (e.g., Schmitt and Reid, 1985). It changes over several orders of magnitude as the fluid is compressed isothermally. In sharp contrast, the denominator (Poynting correction) is quite insensitive to pressure (e.g., at 300 K, and for Urnlid = 1.66 X m3, the Poynting correction changes only from 1 to 1.49 in going from 1to 100 bar). Thus, it is the solute’s fugacity coefficient (and not its chemical potential) that should be studied when trying to relate solubility to microstructure. In contrast with Munoz and Chimowitz (1992),we do not seek to explain the relative constancy of the solute’s chemical potential. In real situations, this constancy is imposed by virtue of the presence of a solidphase in contact with a saturated supercritical solution. We in fact use this relative constancy as a test of the soundness of our calculations. In what follows we study $1 and its dependence on both local and bulk conditions. Of course, we do not model a saturated solution, but one that approaches infinite dilution. However,we know from the work of Debenedetti and Kumar (1986) that, for dilute mixtures, 1> &/@lo> e-l, with
41(TP,~1) = dl” exp[-K(TQ)y,I and where
K >0
(19)
p?
io3
@ 0
lo2
8 B 0
6’
z
10’ 0.1
0.2
0.3
0.4
0.5
-
v;
P‘ Figure 7. Density dependence of the solute’s fugacity coefficient and its reciprocal for the naphthalene-COz system at !P = 1.4 and y1 = 1.0 x 1 ~ .
dl”(TP) = lim d1(TQ,y1)
(20)
Y1 - 0
Thus, to within a factor of order 1(and recall that we are studying a quantity, 41, that changes over several orders of magnitude across the range of conditions of interest here), we can understand 41 by computing 41”. From eq 6, we calculate the solute’s chemical potential using the pair correlation functions, and we then relate this quantity to the fugacity coefficient. As shown in the Appendix, the relationship between 1.11 and 41 is
41= ~ X P ( P C- 1n L ~P
~ A ~ ~ I ~ $ ~ (21)
where P is the calculated pressure a t the given T, p , y1. The calculated solute chemical potential is shown in Figure 6. As the bulk density changes from p* = 0.1 to 0.5, the solute’s chemical potential changes by only 9%. Chemical potential values are within 10% of those calculated by Munoz and Chimowitz, who used slightly different Lennard-Jones parameters for the naphthaleneC02 system. This insensitivity is in agreement with the previous discussion and is taken as a confirmation of the surprising extent to which calculations on simplified models can shed light on important and subtle aspects of supercritical solution behavior. The trend shown in Figure 6 (though not our interpretation of it) is in agreement with the calculations of Munoz and Chimowitz (1992). Figure 7 shows the density dependence of 41 and 41-l (--solubility enhancement). The extreme sensitivity to density is obvious. We now seek to relate this trend to changes in the microstructure that occur as the solvent is compressed isothermally. The relationship between microstructure and fugacity coefficients (or, equivalently,between microstructure and thermodynamics) is examined by computing the fugacity coefficient using a radial cutoff, R, in the integrals for the solute’s chemical potential (eq 6). Beyond this cutoff, the
2124 Ind. Eng. Chem. Res., Vol. 32, No. 9,1993 '_
E
P*
I
k
-
1.5 1
N
0.5
0
Figure 8. Dependence of the computed solute's fugacity coefficient upon thecutoffdistance beyondwhich themixtureisassumed ta be homogeneous. NsphthaleneCO~system at 'P = 1.4 and 3 I = 1.0 x 10-9.
r/02 F-re 10. Densitydependenceofthesolutesolventpaircorrelation function for the neon-xenon system at T* = 1.4 and y ~ 1=.0 X l e. 1.05
6 5.5
. a"
5
4.5
c r : 4
3.5 3
2.5 0
0.1
0.2
0.3
0.4
0.5
0.6
P* Figure 9. Density dependenceof the cutoff distance (beyond which the mixture is assumed to he homogeneous) needed to obtain 98%, 99%, and 99.5% of the asymptotic salute fugacity coefficient. NaphthaleneCOl system at P = 1.4 and y1 = 1.0 X 10-9.
fluid is assumed to be random or unstructured @;j(r)= 1, for r > R and i, j = 1,2). Within this cutoff (r < R),pair correlation values from integral equations are used. This calculation provides a direct measure of the extent of the fluid region surrounding the solute molecules that contributes to that species' fugacity coefficient (and therefore to its solubility). Figure 8 shows the fugacity coefficient at p* = 0.3 as a function of the cutoff distance used in the evaluation ofthe integrals ineq 6. The fugacity coefficient approaches a constant asymptote as R increases. Similar trends are found at the other densities investigated. The radial cutoff required in the evaluation of the integrals in order to obtain a quantity equal to 99.5,99.0,and 98% of the true (R -) fugacity coefficient is plotted in Figure 9 as a function of bulk density. There is a noticeable dependence of the cutoff on bulk density. As the critical density is approached, larger cutoffs are needed. In contrast, Munoz and Chimowitz (1992)found the radial cutoff to obtain 90% of the solute's chemical potential to be relatively invariant with bulk density (0.1< p* < 0.4), withR* ranging between 2.9 and 3.3 at T* = 1.32 and 1.37 for a similar attractive mixture. We find that it is the local environment (microstructure),Le., 3 < R* < 5.5, or 11A > m;hence e22
Ind. Eng. Chem. Res., Vol. 32, No. 9, 1993 2125 CI -4 (0
5 5 2 0
>
1.05 1
0.95 0.9 0.85
m
Figure 12. Three-dimensional representation of the radial dependence of the ratio of local hulk density for the neon-xenon system at I* = 1.4 and y, = 1.0 X l V . I
- 1
0
0.1
0.2
0.3,
0.4
0.5
0.6
P Figure 13. Densitydependenceofthecutoffdistance(heyondwhich themixtweisassumedthe homogeneous)neededtoohbin 100.01% of the asymptotic solute fugacity coefficient. Neon-xenon at l“ = 1.4 and y1 = 1.0 X 1W9.
> €12) at moderate density. Note that rapid disappearance of repulsive behavior upon compression. Figure 13 shows the density dependence of the radial cutoff required in order to obtain a quantity equal (to -) solute fugacity within 0.01%) to the true (R coefficient. The calculation is identical to that for the attractive case except that the cutoff corresponds to 100.01% of the true solute fugacity coefficient because this quantity, when plotted as a function of the radial cutoff, exhibits a maximum beyond which it decays to its bulk value from above. Once again, larger cutoffs are needed close to criticality, but it is local environment (microstructure), Le., 1.3 < R* < 2.8, or 6 A < R < 12 A, which determines the fugacity coefficient. Note the asymmetryofthecurveshowninFigure13whencompared to the corresponding attractive curves of Figure 9.
-
IV. Conclusions In this work, we have used integral equations to investigate the mutual spatial distribution of solvent and solute molecules in two binary Lennard-Jones mixtures: one, attractive, with potential parameters chosen to represent naphthalene in carbon dioxide; the other, repulsive, with potential parameters chosen to represent neon in xenon. We have placed emphasis on understanding differencesbetween local conditions surrounding solute molecules (microstructure) and bulk conditions. As with other types of mixtures, the microstructure is important in determining the solute’s solubility in a given solvent, and hence we also investigated the relationship between microstructure and solute fugacity coefficients. What is unique to supercriticalmixturesis thecoexistenceofwidely different relevant length scales, with microstructure controlling solubility, and long-ranged critically-driven fluctuations controlling the temperature and pressure derivatives of solubility.
In the attractive mixture, the microstructure is found to be solvent-rich with respect to the bulk. This enrichment persists even when the solvent density is averaged over three solvation shells, and is more pronounced in the reduceddensityrangeO.5 Rm, = 20.48~2. To calculate the chemical potential without using coupling parameters, the explicit closed formula initially developed by Lee (1974) and recently reformulated by Kjellander and Sarman (1989) was used:
B. Chemical Potential Calculations. In order to verify the accuracy of eq 6, we also calculated the solute’s chemical potential using Kirkwood’s coupling parameter equation
where A1 = (h2/(2arn1kT))1/2;gll(r,t) and g12(r,5) are calculated with the solute-solute and solute-solvent interactions coupled to the extent 5, i.e., ell([) = E e l 1 and e&) = &12. The Lennard-Jones potentials ull(r) and u12(r) are the full potentials, ell(1) and elz(1). The integration over the coupling parameter is then divided into three parts:
(B9) whereElk(Ylk(r)) is afunction of Ylk(r),and ylk(r) = hlk(r) - clk(r). The exact form of E(y(r))depends on the closure used for the Ornstein-Zernike equation. For the PercusYevick equation, E(y(r)) = ln[l + y(r)l - y(r). The final form of the explicit expression for the solute’s chemical potential is
where Z is the integrand J4a[plull(r) gll(r,D + pzulz(r) g1z(r,5)lr2 dr and 7 is the smallest value of coupling parameter used to calculate pair correlation functions.The first integral is evaluated by using the approximation of Lee and Hubert (1973)
resulting in the following generic integrand
(BIOI The chemical potentials of naphthalene in the naphthalene-C0z mixture using both methods of calculation
Ind. Eng.Chem. Res.,Vol. 32, No. 9,1993 2127 were compared a t p* = 0.1 and 0.3. At these two densities, the values were within 0.2% of each other. C. Derivation of Fugacity Coefficient from Solute Chemical Potential. The left-hand side of eq B1 is the dimensionless residual chemical potential of the solute (BPI'). Denoting the right-hand side of that equation by the symbol A, we have /3pl - In plA13 = A
(C1) For a real mixture, A # 0 and we have j3P = pZ where Z is the compressibility factor. Substituting p1 = y1p and p = /3PZ-1 into eq C1, we have
+
j3pl = ln[j3P~'ylAl31 A
(C2)
Rearranging, j3pl = ln[P(exp A)Z-'yll
+ ln[j3A:]
(C3)
Comparing eq C3 with eqs 14 and 15, we obtain, finally,
Literature Cited Betta, T. A.; Zagrobelny,J.; Bright, F. V. SpectroscopicDetermination of Solute-Fluid Cluster Size in Supercritical NzO. J. Am. Chem. SOC. 1992a, 114,8163. Betta, T. A.; Zagrobelny,J.; Bright, F. V. Elucidation of Solute-Fluid Interactions in Supercritical CFsH by Steady-State and TimeResolved Fluorescence Spectroscopy. In Supercritical Fluid Technology;Bright, F. V., McNally, M. E., Eds.; ACS Symposium Series 488;American Chemical Society: Washington, DC, 1992b. Biggerstaff, D. R.; Wood, R. H. Apparent Molar Volumes of Aqueous Argon, Ethylene, and Xenon from 300 to 716K. J. Phys. Chem. 1988a, 92,1988. Biggerstaff, D. R.; Wood, R.H. Apparent Molar Heat Capacities of AqueousArgon, Ethylene, and Xenon at Temperatures up to 720K and Pressures to 33 MPa. J. Phys. Chem. 1988b, 92,1994. Brennecke, J. F.; Eckert, C. A. Molecular Interactions from Fluorescence Spectroscopy. Proc. Znt. Symp. Supercrit. Fluids (Z) 1988,263. Brennecke, J. F.; Eckert, C. A. Fluorescence Spectroscopy Studies of Intermolecular Interactions in Supercritical Fluids. In Supercritical Fluid Science and Technology; Johnston, K. p., Penninger, J. M. L., Eds.; ACS Symposium Series 406; American Chemical Society: Washington, DC, 1989. Brennecke, J. F.; Tomasko, D. L.; Peshkin, J.; Eckert, C. A. Fluorescence Spectroscopy Studies in Dilute Supercritical Solutions. Znd. Eng. Chem. Res. 1990a, 29, 1682. Brennecke, J. F.; Tomasko, D. L.; Eckert, C. A. Naphthalene/ Triethylamine Exciplex and Pyrene Excimer Formation in Supercritical Fluid Solutions. J. Phys. Chem. 1990b, 94, 7692. Brennecke, J. F.; Debenedetti, P. G.; Eckert, C. A,; Johnston, K. P. AZChE J. 199Oc, 36, 1927. Brey, J. J.; Santos, A. Absence of Criticality in the Hypernetted Chain Equation for a Truncated Potential. Mol. Phys. 1986,57, 149. Bright, F. V., McNally, M. E., Eds. Supercritical Fluid Technology; ACS Symposium Series 488: American Chemical Societv: Washington; DC, 1992. Carlier, C.; Randolph, T. R. Dense-Gas Solvent-Solute Clusters at Near-Infinite Dilution: EPR Spectroscopic Evidence. AIChE J. 1993, 39, 876. Chandler, D. In The Liquid State and Its Electrical Properties; Kunhardt, E.E., Chrietophorou,L. G.,Lueasen,L. H.,Eds.;Plenum Press: New York, 1988. Chandler, D.; Andersen, H. C. Optimized Cluster Expansions for ClassicalFluids 11. Theory of Molecular Liquids. J.Chem. Phys. 1972,57, 1930. Chappelear, D. C.; Elgin, J. C. Phase Equilibria in the Critical Region: Binary Systems with Chlorotrifluoromethane. J. Chem. Eng. Data 1961,6,415. Chialvo, A. A.; Debenedetti, P. G Molecular Dynamics Study of Solute-Solute Microstructure in Attractive and Repulsive Supercritical Mixtures. Znd. Eng. Chem. Res. 1992a, 31, 1391.
Chialvo, A. A.; Debenedetti, P. G. An Automated Verlet Neighbor List Algorithm with a Multiple Time-Step Approach for the Simulation of Large System. Comput. Phys. Commun. 1992b, 70,467. Cochran, H. D.; Lee, L. L. Solvation Structure in Supercritical Fluid Mixtures based on Molecular Distribution Functions. In Supercritical Fluid Science and Technology; Johnston, K. p., Penninger, J. M. L., Eds.;ACS Symposium Series 406;American Chemical Society: Washington, DC, 1989. Cochran, H. D.; Johnson, E.; Lee, L. L. Molecular Theory of Ternary Supercritical Solutions. I. Theory and Methods of Calculation. J. Supercrit. Fluids. 1990,3, 157. Combes,J. R.;Johnston, K. P.; O'Shea, K. E.;Fox, M. A. The Influence of Solvent-Solute and Solute-Solute Clustering on Chemical Reactions in Supercritical Fluids. In Supercritical Fluid Technology;Bright, F. V.; McNally,M. E., Eds.; ACS SymposiumSeries 488; American Chemical Society: Washington, DC, 1992. Cummings, P. T.; Cochran, H. D.; Simonson, J. M.; Mesmer, R. E.; Karaboni, S. Simulation of Supercritical Water and of Supercritical Aqueous Solutions. J. Chem. Phys. 1991,94,5606. Curro, J. G.; Schwiezer, K. S. Theory of Polymer Melts: An Integral Equation Approach. Macromolecules 1987,20,1928. Curro, J. G.; Schweizer,K. S. An Integral Equation Theory of Polymer Blends: Athermal Mixtures. Macromolecules 1990,23,1402. Debenedetti, P. G.; Kumar, S. K. Infinite Dilution Fugacity Coefficienta and the General Behavior of Dilute Binary Systems. AZChE J. 1986.32.1253. Debenedetti, P. G.; Mohammed, R. S. Attractive, Weakly Attractive, and Repulsive Near-Critical Systems. J. Chem. Phys. 1989,90, 4528. Diepen, G. A. M.; Scheffer, F. E. C. On Critical Phenomena of Saturated Solutions in Binary Systems. J.Am. Chem. SOC. 1948a, 70, 4081. Diepen, G. A. M.; Scheffer, F. E. C. The Solubility of Naphthalene in Supercritical Ethylene. J. Am. Chem. SOC.1948b, 70, 4085. Gale, R. W.; Fulton, J. L.; Smith, R. D. Organized Molecular Assembliesin the Gas Phase: Reverse Micellesand Microemulsions in Supercritical Fluids. J. Am. Chem. SOC. 1987, 109, 920. Gillan, M. J. A New Method of Solving the Liquid Structure Integral Equations. Mol. Phys. 1979, 38, 1781. Hannay, J. B.; Hogarth, J. On the Solubility of Solids in Gases. Proc. R. SOC. London 1879,29,324. Johnston, K. P. New Directions in Supercritical Fluid Science and Technology. In Supercritical Fluid Science and Technology; Johnston, K. P., Penninger, J. M. L., Eds.; ACS Symposium Series 406; American Chemical Society: Washington, DC, 1989. Johnston, K. P.; Kim, S.; Combes, J. Spectroscopic Determination of Solvent Strength and Structure in Supercritical Fluid Mixtures: A Review. In SupercriticalFluid Science and Technology; Johnston, K. P., Penninger, J. M. L., Eds.; ACS Symposium Series 406; American Chemical Society: Washington, DC, 1989. Kim, S.; Johnston, K. P. Clustering in Supercritical Fluid Mixtures. AZChE J. 1987a, 33,1603. Kim, S.;Johnston, K. P. Molecular Interactions in Dilute Supercritical Fluid Solutions. Znd. Eng. Chem. Res. 1987b, 26, 1207. Kirkwood, J. G. Statistical Mechanics of Fluid Mixtures. J . Chem. Phys. 1935,3,300. Kirkwood, J. G. Statistical Mechanics of Liquid Solutions. Chem. Rev. 1936, 19, 275. Kirkwood, J. G.; Buff, F. P. The Statistical Mechanical Theory of Solutions I. J. Chem. Phys. 1951, 19, 774. Kiselyov,0.E.; Martynov, G. A. Free Energy in the StatisticalTheory of Fluids. J. Chem. Phys. 1990,93, 1942. Kjellander, R.; Sarman, S. The Chemical Potential of Simple Fluids in a Common Class of Integral Equation Closures. J. Chem. Phys. 1989,90,2768. Knutaon, B. L.; Tomasko, D. L.; Eckert, C. A.; Chialvo, A. A.; Debenedetti, P. G. Local Density Augmentation in Supercritical Solutions. In Supercritical Fluid Technology; Bright, F. V., McNally, M. E., Eds.; ACS Sympoaium Series 488; American Chemical Society: Washington, DC, 1992. Labik, S.; Malijevsky, A.; Vonka, P. A Rapidly Converging Method of Solving the OZ Equation. Mol. Phys. 1985,56, 709. Lee, L. L. Correlation Functions of C h i d Fluids. 111. The Method of Partition Function Variation Applied to the Chemical Potential: Cases of PY and HNC2. J. Chem. Phys. 1974,60, 1197. Lee, L. L. Molecular Thermodynamics of Nonideal Fluids; Butterworths: Boston, 1988. Lee, L. L.; Hubert, H. M. Correlation Functions of Classical Fluids. I. The Radial Distribution Functions of Mixtures of Kihara
2128 Ind. Eng. Chem. Res., Vol. 32,No. 9,1993 Molecules in the Percus-Yevick Approximation and Their Thermodynamic Functions. J. Chem. Phys. 1973,58,44. Lele, A.; Shine, A. D. Morphology of Polymers Precipitated from a Supercritical Solvent. AIChE J. 1992,38, 742. McGuigan, D. B.; Monson, P. A. Analysis of Infinite Dilution Partial Molar Volumes Using a Distribution Function Theory. Fluid Phase Equilib. 1990,57,227. McHugh, M. A.; Krukonis, V. J. Supercritical Fluid Extraction Principles and Practice; Butterworths: Boston, MA, 1986. McQuarrie, D. A.StatisticalMechanics; Harper & Row: New York, 1976. Modell, M.; Reid, R. C. Thermodynamics and its Applications; Prentice Ha& Englewood Cliffs, NJ, 1983. Munoz, F.; Chimowitz, E. H. Integral Equation Calculations of the Solute Chemical Potential in a Near-Critical Fluid Environment. Fluid Phase Equilib. 1992,71, 237. OBrien, J. A.; Randolph, T. W.; Carlier, C.; Shankar, G. Quasicritical Behavior of Dense-Gas Solvent-Solute Clusters at Near-Infinite Dilutions. AZChE J. 1993,39,1061. Percus, J. K.; Yevick, G. J. Analysie of ClassicalStatistical Mechanics by Means of Collective Coordinates. Phys. Reo. 1968,110,1. Petache, I. B.; Debenedetti, P. G. Solute-Solvent Interactions in Infinitely Dilute Supercritical Mixtures: A Molecular Dynamics Investigation. J. Chem. Phys. 1989,91,7075. Petache, I. B.; Debenedetti, P. G. Influence of Solute-Solvent Asymmetry upon the Behavior of Dilute Supercritical Mixtures. J. Phys. Chem. 1991,95,386. Randolph, T. W.; Carlier, C. Free-Radical Reactions in Supercritical Ethane: A Probe of SuDercriticalFluid Structure. J.Phys. Chem. 1992,96,5146. Randolph, T. W.; Blanch, H. W.; Prausnitz, J. M. Enzyme-Catalyzed Oxidation of Cholesterol in Supercritical Carbon Dioxide. AZChE J. 1988a,34,1354. Randolph, T. W.; Clark, D. S.; Blanch, H. W.; Prausnitz, J. M. Enzymatic Oxidation of Cholesterol Aggregates in Supercritical Carbon Dioxide. Science 1988b,239,387.-
Rowlinson, J. S.; Swinton, F. L. Liquids and Liquid Mixtures; Butterworths: Boston, 1982. Schmitt, W. J.; Reid, R. C. The Influence of the Solvent Gas on Solubility and Selectivity in Supercritical Extraction. In Supercritical Fluid Technology; Penninger J. M. L., Radoez, M., McHugh, M. A., Krukonis, V. J., Eds.; Elsevier Science Publishers: Amsterdam, Netherlands, 1985. Smit, B.; De Smedt, P.; Frenkel, D. Computer Simulations in the Gibbs Ensemble. Mol. Phys. 1989,68,931. Smith, R. D.; Yonker, C. R.; Fulton, J. L.; Tingey, J. M. Organized Molecular Assembliesin Supercritical Fluids; Metal Ion Chelates and Reverse Micelles. J. Supercrit. Fluids. 1988,1,7. Sun, Y.; Fox, M. A.; Johnston, K. P. Spectroscopic Studies of p-(N,NDimethy1amino)benzonitrileand Ethyl p-(N,NDimethylarnino)benzoate in Supercritical Trifluoromethane, Carbon Dioxide,and 1992,114,1187. Ethane. J. Am. Chem. SOC. Todd, D. B.; Elgin, J. D. Phase Equilibria in Systems with Ethylene above ita Critical Temperature. AIChE J. 1966, I, 20. Tom, J. W.; Debenedetti, P. G. Formation of Bioerodible Polymeric Microapheres and Microparticles by Rapid Expansion of Supercritical Solutions. Biotechnol. Prog. 1991,7,403. Tsekhanskaya, Y.; Iomtev, M. D.; Mushkina, E. V. Solubility of Diphenylamine and Naphthalene in Carbon Dioxide under Pressure. Russ. J. Phys. Chem. (Engl. Transl.) 1962,36,1177. Wong,J. M.; Johnston, K. P. SolubilizationofBiomoleculeein Carbon Dioxide Based Supercritical Fluids. Biotechnol. Prog. 1986,2,29. Wu, R.S.;Lee, L. L.; Cochran, H. D. Structure of Dilute Supercritical Solutions: Clustering of Solvent and Solute Molecules and the Thermodynamic Effects. Ind. Eng. Chem. Res. 1990,29,977. Received for review January 22, 1993 Revised manuscript received May 18,1993 Accepted June 11, 1993. Abstract published in Advance ACS Abstracts, August 15, 1993. Q