The Boundary Integral Equation for Kinetically Limited Dendrite Growth

Abstract

The boundary integral equation defining the interface function for a curved solid/liquid phase transition boundary is analytically solved in steady-state growth conditions. This solution describes dendrite tips evolving in undercooled melts with a constant crystallization velocity, which is the sum of the steady-state and translational velocities. The dendrite tips in the form of a parabola, paraboloid, and elliptic paraboloid are considered. Taking this solution into account, we obtain the modified boundary integral equation describing the evolution of the patterns and dendrites in undercooled binary melts. Our analysis shows that dendritic tips always evolve in a steady-state manner when considering a kinetically controlled crystallization scenario. The steady-state growth velocity as a factor that is dependent on the melt undercooling, solute concentration, atomic kinetics, and other system parameters is derived. This expression can be used for determining the selection constant of the stable dendrite growth mode in the case of kinetically controlled crystallization.

1. Introduction

Phase transformations from nonequilibrium and metastable states to the solid form are widespread in nature (the freezing of water, the solidification of magma, etc.) and are often used in technological processes and laboratory facilities to obtain alloys with certain structures and properties, as well as for the synthesis of various compounds in the chemical and pharmacological industries. This explains extensive studies of directional and bulk crystallization in the phase transformation region consisting of simultaneous liquid and solid phases [1,2,3,4,5,6]. As this takes place, the evolution of the solid/liquid phase interface in an undercooled or supersaturated liquid is responsible for the processes of impurity redistribution between the growing solid phase elements, pattern development, and microstructure formation. Therefore, the problem of the mathematical description of a phase transformation can be reduced to the problem of the evolution of a solid/liquid phase interface.

A single integral differential equation defining its evolution in an undercooled one-component liquid was derived for the first time by Nash and Glicksman [7,8]. This equation for the interface function was then extended for binary systems solidifying in local-equilibrium and local-nonequilibrium conditions [9,10], as well as for convective thermodiffusion problems [11,12,13]. An important point is that the boundary integral equation (BIE) represents the basis for the morphological stability analysis of dendritic crystals and the selection of their stable growth modes [14,15,16,17,18,19,20,21,22]. Moreover, the BIE is a formal solution of various Stefan-like problems and can be used for the numerical modeling of pattern formation and evolution [23]. In addition, the BIE can also be used for the verification of crystal growth theories and the comparison of model conclusions with computations carried out, for example, by the phase-field and enthalpy-based methods or by using direct simulations of moving-boundary problems [24,25,26,27,28].

It is well known that dendritic crystals reach steady-state growth very quickly [29,30,31,32,33]. Moreover, the solid–liquid interface may reach the kinetic regime of growth when approaching the steady state. Only the attachment/detachment of particles (atoms and molecules) at the interface controls the intensity of the interface migration. Special interest in the analysis of the kinetic mode of growth is also seen in the following:

-

the explanations of the sharp transition from thermally controlled to kinetically limited regimes of growth [34,35,36,37];

-

the interpretation of the data of molecular dynamics where the interface is governed by interfacial undercooling, i.e., by “kinetic undercooling” [38,39,40,41].

Taking a special interest in the study of stationary kinetic modes into account, we analyze the BIE for a unidirectional stationary crystallization scenario with a steady-state velocity $V_{ss}$ and show that the crystal tip evolves with a constant velocity $V=V_{ss}+V_{tr}$ in the same reference frame (here, $V_{tr}$ is a constant translational velocity in the moving reference frame). The theory under question is formulated for two- and three-dimensional problems describing a stable growth mode of dendritic crystals whose tips represent a parabola, paraboloid and elliptic paraboloid.

2. Two-Dimensional BIE for Solidifying Pure Liquid

Let us first consider the simplest case of dendritic growth in a single-component undercooled liquid with a constant (steady-state) velocity $V_{ss}$ in a laboratory reference frame (see Figure 1). We introduce a characteristic spatial scale $\rho$ and time scale $\rho /V_{ss}$, where $\rho$ can be the diameter of the crystal tip or the doubled radius of the curvature. Let us write down the two-dimensional BIE for the solidification of the pure (chemically one-component) liquid with the interface function $\zeta (x,t)$ by assuming these scales as follows [9] (see also the Appendix A):

$$\begin{array}{c}\hfill \begin{array}{c}\hfill \Delta -\frac{d_c}{\rho }K(x,t)-\beta V_{ss}-\beta V_{ss}\frac{\partial \zeta (x,t)}{\partial t}=I_{\zeta }^T(x,t),I_{\zeta }^T(x,t)=\frac{P_T}{2\pi }\underset{0}{\overset{\infty }{\int }}\frac{d\tau }{\tau }\underset{-\infty }{\overset{\infty }{\int }}d{x}_1\\ \hfill \times \left(1+\frac{\partial \zeta (x_1,t-\tau )}{\partial t}\right)exp\left[-\frac{P_T}{2\tau }\left({(x-x_1)}^2+{(\zeta (x,t)-\zeta (x_1,t-\tau )+\tau )}^2\right)\right].\end{array}\end{array}$$

(1)

Here, $\Delta =c_p(T_0-T_{\infty })/Q$ is the dimensionless undercooling, $c_p$ is the thermal capacity, Q is the latent crystallization heat, $T_0$ is the phase transition temperature, $T_{\infty }$ is the liquid temperature far from the solid/liquid interface $\zeta$, K is the interface curvature, $d_c$ is the capillary length, $\beta$ is the kinetic coefficient, $P_T=\rho V_{ss}/\left(2D_T\right)$ is the Péclet number, $D_T$ is the thermal diffusivity, and x and t are the dimensionless spatial and time variables, respectively. Equation (1) describes the interface motion, with velocity $V_{ss}$ defining the driving force $\Delta$, interface curvature $K\propto {\rho }^{-1}$, and the heat transport outgoing the interface; see the r.h.s. of Equation (1). Note that $\Delta$, $d_c$, $\rho$, $\beta$, $V_{ss}$, and $P_T$ in Equation (1) are assumed to be constant. The interface function $\zeta (x,t)$ should be found as a solution of this equation. Equation (1) must be supplemented with appropriate initial and boundary conditions. As an initial condition, one should choose a certain shape of the interface function $\zeta (x,t)$, e.g., a parabola, as is done below. It is also necessary to set boundary conditions reflecting the geometry of solidification region.

In the present work, we consider the special case of negligible heat transport to govern the interface motion. This situation is well known in atomistic simulations where the interface motion in melting and crystallization is dictated by the atomistic kinetics at the interface, and the computational domain has (almost) uniform temperature [38,39,40,41]. In this case, the r.h.s. of the BIE (1) can be omitted. Additionally, we only consider the motion of the crystal tip, where the interface curvature is $K\approx -{\partial }^2\zeta /\partial x^2$.

Taking this into account, we rewrite Equation (1) as

$$\left(\Delta -\beta V_{ss}\right)\frac{\rho }{d_c}+\frac{{\partial }^2\zeta }{\partial x^2}-\frac{\beta V_{ss}\rho }{d_c}\frac{\partial \zeta }{\partial t}=0.$$

We solve this equation using the Laplace integral transform with respect to the time variable t. By introducing the Laplace variable p and Laplace image ${\zeta }^{*}(x,p)={\zeta }^{*}\left(x\right)$ of the interface function $\zeta (x,t)$ as

$${\zeta }^{*}\left(x\right)=\underset{0}{\overset{\infty }{\int }}exp\left(-pt\right)\zeta (x,t)dt,$$

we obtain the following equation in the Laplace image space:

$$\begin{array}{c}\hfill \frac{\Delta -\beta V_{ss}}{p}\frac{\rho }{d_c}+\frac{d^2{\zeta }^{*}}{d{x}^2}-\frac{\beta V_{ss}p\rho }{d_c}{\zeta }^{*}+\frac{\beta V_{ss}\rho }{d_c}\zeta (x,0)=0.\end{array}$$

(2)

Here, we have taken into account the following Laplace transform formula:

$$\zeta (x,t)\to p{\zeta }^{*}\left(x\right)-\zeta (x,0),\mathrm{const}.\to \frac{\mathrm{const}.}{p}.$$

The solution of Equation (2) reads as

$$\begin{array}{c}\hfill \begin{array}{c}\hfill {\zeta }^{*}\left(x\right)=\frac{\sqrt{A}}{2\sqrt{p}}exp\left(-\sqrt{pA}x\right)\underset{0}{\overset{\infty }{\int }}\left(\zeta (x_1,0)-\zeta (-x_1,0)\right)d{x}_1+\frac{\Delta \rho -A{d}_c}{p^{2}A{d}_c}\\ \hfill +\frac{\sqrt{A}}{2}\left[\underset{x}{\overset{\infty }{\int }}\frac{exp\left[-\sqrt{pA}(x_1-x)\right]}{\sqrt{p}}\zeta (x_1,0)d{x}_1+\underset{-\infty }{\overset{x}{\int }}\frac{exp\left[\sqrt{pA}(x_1-x)\right]}{\sqrt{p}}\zeta (x_1,0)d{x}_1\right],\end{array}\end{array}$$

(3)

where $A=\beta V_{ss}\rho /d_c$.

Let us assume $\zeta (x,0)=-x^2/2$ as the initial condition for the interface function in the form of a parabola. By substituting $\zeta (x,0)$ into (3), we arrive at

$$\begin{array}{c}\hfill \begin{array}{c}\hfill {\zeta }^{*}\left(x\right)=\frac{\Delta \rho -A{d}_c}{p^{2}A{d}_c}\\ \hfill -\frac{\sqrt{A}}{4\sqrt{p}}\left[\underset{x}{\overset{+\infty }{\int }}exp\left(\sqrt{pA}(x-x_1)\right)x_1^{2}d{x}_1+\underset{-\infty }{\overset{x}{\int }}exp\left(\sqrt{pA}(x_1-x)\right)x_1^{2}d{x}_1\right].\end{array}\end{array}$$

(4)

By integrating the r.h.s. of Equation (4), we obtain

$$\begin{array}{c}\hfill {\zeta }^{*}\left(x\right)=\left(\frac{\Delta \rho -A{d}_c}{A{d}_c}-\frac{1}{A}\right)\frac{1}{p^2}-\frac{x^2}{2p}.\end{array}$$

(5)

By calculating the inverse Laplace transform of Equation (5), we come to [42,43]

$$\begin{array}{c}\hfill \zeta (x,t)=\left(\frac{\Delta \rho -A{d}_c}{A{d}_c}-\frac{1}{A}\right)t-\frac{x^2}{2}.\end{array}$$

(6)

Expression (6) shows that the dendrite grows in a moving reference frame with the following constant translational velocity:

$$\begin{array}{c}\hfill V_{tr}=V_{ss}\left(\frac{\Delta -\beta V_{ss}}{\beta V_{ss}}-\frac{d_c}{\beta V_{ss}\rho }\right).\end{array}$$

(7)

This means that the dendrite evolves with the velocity $V=V_{ss}+V_{tr}$ in the laboratory reference frame. As is easily seen, the contribution of $V_{tr}$ is caused by atomic kinetics (parameter $\beta$). Figure 2 shows the interface function $\zeta (x,t)$ versus the rescaled undercooling $\Delta /\left(\beta V\right)$ (panel (a)) and spatial coordinate x (panel (b)). As is easily seen, this function increases with an increasing driving force $\Delta /\left(\beta V\right)$ and time t. This is explained by the growth of the dendritic tip in the reference frame moving with the steady-state velocity $V_{ss}$. As can be easily understood, the velocity of this growth is constant and equal to $V_{tr}$.

3. Two-Dimensional BIE for Solidifying Pure and/or Binary Liquid

Let us now consider the case of crystal growth in a binary undercooled liquid with a solute concentration $C_{l\infty }$ in the binary liquid or solution far from the solid/liquid interface. In this case, the BIE takes the following form [9] (see also the Appendix A):

$$\begin{array}{c}\hfill \begin{array}{c}\hfill -\frac{Q}{m{c}_p}\left[\Delta -\frac{d_c}{\rho }K-\beta V_{ss}-\beta V_{ss}\frac{\partial \zeta }{\partial t}-I_{\zeta }^T\right]-C_{l\infty }=I_{\zeta }^C,\\ \hfill I_{\zeta }^C=\frac{(1-k_0)P_C}{2\pi }\underset{0}{\overset{\infty }{\int }}\frac{d\tau }{\tau }\underset{-\infty }{\overset{\infty }{\int }}d{x}_1{C}_i(x_1,t-\tau )\left(1+\frac{\partial \zeta (x_1,t-\tau )}{\partial t}\right)\\ \hfill \times exp\left[-\frac{P_C}{2\tau }\left({(x-x_1)}^2+{(\zeta (x,t)-\zeta (x_1,t-\tau )+\tau )}^2\right)\right].\end{array}\end{array}$$

(8)

Here, m and $k_0$ are the equilibrium liquid’s slope and partition coefficient, respectively, $P_C=\rho V_{ss}/\left(2D_C\right)$, $D_C$ is the diffusion coefficient, and $C_i$ is the interfacial concentration.

Again, let us assume kinetically limited crystal growth with $K\approx -{\partial }^2\zeta /\partial x^2$. This leads to $I_{\zeta }^T\to 0$ and $I_{\zeta }^C\to 0$. In this case, Equation (8) transforms to

$$\begin{array}{c}\hfill \Delta +\frac{m{c}_p{C}_{l\infty }}{Q}+\frac{d_c}{\rho }\frac{{\partial }^2\zeta }{\partial x^2}-\beta V_{ss}-\beta V_{ss}\frac{\partial \zeta }{\partial t}=0.\end{array}$$

(9)

Equation (9) can be solved using the Laplace transform by analogy with the aforementioned problem for a pure liquid. By omitting all mathematical manipulations, let us write out the final result as follows:

$$\begin{array}{c}\hfill \zeta (x,t)=\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{Q\beta V_{ss}}-1-\frac{d_c}{\beta V_{ss}\rho }\right)t-\frac{x^2}{2}.\end{array}$$

(10)

This means that the translational velocity for a binary system becomes

$$\begin{array}{c}\hfill V_{tr}=V_{ss}\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{Q\beta V_{ss}}-1-\frac{d_c}{\beta V_{ss}\rho }\right).\end{array}$$

(11)

As is easily seen, expression (11) transforms into expression (7) in the absence of a solute concentration, i.e., $C_{l\infty }=0$.

4. Three-Dimensional BIE for a Paraboloid of Revolution and Elliptic Paraboloid Growing from Undercooled Binary Liquid

We now consider three-dimensional crystal growth in an undercooled binary melt or solution. Again, considering kinetically controlled interface motion (neglecting heat and mass transport), we write out the corresponding BIE in the form of

$$\begin{array}{c}\hfill \Delta +\frac{m{c}_p{C}_{l\infty }}{Q}+\frac{d_c}{\rho }\left(\frac{{\partial }^2\zeta }{\partial x^2}+\frac{{\partial }^2\zeta }{\partial y^2}\right)-\beta V_{ss}-\beta V_{ss}\frac{\partial \zeta }{\partial t}=0,\end{array}$$

(12)

where the interface curvature is chosen as

$$K\approx -\frac{{\partial }^2\zeta }{\partial x^2}-\frac{{\partial }^2\zeta }{\partial y^2}$$

in the vicinity of a dendritic tip.

We solve this equation by analogy with the two-dimensional case. Considering the two shapes, a paraboloid of revolution and an elliptic paraboloid, we obtain the interface functions in the forms of

$$\begin{array}{c}\hfill \zeta (x,y,t)=\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{\beta V_{ss}Q}-\frac{2d_c}{\beta V_{ss}\rho }-1\right)t-\frac{x^2}{2}-\frac{y^2}{2},\end{array}$$

(13)

$$\begin{array}{c}\hfill \zeta (x,y,t)=\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{\beta V_{ss}Q}-\frac{2d_c}{\beta V_{ss}\rho (1-a^2)}-1\right)t-\frac{x^2}{2(1-a)}-\frac{y^2}{2(1+a)},\end{array}$$

(14)

where a is the ellipticity parameter. The translation velocity for these two shapes has the respective forms of

$$\begin{array}{c}\hfill \begin{array}{c}\hfill V_{tr}=\left\{\begin{array}{c}{\displaystyle V_{ss}\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{\beta V_{ss}Q}-\frac{2d_c}{\beta V_{ss}\rho }-1\right),\mathrm{paraboloid}\mathrm{of}\mathrm{revolution},}\hfill \\ {\displaystyle V_{ss}\left(\frac{\Delta Q+m{c}_p{C}_{l\infty }}{\beta V_{ss}Q}-\frac{2d_c}{\beta V_{ss}\rho (1-a^2)}-1\right),\mathrm{elliptic}\mathrm{paraboloid}.}\hfill \end{array}\right.\end{array}\end{array}$$

(15)

Note that the second line of (15) transforms into the first line at $a=0$ when an elliptical paraboloid becomes a paraboloid of revolution. In addition, the first line of (15) transforms into a two-dimensional case that accounts for the double curvature in three-dimensional geometry, i.e., $K_{3D}=2K_{2D}$.

5. Conclusions

In summary, this paper addresses the BIE for the steady-state growth mode of dendritic crystals. First, we show that a dendrite whose tip is connected with the steady-state crystallization velocity $V_{ss}$ grows in this reference frame with a constant translational velocity $V_{tr}$. As a result, the dendrite tip moves with a constant velocity $V=V_{ss}+V_{tr}$ in a laboratory reference frame. Accordingly, the same velocity V has been observed experimentally as the velocity of the recalescence front during the solidification of undercooled liquid droplets in electromagnetic levitators and other crystallization facilities (see, among others, [44,45]).

The constant translational crystallization velocity $V_{tr}$ has been found to be a function of the established liquid undercooling $\Delta$, dendrite tip diameter $\rho$, attachment kinetics of the atoms to the surface of the growing crystal $\beta$, steady-state velocity $V_{ss}$, and capillary constant $d_c$.

The main outcomes following from the aforementioned analysis are as follows:

(i) Considering kinetically limited crystal growth (neglecting heat and mass transport), we conclude that this process always occurs in a steady-state manner. In other words, there is no time to reach stationary growth, because the crystallization process is always in a steady state.

(ii) Moving to the frame of reference associated with the crystal surface (a system moving with the steady-state velocity $V_{ss}$), we see that a dendrite grows with a constant velocity $V_{tr}$ in this system. Its growth is caused by the driving force (undercooling $\Delta$) and the attachment of the atoms to the interphase boundary (kinetic parameter $\beta$). As this takes place, we can choose a reference frame that is always connected with the dendritic tip, i.e., $V_{tr}=0$, if $V_{ss}=V=\Delta /\beta -d_c/\left(\beta \rho \right)$ for the two-dimensional solidification of a pure liquid (the reference frame following the growth of the crystal tip). In a more general case, we find this velocity from the second line of expression (15) as $V_{ss}=V=\Delta /\beta +m{c}_p{C}_{l\infty }/\left(\beta Q\right)-2d_c/\left(\beta \rho (1-a^2)\right)$ for the three-dimensional thermodiffusion case.

(iii) Considering kinetically limited crystal growth (neglecting heat and mass transport), we can use the aforementioned expressions for $V_{ss}=V$ (at $V_{tr}=0$) as auxiliary relations to find the selection constant entering the criterion of stable dendritic growth [10,20,46].

(iv) The translational velocity $V_{tr}$ is defined by expressions (7), (11), or (15) in the case of kinetically limited crystallization when the integral contributions $I_{\zeta }^T$ and $I_{\zeta }^C$ are small enough. In a more general case, and especially in the case of rapid crystallization, $V_{tr}$ should be numerically found from the BIE [10] that contains the integral contributions $I_{\zeta }^T$ and $I_{\zeta }^C$.

(v) To move to a frame of reference moving together with the dendrite tip (moving with the constant velocity $V=V_{tr}$ relative to the laboratory coordinate system), we can choose $V_{ss}=0$ (since this velocity is arbitrary). In this case, the BIE becomes

$$\begin{array}{c}\hfill \begin{array}{c}\hfill -\frac{Q}{m{c}_p}\left[\Delta -\frac{d_c}{\rho }K-\beta V\left(1+\frac{\partial \zeta }{\partial t}\right)-I_{\zeta }^T\right]-C_{l\infty }=I_{\zeta }^C,\\ \hfill V=V_{tr}=\frac{\Delta Q+m{c}_p{C}_{l\infty }}{\beta Q}-\frac{2d_c}{\beta \rho (1-a^2)},\end{array}\end{array}$$

(16)

where $I_{\zeta }^T$ and $I_{\zeta }^C$ are defined by expressions (1) and (8), respectively, with $P_T=\rho V/\left(2D_T\right)$ and $P_C=\rho V/\left(2D_C\right)$. Here, the characteristic spatial and time scales are chosen in terms of $\rho$ and $\rho /V$, respectively. As this takes place, the Ivantsov solutions and selection criterion should be written in terms of V. Note that the second line of Equation (16) takes place only when $\beta \ne 0$.

The present theory can be extended to high-speed crystallization processes, which are described by a more complex BIE that takes into account the rate of atomic jumps to neighboring vacant positions [10]. From a mathematical point of view, such processes are described by the hyperbolic equation of impurity diffusion changing the concentration part of the BIE. In addition, the presence of a constant translational velocity $V_{tr}$ (the growth velocity of the dendritic tip in any reference frame moving with the steady-state velocity $V_{ss}$ with respect to the laboratory reference frame) should be taken into account when considering the evolution of various patterns and solid/liquid interfaces [37,47,48,49,50,51,52,53,54].