A Full Predictor-Corrector Finite Element Method for the One-Dimensional Heat Equation with Time-Dependent Singularities

Jake L. Nkeck

Department of Mathematics, University of Buea, Buea, Cameroon.

DOI: 10.4236/jamp.2024.124084                

Abstract

The energy norm convergence rate of the finite element solution of the heat equation is reduced by the time-regularity of the exact solution. This paper presents an adaptive finite element treatment of time-dependent singularities on the one-dimensional heat equation. The method is based on a Fourier decomposition of the solution and an extraction formula of the coefficients of the singularities coupled with a predictor-corrector algorithm. The method recovers the optimal convergence rate of the finite element method on a quasi-uniform mesh refinement. Numerical results are carried out to show the efficiency of the method.

Keywords

Singularities, Finite Element Methods, Heat Equation, Predictor-Corrector Algorithm

1. Introduction

The heat equation is a partial differential equation known to model problems in thermodynamics, Brownian motions and particle diffusion. In the last seven decades, the increasing availability of powerful computers has been a major propulsion to the implementation of numerical methods for partial differential equations. Concerning the numerical solution of the heat equation, the usual methods are the finite difference methods (see [1] [2] [3] ), the methods of lines (see [4] [5] [6] ) and more recently the space-time finite element methods (FEMs) (see [7] - [12] ).

The singular behavior of solutions of a partial differential equation is known to affect the convergence rate of the numerical approximation. The singularity can be geometric (due to the existence of reentrant corners or edges, the change of boundary conditions), or can be directly related to a singularity of a component of the equation (sign-changing coefficient, component singular at a point...). The goal of this paper is to use classical FEM to solve the one-dimensional heat problem: Find u such that:

$\begin{array}{lll}\frac{\partial u}{\partial t}-k_w^2\frac{{\partial }^{2}u}{\partial x^2}=f\hfill & \hfill & \text{in}\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right]\hfill \\ u\left(\mathrm{0,}t\right)=u\left(L\mathrm{,}t\right)=0\hfill & \mathrm{<mi>\; </mi>}\hfill & t\in \left(\mathrm{0,}T\right)\hfill \\ u\left(x\mathrm{,0}\right)=0\hfill & \mathrm{<mi>\; </mi>}\hfill & x\in \left(\mathrm{0,}L\right)\hfill \end{array}$ (1)

where $L>0$ , $T>0$ , $k_w>0$ and the source term $f\left(x\mathrm{,}t\right)$ is some function singular at a point $t_0\in \left[\mathrm{0,}T\right]$ .

The space-time FEM is known to give an approximation of the solution of Problem (1) with optimal convergence rates, but the method either uses discontinuous Galerkin methods (see [12] ) that increase the number of degrees of freedom, or a particular interpolation error estimate in order to guarantee the optimal convergence (see [7] ). Those techniques are efficient but involve particular adaptations’ techniques.

Singular complement methods and predictor-corrector FEM (PCFEM) have been introduced to solve the problem of geometric singularities of elliptic boundary value problems (see [13] - [23] ) and more recently on geometric singularities of hyperbolic problems (see [24] ). The treatment of time-dependent singularities by means of PCFEM is also introduced in [25] using the method of lines. The main advantage of PCFEM is its application on any standard type of mesh (structured or unstructured) and on usual interpolation operators. This paper presents a PCFEM for the one-dimensional heat Problem (1). The method relies on a Fourier decomposition of the solution where computable formulas of the singular functions are derived and the solution is recomputed with a corrector algorithm in order to cancel the effect of the singularity caused by the source term f.

This paper is organized as follows: Section 2 presents the initial boundary value problem with the function spaces, the existence and uniqueness of the solution, a stability result and the Fourier decomposition of the solution with an extraction formula of the coefficients of the singularity. Section 3 presents the FEM, the predictor-corrector algorithm and some error estimates. Section 4 gives some numerical results that show the efficiency of the method and Section 5 is the conclusion of the paper.

2. Initial Boundary Value Problem

2.1. Function Spaces and Stability Results

In order to study Problem (1), it is crucial to introduce the following spaces:

1) Given a domain $\Omega \subset {ℝ}^n$ , $n=\mathrm{1,2}$ , the space $L^2\left(\Omega \right)$ of classes of Lebesgue-measurable and square-integrable functions over Ω with the norm:

${\Vert u\Vert }_{L^2\left(\Omega \right)}={\left({\displaystyle {\int }_{\Omega }}{\left|u\left(x\right)\right|}^2\text{d}x\right)}^{1/2}$

where $x$ can be a vector or a scalar depending on n.

2) Given an interval $I=\left(a\mathrm{,}b\right)$ and a Hilbert space V, ${\mathcal{C}}^k\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}V\right)$ will denote the space of bounded k-continuously differentiable functions of the form $t\mapsto u\left(\cdot \mathrm{,}t\right)\in V$ with the norm:

${\Vert u\Vert }_{{\mathcal{C}}^k}=\underset{t\in I\mathrm{,}\left|\beta \right|\le k}{\sup }{\Vert D^{\beta }u\left(\cdot \mathrm{,}t\right)\Vert }_V$

where $\beta =\left({\beta }_1\mathrm{,}\cdots \mathrm{,}{\beta }_n\right)\in {\mathbb{N}}^n$ is a multi-index, $\left|\beta \right|={\beta }_1+\cdots +{\beta }_n$ and $D^{\beta }\phi =\frac{{\partial }^{\left|\beta \right|}\phi }{\partial x_1^{{\beta }_1}\cdots \partial x_n^{{\beta }_n}}$ .

3) Given an interval $I=\left(a\mathrm{,}b\right)$ , a Hilbert space V and $0<\alpha <1$ , define the space ${\mathcal{C}}^{k\mathrm{,}\alpha }\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}V\right)$ be the space of k-continuously differentiable functions of the form $t\mapsto u\left(\cdot \mathrm{,}t\right)\in V$ such that $D^{\beta }u\left(\cdot \mathrm{,}t\right)$ , $\left|\beta \right|=k$ , is α-Hölder. The space ${\mathcal{C}}^{k\mathrm{,}\alpha }\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}V\right)$ is equipped with the norm:

${\Vert u\Vert }_{{\mathcal{C}}^{k\mathrm{,}\alpha }}={\Vert u\Vert }_{{\mathcal{C}}^k}+\underset{\left|\beta \right|=k}{\max }\underset{a\le t_1<t_2\le b}{\sup }\frac{{\Vert D^{\beta }u\left(\cdot \mathrm{,}{t}_1\right)-D^{\beta }u\left(\cdot \mathrm{,}{t}_2\right)\Vert }_V}{{\left|t_1-t_2\right|}^{\alpha }}\mathrm{.}$

4) For a given domain Ω, the space $H^m\left(\Omega \right)$ , $m\in \mathbb{N}$ , $m\ge 1$ denotes the space of Lebesgue-measurable classes of functions u such that $u\in L^2\left(\Omega \right)$ and $D^{\beta }u\in L^2\left(\mathrm{0,}L\right)$ for $0<\left|\beta \right|\le m$ , equipped with the norm:

${\Vert u\Vert }_{H^m\left(\Omega \right)}={\left({\displaystyle {\int }_{\Omega }}{\displaystyle \underset{\left|\beta \right|\le m}{\sum }}{\left|D^{\beta }u\right|}^2\text{d}x\right)}^{\frac{1}{2}}\mathrm{.}$

The space $H_0^1\left(\Omega \right)=\left\{u\in H^1\left(\Omega \right)\mathrm{<mo>\; :\; </mo>}u\left(x\right)=0\text{on}\partial \Omega \right\}$ where $\partial \Omega$ represents the boundary of Ω. The space $H_0^1\left(\Omega \right)$ is equipped with the norm of $H^1\left(\Omega \right)$ .

The space $H^{-1}\left(\Omega \right)$ is the dual space of $H_0^1\left(\Omega \right)$ in the sense of distributions.

5) For an interval $I=\left(a\mathrm{,}b\right)$ and a domain Ω, the space $L^2\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\Omega \right)\right)$ is the completion of $\mathcal{C}\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\Omega \right)\right)$ with respect to the norm:

${\Vert u\Vert }_{L^2\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\Omega \right)\right)}={\left({\displaystyle {\int }_a^b}\left({\displaystyle {\int }_{\Omega }}{\displaystyle \underset{\left|\beta \right|\le 1}{\sum }}{\left|D^{\beta }u\left(x\mathrm{,}t\right)\right|}^2\text{d}x\right)\text{d}t\right)}^{\frac{1}{2}}\mathrm{.}$

The space $L^2\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\Omega \right)\right)$ with the norm:

${\Vert u\Vert }_{L^2\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\Omega \right)\right)}={\left({\displaystyle {\int }_a^b}{\Vert u\Vert }_{H^{-1}\left(\Omega \right)}^2\text{d}t\right)}^{\frac{1}{2}}\mathrm{,}$

where

${\Vert u\Vert }_{H^{-1}\mathrm{<mo\; stretchy="false">\; (\; </mo>}\Omega \mathrm{<mo\; stretchy="false">\; )\; </mo>}}=\underset{v\in H_0^1\left(\Omega \right)}{\sup }\left(\frac{\left|u\left(v\right)\right|}{{\Vert v\Vert }_{H_0^1\left(\Omega \right)}}\right)\mathrm{.}$

The space $H^1\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\Omega \right)\right)$ with the norm:

${\Vert u\Vert }_{H^1\left(a\mathrm{,}b\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\Omega \right)\right)}={\left({\displaystyle {\int }_a^b}\left({\Vert u\Vert }_{H^{-1}\left(\Omega \right)}^2+{\Vert \frac{\partial u}{\partial t}\Vert }_{H^{-1}\left(\Omega \right)}^2\right)\text{d}t\right)}^{\frac{1}{2}}\mathrm{.}$

6) The space

$X=\left\{v\in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)\mathrm{<mo>\; :\; </mo>}v\left(x\mathrm{,0}\right)=0\forall x\in \left(\mathrm{0,}L\right)\right\}\mathrm{.}$

7) The space $Y=\left\{v\in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\mathrm{<mo>\; :\; </mo>}v\left(x\mathrm{,0}\right)=0\forall x\in \left(\mathrm{0,}L\right)\right\}$ .

The variational formulation of Problem (1) is: Find $u\in X$ such that:

$a\left(u\mathrm{,}v\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}fv\text{d}x\forall v\in H_0^1\left(\mathrm{0,}L\right)$ (2)

where $a\left(u\mathrm{,}v\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}\frac{\partial u}{\partial t}v+k_w^2\frac{\partial u}{\partial x}\frac{\partial v}{\partial x}\text{d}x\text{d}t$ .

Theorem 1. ( [26] p. 382 (Theorem 11.3))

Assume that $f\in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ . Then, Problem (1) has a unique solution $u\in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ .

Remark. From the Sobolev embedding theorem (see [27] , p. 98]), the space $H^1\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\subset {\mathcal{C}}^{\mathrm{0,}1/2}\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ , so the solution u of Problem (1) is at least $\frac{1}{2}$ -Hölder continuous in time.

Define the linear continuous map named as the Newton potential (see [7] , p. 5]) $\mathcal{N}\mathrm{<mo>\; :\; </mo>}{L}^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)\to L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)$ by $\mathcal{N}\phi$ such that ${\Vert \phi \Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}=k_w^2{\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}\phi \mathcal{N}\phi \text{d}x\text{d}t$ and ${\Vert \mathcal{N}\phi \Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}={\Vert \phi \Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}$ for any $\phi \in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ , if one uses the Cauchy-Schwartz inequality and the Hölder inequality,

$\begin{array}{c}a\left(u,v\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}\frac{\partial u}{\partial t}v+k_w^2\frac{\partial u}{\partial x}\frac{\partial v}{\partial x}\text{d}x\text{d}t\\ \le \left({\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(0,T;H^{-1}\left(0,L\right)\right)}+k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\right){\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\\ \le \left(\sqrt{2}+k_w^2\right){\left({\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(0,T;H^{-1}\left(0,L\right)\right)}^2+{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\right)}^{1/2}\mathrm{<mi>\; </mi>}{\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\\ \le \left(\sqrt{2}+k_w^2\right){\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}{\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}.\end{array}$

Hence, $a\left(\cdot \mathrm{,}\cdot \right)$ is bounded and

$\begin{array}{c}{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}^2={\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2+{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\\ \mathrm{<mo>\; =\; </mo>}{\Vert \mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2+{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\end{array}$ (3)

so if $u\in X$ , one has:

$\begin{array}{c}a\left(u\mathrm{,}u\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}\frac{\partial u}{\partial t}u+k_w^2\frac{\partial u}{\partial x}\frac{\partial u}{\partial x}\text{d}x\text{d}t\\ =\frac{1}{2}{\displaystyle {\int }_0^T}\frac{\text{d}}{\text{d}t}{\Vert u\Vert }_{L^2\left(\mathrm{0,}L\right)}^2\text{d}t+k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\\ =\frac{1}{2}{\Vert u\left(\cdot \mathrm{,}T\right)\Vert }_{L^2\left(\mathrm{0,}L\right)}^2+k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\ge k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\end{array}$

and using the definition of $\mathcal{N}$ , the Cauchy-Schwartz inequality and (3),

$\begin{array}{c}a\left(u\mathrm{,}\mathcal{N}\frac{\partial u}{\partial t}\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}\frac{\partial u}{\partial t}\mathcal{N}\frac{\partial u}{\partial t}+k_w^2\frac{\partial u}{\partial x}\frac{\partial }{\partial x}\left(\mathcal{N}\frac{\partial u}{\partial t}\right)\text{d}x\text{d}t\\ \ge \frac{1}{k_w^2}{\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2\mathrm{<mi>\; </mi>}-k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}{\Vert \mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\\ \ge \frac{1}{k_w^2}{\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2\mathrm{<mi>\; </mi>}-\frac{k_w^2}{2}\left({\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2+{\Vert \mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\right)\\ \mathrm{<mo>\; =\; </mo>}\frac{1}{k_w^2}{\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2-\frac{k_w^2}{2}{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\end{array}$ (4)

Hence, from (3) and (4),

$\begin{array}{c}a\left(u\mathrm{,}u+\mathcal{N}\frac{\partial u}{\partial t}\right)=a\left(u\mathrm{,}u\right)+a\left(u\mathrm{,}\mathcal{N}\frac{\partial u}{\partial t}\right)\\ \ge k_w^2{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2+\frac{1}{k_w^2}{\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2\mathrm{<mi>\; </mi>}-\frac{k_w^2}{2}{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\\ \ge \frac{k_w^2}{2}{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2+\frac{1}{k_w^2}{\Vert \frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)}^2\\ \ge \max \left(\frac{k_w^2}{2}\mathrm{,}\frac{1}{k_w^2}\right){\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}^2\end{array}$

and from (3),

$\begin{array}{c}{\Vert u+\mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\le 2\left[{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2+{\Vert \mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\right]\\ \mathrm{<mo>\; =\; </mo>2}{\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}^2\mathrm{.}\end{array}$

Then,

$\begin{array}{l}a\left(u\mathrm{,}u+\mathcal{N}\frac{\partial u}{\partial t}\right)\\ \ge \frac{1}{\sqrt{2}}\max \left(\frac{k_w^2}{2}\mathrm{,}\frac{1}{k_w^2}\right){\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}{\Vert u+\mathcal{N}\frac{\partial u}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\end{array}$ (5)

and one has the stability condition:

${\Vert u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}\le \underset{v\in L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\mathrm{,}v\ne 0}{\sup }\frac{a\left(u\mathrm{,}v\right)}{{\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\mathrm{.}$ (6)

2.2. Fourier Decomposition of the Solution

Assuming that the solution u of (1) is an element of ${\mathcal{C}}^{\mathrm{0,}\alpha }\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)$ , $1/2\le \alpha <1$ , and has a line singularity at a point $t_0\in \left[\mathrm{0,}T\right]$ , based on the boundary conditions one can set $u\left(x\mathrm{,}t\right)={\displaystyle {\sum }_{k=1}^{\infty }}{v}_k\left(t^{\alpha }\right)\sin \left(k\pi x/L\right)\mathrm{<mo>\; :\; </mo>}={\displaystyle {\sum }_{k=1}^{\infty }}{\gamma }_k\left(t\right)\sin \left(k\pi x/L\right)$ and obtain use separation of variables to obtain the following result.

Theorem 2. The solution u of Problem (1) can be written as $u\left(x,t\right)={\displaystyle {\sum }_{k=1}^{\infty }}{\gamma }_k\left(t\right)s_k\left(x\right)$ where $s_k\left(x\right)=\sin \left(k\pi x/L\right)$ , and setting ${\beta }_k=\frac{k\pi }{L}{k}_w$ ,

1) If $t_0=0$ ,

${\gamma }_k\left(t\right)=\left[{\displaystyle {\int }_0^{t^{\alpha }}}\frac{1}{\alpha }{f}_k\left({\tau }^{1/\alpha }\right){\tau }^{\frac{1-\alpha }{\alpha }}{\text{e}}^{{\beta }_k^2{\tau }^{1/\alpha }}\text{d}\tau \right]{\text{e}}^{-{\beta }_k^{2}t}$ (7)

and $\frac{\text{d}{\gamma }_k}{\text{d}t}\left(t\right)=f_k\left(t\right)-{\beta }_k^2{\gamma }_k\left(t\right)$ .

2) If $t_0>0$ ,

${\gamma }_k\left(t\right)=\left[{\displaystyle {\int }_{{\left(t_0-t\right)}^{\alpha }}^{t^{\alpha }}}\frac{1}{\alpha }{f}_k\left(t_0-{\tau }^{1/\alpha }\right){\tau }^{\frac{1-\alpha }{\alpha }}{\text{e}}^{-{\beta }_k^2{\tau }^{\mathrm{1\; <mo>\; /\; </mo>}\alpha }}\text{d}\tau \right]{\text{e}}^{{\beta }_k^2\left(t_0-t\right)}$ (8)

and $\frac{\text{d}{\gamma }_k}{\text{d}t}\left(t\right)=f_k\left(t\right)-{\beta }_k^2{\gamma }_k\left(t\right)$ .

3. Finite Element Discretization

3.1. Finite ELement Method

Let ${\mathcal{T}}_h$ be a triangulation of the domain $\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)$ with quasi-uniform meshes of diameter h and define:

$X_h=\left\{v\in X\mathrm{<mo>\; :\; </mo>}{v}_{\mathrm{<mo>\; |\; </mo>}K}\in {\mathcal{P}}^1\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)\forall K\in {\mathcal{T}}_h\right\}$

with the norm of $L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ .

The Galerkin approximation consists to find $u_h\in X_h$ such that:

$a\left(u_h\mathrm{,}{v}_h\right)={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}f{v}_h\text{d}x\text{d}t\forall v\in X_h\mathrm{.}$ (9)

The variational Problem (9) has a unique solution according to Theorem 1 and since $X_h\subset X\subset L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^{-1}\left(\mathrm{0,}L\right)\right)$ , for any $u_h\in X_h$ , from the inequality (5),

$\begin{array}{c}{\Vert u_h\Vert }_{X_h}={\Vert u_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)\cap H^1\left(0,T;H^{-1}\left(0,L\right)\right)}\\ \le C\frac{a\left(u_h,u_h+\mathcal{N}\frac{\partial u_h}{\partial t}\right)}{{\Vert u_h+\mathcal{N}\frac{\partial u_h}{\partial t}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\\ =C\frac{a\left(u_h,u_h+\mathcal{N}\frac{\partial u_h}{\partial t}\right)}{{\Vert u_h+\mathcal{N}\frac{\partial u_h}{\partial t}\Vert }_{X_h}}\\ \le C\underset{v_h\in X_h,v_h\ne 0}{\sup }\frac{a\left(u_h,v_h\right)}{{\Vert v_h\Vert }_{X_h}}\end{array}$ (10)

where $C>0$ is a constant depending on $k_w$ .

Furthermore, one has the Galerkin orthogonality:

$a\left(u-u_h\mathrm{,}{v}_h\right)=0\forall v_h\in X_h$ (11)

The following theorem is a lemma equivalent to the Céa’s Lemma known for elliptic boundary value problems. This theorem has been proved in [7] in a more general case.

Theorem 3. Let u be the solution of the variational Problem (2) and $u_h\in X_h$ be the solution of the discretized variational Problem (9). Then, there exists $C>0$ depending on $k_w$ such that:

${\Vert u-u_h\Vert }_{X_h}\le C\underset{w_h\in X_h}{\inf }{\Vert u-w_h\Vert }_X\mathrm{.}$

Proof. From the discrete stability condition (10) and the Galerkin orthogonality (11), for any $w_h\in X_h$ ,

$\begin{array}{c}{\Vert u_h-w_h\Vert }_{X_h}\le C\underset{v_h\in X_h}{\sup }\frac{a\left(u_h-w_h,v_h\right)}{{\Vert v_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\\ \le C\left(\underset{v_h\in X_h}{\sup }\frac{a\left(u_h-u,v_h\right)}{{\Vert v_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}+\underset{v_h\in X_h}{\sup }\frac{a\left(u-w_h,v_h\right)}{{\Vert v_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\right)\\ =C\underset{v_h\in X_h}{\sup }\frac{a\left(u-w_h,v_h\right)}{{\Vert v_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\\ \le C\underset{v\in X}{\sup }\frac{a\left(u-w_h,v\right)}{{\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\\ \le C{\Vert u-w_h\Vert }_X.\end{array}$

Then, ${\Vert u-u_h\Vert }_{X_h}\le {\Vert u-w_h\Vert }_{X_h}+{\Vert w_h-u_h\Vert }_{X_h}\le \left(1+C\right){\Vert u-w_h\Vert }_X$ for any $w_h\in X_h$ where $C>0$ is an arbitrary constant depending on $k_w$ .

Hence, ${\Vert u-u_h\Vert }_{X_h}\le C\underset{w_h\in X_h}{\inf }{\Vert u-w_h\Vert }_X$ .£

The following theorem proves the energy error estimate when the solution is smooth. This theorem has been proved in [7] in a more general case.

Theorem 4. Let u be the unique solution of Problem (1) and $u_h$ be the unique solution of the equivalent discretized variational Problem (9). Assume $u\in H^2\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)$ , then

${\Vert u-u_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\le Ch{\Vert u\Vert }_{H^2\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)}\mathrm{.}$

Proof.

${\Vert u-u_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}={\Vert u-u_h\Vert }_X\le C\underset{w_h\in X_h}{\inf }{\Vert u-w_h\Vert }_X\le C{\Vert u-{\Pi }_{h}u\Vert }_X$

where ${\Pi }_{h}u\in X_h$ is the unique solution of the variational problem,

$a\left({\Pi }_{h}u\mathrm{,}{v}_h\right)=a\left(u\mathrm{,}{v}_h\right)\forall v_h\in X_h\mathrm{.}$

Note that for any $v\in H^1\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)$ such that $v\left(x\mathrm{,0}\right)=0$ and $v\left(x\mathrm{,}t\right)=0$ for $x=0$ or $x=L$ ,

${\Vert v\Vert }_{H^1\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)}^2={\displaystyle {\int }_0^T}{\displaystyle {\int }_0^L}{\left|v\right|}^2+{\left|\frac{\partial v}{\partial t}\right|}^2+{\left|\frac{\partial v}{\partial x}\right|}^2\text{d}x\text{d}t\ge {\Vert v\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}^2\mathrm{.}$

Then,

$\begin{array}{c}{\Vert u-u_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\le C{\Vert u-{\Pi }_{h}u\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\\ \le C{\Vert u-u_h\Vert }_{H^1\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)}\\ \le C\underset{w_h\in X_h}{\inf }{\Vert u-v_h\Vert }_{H^1\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)}\\ \le Ch{\Vert u\Vert }_{H^2\left(\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)\right)}.\end{array}$

£

3.2. Predictor-Corrector Algorithm

The predictor-corrector method presented is based on the decomposition of the approximate solution $u_{hN}\left(x,t\right)=w_{hN}\left(x,t\right)+{\displaystyle {\sum }_{k=1}^N}{\gamma }_k\left(t\right)s_k\left(x\right)$ where $u_{hN}$ is the solution of (9), $s_k$ and ${\gamma }_k$ are defined as in Theorem 2.

Predictor-Corrector Algorithm

1) Solve the problem: Find $u_h$ such that:

$\begin{array}{lll}\frac{\partial u_h}{\partial t}-k_w^2\frac{{\partial }^2{u}_h}{\partial x^2}=f\hfill & \hfill & \text{in}\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right]\hfill \\ u_h\left(\mathrm{0,}t\right)=u_h\left(L\mathrm{,}t\right)=0\hfill & \mathrm{<mi>\; </mi>}\hfill & t\in \left[\mathrm{0,}T\right]\hfill \\ u_h\left(x\mathrm{,0}\right)=0\hfill & \mathrm{<mi>\; </mi>}\hfill & x\in (\; 0,\; L\; )\hfill \end{array}$

2) For $N\in \mathbb{N}$ fixed, do:

a) For $k=\mathrm{1,}\cdots \mathrm{,}N$ , compute ${\gamma }_k\left(t\right)$ and $\frac{\text{d}{\gamma }_k}{\text{d}t}\left(t\right)$ using Formula (7) or (8).

b) If

$X_{1h}\mathrm{<mo>\; :\; </mo>}=X_h+\text{Span}\left\{s_k\mathrm{,\; <mi>\; </mi>}k=\mathrm{1,}\cdots \mathrm{,}N\right\}\mathrm{,}$

solve the problem: Find $w_{hN}\in X_{1h}$ such that:

$\begin{array}{ll}\frac{\partial w_{hN}}{\partial t}-k_w^2\frac{{\partial }^2{w}_{hN}}{\partial x^2}=f-{\displaystyle \underset{k=1}{\overset{N}{\sum }}}\frac{\text{d}{\gamma }_k}{\text{d}t}{s}_k-k_w^2{\gamma }_k\frac{{\text{d}}^2{s}_k}{\text{d}{x}^2}\hfill & \text{in}\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right]\hfill \\ w_{hN}\left(\mathrm{0,}t\right)=-{\displaystyle \underset{k\mathrm{<mo>\; =\; </mo>1}}{\overset{N}{\sum }}}{\gamma }_{k\tau }\left(t\right)s_k\left(0\right)\hfill & t\in \left[\mathrm{0,}T\right]\mathrm{<mi>\; </mi>}\hfill \\ w_{hN}\left(L\mathrm{,}t\right)=-{\displaystyle \underset{k\mathrm{<mo>\; =\; </mo>1}}{\overset{N}{\sum }}}{\gamma }_k\left(t\right)s_k\left(L\right)\hfill & t\in \left[\mathrm{0,}T\right]\mathrm{<mi>\; </mi>}\hfill \\ w_{hN}\left(x\mathrm{,0}\right)=-{\displaystyle \underset{k=1}{\overset{N}{\sum }}}{\gamma }_k\left(0\right)s_k\left(x\right)\hfill & x\in \left[\mathrm{0,}L\right]\mathrm{<mi>\; </mi>}\hfill \end{array}$

c) Compute $u_{hN}\left(x\mathrm{,}t\right)=w_{hN}\left(x\mathrm{,}t\right)+{\displaystyle {\sum }_{k=1}^N}{\gamma }_k\left(t\right)s_k\left(x\right)$ .

3.3. Error Estimates

From [28] , if $u\in {\mathcal{C}}^{\mathrm{0,}\alpha }\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)$ , $0<\alpha <1$ , then for any $N\in \mathbb{N}$ , $N\ge 1$ ,

${\Vert u-{\displaystyle \underset{k=1}{\overset{N}{\sum }}}{\gamma }_k{S}_k\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\le C{N}^{-\alpha }\ln N$ (12)

where ${\gamma }_k$ and $s_k$ are defined in Theorem 2 and $C>0$ is an arbitrary constant.

The convergence result in the energy norm is given by the following theorem.

Theorem 5. Let $u\in {\mathcal{C}}^{\mathrm{0,}\alpha }\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}^1\left(\mathrm{0,}L\right)\right)$ , $1/2\le \alpha <1$ be the solution of Problem (1) and let $u_{hN}$ , $N\in \mathbb{N}$ , $N\ge 1$ be the approximate solution of u using the Predictor-Corrector Algorithm 3.2 on a quasi-uniform and triangulation ${\mathcal{T}}_h$ of the domain $\left(\mathrm{0,}L\right)\times \left(\mathrm{0,}T\right)$ . Assume further that $N^{-\alpha }\ln N\le h$ , then there is $C>0$ such that:

${\Vert u-u_h\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\le Ch$

Proof. Suppose $u\left(x,t\right)=w_N\left(x,t\right)+{\displaystyle {\sum }_{k=1}^N}{\gamma }_k\left(t\right)s_k\left(x\right)$ is the decomposition of the solution u as in Theorem 2. Set $u_{hN}\left(x,t\right)=w_{hN}+{\displaystyle {\sum }_{k=1}^N}{\gamma }_k\left(t\right)s_k\left(x\right)$ be the approximation of u using the Predictor-Corrector Algorithm 3.2.

From Theorem 4, (12) and the assumption $N^{-\alpha }\ln N\le h$ , one has:

${\Vert u-u_{hN}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}={\Vert w_N-w_{hN}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}\le C{N}^{-\alpha }\ln N\le Ch\mathrm{.}$

£

4. Numerical Experiments

The main purpose of this section is to present some computational results carried out with the Predictor-Corrector Algorithm 3.2. Even though the examples considered may not represent any real world situation, the difficulties due to singularities in the solution are well reflected. We consider two examples with known solutions. All computations are carried out using a code developped in python. The resulting linear systems are solved iteratively with a generalized minimal residual method (GMRES) and all the integrals are computed using a 15-point Legendre Gauss formula on the stifness, mass matrices and load vectors, and a 150-points Legendre Gauss formula for integrals used to compute the Fourier coefficients.

In all the computations, starting from a triangulation ${\mathcal{T}}_{h_0}$ , the refinement are perfomed by dividing the concerned lentgh by two ( $h_{i+1}=h_i/2$ ). The order or rate of convergence $\nu$ is computed using the formula:

$\nu ={\log }_2\left(\frac{{\Vert u-u_{h_i}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}{{\Vert u-u_{h_{i+1}}\Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}}\right)\mathrm{.}$

In all the following examples, ${\nu }_{exp}$ is the expected order of convergence. Figure 1 presents the type of mesh refinements used in the examples.

4.1. Example 1

In this example, we consider Problem (9) with the exact solution $u\left(x\mathrm{,}t\right)=t^{1/2}x\left(x-1\right)$ and $f=\frac{\partial u}{\partial t}-\frac{1}{144}\frac{{\partial }^{2}u}{\partial x^2}$ , $T=1$ , $L=1$ .

The function has a line singularity at $t=0$ and belongs to ${\mathcal{C}}^{\mathrm{0,}1/2}\left(\mathrm{0,1,}{H}^1\left(\mathrm{0,1}\right)\right)$ . Table 1 presents the error estimates in the energy norm ${\Vert \cdot \Vert }_{L^2\left(\mathrm{0,}T\mathrm{<mo>\; ;\; </mo>}{H}_0^1\left(\mathrm{0,}L\right)\right)}$ and the required number of Fourier coefficients N in order to obtain an optimal convergence

Figure 1. The uniform meshes: Refinements 1, 2 and 3

Table 1. Error estimates in the energy norm.

rate. Note that we need $N^{-1/2}\ln N\le h_j$ at each refinement $j=\mathrm{1,}\cdots \mathrm{,5}$ . A convergence of the classical FEM is observed and the convergence rate that is greater than 0.5, which is better than the expected convergence rate 0.25. However the convergence rate of the classical FEM is still not optimal. From the results obtained with the PCFEM, the optimal convergence rate is recovered.

Figures 2-4 present the plots of the exact solution, the approximate solution before adaptation and the aproximate solution using the Predictor-Corrector Algorithm 3.2 respectively at the fixed times $t=1/4$ , $t=1/2$ , $t=3/4$ and $t=1$ . Figures 5-7 present the 3-dimensional line countour plots of the exact solution, the approximate solution before adaptation and the approximate solution after using the predictor-corrector adaptation (PCFEM) with 75 lines. It is observed that the predictor-corrector solution has a better approximation with a least number of nodes in the refinement.

4.2. Example 2

In this example, we consider Problem (9) with the exact solution $u\left(x\mathrm{,}t\right)={\left(t-1\right)}^{3/4}x\left(x-1\right)$ and $f=\frac{\partial u}{\partial t}-\frac{1}{400}\frac{{\partial }^{2}u}{\partial x^2}$ , $T=1$ , $L=1$ .

The solution has a line singularity at $t_0=1$ and is an element of ${\mathcal{C}}^{\mathrm{0,}3/4}\left(\mathrm{0,1,}{H}^1\left(\mathrm{0,1}\right)\right)$ so the number of Fourier coefficients at each refinement on a triangulation of diameter h should satisfy $N^{-3/4}\ln N\le h$ . Figure 8 shows the

Figure 2. The exact solution, the FEM solution and the PCFEM solution at Refinement 3 ( $h=2^{-3}$ ).

Figure 3. The exact solution, the FEM solution and the PCFEM solution at Refinement 4 ( $h=2^{-4}$ ).

Figure 4. The exact solution, the FEM solution and the PCFEM solution at Refinement 5 ( $h=2^{-5}$ ).

Figure 5. 3D plots of the exact solution, the FEM solution and the PCFEM solution at Refinement 3.

Figure 6. 3D plots of the exact solution, the FEM solution and the PCFEM solution at Refinement 4.

Figure 7. 3D plots of the exact solution, the FEM solution and the PCFEM solution at Refinement 5.

Figure 8. The graph of $\ln {\Vert u-u_{h_j}\Vert }_{X_h}$ against $-\ln h_j$ , $h_j=2^{-j}$ , $j=\mathrm{1,}\cdots \mathrm{,5}$ .

error estimates of the FEM and the PCFEM. The convergence rate of the classical FEM is not optimal while the one of the PCFEM is optimal.

5. Conclusions

This paper presents an adaptive predictor-corrector finite element method (PCFEM) for the one-dimensional heat equation with zero initial and boundary conditions. The method is based on a Fourier decomposition with a correction of the solution by a predictor-corrector algorithm. The PCFEM is shown to be efficient on quasi-uniform meshes and then can be used without a particular mesh treatment or a special interpolation error in the energy norm. In the predictor-corrector algorithm, the stiffness matrices are the same, only the load vector is modified to correct the solution, and then the complexity of the PCFEM algorithm and the classical FEM algorithm are of the same order.

The PCFEM developed can be applied to the heat equation with non-zero initial and boundary conditions, and the proof of the stability conditions and the Homogenization methods can be adapted using the proofs developed in [7] provided that the initial and boundary conditions are regular enough. It is worth saying that the method developed is another motivation to continue using standard FEM to solve partial differential equations.

List of Notations

The following are some variables used in the paper.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Fletcher, C. (1988) Computational Techniques for Fluid Dynamics. Springer-Verlag, Berlin. https://doi.org/10.1007/978-3-642-97071-9
[2]	Maturi, D., Alsulami, N. and Alaidarous, E. (2020) Finite Difference Approximation for Solving Transient Heat Conduction Equation of Copper. Advances in Pure Mathematics, 10, 350-358. https://doi.org/10.4236/apm.2020.105021
[3]	Ali, M., Loskor, W., Taher, S. and Bilkis, F. (2022) Solution of a One-Dimension Heat Equation Using Higher-Order Finite Difference Methods and Their Stability. Journal of Applied Mathematics and Physics, 10, 877-886. https://doi.org/10.4236/jamp.2022.103060
[4]	Thomée, V. (2006) Galerkin Finite Element Methods for Parabolic Problems. Springer-Verlag, Berlin.
[5]	Jung, M. and Langer, U. (2013) Methode der Finiten Elemente für Ingenieure. Springer Science Business Media, Wiesbaden. https://doi.org/10.1007/978-3-658-01101-7
[6]	Dautray, R. and Lions, J.-L. (1999) Mathematical Analysis and Numerical Methods for Science and Technology: Evolution Problems II. Vol. 6, Springer-Verlag, Berlin.
[7]	Steinbach, O. (2015) Space-Time Finite Element Methods for Parabolic Problems. Computational Methods in Applied Mathematics, 15, 551-566. https://doi.org/10.1515/cmam-2015-0026
[8]	Schwab, C. and Stevenson, R. (2009) Space-Time Adaptive Wavelet Methods for Parabolic Evolution Problems. Mathematics of Computation, 78, 1293-1318. https://doi.org/10.1090/S0025-5718-08-02205-9
[9]	Neumüller, M. (2013) Space-Time Methods: Fast Solvers and Applications. Ph.D. Thesis, Graz University of Technology, Graz.
[10]	Mollet, C. (2014) Stability of Petrov-Galerkin Discretizations: Application to the Space-Time Weak Formulation for Parabolic Evolution Problems. Computational Methods in Applied Mathematics, 14, 231-255. https://doi.org/10.1515/cmam-2014-0001
[11]	Gander, J. and Neumller, M. (2014) Analysis of a New Space-Time Parallel Multigrid Algorithm for Parabolic Problems. arXiv: 1411.0519. http://arxiv.org/abs/1411.0519
[12]	Hughes, T. and Hulbert, G. (1988) Space-Time Finite Element Methods for Elastodynamics: Formulations and Error Estimates. Computer Methods in Applied Mechanics and Engineering, 66, 339-363. https://doi.org/10.1016/0045-7825(88)90006-0
[13]	Hazard, C. and Lenoir, M. (1996) On the Solution of Time-Harmonic Scattering Problems for Maxwell’s Equations. SIAM Journal of Mathematical Analysis, 27, 1597-1630. https://doi.org/10.1137/S0036141094271259
[14]	Bonnet-Ben Dhia, A., Hazard, C. and Lorenghel, S. (1999) A Singular Field Method for the Solution of Maxwell’s Equations in Polyhedral Domains. SIAM Journal of Applied Mathematics, 59, 2028-2044. https://doi.org/10.1137/S0036139997323383
[15]	Hazard, C. and Lorenghel, S. (2002) A Singular Field Method for the Solution of Maxwell’s Equations: Numerical Aspects for 2D Magnetostatics. SIAM Journal of Applied Mathematics, 40, 1021-1040. https://doi.org/10.1137/S0036142900375761
[16]	Assous, F., Ciarlet, P., Labrunie, S. and Lorenghel, S. (2001) The Singular Complement Method. In: Debit, N., Garbey, M., Hoppe, R., Périaux, J., Keyes, D. and Kuznetsov, Y., Eds., Thirteenth International Conference on Domain Decomposition Methods, Institut Elie Cartan, Villers-lès-Nancy, 159-187.
[17]	Nkemzi, B. and Jung, M. (2013) A Postprocessing Finite Element Strategy for Poisson’s Equation in Polygonal Domains: Computing the Stress Intensity Factors. In: Apel, T. and Steinbach, O., Eds., Advanced Finite Element Methods and Applications, Springer, Berlin, 153-173. https://doi.org/10.1007/978-3-642-30316-6_7
[18]	Nkemzi, B. (2006) On the Solution of Maxwell’s Equations in Polygonal Domains. Mathematical Methods in Applied Science, 29, 1053-1080. https://doi.org/10.1002/mma.717
[19]	Nkemzi, B. (2007) On Singularities of Solution of Maxwell’s Equations in Axisymmetric Domains with Conical Points. Mathematical Methods in Applied Science, 30, 877-888. https://doi.org/10.1002/mma.812
[20]	Nkemzi, B. (2005) On the Solution of Maxwell’s Equations in Axisymmetric Domains with Edges. Zeitschrift für Angewandte Mathematik und Mechanik, 85, 570-592. https://doi.org/10.1002/zamm.200310188
[21]	Nkemzi, B. (2016) On the Coefficients of the Singularities of the Solution of Maxwell’s Equations near Polyhedral Edges. Mathematical Problems in Engineering, 2016, Article ID: 7965642. https://doi.org/10.1155/2016/7965642
[22]	Nkemzi, B. and Tanekou, S. (2018) Predictor-Corrector p-and hp-Versions of the Finite Element Method for Poisson’s Equation in Polygonal Domains. Computer Methods in Applied Mechanics and Engineering, 333, 74-93. https://doi.org/10.1016/j.cma.2018.01.027
[23]	Nkemzi, B. and Nkeck, J. (2020) A Predictor-Corrector Finite Element Method for Maxwell’s Equations in Polygonal Domains. Mathematical Problems in Engineering, 2020, Article ID: 3502513. https://doi.org/10.1155/2020/3502513
[24]	Nkemzi, B. and Jung, M. (2023) The Coefficients in the Asymptotic Expansion of Solutions of Second-Order Hyperbolic Problems in Polygonal Domains. Mathematical Problems in Engineering, 46, 9576-9588. https://doi.org/10.1002/mma.9075
[25]	Nkeck, J. (2024) A Crank-Nicolson Finite Element Treatment of Time-Dependent Singularities of the One Dimensional Heat Equation. https://doi.org/10.21203/rs.3.rs-3989135/v1
[26]	Renardy, M. and Rogers, R. (2004) An Introduction to Artial Differential Equations. 2nd Edition, Springer-Verlag, New York.
[27]	Adams, R. (1975) Sobolev Spaces. Academic Press, New York.
[28]	Jackson, D. (1930) The Theory of Approximation. Vol. XI, AMS Colloquium Publication, New York.