Batch Cooling Crystallization of Plate-like Crystals : A Simulation Study

In the paper, batch cooling crystallization of plate-like crystals is investigated by numerical experimentation. The population balance approach is used to describe a seeded crystallization process, developing a bivariate morphological population balance model including nucleation and size dependent crystal growth. The population balance equation is reduced to a system of equations for the bivariate mixed moments of internal variables, which are solved by means of the quadrature method of moments. A sensitivity analyses is carried out to investigate the effects of the process parameters on dynamic evolution of particulate properties including the characteristic mean particle size and shape.


Introduction
From practical point of view, crystal shape and size define properties like specific surface area and dissolution rate.They affect also some biological properties and have major influence on the downstream operations during the manufacturing like filtering, washing or drying.As a consequence, one of the recent challenges in particulate science is the simultaneous shape and size manipulation of the crystals during the crystallization process [1].The shape often differs from the simple geometries like cube or sphere because of the different growth rates of specific crystal faces [2].The growth rate in each direction can be an individual function of temperature, supersaturation and even the size of the growing face.In this context, the shape and size of an individual crystal is continuously changing during the crystallization process.Moreover, nucleation occurs simultaneously with the crystal growth that beside of the material specific properties depend on the temperature, supersaturation, concentrations, surface of the existing particles, and even on the stirring energy [3].The time evolution of crystal size distribution is sensitive to the initial supersaturation, temperature and the seed properties, i.e. to the quantity and size distribution of seeds.Although some advanced techniques exist for in situ measurement of the particulate properties which undoubtedly eases the data gathering and the development of on-line control strategies [4], currently these would be too expensive for a detailed system investigation.Therefore, mathematical modelling appears to be a useful tool for analysing these processes.
In the last decade, the population balances were extended to multivariate cases, making possible to model not only the size but the shape evolution in time of the crystal populations.These types of population balances are the so called morphological population balances and are typically able to reconstruct the particle shape at each moment of simulation.They are useful for kinetic estimations followed by solving shape optimization and/or control problems and are especially handful for the investigation of the crystallization process itself.
The morphological population balances were often applied to model and simulate crystallization of rod-like or needleshape crystals [4,5], but only several studies were reported related to simulation of crystallization of plate-like crystals [2,3], despite that crystals of such shape can also be play significant role in crystallization practice [6][7][8].
The aim of the current work is to develop a 2D population balance model for seeded batch crystallization of plate-like crystals, and to present and analyse the results of numerical experimentation carried out for investigating the effects of process parameters on the system behaviour and the product properties.

2D population balance model
The simulated experiments were carried out as follows: a hot temperature (T h ) saturated solution was obtained which was carefully cooled until a given temperature at which the seed crystals were added (T s ).In order to avoid the intensive primary nucleation the seeding temperature was kept in metastable zone, however, the crystal surface induced secondary nucleation may appear in these conditions also.The seed crystals are added with well specified quantity and quality conditions, i.e. with size distribution.After seeding, the solution was cooled with a constant cooling rate (cr) until reaching the final temperature (T f ) which was the end of the experiment.Evolution of the concentration of solute and of the solubility as a function of temperature is shown in Fig. 1, where the temperature of solution was varied linearly in time from 35 o C to 25 o C. It is seen that the supersaturation is very small during the course of the process.The plate-like shape can be approximated as a body which has two dimensions of significantly higher size than the third one as it is shown in Fig. 2. For the sake of simplicity, the thickness of the plate is considered that remains constant.In these conditions, the crystal shape can be specified by two dimensions: length (L 1 ) and width (L 2 ).The crystal population is described by the 2D population density function (L 1 ,L 2 ,t)→n(L 1 ,L 2 ,t) by means of which n(L 1 ,L 2 ,t)dL 1 dL 2 expresses the number of crystals from the size domain (L 1 ,L 1 + dL 1 ) × (L 2 ,L 2 + dL 2 ) in a unit volume of suspension at time t.
As being a seeded crystallization, the secondary nucleation mechanism is considered which depends both on stirring energy and specific surface of existing crystals [3] where σ denotes the absolute supersaturation, σ = (c s -c)/c s = S -1, c is the concentration of solute, c s the saturation concentration and ε is the stirring power in a unit mass of slurry.µ 11 is one of mixed moments of internal coordinates The crystal growth rate is size and supersaturation dependent: Note that the temperature effect on the kinetics is taken into consideration indirectly through the temperature dependency of saturation concentration.
The breakage and agglomeration of crystals are considered negligible but the working volume and the suspension density are considered constant.Under the presented conditions the population balance equation governing the temporal evolution of crystal size distribution takes the form: with the appropriate initial and boundary conditions: Note that the initial condition n 0 (L 1 ,L 2 ) for Eq. ( 4) is the seed distribution.
The mass balance for the solute: where R V denotes the rate of change of the total volume of crystals and ρ c is the crystal density.The R V can be expressed as Note that the initial condition for Eq. ( 8) is the saturation concentration computed at temperature T h .
The linear cooling is described as follows where the cr is expressed in o C/s.Note that the initial condition for the Eq. ( 10) is the T s .
The temperature dependency of the solubility concentration is described by the Apelblat model [9]: where c s denotes the saturation concentration on temperature T expressed in kg/m 3 .The parameters a k were obtained by fitting the model (11) to the experimental data existing in literature [3].

Solution of the model equations
The population balance Eq. ( 4) describes the time evolution of the crystal size distribution which is a state-space equation having two internal coordinates (L 1 and L 2 ), the simultaneous solution with the mass balance Eq. ( 8) and temperature Eq. ( 10) is quite complicated.However, some numerical techniques exist but those generally are computationally expensive.The moment method seems to be an attractive alternative which permits the calculation of the mixed moments of the internal size variables and is defined as: Applying the moment transformation (12) on the population balance (4), the infinite hierarchy of moment equation system can be generated: It is obvious that the moment equation system can be closed and solved only if α d = (0, 1).In practice, the size dependency of growth rate is rarely linear so the moment equation system can be applied only in special situations.
In order to avoid the problem of system closure, the quadrature approximation can be used, numerically evaluating the crystal density function as [10].
Applying the quadrature approximation and the moment transformation on the population balance Eq. ( 4) the quadrature moment equation system is obtained where A 1 (k,m) + A 2 (k,m)= B if k = 0 and m = 0 and 0 else.The set of quadrature moment equations (17) with the mass balance Eq. ( 8) and the temperature Eq. ( 10) forms a closed system which describes the dynamics of the cooling batch crystallization of plate-like crystals.

Simulation results and discussions
Computations were carried out in the MatLab environment solving the resulted differential algebraic equation system using the ode15s function of the program which uses backward differentiation formulas.
The basic values of the kinetic and process parameters are listed in Table 1.
The parameters were chosen based on the paper of Oullion [3] so they show observable effects on process.Each parameter which differs from those presented above is specified.For the sake of a more realistic system investigation, in the sensitivity analysis only the process parameters are varied.
When computing the moments of initial (seed) distribution a log-normal density function with 30% standard deviation of mean values was assumed and the quadrature nodes and weights were computed from moments using the direct optimization technique which are presented in Table 2.The seeds will be characterized with (L 1 , L 2 ) pairs which denote the mean sizes expressed, for the sake of simplicity, in micrometres having the upper mentioned distributions.
Figure 3 presents the dynamic evolution of the mean characteristic particle sizes <L 1 > = µ 10 /µ 00 and <L 2 > = µ 01 /µ 00 with different stirring powers.As it seems, the mean particle sizes are monotonically increasing and the lower stirring power favours for the production of the higher particles as it intensifies the crystal nucleation.The initial values are the sizes of seed crystals.
It is known that the seed characteristics, i.e. mass and size have impact on the product properties.In Fig. 4 we can see the mean crystal sizes obtained using different seeding.The finer seeds with 2 % seed addition results the highest crystals which can be explained with the fact that as the average seed crystal size is decreasing, their number and specific surface increases which favours for the crystal growth more than for the secondary nucleation.
The shape of crystals can be characterized by the length to width ratio L 1 /L 2 , called aspect ratio of an individual crystal, but we can define this ratio also by using the mean length to width <L 1 >/<L 2 >, termed mean aspect ratio.Figure 5 presents the final mean aspect ratio obtained with different seeds.Generally, the increased seed quantity reduces the aspect ratio.
Figure 6 presents the effects of cooling rate and stirring power on the product properties.It seems that the mean crystal sizes are sensitive to both of them but the <L 2 > shows increased sensitivity to the cooling rate compared to <L 1 >.This fact suggests that the aspect ratio should be more sensitive to the cooling rate.
Figure 7 presents the steady-state aspect ratios obtained with different cooling rates and stirring powers.The aspect ratio is generally between 2 and 2.5.The interval is wider than in the case of varying the seed properties and generally produces higher ratios which seem to depend weakly on the stirring power.
9.74 x 10 -5 9.80 x 10 -5 1.138 x 10 -4 6.85 x 10 -5 5.74 x 10 -5 6.05 x 10 -5 Table 1 Basic values of process and kinetic parameters The seeding temperature is also significant parameter of the seeded crystallization.As the seeding temperature decreases, the initial supersaturation is inevitably increased.The lower limit of the seeding temperature is metastable zone limit, where the crystallization should always be conducted.The Fig. 8 illustrates the steady state mean crystal sizes obtained with different seeding temperatures and cooling rates.Both sizes are sensitive to the variations of the seeding temperature which has greater impact than the cooling rate.
In Fig. 9 it is seems clearly that, out of the presented analyses, the variation of the mean aspect ratio is the most sensitive to the seeding temperature.These observations are explained later in this article.
Figure 10 presents the temporal evolution of the nucleation rate if different seeding temperatures are applied.Three order of magnitude differences exist in the initial nucleation rate and in the case of lower seeding temperatures a peak is also observed which is followed by a strong decrease in nucleation.The presence of peak can be explained with the appearance of the intensive crystal growth which increases the surface of crystals so intensifies the nucleation.However, the crystal growth quickly reduces the supersaturation so quickly leads to a dramatic reduction of nucleation rate.
Figure 11 presents a three dimensional projection of the 10 dimensional trajectory on the (µ 11 , S, R V ) subspace.In the case of 35 ºC seeding temperature the initial supersaturation ratio is near to the unity while it is almost 2.5 when 29 ºC seeding temperature is applied which has a strong impact both on nucleation and growth rate.It is seen that the increased growth rate decreases the supersaturation quickly and in all cases slows down only when the supersaturation ratio reaches low levels so the process is governed by the growth, what was expected as being a seeded crystallization.The crystal surface is increasing rapidly during the intensive growth.
The earlier presented anomalies can be explained as follows: based in Fig. 10 the initial supersaturation is relatively high as a consequence the initial growth rate is also high.As the supersaturation dependency of growth rate of characteristic crystal faces differs, they are growing with quite different rates.However, if initially a crystal face becomes higher due to the supersaturation dependency, later is expected to grow faster due to its size dependent character.This chain reaction inevitably leads to the increased length to width ratios as presented in Fig. 9 so to distorted crystal shapes compared to ones obtained in any other conditions.A sensitivity analysis is carried out by numerical experimentation, investigating the effects of process parameters like the seed properties, applied stirring power and cooling rate on the particle size and shape properties.It is found that the seed properties, seeding temperature and the cooling rate has strong impact both on mean size and the length to width ratio, thus the particle shape, but the stirring power have no significant impact on the particle shape; it affect only the mean crystal sizes.It is showed that the seeding temperature has the most significant effect on the crystal shape properties.
Based on the presented case study, it is possible to vary the crystal shape and size simultaneously by choosing the adequate seed and setting up the appropriate process parameters like the cooling rate, seeding temperature and stirring power.However, to optimize a real crystallization process an experimental data based model calibration is always needed.The robustness of the model is sensitive to the quality and quantity of the experimental data so real enhancement on a physical system can be achieved only if having the necessary measurements.In this context the presented model is not able and is not attempted to replace the laboratory experiments.

Fig. 1
Fig. 1 Solute concentration and solubility as a function of temperature in a batch seeded cooling crystallization

Fig. 8 Fig. 9
Fig. 8 Effects of cooling rate and seeding temperature on the product properties a) b)

Fig. 11
Fig. 11 Time evolution of nucleation rate with different seeding temperatures

Table 2
Initial quadrature weights and abscissas The financial support of the Hungarian Scientific Research Fund under Grant K 77955 and the Collegium Talentum is gratefully acknowledged.This research was realized in the frames of TÁMOP 4.2.4.A/2-11-1-2012-0001 "National Excellence Program -Elaborating and operating an inland student and researcher personal support system". - s -1 c concentration of solute, kgm -3 cr cooling rate, ºCs -1 c s saturation concentration, kgm -3 g exponent of crystal growth rate G crystal growth rate, ms -1 k rate coefficient of crystal growth, ms -1 k S rate coefficient of nucleation, # m -3 s -1 k V volume shape factor L linear size of crystals, m L i quadrature abscissas, m <L 1 > µ 10 /µ 00 , m <L 2 > µ 01 /µ 00 , m n population density function, # m -3 approximated population density function, # m -3 R V overall volumetric growth rate m 3 s -1