Sublimation Pressure Calculated from High-Pressure Gas−Solid

The proposed method allows one to calculate sublimation pressures of .... numerical techniques, which iteratively refine a single solution vector as t...
0 downloads 0 Views 152KB Size
Download PDF
4824

Ind. Eng. Chem. Res. 2005, 44, 4824-4833

Sublimation Pressure Calculated from High-Pressure Gas-Solid Equilibrium Data Using Genetic Algorithms Jose´ O. Valderrama*,†,‡ and Jack Zavaleta‡,§ Faculty of Engineering, Mechanical Engineering Department, University of La Serena, Casilla 554, La Serena, Chile, Centro de Informacio´ n Tecnolo´ gica, Casilla 724, La Serena, Chile, and Faculty of Engineering, Chemical and Textile Engineering Department, National Engineering University, Lima, Peru´

Sublimation pressures of pure solids are calculated from high-pressure solubility data using genetic algorithms. The Peng-Robinson equation of state with the mixing rules proposed by Wong and Sandler is used as the thermodynamic model to evaluate the fugacity coefficients in the classical solubility equation and the variables P-y, for given values of T-x. The van Laar model was incorporated to evaluate the excess Gibbs free energy included in the Wong-Sandler mixing rule. The sublimation pressure is considered as a parameter to be determined by regression analysis of experimental solubility data. Thus, an optimization problem, in which the difference between the correlated and experimental data of solubility is to be minimized, is solved using a method based on genetic algorithms. This method uses biologically derived techniques such as inheritance, mutation, natural selection, and recombination to find the optimum solution of the optimization problem. Five gas-solid systems, including 16 isotherms and a total of 344 P-T-y data points, were used for the study. The systems studied were binary mixtures containing supercritical carbon dioxide with naphthalene, biphenyl, phenanthrene, anthracene, and pyrene. The proposed method allows one to calculate sublimation pressures of solids with high accuracy. Introduction Phase equilibrium calculations are of fundamental importance in several separation processes, and different methods have been suggested for phase equilibrium property calculations needed in practical applications such as process design and process simulation and in more fundamental applications such as stability analysis or thermodynamic consistency tests.1,2 One of the common methods used in the literature to correlate and predict the phase equilibrium is the so-called “equation of state method”. The use of this method requires an equation of state (EoS) that well relates the variables temperature, pressure, and volume and appropriate mixing rules to express the dependence of the EoS parameters on the concentration. There are several industrially important EoSs, the most popular ones being the so-called cubic equations derived from the equation first proposed by van der Waals. A recent review by one of the authors gives a detailed picture of the present use and applications of these types of equations.3 All practical cubic equations when applied to mixtures involve the use of mixing rules, which include empirical binary interaction parameters kij’s usually calculated from experimental phase equilibrium data. In modern models such as quadratic mixing rules, mixing rules including Gibbs free energy models, or Mansoori and co-workers’s mixing rules, among others, more than one parameter are included. The mixing rule of Wong and Sandler,4 for instance, includes the interaction param* To whom correspondence should be addressed. E-mail: [email protected]. † University of La Serena. ‡ Centro de Informacio´n Tecnolo´gica. § National Engineering University.

eters for the combining rules and the parameters for the chosen Gibbs free energy model included in the mixing rule.5 Phase equilibrium calculations of a solid (2) dissolved in a compressed gas (1), at a pressure P and a temperature T, also called solubility calculations, can be performed using the fundamental equation of phase equilibrium, which leads to a simple equation: y2 ) (Ps2/Pφs2){exp[Vs2(P - Ps2)/RT]}, with φs being the fugacity coefficient of the solid component in the highpressure gas, Ps2 the sublimation pressure of the pure solid, Vs the molar volume of the solid, and R the ideal gas constant.6 Of all of these properties involved in the calculation of the solubility of the solid in the highpressure gas, the sublimation pressure has received less attention in the literature, although it is directly related to the solubility, as seen in the preceding equation. These aspects are discussed later in the paper. The molar volume has a relative influence on the calculations, and the fugacity coefficient φs can be estimated from an appropriate EoS and mixing rules. The sublimation pressure of the pure solid is usually small for common industrially important solids, and experimental techniques cannot be, in many cases, used to obtain accurate values. Therefore, it is proposed here to evaluate the sublimation pressure of a pure solid from the solubility data of the solid in a high-pressure gas. For this, the Peng-Robinson (PR) EoS7 was incorporated into the solubility equation. The mixing rule proposed by Wong and Sandler2 was used, and the van Laar model was included to evaluate the excess Gibbs free energy that appears in this mixing rule. The model parameters and the sublimation pressures are then determined by regression analysis of the solubility data. When the solubility data are selected for analysis, at least these three conditions must be fulfilled to ac-

10.1021/ie0501529 CCC: $30.25 © 2005 American Chemical Society Published on Web 05/19/2005

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005 4825

curately estimate the sublimation pressure of a pure solid using the solubility data of that solid in a compressed gas: (i) the data must be accurate enough and hopefully tested for thermodynamic consistency; (ii) the solubility data must correspond to “true” solid-gas data, hopefully tested for phase stability; (iii) a minimum number of data points, of at least double the model parameters to be calculated, must be available. The majority of the existing methods for solving the phase equilibrium are local in nature and at best yield only local solutions. It has been demonstrated that, for cases such as those containing supercritical carbon dioxide, multiple solutions (local optimum values) are found in a range that seems to be acceptable for correlation purposes.8,9 The optimum value of the interaction parameters depends on the searching interval and on the initial value of the interaction parameters used to start the iterative procedure. In this sense, genetic algorithms (GAs) represent a good alternative to reduce computer time and to find the model parameters for phase equilibrium calculations such as the one analyzed in this paper. Similar to other optimization methods, an acceptable searching interval must be defined to make efficient use of the GA method proposed here and to find values of the model parameters with reasonable physical meaning. Although several applications of GAs have been presented in the literature, few applications to phase equilibrium properties have been done. Rangaiah10 evaluated the application of GAs for phase equilibrium and stability problems. This author studied two stochastic global optimization techniques (GA and simulated annealing) and evaluated their capabilities to analyze phase equilibrium and stability problems. Li and Yang11 used GA to generate the nonlinear binary interaction parameter of the PR equation and correlate phase equilibrium data of binary systems benzoic acid + ethane and benzoic acid + carbon dioxide. Vijande et al.12 used GAs to evaluate the optimum parameters for the interaction of the carbonate-hydroxyl pair using the Nitta-Chao model and experimental data of the vaporliquid equilibrium, activity coefficients at infinite dilution, excess molar enthalpies, and excess molar volumes of carbonate + 1-alkanol binary mixtures. To the best of our knowledge, GAs have not been applied to determine the sublimation pressure of the solid using the solubility data, as proposed here. The model parameters that appear in the mixing rules and the sublimation pressure that appears in the solubility equation are used as optimization variables, and an appropriate objective function is defined. It is also shown that, because the sublimation pressure is a property of the pure solid, its value can be determined using the solubility data of the solid in any compressed gas.

If the solid phase is considered to be a pure substance, then

Ps2φs2 ) y2φg2P

Because the sublimation pressure is normally low, an ideal gas behavior for the gas phase over the pure solid can be assumed, and φs2 ≈ 1. Also, the volume of the solid is considered to be pressure-independent, and the mole fraction of the solute in the gas phase, or solubility, at temperature T and pressure P for component 2 is6 s

y2 )

fs2 ) fg2

(1)

s

Ps2e(V2/RT)(P-P2)

(3)

Pφg2

Here, Ps2 is the sublimation pressure of the pure solid, Vs2 is the solid molar volume, all at temperature T, and φ2 is the fugacity coefficient of the solid (2) at pressure P. The fugacity coefficient is calculated from standard thermodynamic relations as13

RT ln(φi) )

[( )

∫V∞

∂P ∂ni

-

T,V,nj

]

RT dV - RT ln Z (4) V

The PR EoS,7 with the mixing rules proposed by Wong and Sandler,6 is used as the thermodynamic model to evaluate the fugacity coefficient φ2. The PR equation can be expressed as follows:

P)

RT a + V - b V(V + b) + b(V - b)

(5)

a ) 0.457235(R2Tc2/Pc)R(Tr) b ) 0.077796(RTc/Pc) R(Tr)0.5 ) 1 + FPR(1 - Tr0.5) FPR ) 0.37646 + 1.54226ω - 0.26992ω2

(6)

For mixtures

P)

am RT + V - bm V(V + bm) + bm(V - bm)

(7)

In this equation, am and bm are the EoS constants to be calculated using defined mixing rules. The Wong-Sandler (WS) mixing rules for the PR EoS that are used in this work can be summarized as follows:5 N

Solubility Calculations The theory of solid solubility in a compressed gas is found in standard books,6,13 so a summary only is given in what follows. The fundamental equation of phase equilibria establishes that, at a given temperature and pressure, the fugacity of a component, for instance, the solid solute, in the gas phase must be equal to the fugacity of the same component in the solid phase. If subscript 2 stands for the solid component, then

(2)

bm )

(

N

a

)

∑i ∑j yiyj b - RT ij N

1-

yiai

AE∞(y)

(8)

∑i b RT - ΩRT i

(b - RTa )

ij

xaiaj(1 - k ) 1 ) [bi + bj] ij 2 RT

am ) bm

(∑ N

yi a i

i

bi

+

)

AE∞(y) Ω

(9)

(10)

4826

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005

In these equations, am and bm are the EoS constants with kij as an adjustable parameter, Ω ) 0.346 57 for the PR EoS, and AE∞(y) is calculated assuming that AE∞(y) ≈ AE0 (y) ≈ GE0 (y). For the excess Gibbs free energy, GE0 (y) is calculated using an appropriate liquidphase model. In this work, GE0 (y) has been calculated using the van Laar model that has been shown to perform well in high-pressure phase equilibrium calculations.14 For a binary mixture

bm ) y12(b - a/RT)1 + 2y1y2(b - a/RT)12 + y22(b - a/RT)2 y2a2 AE∞(y) y1a1 b1RT b2RT ΩRT

1-

(11)

(b - RTa )

12

xa1a2(1 - k ) 1 ) [b1 + b2] 12 2 RT

(

)

y2a2 AE∞(y) y1a1 + + b1RT b2RT Ω

a m ) bm

(12) (13)

The van Laar model for GE0 (x) is described by the following equation: N

GE0 RT

N

)

∑j

[

N

yi

yjAij

∑i yi 1 - y

i

1-

∑j

N

yi

∑j

]

2

yjAij N

yjAij + (1 - yi)

∑j yjAji

(14)

For a binary mixture, the model reduces to

(A12/RT)y1y2 GE0 ) RT y1(A12/A21) + y2

(15)

The expressions for the fugacity coefficient using the PR equation with then above-described WS mixing and combination rules can be found elsewhere.5 The problem is then reduced here to determine the parameters A12 and A21 in the van Laar model, the k12 parameter included in the combining rule for (b - a/RT)12, and the sublimation pressure Ps2 that appears in the solubility equation (3), using available high-pressure T-P-y data for gas-solid systems. The optimization procedure in the GA scheme requires an objective function that is defined here as

F)

100 ND

ycalc - yexp 2 2

∑ | i)1 ND

yexp 2

|

(16)

i

In this equation, ND is the number of points in the experimental data set and y2 is the solid solute concentration in the gas phase, both experimental (exp) and calculated (calc) values.

Genetic Algorithms, GAs Among the global stochastic optimization techniques, the evolutionary algorithms known as GAs have found many applications in several fields in science and engineering.15 In contrast to more traditional numerical techniques, which iteratively refine a single solution vector as they search for the optimum solution in a multidimensional landscape, GAs operate on entire populations of candidate solutions in parallel. In fact, the parallel nature of a GA stochastic search is one of the main strengths of the genetic approach. This parallel nature implies that GAs are much more likely to locate a global optimum than traditional techniques because they are much less likely to get stuck at local optima. Therefore, the problem of finding a local optimum is greatly minimized because GAs make hundreds, or even thousands, of initial guesses. The most significant differences between GAs and more traditional search and optimization methods are as follows: (i) GAs search a population of points in parallel, not a single point. (ii) GAs do not require derivative information or other auxiliary knowledge; only the objective function and the corresponding fitness levels influence the directions of search. (iii) GAs use probabilistic transition rules, not deterministic ones. (iv) GAs work on an encoding environment of the parameter set rather than the parameter set itself.16 Most applications of GAs correspond to optimization problems. These are problems in which a defined objective function must be fulfilled and the values of some variables and parameters included in this objective function must be searched for. GAs use biologically derived techniques such as inheritance, mutation, natural selection, and recombination to evolve toward better solutions so the objective function is best satisfied. At the beginning of the computation, a number of individuals represented by chromosomes are randomly created, forming a set known as the population. Each chromosome consists of a number of “zeros” and “ones” and represents a value of the parameter to be calculated. A set of individuals forms a matrix that represents the first/initial generation. After this, the objective function is evaluated for each of the individuals. If the optimization criteria are not met, the creation of a new generation starts. Individuals are selected according to their fitness (how close they fulfill the objective function) for the production of offspring (new individuals). All offspring are mutated with a certain probability, and the fitness of the offspring is then computed. The offspring are inserted into the population, replacing the parents and producing a new generation. This cycle is performed until the optimization criteria are reached. During the process of reproduction, mutation, and generation of new individuals, the so-called Gray encoding is incorporated. Gray codes, named after Frank Gray who patented their use for shaft encoders in 1953, constitute a form of representing or encoding numbers such that to go from one number to the next, always entails only one bit change (a characteristic that the usual binary coding does not have). There are many possible Gray codes, but of these many options, one particular form of the Gray code criterion has become popular, the so-called “binary reflected gray code” (BRGC), simply known as “the Gray code”. The word “binary” in the BRGC suggests the fact that the BRGC is a system of representing numbers using ones and zeros and also that it is similar to the binary system of

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005 4827 Table 1. Properties of the Solids Involved in This Studya substance naphthalene biphenyl phenanthrene anthracene pyrene carbon dioxide a

formula Tc (K) Pc (MPa) C10H8 C12H10 C14H10 C14H10 C16H10 CO2

748 773 869 873 936 304.2

4.05 3.38 2.90 2.90 2.61 7.383

ω

V (cm3/mol)

0.3020 0.4028 0.4707 0.4857 0.5074 0.2236

0.1100 0.1312 0.1528 0.1430 0.1591

All data are from Daubert et al.21

representing numbers using ones and zeros that computers regularly utilize. It is so similar that one can easily switch back and forth between regular bit patterns and BRGC bit patterns for any given integer. For further details on GAs and Gray encoding, the works of Davis,15 Pohlheim,16 Escudero,17 Obitko,18 Whitley,19 and Chakraborty and Janikow20 are recommended. The particular sequence of the GA mentioned above corresponds to a single population evolutionary algorithm that has been demonstrated to be powerful and to perform well on a broad class of problems. For other more complex cases, better results have been obtained by introducing many populations, called subpopulations. The multipopulation evolutionary algorithm models the evolution of a species in a way more similar to nature than the single population evolutionary algorithm. Although the problem of finding the sublimation pressure of solids from the solubility data seems to be complex, it is not so for the GA method. Thus, a single population evolutionary algorithm is used in this work, as detailed in a later section.

Table 2. Details on the Phase Equilibrium Data for the Eight Systems Considered in This Studya range of data system CO2 + ND T (K) P (MPa) y2 (×104)

reference

9 13 19 7

308 328 333 338

9-25 8-29 10-29 15-23

75-192 McHugh et al.22 13-538 52-980 247-790

8 8 8

308 318 322

10-44 15-44 15-47

104-154 McHugh et al.22 155-272 178-360

30

313

11-33

0.2-0.4

30 30

323 333

12-35 13-35

0.2-1.4 0.2-1.8

phenanthrene 26

313

11-35

4-23

30 30

323 333

11-34 11-35

3-31 2-40

32

313

10-34

0.5-3.6

32 32

323 333

12-35 11-34

0.4-5 0.2-6

naphthalene

biphenyl

anthracene

pyrene

Anitescu and Tavlarides23

Anitescu and Tavlarides23

Anitescu and Tavlarides23

a The temperature and pressure values have been rounded to the closest integer.

all of those conditions and have been used by several authors for phase equilibrium modeling,24-26 for stability analysis,27 and for thermodynamic consistency.28 Application of the GA

Cases Studied Five binary gas-solid systems containing supercritical carbon dioxide have been considered to show the applicability of GAs to the optimization problem of finding the sublimation pressure of solid substances. The systems studied were binary mixtures containing supercritical carbon dioxide with naphthalene, phenanthrene, anthracene, biphenyl, and pyrene. Table 1 shows the basic properties of the compounds involved in the study. In the table, Tc is the critical temperature, Pc is the critical pressure, ω is the acentric factor, and V is the solid molar volume. Values of the sublimation pressure of each of the solids were obtained from the equations provided in the DIPPR database.15 These values are used to compare the values calculated using the proposed method. The experimental phase equilibrium data taken from the literature are presented in Table 2 (16 isotherms and a total of 344 P-T-y data points). As seen in the table, the temperature and pressure ranges are narrow and go from 308 to 338 K and from 8 to 49 MPa, respectively. The solubility of the solid in the pressurized gas, however, covers a much wider range and goes from 2 × 10-5 to 1.8 × 10-2. This aspect is of special importance because the model used must be able to correlate the data in this wide range of solubility. The problem analyzed in this study is obtaining the optimum value of the sublimation pressure Ps2 using the available solubility data of solid solutes in a supercritical gas. In doing this, the interaction parameter k12 included in the combining rule for (b - a/RT)12 and the parameters A12 and A21 included in the van Laar model are also calculated. As stated before, the data used for the calculation must fulfill some basic requirements. In particular, the data selected for the present study fulfill

To explain in more detail the implementation of the GA, let us consider a set of nine P-T-y data for the system CO2 + naphthalene at 308 K and pressures ranging from 9 to 25 MPa, the first data set in Table 2. To start the GA, the parameters shown in Table 3 must be defined. After these numbers are defined, the steps, part of the GA, to evaluate the NV ) 4 parameters are as follows: (1) The initial population is created. This consists of NI ) 30 individuals, each represented by NV ) 4 chromosomes of L ) 20 genes each. Each block represents the possible values of the parameters k12, A12, A21, and Ps2. There are then 80 binary digits (0 or 1) randomly determined, so the total length is LT ) 80. This forms the chromosome matrix (MatChrom) containing NI × LT elements. (2) The elements of the matrix MatChrom (each of the NV chromosomes) are converted into real numbers in the established range (search space) for k12, A12, A21, and Ps2. This is done by using Gray’s scheme, and a new phenotype matrix MatPhen, containing NI × NV elements, is formed. (3) The objective function (which we want to be a minimum) is evaluated for each of the N ) 30 sets of values k12, A12, A21, and Ps2, and the value of the objective function is called the fitness. The evaluation of the objective function is done by calculating the solubility with eq 3 using the PR/WS/VL model for each value of the set of parameters. The calculated and experimental solubility values are introduced into the objective function, resulting in a vector of NI ) 30 elements, called FunObj. The elements are in this case the deviations between the experimental and calculated values of the solubility, defined by eq 16. Thus, it is

4828

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005

Table 3. Variables and Parameters for the GA Method Developed in This Work To Evaluate the Sublimation Pressure of Naphthalene at 308 K from Experimental Solubility Data variable

notation/method/formula

no. of experimental data no. of variables no. of generations length of a chromosome total length of an individual no. of individuals aptitude definition selection scheme crossover operator crossover probability mutation operator mutation probability search space (range of feasible solutions for the parameters k12, A12, A21, and Ps2)

ND NV Ngen L LT ) NV × L NI lineal Universal Stochastic Sampling multipoint Pcros binary Pmut SS1 for kij

9 4 30 20 80 20

SS2 for A12 SS3 for A21 SS4 for -log Ps2

[1, 8] [1, 2] [1, 5]

objective function

F)

100 ND

expected that the objective function be the lowest positive value. (4) The elements of the vector FunObj are arranged according to their decreasing fitness (the value of the objective function). The element whose value is the closest to zero is assigned an arbitrary value of 2.0. The element whose value is the highest of the vector FunObj is assigned an arbitrary value of 0.0. The other NI 2 ) 28 elements are linearly distributed between 0 and 2 according to their values (fitness). Thus, a new vector named FunRank is formed. (5) The vector FunRank and the MatChrom matrix determine the selection of the individuals for crossover, according to the defined crossover probability, Pcros. The operator used to select the elements for crossover is known as Universal Stochastic Sampling. This consists of a random sampling in which the individuals in the MatChrom matrix whose corresponding values in the FunRank vector are higher have a higher probability of being selected as parents and reproduce offspring. It could happen that some individuals are selected more than once for reproduction. The new arrangement of the individuals forms a new matrix named MatSel. (6) Reproduction of the elements of the MatSel matrix is performed by multipoint crossover. Crossover is done by interchanging some of the genes of one parent with some of the other parent. Crossover points are randomly defined, and genes are interchanged. After this is done, a new population of offspring is created, a population represented by the matrix MatRecomb. (7) Mutation is performed on the new population MatRecomb, according to the defined mutation probability (Pmut ) 0.035). To do this, each gene is assigned a random number that represents the probability of being mutated. If the random number is less than Pmut, the gene is selected for mutation. Mutation is done by changing the binary number (if it is 0, it is changed to 1, and if it is 1, it is changed to 0). The elements of the new matrix formed after mutation are randomly introduced into the chromosome matrix (MatChrom). (8) The new MatChrom represents a new generation of chromosomes (new possible values for the parameters k12, A12, A21, and Ps2), so steps 2-7 are repeated until all defined generations (Ngen ) 30) have been created.

∑|

|

ND

ycalc - yexp 2 2

i)1

yexp 2

value

0.8 0.035 [0.02, 0.20]

minimum positive value

i

(9) The matrix MatPhen formed after Ngen generations contains the individuals representing the optimum value of the parameters k12, A12, A21, and Ps2. The elements of the matrix MatPhen are converted so that all numbers are in the search space for each of the variables (SS1, SS2, SS3, and SS4). The objective function is evaluated for each of the 30 sets (k12, A12, A21, and Ps2]. The set of parameters that gives the minimum value of the objective function represents the solution to the problem. A schematic diagram that summarized the above steps of the proposed numerical method is presented in Figure 1. The diagram includes simplified blocks of each of the main steps. A more detailed explanation of how each step is actually done can be found in books and manuals such as those of Davis,15 Polheim,16 and Goldberg.29 In the above detailed sequence, summarized in Figure 1, the definitions of the parameters required to start the GA method (the number of variables or chromosomes, the number of generations, the length of a chromosome, the search space, the crossover probability, and the mutation probability) were defined according to literature information. Also, the aptitude definition, the selection scheme, the crossover operator, and the mutation operator were defined according to recommendations from the literature.15,16 According to Davis,15 these parameters must be carefully defined to avoid divergence or slownesss in reaching convergence. Of all of these values and definitions, the GA method has shown to be very sensitive to the search space, that is, the range of feasible solutions for the parameters to be calculated: k12, A12, A21, and Ps2. To accelerate the search and to avoid divergence of the method, one can define a large search space and a large value for the mutation probability Pmut (0.1 is a good start). After a solution is reached, the search space is narrowed and Pmut is changed to a smaller number. A final value of Pmut ) 0.05 is recommended for the problem solved here. This exploratory step is also necessary in other optimization methods when applied to these types of phase equilibrium problems.8,9

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005 4829

Figure 1. Flow diagram of the GA used in this work for the calculation of the sublimation pressure.

Results and Discussion The accuracy of the model and the method used to calculate the sublimation pressure is evaluated by considering the deviations between the experimental and calculated values of the solubility of the solid in the high-pressure gas. These deviations are expressed in relative form and absolute form as follows:

∆y2 % )

|∆y2 %| )

100 N 100 N

N

∑ i)1 N

∑ i)1

[ [

]

ycalc - yexp 2 2 yexp 2

(17) i

]

|ycalc - yexp 2 2 | y2exp

(18) i

If these two deviations are known, then a clear conclusion about the accuracy of the model used can be obtained.9 Table 4 shows the results for the five gas-solid systems considered in this study, using the PR/WS/VL model. Both the relative and absolute deviations defined in eqs 17 and 18 are shown in the table. Because the

solubility data were used to evaluate the sublimation pressure at each temperature, accurate modeling of the solubility is essential to obtain accurate values of the sublimation pressure. The results found using the proposed GA method show that the solubility is correlated with average deviations of -1.57% and absolute average deviations of 3.56%. Relative maximum deviations are -4.54% for phenanthrene + CO2, and absolute maximum deviations are 6.08% for naphthalene + CO2. These values are lower than other results presented in the literature for the same systems in the same range of temperature, as shown in Table 5. The literature values shown in Table 5 were taken from Escobedo-Alvarado et al.30 and from Valderrama and Alvarez.9 In both cases, three optimum temperature-independent interaction parameters were determined. Escobedo-Alvarado et al. used the PR EoS with the WS mixing rules, including the NRTL and UNIQUAC models for the excess Gibbs free energy (WS/ NRTL and WS/UNIQUAC). Valderrama and Alvarez used the modified PR EoS with different temperatureindependent mixing rules of Kwak and Mansoori31 (PR/ KM and PR/KM3). As seen in the table, the proposed method using GAs gives low deviations for the solubility,

4830

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005

Table 4. Parameters for the Mixing Rules and Deviations between Calculated and Experimental Solubility Values system CO2 +

ND

T (K)

P (MPa)

range y2 (×103)

k12

A12

A21

∆y %

|∆y| %

naphthalene

9 13 19 7 8 8 8 26 30 30 30 30 30 32 32 32

308 328 333 338 309 318 323 313 323 333 313 323 333 313 323 333

9-25 8-28 11-29 15-23 10-44 15-44 15-47 11-33 12-35 13-35 11-35 11-34 11-35 10-34 11-35 11-34

75-192 13-538 52-980 247-790 104-154 155-272 178-360 0.2-0.4 0.2-1.4 0.2-1.8 4.2-23.4 3.0-31 2.0-40 0.5-3.6 0.4-4.5 0.2-5.5

0.0768 0.0336 0.0949 0.0823 0.0650 0.0875 0.0942 0.2063 0.1852 0.1580 0.1581 0.1557 0.1380 0.2479 0.2082 0.1882

2.8569 3.3386 2.2897 2.5744 3.1149 2.5287 2.3596 5.5855 5.2429 5.4934 4.9749 4.5326 5.1154 6.9635 6.3920 6.0294

2.7653 0.6618 0.8385 0.8333 3.9735 2.3892 1.6576 0.1223 0.5624 2.0741 0.2825 0.7596 0.2586 2.0779 1.9664 0.4953

-0.60 -3.98 -3.49 -3.15 -3.15 0.37 1.06 -1.07 -0.67 -0.97 -0.15 -3.59 -4.54 -0.21 -1.44 -2.59

3.01 5.67 6.08 3.83 5.05 0.62 1.66 2.97 1.58 4.97 2.29 4.79 5.60 2.00 3.49 3.73

biphenyl anthracene phenanthrene pyrene

Table 5. Deviations between Literature and Calculated Values of the Sublimation Pressure, Ps2 a system CO2 + naphthalene

biphenyl anthracene phenanthrene pyrene

T (K)

Ps2(lit.) (bar)

Ps2(calc) (bar)

% ∆Ps

308 328 333 338 309 318 323 313 323 333 313 323 333 313 323 333

2.80 × 10-4 1.57 × 10-3 2.40 × 10-3 3.38 × 10-3 4.17 × 10-5 1.10 × 10-4 1.64 × 10-4 8.60 × 10-8 2.61 × 10-7 7.53 × 10-7 1.10 × 10-6 3.21 × 10-6 9.45 × 10-6 1.03 × 10-7 2.04 × 10-7 4.86 × 10-7

2.89 × 10-4 1.59 × 10-3 2.44 × 10-3 3.37 × 10-3 4.18 × 10-5 1.09 × 10-4 1.65 × 10-4 8.61 × 10-8 2.60 × 10-7 7.53 × 10-7 1.09 × 10-6 3.27 × 10-6 9.46 × 10-6 1.04 × 10-7 2.05 × 10-7 4.93 × 10-7

3.0 1.3 1.7 -0.3 0.2 -0.9 0.6 -0.4 0.0 -0.9 1.9 0.1 1.0 0.5 1.4 -0.4

a % ∆Ps is the percent deviation defined by eq 19. The values Ps2(lit.) were obtained from Daubert et al.21

a necessary characteristic of the method to obtain reliable values of the sublimation pressure. Table 6 shows relative deviations between the estimated and experimental sublimation pressures for the five solids considered in this study. The relative deviation is defined as

% ∆Ps )

Ps2(calc) - Ps2(lit.) Ps2(lit.)

(19)

Comparisons of calculated sublimation pressures are done only with the available data in the DIPPR database. As seen in the table, maximum deviations are on the order of 3% while global averages are less than 1% for both the relative and absolute deviations. Also, Figure 2 shows the sublimation pressure for the five solids considered in this study. These results, shown in Table 5 and Figure 2, indicate that the proposed method is reliable enough to estimate the sublimation pressure of pure solids using the solubility data of the solid in high-pressure gases. In addition to the results shown and discussed above, three points need further discussion: (i) the effect of the sublimation pressure of the solid on solubility; (ii) the effect of the volume on the calculated sublimation pressure of the pure solid; (iii) the solvent effect on the calculation of the sublimation pressure of the pure solid. (i) Effect of the sublimation pressure on the solid solubility. For this, calculation of the solubility of

Table 6. Effect of Different Variables and Properties on Sublimation Pressure and Solubility Calculations Effect of the Sublimation Pressure on the Solid Solubility Ps2(lit.) (×104 bar)

Ps2(calc) (×104 bar)

∆y %

|∆y| %

2.8

2.7 2.8 3.0

6.5 -0.6 5.0

19.5 3.0 5.8

Effect of the Solid Volume on the Calculated Sublimation Pressure Vs2

Ps2(calc) (×104 bar)

∆y %

|∆y| %

0.0943 0.1095 0.1249

2.99 2.89 3.09

1.7 -0.6 0.8

3.8 3.0 3.2

Solvent Effect on the Calculation of the Sublimation Pressure solute

Ps2 (×104 bar)

∆y %

|∆y| %

CO2 ethylene32 ethane33

2.89 2.81 2.84

-0.6 -0.1 1.5

3.0 5.5 4.4

naphthalene on compressed carbon dioxide at 308.15 K has been done using values of the sublimation pressure higher and lower than that reported in the literature (Ps2 ) 2.80 × 10-4 bar). Table 6 shows the results. As seen in the table, a reasonably accurate value of the sublimation pressure of the pure solid must be available if accurate values of the solubility are needed. The average deviation in solubilities is less than 1% using Ps2 ) 2.80 × 10-4 bar. Variations of (10% in Ps2 increase this average solubility deviation to 5% for +10% variation and 19% for -10% variation. (ii) Effect of the volume on the calculated sublimation pressure. For this, calculation of the sublimation pressure of naphthalene and of the solubility of naphthalene on compressed carbon dioxide at 308.15 K has been done using values of the solid volume higher and lower than that reported in the literature (Vs2 ) 0.1095). Table 6 shows the results. As observed in the table, the sublimation pressure is not much affected by the value of the solid volume, while this is between acceptable ranges of (10%. However, deviations in the solid solubility increase. (iii) Solvent effect on the calculation of the sublimation pressure. For this, calculation of the sublimation pressure of naphthalene at 308.15 K has been done using values of the solid solubility of naphthalene in three compressed gases: carbon dioxide, ethylene,32 and ethane.33 Table 6 shows the results. As observed in the

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005 4831

Figure 2. Sublimation pressure for the five solids considered in this study. In these figures, the squared dots (0) are literature data and lines are calculated values.

table, the sublimation pressure is not much affected by the solvent in which the solid is dissolved (less than 3% difference between the calculated sublimation pressure). This is an expected result that, however, needs to be demonstrated, as is done in Table 6. The results shown indicate that the proposed method can be considered solvent-independent, so it can be used to determine the solubility of a pure solid using data of the solubility of that solid in any compressed gas.

estimate the sublimation pressure of the solids for which no data are available (Table 4). (ii) The GA method has been shown to be a good tool to solve the optimization problem studied here, providing an accurate global optimum (Table 4). (iii) The PR/WS/VL model has been shown to be an accurate model for the solubility of the solid solute in the compressed gas phase. Deviations found are lower than those reported in the literature using other models (Table 4).

Conclusions Sublimation pressures of solids from high-pressure solubility data are calculated using GAs. On the basis of the Results and Discussion presented in this study, the following main conclusions are obtained: (i) The low deviations between the experimental and calculated values of the sublimation pressure show that the thermodynamic model PR/WS/VL is appropriate to

Acknowledgment The authors are thankful for the support of the National Commission for Scientific and Technological Research (CONICYT-Chile), through the research grant FONDECYT 1040285, the Direction of Research of the University of La Serena, La Serena, Chile, for permanent support through several research grants, and the

4832

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005

Center for Technological Information (CIT, La Serena, Chile), for computer and library support. Notation Symbols ai, aj, bi, bj ) constants for pure components in the equation of state am, bm ) constants for mixtures in the equation of state aij, bij ) interaction constants in the equation of state for a mixture F ) objective function FPR ) parameter in the R(Tr) function fs2 ) fugacity of the solute in the solid phase fg2 ) fugacity of the solute in the gas phase kij ) binary interaction parameters in the Peng-Robinson equation of state M ) molecular weight N ) total number of components NI ) number of individuals ND ) number of points in a data set NV ) number of variables Ngen ) total generation number Npop ) population size P ) pressure Pc ) critical pressure Ps2 ) sublimation pressure of the pure solid Pcros ) crossover probability Pmut ) mutation probability R ) universal ideal gas constant Tc ) critical temperature Tr ) reduced temperature T ) system temperature V ) molar volume Vs2 ) solid molar volume yi ) mole fraction of component i in the gas phase Abbreviations EoS ) equation of state GAs ) genetic algorithms PR ) Peng-Robinson WS ) Wong-Sandler VL ) van Laar model for Gex Greek Letters R(Tr) ) temperature function in an equation of state ∆ ) deviation ω ) acentric factor φ2 ) fugacity coefficient of the solid solute Super/subscripts aver ) average value exp ) experimental calc ) calculated s ) sublimation sol ) solid phase

Literature Cited (1) Bertucco, A.; Barolo, M.; Elvassore, N. Thermodynamic Consistency of Vapor-Liquid Equilibrium Data at High Pressure. AIChE J. 1997, 43 (2), 547-554. (2) Hua, J. Z.; Brennecke, J. F.; Stadtherr, M. A. Enhanced Interval Analysis for Phase Stability: Cubic Equation of State Models. Ind. Eng. Chem. Res. 1998, 37, 1519-1527. (3) Valderrama, J. O. The State of the Cubic Equations of State. Ind. Eng. Chem. Res. 2003, 42 (7), 1603-1618. (4) Wong, D. S.; Sandler, S. I. A Theoretically Correct Mixing Rule for Cubic Equations of State. AIChE J. 1992, 38, 671-680.

(5) Orbey, H.; Sandler, S. I. Modeling Vapor-Liquid Equilibria. Cubic Equations of State and Their Mixing Rules; Cambridge University Press: New York, 1998. (6) Prausnitz, J. M.; Lichtenthaler, R. N.; Gomes de Azevedo, E. Molecular Thermodynamics of Fluid-Phase Equilibria, 3rd ed.; Prentice Hall Inc.: Upper Saddle River, NJ, 1999. (7) Peng, D. Y.; Robinson, D. B. A New Two-Constant Equation of State. Ind. Eng. Chem. Fundam. 1976, 15 (1), 59-64. (8) Valderrama, J. O.; Alvarez, V. H. Vapor-Liquid Equilibrium in Mixtures Containing Supercritical CO2 using a New Modified Kwak-Mansoori Mixing Rule. AIChE J. 2004, 50 (2), 480-488. (9) Valderrama, J. O.; Alvarez, V. H. Temperature Independent Mixing Rules to Correlate the Solubility of Solids in Supercritical Carbon Dioxide. J. Supercrit. Fluids 2004, 32 (1-3), 37-46. (10) Rangaiah, G. P. Evaluation of Genetic Algorithms and Simulated Annealing for Phase Equilibrium and Stability Problems. Fluid Phase Equilib. 2001, 187-188, 83-109. (11) Li, H.; Yang, S. X. Modeling of Supercritical Fluid Extraction by Hybrid Peng-Robinson Equation of State and Genetic Algorithms. Biosyst. Eng. 2003, 86 (1), 17-25. (12) Vijande, J.; Pieiro, M. M.; Mosteiro, L.; Legido, J. L. Estimation of Carbonate-Alcohol Interaction Parameters for Nitta-Chao Group Contribution Model: Application of a Genetic Algorithm. Fluid Phase Equilib. 2003, 212 (1), 165-174. (13) Walas, S. M. Phase Equilibria in Chemical Engineering; Butterworth-Heinemann: Boston, 1985. (14) Valderrama, J. O.; Zavaleta, J. Generalized Binary Interaction Parameters in the Wong-Sandler Mixing Rules for Mixtures Containing Alkanols and Carbon Dioxide. Fluid Phase Equilib. 2004, submitted for publication. (15) Davis, L., Ed. Handbook of Genetic Algorithms; Van Nostrand Reinhold: New York, 1996. (16) Pohlheim, H. Genetic and Evolutionary Algorithm Toolbox for use with Matlab (GEATbx). http://www.systemtechnik. tu-ilmenau.de/∼pohlheim/GA_Toolbox/index.html, 1997. (17) Escudero, L. A. Algoritmos geneticos. http://www. tierradelazaro.com/mates/alggen.htm, Oct 26, 2004. (18) Obitko, M. Introduction to Genetic Algorithms, http:// cs.felk.cvut.cz/∼xobitko/ga, october 26, 2004. (19) Whitley, D. A Free Lunch Proof for Gray versus Binary Encodings. Proceedings of the Genetic and Evolutionary Computation Conference, GECCO99, Orlando, FL, July 13-17, 1999. (20) Chakraborty, U. K.; Janikow, C. Z. Binary and Gray Encoding in Univariate Marginal Distribution Algorithm, Genetic Algorithm, and Stochastic Hillclimbing. Proceedings of the 7th International Conference on Genetic Algorithms, East Lansing, MI, July 19-23, 1997. (21) Daubert, T. E.; Danner, R. P.; Sibul, H. M.; Stebbins, C. C. Physical and Thermodynamic Properties of Pure Chemicals. Data Compilation; Taylor & Francis: London, U.K., 1996. (22) McHugh, M.; Paulaitis, M. E. Solid Solubilities of Naphthalene and Biphenyl in Supercritical Carbon Dioxide. J. Chem. Eng. Data 1980, 25, 326-329. (23) Anitescu, G.; Tavlarides, L. L. Solubilities of Solids Supercritical FluidssI. New Quasistatic Experimental Method for Polycyclic Aromatic Hydrocarbons (PAHs) + Pure Fluids. J. Supercrit. Fluids 1997, 10, 175-189. (24) Zhong, M.; Han, B.; Ke, J.; Yan, H.; Peng, D. Y. A model for Correlating the Solubility of Solids in Supercritical CO2. Fluid Phase Equilib. 1998, 146, 93-102. (25) Ruckenstein, E.; Shulgin, I. Cubic Equation of State and Local Composition Mixing Rules: Correlations and Predictions. Application to the Solubility of Solids in Supercritical Solvents. Ind. Eng. Chem. Res. 2001, 40, 2544-2549. (26) Kurnik, R. T.; Reid, R. C. Solubility of Solid Mixtures in Supercritical Fluids. Fluid Phase Equilib. 1982, 8, 93-105. (27) Xu, G.; Scurto, A. M.; Castier, M.; Brennecke, J. F.; Stadtherr, M. A. Reliable Computation of High-Pressure SolidFluid Equilibrium. Ind. Eng. Chem. Res. 2000, 39, 1624-1636. (28) Valderrama, J. O.; Zavaleta, J. Advances on Thermodynamic Consistency of High Pressure Phase Equilibrium Data. XXI Interamerican Congress of Chemical Engineering, Lima-Peru´, April 24-27, 2005. (29) Goldberg, D. Genetic Algorithms in Searching, Optimization and Machine Learning; Addison-Wesley: Reading, MA, 1989.

Ind. Eng. Chem. Res., Vol. 44, No. 13, 2005 4833 (30) Escobedo-Alvarado, G. N.; Sandler, S. I.; Scurto, A. M. Modeling of Solid-Supercritical Fluid Phase Equilibria with a Cubic Equation of State-Gex model. J. Supercrit. Fluids 2001, 21, 123-134. (31) Kwak, T. Y.; Mansoori, G. A. Van der Waals Mixing Rules for Cubic Equations of State. Applications for Supercritical Fluid Extraction Modeling. Chem. Eng. Sci. 1986, 41 (5), 1303-1309. (32) Masuoka, H.; Yorizane, M. Solid Solubility of Naphthalene in Gaseous Ethylene at High Pressures. J. Chem. Eng. Jpn. 1982, 15 (1), 5-10.

(33) Johnston, K. P.; Ziger, D. H.; Eckert, C. A. Solubilities of Hydrocarbon Solids in Supercritical Fluids. The Augmented van der Waals Treatment. Ind. Eng. Chem. Fundam. 1982, 21, 191197.

Received for review February 7, 2005 Revised manuscript received April 16, 2005 Accepted April 19, 2005 IE0501529