Energy Stable Numerical Schemes for Ternary Cahn-Hilliard System

Abstract

We present and analyze a uniquely solvable and unconditionally energy stable numerical scheme for the ternary Cahn-Hilliard system, with a polynomial pattern nonlinear free energy expansion. One key difficulty is associated with presence of the three mass components, though a total mass constraint reduces this to two components. Another numerical challenge is to ensure the energy stability for the nonlinear energy functional in the mixed product form, which turns out to be non-convex, non-concave in the three-phase space. To overcome this subtle difficulty, we add a few auxiliary terms to make the combined energy functional convex in the three-phase space, and this, in turn, yields a convex-concave decomposition of the physical energy in the ternary system. Consequently, both the unique solvability and the unconditional energy stability of the proposed numerical scheme are established at a theoretical level. In addition, an optimal rate convergence analysis in the $\ell ^\infty (0,T; H_N^{-1}) \cap \ell ^2 (0,T; H_N^1)$ norm is provided, with Fourier pseudo-spectral discretization in space, which is the first such result in this field. To deal with the nonlinear implicit equations at each time step, we apply an efficient preconditioned steepest descent (PSD) algorithm. A second order accurate, modified BDF scheme is also discussed. A few numerical results are presented, which confirm the stability and accuracy of the proposed numerical scheme.

A Proof of Lemma 6

To prove the inverse inequality (3.16), we begin with the following observation:

$$\begin{aligned} \left( \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^p \right) ^{q/p} \ge \sum _{i,j,k=0}^{2K} ( | f_{i,j,k} |^p )^{q/p} = \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^q , \quad \forall 1\le p \le q < \infty ,\nonumber \\ \end{aligned}$$

(A.1)

due to the fact that $\frac{q}{p} \ge 1$. This is equivalent to

$$\begin{aligned} \left( \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^q \right) ^{1/q} \le \left( \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^p \right) ^{1/p} , \quad \forall 1\le p \le q < \infty . \end{aligned}$$

(A.2)

In turn, we get

$$\begin{aligned} \left( h^3 \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^q \right) ^{1/q} \le h^{3/q - 3/p} \left( h^3 \sum _{i,j,k=0}^{2K} | f_{i,j,k} |^p \right) ^{1/p} , \end{aligned}$$

(A.3)

which is exactly (3.16), because of the fact that $h = \frac{L}{2K+1}$.

For the rest three inequalities, we assume the periodic grid f has the discrete Fourier transformation as given by (3.1). For simplicity of presentation, we also assume that $L=1$. In addition, we denote its extension to a continuous function as

$$\begin{aligned} f_{\mathbf{F}} (x,y,z) = \sum _{\ell , m, n=-K}^{K} {\hat{f}}_{\ell , m, n}^{N} \exp \left( 2\pi \mathrm {i} (\ell x +m y +n z )\right) . \end{aligned}$$

(A.4)

The following estimates could be derived with the help of the Parseval’s identity; also see the related analysis in  [30, 42, 43, 52]:

$$\begin{aligned} \Vert f \Vert _2^2= & {} \Vert f_{\mathbf{F}} \Vert _{L^2}^2 = \sum _{\ell , m, n=-K}^{K} | {\hat{f}}_{\ell , m, n}^{N} |^2 , \end{aligned}$$

(A.5)

$$\begin{aligned} \overline{f}= & {} \int _\Omega \, f_{\mathbf{F}} d \mathbf{x} = {\hat{f}}_{0,0,0}^{N} , \end{aligned}$$

(A.6)

$$\begin{aligned} \Vert \nabla _N f \Vert _2^2= & {} \Vert \nabla f_{\mathbf{F}} \Vert ^2 . \end{aligned}$$

(A.7)

As a result, the discrete Poincaré inequality (3.17) is a direct consequence of its continuous version:

$$\begin{aligned} \Vert f_{\mathbf{F}} \Vert _{H^1}&\le C_2 ( | \overline{ f_{\mathbf{F}}} | + \Vert \nabla f_{\mathbf{F}} \Vert ) , \end{aligned}$$

(A.8)

combined with the identities (A.5)–(A.7).

On the other hand, we recall a key estimate given by Lemma A.2 in an existing work  [43]:

$$\begin{aligned} \Vert f \Vert _p \le \sqrt{\frac{p}{2}} \Vert f_{\mathbf{F}} \Vert _{L^p} , \quad \text{ in } \text{2-D, } p = 4, 6, ... , \end{aligned}$$

(A.9)

In addition, such an analysis could also be extended to the 3-D case

$$\begin{aligned} \Vert f \Vert _p \le ( \frac{p}{2} )^{3/4} \Vert f_{\mathbf{F}} \Vert _{L^p} , \quad \text{ in } \text{3-D, } p = 4, 6, ... . \end{aligned}$$

(A.10)

For the 3-D discrete Sobolev inequality (3.18), we begin with the observation that

$$\begin{aligned} f \in \mathring{{\mathcal {G}}} \, \, \, \text{ is } \text{ equivalent } \text{ to } {\hat{f}}_{0,0,0}^{N}=0\text{, } \text{ so } \text{ that } \int _\Omega \, f_{\mathbf{F}} d \mathbf{x} =0 . \end{aligned}$$

(A.11)

Subsequently, an application of the 3-D estimate (A.10) results in

$$\begin{aligned} \Vert f \Vert _6 \le 3^{3/4} \Vert f_{\mathbf{F}} \Vert _{L^6} \le C_4 \Vert \nabla f_{\mathbf{F}} \Vert = C_4 \Vert \nabla _N f \Vert _2 , \end{aligned}$$

(A.12)

in which a continuous Sobolev embedding has been applied in the second step, due to the fact that $\int _\Omega \, f_{\mathbf{F}} d \mathbf{x} =0$.

Similarly, for the 2-D discrete Sobolev inequality (3.19), we apply the estimate (A.9) and obtain

$$\begin{aligned} \Vert f \Vert _8 \le 2 \Vert f_{\mathbf{F}} \Vert _{L^8} \le C_5 \Vert f_{\mathbf{F}} \Vert _{H^1} = C_5 \Vert f \Vert _{H_N^1} . \end{aligned}$$

(A.13)

Again a continuous Sobolev embedding has been applied in the second step, and the last step is based on the identities (A.5), (A.7).