Organic Process Research & Development

Article

pubs.acs.org/OPRD

Mathematical Modeling, Design, and Optimization of a Multisegment Multiaddition Plug-Flow Crystallizer for Antisolvent Crystallizations

Qinglin Su,*^ Brahim Benyahia^ Zoltan K. Nagy,^* and Chris D. Rielly*^

^Department of Chemical Engineering, Loughborough University, Loughborough, LE11 3TU, U.K. ^School of Chemical Engineering, Purdue University, West Lafayette, Indiana 47907-2100, United States

ABSTRACT: In the pharmaceutical industries, the requirements of rapid process development and scalable design have made the tubular crystallizer a promising platform for continuous manufacturing and crystallization processes, capable of replacing conventional capital- and labor-intensive batch operations. This paper takes a process systems engineering (PSE) approach to the optimal design of a continuous antisolvent addition crystallizer to deliver the most promising product qualities, such as the crystal size distribution. A multisegment multiaddition plug-flow crystallizer (MSMA-PFC) is considered as an example of a continuous antisolvent crystallization process, in which the total number, location, and distribution of antisolvent additions are to be optimized. First-principles dynamic and steady-state mathematical models for the MSMA-PFC are presented, based on example kinetic models for nucleation and growth of paracetamol crystallizing in acetone, with water as the antisolvent. The performances of different crystallizer configurations operated under optimal design conditions are then compared. The configuration in which antisolvent could be added at a variety of different locations along the tube length and at optimal flow rates was able to outperform previous designs in the literature which considered equally spaced antisolvent additions. The use of dynamic models to detect problems during startup of an MSMA-PFC was also highlighted.

■ INTRODUCTION

Crystallization is a common unit operation for separation and purification and is used for approximately 90% of organic molecules in the pharmaceutical and fine chemical sectors.1-3 Traditionally, batch operations have been used for the crystallization stage and for downstream secondary manufacturing processes, as they provide a flexible way to meet the stringent regulations for product quality and the variable demands of the market. However, rising market competition and the need to reduce manufacturing costs are driving the future of pharmaceutical and fine chemical industries toward continuous processes, which have potential for improvements in product quality, through online process monitoring and control, and reduced equipment footprint, energy, and labor

costs.

Over the past decade, the development of continuous manufacturing and purification processes, particularly crystallization, has mainly focused on the modification of existing batch units,8-10 in addition to studying innovative equipment design. Batch crystallizers are often based around stirred-tank technologies, which can be converted to continuous mode as multistage mixed-suspension and mixed product removal (MSMPR) operations,11-16 but they suffer from broad residence time distributions, leading to broad crystal size distributions and problems in downstream processes, such as filtration, isolation, drying, and solids mixing. Consequently, interest in tubular designs of continuous crystallizer has increased in recent years, because of the potential benefits of a much narrower residence time distribution, enhanced control of the supersupersaturation, and the ease of scaling-up from pilot to full-scale operation.

Recent research efforts have been devoted to the experimental investigation of tubular crystallizers. For example, Lawton et al.17 reported the application of a continuous oscillatory baffled crystallizer (COBC) for a cooling crystallization which offered significant advantages in operating cost and processing time. Eder18,19 and his colleagues investigated the impacts of flow rate and seed loading in continuously seeded tubular cooling crystallizers for the production of active pharmaceutical ingredients (APIs). Ferguson et al.20 reported the employment of a variety of process analytical technology (PAT) tools for particle size and solute concentration measurement in a plug flow crystallizer. Such PAT tools provide time-varying measurements at a single spatial location, but do not reveal the evolution of the CSD or solute concentration profile with increasing residence time in the flow device. These experimental studies provide valuable insight into the performance of tubular crystallizers, but with limited information about spatial evolution; moreover they show that there is a complex interaction of effects involving the flow configuration, the use of seeding, the supersaturation profile, and the nucleation and growth kinetics on the final product qualities, such as the crystal size distribution. Without a process systems model, they could only provide qualitative indications of how to optimize the crystallizer design to produce crystals with a target size distribution, or of particular shape or polymorph type. Thus, the design of such continuous processes is complicated by the interaction of numerous operating

Special Issue: Polymorphism & Crystallisation 2015 Received: April 8, 2015

ACS Publications © XXXX American Chemical Society

variables, which affect the supersaturation profile and the product quality, in terms of its CSD. A purely experimental investigation of the crystallizer operating space would be time-consuming and difficult to interpret; hence the operation would be almost impossible to design and optimize without a model framework.

Other studies make use of first-principles models, alongside experimental studies, to deduce the optimal operating conditions to produce crystals of a specified size. For example, Alvarez and Myerson1 modeled and conducted experiments on a multisegment, multiaddition plug-flow crystallizer (MSMA-PFC) with antisolvent addition empirically distributed along the tube and showed how the device could be operated to produce crystals of a small mean size and narrow CSD. Zhou et al.21 studied experimentally a novel concentric annulus to achieve layer crystallization and formulated a model which predicted the profile of the temperature, concentration, and crystal thickness along the pipe length; their model was used to optimize the crystal yield and the layer thickness.

Process systems engineering (PSE) approaches using mathematical modeling, intensification, and optimization are now being applied to the design of tubular crystallizers, e.g., by Lakerveld et al. Majumder and Nagy optimized the temperature profile along a multisegment PFC, which included cooling and heating segments, with the objective of removing crystal fines by controlled dissolution. Vetter et al.23 recently investigated the attainable regions of particle sizes for a singlestage ideal plug-flow crystallizer, which assumed continuous addition of antisolvent along the tube; they compared the achievable product qualities to those manufactured in MSMPR and batch crystallizers. Ridder et al.24 proposed a simultaneous design and control (SDC) approach to optimize the number of segments and antisolvent distributions along an MSMA-PFC for an antisolvent crystallization process. More specifically, recent studies by Kwon et al.25 also focused on the multisegment plug-flow crystallizer to control the crystal shape and size distribution by optimizing the jacket temperature for each segment. These simulations are based on first-principles process systems models using spatially distributed mass, heat, and population balances and include the various complexities in crystallization, such as the effect of antisolvent addition, which may increase or decrease the supersaturation depending on the slope of the solubility curve. Moreover, these flow sheet simulations of continuous crystallizers are often based on bespoke code, since plug flow modules are not yet available in commercial software packages, such as the gCRYSTAL 4.0 software offered by Process System Enterprises Ltd.

The study reported here is part of a larger project to build dynamic simulations of continuous crystallization and downstream secondary manufacturing operations, to optimize performance and test out control strategies for an integrated pharmaceutical production plant flow sheet. This paper presents an example of the first stage, which considers the mathematical modeling, design, and optimization of tubular crystallizers for antisolvent crystallization, to provide a better understanding of flexible design, and to improve the process performance accordingly. Dynamic models are used to examine the nonlinear system response to external manipulations or disturbances for the purpose of process control and can be applied to design start-up and shut-down procedures, where changes of operating states are required with minimal production of waste or off-specification product. The optimal

design of such transient operations is critical for a cost-effective continuous process, particularly when the production capacity is scaled up and/or the manufacturing time is short.26 In contrast, a steady-state model is computationally more efficient when dealing with process design and optimization. However, dynamic models are still required to confirm that the optimal operating conditions from a steady-state model are obtained reliably and lead to stable operation.

The multisegment multiaddition design of a plug-flow crystallizer is a generalized type of antisolvent crystallization which aims to provide accurate control of supersaturation and antisolvent mass fraction profiles along the length of the tube, as shown in Figure 1. The modular design of the PFC

Tube length, Z

Figure 1. Schematic of a multisegment multiaddition plug-flow crystallizer (A-: antisolvent addition for the 1th segment; S,-: fresh solution addition; Z,-: tube length; C,-: main solute concentration; C*,-: main solute solubility in 1th segment; Ttemperature; N: total number of segments).

segments,1 where flexibility is a key aspect, is intended for campaign production of different active pharmaceutical ingredients using the same processing line. This is an important consideration for the design of a continuous manufacturing process.26 For example, varying lengths of tube and longer mean residence times can be obtained by assembling a greater number of PFC segments and the distribution of antisolvent additions along the segments can be optimized to create different supersaturation profiles along the tube,24 as indicated schematically in Figure 1.

This paper is organized as follows. Dynamic and steady-state mathematical models of a general plug-flow crystallizer are presented, followed by their extensions to a multisegment multiaddition design for antisolvent crystallization. Unlike the equally spaced antisolvent addition proposed by Ridder et al.,24 the configuration studied here allows for optimization of both the locations and distributions of antisolvent additions along the tubular crystallizer; a mixed integer nonlinear programming (MINLP) optimization is used for this purpose. The Results and Discussion section considers the application of the proposed design and optimization framework to an antisolvent crystallization of paracetamol in acetone with water as the antisolvent and compares its performance to other configurations reported in the literature. Finally, concluding remarks are made.

■ MATHEMATICAL MODELING OF PLUG-FLOW CRYSTALLIZER

In the pharmaceutical industries, a typical bulk API is usually produced at a rate of >1 tons per day and thus the diameter of a tubular crystallizer could be of the order of a few centimeters,

which is small when compared to the typical tube length (tens of meters). As a result, there is often near perfect mixing in the radial direction and limited dispersion in the axial direction, which means that the assumptions of a plug flow device are likely to apply. A relatively high slurry mean velocity along the tube is required to avoid the sedimentation of crystals (or, in the case of an oscillatory baffled crystallizer,17 high-frequency oscillations are introduced to suspend the crystals and produce a narrow residence time distribution, but not perfect plug flow, as is discussed at the end of this section). Often it is reasonable to assume that the tubular crystallizer is a plug-flow crystallizer and therefore requires only one spatial dimension, i.e., the tube axial length, in a mathematical model.

The governing equations for a general segment of a plug-flow crystallizer are based on a population balance equation for the crystal size distribution and a mass balance equation for solute concentrations, as shown below:

C(z, 0) = C0(z)

dn + d(vzn) + d(Gn) dt dz dL

dC d(vzC) typKv — + +-—

dt dz Mw

S(n, L, z, t)

GLndL = 0

In eqs 1 and 2, n is the number probability density function of crystals in the slurry, #/m3/m; t is the time, s; vz is the slurry superficial velocity along the tube axial, m/s; z is the tube axial coordinate, m; G is the crystal growth rate, m/s; L is the characteristic length of crystal, m; S represents the source and sink terms accounting for crystal breakage and agglomeration, #/m3/s; C is a vector representing molar concentrations of the various chemical species (including the solute, the antisolvent, impurities if any, etc.), kmol/m3; ф is the chemometrics vector for the chemical species mole fractions in the crystal phase; ps is the crystal density, kg/m3; Kv is the volumetric shape factor of crystal; Mw is the mean molecular weight of the crystal, kg/ kmol. In this study, a single crystal phase is considered, but extension of the balance equations to multiple crystal phases or to multiple polymorphs is straightforward.

The plug flow crystallizer is an unsteady and spatially distributed system, and hence the number density n and solute concentrations C are both functions of position along the tube z and time t; their boundary and initial conditions are listed below.

Boundary conditions for n(L, z, t) and C(z, t) at L = L0 or z = 0:

n(Lo, z, t)

n(L, 0, t) = nfeed(L, t)

C(0, t)

■feed

where B is the nucleation rate, #/m3/s; L0 is the nuclei size, m; nfeed is the number probability density function of the feed CSD at the crystallizer inlet, #/m3/m, which for example, is required for a seeded crystallization, or where a feed slurry is supplied from an upstream unit; and Cfeed represents the feed species concentrations at the crystallizer inlet, kmol/m3.

The initial conditions for n (L, z, t) and C (z, t) at t = 0 are as follows:

where n0 (#/m3/m) and C0 (kmol/m3) are the initial number probability density function for the CSD and species concentrations, respectively, inside the crystallizer.

A steady-state model of a plug-flow crystallizer can be obtained straightforwardly by removing the time derivative terms of the crystal number density n and solute concentrations C from the dynamic eqs 1 and 2 and retaining the boundary conditions of eqs 3 to 5 only. The slurry superficial velocity vz along the crystallizer is not necessarily constant; e.g., if the volume changes due to crystallization (formation of solids) along the crystallizer it cannot be neglected.16 Hence, vz will depend on the yield of crystals at each point along the crystallizer and should be rigorously calculated. For simplicity, we assume a constant slurry mean velocity vz along the tube1 and there is no breakage or agglomeration in the tube. Hence the steady-state model can be simplified as follows:

dn 1 d(Gn) — + —- = 0

dz vz dL

(2) dC

dC . 3фрKv , ^t2

GLndL = 0

The above partial differential equations (PDEs) 1 and 2 for the dynamic model can be solved in MATLAB with a highresolution finite-volume method (FVM)28 by discretizing both the tube axial length z and the crystal characteristic length L and then integrating the resulting set of ordinary differential equations (ODEs), together with eq 2, with respect to the system time t (i.e., using the method of lines). The FVM scheme is of second-order accuracy, combining a robust upwind discretization method with the novel k = 1/3 flux limiter, capturing the sharp front of the nucleation boundary in eq 3 without numerical oscillations, and providing a smooth solution. In terms of steady-state modeling of the population balance eq 8, only discretization in crystal length L and integration, together with the mass balance of eq 9 for solute concentrations, along the tube length z are necessary. All the integrations of ODEs were executed with the built-in "ode45" or "ode23" function for nonstiff problems in MATLAB.

Depending on the fineness of discretization in the FVM, the PDEs often result in tens or hundreds of ODEs which can be computationally expensive to solve. Alternatively, when the crystal growth rate G is size independent, the classical method of moments (MOM) can also be applied to integrate the PDE (8) along the crystal length L by conversion into a set of moment ODEs, which further reduces the computational burden for steady-state modeling. This is critically important when a large optimization problem is considered. In the MOM, the integration and moment transformation of the PDE into the first six moments are given below.

= f Lkn(L, z)dL, k = 0 ••• 5 J 0

n(L, z, 0) = n0(L, z)

d^0 dz

klf к- +1 vr

where m is the kth moment of the number CSD probability density function n, mk/m3. Although the full CSD development along the tube z is lost by the integration of eq 10, the physical meanings of the first few moments provide useful information: e.g., m0 is the total particle number per volume of the slurry; is related to the total volume of crystals per volume of the slurry. Therefore, the volume-based mean crystal size (L43) and the coefficient of variation (CV) of the full CSD still can be captured by the MOM, as shown in eqs 13 and 14. Furthermore, after the optimization problem is solved based on the MOM model, the full CSD can be recovered by solving the full model of eqs 8 and 9. In such a way, less-expensive function evaluations are possible for the optimization algorithms.24 Identical results were obtained with MOM and FVM when a suitably fine level of discretization was employed.29,30

CV ■

M5M3 2 M4

The mathematical model could easily be extended to a COBC, which is not a perfect plug flow device. The COBC could be considered as a cascade of MSMPR crystallizers in series or a plug-flow type crystallizer with dispersion terms Dndn2/dz2 and DcdnC2/dz2 added to eqs 1 and 2 for dynamic simulation, or to eqs 8 and 9 for steady-state modeling, respectively. Here Dn and Dc are dispersion coefficients for crystals and liquid solution, respectively. Both the number of stages in a cascaded MSMPR crystallizer and the dispersion coefficients could be estimated from experimental residence time distributions of a COBC, which will be a future development in our work.

■ OPTIMIZATION OF MSMA PLUG-FLOW CRYSTALLIZER

In a similar way to a cascaded multistage continuous MSMPR crystallizer,11'12'15 several segments of tubular crystallizer can be joined together with temperature control applied to each segment and/or additions of antisolvent or fresh solution injected between any two consecutive segments in order to regulate the main solute concentration and supersaturation along the tube, as shown in Figure 1. In the present study, each segment is regarded as an individual plug-flow crystallizer, which forms a part of the multisegment, multiaddition plug-flow crystallizer. The current study considers the addition of antisolvent to create the supersaturation profile.

In terms of the dynamic simulation of the MSMA-PFC, the outlet slurry of one PFC segment is assumed to mix instantly with the added fresh solution or antisolvent stream under an ideal mixing rule. Both solvent and antisolvent have a low viscosity and are fully miscible with each other; moreover, the system does not nucleate very rapidly and hence the assumption of instantaneous mixing may be justified. In practice, if mixing of the antisolvent stream with the slurry is limiting, then a high addition velocity or static mixers inside the tube can be applied to intensify the mixing conditions.1 Hence the crystal number probability density function n and the solute concentrations C at the outlet are diluted accordingly and thus need to be updated and set as inlet boundary conditions for the next PFC segment, viz., eqs 4 and 5. Similar calculations also

apply to the steady-state model. For more information on this, readers are referred to the work by Ridder et al.24

Previously, Alvarez and Myerson1 reported a combination of four equal-length tubular units (length: 4 m X 0.6 m; diameter: 1.27 cm), with empirically designed distribution of antisolvent among the four units, for example, 50%, 50%, 0%, 0%, and with a fixed amount of total antisolvent addition. Ridder et al.24 optimized the distribution of antisolvent for a number of equally spaced injection points, where up to 15 injection points were considered for a 50 m long tube, requiring a unit segment length of 50/15 m. However, the segmentation of their design method is only optimal for a specific crystallization system; it may not be optimal for other systems, and thus their multisegment design is not as versatile as generally demanded by the pharmaceutical industries. Furthermore, depending on the relative competitiveness of crystal growth rate and nucleation rate, varying ranges of supersaturation along the tubular crystallizer may be required to achieve the best product attributes. Hence, without resorting to optimize too many equally spaced injections, it is more desirable and useful to optimize both the locations and distribution of a limited number of antisolvent or fresh solution injections to achieve better control of supersaturation. This new method of optimization of the continuous MSMA-PFC provides the same level of operating flexibility as a batch crystallizer in which antisolvent can be added at any time and hence should result in better control or product quality attributes.

To sum up, a practical multisegment multiaddition tubular crystallizer should have the concept of a modular design, which comprises standardized unit segments, for example, a module of 0.6 m long, as studied by Alvarez and Myerson.1 Then for a tubular crystallizer made up of N segments (hence with N possible addition points), the locations and distribution of antisolvent/fresh solution additions for a total number of m (m < N) injections points could be optimized to obtain better product qualities. At the current stage, the optimization of only the antisolvent addition is considered, as shown in eqs 15 to 18 below; optimization of both antisolvent and fresh solution additions (and temperature) would be straightforward and can be considered in future work. The optimization problem can be stated as to maximize the product qualities (e.g., mean crystal size, product yield, etc.) P at the outlet of the final segment ZN

MAX(P)

subject to a number of linear and nonlinear constraints on the product qualities P, represented by eq 16

f (P)Zn ^ 0

The optimization represented by eq 15 involves

(i) changing the location of m antisolvent additions; in eq 17, U represents an integer list of index numbers of the addition points; there is always an antisolvent addition in segment 1 and the remaining m — 1 additions may occur at the start of segments between 2 and N;

(ii) changing the mass fraction x, of the ith antisolvent addition, where i = 1---m, as represented by the list in A. The total flow rate of antisolvent added in each case is treated as fixed.

U = [u1, u2, u3, ..., um] where u1 = 1

Table 1. Summary of the three case studies

case m U A l43 (^m) CV czn (kg/kg)

1 1 [1] [1.000] 278.5 0.182 0.0795

2 [161] [0.500, 0.500] 273.3 0.160 0.0798

3 [1 41 81] [0.333, 0.333, 0.333] 374.7 0.150 0.0804

4 [1 31 61 91] [0.250, 0.250, 0.250, 0.250] 369.7 0.301 0.0806

5 [1 25 49 73 97] [0.200, 0.200, 0.200, 0.200, 0.200] 327.1 0.196 0.0808

6 [1 21 41 61 81 101] [0.166, 0.166, 0.166, 0.166, 0.166, 0.166] 361.5 0.188 0.0812

2 1 [1] [1.000] 278.5 0.182 0.0795

2 [1 61] [0.245, 0.755] 494.9 0.251 0.0805

3 [1 41 81] [0.272, 0.376, 0.353] 439.5 0.226 0.0806

4 [1 31 61 91] [0.242, 0.025, 0.733, 0.001] 499.0 0.247 0.0805

5 [1 25 49 73 97] [0.256, 0.018, 0.369, 0.341, 0.016] 470.3 0.225 0.0806

6 [1 21 41 61 81 101] [0.233, 0.035, 0.079, 0.546, 0.106, 0.001] 483.5 0.229 0.0806

3 1 [1] [1.000] 278.5 0.182 0.0795

2 [1 62] [0.246, 0.754] 497.2 0.248 0.0806

3 [1 48 62] [0.238, 0.135, 0.627] 511.7 0.248 0.0806

4 [1 46 54 63] [0.241, 0.074, 0.135, 0.551] 519.6 0.223 0.0806

5 [1 53 61 62 64] [0.241, 0.135, 0.114, 0.156, 0.356] 520.0 0.230 0.0806

6 [1 19 50 56 63 67] [0.237, 0.003, 0.071, 0.111, 0.516, 0.063] 521.1 0.238 0.0806

A = [x1, X2, ..., xm] where ^ Xi

Under the proposed design and optimization framework, it is convenient to consider an MSMA-PFC with a desired total length assembled from a fixed number of unit segments1 and also with an optimized antisolvent addition. Thus, it is possible to make the MSMA-PFC flexible and efficient for a variety of crystallization systems. As expected, when a large number of injections are chosen (m « N), there would be only marginal benefits over the equally spaced injections; hence, the list U of injection positions can be fixed and only the distribution vector A is allowed to vary in the optimization problem.

The above optimization problem of eqs 15 to 18 is a mixed integer nonlinear programming problem (MINLP) and can be solved by the genetic algorithm with an integer constraint "intcon" in MATLAB. The genetic algorithm is an adaptive heuristic search method based on the evolutionary ideas of natural selection and genetics and is capable of searching a large or multimodal state-space, offering significant benefits over more deterministic optimization techniques.

Multiobjective optimizations have been widely applied to the crystallization process, e.g. to observe the Pareto front among several key performance indices, such as the volume-based mean crystal size, coefficient of variation (CV), yield, or the ratio of seeded crystals over nucleated crystals.24'31'32 In the current application, the residence time distribution is narrow in a plug-flow crystallizer and should not contribute to a broadening of the CSD through back-mixing effects. Instead an increase in the CV is more likely to result from multiple nucleation events along the tube, triggered by antisolvent additions and potentially leading to a multimodal CSD. Multiobjective optimization of the mean crystal size and the CV is possible; however, here the objective function was simply written in terms of maximizing the mean particle size, but adding the constraints that (i) a minimum yield had to be obtained and (ii) the CV had to be less than 0.30, ensuring a relatively narrow CSD. This approach is intended to minimize the effects of multiple nucleation events and to provide a narrow distribution of large crystals suitable for secondary manufacturing processes, e.g. filtration or continuous mixing.33

Extensions of the approach to other (multi-) objective functions would be straightforward to implement, but the formulation proposed here is sufficient to demonstrate the ability to design optimized processes. Therefore, in this study, only the mean crystal size L43 is maximized with constraints imposed on both the CV and product yield (as shown in Table 1).

■ RESULTS AND DISCUSSION

Tubular Crystallizer Design. To demonstrate the proposed design and optimization framework for the MSMA-PFC, an unseeded antisolvent crystallization system of paracetamol in an acetone (solvent) and water (antisolvent) mixture at a constant temperature of 16 °C is investigated using numerical simulation. The feed solution of paracetamol is first saturated at a water mass fraction of 60% (0.1917 g solute/g solvents) and then injected into a 72 m long tube, which consists of 120 segments. Each modular segment is 0.6 m in length and 1.27 cm in diameter.1 Addition of antisolvent is also assumed to be possible only in the inlet of each modular segment and as described above is assumed to mix instantly with the slurry inside the tube. A total flow rate of 50 mL/min of fresh saturated solution is fed to the first segment, and a total flow rate of 25 mL/min of water as antisolvent is injected along the MSMA-PFC. The mean residence time is assumed to be fixed at 120 min (the mean flow rate vz in each modular segment is different due to the antisolvent addition, which affects the mean residence time in each segment), similar to the typical semibatch crystallization process reported by Woo et al.,34 in which the solubility model and crystallization kinetics equations (minor typo errors therein are corrected) are presented as follows for easy reference.

C*(kg solute/kg solvents)

= 1.302 X 10-V - 1.882 X 10-V + 2.237 X 10-—

+ 5.746 X 10-1 AC(kg solute/kg solvents) = C - C*

B(#/m3/s) = kbACb

(20) (21)

kb = 4.338 X 1058 exp(-1.374w) b = 1.997 X 10-3w2 - 6.237 X 10-1w + 4.042 X 101

G(m/s) = kgACg

kg = -9.6300 X 10-11w3 + 3.3558 X 10-8w2 - 1.2606 X 10-6w + 3.6852 X 10-5

g = -1.108 X 10-4w2 + 1.024 X 10-2w + 1.427 (26)

where w is the mass percentage of antisolvent in the solvent mixture; AC is the absolute supersaturation. Both the crystal growth rate G and the nucleation rate B are dependent on supersaturation and the antisolvent mass fraction; as previously mentioned the antisolvent mass fraction has a complex effect on supersaturation, which makes the control and optimization of product qualities very challenging. For example, in a semibatch crystallizer system,34 the antisolvent flow rate was required to increase exponentially to maintain a constant tradeoff of crystal growth and nucleation (so-called C-control34 for a semibatch process).

The selected crystallization kinetics was taken from batch crystallization experiments by Woo et al.,34 obtained through parameter estimation. They are used for illustrative purposes only and are not intended to perfectly represent the crystallization behavior in a continuous plug-flow crystallizer. The fluid mechanics in these two devices are not comparable, and hence phenomena such as the secondary nucleation processes induced by a moving impeller in a batch crystallizer are unlikely to be the same as in plug-flow crystallizer, through particle—wall, or particle—particle collisions. Realistic nuclea-tion models should be included in design calculations for continuous crystallizers, particularly given the sensitivity of the optimization to the detail of the nucleation kinetics.24 For the design of real crystallization processes, parameter estimation to obtain the rate laws relevant to the fluid flow in the process equipment should be conducted as part of a preliminary experimental study of a PFC design, but is not considered as part of the current work.

Optimization of Antisolvent Addition. As discussed previously, a computationally efficient steady-state model of MSMA-PFC using MOM is employed here for the optimal design problem of the antisolvent additions. The optimal results are then resimulated by the FVM to capture the full CSD development in the tubular crystallizer. This avoids the difficulty of trying to reconstruct the full CSD from its moments and reduces the computational time for the optimization; it is important here, because there is a possibility to obtain multimodal size distributions, which may not be evident from the moments alone.

Practically it would be infeasible to use a large number of injection points. The performance of the MSMA-PFC was studied for up to 6 injection points; the first is always at the entry to the first segment, and the remaining points are selected from the following 119 segments. Three case studies are considered:

• Case 1 serves as an unoptimized bench mark and considers equally spaced injection points with equally distributed antisolvent addition, as studied by Alvarez and Myerson;1

• Case 2 considers also the equally spaced injection points, but with optimized distribution of antisolvent addition, as studied by Ridder et al.;24

• Case 3 optimizes both the location of the addition points and the distribution of antisolvent

For Cases 2 and 3, the objective was to maximize L43, while maintaining CV < 0.30 and C(ZN) < 0.806 kg/kg, i.e. without broadening the crystal size distribution too much and with an acceptable yield. For an unseeded antisolvent crystallization in this MSMA-PFC, the width of the CSD is mainly determined by the magnitudes of multiple nucleation events. Paracetamol in acetone is a relatively slow growing system and is dominated by nucleation; hence, maximizing the mean crystal size under a minimum yield constraint will automatically limit these multiple nucleation events.

The genetic algorithm for the MINLP problem in MATLAB 2013b was implemented for the optimization problems of Cases 2 and 3. To cope with the stochastic nature of the genetic algorithm, population sizes of 30, 60, 100, and 200 were implemented together with maximum generations of 150 for each optimization scenario of Cases 2 and 3; the best solution of each was then chosen as the final optimal result.

The final optimization results of the case studies are summarized in Table 1. For example, in Case 2 of the MSMA-PFC with totally four injection points (m = 4), the location index U = [1 31 61 91] means the third antisolvent addition is located at the 61st segment, or axial position of z = (61 — 1) X 0.60 m = 36.0 m. The corresponding antisolvent feed distribution is given by A = [0.242, 0.025, 0.733, 0.001] and shows, for example, that the flow rate of the third antisolvent addition was 0.733 X 25 mL/min = 18.3 mL/min. Comparisons of the three case studies for four injection points in the MSMA-PFC are shown in Figures 2 to 4 for paracetamol concentration, evolution of crystal number density, and the volume-based CSD profiles at specific axial positions, respectively.

Figure 2 shows that all three cases maintain a certain level of supersaturation, which reduces between each addition; the first antisolvent addition at the inlet of the crystallizer produces the first burst of nuclei, which passes through the following PFC segments for further growth and finally contributes to the generation of the majority of the large particles at the outlet, as is illustrated in Figures 3 and 4. After the second addition, the equally distributed addition in Case 1 creates a large supersaturation, which results in the second peak in the broadened CSD (see Figure 4: case 1 for z > 18 m). Case 2 reduces the amount of antisolvent addition at the same position compared to Case 1 at z = 18 m; this extends the supersaturation level achieved in the previous segment so that before the third addition (where a large amount of antisolvent is added) there are enough medium-size crystals for growth to compete with the nucleation effect and thereby consume the supersaturation. In such a way, there is a significant increase of the final mean crystal size from 370 j«m of Case 1 to 499 j«m of Case 2 as shown in Table 1; the CV was constrained to be <0.30. Interestingly, due to the limited residence time for crystal growth in the final segment, there are only small amounts of antisolvent added from the last injection points close to the end of PFC for 4, 5, and 6 addition points; see Case 2 (see Table 1). Figure 3 also shows that there is very little change in the CSD in the last few segments of the MSMA-PFC. Hence, when the total number of injection points is limited by

Figure 2. Optimization results of concentrations for four injection points in MSMA-PFC (solid line: solute concentration; dashed line: solute solubility).

practical constraints, the optimal locations of the injection points are rather important. With this in mind, in Case 3 both the locations and the distribution of the antisolvent additions are optimized and, as a result, a nearly constant supersaturation level is maintained up to the fourth addition at the 63rd segment, after which a large supersaturation level is also generated, analogous to that in Case 2. A further increase in the mean crystal size with a lower CV is obtained for Case 3 compared to Case 2 as shown in Table 1 and Figures 3 and 4. The latter indicates the strategy employed in Case 3 is more successful in controlling the nucleation rate throughout the whole tube, giving the largest mean crystal sizes and maintaining a CV well below the constraint of 0.30.

The effect of the total number of injection points on the mean crystal size at the exit of the MSMA-PFC is summarized in Table 1 and also depicted in Figure 5; Case 3 always produces the largest mean crystal size, and increasing the number of injection points beyond four only contributes marginally to improve the performance of the MSMA-PFC. Most importantly, consistent optimal locations for additions were found for Case 3, where nearly all of the antisolvent additions are added in the first half of the crystallizer tube. In

Figure 3. Optimization results of number-based CSD for four injection points in MSMA-PFC.

contrast, variations in the optimization results were obtained for equally spaced injection points of Cases 1 and 2, which were also reported by Ridder et al.24 In Cases 1 and 2 the antisolvent is added more uniformly along the length of the crystallizer, as the number of injection points increases; however, this strategy fails to provide enough antisolvent addition in the first half of the tube. Thus, the results of Figures 3 to 5 show that, for an MSMA-PFC, an optimization of both the location and amount of antisolvent additions should be taken into account to achieve the best flexible, efficient, and cost-effective design.

Figure 4. Optimization results of volume-based CSD for four injection points in MSMA-PFC (upper four: Case 1; center four: Case 2; lower four: Case 3).

For further comparison, the batch crystallization process with an optimal C-control strategy (a methodology to manipulate the supersaturation to trade-off the nucleation and crystal growth rates), using a seed (mean size 220 ^m), achieved a final crystal mean size of 556 ^m for a batch time of 120 min in ref 34. Therefore, the proposed design and optimization framework of MSMA-PFC, producing final crystals with good quality, shows the potential to implement an innovative continuous crystallizer design to replace existing batch crystallization processes.

In terms of the robustness analysis of the optimization results, Ridder et al.24,35 have thoroughly illustrated that, for an

unseeded antisolvent crystallization in an MSMA-PFC, the final product qualities are rather sensitive to uncertainties in the nucleation kinetics and the first antisolvent addition rate. Obviously, the nucleation kinetics plays a decisive role in an unseeded crystallization, as they determine the total number of crystals generated. Therefore, the first antisolvent addition is important since it generates the supersaturation which is responsible for the largest nucleation rates. For sufficiently long tubes, the nascent crystals continue to grow until the supersaturation is exhausted, which leads to a reduced sensitivity of the optimization results to crystal growth kinetics and the subsequent antisolvent additions. Similar trends are

Figure 5. Effect of the number of injection points on optimization results of MSMA-PFC.

also found in the current work. For example, comparisons of the robustness of the optimal results (m = 4) in Table 1 under uncertainties in the nucleation kinetic parameter kb are shown in Figure 6, where kb is multiplied by a parameter kuc to represent the uncertainty. It is found that similar sensitivities are obtained for Cases 2 and 3, as mentioned above, although the proposed optimization framework maintains a better performance over its counterparts. Therefore, it is much more desirable to reduce these sensitivities by introducing

proper seeding techniques or an elaborate design of a nucleator in the first PFC segment, which is an ongoing project in our group.

Dynamic Simulation. A dynamic simulation of the above MSMA-PFC design was conducted with the optimized four injection locations and antisolvent distributions of Case 3, as shown in Table 1. First, the MSMA-PFC is assumed to start up by injections of antisolvent into the crystallizer according to the optimal results of antisolvent addition, by which the crystallizer can be first filled and purged using the less expensive antisolvent (water). Then the fresh feed saturated solution, with a constant flow rate of 50 mL/min, is injected and mixed with the first antisolvent addition at the inlet of first PFC segment.

The evolution of the paracetamol concentration along the crystallizer axial is depicted in Figure 7. The discontinuities in the concentration are due to the dilution effect of the antisolvent additions. Apart from the transient response at the inlet for the first few minutes, the start-up of the MSMA-PFC is relatively smooth compared to a stirred-tank design (which typically requires more than four mean residence times to reach steady-state), approaching the steady-state just after the first mean residence time. The observed transient response

Figure 6. Robustness analyses of optimization results (m = 4) of MSMA-PFC with regard to nucleation uncertainties.

Figure 7. Dynamic simulation of main solute concentration in MSMA-PFC for Case 3 (m =4).

at the inlet is due to the high supersaturation generated by the instant mixing of fresh solution with the first antisolvent addition of 0.241 X 25 mL/min, leading to a relatively large nucleation effect and its propagation into the following segments, as is more clearly shown in Figures 8 and 9. The

Figure 8. Dynamic simulation of main solute concentration at different tube locations of Case 3 (m = 4) (+: after antisolvent addition; —: before antisolvent addition at addition locations).

initial peak at the solute concentration profile for z = 0.0+ m is due to the fast depletion of solute by spontaneous nucleation after the generation of a large amount of supersaturation at the inlet. After that event, the nucleation effect is reduced (see Figure 9 for z = 0.0+) as the concentration remains nearly constant and reaches its steady-state; the already generated nuclei in the inlet then serve as the seeding crystals in the remaining sections of the crystallizer.

The large numbers of nuclei generated at the entry pass through the following PFC segments, consuming supersaturation through growth and resulting in an initially slightly lower concentration than that at the steady-state condition, as shown in Figure 8 for different tube locations. Figure 9 shows that the effects of this initial burst of nuclei quickly vanish from the dynamic evolution of the CSD inside the tube. Although this high supersaturation only occurs at the beginning of the

Figure 9. Dynamic simulation of CSD in MSMA-PFC for Case 3 (m =

start-up process, its effect on the crystal size distribution could be important in practice, causing crash nucleation and fouling at the inlet.35 The use of a dynamic simulation illustrates the problems that might occur on start-up and allows design of a mitigation strategy, e.g. to initially provide some seeding into the fresh feed saturated solution at the beginning of the start-up process and then move more slowly back to an unseeded operation. For example, a 2.0% seeded fresh slurry was used to start up the MSMA-PFC with an optimal antisolvent distribution of Case 3 (m = 4), as shown in Figures 10 and

Tube length, z (m)

Figure 10. Main solute concentration at the inlet during seeded startup of MSMA-PFC for Case 3 (m = 4).

11 for the paracetamol concentration and CSD evolution at the inlet, respectively. Here, a log-normal seed crystal size

1000 0 Time, f (min)

Crystal size, L (цт)

Figure 11. Crystal size distribution at the inlet during seeded start-up of MSMA-PFC for Case 3 (m = 4).

distribution with a mean size of 60 ^m and a standard deviation of 1.5 was continuously seeded for 60 min. When compared to the unseeded start-up in Figures 7 and 9 (for z = 0.0+), the seeded crystallization over the first 60 min decreased the paracetamol concentration significantly, while initializing the crystallization at the inlet. After stopping the seed feed addition, the process slowly resumed to a high concentration at the inlet and also reached the same optimal steady state as Case 3 (m = 4), without causing an initial large burst of nucleation.

In summary, it is obvious that the MSMA-PFC, by optimized process design, shows the advantages of a quick start-up and potentially simplified process control needed to reach a steady-state operation, compared to that observed in cascaded multistage MSMPR crystallizers.16,36

■ CONCLUSIONS AND FUTURE WORK

The multisegment, multiaddition plug-flow crystallizer (MSMA-PFC) has shown potential benefits for the design of optimized continuous crystallization processes, which can replace existing batch operations. The current work has

extended the previous work in mathematical modeling, design, and optimization of an MSMA-PFC, proposing a conceptual design based on a number of standardized modular units, an optimization framework for finding the best locations and amounts of antisolvent additions, and a dynamic simulation to study its start-up. Improvements to the previous optimization frameworks reported in the literature by Alvarez and Myerson1 and Ridder et al.24 were illustrated, showing that larger mean crystal sizes could be obtained, with the CV maintained below a target level. The proposed design framework avoids the formation of multiple large nucleation events, which give rise to multimodal and broad size distributions of the crystal product. The method is quite general and can be adapted to take into account different definitions of the process objective and to target a variety of product quality attributes. Nevertheless, the simulations require input in the form of kinetic rate laws, which should be obtained from experiments conducted under the relevant flow conditions. Future work will consider the effects of initial seeding at the inlet, additions of fresh solution along the tube, and a temperature profile for a combined cooling and antisolvent crystallization process. Further extensions of the plug-flow crystallizer to a continuous oscillatory baffled crystallizer (COBC) with some axial dispersion are also under investigation.

■ AUTHOR INFORMATION Corresponding Authors

*E-mail: Q.Su@lboro.ac.uk.

*E-mail: C.D.Rielly@lboro.ac.uk. Tel: +44 (0) 1509 222504. Funding

Engineering and Physics Science Research Council of United Kingdom (UK EPSRC EP/K014250/1)

The authors declare no competing financial interest.

■ ACKNOWLEDGMENTS

This work was performed within the UK EPSRC funded project (EP/K014250/1) 'Intelligent Decision Support and Control Technologies for Continuous Manufacturing and Crystallization of Pharmaceuticals and Fine Chemicals' (ICT-CMAC). The authors would like to acknowledge financial support from EPSRC, AstraZeneca, and GSK. The authors are also grateful for useful discussions with industrial partners from AstraZeneca, GlaxoSmithKline, Mettler-Toledo, Perceptive Engineering, and Process Systems Enterprise.

■ ABBREVIATIONS

ATR-FTIR, attenuated total reflectance Fourier-transform infrared spectroscopy

COBC, continuous oscillatory baffled crystallizer CSD, crystal size distribution CV, coefficient of variation FBRM, focused-beam reflectance measurement FVM, finite volume method

MINLP, mixed integer nonlinear programming problem

MOM, method of moments

MSMA, multisegment multiaddition

MSMPR, mixed-suspension mixed-product-removal

ODE, ordinary differential equation

PAT, process analytical technology

PDE, partial differential equation

PFC, plug-flow crystallizer

PSE, process systems engineering PVM, particle vision measurement SDC, simultaneous design and control

■ REFERENCES

(1) Alvarez, A. J.; Myerson, A. S. Cryst. Growth Des. 2010, 10, 22192228.

(2) Chen, J.; Sarma, B.; Evans, J. M. B.; Myerson, A. S. Cryst. Growth Des. 2011, 11, 887-895.

(3) Majumder, A.; Nagy, Z. K. ATChE J. 2013, 59, 4582-4594.

(4) Plumb, K. Chem. Eng. Res. Des. 2005, 83, 730-738.

(5) Reklaitis, G. V.; Khinast, J.; Muzzio, F. Chem. Eng. Sci. 2010, 65, iv-vii.

(6) Kessel, M. Nat. Biotechnol. 2011, 29, 27-33.

(7) Mascia, S.; Heider, P. L.; Zhang, H.; Lakerveld, R.; Benyahia, B.; Barton, P. I.; Braatz, R. D.; Cooney, C. L.; Evans, J.; Jamison, T. F.; Jensen, K. F.; Myerson, A. S.; Trout, B. L. Angew. Chem., Int. Ed. 2013, 52, 12359-12363.

(8) Buchholz, S. Chem. Eng. Process. 2010, 49, 993-995.

(9) Aksu, B.; De Beer, T.; Folestad, S.; Ketolainen, J.; Linden, H.; Lopes, J. A.; de Matas, M.; Oostra, W.; Rantanen, J.; Weimer, M. Eur. J. Pharm. Sci. 2012, 47, 402-405.

(10) Banholzer, W. F.; Jones, M. E. ATChE J. 2013, 59, 2708-2720.

(11) Griffin, D. W.; Mellichamp, D. A.; Doherty, M. F. Chem. Eng. Sci. 2010, 65, 5770-5780.

(12) Alvarez, A. J.; Singh, A.; Myerson, A. S. Cryst. Growth Des. 2011, 11, 4392-4400.

(13) Quon, J. L.; Zhang, H.; Alvarez, A.; Evans, J.; Myerson, A. S.; Trout, B. L. Cryst. Growth Des. 2012, 12, 3036-3044.

(14) Zhang, H.; Quon, J.; Alvarez, A. J.; Evans, J.; Myerson, A. S.; Trout, B. Org. Process Res. Dev. 2012, 16, 915-924.

(15) Sen, M.; Rogers, A.; Singh, R.; Chaudhury, A.; John, J.; Ierapetritou, M. G.; Ramachandran, R Chem. Eng. Sci. 2013, 102, 5666.

(16) Su, Q.; Nagy, Z. K.; Rielly, C. D. Chem. Eng. Process. 2015, 89, 41-53.

(17) Lawton, S.; Steele, G.; Shering, P.; Zhao, L.; Laird, I.; Ni, X. W. Org. Process Res. Dev. 2009, 13, 1357-1363.

(18) Eder, R J. P.; Radl, S.; Schmitt, E.; Innerhofer, S.; Maier, M.; Gruber-Woelfler, H.; Khinast, J. G. Cryst. Growth Des. 2010, 10, 22472257.

(19) Eder, R. J. P.; Schmitt, E. K.; Grill, J.; Radl, S.; Gruber-Woelfler, H.; Khinast, J. G. Cryst. Res. Technol. 2011, 46, 227-237.

(20) Ferguson, S.; Morris, G.; Hao, H.; Barrett, M.; Glennon, B. Chem. Eng. Sci. 2012, 77, 105-111.

(21) Zhou, L.; Su, M.; Benyahia, B.; Singh, A.; Barton, P. I.; Trout, B. L.; Myerson, A. S.; Braatz, R D. ATChE J. 2013, 59, 1308-1321.

(22) Lakerveld, R.; Kramer, H. J. M.; Stankiewicz, A. I.; Grievink, J. Chem. Eng. Process. 2010, 49, 979-991.

(23) Vetter, T.; Burcham, C. L.; Doherty, M. F. Chem. Eng. Sci. 2014, 106, 167-180.

(24) Ridder, B. J.; Majumder, A.; Nagy, Z. K. Tnd. Eng. Chem. Res. 2014, 53, 4387-4397.

(25) Kwon, J. S.; Nayhouse, M.; Orkoulas, G.; Christofides, P. D. Chem. Eng. Sci. 2014, 119, 30-39.

(26) Byrn S., Futran M., Thomas H., Jayjock E., Maron N., Meyer R F., Myerson A. S., Thien M. P., Trout B. L. International Symposium on Continuous Manufacturing of Pharmaceuticals. May 20-21, 2014.

(27) Leuenberger, H. Eur. J. Pharm. Biopharm. 2001, 52, 289-296.

(28) Su, Q. L.; Chiu, M. S.; Braatz, R D. ATChE J. 2014, 60, 28282838.

(29) Qamar, S.; Elsner, M. P.; Angelov, I. A.; Warnecke, G.; SeidelMorgenstern, A. Comput. Chem. Eng. 2006, 30, 1119-1131.

(30) Mesbah, A.; Kramer, H. J. M.; Huesman, A. E. M.; Van den Hof, P. M. J. Chem. Eng. Sci. 2009, 64, 4262-4277.

(31) Deb, K.; Pratap, A.; Agarwal, S.; Meyarivan, T. TEEE Trans. Evol. Comput. 2002, 6, 182-197.

(32) Su, Q. L.; Braatz, R D.; Chiu, M. S. J. Process Control 2014, 24, 415-421.

(33) Abel M. J. PhD Thesis, Massachusetts Institute of Technology, 2009.

(34) Woo, X. Y.; Nagy, Z. K.; Tan, R B. H.; Braatz, R D. Cryst. Growth Des. 2009, 9, 182-191.

(35) Ridder, B. J., Majumder, A, Nagy, Z. K. Proceedings of the American Control Conference (ACC 2014); IEEE Press: Piscataway, NY, Seattle, USA, 2014.

(36) Yang, Y.; Nagy, Z. K. Chem. Eng. Sci. 2015, 127, 362-373.