Thermodynamic Property Modeling for Chemical ... - ACS Publications

Ind. Eng. Chem. Res. 2009, 48, 4619–4637

4619

Thermodynamic Property Modeling for Chemical Process and Product Engineering: Some Perspectives John P. O’Connell,*,† Rafiqul Gani,‡ Paul M. Mathias,§ Gerd Maurer,| James D. Olson,⊥ and Peter A. Crafts# Department of Chemical Engineering, UniVersity of Virginia, CharlottesVille, Virginia 22904-4741, CAPEC, Technical UniVersity of Denmark, DK-2800, Lyngby, Denmark, Fluor Corporation, 3 Polaris Way, Aliso Viejo, California 92698, Applied Thermodynamics, UniVersity of Kaiserslautern, D-67653, Kaiserslautern, Germany, Dow Chemical Company, South Charleston, West Virginia 25303, AstraZeneca Pharmaceutical Ltd., Macclesfield, Cheshire, SK10 2NA, U. K.

Thermodynamic properties have always played essential roles in the engineering of chemical products and in the processes that manufacture them. Further, contemporary and future chemical technologies depend more than ever on property model formulation and application. This work explores how properties are utilized in process and product engineering, including opportunities and constraints of current property models, the current status of data availability and needs, and the interplay of data and models. Several case studies are given to illustrate underlying concepts, strategies for development, and methods of application to some industrial systems. 1. Introduction Over the past century, properties have been important in industrial and engineering chemistry. In 1909, the first issue of the Journal of Industrial and Engineering Chemistry1 contained 6 editorials, 16 articles and notes, several book reviews and new books notices, trade and industrial notes, and official regulations and rulings. About half of these pieces used the word “properties” to identify standards and characteristics of products, while several book titles contained the word “properties”. Other articles, such as “A New Bomb Calorimeter”, by C. J. Emerson and “Report of the Committee on the Analysis of Phosphate Rock”, were also oriented toward describing the attributes of chemical substances in commercial and personal use at the time, ranging from dyes to castor oil to whiskey. As chemical technology has broadened and deepened, the uses of properties have become much more sophisticated. Process designs are now developed via computation based on accurate data and complex models to reveal the conditions needed to attain desired product content and quality, to optimize efficiencies for sustainability, and to suggest alternative molecular structures for novel applications in health, comfort, and defense. This work describes a view, and gives examples, of contemporary and future applications for properties in chemical process and product engineering. We seek to give some perspectives to both property model developers and users who might benefit from seeing a conceptual framework along with some specific cases, to enhance their efficiency in performing process and product simulation and simulation of complex chemical systems. * To whom correspondence should be addressed. Tel.: 1 (434) 9243428. Fax: 1 (434) 982-2658. E-mail: [email protected]. † University of Virginia. ‡ Technical University of Denmark. § Fluor Corporation. | University of Kaiserslautern. ⊥ Dow Chemical Company. # AstraZeneca Pharmaceutical Ltd.

We start with the nature of process and product engineering and of property models as well as their possible roles in process and product design. This is followed by a brief discussion of information sources for property descriptions and predictions. We then give a set of cases to provide current context and a sense of future developments. Of necessity, important situations have had to be omitted; our focus in this paper is on chemical and pharmaceutical systems. As a result, areas such as advanced biotechnology, nanotechnology, interfaces, polymers, near-critical systems, petroleum fractions, and ionic liquids have not been included. Most of these are of great practical importance and represent many challenges. However, the nature of modeling their properties is either an extension of techniques we give here or involve concepts and details that require more discussion than we can provide. In any case, we hope that the approaches described here will stimulate workers in such areas to find value in our suggested modeling framework and methodologies, such as careful data collection for maximum information content, shrewd model development for efficient application, and appropriate computational schemes for obtaining optimal results. 2. Property and Process Modeling in Chemical Engineering The principal goals of product and process engineering are to bring to commercial reality both innovations and enhancements for the myriad of items, substances, and systems that nature allows by manipulation at the molecular and macroscopic levels. Contemporary and future approaches build upon traditional analysis and design, while using greater understanding of natural laws and relationships plus advanced tools, especially modern instruments and computers. Incredible advances have been made in chemical technology during the last century. A comparison of article titles in a current Industrial and Engineering Chemistry Research issue with those of the early issues of the Journal of Industrial and Engineering

10.1021/ie801535a CCC: $40.75  2009 American Chemical Society Published on Web 04/16/2009

4620

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 1. Process and property model relationships.

Chemistry leads to an appreciation of how the issues and techniques have advanced. Several valuable perspectives on applied thermodynamics in process modeling are available.2-5 Another illustration of change is in the kinds of failures and systems investigated by Zudkevitch in his 1980 paper on “Forensic Thermodynamics”6 compared to those of today. Earlier errors about phase boundaries and incorrect heat effects are much less frequent for well-defined systems, though they are still found in many current and complex systems. These developments have come about because there is much fuller recognition about what product and process engineering, especially design, require in order to be accurate and reliable. It is important to recognize the enormous number and manifold varieties of systems and substances for which models of many different properties are now expected to be found. As of March 2009, there were over 45 000 000 organic and inorganic substances and nearly 61 000 000 chemical sequences in the CAS registry.7 Dealing with the infinity of mixtures that can be generated from these compounds is daunting, but classification and careful strategies, combined with computational power, have allowed significant progress in treating more and more complex situations in recent years. Gani and O’Connell8 provide a list of some of the common types of chemical systems for which property models may be found in the literature, in process simulators, and specialty software. In general, for pure component and mixture thermodynamic properties, the principal model forms are PVTx equations of state and models for excess Gibbs energy, though these are generally limited to fluids. Additionally standard state properties of formation of substances are needed for reaction equilibria. They lead to bulk properties such as volumes, enthalpies, and entropies for equipment sizing and energy analyses, as well as to component properties, such as chemical potentials and fugacities for phase and reaction equilibria. In our work, we distinguish process models and property models (see Figure 1). Process models are the sets of mass and energy balance equations, as well as imposed physical, chemical, and economic constraints, of a process situation. These contain two kinds of properties: measurables (y), such as temperature (T), pressure (p), composition (x), and mass flows (m) of the streams, as well as conceptuals (θ), such as enthalpy, entropy, fugacity or chemical potential, etc. To obtain quantitative equilibrium (no time dependence) or dynamic behavior (t varies), process model relations require values of both property types, either directly by empirical correlation of measured data, or by calculations with property models. Product engineering first requires specification of the desired behavior, commonly expressed in both properties and microscopic structure, and then establishment of the process steps to make the desired product. Evaluations of process behavior or product performance often use conservation of mass and/or energy; in such cases, the property variables are internal variables and are part of the set of constitutive equations used to find the system variables. On the other hand, when process states or product qualities need to be determined, balance relations are not usually involved, so

property variables become the unknowns and are found on a stand-alone basis. This means that product design problems are distinct from process design problems. In processing, the chemicals are known and the behavior is to be solved for. With products, the desired properties are known, but the chemical identities (molecular or atomic structure) or their mixture compositions are unknown. For products, property models are typically used in an iterative mannersthe generate-and-test paradigm. That is, for a generated molecular structure or mixture, the target properties are calculated or tested. If they are not within the desired range, another alternative is generated, and the calculation cycle is repeated until a generated set of structures meets the specifications. Note that the property models of a product design problem also become the constitutive equations of the process models for designing the process to manufacture the product. While process relations can be viewed as always true for the system under consideration, property relations are expected to have limited accuracy and reliability, thereby introducing uncertainties to the design. 2.1. Roles of Property Models in Process and Product Engineering. Gani and O’Connell8 and O’Connell and Neurock9 have described the roles of property models in computer aided process and product engineering. Property models play three distinctive roles, particularly in process and product design. There is the common serVice role, where a specified set of property values is provided when requested. Additionally, there is a serVice plus adVice role, where models provide information about feasibility, in addition to property values. Finally, there is the integration role where property models further, and directly, contribute to the technique of problem solution. Typically, properties play the serVice role in process simulation, the adVice role in process and product design, and the integration role in developing efficient and flexible integrated simulationdesign strategies. Though property models are generally recognized for the traditional serVice role, the other roles are actually more powerful. The adVice and integration roles often improve designs, widening or narrowing search spaces and increasing solution method efficiency. Simulation-based forward design approaches can become bogged down with complex models. In fact, in many design and/or simulation problems, the target conditions (which are usually dependent on properties) are known but are only used to verify/analyze such simulation results and as a basis for initiating a simulation for an alternative design. The “reverse approach” of searching with property models for the substances and conditions that take the known inputs to the desired outputs can be much more efficient than forward trial methods. The roles of property models need to be considered for determining strategies for their selection, use, and ultimate problem solution as described in the work of Kontogeorgis and Gani3 and Gani and O’Connell.8 2.2. Property Models: What and How? We complete our introduction by describing the nature of property model relations, suggest some strategies for model selection, list information sources for quantitative use of models, mention calculational aspects, and end with techniques for model development. This is followed by several examples to illustrate how these elements can be put together in some contemporary, and complex, systems. 2.2.1. Model Relations. Property models inevitably involve a mathematical relation from which properties can be computed when the system is fully defined. This also means a property model is a general algebraic relation containing parameters

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

whose values characterize the system conditions and components. Property relations may be linear or nonlinear, and the parameters may be constants or functions of variables such as temperature, pressure or density, and composition. As indicated below, there can be several possible sources of parameter values. The most common model formulation is a generalized expression, often in the form of corresponding states where the variables are scaled by with pure component critical properties and combining rules are used for mixtures, as determined from a few substances and mixtures, and then applied to other systems. The issue is whether the corresponding states similarities among substances and mixtures will hold for predicting the property for new systems. Thus, the acentric factor formulation was established for “normal” fluidssnonpolar and weakly polar substances whose deviations from noble gases are principally from the intermolecular forces of globularity and nonsphericity. One should not expect to accurately describe systems with large dipoles or, especially with association and solvation interactions such as hydrogen bonding, with simple rules. 2.2.2. Model Selection. The process or product situation determines which property model descriptions are needed. Consider specification of the volume of a tank to hold a liquid mixture. If the conditions will be for a single phase at a specified temperature and pressure, only a single density would be needed. However, if the temperature and composition of the phase varied, a model for calculating liquid density as a function of these variables is necessary. When multiple phases could occur in the system, models for fugacity or activity to obtain the equilibrium phase compositions are required to determine the number, relative amounts, and density of all possible phases of the mixed system. If heat is to be added or removed to manipulate the number of phases, values for mixture enthalpies are needed at different conditions of state (temperature, pressure, composition, and phase). Finally, if reactions occur in the tank, the calculations involve simultaneous chemical and physical equilibrium and require properties of formation for the species of the system. This apparently simple case demonstrates the possible range of demands for adequate descriptions of chemical systems. For the more complicated situations, model selection may not be straightforward. In all property roles, important decisions must be made about model form and content to balance effort with accuracy, reliability, and ease of prediction. Dependability, accessibility, generality, and effort usually need to be weighed carefully. Even when expectations for outcomes are clear, the optimal procedure is often not apparent. For example, if vapor pressure is needed, should data be taken and correlated, Antoine parameters from the literature be accepted, corresponding states estimation be used, group or atom contribution methods be implemented, or descriptors from molecular calculations be tried? Is the accuracy adequate for the service role in detailed process design? Is the generality suitable for the advising role when creating viable alternatives? Can the information content be sufficient for the integration role for the elements of a process simulation? Different answers to these questions lead to different paths as well as variations in design efficiency besides different outcomes. Often it seems that no new model is needed, and the only decision is to choose the best model from among those available. But we note that, even for experienced users, the multitude of model possibilities, and their complexities to deal with the variations of pure components and mixtures, can make model and data selection quite difficult. That means a property-model user with limited knowledge and experience faces real challenges to find a model for maximum effectiveness rather than

4621

going with familiarity, hearsay, or ease of accessibility. Our warning is that there can be unpleasant consequences of using an inappropriate property description (e.g., in either the model or the parameters) with wrong results. The design can have bottlenecking or be oversized, as well as give too simplistic or overly complex, process configurations.10,11 2.2.3. Model Information Sources. Property models are intended to give quantitative descriptions of nature. As a result, the information used must be reliable and appropriately accurate. Three principal sources for property or parameter values are used: • Retrieval from literature or computerized databases • Estimation via correlation, prediction, or computation • New experimental measurement These sources are given in increasing order of time and resource requirements. Below are some details of their uses. Literature and Databases. A part of the contemporary revolution and explosion of information availability is storage and computerized access to property data. Electronic databases are preferred because, in addition to easy searching, they can be more readily corrected and kept up to date. Currently available electronic property and phase equilibria databases include: Thermodynamics Research Center (TRC) now at NIST,12 NIST WebBook,13 AICHE Design Institute for Physical Properties (DIPPR),14 Physical Properties Data Service (PPDS),15 National Standard Reference Data Series (NSRDS),16 Dortmund Data Bank,17 DECHEMA Chemistry Data Series,18 and the Engineering Science Data Unit-Glasgow.19 In addition, there are many handbooks and other printed sources of data, as recently reviewed by Harvey.20 A key issue in retrieved data is data-quality review and information management. This issue has been explored particularly well by Frenkel et al. in their discussion of the THERMO_XML initiative.21 A great variety of errors of inconsistency, tabulation, and omission are known to occur, and users should always verify the reliability of vital data by comparisons among sources and by using fundamental variations with state variables as given by thermodynamic equations and derivatives. For example, the expected signs of property derivatives can often be checked for consistency with respect to temperature, pressure, and composition. Some journals publishing properties work are now beginning a cooperation with the Thermodynamics Research Center12 to screen manuscript submissions with new data for internal consistency, as well as agreement, duplication, and plagiarism related to the literature. Estimation. Engineers may identify the estimation source of information as prediction, calculation, correlation, data-extension, or, perhaps somewhat misleadingly, “modeling”. Estimation of values is required because, particularly in process-design engineering, all the needed data at all the states, and their variations with conditions, cannot possibly be measured. An extensive compilation and review of estimation methods is given by Poling, Prausnitz, and O’Connell.22 The “how to” of estimation is most obviously related to the level of empiricism found in the different approaches highlighted in Figure 2. At the empirical end of the scale, the fundamental behavior of the needed quantity is unknown or excessively complex. Obtaining the required continuous variables is by curve-fitting sparse experimental data using, for example, polynomials, log-log plots, analysis of variable statistical methods (ANOVA), and time-series analysis. This technique should be restricted to data interpolation, and the resulting

4622

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 2. Iterative steps of property model development.

equations with fitted parameters should not be extrapolated outside the space of measured variables. At the theoretical end of the scale lie the methods of computational chemistry which are calculations from first principles where it may be claimed that “no data are needed.” These include quantum mechanics and statistical mechanics. The applied statistical mechanics technique of molecular simulation, and its application to process and product design, is discussed below. In between a full empiricism/and basic theory is the approach of “enlightened empiricism”, which uses actual or adapted forms of rigorous equations or models from chemical theory. Here, carefully selected collections of quantities include parameters adjusted to fit data. Examples include the following: activity coefficients, compressibility factors, residual enthalpies, and entropies and fugacity coefficients. Relations include equations-of-state, group-contribution methods, and corresponding-states formulations, which may be predictive, or at least allow cautious extrapolations outside the range of the data regressed. Molecular simulation computes properties using statistical mechanical relationships for properties from assumed inter- and intramolecular potential functions, or force fields. Computational chemistry obtains chemical energies and reaction paths from approximate solutions to quantum relations such as the Schrodinger equation. In addition to potentially costing less time and money than real experiments, “data” values can be obtained for states that would be difficult, dangerous, or impossible to measure in the laboratory. For example, the solubility of oxygen in flammable solvents and the formation properties of highly toxic substances can be estimated from molecular simulation and quantum calculation. Molecular simulation can also validate theoretically based models for properties,23-27 provide detailed molecular structures for connections to properties, and confirm or debunk molecular speculations associated with process or product design.28 Our view is that, over the last twenty-five years, molecular simulation methods have mainly provided qualitative molecular insight and identify erroneous assumptions about molecular structure and behavior with particular force field models, while fewer have provided quantitative estimation of properties and phase equilibria by comparisons with experimental data. In addition to typical petrochemicals,25 molecular simulation methods have been applied to polymers, ordered materials, electronic materials,29 and to nanostructures.30 The American Institute of Chemical Engineers has a forum dedicated to Computational Molecular Science and Engineering (CoMSEF).31 Experiment. Finally, laboratory measurement is the slowest and most expensive route for getting information about product and process design. Measurements should be reserved, but strongly considered, when uncertain quantities sufficiently influence process condition and configuration outcomes to justify the time and expense of reliable determination of properties. Experiments can be used to answer immediate property or design questions, to provide the “raw material” from which estimation methods are developed, or to provide validation of theoretical

methods. The importance of, and trends in, experimental measurement have recently been reviewed by Gupta and Olson.32,33 2.2.4. Model Calculation. While there are relatively few calculational issues for property models, they can be significant. Algorithm efficiency is important only for extensive process simulations and optimizations. However, convergence can be problematic if models do not give consistent results, such as proper phase or partitioning, at conditions encountered away from final states. Additional convergence difficulties can occur if reliable derivatives are not obtained from a model. Regression to reliable parameter values may not occur if calculated values are too sensitive, or not sensitive enough, to parameter variation, or if some of the parameters are correlated so that uniqueness cannot be established. Thus, care about computational aspects also needs to be exercised in model formulation. Gani et al.34 have proposed a systematic property model analysis to identify the complete and consistent set of property model equations to be solved and classifications of the variables to be specified (model parameters, chemical system characteristics, and problem specifications) to solve for the unknown properties. 2.2.5. Model Development. When no, or only inadequate, models exist, the level of investment to establish a reliable properties model, from basis to operating code, demands an efficient development procedure. An essential element is appropriate specification of the problem type and the expected outcome. If one formulates the type precisely such as, “develop a model for the estimation of the average density of polymers for the pressures and temperatures encountered in extruding the product”, there is a clear delineation of model type, application range, and expected users. The task is fairly clear. On the other hand, if the statement is, “develop a model to predict the activity coefficients of liquid solutions”, either constraints need to be made on the types of systems and conditions for which the model will be applicable or the model must treat all types of systems under all conditions. The latter may be impossible to achieve, even with current knowledge and techniques. One development strategy for property models involves the two-part, iterative process shown in Figure 3. Problem specification, involving important decisions, is first, where problem type, system classification, information sources, and expertise of intended users are chosen. These must include accuracy requirements on the properties to be obtained and the amount of effort to be expended in application. Then, the model development is done via construction, solution, and verification for internal consistency. Gani et al.34 describe structuring models for rapid coding as well as for recognizing the many relation-

Figure 3. Iterative steps of property model development.

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4623

ships, the independent and dependent variables, and the parametrizations embedded in models for simple and complex systems. It has also become common to formulate models in terms of contributions from groups of atoms, such as UNIFAC,35 as well as with dielectric continuum solvation models with charge screening, such as COSMO-RS36,37 and COSMO-SAC.38 3. Contemporary Examples of Property Modeling in Process and Product Design The above introduction describes our perspectives on property modeling, for process and product engineering, especially for design. The following is a set of specific examples which illustrate aspects of property model implementation so that readers can appreciate current approaches to a variety of systems and complexities. The first four cases involve different chemical processes with complex chemistry, or speciation. They require careful analysis, somewhat elaborate development, and extensive data. These models are principally used in the service role for design. Their success comes from ultimately achieving a robust and reliable functional form for both the chemical speciation and physical effect parameters as found from an array of data. The focus of our discussion is the final model form and content; the process of model development usually took considerable time and persistence. The next set of three cases, drawn from the pharmaceutical industry, gives insights about the service and advice role of properties in solvent selection. Finally, two emerging methodologies, estimation of group parameters and molecular calculations, illustrate extensions of current approaches to predicting behavior in broader applications. 3.1. Models for Complex Solutions. The design of staged separation processes such as distillation, extraction or crystallization typically uses phase equilibrium as a basis for connecting to real systems. Thus, property models for component properties in vapors and liquids are heavily involved. Many important commercial chemicals strongly interact or dissociate to form detectable or hypothesized species, especially in aqueous and alcoholic solutions. Thus, even a system made from only two components can become a complex, multispecies mixture. These effects dramatically influence phase behavior and component partitioning. Describing such systems generally uses the same properties as for simpler solutions, but the treatment must deal with phase and reaction equilibria simultaneously. 3.1.1. Vapor-Liquid Equilibria for Formaldehyde with Water. A well-studied, and commercially important, system is the aqueous binary of formaldehyde plus water. Figure 4 shows a schematic isothermal vapor-liquid equilibria (VLE) diagram based on data referenced by Kuhnert et al.39 The liquid-vapor phase region is restricted to low formaldehyde concentrations since solids characterized as oligomers of poly(oxymethylene) glycols, formed from varying numbers of formaldehyde and water molecules, precipitate at higher concentrations. However, even at the low formaldehyde concentrations, complicated phase behavior arises, as shown in the insert of Figure 4. This is attributed to the formation of oligomers that remain in solution, preventing formaldehyde from volatilizing as well as complexing to form methylene gycol which can appear in the vapor phase. The selection of species for this case was confirmed by NMR spectroscopic measurements. Thermodynamic properties are usually insufficient to determine speciation, because there are too many different options, with too many parameters, that may correlate the data satisfactorily. Many of the models based on property data would be unreliable for extrapolation. The model was established some years ago when computational chemistry would not have been adequate to explore the most stable species.

Figure 4. Schematic phase equilibrium diagram of the binary system formaldehyde + water. V ) vapor; L ) liquid; S ) solid. (insert) Lower temperatures where an azeotrope exists at dilute formaldehyde.

Figure 5. Vapor and liquid species equilibrium model for the system formaldehyde-water.

However, even now, it is essential to validate calculations with appropriate measurements. Figure 5 shows a schematic of the vapor and liquid species equilibria with formaldehyde (FA), water (W), methylene glycol (MG ) HO(CH2O)H), and its oligomers (MGi ) HO(CH2O)iH, i > 1). The liquid phase, instead of being a binary, contains at least four species. The vapor phase is considered as a ternary mixture of the volatile species FA, W, and MG, since the oligomers should have vapor pressures that are extremely low. Model Relations and Selection. The thermodynamic problem is simultaneous vapor-liquid equilibrium for FA, W, and MG and chemical-reaction equilibrium for the formation of MG in both phases and the formation of MGi in the liquid phase. For the expected low-pressure distillation, the phase equilibrium relation chosen here assumes ideal gas vapor and nonideal liquid solution: psi xiγi ) pyi

(1)

where psi is the saturation pressure of the three volatilizing species i (i ) FA, W, and MG), with xi and yi being the mole fractions of volatilizing species i in the liquid phase and in the vapor phase, respectively, and γi being the activity coefficient of species i in the liquid phase. The total pressure is indicated by p.

4624

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 6. Equilibrium concentrations of methylene glycol (MG) and poly(oxymethylene) glycols MG2 and MG3 in aqueous solutions of formaldehyde at 338 and 368 K: (experiment) 0 Hahnenstein et al.,41 O Balashov et al.;42 (calculated) s Albert et al.43

The chemical reaction in the vapor phase is the formation of methylene glycol from formaldehyde and water: CH2O + H2O h HOCH2OH

(I)

The species concentrations in the ideal vapor phase are related to the equilibrium constant, a function only of temperature, through Kgas 1 (T) )

yMG p(0) yFAyW p

(2)

where p(0) ) 0.1 MPa is the standard state pressure. For reaction I, the equilibrium constant in the liquid phase KI can be expressed with K1gas and the saturation pressures of the pure components and is related to the liquid phase species activities: K1(T) ≡ Kgas 1 (T)

s psW pFA s pMG p(0)

)

xMG γMG xFAxW γFAγW

(3)

methylene glycol, as it does not exist as a pure substance. Thus, psMG(T) for MG was estimated or treated as an adjustable parameter when regressing VLE data. The parameters for the vapor reaction equilibrium constant, K1gas(T), were determined by the correlation of experimental gas phase density data.40 The parameters for the other Kj(T) were found by fitting experimental NMR data.41,42 Some of the UNIFAC parameters for the species were adjusted from those in the literature by regression of new VLE data. Model Calculation and Development. The calculations for the simultaneous phase and reaction equilibria were straightforward, and no new models needed to be developed. Results. Figure 6 shows the MG species concentrations in aqueous formaldehyde mixtures as calculated from the model43 in comparison with the experimental NMR data used to obtain the equilibrium constants.41,42 Since in the liquid phase more than 99% of the formaldehyde is converted to methylene glycol and poly(oxymethylene) glycols under the conditions shown in Figure 6, the mole fraction of (monomeric) formaldehyde is too small to show. Figure 7 shows a typical comparison between experimental data43-47 and correlation results44 for the vapor-

The formation of poly(oxymethylene) glycols involves other equilibria: HO(CH2O)n-1H + HOCH2OH h HO(CH2O)nH + H2O n > 1 (II) The equilibrium constant for the oligomer of degree n is Kn(T) )

xMGnxW

γMGnγW

xMGn-1xMG γMGn-1γMG

(n ) 2, 3, 4, 5,...)

(4)

The variation with T of psi (T) was the Antoine form, while that of various equilibrium constants was the usual parametrized form of ln Kj ) aj - bj /T

(5)

These were selected because the range of T was limited and no more complicated relations, i.e., added parameters, could be justified. The form for γi(T,x)was chosen to be the UNIFAC group contribution approach because it is predictive and the number of groups in this system is limited. Information Sources. The vapor pressures of formaldehyde and water are available in the literature, but none are found for

Figure 7. vapor-liquid partitioning of formaldehyde in the system formaldehyde-water at 363 and 413 K: (experiment) 3 Credali et al.,45 4 Kogan,46 ] Maurer,47 9 Albert et al.,43 b Albert et al.;44 (calculated) s Albert et al.44

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

liquid equilibrium of the binary system (formaldehyde + water) at 363 and 413 K. Extensions. Formaldehyde also reacts with alcohols. For example, it forms hemiformal and poly(oxymethylene) hemiformals with methanol. The above model was extended in a straightforward manner to the binary formaldehyde + methanol, the ternary formaldehyde + water + methanol, and to multicomponent systems containing trioxane and some reaction side products.39 The model correctly predicts that at low temperatures the presence of methanol results (at very low methanol concentrations) in a higher volatility of formaldehyde, whereas at higher methanol concentrations the volatility of formaldehyde is lowered. The thermodynamic model was also extended to describe caloric properties of these mixtures.48,49 This whole framework has been successfully applied by many companies to do basic engineering of processes involving aqueous solutions of formaldehyde. The important feature is that the species are established by, and their predicted amounts under some conditions are compared with, molecular measurements such as NMR and/ or UV-VIS spectroscopy. 3.1.2. VLE for CO2 and H2S in Aqueous Amine Solutions over Wide Ranges of Conditions. Sour gases, e.g., carbon dioxide and hydrogen sulfide, are commonly removed from natural or synthesis gas by “chemical” absorption in aqueous solutions of amines (such as, 2,2′-methyliminodiethanol ) N-methyldiethanolamine ) MDEA) or amine mixtures (e.g., MDEA + piperazine). While the competitive chemical absorption of CO2 and H2S is kinetically controlled, departures from equilibrium are the driving forces of such processes. Hence, the reliable design and optimization of the separation equipment requires knowing the chemical reaction thermodynamics and the vapor-liquid equilibria, along with information about the energy to vaporize/condense the mixtures. Gas absorption plants are run at ambient temperatures and a range of pressures as high as 4 MPa, whereas solvent regeneration is in a stripper (i.e., gas desorption) at elevated temperatures (over 390 K) and low pressures. Composition measurements show that the liquids leaving the absorption tower contain nearly no neutral amine and very small amounts of neutral sour gases, though there are large amounts of electrolyte reaction products such as protonated amines, bicarbonate, carbonate, and carbamate. In contrast, the liquids leaving the regeneration unit contain nearly no electrolytes, the sour gases have been stripped off, and the amines are mostly neutral. The latest references addressing this approach and earlier works are by Maurer and co-workers.50,51 A recent similar analysis related to postcombustion carbon dioxide capture in aqueous ammonia is given by Mathias et al.52 As a typical example, Figure 8 shows a speciation scheme for the solubility of CO2 in aqueous solutions of MDEA and piperazine (PIP). The vapor phase is considered to have only CO2 (C) and water (W), though solvent volatilization might need to be treated in full process design. As can be seen, the liquid phase is extremely complicated, containing more than a dozen species, neutral and ionic. Model Relations and Selection. This broad range of species and compositions requires a model that is able to describe phase behavior over a very wide range of temperature and pressure, as well as high loading of amine and CO2. The vapor-liquid equilibrium relation includes vapor-phase nonideality and the effect of pressure on the liquid phase. For water, the phase equilibrium relation is

(

psWφsW exp

)

VW(p - psW) aW ) yWpφ′′W RT

(6)

4625

Figure 8. VLE and chemical reactions in the CO2/MDEA/piperazine/H2O system.

while the extended Henry’s law standard state on the molality scale, (m) s (T; pW ) is used for carbon dioxide because it is supercritical HCW at most conditions of interest here. Consistent treatment of pressure and nonideality with eq 6 must be implemented.

(

s H(m) CW(T, pW) exp

)

V∞CW(p - psW) mCγ*C ) yCpφC′′ RT

(7)

Here the solute composition is given in molality (moles solute per kilogram of water), mC, and its activity coefficient (normalized according to Henry’s law on the molality scale) in the liquid phase is γ*C. The fugacity coefficients of saturated water vapor, water in the vapor mixture, and CO2 in the vapor are φHs 2O, φH′′2O, and φCO ′′ 2, respectively. These were predicted by the second virial equation of state since the pressures were not extremely high and the coefficients are generally more reliable for aqueous systems than cubic equations of state based on corresponding states. The Poynting factors for liquid phase pressure effects use the pure liquid molar volume for water, VW, and the partial molar ∞ volume at infinite dilution for CO2 in water, VCW , since these are good estimates and the effects on them of composition and pressure can be ignored. Chemical reactions dominate the liquid phase properties as shown in Figure 8. As in eqs 2-4, chemical equilibrium for reaction Rk is expressed using activities: KRk(T) )

∏a

νi,Rk i

(7a)

i

where the activity of a solute species i (i.e., all species except water) is the product of its stoichiometric molality and its activity coefficient appropriate for the Henry’s Law standard state, denoted with *: ai ) miγ*i

(8)

For water, the activity, aw, is calculated via integration of the Gibbs-Duhem equation using the activities of all the solutes. The temperature-dependence of the chemical equilibrium constant for the autoprotolysis of water (KR1) was adopted from the work of Edwards et al.,53 while KR2 and KR3 for the formation and dissociation of bicarbonate were taken from the work of Patterson et al.54,55 For the protonation of methyldiethanolamine,

4626

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 9. Partial pressure of carbon dioxide (left diagram) and total pressure (right diagram) above liquid mixtures of (CO2 + MDEA + H2O), m j MDEA ≈ 2 mol · kg-1: (experiment)60 4 313, O 353, 0 393 K; (experiment)61 2 313, [ 333, * 373, 9 393, × 413 K; (calculation) s correlation from all data, - - - from only high pressure data.60

Figure 10. Carbon dioxide partial pressure (left) and total pressure (right) above liquid mixtures of (CO2 + MDEA + H2O), m j MDEA ≈ 8 mol · kg-1: (experiment)59 4 313,O 353, 0 393 K; (experiment)60,61 2 313.7, b 354.4, 9 395 K; (calculations) s correlation from all data, - - - from only high pressure data.60

KR4 was derived from electrochemical investigations by Pe´rezSalado Kamps and Maurer.56 The corresponding constants for the protonation and diprotonation of piperazine (KR5 and KR6) were taken from the work of Hetzer et al.,57 while those for formation of piperazine carbamate, piperazine dicarbamate, and protonated piperazine carbamate (KR7, KR8, and KR9) were determined from the results of NMR-spectroscopic investigations by Ermatchkov et al.58 These parametrizations were chosen because the correlations were carefully done with precise data from spectroscopy and other methods, which were directly related to the reactions and species, rather than buried in combined phase and reaction equilibria. Activity coefficients of both molecular and ionic solute species were calculated from a modification of Pitzer′s equation for the excess Gibbs energy of aqueous electrolyte solutions.59 This model has been shown to accurately describe dilute and concentrated aqueous electrolytes. Information Sources. The scheme shown in Figure 8 involves extensive information. Parameters for the equilibrium constants were provided from the literature cited above. The Henry’s law constant is that for unreacted carbon dioxide in s , VW, and the second virial coefwater, which along with pW ∞ were ficients, which, along with the estimation method for VCW taken from in the literature. The Pitzer model for Gibbs excess energy requires binary and ternary parameters to describe the interactions between solute species from low gas loadings (i.e., at low partial pressures of carbon dioxide) to high gas loadings (i.e., at high partial pressures of carbon dioxide). The only source for obtaining reliably the most important parameter values is experimental data of the solubility of CO2 in aqueous solutions of MDEA and of piperazine at low, as well as at high, CO2 partial pressures. Such investigations were performed with two different types of experimental equipment.50,51,60 Results. The solid lines in Figure 9 and 10 show the excellent comparisons between experimental60,61 and model results60 for CO2 solubility in aqueous solutions at low (Figure 9) and high (Figure 10) MDEA molalities. In addition, the broken lines in both figures show predictions at low partial pressures of carbon dioxide when all interaction parameters were estimated by using only high-pressure gas solubility data. These predictions agree well with low pressure experimental data at low and moderate j MDEA values) but are less amine concentrations (high m j C/m accurate at higher amine concentrations where parameters characterizing important interactions between molecular MDEA and other solute species were not included, as they cannot be determined from high pressure gas solubility data.

This example illustrates the very great range of data that must be assembled in order to reliably develop a full model for a truly complex system. Note that every aspect was compared with data and sensitivity to parameters was tested. Compared to the first case, fewer molecular measurements were available, though the speciation was generally known. It is possible that computational chemistry methods could apply here, but this work was done before they were ready for implementation. 3.1.3. LLE Extraction of Carboxylic Acids from Aqueous Solutions. Many carboxylic acid products are produced by fermentation. Product recovery from the dilute aqueous solutions is achieved by reaction with a hydrophobic component, e.g. tri-n-octylamine (TnOA), with the resulting complexes extracted into an organic phase. Subsequently, the carboxylic acid must be separated and the auxiliary component regenerated. The design of such extraction and recovery processes requires a thermodynamic model for liquid-liquid equilibria that accounts for electrolytes in both aqueous and organic liquid solutions. Our example is for citric acid partitioning between water and different organic solvents in the presence of TnOA, including the effect of salts on the TnOA partitioning for recovery and regeneration. The latest references summarizing this approach are by Maurer and co-workers.62,63 At low aqueous-phase molalities of citric acid and fixed TnOA, the ratio of organic to aqueous phase concentrations of citric acid increases with increasing acid concentration. It then passes through a maximum and decreases at the highest concentrations. This behavior is attributed to two competing effects. At low acid concentrations in the aqueous phase, the dissociation equilibrium for citric acid is shifted to its ionic species. As only neutral acid molecules can be extracted into the organic solvent, the partition coefficient increases when the amount of dissolved neutral citric acid increases. The decrease of the partition coefficient observed at high acid concentrations results from complete complexation of the TnOA, so additional citric acid cannot be bound, and the acid remains in the aqueous phase. In such processes the aqueous phase may also contain strong electrolytes. While most strong electrolytes reduce the solubility of an organic compound in an aqueous phase, i.e., “salt-out” the organic compound, the presence of salt actually reduces the partitioning of a carboxylic acid to the organic phase when a complexing agent, such as TnOA, is present. Thus, Figure 11 shows the influence of NaNO3 on the partition coefficient of

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4627

(aq),(0) ) 0.05 volumes of aqueous and organic feed solutions (m ˜ NaNO 3 (org),(0) (aq) -1 -1 ˜ TnOA ) 1.24 mol kg ) giving m ˜ Cit.tot ) 0.02 mol mol kg ; m -1

(m) (org) (aq) Figure 11. Partition coefficient of citric acid (PCit ) m ˜ Cit /m ˜ Cit ) in the aqueous/organic two-phase system (citric acid + water + MIBK + TnOA + NaNO3) at 298.15 K for equal volumes of the aqueous and the organic (org),(0) feed solutions at constant TnOA concentration (m ˜ TnOA ) 1.24 molal) in the organic feed and several salt concentrations in the aqueous feed (aq),(0) solution:63 (experiment) s m ) 0, ] ) 0.01, O ) 0.05, 0 ) 0.1 ˜ NaNO 3 molal; - b - prediction.

Figure 12. Speciation in liquid-liquid phase equilibria for the system (citric acid + water + toluene + TnOA).

citric acid in the organic phase relative to the aqueous phase for the system (water + methyl isobutyl ketone (MIBK) + TnOA).64 Model Relations, Selection, and Information Sources. For this complex system, the fugacities of the coexisting liquid phases are treated with the aqueous phase as an electrolyte solution as in eqs 6 and 8, while the organic phase has only neutral species as in eq 1. Pitzer’s excess Gibbs energy model59 is used in both phases. It has terms associated with electrostatics and ions in the aqueous phase, while a power-law equation is used for the organic phase. The organic-phase complexes of citric acid and TnOA are in chemical reaction equilibrium; most contain water, as verified by IR-spectroscopy. The stoichiometry of the complexes depends on the organic solvent and can be complicated. For example, two complexes (citric acid:TnOA:water ) 2:3:2 and 1:1:1, respectively) were found for toluene, whereas four complexes (1:0:3, 1:2:3, 1:1:3, and 2:1:6) were in methyl isobutyl ketone (MIBK). The equilibrium constants for the various reactions were expressed as in eqs 4, 5, 7, and 8. Figure 12 shows the resulting speciation for liquid-liquid equilibria with toluene. Results. Comparisons between predictions and experimental data are also shown in Figure 11; the agreement is within experimental uncertainty. When NaNO3 is added to equal

kg in the equilibrated aqueous solution, the partition coefficient of citric acid to the organic phase is about 0.4 compared to about 40 in the salt-free system. When NaCl is the salt under the same conditions, the value is about 10. Extensions. A number of applications have been made by adding other reactive components to the scheme of Figure 12. The partitioning of inorganic acids (HCl, HNO3, and H2SO4) in TnOA-containing two-phase systems with toluene65 and MIBK64 and chloride partitioning in the system (citric acid + water + MIBK + TnOA + NaCl) have been studied. The addition of the salt of the carboxylic acid was predicted not to affect the partitioning of the acid, as there is no competition of different acids for the amine; this is found to be true for monocarboxylic acids like acetic acid, though not for acids with more than one carboxylic group. The model predicts all of the above behavior quantitatively, verifying the thermodynamic framework. Again, comprehensive use of, and comparisons with, many different data, along with careful speciation, allows quantitative description of many variations of these complex solutions. The results can be used for testing, for exploring options in process synthesis, and for process optimization. 3.1.4. Vapor-Liquid Equilibria for Oleum. Oleum, also known as “fuming” sulfuric acid, consists of SO3 dissolved in 100% H2SO4. Thus, for example, 20 mol % oleum consists of 20% SO3 and 80% H2SO4 by moles. The modeling objective is to provide VLE and heat of vaporization information for the system components in a case where corrosion and toxicity make extensive measurements extremely challenging. There are known to be many different species complexes in the liquid phase of oleum in addition to the H2SO4 and SO3 components. However, the dominant complex is H2S2O7,66 which is formed from 1 mol each of H2SO4 and SO3 by the complexation reaction proposed by Nilges and Schrage67 and by Mathias et al.:68 H2SO4+SO3 h H2S2O7

(III)

No ions are considered to exist in the solution. Model Relations, Selection, and Information Source. The phase equilibrium is obtained with liquid fugacity expressions of the form of eq 1 and vapor fugacity expressions of the form of eqs 6 and 7. The activity coefficients were obtained from the NRTL excess Gibbs energy model69 since the nonideality is strong. For this system where SO3 is the principal volatile component, the effect on VLE of parameters for the H2SO4-H2S2O7 pair is very small and is therefore ignored. However, parameters for the SO3-H2SO4 and SO3-H2S2O7 pairs must be valid over a range of temperature both for VLE and for calorimetric properties. Since the pressure is elevated, there is vapor-phase nonideality that was described by the Redlich-Kwong equation69 with standard parameters based on critical properties. This is not expected to be rigorous, but the conditions are such that the limited composition dependence of the vapor nonideality minimally affects the predicted behavior. The reaction equilibria of the system were defined by the equilibrium constant for reaction III plus the nonidealities of the three species in their liquid mixture: Koleum(T) )

xH2S2O7γH2S2O7 xH2SO4xSO3γH2SO4γSO3

(9)

4628

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

This is another example where success is obtained with an appropriate conceptual model and multiproperty fitting of quality data, though no molecular measurements were involved. In particular, calorimetric data proved quite valuable in validating a proposed “chemistry model”, given by the “signature” in the enthalpy of vaporization vs concentration curve. This case shows how a speciation can be validated by careful use of both phase equilibrium and calorimetric data, giving a better chance of finding the most appropriate description. 4. Solvent Selection for Pharmaceutical Production Figure 13. Comparisons of experimental data71 and model calculations (s) for vapor pressures of oleum mixtures.

Figure 14. Comparison between experimental data70 and model calculations for enthalpy of vaporization of SO3 from oleum mixtures of various concentrations at 30 °C.

where the temperature dependence of Koleum(T) is that of eq 5 plus polynomials in T. Miles et al.70 used two methods to measure the enthalpy of vaporization: (1) evaporation of SO3 from oleum under reduced pressure and (2) heat of solution (three sets of data) of addition of SO3 to oleum. The enthalpy of mixing used for matching the calorimetry data is obtained from the temperature dependence of Koleum(T) and of GE via the Gibbs-Helmholtz relation

(

∂ln Koleum ∂T

[

)

)-

]

)-

∂(GE /RT) ∂T

Px

Px

∆H0 RT2

(10)

hE RT2

(11)

Results. Comparisons between VLE data71 and the model, as shown in Figure 13, demonstrate that the model provides essentially quantitative agreement with all of these data. Figure 14 shows experimental data70 and model calculations for the enthalpy of vaporization of SO3 from oleum mixtures of various concentrations at 30 °C. The “enthalpy of vaporization” here is the negative of the enthalpy change that occurs at 30 °C when 1 kg of gaseous SO3 is dissolved in a large quantity of an oleum mixture of the given concentration. The sigmoidal shape of the enthalpy of vaporization curve is a clear “fingerprint” that strong chemical reactions are involved. At low SO3 concentrations, essentially all the added SO3 combines with H2SO4 to form H2S2O7. Hence, the total enthalpy of vaporization is approximately equal to the enthalpy of vaporization of pure SO3 (≈ 540 kJ kg-1) plus the heat of reaction (≈ 210 kJ kg-1). At about 50 mol % SO3 (45 wt % SO3), the amount of free H2SO4 has substantially decreased, so reaction does not occur and the enthalpy of vaporization rapidly decreases toward the enthalpy of vaporization of pure SO3.

We give here three examples of the service and advice role, where more effective solvents were selected for solution separation by using group contribution and other molecular thermodynamic methods. The pharmaceutical industry is increasingly challenged to develop more efficient and environmentally friendly processes for API manufacture, in a modern development environment that is characterized by high attrition. On average, only 1 in 10 new drug candidates survives through clinical trials to enter the market. In early process development, there are usually few physical property data available. Further, it is not economical to collect large quantities of thermophysical data since the systems frequently change as new chemistry is explored and separation processes are adapted. Thus, pharmaceutical process design needs predictive thermodynamic models to reduce the experimental search space and focus laboratory effort in the areas of greatest potential success. Modern drugs are functionally complex72 and often fall beyond the capabilities of traditional predictive models like UNIFAC. The reliable prediction of crystal structure and solid-state properties is computationally demanding, and still years away from mainstream application.73 Typically, the industry deals with complex chemistry,74 phase equilibria involving organic salts, and aqueous electrolytes. These factors make it clear that pharmaceutical systems challenge the capability of modern property tools. The following examples demonstrate how current methods can still usefully applied, as long as the problems are broken down into tractable sub-parts and appropriate modeling tools are used at each step. 4.1. Anisole Removal during Washing and Drying Operations. In this problem, an environmentally friendly solvent is needed to wash and dry a crystalline pharmaceutical intermediate. The solvent must efficiently remove anisole residues from upstream chlorination and coupling reactions, where the product is precipitated as an organic HCl salt and separated by pressure filtration to yield a 40% w/w anisole wet cake. Due to anisole’s low volatility, removing the residual anisole via inert gas drying is very slow and not commercially viable. Washing first with a more volatile solvent can increase the drying rate, and MTBE was used in early process development. However, at production scale, MTBE would possibly have required specific VOC abatement equipment, so a search for a better wash solvent was initiated. The final wash solvent was to be environmentally preferred, fully miscible with anisole, and promote vaporization of residual anisole. Model Relations, Selection, and Information Sources. To solve the above solvent selection problem, the solvent needs are translated into target properties. Table 1 lists an appropriate set of pure component and mixture property values. The boiling point and melting point indicate the liquid range. The key property to characterize the ease of anisole removal is the partial pressure, yip. As computed via eq 1, this quantity is increased by selecting a wash solvent that gives anisole activity coefficients as large as possible. Since the anisole will be dilute after washing, the relevant activity coefficient for screening is its limiting value at infinite

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4629

Table 1. Selection Criteria for Wash Solvent with Anisole property

minimum

maximum

comments reasonable volatility for drying and allows recovery by simple condensation liquid at ambient temperatures see eq 10; lower values mean lower residual anisole after drying should be less polar than IPA (δp ∼ 8) should have lower hydrogen bonding potential than ethyl acetate (δh ∼ 7) completely miscible with anisole CH3, CH2, CH, C, OH, CH3CO, CH2CO, CHO, CH3COO, CH2COO, HCOO, CH3O, CH2O, CH-O, COOH, CH2Cl, CHCl, CCl, CH2Cl2, CHCl2, CCl2, CCl3, CF3, CF2, CF, COO, CCl2F, HCCl2F, HCClF, CClF2, HCClF2, CClF3, CCl2F2, F

Tboil, °C

50

100

Tmelt, °C Sp

none 0.0

0.0 0.3

Hansen δp, MPa1/2 Hansen δh, MPa1/2

0.0 0.0

8.0 7.0

miscibility permissible solvent types and functional groups

alcohols, ketones, aldehydes, acids, phenols, esters, ethers, Cl, and F

Table 2. A Selection of Wash Solvent Alternatives solvent 2,2-dimethylpentane heptane MTBE CH2(OC2H5)2 (ethylal)

δh, MPa1/2 δp, MPa1/2 Tmelt,°C Tboil,°C 0 0 5 6.9

0 0 3.4 4.5

-124 -91 -109 -66.5

79 98 55.2 88

Sp 0.081 0.088 7.3 13.8

omers of chiral molecules. For esters, it is sometimes possible to use a catalytic lipase enzyme in an aqueous alcohol mixture to selectively dissociate the undesirable enantiomer into its acid via the reaction:

∞ dilution, γA,S . Thus, the solvent power, Sp, is used to rank prospective wash solvents for anisole volatilization

Sp )

1 MwA ∞ Mw γΑ,S S

(12)

The miscibility criterion is explicit. In addition to Sp, the Hansen solubility parameters, δh and δp,75,76 were used to reduce the solvent search space since the feasible solvent candidates are to have low affinity for ionic solutes; i.e., small terms for hydrogen bonding and polarity. Finally, the functional groups were limited to those with good environmental profiles. The ProCAMD solvent search software77 was used to find solvent candidates matching the target properties listed in Table 1. The pure component target properties were first estimated from generated molecular structural information via the Constantinou-Gani or Marerro-Gani78 method. Those structures (molecules) satisfying the target values were then examined for Sp and miscibility with the UNIFAC-LLE model.79 Results. A set of solvent candidates was identified among the 26 837 molecular structures generated by excluding 13 698 substances by the δh parameter; 3533 by δp; 2925 by Tmelt; 6484 by Tboil; 62 by Sp; and 7 by miscibility. That left 47 acceptable candidates as possible solvents; Table 2 shows representative properties of four candidates. Of these, heptane was selected as the wash solvent for its Sp value, commercial availability, volatility, and environmental profile (fugitive releases to the atmosphere). Figures 15 and 16 show activity coefficients for anisole with the original MTBE wash solvent and the heptane replacement as predicted by the ∞ in heptane is about a UNIFAC method.79 The value of ln γA,S factor of 3 greater than that for MTBE. In pilot plant trials with heptane, the drying step was considered rapid at about 7 h, and residual anisole levels were significantly reduced. Perhaps most importantly, the VOC emissions were reduced from 0.15 mol fraction for MTBE to 0.02 mol fraction for heptane. This example of solvent substitution demonstrates how group contribution methods may be applied in reverse and thus narrow the size of a search space and minimize time-consuming laboratory experiments. 4.2. Solvent Selection for an Enantiomeric Pharmaceutical. The synthesis of drugs often results in intermediates containing racemic mixtures of left- and right-handed enanti-

After dissociation, the chirally resolved ester is easily separated from the acid and alcohol using pH-buffered liquid-liquid extraction. A final crystallization yields the desired product. Finding the optimal conditions for the dissociation in a pH-buffered liquid extraction depends on estimating accurate values of the acid and ester dissociation equilibrium constants

Figure 15. Binary solution activity coefficients for MTBE and anisole predicted by UNIFAC.79

Figure 16. Binary solution activity coefficients for heptane and anisole predicted by UNIFAC.79

4630

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Table 3. Estimated pKa functional group

pKa’s for ester

pKa’s for acid

pyrazine N1 pyrazine N2 carboxylic H

7.8 3.4

7.8 3.9 2.9

or pKa, while selection of optimal solvent(s) is needed for both partitioning and for crystallization. Model Relations, Selection, and Information Sources. Several properties must be estimated for the process steps: the pKa values for the ester and acid as well as LLE for extraction and SLE for crystallization of the product with the proposed solvents. A review of pKa predicting methods is presented in ref 80. Note that there are two pKa’s for the ester and three for the acid. Software from ACD Laboratories81 was used to predict the pKa’s for these organic molecules with the values given in Table 3. The OH groups on the reaction product alcohol and on the t-butanol are relatively stable and do not need to be considered. Using straightforward calculations, the distributions of the species can be found as functions of pH. These are shown for the ester and acid in Figures 17 and 18. Above pH 10, the ester is present in the uncharged state, while the acid is fully deprotonated with a charge of negative one. Under these conditions, the partitioning gives maximum separation efficiency, since charged species prefer the aqueous phase, while the neutral ester prefers the organic phase. The pH is adjusted by adding a bicarbonate salt, since this acts as a sufficiently strong inorganic base for buffering during extraction but is not strong enough to hydrolyze the ester. After a water wash to remove the bicarbonate and subsequent decanting, the purified ester product is in a stream rich in t-butanol and saturated with approximately 30% w/w water. However, direct crystallization of the product by cooling this solution gave a poor yield; the product solubility in aqueous t-butanol was too high. Thus, a solvent was sought where the product could be extracted for crystallization. Results. Experimental screening of product solubility suggested toluene would be a good crystallization solvent, but the yield from the actual process solution turned out to be poor. The reason is clear from the ternary LLE diagram of Figure 19, predicted with the original UNIFAC LLE parameters.79 It

Figure 17. Species distributions for racemic ester dissociation as functions of pH.

Figure 18. Species distributions for racemic acid formation as functions of pH.

Figure 19. Liquid phase compositions for aqueous t-butanol with toluene at T ) 25 °C, 1 atm from UNIFAC-LLE.79

shows tie lines and a binodal curve with t-butanol favoring the toluene-rich organic phase rather than the aqueous phase. Alternatives to water, while keeping toluene, were sought using the search criteria of Table 4, similar to the process above. In this case, the normal melting point and boiling point were estimated from group-contribution methods for pure component properties; for the liquid density, the Rackett equation22 was used, while for the selectivity of the product for the organic phase and miscibility calculations, the UNIFAC-LLE model79 with the associated parameters was used. The search led to three substances expected to be commercially available, with the properties shown in Table 5. With the best solvent, 1,3-propylene glycol, the LLE phase diagram appeared as in Figure 20. The tie lines between the toluenerich and glycol-rich phases show the desired selectivity for t-butanol, but the two-liquid region extended too little toward the t-butanol apex. It was concluded that toluene would not provide a commercially viable batch extraction process. Further solubility screening identified cyclohexane as a potential crystallization solvent, so it was examined as an extraction solvent, by modifying the polar solvent search criteria given in Table 4 to include selectivity with cyclohexane and an updated density limit of 0.8 g/mL. Again, the propylene glycols were identified as the topranking candidates. The phase diagram of Figure 21 for 1,3propylene glycol with cyclohexane and t-butanol shows fully desirable characteristics. The two-liquid region extends to a 50: 50 volume ratio of cyclohexane to t-butanol, and the t-butanol

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4631

Table 4. Selection Criteria for Polar Solvent for Extraction of Pharmaceutical Ester property

min

Tboil, °C Tmelt, °C density, g cm-3

60

selectivity at feed composition miscibility permissible solvent types and functional groups

5.0

max 0.0

0.9

alcohols, ketones, aldehydes, esters, ethers

comments relatively low volatility to minimize VOC emissions be liquid at ambient temperature for ease of handling more dense than toluene (0.86 g cc-1) to separate by gravity at the base of the reactor. lower values mean higher residual t-butanol in polar phase. The feed mole fraction of t-butanol is 0.3. immiscible with feed over a significant composition range CH3, CH2, CH, C, OH, CH3CO, CH2CO, CHO, CH3COO, CH2COO, HCOO, CH3O, CH2O, CH-O, COO

Table 5. Polar Solvents for Pharmaceutical Ester Crystallization with Toluene solvent

Tboil, °C

Tmelt, °C

density, g cm-3

selectivity

1,3-propylene glycol 1,2-propylene glycol 1-hydroxyacetone (acetol)

215 187 146

-27 -60 -17

1.05 1.03 1.08

12.4 12.3 7.8

partition coefficients are more appropriate for productivity. In fact, just two washes with propylene glycol reduces the t-butanol content to about 3% w/w with a residual level of propylene glycol in the organic phase of only 0.5% w/w. This scheme gives good productivity in generic batch manufacturing equipment, which is common to the pharmaceuticals industry. Adequate success in this conceptual design project was obtained by flexibly searching for alternative solvents using phase equilibrium representations of multicomponent systems as computed from group contribution methods. Identifying the desired partitioning in the ternary system was particularly crucial. 4.3. Selection of Binary Solvent Mixtures for a Crystallization Process. This example concerns a pharmaceutical intermediate produced by reaction in tetrahydrofuran (THF) with the desired product obtained by crystallization. The yield from THF alone was found to be poor, and water was introduced as an

Figure 20. Ternary LLE for t-butanol, 1,3-propylene glycol, and toluene at T ) 25 °C, 1 atm from UNIFAC-LLE.79

Figure 21. Ternary LLE for t-butanol, 1,3-propylene glycol, and cyclohexane at T ) 25 °C, 1 atm from UNIFAC-LLE.79

Figure 22. Regression results of NRTL-SAC model83 parameters for 25 pure solvents.

antisolvent to increase the yield. Laboratory test results showed inconsistencies, and it was suggested that the actual amount of water in the crystallization solvent was varying because of carryover from an upstream washing step. Predictions of the solubility in the aqueous THF solution were made to determine how water content could affect the product solubility. Details of this application can be found in ref 82. This example highlights the use of techniques other than group contribution (GC) when the necessary parameters for a GC-based method are not available. Model Relations, Selection, and Information Sources. The structure of many heterocyclic pharmaceutical molecules cannot be treated with GC-based models like UNIFAC, due to missing functional groups or interaction parameters. The NRTL-SAC model83 is an alternative approach which uses characteristic surface segments to describe intermolecular interactions from surface charge density. The NRTL binary interaction parameters between the segments are fixed, and there are adjustable characteristic segment values for each molecular species. The NRTL-SAC database contains segment profiles for 130 solvents derived from available literature VLE and LLE data. For a particular solvent, the 4 solute segment parameters are regressed from solubilities in at least 4, and up to 10, pure solvents spanning the expected surface segment values. This will allow prediction of solid solubility of a solute, in pure or binary solvents of the components in the database, with sufficient accuracy for solvent-ranking and trends in ternary systems. For the present system, an NRTL-SAC model was established from existing solubility data in 25 pure solvents over the temperature range of 10-80 °C. Results. The solubility of the pharmaceutical intermediate in mixtures of THF and water with the regression and prediction results shown in Figure 22, along with solute solubilities in the mixed THF-water solvent shown in detail in Figure 23. While

4632

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009 Table 6. Calculated Group Interaction Parameters through the CI Method for the 1,2-Dichloroethane-DMSO System

CCl DMSO

Figure 23. Experimental and predicted solubilities of product in THF-water mixtures.

the calculated values were sometimes far from experiment, they were adequate for the purposes of the problem. Figure 23 shows that water was acting unexpectedly. At low concentrations, it is a cosolvent, increasing the solubility. At higher amounts, the water depresses solubility. The strong variation with the water fraction suggested why the laboratory tests on the process stream were not consistent: fluctuating carryover combined with extreme sensitivity to composition. From these results, a different and more robust crystallization process was developed. The key was the estimation of solubilities and careful scrutiny of the sensitivity of properties to variations in conditions. 5. Emerging Methods for Property Estimation The examples given above demonstrate current approaches and capabilities for property modeling when data sources are available or the opportunity for new measurements exists. We now give two examples using advanced techniques to overcome the limitations of current approaches such as group-contribution methods. They use contemporary computational techniques either directly or indirectly to predict phase equilibria. The first describes how unavailable group-contribution parameters can be obtained with a new method based only on chemical structure, with application to VLE. The concept is appealing, and it suggests an avenue for future developments. The other is for a separations process and compares different modeling techniques, including molecular methods. The results have implications for efforts to improve predictive modeling capabilities. 5.1. Estimation of Group-Contribution Parameters. Recently, Gani et al.84 have suggested how already available experimental data might be used to predict group-contribution model parameters that are missing in a host tabulation, such as the Marrero-Gani78 group-contribution method for pure component properties. The basis is an atom-connectivity index, developed under the principle of additivity of contributions of different descriptors for a specific property that gives contributions to molecular properties by atoms and their connectivities. With atoms, many fewer parameters are needed to represent groups of atoms. Further, index parameter values for connectivity indices can be found from the same available experimental data as for regressing group-contribution parameters. Combining known group contributions (GC) with estimated group contributions from atom-connectivity indexes (CI) results in an approach called GCplus. The method can be applied to any host group contribution model.

CCl

DMSO

0.0 259.82

-198.22 0.0

Extension to a wide range of property models for pure component properties and to average properties of polymer repeat units has been made.85 Gonzalez et al.86 have applied the GCplus approach to predict missing group-interaction parameters when the host method is the UNIFAC model for activity coefficients (GC). The available experimental data used for UNIFAC group contributions were employed in regressing the interaction parameters for the atomconnectivity indices (CI). Then, the CI values were used to estimate missing group-interaction parameters. Examples applying GCplus to pure component properties are given in refs 84-86. Here, we illustrate GCplus for mixture properties. In each case, the chemicals and the phases of interest are given along with the host model and the group(s) with missing values. Then predicted UNIFAC group parameters are given, along with comparison of predicted and measured phase behavior results. VLE for 1,2-Dichloroethane-DMSO. If the original UNIFAC-VLE model87 is the base method, the missing group interaction parameters are for the pair CCl and DMSO. Using the CI method, estimates of the missing group interaction parameters are listed in Table 6. Figure 24 shows Txy VLE comparisons from using these parameters along with the GC parameters in the UNIFAC table with measured VLE.88 These data were not used for the CI-model parameter estimation, but the agreement is excellent. SLE for Acetaminophen (Paracetamol) with 1-Butanol. If a later revision of the original UNIFAC-VLE model86 is the host for this system, only the ACNH2-CH2CO interactions are missing. However, for illustration, we also present results when all the group interaction parameters are obtained via CI. Table 7 lists GCplus group interaction parameters for the system where those for ACNH2/CH2CO (in bold) are estimated with CI. Table 8 lists the parameters when all are estimated from CI. The values are different. Figure 25 shows the comparisons of the results for both cases with data.89 Over the limited range of paracetamol compositions, both sets of parameters describe the data well. Thus, while parameter values may differ, the data can be adequately predicted with both methods. VLE for Methylethylketone-n-Heptane. This system involves only the CH2 and CH2CO groups, but there is a significant

Figure 24. Txy VLE diagram at 0.953 bar for 1,2-dichloroethane with DMSO from the UNFAC-CI method with parameters not used for the CI model regression compared to measured values: (experiment).88

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009 87

Table 7. Original UNIFAC and CI-Generated ACNH2/CH2CO UNIFAC Interaction Parameters for the System Acetaminophen (Paracetamol) with 1-Butanol

CH2 ACH OH ACOH CH2CO ACNH2

CH2

ACH

OH

ACOH

CH2CO

ACNH2

0 -11.12 156.04 275.8 26.76 1139

61.13 0 89.6 25.34 140.1 247.5

932.65 636.100 0 -451.60 164.5 -17.4

1333 1329 -259.7 0 -133.1 -253.1

476.39 25.77 84 -356.1 0 4648

920.7 648.2 -52.39 119.9 -44.56 0

Table 8. CI-Generated UNIFAC Interaction Parameters for the System Acetaminophen (Paracetamol) with 1-Butanol

CH2 ACH OH ACOH CH2CO ACNH2

CH2

ACH

OH

ACOH

CH2CO

ACNH2

0 -46.12 52.07 29.3 -24.83 35.65

112.34 0 98.03 777.92 -60.15 371.2

932 571.97 0 135.08 432.62 -101.82

327.64 -207.56 -379.27 0 423.34 -1193.02

531.84 195.38 -135.03 567.88 0 4648.0

1232.09 966.75 -731.19 -318.1 -44.56 0

Table 9. CI-Computed Parameter Matrix for Dortmund Group Contributions90 aij

CH2 CH2CO

bij

cij

CH2

CH2CO

CH2

CH2CO

CH2

CH2CO

0 -1067.87

7824.76 0

0 3.82

-35.25 0

0 0

0 0

temperature variation to be dealt with. The later UNIFAC model90 has parameter values, so a comparison can be made among data, GC prediction, and CI prediction. The CI-computed parameter matrix with temperature dependence is given in Table 9. Figure 26 shows the calculated Pxy diagram with the UNIFAC Dortmund parameters and the CI-computed parameters, along with data from ref 91. The agreement for the CI method is not as good as with the GC method over the whole data range, but the pressure and composition of the azeotrope are given reasonably accurately. These examples illustrate the possibilities of using models based on limited information, such as connectivity indices, as well as a level of compromise encountered when they are used in place of more elaborate methods such as UNIFAC. In general, CI may be a reliable expedient to determine unavailable, and perhaps less sensitive, parameters for use with incomplete group contribution methods. It is not proposed as a replacement for experiments, but rather to focus on a few experiments through which the extension can be verified. 5.2. Molecular Calculations. While methods such as CI are easily used to obtain group-contribution parameters, their accuracy

Figure 25. Solubility of paracetamol in 1-butanol as a function of temperature estimated with CI-generated values: s only for ACNH2/CH2CO groups (original UNIFAC87 parameters for other groups), - - - all parameters estimated from CI; (experiment)89 [.

4633

and generality may be limited. An alternative which does not, in principle, require data for model parameter regression is quantum chemistry calculations for inter- and intramolecular force fields followed by molecular simulation or statistical thermodynamic methods to obtain properties. Solvents for Extractive Distillation of 1,3-Butadiene. Mathias et al.92 describe an investigation using quantum mechanics and molecular simulation to improve process simulation for the classical problem of 1,3-butadiene recovery from steam cracker C4 hydrocarbons by determining the relative effectiveness of n,ndimethylformamide (DMF) and acetonitrile (ACN) as extractivedistillation solvents. The principal properties obtained were the activity coefficients of the hydrocarbon components in the presence of the extractive solvent for use in eq 1. Comparisons were made among a quantum mechanical and statistical mechanical method, COSMO-RS,36,37 a molecular dynamics simulation approach, SPEADMD,93 group contributions from UNIFAC,94 and “thermodynamic intuition”. Mathias et al.92 describe the methods and results in some detail; only a brief summary is given here to indicate the findings. The COSMO-RS method reliably predicted the trends of infinitedilution activity coefficients with accuracy comparable to UNIFAC, but only with systematic empirical corrections. This limited the true predictive capability of the method. The SPEADMD molecular simulation used a force field from the principle of transferability,95,96 which assumes that forces inferred from experimental data for one set of mixtures can be applied to other compounds and mixtures. The computed results provided unique qualitative structural and orientational insights at the molecular scale about the solvation interactions between the polar solvents and the olefinic moieties in the hydrocarbon compounds. The differences in accessibility for DMF and ACN and the sizes and shapes that affect intermolecular contacts were reliably characterized. However, to achieve accuracy for activity coefficients, the molecular simulations required refinement of the interaction potentials by regression to data, similar to finding UNIFAC parameters. Extensions of Molecular Calculations. The experience of Mathias et al.92 suggests some of the limitations and future prospects of molecular simulation, as do the International Fluid Property Simulation Challenges (IFPSC).97 The present important question about the potential for molecular simulation as a routine tool to provide quantitative property data for process and product design is, “Are we there yet?”. In our opinion, the answer is a qualified no. While progress is being made, the results are often like the butadiene example: good, perhaps adequate for the advice role without high accuracy, but not sufficient for the service role. A common shortcoming of molecular simulation methods is the lack of easily available and suitable force fields98 to solve the wide variety of problems under consideration.99 In particular, what should be done when no experimental data exist for empirically fitting a force field? As an example, consider the molecule whose chemical formula is C10H19N and structure is the following:

There are no experimental data and not even a CAS number has been assigned. If this molecule is of interest in a product design for a particular application or is an impurity that must be effectively removed in a process design, molecular simulation could not be used for property estimation unless ab inito quantum methods alone could produce a force field. Figure 27

4634

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Figure 26. VLE for the methylethylketone-n-heptane system at 318.15 K: s calculated with UNIFAC90 group-contribution parameters, - - - calculated with CI-estimated group-contribution parameters; (experiment)90 [ y, 9 x.

erroneous results. In such cases, what would nonexperts find? We mention two efforts to produce standard molecular simulation tools which are the TOWHEE project102 and the LAMMPS project;103 others should appear in the future. 5.3. Other Aspects of Emerging Methods. Recent years have shown two particularly important developments in properties modeling that combine molecular theory and thermodynamic relations. One is the SAFT equation of state104 whose original formulation was computationally complicated and initially had restricted applications. SAFT and PC-SAFT105 are now appearing in many variations with improved efficiency and reliability. In particular, formulations for polymers,106 electrolytes,107 ionic liquids,108 and group contributions109 have appeared. Another development of the significance is the continued extensions and refinement of COSMO methodology to a variety of systems.110,111 As time goes along, we expect both approaches to provide practical results for challenging systems and be worthy of consideration by model developers and users into the future. Figure 27. Schematic for developing molecular simulation force fields for cases with different amounts of available data.

6. Conclusions

suggests a strategy to obtain force field parameters of molecules of industrial interest. Other limitations in simulations occur at lower temperatures and for larger molecules (high density) and for transport properties. For quantum calculations, the system size may be limiting, the fundamental basis sets might not be accurate enough, and the way to improve results may not be clear. Also, while only a few molecular-simulation researchers and reviewers now list computing machinery and computing resources as a major limitation in the extension of MC and MD to new fluid property applications, computational capabilities beyond those currently available are needed for computational chemistry to directly treat many practical systems. Finally, a key issue for practical application molecular simulation is the lack of availability of tools for nonexpert users. There is no standard toolbox for molecular simulation as pointed out by Wei.100,101 Our experience is that even experts can have problems, giving one pause about very widespread application of computation. For example, in at least two entries during IFPSC contests, expert researchers incorrectly transcribed molecular parameters that then produced nonrepresentative and

For 100 years, the Industrial and Engineering Chemistry journals have published indispensable information about properties: data, models, behaviors, designs, and analyses. We have presented a contemporary perspective of the importance, stateof-the-art, and the potential for continued improvement in accuracy, reliability, and efficiency for thermodynamic properties in the design and optimization of processes and products. Emphasis has been given to skillful modeling strategies, intelligent resource utilization, and adequate validation via experiment and consistency checks. Several examples focused on complexing solutions for process simulation, where the role of property models is in service; model development for solvent substitution, where the role is service plus advice; and emerging techniques for prediction of properties from molecular structure. The ranges of these systems illustrate the kinds of systems that can now be directly addressed when appropriate tools, experimental information, and computational methods are applied. As the sophistication of chemical products and processes increases, along with greater societal demands for sustainability, health, safety, and economy, the ability to suitably estimate property values, and implement property models, must grow in effective-

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

ness and comprehensiveness. We believe that the concepts, structure, and information described here can provide a sound basis for developing and implementing such capability. Nomenclature Symbols ai ) activity of species i GE ) excess Gibbs energy, kJ mol-1 H ) Henry’s constant, MPa ∆H0 ) standard state enthalpy of reaction, kJ mol-1 k ) equilibrium constant for reaction m ) mass flow rate of process stream, kg s-1 mi ) molal composition, mol kg-1 Mw ) molecular weight p ) pressure, MPa p(0) ) standard state pressure, 0.1 MPa ps ) vapor pressure, MPa R ) gas constant, 8.314 J K-1 mol-1 Sp ) solvent power, eq 10 t ) time, s T ) temperature, K V ) set of measurable variables Vi ) molar volume, cm3 mol-1 x, xi ) liquid mole fraction yi ) vapor mole fraction of species i Greek Symbols γt, γ*, t activity coefficient of species i relative to pure component standard state, Henry’s law standard state δh, δp ) Hansen solubility parameters for special interactions θ ) conceptual variable νt ) stoichiometric coefficient of species i in a reaction φi, φis ) mixture fugacity coefficient and saturation fugacity coefficient of species i

Literature Cited (1) Ind. Eng. Chem. 1909, 1 (1). (2) Chen, C.-C.; Mathias, P. M. Applied Thermodynamics for Process Modeling. AIChE J. 2004, 48 (2), 200. (3) Kontogeorgis, G. M.; Gani, R. Computer Aided Property Estimation for Process and Product Design. Computer-Aided Chemical Engineering; Elsevier; Amsterdam, 2004; Vol. 19. (4) Mathias, P. M. Applied thermodynamics in chemical technology: current practice and future challenges. Fluid Phase Equilib. 2005, 228229, 49–57. (5) Carlson, E. C. Don’t gamble with physical properties for simulation. Chem. Eng. Progress 1996, (October), 35. (6) Zudkevitch, D. Forensic thermodynamics - Erroneous decisions on thermodynamic data can cause plant failures; EFCE Publication Series; Taylor & Francis: U.K., 1980; Issue 11, pp 885-905. (7) http://www.cas.org/expertise/cascontent/registry, accessed March 28, 2009. (8) Gani, R.; O’Connell, J. P. Properties and CAPE: From Present Uses to Future Challenges. Comput. Chem. Eng. 2001, 25 (1), 3–14. (9) O’Connell, J. P.; Neurock, M. Trends in property estimation for process and product design. AIChE Symp. Ser. 2000, 96 (323), 5–23. (10) Xin, Y.; Whiting, W. B. Case Studies of Computer-Aided Design Sensitivity to Thermodynamic Data and Models. Ind. Eng. Chem. Res. 2000, 39, 2998–3006. (11) Brignole, E. A.; Gani, R.; Romagnoli, J. A. A simple algorithm for sensitivity and operability analysis of separation processes. Ind. Eng. Chem. Process. Des. DeVel. 1985, 24, 42–48. (12) http://trc.nist.gov/, accessed March 28, 2009. (13) http://webbook.nist.gov/chemistry/, accessed March 28, 2009. (14) http://dippr.byu.edu/, accessed March 28, 2009. (15) http://www.ppds.co.uk/default.asp, accessed March 28, 2009. (16) http://www.nist.gov/srd/nsrds.htm, accessed March 28, 2009. (17) http://www.ddbst.de/new/default.htm, accessed March 28, 2009. (18) http://i-systems.dechema.de/detherm/mixture.php?, accessed March 28, 2009.

4635

(19) http://www.lib.gla.ac.uk/Resources/Databases/esdu.shtml, accessed March 28, 2009. (20) Harvey, A. H. Physical and Chemical Properties. In Albright’s Chemical Engineering Handbook; Albright, L., Ed.; CRC Press: Boca Raton, FL, 2008; Chapter 1. (21) Frenkel, M.; Chirico, R. D.; Diky, V.; Dong, Q.; Marsh, K. N.; Dymond, J. H.; Wakeham, W. A.; Stein, S. E.; Konigsberger, E.; Goodwin, A. R. H. XML-based IUPAC standard for experimental, predicted, and critically evaluated thermodynamic property data storage and capture (ThermoML). Pure App. Chem. 2006, 78, 541. (22) Poling, B. E.; Prausnitz, J. M.; O’Connell, J. P. The Properties of Gases and Liquids, fifth ed.; McGraw-Hill: New York; 2001; errata at http:// www.che.virginia.edu/PGL5/. (23) McQuarrie, D. A. Statistical Mechanics; Harper & Row: New York, 1976; pp 255-288. (24) Haile, J. M. Molecular Dynamics Simulation; Wiley-Interscience: New York, 1992. (25) Ungerer, P.; Nieto-Draghi, C.; Rousseau, B.; Ahunbay, G.; Lachet, V. Molecular simulation of the thermophysical properties of fluids: From understanding toward quantitative predictions. J. Mol. Liq. 2007, 134, 71. (26) Panagiotopoulos, A. Z. Direct Determination of Fluid Phase Equilibria by Simulation in the Gibbs Ensemble: A Review. Mol. Simul. 1992, 9, 1. (27) Cummings, P. T. Molecular Dynamics Simulation of Realistic Systems. Fluid Phase Equilib. 1996, 116, 237. (28) Prausnitz, J. M. Athena, Hercules and Nausica: Three dimensions of chemical engineering in the twenty-first century. Fluid Phase Equilib. 2007, 261, 3. (29) de Pablo, J. J. Molecular Simulations in Chemical Engineering: Present and Future. AIChE J. 2002, 48, 2716. (30) Ghorai, P. K.; Glotzer, S. C. Molecular Dynamics Simulation Study of Self-Assembled Monolayers of Alkanethiol Surfactants on Spherical Gold Nanoparticles. J. Phys. Chem. C 2007, 111, 15857. (31) http://www.aiche.org/DivisionsForums/ViewAll/CoMSEF.aspx, accessed March 28, 2009. (32) Olson, J. D. Importance of Experimental Measurements. In Report on Forum 2000: Fluid Properties for New Technologies-Connecting Virtual Design with Physical Reality; NIST Special Publication 975; U.S. Government Printing Office: Washington, DC, 2001; pp 102-104. (33) Gupta, S.; Olson, J. D. Industrial Needs in Physical Properties. Ind. Eng. Chem. Res. 2003, 42, 6359. (34) Gani, R.; Muro-Sune, N.; Sales-Cruz, M.; Leibovici, C.; O’Connell, J. P. Mathematical and numerical analysis of classes of property models. Fluid Phase Equilib. 2006, 250, 1–32. (35) Fredenslund, A.; Gmehling, J.; Rasmussen, P. Vapor-Liquid Equilibria Using UNIFAC; Elsevier: Amsterdam, 1977. (36) (a) Klamt, A.; Eckert, F. COSMO-RS: A novel and efficient method for the apriori prediction of thermophysical data of liquids. Fluid Phase Equilib. 2000, 172, 43. (b) Klamt, A.; Eckert, F. Erratum. Fluid Phase Equilib. 2003, 205 (2), 357. (37) Zilnik, L. F.; Jazbinsek, A.; Hvala, A.; Vrecer, F.; Klamt, A. Solubility of sodium diclofenac in different solvents. Fluid Phase Equilib. 2007, 261 (1-2), 140–145. (38) Wang, S.; Sandler, S. I.; Chen, C.-C. Refinement of COSMO-SAC and the Applications. Ind. Eng. Chem. Res. 2007, 46, 7275. (39) (a) Kuhnert, C.; Albert, M.; Breyer, S.; Hahnenstein, I.; Hasse, H.; Maurer, G. Phase Equilibrium in formaldehyde containing multicomponent mixtures: Experimental results for the vapor-liquid equilibrium of (formaldehyde + methylal + (water or methanol)) and (formaldehyde + water + methanol + methylal) and comparison with predictions. Ind. Eng. Chem. Res. 2006, 45 (14), 5155–5164. (b) Kuhnert, C.; Albert, M.; Breyer, S.; Hahnenstein, I.; Hasse, H.; Maurer, G. Ind. Eng. Chem. Res. 2006, 45 (17), 6093–6094. (40) Kogan, L. V. State of the Vapor Phase Above Solutions of Formaldehyde in Water and Methanol. Zhur. Prikl. Khim. 1979, 52, 2722. (41) Hahnenstein, I.; Hasse, H.; Kreiter, C. G.; Maurer, G. 1H- and 13C NMR-spectroscopic study of chemical equilibria in solutions of formaldehyde in water, deuterium oxide and methanol. Ind. Eng. Chem. Res. 1994, 33, 1022–1029. (42) Balashov, A. L.; Danov, S. M.; Golovkin, A. Yu.; Krasnov, V. L.; Ponomarev, A. N.; Borisova, I. A. Equilibrium Mixture of Polyoxymethylene Glycols in Concentrated Aqueous Formaldehyde Solutions. Russ J. App. Chem 1996, 2, 190–192. (43) Albert, M.; Coto Garcı´a, B. C.; Kreiter, C. G.; Maurer, G. Vaporliquid and chemical equilibria of formaldehyde-water mixtures. AIChE J. 1999, 45, 2024–2033. (44) Albert, M.; Hahnenstein, I.; Hasse, H.; Maurer, G. Vapor-Liquid Equilibrium of Formaldehyde Mixtures: New Experimental Data and Model Revision. AIChE J. 1996, 42 (6), 1741–1752.

4636

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

(45) Credali, L.; Mortillaro, L.; Galiazzo, G.; Russo, M.; De Checci, C. Pressione di vapore sul sistema H2O-CH2O liquido e solido. La Chim. e L’Ind. 1965, 47 (7), 732–736. (46) Kogan, L. V. NMR-Study of the State of Aqueous Methanol Solutions of Formaldehyde. Zhur. Prikl. Khim. 1979, 52 (12), 2725–2730. (47) Maurer, G. Vapor-Liquid Equilibrium of Formaldehyde- and WaterContaining Multicomponent Mixtures. AIChE J. 1986, 32, 932–948. (48) Hasse, H.; Maurer, G. Heat of Dilution in Aqueous and Methanolic Formaldehyde Solutions. Ber. Bunsen-Ges. Phys. Chem. 1992, 96, 83–96. (49) Liu, Y.-Q.; Hasse, H.; Maurer, G. Enthalpy Change on Vaporization of Aqueous and Methanolic Formaldehyde Solutions. AIChE J. 1992, 38, 1693–1702. ´ .; Maurer, G. Solubility (50) Ermatchkov, V.; Pe´rez-Salado Kamps, A of Carbon Dioxide in Aqueous Solutions of N-Methyldiethanolamine in the Low Gas Loading Region. Ind. Eng. Chem. Res. 2006, 45, 6081–6091. ´ .; Speyer, D.; Maurer, G. (51) Ermatchkov, V.; Pe´rez-Salado Kamps, A Solubility of Carbon Dioxide in Aqueous Solutions of Piperazine in the Low Gas Loading Region. J. Chem. Eng. Data 2006, 51, 1788–1796. (52) Mathias, P. M.; Reddy, S.; O’Connell, J. P. Quantitative Evaluation of the Aqueous-Ammonia Process for CO2 Capture Using Fundamental Data and Thermodynamic Analysis. Energy Procedia, in press. (53) Edwards, T. J.; Maurer, G.; Newman, J.; Prausnitz, J. M. Vaporliquid equilibria in multicomponent aqueous solutions of volatile weak electrolytes. AIChE J. 1978, 24, 966–976. (54) Patterson, C. S.; Slocum, G. H.; Busey, R. H.; Mesmer, R. E. Carbonate equilibria in hydrothermal systems: First ionization of carbonic acid in NaCl media to 300 °C. Geochim. Cosmochim. Acta 1982, 46, 1653–1663. (55) Patterson, C. S.; Busey, R. H.; Mesmer, R. E. Second ionization of carbonic acid in NaCl media to 250 °C. J. Sol. Chem 1984, 13, 647–661. ´ .; Maurer, G. Dissociation constant of (56) Pe´rez-Salado Kamps, A N-methyldiethanolamine in aqueous solution at temperatures from 278 to 368 K. J. Chem. Eng. Data 1996, 41, 1505–1513. (57) Hetzer, H. B.; Robinson, R. A.; Bates, R. G. Dissociation constants of piperazinium ion and related thermodynamic quantities from 0 to 50°. J. Phys. Chem. 1968, 72, 2081–2086. ´ .; Maurer, G. Chemical (58) Ermatchkov, V.; Pe´rez-Salado Kamps, A Equilibrium Constants for the Formation of Carbamates in the System CO2 + Piperazine + Water from 1H-NMR-spectroscopy. J. Chem. Thermodyn 2003, 35, 1277–1289. (59) Pitzer, K. S. Thermodynamics of electrolytes. 1. Theoretical basis and general equations. J. Phys. Chem. 1973, 77, 268–277. ´ .; Balaban, A.; Jo¨decke, M.; Kuranov, G.; (60) Pe´rez-Salado Kamps, A Smirnova, N. A.; Maurer, G. Solubility of Single Gases Carbon Dioxide and Hydrogen Sulfide in Aqueous Solutions of N-Methyldiethanolamine at Temperatures from 313 to 413 K and Pressures up to 7.6 MPa: New Experimental Data and Model Extension. Ind. Eng. Chem. Res. 2001, 40, 696–706. (61) Kuranov, G.; Rumpf, B.; Maurer, G.; Smirnova, N. A. Vapor-liquid equilibrium in aqueous systems containing methyldiethanolamine, carbon dioxide and hydrogen sulfide. Fluid Phase Equilib. 1997, 136, 147–162. (62) Maurer, G. Modeling the Liquid-Liquid Equilibrium for the Recovery of Carboxylic Acids from Aqueous Solutions: A Review. Fluid Phase Equilib. 2006, 241, 86–95. (63) Schunk, A.; Maurer, G. On the Influence of some Strong Electrolytes on the Partitioning of Acetic Acid to Aqueous/Organic Two-Phase Systems in the Presence of Tri-n-octylamine. Part I: Methyl isobutyl ketone as Organic Solvent. Fluid Phase Equilib. 2006, 239, 223–239. (64) Schunk, A.; Maurer, G. Distribution of Hydrochloric, Nitric, and Sulfuric Acid between Water and Organic Solutions of Tri-n-octylamine: Part II. Methylisobutylketone as Organic Solvent. Fluid Phase Equib 2003, 211, 189–209. (65) Schunk, A.; Maurer, G. Distribution of Hydrochloric, Nitric, and Sulfuric Acid between Water and Organic Solutions of Tri-n-octylamine. Part I. Organic Solvent Toluene. Fluid Phase Equilib. 2003, 207, 1–21. (66) Stopperka, K.; Kilz, F. Die Zusammensetzung der Gasphase u¨ber dem Flu¨ssigen System H2SO4-SO3 in Abha¨ngigkeit von der Temperatur. Z. Annorg. Allg. Chem. 1969, 370, 59–66. (67) Nilges, J.; Schrage, J. Vapor-Liquid Equilibrium of the System H2SO4-SO3. Part II. Thermodynamic Description with Regard to the Formation of H2S2O7. Fluid Phase Equilib. 1991, 68, 247–261. (68) Mathias, P. M.; Chen, C.-C.; Zou, B.; Randolph, D., III; Doering, F. A Representation of the Thermodynamic Properties of Sulfuric Acid and Oleum. Presented at the AIChE Annual Meeting, Reno, Nevada, November 4-9, 2001. (69) Prausnitz, J. M.; Lichtenthaler, R.; de Azevedo, E. G. Molecular Thermodynamics of Fluid Phase Equilibria, 3rd ed.; Prentice-Hall: Englewood Cliffs, NJ, 2001. (70) Miles, F. D.; Niblock, H.; Smith, D. The Heat of Formation of Oleum. Trans. Faraday Soc. 1944, 40, 281–295.

(71) Schrage, J. Vapor-Liquid Equilibrium of the System H2SO4-SO3. Part I. Vapor Pressure Measurements with a New Static Apparatus. Fluid Phase Equilib. 1991, 68, 229–245. (72) Llinas, A.; Glen, R. C.; Goodman, J. M. Solubility Challenge: Can you predict solubilities of 32 molecules using a database of 100 reliable measurements? J. Chem. Inf. Model. 2008, 48, 1289. (73) Price, S. L. Computational prediction of organic crystal structures and polymorphism. Int. ReV. Phys. Chem. 2008, 27 (3), 541. (74) Abbas, S.; Leigh, F.; Norton, A. K.; Powell, L.; Robinson, G. E.; Siedlecki, P.; Southworth, R. J.; Stark, A.; Williams, E. G. Application of an Enantiomerically Pure Bicyclic Thiolactone in the Synthesis of a Farnesyl Transferase Inhibitor. Org. Process Res. DeV. 2008, 12 (2), 202. (75) Barton, A. F. M. CRC handbook of solubility parameters and other cohesion parameters, 2nd ed.; CRC Press; Boca Raton, FL, 1991. (76) Hansen, C. M. Hansen solubility parameters, a users handbook; CRC Press: Boca Raton, FL, 2000. (77) Gani, R. ProCAMD Manual; Technical University of Denmark: Lyngby, Denmark, 2003. (78) Marrero, J.; Gani, R. Group Contribution Based Estimation of Pure Component Properties. Fluid Phase Equilib. 2001, 183-184, 183–208. (79) Magnusson, T.; Rasmussen, P.; Fredenslund, A. UNIFAC Parameter Table for Prediction of Liquid-Liquid Equilibria. Ind. Eng. Chem. Process Des. DeV. 1981, 20 (2), 331–339. (80) Meloun, M.; Bordovska, S. Benchmarking and validating algorithms that estimate pKa values of drugs based on their molecular structures. Anal. Bioannal. Chem. 2007, 389, 1267. (81) www.acdlabs.com/products/phys_chem_lab/pka, accessed March 28, 2009. Advanced Chemistry Development Inc.: Toronto, Canada. (82) Crafts, P. A. The Role of Crystallisation and Solubility Modelling in the Design of Active Pharmaceutical Ingredients. In Computer-Aided Chemical Engineering; Ng, K. M., Gani, R., Dam-Johansen, K., Eds.; Elsevier BV: Amsterdam, 2006; Vol. 23. (83) Chen, C.-C.; Crafts, P. A. Correlation and Prediction of Drug Molecule Solubility in Mixed Solvent Systems with the Nonrandom TwoLiquid Segment Activity Coefficient (NRTL-SAC) Model. Ind. Eng. Chem. Res. 2006, 45, 4816–4824. (84) Gani, R.; Harper, P. M.; Hostrup, M. Automatic creation of missing groups through connectivity index for pure-component property prediction. Ind. Eng. Chem. Res. 2005, 44 (18), 7262. (85) Satyanarayana, K. C.; Gani, R.; Abildskov, J. Polymer Property Modeling Using Grid Technology for Design of Structured Products. Fluid Phase Equilib. 2007, 261, 58–63. (86) Gonza´lez, H. E.; Abildskov, J.; Gani, R.; Rousseauz, P.; Le Bert, B. Prediction of Missing UNIFAC Group-Interaction Parameters through Connectivity Indices. AIChE J. 2007, 53 (6), 1620–1632. (87) Hansen, H. K.; Rasmussen, P.; Fredenslund, A.; Schiller, M.; Gmehling, J. Vapor-liquid equilibria by UNIFAC group-contribution. 5. Revision and extension. Ind. Eng. Chem. Res. 1991, 30, 2352–2355. (88) Radhamma, M.; Hsieh, C.-T.; Venkatesu, P.; Rao, M. V. P.; Lee, M.J.; Lin, H.-m. Isobaric Vapor-Liquid Equilibrium for Dimethylsulfoxide with Chloroethanes and Chloroethenes. J. Chem. Eng. Data 2008, 53, 374–377. (89) Granberg, R. A.; Rasmuson, Aa. C. Solubility of Paracetamol in Binary and Ternary Mixtures of Water + Acetone + Toluene. J. Chem. Eng, Data 2000, 45 (3), 478–483. (90) Wittig, R.; Lohmann, J.; Gmehling, J. Vapor-liquid equilibria by UNIFAC Group Contribution. 6. Revision and Extension. Ind. Eng. Chem. Res. 2003, 42, 183–188. (91) Takeo, M.; Nishi, K.; Nitta, T.; Katayama, T. Isothermal vaporliquid equilibriums for two binary mixtures of heptane with 2-butanone and 4-methyl-2-pentanone measured by a dynamic still with a pressure regulation. Fluid Phase Equilib. 1979, 3, 123–131. (92) Mathias, P. M.; Elliott, J. R., Jr.; Klamt, A. Butadiene Purification Using Polar Solvents. Analysis of Solution Nonideality Using Data and Estimation Methods. Ind. Eng. Chem. Res. 2008, 47 (15), 4996–5004. (93) Unlu, O.; Gray, N. H.; Gerek, Z. N.; Elliott, J. R., Jr. Transferable Step Potentials for the Straight-Chain Alkanes, Alkenes, Alkynes, Ethers, and Alcohols. Ind. Eng. Chem. Res. 2004, 43, 1788. (94) Lohmann, J.; Gmehling, J. From UNIFAC to Modified UNIFAC (Dortmund). Ind. Eng. Chem. Res. 2001, 40, 957. (95) Martin, M. G.; Siepmann, J. I. Transferable potentials for phase equilibria. 1. United-atom description of n-alkanes. J. Phys. Chem. B 1998, 102, 2569. (96) Elliott, J. R., Jr.; Gray, N. H. Asymptotic trends in thermodynamic perturbation theory. J. Chem. Phys. 2005, 123, 184902. (97) http://fluidproperties.org, accessed March 28, 2009. (98) http://www.ctcms.nist.gov/∼fstarr/ptpfms/home.html, accessed March 28, 2009. (99) Spiegel, K.; De Grado, W. F.; Klein, M. L. Structural and dynamical properties of manganese catalase and the synthetic protein DF1 and their

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009 implication for reactivity from classical molecular dynamics calculations. Proteins: Struct., Funct., Bioinf. 2006, 65, 317. (100) Wei, J. Paradigms of Chemical Engineering. 10th Internationl Conference on Properties and Phase Equilibria for Product and Process Design, Snowbird, UT, May 16-21, 2004. (101) Wei, J. Molecular Structure and Property: Product Engineering. Ind. Eng. Chem. Res. 2002, 41, 1917. (102) http://towhee.sourceforge.net/. (103) http://lammps.sandia.gov/. (104) Chapman, W. G.; Gubbins, K. E.; Jackson, G.; Radosz, M. SAFT - Equation-of-State Solution Model for Associating Fluids. Fluid Phase Equilib. 1989, 52, 31–38. (105) Gross, J.; Sadowski, G. Perturbed-Chain SAFT: An Equation of State Based on a Perturbation Theory for Chain Molecules. Ind. Eng. Chem. Res. 2001, 40, 1244–1260. (106) Kouskoumvekaki, I.; von Solms, N.; Michelsen, M. L.; Kontogeorgis, G. Applied perturbed chain SAFT equation of state to complex polymer systems using simple mixing rules. Fluid Phase Equilib. 2004, 215, 71–78.

4637

(107) Zhao, H.; dos Ramos, M. C.; McCabe, C. Development of an equation of state for electrolyte solutions by combining the statistical associating fluid theory and the mean spherical approximation for the nonprimitive model. J. Chem. Phys. 2007, 126 (24), 244503. (108) Kroon, M. C.; Karatsani, E. K.; Economou, I. G.; Witkamp, G. J.; Peters, C. J. Modeling the carbon dioxide solubility in imidazolium-based ionic liquids with the tPC-SAFT equation of state. J. Phys. Chem. B 2006, 110, 9262. (109) Tihic, A.; Kontogeorgis, G. M.; von Solms, N.; Michelsen, M. L.; Constantinou, L. A Predictive Group-Contribution Simplified PC-SAFT Equation of State: Application to Polymer Systems. Ind. Eng. Chem. Res. 2008, 47 (15), 5092–5101. (110) http://www.cosmologic.de/Symposium/symposium.html. (111) http://www.design.che.vt.edu/VT-Databases.html.

ReceiVed for reView October 11, 2008 ReVised manuscript receiVed March 11, 2009 Accepted March 13, 2009 IE801535A