A Review On Computational Modelling Of Phase-transition .
rsta.royalsocietypublishing.orgA review on computationalmodelling of phase-transitionproblemsHector Gomez, Miguel Bures and AdrianReviewMoureSchool of Mechanical Engineering, Purdue University,Article submitted to journalSubject Areas:computational mechanics,computational physics, computationalmodelling and simulationKeywords:computational methods, phasetransition, phase field modelling585 Purdue Mall, West Lafayette, 47907 IN, USAPhase-transition problems are ubiquitous in scienceand engineering. They have been widely studiedvia theory, experiments and computations. Thispaper reviews the main challenges associated tocomputational modelling of phase-transition problems,addressing both model development and numericaldiscretization of the resulting equations. We focus onclassical phase-transition problems, including liquidsolid, gas-liquid and solid-solid transformations. Ourreview has a strong emphasis on the treatment ofinterfacial phenomena and the phase-field method.Author for correspondence:Hector Gomeze-mail: hectorgomez@purdue.edu1. IntroductionPerhaps the most classical example of a phase transitionis a change between the gaseous, liquid and solid phasesof a single component, for example, ice melting to liquidwater, steam condensing into liquid water or dry icesubliming to gaseous CO2 . The term phase transitionis also used more broadly to describe phenomena thatare at a first glance very different from the classicalexample of, e.g., liquid-to-gas transformations, but sharemany important properties. For example, we usuallyconsider a phase transition the process whereby amagnet loses its magnetic properties upon a temperatureincrease. For most phase transformations problems, theappearance of interfaces is inherent to the problem. Froma modelling point of view, the presence of interfacestransforms the problem into a free-boundary problem.Free-boundary problems are very difficult to treatmathematically and computationally. Here, we reviewdifferent approaches to computational modelling ofphase-transition problems. By computational modelling,we refer to both model development and numericaldiscretization of the resulting equations. Due to spacelimitations, we focus on specific examples of liquid-solid,cThe Authors. Published by the Royal Society under the terms of theCreative Commons Attribution License http://creativecommons.org/licenses/by/4.0/, which permits unrestricted use, provided the original author andsource are credited.
2. The phase-field approach(a) Mean curvature flowLet us consider a surface Γ (t) that evolves in time. We say that the surface evolves bymean curvature flow if the normal velocity of a point on the surface is proportional to themean curvature of the surface at that point. Computing the time evolution of Γ requires aparameterization of the surface which changes in time. This may initially look like a very simpleproblem to treat, but even in this case, the situation is quite challenging because it is known thatmean curvature flow develops singularities and these are usually unavoidable [1]. For example,for an initial closed surface, the flow contracts the surface until it collapses to a point. Anotherpotential complication is that the topology of the surface can change through the evolution. Thesetwo examples illustrate the challenge of treating mean curvature flow computationally by usingan explicit parametric representation of the surface. Mean curvature flow may be understood asa very simple phase transition problem. More physically relevant phase transition problems havea higher level of complexity and may involve the solution of partial differential equations (PDEs)in the interior and exterior of the surface. Because the surface evolves in time, the numericalapproximation of the PDEs requires geometrically flexible computational methods (e.g., the finiteelement method) and a procedure to update the mesh. Interface motion is usually not onlycontrolled by the geometry of the surface (like in mean curvature flow), but also by the valueof the unknowns to be solved in the interior and exterior of the surface. For these reasons,an alternative formulation that does not require an explicit parameterization of the surface isdesirable. There are multiple approaches that permit to avoid the surface parameterization [2].In our opinion, the two most successful are the level set approach [3–5] and the phase-fieldmethod [6–9]. Both formulations employ an auxiliary field φ(x, t), which is a function of spaceand time. The surface Γ (t) is defined as Γ (t) {x φ(x, t) 0}. Level set methods usually defineφ as φ(x, t) d(x, Γ (t)), where d(x, Γ (t)) represents the signed distance from x to Γ at timet. The phase-field however, defines φ as a hyperbolic tangent function of d, such that method,d(x,Γ (t)) , where is a length scale. The phase field may be understood as aφ(x, t) tanh2 diffuse interface between the phases that replaces the original sharp interface. The length scale isa measure of the thickness of the diffuse interface. As 0, we recover the sharp-interface limit.We favor the phase-field approach because, in addition to the purely geometric interpretation thatwe have just described, it may also be understood as a generalized approach to thermomechanicsthat permits to derive mathematical models for interface problems by using the Coleman-Nollapproach and classical balance laws for mass, linear momentum, angular momentum and energy[6]. This has led to an enormous number of applications of the phase-field method to materialscience [10,11], solid mechanics [12–15], fluid mechanics [7,16–19], biomechanics [20–24] andinterface problems in general [25–28].Next, we show how to derive a phase-field formulation of mean curvature flow [29,30]. Themathematical formulation of mean curvature flow can be stated as follows: Let us assume that thetime evolving surface Γ (t) R3 is defined parametrically as Γ (t) {xΓ (p, t) p Rp } where prepresents two parametric coordinates defined in suitable ranges Rp . Mean curvature flow can bedefined as xΓ(p, t) γκ,(2.1) twhere κ is the curvature vector of the surface and γ is a positive constant. To derive a phase-fieldformulation of mean curvature flow, we introduce the phase-field d(x, Γ (t)) φ(x, t) f, where f (z) tanh(z).(2.2)2 2rsta.royalsocietypublishing.org Phil. Trans. R. Soc. A 0000000.liquid-gas and solid-solid transformations. This allows us to cover classical examples such assolidification of a pure material, vaporization and condensation, martensitic transformations andgrain growth. Our review has a strong emphasis on the phase-field method.
Next, by using the identities f φ, f 0 1 f 2 and f 00 2f f 0 , we obtain, from Eq. (2.3), 2 12φ2 d κ d n φ,and φ φ .m21 φ21 φ22(1 φ2 )(2.4)Let us define the function φ̃(t) φ(xΓ (·, t), t), which is identically zero according to the definitionof the phase field; see Eq. (2.2). Using the chain rule we can write0 x φdφ̃ φ Γ .dt t t(2.5)The motion of an interface is determined by its normal velocity, which in this case, can becomputed as xΓ2 xΓ2 φ·n · φ ,(2.6) t1 φ2 t1 φ2 twhere we used Eqs. (2.4) and (2.5). Then, substituting Eqs. (2.4) and (2.6) into (2.1), we obtain thephase-field formulation of mean curvature flow φ γ2φ2 φ φ .(2.7) t21 φ2To write Eq. (2.7) as the classical Allen-Cahn equation, we need to introduce the function2W (φ) 41 (1 φ2 )2 . From the relation f 0 1 f 2 , it may be shown that 2 φ 2 14 (1 φ2 )2 .This implies thatbecomes2φ φ 21 φ2 W 0 (φ) 2and the phase-field formulation of mean curvature flow φ γ t2 W 0 (φ) φ , 2(2.8)which corresponds to the classical Allen-Cahn equation [31,32].(b) Thermomechanics approach to Allen-CahnThe Allen-Cahn equation can be derived as an approximation to mean curvature flow, butalso from fundamental balance laws and the Coleman-Noll approach [33,34]. The ColemanNoll method starts by postulating a free energy functional that is enforced to decrease as thesolution evolves in time. This is a rational approach to impose restrictions to the constitutiverelations. The energy functional may be expressed in terms of the free energy density Ψ asRF[φ] Ω Ψ dx, where Ω denotes the spatial domain. The free-energy density Ψ is assumed tobelong to the constitutive class Ψ Ψb(φ, φ). The energy-dissipation property that we postulatecan be expressed as Ḟ W(Ω) D(Ω), where Ḟ denotes the time derivative of F. In addition,W(Ω) is referred to as the working term and accounts for external forces or energy supplies comingthrough the boundary of Ω. The dissipation term, D(Ω) must be non-negative for all conceivableprocesses. The Allen-Cahn equation represents non-conservative phase dynamics. Therefore, themass balance equation may be expressed as φ R, t(2.9)where the function R is determined to ensure energydissipation. The constitutive class of R is b φ, µ), where µ Ψb · Ψb is the variational derivative of F withgiven by R R(φ, φ φ3rsta.royalsocietypublishing.org Phil. Trans. R. Soc. A 0000000.We adopt the sign convention that d is positive outside Γ (t) and negative inside. The unit outwardnormal vector to Γ (t) can be defined as n d and the mean curvature as κm d2 . The spatialderivatives of the phase field can be computed as d d 1 φ 2φ1 d d1 2d , f0 2 f 00 f0.(2.3) xi xi xj xi xj2 2 xi 2 2 xi xj
4(2.10)where Ω is the boundary of Ω and n denotes the unit outward normal to Ω. In Eq. (2.10),Rwe identify D(Ω) Ω µR dx and the remaining term on the right-hand side as W(Ω). Toachieve energy dissipation, we take R m(φ)µ, with m(φ) 0. If we define Ψ using the classicalGinzburg-Landau energy, that is, Ψ W (φ) 2 φ 2 , then substituting R m(φ)µ in Eq. (2.9),we get the canonical Allen-Cahn equation φ m(φ) W 0 (φ) 2 φ .(2.11) tTaking m(φ) γ,2 2we recover the phase-field form of mean curvature flow; see Eq. (2.8).3. Liquid-solid phase transformationsLiquid-solid phase transformations usually refer to solidification and melting. This type of phasetransformations have been widely studied in the literature [28,35,36].(a) Generalized Stefan problemThe generalized Stefan problem considers a solid-liquid system in the open, spatial domain Ω.The spatial domain Ω can be decomposed as Ω Ωs Ω , where Ω and Ωs denote, respectively,the liquid and solid subdomains. Due to the phase transformations, both Ω and Ωs change withtime and their time evolution is actually one of the unknowns of the problem. The solid-liquidinterface Γ s is defined as Γ s Ω Ωs . The generalized Stefan problem can be expressed as e ·q 0 tρ Vn [[q]] · n sinΩs Ω ,(3.1)onΓ s ,(3.2)[[θ]] 0onΓ s ,(3.3)θm θρH σ(2κm ωVn )θmonΓ s ,(3.4)ρwhere ρ denotes density, e Cv θ χ is internal energy per unit mass, θ is the temperature,Cv is heat capacity per unit mass, is specific latent heat (energy per unit mass) and χ is acharacteristic function of the liquid phase, i.e., a discontinuous function that takes the values 1in the liquid phase and 0 in the solid phase. We use Fourier equation, so that q k θ where kis the thermal conductivity (energy per unit length, per unit time and per unit temperature). Vnis the normal velocity of the interface (positive if directed towards the liquid), [[q]] is the heat fluxjump across the interface with the sign convention that [[f ]] fs f for any function quantitythat is discontinuous across the interface, n s is the unit normal to Γ s pointing towards theliquid, θm is the melting temperature, H represents the interfacial enthalpy per unit mass, ω is thekinetic undercooling coefficient, κm is the mean curvature of the interface (positive for sphericalsolid crystals) and σ is the surface tension. By following the same procedure employed to derivea phase-field formulation of mean curvature flow we can obtain a phase-field approximation tothe generalized Stefan problem that reads asρ Cv θ φ ρh0 (φ) · k(φ) θ ; t tωW 0 (φ)ρH φθ θm φ G(φ). tθm 22 σ(3.5)Here, φ is a phase field that transitions smoothly from 1 in the solid to 1 in the liquid, h isa monotone interpolatory function that verifies h( 1) 1, h( 1) 0, e.g., h(φ) 12 (1 φ) andW is the double-well potential introduced in Sect. 2. The thermal conductivity is a functionof the phase-field to account for a potentially different thermal conductivity of the liquid andsolid phases. We take a function that satisfies k( 1) k and k( 1) ks , for example, k(φ) 112 (1 φ)k 2 (1 φ)ks . The function G(φ) vanishes in the pure phases. There are multiplersta.royalsocietypublishing.org Phil. Trans. R. Soc. A 0000000.respect to φ. By using the definition of the free energy F and the chain rule, we obtainZZZd Ψb φḞ Ψ dx µR dx ·nda,dt Ω φ tΩ Ω
1.050.5-0.5-1.0Figure 1: Crystal growth on the two-dimensional domain Ω [0, L]2 with L 0.1 cm. We used auniform mesh with 5122 C 1 -quadratic elements. We employed the parameters ρ 8.91 g cm 3 ,Cv 0.6083 J K 1 g 1 , 274.98 J g 1 , ks k 0.8401 J K 1 cm 1 s 1 , ω 130 s cm 2 , 2 · 10 4 cm, H 4.751 J g 1 , σ 3.7 · 10 5 J cm 2 , θm 1728 K, δ 0.05, α0 π/4 and q 4.possibilities for G(φ). Depending on the functional form that we choose for G(φ), the phase-fieldformulation will converge faster or slower to the generalized Stefan problem. Common choices inthe literature are G(φ) 1 φ2 and G(φ) (1 φ2 )2 . Here, we will use the second one. To attaingood agreement between simulations and experiments, it is common to introduce anisotropy inthe material surface tension by assuming that σ depends on the unit normal to the liquid-solidinterface. We use the expression σ/σ 1 δ cos[q(α α0 )], where σ is the mean value of σ, δ isthe strength of the anisotropy, q is the mode number and α0 is the preferred angle. The angle ofthe normal to the surface, α, can be easily defined from the phase field. Fig. 1 shows the results ofa simulation performed using isogeometric analysis [37]. The computation shows how a crystallocated in an undercooled region grows creating a dendritic pattern.(b) Free-energy approach to solidification(i) Wang’s solidification modelPhase-field models of solidification can also be derived using thermomechanics [38,39]. Theinternal energy and entropy of any subvolume V of the system of interest can be defined as ZZ 2E ρ e dx, S s φ 2 dx,(3.6)2VVwhere e is the specific internal energy and s(e, φ) is an entropy density. The key difference withclassical thermodynamics is that S depends not only on φ, but also on its gradient. By applyingthe first and second laws of thermodynamics we obtainZZ q 2 φ̇ φ n da 0.(3.7)Ė q · n da 0, Ṡ V V θwhere V is the boundary of V. By defining a suitable expression for the specific internal energye(θ, φ), assuming that the latent heat is a constant and using the Coleman-Noll approach, oneobtains the following model W 0 (φ) θθm 0 φ φ11ρ Cv ρh0 (φ) · k(φ) θ ; ω φ h(φ) . (3.8) t t tθmθ 2H 2This model is identical to Eqs. (3.5), except for the last term in the phase-field equation. To ourknowledge, this model converges to the sharp interface solidification theory as 0 only whenθ θm . However, the functional form of the last term in the phase-field equation depends on theexpression of e(θ, φ). Therefore, a judicious choice of e(θ, φ) may lead to a model that converges,for all temperatures, to the expected sharp-interface limit as 0. More details on this may befound in [40]. To our knowledge, deriving thermodynamically consistent solidification modelsthat converge to the generalized Stefan problem as 0 in the anisotropic case is an openproblem.rsta.royalsocietypublishing.org Phil. Trans. R. Soc. A 0000000.0.0
(ii) Phase-field crystal approach to atomic scale solidification6where ρ represents a local atomistic density field and Φ(ρ) 2r ρ2 14 ρ4 . Here, r is the so-calledundercooling coefficient [42]. By impossing mass conservation and free energy dissipation, onecan derive the equation (see [43] for more details) ρ Φ0 (ρ) ρ 2 ρ 2 ρ .(3.10) tThe phase-field crystal equation has been widely studied in the computational physics [44–47]and condensed matter physics literature [48–50]. Recently, the phase-field crystal model has beengeneralized to account for faster dynamics [51–54].For the computational study of the phase-field crystal equation, instead of resorting to classicalfinite element methods, we will illustrate here how very fast and simple algorithms can be derivedusing first-order accurate semi-implicit time integration schemes and fast Poisson solvers forthe spatial dicretization. We begin by performing time discretization. Rearranging Eq. (3.10) andapplying standard finite differences to the time derivative, we getρn 1 ρn (1 r) ρn 1 (ρ3n ) 2 2 ρn 3 ρn 1 , t(3.11)where ρn is the time discrete approximation to ρ(x, tn ) and tn n t for n 0. Note that some ofthe terms on the right-hand side of Eq. (3.11) have been approximated at tn and some others attn 1 . We decided whether a linear term is approximated at tn (explicitly) or at tn 1 (implicitly)based on stability considerations. Linear terms that contribute to increase the L2 norm wereapproximated explicitly. This is a common approach to attain stability in time discretizationschemes; see [6]. The cubic, nonlinear term was approximated explicitly to keep the algorithmlinear. Basic manipulations allow us to rewrite Eq. (3.11) ashiI d t(1 r) t 3 ρn 1 ρn 2 t 2 ρn t (ρ3n ),(3.12)where I d is the identity operator. What remains to be done at this point is to discretize inspace, which in view of Eq. (3.12) boils down to constructing a discrete approximation ofthe Laplace operator. Given this situation, spatial discretization can be accomplished usinga fast Poisson solver [55], which may be simply understood as a very fast implementationof a second-order finite difference method on a uniform mesh. For example, for a Poissonequation u f 0 in two dimensions, the discretized equation at an interior node (i, j) is4ui,j ui,j 1 ui 1,j ui 1,j ui,j 1 bij , where ui,j u(ih, jh), h is the mesh size andbij h2 f (ih, jh). The indices i and j go from 1 to N . After collecting the equations for allthe interior nodes and boundary conditions, we obtain a linear system of equations given byKu b. There are multiple ways to solve linear systems of equations, but they very rarelyinvolve the use of the matrix eigenvalues. However, for the matrix K that arises from centralrsta.royalsocietypublishing.org Phil. Trans. R. Soc. A 0000000.The solidification models described up to now represent continuum descriptions of the physicalprocess. Continuum models permit simulating long times and large systems, but sometimesfall short of incorporating fundamental physical phenomena that occurs at very small scales.Molecular dynamics [41] is an alternative approach to continuum modelling that may accountfor phase-transition physics at atomic length scales. The disadvantage of molecular dynamicsis that simulations are restricted to very small systems and microsecond time periods at most.Recently, the phase-field crystal equation has been proposed as a model to describe two-phasesystems on atomic length scales and diffusive time scales. The use of diffusive time scales allowsto study much longer time intervals than those reachable with molecular dynamics simulations.The phase-field crystal equation has been employed to model multiple physical phenomena,including c
computational physics, computational modelling and simulation Keywords: computational methods, phase transition, phase field modelling Author for correspondence: Hector Gomez . approach and classical balance laws for mass, linear momentum, angular momentum and energy [6]. This has led to an enormous number of applications of the phase-field .
Computational fluid dynamics modelling of the air movement in an environmental test chamber with a respiring manikin Article Accepted Version Mahyuddin, N., Awbi, H. B. and Essah, E. A. (2015) Computational fluid dynamics modelling of the air movement in
Technology has been looking at UK computational modelling capability and how it could be better leveraged in both the public and private sector. Our aim for this report is to demystify computational modelling, to demonstrate our capabilities, and to consider steps which could be taken to f
theoretical framework for computational dynamics. It allows applications to meet the broad range of computational modeling needs coherently and with fast, structure-based computational algorithms. The paper describes the SOA computational ar-chitecture, the DARTS computational dynamics software, and appl
and simplified method to describe masonry vaults in global seismic analyses of buildings. Fig. 1 summarizes three different modelling techniques for ma sonry modelling, respectively, mi cro- , macro- and simplified micro modelling. In the case a micro modelling approach is take n, the challenge is to describe the complex behavior of the
Agile Modelling is a concept invented in 1999 by Scott Ambler as a supplement to Extreme Pro-gramming (XP) [Source: Agile Modelling Values]. Strictly defined, Agile Modelling (AM) is a chaordic, practices-based methodology for effective modelling and documentation [Source: Interview with SA by Clay Shannon].
equately support part modelling, i.e. modelling of product elements that are manufactured in one piece. Modelling is here based on requirements from part-oriented applica-tions, such as a minimal width for a slot in order to be able to manufacture it. Part modelling systems have evolved for some time now, and different modelling concepts have
5. Who can grow the largest crystal from solution? Modelling crystals 15 . 1. Modelling a salt crystal using marshmallows 2. Modelling crystals using cardboard shapes 3. Modelling diamond and graphite 4. Modelling crystal growth using people. More about crystals 21 . 1. Crystalline or plastic? 2. Make a crystal garden. Putting crystals to use .
I am My Brother’s Keeper (2004) As our New Year’s celebration draws near, I once again find myself pondering the enigmatic story that our tradition places before us at this time—the story of the Binding of Isaac. Once again, I walk for those three long days with father Abraham and ponder the meaning of his journey with his son to the mountain. And once again, I find fresh meaning in the .
