**Modeling Ice Crystal Formation of Water in Biological System**

Zhi Zhu He and Jing Liu

cubation method: a correlated study on the hypertrophic zone of the rat epiphyseal

[228] Moradian-Oldak J. Protein-mediated enamel mineralization. *Front Biosci* 2012; 17,

[229] Prince CW, Rahemtulla F, Butler WT. Metabolism of rat bone proteoglycans *in vivo*.

[230] Amizuka N, Hasegawa T, Oda K, Luiz de Freitas PH, Hoshi K, Li M, Ozawa H. His‐ tology of epiphyseal cartilage calcification and endochondral ossification. *Front Biosci*

[231] Marin F, Le Roy N, Marie B. The formation and mineralization of mollusk shell. *Front*

[232] Veis A. Organic matrix-related mineralization of sea urchin spicules, spines, test and

[233] Nagasawa H. The crustacean cuticle: structure, composition and mineralization.

[234] Hincke MT, Nys Y, Gautron J, Mann K, Rodriguez-Navarro AB, McKee MD. The eggshell: structure, composition and mineralization. *Front Biosci* 2012; 17, 1266-1280

[235] Smales FC. Structural subunit in prisms of immature rat enamel. *Nature* 1975; 258,

[236] Bonucci E. Calcification and silicification: a comparative survey of the early stages of

[237] Wang X, Müller WE. Complex structures - smart solutions: Formation of siliceous

[238] Wang X, Schlossmacher U, Wiens M, Batel R, Schröder HC, Müller WE. Silicateins, silicatein interactors, and cellular interplay in sponge skeletogenesis: Formation of glass fiber-like spicules. *FEBS J* 2012; 279, 1721-1736 (doi: 10.1111/j.

[239] Andre R, Tahir MN, Natalio F, Tremel W. Bioinspired synthesis of multifunctional inorganic and bio-organic hybrid materials. *FEBS J* 2012; 279, 1737-1749 (doi:

biomineralization. *J Bone Miner Metab* 2009; 27, 255-264

spicules. *Commun Integr Biol* 2011; 4, 684-688

1742-4658.2012.08533.x)

10.1111/j.1742-4658.2012.08584.x)

cartilage. *Calcif Tissue Int* 1996; 58, 273-282

1996-2023

184 Advanced Topics on Crystal Growth

772-774

*Biochem J* 1983; 216, 589-596

*(Elite Ed.)* 2012; 4, 2085-2100

*Biosci* (Schol Ed) 2012; 4, 1099-1125

teeth. *Front Biosci* 2011; 17, 2540-2560

*Front Biosci (Elite Ed)* 2012; 4, 711-720

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/54098

### **1. Introduction**

Ice crystal formation of water plays a very important role in broadest scientific and engi‐ neering fields such as cloud evolution in atmosphere and brine pockets in marine ice. One of the most attractive attentions of this issue is about the mechanisms of cell injury induced by freezing in cryobiology [1]. During the past few decades, a series of cryogen‐ ic techniques have been developed for clinical applications in cryosurgery and cryopre‐ servation. The former is aimed to destroy the target diseased tissues (such as tumor) while minimize pain and bleeding of the tissue as much as possible. It has been demon‐ strated to be a safe, minimally invasive and cost-effective treatment approach. In contrast to the cryosurgery, the task of cryopreservation is however to preserve cell, tissues or or‐ gans at a super-low temperature condition with the addition of protective agents, and ensure that they can be normally used after several re-warming stages. Clearly, no mat‐ ter which process it involves, the formation of the ice crystal is ineluctable and decisive to the final fate of those cells and tissues subject to treatment.

Cryoinjury induced by ice formation is usually considered as a main killer of cells during cryosurgery or cryopreservation. A two-factor hypothesis [2], "intracellular ice formation" (IIF) at rapid cooling rates and "solution effects" at slow cooling rates, proposed by Mazur, has been serving well to explain the mechanism of freeze injury. The biological materials would subject to damage due to extracellular and intracellular ice formation during freez‐ ing. For the cryopreservation at super-cooling temperature, the occurrence of IIF is the sin‐ gle most important factor determining whether or not cells will survive cryopreservation [3]. So far, the ice crystal formation problems have been partially addressed in material and at‐ mospheric sciences. However, similar researches are hard to find in cryobiology area due to complexity of the water structure in a biological system. This chapter is dedicated to sum‐

© 2013 He and Liu; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. © 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

marize the typical physical and mathematical model for ice nucleation and growth in cryobi‐ ology area. In addition, the heat and solute transport, and interaction between the cell and ice during freezing are also discussed.

ume of phase *β*, and *σ αβ*

spectively modified as

ticle could read as

( )

where *g* =(1 + *x* <sup>2</sup> −2*ηx*)1/2

*c*

h

barrier*ΔGi*

2

could be given by

h

is the interfacial free energy per unit area. According to thermody‐

( ) ( )

s

 p ab 3

2

(4)

 h

 h

<sup>+</sup> (2)

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

(3)

187

(5)

(6)

and nucleation

*c* is re‐

namics, Gibbs free-energy difference between *α* and *β* could be given as

*t*

D =

*<sup>f</sup>* 2

*H T <sup>T</sup> <sup>G</sup>*

water molecules, the following relationships for the critical cluster could be derived:

<sup>2</sup> <sup>16</sup> and = 3 *<sup>c</sup> c i*

 =- D D D

*r G*

ab

s

D D

*f f*

where *ΔHf* is the heat of fusion, *Tf* is the equilibrium freezing temperature, and *ΔT* =*Tf* −*T* . The free energy of form of a cluster of ice is initially positive and goes through a maximum as the cluster size is increased. Maximizing with respect to *i* and assuming closely packed

*<sup>t</sup> <sup>t</sup>*

*G G*

where *rc* is the radius of the critical cluster. For heterogeneous nucleation theory *ΔGi*

( ) , = , *c c G f xG i h <sup>i</sup>* D D h

3 3 1 1 3 2 , 1 2 3 3 1

> ' exp exp *<sup>c</sup> <sup>i</sup> <sup>i</sup> nkT <sup>G</sup> <sup>G</sup> <sup>J</sup> h kT kT* æ ö <sup>D</sup> æ ö <sup>D</sup> =- -ç ÷ ç ÷ ç ÷ è ø è ø

ì ü é ù ï ï æ ö æ öæ ö æ ö - -- - =+ + - + + - ê ú í ý ç ÷ ç ÷ç ÷ ç ÷ ï ï è ø è øè ø è ø î þ ë û

*g gg g*

hh

*x xx x f x <sup>x</sup> <sup>x</sup>*

The homogeneous nucleation rate *J* considering both diffusion barrier*ΔGi*'

where *x* is defined as *x*=*Rs* /*rp*, *Rs* is the radius of the foreign particle, and *rp* is the radius of curvature of the critical nuclei. *η* is approximately regarded as cos*θ* where *θ* is the contact angle between the nucleation phase and the substrate. The factor f for convex surface of par‐

*T TT*

### **2. Ice nucleation**

Freezing often denotes a process of new phase defined by ice creation via nucleation and subsequent crystal growth from supercooled water. This process undergoes three stages including the chemical diffusion stage for water molecule transportation, the surface kinetic stage for hydrogen bonds formation, and the final heat conduction stage for heat releasing and diffusing. It implies that the ice nucleation—whose forma‐ tion determines the growth of ice crystal to a large extent—turns out to be the initial and crucial step during crystallization. In terms of whether there exist foreign bodies or not, the nucleation is usually classified as homogeneous or heterogeneous cases. By definition, homogeneous ice nucleation requests a severe condition that the liquid wa‐ ter is absolutely dust-free, or there exist foreign particles indeed, but whose sizes are small enough to be negligible. Since this rigid demand is hard to be satisfied, the nu‐ cleation is in fact regarded as heterogeneous in most cases and may be influenced by many external factors, such as the size of foreign particles, surface roughness and the effective size of the nucleators at different supercoolings. There exist a great many bar‐ riers to build up an integrated theory on nucleation, because of the dissimilar proper‐ ties might exhibit between bulk fluid and the micro/nanoscale molecules. Up to now, the two major theoretical approaches to characterize the formation of ice can be clas‐ sified as thermodynamic approximation and molecular simulations [4].

#### **2.1. Thermodynamic model of ice nucleation**

The ice nucleation does not happen necessarily if only temperature is below the freezing point. This process is determined by the free energy of system, which implies that the ice germ growth or disappearance depends on the free energy decrease or increase. For constant rate cooling protocols, the temperature of system is given by *T* (*t*)=*T*<sup>0</sup> − *Bt*, where *B* is the cooling rate. Considering a cluster of ice denoted by *β* with *i* molecules from its mother water phase *α*, the driving force of nucleation is derived from the Gibbs free energy difference*ΔGi* [5]

$$
\Delta G\_{\rm i} = \dot{\nu}^{\rho} \Delta G\_{\rm t} + \text{(36\pi)}^{1/3} \dot{\nu}^{2/3} \sigma^{a\beta} \{\nu^{\beta}\}^{2/3} \tag{1}
$$

where *ΔGt* is negative and denotes the Gibbs free-energy difference between *α* and *β* phase per volume, according to classical homogeneous nucleation theory; *ν <sup>β</sup>* is the molecular vol‐ ume of phase *β*, and *σ αβ* is the interfacial free energy per unit area. According to thermody‐ namics, Gibbs free-energy difference between *α* and *β* could be given as

$$
\Delta \mathbf{G}\_t = \frac{\Delta H\_f \Delta T}{T\_f} \frac{2T}{T\_f + T} \tag{2}
$$

where *ΔHf* is the heat of fusion, *Tf* is the equilibrium freezing temperature, and *ΔT* =*Tf* −*T* . The free energy of form of a cluster of ice is initially positive and goes through a maximum as the cluster size is increased. Maximizing with respect to *i* and assuming closely packed water molecules, the following relationships for the critical cluster could be derived:

$$r\_c = -\frac{2\sigma^{\alpha\beta}}{\Delta G\_t} \text{ and } \Delta G\_{i\_c} = \frac{16\pi}{3} \frac{\left(\sigma^{\alpha\beta}\right)^3}{\left(\Delta G\_t\right)^2} \tag{3}$$

where *rc* is the radius of the critical cluster. For heterogeneous nucleation theory *ΔGi c* is re‐ spectively modified as

$$\Delta G\_{i\_c, h} \triangleq f\left(\eta, \mathbf{x}\right) \Delta G\_{i\_c} \tag{4}$$

where *x* is defined as *x*=*Rs* /*rp*, *Rs* is the radius of the foreign particle, and *rp* is the radius of curvature of the critical nuclei. *η* is approximately regarded as cos*θ* where *θ* is the contact angle between the nucleation phase and the substrate. The factor f for convex surface of par‐ ticle could read as

$$f(\eta, \mathbf{x}) = \frac{1}{2} \left\{ 1 + \left( \frac{1 - \eta \mathbf{x}}{g} \right)^3 + \mathbf{x}^3 \left[ 2 - 3 \left( \frac{\mathbf{x} - \eta}{g} \right) + \left( \frac{\mathbf{x} - \eta}{g} \right)^3 \right] + 3\eta \mathbf{x}^2 \left( \frac{\mathbf{x} - \eta}{g} - 1 \right) \right\} \tag{5}$$

where *g* =(1 + *x* <sup>2</sup> −2*ηx*)1/2

marize the typical physical and mathematical model for ice nucleation and growth in cryobi‐ ology area. In addition, the heat and solute transport, and interaction between the cell and

Freezing often denotes a process of new phase defined by ice creation via nucleation and subsequent crystal growth from supercooled water. This process undergoes three stages including the chemical diffusion stage for water molecule transportation, the surface kinetic stage for hydrogen bonds formation, and the final heat conduction stage for heat releasing and diffusing. It implies that the ice nucleation—whose forma‐ tion determines the growth of ice crystal to a large extent—turns out to be the initial and crucial step during crystallization. In terms of whether there exist foreign bodies or not, the nucleation is usually classified as homogeneous or heterogeneous cases. By definition, homogeneous ice nucleation requests a severe condition that the liquid wa‐ ter is absolutely dust-free, or there exist foreign particles indeed, but whose sizes are small enough to be negligible. Since this rigid demand is hard to be satisfied, the nu‐ cleation is in fact regarded as heterogeneous in most cases and may be influenced by many external factors, such as the size of foreign particles, surface roughness and the effective size of the nucleators at different supercoolings. There exist a great many bar‐ riers to build up an integrated theory on nucleation, because of the dissimilar proper‐ ties might exhibit between bulk fluid and the micro/nanoscale molecules. Up to now, the two major theoretical approaches to characterize the formation of ice can be clas‐

sified as thermodynamic approximation and molecular simulations [4].

The ice nucleation does not happen necessarily if only temperature is below the freezing point. This process is determined by the free energy of system, which implies that the ice germ growth or disappearance depends on the free energy decrease or increase. For constant rate cooling protocols, the temperature of system is given by *T* (*t*)=*T*<sup>0</sup> − *Bt*, where *B* is the cooling rate. Considering a cluster of ice denoted by *β* with *i* molecules from its mother water phase *α*, the driving force of nucleation is derived from the Gibbs

1/3 2/3 2/3 (36 ) ( ) *Gi G i i t*

where *ΔGt* is negative and denotes the Gibbs free-energy difference between *α* and *β* phase

 p

D= D+

ab b

(1)

is the molecular vol‐

 sn

b

per volume, according to classical homogeneous nucleation theory; *ν <sup>β</sup>*

n

**2.1. Thermodynamic model of ice nucleation**

free energy difference*ΔGi* [5]

ice during freezing are also discussed.

**2. Ice nucleation**

186 Advanced Topics on Crystal Growth

The homogeneous nucleation rate *J* considering both diffusion barrier*ΔGi*' and nucleation barrier*ΔGi c* could be given by

$$J = \frac{nkT}{h} \exp\left(-\frac{\Delta G\_{i^\*}}{kT}\right) \exp\left(-\frac{\Delta G\_{i\_c}}{kT}\right) \tag{6}$$

where *h* is Planck constant, *k* the Boltzmann constant, and *n* denotes molecules per unit vol‐ ume of mother phase. Thus the heterogeneous nucleation rate *Jhet* could be obtained by modifying *ΔGi c* as *ΔGi c*,*h* in the expression of *J*.

#### **2.2. Diffusion limited model of ice growth inside cell**

Karlssonet al. [6] has developed a physiochemical theory of ice growth inside biological cells by coupling the water transport models, theory of ice nucleation, and theory of crystal growth. Considering a given cell with the control volume of *Vcv*(*t*), which is composed of water-NaCl-glycerol solution. Evolution of the number of ice nuclei denoted by *N* (*t*) in a given cell is a stochastic process given as

$$N(t) = \text{int}\left[\int\_0^t f(t)V\_{cv}(t)dt\right] \tag{7}$$

According the water-transport model developed by Mazur [7] to cells containing a ternary

<sup>+</sup> 1 1 ln

*w g cv s s g g w s s g*

*R TT V nv nv v n n*

where *L <sup>p</sup>* is the water permeability; *A* is the cell membrane surface area; *Rg* gas constant; *vw*, *vs* and *vg* denote specific volumes of water, NaCl and the glycerol, respectively; *ns* and *ng* are the number of moles of NaCl and glycerol in cell. Coupling the Eq. (10) and Eq. (11), one could calculate the crystallized volume fraction during freezing. Recently more sophisticat‐ ed models [8-10] of intracellular ice formation and its growth have been developed to obtain

Molecular dynamics (MD) is a powerful computer simulation technique to track time evolu‐ tion of a system of interacting particles, such as atoms, molecules and coarse-grained parti‐ cles. It has been widely used in the study of the structure, dynamics, and thermodynamics of water phase change and their complexes addressing a variety of issues including ice nu‐ cleation. The basic idea of the MD simulations is to generate the atomic trajectories of a sys‐ tem of finite particles by numerical integration, of a set of Newton's equations of motion for all particles in the system with certain given boundary conditions, an initial set of positions and velocities. In addition, the potential energy of the particles system needs to be specified in order to describe the interaction between atoms or molecules. The detailed potential ener‐ gy expression and its parameter are derived from both experimental work and high-level

The dynamics of the ice–water interface determines the ice nucleation and growth such that considerable MD methods have been developed. Karim and Haymet [11, 12] have used the TIP4P water model to simulate the contact of the basal plane of ice with water and investigate dynamics of the ice–water interface. Their results show that the interface was found to be stable at 240 K. In addition, there appears a gradual change in the ori‐ entation and mobility of the molecules in the interfacial region. A very extensive MD simulation by Matsumoto et al. [13] predicted the freezing of pure water at different su‐ percooling. They found that the occurrence of the ice nucleation was impressed exten‐ sively by the number of spontaneously developing long-lived hydrogen bonds at the same location. Handel et al. [14] had performed direct calculations of the ice Ih-water in‐ terfacial free energy for the TIP4P model by extending the cleaving method to molecular

MD simulations also could be used to investigate the effects of salt on ice nucleation, which is more close to biological environments. Smith et al. [15] had studied the poten‐ tial of mean force of single ions as a function of the distance from the ice-water inter‐

é ù <sup>D</sup> æ ö æ ö - ê ú ç ÷ = - ç ÷ - - ê ú ç ÷ è ø - ë û è ø

( )

++ +

n

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

(11)

189

0 ( ) ( )

water-NaCl-glycerol solution, one could derive the time evolution of cell volume

*pg f cv s s g g cv*

*L AR T H V nv nv dV*

*dt*

n

more close to the real cell environments.

quantum mechanical calculations.

systems, which agrees with the experiments results.

**2.3. Molecular simulation of ice nucleation**

where int denote the nearest integer truncated from the ensemble average. When the crystal‐ lizing phase is of different composition than the bulk solution, such as the case with ice crys‐ tallization from aqueous solutions, the rate of crystal growth is believed to be limited by the diffusion water molecules to the ice-solution interface. Considering a spherically symmetric ice crystal under isothermal conditions, its radius at time t could be given as [6]

$$r\left(t, t\_i\right) = \left[\int\_{T\left(t\_i\right)}^{T\left(t\right)} \rho^2 \overline{D}\left(-\frac{dT}{B}\right) dt\right]^{1/2} \tag{8}$$

where *<sup>D</sup>*¯ is the effective diffusivity of water and is calculated according to the viscosity of the intracellular solution. *φ* is the nondimensional crystal growth parameter and is de‐ rived from the equation obtained by curve fitting in advance based on the nondimen‐ sional supersaturation [6]. If the crystals within a given cell are small enough such that neighboring crystals and their depletion regions do not overlap, the total crystallized vol‐ ume could be given as

$$V\_{ice}\left(t\right) = \sum\_{i=1}^{N(t)} \frac{4\pi}{3} r^3\left(t, t\_i\right) \tag{9}$$

Thus the crystallized volume fraction considering the impingement appearance could be given as

$$X\_{ice}\left(t\right) = 1\cdot \exp\left(\frac{V\_{ice}}{V\_{cv}}\right) \tag{10}$$

According the water-transport model developed by Mazur [7] to cells containing a ternary water-NaCl-glycerol solution, one could derive the time evolution of cell volume

$$\frac{dV\_{cv}}{dt} = -\frac{L\_p A R\_g T}{\nu\_w} \left[ \frac{\Delta H\_f}{R\_g} \left( \frac{1}{T\_0} - \frac{1}{T} \right) - \ln \left( \frac{V\_{cv} - \left( n\_s \upsilon\_s \star n\_g \upsilon\_g \right)}{V\_{cv} - \left( n\_s \upsilon\_s \star n\_g \upsilon\_g \right) + \upsilon\_w \left( n\_s \upsilon\_s \star n\_g \right)} \right) \right] \tag{11}$$

where *L <sup>p</sup>* is the water permeability; *A* is the cell membrane surface area; *Rg* gas constant; *vw*, *vs* and *vg* denote specific volumes of water, NaCl and the glycerol, respectively; *ns* and *ng* are the number of moles of NaCl and glycerol in cell. Coupling the Eq. (10) and Eq. (11), one could calculate the crystallized volume fraction during freezing. Recently more sophisticat‐ ed models [8-10] of intracellular ice formation and its growth have been developed to obtain more close to the real cell environments.

### **2.3. Molecular simulation of ice nucleation**

where *h* is Planck constant, *k* the Boltzmann constant, and *n* denotes molecules per unit vol‐ ume of mother phase. Thus the heterogeneous nucleation rate *Jhet* could be obtained by

Karlssonet al. [6] has developed a physiochemical theory of ice growth inside biological cells by coupling the water transport models, theory of ice nucleation, and theory of crystal growth. Considering a given cell with the control volume of *Vcv*(*t*), which is composed of water-NaCl-glycerol solution. Evolution of the number of ice nuclei denoted by *N* (*t*) in a

> ( ) ( ) ( ) <sup>0</sup> int *<sup>t</sup> N t J t V t dt cv* é ù <sup>=</sup> ê ú ë û

ice crystal under isothermal conditions, its radius at time t could be given as [6]

*dT rtt D dt*

j

é ù æ ö <sup>=</sup> ê ú ç ÷ -

where *<sup>D</sup>*¯ is the effective diffusivity of water and is calculated according to the viscosity of the intracellular solution. *φ* is the nondimensional crystal growth parameter and is de‐ rived from the equation obtained by curve fitting in advance based on the nondimen‐ sional supersaturation [6]. If the crystals within a given cell are small enough such that neighboring crystals and their depletion regions do not overlap, the total crystallized vol‐

> ( ) ( ) ( )

> > =1

*N t ice i i V t r tt* p

3

<sup>4</sup> , <sup>3</sup>

Thus the crystallized volume fraction considering the impingement appearance could be

*cv V*

*V* æ ö <sup>=</sup> ç ÷ ç ÷ è ø

( ) 1-exp *ice*

*ice*

*X t*

( ) ( )

<sup>2</sup> , *i T t i T t*

where int denote the nearest integer truncated from the ensemble average. When the crystal‐ lizing phase is of different composition than the bulk solution, such as the case with ice crys‐ tallization from aqueous solutions, the rate of crystal growth is believed to be limited by the diffusion water molecules to the ice-solution interface. Considering a spherically symmetric

( ) 1 2

*B*

ò (7)

ë û è ø ò (8)

<sup>=</sup> å (9)

(10)

*c*,*h* in the expression of *J*.

**2.2. Diffusion limited model of ice growth inside cell**

modifying *ΔGi*

188 Advanced Topics on Crystal Growth

*c* as *ΔGi*

ume could be given as

given as

given cell is a stochastic process given as

Molecular dynamics (MD) is a powerful computer simulation technique to track time evolu‐ tion of a system of interacting particles, such as atoms, molecules and coarse-grained parti‐ cles. It has been widely used in the study of the structure, dynamics, and thermodynamics of water phase change and their complexes addressing a variety of issues including ice nu‐ cleation. The basic idea of the MD simulations is to generate the atomic trajectories of a sys‐ tem of finite particles by numerical integration, of a set of Newton's equations of motion for all particles in the system with certain given boundary conditions, an initial set of positions and velocities. In addition, the potential energy of the particles system needs to be specified in order to describe the interaction between atoms or molecules. The detailed potential ener‐ gy expression and its parameter are derived from both experimental work and high-level quantum mechanical calculations.

The dynamics of the ice–water interface determines the ice nucleation and growth such that considerable MD methods have been developed. Karim and Haymet [11, 12] have used the TIP4P water model to simulate the contact of the basal plane of ice with water and investigate dynamics of the ice–water interface. Their results show that the interface was found to be stable at 240 K. In addition, there appears a gradual change in the ori‐ entation and mobility of the molecules in the interfacial region. A very extensive MD simulation by Matsumoto et al. [13] predicted the freezing of pure water at different su‐ percooling. They found that the occurrence of the ice nucleation was impressed exten‐ sively by the number of spontaneously developing long-lived hydrogen bonds at the same location. Handel et al. [14] had performed direct calculations of the ice Ih-water in‐ terfacial free energy for the TIP4P model by extending the cleaving method to molecular systems, which agrees with the experiments results.

MD simulations also could be used to investigate the effects of salt on ice nucleation, which is more close to biological environments. Smith et al. [15] had studied the poten‐ tial of mean force of single ions as a function of the distance from the ice-water inter‐ face, and the free energy of transfer of Na<sup>+</sup> and Cl ions from bulk water into the ice and to the interface. Their calculations suggested that the interface may have a potential that is attractive for Cl and repulsive for Na<sup>+</sup> . The investigation from Vrbka and Jungwirth [16], and Carignano et al. [17] indicated that the growth rate of ice was much lower in solutions with salt than in pure water.

retical strategy to characterize the ice crystal growth under electrostatic field and describe

A dimensionless phase field model [33] is adopted to characterize the ice growth. A nonconserved phase field *ϕ* (*ϕ* =0 is corresponding to ice and *ϕ* =1 to water) describes the transi‐ tion from water to ice, which is driven by undercooling. A conserved dimensionless temperature field *u* is used to characterize the energy variation due to temperature diffusion

( ( ) ) <sup>2</sup> (1 ) 0.5 30 (1 ) ,

 e f f

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

191

*xyt* - =Ñ - - + (13)

<sup>D</sup> = == = (15)

, *TM* =273.15*K*, *cp* =4.212*<sup>J</sup>* / *cm*<sup>3</sup> <sup>⋅</sup> *<sup>K</sup>* and

*cm*. Thus the phase field

) under the

and *Δt* =10−<sup>4</sup>

(14)

*<sup>t</sup> mB m aSu* - = Ñ× Ñ + - - + - é ù ë û (12)

 y

2

 g q

 ff

( ) ( ) <sup>2</sup> 2 12 30 1 ,, , *t t u uS* f

( ) ( ) ( ) ( )

2

g q

g qgq

direction and *x* axis. Four dimensionless parameters are defined as follow

*m aS*

e

d

( ) ( ) ( )

In the above equations, the term *ψ* represents the thermal noise, which indicates the fluctua‐ tion of temperature, and *u* =(*T* −*TM* ) / *ΔT* . *γ*(*θ*)=1 + *ς*cos*n*(*θ* −*θ*0) denotes dimensionless ani‐ sotropic surface tension and *θ* =arctan(*ϕ<sup>y</sup>* / *ϕx*), where *θ*0 is the angle between reference

> <sup>2</sup> ,, , , <sup>12</sup> *<sup>m</sup> <sup>p</sup> <sup>T</sup> <sup>w</sup> c T*

*Lw d L*

The physical properties of pure water [30] used in the simulation are listed as:

model parameters for all the simulations, according to Ref. [33], can be fixed with *ε* =0.005, *m*=0.035, *a* =400 and *S* =0.5, which is corresponding to the cooling temperature 233.56*K*.

An explicit time-differencing scheme is adopted to solve *ϕ*-equation, and implicit Crank-Nicolson algorithm is chosen for solving *u*-equation. All the simulations are performed in a

. The length scale *w* is considered as 2.7×10−<sup>4</sup>

 g qgq

æ ö - ¢ <sup>=</sup> ç ÷ ¢ è ø

 ff

multiple ice nucleuses' competitive growth behavior.

 q f

*B*

ns

[33].

*cm*<sup>2</sup> /*s*, *L* =333.5*J* / *cm*<sup>3</sup>

two dimensional domain with 1500×1500 grid points (*Δx* =10−<sup>2</sup>

**3.1. Ice crystal growth from pure water**

*J* / *cm*<sup>2</sup>

k

q

 ef

and phase transition. It reads as

f

where *<sup>d</sup>* <sup>=</sup>*cpσTm* / *<sup>ρ</sup><sup>L</sup>* <sup>2</sup>

*κ* =1.31×10−<sup>3</sup>

*σ* =7.654×10−<sup>6</sup>

where

The mechanism of ice growth inhibition by antifreeze proteins (AFP) is still poorly un‐ derstood, owing to the difficulty in elucidating the molecular-scale structure of an ice– water interface, and in experimentally observing the structure and dynamics of AFP at an ice–water interface. MD is one of most popular methods that currently have the po‐ tential to clarify the mechanism of ice growth inhibition by AFP at the molecular level [18-21]. Wathen et al. [19] developed a MD which combines molecular representation and detailed energetics calculations of molecular mechanics techniques with the less-so‐ phisticated probabilistic approach to study systems containing millions of water mole‐ cules undergoing billions of interactions. Such model could enable the 3-D shape and surface properties of the molecules to directly affect crystal formation. They have applied this technique to study the inhibitory effects of antifreeze proteins (AFPs) (from both fish and insect) on ice-crystal formation, including the replication of ice-etching patterns, icegrowth inhibition, and specific AFP-induced ice morphologies. Their work suggests that the degree of AFP activity results more from AFP ice-binding orientation than from AFP ice-binding strength. The simulation results were consistent with experimental observa‐ tions very well. The effects of cryoprotectants such as ethylene glycol and glycerol solu‐ tions on ice nucleation have also been investigated [22, 23].

The electrostatic field was recently proposed to improve the output of a cryopreservation process [24-26]. As the molecular dynamics simulated [27], a DC field with magnitude of *E* =5.0×10<sup>9</sup> V / m could induce crystallization of supercooled water. Recently, experiments [24-26, 28] have demonstrated that the external electric field can effectively affect the water nucleation and ice growth. A lower electrostatic field strength in the range of *E* =103∼10<sup>5</sup> V / m, which is easily obtained in experiments, can induce increase in the nuclea‐ tion temperature of water [25] and the asymmetric growth of ice [24].

### **3. Ice crystal growth using phase field method**

To better understand the detailed events in cell damage due to extracellular and intracellu‐ lar ice formation during freezing, it is rather important to model and predict the micro-scale ice crystal formation and growth behavior, and thus characterize the effects of several core factors, such as undercooling, anisotropy, thermal noise, and external electric field etc. Among various ways ever developed, the phase field method is rather flexible in character‐ izing the crystal formation and has been extensively adopted to model solidification process in material science [29]. Recently some authors [30, 31] begun to pay attention to possible application of this method in cryopreservation. He and Liu [32] extended deeply this theo‐ retical strategy to characterize the ice crystal growth under electrostatic field and describe multiple ice nucleuses' competitive growth behavior.

A dimensionless phase field model [33] is adopted to characterize the ice growth. A nonconserved phase field *ϕ* (*ϕ* =0 is corresponding to ice and *ϕ* =1 to water) describes the transi‐ tion from water to ice, which is driven by undercooling. A conserved dimensionless temperature field *u* is used to characterize the energy variation due to temperature diffusion and phase transition. It reads as

$$\phi\_t = m\nabla \cdot \left( B\left(\theta\right) \nabla \phi \right) + m\omega^{-2} \phi (1-\phi) \left[ \phi - 0.5 + 30\alpha a S u \phi (1-\phi) \right],\tag{12}$$

$$\ln \mu\_t = \nabla^2 \mu - \mathbf{S}^{-1} \Im 0 \phi^2 \left( \mathbf{1} - \phi \right)^2 \phi\_t + \nu \left( \mathbf{x}, y, t \right), \tag{13}$$

where

face, and the free energy of transfer of Na<sup>+</sup>

solutions with salt than in pure water.

and repulsive for Na<sup>+</sup>

tions on ice nucleation have also been investigated [22, 23].

tion temperature of water [25] and the asymmetric growth of ice [24].

**3. Ice crystal growth using phase field method**

is attractive for Cl-

190 Advanced Topics on Crystal Growth

*E* =5.0×10<sup>9</sup>

*E* =103∼10<sup>5</sup>

and Cl-

to the interface. Their calculations suggested that the interface may have a potential that

[16], and Carignano et al. [17] indicated that the growth rate of ice was much lower in

The mechanism of ice growth inhibition by antifreeze proteins (AFP) is still poorly un‐ derstood, owing to the difficulty in elucidating the molecular-scale structure of an ice– water interface, and in experimentally observing the structure and dynamics of AFP at an ice–water interface. MD is one of most popular methods that currently have the po‐ tential to clarify the mechanism of ice growth inhibition by AFP at the molecular level [18-21]. Wathen et al. [19] developed a MD which combines molecular representation and detailed energetics calculations of molecular mechanics techniques with the less-so‐ phisticated probabilistic approach to study systems containing millions of water mole‐ cules undergoing billions of interactions. Such model could enable the 3-D shape and surface properties of the molecules to directly affect crystal formation. They have applied this technique to study the inhibitory effects of antifreeze proteins (AFPs) (from both fish and insect) on ice-crystal formation, including the replication of ice-etching patterns, icegrowth inhibition, and specific AFP-induced ice morphologies. Their work suggests that the degree of AFP activity results more from AFP ice-binding orientation than from AFP ice-binding strength. The simulation results were consistent with experimental observa‐ tions very well. The effects of cryoprotectants such as ethylene glycol and glycerol solu‐

The electrostatic field was recently proposed to improve the output of a cryopreservation process [24-26]. As the molecular dynamics simulated [27], a DC field with magnitude of

[24-26, 28] have demonstrated that the external electric field can effectively affect the water nucleation and ice growth. A lower electrostatic field strength in the range of

To better understand the detailed events in cell damage due to extracellular and intracellu‐ lar ice formation during freezing, it is rather important to model and predict the micro-scale ice crystal formation and growth behavior, and thus characterize the effects of several core factors, such as undercooling, anisotropy, thermal noise, and external electric field etc. Among various ways ever developed, the phase field method is rather flexible in character‐ izing the crystal formation and has been extensively adopted to model solidification process in material science [29]. Recently some authors [30, 31] begun to pay attention to possible application of this method in cryopreservation. He and Liu [32] extended deeply this theo‐

V / m could induce crystallization of supercooled water. Recently, experiments

V / m, which is easily obtained in experiments, can induce increase in the nuclea‐

ions from bulk water into the ice and

. The investigation from Vrbka and Jungwirth

$$B(\theta) = \begin{pmatrix} \gamma^2(\theta) & -\gamma'(\theta)\gamma(\theta) \\ \gamma'(\theta)\gamma(\theta) & \gamma^2(\theta) \end{pmatrix} \tag{14}$$

In the above equations, the term *ψ* represents the thermal noise, which indicates the fluctua‐ tion of temperature, and *u* =(*T* −*TM* ) / *ΔT* . *γ*(*θ*)=1 + *ς*cos*n*(*θ* −*θ*0) denotes dimensionless ani‐ sotropic surface tension and *θ* =arctan(*ϕ<sup>y</sup>* / *ϕx*), where *θ*0 is the angle between reference direction and *x* axis. Four dimensionless parameters are defined as follow

$$m = \frac{\nu \sigma T\_m}{\kappa L}, \; \varepsilon = \frac{\delta}{w}, \; a = \frac{\sqrt{2}w}{12d}, \; S = \frac{c\_p \Delta T}{L}, \tag{15}$$

where *<sup>d</sup>* <sup>=</sup>*cpσTm* / *<sup>ρ</sup><sup>L</sup>* <sup>2</sup> [33].

#### **3.1. Ice crystal growth from pure water**

The physical properties of pure water [30] used in the simulation are listed as: *κ* =1.31×10−<sup>3</sup> *cm*<sup>2</sup> /*s*, *L* =333.5*J* / *cm*<sup>3</sup> , *TM* =273.15*K*, *cp* =4.212*<sup>J</sup>* / *cm*<sup>3</sup> <sup>⋅</sup> *<sup>K</sup>* and *σ* =7.654×10−<sup>6</sup> *J* / *cm*<sup>2</sup> . The length scale *w* is considered as 2.7×10−<sup>4</sup> *cm*. Thus the phase field model parameters for all the simulations, according to Ref. [33], can be fixed with *ε* =0.005, *m*=0.035, *a* =400 and *S* =0.5, which is corresponding to the cooling temperature 233.56*K*.

An explicit time-differencing scheme is adopted to solve *ϕ*-equation, and implicit Crank-Nicolson algorithm is chosen for solving *u*-equation. All the simulations are performed in a two dimensional domain with 1500×1500 grid points (*Δx* =10−<sup>2</sup> and *Δt* =10−<sup>4</sup> ) under the Neumann boundary conditions. The initial condition of the phase field for single ice nucleus is set for *<sup>ϕ</sup>* <sup>=</sup> <sup>1</sup> <sup>2</sup> tanh( *<sup>r</sup>* <sup>−</sup>*rc* 2 2*ε* ) + 1 . The initialization of the dimensionless temperature is *u* = −*ϕ* for all the simulations, which is corresponding to initial temperature *TM* in solid and *TM* −*ΔT* in liquid.

search. The influence of low electrostatic field strength on water nucleation and ice growth may be attributed to dipole polarization of the water molecules by the electrostatic field. Un‐ der the electrostatic field, the phase transition must overcome excess energy *ΔG <sup>E</sup>* =*EΔμ*,

and can be negligible. In other words, the water molecules with dipole moments may follow the Boltzmann distribution function *p* = *A*exp(*μwE*cos*φ* / *kT* ) [25], and get in the most stable state with the maximum energy along the direction of the electrostatic field (*φ* =0). Although for the low electrostatic field magnitude, the value variation of *p* with respect to *φ* is very small in bulk phase, the effects of electrostatic field may become largely strong at the inter‐ face and affect dynamical state of interface [24]. The polar water molecules (or clusters) may be torqued, rearranged, and forced to join the ice lattice via a special orientation and posi‐ tion under action of the electrostatic field. Thus the rate of attachment (or detachment) of water molecules on the ice interface is affected by external electrostatic field. In order to characterize such influence of the electrostatic field, the surface tension should be modified

From the above dissection, the effects of electrostatic field with low magnitude can induce change of the surface tension with the form *σ<sup>E</sup>* =exp(*λ*cos(*θ* −*φ*))*σ*. Thus in phase field model, the dimensionless anisotropic surface tension should be modified as *γE*(*θ*)=exp(*λ*cos(*θ* −*φ*))*γ*(*θ*). The value of *λ* is proportional to *μwE* / *kT* , which characterizes the strength of electric field. In our simulations, the value of *λ* is chosen in the range 0.02~0.2

**Figure 2.** The influence of electric field on ice crystal growth at (a) γ =0.04, *n* =4 and λ =0.1, at (b) γ =0.02, *n* =6 and

λ =0.1, at (c), (d) γ =0.02, *n* =6 and λ =0.2. The angle θ<sup>0</sup> =0 for (a-c) and θ<sup>0</sup> =π / 9 for (d). [32]

as exp(*λ*cos(*θ* −*φ*))*σ* during the ice growth, where *λ* is proportional to *μwE* / *kT* .

and *φ* =0, which indicates that the orientation of electric field is along positive *x* axis.

is the difference of electrical dipole moments between water and ice.

V/m), *ΔG <sup>E</sup>* is much smaller than *L*

http://dx.doi.org/10.5772/54098

193

Modeling Ice Crystal Formation of Water in Biological System

where *Δμ* =*μw* −*μi*

However, for the low electrostatic field strength (*E* <10<sup>5</sup>

**Figure 1.** Micro-scale dendritic ice crystal patterns with (a) 4-fold structures without thermal noise, (b) 6-fold struc‐ tures with thermal noise, (c) two seeds and (d) four seeds competitive growth. The parameters γ =0.04 and *n* =4 in (a), (d), and γ =0.02 and *n* =6 in (b), (c). The angle θ<sup>0</sup> =0 for (a-d). [32]

Simulations indicate that the influence of surface tension anisotropy on the shape of ice is more evident compared with that on the growth rate. The variation of the anisotropy strength does not significantly change the velocity at the dendritic tip, which is strongly af‐ fected by the local cooling strength. Four-fold and six-fold dendritic ice crystals are often ob‐ served in experiments. Fig. 1(a) shows a four-fold dendritic ice crystal evolving from an initial circle shaped seed [32]. The thermal noise, which often induces the temperature fluc‐ tuation on the interface, is considered as the main reason of sidebranching growth [34]. A six-fold dendritic ice with sidebranching is shown in Fig. 1(b). It is necessary to consider the influence of the other crystal seeds on the ice growth. Fig. 1(c), (d) represent two seeds and four seeds competitive growth, which results in the crease of growth for both collided den‐ drites and preferred growth in directions of unhindered grains.

### **3.2. The effects of external electric on Ice crystal growth**

The growth rate depends on the dynamical state of the ice-water interface, such as the sur‐ face tension and the temperature at interface. In addition, the effects of the external electric field on the ice formation have been demonstrated in both theoretical and experimental re‐ search. The influence of low electrostatic field strength on water nucleation and ice growth may be attributed to dipole polarization of the water molecules by the electrostatic field. Un‐ der the electrostatic field, the phase transition must overcome excess energy *ΔG <sup>E</sup>* =*EΔμ*, where *Δμ* =*μw* −*μi* is the difference of electrical dipole moments between water and ice. However, for the low electrostatic field strength (*E* <10<sup>5</sup> V/m), *ΔG <sup>E</sup>* is much smaller than *L* and can be negligible. In other words, the water molecules with dipole moments may follow the Boltzmann distribution function *p* = *A*exp(*μwE*cos*φ* / *kT* ) [25], and get in the most stable state with the maximum energy along the direction of the electrostatic field (*φ* =0). Although for the low electrostatic field magnitude, the value variation of *p* with respect to *φ* is very small in bulk phase, the effects of electrostatic field may become largely strong at the inter‐ face and affect dynamical state of interface [24]. The polar water molecules (or clusters) may be torqued, rearranged, and forced to join the ice lattice via a special orientation and posi‐ tion under action of the electrostatic field. Thus the rate of attachment (or detachment) of water molecules on the ice interface is affected by external electrostatic field. In order to characterize such influence of the electrostatic field, the surface tension should be modified as exp(*λ*cos(*θ* −*φ*))*σ* during the ice growth, where *λ* is proportional to *μwE* / *kT* .

Neumann boundary conditions. The initial condition of the phase field for single ice nucleus

for all the simulations, which is corresponding to initial temperature *TM* in solid and

**Figure 1.** Micro-scale dendritic ice crystal patterns with (a) 4-fold structures without thermal noise, (b) 6-fold struc‐ tures with thermal noise, (c) two seeds and (d) four seeds competitive growth. The parameters γ =0.04 and *n* =4 in (a),

Simulations indicate that the influence of surface tension anisotropy on the shape of ice is more evident compared with that on the growth rate. The variation of the anisotropy strength does not significantly change the velocity at the dendritic tip, which is strongly af‐ fected by the local cooling strength. Four-fold and six-fold dendritic ice crystals are often ob‐ served in experiments. Fig. 1(a) shows a four-fold dendritic ice crystal evolving from an initial circle shaped seed [32]. The thermal noise, which often induces the temperature fluc‐ tuation on the interface, is considered as the main reason of sidebranching growth [34]. A six-fold dendritic ice with sidebranching is shown in Fig. 1(b). It is necessary to consider the influence of the other crystal seeds on the ice growth. Fig. 1(c), (d) represent two seeds and four seeds competitive growth, which results in the crease of growth for both collided den‐

The growth rate depends on the dynamical state of the ice-water interface, such as the sur‐ face tension and the temperature at interface. In addition, the effects of the external electric field on the ice formation have been demonstrated in both theoretical and experimental re‐

) + 1 . The initialization of the dimensionless temperature is *u* = −*ϕ*

is set for *<sup>ϕ</sup>* <sup>=</sup> <sup>1</sup>

*TM* −*ΔT* in liquid.

192 Advanced Topics on Crystal Growth

<sup>2</sup> tanh( *<sup>r</sup>* <sup>−</sup>*rc*

2 2*ε*

(d), and γ =0.02 and *n* =6 in (b), (c). The angle θ<sup>0</sup> =0 for (a-d). [32]

drites and preferred growth in directions of unhindered grains.

**3.2. The effects of external electric on Ice crystal growth**

From the above dissection, the effects of electrostatic field with low magnitude can induce change of the surface tension with the form *σ<sup>E</sup>* =exp(*λ*cos(*θ* −*φ*))*σ*. Thus in phase field model, the dimensionless anisotropic surface tension should be modified as *γE*(*θ*)=exp(*λ*cos(*θ* −*φ*))*γ*(*θ*). The value of *λ* is proportional to *μwE* / *kT* , which characterizes the strength of electric field. In our simulations, the value of *λ* is chosen in the range 0.02~0.2 and *φ* =0, which indicates that the orientation of electric field is along positive *x* axis.

**Figure 2.** The influence of electric field on ice crystal growth at (a) γ =0.04, *n* =4 and λ =0.1, at (b) γ =0.02, *n* =6 and λ =0.1, at (c), (d) γ =0.02, *n* =6 and λ =0.2. The angle θ<sup>0</sup> =0 for (a-c) and θ<sup>0</sup> =π / 9 for (d). [32]

Under the electrostatic field, the growth of dendritic ice crystals represents asymmetry be‐ havior. Fig 2 shows that the main branches parallel to the electrostatic field grow faster than the other branches [32]. One can also find that the growth of sidebranching is strengthened along the direction of electric field and weakened along the inverse direction. The above in‐ fluence becomes clearer by increasing the value of *λ* (comparing Fig. 2(b) with (c)), which indicates that a more strong magnitude of electrostatic field would strengthen the asymmet‐ ric growth of ice. Fig. 2(d) shows how the reference direction affects the ice growth under electric field. The above results concerning the effects of electric field are in good accordance with the qualitative analysis in Ref. [25] and experimental observation in Ref. [24].

is also considered and shown to yield significant differences in the cell response in com‐

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

195

Another important process is from Chang et al. [36]. They developed a modified level set method for modeling the interaction of biological cells with ice front during a directional solidification. The simulation result has shown that different types of cell interaction with the ice interface can lead to significant differences in the final volume of individual cells. The cell is easily exposed to high concentrations of electrolytes in the interdendritic regions due to cell entrapment. Cell pushing can be used to control the dehydration of the cell more easi‐ ly, which avoids exposing the cell to high concentrations that develop in the interdendritic spaces. Cell engulfment protects the cell from high concentration regions and leaves the cell vulnerable to intracellular ice formation. However, it could be used to control the dehydra‐ tion of the cell and provide limited interaction with high electrolyte concentrations through combination of cell pushing and cell engulfment. The numerical results also indicate that the interactions between cells are very important for cryopreservation conditions, where sus‐ pensions have higher volume fractions of cells. Level set method could deal with the interac‐ tion of clusters of cells with ice by assuming that the clusters act like a single particle with a radius that represents the overall size of the cell cluster. Increasing the radius of the repre‐ sentative particle would increase the probability of engulfment and pushing by the inter‐ face. Although the intermolecular forces between cells and the fluid flow around various arrays of single cells and cell clusters have not been included in the current model, these phenomena could be easily included in the simulations in order to study their effects on the

The present chapter has presented an overview on the theoretical modeling of the micro‐ scale phase change of biological system subject to freezing with special focus on the ice crys‐ tal formation and growth. In addition, the heat and solute transport, and interaction between the cell and ice during freezing have also been discussed. Although considerable investigation focus on ice crystal formation and its growth, there still exists some outstand‐ ing questions [37-40] including ice structure and ice phase diagram. Here we discuss some important issues and challenges about ice formation in biological system in brief. Firstly, the effects of cell membrane on ice nucleation and growth inside a cell with restricted space dur‐ ing freezing have not been investigated deeply, which is in fact very important to clarify the puzzle about the ice propagation between cells. In addition, the improved understanding on the effects of cell–cell interactions suggested that the mechanism of ice propagation should be mediated to some extent, but not at all, by the gap junction, if uncontrolled factors are eliminated. Secondly, a detailed three dimensional model needs to be developed to charac‐ terize the complicated freezing process that coupling micro heat and mass transfer, and ice nucleation and growth, which involves the moving water-ice interface and deformed cell membrane. Thirdly, the advanced experimental measures need to be designed to determine

parison to the isothermal case.

partitioning of cells by the ice interface.

**5. Summary**

It is noteworthy that the biological media is saline solution including various ions, such as Na+ and Cl<sup>−</sup> . Under the electrostatic field, these positive and negative ions have com‐ pletely different movement and gather near the corresponding electrodes, which induce non-uniform distribution of ions. In addition, the heterogeneous structure of ions, such as their different size, also leads to different ions' behaviors. These behaviors induced by the electrostatic field will affect the dynamical behavior of ice interface and strengthen the asymmetric growth of ice [4]. Thus, the effects of the ions under the electrostatic field are also indicated via modifying the surface tension. The effects of ions coupling with the electrostatic field on the dynamical behavior of ice interface should be incorpo‐ rated to a general model. The present investigation focuses on the effects of the electro‐ static field on the ice interface but is not to completely characterize the effects of ions, which will be addressed in further work.

### **4. The interaction between ice growths with cell**

The interaction between the cell and ice during freezing is very important to investigate the cell cryoinjury during freezing. However, it is still a rather difficult task to track the evolu‐ tion of the solidification front and characterize the interaction between the cell membrane and ice in detail. Recently, the progress [35, 36] gave some promising research directions to solve this issue.

Mao et al. [35] developed a sharp interface method to simulate the response of a biologi‐ cal cell during freezing. The cell is modeled as an aqueous salt solution surrounded by a semi-permeable membrane. The concentration and temperature fields both inside and outside a single cell are computed taking into account heat transfer, mass diffusion, membrane transport, and evolution of the solidification front. The external ice front is computed for both stable and unstable growth modes. It is shown that for the particular geometry chosen in this study, the instabilities on the front and the diffusional transport have only modest effects on the cell response. For the cooling conditions, solute and cell property parameters used, the low Peclet regime applies. The computational results are therefore validated against the conventional membrane-limited transport (Mazur) model. Good agreement of the simulation results with the Mazur model are obtained for a wide range of cooling rates and membrane permeabilities. A spatially non-isothermal situation is also considered and shown to yield significant differences in the cell response in com‐ parison to the isothermal case.

Another important process is from Chang et al. [36]. They developed a modified level set method for modeling the interaction of biological cells with ice front during a directional solidification. The simulation result has shown that different types of cell interaction with the ice interface can lead to significant differences in the final volume of individual cells. The cell is easily exposed to high concentrations of electrolytes in the interdendritic regions due to cell entrapment. Cell pushing can be used to control the dehydration of the cell more easi‐ ly, which avoids exposing the cell to high concentrations that develop in the interdendritic spaces. Cell engulfment protects the cell from high concentration regions and leaves the cell vulnerable to intracellular ice formation. However, it could be used to control the dehydra‐ tion of the cell and provide limited interaction with high electrolyte concentrations through combination of cell pushing and cell engulfment. The numerical results also indicate that the interactions between cells are very important for cryopreservation conditions, where sus‐ pensions have higher volume fractions of cells. Level set method could deal with the interac‐ tion of clusters of cells with ice by assuming that the clusters act like a single particle with a radius that represents the overall size of the cell cluster. Increasing the radius of the repre‐ sentative particle would increase the probability of engulfment and pushing by the inter‐ face. Although the intermolecular forces between cells and the fluid flow around various arrays of single cells and cell clusters have not been included in the current model, these phenomena could be easily included in the simulations in order to study their effects on the partitioning of cells by the ice interface.

### **5. Summary**

Under the electrostatic field, the growth of dendritic ice crystals represents asymmetry be‐ havior. Fig 2 shows that the main branches parallel to the electrostatic field grow faster than the other branches [32]. One can also find that the growth of sidebranching is strengthened along the direction of electric field and weakened along the inverse direction. The above in‐ fluence becomes clearer by increasing the value of *λ* (comparing Fig. 2(b) with (c)), which indicates that a more strong magnitude of electrostatic field would strengthen the asymmet‐ ric growth of ice. Fig. 2(d) shows how the reference direction affects the ice growth under electric field. The above results concerning the effects of electric field are in good accordance

It is noteworthy that the biological media is saline solution including various ions, such

pletely different movement and gather near the corresponding electrodes, which induce non-uniform distribution of ions. In addition, the heterogeneous structure of ions, such as their different size, also leads to different ions' behaviors. These behaviors induced by the electrostatic field will affect the dynamical behavior of ice interface and strengthen the asymmetric growth of ice [4]. Thus, the effects of the ions under the electrostatic field are also indicated via modifying the surface tension. The effects of ions coupling with the electrostatic field on the dynamical behavior of ice interface should be incorpo‐ rated to a general model. The present investigation focuses on the effects of the electro‐ static field on the ice interface but is not to completely characterize the effects of ions,

The interaction between the cell and ice during freezing is very important to investigate the cell cryoinjury during freezing. However, it is still a rather difficult task to track the evolu‐ tion of the solidification front and characterize the interaction between the cell membrane and ice in detail. Recently, the progress [35, 36] gave some promising research directions to

Mao et al. [35] developed a sharp interface method to simulate the response of a biologi‐ cal cell during freezing. The cell is modeled as an aqueous salt solution surrounded by a semi-permeable membrane. The concentration and temperature fields both inside and outside a single cell are computed taking into account heat transfer, mass diffusion, membrane transport, and evolution of the solidification front. The external ice front is computed for both stable and unstable growth modes. It is shown that for the particular geometry chosen in this study, the instabilities on the front and the diffusional transport have only modest effects on the cell response. For the cooling conditions, solute and cell property parameters used, the low Peclet regime applies. The computational results are therefore validated against the conventional membrane-limited transport (Mazur) model. Good agreement of the simulation results with the Mazur model are obtained for a wide range of cooling rates and membrane permeabilities. A spatially non-isothermal situation

. Under the electrostatic field, these positive and negative ions have com‐

with the qualitative analysis in Ref. [25] and experimental observation in Ref. [24].

as Na+ and Cl<sup>−</sup>

194 Advanced Topics on Crystal Growth

solve this issue.

which will be addressed in further work.

**4. The interaction between ice growths with cell**

The present chapter has presented an overview on the theoretical modeling of the micro‐ scale phase change of biological system subject to freezing with special focus on the ice crys‐ tal formation and growth. In addition, the heat and solute transport, and interaction between the cell and ice during freezing have also been discussed. Although considerable investigation focus on ice crystal formation and its growth, there still exists some outstand‐ ing questions [37-40] including ice structure and ice phase diagram. Here we discuss some important issues and challenges about ice formation in biological system in brief. Firstly, the effects of cell membrane on ice nucleation and growth inside a cell with restricted space dur‐ ing freezing have not been investigated deeply, which is in fact very important to clarify the puzzle about the ice propagation between cells. In addition, the improved understanding on the effects of cell–cell interactions suggested that the mechanism of ice propagation should be mediated to some extent, but not at all, by the gap junction, if uncontrolled factors are eliminated. Secondly, a detailed three dimensional model needs to be developed to charac‐ terize the complicated freezing process that coupling micro heat and mass transfer, and ice nucleation and growth, which involves the moving water-ice interface and deformed cell membrane. Thirdly, the advanced experimental measures need to be designed to determine the parameter in the physical model and confirm the model validity. All the efforts will war‐ rant a bright future for better understanding and practice of cryobiology.

[7] Mazur P: Kinetics of Water Loss from Cells at Subzero Temperatures and Likelihood

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

197

[8] Zhao G, Luo DW, Gao DY: Universal Model for Intracellular Ice Formation and Its

[9] Yang G, Zhang AL, Xu LX, He XM: Modeling the Cell-type Dependence of Diffusionlimited Intracellular Ice Nucleation and Growth During both Vitrification and Slow

[10] Yan JF, Liu J: Characterization of the Nanocryosurgical Freezing Process through Modifying Mazur's Model. *Journal of Applied Physics* 2008, 103: 084311- 084321.

[11] Karim OA, Haymet ADJ: The Ice Water Interface - a Molecular-Dynamics Simulation

[12] Karim OA, Haymet ADJ: The Ice Water Interface. *Chemical Physics Letters* 1987,

[13] Matsumoto M, Saito S, Ohmine I: Molecular Dynamics Simulation of the Ice Nuclea‐ tion and Growth Process Leading to Water Freezing. *Nature* 2002, 416:409-413.

[14] Handel R, Davidchack RL, Anwar J, Brukhno A: Direct Calculation of Solid-liquid In‐ terfacial Free Energy for Molecular Systems: TIP4P Ice-water Interface. *Physical Re‐*

[15] Smith EJ, Bryk T, Haymet ADJ: Free Energy of Solvation of Simple Ions: Moleculardynamics Study of Solvation of Cl- and Na+ in the Ice/water Interface. *Journal of*

[16] Vrbka L, Jungwirth P: Brine Rejection from Freezing Salt Solutions: A Molecular Dy‐

[17] Carignano MA, Shepson PB, Szleifer I: Ions at the Ice/vapor Interface. *Chemical Phys‐*

[18] Nada H, Furukawa Y: Growth Inhibition at the Ice Prismatic Plane Induced by a Spruce Budworm Antifreeze Protein: A Molecular Dynamics Simulation Study. *Phys‐*

[19] Wathen B, Kuiper M, Walker V, Jia ZC: A New Model for Simulating 3-D Crystal Growth and Its Application to the Study of Antifreeze Proteins. *Journal of the Ameri‐*

[20] Yang C, Sharp KA: The Mechanism of the Type III Antifreeze Protein Action: A

[21] Nada H, Furukawa Y: Growth Inhibition Mechanism of An Ice-water Interface by A Mutant of Winter Flounder Antifreeze Protein: A Molecular Dynamics Study. *Journal*

of Intracellular Freezing. *Journal of General Physiology* 1963, 47:347-369.

Freezing. *Journal of Applied Physics* 2009, 105:114701-114711.

Study. *Journal of Chemical Physics* 1988, 89:6889-6896.

*view Letters* 2008, 100:036104-036107.

*ics Letters* 2007, 436:99-103.

*Chemical Physics* 2005, 123:03470601-03470616.

namics Study. *Physical Review Letters* 2005, 95:148501-148504.

Computational Study. *Biophysical Chemistry* 2004, 109:137-148.

*ical Chemistry Chemical Physics* 2011, 13:19936-19942.

*can Chemical Society* 2003, 125:729-737.

*of Physical Chemistry B* 2008, 112:7111-7119.

138:531-534.

Growth. *AICHE Journal* 2006, 52:2596-2606.

### **Acknowledgements**

This work is supported by the NSFC under Grant No. 51006114 & 81071255, the Specialized Research Fund for the Doctoral Program of Higher Education, and Research Fund from Tsinghua University under Grant 523003001.

### **Author details**

Zhi Zhu He1 and Jing Liu1,2

\*Address all correspondence to: jliu@mail.ipc.ac.cn

1 Beijing Key Laboratory of Cryo-Biomedical Engineering, and Key Laboratory of Cryogen‐ ics, Technical Institute of Physics and Chemistry, Chinese Academy of Sciences, Beijing, China

2 Department of Biomedical Engineering, School of Medicine, Tsinghua University, Beijing, China

### **References**


[7] Mazur P: Kinetics of Water Loss from Cells at Subzero Temperatures and Likelihood of Intracellular Freezing. *Journal of General Physiology* 1963, 47:347-369.

the parameter in the physical model and confirm the model validity. All the efforts will war‐

This work is supported by the NSFC under Grant No. 51006114 & 81071255, the Specialized Research Fund for the Doctoral Program of Higher Education, and Research Fund from

1 Beijing Key Laboratory of Cryo-Biomedical Engineering, and Key Laboratory of Cryogen‐ ics, Technical Institute of Physics and Chemistry, Chinese Academy of Sciences, Beijing, China

2 Department of Biomedical Engineering, School of Medicine, Tsinghua University, Beijing,

[1] Bischof JC: Quantitative Measurement and Prediction of Biophysical Response Dur‐ ing Freezing in Tissues. *Annual Review of Biomedical Engineering* 2000, 2:257-288.

[2] Mazur P, Leibo SP, Chu EH: A Two-factor Hypothesis of Freezing Injury. Evidence

[3] Mazur P, Seki S, Pinn IL, Kleinhans FW, Edashige K: Extra- and Intracellular Ice For‐

[4] Li FF, Liu J: Characterization of Micro-/Nano-Scale Ice Crystal Formation in Cryo-Bi‐ omedical Engineering: A Review. *Journal of Computational and Theoretical Nanoscience*

[5] Toner M, Cravalho EG, Karel M: Thermodynamics and Kinetics of Intracellular Ice Formation during Freezing of Biological Cells. *Journal of Applied Physics* 1990,

[6] Karlsson JOM, Cravalho EG, Toner M: A Model of Diffusion-Limited Ice Growth in‐ side Biological Cells during Freezing. *Journal of Applied Physics* 1994, 75:4442-4445.

from Chinese Hamster Tissue-culture Cells. *Exp Cell Res* 1972, 71:345-355.

mation in Mouse Oocytes. *Cryobiology* 2005, 51:29-53.

rant a bright future for better understanding and practice of cryobiology.

**Acknowledgements**

196 Advanced Topics on Crystal Growth

**Author details**

Zhi Zhu He1

China

**References**

2010, 7:85-96.

67:1582-1593.

Tsinghua University under Grant 523003001.

and Jing Liu1,2

\*Address all correspondence to: jliu@mail.ipc.ac.cn


[22] Chen C, Li WZ, Song YC, Yang J: Molecular Dynamics Simulation Studies of Cryo‐ protective Agent Solutions: The Relation Between Melting Temperature and the Ra‐ tio of Hydrogen Bonding Acceptor to Donor Number. *Molecular Physics* 2009, 107:673-684.

[36] Chang A, Dantzig JA, Darr BT, Hubel A: Modeling the Interaction of Biological Cells with a Solidifying Interface. *Journal of Computational Physics* 2007, 226:1808-1829.

Modeling Ice Crystal Formation of Water in Biological System

http://dx.doi.org/10.5772/54098

199

[37] Bartels-Rausch T, Bergeron V, Cartwright JHE, Escribano R, Finney JL, Grothe H, Gutierrez PJ, Haapala J, Kuhs WF, Pettersson JBC, et al: Ice Structures, Patterns, and Processes: A View Across the Icefields. *Reviews of Modern Physics* 2012, 84:885-944.

[38] Libbrecht KG: The Physics of Snow Crystals. *Reports on Progress in Physics* 2005,

[39] Blackford JR: Sintering and Microstructure of Ice: A Review. *Journal of Physics D-Ap‐*

[40] Petzold G, Aguilera JM: Ice Morphology: Fundamentals and Technological Applica‐

68:855-895.

*plied Physics* 2007, 40:355-385.

tions in Foods. *Food Biophysics* 2009, 4:378-396.


[36] Chang A, Dantzig JA, Darr BT, Hubel A: Modeling the Interaction of Biological Cells with a Solidifying Interface. *Journal of Computational Physics* 2007, 226:1808-1829.

[22] Chen C, Li WZ, Song YC, Yang J: Molecular Dynamics Simulation Studies of Cryo‐ protective Agent Solutions: The Relation Between Melting Temperature and the Ra‐ tio of Hydrogen Bonding Acceptor to Donor Number. *Molecular Physics* 2009,

[23] Chen C, Li WZ, Song YC, Yang J: A Molecular Dynamics Study of Cryoprotective Agent - Water-sodium Chloride Ternary Solutions. *Journal of Molecular Structure-The‐*

[24] Tao LR: Microscopic Experimental Study of the Freezing Process in Aqueous Solu‐ tion. Ph.D. Thesis, University of Shanghai for Science and Technology, 2000.

[25] Sun W, Xu XB, Zhang H, Xu CX: Effects of Dipole Polarization of Water Molecules

[26] Petersen A, Schneider H, Rau G, Glasmacher B: A New Approach for Freezing of Aqueous Solutions under Active Control of the Nucleation Temperature. *Cryobiology*

[27] Svishchev IM, Kusalik PG: Crystallization of Liquid Water in a Molecular-Dynamics

[28] Braslavsky I, Lipson SG: Electrofreezing Effect and Nucleation of Ice Crystals in Free

[29] Boettinger WJ, Warren JA, Beckermann C, Karma A: Phase-field Simulation of Solidi‐

[30] Li FF, Liu J, Yue K: Numerical Simulation on Ice Crystal Formulation in Cellular Lev‐ el Based on Phase Field Theory. *Chinese Journal of Low Temperature Physics* 2008,

[31] Xu Y, Mcdonough JM, Tagavi KA, Gao DY: Two-Dimensional Phase-Field Model Applied to Freezing into Supercooled Melt. *Cell Preservation Technology* 2004,

[32] He ZZ, Liu J: Characterizing Ice Crystal Growth Behavior Under Electric Field Using Phase Field Method. *ASME Journal of Biomechanical Engineering* 2009, 131:

[33] Wang SL, Sekerka RF, Wheeler AA, Murray BT, Coriell SR, Braun RJ, Mcfadden GB: Thermodynamically-Consistent Phase-Field Models for Solidification. *Physica D*

[34] Karma A, Rappel WJ: Phase-field Model of Dendritic Sidebranching with Thermal

[35] Mao L, Udaykumar HS, Karlsson JOM: Simulation of Micro-scale Interaction be‐ tween Ice and Biological Cells. *International Journal of Heat and Mass Transfer* 2003,

on Ice Formation under An Electrostatic Field. *Cryobiology* 2008, 56:93-99.

Simulation. *Physical Review Letters* 1994, 73:975-978.

Growth Experiments. *Applied Physics Letters* 1998, 72:264-266.

fication. *Annual Review of Materials Research* 2002, 32:163-194.

107:673-684.

198 Advanced Topics on Crystal Growth

*ochem* 2009, 916:37-46.

2006, 53:248-257.

30:84-92.

2:113-124.

0745021-0745023.

1993, 69:189-200.

46:5123-5136.

Noise. *Physical Review E* 1999, 60:3614-3625.


**Chapter 8**

**Crystallization: From the Conformer to the Crystal**

Crystallization, commonly defined as a process of formation of a crystalline solid from a supersaturated solution, melt or vapor phase, is an old technique widely used in laboratory and in industrial processes to separate and purify substances. In various modern industries, crystalline forms with a certain habit, size and structure, constitute the basic materials for the production of highly sophisticated materials [1, 2]. Integrat‐ ed circuits as well as piezoelectric and optical materials are just a few examples of de‐ vices whose properties are dependent on the crystal structure. Also, in organic chemistry, molecular crystals with determined characteristics are now-a-days of utmost importance for the production of pharmaceuticals, dyestuffs, pigments, foodstuffs, chemicals, cosmetics, etc. For all these reasons, crystal growth has become an impor‐

Crystallization from solution is one of the preferred methods to obtain a crystal since it can be carried out under different experimental conditions and provides a wide variety of prod‐ ucts. In fact, it can occur by lowering the temperature of a supersaturated solution, partial evaporation of the solvent, precipitation by adding an anti-solvent or vapor diffusion of a gas into the solution. Furthermore, solvents with different properties can be used. Such a di‐ versity of experimental conditions has great influence in the output, *i.e.*, in the type of crys‐ talline phase that can be obtained. Therefore, any aprioristic selection of the conditions leading to a desired final product is a fundamental challenge but it is also an almost unfeasi‐ ble task [3-5]. On respect to the harvest, crystallization is, in a certain extension, a trial-and-

The various steps of crystallization from solvents are summarized in the following scheme:

and reproduction in any medium, provided the original work is properly cited.

© 2013 Redinha et al.; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

© 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution,

J.S. Redinha, A.J. Lopes Jesus, A.A.C.C. Pais and

Additional information is available at the end of the chapter

J. A. S. Almeida

**1. Introduction**

error operation.

http://dx.doi.org/10.5772/54447

tant and attractive research field.

## **Crystallization: From the Conformer to the Crystal**

J.S. Redinha, A.J. Lopes Jesus, A.A.C.C. Pais and J. A. S. Almeida

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/54447

### **1. Introduction**

Crystallization, commonly defined as a process of formation of a crystalline solid from a supersaturated solution, melt or vapor phase, is an old technique widely used in laboratory and in industrial processes to separate and purify substances. In various modern industries, crystalline forms with a certain habit, size and structure, constitute the basic materials for the production of highly sophisticated materials [1, 2]. Integrat‐ ed circuits as well as piezoelectric and optical materials are just a few examples of de‐ vices whose properties are dependent on the crystal structure. Also, in organic chemistry, molecular crystals with determined characteristics are now-a-days of utmost importance for the production of pharmaceuticals, dyestuffs, pigments, foodstuffs, chemicals, cosmetics, etc. For all these reasons, crystal growth has become an impor‐ tant and attractive research field.

Crystallization from solution is one of the preferred methods to obtain a crystal since it can be carried out under different experimental conditions and provides a wide variety of prod‐ ucts. In fact, it can occur by lowering the temperature of a supersaturated solution, partial evaporation of the solvent, precipitation by adding an anti-solvent or vapor diffusion of a gas into the solution. Furthermore, solvents with different properties can be used. Such a di‐ versity of experimental conditions has great influence in the output, *i.e.*, in the type of crys‐ talline phase that can be obtained. Therefore, any aprioristic selection of the conditions leading to a desired final product is a fundamental challenge but it is also an almost unfeasi‐ ble task [3-5]. On respect to the harvest, crystallization is, in a certain extension, a trial-anderror operation.

The various steps of crystallization from solvents are summarized in the following scheme:

© 2013 Redinha et al.; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. © 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

As stated before, depending on the experimental conditions and intrinsic molecular features of the molecules, different crystalline phases or polymorphs can be obtained by crystalliza‐ tion. The polymorphs differ from one another in the arrangement (packing) and/or confor‐ mation of the molecules in the crystal lattice [26-29]. They are commonly called packing polymorphs in the first case and conformational polymorphs in the second one. Polymor‐ phism has significant implications in a wide range of areas and for this reason it has been matter of extensive investigation. Some typical examples are here reported with the objec‐

Fluctuations of the order of a solution, for example density, are accompanied by the forma‐ tion of solute molecular clusters. The driving force to bring the molecules close together is

> \* ln *<sup>C</sup> kT C*

In this equation, *μ* and *μ*\* are the chemical potentials of the solute at the actual (*C*) and satu‐

<sup>3</sup> 4/3 <sup>2</sup> <sup>4</sup> *<sup>r</sup> g r v* pD= D+

m ps

where *v* is the molecular volume and *σ* the interfacial energy per surface. The first term of the right hand equation gives the work produced by the intermolecular attractive forces (Δ*g*b) while the second one corresponds to the work required to increase the surface area in

With the increase of concentration, the amplitude of the fluctuations also increase and so do the size of the clusters. For smaller size clusters, Δ*g*s is the dominant contribution and Δ*g*/d*r* > 0, *i.e.*, the clusters are unstable. Conversely, for larger size clusters, Δ*g*<sup>b</sup> is dominant and Δ*g*/d*r* < 0, meaning that the clusters are stable solid particles. The curve representing the var‐ iation of Δ*g* with the radius, shown in Figure 1, passes through a maximum at *r* = *r*c (dΔ*g*/d*r* = 0), with *r*c representing the limiting size for stable clusters. A cluster of this size is the em‐ bryo of the solid form nucleus. The following relationships involving *r*c can be established:

(Δ*g*s). Owing to the meaning of both terms, equation (2) can be rewritten as

) concentrations, respectively. The Gibbs energy change per molecule (Δ*g*) corre‐

\* - = (1)

*b s* D =D +D *gg g* (3)

(2)

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

203

tive of pointing out the diversity of structures that can obtained by crystallization.

**2. Nucleation and crystal growth**

ration (*C*\*

1cm2

given by the following chemical potential difference:

m m

sponding to the formation of a spherical cluster of radius *r* is [30-32]

**Scheme 1.** Main crystallization pathways from solution

The competition between the solute-solvent and solute-solute interactions to form the nu‐ clei, and between the specific aggregation forces and packing in order to minimize repulsive interactions, determine the structure of the new crystalline phase. During this process the molecular conformation can change drastically in an unpredicted way. Various computa‐ tional approaches have been employed to predict the crystalline structure (unit cell, space group and atomic coordinates) formed by a given molecule [6-8], most of them relying on the assumption that the most probable forms are those with lowest lattice energy. Despite the recent success of these methods in the prediction of crystalline structures of small rigid molecules [9], their applicability to high size flexible molecules is not yet satisfactory [10, 11]. The main drawbacks lie in the difficultly of describing accurately the high complexity of inter and intramolecular interactions (e.g. hydrogen bonds and van der Waals interactions) and of selecting, among the energetically feasible polymorphs, that or those that really exist [3]. In truth, crystal prediction is still a very hard task so far.

Nevertheless, computational calculation can be a powerful tool to get information about the different steps of crystallization. The main goal of the present chapter is included in this re‐ search field. Quantum-chemical calculations or molecular dynamics simulations, adapted to the molecular complexity can be used for this purpose. In molecular terms, a compound with a certain conformational flexibility exists in a supersaturated aqueous solution as a mixture of conformers with different populations. The knowledge of the conformational equilibrium in this medium is a way of identifying the conformers that are likely to be involved in the early stages of nucleation and to interpret the mechanism of the crystallization process. In this chap‐ ter, these topics are addressed by presenting the results obtained for erythritol and *L*-threitol, two isomeric alditols with four carbon atoms, and glutamic acid. These compounds were chos‐ en because of the role played by their conformational features on crystallization. Furthermore, the two polyols exhibit a different crystallization behavior despite being chemically similar. The energies of the different conformers were obtained through theoretical calculations at the DFT level of theory [12, 13] using the popular B3LYP (Becke, 3-parameter, Lee-Yang-Parr) ex‐ change-correlation functional [14-16] combined with the aug-cc-pVDZ [17, 18] or 6-311+ +G(d,p) [19] basis set. Solvent effects have been accounted by using the well-established con‐ ductor-like polarizable continuum model (CPCM) [20-22].

As crystallization proceeds, the self-association of the solute molecules gives rise to molecu‐ lar aggregates which play an important role in the resulting crystalline structure [23-25]. Da‐ ta from the application of molecular dynamics simulations in the study of aqueous solutions of erythitol and *L*-threitol near the saturation concentration are reported in an attempt to give an insight into the type of supramolecular species formed in solution, and hence a step forward to elucidate about their crystallization behavior.

As stated before, depending on the experimental conditions and intrinsic molecular features of the molecules, different crystalline phases or polymorphs can be obtained by crystalliza‐ tion. The polymorphs differ from one another in the arrangement (packing) and/or confor‐ mation of the molecules in the crystal lattice [26-29]. They are commonly called packing polymorphs in the first case and conformational polymorphs in the second one. Polymor‐ phism has significant implications in a wide range of areas and for this reason it has been matter of extensive investigation. Some typical examples are here reported with the objec‐ tive of pointing out the diversity of structures that can obtained by crystallization.

### **2. Nucleation and crystal growth**

Molecule Discrete aggregates Nuclei Crystals

The competition between the solute-solvent and solute-solute interactions to form the nu‐ clei, and between the specific aggregation forces and packing in order to minimize repulsive interactions, determine the structure of the new crystalline phase. During this process the molecular conformation can change drastically in an unpredicted way. Various computa‐ tional approaches have been employed to predict the crystalline structure (unit cell, space group and atomic coordinates) formed by a given molecule [6-8], most of them relying on the assumption that the most probable forms are those with lowest lattice energy. Despite the recent success of these methods in the prediction of crystalline structures of small rigid molecules [9], their applicability to high size flexible molecules is not yet satisfactory [10, 11]. The main drawbacks lie in the difficultly of describing accurately the high complexity of inter and intramolecular interactions (e.g. hydrogen bonds and van der Waals interactions) and of selecting, among the energetically feasible polymorphs, that or those that really exist

Nevertheless, computational calculation can be a powerful tool to get information about the different steps of crystallization. The main goal of the present chapter is included in this re‐ search field. Quantum-chemical calculations or molecular dynamics simulations, adapted to the molecular complexity can be used for this purpose. In molecular terms, a compound with a certain conformational flexibility exists in a supersaturated aqueous solution as a mixture of conformers with different populations. The knowledge of the conformational equilibrium in this medium is a way of identifying the conformers that are likely to be involved in the early stages of nucleation and to interpret the mechanism of the crystallization process. In this chap‐ ter, these topics are addressed by presenting the results obtained for erythritol and *L*-threitol, two isomeric alditols with four carbon atoms, and glutamic acid. These compounds were chos‐ en because of the role played by their conformational features on crystallization. Furthermore, the two polyols exhibit a different crystallization behavior despite being chemically similar. The energies of the different conformers were obtained through theoretical calculations at the DFT level of theory [12, 13] using the popular B3LYP (Becke, 3-parameter, Lee-Yang-Parr) ex‐ change-correlation functional [14-16] combined with the aug-cc-pVDZ [17, 18] or 6-311+ +G(d,p) [19] basis set. Solvent effects have been accounted by using the well-established con‐

As crystallization proceeds, the self-association of the solute molecules gives rise to molecu‐ lar aggregates which play an important role in the resulting crystalline structure [23-25]. Da‐ ta from the application of molecular dynamics simulations in the study of aqueous solutions of erythitol and *L*-threitol near the saturation concentration are reported in an attempt to give an insight into the type of supramolecular species formed in solution, and hence a step

**Scheme 1.** Main crystallization pathways from solution

202 Advanced Topics on Crystal Growth

[3]. In truth, crystal prediction is still a very hard task so far.

ductor-like polarizable continuum model (CPCM) [20-22].

forward to elucidate about their crystallization behavior.

Nucleation Crystal growth

Fluctuations of the order of a solution, for example density, are accompanied by the forma‐ tion of solute molecular clusters. The driving force to bring the molecules close together is given by the following chemical potential difference:

$$
\mu - \mu^\* = kT \ln \frac{\mathbb{C}}{\mathbb{C}^\*} \tag{1}
$$

In this equation, *μ* and *μ*\* are the chemical potentials of the solute at the actual (*C*) and satu‐ ration (*C*\* ) concentrations, respectively. The Gibbs energy change per molecule (Δ*g*) corre‐ sponding to the formation of a spherical cluster of radius *r* is [30-32]

$$
\Delta \mathbf{g} = \frac{4 \times 3\pi r^3}{v} \Delta \mu + 4\pi r^2 \sigma \tag{2}
$$

where *v* is the molecular volume and *σ* the interfacial energy per surface. The first term of the right hand equation gives the work produced by the intermolecular attractive forces (Δ*g*b) while the second one corresponds to the work required to increase the surface area in 1cm2 (Δ*g*s). Owing to the meaning of both terms, equation (2) can be rewritten as

$$
\Delta \mathbf{g} = \Delta \mathbf{g}\_b + \Delta \mathbf{g}\_s \tag{3}
$$

With the increase of concentration, the amplitude of the fluctuations also increase and so do the size of the clusters. For smaller size clusters, Δ*g*s is the dominant contribution and Δ*g*/d*r* > 0, *i.e.*, the clusters are unstable. Conversely, for larger size clusters, Δ*g*<sup>b</sup> is dominant and Δ*g*/d*r* < 0, meaning that the clusters are stable solid particles. The curve representing the var‐ iation of Δ*g* with the radius, shown in Figure 1, passes through a maximum at *r* = *r*c (dΔ*g*/d*r* = 0), with *r*c representing the limiting size for stable clusters. A cluster of this size is the em‐ bryo of the solid form nucleus. The following relationships involving *r*c can be established:

$$\begin{aligned} r\_c &= -\frac{2v\sigma}{\Delta\mu} & \mathbf{a} \\ r\_c &= -\frac{2\sigma}{\Delta g\_b} & \mathbf{b} \\ \Delta g\_c &= \frac{4\pi r\_c^2 \sigma}{3} & \mathbf{c} \end{aligned} \tag{4}$$

tion between a group and its complementar in nature. For example, to a positive charge there would be a negative one, or to a H-bond donor a H-bond acceptor. The molecular packing bringing the molecules close together also gives rise to steric repulsion. Hence, the equilibrium structure is the result between the attraction force and the steric repulsion. This compromise

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

205

*<sup>r</sup>*<sup>c</sup> *r*

Δ*g*<sup>s</sup> + Δ*g*<sup>b</sup>

Let us consider the conformational variation from the gas phase to the solid for a simple molecule as paracetamol (acetaminophen). This is a drug of great commercial interest in the pharmaceutical industry owing to its wide use as antipyretic and analgesic agent. This mole‐ cule exhibits two relevant conformers differing from one another in the orientation of the OH relatively to the carbonyl group [35]. In one of them the dihedral angle formed by the H-O(1)-C-O(2) atoms is close to 180° (*trans* conformation, Figure 2), while in the other the angle is approximately 0° (*cis* conformation). Geometry optimization of the two conformers fol‐ lowed by vibrational frequencies calculation at the B3LYP/aug-cc-pVDZ level shows that in the gas phase the population (%) of *cis* relative to *trans* is 57:43. Including water effects by applying the CPCM model, both conformers are practically isoenergetic and so their relative population is 50:50. However, in the crystalline phase only the *trans* conformer is presented. This means that the lattice energy yielded by the *trans* conformer is more negative than that of the *cis* conformer. Along this chapter some examples will be given concerning the rela‐

makes the molecule to adopt a conformation corresponding to the lowest Gibbs energy.


D*g*

Δ*g*<sup>c</sup>

**Figure 1.** Variation of ∆g and its components with the radius of the molecular clusters.

tionship between the conformers' population in the three states of matter.

**Figure 2.** Optimized geometry of the *trans* conformer of paracetamol at the B3LYP/aug-cc-pVDZ level.

where Δg<sup>c</sup> corresponds to the value of Δg at *r*c. The supersaturation concentration (*S*) can be related with the clusters size by the Gibbs-Thomson equation:

$$\ln\frac{C}{C} = \ln S = \frac{2v\sigma}{kTr\_c} \tag{5}$$

Combining (4c) and (5), the following expression is obtained for Δ*g*c as function of the super‐ saturation:

$$
\Delta \mathbf{g}\_c = \frac{16\pi \sigma^3 v^2}{3kT \ln \mathbf{S}^2} \tag{6}
$$

Δ*g*<sup>c</sup> represents the energy barrier for the nucleation process and the rate of nucleation (J) will be given by:

$$J = Ae^{\left(\frac{-\Delta \mathbf{g}\_c}{kT}\right)}\tag{7}$$

This equation shows that the rate of nucleation is strongly dependent on the supersaturation and tends to a finite value for a given supersaturation.

Organic crystals are supramolecules, *i.e*, an ensemble of molecules, bonded by non-covalent forces involving ions, polar or non-polar groups. Such interactions, whatever their strength, do not change the molecular individuality. Hydrogen bonding occupies a special place among the intermolecular interactions occurring in molecular crystals. Its energy ranges from 160 kJ mol-1 (strongly covalent) down to 16 kJ mol-1 (mostly electrostatic) [33] and it is highly dependent on the orientation of the hydrogen donor group (D–H) relative to the ac‐ ceptor atom (A). The maximum strength occurs when the D–H‧‧‧A angle is approximately 180°, though values greater than 110° are acceptable for a hydrogen bond [34].

It is obvious that a crystal packing in which the specific interactions between the functional groups (molecular recognition) are more favorable is, in principle, that leading to the lowest energy. Thus, the lowest energy crystalline structure would correspond to the stronger interac‐ tion between a group and its complementar in nature. For example, to a positive charge there would be a negative one, or to a H-bond donor a H-bond acceptor. The molecular packing bringing the molecules close together also gives rise to steric repulsion. Hence, the equilibrium structure is the result between the attraction force and the steric repulsion. This compromise makes the molecule to adopt a conformation corresponding to the lowest Gibbs energy.

**Figure 1.** Variation of ∆g and its components with the radius of the molecular clusters.

c

*r*

2

= - <sup>D</sup>

= - <sup>D</sup>

*v*

s

m

4

c

*r*

related with the clusters size by the Gibbs-Thomson equation:

and tends to a finite value for a given supersaturation.

saturation:

204 Advanced Topics on Crystal Growth

be given by:

c

D =

*g*

2

3

*b c*

*g r*

p s

s

<sup>2</sup> <sup>b</sup>

where Δg<sup>c</sup> corresponds to the value of Δg at *r*c. The supersaturation concentration (*S*) can be

\* <sup>c</sup> <sup>2</sup> ln ln *C v <sup>S</sup> C kTr*

Combining (4c) and (5), the following expression is obtained for Δ*g*c as function of the super‐

3 2

*v*

Δ*g*<sup>c</sup> represents the energy barrier for the nucleation process and the rate of nucleation (J) will

*c g*

This equation shows that the rate of nucleation is strongly dependent on the supersaturation

Organic crystals are supramolecules, *i.e*, an ensemble of molecules, bonded by non-covalent forces involving ions, polar or non-polar groups. Such interactions, whatever their strength, do not change the molecular individuality. Hydrogen bonding occupies a special place among the intermolecular interactions occurring in molecular crystals. Its energy ranges from 160 kJ mol-1 (strongly covalent) down to 16 kJ mol-1 (mostly electrostatic) [33] and it is highly dependent on the orientation of the hydrogen donor group (D–H) relative to the ac‐ ceptor atom (A). The maximum strength occurs when the D–H‧‧‧A angle is approximately

It is obvious that a crystal packing in which the specific interactions between the functional groups (molecular recognition) are more favorable is, in principle, that leading to the lowest energy. Thus, the lowest energy crystalline structure would correspond to the stronger interac‐

æ -D ö ç ÷

c 2 16 3 ln

*<sup>g</sup> kT S* ps

*kT J Ae*

180°, though values greater than 110° are acceptable for a hydrogen bond [34].

a

(4)

c

s

= = (5)

D = (6)

è ø <sup>=</sup> (7)

Let us consider the conformational variation from the gas phase to the solid for a simple molecule as paracetamol (acetaminophen). This is a drug of great commercial interest in the pharmaceutical industry owing to its wide use as antipyretic and analgesic agent. This mole‐ cule exhibits two relevant conformers differing from one another in the orientation of the OH relatively to the carbonyl group [35]. In one of them the dihedral angle formed by the H-O(1)-C-O(2) atoms is close to 180° (*trans* conformation, Figure 2), while in the other the angle is approximately 0° (*cis* conformation). Geometry optimization of the two conformers fol‐ lowed by vibrational frequencies calculation at the B3LYP/aug-cc-pVDZ level shows that in the gas phase the population (%) of *cis* relative to *trans* is 57:43. Including water effects by applying the CPCM model, both conformers are practically isoenergetic and so their relative population is 50:50. However, in the crystalline phase only the *trans* conformer is presented. This means that the lattice energy yielded by the *trans* conformer is more negative than that of the *cis* conformer. Along this chapter some examples will be given concerning the rela‐ tionship between the conformers' population in the three states of matter.

**Figure 2.** Optimized geometry of the *trans* conformer of paracetamol at the B3LYP/aug-cc-pVDZ level.

### **3. Polymorphism**

According to McCrone [36], polymorphism of a compound corresponds to its ability to crys‐ tallize into more than one crystalline structure. From the thermodynamic point of view, only the lowest Gibbs energy form exists. However, higher energy polymorphs can remain as metastable forms for a period of time long enough to be used for practical purposes, provid‐ ing the height of the energy barrier separating these polymorphs and the most stable one is sufficiently high [28]. Since different solid-state modifications exhibit distinct physicochemi‐ cal properties, such for example the melting point, solubility, dissolution rate and density, polymorphism has a great impact in the pharmaceutical, food and dyes technologies, among others. For example, in the pharmaceutical industry, the desired polymorph for a given ac‐ tive pharmaceutical ingredient (API) would be the one with highest bioavailability and structural stability during shelf life [37]. A typical example of a drug exhibiting polymorphism is paracetamol. Three polymorphs, labeled as I, II and III, have been identified for this compound [38-40]. The first (form I) is the thermodynamically stable form while the second (form II) is metastable but exists long enough to be experimentally studied. It is usually prepared from cooled melt [41]. Form III is very unstable [41] and thus it is not possible to investigate its structure and properties. Form I crystallizes in the monoclinic system (space group *P21/a*) [42, 43]. As shown in Figure 3, the molecules in the crystalline structure are linked through O–H···O=C and N– H···O–H hydrogen bonds forming chains which in turn give rise to pleated layers [frames (a) and (a')]. Regarding form II, it crystallizes in the orthorhombic system (space group *Pcab*) [44, 45]. The H-bonding system is identical to that existing in form I but

the layers are plane interconnected by van der waals forces [frames (b) and (b')].

Figure 3. Views of the crystalline structures of forms I (a) and II (b) of paracetamol and detail of the H-bonds formed in both polymorphs [(a') and (b')]. **Figure 3.** Views of the crystalline structures of forms I (a) and II (b) of paracetamol and detail of the H-bonds formed in both polymorphs [(c) and (d)].

The geometric parameters of the H-bonds in forms I and II are summarized in Table 1. In both forms the O–H···O=C hydrogen bond is stronger than the N–H···O–H one. In addition, the hydrogen bonds in the stable polymorph are stronger than in the metastable one. The two polymorphs were also characterized by infrared spectroscopy [46, 47]. The values of the stretching vibration frequencies for the groups involved in hydrogen bonds are given in Table 2. To work as reference, the frequencies of the "free" groups, obtained from a spectrum of paracetamol in a CDCl3 solution were also included in the Table. From the frequency shift (Δν) of the NH and OH groups one can estimate the enthalpy involved in the two hydrogen bonds by applying the empirical relationship proposed by Iogansen: Δ*H*2 = 1.92(Δν – 40) [48]. The value of Δ*H* estimated for the N–H···O–H H-bond was found to be 12 kJ mol-1 A typical example of a drug exhibiting polymorphism is paracetamol. Three polymorphs, la‐ beled as I, II and III, have been identified for this compound [38-40]. The first (form I) is the thermodynamically stable form while the second (form II) is metastable but exists long enough to be experimentally studied. It is usually prepared from cooled melt [41]. Form III is very unstable [41] and thus it is not possible to investigate its structure and properties. Form I crystallizes in the monoclinic system (space group *P21/a*) [42, 43]. As shown in Figure 3, the molecules in the crystalline structure are linked through O–H O=C and N–H O–H hy‐ drogen bonds forming chains which in turn give rise to pleated layers [frames (a) and (a')].

Table 1. Distances and angles of the hydrogen bonds in the paracetamol polymorphs.*<sup>a</sup>*

Form

*<sup>a</sup>* Values taken from ref. [43, 45]

and in a CDCl3 solution.*<sup>a</sup>*

*<sup>a</sup>* Values taken from ref. [46, 47]

for both polymorphs, while that estimated for the O–H···O=C one was found to be 28 kJ mol-1 for form I and 19 kJ mol-1 for form II.

I 1.719 2.653 167.5 2.032 2.904 161.6 II 1.838 2.709 170.6 2.069 2.943 163.8

form I 3324 3160 1656 form II 3326 3205 1666 ; 1656 CDCl3 solution 3438 3600 1683 Table 2. Vibrational frequencies (ν/cm-1) corresponding to the NH, OH and CO stretching vibrations of paracetamol in the two crystalline phases

**O–H···O=C N–H···O–H** H···O/Å O···O/Å O–H···O/º H···O/Å N···O/Å N–H···O/º

**νNH νOH νCO** 

Regarding form II, it crystallizes in the orthorhombic system (space group *Pcab*) [44, 45]. The H-bonding system is identical to that existing in form I but the layers are plane interconnect‐

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

207

The geometric parameters of the H-bonds in forms I and II are summarized in Table 1. In both forms the O–H‧‧‧O=C hydrogen bond is stronger than the N–H‧‧‧O–H one. In addition,

The two polymorphs were also characterized by infrared spectroscopy [46, 47]. The values of the stretching vibration frequencies for the groups involved in hydrogen bonds are given in Table 2. To work as reference, the frequencies of the "free" groups, obtained from a spec‐ trum of paracetamol in a CDCl3 solution were also included in the Table. From the frequen‐ cy shift (Δν) of the NH and OH groups one can estimate the enthalpy involved in the two hydrogen bonds by applying the empirical relationship proposed by Iogansen: Δ*H*<sup>2</sup> = 1.92(Δν – 40) [48]. The value of Δ*H* estimated for the N–H‧‧‧O–H H-bond was found to be 12 kJ mol-1 for both polymorphs, while that estimated for the O–H‧‧‧O=C one was found to be

> **O–H···O=C N–H···O–H H···O/Å O···O/Å O–H···O/º H···O/Å N···O/Å N–H···O/º**

> > **νNH νOH νCO**

I 1.719 2.653 167.5 2.032 2.904 161.6 II 1.838 2.709 170.6 2.069 2.943 163.8

form I 3324 3160 1656 form II 3326 3205 1666 ; 1656 CDCl3 solution 3438 3600 1683

**Table 2.** Vibrational frequencies (ν/cm-1) corresponding to the NH, OH and CO stretching vibrations of paracetamol in

The higher stability of form I is also confirmed by the thermodynamic properties obtained for sublimation, fusion and vaporization of the two polymorphic forms [38]. The value of the packing density, determined by flotation or X-ray (or neutron) diffraction [49] is higher in form II than in form I. In fact, at 298K the density of the former is 1.336 g cm-3 while that of the latter is 1.297 g cm-3. That is, the higher stability of form I results from a stronger hy‐ drogen bond interaction rather than from a more favorable packing. In this compound, the differences between specific interactions resulting from the conformational recognition over‐

**Table 1.** Distances and angles of the hydrogen bonds in the paracetamol polymorphs.*<sup>a</sup>*

the hydrogen bonds in the stable polymorph are stronger than in the metastable one.

ed by van der waals forces [frames (b) and (b')].

28 kJ mol-1 for form I and 19 kJ mol-1 for form II.

**Form**

Values taken from ref. [43, 45]

Values taken from ref. [46, 47]

the two crystalline phases and in a CDCl3 solution.*<sup>a</sup>*

come that due to the crystal packing.

*a*

*a*

Regarding form II, it crystallizes in the orthorhombic system (space group *Pcab*) [44, 45]. The H-bonding system is identical to that existing in form I but the layers are plane interconnect‐ ed by van der waals forces [frames (b) and (b')].

The geometric parameters of the H-bonds in forms I and II are summarized in Table 1. In both forms the O–H‧‧‧O=C hydrogen bond is stronger than the N–H‧‧‧O–H one. In addition, the hydrogen bonds in the stable polymorph are stronger than in the metastable one.

A typical example of a drug exhibiting polymorphism is paracetamol. Three polymorphs, labeled as I, II and III, have been identified for this compound [38-40]. The first (form I) is the thermodynamically stable form while the second (form II) is metastable but exists long enough to be experimentally studied. It is usually prepared from cooled melt [41]. Form III is very unstable [41] and thus it is not possible to investigate its structure and properties. Form I crystallizes in the monoclinic system (space group *P21/a*) [42, 43]. As shown in Figure 3, the molecules in the crystalline structure are linked through O–H···O=C and N– H···O–H hydrogen bonds forming chains which in turn give rise to pleated layers [frames (a) and (a')]. Regarding form II, it crystallizes in the orthorhombic system (space group *Pcab*) [44, 45]. The H-bonding system is identical to that existing in form I but The two polymorphs were also characterized by infrared spectroscopy [46, 47]. The values of the stretching vibration frequencies for the groups involved in hydrogen bonds are given in Table 2. To work as reference, the frequencies of the "free" groups, obtained from a spec‐ trum of paracetamol in a CDCl3 solution were also included in the Table. From the frequen‐ cy shift (Δν) of the NH and OH groups one can estimate the enthalpy involved in the two hydrogen bonds by applying the empirical relationship proposed by Iogansen: Δ*H*<sup>2</sup> = 1.92(Δν – 40) [48]. The value of Δ*H* estimated for the N–H‧‧‧O–H H-bond was found to be 12 kJ mol-1 for both polymorphs, while that estimated for the O–H‧‧‧O=C one was found to be 28 kJ mol-1 for form I and 19 kJ mol-1 for form II.


*a* Values taken from ref. [43, 45]

**3. Polymorphism**

206 Advanced Topics on Crystal Growth

structural stability during shelf life [37].

(b')].

both polymorphs [(c) and (d)].

metastable one.

Form

**Figure 3.** Views of the crystalline structures of forms I (a) and II (b) of paracetamol and detail of the H-bonds formed in

A typical example of a drug exhibiting polymorphism is paracetamol. Three polymorphs, la‐ beled as I, II and III, have been identified for this compound [38-40]. The first (form I) is the thermodynamically stable form while the second (form II) is metastable but exists long enough to be experimentally studied. It is usually prepared from cooled melt [41]. Form III is very unstable [41] and thus it is not possible to investigate its structure and properties. Form I crystallizes in the monoclinic system (space group *P21/a*) [42, 43]. As shown in Figure 3, the molecules in the crystalline structure are linked through O–H O=C and N–H O–H hy‐ drogen bonds forming chains which in turn give rise to pleated layers [frames (a) and (a')].

*<sup>a</sup>* Values taken from ref. [43, 45]

and in a CDCl3 solution.*<sup>a</sup>*

*<sup>a</sup>* Values taken from ref. [46, 47]

Table 1. Distances and angles of the hydrogen bonds in the paracetamol polymorphs.*<sup>a</sup>*

According to McCrone [36], polymorphism of a compound corresponds to its ability to crys‐ tallize into more than one crystalline structure. From the thermodynamic point of view, only the lowest Gibbs energy form exists. However, higher energy polymorphs can remain as metastable forms for a period of time long enough to be used for practical purposes, provid‐ ing the height of the energy barrier separating these polymorphs and the most stable one is sufficiently high [28]. Since different solid-state modifications exhibit distinct physicochemi‐ cal properties, such for example the melting point, solubility, dissolution rate and density, polymorphism has a great impact in the pharmaceutical, food and dyes technologies, among others. For example, in the pharmaceutical industry, the desired polymorph for a given ac‐ tive pharmaceutical ingredient (API) would be the one with highest bioavailability and

the layers are plane interconnected by van der waals forces [frames (b) and (b')].

for both polymorphs, while that estimated for the O–H···O=C one was found to be 28 kJ mol-1 for form I and 19 kJ mol-1 for form II.

I 1.719 2.653 167.5 2.032 2.904 161.6 II 1.838 2.709 170.6 2.069 2.943 163.8

form I 3324 3160 1656 form II 3326 3205 1666 ; 1656 CDCl3 solution 3438 3600 1683 Table 2. Vibrational frequencies (ν/cm-1) corresponding to the NH, OH and CO stretching vibrations of paracetamol in the two crystalline phases

**O–H···O=C N–H···O–H** H···O/Å O···O/Å O–H···O/º H···O/Å N···O/Å N–H···O/º

**νNH νOH νCO** 

**Table 1.** Distances and angles of the hydrogen bonds in the paracetamol polymorphs.*<sup>a</sup>*


Figure 3. Views of the crystalline structures of forms I (a) and II (b) of paracetamol and detail of the H-bonds formed in both polymorphs [(a') and **Table 2.** Vibrational frequencies (ν/cm-1) corresponding to the NH, OH and CO stretching vibrations of paracetamol in the two crystalline phases and in a CDCl3 solution.*<sup>a</sup>*

The geometric parameters of the H-bonds in forms I and II are summarized in Table 1. In both forms the O–H···O=C hydrogen bond is stronger than the N–H···O–H one. In addition, the hydrogen bonds in the stable polymorph are stronger than in the The two polymorphs were also characterized by infrared spectroscopy [46, 47]. The values of the stretching vibration frequencies for the groups involved in hydrogen bonds are given in Table 2. To work as reference, the frequencies of the "free" groups, obtained from a spectrum of paracetamol in a CDCl3 solution were also included in the Table. From the frequency shift (Δν) of the NH and OH groups one can estimate the enthalpy involved in the two hydrogen bonds by applying the empirical relationship proposed by Iogansen: Δ*H*2 = 1.92(Δν – 40) [48]. The value of Δ*H* estimated for the N–H···O–H H-bond was found to be 12 kJ mol-1 The higher stability of form I is also confirmed by the thermodynamic properties obtained for sublimation, fusion and vaporization of the two polymorphic forms [38]. The value of the packing density, determined by flotation or X-ray (or neutron) diffraction [49] is higher in form II than in form I. In fact, at 298K the density of the former is 1.336 g cm-3 while that of the latter is 1.297 g cm-3. That is, the higher stability of form I results from a stronger hy‐ drogen bond interaction rather than from a more favorable packing. In this compound, the differences between specific interactions resulting from the conformational recognition over‐ come that due to the crystal packing.

The polymorphism of paracetamol assumes great importance from the practical point of view in so far the commercialized polymorph is not the most suitable solid for formulation [45, 50, 51]. Despite it is easily obtained from various solvents, it has the inconvenient to re‐ quire a binding agent to make tablets for compression. During this process it gives a brittle solid, consequence of its rough molecular layers. Conversely, form II is constituted by thin plates that glide on pressing. Its plasticity allows the formulation into tablets by direct com‐ pression without the need of incorporating any binder agent [45]. The main difficulty to ob‐ tain this form comes up against the existence of an adequate method for their preparation at an industrial scale [41].

An interesting feature of aspirin single crystals is that they exhibit simultaneously domains of forms I and II of variable relative size. Deran *et al*. [58] explain this type of behavior as an accidental degeneracy. According to these authors the two polymorphs are isoenergetic and during crystallization an intergrowth of both forms occurs. The more favorable internal con‐ formation in form I is compensated by an enhanced cooperativity of the catameric H-bond‐ ing network in form II. These results raise the question of whether or not these two forms of aspirin can be called polymorphs. According to the definition given by McCrone [36], they

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

209

Besides polymorphism, further crystalline structure complexity arises with the existence of disorder in crystals, situation quite common in many organic crystals [60-64]. This occurs when portions of the structure under analysis occupy two or more positions. For example, molecular fragments of organic molecules not involved in intermolecular interactions can be free enough to acquire different conformations. Such solids are designated by the conflicting denomination of disordered crystals. An example of a compound exhibiting crystal disorder is betaxolol hydrochloride, 1-{4-[2-(cyclopropylmethoxy)ethyl]phenoxy}-3-isopropyl-ami‐ no-2-propanol hydrochloride (Figure 5). It is a cardio selective β adrenergic drug that crys‐

Common ordered head Disordered tail

**N O3**

**Figure 5.** Conformations of betaxolol hydrochloride in the crystalline structure. Hydrogen bonds between the amino and hydroxyl groups with the chloride ion of the same molecule are represented by dotted lines. For better visualiza‐

OH groups of the isopropylaminoethanol moiety form H-bonds with the chloride ion of the

also H-bonded to oxygen atom of the OH group of a second molecule while the other hydro‐ gen of the same group is linked to the chloride ion of a third molecule [63, 64]. That is, the betaxolol molecules are connected each other just through the isopropylaminoethanol frag‐ ment, leaving a long molecular chain free to move. The X-ray diffraction data shows the ex‐ istence of an ordered head from the isopropyl group up to O2 and a disordered tail from this atom to the end. By refinement, the crystalline structure can be satisfactory interpreted as being constituted by two conformations which are displayed in Figure 5. Since we have a no well defined unit cell, the structure of betaxolol hydrochloride can not be included in the

+ and the

+

group is

The intermolecular arrangement in the crystalline structure is as follows: the NH2

same molecule as shown in Figure 5. In addition, one of the hydrogens of the NH2

**O1 O2**

do not fit into the concept.

tion the hydrogen atoms are not shown.

polymorphism or isomorphism definitions.

tallizes into the triclinic system with *P*1 symmetry [65].

**Cl<sup>−</sup>**

Aspirin (acetylsalicylic acid) is a drug commonly used as analgesic, antipyretic and anti-in‐ flammatory. It is also known to act as anticoagulant by reducing the blood platelet aggrega‐ tion [52]. About 35,000 tons of aspirin are taken a day throughout the world.

Using X-ray diffraction, the structure of aspirin has been characterized for the first time in 1964 as a monoclinical crystal, space group P21/c [53]. In 2004, the existence of a second pol‐ ymorph (form II) was predicted computationally [54] and observed experimentally one year later [55]. However, in view of the slight differences between both structural modifications, it was speculated that this new form might result from uncertainties of the methods used in the structural characterization [56, 57]. Just recently it has been settled the existence of a sec‐ ond polymorph for this compound [58]. In both forms, the carboxyl groups of two neighbor molecules are connected by O–H‧‧‧O=C H-bonds giving rise to centrosymmetric dimers as illustrated in Figure 4. These dimers are arranged into stacked layers, which are identical in the two polymorphs. The acetyl groups of molecules belonging to different layers are in turns interconnected through C–H‧‧‧O interactions. In Form I, these interactions have cen‐ trosymmetric geometry, *i.e.*, the CH3 and CO groups of one molecule are bonded to the com‐ plementary groups of a neighbor molecule. The structure is a sequence of AAAA layers. Unlike, the structure of form II is a sequence of ABAB layers [59]. The difference between the B and A layers lies in a relative translation of the lattice repetition. This means that in the latter form the CH3 and CO groups of one molecule are connected to the complementary groups of two neighbor molecules, yielding catemers.

**Figure 4.** H-bonds in the crystalline structures of forms I and II of aspirin.

An interesting feature of aspirin single crystals is that they exhibit simultaneously domains of forms I and II of variable relative size. Deran *et al*. [58] explain this type of behavior as an accidental degeneracy. According to these authors the two polymorphs are isoenergetic and during crystallization an intergrowth of both forms occurs. The more favorable internal con‐ formation in form I is compensated by an enhanced cooperativity of the catameric H-bond‐ ing network in form II. These results raise the question of whether or not these two forms of aspirin can be called polymorphs. According to the definition given by McCrone [36], they do not fit into the concept.

The polymorphism of paracetamol assumes great importance from the practical point of view in so far the commercialized polymorph is not the most suitable solid for formulation [45, 50, 51]. Despite it is easily obtained from various solvents, it has the inconvenient to re‐ quire a binding agent to make tablets for compression. During this process it gives a brittle solid, consequence of its rough molecular layers. Conversely, form II is constituted by thin plates that glide on pressing. Its plasticity allows the formulation into tablets by direct com‐ pression without the need of incorporating any binder agent [45]. The main difficulty to ob‐ tain this form comes up against the existence of an adequate method for their preparation at

Aspirin (acetylsalicylic acid) is a drug commonly used as analgesic, antipyretic and anti-in‐ flammatory. It is also known to act as anticoagulant by reducing the blood platelet aggrega‐

Using X-ray diffraction, the structure of aspirin has been characterized for the first time in 1964 as a monoclinical crystal, space group P21/c [53]. In 2004, the existence of a second pol‐ ymorph (form II) was predicted computationally [54] and observed experimentally one year later [55]. However, in view of the slight differences between both structural modifications, it was speculated that this new form might result from uncertainties of the methods used in the structural characterization [56, 57]. Just recently it has been settled the existence of a sec‐ ond polymorph for this compound [58]. In both forms, the carboxyl groups of two neighbor molecules are connected by O–H‧‧‧O=C H-bonds giving rise to centrosymmetric dimers as illustrated in Figure 4. These dimers are arranged into stacked layers, which are identical in the two polymorphs. The acetyl groups of molecules belonging to different layers are in turns interconnected through C–H‧‧‧O interactions. In Form I, these interactions have cen‐ trosymmetric geometry, *i.e.*, the CH3 and CO groups of one molecule are bonded to the com‐ plementary groups of a neighbor molecule. The structure is a sequence of AAAA layers. Unlike, the structure of form II is a sequence of ABAB layers [59]. The difference between the B and A layers lies in a relative translation of the lattice repetition. This means that in the latter form the CH3 and CO groups of one molecule are connected to the complementary

(form I) (form II)

tion [52]. About 35,000 tons of aspirin are taken a day throughout the world.

groups of two neighbor molecules, yielding catemers.

**Figure 4.** H-bonds in the crystalline structures of forms I and II of aspirin.

an industrial scale [41].

208 Advanced Topics on Crystal Growth

Besides polymorphism, further crystalline structure complexity arises with the existence of disorder in crystals, situation quite common in many organic crystals [60-64]. This occurs when portions of the structure under analysis occupy two or more positions. For example, molecular fragments of organic molecules not involved in intermolecular interactions can be free enough to acquire different conformations. Such solids are designated by the conflicting denomination of disordered crystals. An example of a compound exhibiting crystal disorder is betaxolol hydrochloride, 1-{4-[2-(cyclopropylmethoxy)ethyl]phenoxy}-3-isopropyl-ami‐ no-2-propanol hydrochloride (Figure 5). It is a cardio selective β adrenergic drug that crys‐ tallizes into the triclinic system with *P*1 symmetry [65].

**Figure 5.** Conformations of betaxolol hydrochloride in the crystalline structure. Hydrogen bonds between the amino and hydroxyl groups with the chloride ion of the same molecule are represented by dotted lines. For better visualiza‐ tion the hydrogen atoms are not shown.

The intermolecular arrangement in the crystalline structure is as follows: the NH2 + and the OH groups of the isopropylaminoethanol moiety form H-bonds with the chloride ion of the same molecule as shown in Figure 5. In addition, one of the hydrogens of the NH2 + group is also H-bonded to oxygen atom of the OH group of a second molecule while the other hydro‐ gen of the same group is linked to the chloride ion of a third molecule [63, 64]. That is, the betaxolol molecules are connected each other just through the isopropylaminoethanol frag‐ ment, leaving a long molecular chain free to move. The X-ray diffraction data shows the ex‐ istence of an ordered head from the isopropyl group up to O2 and a disordered tail from this atom to the end. By refinement, the crystalline structure can be satisfactory interpreted as being constituted by two conformations which are displayed in Figure 5. Since we have a no well defined unit cell, the structure of betaxolol hydrochloride can not be included in the polymorphism or isomorphism definitions.

### **4. Erythritol and threitol: identical chemical structure, different crystalline assembling**

Erythritol and *L*-threitol (see Figure 8) are two diastereomers of 1,2,3,4-butanetetrol. The first is a *meso* compound while the second is optically active, existing as *D*- or *L*-threitol. An im‐ portant difference between both alditols lies in their crystallization ability. Whilst erythritol is easily crystallized from aqueous solution or melt, *L*-threitol is much more difficult to crys‐ tallize [66, 67]. As illustrated in the DSC curves shown in Figure 6, molten erythritol trans‐ forms readily into a crystalline phase providing it is cooled at scanning rates lower than or equal to 10°C/min. On the contrary, no crystallization peak is observed for *L*-threitol, even when the melt is cooled very slowly, say 2°C/min. Since the two compounds are chemically similar, their different crystallization tendency is likely to have a conformational origin [68].

The Boltzmann populations of the erythritol and *L*-threitol conformers in gas phase and aque‐ ous solution are displayed in Figure 7. The geometries of the relevant conformers in both media, as well as the conformations exiting in the crystalline structures, are depicted in Figure 8. For the sake of simplicity the conformers were grouped according to their backbone family, being la‐

), *g*' (gauche–

the O(1)-C(1)-C(2)-O(2), C(1)-C(2)-C(3)-C(4) and O(3)-C(3)-C(4)-O(4) dihedrals, respectively.

As shown in Figure 7, isolated erythritol exists preferentially in the *tgg* and *ggg* bent confor‐ mations (both make up 38% of the conformational mixture at 298.15K). Their stability results from the formation of intramolecular hydrogen bonds which are particularly strong in the *tgg* conformers due to the participation of the terminal OH groups. In the presence of the solvent, the population of these conformers decreases significantly (11%). Conversely, two higher energy conformers in gas phase (*g'tg* and *gtg*), both with a distended carbon back‐ bone, emerge as the most stable forms in solution (49%). This population change can be ex‐ plained by the stronger hydration of conformers with the OH groups less internally Hbonded and therefore more available to interact with water. In the crystal, the neutron diffraction data available for this compound indicates that the molecule exhibits two confor‐ mations, labeled as A and B, differing each other by the positions of the H[O(2)] and H[O(3)] hydrogen atoms. At 22.6 K the occupancy of conformations A and B are 85% and 15%, re‐ spectively [71]. In terms of backbone, both have a *g'tg* conformation which is the same as that exhibited by the preferred conformers of erythritol in aqueous solution. As shown in Figure 8, the geometry of the most stable conformer in solution matches that of the crystal

*gg'g' g'tg' g'tg ttt gtt ggg tg't gtg tgt*

**Figure 7.** Conformational distribution of erythritol and *L*-threitol in gas phase and aqueous solution. Boltzmann popu‐

lations were estimated from the Gibbs energies of the conformers [B3LYP/6-311++G(d,p)] in both media.

*L*-threitol

*g'tg gtg gg'g' ttg' ggg g'tt gtg' gg'g tgg ggt ttt g'gt g'gg*

Gas phase

Aqueous solution Erythritol

) and *t* (trans), to specify the orientation of

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

211

beled with three italic letters: *g* (gauche+

conformation having the highest occupancy percentage.

Boltzmann population (%)

**Figure 6.** DSC cooling curves of molten erythritol (m.p = 118°C) and *L*-threitol (m.p= 89°C) obtained at different cool‐ ing rates

Conformationally, both diastereomers are highly flexible. In fact, the existence of seven inde‐ pendent dihedral angles, each of them assuming three standard orientations: gauche+ (60°), gauche– (-60°) and trans (180°), originate a total of 2187 possible conformers. Such a large num‐ ber of structures make their systematic exploration by means of *ab initio* calculations unfeasi‐ ble. Therefore, a previous conformational search method using random generation and subsequent molecular mechanics energy minimization has been employed to sampling the most relevant conformations. All generated structures falling within an energy window of 20 kJ mol-1 above the global minimum have been then optimized at DFT level using the B3LYP functional and the 6-311++G(d,p) basis set. In order to simulate solvent effects, their energy was also evaluated by single-point calculations in aqueous solution using the CPCM continuum model. Details of the computational calculations are given elsewhere [69, 70].

The Boltzmann populations of the erythritol and *L*-threitol conformers in gas phase and aque‐ ous solution are displayed in Figure 7. The geometries of the relevant conformers in both media, as well as the conformations exiting in the crystalline structures, are depicted in Figure 8. For the sake of simplicity the conformers were grouped according to their backbone family, being la‐ beled with three italic letters: *g* (gauche+ ), *g*' (gauche– ) and *t* (trans), to specify the orientation of the O(1)-C(1)-C(2)-O(2), C(1)-C(2)-C(3)-C(4) and O(3)-C(3)-C(4)-O(4) dihedrals, respectively.

**4. Erythritol and threitol: identical chemical structure, different**

Erythritol and *L*-threitol (see Figure 8) are two diastereomers of 1,2,3,4-butanetetrol. The first is a *meso* compound while the second is optically active, existing as *D*- or *L*-threitol. An im‐ portant difference between both alditols lies in their crystallization ability. Whilst erythritol is easily crystallized from aqueous solution or melt, *L*-threitol is much more difficult to crys‐ tallize [66, 67]. As illustrated in the DSC curves shown in Figure 6, molten erythritol trans‐ forms readily into a crystalline phase providing it is cooled at scanning rates lower than or equal to 10°C/min. On the contrary, no crystallization peak is observed for *L*-threitol, even when the melt is cooled very slowly, say 2°C/min. Since the two compounds are chemically similar, their different crystallization tendency is likely to have a conformational origin [68].


Erythritol 2ºC/min

Temperature /ºC

**Figure 6.** DSC cooling curves of molten erythritol (m.p = 118°C) and *L*-threitol (m.p= 89°C) obtained at different cool‐

Conformationally, both diastereomers are highly flexible. In fact, the existence of seven inde‐ pendent dihedral angles, each of them assuming three standard orientations: gauche+ (60°),

ber of structures make their systematic exploration by means of *ab initio* calculations unfeasi‐ ble. Therefore, a previous conformational search method using random generation and subsequent molecular mechanics energy minimization has been employed to sampling the most relevant conformations. All generated structures falling within an energy window of 20 kJ mol-1 above the global minimum have been then optimized at DFT level using the B3LYP functional and the 6-311++G(d,p) basis set. In order to simulate solvent effects, their energy was also evaluated by single-point calculations in aqueous solution using the CPCM continuum

model. Details of the computational calculations are given elsewhere [69, 70].

(-60°) and trans (180°), originate a total of 2187 possible conformers. Such a large num‐

*L*-threitol

Endo

dQ/dt

ing rates

gauche–

10ºC/min

2ºC/min

10ºC/min

**crystalline assembling**

210 Advanced Topics on Crystal Growth

As shown in Figure 7, isolated erythritol exists preferentially in the *tgg* and *ggg* bent confor‐ mations (both make up 38% of the conformational mixture at 298.15K). Their stability results from the formation of intramolecular hydrogen bonds which are particularly strong in the *tgg* conformers due to the participation of the terminal OH groups. In the presence of the solvent, the population of these conformers decreases significantly (11%). Conversely, two higher energy conformers in gas phase (*g'tg* and *gtg*), both with a distended carbon back‐ bone, emerge as the most stable forms in solution (49%). This population change can be ex‐ plained by the stronger hydration of conformers with the OH groups less internally Hbonded and therefore more available to interact with water. In the crystal, the neutron diffraction data available for this compound indicates that the molecule exhibits two confor‐ mations, labeled as A and B, differing each other by the positions of the H[O(2)] and H[O(3)] hydrogen atoms. At 22.6 K the occupancy of conformations A and B are 85% and 15%, re‐ spectively [71]. In terms of backbone, both have a *g'tg* conformation which is the same as that exhibited by the preferred conformers of erythritol in aqueous solution. As shown in Figure 8, the geometry of the most stable conformer in solution matches that of the crystal conformation having the highest occupancy percentage.

**Figure 7.** Conformational distribution of erythritol and *L*-threitol in gas phase and aqueous solution. Boltzmann popu‐ lations were estimated from the Gibbs energies of the conformers [B3LYP/6-311++G(d,p)] in both media.

Regarding *L*-threitol, its gas phase conformational mixture is largely dominated by the *gtg* conformers which are characterized by a cyclic, symmetric and cooperative hydrogen bond‐ ing network. Since this OH groups' arrangement is very unfavorable to interact with water, they practically disappear in solution and are replaced by the *gg'g'*, *g'tg'* or *g'tg* conforma‐ tions which are not so highly intramolecularly bonded. Unlike erythritol, none of these con‐ formations has been identified in the crystalline *L*-threitol [72]. Here, the molecules adopt a distended conformation (*ttt*) which has a relatively low weight in solution (8%) and does not exist in the vacuum.

**5. Conformational variation during the molecular incorporation into the**

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

213

Glutamic acid, 2-aminopentanedioic acid (C5H9NO4), is another example illustrating the role played by conformation on crystallization. Two conformational polymorphs have been identified for this compound, labeled as α and β, with the latter being the thermodynamical‐ ly stable form over all temperature range [74-77]. The crystals have different morphologies: the metastable form exhibits a prismatic shape while the stable one has a needle-like shape. Both belong to the ortorrombic space group with four molecules in the unit cell (*Z*=4) [78, 79]. In the crystal lattice the molecules are in the zwitterionic state and all functional groups

The conformations adopted by the *L*-glutamic acid molecules in the two polymorphs are displayed in Figure 9. Their main difference lies in the orientation of the two C-C-C-C dihe‐ dral angles. In the crystal lattice the C(1)-C(2)-C(3)-C(4) and C(2)-C(3)-C(4)-C(5) dihedrals assume, respectively, values of 59.2°, 68.3° in α and -171.1°, -73.1° in β. Using the same dihe‐ dral labeling scheme previously adopted for the polyols, the conformation of glutamic acid in the α and β-crystals is *gg* and *tg'*, respectively. Depending on the experimental conditions, namely the temperature, both polymorphs can be crystallized from aqueous solution. Cool‐ ing a supersaturated solution to temperatures below 25°C originates almost pure α-crystals, whereas as the temperature raises the proportion of form β in the precipitated crystals in‐

The conformational behavior of zwitterionic glutamic acid is aqueous solution is a useful starting point to understand the crystallization of this compound in molecular terms. This has been done theoretically by performing full geometry optimizations using the CPCM continuum solvation model and the B3LYP/aug-cc-pVDZ model chemistry, both imple‐ mented in the Gaussian 03 program. The cavity was built with the Bondii radii which have been found to yield accurate results for the hydration of similar molecules in the zwitterion‐ ic state, such for example glycine [80]. Nine starting geometries were built by assuming the three standard orientations (*g*, *g*' and *t*) for each one of the C-C-C-C dihedrals. The remain‐ ing dihedrals were kept in their preferred orientations. The optimized structures were fur‐ ther submitted to a vibrational frequency calculation at the same level to ensure that they correspond to minima on the aqueous potential energy surface and also to calculate the ther‐

The Gibbs energies of the conformers in aqueous solution at 298.15 K (*G*sol) can be expressed as: *G*sol = *E*int + Δ*G*thermal + Δ*G*hyd, with *E*int representing the intrinsic (abbreviated as "int") electronic energy of the conformer at 0 K, *i.e.*, excluding solvent effects, Δ*G*thermal the thermal correction to the Gibbs energy from 0 to 298K and Δ*G*hyd the Gibbs energy of hydration at 298.15K. The val‐ ues of *E*int, *G*sol, Δ*G*hyd, as well as the equilibrium populations at 298.15 K, are given in Table 3. In this table are also included the values of the C(1)-C(2)-C(3)-C(3) and C(2)-C(3)-C(4)-C(5) dihe‐ drals, abbreviated as D(1) and D(2), respectively. The relative stability of the conformers in sol‐ ution is governed by a balance between two effects: interaction with the solvent (quantified by Δ*G*hyd) and intramolecular interaction (quantified by *E*int). In general, conformers with a favora‐

participate in a complex intermolecular hydrogen bonding network [63].

creases. For example at 45°C, the fraction of β-crystals is 45% [74, 75].

**crystal: Glutamic acid**

mal corrections at 298.15K.

The results just presented provide an important basis to understand the different crystalliza‐ tion ability of erythritol and *L*-threitol. In the *meso* form, this process is facilitated due to the resemblance between the backbone conformation of the molecules existing in the crystal lat‐ tice and that characteristic of the predominant forms in solution. The fact that such resem‐ blance is not found to exist in *L*-threitol may contribute to its lower crystallization tendency. This behavior is in agreement with the results obtained for other alditols [68, 73]

**Figure 8.** Most stable conformers of erythritol and *L*-threitol in gas phase and aqueous solution, as well as the respec‐ tive crystal conformations. The crystal conformation of erythritol corresponds to that having highest occupancy in the crystal lattice (Conformation A).

### **5. Conformational variation during the molecular incorporation into the crystal: Glutamic acid**

Regarding *L*-threitol, its gas phase conformational mixture is largely dominated by the *gtg* conformers which are characterized by a cyclic, symmetric and cooperative hydrogen bond‐ ing network. Since this OH groups' arrangement is very unfavorable to interact with water, they practically disappear in solution and are replaced by the *gg'g'*, *g'tg'* or *g'tg* conforma‐ tions which are not so highly intramolecularly bonded. Unlike erythritol, none of these con‐ formations has been identified in the crystalline *L*-threitol [72]. Here, the molecules adopt a distended conformation (*ttt*) which has a relatively low weight in solution (8%) and does not

The results just presented provide an important basis to understand the different crystalliza‐ tion ability of erythritol and *L*-threitol. In the *meso* form, this process is facilitated due to the resemblance between the backbone conformation of the molecules existing in the crystal lat‐ tice and that characteristic of the predominant forms in solution. The fact that such resem‐ blance is not found to exist in *L*-threitol may contribute to its lower crystallization tendency.

Gas

Solution

*gtg*

gg'g'

*ttt*

Crystal

**Erythritol** *L***-Threitol**

**Figure 8.** Most stable conformers of erythritol and *L*-threitol in gas phase and aqueous solution, as well as the respec‐ tive crystal conformations. The crystal conformation of erythritol corresponds to that having highest occupancy in the

This behavior is in agreement with the results obtained for other alditols [68, 73]

*tgg*

C(3) O(3)

C(4)

O(4)

C(2)

O(2)

O(1)

C(1)

g'tg

g'tg

crystal lattice (Conformation A).

exist in the vacuum.

212 Advanced Topics on Crystal Growth

Glutamic acid, 2-aminopentanedioic acid (C5H9NO4), is another example illustrating the role played by conformation on crystallization. Two conformational polymorphs have been identified for this compound, labeled as α and β, with the latter being the thermodynamical‐ ly stable form over all temperature range [74-77]. The crystals have different morphologies: the metastable form exhibits a prismatic shape while the stable one has a needle-like shape. Both belong to the ortorrombic space group with four molecules in the unit cell (*Z*=4) [78, 79]. In the crystal lattice the molecules are in the zwitterionic state and all functional groups participate in a complex intermolecular hydrogen bonding network [63].

The conformations adopted by the *L*-glutamic acid molecules in the two polymorphs are displayed in Figure 9. Their main difference lies in the orientation of the two C-C-C-C dihe‐ dral angles. In the crystal lattice the C(1)-C(2)-C(3)-C(4) and C(2)-C(3)-C(4)-C(5) dihedrals assume, respectively, values of 59.2°, 68.3° in α and -171.1°, -73.1° in β. Using the same dihe‐ dral labeling scheme previously adopted for the polyols, the conformation of glutamic acid in the α and β-crystals is *gg* and *tg'*, respectively. Depending on the experimental conditions, namely the temperature, both polymorphs can be crystallized from aqueous solution. Cool‐ ing a supersaturated solution to temperatures below 25°C originates almost pure α-crystals, whereas as the temperature raises the proportion of form β in the precipitated crystals in‐ creases. For example at 45°C, the fraction of β-crystals is 45% [74, 75].

The conformational behavior of zwitterionic glutamic acid is aqueous solution is a useful starting point to understand the crystallization of this compound in molecular terms. This has been done theoretically by performing full geometry optimizations using the CPCM continuum solvation model and the B3LYP/aug-cc-pVDZ model chemistry, both imple‐ mented in the Gaussian 03 program. The cavity was built with the Bondii radii which have been found to yield accurate results for the hydration of similar molecules in the zwitterion‐ ic state, such for example glycine [80]. Nine starting geometries were built by assuming the three standard orientations (*g*, *g*' and *t*) for each one of the C-C-C-C dihedrals. The remain‐ ing dihedrals were kept in their preferred orientations. The optimized structures were fur‐ ther submitted to a vibrational frequency calculation at the same level to ensure that they correspond to minima on the aqueous potential energy surface and also to calculate the ther‐ mal corrections at 298.15K.

The Gibbs energies of the conformers in aqueous solution at 298.15 K (*G*sol) can be expressed as: *G*sol = *E*int + Δ*G*thermal + Δ*G*hyd, with *E*int representing the intrinsic (abbreviated as "int") electronic energy of the conformer at 0 K, *i.e.*, excluding solvent effects, Δ*G*thermal the thermal correction to the Gibbs energy from 0 to 298K and Δ*G*hyd the Gibbs energy of hydration at 298.15K. The val‐ ues of *E*int, *G*sol, Δ*G*hyd, as well as the equilibrium populations at 298.15 K, are given in Table 3. In this table are also included the values of the C(1)-C(2)-C(3)-C(3) and C(2)-C(3)-C(4)-C(5) dihe‐ drals, abbreviated as D(1) and D(2), respectively. The relative stability of the conformers in sol‐ ution is governed by a balance between two effects: interaction with the solvent (quantified by Δ*G*hyd) and intramolecular interaction (quantified by *E*int). In general, conformers with a favora‐ ble intramolecular interaction are weakly hydrated and *vice-versa*. The most abundant con‐ former in solution, *tg* (Pop. = 67%), is the one having a less negative ∆*G*hyd; however, the intrinsic energy decrease arising from the formation of an internal N+ –H‧‧‧O hydrogen bond (see Figure 9), turns it into the lowest energy conformer. Unlike, the second most stable con‐ former (*tt*) has a fully distended backbone that hampers the formation of the just referred intra‐ molecular hydrogen bond and so it has a relatively high intrinsic energy. Its stabilization results essentially from a favorable interaction with the solvent. Both conformers represent *ca*. 81% of the conformational composition in aqueous solution at 298.15K.

this purpose, we have used the synchronous transit-guided quasi-Newton method (STQN) [81] in its QST2 variety at the CPCM/B3LYP/aug-cc-pVDZ level. The activation energies corre‐ sponding to the interconversion of the most stable conformer into conformers α and β were cal‐

explain the preferential crystallization of form α at temperatures below 25°C and the increase of the percentage of form β in the crystals obtained at higher temperatures. Apparently, the de‐ pendence of the crystallization behavior of glutamic acid on temperature is not due to a change on the relative abundance of the crystallizing conformers in solution with temperature [75], but rather to the value of the energy barriers for the interconversion of the most stable conformer in

**kJ mol-1**

*tg* 178.5 75.3 13.1 184.8 0.0 67.1 *tt* 168.9 -177.2 56.9 223.7 3.9 13.6 *g'g'* -52.5 -79.9 0.0 168.2 4.8 9.8 *gt* 63.1 178.8 36.4 196.1 7.4 3.5 *g´t* -56.5 178.7 52.1 213.7 7.7 3.0 *tg'*(β) 172.5 -69.8 59.6 221.9 8.4 2.3 *gg'* 62.4 -85.1 49.8 208.1 12.5 0.4 *gg*(α) 66.1 61.1 54.7 210.5 14.0 0.2 *g´g* -61.7 81.4 71.1 225.3 18.3 0.0

**Table 3.** Dihedral angles, intrinsic electronic energies at 0K (*E*int), Gibbs energies in aqueous solution (*G*sol), Gibbs energies of hydration (Δ*G*hyd) and Boltzmann populations (Pop.) at 298.15K of the zwitterionic glutamic acid

The molecular aggregates individualized in the early stages of nucleation (pre-nucleation) are windows to follow the crystallization process. Data on the self-assembly of erythritol and *L*-threitol in supersaturated aqueous solutions are reported through molecular dynam‐ ics (MD) simulations. MD is one of the most comprehensive computer simulation tools, with enormous application in the study of large systems and molecular phenomena requiring longer observation times. Relevant molecular phenomena are often adequately described by classical mechanics, providing reliable force-fields are available for a variety of systems. The results of the MD simulations presented in this chapter were carried out in the NpT ensem‐ ble and under periodic boundary conditions, resorting to the GROMACS package, version 4.5.4 [82] and GROMOS 96 43a1 force field [83]. Long-range electrostatics was computed us‐ ing the particle mesh Ewald (PME) method, while for Lennard-Jones energies a cut-off of 1.4

**D(1) D(2) (%)**

**Dihedral angles/º** *E***int/**

= 18 and 29 kJ mol-1, respectively. This barrier height difference may help to

**-Δ***G***hyd/ kJ mol-1**

*G***sol/ kJ mol-1**

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

215

**Pop.**

culated to be ∆*G*‡

**Conf.**

solution into the α and β conformers.

conformers, calculated at the CPCM/B3LYP/aug-cc-pVDZ level.

**6. Molecular dynamics study of pre-nucleation clusters**

**Figure 9.** CPCM/B3LYP/aug-cc-pVDZ optimized geometries of relevant conformers of the zwitterionic glutamic acid in aqueous solution. Dashed line represents an intramolecular hydrogen bond. Atom numbering scheme is shown for the most stable conformer.

Regarding the two crystal conformations, *tg*' (β) and *gg* (α), they are stable forms in aqueous solution with the values of the C-C-C-C dihedral angles in solution not differing significantly from those in the crystal. However, somehow surprisingly, their estimated equilibrium popu‐ lations in solution are too low: 2.3 and 0.2% for conformers β and α, respectively. Glutamic acid thus constitutes an interesting example of a molecule where the conformers experimentally ob‐ served in the crystalline structure are much diluted in aqueous solution. Two immediate con‐ clusions can be drawn from this result: (1) the conformational composition in aqueous solution does not determine the conformations exhibited by this compound in the final crystalline state; (2) crystallization is accompanied by a conformational change. The second conclusion implies that conformers *tg* and *tt*, namely the first, are interconverted into the *gg* and *tg'*. These inter‐ conversions can be characterized by the energy barriers (∆*G*‡ ) separating these conformers. For this purpose, we have used the synchronous transit-guided quasi-Newton method (STQN) [81] in its QST2 variety at the CPCM/B3LYP/aug-cc-pVDZ level. The activation energies corre‐ sponding to the interconversion of the most stable conformer into conformers α and β were cal‐ culated to be ∆*G*‡ = 18 and 29 kJ mol-1, respectively. This barrier height difference may help to explain the preferential crystallization of form α at temperatures below 25°C and the increase of the percentage of form β in the crystals obtained at higher temperatures. Apparently, the de‐ pendence of the crystallization behavior of glutamic acid on temperature is not due to a change on the relative abundance of the crystallizing conformers in solution with temperature [75], but rather to the value of the energy barriers for the interconversion of the most stable conformer in solution into the α and β conformers.

ble intramolecular interaction are weakly hydrated and *vice-versa*. The most abundant con‐ former in solution, *tg* (Pop. = 67%), is the one having a less negative ∆*G*hyd; however, the

(see Figure 9), turns it into the lowest energy conformer. Unlike, the second most stable con‐ former (*tt*) has a fully distended backbone that hampers the formation of the just referred intra‐ molecular hydrogen bond and so it has a relatively high intrinsic energy. Its stabilization results essentially from a favorable interaction with the solvent. Both conformers represent *ca*.

*tg tt*

*tg'* (β) *gg* (α)

**Figure 9.** CPCM/B3LYP/aug-cc-pVDZ optimized geometries of relevant conformers of the zwitterionic glutamic acid in aqueous solution. Dashed line represents an intramolecular hydrogen bond. Atom numbering scheme is shown for

Regarding the two crystal conformations, *tg*' (β) and *gg* (α), they are stable forms in aqueous solution with the values of the C-C-C-C dihedral angles in solution not differing significantly from those in the crystal. However, somehow surprisingly, their estimated equilibrium popu‐ lations in solution are too low: 2.3 and 0.2% for conformers β and α, respectively. Glutamic acid thus constitutes an interesting example of a molecule where the conformers experimentally ob‐ served in the crystalline structure are much diluted in aqueous solution. Two immediate con‐ clusions can be drawn from this result: (1) the conformational composition in aqueous solution does not determine the conformations exhibited by this compound in the final crystalline state; (2) crystallization is accompanied by a conformational change. The second conclusion implies that conformers *tg* and *tt*, namely the first, are interconverted into the *gg* and *tg'*. These inter‐

–H‧‧‧O hydrogen bond

) separating these conformers. For

intrinsic energy decrease arising from the formation of an internal N+

81% of the conformational composition in aqueous solution at 298.15K.

C(4)

C(5)

C(1) C(2)

214 Advanced Topics on Crystal Growth

the most stable conformer.

C(3)

conversions can be characterized by the energy barriers (∆*G*‡


**Table 3.** Dihedral angles, intrinsic electronic energies at 0K (*E*int), Gibbs energies in aqueous solution (*G*sol), Gibbs energies of hydration (Δ*G*hyd) and Boltzmann populations (Pop.) at 298.15K of the zwitterionic glutamic acid conformers, calculated at the CPCM/B3LYP/aug-cc-pVDZ level.

### **6. Molecular dynamics study of pre-nucleation clusters**

The molecular aggregates individualized in the early stages of nucleation (pre-nucleation) are windows to follow the crystallization process. Data on the self-assembly of erythritol and *L*-threitol in supersaturated aqueous solutions are reported through molecular dynam‐ ics (MD) simulations. MD is one of the most comprehensive computer simulation tools, with enormous application in the study of large systems and molecular phenomena requiring longer observation times. Relevant molecular phenomena are often adequately described by classical mechanics, providing reliable force-fields are available for a variety of systems. The results of the MD simulations presented in this chapter were carried out in the NpT ensem‐ ble and under periodic boundary conditions, resorting to the GROMACS package, version 4.5.4 [82] and GROMOS 96 43a1 force field [83]. Long-range electrostatics was computed us‐ ing the particle mesh Ewald (PME) method, while for Lennard-Jones energies a cut-off of 1.4 nm was applied. Temperature (298K) and pressure (1bar) were coupled to Berendsen exter‐ nal baths, with coupling constants of 0.1 and 0.5 ps, respectively. Each system was firstly subjected to an energy minimization step, and then left to evolve up to 80ns, using in both parts a standard time step of 2 fs. The last 40 ns of production runs were subsequently sub‐ jected to standard analysis, such as radial distribution functions [RDF, g(*r*)].

of all RDFs up to 0.28 nm (cut-off distance for the formation of an H-bond [85]) reveals that

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

217

Regarding *L*-threitol, the RDF profiles are substantially different from those obtained for eryth‐ ritol. In fact, the most relevant RDFs, labeled in Figure 10 as c1 to c5 and corresponding to the O(3)H–O(3)H, O(2)H–O(3)H, O(2)H–O(4)H, O(1)H–O(3)H and O(2)H–O(2)H interactions, show relatively well defined first peaks centered at 0.15 nm. The first peaks heights are much lower than in erythritol which means that the organized part of the *L*-threitol clusters is much smaller as compared to that in erythitol. The remaining part of the clusters has a disordered structure as evidenced by the broad peaks located at larger distances. This result is an argu‐ ment in favor of the lower crystallization tendency of this compound, in addition to the confor‐ mational one discussed in section 4. The snapshots of *L*-threitol clusters during the simulation

show the existence of various types of disordered oligomers (Figure 12).

**Figure 11.** Snapshots of the erythritol clusters formed in aqueous solution.

**Figure 12.** Snapshots of the *L*-threitol clusters formed in aqueous solution.

their probability of formation is lower than 15%.

Useful data about the type of aggregates formed in solution can be given by the RDFs for the different pairs of solute OH groups (OH–OH) which are depicted in Figure 10. They were obtained from the simulation of an aqueous solution containing 11 solute molecules and 4000 solvent molecules, corresponding to a concentration of 75g /100 cm3 , which is near to the saturation concentration [84].

**Figure 10.** OH–OH RDFs for erythritol and *L*-threitol in aqueous solution. For better visualization the most relevant curves are numbered as c1, c2, etc.

In erythritol, the most prominent RDFs are found for the O(2)H–O(3)H and O(1)H–O(4)H pairs. Both RDFs, labeled in Figure 10 as c1 and c2, respectively, exhibit sharp and intense peaks at *r* (nearest distance between groups) = 0.16 and 0.25 nm. This indicates that the selfassembly of erythritol in solution yields well organized clusters formed by hydrogen bonds involving preferentially the O(1)H and O(4)H groups or the O(2)H and O(3)H groups. Anal‐ ysis of the MD trajectory shows that these clusters are mostly dimers with some of them in‐ volving simultaneously the two just referred H-bonds (cyclic dimers). Interestingly, the cyclic dimer formed by these two H-bonds is the one existing in the crystalline structure. Some illustrative snapshots of these dimers are displayed in Figure 11. From the other OH– OH RDFs shown in Figure 10 one can conclude that the formation of intermolecular hydro‐ gen bonds involving the remaining groups is very unlikely to occur. In fact, the integration of all RDFs up to 0.28 nm (cut-off distance for the formation of an H-bond [85]) reveals that their probability of formation is lower than 15%.

Regarding *L*-threitol, the RDF profiles are substantially different from those obtained for eryth‐ ritol. In fact, the most relevant RDFs, labeled in Figure 10 as c1 to c5 and corresponding to the O(3)H–O(3)H, O(2)H–O(3)H, O(2)H–O(4)H, O(1)H–O(3)H and O(2)H–O(2)H interactions, show relatively well defined first peaks centered at 0.15 nm. The first peaks heights are much lower than in erythritol which means that the organized part of the *L*-threitol clusters is much smaller as compared to that in erythitol. The remaining part of the clusters has a disordered structure as evidenced by the broad peaks located at larger distances. This result is an argu‐ ment in favor of the lower crystallization tendency of this compound, in addition to the confor‐ mational one discussed in section 4. The snapshots of *L*-threitol clusters during the simulation show the existence of various types of disordered oligomers (Figure 12).

**Figure 11.** Snapshots of the erythritol clusters formed in aqueous solution.

nm was applied. Temperature (298K) and pressure (1bar) were coupled to Berendsen exter‐ nal baths, with coupling constants of 0.1 and 0.5 ps, respectively. Each system was firstly subjected to an energy minimization step, and then left to evolve up to 80ns, using in both parts a standard time step of 2 fs. The last 40 ns of production runs were subsequently sub‐

Useful data about the type of aggregates formed in solution can be given by the RDFs for the different pairs of solute OH groups (OH–OH) which are depicted in Figure 10. They were obtained from the simulation of an aqueous solution containing 11 solute molecules

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

*L*-Threitol

Erythritol

*r* / nm

**Figure 10.** OH–OH RDFs for erythritol and *L*-threitol in aqueous solution. For better visualization the most relevant

In erythritol, the most prominent RDFs are found for the O(2)H–O(3)H and O(1)H–O(4)H pairs. Both RDFs, labeled in Figure 10 as c1 and c2, respectively, exhibit sharp and intense peaks at *r* (nearest distance between groups) = 0.16 and 0.25 nm. This indicates that the selfassembly of erythritol in solution yields well organized clusters formed by hydrogen bonds involving preferentially the O(1)H and O(4)H groups or the O(2)H and O(3)H groups. Anal‐ ysis of the MD trajectory shows that these clusters are mostly dimers with some of them in‐ volving simultaneously the two just referred H-bonds (cyclic dimers). Interestingly, the cyclic dimer formed by these two H-bonds is the one existing in the crystalline structure. Some illustrative snapshots of these dimers are displayed in Figure 11. From the other OH– OH RDFs shown in Figure 10 one can conclude that the formation of intermolecular hydro‐ gen bonds involving the remaining groups is very unlikely to occur. In fact, the integration

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

c6 c5 c4 c3

c1: O(2)-O(3) c2: O(1)-O(4) c3: O(1)-O(3) c4: O(1)-O(2) c5: O(3)-O(4) c6: O(2)-O(4)

c1: O(3)-O(3) c2: O(2)-O(3) c3: O(2)-O(4) c4: O(1)-O(3) c5: O(2)-O(2) , which is near

jected to standard analysis, such as radial distribution functions [RDF, g(*r*)].

and 4000 solvent molecules, corresponding to a concentration of 75g /100 cm3

to the saturation concentration [84].

216 Advanced Topics on Crystal Growth

c5 c4 c3 c2 c1

g(*r*)

curves are numbered as c1, c2, etc.

c2

c1

**Figure 12.** Snapshots of the *L*-threitol clusters formed in aqueous solution.

### **7. Conclusion**

Structural aspects of crystallization from solvents were pointed out throughout this chapter. Particular attention was paid to the conformational variation of flexible molecules during this process. Three compounds were taken as examples: erythritol, *L*-threitol and glutamic acid. The different crystallization behavior shown by the two alditols was understood in terms of the resemblance between the solution and crystal conformers. Regarding glutamic acid, it is a quite peculiar example since the conformers existing in the two identified confor‐ mational polymorphs have a negligible weight in aqueous solution. Their selective crystalli‐ zation has been interpreted on the grounds of the energy barriers separating the dominant conformer in solution and those found in both crystalline forms.

[5] Lommerse JPM, Motherwell WDS, Ammon HL, Dunitz JD, Gavezzotti A, Hofmann DWM, et al. A test of crystal structure prediction of small organic molecules. Acta

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

219

[6] Oganov AR, editor. Modern Methods of Crystal Structure Predictions. Weinheim:

[7] Woodley SM, Catlow R. Crystal structure prediction from first principles. Nature

[8] Oganov AR, Glass CW. Crystal structure prediction using ab initio evolutionary techniques: Principles and applications. Journal of Chemical Physics. 2006;124(24):

[9] Day GM, Cooper TG, Cruz-Cabeza AJ, Hejczyk KE, Ammon HL, Boerrigter SXM, et al. Significant progress in predicting the crystal structures of small organic molecules - a report on the fourth blind test. Acta Crystallographica Section B. 2009;65(2):

[10] Kim S, Orendt AM, Ferraro MB, Facelli JC. Crystal structure prediction of flexible molecules using parallel genetic algorithms with a standard force field. Journal of

[11] Day GM, S. Motherwell WD, Jones W. A strategy for predicting the crystal structures of flexible molecules: the polymorphism of phenobarbital. Physical Chemistry Chem‐

[12] Hohenberg P, Kohn W. Inhomogeneous Electron Gas. Physical Review.

[13] Kohn W, Sham LJ. Self-Consistent Equations Including Exchange and Correlation Ef‐

[14] Becke AD. Density-Functional Exchange-Energy Approximation with Correct

[15] Becke AD. Density-Functional Thermochemistry.3. The Role of Exact Exchange. Jour‐

[16] Lee CT, Yang WT, Parr RG. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron-Density. Physical Review B. 1988 Jan

[17] Dunning J, Thom H. Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen. Journal of Chemical Physics. 1989

[18] Woon DE, Dunning J, Thom H. Gaussian basis sets for use in correlated molecular calculations. III. The atoms aluminum through argon. Journal of Chemical Physics.

Asymptotic-Behavior. Physical Review A. 1988 Sep 15;38(6):3098-100.

Crystallographica Section B. 2000;56(4):697-714.

Computational Chemistry. 2009;30(13):1973-85.

fects. Physical Review. 1965;140(4A):A1133-A8.

nal of Chemical Physics. 1993 Apr 1;98(7):5648-52.

ical Physics. 2007;9(14):1693-704.

1964;136(3B):B864-B71.

15;37(2):785-9.

1989/01/15/;90(2):1007-23.

1993 1993/01/15/;98(2):1358-71.

Wiley-VCH; 2011.

244704-15.

107-25.

Materials. 2008;7(12):937-46.

The investigation of the structure of the molecular aggregates formed in solution by molecu‐ lar dynamics, here exemplified for erythritol and *L*-threitol, has proven to be a valuable con‐ tribution to better understand crystallization.

Polymorphism, an important property of the solid state structure with various practical im‐ plications, was called in the present chapter as a tangly pathway in the crystallization proc‐ ess. Additional crystalline structure complexity may result from crystal disorder, hardly to be included in the polymorphism concept.

### **Author details**


### **References**


[5] Lommerse JPM, Motherwell WDS, Ammon HL, Dunitz JD, Gavezzotti A, Hofmann DWM, et al. A test of crystal structure prediction of small organic molecules. Acta Crystallographica Section B. 2000;56(4):697-714.

**7. Conclusion**

218 Advanced Topics on Crystal Growth

**Author details**

J.S. Redinha1

**References**

B. 2005;61(5):511-27.

Structural aspects of crystallization from solvents were pointed out throughout this chapter. Particular attention was paid to the conformational variation of flexible molecules during this process. Three compounds were taken as examples: erythritol, *L*-threitol and glutamic acid. The different crystallization behavior shown by the two alditols was understood in terms of the resemblance between the solution and crystal conformers. Regarding glutamic acid, it is a quite peculiar example since the conformers existing in the two identified confor‐ mational polymorphs have a negligible weight in aqueous solution. Their selective crystalli‐ zation has been interpreted on the grounds of the energy barriers separating the dominant

The investigation of the structure of the molecular aggregates formed in solution by molecu‐ lar dynamics, here exemplified for erythritol and *L*-threitol, has proven to be a valuable con‐

Polymorphism, an important property of the solid state structure with various practical im‐ plications, was called in the present chapter as a tangly pathway in the crystallization proc‐ ess. Additional crystalline structure complexity may result from crystal disorder, hardly to

[1] Brian Henderson and Ralph HB. Crystal-Field Engineering of Solid-State Laser Mate‐

[2] Brooks JS. Organic crystals: properties, devices, functionalization and bridges to bio-

[3] Day GM, Motherwell WDS, Ammon HL, Boerrigter SXM, Della Valle RG, Venuti E, et al. A third blind test of crystal structure prediction. Acta Crystallographica Section

[4] Motherwell WDS, Ammon HL, Dunitz JD, Dzyabchenko A, Erk P, Gavezzotti A, et al. Crystal structure prediction of small organic molecules: a second blind test. Acta

and J. A. S. Almeida1

conformer in solution and those found in both crystalline forms.

, A.J. Lopes Jesus1,2, A.A.C.C. Pais1

1 University of Coimbra / Department of Chemistry, Portugal

rials. New York: Cambridge University Press; 2000.

Crystallographica Section B. 2002;58(4):647-61.

molecules. Chemical Society Reviews. 2010;39(7):2667-94.

2 University of Coimbra / Faculty of Pharmacy, Portugal

tribution to better understand crystallization.

be included in the polymorphism concept.


[19] Krishnan R, Binkley JS, Seeger R, Pople JA. Self-consistent molecular orbital meth‐ ods. XX. A basis set for correlated wave functions. Journal of Chemical Physics. 1980;72(1):650-4.

[34] Elangannan Arunan, Gautam R. Desiraju, Roger A. Klein, Joanna Sadlej, Steve Scheiner, Ibon Alkorta, et al. Definition of the hydrogen bond (IUPAC Recommenda‐

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

221

[35] Lee SJ, Min A, Kim Y, Ahn A, Chang J, Lee SH, et al. Conformationally resolved structures of jet-cooled acetaminophen by UV-UV hole-burning spectroscopy. Physi‐

[36] McCrone WC. In: Fox D, Labes MM, Weissberger A, editors. Physics and Chemistry of the Organic Solid State. London: Interscience Publishers; 1965. p. 725-67.

[37] Singhal D, Curatolo W. Drug polymorphism and dosage form design: a practical per‐

[38] Perlovich G, Volkova T, Bauer-Brandl A. Polymorphism of paracetamol. Journal of

[39] Kolesov BA, Mikhailenko MA, Boldyreva EV. Dynamics of the intermolecular hydro‐ gen bonds in the polymorphs of paracetamol in relation to crystal packing and con‐ formational transitions: a variable-temperature polarized Raman spectroscopy study.

[40] Boldyreva E, Drebushchak V, Paukov I, Kovalevskaya Y, Drebushchak T. DSC and adiabatic calorimetry study of the polymorphs of paracetamol. Journal of Thermal

[41] Di Martino P, Conflant P, Drache M, Huvenne JP, Guyot-Hermann AM. Preparation and physical characterization of forms II and III of paracetamol. Journal of Thermal

[42] Haisa M, Kashino S, Kawai R, Maeda H. The Monoclinic Form of p-Hydroxyacetani‐ lide. Acta Crystallographica Section B-Structural Science. 1976;32(4):1283-5.

[43] Bouhmaida N, Bonhomme F, Guillot B, Jelsch C, Ghermani NE. Charge density and electrostatic potential analyses in paracetamol. Acta Crystallographica Section B.

[44] Haisa M, Kashino S, Maeda H. The orthorhombic form of p-hydroxyacetanilide. Acta

[45] Nichols G, Frampton CS. Physicochemical characterization of the orthorhombic poly‐ morph of paracetamol crystallized from solution. Journal of Pharmaceutical Sciences.

[46] Burgina EB, Baltakhinov VP, Boldyreva EV, Shakhtschneider TP. IR Spectra of Para‐ cetamol and Phenacetin. 1. Theoretical and Experimental Studies. Journal of Structur‐

[47] Ivanova BB. Monoclinic and orthorhombic polymorphs of paracetamol - solid state linear dichroic infrared spectral analysis. Journal of Molecular Structure.

tions 2011). Pure and applied chemistry. 2011;83(8):1637-41.

spective. Advanced Drug Delivery Reviews. 2004;56(3):335-47.

Physical Chemistry Chemical Physics. 2011;13(31):14243-53.

cal Chemistry Chemical Physics. 2011;13(37):16537-41.

Thermal Analysis and Calorimetry. 2007;89(3):767-74.

Analysis and Calorimetry. 2004;77(2):607-23.

Analysis and Calorimetry. 1997;48(3):447-58.

Crystallographica Section B. 1974;30(10):2510-2.

2009;65(3):363-74.

1998;87(6):684-93.

2005;738(1-3):233-8.

al Chemistry. 2004;45(1):64-73.


[34] Elangannan Arunan, Gautam R. Desiraju, Roger A. Klein, Joanna Sadlej, Steve Scheiner, Ibon Alkorta, et al. Definition of the hydrogen bond (IUPAC Recommenda‐ tions 2011). Pure and applied chemistry. 2011;83(8):1637-41.

[19] Krishnan R, Binkley JS, Seeger R, Pople JA. Self-consistent molecular orbital meth‐ ods. XX. A basis set for correlated wave functions. Journal of Chemical Physics.

[20] Cossi M, Rega N, Scalmani G, Barone V. Energies, structures, and electronic proper‐ ties of molecules in solution with the C-PCM solvation model. Journal of Computa‐

[21] Cossi M, Barone V, Cammi R, Tomasi J. Ab initio study of solvated molecules: A new implementation of the polarizable continuum model. Chemical Physics Letters. 1996

[22] Klamt A, Schuurmann G. Cosmo - a New Approach to Dielectric Screening in Sol‐ vents with Explicit Expressions for the Screening Energy and Its Gradient. Journal of

[23] Weissbuch I, Kuzmenko I, Vaida M, Zait S, Leiserowitz L, Lahav M. Twinned Crys‐ tals of Enantiomorphous Morphology of Racemic Alanine Induced by Optically Re‐ solved.alpha.-Amino Acids; A Stereochemical Probe for the Early Stages of Crystal

[24] Etter MC. Hydrogen bonds as design elements in organic chemistry. Journal of Phys‐

[25] Weissbuch I, Lahav M, Leiserowitz L. Toward Stereochemical Control, Monitoring, and Understanding of Crystal Nucleation. Crystal Growth & Design. 2003

[26] Bernstein J. Polymorphism in Molecular Crystals. Oxford: Oxford University Press;

[27] Hilfiker R, editor. Polymorphism in the Pharmaceutical Industry. Weinheim: Willey-

[28] Brittain HG, editor. Polymorphism in Pharmaceutical Solids. 2 ed. New York: Infor‐

[29] Nangia A. Conformational Polymorphism in Organic Crystals. Accounts of Chemical

[31] Schmelzer JWP, editor. Nucleation Theory and Applications: WILEY-VCH; 2005.

[32] Laaksonen A, Talanquer V, Oxtoby DW. Nucleation: Measurements, Theory, and At‐ mospheric Applications. Annual Review of Physical Chemistry. 1995;46(1):489-524.

[33] Steiner T. The hydrogen bond in the solid state. Angewandte Chemie, International

[30] J.W.Mullin. Crystallization. 4 ed. Oxford: Butterworth-Heinemann; 2001.

the Chemical Society- Perkin Transactions 2. 1993(5):799-805.

Nucleation. Chemistry of Materials. 1994 2012/09/14;6(8):1258-68.

1980;72(1):650-4.

220 Advanced Topics on Crystal Growth

Jun 14;255(4-6):327-35.

2012/09/14;3(2):125-50.

ma Healthcare USA, Inc; 2009.

Research. 2008;41(5):595-604.

Edition in English. 2002;41(1):48-76.

2002.

VCH; 2006.

tional Chemistry. 2003 Apr 30;24(6):669-81.

ical Chemistry. 1991 2012/09/14;95(12):4601-10.


[48] Iogansen AV. Direct proportionality of the hydrogen bonding energy and the intensi‐ fication of the stretching v(XH) vibration in infrared spectra. Spectrochimica Acta, Part A: Molecular and Biomolecular Spectroscopy. 1999;55(7-8):1585-612.

[62] Habgood M, Grau-Crespo R, Price SL. Substitutional and orientational disorder in organic crystals: a symmetry-adapted ensemble model. Physical Chemistry Chemical

Crystallization: From the Conformer to the Crystal

http://dx.doi.org/10.5772/54447

223

[63] Redinha JS, Lopes Jesus AJ. Molecular Recognition and Crystal Growth. In: McEvoy JA, editor. Molecular Recognition: Biotechnology, Chemical Engineering and Materi‐

[64] Redinha JS, Lopes Jesus AJ. Crystal Growth of Pharmaceuticals from Melt. In: Bori‐ senko E, editor. Crystallization and material science of modern artificial and natural

[65] Pascard C, Tran Huu Dau E, Manoury P, Mompon B. Betaxolol hydrochloride: 1-{4- [2-(cyclopropylmethoxy)ethyl]phenoxy}-3-isopropylamino-2-propanol hydrochlor‐ ide, C18H30NO3+.Cl. Acta Crystallographica, Section C: Crystal Structure

[66] Lopes Jesus AJ, Nunes SCC, Ramos Silva M, Matos Beja A, Redinha JS. Erythritol: Crystal growth from the melt. International Journal of Pharmaceutics. 2010;388(1-2):

[67] Lopes Jesus AJ, Redinha JS. On the structure of erythritol and L-threitol in the solid state: An infrared spectroscopic study. Journal of Molecular Structure. 2009;938(1-3):

[68] Yu L, Reutzel-Edens SM, Mitchell CA. Crystallization and Polymorphism of Confor‐ mationally Flexible Molecules: Problems, Patterns, and Strategies. Organic Process

[69] Lopes Jesus AJ, Tomé LIN, Rosado MTS, Leitão MLP, Redinha JS. Conformational study of erythritol and threitol in the gas state by density functional theory calcula‐

[70] Lopes Jesus AL, Tomé LIN, Eusébio MES, Redinha JS. Determination of the Enthalpy of Solute-Solvent Interaction from the Enthalpy of Solution: Aqueous Solutions of Er‐ ythritol and L-Threitol. Journal of Physical Chemistry B. 2006 May 11, 2006;110(18):

[71] Ceccarelli C, Jeffrey GA, McMullan RK. A Neutron-Diffraction Refinement of the Crystal-Structure of Erythritol at 22.6 K. Acta Crystallographica Section B.

[72] Kopf J, Morf M, Zimmer B, Haupt ETK, Jarchow O, Köll P. The crystal and molecular

[73] Jeffrey GA, Kim HS. Conformations of the alditols. Carbohydrate Research.

[74] Kitamura M. Polymorphism in the crystallization of L-glutamic acid. Journal of Crys‐

structure of threitol. Carbohydrate Research. 1993;247:119-28.

Research & Development. 2000 2012/09/14;4(5):396-402.

tions. Carbohydrate Research. 2005;340(2):283-91.

Physics. 2011;13(20):9590-600.

crystals. Rijeka: InTech; 2012.

129-35.

156-64.

9280-5.

1980;36:3079-83.

1970;14(2):207-16.

tal Growth. 1989;96(3):541-6.

Communications. 1984;40(8):1430-2.

als Applications: Novapublishers; 2011.


[62] Habgood M, Grau-Crespo R, Price SL. Substitutional and orientational disorder in organic crystals: a symmetry-adapted ensemble model. Physical Chemistry Chemical Physics. 2011;13(20):9590-600.

[48] Iogansen AV. Direct proportionality of the hydrogen bonding energy and the intensi‐ fication of the stretching v(XH) vibration in infrared spectra. Spectrochimica Acta,

[49] Wendy C. Duncan-Hewitt, Grant DJW. True density and thermal expansivity of pharmaceutical solids: comparison of methods and assessment of crystallinity. Inter‐

[50] Joiris E, Martino PD, Berneron C, Guyot-Hermann A-M, Guyot J-C. Compression Be‐ havior of Orthorhombic Paracetamol. Pharmaceutical Research. 1998;15(7):1122-30.

[51] Di Martino P, Guyot-Hermann AM, Conflant P, Drache M, Guyot JC. A new pure paracetamol for direct compression: The orthorhombic form. International Journal of

[52] Roth GJ, Calverley DC. Aspirin, Platelets, and Thrombosis: Theory and Practice.

[53] Wheatley PJ. The crystal and molecular structure of aspirin. Journal of the Chemical

[54] Ouvrard C, Price SL. Toward Crystal Structure Prediction for Conformationally Flex‐ ible Molecules: The Headaches Illustrated by Aspirin. Crystal Growth & Design.

[55] Vishweshwar P, McMahon JA, Oliveira M, Peterson ML, Zaworotko MJ. The Predict‐ ably Elusive Form II of Aspirin. Journal of the American Chemical Society. 2005

[56] Bond AD, Boese R, Desiraju GR. On the Polymorphism of Aspirin. Angewandte

[57] Chang C-J, Díaz LE, Morin F, Grant DM. Solid-state 13C NMR study of drugs: Aspir‐

[58] Wen S, Beran GJO. Accidental Degeneracy in Crystalline Aspirin: New Insights from High-Level ab Initio Calculations. Crystal Growth & Design. 2012 2012/09/10;12(5):

[59] Bond AD, Boese R, Desiraju GR. What Is A Polymorph? Aspirin As A Case Study.

[60] Glusker JP, Lewis M, Rossi M. Crystal Structure Analysis for Chemists and Biologists

[61] Wilson CC. A basic introduction to thermal motions of atoms in crystal structures, the underlying potentials and the physical information available from their analysis.

(Methods in Stereochemical Analysis). New York: Wiley-VCH; 1994.

Part A: Molecular and Biomolecular Spectroscopy. 1999;55(7-8):1585-612.

national Journal of Pharmaceutics. 1986;28(1):75-84.

Pharmaceutics. 1996;128:1-8.

Society (Resumed). 1964:6036-48.

Blood. 1994;83(4):885-98.

222 Advanced Topics on Crystal Growth

2004;4(6):1119-27.

2169-72.

2012/09/10;127(48):16802-3.

Chemie International Edition. 2007;46(4):615-7.

American Pharmaceutical Review. 2007:1-4.

Crystallography Reviews. 2009 2012/09/11;15(1):3-56.

in. Magnetic Resonance in Chemistry. 1986;24(9):768-71.


[75] Kitamura M. Controlling factor of polymorphism in crystallization process. Journal of Crystal Growth. 2002;237-239:2205-14.

**Section 2**

**Inorganic Systems**


**Section 2**

**Inorganic Systems**

[75] Kitamura M. Controlling factor of polymorphism in crystallization process. Journal

[76] Davey RJ, Blagden N, Potts GD, Docherty R. Polymorphism in Molecular Crystals: Stabilization of a Metastable Form by Conformational Mimicry. Journal of the Ameri‐

[77] Ono T, ter Horst JH, Jansens PJ. Quantitative Measurement of the Polymorphic Transformation of l-Glutamic Acid Using In-Situ Raman Spectroscopy. Crystal

[78] Lehmann M, Koetzle T, Hamilton W. Precision neutron diffraction structure determi‐ nation of protein and nucleic acid components. VII. The crystal and molecular struc‐ ture of the amino acid L-lysine monohydrochloride dihydrate. Journal of Chemical

[79] Lehmann MS, Nunes AC. A short hydrogen bond between near identical carboxyl groups in the [alpha]-modification of l-glutamic acid. Acta Crystallographica, Section

[80] Kim CK, Park B-H, Lee HW, Kim CK. Comprehensive Studies on the Free Energies of Solvation and Conformers of Glycine: A Theoretical Study. Bulletin of the Korean

[81] Peng CY, Ayala PY, Schlegel HB, Frisch MJ. Using redundant internal coordinates to optimize equilibrium geometries and transition states. Journal of Computational

[82] Hess B, Kutzner C, van der Spoel D, Lindahl E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. Journal of Chemical

[83] Schuttelkopf AW, van Aalten DMF. PRODRG: a tool for high-throughput crystallog‐ raphy of protein-ligand complexes. Acta Crystallographica Section D. 2004;60(8):

[84] Cohen S, Marcus Y, Migron Y, Dikstein S, Shafran A. Water sorption, binding and solubility of polyols. Journal of the Chemical Society, Faraday Transactions.

[85] Gilli G, Gilli P. Towards an unified hydrogen-bond theory. Journal of Molecular

of Crystal Growth. 2002;237-239:2205-14.

224 Advanced Topics on Crystal Growth

can Chemical Society. 1997;119(7):1767-72.

Growth & Design. 2004;4(3):465-9.

Crystallography. 1972;2(5):225-33.

B: Structural Science. 1980;36(7):1621-5.

Chemical Society. 2011;32(6):1985 - 92.

Theory and Computation. 2008;4(3):435-47.

Chemistry. 1996;17(1):49-56.

1355-63.

1993;89(17):3271-5.

Structure. 2000;552(1-3):1-15.

**Chapter 9**

**Chemical Vapor Transport Reactions–Methods,**

A variety of processes of crystal growth proceeds via the gas phase. A short comparative overview on gas phase transports is given here. However, in the main we deal with the concept of Chemical Vapor Transport Reactions [1, 2]. The term "*Chemical Vapor Transport*" (CVT) summarizes heterogeneous reactions which show one shared feature: a condensed phase, typically a solid, has an insufficient pressure for its own volatilization. But the pertinent phase can be volatilized in the presence of a gaseous reactant, the *transport agent*, and deposits elsewhere, usually in the form of crystals. The deposition will take place if there are different external conditions for the chemical equilibrium at the position of crystallization than at the position of volatilization. Usually, different temperatures are applied for volatilization and

**Figure 1.** Scheme of CVT experiments for crystallization of solids in a temperature gradient.

© 2013 Schmidt et al.; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

and reproduction in any medium, provided the original work is properly cited.

© 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution,

**Materials, Modeling**

http://dx.doi.org/10.5772/55547

**1. Introduction**

crystallization, Figure 1.

Peer Schmidt, Michael Binnewies, Robert Glaum and Marcus Schmidt

**1.1. Chemical vapor transport reactions**

Additional information is available at the end of the chapter
