Phase-field Equations Research Articles

Consistent, energy-conserving momentum transport for simulations of two-phase flows using the phase field equations

Realistic two-phase flow problems of interest often involve high Re flows with high density ratios. Accurate and robust simulation of such problems requires special treatments. In this work, we present a consistent, kinetic energy conserving momentum transport scheme in the context of a second order conservative phase field method. This is achieved by (1) accounting for the mass flux associated with the right-hand-side of the phase field equation in the convective flux of the conservative form of the momentum transport equation—a correction absent in previous phase field simulations (2) utilization of non-dissipative spatial discretization. We demonstrate accuracy and robustness improvements from our proposed scheme using numerical tests, including a realistic jet in cross-flow attaining turbulent conditions. Our proposed modifications to the momentum transport equation can be readily extended to other conservative phase field methods.

Error estimates for fully discrete BDF finite element approximations of the Allen–Cahn equation

Abstract For a class of compatible profiles of initial data describing bulk phase regions separated by transition zones, we approximate the Cauchy problem of the Allen–Cahn (AC) phase field equation by an initial-boundary value problem in a bounded domain with the Dirichlet boundary condition. The initial-boundary value problem is discretized in time by the backward difference formulae (BDF) of order $1\leqslant q\leqslant 5$ and in space by the Galerkin finite element method of polynomial degree $r-1$, with $r\geqslant 2$. We establish an error estimate of $O(\tau ^q\varepsilon ^{-q-\frac 12}+h^{r}\varepsilon ^{-r-\frac 12}+{e}^{-c/\varepsilon })$ with explicit dependence on the small parameter $\varepsilon$ describing the thickness of the phase transition layer. The analysis utilizes the maximum-norm stability of BDF and finite element methods with respect to the boundary data, the discrete maximal $L^p$-regularity of BDF methods for parabolic equations and the Nevanlinna–Odeh multiplier technique combined with a time-dependent inner product motivated by a spectrum estimate of the linearized AC operator.

A combined XFEM phase-field computational model for crack growth without remeshing

This paper presents an adaptive strategy for phase-field simulations with transition to fracture. The phase-field equations are solved only in small subdomains around crack tips to determine propagation, while an extended finite element method (XFEM) discretization is used in the rest of the domain to represent sharp cracks, enabling to use a coarser discretization and therefore reducing the computational cost. Crack-tip subdomains move as cracks propagate in a fully automatic process. The same mesh is used during all the simulation, with an h-refined approximation in the elements in the crack-tip subdomains. Continuity of the displacement between the refined subdomains and the XFEM region is imposed in weak form via Nitsche’s method. The robustness of the strategy is shown for some numerical examples in 2D and 3D, including branching and coalescence tests.

Adaptive higher-order phase-field modeling of anisotropic brittle fracture in 3D polycrystalline materials

Due to its ability to simulate complex microstructure evolution, the phase-field modeling has been extensively developed to investigate brittle fracture in recent years. However, low computational efficiency still imposes substantial difficulties in the development of phase-field models. In this work, we develop a novel adaptive phase-field approach based on the isogeometric meshfree collocation method (IMCM) to simulate the crack propagation in 2D and 3D polycrystalline materials. The concept of IMCM is based upon the correspondence between the isogeometric collocation and reproducing kernel meshfree method to facilitate a robust mesh adaptivity in isogeometric collocation. The strong form collocation formulation further enhances the computational efficiency of phase-field modeling by reducing the number of point evaluations. The present numerical framework is utilized for the adaptive phase-field modeling which introduces the anisotropy of fracture resistance for each grain in polycrystals. Furthermore, the discrete displacement and phase-field equations are generalized to enable the calculation of both second- and fourth-order gradients, which are required to solve the phase-field models using IMCM. The smoothness and higher-order continuity of IMCM enable the fourth-order phase-field equation to be solved directly without splitting it into two second-order differential equations. The fourth-order model can capture the crack surface accurately with fewer nodes than the second-order model. Several numerical examples of isotropic and anisotropic brittle fracture in polycrystalline materials are investigated to demonstrate the effectiveness and robustness of the proposed approach.

Unified simplified multiphase lattice Boltzmann method for ferrofluid flows and its application

A unified simplified multiphase lattice Boltzmann method (USMLBM) is constructed in this work for simulating complex multiphase ferrofluid flows with large density and viscosity ratios. In USMLBM, the Navier–Stokes equations, the Poisson equation of the magnetic potential, and the phase-field equation are utilized as the ferrohydrodynamics behavior modeling and interface tracking algorithm. Solutions of the macroscopic governing equations are reconstructed with the lattice Boltzmann framework and resolved in a predictor–corrector scheme. Various benchmark tests demonstrate the efficiency and accuracy of USMLBM in simulating multiphase ferrofluid flows. We further adopt USMLBM to analyze in detail the mechanisms of bubble merging inside a ferrofluid under a uniform external magnetic field. The numerical results indicate that the bubbles tend to move toward each other and further merge together, even for a large initial separation between the bubbles. Due to complex interaction between the bubbles and the ferrofluid during the magnetophoretic acceleration process, the nonlinear effect on bubble merging is observed when the initial separation increases. Moreover, at a larger initial separation, the shape of bubbles seems to be not sensitive to the initial separation.

Development of phase-field model based on balance laws and thermodynamic discussion

In this work, a phase-field model for recrystallization is developed based on the conservation laws. There has been no attempt to develop a phase-field model of recrystallization based on the conservation laws, even though various phase-field simulation models to reproduce the recrystallization phenomenon have been proposed. However, it is unclear what conservation laws are required for such a model. In the previous paper, toward solving this problem, we developed conservation laws of mass, momentum, angular momentum, and energy and a law of entropy at the lattice scale for the process of recrystallization. In this paper, first, two continuous variables, i.e., the order parameter and crystal orientation, are introduced into the balance equation of mass for a single phase and that of angular momentum for the lattice, respectively. Next, the fluxes of the order parameter and crystal orientation are derived from the law of entropy by the use of rational thermodynamics. Moreover, the diffusion coefficient and mass source are modeled to derive the evolution equations, i.e., phase-field equations of the order parameter and crystal orientation. Finally, for the phase-field equation of the crystal orientation, neglecting the conservative part and integrating the equation with respect to time under the first-order approximation, a phase-field model that is used for stable calculations is developed. This work aims to develop a phase-field theory on the basis of the change in crystal lattice during recrystallization. This paper gives a physical background to the methodological phase-field approach in the case of recrystallization.

Explicit nonlinear finite element approach to the Lagrangian-based coupled phase field and elasticity equations for nanoscale thermal- and stress-induced martensitic transformations

In this paper, a nonlinear finite element procedure is developed to solve the coupled system of phase field and elasticity equations at large strains for martensitic phase transformations at the nanoscale. The transformation is defined based on an order parameter which varies between 0 for austenite to one for martensite. The phase field equation relates the rate of change of the order parameter to the Helmholtz free energy. Both the Helmholtz free energy and the transformational strain tensor depend on the order parameter, and the coupling between the phase field and elasticity equations occurs due to the presence of elastic energy in the Helmholtz free energy and the transformational strain in the total strain. The deformation gradient tensor is introduced as a multiplicative decomposition of elastic and transformational gradient tensors. A staggered strategy is used to solve the coupled system of equations so that first, the principle of virtual work is utilized to obtain the integral form of the Lagrangian equation of motion. It is then discretized and solved using the Newton–Raphson method which gives the displacements and consequently, the deformation gradient tensor. Next, since the order parameter is obtained from the phase field equation from the previous time step, the elastic deformation gradient tensor and the first Piola–Kirchhoff stress can be calculated and substituted in the phase field equation to find the order parameter for the current time step. The weighted residuals method is used to derive the finite element form of the phase field equation and the explicit method is used for its time discretization. The finite element procedure is well verified by comparing the obtained results with those from other simulations. Examples of thermal- and stress-induced phase transformations are presented, including austenite–martensite interface propagation and growth of preexisting martensitic regions. The developed algorithm and code can be advanced to solve kinetic problems coupled with mechanics at large strains for phase transformations and similar phenomena at the nanoscale.

A framework for efficient isogeometric computations of phase-field brittle fracture in multipatch shell structures

We present a computational framework for applying the phase-field approach to brittle fracture efficiently to complex shell structures. The momentum and phase-field equations are solved in a staggered scheme using isogeometric Kirchhoff–Love shell analysis for the structural part and isogeometric second- and fourth-order phase-field formulations for the brittle fracture part. For the application to complex multipatch structures, we propose penalty formulations for imposing all the required interface constraints, i.e., displacement (C0) and rotational (C1) continuity for the structure as well as C0 and C1 continuity for the phase field, where the latter is required only in the case of the fourth-order phase-field model. All involved penalty terms are scaled with the corresponding problem parameters to ensure a consistent scaling of the penalty contributions to the global system of equations. As a consequence, all coupling terms are controlled by one global penalty parameter, which can be set to 103 independent of the problem parameters. Furthermore, we present a multistep predictor–corrector algorithm for adaptive local refinement with LR NURBS, which can accurately predict and refine the region around the crack even in cases where fracture fully develops in a single load step, such that rather coarse initial meshes can be used, which is essential especially for the application to large structures. Finally, we investigate and compare the numerical efficiency of loosely vs. strongly staggered solution schemes and of the second- vs. fourth-order phase-field models.

Machine learning materials physics: Multi-resolution neural networks learn the free energy and nonlinear elastic response of evolving microstructures

Many important multi-component crystalline solids undergo mechanochemical spinodal decomposition: a phase transformation in which the compositional redistribution is coupled with structural changes of the crystal, resulting in dynamically evolving microstructures. The ability to rapidly compute the macroscopic behavior based on these detailed microstructures is of paramount importance for accelerating material discovery and design. Here, our focus is on the macroscopic, nonlinear elastic response of materials harboring microstructure. Because of the diversity of microstructural patterns that can form, there is interest in taking a purely computational approach to predicting their macroscopic response. However, the evaluation of macroscopic, nonlinear elastic properties purely based on direct numerical simulations (DNS) is computationally very expensive, and hence impractical for material design when a large number of microstructures need to be tested. A further complexity of a hierarchical nature arises if the elastic free energy and its variation with strain is a small-scale fluctuation on the dominant trajectory of the total free energy driven by microstructural dynamics. To address these challenges, we present a data-driven approach, which combines advanced neural network (NN) models with DNS to predict the homogenized, macroscopic, mechanical free energy and stress fields arising in a family of multi-component crystalline solids that develop microstructure. The microstructures are numerically generated by solving a coupled, mechanochemical spinodal decomposition problem governed by nonlinear strain gradient elasticity and the Cahn–Hilliard phase field equation. The hierarchical structure of the free energy’s evolution induces a multi-resolution character to the machine learning paradigm: We construct knowledge-based neural networks (KBNNs) with either pre-trained fully connected deep neural networks (DNNs), or pre-trained convolutional neural networks (CNNs) that describe the dominant characteristic of the data to fully represent the hierarchically evolving free energy. We demonstrate multi-resolution learning of the materials physics; specifically of the nonlinear elastic response for both fixed and evolving microstructures.

An efficient Pseudo-spectral based phase field method for dendritic solidification

In the current study, Fourier Pseudo-spectral based dealiasing scheme Random Phase Shift Method (RPSM) using the RK4 temporal scheme is implemented for simulating some benchmark problems of Phase Field based pure metal solidification. This high-resolution scheme is capable of resolving the thin solid-liquid interface accurately without using Adaptive Mesh Refinement (AMR) algorithm even using just adequate number of grid points. Both the conservative and non-conservative forms of the Phase Field equations are studied. It is found that the conservative RPSM scheme gives better computational performances as indicated by the lower values of average CPU time and different error norms calculated for the Phase Field variable. The RPSM based Phase Field model is applied for simulating dendritic solidification and it is found that dendrite shape and dimensionless tip velocity agree well with the existing literature. The RK4 based RPSM scheme is also compared with a 10th order exponential Fourier smoothing Filter in the evolution of four arm dendritic structure and it is observed that the conservative Phase Field equations with filtering scheme display numerical instabilities at later simulation time, while both forms of the RPSM scheme show stable and accurate solutions. Subsequently, the capability of the RPSM based Phase Field model in simulating microstructure evolution consisting of multiple dendrites with random location and orientation is demonstrated. The present study shows that the RPSM based Phase Field model computes dendritic growth accurately and efficiently.

Study of phase-field lattice Boltzmann models based on the conservative Allen-Cahn equation.

Conservative phase-field (CPF) equations based on the Allen-Cahn model for interface tracking in multiphase flows have become more popular in recent years, especially in the lattice-Boltzmann (LB) community. This is largely due to their simplicity and improved efficiency and accuracy over their Cahn-Hilliard-based counterparts. Additionally, the improved locality of the resulting LB equation (LBE) for the CPF models makes them more ideal candidates for LB simulation of multiphase flows on nonuniform grids, particularly within an adaptive-mesh refinement framework and massively parallel implementation. In this regard, some modifications-intended as improvements-have been made to the original CPF-LBE proposed by Geier etal. [Phys. Rev. E 91, 063309 (2015)PLEEE81539-375510.1103/PhysRevE.91.063309] which require further examination. The goal of the present study is to conduct a comparative investigation into the differences between the original CPF model proposed by Geier etal. [Phys.Rev. E 91, 063309 (2015)PLEEE81539-375510.1103/PhysRevE.91.063309] and the so-called improvements proposed by Ren etal. [Phys. Rev. E 94, 023311 (2016)2470-004510.1103/PhysRevE.94.023311] and Wang etal. [Phys. Rev.E 94, 033304 (2016)2470-004510.1103/PhysRevE.94.033304]. Using the Chapman-Enskog analysis, we provide a detailed derivation of the governing equations in each model and then examine the efficacy of the above-mentioned models for some benchmark problems. Several test cases have been designed to study different configurations ranging from basic yet informative flows to more complex flow fields, and the results are compared with finite-difference simulations. Furthermore, as a development of the previously proposed CPF-LBE model, axisymmetric formulations for the proposed model by Geier etal. [Phys. Rev. E 91, 063309 (2015)PLEEE81539-375510.1103/PhysRevE.91.063309] are derived and presented. Finally, two benchmark problems are designed to compare the proposed axisymmetric model with the analytical solution and previous work. We find that the accuracy of the model for interface tracking is roughly similar for different models at high viscosity ratios, high density ratios, and relatively high Reynolds numbers, while the original CFP-LBE without the additional time-dependent terms outperforms the so-called improved models in terms of efficiency, particularly on distributed parallel machines.

Generalized conservative phase field model and its lattice Boltzmann scheme for multicomponent multiphase flows

In this paper a generalized conservative phase field-lattice Boltzmann model is proposed for immiscible incompressible flows with any number of components. In this phase field model, the each phase field variable satisfies a generalized eikonal equation and the volume fraction constraint. By introducing a quadratic functional, the dynamics of phase field variables can be treated as a gradient flow method. The phase field equations are reformulated based on the mass conservative law and the reduction-consistency conditions. Moreover, based on the generalized free energy functional, a family of generalized continuous surface tension force formulations is deduced by using the virtual work method. A lattice Boltzmann model (LBM) is developed to solve the generalized conservative phase field equations and the hydrodynamics equations. In the numerical simulation section, the motion of circular interfaces and Poiseuille flow with four phases are simulated first. The convergence rate and reduction-consistent of present model for the phase field equations are validated. The computational accuracy of present model for the hydrodynamics equations is also validated. Then serval two-dimensional (2D) examples with pairwise surface tension effect, including the droplets immersed in another fluid with four fluid phases, spreading of liquid lenses and spinodal decomposition with four fluid phases, are simulated. The obtained numerical results are in good agreement with the analytical solutions. Four continuous surface tension force formulations are also comparative studied. Finally, as a practical application, a three-dimensional (3D) bubble rising in a three-layer liquid is investigated.

Influence of Interfacial Stress on Microstructural Evolution in NiAl Alloys

A phase-field model for the phase transition between austenite and martensite and twinning between two martensitic variants is presented from our previous theory [1] with the main focus on the influence of interfacial stress that is consistent with the sharp interface limit. Each variant-variant transformation can be represented by only one order parameter. Thus, it allows us to get the analytical solution of interface energy and width. Coupled phase-field and elasticity equations are solved for cubic-to-tetragonal phase transformation in NiAl shape memory alloy. The effects of interfacial stress are studied for martensite-martensite interfaces in detail, which was absent in [1]. Additionally, stress and temperature-induced growth of the martensitic phase inside austenite and twining are simulated. Some of the nontrivial experimentally observed microstructures reproduced in the simulations [1] are analyzed in detail. It includes tip splitting and bending, and twins crossing. This theory can be extended for electric, reconstructive, and magnetic phase transformations.

Isogeometric continuity constraints for multi-patch shells governed by fourth-order deformation and phase field models

This work presents numerical techniques to enforce continuity constraints on multi-patch surfaces for three distinct problem classes. The first involves structural analysis of thin shells that are described by general Kirchhoff–Love kinematics. Their governing equation is a vector-valued, fourth-order, nonlinear, partial differential equation (PDE) that requires at least C1-continuity within a displacement-based finite element formulation. The second class are surface phase separations modeled by a phase field. Their governing equation is the Cahn–Hilliard equation – a scalar, fourth-order, nonlinear PDE – that can be coupled to the thin shell PDE. The third class are brittle fracture processes modeled by a phase field approach. In this work, these are described by a scalar, fourth-order, nonlinear PDE that is similar to the Cahn–Hilliard equation and is also coupled to the thin shell PDE. Using a direct finite element discretization, the two phase field equations also require at least a C1-continuous formulation. Isogeometric surface discretizations – often composed of multiple patches – thus require constraints that enforce the C1-continuity of displacement and phase field. For this, two numerical strategies are presented: A Lagrange multiplier formulation and a penalty method. The curvilinear shell model including the geometrical constraints is taken from Duong et al. (2017) and it is extended to model the coupled phase field problems on thin shells of Zimmermann et al. (2019) and Paul et al. (2020) on multi-patches. Their accuracy and convergence are illustrated by several numerical examples considering deforming shells, phase separations on evolving surfaces, and dynamic brittle fracture of thin shells.

An Efficient Implementation of Phase Field Method with Explicit Time Integration

The phase field method integrates the Griffith theory and damage mechanics approach to predict crack initiation, propagation, and branching within one framework. No crack tracking topology is needed, and complex crack shapes can be captures without user intervention. In this paper, a detailed description of how the phase field method is implemented with explicit dynamics into LS-DYNA is provided. The displacement field and the damage field are solved in a staggered approach and the phase field equation is solved every Nth time step (N is refered to as calculation cycle) to save computational time. An N value smaller than 1/400 of the total time step numbers is suggested. Several simulations are presented to demonstrate the feasibility of this solving scheme.

Effects of interfacial stress in phase field approach for martensitic phase transformation in NiAl shape memory alloys

A phase field approach for phase transition between austenite to martensitic variants and twinning between martensitic variants is presented with the main focus on the effects of interfacial stress on phase transformations. In this theory, each variant-variant phase transformations and twinnings within each martensitic variant can be represented by only one-order parameter. Thus, it allows us to get the analytical solution of the martensite-martensite interface profile, energy, and width. Moreover, this model allows us to include interface stress which is consistent with the sharp interface limit. The finite element method is utilized to solve the coupled phase field and elasticity equations for a cubic to tetragonal phase transformation in NiAl shape memory alloy. The stress fields are obtained and the effects of interfacial stress on the stress field at the interfaces are studied for both austenite–martensite and martensite–martensite interfaces in detail. Additionally, the temperature-induced growth of the martensitic phase inside austenite, martensitic phase transformation, and twining are solved. The evolution of the microstructure and stress fields are obtained and the effects of interfacial stress on the morphological evolution of martensitic nanostructures are examined. It is found that the interfacial stress is an important factor, influencing the stress distribution at interfaces and the phase field solution significantly. This theory can be extended for electric, reconstructive, and magnetic phase transformations.

Efficient numerical methods for phase-field equations

In this article, we overview recent developments of numerical methodsfor phase-field equations. The main difficulty fornumerically solving phase-field equations is about a severe restrictionon the time step due to nonlinearity and high order differential terms, while it usually requires a very long timesimulation to reach the steady state. It is knownthat phase-field models satisfy a nonlinear stability relationship, called energy stability,which means that the free energy functional decays in time.It has attracted more and more attention to design numerical schemes inheriting the energy stability so thatthe numerical simulation may use large time steps and keep the accuracy.For some popularly studied phase-field equations, this article will present several widely used highly efficient numerical schemes andshow an adaptive time-stepping strategy based on the changing rate in time of the energy functional, which could guarantee the accuracy and stability of the numerical solution and improves the computational efficiency significantly.

The numerical solutions for the energy-dissipative and mass-conservative Allen–Cahn equation

There are three criteria when we model with phase-field equations. The first one is the evolution of phase separation. The second one is energy dissipation law. The last one is whether we have the law of conservation of mass or not. If we minimize the Ginzburg–Landau free-energy functional, the kinetics of phase separation inherently follows. However, the other properties such as energy dissipation and mass conservation are more difficult to model. Hitherto, there have been a few proposals to conserve mass using the Allen–Cahn (AC) equation, but we do not have the energy dissipation law in these models when considered in light of local effects of mass-correction. For the reasons, certain researchers still prefer to use the Cahn–Hilliard equation, e.g. biharmonic equation instead of using the Allen–Cahn equation when we consider a mass conserving case. This paper is concerned with numerical solutions of the AC equation in the mass-conserving space, where the solution is corrected at each iteration for mass conservation. We introduce two different methods for solving the conservative AC equation with periodic boundary conditions. In this work, we have the AC equation in common, but use two different mass correction operators. Note that both of mass correction operators satisfy the discrete energy dissipation property. First, we project the solution of the AC equation to the mass conserving space, and then recover the linear convex splitting scheme for the conventional mass-conserving equation. The linear convex splitting scheme is well-known but the convergence of numerical solutions is elaborately studied in terms of gradient-projection. Next, we propose the other energy-minimizing method on the mass-conserving space. It translates the solution of the AC equation to the mass conserving space via an energy-dissipation operator. The dissipation of discrete energy functional for the proposed method is also established. As far as the author knows, no report has been presented for the numerical scheme having discrete energy dissipation property of the mass-conserving AC equation without using the nonlocal and constant multiplier. To verify the common properties of mass-projection and proposed methods numerically, we simulate the evolutions for phase separation, dissipation of energy, and mass conservation. We moreover consider the short and long time stages of evolutions to see the difference between two methods. These several tests are provided for its applications of the proposed method for modeling natural phenomena.

Thin interface limit of the double-sided phase-field model with convection.

The thin interface limit of the phase-field model is extended to include transport via melt convection. A double-sided model (equal diffusivity in liquid and solid phases) is considered for the present analysis. For the coupling between phase-field and Navier-Stokes equations, two commonly used schemes are investigated using a matched asymptotic analysis: (i) variable viscosity and (ii) drag force model. While for the variable viscosity model, the existence of a thin interface limit can be shown up to the second order in the expansion parameter, difficulties arise in satisfying no-slip boundary condition at this order for the drag force model. Nevertheless, detailed numerical simulations in two dimensions show practically no difference in dendritic growth profiles in the presence of forced melt flow obtained for the two coupling schemes. This suggests that both approaches can be used for the purpose of numerical simulations. Simulation results are also compared to analytic theory, showing excellent agreement for weak flow. Deviations at higher fluid velocities are discussed in terms of the underlying theoretical assumptions. This article is part of the theme issue 'Patterns in soft and biological matters'.

Long-time simulation of the phase-field crystal equation using high-order energy-stable CSRK methods

The phase-field crystal (PFC) equation is derived by the gradient flow for the Swift–Hohenberg free energy functional; thus, the numerical method requires the energy of the functional to decrease. Convex Splitting Runge–Kutta (CSRK) methods can be suitably applied to achieve high-order temporal accuracy as well as unconditional energy stability and unique solvability. For the PFC equation, we prove the unconditional energy stability and unique solvability of the CSRK methods and provide one family of parameters of the second-order CSRK methods and possible examples of third-order CSRK methods. We present numerical experiments to demonstrate the accuracy and energy stability of the methods. Specifically, based on the high-order accuracy and energy stability of the CSRK method, we propose an indicator function capable of characterizing the pattern formation of the phase-field crystal model for long-time simulation.