Modelling of the Polymorph Nucleation Based on Classical Nucleation Theory

To elucidate the relative nucleation rates of different polymorphs, a competitive kinetic model is developed based on classical nucleation theory to describe the time evolution of two different polymorphic cluster size distributions controlled by the association and dissociation of the solute molecules during polymorph nucleation. Although there is only one type of the solute molecules, the agglomerated solute clusters are divided into two types–A form and B form, which resemble the structures and morphologies of the different mature polymorphs and eventually lead to the formation of two polymorphic crystals. A dissociation kernel is incorporated into the proposed model to account for gradual dissolution of the solute clusters smaller than a critical nucleus size due to the thermodynamic instability. By fitting the experimental induction period data and the final measured weight fractions of eflucimibe polymorphs with the proposed model, the association and dissociation rate constants for two polymorphs are determined. The developed model is satisfactory to explain the competitive mechanism of polymorph nucleation for eflucimibe that B form dominates at higher supersaturation while A form dominates at lower supersaturation. The results also indicate that A form is more stable than B form with a transition energy of 3.1 kJ/mole at 35 ◦C.


Introduction
A substance capable of crystallizing into different crystalline forms is said to exhibit polymorphism.Polymorphism plays an important role for pharmaceutical products.Different polymorphs of the same substance might lead to different biological activity due to the possibility of conversion among polymorphic forms affecting the dissolution rate [1,2].
Solute molecules in the supersaturated solution can aggregate and form clusters.Clusters can further aggregate to a bigger size or break into some smaller particles due to the thermodynamic instability.When the size of a cluster exceeds a critical size, it becomes stable and the subsequent growth leads to a new crystal [3][4][5].The existence of solute clusters in the supersaturated solution for the single polymorphic system has been reported by many researchers [6][7][8][9][10].For the multiple polymorphic system, a general hypothesis is accepted in the literature that the solute molecules in the supersaturated solution assemble to form different polymorphic clusters, which resemble the structures and morphologies of the various mature polymorphs and eventually lead to the formation of multiple polymorphic crystals [11,12].Recently, Van Driessche et al. [13] observed polymorph nucleation events that are driven by oriented attachments between subcritical clusters of the protein glucose isomerase that already exhibit a degree of crystallinity.Consequently, the relative nucleation rates of different polymorphs should be closely related to the time evolution of multiple polymorphic clusters in nucleation.
Understanding the nucleation and growth mechanisms of polymorphism is crucial to better control of the desired forms during polymorph crystallization [14][15][16][17][18][19][20][21][22][23][24][25][26][27].Although various nucleation theories for the single polymorphic system have been developed based on either the thermodynamics of the process [4,5] or the kinetics of the process [28][29][30], nucleation theory for the multiple polymorphic system is comparatively less studied in the literature.Hammond et al. [31] used a molecular modeling approach to study the stability of different polymorphic forms of L-glutamic acid through building and optimizing molecular clusters of different sizes and shapes.ter Horst et al. [32] adopted a combination method of molecular simulations and process modeling to predict the polymorphic fraction and crystal size distribution during polymorph crystallization.Deij et al. [33] applied the growth probability method combined with a Monte Carlo routine to simulate polymorph formation for some dimorphic systems.
When multiple polymorphs of the same substance crystallize together out of a solution, the relative productivity of different polymorphs depends on supersaturation, temperature, cooling rate, solvent, agitation, PH, additive, impurity, seeding, etc.Among these factors, supersaturation is of utmost significance [34].For example, Ni and Liao [35] and Qu et al. [36] indicated for L-glutamic acid that the metastable α crystals are favored at lower supersaturations while the stable β crystals are favored at higher supersaturations.Sun et al. [37] reported that either spontaneously nucleating quiescent aqueous L-glycine solution or nonphotochemical laser-induced nucleation tends to produce the intermediate stable α glycine at lower supersaturations and the most stable γ glycine at higher supersaturations.Gracin and Rasmuson [38] observed for the polymorph nucleation of p-aminobenzoic acid that cooling crystallization can be performed to produce the most stable pure β form in water and in ethyl acetate by careful control of supersaturation and temperature.Teychene and Biscans [39] investigated the nucleation kinetics of two eflucimibe polymorphs by induction time measurements and found that the stable A-form crystals are favored at lower supersaturations while the metastable stable B-form crystals are favored at higher supersaturations.Zhu et al. [40] indicated that the polymorph nucleation of gestodene in ethanol depends on both the supersaturation and crystallization temperature.Roelands et al. [41] and Wantha et al. [42] reported that the metastable polymorph increases with increasing supersaturation for the crystallization of L-histidine from aqueous solution with the antisolvent ethanol.
Stranki and Totomanov [43] attempted an explanation of the Ostwald rule of stages from the corresponding rates of crystal nucleation and argued that the first nucleated phase is the phase that has the lowest free-energy barrier of formation, i.e., the one which has the fastest nucleation rate, rather than the most stable phase; afterwards, the system may undergo a polymorphic form transition toward another metastable phase, or directly to the stable phase.Tahri et al. [44] evaluated the competition between the nucleation, the growth, and the Ostwald ripening of the different phases by means of the kinetic equation model and concluded that Ostwald ripening can induce the total dissolution of the slow growing stable polymorph nuclei, leading to a result in agreement with the Ostwald rule of stages.
As supersaturation plays an important role during polymorph crystallization, understanding the competition between nucleation rates of different polymorphs is crucial for polymorphism control in a supersaturated solution.In the current paper a competitive kinetic model for polymorph nucleation is proposed to elucidate the effects of supersaturation on polymorph selection.

Theory
A competitive kinetic model is developed to describe the time evolution of two different polymorphic cluster size distributions in a supersaturated solution, where the simultaneous nucleation of the two polymorphs is controlled by the association and dissociation mechanisms during the induction time period.Although the Ostwald rule of stages postulates that a crystallization system progresses from the supersaturated state to equilibrium in stages, each stage representing the smallest possible change in free energy [4].In the dimorphic system this means the initial appearance of the metastable crystals, followed by their transformation to the stable crystals.However, in the actual crystallization system kinetics are often more important than thermodynamics.It is in fact, if more than one phase is thermodynamically possible, the resulting phase is not just the one that is thermodynamically most likely.Instead, the resulting crystals are determined by the relative rates of crystal nucleation and growth of the metastable and stable forms, which usually depends on the process conditions, e.g., supersaturation, crystallization temperature, and the solvent type [4,[45][46][47].
In derivation of the competitive kinetic model for polymorph nucleation, it is assumed that: (1) Although there is only one type of the solute molecules (n 1 ), the solute clusters are divided into two types-A form (n i,A , i ≥ 2) and B form (n i,B , i ≥ 2); (2) Depending on the steric colliding direction between two solute molecules, two solute molecules can collide to form either an A-type or B-type cluster; (3) Depending on the steric colliding direction, a solute molecule can collide with an A-type cluster to form a larger A-type cluster or collide with a B-type cluster to form a larger B-type cluster; (4) Only two solute clusters of the same type can collide to form a larger solute cluster of the same type, i.e., the collision between A-type clusters and B-type clusters results in no formation of a larger solute cluster due to the structure difference; (5) A solute cluster smaller than the critical nucleus size for n i,A (2 ≤ i ≤ g A − 1) or n i,B (2 ≤ i ≤ g B − 1) might dissociate to form one primary particle and one smaller solute cluster of the same type due to the thermodynamic instability; (6) As the size of a cluster exceeds the critical nucleus size for n i,A (i ≥ g A ) or n i,B (i ≥ g B ), it becomes stable and no dissociation occurs.The subsequent growth leads to a new crystal of the same type.
Various association and dissociation processes are depicted in simplified form in Table 1.The time evolution of n 1 , n i,A (i ≥ 2) and n i,B (i ≥ 2) for polymorph nucleation in a super saturation solution can be derived as follows.Note that n 1,A = n 1,B = n 1 due to only one type of the solute molecules.Based on Smoluchowski's agglomeration theory [48,49], the net formation rate of n 1 , n i,A (i ≥ 2), and n i,B (i ≥ 2), by association can be described respectively as: RA where k A and k B represents the association rate constant of A-type clusters and B-type clusters, respectively.Note that n t,A = ∑ ∞ i=2 n i,A and n t,B = ∑ ∞ i=2 n i,B .
Table 1.Various association and dissociation processes for polymorph nucleation.

Association and Dissociation
Processes According to classical nucleation theory (CNT) [4,5], a cluster becomes stable when its size reaches the critical nucleus size.Thus, dissociation only occurs for n i,A (2 ≤ i ≤ g A − 1) or n i,B (2 ≤ i ≤ g B − 1).As similar to the solute clustering process in a supersaturation solution [50], gradual dissociation kernel is proposed in this work.For gradual dissociation kernel, one large aggregate, n i,A (2 ≤ i ≤ g A − 1), dissociates into one primary particle, n 1 , and one smaller aggregate, n i−1,A .Similar dissociation process is applied for n i,B (2 ≤ i ≤ g B − 1).Thus, gradual dissociation kernel corresponds to gradual dissolution of the solute clusters smaller than a critical nucleus size due to the thermodynamic instability.
The net formation rate of n 1 by dissociation for gradual disruption kernel is given by where the first and second terms on the right-hand side of Equation ( 4) represent the birth term due to dissociation of n 2,A and n 2,B while the third and fourth terms represent the birth term due to dissociation of n j,A (3 ≤ j ≤ g A − 1) and n j,B (3 ≤ j ≤ g B − 1).In this work, k dA,j = k dA0 (j − 1) and k dB,j = k dB0 (j − 1) are assumed due to more molecules available for dissociation for larger clusters.As dissociation occurs for n i,A (2 ≤ i ≤ g A − 1), the net formation rate of n i,A (i ≥ 2) by dissociation can be described as Note that no dissociation occurs for n i,A (i ≥ g A ). Similarly, as dissociation occurs for n i,B (2 ≤ i ≤ g B − 1), RD i,B (2 ≤ i ≤ g B − 1) can be derived.Note that no dissociation occurs for n i,B (i ≥ g B ).
To determine the time evolution of n 1 , n i,A (i ≥ 2), and n i,A (i ≥ 2), the net formation rate by both association and dissociation should be considered.Thus, one obtains Substituting Equations (1,4) into Equation ( 8) yields dn 1 Substituting Equations (2,5-7) into Equation ( 9) yields Summing Equations (12)(13)(14) where ∑ ∞ i=2 ∑ i−1 j=1 n j,A n i−j,A = (n 1 + n t,A ) 2 (see Equation (A1) in Appendix A).Multiplying Equations ( 14) by i and summing the resulting equation from i where Summing Equations (17)(18)(19) Multiplying Equations ( 19) by i and summing the resulting where When solutes exist as free molecules initially, n 1 (0) = C 0 , n i,A (0) = n i,B (0) = 0 (i ≥ 2), n t,A (0) = n t,B (0) = 0, and M C,A (0) = M C,B (0) = 0. Thus, Equations (11-13, 15-18, 20, 21) constitute of a set of g A + g B + 1 differential equations that can be simultaneously solved for the time evolution of The induction time is defined as the time at which the volume fraction of the newly formed solid phase in solutions becomes detectable [4].For simplicity, the minimum detectable volume fraction of the newly formed solid phase, f V , at the induction time is assumed corresponding to the total volume fraction of stable A-type and B-type clusters.As n i,A (i ≥ g A ) and n i,B (i ≥ g B ) represents the number concentration of stable A-type and B-type clusters, respectively, one obtains At where V m denotes the volume of a solute molecule.As V m M C,A and V m M C,B represents the volume fraction of stable A-type and B-type clusters in solutions, respectively, the weight fraction of the stable A-type crystals at the induction time is defined as Note that f V depends on the measurement device and on the substance.Based on the study of 28 inorganic systems, Mersmann and Bartosch [51] estimated f V = 10 −4 − 10 −3 .As the intermediate value, f V = 4 × 10 −4 , was adopted at the detection of the nucleation point for the Lasentec focus beam reflectance measurements reported by Lindenberg and Mazzotti [52] and for the turbidity measurements reported by Shiau and his coworkers [53][54][55], this value is also adopted in this study.
The proposed model is applied to the polymorph nucleation based on CNT.However, for the two-step mechanism of nucleation, nucleation is thought to consist of two steps in series, i.e., the formation of the dense liquid clusters in solutions followed by the formation of the crystalline nucleus inside the dense liquid clusters [56].Thus, the derivation above is not applicable to the two-step mechanism of nucleation.

Results and Discussion
The experimental results for eflucimibe polymorph nucleation reported by Teychene and Biscans [39] are illustrated to verify the developed model.Eflucimibe is a new drug inhibiting acyl-coenzyme A: cholesterol acyltransferase (ACAT), an enzyme which inhibition may lead to lower serum cholesterol concentration.Eflucimibe crystallizes from a mixture of ethanol and n-heptane into two polymorphic forms-A form and B form.A form is the stable form while B form is the metastable form.These two forms have different solubilities and interfacial energies.The experimental induction time data and the final measured weight fractions of eflucimibe polymorphs for various supersaturation at 35 • C are listed in Table 2. Based on CNT, critical nucleus sizes for A-type and B-type clusters can be respectively estimated as [4,5,57,58] In Table 2, g A and g B are the closest integers calculated for each S A based on γ A = 5.17erg/cm 2 and γ B = 4.23erg/cm 2 , respectively.Note that both g A and g B decrease with increasing S A .
The following values at 35 • C are used in the calculations: ρ c = 1.3g/cm 3 , C eq,A = 2.12 × By fitting the experimental induction time data and measured weight fractions of eflucimibe polymorphs with the proposed model, the following regression procedure is adopted to determine four parameters-k A , k B , k dA0 and k dB0 : (1) guess a set of k A , k B , k dA0 and k dB0 ; (2) determine the time evolution of n 1 , n i,A (2 ≤ i ≤ g A − 1), n t,A , n i,B (2 ≤ i ≤ g B − 1), n t,B , M C,A and M C,B for each S A by solving Equations (11-13, 15-18, 20, 21) simultaneously; (3) calculate t ind,the , at which it leads to V m (M C,A + M C,B ) = f V , based on Equations ( 22) for each S A ; (4) calculate W A,the at t ind,exp based on Equation ( 23) for each S A ; (5) calculate ε t from Equation ( 26); ( 6) calculate ε w from Equation (27).
The average relative deviation between the experimental and theoretical induction time is defined as where H is the number of the experimental runs, t ind,exp is the induction time measured during the experiments, and t ind,the is the induction time calculated from Equation (22).As shown in Table 2, pure A form is crystallized at low supersaturation and pure B form is crystallized at high supersaturation while a mixture of A form and B form is crystallized at intermediate supersaturation.For simplicity, the weight fractions of A-form eflucimibe measured at the end of the experiments are assumed close to the weight fractions of A-form eflucimibe nucleated at the induction time.The average relative deviation between the experimental and theoretical weight fraction of A-type crystals at the induction time is defined as where w A,exp is the weight fractions of A-form eflucimibe measured at the end of the experiments, and w A,the is the weight fractions of A-form eflucimibe calculated from Equation (23).
It should be noted in the regression procedure that, for any guessed values of k A , k B , k dA0 and k dB0 , a set of 184 differential equations for S A = 1.86 (g A = 91, g B = 92) need to be simultaneously solved for the time evolution of n 1 , n i,A (2 ≤ i ≤ 90), n t,A , n i,B (2 ≤ i ≤ 91), n t,B , M C,A and M C,B .However, as supersaturation is increased to S A = 3.4 (g A = 12, g B = 9), only a set of 22 differential equations need to be simultaneously solved for the time evolution of n 1 , n i,A (2 , n t,B , M C,A and M C,B . By repeating the regression procedure from ( 1) to ( 6), the optimal values of k A , k B , k dA0 and k dB0 with the smallest sum of ε t and ε w are determined based on f V = 4 × 10 −4 in Table 3, which indicates k A = 7.7 × 10 −22 cm 3 /s, k B = 1.4 × 10 −21 cm 3 /s, k dA0 = 5.3 × 10 −4 s −1 and k dB0 = 3.2 × 10 −3 s −1 , leading to k B k A = 2 and k dB0 k dA0 = 6.Thus, one obtains k B > k A and k dB0 > k dA0 .
Table 3.The optimal values of k A , k B , k dA0 and k dB0 with the smallest sum of ε t and ε w for Due to the different activation energy for bonding between solute molecules, the association rate constants for A-type and B-type clusters can be expressed respectively as where exp − represents the fraction of the successful collisions which overcomes the activation energy of association for A-type and B-type clusters, respectively.The collision rate constant, k 0 , is proportional to the diffusivity for solute clusters in a supersaturated solution [4].Combining Equations (28)(29) Substituting the optimal values of k A and k B at 35 • C into Equation ( 30) yields E * A − E * B = 1.5kJ/mole.Thus, we conclude E * A > E * B for association.The dissociation rate constants for A-type and B-type clusters can be expressed respectively as where the dissociation frequency factor, k d0 , mainly depends on the solution property.Combining Equations (32)(33) yields As E * A − E * B = 1.5kJ/mole, substituting the optimal values of k dA0 and k dB0 at 35 • C into Equation ( 33) yields E B − E A = 3.1kJ/mole.Thus, it yields E B > E A , which is consistent with the literature that form A is more stable than form B [39].Based on the results using Equations (30,33), the change in potential energy during polymorph crystallization is plotted in Figure 1.Based on the optimal values of k A , k B , k dA0 and k dB0 in Table 3, Figure 2 shows comparison of t ind,the and t ind,exp for various S A , where t ind,the is calculated from Equation ( 22) using the time evolution M C,A and M C,B .The dashed line represents the calculated t ind,the while the solid circle represents the experimental t ind,exp .It is found that t ind,the and t ind,exp decreases with increasing S A as both g A and g B decrease with increasing S A .Comparison of t ind,the and t ind,exp for various S A , where t ind,the is calculated from Equation ( 22) using the optimal values of k A , k B , k dA0 and k dB0 .The dashed line represents the calculated t ind,the while the solid circle represents the experimental t ind,exp .
Figure 3 shows comparison of W A,exp and W A,the at the end of the experiments for various S A , where W A,the is calculated from Equation ( 23) using the time evolution M C,A and M C,B .The dashed line represents the calculated W A,the while the solid circle represents the experimental W A,exp .It is found that W A,exp and W A,the decreases with increasing S A .These results can be attributed to k B > k A and k dB0 > k dA0 .At higher S A due to smaller critical nucleus size and shorter induction time, association rate among n 1 , n i,A (2 ≤ i ≤ g A − 1), and n i,B (2 ≤ i ≤ g B − 1) plays a more important role than dissociation rate for n i,A (2 ≤ i ≤ g A − 1) and n i,B (2 ≤ i ≤ g B − 1) in nucleation.As k B > k A , it is easier for the metastable B-form to grow to smaller critical nucleus size at higher S A .Thus, nucleation of the metastable B-form dominates at higher S A , leading to M C,B > M C,A and w B > w A at the end of the experiments.However, at lower S A due to greater critical nucleus size and longer induction time, dissociation rate becomes important for n i,A (2 ≤ i ≤ g A − 1) and n i,B (2 ≤ i ≤ g B − 1) in nucleation.As k dB0 > k dA0 , it becomes more difficult for the metastable B-form to grow to larger critical nucleus size at lower S A .Subsequently, nucleation of the stable A-form dominates at lower S A , leading to M C,A > M C,B and w A > w B at the end of the experiments.As displayed in Figures 2 and 3, the fourteen experimental data points are fitted well to the corresponding calculated results by the developed model using the optimal values of k A , k B , k dA0 and k dB0 .The following time evolution of n 1 , n i,A (2 ≤ i ≤ g A − 1), n t,A , n i,B (2 ≤ i ≤ g B − 1), n t,B , M C,A and M C,B are calculated based on the optimal values of k A , k B , k dA0 and k dB0 in Table 3.For example, a set of 184 differential equations for S A = 1.86 (g A = 91, g B = 92) are simultaneously solved for the time evolution of n 1 , n i,A (2 ≤ i ≤ 90), n t,A , n i,B (2 ≤ i ≤ 91), n t,B , M C,A and M C,B .Similarly, a set of 22 differential equations for S A = 3.4 (g A = 12, g B = 9) are simultaneously solved for the time evolution of n 1 , n i,A (2 ≤ i ≤ 11), n t,A , n i,B (2 ≤ i ≤ 8), n t,B , M C,A and M C,B .
The variations of n 1 , n 2,A , n 2,B , n 3,A , n 3,B , n t,A and n t,B with increasing time are displayed for various S A in Figure 4, where n 1 decreases monotonically with increasing time; however, n 2,A , n 2,B , n 3,A , n 3,B , n t,A and n t,B passes through a maximum at certain time and then declines slowly for each S A .A general trend is observed for various S A that, although n 1 is significantly greater than n 2,A , n 2,B , n 3,A and n 3,B for 0 < t < 0.1 t ind , n 2,A , n 2,B , n 3,A and n 3,B become not negligible compared to n 1 for t > 0.1 t ind .Thus, as n 1 dominates in the earlier stage of nucleation, association between two solute clusters of the same type, n j,A (j ≥ 2) or n j,B (j ≥ 2), is nearly negligible for 0 < t < 0.1 t ind .However, such association between two solute clusters of the same type becomes significant in the later stage of nucleation for t > 0.1 t ind .In Figure 4, n 1 /C 0 is decreased to 0.04 at t ind = 80 s for S A = 3.4; n 1 /C 0 is decreased to 0.02 at t ind = 1614s for S A = 2.5; n 1 /C 0 is decreased to 0.01 at t ind = 28000 s for S A = 1.86.Thus, n 1 /C 0 at t ind becomes smaller for a lower S A due to longer induction time available for association.The variations of M C,A and M C,B with increasing time are displayed for various S A in Figure 5, where M C,A and M C,B increases monotonically with increasing time for each S A .For S A = 1.86 − 1.97, M C,A is significantly greater than M C,B at t ind , leading to w A,the = 1 at t ind based on Equation (23), which is consistent with w A,exp = 0.96 − 0.97 in Table 2.Note that M C,B /C 0 remains nearly zero in the range 0 − t ind .Thus, nearly only form A is obtained at the end of experiments.For S A ≥ 2.3, M C,B becomes significant compared to M C,A at t ind .For example, M C,A is only slightly than M C,B at t ind for S A = 2.5, leading to w A,the = 0.55 at t ind based on Equation (23), which is close to w A,exp = 0.64 in Table 2. Thus, form A is slightly more than form B at the end of experiments.For S A = 2.7, M C,B becomes greater than M C,A at t ind , leading to w A,the = 0.33 at t ind based on Equation (23), which is consistent with w A,exp = 0.33 in Table 2. Thus, form B is more than form A at the end of experiments.For S A = 3.4, M C,B is significantly greater than M C,A at t ind , leading to w A,the = 0 at t ind based on Equation (23), which is consistent with w A,exp = 0 in Table 2.Note that M C,A /C 0 remains nearly zero in the range 0 − t ind .Thus, nearly only form B is obtained at the end of experiments.A − E * B = 1.5kJ/mole and E B − E A = 3.1kJ/mole.This is consistent with Ostwald's rule of stages that the metastable B-form tends to crystallize out more easily than the stable A-form due to E * A > E * B while the stable A-form has a lower potential energy than the metastable B-form due to E B > E A .Supersaturation is crucial in the polymorph selection for nucleation of eflucimibe.Association rate plays a more important role than dissociation rate at higher supersaturation due to smaller critical nucleus size and shorter induction period.As k B > k A , it is easier for the metastable B-form to grow to smaller critical nucleus size at higher supersaturation.Thus, nucleation of the metastable B-form dominates at higher supersaturation, leading to w B > w A at the end of the experiments.However, dissociation rate becomes important at lower supersaturation due to larger critical nucleus size and longer induction period.As k dB0 > k dA0 , it becomes more difficult for the metastable B-form to grow to larger critical nucleus size at lower supersaturation.Subsequently, nucleation of the stable A-form dominates at lower supersaturation, leading to w A > w B at the end of the experiments.These findings are consistent with the experimental results obtained by Teychene and Biscans [39], indicating that the metastable B-form dominates at higher supersaturation while the stable A-form dominates at lower supersaturation.

Conflicts of Interest:
The authors declare no conflict of interest.

Notation
C 0 = initial concentration of solutes (molecules/cm 3 ) C eq = saturated concentration of solutes (molecules/cm 3 ) E 0 = potential energy of the solute molecules (kJ/mole) E A = potential energy of A-type crystals (kJ/mole) E B = potential energy of B-type crystals (kJ/mole) E * A = activation energy for the association between A-type solute clusters (kJ/mole) E * B = activation energy for the association between B-type solute clusters (kJ/mole) g A = critical nucleus size of A-type crystals (dimensionless) g B = critical nucleus size of B-type crystals (dimensionless) k A = association rate constant of A-type clusters (cm 3 /s) k B = association rate constant of B-type clusters (cm 3 /s) k 0 = collision rate constant (cm 3 /s) k d0 = dissociation frequency factor (cm 3 /s) k dA,j = dissociation rate constant of n j,A (1/s) k dA0 = dissociation rate coefficient of A-type clusters (1/s) k dB,j = dissociation rate constant of n j,B (1/s) k dB0 = dissociation rate coefficient of B-type clusters (1/s) M C,A = total number of molecules for n i,A (i ≥ g A ) per unit volume, (#/cm 3 ) M C,B = total number of molecules for n i,B (i ≥ g B ) per unit volume, (#/cm 3 ) M W = molecular weight (g/mol) N A = Avogadro number (1/mol) n i,A = number concentration of A-type clusters with i solute molecules (#/cm 3 ) n i,B = number concentration of B-type clusters with i solute molecules (#/cm 3 ) n t,A = total number concentration of A-type clusters (#/cm 3 ) n t,B = total number concentration of B-type clusters (#/cm 3 ) RA i,A = net formation rate of n i,A due to association (#/cm 3 -s) RD i,A = net formation rate of n i,A due to dissociation (#/cm 3 -s) RA i,B = net formation rate of n i,B due to association (#/cm 3 -s) RD i,B = net formation rate of n i,B due to dissociation (#/cm 3 -s) S A = supersaturation based on A-type crystals, S A = C 0 /C eq,A (dimensionless) S B = supersaturation based on B-type crystals, S B = C 0 /C eq,B (dimensionless) T = temperature (K) t = time (s) t ind,exp = experimental induction time (s) t ind,the = theoretical induction time (s) V m = volume of a solute molecule (cm 3 ) w A,exp = experimental weight fraction of A-type crystals (dimensionless) w A,the = theoretical weight fraction of A-type crystals (dimensionless) Greek letters γ A = interfacial energy for form A (erg/cm 2 ) γ B = interfacial energy for form B (erg/cm 2 ) ρ C = crystal density (g/cm 3 ) ε t = average relative deviation between the experimental and theoretical induction time (dimensionless) ε w = average relative deviation between the experimental and theoretical weight fraction of A-type or B-type crystals (dimensionless).

Figure 1 .
Figure 1.Schematic diagram for the change in potential energy during polymorph crystallization from the initial state to two different polymorphs.

Figure 2 .
Figure 2. Comparison of t ind,the and t ind,exp for various S A , where t ind,the is calculated from Equation (22) using the optimal values of k A , k B , k dA0 and k dB0 .The dashed line represents the calculated t ind,the while the solid circle represents the experimental t ind,exp .

Figure 3 .
Figure 3.Comparison of w A,exp and w A,the at the end of the experiment for various S A , where w A,the is calculated from Equation (23) using the optimal values of k A , k B , k dA0 and k dB0 .The dashed line represents the calculated w A,the while the solid circle represents the experimental w A,exp .

Figure 4 .
Figure 4. Variation of n 1 , n 2,A , n 2,B , n 3,A , n 3,B , n t,A and n t,B with increasing time for various S A in the range t = 0 − t ind .

Figure 5 .
Figure 5. Variation of M C,A and M C,B with increasing time for various S A in the range t = 0 − t ind .4. Conclusions A competitive kinetic model for polymorph nucleation is developed in this work to describe the time evolution of two different polymorphic cluster size distributions in a supersaturated solution.By fitting the experimental induction time data and measured weight fractions of eflucimibe polymorphs with the proposed model, the association and dissociation rate constants for two polymorphs are determined, leading to k A = 7.7 × 10 −22 cm 3 /s, k B = 1.4 × 10 −21 cm 3 /s, k dA0 = 5.3 × 10 −4 s −1 and k dB0 = 3.2 × 10 −3 s −1 .Thus, one obtains k B k A = 2 and k dB0 k dA0 = 6, leading to E *A − E * B = 1.5kJ/mole and E B − E A = 3.1kJ/mole.This is consistent with Ostwald's rule of stages that the metastable B-form tends to crystallize out more easily than the stable A-form due to E * A > E *

Funding:
This research was funded by Chang Gung Memorial Hospital (CMRPD2G0242) and Ministry of Science and Technology of Taiwan (MOST107-2221-E-182-031).Acknowledgments: The author would like to thank Chang Gung Memorial Hospital (CMRPD2G0242) and Ministry of Science and Technology of Taiwan (MOST107-2221-E-182-031) for financial support of this research.

Table 2 .
35e experimental induction time data and the final measured weight fractions of eflucimibe polymorphs for various supersaturation at35 [39]reported by Teychene and Biscans[39].S