Antisolvent crystallization: Model identification, experimental validation and dynamic simulation

Chemical Engineering Science 63 (2008) 5457 -- 5467

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: w w w . e l s e v i e r . c o m / l o c a t e / c e s

Antisolvent crystallization: Model identification, experimental validation and dynamic simulation S. Mostafa Nowee a , Ali Abbas a,∗ , Jose A. Romagnoli b a b

School of Chemical and Biomolecular Engineering, University of Sydney, Australia Department of Chemical Engineering, Louisiana State University, USA

A R T I C L E

I N F O

Article history: Received 31 July 2008 Accepted 1 August 2008 Available online 7 August 2008 Keywords: Antisolvent crystallization Population balance Kinetics Parameter identification Particle size Pharmaceuticals

A B S T R A C T

This paper is concerned with the development, simulation and experimental validation of a detailed antisolvent crystallization model. A population balance approach is adopted to describe the dynamic change of particle size in crystallization processes under the effect of antisolvent addition. Maximum likelihood method is used to identify the nucleation and growth kinetic models using data derived from controlled experiments. The model is then validated experimentally under a new solvent feedrate profile and showed to be in good agreement. The resulting model is directly exploited to understand antisolvent crystallization behavior under varying antisolvent feeding profiles. More significantly, the model is proposed for the subsequent step of model-based optimization to readily develop optimal antisolvent feeding recipes attractive for pharmaceutical and chemicals crystallization operations. © 2008 Elsevier Ltd. All rights reserved.

1. Introduction Crystallization is a widely used technique in solid–liquid separation processes. The driving force in crystal formation is supersaturation. The trend of supersaturation generation during the process has a direct and significant role on crystal characteristics such as size, morphology, and purity. Keeping the supersaturation constant during a crystallization operation has been a long-standing technique for arguably optimal operation (Jones and Mullin, 1974; Mullin and Nyvlt, 1971). Cooling and solvent evaporation are two ways of achieving supersaturation control. In the last decade, salting-out as a means to induce supersaturation has been drawing more attention. In this method, a secondary solvent known as antisolvent or precipitant is added to the solution resulting in the reduction of the solubility of the solute in the original solvent and consequently generating a supersaturation driving force. The rate of supersaturation generation in antisolvent crystallization is highly dependent on the antisolvent addition rate. Antisolvent crystallization is an advantageous method where the substance to be crystallized (solute) is highly soluble, has solubility that is a weak function of temperature, is heat sensitive, or unstable in high temperatures (Doki et al., 2002). This technique is an energy-saving alternative to evaporative crystallization, if the antisolvent can be separated at low (energy) costs. Disadvantages of antisolvent crystallization include its high-dependence on mixing,

* Corresponding author. Fax: +61 2 9351 2854.

E-mail address: [email protected] (A. Abbas). 0009-2509/$ - see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2008.08.003

for instance, where poor mixing regimes exist, high local supersaturation at antisolvent induction zones exists leading to excess primary nucleation consequently resulting in fine crystal particle formation that easily tend to agglomerate (Takiyama et al., 1998). The strive to understand antisolvent crystallization phenomena and to develop systematic operational policies towards crystal size control have been the interest of many recent and overly experimental investigations. Studies looking at the effect of various operating conditions are numerous. The process of the salting-out precipitation of cocarboxylase hydrochloride from its aqueous solution by addition of acetone was studied by Budz et al. (1986). They found that the crystallization may be remarkably improved after careful selection of operating parameters viz. seeding, agitation and acetone dosage rate. Takiyama et al. (1998) investigated the effect of different concentrations of aqueous and antisolvent solutions on crystal shape and distribution. The main variable they studied was supersaturation. They proposed a mechanistic formulation in finding the relation between nucleation and supersaturation. Other works considered antisolvent concentration and feed rate effects on final crystal habit (Doki et al., 2002; Holmback and Rasmuson, 1999; Kaneko et al., 2002; Kim et al., 2003; Kitamura and Nakamura, 2002; Kitamura and Sugimoto, 2003; Oosterhof et al., 1999a; Taboada et al., 2000; Yu et al., 2005). A study of the salting-out precipitation process was made to obtain fine particle size by varying initial supersaturation (Toth et al., 2005). It was shown that the process was heavily dependent on initial supersaturation and it was possible to produce particles with mean sizes that varied within a wide range. The effect of seeding at constant supersaturation was investigated

5458

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

in antisolvent crystallization of paracetamol from a water–acetone mixture (Yu et al., 2006). Since these kinds of systems are complex due to three (or sometimes more) components interacting with each other in two phases, some thermodynamic and mass transfer studies are essential. Ho-Gutierrez et al. (1994) plotted solubility and phase diagrams for the two ternary systems of polyethylene glycol and aqueous solutions of sodium sulphate and sodium chloride. Oosterhof et al. (1999b) also investigated different antisolvents for sodium carbonate–water system and they found a correlation between salt solubility and water–antisolvent mixture decomposition at constant temperature. Unlike most works that study inorganic solute crystallization by organic antisolvents, Holmback and Rasmuson (1999) and Yu et al. (2005) considered organic solutions and water as second solvent. Size and morphology of benzoic acid crystals from ethanol–water system were investigated by Holmback and Rasmuson (1999) who varied the feed and bulk solution composition and feed rate. They observed that supersaturation mainly governs the mean size while solvent composition has significant effect on the crystal shape. Agglomeration and habit of paracetamol crystals were the focus of the study by Yu et al. (2005) who varied the agitation speed and feedrate. In their study they concluded that low agitation speed and high feed rate will result in excessive nucleation due to high local supersaturation leading to highly agglomerated product with lower mean size. In some other works antisolvent crystallization was applied in the production of anhydrous crystals. For instance, Oosterhof et al. (2001) studied production of anhydrous sodium carbonate at room temperature using ethylene glycol and diethylene glycol as antisolvents. Unlike conventional crystallization methods that cannot produce anhydrous sodium carbonate (due to its instability in aqueous solutions at temperatures lower than boiling point of saturated soda ash solution), their method incorporating antisolvent showed success in doing so. Zhou et al. (2006) have carried out concentration controlled seeded antisolvent crystallization of a pharmaceutical compound using an algebraic equation for the solubility as a function of solvent. The main objective of their feedback concentration control system was to keep the supersaturation low and constant. For this, different constant supersaturation values were investigated and their influence over nucleation discussed. In addition, they present simulation results investigating the antisolvent crystallization of paracetamol. Nonoyama et al. (2006) also present a simulation study on seeded solvent crystallization of an active pharmaceutical ingredient (API) by water addition to original solution. The model they incorporate neglects nucleation, breakage and agglomeration, and only considers a size independent growth kinetic derived from experimental studies. Their results show a prevention of nucleation of undesired polymorph by keeping supersaturation at a certain level. Few crystallization works with emphasis on particle size control have been carried out using gas as the antisolvent. Such processes are becoming more and more attractive because of economical and more importantly environmental objectives. Tai and Cheng (1998) investigate the saturation behavior of solutes, and the growth of crystals in the gas antisolvent crystallization (GAS) process. It is possible to tune and regulate the high pressures (needed to solubilize the gas antisolvent) or the gas antisolvent addition rate in order to achieve the desired particle size and/or size distribution. Muhrer et al. (2002) provide a simplified model of the GAS crystallization process that describes the vapour–liquid equilibrium as well as particle formation and growth. From their experiments, they identify different patterns of behavior in terms of antisolvent addition rates and they are also able, through these findings, to tune the product average particle size and size distribution to a large extent. They validate

the developed model and suggest its use as a tool to design and optimize GAS recrystallization processes. In spite of all the previous works, there still exists a big gap in the systematic control of the product crystal size and other properties in antisolvent mediated crystallization. Model-based optimization to determine optimal antisolvent crystallization strategies are poised to fill this gap and as far as the authors are aware, are not currently found in the literature. The first step towards model-based optimal strategies is model development. This is the concern of this paper, in which we present extensions to previously developed models (Abbas, 2003; Nowee et al., 2007), and address antisolvent effects and interactions. The model which is based on the population balance theory is used in optimization-based parameter estimation to arrive at solubility and kinetic sub-models. Finally, the developed model is experimentally validated and exploited in simulation studies with the aim to understand antisolvent crystallization behavior and the control of particle size. The use of the model for optimizing the feedrate and thus operation of antisolvent crystallization is presented in a sequel to this paper (Nowee et al., 2008). 2. Model development An appropriate approach for modeling the crystallization process is the statistical mechanical formulation proposed by Hulburt and Katz (1964). A population balance is employed and accounts for the evolution of crystal particles across temporal and size domains. For a crystallization system with crystal growth assumed to be nondispersed and independent of crystal size and where agglomeration and attrition are considered negligible, the population balance equation (PBE) simplifies to n(L, t) dV jn(L, t) jn(L, t) = −G +B− jt jL V dt

(1)

where n(L, t) is the number density of crystals, t is time, L is the characteristic crystal size, V is the suspension volume, G is the growth rate of the crystals and B is the nucleation rate, which is equal to zero for all sizes but the first. It follows that B = B0 (L)

(2)

where B0 is the nucleation rate and (L) is the Dirac delta function. A set of operating equations for the ternary solvent–solute–antisolvent crystallization is needed. We firstly define the mass fraction of each component as mi xi =  mi

(3)

and without loss of generality, xi and mi are the mass fraction and mass of the ith component respectively. For our case, the number of components, i, is 3 and so x1 , x2 and x3 are the primary solvent, solute and antisolvent mass fractions, respectively, and similarly m1 , m2 and m3 are the primary solvent, solute and antisolvent masses, respectively. Changes in suspension volume are evaluated by dV dVL dVS = + dt dt dt In Eq. (4), the liquid phase volume, VL , is given by  mi VL =



(4)

(5)

where  is the solution density, while the solute phase volume is given by  ∞ nL3 dL (6) VS = V v 0

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

where v is the crystal volume shape factor. Accordingly the numerator of Eq. (4) bears the antisolvent mass addition term as well as total solute mass in liquid phase. The mass balance, after Rawlings et al. (1993), for the solute mass accumulation, is incorporated as  ∞ dm2 = −3V c v G nL2 dL (7) dt 0

5459

The density–concentration correlation (Eq. (10)) is adopted from Galleguillos et al. (2003) who experimentally determined the set of constant coefficients (Ai ) for the NaCl–water–ethanol mixture which is used in this study as the model antisolvent crystallization system. The set of equations above representing the antisolvent crystallization model was solved in gPROMS environment (Process Systems Enterprise Ltd., UK) using backward finite discretization across the size domain, transforming the PBE into a set of ordinary differential equations. The solubility and kinetic sub-models for the NaCl–water–ethanol system are discussed in the next sections.

they develop kinetic models for mixed suspension mixed product removal (MSMPR) crystallization. Their growth rate model is size dependent while its kinetic parameters are functions of antisolvent to water weight ratio. In yet another work, the same authors (Mydlarz and Jones, 1991) consider agglomeration kinetics in addition to nucleation and growth. Later, Nyvlt (1992) mathematically investigated a kinetic model for nucleation rate. His model is a function of growth rate and supersaturation that depends on antisolvent addition rate. Nyvelt's conclusion suggests that the antisolvent addition should be slow in the early stages of the batch and increase proportionally with the crystal surface area. Granberg et al. (1999) substantiate Nyvelt's results in their study of paracetamol antisolvent crystallization. They present results of nucleation experiments showing that the higher the initial supersaturation is, the higher is the rate of supersaturation decay, suggesting that conditions of higher initial supersaturation generate larger total growing crystal surface area. Granberg et al. (1999) propose a relation between the surface energy and the solubility to describe the nucleation kinetics. Such kinetic model developments are rare in the literature and present two opportunities; the first is to repostulate the mechanisms behind such complex kinetic phenomena while the second is to develop systematic methods for the determination of the kinetic parameters. In this study, we start by incorporating a growth kinetic formulation proposed by Linnikov (2006) who neglected the effect of secondary nucleation and implemented seed addition. Seeding is not considered in the current study, and so a nucleation kinetic model is needed. Such a nucleation model is not available in the literature, however kinetics for cooling crystallization of sodium chloride as reported by Akal et al. (1986) are initially used. Together, the nucleation kinetics of Akal et al. (1986) and growth kinetics of Linnikov (2006) are referred to as `Literature' kinetics in the rest of the paper and are represented by the following equations:

3. Solubility model identification

B = kb C b MT

(12)

G = kg Ce−E/RT

(13)

where c is the crystal density. The concentration of the solute, C, is then defined as C=

m2 VL

(8)

¯ both being Solution density, , and volume mean size of particles, L, measured experimental variables are introduced to the model as:  n(L, t)L4 dL L¯ =  n(L, t)L3 dL

(9)

 = (A0 + A1 x2 + A2 x3 + A3 x2 x3 + A4 x2 x23 ) 2 × e(A5 (x3 /x1 )+A6 (x3 /x1 ) )

(10)

Knowledge of the equilibrium condition (solubility) of the crystallization system is crucial to the control of the particle size. We develop a model describing the solubility of the NaCl in the solvent–antisolvent mixture based on experimental solubility data (Yeo et al., 2007). This data corroborates with other data presented elsewhere for the same ternary system (Farelo et al., 2004). This model describing the saturated concentration dependence on solvent concentration in solute-free mixture is shown in Eq. (11): C∗ =

(1 − z3 ) ∗ C (K − z3 ) water

(11)

∗ is molalwhere C ∗ is the molality of solute in the mixture and Cwater ity of solute in pure water and has a value of 33.80 g/100 g of water ∗ at constant temperature 25 ◦ C. K, like Cwater is constant at constant temperature and is equal to 0.83094. z3 is the mass fraction of antisolvent in solute-free mixture.

where C is the supersaturation, kb and b are the nucleation rate parameters, MT is the magma density, kg and E are growth rate parameters. However, the growth model structure used in this study differs subtly from Eq. (13) and is proposed as G = kg S g

(14)

where the term kg e−E/RT in Eq. (13) is lumped into a new parameter kg . We can do this since the temperature of the process remains constant throughout. S is the supersaturation term which could be chosen to be absolute or relative supersaturation. Finally, the exponent g is introduced to the model. In this growth rate power law model structure, the parameters kg and g may have a constant value or may be defined as functions of antisolvent mass fraction in solutefree mixture (z3 ):

4. Kinetic model identification

kg = k0 + k1 z3 + k2 z32

(15)

Limited work has been presented in the literature on the modeling and identification of antisolvent crystallization kinetics. No comprehensive studies are available that address model-based identification in this field. The most relevant research addressing the antisolvent kinetics of nucleation and growth seems to be that of Mydlarz and Jones (1989, 1990). In Mydlarz and Jones (1989), an extensive experimental study is carried out with over 130 experiments conducted leading to the determination and rigorous analysis of the growth kinetics of potassium sulphate antisolvent crystallization in aqueous 2-propanol solutions. A power law model is shown to represent well the overall growth rate. In Mydlarz and Jones (1990),

g = g0 + g1 z3

(16)

The kinetic sub-model (Eqs. (12) and (14)–(16)) contains a set of parameters that needs to be identified via a parameter estimation step. Essentially the parameter estimation problem is an optimization one where the model is fitted to the experimental data. Once parameters are estimated, the identified model may be used for predictive simulation analysis where the various operational variables may be studied, and for optimization and control scheme development. Initial works in crystallization parameter estimation relied on population density data fitting (Chen and Larson, 1985). A common

5460

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

method presented throughout literature is the use of the analytical solution of the PBE along with experimental data from MSMPR crystallizers. This method is limited because (1) samples need to be drawn from the MSMPR at steady-state conditions which may be difficult to achieve, (2) estimates of parameters made from steadystate data would not be suitable for the description of the varying dynamics occurring during the batch time, and (3) nonlinearities in the size data may appear in some systems which do not obey McCabe's L law. In recent years, rigorous parameter estimation methods have been employed in crystallization kinetics identification but of all techniques, posing the estimation as an optimization problem has been regarded as most relevant (Rawlings et al., 1993; Mignon et al., 1996). Non-linear optimization techniques such as maximum likelihood are now available paving the way for more efficient estimation and with the advances in power of modern computers, optimization methods are becoming more prevalent in the literature. The mathematical background behind these techniques is intricate and there are several text available detailing the mathematical derivations (Bard, 1974; Beck and Arnold, 1977; Leonard and Hsu, 1999) for the reader to refer to. The crystallization model described in Sections 2–4 consists of a coupled set of non-linear integro-differential equations and can ¯ be represented by Y(t) = [L(t), (t)] being the differential output variables—mean size and solution density, respectively, U(t) = [F(t)] being the time-varying input variable which is the feed rate of antisolvent, and  being the vector of parameters to be identified which in this case is the set [kb , kg , b, g] when considering kg and g to hold constant values, or the set [kb , k0 , k1 , k2 , b, g0 , g1 ] when considering kg and g to be functions of antisolvent mass fraction in solute-free media. We can then represent the model as follows: f (Y(t), U(t), ) = 0

5. Experimental For the estimation step and model refinement, three experiments were conducted under antisolvent feeding conditions. Fig. 1 shows the schematic of the experimental apparatus, instrumentation and control system used. In all experiments, purified water by a Milli-Q system was used. The purities of sodium chloride (NaCl) salt (Merck) and ethanol (Merck) used in the experiments were 99.5% and 99.9%, respectively. Ethanol was added to the aqueous NaCl solution using a calibrated digital dosing pump (Grundfos, Denmark). The ethanol addition profile was implemented and controlled from within a distributed control system (DCS) environment (Honeywell, USA). Temperature was controlled, at 25 ◦ C for all experiments, using a Pt100 thermocouple connected to heating/cooling circulator (Lauda, Germany). Particle chord length was measured in situ every 2 s using focused beam reflectance measurement (FBRM) probe (MettlerToledo Lasentec Products, USA). Infrequent samples were removed iso-kinetically from the crystallizer for offline particle size and density measurements. Particle size was measured using Mastersizer 2000 particle size analyzer (Malvern Instruments, UK) while density was measured using a density meter (Anton Paar DMA 5000, Austria). Samples were filtered using a syringe microfilter prior to density measurements. The concentration of the solution was inferred from the density reading using the correlation of Galleguillos et al. (2003). The FBRM signal was converted to 4–20 mA signals using an eight-channel analog output PCI card inserted in the FBRM

(17)

Eq. (17) is accompanied by a set of initial conditions. Measured variables from experiments are essential in the parameter estimation ˆ ˆ exercise and are denoted by the vector Y(t) = [L(t), ˆ (t)]. The maximum likelihood criterion that describes the highest probability of the model predicting the real data is given by the following objective function: ⎛

⎤⎞ ij ⎡ i    ˆ  −Y Y M 1 ⎣ln( 2 )+ ijk ijk ⎦ ⎟ (k, )= ln(2)+ min ⎜ ⎝ ⎠ ijk 2 2 k,

2ijk i=1 j=1 k=1

(18)

where M is the total number of measurements taken during all experiments, , i and ij are, respectively, the number of experiments, the number of variables measured in the ith experiment and the number of measurements of the jth variable in the ith experiment,

2ijk is the variance of the kth measurement of variable j in experiment i. When the error structure of the data is known, 2ijk can be described by any of a number of models available (Process System Enterprise Ltd., 2001). In previous works, we successfully applied this method to the identification of ammonium sulfate kinetics for crystallizations under cooling mode of operation (Abbas and Romagnoli, 2007; Nowee et al., 2007). The same method is applied here with differences being modifications in the current model as well as in the experimental conditions to describe antisolvent crystallization under constant temperature. Here again the implementation of the maximum likelihood optimization is carried out using the gEST facility in gPROMS. The kinetic parameter set that provides the highest probability of the model predicting the real data is identified.

Fig. 1. Experimental setup. A: crystallization vessel, B: temperature control system, C: Pt100 thermocouple, D: ethanol addition line, E: pump, F: ethanol reservoir, G: FBRM probe, H: analog output card, I: FBRM control computer, J: DCS station, K: controller I/O terminals, L: DCS server, M: feed profile to pump.

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

computer. This was then connected to the DCS for data monitoring and archiving.

5461

• Model 4: Growth kinetic is a function of relative supersaturation (C/C ∗ ) and the parameters of growth are functions (Eqs. (15) and (16)) of antisolvent mass fraction in solute-free mixture.

6. Results and discussion 6.1. Model identification and kinetic parameter analysis A series of experiments were conducted with the aim to generate data to be used in the estimation of the parameters of nucleation and growth kinetic models. Either of the parameter sets 1 = [kb , kg , b, g] and 2 = [kb , k0 , k1 , k2 , b, g0 , g1 ] was determined during an estimation calculation depending on which of the kinetic model structures presented in Section 4 was being used. This was possible from data gathered from four experiments under three constant antisolvent feeding profiles—one experiment at a feedrate of 194 ml h−1 , two at a medium rate of 98.4 ml h−1 and one at a lower rate of 50 ml h−1 . We took a further step here beyond the literature kinetics and investigated a number of structural changes to the growth kinetic model (Eq. (14)). The following four cases were examined: • Model 1: The literature kinetic model, where, growth kinetic model is a function of absolute supersaturation (C) and the parameters of growth (kg and g) are constant. • Model 2: Growth kinetic model is a function of relative supersaturation (C/C ∗ ) and the parameters of growth (kg and g) are constant. • Model 3: Growth kinetic model is a function of absolute supersaturation (C) and the parameters of growth are functions (Eqs. (10) and (11)) of antisolvent mass fraction in solute-free mixture.

Temporal mean size data from the Mastersizer and concentration data inferred from the density meter were necessary in this estimation step. The results of each of the four experiments are displayed in Figs. 2–5 against the data from the identified kinetic models of each of the four models above. The estimation activities for all models were executed successfully and parameter values obtained are displayed in Table 1. The goodness of the fits varies across the four models, yet all can be deemed satisfactory. This is more so for the concentration profile data. As for the mean size profiles, it can be visually observed from the overlay plots that Model 4 offers, in general a better fit, closely followed by Model 3. This is corroborated by the estimation residual values (Table 1) showing a relatively lower residual error in the case of Model 4 in comparison to the other three models. Model 4 and its corresponding estimated parameters is favored and is adopted in a subsequent experimental validation exercise presented in Section 6.2. Both Models 3 and 4, in contrast to Models 1 and 2, carry growth rate coefficients with functionality of antisolvent mass fraction. This suggests (since these two models offer better fits) the significance of including antisolvent mass fraction in the growth kinetics. Models 3 and 4 differ in their formulation of the supersaturation term, but no profound effect of this difference can be observed in the fits, except that in the early stages of the process, Model 4 mean size profile always rises more rapidly. This suggests that the two supersaturation formulations employed differ with respect to their effect on

140

300 Model 1 Model 2 Model 3 Model 4 Experiment

250 100 80 60 Model 1 Model 2 Model 3 Model 4 Experiment

40 20

Csolute (gr/Ltr)

Mean Size (micron)

120 200 150 100 50

0

0 0

10

20

30 40 Time (min)

50

60

0

10

20

30

40

50

60

Time (min)

Fig. 2. Results of Experiment 1 at antisolvent addition of 194 ml h−1 . Overlay plots of experimental mean size (left) and solution concentration (right) against model predictions using the different kinetic models.

140

300 Model 1 Model 2 Model 3 Model 4 Experiment

250

100 80 60 Model 1 Model 2 Model 3 Model 4 Experiment

40 20 0

Csolute (gr/Ltr)

Mean Size (micron)

120

200 150 100 50 0

0

20

40

60 80 Time (min)

100

120

0

20

40

60 80 Time (min)

100

120

Fig. 3. Results of Experiment 2 at antisolvent addition of 98.4 ml h−1 . Overlay plots of experimental mean size (left) and solution concentration (right) against model predictions using the different kinetic models.

5462

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

140

300 Model 1 Model 2 Model 3 Model 4 Experiment

250

100 80 60 Model 1 Model 2 Model 3 Model 4 Experiment

40 20

Csolute (gr/Ltr)

Mean Size (micron)

120

200 150 100 50

0

0 0

20

40

60 80 Time (min)

100

120

0

20

40

60 80 Time (min)

100

120

Fig. 4. Results of Experiment 3 at antisolvent addition of 98.4 ml h−1 . Overlay plots of experimental mean size (left) and solution concentration (right) against model predictions using the different kinetic models.

180

300 Model 1 Model 2 Model 3 Model 4 Experiment

250

140 120 100 80 60

Model 1 Model 2 Model 3 Model 4 Experiment

40 20

Csolute (gr/Ltr)

Mean Size (micron)

160

200 150 100 50

0

0 0

50

100 150 Time (min)

200

0

50

100

150

200

Time (min)

Fig. 5. Results of Experiment 4 at antisolvent addition of 50 ml h−1 . Overlay plots of experimental mean size (left) and solution concentration (right) against model predictions using the different kinetic models.

Table 1 Values of the estimated parameters for the different growth kinetic model structures along with the estimation residual values for each model Model 1

Model 2

Model 3

Model 4

bn kb g kg g0 g1 k0 k1 k2

3.464 2.11E + 07 2.576 3.10E − 05

1.494 2.12E + 07 0.2302 7.67E − 07

1.937 9.66E + 06

0.9655 9.68E + 06

1.051 4.875 1.86E − 05 0 0

0.8006 0.7379 9.88E − 05 −0.00241 0.01681

Residual value

61971

2681

100

2113

Volume Percentage

Parameter

Model 1 Model 2 Model 3 Model 4 Experimental

10 8 6 4 2 0

nucleation activities in the early stages of the process. The exact nature of this is worthy of further investigation. The weaker models, i.e. Models 1 and 2, deviate more from the experimental data. Model 1 almost always overpredicts the mean size in the later phases of the process, but there is some consistency in this, that is the absolute overprediction decreases for decreasing constant feedrate of the antisolvent. Similar consistency is observed for Model 2, but this time there is underprediction at the higher constant antisolvent feedrate and this absolute underprediction moves to become an overprediction for the lower feedrate experimental data. In any case, the exclusion of the antisolvent mass fraction from the growth kinetics as is done in Models 1 and 2, reduces the goodness of their fits. Out of all the models, the fit of Model 1 behaves the worst at the early times of the process in the form of severe underprediction to the data

0.1

1

10 Size (μ)

100

1000

Fig. 6. End of run product PSD overlay plot for Experiment 1 at antisolvent addition rate of 194 ml h−1 against model predictions using the different kinetic models.

under all feedrates tested. Further a sigmoidal shape to this model is observed in the beginning which may be explained by the higher orders estimated for the kinetic exponents (see Table 1). The comparison between the model predictions and experimental data is extended to the end of process crystal size distribution (CSD) (Figs. 6–9). The experimentally obtained CSD is always broader. It is difficult to conclude from this analysis which of the four models does a better job at predicting the CSD since all models are close

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

5463

5000

8

4000

3000

6

2000

4 1000

2

0

0 0.1

1

10

100

1000

Size (μ) Fig. 7. End of run product PSD overlay plot for Experiment 2 at antisolvent addition rate of 98.4 ml h−1 against model predictions using the different kinetic models.

0

50

100 150 Time (min)

200

250

Fig. 10. FBRM counts in the 1–5 m range over time for the four antisolvent crystallization experiments at the following antisolvent feedrates (in ml h−1 ); Experiment 1 at 194 ml h−1 , Experiments 2 and 3 at 98.4 ml h−1 and Experiment 4 at 50 ml h−1 .

Model 1 Model 2 Model 3 Model 4 Experimental

10

Volume Percentage

Experiment 1 Experiment 2 Experiment 3 Experiment 4

#/s

Volume Percentage

10

Model 1 Model 2 Model 3 Model 4 Experimental

8 6 4 2 0 0.1

1

10 Size (μ)

100

1000

Fig. 8. End of run product PSD overlay plot for Experiment 3 at antisolvent addition rate of 98.4 ml h−1 against model predictions using the different kinetic models.

12 Model 1 Model 2 Model 3 Model 4 Experimental

Volume Percentage

10 8 6 4 2 0 0.1

1

10 Size (μ)

100

1000

Fig. 9. End of run product PSD overlay plot for Experiment 4 at antisolvent addition rate of 50 ml h−1 against model predictions using the different kinetic models.

Fig. 11. End of run product SEM micrograph from Experiment 1 at antisolvent addition rate of 194 ml h−1 .

to the experimental CSD. An exception is Model 3, which positively shows some attempt to predict the lower tail of the CSD. More work is needed to further improve the prediction of CSD which remains an open challenge due to the multi-dimensionality of the problem when considering the distribution and not only the mean size profile. It may be beneficial to revisit the PBE as well as undertake a detailed examination of the prevalent kinetics, since for instance neglecting agglomeration and attrition kinetics could be a source of this CSD model–experiment mismatch. In this same regard, well-designed and controlled experiment as well as the selection of strategic manipulations and combinations thereof could lead to better handling over the crystallization process and consequently the CSD. The FBRM chord-length counts for the size range 1–5 m for the four experiments are shown in Fig. 10. There is a clear trend in the data which is consistent with previous knowledge. Under the higher constant feedrate, there exists the highest rate of generation of nuclei as is observed in Experiment 1 (Feedrate = 194 ml h−1 ). The generation of crystal nuclei decreases as the feedrate decreases. This is best explained in terms of supersaturation. At higher feedrates of antisolvent, there are higher levels of supersaturation generation, and so the driving force is larger for nucleation. One consistent trend is that the FBRM chord-length count profiles go through a maximum

5464

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

in the very first minutes after antisolvent feeding starts. This maximum is achieved sharply for all feedrates. The counts start to decline after the maximum and keep doing so till the end of the process. The

supersaturation profile generated under the different feeding rates has a big role to play here. Supersaturation was not directly measured here but is expected to also be going through a maximum early on in the process. More nuclei form in the beginning due to primary nucleation. As supersaturation is depleted steadily, primary nucleation rate reduces and the nuclei initially formed grow out of the size range 1–5 m. This and the fact that secondary nucleation becomes more prevalent as the mass of crystals increases in the process are responsible for the decline in the nuclei counts. One important note here is that our model currently does not consider secondary nucleation and so has the drawback that it is incapable of predicting the process-time dependent nucleation events (primary nucleation dominant in the beginning while secondary nucleation dominant later). As a consequence all the models tested in the estimation activity fail to predict well the rise (or hump) in the mean size profile early in the process. All models tested underpredict this and this suggests the importance of including detailed nucleation kinetics to describe primary and secondary nucleation events. Having said that, we also point to the relationship between the FBRM signal and the size of the particles detected. Bigger particles are easier seen by the laser beam of the FBRM, i.e. bigger particles have a higher probability of being detected by the laser. This fact could contribute to the level of decline of the counts observed in Fig. 10 as the population of bigger particles increases at the later stages of the process. Microscopic analysis was conducted on the end of process samples under the different constant feedrates (Figs. 11–13). At the higher feedrate, the population of crystals is confirmed to be dominated by smaller sized crystals corroborating the above discussion. Bigger sized crystals are observed in images corresponding to lower feedrates, but small crystals still exist, suggesting that nucleation persists throughout the process, but is less under the condition of lower feedrate. This also matches with the experimental data for the CSD of Figs. 6–9, which show in essence bimodal distributions where each mode represents the small and the large population, respectively. The bigger particle mode is seen to dominate because the CSDs are presented as volume-based.

Fig. 12. End of run product SEM micrograph from Experiment 3 at antisolvent addition rate of 98.4 ml h−1 .

6.2. Experimental model validation under a new antisolvent feedrate profile In a further exercise, we undertook the testing of the antisolvent crystallization model with the chosen kinetic model (Model 4). A totally new and non-constant feedrate of antisolvent was defined arbitrarily and implemented. Fig. 14 shows this feedrate profile along with the mean size profiles of experiment versus model prediction. A good agreement is found, but yet again the model underpredicts the `hump' region. As discussed earlier, this is attributed to the model's lack of a description for secondary nucleation events. The experimental and model predicted CSDs are also shown in Fig. 14.

6 5 Mean Size (Expt) Mean Size (Model) Ethanol Feed Rate

4 3 2 1

0

0 50 100 150 200 250 300 Time (min)

10 Volume Percentage

180 160 140 120 100 80 60 40 20 0

Feed Rate (ml/min)

Mean Size (micron)

Fig. 13. End of run product SEM micrograph from Experiment 4 at antisolvent addition rate of 50 ml h−1 .

Model 4 Experiment

8 6 4 2 0 0.01

0.1

1 10 Size (μ)

100

1000

Fig. 14. Experimental validation results showing the model prediction versus experimental data for mean size profile (left) and CSD (right) under a new antisolvent feedrate profile.

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

1 ml.min 5 ml/min Increasing 1-4 ml/min Decreasing 4-1 ml/min

Feed Rate (ml/min)

4

3

2

1

0.025 FAntisolvent = 1 ml/min FAntisolvent = 5 ml/min

Relative Supersaturation σ

5

5465

0

FAntisolvent = 1-4 ml/min

0.020

FAntisolvent = 4-1 ml/min

0.015

0.010

0.005

0.000 50

0

100 150 Time (min)

200

250

300

50

0

200

100 150 Time (min)

250

140 5e+10 120 4e+10 m0 (#/m3)

Mean Size (μ)

100 80 60 40

FAntisolvent = 1 ml/min

2e+10 FAntisolvent = 1 ml/min

1e+10

FAntisolvent = 5 ml/min

20

3e+10

FAntisolvent = 5 ml/min

FAntisolvent = 1-4 ml/min

FAntisolvent = 1-4 ml/min

FAntisolvent = 4-1 ml/min

0 0

50

100 150 Time (min)

200

FAntisolvent = 4-1 ml/min

0 250

0

50

100 150 Time (min)

200

250

Fig. 15. Simulated antisolvent crystallization results. The relative supersaturation, the mean size and the zeroth moment profiles under different antisolvent feedrate profiles.

The predicted CSD falls in the same size range but has some room for improvement. Overall, this is a satisfactory validation outcome considering the extensive variation of the new antisolvent feed profile.

12 10

6.3. Model simulation analysis

8 Volume %

The validated model is now used in a simulation analysis to elucidate the behavior of the antisolvent crystallization process. Four simulated experiments were conducted with a fixed total antisolvent addition of 250 ml, each under a different antisolvent feedrate profile. A low feedrate of 1 ml min−1 was used in the first simulation and a higher rate of 5 ml min−1 was used in a second simulation. The third simulated experiment had the feedrate increased from the low (1 ml min−1 ) feedrate to the high rate in stepwise manner. Finally, the fourth and last experiment had a feedrate that decreased in exactly the opposite direction to the third experiment. The results of these simulations showing key variables are shown in Figs. 15 and 16. For all feedrates there is a high peak of supersaturation at the early stages of the process when primary nucleation occurs, though there are large differences in the magnitude. As expected for high feedrates (4 and 5 ml min−1 ) this peak is much higher resulting in excessive nuclei formation. This is tracked by the zeroth moment of particle size distribution representing the number count of particles per unit volume. In addition, in this figure it is worth pointing out that the slope of growth in the number of particles also highlights the high rate of nucleation for cases 2 and 4 at the early stages. It

FAntisolvent = 1 ml/min FAntisolvent = 5 ml/min FAntisolvent = 1-4 ml/min FAntisolvent = 4-1 ml/min

6 4 2 0 1

10

100

1000

Size (μ) Fig. 16. Simulated antisolvent crystallization results. The end of process crystal size distributions under different antisolvent feedrate profiles.

should be noted that despite the different antisolvent feeding profiles for all four cases, the feedrate at the early stages of the process is critical enough that it defines to a large extent the outcomes of the final particle size characteristics. This fact can be deduced from the mean size and zeroth moment profiles of Fig. 15. These profiles

5466

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

show that for cases 1 and 3, both having 1 ml min−1 as the initial step of feeding, the final mean size and number counts are almost identical. Similarly, for cases 2 and 4, both having similar high initial step of feeding, have very identical values of final mean size and number count of particles. The same behavior is observed in the CSDs. Fig. 16 also demonstrates this (note that final CSDs of cases 1 and 3 coincide with each other).

7. Conclusions In this paper, we have presented a detailed model for antisolvent crystallization. Optimization based parameter estimated paved the way for the systematic estimation of unknown kinetic parameters. The structure of the kinetic model was investigated and it was found that a growth rate relation that is a function of relative supersaturation results in sharper mean size increases in the early stages of the feeding. It was also found that a growth rate model that considers its parameters to be functions of antisolvent mass fraction is superior to a model without such functionality. The absence of secondary nucleation relation was recognized to be a cause of the model–experiment mismatch in the early stages of the process. The model developed was validated with satisfactory outcomes under a new and totally different profile of antisolvent feeding. The model was subsequently shown to be a useful tool in analyzing, via simulation, the effect of antisolvent behavior of the crystallization process. The model is thus proposed for use in model-based optimal profile generation as well as model-based control applications.

Notation Ai

constant coefficients

B

nucleation rate, no. (m3 s)−1

B0 b

nucleation rate of zero-sized crystals, no. (m3 s)−1 nucleation rate order

C

solution concentration, g (100 g soln)−1

C∗

saturated solution concentration, g (100 g soln)−1

∗ Cwater

molality of solute in pure water, g (100 g water)−1

C

E F G g g0 , g1 k0 , k1 , k2 kb kg kg L Lˆ L¯ M MT mi n S t U V VL VS

supersaturation, g (100 g soln)−1 growth rate coefficient feed rate of antisolvent, ml min−1 growth rate, m s−1 growth rate index constant coefficients constant coefficients nucleation rate coefficient, no. (kg s)−1 lumped growth rate coefficient, m s−1 growth rate coefficient, m s−1 crystal size, m measured mean size, m mean size, m total number of measurements magma density, kg m−3 massof ith component, kg particle number density, no. (m3 m)−1 relative supersaturation time, s time varying input variables crystallizer suspension volume, m3 liquid phase volume, m3 solid phase volume, m3

xi Y Yˆ z3

mass fraction of ith component differential output variables vector of measured variables mass fraction of antisolvent in solute-free mixture

Greek letters

 v i ij    ˆ c

2ijk 

number of experiments particle volume shape factor number of variables measured in the ith experiment number of measurements of the jth variable in the ith experiment dirac delta function vector of kinetic parameters solution density, kg m−3 measure solution density, kg m−3 crystal density, kg m−3 variance of the kth measurement of jth variable in the ith experiment objective function

Acknowledgment Experimental work was carried out at Nanyang Technological University (Singapore) with the support of the School of Chemical and Biomedical Engineering. References Abbas, A., 2003. A model-based framework for advanced operation of crystallisation processes: a pilot-scale study. Ph.D. Thesis, University of Sydney, 2003. Abbas, A., Romagnoli, J.A., 2007. Multiscale modeling, simulation and validation of batch cooling crystallization. Separation and Purification Technology 53 (2), 153 –163. Akal, M.M., Zakaria, M., Ebrahim, A., Nassar, M.M., 1986. Secondary nucleation rate of sodium chloride under different stirring conditions. Journal of Crystal Growth 78 (3), 528–532. Bard, Y., 1974. Nonlinear Parameter Estimation. Academic, New York. Beck, J.V., Arnold, K.J., 1977. Parameter Estimation in Engineering and Science. Wiley, New York. Budz, J., Karplirskl, P.H., Mydlarz, J., Nyvlt, J., 1986. Salting-out precipitation of cocarboxylase hydrochloride from aqueous solution by addition of acetone. Industrial Engineering and Chemical Product Research and Development 25 (4), 657–664. Chen, M.R., Larson, M.A., 1985. Crystallization kinetics of calcium nitrate tetrahydrate from MSMPR crystallizer. Chemical Engineering Science 40 (7), 1287–1294. Doki, N., Kubota, N., Yokota, M., Kimura, S., Sasaki, S., 2002. Production of sodium chloride crystals of uni-modal size distribution by batch dilution crystallization. Journal of Chemical Engineering of Japan 35 (11), 1099–1104. Farelo, F., Lopes, A.M.C., Ferra, M.I., 2004. Solubilities of sodium chloride and potassium chloride in water + ethanol mixtures from 298 to 323 K. Journal of Chemical and Engineering Data 49 (6), 1782–1788. Galleguillos, H.R., Taboada, M.E., Graber, T.A., Bolado, S., 2003. Compositions, densities, and refractive indices of potassium chloride plus ethanol plus water and sodium chloride plus ethanol plus water solutions at 298.15 and 313.15 K. Journal of Chemical and Engineering Data 48 (2), 405–410. Granberg, R.A., Bloch, D.G., Rasmuson, A.C., 1999. Crystallization of paracetamol in acetone–water mixtures. Journal of Crystal Growth 198–199 (Part 2), 1287– 1293. Ho-Gutierrez, I.V., Cheluget, E.L., Vera, J.H., Weber, M.E., 1994. Liquid–liquid equilibrium of aqueous mixtures of poly(ethylene glycol) with Na2 SO4 or NaCl. Journal of Chemical and Engineering Data 39 (2), 245–248. Holmback, X., Rasmuson, A.C., 1999. Size and morphology of benzoic acid crystals produced by drowning-out crystallization. Journal of Crystal Growth 198–199 (Part 1), 780–788. Hulburt, H.M., Katz, S., 1964. Problems in particle technology. A statistical mechanical formulation. Chemical Engineering Science 19, 555–574. Jones, A.G., Mullin, J.W., 1974. Programmed cooling crystallization of potassium sulphate solutions. Chemical Engineering Science 29 (1), 105–118. Kaneko, S., Yamagami, Y., Tochihara, H., Hirasawa, I., 2002. Effect of supersaturation on crystal size and number of crystals produced in antisolvent crystallization. Journal of Chemical Engineering of Japan 35 (11), 1219–1223. Kim, Y., Haam, S., Shul, Y.G., Kim, W.S., Jung, J.K., Eun, H.C., Koo, K.K., 2003. Pseudopolymorphic crystallization of L-ornithine-L-aspartate by drowning out. Industrial & Engineering Chemistry Research 42 (4), 883–889. Kitamura, M., Nakamura, K., 2002. Dependence of polymorphic transformation on antisolvent composition and crystallization behavior of thiazole-derivative pharmaceutical. Journal of Chemical Engineering of Japan 35 (11), 1116–1122.

S. Mostafa Nowee et al. / Chemical Engineering Science 63 (2008) 5457 -- 5467

Kitamura, M., Sugimoto, M., 2003. Antisolvent crystallization and transformation of thiazole-derivative polymorphs—I: effect of addition rate and initial concentrations. Journal of Crystal Growth 257 (1–2), 177–184. Leonard, T., Hsu, J.S.J., 1999. Bayesian Methods: An Analysis for Statisticians and Interdisciplinary Researchers. Cambridge University Press, Cambridge, U.K. Linnikov, O.D., 2006. Spontaneous crystallization of sodium chloride from aqueous–ethanol solutions. Part 1. Kinetics and mechanism of the crystallization process. Crystal Research and Technology 41 (1), 10–17. Mignon, D., Manth, T., Offermann, H., 1996. Kinetic modelling of batch precipitation reactions. Chemical Engineering Science 51 (11), 2565–2570. Muhrer, G., Lin, C., Mazzotti, M., 2002. Modeling the gas antisolvent recrystallization process. Industrial & Engineering Chemistry Research 41 (15), 3566– 3579. Mullin, J.W., Nyvlt, J., 1971. Programmed cooling of batch crystallizers. Chemical Engineering Science 26 (3), 369–377. Mydlarz, J., Jones, A.G., 1989. Growth and dissolution kinetics of potassium-sulfate crystals in aqueous 2-propanol solutions. Chemical Engineering Science 44 (6), 1391–1402. Mydlarz, J., Jones, A.G., 1990. Potassium sulfate water–alcohols systems: composition and density of saturated solutions. Journal of Chemical and Engineering Data 35 (2), 214–216. Mydlarz, J., Jones, A.G., 1991. Crystallization and agglomeration kinetics during the batch drowning-out precipitation of potash alum with aqueous acetone. Powder Technology 65 (1–3), 187–194. Nonoyama, N., Hanaki, K., Yabuki, Y., 2006. Constant supersaturation control of antisolvent-addition batch crystallization. Organic Process Research & Development 10 (4), 727–732. Nowee, S.M., Abbas, A., Romagnoli, J.A., 2007. Optimization in seeded cooling crystallization: a parameter estimation and dynamic optimization study. Chemical Engineering and Processing 46 (11), 1096–1106. Nowee, S.M., Abbas, A., Romagnoli, J.A., 2008. Model-based optimal strategies for controlling particle size in antisolvent crystallization operations. Crystal Growth & Design 8 (8), 2698–2706. Nyvlt, J., 1992. Batch salting-out crystallization. Chemical Engineering and Processing 31 (1), 39–42.

5467

Oosterhof, H., Geertman, R.M., Witkamp, G.J., van Rosmalen, G.M., 1999a. The growth of sodium nitrate from mixtures of water and isopropoxyethanol. Journal of Crystal Growth 198–199 (Part 1), 754–759. Oosterhof, H., Witkamp, G.J., van Rosmalen, G.M., 1999b. Some antisolvents for crystallisation of sodium carbonate. Fluid Phase Equilibria 155 (2), 219–227. Oosterhof, H., Witkamp, G.J., Van Rosmalen, G.M., 2001. Antisolvent crystallization of anhydrous sodium carbonate at atmospheric conditions. A.I.Ch.E. Journal 47 (3), 602–608. Process System Enterprise Ltd., 2001. gPROMS Advanced User Guide, UK. Rawlings, J.B., Miller, S.M., Witkowski, W.R., 1993. Model identification and control of solution crystallization processes: a review. Industrial & Engineering Chemistry Research 32 (7), 1275–1296. Taboada, M.E., Graber, T.A., Asenjo, J.A., Andrews, B.A., 2000. Drowning-out crystallisation of sodium sulphate using aqueous two-phase systems. Journal of Chromatography, B: Biomedical Sciences and Applications 743 (1+2), 101–105. Tai, C.Y., Cheng, C.S., 1998. Supersaturation and crystal growth in gas antisolvent crystallization. Journal of Crystal Growth 183 (4), 622–628. Takiyama, H., Otsuhata, T., Matsuoka, M., 1998. Morphology of NaCl crystals in drowning-out precipitation operation. Chemical Engineering Research and Design 76 (A7), 809–814. Toth, J., Kardos-Fodor, A., Halasz-Peterfi, S., 2005. The formation of fine particles by salting-out precipitation. Chemical Engineering and Processing; Pneumatic Conveying and Handling of Particulate Solids 44 (2), 193–200. Yeo, P., Abbas, A., Foddi, O., Romagnoli, J.A., 2007. In-situ measurement of NaCl solubility in mixed solvent system using laser-based measurement, submitted for publication. Yu, Z.Q., Tan, R.B.H., Chow, P.S., 2005. Effects of operating conditions on agglomeration and habit of paracetamol crystals in antisolvent crystallization. Journal of Crystal Growth 279 (3–4), 477–488. Yu, Z.Q., Chow, P.S., Tan, R.B.H., 2006. Seeding and constant-supersaturation control by ATR-FTIR in antisolvent crystallization. Organic Process Research & Development 10 (4), 717–722. Zhou, G.X., Fujiwara, M., Woo, X.Y., Rusli, E., Tung, H.H., Starbuck, C., Davidson, O., Ge, Z., Braatz, R.B., 2006. Direct design of pharmaceutical antisolvent crystallization through concentration control. Crystal Growth & Design 6 (4), 892–898.