Article pubs.acs.org/IECR
Combined Cooling and Antisolvent Crystallization in Continuous Mixed Suspension, Mixed Product Removal Cascade Crystallizers: Steady-State and Startup Optimization Yang Yang and Zoltan K. Nagy* School of Chemical Engineering, Purdue University, West Lafayette, Indiana 47907, United States ABSTRACT: A continuous multistage mixed suspension, mixed product removal (MSMPR) cascade crystallization system has been modeled and optimized for both the steady-state and the startup behaviors when both cooling and antisolvent addition are applied as the crystallization driving forces. The optimal steady-state operating profiles of MSMPR cascade crystallizers have been obtained and found to be similar to the optimal batch trajectory. In addition, the minimum startup duration time is used as the objective function to optimize the startup behavior. Methods, including changing the initial solution compositions, applying dynamic antisolvent profiles, applying dynamic temperature profiles, and applying dynamic combined antisolvent and temperature profiles, are discussed and used to formulate the optimization problems. The optimal results indicate that for a cascade of two MSMPR crystallizers, it is possible to achieve a reduction of approximately 50% of the startup duration time and amount of waste. A model-free or direct design approach via closed-loop control to improve the startup performance is also proposed.
1. INTRODUCTION Crystallization is an essential step in the separation and purification of fine chemical and pharmaceutical products.1 Batch crystallization is widely used in industry because of its great operation flexibility, easy design, and well-established standard operation procedure (SOP).2 However, in batch crystallization, the operating conditions change with time, resulting in crystal properties that are difficult to control and inconsistent from batch to batch.3 On the other hand, continuous crystallization is promising because of its enhanced reproducibility, reduced operating cost, simplified scale up, and improved efficiency of process and equipment.4 As a result, the development of robust continuous crystallization systems is of great interest. The most popular approach of continuous crystallization system is to use a mixed suspension, mixed product removal (MSMPR) crystallizer. It is typically treated as an ideally wellmixed vessel. During the operation, the feed solution continuously enters into the vessel and the crystal slurry is continuously withdrawn from the vessel. The output slurry is assumed to have exactly the same composition as the vessel content.4 It is operated at steady-state, which can be close to batch equilibrium conditions by increasing the number of stages. Such a MSMPR cascade system shows a possibility to easily convert batch knowledge to a continuous system.5 Similar to the batch crystallization, either cooling or antisolvent addition can be used as the driving force for crystallization in MSMPR cascade crystallizers. Alternatively, one could consider combined cooling and antisolvent crystallization (CCAC) in the MSMPR cascade crystallizers because of the higher yield and the larger control freedom on multiple crystal properties like crystal shape and polymorphic form of CCAC when compared with cooling or antisolvent addition only crystallization.6 © 2015 American Chemical Society
The control of cooling or antisolvent addition only crystallization in a batch system has been studied extensively in the literature,7−11 whereas the control of CCAC was not of great interest until Nagy et al.12 first studied the control problem of CCAC via both model-free and model-based approaches in a batch system systematically. They proved that better crystal quality (e.g., larger crystals) could be obtained under the CCAC approach than under the cooling or antisolvent only approach. CCAC was also applied to produce consistent crystal size distribution (CSD) or to improve the yield in batch crystallization systems.13,14 On the other hand, both academics and industry professionals have been working hard on the development of continuous MSMPR cascade systems in recent years. Myerson, Wong and co-workers developed a single-stage MSMPR crystallizer with recycling which could obtain high yield and purity.5 Multistage continuous MSMPR crystallizers were also developed by the same group.3,4 These cascade systems were proven experimentally to be able to improve yield and purity for crystallization of compounds like cyclosporine and aliskiren hemifumarate. In addition, Myerson, Zhang, and co-workers experimentally demonstrated the capability of performing CCAC in MSMPR cascade crystallizers to crystallize an organic compound. 15 Vetter et al. modeled MSMPR cascade crystallizers using population balance equation (PBE) models and obtained the attainable regions of particle sizes of aspirin under various residence times.16 In addition, recently Su et al. developed a rigorous mathematical model and a concentration control strategy for continuous MSMPR cascade crystallizers.17 Received: Revised: Accepted: Published: 5673
August 29, 2014 February 8, 2015 May 8, 2015 May 8, 2015 DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research It was extended to facilitate the convenient design of startup procedure and online process control. The startup problem is crucial in continuous crystallization because the crystal properties are inconsistent before the continuous system reaches the steady-state. As a result, the time and the products are wasted during startup. Hence, feasible methods to shorten the startup duration time and reduce wasted material in the continuous crystallization systems are highly desirable. However, currently there are very few works reporting the control and optimization of the startup behavior of continuous crystallization system such as the MSMPR cascade crystallizers.17 In this work, both the steady-state and the startup behaviors of CCAC in multistage MSMPR cascade crystallizers are evaluated and optimized through a model-based framework. The PBE models are used and solved using the method of moments. First, the optimal steady-state operating trajectories are obtained and compared with the optimal batch trajectory. Then, different methods to shorten the startup duration and reduce the amount of waste are discussed. The methods include changing the initial solution compositions, applying dynamic antisolvent profiles, applying dynamic temperature profiles, and applying dynamic antisolvent and temperature combined profiles. Lastly, a concept of implementing a model-free approach to optimize the startup behavior is proposed. Aspirin is used as the model compound in the case study, and ethanol and water are used as solvent and antisolvent, respectively.
dμj(i) dt
τi = αi =
dμj(i) dt
L jni(L , t ) dL ,
dμ0(i) dt
dt
(i − 1) (i − 1) (i) (i) (i) = Fout μj − Fout μj + VG i ijμj − 1 ,
(6)
=
1 1 μ(i − 1) − μj(i) + jGiμj(−i)1 , αi − 1τi − 1 j τi
j≥1
(7)
=
1 1 μ0(i − 1) − μ0(i) + Bi αi − 1τi − 1 τi
(8) −1
where Bi represents the nucleation rate (no./cm min ) at stage i. The mass balance equation for the solute at stage i is derived as (i)
Vi
dμ dci (i − 1) (i) = Fout ci − 1 − Fout ci − kvρc Vi 3 dt dt
(9)
where ci is the solute concentration (g solute/g solvent mixture) at stage i. According to eqs 2, 5, and 6, eq 9 can be simplified as18 dci ci − 1 c = − i − 3kvρc Giμ2(i) αi − 1τi − 1 τi dt
(10)
where kv is the volume shape factor (π/6), ρc the crystal density (1.4 g/cm3).6,18,22 The mass balance equation for the solvent ethanol at stage i is derived as Vi
dvi (i − 1) (i) = Fout vi − 1 − Fout vi dt
(11)
where vi is the volume fraction of solvent in the solvent mixture. Similarly, it could be simplified as18 dvi vi − 1 v = − i αi − 1τi − 1 τi dt
(2)
(12)
In this study, the following empirical equations are used to express the growth rate and the nucleation rate:6,18,22
Apply j ≥ 1 on both sides of eq 1. To make the equations mathematically simple, assume the growth rate, Gi, is size-independent and the solution volume, Vi, is constant. The constant volume assumption can be easily implemented in practice by using a level controller in the vessel. Then, eq 3 can be obtained:18 Vi
Vi + 1 Vi
3
j ∫∞ 0 L (eq.1) dL,
dμj(i)
(5)
Similarly, the zeroth moment can be written as
(1)
j = 0, 1, ..., ∞
(i) Fout
18
where ni(L,t) is the crystal number density (no./cm ) at stage i, L is the crystal size (cm), Gi is the growth rate (cm min−1) at stage i, Vi is the volume of slurry (cm3) at stage i, F(i) out is the outlet flow rate (cm3 min−1) of stage i, and t is time (min). Define the jth moment at stage i by21
∫0
Vi
Then eq 4 can be written as the following ordinary differential equation (ODE):18
3
μj(i)(t ) =
j≥1 (4)
∂Vn ∂VG n (L , t ) i i(L , t ) (i − 1) (i) = Fout ni − 1(L , t ) − Fout ni(L , t ) − i i i , ∂t ∂L
∞
(i − 1) Fout Vi − 1 (i − 1) F (i) μj − out μj(i) + jGiμj(−i)1 , Vi − 1 Vi Vi
Define the residence time and dilution factor at stage i as τi and αi, respectively:18
2. MATHEMATICAL MODELING In this work, a continuous N-stage MSMPR cascade system is modeled via population balance model following the literature.18 All the vessels are assumed to be well-mixed. The compositions of the outlet flow are identical to the slurry inside the crystallizer. The solvent−antisolvent mixing and the crystallization process do not affect the solution volume. The crystal breakage and agglomeration are assumed to be negligible. Then, the PBE for the ith stage can be written as18−20
i = 1, 2, ..., N
=
j≥1
⎛ k ⎞ G = kG1 exp⎜ − G2 ⎟(c*(S − 1))kG3 ⎝ RT ⎠
(13)
⎛ k ⎞ ⎛ k ⎞ B = kB1 exp⎜ − B2 ⎟ exp⎜ − B2 3 ⎟ ⎝ RT ⎠ ⎝ ln S ⎠
(14)
S = c / c*
(15)
where c is the solute concentration (g solute/g solvent mixture), c* is the solubility (g solute/g solvent mixture) at a certain temperature and solvent fraction, S is the super-
(3)
Rearrange eq 3: 5674
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
Figure 1. Schematic of performing CCAC in MSMPR cascade crystallizers.
19 fix the feed solution temperature and the final temperature at stage N, while constraint 20 fixes the final solvent composition at stage N. The last equality constraint 21 fixes the total residence time as 30 min. The objective function, J, is some desired crystal properties at stage N. In this study the number-average mean size Ln is considered. Width of CSD can be represented by number-based coefficient of variation (CVn). However, because there is often a trade-off between crystal mean size and CV, only a single objective (crystal mean size) is included here in order to simplify the study. The optimization problem is solved using a sequential quadratic programming (SQP) algorithm. Matlab function f mincon and default specifications are used. The above solution approach is used throughout the entire study. Because it is not practical to set up more than 3 stages in the cascade in the pharmaceutical industry, the scenarios when N = 1, 2, 3 are computed. In addition, the batch system that has the same initial and final conditions as the feed solution and the Nth stage of the cascade system is compared as well. The batch time is 30 min, which is exactly the same as the total residence time of the cascade system. The initial solutions at each stage are set to have the same ethanol volume fraction (79.2%) as the feed solution. The solute concentration of the initial solutions is at saturation of 25 °C, which equals 0.222 g/g solvent mixture. The optimal operating profiles are presented in Table 1. In a certain MSMPR system, it can be observed that the optimal residence time τi would increase as the stage number i increases. This phenomenon does make sense, because in a typical crystallization operation, the time for growth should be much longer than that for nucleation in order to obtain large crystals. In a MSMPR cascade system, the first stage generally acts as a nucleation stage, whereas the rest of the stages are designed for growth. Also one can note that when more stages are used in the cascade system, the obtained crystals will be larger with smaller size variations. The batch system can achieve better products as compared to the cascade system because it corresponds to the limiting optimal operating trajectory that could be achieved using an infinite number of MSMPR crystallizers in the cascade. The batch crystallization is operated in equilibrium state whereas the continuous crystallization is operated in steady-state.23 The optimal trajectories of the four scenarios in Table 1 are compared in the three-dimensional (3D) phase diagram, as shown in Figure 2a. It can be noticed that the optimal cascade trajectory becomes more and more similar to the optimal batch trajectory when N increases, especially in the x and y directions, as plotted in Figure 2b. This phenomena demonstrates that the
saturation ratio. kG1, kG2, kG3, kB1, kB2, kB3 and solubility data used in this work are from literature.6,18,22
3. RESULTS AND DISCUSSION 3.1. Steady-State Optimization. Figure 1 presents the schematic of using both cooling and antisolvent addition to generate supersaturation in multistage MSMPR cascade crystallizers. The feed solution is free of crystals so that the first stage is a nucleation-dominated stage while the remaining stages are growth-dominated stages. Because there are two means to manipulate the supersaturation driving force, one could consider using temperature only, antisolvent only, or both to control the supersaturation level for a certain stage. The continuous MSMPR cascade system can finally reach steady-state and produce consistent products under certain fixed operating conditions. In this study, the optimal steadystate operating profiles of MSMPR cascade crystallizers of N stages are obtained, using the maximum mean size as the objective. The temperature and the ethanol volume fraction of the feed solution are fixed as (T0 = 40 °C, v0 = 79.2%) (79.2% in volume = 75% in mass), respectively. The temperature and the ethanol volume fraction of the final stage, or Nth stage, are fixed as (TN = 25 °C, vN = 45.8%), (45.8% in volume = 40% in mass), respectively. The volumetric flow rate of the feed solution is a constant of F0 = 5 mL/min. The optimization problem can be written as max{Jn . m . s . = Ln(N ) = μ1(N ) /μ0(N )},
Ti , Fi , τi
1≤i≤N
(16)
subject to model eqs 7, 8, 10, and 12−15 and T0 ≥ T1 ≥ T2 ≥ ···≥TN
(17)
T0 = 40 °C
(18)
TN = 25 °C
(19)
N
∑ Fi = 3.65 mL/min i=1
(20)
N
∑ τi = 30 min i=1
(21)
where Ti is the temperature (°C) at stage i and Fi is the antisolvent addition rate (mL min−1) at stage i. They are the optimization variables. TN is the temperature (°C) of the final stage or Nth stage. The inequality constraint 17 limits the cooling rate to be nonnegative. The equality constraints 18 and 5675
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
crystallization), as well as for the 3D crystallization (i.e., CCAC). In other words, the knowledge of the optimal operating procedure for a CCAC process in a batch system can be transferred to a nearly equivalent CCAC operation in a MSMPR cascade system. For instance, if it is already known that cooling first is desired in CCAC in a batch system to produce large crystals, then it is very likely that cooling first should also be preferred in CCAC in a MSMPR cascade system. This would be implemented by operating the first stage as a cooling crystallization and the second as antisolvent. If the batch operation requires simultaneous cooling and antisolvent addition, then most likely both stages of the MSMPR cascade crystallizers would be operated using cooling as well as antisolvent addition. 3.2. Startup Optimization. The startup problem is very important in continuous manufacturing. The crystal properties are inconsistent before the continuous system can reach steady state. As a result, the time and the products are wasted during startup. Therefore, the development of feasible methods to shorten the startup duration time and reduce the waste in the continuous crystallization systems can have high impact in continuous crystallization development and operation. 3.2.1. Startup Optimization via Using Different Initial Solutions. Typically, it is necessary to prepare some initial solutions for each stage before one turns on a continuous MSMPR cascade system. Theoretically, the compositions of the initial solution at each stage should have some impact on the startup dynamic behavior and perhaps even the final steadystate product quality if multiplicity exists in the model. Hence, the multiplicity of the MSMPR cascade crystallizers when N equals 1, 2, and 3 is studied first. To simplify the analysis, it is assumed that the initial solutions for each stage in a cascade system are the same when N equals 2 or 3. Then the process is simulated via eqs 7, 8, 10, and 12−15, with different initial solvent volume fractions and solute concentrations. The results are plotted in Figure 3, which indicates that the cascade system will always reach the same steady-state independently of the composition of the initial solutions in the vessels. Hence, there is no multiplicity in the MSMPR cascade system in our case. In other words, the cascade system is robust to the variation of initial conditions. It also indicates that the results of the steadystate optimization are independent of how the initial solution compositions are chosen. However, the initial conditions can provide a way to shorten the startup duration time and reduce the waste by optimizing the compositions of the initial solutions. Then, a 2-stage MSMPR cascade system is modeled here because N = 2 is the most often used configuration in the pharmaceutical industry considering both performance and the operation complexity. The solvent volume fraction and solute concentration of the initial solution at each stage are the optimization variables. The objective is the minimum t99 or t95, which can be defined as the time required for the system to reach a controlled state of operation (CSO) that has less than 1% or 5% variation, respectively. The optimization problem can be written as
Table 1. Optimal Operating Temperatures and Antisolvent Addition Rates for Cascade of N = 1, 2, 3 and for the Batch System 1 stage
2 stages
3 stages
batch
T1 (°C) T2 (°C) T3 (°C)
25
25 25
26.09 25 25
− − −
F1 (mL min−1) F2 (mL min−1) F3 (mL min−1)
3.65
1.05 2.60
0.00 1.60 2.05
− − −
τ1 (min) τ2 (min) τ3 (min)
30
10.74 19.26
6.61 11.37 12.03
− − −
V1 (mL) V2 (mL) V3 (mL)
259.5
64.98 166.60
33.05 75.04 104.06
− − −
Ln (μm) CVn (−)
257.39 1.000
369.55 0.790
449.90 0.683
704.69 0.250
Figure 2. Optimal trajectories of 1-stage (dotted line), 2-stage (dashdotted line), 3-stage (dashed line) cascade system and batch (solid line) system, shown in (a) 3D phase diagram and (b) its 2D projection.
min t 99 , ci , vi
1≤i≤N
(22)
subject to model eqs 7, 8, 10, and 12−15 and
batch trajectory can be regarded as a trajectory of a MSMPR cascade system of an infinite number of stages for both the simple 2D crystallizations (i.e., cooling or antisolvent only 5676
0 ≤ v1 , v2 ≤ 1
(23)
0 ≤ c1 ≤ solubility(T1, v1)
(24) DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
compositions discussed in section 3.1 is also used here as base case for comparison (v = 79.2%, c = 0.222 g/g solvent mixture or 25 °C saturated solution). This solution is used as the initial solution for both stages. The optimization results are shown in Table 2 and Figure 4. The waste material produced in each scenario (amount of Table 2. Comparison between the Compositions, t99 and Waste(t99) of the Base Case and the Optimal Case, Which Uses min{t99} as the Objective
Figure 3. Final crystal mean size under various solvent volume fractions and solute concentrations of the initial solution for cascades that have different number of stages: (a) N = 1, (b) N = 2, and (c) N = 3.
0 ≤ c 2 ≤ solubility(T2 , v2)
base case (t99)
optimal case (t99)
c1 (g/g solvent mixture) v1 (−) c2 (g/g solvent mixture) v2 (−)
0.222 79.2% 0.222 79.2%
0.201 74.7% 0.026 81.8%
t99 (min) Waste(t99) (g)
143.1 156.6
106.7 94.0
Figure 4. Mean size Ln(t) and Waste(t99) (a) before and (b) after optimizing the compositions of the initial solutions when using min{t99} as the objective.
(25)
where v1, v2, c1, and c2 are the ethanol volume fractions and solute concentrations of stage 1 and 2, respectively. The inequality constraints 23−25 can ensure that the optimization variables are in a proper range and that the solutions can be easily prepared in practice (no supersaturated initial solution is allowed). A constant operating profile obtained from Table 1, when N = 2, is applied in this section. A scenario with the initial solution
solute crystallized until the controlled state of operation is achieved) is also computed using eq 26: Waste(t 99) =
∫0
t 99
(N ) αvρc Fout (t )μ3(N )(t ) dt
(26)
The optimal case in Table 2 and Figure 4 indicates an improvement on the startup duration time t99 as well as 5677
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research Waste(t99). However, the dynamic behavior of Ln(t) did not change too much, in terms of the overshoots and the oscillations. Additionally, it can be noticed that in the optimal case, the initial concentration of the second stage is much smaller than that of the first stage. This is because the second stage, or the growth stage, needs relatively smaller supersaturation as compared with the first stage or the nucleation stage. In practice, instead of 1% deviation, sometimes 5% deviation is allowable. Hence, min{t95} is also considered here as the optimization objective, which can be written as min t 95 , ci , vi
1≤i≤N
(27)
where the constraints remain the same as in eqs 23-25. Then the waste that uses t95 as the criterion should be defined as Waste(t 95) =
∫0
t 95
(N ) αvρc Fout (t )μ3(N )(t ) dt
(28)
The optimization results using min{t95} as the objective are listed in Table 3 and Figure 5. A similar improvement could be Table 3. Comparison between the Compositions, t95 and Waste Amount of the Base Case and the Optimal Case, Which Uses min{t95} as the Objective base case (t95)
optimal case (t95)
0.222 79.2% 0.222 79.2% 100.1 97.2
0.207 75.0% 0 73.4% 59.9 32.6
c1 (g/g solvent mixture) v1 (−) c2 (g/g solvent mixture) v2 (−) t95 (min) Waste(t95) (g)
Figure 5. Mean size Ln(t) and Waste(t95) (a) before and (b) after optimizing the compositions of the initial solutions when using min{t95} as the objective.
observed as using min{t99} as the objective. The difference is that when using t95, the improvement is more significant. This can be explained by the larger flexibility for optimization offered by the looser limits that define the controlled state of operation with 5% deviation as compared to the case with 1% deviation. 3.2.2. Startup Optimization via Applying Dynamic Operating Profiles. Typically, the continuous MSMPR cascade system should be operated under some constant operating conditions to maintain steady-state and consistent product properties. However, during the startup it might be helpful to apply dynamic operating profiles to shorten the startup duration and reduce the waste. In the dynamic startup study, a 2-stage MSMPR cascade system is investigated. The computation of the optimal dynamic temperature or/and antisolvent addition profiles are formulated as a nonlinear optimization problem solved using a SQP approach. Because the startup duration time is typically around 100 min, as shown before, a time period of 300 min is prepared for observation. The 300 min time period is then divided into 30 intervals, and the temperatures and the antisolvent addition rates at every discrete time are considered as the optimization variables. 3.2.2.1. Scenario 1: Dynamic Antisolvent Addition Profiles with min{t99} as the Objective. If using t99 as the criterion and dynamic antisolvent flow rates with constant temperatures as operating profiles, the optimization problem can be written as min
F1(k), F2(k)
t 99 ,
0 ≤ k ≤ 30
0 ≤ F1(k) ≤ 10,
k = 0, 1, ..., 15
(30)
0 ≤ F2(k) ≤ 10,
k = 0, 1, ..., 15
(31)
F1(k) = F1,ss ,
k = 16, 17, ..., 30
(32)
F2(k) = F2,ss ,
k = 16, 17, ..., 30
(33)
where F1(k) and F2(k) are the antisolvent addition rates at the kth discrete time for stage 1 and 2, respectively. F1,ss and F2,ss are the steady-state operating values obtained from Table 1. Inequality constraints 30 and 31 limit the antisolvent addition rate in the proper operating range. Equality constraints 32 and 33 force the antisolvent addition rates back to the desired steady-state values to ensure that the dynamic antisolvent profiles will not affect the product properties at the final steadystate. The initial solution compositions are the same as those in section 3.1, or the base case in section 3.2.1. All the other variables, such as the temperature, volume, and feed solution flow rate and composition remain the same as the steady-state values listed in Table 1. Thus, the base case in Table 2 and Figure 4a can also be used as the base case here for comparison. The optimal dynamic antisolvent profiles when using t99 as the criterion are plotted in Figure 6. The antisolvent flow rates during startup are larger than the steady-state values, which makes sense because high supersaturation is necessary to make the system crystallize more quickly and strongly during startup than in the steady-state. The dynamic behavior of the mean
(29)
subject to model eqs 7, 8, 10, and 12−15 and 5678
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
where T1(k) and T2(k) are the temperatures at the kth discrete time for stage 1 and 2, respectively. Equations 37 and 38 contain more equality constraints than 32 and 33. In 37 and 38, the equality constraints include T1(0) = T1,ss and T2(0) = T2,ss, whereas in 32 and 33 of scenario 1, F1(0) = F1,ss and F1(0) = F1,ss are not listed as constraints. This is because in practice the antisolvent flow rates can achieve a large step change within a negligible time period, whereas the temperatures cannot. The optimization results of scenario 2 are plotted in Figure 7, in which some phenomena similar to those in scenario 1 can be
Figure 6. (a) Optimal dynamic antisolvent addition profiles and (b) corresponding mean size Ln(t) and Waste(t99) when using min{t99} as the objective.
size, Ln(t), shows excellent trajectory, with almost no overshoots or oscillations anymore. The t99 and Waste(t99) in Figure 6 are also reduced significantly, as compared to those of the base case (Figure 4a). Thus, unlike changing the initial solution composition, the dynamic profiles can have much stronger influence on the dynamic behavior during startup. In addition, another scenario that uses minWaste(t99) instead of min{t99} as the objective is also computed. Finally, an optimal result similar to that in scenario 1 is obtained. This is because overshoot and oscillations are eliminated; hence, the lower t99 is, the lower Waste(t99) would be. Thus, it is reasonable to consider min{t99} as the only objective. 3.2.2.2. Scenario 2: Dynamic Temperature Profiles with min{t99} as the Objective. Alternatively, dynamic temperature profiles could be applied instead of dynamic antisolvent addition profiles. The optimization problem could be formulated in a similar way: min t 99 ,
T1(k), T2(k)
0 ≤ k ≤ 30
Figure 7. (a) Optimal dynamic temperature profiles and (b) corresponding mean size Ln(t) and Waste(t99) when using min{t99} as the objective.
observed. The temperatures during startup are lower than the steady-state values, resulting in relatively higher supersaturation and quicker crystallization during startup. The optimal performance in scenario 2 is also very close to the performance in scenario 1, and both of these are much better than the base case (Figure 4a). 3.2.2.3. Scenario 3: Dynamic Combined Profiles with min{t99} as the Objective. In this case, startup optimization is considered using dynamic antisolvent profiles and dynamic temperature profiles simultaneously. The optimization problem using the dynamic combined profiles can be written as
(34)
subject to model eqs 7, 8, 10, and 12−15 and 15 ≤ T1(k) ≤ 40,
k = 0, 1, ..., 15
(35)
15 ≤ T2(k) ≤ 40,
k = 0, 1, ..., 15
(36)
T1(k) = T1,ss ,
k = 0, 16, 17, ..., 30
(37)
T2(k) = T2,ss ,
k = 0, 16, 17, ..., 30
(38)
min
F1(k), F2(k), T1(k), T2(k)
t 99 ,
0 ≤ k ≤ 30
(39)
subject to model eqs 7, 8, 10, and 12−15 and constraints 30−33 and 35−38. The optimization results of scenario 3 are shown in Figure 8. The performance and conclusions are generally similar to those of the previous two scenarios. The dynamic antisolvent flow 5679
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
dynamic antisolvent profiles only and min{t95} as the objective is considered. The only difference between scenario 1 and scenario 4 is the objective function. All the other optimization formulations are the same as in scenario 1. The optimal results are plotted in Figure 9, which indicates a great improvement on
Figure 9. (a) Optimal dynamic antisolvent addition profiles and (b) corresponding mean size Ln(t) and Waste(t95) when using min{t95} as the objective.
startup duration time and waste when compared with the base case (Figure 5a). However, because 5% deviation still allows room for overshoot and oscillation, they can hardly be eliminated when setting min{t95} as the objective. Finally, the performance of all scenarios and base cases are summarized in Table 4 for comparison, from which at least about 50% improvement can be observed when using dynamic operating profiles during startup.
Figure 8. (a) Optimal dynamic antisolvent addition profiles, (b) optimal dynamic temperature profiles, and (c) corresponding mean size Ln(t) and Waste(t99) when using min{t99} as the objective.
Table 4. Comparison among the Performances of All Scenarios and Base Cases
rates here have much stronger effect on decreasing startup duration than the dynamic temperatures. This is because the antisolvent flow rates can also affect the residence time whereas the temperatures cannot. Considering the dynamic operation complexity and cost, a practical way to implement this design is to apply dynamic antisolvent flow rates only while keeping the temperature profiles constant, as in scenario 1. 3.2.2.4. Scenario 4: Dynamic Antisolvent Profiles with min{t95} as the Objective. As shown in the previous section, the impact on the system of antisolvent flow rates is stronger than that of the operating temperatures. Thus, a scenario using
without dynamic profiles (base case) dynamic antisolvent profiles (scenarios 1 and 4) dynamic temperature profiles (scenario 2) dynamic combined profiles (scenario 3) 5680
t99 (min)
Waste(t99) (g)
t95 (min)
Waste(t95) (g)
143.1
156.6
100.1
97.2
83.0
80.4
27.5
6.6
83.1
81.0
−
−
80.1
77.4
−
−
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research Notes
3.2.3. Startup Optimization via a Model-Free Approach. Typically, the methods to control crystallization can be distinguished as model-based or model-free (direct design) approaches.24,25 In the previous discussion, a model-based approach is used, which is an offline control method, and it may not be robust to the outside disturbances as well as the model uncertainties. Thus, a model-free, or direct design approach can also be considered. In the model-free approach of the batch system, typically two methods are used via closed-loop control. One is the direct nucleation control (DNC), which uses focused beam reflectance measurement (FBRM) as the sensor to control the crystal counts in the slurry in a certain range.26,27 The other method is the supersaturation control (SSC), which uses attenuated total reflectance (ATR) spectroscopy to measure the solute concentration and controls the supersaturation level to follow a certain trajectory.28,29 These two methods have been demonstrated to be very efficient in the control of batch operation. Two similar methodologies can be proposed in the MSMPR cascade system to improve the startup behavior. By using FBRM as the sensor in the closed-loop, the crystal counts are measured in real time. When the counts are lower than the desired steady-state value, the control loop can result in a higher antisolvent addition rate or lower temperature to increase the supersaturation level, and vice versa. By using ATR spectroscopy, a higher antisolvent addition rate or lower temperature should be actuated when the supersaturation is measured lower than the desired steady-state value, and vice versa. These two model-free approaches are promising in the automated control of the startup behavior of continuous MSMPR cascade crystallizers and are currently being evaluated in practice in our group.
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS The authors thank Eli Lilly, IN, United States for partial financial support. The authors are grateful to Daniel Jarmer, Christopher Burcham, and Christopher Polster for the industrial input through discussions and comments provided.
■
4. CONCLUSIONS The CCAC process in continuous MSMPR cascade crystallizers is modeled using PBE in this work. Both the steady-state and the startup behaviors are investigated and optimized. Aspirin− ethanol−water is used as the model system in the study. In the steady-state optimization, the optimal steady-state operating profiles have been obtained. In addition, the optimal trajectories for both MSMPR cascade and batch systems have been plotted together in a 3D phase diagram from which the fact that the batch trajectory can be regarded as a trajectory of MSMPR cascade system of an infinite number of stages has been observed in 3D crystallization (i.e., CCAC) for the first time. In addition, different methods to shorten the startup duration time and reduce the waste have been discussed and proven to be feasible. These methods include changing the initial solution compositions, applying dynamic antisolvent profiles, applying dynamic temperature profiles, and applying dynamic antisolvent and temperature combined profiles. The most efficient and practical way to minimize startup duration and waste is to use dynamic antisolvent profiles only during startup. It is shown that this approach can reduce the startup duration and the amount of waste by approximately 50% while having no impact on the final steady-state product properties. Finally, a model-free approach to optimize the startup behavior is proposed using feedback control through FBRM or ATR spectroscopy techniques.
■
■
NOTATION αi = dilution factor at stage i (−) Bi = nucleation rate at stage i (no./cm3 min−1) c* = solubility (g solute/g solvent mixture) c0 = solute concentration of feed solution (g solute/g solvent mixture) ci = solute concentration at stage i (g solute/g solvent mixture) F0 = volumetric flow rate of feed solution (mL min−1) Fi = volumetric antisolvent addition rate at stage i (mL min−1) −1 F(i) out = volumetric outlet flow rate at stage i (mL min ) −1 Gi = growth rate at stage i (cm min ) kG1 = growth rate parameter, (3.21 ± 0.18) × 10−4 m s−1 kG2 = growth rate parameter, (2.58 ± 0.14) × 104 J mol−1 kG3 = growth rate parameter, (1.00 ± 0.01) kB1 = nucleation rate parameter, (1.15 ± 0.51) × 1021 m−3 s−1 kB2 = nucleation rate parameter, (7.67 ± 0.11) × 104 J mol−1 kB3 = nucleation rate parameter, (0.16 ± 0.01) kv = volume shape factor (−) L = particle characteristic size (cm) Ln = number-average mean size (μm) v0 = ethanol volume fraction of feed solution (−) vi = ethanol volume fraction at stage i (−) Vi = volume of mother liquor/solvent mixture at stage i (cm3) ni = number density at stage i (no./cm3) R = ideal gas constant (J mol−1 K−1) S = supersaturation ratio (−) t = time (min) Ti = temperature at stage i (°C) τi = residence time at stage i (min) ρc = crystal density (g cm−3) μ(i) j = jth moment at stage i REFERENCES
(1) Myerson, A. S. Handbook of Industrial Crystallization, 2nd ed.; Butterworth-Heinemann: Boston, 2002. (2) Nagy, Z. K.; Aamir, E. Systematic design of supersaturation controlled crystallization processes for shaping the crystal size distribution using an analytical estimator. Chem. Eng. Sci. 2012, 84, 656−670. (3) Quon, J. L.; Zhang, H.; Alvarez, A.; Evans, J.; Myerson, A. S.; Trout, B. L. Continuous crystallization of aliskiren hemifumarate. Cryst. Growth Des. 2012, 12, 3036−3044. (4) Alvarez, A. J.; Singh, A.; Myerson, A. S. Crystallization of cyclosporine in a multistage continuous MSMPR crystallizer. Cryst. Growth Des. 2011, 11, 4392−4400. (5) Wong, S. Y.; Tatusko, A. P.; Trout, B. L.; Myerson, A. S. Development of continuous crystallization processes using a singlestage mixed-suspension, mixed-product removal crystallizer with recycle. Cryst. Growth Des. 2012, 12, 5701−5707.
AUTHOR INFORMATION
Corresponding Author
*Tel.: +1-765-494-0734. Fax: +1-765-494-0805. E-mail:
[email protected]. 5681
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
Article
Industrial & Engineering Chemistry Research
(28) Woo, X. Y.; Nagy, Z. K.; Tan, R. B. H.; Braatz, R. D. Adaptive concentration control of cooling and antisolvent crystallization with laser backscattering measurement. Cryst. Growth Des. 2009, 9 (1), 182−191. (29) Saleemi, A. N.; Rielly, C. D.; Nagy, Z. K. Comparative investigation of supersaturation and automated direct nucleation control of crystal size distributions using ATR-UV/vis spectroscopy and FBRM. Cryst. Growth Des. 2012, 12, 1792−1807.
(6) Yang, Y.; Nagy, Z. K. Model-based systematic design and analysis approach for unseeded combined cooling and antisolvent crystallization (CCAC) systems. Cryst. Growth Des. 2014, 14, 687−698. (7) Aamir, E.; Nagy, Z. K.; Rielly, C. D. Optimal seed recipe design for crystal size distribution control for batch cooling crystallisation processes. Chem. Eng. Sci. 2010, 65, 3602−3614. (8) Nowee, S. M.; Abbas, A.; Romagnoli, J. A. Antisolvent crystallization: Model identification, experimental validation and dynamic simulation. Chem. Eng. Sci. 2008, 63, 5457−5467. (9) Acevedo, D.; Nagy, Z. K. Systematic classification of unseeded batch crystallization systems for achievable shape and size analysis. J. Cryst. Growth 2014, 394, 97−105. (10) Simone, E.; Saleemi, A. N.; Nagy, Z. K. Raman, UV, NIR, and mid-IR spectroscopy with focused beam reflectance measurement in monitoring polymorphic transformations. Chem. Eng. Technol. 2014, 37, 1305−1313. (11) Simone, E.; Saleemi, A. N.; Tonnon, N.; Nagy, Z. K. Active polymorphic feedback control of crystallization processes using a combined Raman and ATR-UV/Vis spectroscopy approach. Cryst. Growth Des. 2014, 14, 1839−1850. (12) Nagy, Z. K.; Fujiwara, M.; Braatz, R. D. Modelling and control of combined cooling and antisolvent crystallization processes. J. Process Control 2008, 18, 856−864. (13) Hermanto, M. W.; Chow, P. S.; Tan, R. B. H. Operating strategy to produce consistent CSD in combined antisolvent-cooling crystallization using FBRM. Ind. Eng. Chem. Res. 2012, 51, 13773− 13783. (14) Sheikhzadeh, M.; Trifkovic, M.; Rohani, S. Adaptive MIMO neuro-fuzzy logic control of a seeded and an unseeded anti-solvent semi-batch crystallizer. Chem. Eng. Sci. 2008, 63, 1261−1272. (15) Zhang, H.; Quon, J.; Alvarez, A. J.; Evans, J.; Myerson, A. S.; Trout, B. Development of continuous anti-solvent/cooling crystallization process using cascaded mixed suspension, mixed product removal crystallizers. Org. Process Res. Dev. 2012, 16, 915−924. (16) Vetter, T.; Burcham, C. L.; Doherty, M. F. Regions of attainable particle size in continuous and batch crystallization processes. Chem. Eng. Sci. 2014, 106, 167−180. (17) Su, Q.; Nagy, Z. K.; Rielly, C. D. Pharmaceutical crystallization processes from batch to continuous operation using MSMPR stages: Modelling, design, and control. Chem. Eng. Process. 2015, 89, 41−53. (18) Yang, Y.; Nagy, Z. K. Advanced control approaches for combined cooling/antisolvent crystallization in continuous mixed suspension mixed product removal cascade crystallizers. Chem. Eng. Sci. 2015, 127, 362−373. (19) Ramkrishna, D. Population balances: theory and applications to particulate systems in engineering. Academic Press: San Diego, 2000. (20) Rawlings, J. B.; Miller, S. M.; Witkowski, W. R. Model identification and control of solution crystallization processes: A review. Ind. Eng. Chem. Res. 1993, 32, 1275−1296. (21) Randolph, A.; Larson, M. A. Theory of Particulate Processes, 2nd ed.; Academic Press: San Diego, 1988. (22) Lindenberg, C.; Krattli, M.; Cornel, J.; Mazzotti, M. Design and optimization of a combined cooling/antisolvent crystallization process. Cryst. Growth Des. 2009, 9 (2), 1124−1136. (23) Chen, J.; Sarma, B.; Evans, J. M. B.; Myerson, A. S. Pharmaceutical crystallization. Cryst. Growth Des. 2011, 11, 887−895. (24) Nagy, Z. K.; Fevotte, G.; Kramer, H.; Simon, L. L. Recent advances in the monitoring, modelling and control of crystallization systems. Chem. Eng. Res. Des. 2013, 91 (10), 1903−1922. (25) Nagy, Z. K.; Braatz, R. D. Advances and new directions in crystallization control. Annu. Rev. Chem. Biomol. Eng. 2012, 3, 55−75. (26) Saleemi, A. N.; Rielly, C. D.; Nagy, Z. K. Automated direct nucleation control for in situ dynamic fines removal in batch cooling crystallization. CrystEngComm 2012, 14, 2196−2203. (27) Bakar, M. R. A.; Nagy, Z. K.; Saleemi, A. N.; Rielly, C. D. The impact of direct nucleation control on crystal size distribution in pharmaceutical crystallization processes. Cryst. Growth Des. 2009, 9 (3), 1378−1384. 5682
DOI: 10.1021/ie5034254 Ind. Eng. Chem. Res. 2015, 54, 5673−5682
