Isogeometric analysis of the Cahn–Hilliard equation – a convergence study
Introduction
Phase-field models have become a powerful tool for the modelling of phase transformations and morphological changes in different fields of physics as well as materials and engineering science. Compared to sharp-interface models their advantage is that topological changes are avoided, since interfaces are treated in a diffuse manner which is achieved through a parameter which varies continuously. This phase-field parameter accounts for the different material phases and/or the concentration of the different components. Such an approach allows to fully capture the physics of the individual interfaces without the need to explicitly track them. Typical applications include the modelling and simulation of solidification processes, spinodal decomposition, coarsening of precipitate phases, shape memory effects, re-crystallisation, and dislocation dynamics [1], [2], [3], [4]. Phase-field models have been successfully applied to predict microstructural changes under external fields [5], [6], model tumor growth [7], [8], [9], and image impainting [10]. The phase-field approach has also been used to model crack propagation [11], [12], [13], [14], [15].
Different numerical techniques have been utilised to solve the governing equations of phase-field models [16], including finite difference, finite element, and spectral methods. A difficulty regarding the solution of phase-field models is that they typically involve spatial differential operators that are higher than second-order. Therefore, standard finite elements based on -continuous Lagrangian polynomial shape functions do not provide converging solutions when directly applied to the phase-field equation. Instead, mixed formulations [17], [18], where the governing equations are split into a coupled system of lower-order differential equations, can be used. However, this approach introduces additional unknowns and is prone to stability issues [19], [20]. Alternatively, Wells et al. [21] have combined features of continuous and discontinuous Galerkin methods, cf. [22], to use a lower-order continuity basis than is actually required by the weak form. Hence, -continuous Lagrangian basis functions can be used, while the first-order inter-element continuity is weakly enforced by Nitsche's method. Although efficient as no additional unknowns are introduced, the approach changes the assembly procedures for the global system of equations and requires a careful choice of the penalty parameter. Starting with the work of Gómez [23] a new approach to the numerical modelling of structural evolution processes has emerged, which combines phase-field models with spline-based approximations. Such an isogeometric analysis framework allows for an accurate and efficient resolution of steep gradients that can occur when higher-order derivatives are introduced. Moreover, the higher-order continuity that is provided by these approximations allows for a direct discretisation without additional degrees of freedom. Differently from earlier approaches based on Hermite elements [24], [25] it can straightforwardly be applied to two- and three-dimensional problems.
Analytic solutions for multidimensional phase-field problems are generally not available. For this reason, the comparison of different discretisations in terms of a quantitative convergence analysis is cumbersome. Approaches commonly use linearised or one-dimensional models in order to carry out accuracy tests. Alternatively, the numerical result that stems from the finest discretisation is considered as the reference solution [26]. Zhang et al. [27] have presented a convergence analysis based on the method of manufactured solutions [28]. By assuming an appropriate solution and extending the underlying differential equation by a source term, different discretisations of the same phase-field model can be compared, see also Aristotelous et al. [29] and Kay et al. [30]. To the authors' knowledge a quantitative convergence study for diffuse interface models has never been carried out within the framework of isogeometric analysis. For this reason, the Cahn–Hilliard equation is used here as a prototypical higher-order equation, and the numerical properties of different discretisations are analysed by a manufactured solution approach.
Originally derived to model spinodal decomposition of binary mixtures [31], [32], the Cahn–Hilliard equation is one of the most commonly used diffuse interface models. Taking the concentration c of one of the mixtures' components as an appropriate phase-field parameter, the governing equation can be stated as Herein, λ is the interface parameter that governs the width of the interface which separates two adjacent phases, is the concentration dependant mobility, while and designate the derivative and the spatial derivative , respectively. Eventually, represents the configurational or bulk free energy, which can be described by the generalised logarithmic potential [33] with constants A and B. The bulk free energy accounts for phase separation which dominates the early stages of the evolution process. A second term, called the gradient or the interfacial free energy, contributes to the total Ginsburg–Landau free energy and accounts for the influence of the individual interfaces on the system. While the total free energy is minimised in the course of the structural evolution, different phenomena can occur. Both quantities, and , can be used to explain these phenomena from an energetic point of view.
This paper is organised as follows: In Section 2 different finite element formulations are presented for the numerical solution of the Cahn–Hilliard problem, i.e., a mixed form and a direct approach. For the latter, approximations with at least inter-element continuity are required. In Section 3, these discretisations are compared. The numerical properties of both approaches are analysed in terms of the error norms and the convergence rates by making use of a manufactured solution. In Section 4 numerical examples for the spinodal decomposition are presented.
Section snippets
Finite element method
In this contribution, isogeometric analysis is applied to the spatial discretisation of the Cahn–Hilliard equation. The spline-based discretisations offer higher-order continuity properties and are therefore suitable for the direct discretisation of weak forms involving higher-order spatial derivatives. For the temporal discretisation, the generalised-α method [34], [35] is used. Following Gómez et al. [23], the method allows for the development of an efficient time stepping algorithm based on
Convergence analysis based on a manufactured solution
In this section the properties of the presented FE formulations are compared for different approximations in terms of a quantitative analysis of error levels and convergence rates. As no analytical solution is available, a manufactured solution [27], [30] is utilised. The general idea is to use an arbitrary function as analytic reference solution. Since this function will normally not fulfil the governing differential equations exactly, a residual has to be added to the right hand side of the
Numerical examples
In the following examples the predictions of the mixed formulation and the discretisation for the spinodal decomposition of binary systems will be compared. As the different stages of the structural evolution can be identified by the free energy and its individual parts, special attention is paid to them. Starting from an initial concentration distribution in which small perturbations promote the evolution of the system, each problem is considered until a steady state is reached.
The course of
Conclusion
The present work provides a detailed convergence analysis of the Cahn–Hilliard phase-field model. Two different formulations of the model have been implemented in a Matlab-based finite element code that provides Lagrange polynomials and B-splines. Due to the lack of an analytic reference solution, the method of manufactured solutions has been adopted. With respect to the spatial discretisation, convergence rates are obtained that match analytical error estimates available for linear problems.
Acknowledgements
The present study is funded by the German Research Foundation (DFG), Priority Programme (SPP) 1713, grant KA 3309/5-1. This support is gratefully acknowledged. M. Kästner acknowledges financial support of the German Academic Exchange Service (DAAD) during his research visit to the University of Glasgow.
References (44)
- et al.
An introduction to phase-field modeling of microstructure evolution
Calphad
(2008) - et al.
The influence of lattice strain on pearlite formation in Fe-C
Acta Mater.
(2007) - et al.
Three-dimensional multispecies nonlinear tumor growth – I: model and numerical method
J. Theor. Biol.
(2008) - et al.
A continuum phase field model for fracture
Eng. Fract. Mech.
(2010) - et al.
Physical, mathematical, and numerical derivations of the Cahn–Hilliard equation
Compos. Mater. Sci.
(2014) - et al.
A discontinuous Galerkin method for the Cahn–Hilliard equation
J. Comput. Phys.
(2006) - et al.
Continuous/discontinuous finite element approximations of fourth-order elliptic problems in structural and continuum mechanics with applications to thin beams and plates, and strain gradient elasticity
Comput. Methods Appl. Mech. Eng.
(2002) - et al.
Isogeometric analysis of the Cahn–Hilliard phase-field model
Comput. Methods Appl. Mech. Eng.
(2008) - et al.
A numerical method for the nonlinear Cahn–Hilliard equation with nonperiodic boundary conditions
Compos. Mater. Sci.
(1995) - et al.
Accurate, efficient, and (iso)geometrically flexible collocation methods for phase-field models
J. Comput. Phys.
(2014)
A quantitative comparison between and elements for solving the Cahn–Hilliard equation
J. Comput. Phys.
Computationally efficient solution to the Cahn–Hilliard equation: adaptive implicit time schemes, mesh sensitivity analysis and the 3D isoperimetric problem
J. Comput. Phys.
A generalized-α method for integrating the filtered Navier–Stokes equations with a stabilized finite element method
Comput. Methods Appl. Mech. Eng.
An improved NURBS-based isogeometric analysis with enhanced treatment of essential boundary conditions
Comput. Methods Appl. Mech. Eng.
Isogeometric shell analysis with Kirchhoff–Love elements
Comput. Methods Appl. Mech. Eng.
Isogeometric analysis and error estimates for high order partial differential equations in fluid dynamics
Comput. Fluids
Isogeometric analysis of the advective Cahn–Hilliard equation: spinodal decomposition under shear flow
J. Comput. Phys.
Phase-field models for microstructure evolution
Annu. Rev. Mater. Res.
Advances of and by phase-field modelling in condensed-matter physics
Adv. Phys.
Phase-field models in materials science
Model. Simul. Mater. Sci. Eng.
Phase-field modeling of microstructure evolutions in magnetic materials
Sci. Technol. Adv. Mater.
Numerical simulation of a thermodynamically consistent four-species tumor growth model
Int. J. Numer. Methods Biomed. Eng.
Cited by (59)
Local volume-conservation-improved diffuse interface model for simulation of Rayleigh–Plateau fluid instability
2024, Computer Physics CommunicationsAdaptive isogeometric phase-field modeling of the Cahn–Hilliard equation: Suitably graded hierarchical refinement and coarsening on multi-patch geometries
2023, Computer Methods in Applied Mechanics and EngineeringThe study of diffuse interface propagation of dynamic failure in advanced ceramics using the phase-field approach
2023, Computer Methods in Applied Mechanics and EngineeringA local meshless method for transient nonlinear problems: Preliminary investigation and application to phase-field models
2022, Computers and Mathematics with ApplicationsCitation Excerpt :The second term in the integrand of the energy functional (1), penalizes the large gradients and hence accounts for the interfacial energy [89–91]. The interfacial parameter λ governs the thickness of the interface separating adjacent phases [46,92]. The Allen-Cahn equation can be applied to describe the evolution of a nonconserved field variable, resulting in formation of a two-phase configuration in a system [93].
A spatio-temporal adaptive phase-field fracture method
2022, Computer Methods in Applied Mechanics and EngineeringA novel method to impose boundary conditions for higher-order partial differential equations
2022, Computer Methods in Applied Mechanics and Engineering