A maximum principle of the Fourier spectral method for diffusion equations

1.   Introduction

In this study, we investigate the discrete maximum principle of the Fourier spectral method (FSM) for the following diffusion equation:

where $u({\bf x}, t)$ is the density of the diffusing material at location ${\bf x}$ and time $t$; ${\bf n}$ is the unit normal vector to the domain boundary $\partial \Omega$. Equation (1.1) can be derived from a total free energy functional as a gradient flow:

i.e.,

where $\delta \mathcal{E} (u) / \delta u$ is the variational derivative of $\mathcal{E}(u)$ with respect to $u$. In ^[1], the FSM was introduced for fractional diffusion with Dirichlet and Neumann boundary conditions. This approach can be extended to other equations. The study in ^[2] investigated spectral methods and the maximum principle for boundary problems involving the second derivative. Li ^[3] proposed effective maximum principles that allow numerical solutions to deviate from the sharp bound through the use of a controllable discretization error. The author analyzed FSMs for various partial difference equations (refer to ^[3] for further details).

Lee ^[4] showed that the discrete maximum principle holds for the implicit finite difference scheme of the diffusion equation in the operator splitting method for the following Allen–Cahn (AC) equation:

where $\epsilon$ is a small positive parameter. The AC equation models the process of phase separation in binary alloy mixtures. Ayub et al. ^[5] solved the AC equation by using unconditionally stable operator splitting schemes. Ham et al. ^[6] proposed an unconditionally stable scheme for a parabolic sine-Gordon equation based on an operator splitting method. In ^[5,6], the authors solved the split diffusion equation by using the FSM, and in this case, the maximum principle of the diffusion equation was used to demonstrate the unconditional stability of the proposed method. Sun et al. ^[7] developed a maximum-principle-preserving method for the AC equation. The authors extended the method by incorporating the integrating factor two-step Runge–Kutta method and proposed sufficient conditions for the preservation of the discrete maximum principle. Lee ^[8] used the FSM for space to propose high-order and mass conservative methods for the conservative AC equation. In ^[9], the authors proposed several types of methods for solving the non-local AC equation, including the Fourier spectral operator splitting method. Following this study, Weng and Tang showed that the proposed method is unconditionally stable and has second-order accuracy with respect to time by using the FSM for the diffusion term in ^[10]. Alzahrani and Khaliq ^[11] investigated solving the space-fractional reaction-diffusion equation based on the FSM. Their proposed method has spectral convergence, represents fractional operators diagonally and can be applied to multi-dimensional space.

Regarding convection-diffusion equations, Chertock et al. ^[12] proposed a fast explicit operator splitting method. In the proposed method, the parabolic equation was solved by using the FSM. Abbaszadeh and Amjadian ^[13] analyzed the spectral element method for the fractional advection-diffusion equation. Verrall and Read ^[14] solved the combined advection-diffusion-reaction equation by using analytical solutions for decoupled equations. The diffusion problem was solved by using a pseudo-spectral approach. Considering the above-mentioned equations, the diffusion equation is an important building block equation for modeling scientific phenomena.

Due to its efficiency and accuracy, the FSM is extensively employed to address various problems, including those equations mentioned above. In particular, when the linear term takes the form of the diffusion equation, the operator splitting method is used to solve these problems. In this process, the maximum principle of the diffusion equation is used to demonstrate the stability of the numerical schemes ^[6,15]. However, due to the oscillatory profile of the non-smooth initial conditions of a continuous representation in the FSM, it is possible for the discrete maximum principle to be violated. The primary objective of this study is to present counterexamples demonstrating that the discrete maximum principle may not hold for the FSM when applied to the diffusion equation.

The contents of this paper are as follows. In Section 2, we briefly describe the FSM for the diffusion equation. In Section 3, we present numerical experiments to demonstrate that the discrete maximum principle may not hold for the FSM. Conclusions are given in Section 4.

2.   Fourier spectral method

Now, we describe the FSM for the diffusion equation (1.1). In this paper, to apply the homogeneous Neumann boundary condition (1.2), we consider the discrete cosine transform (DCT) and inverse DCT (IDCT). For simplicity of exposition, we consider the one-dimensional space $\Omega = (0, L)$. Let $x_i = (i-0.5)h$ for $i = 1, \dots, N$, where $h = L/N$ and $N$ is a positive integer. Let $t_n = n \Delta t$ for $n = 1, \dots$, where $\Delta t$ is the time step size. For the given values $u(x_i, t_n)$ for $i = 1, \dots, N$ and some $n$, the DCT and IDCT are defined as follows:

respectively, where $\alpha_1 = \sqrt{1 / N}$, $\alpha_p = \sqrt{2 / N}$ for $p\geq2$, and $\xi_p = (p-1)/L$ for $p = 1, \dots, N$ ^[5,9,15,16]. In Eq (2.2), let us replace $x_i$ and $t_n$ by $x$ and $t$, respectively; then, we have a continuous representation function, as follows:

We can use the partial derivatives of $U(x, t)$ in Eq (2.2) with respect to $t$ and $x$ as follows:

Therefore, Eq (1.1) can be written as follows:

which results in

Equation (2.4) is an ordinary differential equation and we can solve it analytically as follows:

Finally, we have the numerical solution at time $t = t_{n+1}$ by using Eqs (2.1), (2.3) and (2.5) as follows:

3.   Numerical experiments

The FSM is a very fast, efficient and accurate method for solving diffusion equations. However, the discrete maximum principle may not hold for the FSM and we present numerical experiments to demonstrate that fact. Let us consider the following smooth (Eq (3.1)) and stiff (Eq (3.2)) initial conditions on the domain $\Omega = (0, \; 2)$ as shown in Figure 1(a), (b), respectively.

Figure 1.  Discrete (circles) and continuous (solid line) representations of the following initial conditions: (a) $u(x, 0) = \cos(0.5\pi x^2)$ and (b) $u(x, 0) = 1$ if $x > 1$; $u(x, 0) = 0.5$ if $x = 1$; $u(x, 0) = 0$ if $x < 1$. Here, the dashed lines are the given initial conditions.

DownLoad: Full-Size Img PowerPoint

As shown in Figure 1(a), if the initial condition is smooth, then the continuous representation of the initial condition using the FSM is almost identical to the given initial condition, Eq (3.1). In the case of the stiff step function, the discrete initial values at the discrete domain coincide with the initial condition. From the continuous representation function given by Eq (2.3), we have the following discrete initial values at the discrete domain, which are the same as for the initial condition (see the open circles in Figure 1(b)):

However, the following continuous function has different values from the initial condition at the points which are not on the discrete domain (see the solid line in Figure 1(b)):

This fact may lead to the violation of the maximum principle for the discrete numerical solution for the diffusion equation by using the FSM.

We investigate the reason why the discrete solution of the diffusion equation by using the FSM may not satisfy the discrete maximum principle. We consider the step function (3.2) in the one-dimensional computational domain $\Omega = (0, 2)$ as an initial condition for the diffusion equation with the homogenous Neumann boundary condition (1.2). For better visualization, we use only six uniform discrete spatial points, i.e., $x = 1/6, \; 3/6, \; 5/6, \; 7/6, \; 9/6, \; 11/6$ and $N = 6$. Let us decompose the continuous solution $U(x, t)$ into $\sum_{p = 1}^6 v_p(x, t)$, where $v_p(x, t)$ $= \alpha_p \hat{u}(p, t) \cos(\xi_p \pi x)$. Figure 2(a) shows the discrete initial condition $u(x_i, 0)$ for $i = 1, \ldots, 6$, continuous representation of the initial condition $U(x, 0)$ and some of its decomposed components $v_p(x, 0)$ for $p = 1, 2, 4, 6$. Figure 2(b) shows the temporal evolution of the discrete and continuous solutions. We can observe that the discrete maximum principle is violated, i.e., the values at $x = 1.5$ are greater than one until $t = 0.0215$. Figure 2(c) shows the temporal evolution of the discrete solutions at $x = 1.5$ of the decomposed components $v_p(x, 0)$ for $p = 1, 2, 3, 4, 5, 6$. Figure 2(d) plots the temporal evolution of the continuous solution $U(1.5, t)$ from $t = 0$ to $t = 0.025$. The discrete solutions on $x = 1.5$ at $t = 0.01$, $0.015$ and $0.02$ do not satisfy the discrete maximum principle, i.e., the solutions are greater than one. As shown in Figure 2(c), the high frequency cosine function with a negative value ($p = 6$) at $x = 1.5$ is damped faster than the low frequency cosine functions. Therefore, the discrete composed solution $u(1.5, t)$ is greater than one on $x = 1.5$ as shown in Figure 2(d).

Figure 2.  (a) Discrete initial condition $u(x_i, 0)$ for $i = 1, \ldots, 6$, continuous representation of the initial condition $U(x, 0)$ and some of its decomposed components $v_p(x, 0)$ for $p = 1, 2, 4, 6$. (b) Temporal evolution of the discrete and continuous solutions. (c) Temporal evolution of the continuous solutions at $x = 1.5$ of the decomposed components $v_p(x, 0)$ for $p = 1, 2, 3, 4, 5, 6$. (d) Temporal evolution of the composed solution $u(1.5, t)$ from $t = 0$ to $t = 0.025$.

DownLoad: Full-Size Img PowerPoint

Figure 3 shows solutions of continuous and discrete functions from the FSM. Here, we use uniform discrete points: $N = 50$ and $h = 2/N$. Figure 3(a) shows the initial values of continuous and discrete functions. Figure 3(b) shows the temporal evolution of solutions with the initial condition. We can observe that the solutions at times $t = T_1, \; T_2, \; T_3$ are greater than one. Figure 3(c) shows the maximum values of the continuous (solid line) and discrete (dashed line) solutions over time. As shown in Figures 3(b), (c), the discrete solution is not bounded by the discrete initial values 0 and 1. However, it is bounded by the maximum and the minimum values of the continuous initial values. These results demonstrate that the continuous solution satisfies the maximum principle, although the discrete solution may not satisfy that.

Figure 3.  (a) Initial values of continuous and discrete functions. (b) Temporal evolution of solutions with the initial condition. (c) Maximum values of the continuous (solid line) and discrete (dashed line) solutions over time. Here, $T_0 = 0$, $T_1 = 0.25e$-4, $T_2 = 1.43e$-4, $T_3 = 2.5e$-4 and $T_4 = 04e$-4.

DownLoad: Full-Size Img PowerPoint

Now, we investigate the effect of the number of grid points $N$ on the evolution dynamics of the numerical solutions for the diffusion equation. We solve the diffusion equation by using the FSM with the step function (3.2) as an initial condition for the different number of grid points $N = 12$, $16$ and $20$. Figure 4(a)–(c) show the temporal evolutions of the continuous and discrete solutions with $N = 12$, $16$ and $20$, respectively. Figure 4(d) shows the temporal evolutions of the maximum values of the discrete solutions with $N = 12$, $16$ and $20$. As shown in Figure 4(d), we can observe that a larger value of $N$ has the effect that the artificial maximum of the discrete solution decreases faster.

Figure 4.  (a)–(c) Temporal evolutions of the continuous and discrete solutions with $N = 12$, $16$ and $20$, respectively. (d) Temporal evolutions of the maximum values of the continuous solutions with $N = 12$, $16$ and $20$.

DownLoad: Full-Size Img PowerPoint

Let us consider the following two smooth and discontinuous functions as the initial condition of the diffusion equation:

which are both bounded by zero and one. For two initial conditions, in Figure 5, we use the uniform discrete points $N = 40$ and grid size $h = 2/N$ on the computational domain $\Omega = (0, 2)$. Figure 5(a), (b) show continuous (solid line) and discrete (circles) representations of the smooth and discontinuous initial conditions, respectively. Figures 5(c), (d) show the temporal evolution of the maximum values of the continuous and discrete representation solutions for the diffusion equation for the smooth initial condition and discontinuous initial condition, respectively. From the results of Figure 5, in the case of the smooth initial condition, the continuous representation is consistent with the discrete representation and the solution satisfies the maximum principle. Conversely, in the case of the discontinuous initial condition, the continuous representation is inconsistent with the discrete representation and the discrete solution does not satisfy the discrete maximum principle.

Figure 5.  Continuous (solid line) and discrete (circles) representations of the following initial conditions: (a) $u(x, 0) = (\cos(0.5\pi x^2)+1)/2$ and (b) $u(x, 0) = 1$ if $x > 1$; $u(x, 0) = 0.5$ if $x = 1$; $u(x, 0) = 0$ if $x < 1$. Temporal evolution for the maximum values of continuous (solid line) and discrete (dashed line) representation solutions of the diffusion equation for the (c) smooth initial condition and (d) discontinuous initial condition.

DownLoad: Full-Size Img PowerPoint

From the results of several numerical experiments, we find that the numerical solution of the diffusion equation obtained by using the FSM satisfies the maximum principle if the initial conditions are smooth enough, but otherwise, they may not. Initial conditions that are not smooth enough cause a phenomenon in which continuous Fourier representation does not correspond to discrete representation. That is, the numerical solution of the diffusion equation by using the FSM is bounded to the continuous Fourier representation of the initial condition. Hence, the maximum principle for the numerical solution of the diffusion equation by using the FSM should be used as follows:

for all $t_n$, where $\Vert U(\cdot, t)\Vert_\infty = \sup_{x \in \Omega} \vert U(x, t) \vert$.

4.   Conclusions

In this study, we investigated the discrete maximum principle of the FSM for the diffusion equation. It is well known that the FSM is fast, efficient and accurate. The maximum principle holds for the diffusion equation: A solution satisfying the diffusion equation has the maximum value under the initial condition or on the boundary points. The same result holds for the discrete numerical solution by using the FSM when the initial condition is smoother. However, if the initial condition is not smooth, then we may have an oscillatory profile of a continuous representation of the initial condition in the FSM, which can cause a violation of the discrete maximum principle. We have presented counterexamples of this phenomenon based on computational experiments in which the numerical solutions are bounded by the supremum of the continuous representation of the initial condition. To explain this phenomenon, we decomposed and analyzed the numerical solution into each Fourier mode. We also proposed a maximum principle formula for obtaining the discrete numerical solution of the diffusion equation by using the FSM. Furthermore, investigating the discrete maximum bound principle of the numerical solution of the AC equation by using the FSM will also be interesting future work.

Use of AI tools declaration

The authors have not used artificial intelligence tools in the creation of this article.

Acknowledgments

The first author (J. S. Kim) expresses thanks for the support from the BK21 FOUR program. The authors would like to express sincere thanks to the reviewers for their valuable suggestions and insightful comments to improve the paper.

Conflict of interest

The authors declare that there is no conflict of interest.