Negative time splitting is stable

by   Dong Li, et al.

For high order (than two) in time operator-splitting methods applied to dissipative systems, a folklore issue is the appearance of negative-time/backward-in-time linear evolution operators such as backward heat operators interwoven with nonlinear evolutions. The stability of such methods has remained an ensuing difficult open problem. In this work we consider a fourth order operator splitting discretization for the Allen-Cahn equation which is a prototypical high order splitting method with negative time-stepping, i.e. backward in time integration for the linear parabolic part. We introduce a new theoretical framework and prove uniform energy stability and higher Sobolev stability. This is the first strong stability result for negative time stepping operator-splitting methods.



There are no comments yet.


page 1

page 2

page 3

page 4


On the energy stability of Strang-splitting for Cahn-Hilliard

We consider a Strang-type second order operator-splitting discretization...

The operator-splitting method for Cahn-Hilliard is stable

We prove energy stability of a standard operator-splitting method for th...

Single-Forward-Step Projective Splitting: Exploiting Cocoercivity

This work describes a new variant of projective splitting in which cocoe...

Splitting Schemes for Some Second-Order Evolution Equations

We consider the Cauchy problem for a second-order evolution equation, in...

Enforcing strong stability of explicit Runge–Kutta methods with superviscosity

A time discretization method is called strongly stable, if the norm of i...

Uniform preconditioners of linear complexity for problems of negative order

We propose a multi-level type operator that can be used in the framework...

A posteriori learning of quasi-geostrophic turbulence parametrization: an experiment on integration steps

Modeling the subgrid-scale dynamics of reduced models is a long standing...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

We consider the Allen-Cahn equation


where the unknown . The parameter is called the mobility coefficient and we fix it as a constant for simplicity. The nonlinear term takes the form , where is the standard double well. To minimize technicality, we take the spatial domain in (1.1) as the -periodic torus . With some additional work our analysis can be extended to physical dimensions . The system (1.1) arises as a -gradient flow of a Ginzburg-Landau type energy functional , where


The basic energy conservation law takes the form


It follows that


Besides the -type conservation law, there is also -type control. Due to special form of the nonlinearity, for smooth solutions we have the maximum principle


As a consequence the long-time wellposedness and regularity is not an issue for (1.1).

The objective of this work is to establish strong stability of a fourth order in time operating splitting algorithm applied to the Allen-Cahn equation (1.1). This is a part of our on-going program to develop a new theory for the rigorous analysis of stability and convergence of operator-splitting methods applied to dissipative-type problems. Due to various subtle technical obstructions, there were very few rigorous results on the analysis of the operator-splitting type algorithms for the Allen-Cahn equation, the Cahn-Hilliard equation and similar models. Prior to our recent series of works [17, 18], most existing results in the literature are conditional one way or another. To put things into perspective, we briefly review a few closely related representative works and more recent developments.

  • The work of Gidey-Reddy. Gidey and Reddy considered in [24] a convective Cahn-Hilliard model of the form


    where . They adopted an operator-splitting of (1.6) into the hyperbolic part, nonlinear diffusion part and diffusion part respectively. Several conditional results concerning certain weak solutions were obtained.

  • The work of Weng-Zhai-Feng. In [25], Weng, Zhai and Feng studied a viscous Cahn-Hilliard model:


    where the parameter . The authors employed a fast explicit Strang-type operator splitting and proved the stability and the convergence (see Theorem 1 on pp. 7 of [25]) under the assumption that , stay bounded, and satisfy a technical condition . Here denotes the numerical solution.

  • The work of Cheng-Kurganov-Qu-Tang. In [23], Cheng, Kurganov, Qu and Tang considered the Cahn-Hilliard equation


    and the MBE equation


    Concerning the Cahn-Hilliard equation, the authors considered a Strang-type splitting approximation of the form




    and is the nonlinear propagator


    Various conditional results were given in [23] but the rigorous analysis of energy stability was a long-standing open problem. This problem and several related open problems were settled in our recent proof [18].

In recent [17], we carried out the first energy-stability analysis of a first order operator-splitting approximation of the Cahn-Hilliard equation. More precisely denote as the exact PDE solution to the Cahn-Hilliard equation . We considered a splitting approximation of the form:


where and solves


We introduced a novel modified energy and rigorously proved monotonic decay of the new modified energy which is coercive in -sense. Moreover we obtained uniform control of higher Sobolev regularity and rigorously justified the first order convergence of the method on any finite time interval.

In [18], we settled the difficult open problem of energy stability of the Strang-type algorithm applied to the Cahn-Hilliard equation which was introduced in the work of Cheng, Kurganov, Qu and Tang [23]. One should note that the new theoretical framework developed in [18] is completely different from the first order case [17]. In the second-order case, for most generic data one no longer has strict energy-monotonicity at disposal and several new ideas such as dichotomy analysis, an absorbing-set approach were introduced in [18] in order to settle long time unconditional (i.e. independent of time-step or time interval) Sobolev bounds of the numerical solution.

In this work we develop further the program initiated in [17, 18] and turn to the analysis of higher order in time splitting methods. A folklore issue is that operator-splitting methods with order higher than two require negative time-stepping, i.e. backward time integration for each time step. In [33], Sheng considered dimensional splitting for the two-dimensional parabolic problem . Using semi-discretization one obtains the ODE system , where the matrices and correspond to the discretization of and respectively. One is then naturally led to the approximation of . In [33], Sheng showed that if the non-commuting matrices and have eigen-values in the left half of the complex-plane, and one employs the approximation


then the highest order of a stable approximation is two even if is chosen to be large. Put it differently, Sheng’s fundamental result states that for an -order () partitioned split-step schemes, at least one of the solution operators must be applied with negative time step, i.e. backwardly. In [34] Suzuki adopted an elegant time-ordering principle from quantum mechanics and proved a general nonexistence theorem of positive decomposition for high order splitting approximations. In [35], Goldman and Kaper strengthened these results further and showed that even with partitioned schemes, each solution operator within a convex partition must be performed with at least one negative/backward fractional time step.

For deterministic Hamiltonian-type systems, backward time stepping in general does not create instabilities and high-order operating splitting methods have shown promising effectiveness [32]. On the other hand, due to the negative time stepping, there are some arguments that higher-than-three operator splitting methods cannot be applied to parabolic equations with diffusive terms ([33, 36, 37]). As we shall explain momentarily, these concerns are not unsubstantiated and even turn up in the well-defined-ness of these algorithms.

Although rigorous analysis of these issues were not available before, recently Cervi and Spiteri [31] considered three third-order operating-splitting methods and demonstrated via extensive numerical simulations the effectiveness of higher order splitting methods for representative cardiac electrophysiology simulations. These compelling numerical evidences propel us to re-examine in detail the stability property of general negative/backward time-stepping methods in parabolic problems. Indeed, the very purpose of this work is to break these aforementioned technical barriers and establish a new stability theory for these problems.

To set the stage and minimize technicality, we consider the Allen-Cahn equation (1.1) posed on the -periodic torus . To build some intuition, let and consider


and let solve the equation


An explicit formula for is readily available, indeed it is not difficult to check that


Define a second-order Strang-type propagator


Following Yoshida [29], we consider a 4th order integrator obtained by a symmetric repetition (product) of the 2nd order integrator:




Somewhat surprisingly, we first show that the propagator is ill-defined if one does not introduce a judiciously chosen spectral cut-off.

Proposition 1.1 (Ill-definedness of with a spectral cut-off).

The propagator is ill-defined in general.

The proof of Proposition 1.1 is given in Section 2. Armed with this important observation, we are led to introduce a spectral cut-off condition for the propagators. Let be an integer and we shall consider the projection operator defined for via the relation


where denotes the Fourier coefficient of . We introduce the following spectral cut-off condition.

Definition (Spectral condition). Let be the spectral truncation parameter. We shall say satisfy the spectral condition if for some .

Remark 1.1.

This constraint in is reminiscent of the CFL condition in hyperbolic problems. Here in the parabolic setting we require that the operator to remain bounded when restricted to the spectral cut-off .

Let . We now consider the following modified propagators:

Theorem 1.1 (Stability of negative-time splitting methods).

Let , and consider the Allen-Cahn equation (1.1) on the one-dimensional -periodic torus . Assume the spectral cut-off condition for some . Assume the initial data ( is an integer). Let and define


There exists a constant depending only on , and , such that if , then


where depends on (, , , ).

Remark 1.2.

One can also show the fourth-order convergence of the operator splitting approximation. Namely if and let be the exact PDE solution corresponding to initial data . Let bet the same as in Theorem 1.1. Then for any , we have


where depends on (, , , ). The regularity assumption on initial data can be lowered. We shall refrain from stating such a pedestrian result here and leave its justification to interested readers as exercises.

The rest of this paper is organized as follows. In Section we set up the notation and collect some preliminary lemmas. In Section we give the proof of Theorem 1.1.

2 Notation and preliminaries

For any two positive quantities and , we shall write or if for some constant whose precise value is unimportant. We shall write if both and hold. We write if the constant depends on some parameter . We shall write if and if .

We shall denote if for some sufficiently small constant . The smallness of the constant is usually clear from the context. The notation is similarly defined. Note that our use of and here is different

from the usual Vinogradov notation in number theory or asymptotic analysis.

For any , we denote , and . Also occasionally we use the Japanese bracket notation:

We denote by the usual -periodic torus. For and any function , we denote the Lebesgue -norm of as

If is a sequence of complex numbers and is the index set, we denote the discrete -norm as


For example, . If

is a vector-valued function, we denote

, and . We use similar convention for the corresponding discrete norms for the vector-valued case.

We use the following convention for the Fourier transform pair:


and denote for ,

Proof of Proposition 1.1.

To prove Proposition 1.1, it suffices for us to examine the following statement. Consider where is the periodic Dirac comb on . Let and consider


Claim: for any .

Proof of Claim. Observe that the Fourier coefficients of are all positive and


Clearly the desired conclusion follows. ∎

3 Proof of Theorem 1.1

Lemma 3.1 (One-step stability).

Let , and assume the spectral condition for some constant . Suppose , and . Let . There exists sufficiently small such that if , then


where depends on (, , , , , ).

Remark 3.1.



Later we shall apply Lemma 3.1 -times. In particular we need to bound the expression


Since the constant depends on the -norm of the iterates, it is of importance to give uniform control of the -norm. To resolve this issue, we first choose , and such that


We choose such that


Clearly in the course of iteration, the -norm of the iterates never exceeds . In particular

Proof of Lemma 3.1.

We begin by noting that for the ODE


we have the estimate (note that

is an algebra for )


It follows that for sufficiently small,


The desired estimates then easily follows from the above using the spectral condition. ∎

Lemma 3.2 (Multi-step stability and higher regularity).

Let , and assume the spectral condition for some constant . Suppose and for some constant . Define


where , . For , define


where . There exist and such that if (we may assume ), then


where depends on (, , ).

Proof of Lemma 3.2.

Observe that


where . Clearly the operators , , fulfill the conditions of Lemma 3.1. The estimates then follow easily from the computations outlined in Remark 3.1 with some necessary adjustment of the constants.

To establish (3.16), we note that for ,


provided . We then rewrite as


where . One can then bootstrap the higher regularity from this. We omit further details. ∎

Remark 3.2.

It is not difficult to check that under the assumption of uniform high Sobolev regularity, we have (below we assume )


It follows that

Lemma 3.3 (Control of the energy flux).

Let and . Suppose satisfies and




where depends only on . Furthermore if and where is a sufficiently small constant depending on (, ), then


where depends only on (, ).


The estimate (3.23) follows from energy estimates using (3.22). Note that the condition is used in the identity . The estimate (3.24) follows from (3.23) and Lemma 3.1. ∎

Lemma 3.4 (One-step strict energy dissipation with nontrivial energy flux).

Let , and . Suppose with and satisfies


where is a given constant. There exists sufficiently small such that if , then

Proof of Lemma 3.4.

By using Remark 3.2, it suffices for us to show


for some . Denote and observe that


It follows that for sufficiently small,


where is a constant. Now we clearly have


The desired result clearly follows. ∎

Theorem 3.1.

Let , and . Assume and . Define and


There exists sufficiently small such that if , then


where depends only on (, , ).

Proof of Theorem 3.1.



where the constants , are the same as in Lemma 3.3. In the argument below we shall assume is sufficiently small. The needed smallness (i.e. the existence of ) can be easily worked out by fulfilling the conditions needed in Lemma 3.2, Lemma 3.3 and Lemma 3.4.

By Lemma 3.2, we can find such that


Claim: We have


To prove the claim we argue by contradiction. Suppose is the first integer such that


By Lemma 3.3 we must have


Since , we have for some integer . By using smoothing estimates we obtain


where depends on (, , ). Since depends on (, , ), we have depends only on (, , ). By (3.37), (3.38) and Lemma 3.4, we obtain for sufficiently small that


which is clearly a contradiction to (3.36). Thus we have proved the claim. ∎

Proof of Theorem 1.1.

The estimate follows from Theorem 3.1. Higher order estimates follow from the smoothing estimates. ∎


  • [1] A. Baskaran, J. S. Lowengrub, C. Wang, and S. M. Wise. Convergence Analysis of a Second Order Convex Splitting Scheme for the Modified Phase Field Crystal Equation. SIAM J. Numer. Anal., 51(2013), 2851šC2873.
  • [2] J.W. Cahn, J.E. Hilliard. Free energy of a nonuniform system. I. Interfacial energy free energy, J. Chem. Phys. 28 (1958) 258–267.
  • [3] L.Q. Chen, J. Shen. Applications of semi-implicit Fourier-spectral method to phase field equations. Comput. Phys. Comm., 108 (1998), pp. 147–158.
  • [4] A. Christlieb, J. Jones, K. Promislow, B. Wetton, M. Willoughby. High accuracy solutions to energy gradient flows from material science models. J. Comput. Phys. 257 (2014), part A, 193–215.
  • [5] C.M. Elliott and A.M. Stuart. The global dynamics of discrete semilinear parabolic equations. SIAM J. Numer. Anal., 30 (1993), pp. 1622–1663.
  • [6] D. J. Eyre. Unconditionally gradient stable time marching the Cahn-Hilliard equation. Computational and mathematical models of microstructural evolution (San Francisco, CA, 1998), 39–46, Mater. Res. Soc. Sympos. Proc., 529, MRS, Warrendale, PA, 1998.
  • [7] X. Feng, T. Tang and J. Yang. Long time numerical simulations for phase-field problems using p-adaptive spectral deferred correction methods. SIAM J. Sci. Comput. 37 (2015), no. 1, A271–A294.
  • [8] H. Gomez and T.J.R. Hughes. Provably unconditionally stable, second-order time-accurate, mixed variational methods for phase-field models. J. Comput. Phys., 230 (2011), pp. 5310-5327
  • [9] Z. Guan, C. Wang and S. M. Wise, A convergent convex splitting scheme for the periodic nonlocal Cahn-Hilliard equation, Numer. Math., 128 (2014), 377–406.
  • [10] J. Guo, C. Wang, S. M. Wise and X. Yue, An convergence of a second-order convex-splitting, finite difference scheme for the three-dimensional Cahn-Hilliard equation, Commu. Math. Sci., 14 (2016), 489–515
  • [11] Y. He, Y. Liu and T. Tang. On large time-stepping methods for the Cahn-Hilliard equation. Appl. Numer. Math., 57 (2007), 616–628.
  • [12] F. Liu and J. Shen. Stabilized semi-implicit spectral deferred correction methods for Allen-Cahn and Cahn-Hilliard equations. Math. Methods Appl. Sci. 38 (2015), no. 18, 4564–4575.
  • [13] D. Li, Z. Qiao and T. Tang. Characterizing the stabilization size for semi-implicit Fourier-spectral method to phase field equations, SIAM J. Numer. Anal., 54 (2016), 1653–1681
  • [14] D. Li, F. Wang and K. Yang. An improved gradient bound for 2D MBE. Journal of Differential Equations, 269(12): 11165-11171, 2020.
  • [15] D. Li, T. Tang. Stability of the Semi-Implicit Method for the Cahn-Hilliard Equation with Logarithmic Potentials. Ann. Appl. Math., 37 (2021), 31–60.
  • [16] D. Li, C. Quan, and T. Tang. Stability analysis for the implicit-explicit discretization of the Cahn-Hilliard equation.
  • [17] D. Li and C. Quan. The operator-splitting method for Cahn-Hilliard is stable. arXiv:2107.01418, 2021.
  • [18] D. Li and C. Quan. On the energy stability of Strang-splitting for Cahn-Hilliard. arXiv:2107.05349, 2021.
  • [19] C.B. Schönlieb and A. Bertozzi. Unconditionally stable schemes for higher order inpainting. Commun. Math. Sci. 9 (2011), no. 2, 413–457.
  • [20] J. Shen and X. Yang. Numerical approximations of Allen-Cahn and Cahn-Hilliard equations. Discrete Contin. Dyn. Syst. A, 28 (2010), 1669–1691.
  • [21] J. Shen, J. Xu and J. Yang. A new class of efficient and robust energy stable schemes for gradient flows. SIAM Rev. 61 (2019), no. 3, 474–506.
  • [22] H. Song and C. Shu. Unconditional energy stability analysis of a second order implicit-explicit local discontinuous Galerkin method for the Cahn-Hilliard equation. J. Sci. Comput. 73 (2017), no. 2-3, 1178–1203.
  • [23] Y. Cheng, A. Kurganov, Z. Qu, and T. Tang. Fast and stable explicit operator splitting methods for phase-field models. Journal of Computational Physics, 303:45–65, 2015.
  • [24] H.H. Gidey and B.D. Reddy. Operator-splitting methods for the 2D convective Cahn-Hilliard equation. Computers & Mathematics with Applications, 2019, 77(12): 3128–3153.
  • [25]

    Z. Weng, S. Zhai and X. Feng. Analysis of the operator splitting scheme for the Cahn-Hilliard equation with a viscosity term. Numer. Meth. for Partial Differential Equations. 35 (2019), no. 6, 1949–1970.

  • [26] C. Xu and T. Tang. Stability analysis of large time-stepping methods for epitaxial growth models. SIAM J. Numer. Anal. 44 (2006), no. 4, 1759–1779.
  • [27] J. Zhu, L.-Q. Chen, J. Shen, and V. Tikare. Coarsening kinetics from a variable-mobility Cahn-Hilliard equation: Application of a semi-implicit Fourier spectral method, Phys. Rev. E (3), 60 (1999), pp. 3564–3572.
  • [28] R. I. Issa. Solution of the implicitly discretised fluid flow equations by operator-splitting. Journal of Computational Physics 62, no. 1 (1986): 40–65.
  • [29] H. Yoshida. Construction of higher order symplectic integrators. Physics letters A 150, no. 5–7 (1990): 262-268.
  • [30] R. I. McLachlan and G. R. W. Quispel. Splitting methods. Acta Numerica 11 (2002): 341–434.
  • [31] J. Cervi and R. J. Spiteri. High-order operator splitting for the bidomain and monodomain models. SIAM Journal on Scientific Computing 40, no. 2 (2018): A769–A786.
  • [32] M. Thalhammer, M. Caliari and C. Neuhauser. High-order time-splitting Hermite and Fourier spectral methods. J. Comput. Phys. 228 (2009), pp. 822–832.
  • [33] Q. Sheng. Solving linear partial differential equations by exponential splitting. IMA J. Numer. Anal. 9 (1989), pp. 199–212.
  • [34] M. Suzuki. General theory of fractal path integrals with applications to many-body theories and statistical physics. J. Math. Phys., 32 (1991), pp. 400–407.
  • [35] G. Goldman and T.J. Kaper. Nth-order operator splitting schemes and nonreversible systems. SIAM J. Numer. Anal. 33 (1996), pp. 349–367.
  • [36] S.A. Chin. Quantum statistical calculations and symplectic corrector algorithms. Phys. Rev. E. 69 (2004), 046118.