Source-linked AI summary
A second order in time, uniquely solvable, unconditionally stable numerical scheme for Cahn-Hilliard-Navier-Stokes equation
Daozhi Han, Xiaoming Wang
TL;DR
The paper addresses stable and uniquely solvable time discretization for the matched-density Cahn-Hilliard-Navier-Stokes model, where diffuse interfaces create stiffness and nonlinear advection complicates analysis. It proposes a second-order convex-splitting and pressure-projection scheme with weak coupling, and reports unconditional stability, unique solvability, conservation, energy dissipation, and second-order L2 accuracy. The method is also decoupled computationally and tested with finite elements and numerical experiments.
Problem
Diffuse-interface CHNS computation must handle severe stiffness and nonlinear advection while establishing unconditional stability and unique solvability.
Method
The paper combines second-order convex splitting for Cahn-Hilliard, pressure projection for Navier-Stokes, weak coupling, Picard decoupling, and mixed finite elements.
Results
The scheme is unconditionally stable and uniquely solvable, while numerical experiments verify conservation, energy dissipation, and second-order accuracy in the L2 norm.
Takeaways & Limitations
The method provides an efficient second-order CHNS discretization that decouples pressure, velocity, and phase-field variables while retaining unconditional stability and unique solvability.
Takeaways & Limitations
The work assumes bounded mobility rather than the potentially more physically relevant degenerate mobility, whose numerical resolution is beyond scope; pressure accuracy is only first order for the splitting method.
Abstract
from arXiv · showhide
We propose a novel second order in time numerical scheme for Cahn-Hilliard-Navier- Stokes phase field model with matched density. The scheme is based on second order convex-splitting for the Cahn-Hilliard equation and pressure-projection for the Navier-Stokes equation. We show that the scheme is mass-conservative, satisfies a modified energy law and is therefore unconditionally stable. Moreover, we prove that the scheme is uncondition- ally uniquely solvable at each time step by exploring the monotonicity associated with the scheme. Thanks to the weak coupling of the scheme, we design an efficient Picard iteration procedure to further decouple the computation of Cahn-Hilliard equation and Navier-Stokes equation. We implement the scheme by the mixed finite element method. Ample numerical experiments are performed to validate the accuracy and efficiency of the numerical scheme.
1 Introduction
The paper develops a numerical method for the matched-density Cahn-Hilliard-Navier-Stokes model, whose diffuse-interface formulation captures two-fluid interface dynamics but presents stiffness, stability, and solvability challenges. Its proposed second-order scheme combines convex splitting and pressure projection, with weak coupling enabling efficient computation.
- 1 Introduction: The CHNS phase-field model describes interface dynamics in binary incompressible, macroscopically immiscible Newtonian fluids with matched density and viscosity.The diffuse-interface formulation represents the interface as a thin transition layer of width ϵ and uses an order parameter varying between the bulk-fluid values.
- 1 Introduction: The model is mass-conservative under the stated boundary conditions, and its total energy combines kinetic and surface-energy contributions.The surface energy measures the fluid system’s interfacial energy.
- 1 Introduction: Degenerate mobility may better preserve the physical bound φ ∈[−1, 1], but its numerical resolution is beyond this paper’s scope.Uniqueness of weak solutions remains open even for the Cahn-Hilliard equation in that case.
- 1 Introduction: Small interfacial width ϵ creates severe stiffness, while nonlinear advection breaks the symmetry used by existing variational unique-solvability arguments.These issues motivate schemes that are unconditionally stable and unconditionally uniquely solvable.
- 1 Introduction: The proposed method is second order in time, uses convex splitting for Cahn-Hilliard and pressure projection for Navier-Stokes, and targets unconditional stability and unique solvability.A weakly coupled formulation also supports Picard-based decoupling and mixed finite-element implementation.
2 A Discrete Time, Continuous Space Scheme
The discrete scheme combines second-order time discretizations, convex splitting, and pressure projection while preserving a weak coupling between phase-field and flow computations. This structure supports stability, unique solvability, and further computational decoupling.
- 2 A Discrete Time, Continuous Space Scheme: The skew-symmetric advection form satisfies b(u,v,v)=0 even when the velocity arguments are not divergence-free, helping preserve stability after spatial discretization.This property follows from the form’s skew symmetry.
- 2 A Discrete Time, Continuous Space Scheme: The overall scheme uses Crank-Nicolson time discretization with second-order Adams-Bashforth extrapolation for the coupled system.The chemical-potential nonlinear term is treated through a second-order convex-splitting construction.
- 2 A Discrete Time, Continuous Space Scheme: Convex splitting treats the concave free-energy part explicitly and the convex part implicitly, enabling unconditional stability and unique solvability.The construction is designed for the overall CHNS scheme despite the Navier-Stokes advection term’s broken symmetry.
- 2 A Discrete Time, Continuous Space Scheme: The Navier-Stokes update uses a second-order incremental pressure projection with linear extrapolation for nonlinear advection.The projection maps an intermediate velocity that is not divergence-free into the divergence-free space via the Leray projection.
- 2 A Discrete Time, Continuous Space Scheme: The pressure-correction splitting is second-order accurate for velocity in l2(0, T; L2(Ω)) but only first-order accurate for pressure in l∞(0, T; L2(Ω)).The reported pressure-accuracy loss is attributed to the artificial pressure boundary condition.
- 2 A Discrete Time, Continuous Space Scheme: The pressure projection is decoupled, while phase-field and velocity equations interact only through advection velocity and elastic forcing chemical potential.This weak coupling enables Picard iteration and a further separation of nonlinear Cahn-Hilliard from linear Navier-Stokes computation.
- 2 A Discrete Time, Continuous Space Scheme: A monotonicity argument reduces the coupled system to a scalar chemical-potential equation and establishes unconditional unique solvability through the Browder-Minty lemma.The associated operator is strictly monotone, with equality only when its two arguments coincide.
3 Properties of the scheme
The scheme is mass-conservative, unconditionally stable through a modified energy law, and unconditionally uniquely solvable at every time step. These properties follow from energy estimates and monotonicity-based solvability analysis.
- Mass conservation and stability: The time-discrete scheme conserves mass and satisfies a modified energy law.The convective term vanishes through skew-symmetry, supporting the discrete energy estimate.
- Mass conservation and stability: The modified energy law establishes unconditional stability, allowing large time steps.The stability proof combines the Cahn-Hilliard and Navier-Stokes estimates.
- Unique solvability: The pressure equation is completely decoupled, so unique solvability reduces to the remaining coupled equations.Once the intermediate velocity is known, pressure and final velocity follow from a Darcy problem or pressure Poisson update.
- Unique solvability: Under the stated regularity assumptions, the full weak system has a unique solution at each time step.The result combines the solvability of the phase-field and velocity subproblems with the decoupled pressure step.
- Unique solvability: Browder-Minty analysis proves existence and uniqueness by establishing boundedness, continuity, coercivity, and strict monotonicity of an operator.The Cahn-Hilliard and velocity subproblems provide the solution operators used in the scalar formulation.
4 Mixed Finite Element Formulation
The scheme is spatially discretized with mixed finite elements under stable approximation assumptions. Its discrete conservation, stability, and unique-solvability properties are preserved, while Picard iteration weakly decouples the nonlinear phase-field and linear flow computations.
- Finite element formulation: Mixed finite element spaces are introduced on a quasi-uniform triangulation for the phase field, velocity, pressure, and chemical potential.The formulation assumes stable approximation spaces and an inf-sup condition for pressure stability.
- Finite element formulation: The fully discrete formulation preserves mass conservation, unconditional stability, and unconditional unique solvability.The pressure projection is formulated as a Darcy problem, yielding an optimal condition number for the pressure operator.
- Iterative solution procedure: Picard iteration decouples the nonlinear Cahn-Hilliard solve from the linear Navier-Stokes solve by iterating on velocity.The Cahn-Hilliard subproblem uses Newton’s method, followed by the linear flow solve until a fixed relative-error tolerance is reached.
- Implementation considerations: The method is initialized with a coupled first-order step because it is a two-step time-discretization method.Numerical simulations suggest that at least four grid elements across the interfacial region are needed for accuracy.
- Implementation considerations: Adaptive mesh refinement is used to improve algorithmic efficiency for the interfacial region.The implementation uses FreeFem++ with variable-metric/Delaunay automatic meshing.
5 Numerical Experiments
Numerical experiments assess convergence, energy dissipation, mass conservation, shape relaxation, adaptive refinement, shear-driven flow, and spinodal-decomposition coarsening.
- Numerical setup: The mixed finite element experiments use P1–P1 and P1b–P1 spaces, which satisfy the relevant inf-sup conditions.Other inf-sup-compatible spaces, including P2–P2 and Taylor–Hood P2–P1, can also be used.
- Convergence, energy dissipation, mass conservation: The scheme is second-order accurate for φ and u in L2 norm, while pressure convergence appears first-order.This conclusion comes from a Cauchy convergence test with grid refinement and δt = 0.2.
- Convergence, energy dissipation, mass conservation: Both discrete energy functionals are non-increasing in time, and their qualitative evolution is virtually identical.The auxiliary energy is a second-order approximation of the discrete energy in δt.
- Shape relaxation: Surface tension relaxes an isolated irregular shape toward a circular shape, while adaptive refinement resolves the diffuse interface.The refinement is shown at t = 0.02 and t = 0.4, with at least four grid cells across the interface.
- Spinodal decomposition: Hydrodynamic effects speed coarsening by promoting droplet coalescence, with stronger surface tension producing a more dramatic coalescence effect.For γ = ϵ, islands merge rapidly and exhibit fewer isolated drops; γ = 0 and γ = 0.1ϵ have nearly identical morphology over the evolution.
- Spinodal decomposition: The long-time simulation suggests a t^1/2 growth law for average domain size, while wall effects become influential near t = 10^4.Surface-energy decay is used as a proxy for phase-coarsening rate because it scales inversely with domain diameter near equilibrium.
6 Conclusions
The paper presents a second-order method for matched-density Cahn-Hilliard-Navier-Stokes flow that is efficient, unconditionally stable, and uniquely solvable. Numerical experiments verify conservation, energy dissipation, second-order accuracy, adaptive-mesh effectiveness, and a late-stage coarsening growth rate consistent with prior work, while several extensions remain open.
- The method decouples pressure from velocity and phase-field variables while retaining unconditional stability and unique solvability.The fully discrete finite-element methods have similar conclusions.
- Numerical experiments verify that the scheme is conservative, energy-dissipative, and second-order accurate in the L2 norm.
- Adaptive mesh refinement is effective for simulating shape relaxation with and without applied shear.
- A long-time simulation suggests a late-stage coarsening growth rate of t^1/2 for a large system, agreeing with prior work.
- Open directions include fully decoupling pressure, velocity, and phase field, extending to unmatched density or Cahn-Hilliard-Stokes-Darcy systems, and rigorously analyzing adaptive-mesh errors.