Stability for discrete time waveform relaxation methods based on Euler schemes

1.   Introduction

Consider the initial problem

Its waveform relaxation (WR) methods based on Euler methods can be constructed as follows:

(ⅰ) Taking the splitting function $F(t, x, x) = f(t, x)$, construct the iterative scheme

with $k\in \mathbb{Z}^+$ and $x^{(0)}\equiv x_0$, where $\mathbb{Z}^+$ denotes all nonnegative integers;

(ⅱ) Discretizing (1.2) by employing the explicit and implicit Euler methods, one arrives at the expected methods

and

where the step size $h$ satisfies the following: $Nh = T$, mesh points $t_n = nh$ and $x^{(0)}_n = x_0$ for any $n$.

The WR method was introduced for the first time by Lelarasmee et al. ^[1] for the time domain analysis of large-scale nonlinear dynamical systems. Its two advantages having a multirate property and inherent parallelism, can be observed, making it quite competitive with the classical methods based, for examples, on discrete variable methods such as the Runge-Kutta, linear multistep or predictor-correction methods for differential systems ^[2,3,4]. The convergence of WR methods has been studied for various types of ordinary differential equations (ODEs). Here we only mention the references ^{[2,5,6,7,8,9,10]}.

However, a convergent WR method may be impractical. To confirm this a convergent and unstable example of WR methods is given as follows:

Applying (1.4) to the equation

one arrives at the WR method given by

This method is convergent, that is $\lim\limits_{k\rightarrow \infty, N\rightarrow \infty} \max\limits_{0\leqslant n\leqslant T/N}|x^{(k)}_n-x(t_n)| = 0$. For simplicity we assume that all calculations associated with (1.5) are exact and only initial values have tiny perturbations. Let $\{\varepsilon^{(k)}_n\}$ denote these perturbations; then, the resulting perturbed solution $\{\tilde{x}^{(k)}_n\}$ satisfies

In Example 3.9 of Section 3 we will show that the differences $x^{(k)}_n-\tilde{x}^{(k)}_n, k, n\in\mathbb{Z}^+$ may be unbounded for any $h$ when $\alpha = 1$. This means that the tiny perturbations of initial values may lead to a huge change in the solution for the above WR method. A method with a similar property is called unstable and is clearly impractical. Hence, one should study the stability property of WR methods. Two simplified approaches have been used to understand partly the stability of WR methods. One approach studies stability of the numerical method generated by letting the iteration index of a WR method go to infinity ^[11,12]. Another approach studies the stability of the time-point relaxation method, a variant of WR in which each window size equals the step size used in numerical integration ^[11,13]. In ^[14] a very different approach is adopted, which investigates whether or not the methods obtained by applying some standard discretized time WR methods to a dissipative system can preserve the contractivity properties of the system.

Now, there exist few studies on the stability properties of WR methods except for the aforementioned references. Thus, further investigation is necessary. In this paper we study the stability of WR methods by using a similar approach as that used in ^[14]. First, apply the WR methods to some stable test equations. Second, study under what conditions the WR methods can preserve the stability of exact solutions. Last, explore the key factors for determining the stability of WR methods. In fact, this approach has been widely used to study the stability of classical numerical methods of ODEs. Being different from ^[14] in this paper we will consider more general stability properties than contractivity, use more splittings than the three special splittings used in ^[14] and derive some new interesting results. Based on these stability results, it is expected that the discrete-time WR methods can be applied to some significant inverse problems of mathematical physics ^{[15,16,17,18]}, and we also refer the reader to ^{[19,20,21,22,23,24,25]} for more related discussions.

The organization of this paper is as follows. In Section 2, we will present the definition of the stability of WR methods and prove two useful lemmas. Section 3 investigates the stability of the WR methods applied to (2.1) and explores the key factors impacting the stability of WR methods. In Section 4, we study the stability of the WR methods applied to (2.2) and try to provide some stable conditions. At last, Section 5 presents some numerical experiments to illustrate the theories obtained.

2.   Preliminary

In this section we provide some stability definitions for the WR methods and two lemmas for the study of stability in the next two sections.

Consider the following two dissipative systems:

with $A\in \mathbb{R}^{{{\text{d}}} \times {{\text{d}}}}$ and $\max\limits_i\Re(\lambda_i(A)) < 0$ ($\Re(\lambda_i(A))$ denoting the real part of eigenvalues $\lambda_i$ of $A$), and

where $f$ satisfies the one-sided Lipschitz condition

Here, $x_1, x_2\in \mathbb{R}^{{{\text{d}}}}$ and $c > 0$. Here and hereafter $\|\cdot\|$ denotes, with the exception of Section 4, a vector norm or its induced matrix norm.

Assume that all calculations in the WR methods are exact with the exception of initial values and let $\{\varepsilon^{(0)}_n, \varepsilon^{(k)}_0, n, k \in \mathbb{Z}^+\}$ be the perturbations of initial values. This is a classical assumption for an investigation into the stability of numerical methods for ODEs. Let $\{\{x^{(k)}_n, n\in \mathbb{Z}^+\}, k\in \mathbb{Z}^+\}$ denote the approximate solution generated by a WR method with the initial values $\{x^{(0)}_n, x^{(k)}_0, n, k\in \mathbb{Z}^+\}$, and let $\{\{\tilde{x}^{(k)}_n, n\in \mathbb{Z}^+\}, k\in \mathbb{Z}^+\}$ denote the perturbed solution generated by the perturbed system obtained by replacing the initial values of the above WR method by $\{x^{(0)}_n+\varepsilon^{(0)}_n, x^{(k)}_0+\varepsilon^{(k)}_0, n, k\in \mathbb{Z}^+\}$.

In general, we say that a numerical method is stable if the differences between the approximate solution and its perturbed solution are controllable when it is applied to a dissipative system.

Definition 2.1. A WR method is called stable if the approximate solution $\{\{x^{(k)}_n, n\in \mathbb{Z}^+\}, k\in \mathbb{Z}^+\}$ and its perturbed solution $\{\{\tilde{x}^{(k)}_n, n\in \mathbb{Z}^+\}, k\in \mathbb{Z}^+\}$ generated by applying it to the dissipative systems (2.1) or (2.2) satisfy that

The WR method is called contractive if it is stable and $C < 1$. The WR method is called asymptotically stable if it is stable and satisfies

that is

Lemma 2.2. Suppose that $a\geqslant 0, b\geqslant 0, c\geqslant 0, \rho = a+b+c < 1$ and the sequence of positive numbers $\{u^{(k)}_n, n, k\in \mathbb{Z}^+\}$ satisfy

Then,

for all $i\geqslant 1$.

Proof. Let $M$ denote $\sup\limits_{j\geqslant 0}\left\{\max\left\{u^{(0)}_j, u^{(j)}_0\right\}\right\}$. We will show that (2.5) holds for all $i\geqslant 1$ by the induction. Our proof is divided into two steps.

Step 1: Show that (2.5) holds given that $i = 1$, that is

for all $j\geqslant 1$. For this the mathematical induction is used again. By virtue of (2.4) one can derive easily that (2.6) holds for $j = 1$. Suppose that (2.6) holds for the index $j$. Then, it is enough to show that (2.6) also holds when $j$ is replaced with $j+1$. By (2.4) and (2.6), we have

and

Hence, (2.6) is true when $j$ is replaced with $j+1$.

Step 2: Assume that (2.5) holds for the index $i$, and we will show that it still holds if $i$ is replaced by $i+1$, that is

for all $j\geqslant i+1$. By (2.4) and (2.5) one can derive easily that

Hence, (2.7) is true for $j = i+1$. Now assume that (2.7) is true for the index $j$. Then, it is enough to show that (2.7) holds when $j$ is replaced by $j+1$ by the induction. By virtue of (2.4), (2.5) and (2.7), we can obtain that

and

which imply that (2.7) is true when $j$ is replaced with $j+1$.

Consequently, Steps 1 and 2 show that (2.5) holds for all $i\geqslant 1$ by the induction. The proof is complete. □

Lemma 2.3. Suppose that $a\geqslant 0, b\geqslant 0, c\geqslant 0, \rho = a+b+c\geqslant 1$ and the sequence of positive numbers $\{u^{(k)}_n, n\in \mathbb{Z}^+, k\in \mathbb{Z}^+\}$ satisfy

Then,

for all $i\geqslant 1$.

Proof. The proof of the lemma is similar to that of Lemma 2.2, so we omit it. □

3.   Stability of WR methods for Eq (2.1)

Choosing the splitting function $F(t, x, x) = A_1x+A_2x = Ax$ and applying the WR methods (1.3) and (1.4) to Eq (2.1) one arrives at

and

It is easy to give the sufficient conditions for the stability of WR methods (3.1) and (3.2) by using Lemma 2.2.

Theorem 3.1. The WR method (3.1) is contractive and asymptotically stable if $\|I+hA_1\|+\|hA_2\| < 1$.

Proof. Let $\{x^{(k)}_n\}$ be the approximate solution generated by (3.1) and $\{\varepsilon^{(0)}_n, \varepsilon^{(k)}_0\}$ be the perturbations of the initial values. Then, the perturbation solution $\{\tilde{x}^{(k)}_n\}$ caused by $\{\varepsilon^{(0)}_n, \varepsilon^{(k)}_0\}$ satisfies

Let $e^{(k)}_n$ denote $\|x^{(k)}_n-\tilde{x}^{(k)}_n\|$. By virtue of (3.1) and (3.3) we have

This together with Lemma 2.2, proves the theorem. □

Theorem 3.2. The WR method (3.2) is contractive and asymptotically stable if $\|(I-hA_1)^{-1}\|+\|(I-hA_1)^{-1}hA_2\| < 1$.

Proof. The proof is similar to that of Theorem 3.1, so we omit it. □

Now we try to derive some interesting results of stability for the WR methods (3.1) and (3.2) by using Theorems 3.1 and 3.2. We call the WR methods (3.1) and (3.2) the block Gauss-Jacobi (BGJ) WR methods if

and the block Gauss-Seidel (BGS) WR methods if

where $A_{ii}, i = 1, 2, \cdots, s$ denotes the square matrices, and $\tilde{a}_{ij}$, the $i$-th line and $j$-th row of $A_1$, is equal to $a_{ij}$ when $\tilde{a}_{ij}\neq 0$. The BGJ and BGS splitting may be the most common splitting. Bellen and his coauthors in ^[14] have studied the stability of WR methods with the Gauss-Jacobi (GJ) and Gauss-Seidel (GS) splittings, and their study shows that the methods are stable when $A$ is diagonally dominant with negative diagonal elements; however, the techniques developed in ^[14] may be unsuitable to study the stability of WR methods with the BGJ and BGS splitting. The following two corollaries show that the WR methods are stable when using the splitting matrix satisfying some conditions including the BGJ and BGS splittings.

Lemma 3.3. (^[26]) Let $A = (a_{ij})_{d\times d}$ and $r_i(A) = \sum\limits_{i\neq j}|a_{ij}|$. If $A$ is diagonally dominant, then,

Here, and throughout, $\|\cdot\|_{\infty}$ denotes the maximum norm in $\mathbb{C}^\mathrm{d}$ or the maximum row sum matrix norm.

Corollary 3.4. Let $A = (a_{ij})_{d\times d}$ be diagonally dominant with negative diagonal elements. Let $D$ denote the diagonal matrix $\mathrm{diag}(a_{11}, a_{22}, \cdots, a_{dd})$, and let $U, V$ be any $d\times d$ matrix satisfying $U+V = A-D$. Suppose that

Then, the method (3.1) with the splitting $A_1 = D+U$ is contractive and asymptotically stable if

Proof. By (3.4), we have

This and (3.5) yield

Multiplying $h$ and adding 1 to the two sides of (3.6) we get

Hence, we have

that is

Noting that (3.5) implies that $1-ah > 0$, we get the following from the above inequality

This means that

Note that $A_2 = V$ if $A_1 = D+U$ for the method (3.1). Thus, the above inequality yields

This together with Theorem 3.1, completes the proof of Corollary 3.4. □

Corollary 3.5. Let $A, D, U$ and $V$ be the same as those in Corollary 3.4 and $D+U$ be diagonally dominant. Suppose that

Then, the method (3.2) with the splitting $A_1 = D+U$ is contractive and asymptotically stable for any step size $h$.

Proof. Noting that $|\Re(a_{ii})| < |a_{ii}|$ we can obtain from (3.7) that

The above condition and (3.7) show that for any step size $h$

This means that

that is

Hence, we get

This together with Lemma 3.3, shows that

Noting that if $A_1 = D+U$ then $A_2 = V$ for the WR method (3.2); the above condition and Theorem 3.2 complete the proof of Corollary 3.5. □

A severe restriction on the matrix is required to satisfy the condition of diagonal dominance in the presence of negative diagonal elements. The following two corollaries show that the condition that all eigenvalues of $A$ in (2.1) have negative real parts is enough to guarantee the stability of the WR methods (3.1) and (3.2) whenever the splitting is suitable.

Corollary 3.6. Let A be the matrix in (2.1) and T be a non-singular matrix satisfying

where $\lambda_i$ is the eigenvalue of $A$ and $\mu_i = 0$ or $1$. Assume that

Then, the WR method (3.1) is asymptotically stable when the step size $h$ satisfies

Proof. Let $y^{(k)}_n = Tx^{(k)}_n, \forall k, n$. By virtue of (3.1) we have

where $\Lambda_1 = \left(\begin{array}{cccc} \lambda_1 & & & 0 \\ & \lambda_2 & & \\ & & \ddots & \\ 0 & & & \lambda_d \\ \end{array} \right), \Lambda_2 = \left(\begin{array}{cccc} 0 & \mu_1 & & 0 \\ & 0 & \ddots & \\ & & \ddots & \mu_{d-1} \\ 0 & & & 0 \\ \end{array} \right).$

Hence, (3.9) is asymptotically stable by Theorem 3.1 if

where $\|\cdot\|$ is an operator norm of the matrix.

It is clear that (3.1) is asymptotically stable if (3.9) is asymptotically stable. Thus, it is enough to show that (3.8) implies (3.10).

Taking $S = \mathrm{diag}(1, t, \cdots, t^d)$ and $\|A\|_S = \|SAS^{-1}\|_{\infty}$, we have

By (3.8) we have

Noting that

we can obtain

Hence,

By (3.11) and (3.12) we have the following for a sufficiently large $t$

This together with Theorem 3.1, completes the proof of Corollary 3.6.

□

Corollary 3.7. If $A_1$ and $A_2$ are the same as those in Corollary 3.6, then the WR method (3.2) is asymptotically stable for any step size $h$.

Proof. Let $y^{(k)}_n = Tx^{(k)}_n, \forall k, n$. By virtue of (3.2) we have

Here, $\Lambda_1$ and $\Lambda_2$ are the same as those in the proof of Corollary 3.6.

Take $S = \mathrm{diag}(1, t, \cdots, t^d)$ and $\|A\|_S = \|SAS^{-1}\|_{\infty}$. By Theorem 3.2, the method (3.13) is asymptotically stable if

Clearly, it is enough to prove that

that is

Note that $\Re(\lambda_i) < 0, \forall i$. Hence,

Let $t$ be large enough such that

Thus, for $t$ large enough, we have

that is (3.14) holds.

Consequently, the proof is complete by the fact that the asymptotically stable property of (3.2) and (3.13) is the same. □

Lemma 2.3 shows that if $\|I+hA_1\|+\|hA_2\| > 1$ and $e^{(k+1)}_{n+1} = \|I+hA_1 \|e^{(k+1)}_{n}+\|hA_2\|e^{(k)}_{n}$, then the WR method (3.1) is unstable.

Example 3.8. For the differential equation

and its WR method

let $\tilde{x}^{(k)}_n$ denote the perturbed value of $x^{(k)}_n, n, k\in \mathbb{Z}^+$, which satisfies

Let $\epsilon^{(k)}_n$ denote $\tilde{x}^{(k)}_n-x^{(k)}_n$. When taking $\epsilon^{(k)}_0 > 0$ for $k\in \mathbb{Z}^+$, $\epsilon^{(0)}_{2l} > 0$ and $\epsilon^{(0)}_{2l+1} < 0$ for $l\in \mathbb{Z}^+$, $\lambda_1, \lambda_2 < 0$ and $1+\lambda_1h < 0$, we have

Suppose that $\inf\limits_{j\geqslant 0}\left\{\min\left\{\left|\epsilon^{(0)}_j\right|, \left|\epsilon^{(j)}_0\right|\right\}\right\} > 0$. Then, the WR method (3.16) is unstable when $h > \max\left\{-\dfrac{1}{\lambda_1}, -\dfrac{2}{\lambda}\right\}$.

Similarly, the WR method (3.2) may also be unstable if $\|(I-hA_1)^{-1}\|+\|(I-hA_1)^{-1}hA_2\| > 1$ by Lemma 2.3.

Example 3.9. Applying (3.2) to Eq (3.15) we get the following WR method

Let $\lambda_2 < \lambda_1 < 0$. Let $\tilde{x}^{(k)}_n$ denote the perturbed value of $x^{(k)}_n$ and $\epsilon^{(k)}_n$ denote $\tilde{x}^{(k)}_n-x^{(k)}_n$ for $n, k\in \mathbb{Z}^+$. When taking $\epsilon^{(0)}_n < 0$ for $n\in \mathbb{Z}^+$ and $\epsilon^{(2l)}_{0} < 0$ and $\epsilon^{(2l+1)}_{0} > 0$ for $l\in \mathbb{Z}^+$, we have

Suppose that $\inf\limits_{j\geqslant 0}\left\{\min\left\{\left|\epsilon^{(0)}_j\right|, \left|\epsilon^{(j)}_0\right|\right\}\right\} > 0$. Then, the WR method (3.17) is unstable for any $h$.

4.   Stability of WR methods for Eq (2.2)

Let $\langle\cdot, \cdot \rangle$ denote an inner product and $\|\cdot\|$ the corresponding inner product norm. Assume that the splitting function $F(t, x, y)$ satisfies the following conditions:

(A1) there exists $C > 0$ such that

and

(A2) there exists $L_1, L_2 > 0$ such that

Consider the WR methods (1.3) and (1.4) of Eq (2.2). Let $\{\varepsilon^{(0)}_n, \varepsilon^{(k)}_0\}$ be the perturbations of the initial values of (1.3) and (1.4). Then, the perturbed solution $\{\tilde{x}^{(k)}_n\}$ corresponding to (1.3) satisfies the following:

and the perturbed solution $\{\tilde{x}^{(k)}_n\}$ corresponding to (1.4) satisfies the following:

Theorem 4.1. Suppose that

and the assumptions (A1) and (A2) hold. Then, the WR method (1.3) is contractive and asymptotically stable if

Proof. Let $\{x^{(k)}_n\}$ and $\{\tilde{x}^{(k)}_n\}$ denote the approximate values respectively generated by (1.3) and (4.1), and let $\delta^{(k)}_n, F_n$ and $\tilde{F}_n$ denote $x^{(k)}_n-\tilde{x}^{(k)}_n, F(t_n, x^{(k+1)}_n, x^{(k)}_n)$ and $F(t_n, \tilde{x}^{(k+1)}_n, \tilde{x}^{(k)}_n)$, respectively. By (1.3) and (4.1) we have

Using the assumptions (A1) and (A2), we can obtain

and

By (4.5), (4.6) and (4.7), we get

It is not difficult to prove that

under the conditions (4.3) and (4.4). This together with Lemma 2.2, completes the proof of Theorem 4.1. □

Theorem 4.2. Suppose that the condition (4.3) and the assumptions (A1) and (A2) hold. Then, the WR method (1.4) is contractive and asymptotically stable for any step size $h$.

Proof. Let $\{x^{(k)}_n\}$ and $\{\tilde{x}^{(k)}_n\}$ denote the approximate values respectively generated by (1.4) and (4.2), and let $\delta^{(k)}_n, F_{n+1}$ and $\tilde{F}_{n+1}$ denote $x^{(k)}_n-\tilde{x}^{(k)}_n, F(t_{n+1}, x^{(k+1)}_{n+1}, x^{(k)}_{n+1})$ and $F(t_{n+1}, \tilde{x}^{(k+1)}_{n+1}, \tilde{x}^{(k)}_{n+1})$, respectively. Using the property of the inner product we can obtain

Note that (1.4) and (4.2) imply that

Thus,

By (A1), (A2) and (4.9), we get

By using (4.8) and (4.10), we can derive that

Note that the condition (4.3) implies the following for any $h$

This together with (4.11) and Lemma 2.2, completes the proof of Theorem 4.2. □

5.   Numerical experiments

In this section we will present some numerical experiments to verify the theories developed in Section 3. Let $\sideset{_i}{_{n}^{(k)}}{\mathop{\xi}}, i = 1, 2, \cdots, d, k\in \mathbb{Z}^+, n\in \mathbb{Z}^+(k\cdot n = 0)$ be independent random variables, each uniformly distributed on the interval $[-0.5, 0.5]$. Take perturbations of the initial values $x^{(k)}_n$ as $\varepsilon^{(k)}_n = (\sideset{_1}{_{n}^{(k)}}{\mathop{\xi}}, \sideset{_2}{_{n}^{(k)}}{\mathop{\xi}}, \cdots, \sideset{_d}{_{n}^{(k)}}{\mathop{\xi}})^{\mathrm{T}}$, where $k, n\in \mathbb{Z}^+$ and $k\cdot n = 0$. Let $\{x^{(k)}_n\}$ and $\{\tilde{x}^{(k)}_n\}$ denote the numerical solutions respectively generated by (3.1) and (3.3), or (3.2) and its perturbed system given by

and let $e^{(k)}_n$ denote $\|x^{(k)}_n-\tilde{x}^{(k)}_n\|$.

In our experiments we have investigated the relation between the errors $e^{(k)}_n$ and factors such as the splitting, the step size and the iteration number. The results are presented as figures, where we have plotted $\log(e^{(k)}_n)$ versus the time $nh$. By Definition 2.1, the WR methods (3.1) and (3.2) are stable if $\lim\limits_{n, k\rightarrow \infty}\log (e^{(k)}_n) = -\infty$ and unstable if $\lim\limits_{n, k\rightarrow \infty}\log (e^{(k)}_n) = \infty$. Here and hereafter, the use of $\log$ denotes the logarithm with base ${\it e}$.

Let $M_1, M_2$ and $M_3$ denote the following matrices:

and

respectively. We have verified Corollaries 3.4 and 3.5 by investigating the stability of (3.1) and (3.2) with the BGJ splitting $A_1 = M_2, A_2 = M_1-M_2$ and the BGS splitting $A_1 = M_3, A_2 = M_1-M_3$. By Corollaries 3.4 and 3.5, the WR method (3.1) is stable for $h < 0.2,$ and (3.2) is stable for any $h$ when using the above BGJ and BGS splittings, which is in agreement with the results of the numerical experiments plotted in Figures 1 and 2.

To verify Corollaries 3.6 and 3.7, we take the matrix $A$ in (2.1) to be

whose eigenvalues have negative parts, but which do not satisfy the conditions of Corollaries 3.4 and 3.5. We have shown that the WR methods (3.1) and (3.2) in this situation are unstable when using the classical GJ and GS splittings, but stable when using the eigenvalue (EV) splitting given in Corollary 3.6 in this paper.

Let

Noting that

and

we can write the splittings of $A = M_4$ as the EV splitting $A_1 = M_6, A_2 = M_4-M_6$, the GJ splitting $A_1 = M_8, A_2 = M_4-M_8$ and the GS splitting $A_1 = M_{10}, A_2 = M_4-M_{10}$, and the splittings of $A = M_5$ as the EV splitting $A_1 = M_7, A_2 = M_5-M_7$, the GJ splitting $A_1 = M_9, A_2 = M_5-M_9$ and the GS splitting $A_1 = M_{11}, A_2 = M_5-M_{11}$.

By Corollary 3.6 the WR method (3.1) with $A = M_4$ is stable for the step size $h < 2/3$ when using the EV splitting, and that with $A = M_5$ is stable for the step size $h < 2$ when using the EV splitting. These are supported by the results of experiments plotted in Figure 3. By Corollary 3.7 the WR method (3.2) with $A = M_4 \text{ or } M_5$ is stable for any $h$ when using the EV splitting, which is consistent with the results of experiments plotted in Figure 4.

One may be interested to know about the stability of the WR methods with $A = M_4$ or $M_5$ when another splitting is used. Here, we provide insight into the problem. Taking the GJ and GS splittings of $A = M_4$ and $M_5$, we have obtained two groups of the numerical solutions of WR methods (3.1) and (3.2) displayed in Figures 5 and 6, which show that the WR (3.1) and (3.2) with $A = M_4 \text{ or } M_5$ are unstable for the GJ and GS splittings.

One may also want to know whether or not the instability of WR methods originates purely from the splitting function or from the underlying explicit Euler method. In Sections 3 and 4 only some sufficient conditions of stability are obtained for the WR methods, which are not enough to deduce the source of instability in WR methods. However, some reasonable conclusions can be reached by performing numerical experiments.

For the explicit Euler method of (2.1) there exists the maximum step size $h_{\max}$ such that the method is stable if and only if the step size $h < h_{\max}$. It seems that there also exists a similar $h_{\max}$ for the WR method (3.1) of (2.1) by Theorem 3.1 and Corollaries 3.4 and 3.6. We have explored whether or not $h_{\max}$ of (3.1) is the same as that of the underlying explicit Euler method by conducting numerical experiments. For comparison in a uniform standard, we take $\varepsilon^{(k)}_0 = \varepsilon^{(0)}_n = (1, 1, \cdots, 1)^T$ for all $k, n\in \mathbb{Z}^+$. The approximate values of $h_{\max}$ for the explicit Euler method and the WR method (3.1) as applied to (2.1) have been obtained, which are listed in Tables 1–3, for the following three cases:

Case 1. $A = -3, A_1 = \alpha, A_2 = A-A_1;$

Case 2. $A = M_1, A_1 = \left(\begin{array}{ccccc} -5 & & & & \\ \alpha & -5 & & & \\ & & -6 & & \\ & & & -5 & \\ & & & & -6 \\ \end{array} \right), A_2 = A-A_1;$

Case 3. $A = M_5, A_1 = \left(\begin{array}{ccc} -1 & & \\ & \alpha & \\ & & -1 \\ \end{array} \right), A_2 = A-A_1.$

From Table 1, we see that the $h_{\max}$ of the WR method (3.1) for Case 1 and of the underlying explicit Euler method are the same; particularly, the $h_{\max}$ of (3.1) does not change with the splitting, i.e., according to the value of $\alpha$. Hence, the instability of (3.1) for Case 1 originates purely from underlying explicit Euler methods and is independent of the splitting. However, from Tables 2 and 3, as well as the comparison of Figures 3 and 5, we can find convincing evidence of a link between the instability of (3.1) and the splitting matrix used. Consequently, the stability of (3.1) for Cases 2 or 3 may change with the splitting.

Hence, we reasonably conclude that the instability of WR method (3.1) originates purely from neither the splitting function used nor the underlying explicit Euler method.

6.   Conclusions

In this paper, we have discussed the stability properties of DWR methods based on Euler schemes. We pointed out that the DWR methods may be unstable by virtue of a numerical example. There is hence a need to study the stability of DWR methods. In this paper, we have applied DWR methods based on Euler schemes to two dissipative systems and analyzed the properties of the numerical solutions generated. We then obtained some conditions under which the DWR methods can preserve the stability of the exact solutions. These conditions show that underlying methods and splitting ways are two key factors in determining the stability of DWR methods.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This work was supported by the Natural Science Foundation of Fujian Province, China (grant number 2021J011031) and Research Project of Fashu Foundation (MFK23013).

Conflict of interest

The authors declare no conflict of interest.