Mathematical modeling and experimental validation of continuous slug-flow tubular crystallization with ultrasonication-induced nucleation and spatially varying temperature

Continuous slug-flow tubular crystallization has been explored by several research groups in academia and industry as a way to produce crystals while having low capital equipment costs. In this crystallization type, slugs of slurry and gas consecutively travel through a tube, with a high degree of mixing and temperature uniformity within each slug. This article presents an experimental system for slug-flow tubular crystallization that employs a spatial temperature profile and directed non-contact ultrasonication to induce primary nucleation to enable the generation of a wide variety of crystal size distributions. The crystal size distributions are compared for data collected from a full-factorial experimental design (27 experiments in total) to predictions from a population balance model that includes the effects of ultrasonication on primary nucleation. This population balance model for tubular crystallization is the first that incorporates the effects of ultrasonication and dissolution on the crystal size distribution. The crystal size distributions are reasonably consistent with the model, within 20% prediction error, for all experiments in which the spatial temperature profile is monotonically decreasing and at low to moderate supersaturation. Potential causes for weaker agreement for other experiments are discussed. © 2021 The Author(s). Published by Elsevier B.V. on behalf of Institution of Chemical Engineers. This is an open access article under the CC BY license (http://creativecommons. org/licenses/by/4.0/).


Introduction
While industrial crystallization is predominately performed under (semi)batch conditions in the pharmaceutical industry, much interest has been in increasing operational flexibil-2017). An early crystallizer design employed a pressure-driven pump to combine a pentaerythritol tetranitrate-acetone solution with an aqueous anti-solvent through a short Kenics-type static mixer positioned directly above a filter screen (Rivera and Randolph, 1978). Myerson and coworkers later demonstrated a multi-segment plug flow crystallizer constructed from four jacketed glass tubes packed with Kenics mixers (Alvarez and Myerson, 2010), which was demonstrated for three commercial active pharmaceutical ingredients (APIs) using multi-point anti-solvent addition to mediate the process. Lemaître and colleagues (Lemaître et al., 1997;Vacassy et al., 2004) developed a segmented flow tubular reactor (SFTR) using 10 m of 4 mm ID PTFE tubing. Peristaltic pumps were used to combine two aqueous solutions with air through a static mixer, forming alternating slugs of liquid and gas that traveled through the tube. Given its minimalist hardware and simplicity of operation, the original SFTR underwent some rapid re-designs, enhancing flexibility as a potential platform technology. Eder and colleagues (Eder et al., 2010;Eder et al., 2012) probed the application of continuous seeding to SFTRs with demonstration to acetylsalicylic acid (ASA) crystallization. These studies demonstrated temperature-gradient SFTRs (TB-SFTR), using a single length of silicone tubing coiled into multiple thermal baths set at predetermined temperatures, and showed the nucleation could be induced in particle-free solutions of ASA using small-scale ultrasonication baths. Industry-led experiments applied a similar mesoscale version of this system to the mitigation of clogging during flow-through pH-swing crystallizations (Furuta et al., 2016). TB-SFTRs were later coupled with flow-through seed generation methods that leveraged rapid temperature shifts and directed spatially localized ultrasonic irradiation (Jiang et al., 2014a;Jiang et al., 2015). Flow conditions and tubing dimensions were identified that supported spontaneous, stable slug formation in the absence of multiple upstream mixers comprising more than two inlets (Jiang et al., 2014b). Building upon these advances, this article presents an experimental system for slug-flow tubular crystallization that employs a spatial temperature profile implemented using modular counter-current heat exchangers (CCHEXs) comprising fully disposable FDA-compliant flow paths of variable length. This study employs L-asparagine monohydrate (LAM) as the model compound (Jiang et al., 2012;Rasche et al., 2016;Jiang and Braatz, 2019), due to its representative solubility and operational difficulty in generating uniformly sized crystals. Primary nucleation is induced using an ultrasonication probe. The crystal nuclei in flow are combined with an air flow at flow rates that spontaneously form alternatively slugs of slurry and air which traverse the heat exchangers. Experimental results are presented for a 3 3 full-factorial experimental design for the crystallization of L-asparagine monohydrate in aqueous solution. The images of product crystals are analyzed using an in-house analysis algorithm for object segmentation, counting, and sizing to quantitatively investigate the effects of spatial temperature profiles and ultrasonication to the crystal size distribution (CSD). The experimentally determined CSDs are used to estimate crystallization kinetics and then compared to predictions from a population balance model (PBM) of the system to experimentally validate the developed model. This population balance model for tubular crystallization is the first to incorporate the effects of ultrasonication and dissolution on the crystal size distribution.

Materials
L-asparagine monohydrate (LAM, A8381; Sigma-Aldrich, St. Louis, MO), the commercially available material in powder form, was used as the model small-molecule solute with deionized (DI) water ( ≥ 18 M cm) for all crystallization experiments. Fig. 1 is a schematic of the system, which comprises four modular shell and tube CCHEXs arranged in series. Crystallization occurs only in the tube portion of each module, and coolant (DI water) is confined to the shell of each CCHEX. Coolant is driven through the shell sides using pairs of peristaltic pumps. One pump in each pair withdraws cold DI water from a common recirculating batch and a second common recirculating batch provides hot DI water. This scheme allows the flow rate and temperature of the coolant within a given CCHEX shell to be maintained independently of all other modules. At the inlet, a 'pulse' of thermally equilibrated LAM solution is supplied to the crystallizer from a well-mixed, threenecked, round-bottom flask. The LAM solution is maintained at temperatures in excess of 70 • C for many hours, while preventing evaporation of the solvent, thereby ensuring that all experiments are performed with constant-concentration, particle-free feedstocks. After entering the crystallizer, the particle-free aqueous phase traverses a short segment of insulated tubing before transiting through the first CCHEX (the 'pre-HEX'; 2 m) and passing beneath a Sonics 13 mm ultrasonication probe (630-0219, ultrasonic processor Sonics VCX750). The purpose of the pre-HEX is to rapidly cool the LAM feed from its undersaturated condition at the system inlet to a temperature amenable to the formation of nuclei. The implementation of temperature control upstream of an ultrasonication-based nucleation stage is an improvement over past designs (Eder et al., 2010;Eder et al., 2012;Jiang et al., 2015;Jiang et al., 2014b;Jiang et al., 2017b;Eder et al., 2011). Following ultrasonication, the LAM solution is mixed with filtered air to form stable liquid slugs and immediately enters the second heat exchanger module ('HEX 1'). Slugs travel with gravity along the length (11 m) of the coil comprising HEX 1 and traverse the third ('HEX 2') and fourth ('HEX 3') shell and tube modules (both 14.8 m long) in a manner identical to their passage through HEX 1. Stable slug flow was induced on the tube side of the crystallizer using a method similar to that described by Jiang et al. (Jiang et al., 2015;Jiang et al., 2014b) Combinations of air and liquid flow rates that yielded stable slugs were analyzed theroetically (Pirkle et al., 2020) and identified empirically (Jiang et al., 2015;Jiang et al., 2014b), and all experiments were conducted under air and liquid flow rates of 6 mL/min.

Experimental apparatus
Seven ungrounded transition-joint thermocouples provide real-time measurements of temperature at critical positions along the length of the crystallizer. Shell-side temperature readings are acquired by placing a thermocouple in direct contact with the flowing coolant. By contrast, to prevent errant nucleation, the temperature of LAM slugs is measured indirectly by fastening thermocouples to the outer wall of the tubing comprising tube-side of the system. An empirically- derived linear correlation implemented as part of the custom software required for overall system operation is used to calculate internal LAM slug temperatures from temperatures measured at the wall at the time of acquisition. Temperature setpoint tracking is achieved via the implementation of proportional-integral (PI) controllers that ratiometrically vary the relative flow rates of each pump in a given set. All system hardware was operated via custom software written and implemented using MATLAB 2016b (v9.1.0.441655; The Math-Works, Inc.). The software includes standard signal processing subroutines. Leveraging the built-in timer class provided by MATLAB, the same series of methods run at 1 s intervals for the duration of an experiment: (1) all thermocouples are probed in a single call to a data acquisition module, (2) temperature signals are filtered, (3) PI controllers specify new shell-side flow rates, (4) shell-side pumps are actuated sequentially, and (5) new thermocouple and pump data are appended to the historian. Details on the control designs and their experimental validation are given in an Appendix. MATLAB timers were used to provide transparent low-level task scheduling functionality.

Full-factorial design
A 3 3 full-factorial design was carried out with the aim of probing the effects of ultrasonication power output (probe vibrational amplitude, %) and two upstream temperatures (pre-HEX and HEX 1, Table 1). The probe vibrational amplitude was 114 microns (0.0045 inch) with the amplitude control set at 100%. Shell-side temperature and flow rate setpoints were selected with the aid of the custom thermodynamic model to ensure that the slug-side temperatures in Table 1, and to mitigate temperature crossovers near the center of each heat exchanger. All LAM slugs entered the crystallizer at a concentration of 0.1058 g LAM / g dH 2 O (T sat = ∼53 • C), and the temperature of the solution is 93-95 • C. The pre-HEX temperature set points were selected to correspond to the saturated, metastable, and supersaturated regions of the LAM solubility diagram (Jiang et al., 2012).

Crystal size distribution (CSD) acquisition
The collected LAM slug was subsequently vortexed to induce settling in the suspension of LAM crystals. Thirty-two transmitted light photographs were taken using a stereomicroscope and color USB camera at maximum intensity and a magnification of 12.22 m/px for each experiment.
The custom image analysis algorithm capable of detecting the size of all crystals is schematized in Fig. 2. The process was implemented using MATLAB 2018a (The MathWorks, Inc.). The color photographs were converted to grayscale images as a weighted sum of the R, G, B components. An alternating sequential filtering operation was employed to estimate and remove variations in background intensity, returning the dilated image. The image of objects in the foreground are obtained by subtracting the intensity of background and scaling the divided image. The label matrix which assigns 0 to the background and different positive integers to each object was then converted into an RGB color image to visualize successfully segmented crystals. In the conversion, pixels labeled 1 form one object; pixels labeled 2 form a second object; and so on. The distance between each adjoining pixel around the border of an object was calculated to compute the size of each object. By applying the same procedure to all experiments, the lumped crystal size distribution was obtained.
Given that all experiments generating LAM CSDs were conducted at steady state and operated under saturated outlet concentrations, the mass balance was used to calibrate CSD in a reasonable range, where the value of crystal density ( c ) is 5.6966 × 10 3 kg/m 3 , the dimensionless shape factor (k v ) is 0.39551, and the difference in solution concentration between inlet and outlet, C in − C out (T out ) is 84 kg/m 3 , which is the experimental conditions (Mozdzierz, 2018). The dimensionless scaling factor ˛ corrects for the sampling bias associated with acquiring crystal images in only one focal plane andf (z N , L) is the observed CSD at slug segment from outlet. The true CSD f (z N , L) is obtained by mul-  tiplying a scaling factor and the observed CSD. For discrete observed data, the mass balance Eq. (1) is modified to Outliers were removed when fitting ˛ to CSD data.

Thermal model
The thermal model in this article is a variant of those reported in the literature (Estel et al., 2000;Arbaoui et al., 2007;Ansari and Mortazavi, 2006). Energy balances of the shell-and-tube crystallizer explicitly account for all significant heat fluxes within the system. All experiments were run in quasi-steady conditions so the time derivative terms are dropped, and separate energy balances are written around the control volume within the shell and tube sides (also respectively 'cold' and 'hot') of the CCHEX. Three additional assumptions are introduced: (1) Fluid heat capacity C p and density are constants for both the shell and tube sides of the crystallizer.
(2) The heat transfer between the slugs can be neglected because the indirect heat transfer coefficient between slugs and coolant is far larger than slugs and air, indicating that cooling water plays a dominant role in the thermal dynamics of the system. (3) The temperature distribution in the radial direction in slug tube is neglected.
Generally, the concentration and temperature are uniformly distributed in a single slug, which is supported by theoretical analysis of the flow within the slugs (Pirkle et al., 2020). The outlet temperature of the shell-side fluid should be predicted and updated in any numerical simulation. When the predicted value at the shell-side inlet matches the observed data, the calculation is stopped. The discretized tube-and shell-side energy balances can be written as by employing finite volume method, where the subscript 'slug' refers to the discretized tube-side slug flow, the subscript 'shell' refers to shell-side flow,ṁ z is the mass flow rate, andQ is the heat of crystallization or conductive heat transfer rate at steady state. All material properties are treated as constants, and a cylindrical coordinate system with the +z axis oriented along the direction of tube-side flow are assumed. The heat fluxes in the above equations are described by standard heat transfer relationships: where UA is the product of the overall heat transfer coefficient and the surface area for heat transfer, m is the mass of a slug within the tube side, is residence time, and H c and M c are the heat of crystallization and molecular weight of LAM which is −35,700 J/mol and 0.150 kg/mol, respectively (Rasche et al., 2016). The convective heat transfer coefficient used to derive UA is determined by the convection correlations for flow in a circular tube (Bergman et al., 2011). The uninsulated connecting regions between heat exchangers and interaction between a heat exchanger and the ambient environment were also considered. The above complete set of linear equations can be solved using any linear solver.

Population balance model
The fact that the liquid slugs can be treated as well-mixed batch crystallizers has been shown by theoretical analysis and in experiments (Jiang et al., 2014b;Pirkle et al., 2020). The population balance model (PBM) (Rasche et al., 2016) ∂f ( , L) ∂ + ∂ G S, T; Â g f ( , L) describes the evolution of the crystal size distribution f ( , L) at residence time , where L is the characteristic crystal length in units of m, ı(L) is the Dirac delta function, G S, T; Â g and B (S, T; Â b ) are the crystal growth and nucleation rates respectively, S = C/C sat is the supersaturation ratio, T is the temperature, and Â g and Â b are vectors of growth and nucleation kinetic parameters respectively. The growth and nucleation are considered as the only dominant phenomena in the continuous crystallizer (i.e., no breakage or agglomeration). Primary nucleation occurs from solution at high supersaturation when no crystals are present. Although the two-step primary nucleation mechanism has been shown to be applicable for some systems, the nucleation rate expression from classical nucleation theory (CNT) (Erdemir et al., 2009;Vekilov, 2010;Zahn, 2015) has been observed to able to fit the observed experimental data for most measured primary nucleation rates, and is used here. CNT takes into account the free energy cost of the nucleus at the top of the nucleation barrier, the rate at which molecules can attach to the nucleus, and the probability that a nucleus at the top of the barrier will go on to form the new phase rather than dissolve (Zahn, 2015;Wirges, 1986). This theory gives the nucleation rate as where is the surface tension of nucleus and is the molecular volume of solute (Zahn, 2015;Wirges, 1986). The pre-exponential factor A is a function of molecular parameters that can be derived by the rates of condensation and dissociation, which results in the expression where Á is the core repulsion range and n 0 denotes the number density of the condensing molecules in solution at an infinite distance from the crystal (Ruckenstein and Nowakowski, 1990). This expression takes into account that the time scales of the formation and dissociation of smaller particles are much shorter than the characteristic times for formation and growth of the critical sized nucleus (Ruckenstein and Nowakowski, 1990). Multiple numerical rate expressions have been proposed for extended classical nucleation theory to model the effects of ultrasonication, taking into account such effects as the surface area of the bubbles, variation in the surface free energy, the adsorption properties of the forming cluster, and the chemical potential for nucleation induction (Kurotani et al., 2009;Miyasaka et al., 2006;Kordylla et al., 2009;Virone et al., 2006). Here we apply a similar nucleation expression to that derived in the most rigorously justified paper (Virone et al., 2006), where the three parameters A ins. , B ins. , and C ins. depend on the particular ultrasonication environment and x is the amplitude set point of ultrasonication. This equation reduces to the nucleation rate expression from classical nucleation theory when the ultrasonication power is zero.
Secondary nucleation is the dominant nucleation mechanism in most industrial crystallizations (Agrawal and Paterson, 2015). As the nucleation occurs under the presence of pre-existing crystals, factors such as the supersaturation, rate of cooling, liquid velocity, and impurity concentration affect the generation of secondary nuclei (Myerson, 2002;Denk and Botsaris, 1972;Omar and Rohani, 2017). Slug-flow crystallization lacks a mixer that would induce secondary nucleation in a mixed tank, and secondary nucleation is only observed at very high supersaturations (Jiang et al., 2015;Jiang et al., 2014b), indicating that the mechanism for slug-flow crystallization is different than in mixed tanks. The secondary nucleation mechanism is not well understood slug-flow systems, so the empirical power law model (Myerson, 2002) with the suspension density are usually used to fit experimental data. Fundamentally, the crystal growth rate is expressed as the linear growth rate of a particular face. The rate of growth is mainly controlled by the cooling rate of the solution, the anti-solvent concentration, pH, and the concentration of any impurities (Garnier et al., 2002;Meyer, 1984;Sangwal, 1996). A popular mechanistic model for crystal growth assumes a screw dislocation with a surface developing growth of layer
Crystal dissolution occurs during undersaturation (C < C sat (T)). Insertion of the lower bound on the Sherwood number of 2 (e.g., see Ref. (Deen, 2012)) into the standard dissolution rate expression (e.g., Refs. (Hintz and Johnson, 1989;Lu et al., 1993;Carstensen and Dali, 1999)) and insertion of that expression into a solute mass balance shows that the dissolution time scale is faster than ∼1 s, which is much less than the residence time of 10 min when the slugs are in HEX 1. As such, the dissolution is complete by the outlet of the tube leaving HEX 1, and the effect of dissolution can be simply calculated by reducing the total mass of crystals and shifting the distribution to the range of smaller sized crystals.
Applying first-order discretization to the characteristic crystal length L and equations in Table 2, the partial differential PBM (Eq. (8)) reduces to the set of ordinary differential equations (ODEs) with the initial condition The finite volume method for characteristic crystal length (L) is used to convert the PBM. Since the mass reduction of total solute induced by nucleation is negligible, the growth rate is solely involved in mass balance equations of LAM solute, where 2 is the second-order moment proportional to the total surface area of crystals, which is L 2 n f n L for large N.
Eqs. (14) and (17) are solved simultaneously with boundary conditions given by Eqs. (16) and (19) over the length of the tubular crystallizer. Based on the spatial trajectory of temperature obtained from a thermodynamic model, the trend of LAM solubility (Jiang et al., 2012) supersaturation ratio, and nucleation rate and growth rate for z-axis direction are derived and discretized in the same way as the population equation and the concentration equation.

Experimental crystal size distributions for varying sonication power and temperature profiles
The CSDs obtained from the image analysis algorithm in Section 2.4 across all 27 experiments are summarized in the box-and-whisker plot in Fig. 3. The experiments, which span a wide range of spatial temperature profiles, result in a wide range of skewness and width of the CSDs (cf. Experiment 10 and 5).
For the strictly monotonically decreasing spatial temperature profiles with equal outlet temperatures of HEX 1 and the most gradual temperature reduction (Exps. 19-21, Table 1), the temperature of pre-HEX determines the level of supersaturation at which primary nucleation is induced. The mean and standard deviation of the particle size distribution in Exp. 21 at high ultrasonication power are 39% and 42% smaller than for Exp. 19 which has no ultrasonication, and 16% and 9% smaller than for Exp. 20 of intermediate ultasonication power, respectively (Table S1). As such, under these conditions, using ultrasonication resulted in smaller median particle sizes due to increased nucleation, as observed in past studies (Jiang et al., 2015;Ruecroft et al., 2005).
In the range of temperature of HEX 1 set points from 42.6 • C to 51.2 • C, some of the experiments have increasing temperature for part of the tube (Table 1). When the slug temperature is increased, some of the nucleated crystals formed upstream can subsequently dissolve until the solution concentration is equal to the solubility, resulting in smaller product crystals. For example, Exps. 13-15 have same operating conditions as Exps. 16-18 except for the HEX 1 temperature rise in the latter. In this case, the mean particle sizes in Exps. 16-18 with temperature rise are 11%, 24%, and 15% smaller than for Exps. 13-15 without temperature rise. Table 3 summarizes the results of kinetic parameter estimation of the mathematical model of the continuous crystallizer. The CSDs from experiments operated in the ultrasonicationfree environment (Exps. 1,4,7,10,13,16,19,22,and 25) were fitted without a term for ultrasound-induced nucleation. For the overall experiments, the estimated cumulative distribution function (ECDF) of observations was established by shape-preserving piecewise cubic interpolation with observed populations of crystals on a regular interval (12.22 m) and   precisely compared with model predicted ECDF point by point. The kinetic parameters that result in minimum error between the ECDFs were numerically determined by using a parallel implementation of a genetic algorithm (GA). Among the parameters, the multipliers of primary nucleation (Â pri.,1 and Â ins.,1 ) show the ultrasonication-induced nucleation is about 500-1800 times higher than the stochastically generated nucleation. The exponent values (Â sec.,2 and Â sec.,3 ) of the sus-pension density of crystals and the concentration difference between bulk and saturated solution in secondary nucleation are reasonable when compared with an earlier correlation (Myerson, 2002). Fig. 4 summarizes the comparison between experimental ECDF and model-predicted ECFD for all of the experiments. Most experiments have a coefficient of determination (Rsquared) near 1, indicating high agreement between model predictions and experiments, whereas other experiments have high R-squared and poor agreement. Experiments in which the temperature increased during the experiment, which can induce dissolution, are not fitted well with observed data when accompanied with ultrasonication-induced nucleation. The predicted mean length of crystals for Exps. 9 and 18 are 41.0% and 31.8% smaller than observed mean length, which means that the effect of ultrasonication on primary nucleation is overestimated when followed by a temperature rise. Potential causes for this shift are that (1) the highest ultrasonication power could produce a strongly localized temperature rise that is not included in the model and (2) the population balance model (Eq. (15)) does not account for change in the shape of crystals that occurs during dissolution (Snyder and Doherty, 2007). The measurements from experiments with monotonically decreasing temperature, in which dissolution cannot occur, and lower supersaturation correspond well to the model predictions. When comparing the predicted and observed data of crystals in a parity plot (Fig. 5), most of the model results are quite consistent with the data, within 20% prediction error.

Comparison of model and experimental crystal size distributions
The spatial temperature profiles and CSD of representative experiments are in Fig. 6 (the plots for all experiments are in the Supplementary Material). The measured temperatures from installed thermocouples are well described by the thermal model for all experiments (Fig. 6ace). The model agreement for the estimated cumulative distribution function is intermediate for Exp. 1, which had the highest supersaturation (Fig. 6b). The largest deviation in the estimated cumulative distribution function occurs for Exp. 9 (Fig. 6d), which has the most extreme ultrasonication power, the temperature drop in pre-HEX, and temperature rise in HEX 1 ( Table 1). The nuclei vigorously induced by the higher power of ultrasonication are dissolved, and the effect of ultrasonication is canceled out. In this case, the model does not sufficiently describe the dissolution process, and the number of smaller size crystals in the distribution are overestimated by the model. In contrast, the model agreement with experiments is very high for Exp. 19 (Fig. 6f), in which the ultrasonication power is 0% and the spatial temperature profile decreases the most gradually.

Conclusions
A continuous tubular crystallizer is presented that exploits modular counter-current heat exchangers and indirect ultrasonication to manipulate the nucleation and growth processes. Twenty-seven experiments that crystallize L-arginine monohydrate are performed under hydrodynamically stable slug flow. Model parameters in a population balance model with mechanistic nucleation and growth kinetics are fitted to the results of a full-factorial design of experiments. From the crystal size distribution obtained from the image analysis algorithm (Section 2.4), the effects of the spatial temperature profile along the tube and the indirect ultrasonication power on the crystal generation are observed.
Primary nucleation is effectively induced by ultrasound for lower supersaturations, which is consistent with past studies. Additionally, the product crystals generated by a spatial profile in which the temperature rises in the first heat exchanger results in the partial dissolution of the crystals nucleated upstream and in smaller product crystals (Table S1).
With the estimated parameters, the mathematical model reasonably predicts the crystal size distribution, within 20% prediction error, for monotonically decreasing spatial temperature profiles. The population balance model incorporates the effect of ultrasonication on primary nucleation kinetics and dissolution, but with weaker agreement between model predictions and data at extreme ultrasonication power and temperature reduction and rise which would maximize dissolution.
In addition to presenting the results of a full factorial design of experiments, this article presents the first population balance model for tubular crystallization to incorporate the effects of ultrasonication and dissolution on the crystal size distribution. The comparisons between model predictions and experiments motivate future development for improving the modeling of the effects of very high ultrasonication and dissolution.

Declaration of interests
None.

Declaration of Competing Interest
The authors report no declarations of interest.   The control parameters were then manually tuned to achieve robust performance over the range of operating flow rates predicted for each module, with the final tuning parameters and operating ranges listed in Table A2. Fig. A1 depicts the results of disturbance rejection studies designed to demonstrate how the control system responds to sudden changes in hot reservoir temperature (Fig. 1). The disturbances are deliberately of substantially greater magnitude than would ever be experience in crystallization operations, as a stress test for the platform (even for a complete coolant bath failure, the temperature of the hot or cold coolant would change only gradually, ∼1-3 • C/min). In spite of operating at a modest shell-side flow rate of 130 mL/min, the aggressively tuned pre-HEX module rejected disturbances of ∼−40% (dark grey dashed line) and ∼115% (light grey dashed line) within 200 s. The PI controllers governing HEX 1 and HEX 2 responded similarly well to the ∼−40% disturbance. As anticipated, the worst overall performance was for HEX 1, which was designed to operate at the lowest shell-side flow rate under the least aggressive control. Specifically, HEX 1 required nearly 10 min to fully reject the ∼115% disturbance. Comparable data and conclusions apply to the cold reservoir disturbance rejection tests (Fig. A2).