## Abstract

An exact analytical solution of an integro-differential model describing the transient nucleation of solid particles (nuclei) and their growth with fluctuating rates at the intermediate stage of bulk phase transitions in metastable systems is constructed. Two important cases of the Weber–Volmer–Frenkel–Zel'dovich and Mier nucleation kinetics are detailed for supercooled melts and supersaturated solutions.

## 1. Introduction

Supersaturation or supercooling of a metastable system is a main driving force for the two dominant phenomena of bulk phase transitions: nucleation and crystal growth. Nucleation is a process taking place on a nano- or mesoscale in various atomic and molecular systems. Nucleation is an important factor for numerous phenomena in fundamental and applied research, ranging from applications to the condensation and freezing of liquids to chemical industry and life science such as production of food, pharmaceuticals and specialty chemicals [1–5].

Dynamics of nucleation and growth of crystals is complex, and it is always difficult to describe the entire phase transition process with a single model. From the theoretical point of view, it is reasonable to distinguish three main stages of the nuclei system evolution after the preliminary period when a supercooled (supersaturated) state establishes (e.g. the liquid is instantaneously cooled to below the crystallization temperature by the constant amount Δ*θ*_{0}). The first stage corresponds to the initiation of critical nuclei of the new phase. The state of a metastable system at this initial stage is mostly not affected by the emerging nuclei, so that each of them might be regarded as evolving under the constant metastability level (the same level Δ*θ*_{0} as at the preliminary stage). The second intermediate stage comprises a combined process of the growth of existing particles and of the initiation of additional nuclei reducing metastability of the parent liquid phase. A theoretical analysis of this stage is greatly complicated owing to the negative feedback between the nucleation process and crystal growth dependent on a transient degree of the system metastability (its supercooling or supersaturation) and owing to the gradual reduction of metastability by the growing solid phase (in the case of a supercooled melt, its supercooling Δ*θ* decreases with time, 0<Δ*θ*<Δ*θ*_{0}). The late stage represents the relaxation process stage such as Ostwald ripening and agglomeration when the origination of new nuclei almost ceases and may be safely overlooked (this stage occurs when Δ*θ*/Δ*θ*_{0}≪1).

This theoretical subdivision into stages is connected with typical relaxation times: *τ*_{1}≪*τ*_{2}≪*τ*_{3} of the process. Here, *τ*_{1} is the time for attaining the steady nucleation rate, *τ*_{2} is the time during which the growing crystals virtually desupercool (desupersaturate) the melt (solution) and *τ*_{3} is the time of entire phase transition. The initial and final stages of the process respectively connected with characteristic times *τ*_{1} and *τ*_{3} have been studied theoretically in detail (see, among others, [6–9]), but the intermediate behaviour (characteristic process time *τ*_{2}) is much less understood. This is connected with unsurmountable mathematical difficulties of solving the nonlinear integro-differential equations for the intermediate stage of bulk phase transitions.

The main task of this study is to develop a new theoretical approach for the construction of exact analytical solutions of the integro-differential model, describing the nucleation kinetics and particle coarsening with fluctuating rates at the intermediate stage of phase transitions for a broad range of process conditions.

## 2. The model

Experimental studies show that the growth rate of particles frequently undergoes random fluctuations [10–12]. Let us analyse a system that does not initially contain any nuclei. We shall assume that the initial supercooling (supersaturation) of a metastable medium is Δ*θ*_{0} (Δ*C*_{0}). Then, solid phase crystals appear continuously, evolve and form a homogeneous mixture of spherical nuclei. With allowance of fluctuations in the rate of their growth, the nucleus radius *r*(*τ*) is a random quantity obeying the stochastic equation [13]
2.1where *τ* is the growth time, *g*(*τ*,*r*) is the growth rate of a particle, which has no fluctuations, *D* is a function, which expresses the rate of fluctuations (the coefficient of mutual Brownian diffusion of the particle) and *W* is the Wiener process.

The particle-radius density distribution function *f*(*τ*,*r*) is described by the corresponding Fokker–Plank type kinetic equation [14]
2.2supplemented by the following initial condition
2.3Here, *r*_{*} stands for the radius of critical nuclei.

The heat (mass) balance equation for the dimensionless supercooling *w*=Δ*θ*/Δ*θ*_{0} (supersaturation *w*=Δ*C*/Δ*C*_{0}) can be written as [15]
2.4where Δ*θ* and Δ*C* are the dimensional supercooling and supersaturation of a metastable medium, *b*=4*πL*_{V}/(3*ρ*_{m}*C*_{m}Δ*θ*_{0}) if the system under consideration is a supercooled melt, and *b*=4*πC*_{p}/(3Δ*C*_{0}) if the system is represented by a supersaturated solution. Here, *L*_{V} is the latent heat parameter, *C*_{p} is the concentration at saturation, and *ρ*_{m} and *C*_{m} are the density and specific heat of the mixture. Equation (2.4) reflects the reduction in supercooling (supersaturation) caused by the crystal growth.

The initial condition for *w* takes the form
2.5

By equating the flux of crystals of minimum size *r*_{*} to the nucleation frequency *I*, we obtain the first boundary condition for the density distribution function [16–18]
2.6

The nucleation frequency can be written as an exponential function of the energy barrier height [14–16,19]
2.7in the case of supercooled melts and
2.8in the case of supersaturated solutions. Here, *p* has the meaning of dimensionless Gibbs number for supercooled melts or supersaturated solutions, *γ*_{i} is the surface tension, *θ*_{p} is the phase transition temperature, *k*_{B} is the Boltzmann constant, *C*_{l} is the solute concentration of a metastable system, *ρ*_{s} is the density of the solid phase, *R*_{g} is the universal gas constant and *θ*_{s} is the temperature of the solution. Final determination of the nucleation frequency *I* depends on calculation of *I*_{*} which has been done in a number of previous works [20–22]. In particular, these investigations showed that *I*_{*} depends slightly on *w*. For this reason, we shall consider *I*_{*} as a known constant.

Let us especially emphasize that the Gibbs–Thomson effect (the dependence of the phase transition temperature *θ*_{p} on the liquid–solid surface tension *σ* and the interface curvature 1/*R*), *θ*_{p}=*θ*_{p0}[1−*σ*/(*L*_{V}*R*)], can be neglected at the intermediate stage of phase transitions (*θ*_{p0} is the phase transition temperature in the absence of surface tension and curvature effects). This is due to the fact that a characteristic radius of growing crystals *r*∼*R*≫*σ*/*L*_{V}. This is reflected in the fact that *σ*/*L*_{V}∼10^{−12} *to* 10^{−10} m, whereas *R* ranges from *r*=*r*_{*}∼10^{−9} m to *r*∼10^{−6} m (see, among others, [23–25]).

Let us also write down the empirical expressions for the nucleation frequency which are frequently used in analysing many industrial processes for supercooled melts 2.9and supersaturated solutions 2.10Note that expressions (2.7)–(2.10) are frequently called the Weber–Volmer–Frenkel–Zel'dovich and Mier nucleation kinetics.

The distribution function vanishes at infinity, i.e. 2.11

To finalize the problem's formulation, one must find the coefficient *D* entering in the Fokker–Plank-type equation. The exact determination of *D* is a difficult problem of statistical physics. However, in a first approximation, it is natural to assume that *D* is proportional to the growth rate of nuclei [14,15]
2.12where *d*_{1} is a pertinent factor.

The rate of particle growth *g* can be found from a Stefan-type problem of evolution of a spherical aggregate in a metastable system [15,26]. The final result can be expressed in the form of
2.13if the system is a supercooled melt and
2.14if the system is a supersaturated solution. Here, *β*_{*} is the kinetic coefficient, *λ*_{l} and *D*_{l} are the temperature conductivity and diffusion coefficients. Expressions (2.13) and (2.14) demonstrate that the growth rate represents a linear function of the supercooling (supersaturation). In the particular case of very fine crystals, *r*≪(*β*_{*}*q*)^{−1}, the rate of growth does not depend on particle size. This growth mode is known as ‘kinetic’, because it is fully defined by surface processes. The growth rate of crystals, whose dimensions exceed markedly the value (*β*_{*}*q*)^{−1}, is controlled by the rate of heat removal. This regime is frequently termed as diffusionally controlled growth.

There exist additional growth theories in which expressions different from (2.13) and (2.14) are used (see, among others, [3,27–29]). Equations (2.13) and (2.14) do not make allowance for hindrance to the transfer of heat and mass through the two-phase mixture to each of the solid particles. Note that they are approximately correct only for systems with low nuclei concentration. When necessary, allowance for this can be made on the basis of theoretical considerations [30]. In addition, expressions (2.13) and (2.14) do not include the dependence of the growth rate on the curvature of the interface between the two phases and, by virtue of this fact, for the transfer of substance from fine to coarser nuclei. The expressions for the rate of crystal growth for this situation were analysed by a number of authors (see, among others, [14,25,31]). Note that expression (2.14) describes the limiting case of the final stage of phase transitions (Ostwald ripening) [25]. It is significant that more simple expressions (2.13) and (2.14) for the growth rate have been checked against experimental data. It was found that these rates adequately describe the corresponding phase transition process at its intermediate stage in many cases (see, among others, [3,27]. On this basis, we shall restrict ourselves to analysing growth rates (2.13) and (2.14), while considering systems with low solids concentrations at the intermediate stage of phase transformations. An important point is that the developed methods can be easily adjusted for use with a large variety of expressions for the growth rate.

Expressions (2.1)–(2.14) comprise the total system of integro-differential equations for describing the intermediate stage of bulk phase transitions in supercooled melts and supersaturated solutions. Below, we demonstrate a theoretical approach for its analytical solution.

## 3. Exact analytical solutions

Let us introduce the following dimensionless variables and parameters (*I*_{0}=*I* at *w*=1)
where
in the case of supercooled melts and
in the case of supersaturated solutions. Here, Δ*θ*=*θ*_{p}−*θ*_{l}>0, Δ*C*=*C*_{l}−*C*_{p}>0, Δ*θ*_{0}=*θ*_{p}−*θ*_{0}, Δ*C*_{0}=*C*_{0}−*C*_{p}, *θ*_{l} is the melt temperature, *θ*_{0} and *C*_{0} are the initial values of *θ*_{l} and *C*_{l}.

The aforementioned equations can be rewritten in dimensionless form valid for supercooled melts and supersaturated solutions as (*s*_{*}=*r*_{*}/*l*_{0})
3.1
3.2
3.3
3.4
and
3.5where function *φ*(*w*) is presented in table 1.

The particle growth rate in dimensionless form is expressed by the following differential equation
Its integration gives the growth rate in the form
3.6at *t*≥*t*_{*} and *g*_{0}(*t*)=0 at *t*<*t*_{*}, where *t*_{*} is the time of initiation of a critical nucleus with radius *s*_{*}.

Now introducing the modified time
3.7we rewrite equations (3.1)–(3.3) and (3.5) as (*F*=*F*(*x*,*s*))
3.8and
3.9where *z*=*s*−*s*_{*}. Note that substitution (3.6) into (3.7) gives us the growth rate in the form
3.10

Equation (3.8) with the boundary conditions (3.9) can be solved by means of the Laplace integral transform with respect to *x* with parameter *λ*. The transforms of the corresponding quantities will be denoted by the index *L*. Then, the boundary-value problem (3.8)–(3.9) has the following solution
3.11The solution in real space (*x*,*z*) is then obtained by inverting the transform given by expression (3.11). Thus, the density distribution function *F*(*x*,*s*) can be written in the form of
3.12where *z*(*s*)=*s*−*s*_{*}.

Substitution of expression (3.12) into (3.4) gives the following integral equation for the dimensionless supercooling (supersaturation) 3.13Note that the last integral can be analytically calculated in an explicit form in terms of special functions.

In order to solve the integral equation (3.13), we change the variable in (3.13) considering the inverse function *y*=*y*(*w*). Taking into account expressions (3.5) and (3.10), we have
where *P*(*y*,*x*)=[1+*Q*(*x*−*y*+*s*_{*})]. Differentiating this equation with respect to *w*, we arrive at a differential equation with separable variables. Its integration gives the implicit function *w*(*x*) of the form
3.14It is significant that the integral in the right-hand side can be simplified analytically or easily calculated numerically.

Let us now calculate the left-hand side of equation (3.14) for different kinetic mechanisms under consideration. So, in the case of the Mier kinetics for supercooled melts and supersaturated solutions, we obtain the following explicit expressions
3.15In the case of the Weber–Volmer–Frenkel–Zel'dovich kinetics, the analytical solution can be written in an implicit form as
3.16for supercooled melts and
3.17for supersaturated solutions. Here, *E*(*w*) and *R*(*w*) are the following integral functions (, )
Note that all analytical solutions (3.14)–(3.17) satisfy the original integral equation (3.13). This result can be easily checked by means of direct substitution of the aforementioned explicit and implicit solutions into (3.13).

Now rewriting *x*(*t*) from expression (3.7) in terms of its inverse function *t*(*x*), we obtain the solution in a closed parametric form (with parameter *x*)
3.18where *w*(*x*_{1}) is determined by expressions (3.14)–(3.17) for different kinetic regimes. This means that we were able to obtain a complete analytic definition of the kinetics of desupercooling (desupersaturation) and the distribution of crystal sizes in a supercooled (supersaturated) system.

## 4. Concluding remarks

Figures 1 and 2 illustrate the exact analytical solutions shown in accordance with expressions (3.12), (3.15) and (3.18) for the Mier kinetics of particle nucleation. At small times, when supercooling (supersaturation) plays an important role, the density distribution is a high and narrow function of the dimensionless radius. This is caused by the active nucleation process at the beginning of the intermediate stage. After the period of continuing initiation of new nuclei, the metastability reduction becomes lower, and the distribution function becomes more and more broad. This means that the density distribution of the nuclei by size is characterized by increasing dispersion. In accordance with this, the maximum value of this function decreases as time increases. Figure 3 shows that the exact analytical solution describes experimental data well over wide range of time scales. It is significant that the theory under consideration takes into account that the dimensionless supercooling (supersaturation) *w*(*τ*) can be an upwards (downwards) convex function at small (large) *τ* (figure 3). This behaviour reflects a smooth transition between the intermediate, initial and final stages of the process.

It should be emphasized in conclusion that expressions (3.12), (3.14)–(3.18) represent an exact analytical solution of the integro-differential model of particle nucleation and coarsening at the intermediate stage of bulk phase transitions. The complete solution is found in a form which is valid for different kinetics (nucleation frequencies) expressed in terms of arbitrary dependence *φ*(*w*_{1}). The analysis presented herein thus makes it possible to solve all the principal problems related to kinetics of removal of supercooling (supersaturation) in a single-component system. These results can also be extended to binary systems in the spirit of works [15,33,34].

An important point of the analytical approach under consideration is that it can be used now to describe similar nucleation phenomena met in metallurgy (solidification of alloys with a two-phase layer [35]), geophysics (freezing processes in terrestrial magma oceans [36] and crystallization at the inner core boundary of the Earth [37]) and applied physics (aggregate growth in colloids, magnetic fluids and other systems of a like nature [18,38,39]). As a final note, the analytical solutions under consideration can be used for determining the initial distribution of nuclei for the concluding stage of phase transition, when there is virtually no supercooling (supersaturation).

## Acknowledgements

Partial support from the Ministry of Education and Science of the Russian Federation (projects nos. 14.A18.21.0858 and 14.A18.21.1126) and from the Russian Foundation for Basic Research (grants nos. 11-01-00137, 13-01-96013-Ural and 13-08-96013-Ural) is acknowledged.

- Received September 28, 2013.
- Accepted November 12, 2013.

- © 2013 The Author(s) Published by the Royal Society. All rights reserved.