A variational approach to the sum splitting scheme

Nonlinear parabolic equations are frequently encountered in applications and efficient approximating techniques for their solution are of great importance. In order to provide an effective scheme for the temporal approximation of such equations, we present a sum splitting scheme that comes with a straight forward parallelization strategy. The convergence analysis is carried out in a variational framework that allows for a general setting and, in particular, nontrivial temporal coefficients. The aim of this work is to illustrate the significant advantages of a variational framework for operator splittings and use this to extend semigroup based theory for this type of scheme.


Introduction
Nonlinear parabolic equations, which we state as abstract evolution equation of the form u (t) + A(t)u(t) = f (t), t ∈ (0, T ) and u(0) = u 0 , (1.1) are frequently encountered in applications appearing in physics, chemistry and biology; see Aronsson et al. (1996) and (Vázquez, 2007 Here, the first and second operator is referred to as the p-Laplacian and the porous medium operator, respectively. Due to the problems' significance, effective techniques for their approximations become crucial. As we consider parabolic equations, for stability reasons the temporal approximation schemes need to be implicit. For equations which in addition are given in several spatial dimensions the resulting spatial and temporal approximation schemes require large scale computations. This typically demands implementations in parallel on a distributed hardware. One possibility to design temporal approximation schemes that can directly be implemented in a parallel fashion is to utilize operator splitting; see, e.g., Hundsdorfer & Verwer (2003) for an introduction. Note that the solutions of nonlinear parabolic problems typically lack high-order spatial and temporal regularity. Thus, there is little use to consider high-order time integrators.
In order to illustrate the splitting concept, consider the simplest implicit scheme, namely, the backward Euler method. For N temporal steps, a step size k = T /N and the starting value U 0 = u 0 the backward Euler approximation U n of u(nk) is given by the recursion 1 k (U n − U n−1 ) + kA n U n = f n , n ∈ {1, . . . , N}, where (A n ) n and (f n ) n are suitable approximations for A and f , respectively. Assuming that the nonlinear resolvent of A n exists, we find the reformulation U n = (I + kA n ) −1 U n−1 + kf n ), n ∈ {1, . . . , N}.
To implement one step of the backward Euler scheme in parallel, we split the Euler step into s independently solvable problems. To this end, we decompose A n and f n as Thus, we obtain s possibly easier subproblems that are solved after each other. For a straightforward parallelization it is more convenient to choose a splitting, where the single steps can be computed at the same time. The sum splitting (I + kA n ) −1 ≈ 1 s s ∑ =1 (I + ksA n ) −1 offers this crucial advantage. The s fractional steps are solved at the same time and their average is used as an approximation. The decomposition (1.3) can be utilized in many different ways. A first possible application is a source term splitting, where the high-order terms are split from the low-order terms. For example a source term splitting of the reaction-diffusion equation governed by A(t)v = −∇ · α(t)∇v + p(t, v) would have the form A n 1 v = −∇ · α(nk)∇v and A n 2 v = p(nk, v).
Here, the actions of (I + kA n 1 ) −1 can be evaluated by a standard fast linear elliptic solver and the actions of the nonlinear resolvent (I + kA n 2 ) −1 can often be expressed in a closed form. Examples of studies dealing with various source term splittings can be found in Arrarás et al. (2017); Hansen & Stillfjord (2013); Koch et al. (2013);Eisenmann (2019).
Another possibility is a dimension splitting, where each spatial derivative is considered as a separate differential operator. For example, the dimension splitting of the nonlinear porous medium operator and the third operator in (1.2) are given by respectively. This splitting yields that the action of each nonlinear resolvent (I + kA n ) −1 can be separated into lower-dimensional subproblems that can be solved on their own. Note that the p-Laplacian lacks a natural dimension splitting. Examples of convergence results for the dimension splitting of the third equation in (1.2) can be found in Temam (1968), where the Lie scheme is used, and in Hansen & Ostermann (2008), where the sum, Lie and Peaceman-Rachford schemes are considered for the autonomous case.
A limitation of the dimension splitting approach is the rather large need of communication between the subproblems, which can impede an effective distributed implementation. Dimension splitting is also quite restrictive in terms of the spatial domains that can be considered. A modern alternative to dimension splitting, which is applicable to a very general class of spatial domains, is the domain decomposition based splitting. Here, the subproblems are given on s spatial subdomains that share a small overlap. As an example consider the three nonlinear diffusion operators (1.2) and introduce a partition of unity (χ ) s =1 , where each weight function χ vanishes outside its corresponding spatial subdomain. The domain decompositions A n v are then respectively. This approach is well suited for a parallel computation, as the actions of (I +kA n ) −1 can be solved independently of each other and the communication required is small, due to the small overlaps between the subdomains. Studies regarding domain decomposition based splittings applied to linear and autonomous parabolic equations include Arrarás et al. (2017); Hansen & Henningsson (2016); Mathew et al. (1998);Vabishchevich (2008). Convergence for the Lie and sum splittings are given in Eisenmann & Hansen (2018) for the autonomous p-Laplace and porous medium equations.
Operator splitting schemes are typically analyzed in a semigroup framework, which yields convergence for a wide range of temporal approximation schemes, including the Lie and sum schemes; see Barbu (1976) for more details on the solution concept. However, there does not seem to be a straightforward way to extend the semigroup based convergence analysis to nonautonomous evolution equations. Furthermore, the semigroup framework requires some additional regularity conditions to relate the intersection of the domains D(A ), ∈ {1, . . . , s}, with the domain D(A) of the full operator. The latter, e.g., implies restrictions on the domain decomposition of the p-Laplace equation Eisenmann & Hansen (2018, Section 6).
In a variational setting this problem is avoided in a natural way while at the same time the analysis of nonautonomous problems is accessible. Also the structure of this approach is well suited to include a Galerkin scheme and therefore, in particular, the finite element method. However, the analysis typically needs to be tailored for each method. The variational setting is a standard tool for existence theories Emmrich (2004); Roubíček (2013); Zeidler (1990) and has been used in several works in the context of "unsplit" time integrators Emmrich (2009,b,c); Emmrich & Thalhammer (2010). However, in the context of temporal splitting schemes for nonlinear parabolic equations the only variational studies that we are aware of is Temam (1968). Here, a variational analysis is employed when proving the convergence of the Lie scheme applied to nonautonomous evolution equations and, as already stated, is applied to the dimension splitting of the third equation in (1.2).
Hence, the aim of this paper is threefold. Firstly, we aim to generalize the previous semigroup based analysis for the sum scheme to nonautonomous evolution equations without any implicit regularity assumptions. The latter generalization will be applicable to splittings of reaction-diffusion, dimension and domain decomposition type. Secondly, we intend to extend the abstract variational convergence results for the Lie scheme to the sum splitting scheme. As this requires a tailored convergence proof, it is not a trivial implication. Thirdly, we also strive to illustrate the advantages of a variational approach in the context of splitting analyses. This paper is organized as follows: In Section 2, we state the exact assumptions that are needed on the abstract variational framework considered in the paper. This section also contains an example that shows that the relevant application of domain decomposition integrators for the p-Laplacian operator fits into our abstract framework. This in mind, we prove the well-posedness of the sum scheme, as well as suitable a priori bounds in Section 3. The main convergence results are proven in Section 4; see Theorem 4.1 and Theorem 4.2.

Abstract setting
In this section, we introduce an abstract setting for the convergence analysis of the sum splitting scheme. We begin by presenting the exact assumptions made on the data and present the temporal discretization of the problem. This at hand, we can state the scheme that we will work with in this paper. The section ends with a more concrete setting that exemplifies the abstract framework.
Assumption 1. Let (H, (·, ·) H , · H ) be a real, separable Hilbert space and let (V, · V ) be a real, separable, reflexive Banach space such that V is continuously and densely embedded into H. Further, there exist a seminorm | · | V on V and c V ∈ (0, ∞) such that · V c V · H + | · | V is fulfilled.
Identifying H with its dual space H * , we obtain the Gelfand triples The next assumption states the properties of the differential operator that are of importance.
Assumption 2. Let H and V be given as stated in Assumption 1. Furthermore, for T > 0 and p > 1, let {A(t)} t∈[0,T ] be a family of operators such that A(t) : V → V * satisfy the following conditions: (1) The mapping Av : , fulfills a monotonicity condition such that there exists η 0 with , fulfills a coercivity condition such that there exist µ > 0 and λ 0 with Now, we can combine Assumption 1 and Assumption 2 to state a decomposition of the operator family {A(t)} t∈[0,T ] that we employ in the analysis of the sum splitting scheme.
. . , s}, of operators do not necessarily have to be the same. For the sake of simplicity, we assume that these coefficients coincide.
We also consider the differential operators of Assumption 3 as Nemytskii operators acting on spaces of Bochner integrable functions. For an introduction to Bochner integrable functions we refer the reader to (Diestel & Uhl, 1977, Chapter II) or (Papageorgiou & Winkert, 2019, Section 4.2). Some properties of such Nemytskii operators are collected in the next lemma. The proofs can be found in (Emmrich, 2004, Lemma 8.4.4).
. This operator is radially continuous, i.e., the mapping τ → A(u + τv), w L q (0,T ;V * )×L p (0,T ;V ) is continuous on [0, 1] for all u, v, w ∈ L p (0, T ;V ). Furthermore, it fulfills a monotonicity, a boundedness and a coercivity condition such that . . , s}, as introduced in Assumption 3 also fulfills the same bounds with V replaced by V . To make our setting complete, it remains to state the assumptions on f .
of Ω ⊂ R 2 , with s = 4 subdomains (left) and s = 2 subdomains that are further decomposed into disjoint sets (right).
We can now state the abstract evolution equation that we want to consider. In the following, let T ] be as stated in Assumption 3, let f fulfill Assumption 4 and let u 0 ∈ H be given. It is our overall goal to find an approximation to the solution u of see (Lions & Strauss, 1965, Section 2.7) and (Roubíček, 2013, Chapter 7-8) for further details. In the following analysis, we employ the sum splitting in order to obtain a temporal discretization of (2.1). To this end, we consider an equidistant grid on [0, T ], where N ∈ N, k = T N and t n = nk for n ∈ {0, . . . , N}. For ∈ {1, . . . , s} and n ∈ {1, . . . , N} we introduce We use this to construct an approximation U n ≈ u(t n ) of the solution u of (2.1) for n ∈ {0, . . . , N}. This approximation is given through a recursion for n ∈ {1, . . . , N} with U 0 = u 0 .
EXAMPLE 2.1 A useful example that fits into our abstract setting is to approximate the solution of the where Ω is a bounded domain and the boundary ∂ Ω is Lipschitz. For p 2 we consider the problem where n denotes outer pointing normal vector. The function α : On these subdomains let the partition of unity {χ } s =1 ⊂ W 1,∞ (Ω ) be given such that we use L 2 (Ω ) the space of square integrable functions on Ω with the usual norm and inner product. The energetic spaces V and V are then given as which are equipped with the norms Together with the partition of unity we define the decomposed energetic operators It is also possible to allow for more general coefficients α : ·, ·) fulfills the condition stated in (Eisenmann & Hansen, 2018, Assumption 3). We assume that for g : . . , d}. These functions are not necessarily unique unless we exchange V = W 1,p (Ω ) by V = W 1,p 0 (Ω ), compare (Leoni, 2009, Theorem 10.41, Corollary 10.49). This in mind, we introduce f (t) for a.e. t ∈ (0, T ) as Note that in this type of setting, we can also consider homogeniuous Direchlet boundary conditions in (2.4). Then an additional condition on the partition of unity becomes necessary. In this case, we have to make the further assumption that for every function χ there exists ε 0 > 0 such that for all ε ∈ (0, ε 0 ) Further examples that fit our framework are a domain decomposition scheme for the porous medium equation as presented in (Eisenmann & Hansen, 2018, Section 7) or a source term splitting as in (Eisenmann, 2019, Section 3.3). An application to the third equation of (1.2) is presented in Temam (1968). Numerical experiments for this equation and the p-Laplace equation can be found in Eisenmann & Hansen (2018) and Hansen & Ostermann (2008), respectively.

Solvability and a priori bounds for the discrete scheme
The abstract setting from the previous section in mind we are now well-prepared to state some properties of the solution of the numerical scheme (2.3). Since the scheme is implicit, we start to verify that (2.3) is uniquely solvable. Once this is at hand, we can provide a priori bounds of the solution. These bounds are a crucial part of the further analysis and allow for the convergence analysis in Section 4 LEMMA 3.1 Let Assumptions 3 and 4 be fulfilled. Then the semidiscrete problem (2.3) is uniquely solvable.
We can now turn our attention to the a priori bounds.
Proof. In the following, let i ∈ {1, . . . , N} and ∈ {1, . . . , s} be arbitrary but fixed. Recall the identity and the inequality · V c 1 · H + | · | V with c 1 = max ∈{1,...,s} c V stated in Assumption 1. Using the weighted Young inequality, see (Evans, 1998, Appendix B.2.d)), we obtain Thus, together with the coercivity condition from Assumption 2 (5) it follows that Employing the specific structure of U i−1 , we obtain for i ∈ {2, . . . , N} due to the Cauchy-Schwarz inequality for sums. Inserting this inequality in (3.5), summing up from = 1 to s as well as dividing by s, yields After a summation from i = 1 to n ∈ {1, . . . , N} and using the telescopic structure, we obtain For the right-hand side we can bound the summands using Assumption 4 and Hölder's inequality Thus, we get that (3.8) As this is fulfilled for every n ∈ {1, . . . , N}, it also follows that max n∈{1,...,N} We abbreviate the terms to obtain x 2 2ax + b 2 . This implies, in particular, that Taking the square root on both sides, this yields |x − a| a 2 + b 2 a + b.
As x − a |x − a| is fulfilled, we obtain x 2a + b after adding a to both sides of the inequality. This shows that max n∈{1,...,N} where M 1 0 is independent of k. Using the norm inequality from Assumption 1, this implies that there exists M 2 0, which does not depend on k, such that Altogether, we have proved the first a priori bound (3.2). In order to prove (3.3), we test (2.3) with v ∈ V and use Assumption 2 (4) to see that where c 3 is the maximal embedding constant of V into V for ∈ {1, . . . , s}. Thus, we can estimate the V * -norm by This bound can be used to see that there exists M 0 such that Due to the first a priori bound (3.2) and (3.7) the constant M is independent of k.

Convergence analysis
In the following, we introduce prolongations of the solution of the discrete problem (2.3) to the interval [0, T ]. The main goal of this section is to prove that the sequence of such prolongations converges to the exact solution u of (1.1). Corresponding to the grid 0 = t 0 < t 1 < · · · < t N = T with k = T N and t n = nk, n ∈ {0, . . . , N}, we construct piecewise constant and piecewise linear functions on the interval [0, T ]. We consider the piecewise constant functions for t ∈ (t n−1 ,t n ], n ∈ {1, . . . , N}, and ∈ {1, . . . , s} given by as well as the piecewise linear functioñ with U k (0) = U k (0) =Ũ k (0) = u 0 , A k (0) = A 1 and f k (0) = f 1 . As we consider step sizes k = T N for N ∈ N, we denote the sequences U T N N∈N as (U k ) k>0 for ∈ {1, . . . , s} in the following to keep the notation more compact. The same simplification in notation is used for the other functions introduced above. Due to the a priori bound (3.2) we see that U k ∈ L p (0, T ;V ) ∩ L ∞ (0, T ; H), U k ,Ũ k ∈ L ∞ (0, T ; H), and f k ∈ L q (0, T ;V * ).
Furthermore, due to Lemma 2.1 the operator A k maps the space L p (0, T ;V ) into L q (0, T ;V * ). Using the prolongations introduced above, we can state a discrete version of the differential equation. We first note that after summing up (2.3) from 1 to s and dividing by s, we obtain 1 ks Thus, we see that where (Ũ k ) is the weak derivative ofŨ k . In the following, we will consider the limiting process of all the appearing terms to connect to the original problem (2.1) with (4.3).
LEMMA 4.1 Let Assumption 3 be fulfilled and let W ∈ L p (0, T ;V ) be given. For ∈ {1, . . . , s} it follows that A k (t)W (t) → A (t)W (t) in V * as k → 0 for a.e. t ∈ (0, T ). Furthermore, it holds true that A k W → A W in L q (0, T ;V * ) as k → 0.
Proof. Let ∈ {1, . . . , s} and ε > 0 be arbitrary. Due to the continuity condition on A , for almost every t ∈ (0, T ) we find δ > 0 such that for all k < δ it follows that where t is within an interval (t n−1 ,t n ], n ∈ {1, . . . , N}. The second assertion of the lemma is a consequence of Lebesgue's theorem of dominated convergence and the boundedness condition (4) from Assumption 2.
Proof. The statement above can easily be verified for a function from the space C([0, T ];V * ). As the space C([0, T ];V * ) is a dense subspace of L q (0, T ;V * ) for ∈ {1, . . . , s} a density argument can be used to verify the claimed statement.
LEMMA 4.3 Let Assumptions 3 and 4 be fulfilled. Then there exists a subsequence (k i ) i∈N of step sizes as well asŨ for every ∈ {1, . . . , s} as i → ∞. Here, U denotes the weak derivative of U.

Proof.
In the following proof, we do not distinguish between subsequences by notation. Using Lemma 3.2, we obtain that as k → 0. In the following, we prove that U 1 = · · · = U s =: U is fulfilled. As s =1 V = V and the norm ∑ s =1 · V is equivalent to · V this implies that U ∈ L p (0, T ;V ). We show that U 1 = U 2 in L p (0, T ;V 1 ∩V 2 ) ∩ L ∞ (0, T ; H), the other equalities follow in an analogous manner. We can write for t ∈ (t n−1 ,t n ], n ∈ {1, . . . , N}, as U k 1 (t n−1 ) = U n−1 = U k 2 (t n−1 ) holds true by the construction of our scheme. Therefore, we obtain for ∈ {1, 2} which is bounded independently of n and k due to the a priori bound (3.2). Thus, it follows that U k 1 (t) −U k 2 (t) V * → 0 as k → 0 for every t ∈ [0, T ] and it is bounded by a constant independent of t. Using Lebesgue's theorem of dominated convergence and the fact that L q (0, T ;V * ) → L 1 (0, T ;V * ), it follows that for c 1 ∈ (0, ∞). This shows that U 1 −U 2 = 0 in L 1 (0, T ;V * ). By assumption the embedding L p (0, where we used the a priori bound (3.3). Making use of the continuous embedding L ∞ (0, T ; H) → L 1 (0, T ;V * ), it follows that U andŨ coincide in L ∞ (0, T ; H). Last, we prove that the limit W ∈ L q (0, T ;V * ) is the weak derivative of U. To this end, let v ∈ V and Applying (Gajewski et al., 1974, Kapitel IV, Lemma 1.7), we obtain W = U in L q (0, T ;V * ). The next lemma is an auxiliary result to identify the limit from the previous lemma with the solution of (2.1).
Combining the prior lemmas, we can now state one of the main results of this section. We prove that the limit of the sequence of prolongations is the solution of (2.1).
Proof. In the following, we will not distinguish between different subsequences by notation. Due to Lemma 2.1 as well as the a priori bound (3.2) there exists a constantM > 0 such that for every k > 0, t ∈ [0, T ], and ∈ {1, . . . , s} M and U k (t) H M is fulfilled. Therefore, we can extract a subsequence of step sizes such that there exits B ∈ L q (0, T ;V * ) and y t ∈ H with as k → 0 for ∈ {1, . . . , s}. In the following, we abbreviate B := ∑ s =1 B . Next we identify the derivative of U with equation (2.1). Due to Lemma 4.2 and Lemma 4.3, the following equality holds true where w-lim denotes the weak limit. The limit U obtained in Lemma 4.3 is an element of W p (0, T ) → C([0, T ]; H). Thus, we can work with the continuous representative of U in the following. This in mind, we prove y t = U(t) and u 0 = U(0) for t ∈ [0, T ]. To this end, let v ∈ V and ϕ ∈ C([0, T ]) ∩C 1 (0, T ) be arbitrary. Recalling the equation for the time discrete values (4.3), we can write for t ∈ (t n−1 ,t n ], n ∈ {1, . . . , N}. Applying integration by parts and the fact that the linear and the constant interpolations always coincide at the grid points then shows that Recall that in Lemma 4.3 we have proved thatŨ k * U in L ∞ (0, T ; H) as k → 0. Therefore, (4.10) and the fact that ( f k − B ) k>0 is a bounded sequence in L q (0, T ;V * ) for every ∈ {1, . . . , s} shows that is fulfilled. This implies U(t) = y t and U(0) = u 0 for every t ∈ [0, T ].
It remains to prove that B = AU is fulfilled. To this end, we use Lemma 4.4. Applying Lemma 4.3, it follows that U k U in L p (0, T ;V ) as k → 0 for ∈ {1, . . . , s}. Further, in (4.10), we have seen that A k U k B in L q (0, T ;V * ) as k → 0 for every ∈ {1, . . . , s}. Therefore, we still have to verify in order to apply Lemma 4.4. To this end, we test the semidiscrete problem (2.3) with U n for ∈ {1, . . . , s} and n ∈ {1, . . . , N} to obtain that U n − U n−1 + skA n U n , U n V * ×V = skf n , U n V * ×V .
Summing up the equation form = 1 to s, dividing by s, and applying the identity from (3.4), it follows that 1 2s After another summation for n = 1 to N and an application of and (3.6), we can rewrite this inequality to 1 2s due to a telescopic sum structure. Inserting the definition of the prolongations from (4.1) and (4.2), we see that 1 2s Together with (3.6) and the weak lower semicontinuity of the norm this yields Applying Lemma 4.4, this verifies that B = AU is fulfilled in L q (0, T ;V * ). Thus, U is a variational solution to the original problem (2.1). As this problem has a unique solution u ∈ W p (0, T ), it follows that U = u. Next, we argue that the original sequence (U k ) k>0 converges weakly to the unique solution u of (2.1) in L p (0, T ;V ) for every ∈ {1, . . . , s}. The arguments above show that every converging subsequence of the bounded sequence (U k ) k>0 has the limit u. Applying the subsequence principle, see, e.g, (Zeidler, 1986, Proposition 10.13), yields that the entire sequence converges to this limit which proves (4.5). An analogous argumentation shows that (4.6)-(4.8) hold true for the original sequence. To prove (4.9), we recall (4.8) and the statement of Lemma 4.2. Inserting these two limiting process in (4.3) yields (4.9). THEOREM 4.2 Let Assumptions 3 and 4 be fulfilled. Then for step sizes k = T N the sequence (U k ) k>0 defined in (4.1) fulfills U k (t) → u(t) in H as k → 0, for t ∈ [0, T ], where u is the solution (2.1). If η in Assumption 2 (3) is strictly positive then the sequence (U k ) k>0 converges strongly to u in L p (0, T ;V ) for ∈ {1, . . . , s}.
Proof. For the analysis we split up the terms as follows 1 s s ∑ =1 u(t) −U k (t) 2 H + 2