3.3 New operator splittings and their analysis
3.3.2 Additive splitting
The split solutions that we have considered before are weakly consistent at the mesh points of the fixed mesh ωτ. When we solve the split sub-problems numerically (which is almost always necessary), for the step-size of the numerical method we can use smaller values than the splitting step-size. Our aim is to use this intermediate values as the approximation to the unknown function. (We discuss this question in Section 3.5 in more details.) Therefore, to investigate the consistency, we have to compare the split solution with the exact solution at any point of the split time-interval, i.e., on the whole [0, τ].
First we should extend the split (discrete) solution to this interval. It is rather plausible to do it, e.g., for the sequential splitting, as follows:
w(Nseq)(t) = exp(tA2) exp(τ A1)w0. (3.3.19) Hence, for the local splitting error at any timet ∈[0, τ] we get
w(t)−w(Nseq)(t) = (t−τ)A1w0+O(t2). (3.3.20) This shows that Errseq(τ) = O(τ), i.e., the method is not weakly consistent for arbitrary t→0 and a fixed τ, only in the case where t=τ →0.
Remark 3.3.7 For the Strang-Marchuk splitting and symmetrically weighted sequential splitting it is quite natural to define the extension of the split solution to the whole interval [0, τ] as
w(NSM)(t) = exp(0.5tA1) exp(τ A2) exp(0.5τ A1)w0. (3.3.21) and
wswss(N)(t) = 0.5 (exp(tA1) exp(τ A2) + exp(tA2) exp(τ A1))w0. (3.3.22) An easy computation shows that they are weakly consistent only at t =τ.
Hence, it is reasonable to introduce the following definition.
Definition 3.3.8 A splitting method with a fixed step-size τ is called continuously weakly consistent when
Errspl(t) =w(t)−w(Nspl)(t) =O(tp+1) (3.3.23) for all t∈(0, τ] with some p >0.
The algorithm of the additive splitting is based on the following idea. Simultaneously (in a parallel way) we solve the Cauchy problems consisting only of one operator Ai and using the same initial value for each sub-problem, namely, the split solution at the previous time level. Then, by special averaging of the results (in order to ensure the weak consistency), we define the split solution at the next time level.
Definition 3.3.9 The operator splitting corresponding to the choice
If we apply the above approximations on the meshωτ, then the iteration (3.1.4) reads as w(Nas )((n+ 1)τ) = ras(τ A)w(Nas )(nτ), n = 0,1, . . . , N, (3.3.25) In order to define the required exponentials in the algorithm, we solve the corresponding Cauchy problems, namely the realization is the following.
1. For each n= 1,2, . . . , N we solve the Cauchy problems:
2. The split solution is defined as w(N)as (nτ) =
Xd
i=1
wni(nτ)−(d−1)w(N)as ((n−1)τ). (3.3.28) Herew(Nas )(0) =w(0) is known from the original continuous problem (3.1.1).
That is, the algorithm is as follows:
A1
In the following we compute the order of the additive splitting. For the split solution we get
Using the formula (3.2.13) for the solution of the un-split problem, for the local splitting error of the additive splitting we obtain
Erras(τ) = τ2 2
Xd
i,j=1 i6=j
AiAjw0+O(τ3). (3.3.29)
Let us notice that (3.3.29) also holds for an arbitrary t∈(0, τ], i.e., Erras(t) = t2
2 Xd
i,j=1 i6=j
AiAjw0+O(t3). (3.3.30) Hence we have
Theorem 3.3.10 For arbitrary operators the additive splitting is continuously consistent, and it has first order accuracy.
As opposed to the sequential splitting, the Strang-Marchuk splitting and the weighted se-quential splitting, the order of the additive splitting is not influenced by the commutativity of the operators. However, if any pair Ai and Aj anticommute, i.e., if AiAj +AjAi = 0, then the local splitting error has second order. We note that in contrast to the previous splittings, where the pairwise commutativity resulted in the exactness of the splitting, for the additive splitting the anticommutativity of the operators does not imply higher than second-order accuracy.
To check the theoretically obtained order, we applied the additive splitting to two systems of simple ordinary differential equations with constant coefficients, where both the original problems and the sub-problems were solved exactly. In the first case the sub-matrices did not anticommute. Figure 3.3.1 shows the results on a log-log diagram.
The points are located along a line with slope close to 2, which means that the method has first order.
In the second case we chose anticommuting matrices (so-called Pauli matrices, well-known in quantummechanics [102]). Here, as Figure 3.3.2 shows, second-order accuracy was achieved.
The subtraction in the second step of (3.3.27) - (3.3.28) causes significant theoretical difficulties in order to show the stability. Therefore we modify the method. We execute the separate splitting steps with the sub-operators multiplied by d, and then we compute the splitting approximation as an average of these results. Namely, the algorithm reads as follows:
1. For each n= 1,2, . . . , N we successively solve the Cauchy problems:
dwni
dt (t) = dAiwni(t), (n−1)τ < t≤nτ, wni((n−1)τ) =w(Nmas)((n−1)τ)
(3.3.31) where i= 1,2, . . . , d.
2. The split solution is defined as
w(Nmas)(nτ) = 1 d
Xd
i=1
win(nτ). (3.3.32)
10−3 10−2 10−1
τ
Theoretical error Fitted line
Figure 3.3.1: The local error of the additive splitting on a log-log diagram, when the sub-problems are solved exactly. The slope of the fitted line is 2.025.
10−3 10−2 10−1
τ
Theoretical error Fitted line
Figure 3.3.2: The local error of the additive splitting on a log-log diagram, when the sub-matrices anticommute, and the sub-problems are solved exactly. The slope of the fitted line is 3.002.
Herew(Nmas)(0) =w(0) is known from the original continuous problem (3.1.1).
We will refer to this method as modified additive splitting (mas).
Remark 3.3.11 It is possible to replace the operators dAi with Ai in (3.3.31). If such a replacement is performed, then the integration should be carried out on the interval
(n−1)τ < t≤(n+d−1)τ, and (3.3.32) should be replaced withw(N)mas(nτ) = 1dPd
i=1win((n+
d−1)τ). If the operators are linear, then the two approaches sketched above are equivalent.
However, for non-linear operators the results obtained by using these two approaches will in general be different.
The investigation of the local splitting error (for bounded linear operators) for the modified additive splitting (3.3.31)-(3.3.32) is similar as it was done in Theorem 3.3.10.
Theorem 3.3.12 For bounded linear operators the modified additive splitting is a first order accurate, continuously consistent method.
Proof. The solution of the modified additive splitting att ∈(0, τ] can be given as w(N)mas(t) = 1
Comparing this expression with the similar Taylor expansion of the exact solution, we get the local splitting error
As we have seen, the sequential splitting is not continuously consistent. The reason for this is the fact that on each sub-problem we passed through the time sub-interval by using one operator only. This observation gives the inspiration to modify the sequential splitting in such a way that each operator is involved in the sub-problems, but only one of them is applied to the unknown function, the other operators are applied to some known (previously defined) approximation. Then, we iterate consecutively on the same time interval m ∈ IN times. The method, which will be called iterated splitting, is defined for two operators, (i.e., d= 2 in (3.1.1)) and was introduced in [39]. On some fixed split time interval [(n−1)τ, nτ] the algorithm of the iterated splitting reads as follows: