# A posteriori error analysis for a space-time parallel discretization of
parabolic partial differential equations

Jehanzeb H. Chaudhry  Department of Mathematics and Statistics, University of
New Mexico. Email<EMAIL_ADDRESS>Don Estep  Department of Statistics and
Actuarial Science, Simon Fraser University. Email<EMAIL_ADDRESS>Simon Tavener  Department of Mathematics, Colorado State University. Email:
<EMAIL_ADDRESS>
###### Abstract

We construct a space-time parallel method for solving parabolic partial
differential equations by coupling the Parareal algorithm in time with
overlapping domain decomposition in space. The goal is to obtain a
discretization consisting of “local” problems that can be solved on parallel
computers efficiently. However, this introduces significant sources of error
that must be evaluated. Reformulating the original Parareal algorithm as a
variational method and implementing a finite element discretization in space
enables an adjoint-based _a posteriori_ error analysis to be performed.
Through an appropriate choice of adjoint problems and residuals the error
analysis distinguishes between errors arising due to the temporal and spatial
discretizations, as well as between the errors arising due to incomplete
Parareal iterations and incomplete iterations of the domain decomposition
solver. We first develop an error analysis for the Parareal method applied to
parabolic partial differential equations, and then refine this analysis to the
case where the associated spatial problems are solved using overlapping domain
decomposition. These constitute our Time Parallel Algorithm (TPA) and Space-
Time Parallel Algorithm (STPA) respectively. Numerical experiments demonstrate
the accuracy of the estimator for both algorithms and the iterations between
distinct components of the error.

Keywords — A posteriori error analysis, adjoint based error estimation,
parabolic partial differential equations, Parareal, domain decomposition,
Schwarz algorithms, time-parallel.

##  1 Introduction

Parallel computing approaches for solving complex problems modeled by partial
differential equations have become increasingly attractive over the last
decade, as Moore’s Law has continued, but transitioned from fabricating more
transistors on a chip to constructing more cores. Spatial parallelization of
partial differential equations is well established [3, 13, 32, 33, 41], but
the improvement available through spatial parallelization alone typically
decays as the number of processors increases for any fixed-size problem (the
so-called strong scaling behavior), see e.g., [20]. Parallel in time methods
provide an additional way to effectively utilize modern high performance
computer architectures [15, 19, 22, 28, 29, 35, 38]. Unfortunately parallel-
in-time methods also introduce significant additional sources of error which
must be estimated for the purposes of validation when seeking to construct
adaptive algorithms, or when performing uncertainty quantification.

We consider a PDE system of the form: Find $u(x,t)\in
L^{2}(0,T;H^{1}_{0}(\Omega))$ such that

$\displaystyle(\dot{u},v)+a(u,v)$ $\displaystyle=l(v),\quad\forall v\in
H^{1}_{0}(\Omega)\text{ and }t\in(0,T],$ (1) $\displaystyle u(x,0)$
$\displaystyle=u_{0}(x),$

where $\dot{u}=\frac{\partial u}{\partial t}$, $a(\cdot,\cdot)$ is a coercive,
positive definite bilinear form, $l(\cdot)$ is a linear form, $\Omega$ is a
convex polygonal domain, and $u_{0}\in H^{1}_{0}$ is the initial condition.
Here $(\cdot,\cdot)$ denotes the $L_{2}(\Omega)$ inner product, so that
$(w,v)=\int_{\Omega}wv\,{\rm d}x$. We assume that the forms $a$ and $l$
satisfy sufficient regularity assumptions so that a weak solution $u(x,t)\in
L^{2}(0,T;H^{1}_{0}(\Omega))$ exists and is unique [31]. Further, we consider
a quantity of interest (QoI) given by

$\mathcal{Q}(u)=\int_{\Omega}\psi(x)u(x,T)\,{\rm d}x$ (2)

where $\psi\in L^{2}(\Omega)$.

The aim of this work is to derive error estimates for the numerically computed
value of $\mathcal{Q}$ when time-parallel and spatial domain decomposition
techniques are employed for the approximation of $u(x,t)$. In particular, we
consider two algorithms which we call the Time Parallel Algorithm (TPA) and
the Space-Time Parallel Algorithm (STPA). The TPA implements the Parareal
algorithm[14, 21, 26, 34, 36], a time-parallel method to solve (1), while the
STPA also employs a parallel domain decomposition strategy in space [30, 37,
39, 40] in addition to the time-parallelism of the Parareal method.

The analysis presented here extends earlier work on the _a posteriori_ error
analysis of the Parareal algorithm in [7], and incorporates the work of [4] on
the _a posteriori_ error analysis of overlapping Schwarz domain decomposition
algorithms to provide an estimate of the spatial discretization error. One
significant difference between the earlier work and the current work is that
we analyze the effects of using an iterative method to solve the discrete
equations arising from an implicit time integration scheme. The analysis
builds on early work of adjoint-based _a posteriori_ error analysis [1, 2, 17,
27] and more recent work addressing iterative and multi-rate methods [6, 10,
18].

###  1.1 Notation

Our analysis requires consideration of both exact and discrete solutions.
Analytical solutions are indicated using a lower case letter, discrete
solutions with an upper case letter. Solutions of the PDE will be designated
$u$ and adjoint solutions $\phi$. A hat indicates a coarse scale entity.
Errors are indicated using $e$ and residuals with $\mathcal{R}$. A dot
indicates a partial differential in time. We use $N$ to indicate the number of
finite elements, $Q$ to indicate the degree of interpolation, $P$ to indicate
the number of temporal or spatial subdomains, and $K$ to indicate the number
of iterations. A subscript $t$ for these symbols indicates a temporal
quantity, a subscript $s$ indicates a spatial quantity. A full list of
variables and discretization choices appears in A.

For a function $w(x,t)$, we denote $w(t):=w(\cdot,t)$, i.e., $w(t)$ refers to
a function of space only for a fixed time $t$. We use the notation
$(w,v)(t):=(w(t),v(t))=\int_{\Omega}w(x,t)v(x,t)\,{\rm d}x$, so for example
the QoI (2) may be expressed as $\mathcal{Q}(u)=(\psi,u)(T)$. Finally, we let
$\langle\cdot\rangle_{I}$ indicate integration over time interval $I$, i.e.
$\langle\cdot\rangle_{I}=\int_{I}\,\cdot\,{\rm d}t$.

###  1.2 Outline

We begin the analysis from the perspective of Parareal integration in time,
presented in [7]. The Parareal approach utilizes both coarse and fine scale
discretizations in time. Extending this analysis to PDEs requires
corresponding coarse and fine discretizations in space. Since we implement a
domain decomposition strategy in space, we employ iterative solution methods
in both time and space.

We introduce two algorithms in section 2: a Time Parallel Algorithm (TPA), and
a Space-Time Parallel Algorithm (SPTA). The TPA employs the time-parallel
algorithm, Parareal, for the solution of PDEs, whereas the SPTA additionally
utilizes additive Schwarz domain decomposition to achieve parallelization in
space. This section also contains precise details of the spatial and temporal
discretizations. We recall some essential mathematical results needed for the
subsequent analysis for the two algorithms in section 3. We adjust the
previous _a posteriori_ error analysis for Parareal time integration to
analyze the TPA in in section 4. We extend this analysis to account for the
errors arising using additive Schwarz domain decomposition in section 5.
Numerical results supporting the accuracy of the _a posteriori_ error
estimates are provided in section 6.

##  2 Discretizations: The “Time Parallel Algorithm” and the “Space-Time
Parallel Algorithm”

We consider two algorithms: (1) An algorithm which is parallel in time only,
referred to as Time Parallel Algorithm (TPA), and (2) An algorithm which is
parallel in both space and time, referred to as Space-Time Parallel Algorithm
(STPA). Both algorithms use the Parareal algorithm for time-parallelism. The
spatial parallelization in STPA is achieved through the use of overlapping
Schwarz domain decomposition method. The discretization also involves using
the finite element method in space and the implicit Euler method in time. The
implicit Euler method is appropriate for dissipative problems modeled by (1).
We describe these algorithms in detail below.

###  2.1 Time Parallel Algorithm (TPA)

The Parareal algorithm is based on a pair of coarse and fine scale solvers,
$\widehat{G}^{p}[\alpha(x)](t)$ and $F^{p}[\alpha(x)](t)$ respectively. Let
$0=T_{0}<\cdots<T_{p-1}<T_{p}<\cdots<T_{P}=T$ be a partition of the time
domain. We view the times $T_{p}$ as “synchronization” times at which point
the coarse and fine scale solutions may exchange information. We call the time
interval $[T_{p-1},T_{p}]$ the temporal subdomain $p$, see Figure 1. Note that
the temporal subdomains are non-overlapping.

The numerical solver $F^{p}[\alpha(x)](t)$ is more accurate than
$\widehat{G}^{p}[\alpha(x)](t)$, e.g., employing a finer discretization or a
higher order method. The super-script $p$ indicates that the solvers produce a
space-time solution for the temporal subdomain $p$. That is, for
$t\in(T_{p-1},T_{p}]$, $\widehat{G}^{p}[\alpha(x)](t)$ (resp.
$F^{p}[\alpha(x)](t)$) indicates the space-time coarse scale (resp. fine
scale) solution at $t$, where $\alpha(x)$ denotes the initial conditions for
the solver at $t=T_{p-1}$. Specific examples of the solvers are presented in
2.1.1.

The standard version of the Parareal algorithm defines the solutions only at
the times $T_{p}$ and is not amenable to adjoint-based _a posteriori_
analysis, which requires solutions to be in variational form. An equivalent,
variational version of the Parareal algorithm, suitable for an adjoint-based
_a posteriori_ error analysis, is provided in Algorithm 1. The standard
Parareal algorithm is provided in B and a proof of its equivalence with the
variational version is given in [7]. Here $\widehat{U}^{p,\\{k_{t}\\}}(x,t)$
and $U^{p,\\{k_{t}\\}}(x,t)$ are the coarse and fine scale solutions for the
variational Parareal algorithm at iteration $k_{t}$ on temporal subdomain $p$.
The fine-scale Parareal solution after $k_{t}$ iterations is defined to be
$U^{\\{k_{t}\\}}=U^{p,\\{k_{t}\\}}$ for $t\in(T_{p-1},T_{p}]$. The initial
condition approximation to $u_{0}$ is referred to as $\widehat{U}_{0}$. Note
that since $\widehat{U}_{0}$ is represented in a finite dimensional space, it
is possible that $u_{0}\neq\widehat{U}_{0}$.

Algorithm 1 Variational form of the Parareal algorithm

procedure VPAR($P_{t},K_{t},\widehat{U}_{0}$ )

$C^{p,\\{0\\}}=0,p=0,\ldots,P_{t}$ $\triangleright$ Initialize corrections

for $k_{t}=1,\ldots,K_{t}$ do

$\widehat{U}^{0,\\{k_{t}\\}}(x,0):=\widehat{U}_{0}$ $\triangleright$ Set
initial conditions

for $p=1,\ldots,P_{t}$ do

$\widehat{U}^{p,\\{k_{t}\\}}(x,t)=\widehat{G}^{p}[\widehat{U}^{p-1,\\{k_{t}\\}}(x,T_{p-1})+C^{p,\\{k_{t}-1\\}}(x)](t)$
$\triangleright$ Serial computation

$U^{p,\\{k_{t}\\}}(x,t)=F^{p}[\widehat{U}^{p,\\{k_{t}\\}}(x,T_{p-1})](t)$
$\triangleright$ Parallel computation

$C^{p,\\{k_{t}\\}}(x)=U^{p,\\{k_{t}\\}}(x,T_{p})-\widehat{G}^{p}[\widehat{U}^{p,\\{k_{t}\\}}(x,T_{p-1})](T_{p})$
$\triangleright$ Update corrections

end for

end for

return $U^{p,\\{K_{t}\\}}(x,t)$

end procedure

####  2.1.1 Temporal and spatial discretizations of the fine scale and coarse
scale solvers

We consider the implicit Euler method for time discretization and finite
elements for the spatial discretization to define coarse and fine scale
solutions on $[T_{p-1},T_{p}]$. We omit the superscript $p$ in this section
when the temporal subdomain is clear.

We discretize the spatial domain $\Omega$ into a quasi-uniform triangulation
$\mathcal{T}_{h}$, where $h$ denotes the maximum diameter of the elements.
This triangulation is chosen so the union of the elements of $\mathcal{T}_{h}$
is $\Omega$ and the intersection of any two elements is either a common edge,
node, or is empty. Let ${V}_{s}^{q}$ be the standard space of continuous
piecewise polynomials of degree $q$ on $\mathcal{T}_{h}$, such that if
$U\in{V}_{s}^{q}$ then $U=0$ on $\partial\Omega$. In particular, we consider a
coarse space ${V}_{s}^{\widehat{q}_{s}}$ and fine space ${V}_{s}^{q_{s}}$ such
that $\widehat{q}_{s}<q_{s}$.

On each time subdomain $p$, the coarse temporal discretization, uses
$\widehat{N}^{p}_{t}$ time steps while the fine temporal discretization uses
$N^{p}_{t}$ time steps. Thus, the total coarse time steps over $[0,T]$ are
$\widehat{N}_{t}=\sum_{p=1}^{P_{t}}\widehat{N}^{p}_{t}$, and the total fine
time steps are $N_{t}=\sum_{p=1}^{P_{t}}N^{p}_{t}$. On the coarse scale, we
partition $[T_{p-1},T_{p}]$ as
$T_{p-1}=\widehat{t}^{p}_{0}<\cdots<\widehat{t}^{p}_{\widehat{N}^{p-1}_{t}}<\widehat{t}^{p}_{\widehat{N}^{p}_{t}}=T_{p}$,
as shown in Figure 1, with
$\widehat{\mathcal{I}}^{p}_{\hat{n}}=[\widehat{t}^{p}_{\hat{n}-1},\widehat{t}^{p}_{\hat{n}}]$.
Let
$\widehat{\mathcal{T}}^{p}=\\{\widehat{\mathcal{I}}^{p}_{1},\cdots,\widehat{\mathcal{I}}^{p}_{\hat{n}},\cdots,\widehat{\mathcal{I}}^{p}_{\widehat{N}^{p}_{t}}\\}$.
On the fine scale we introduce a finer discretization
$T_{p-1}=t^{p}_{0}<\cdots<t^{p}_{N^{p-1}_{t}}<t^{p}_{N^{p}_{t}}=T_{p}$ as
shown in Figure 1. We let $\mathcal{I}^{p}_{n}=[t^{p}_{n-1},t^{p}_{n}]$ and
$\mathcal{T}^{p}=\\{\mathcal{I}^{p}_{1},\cdots,\mathcal{I}^{p}_{n},\cdots,\mathcal{I}^{p}_{N^{p}_{t}}\\}$.

Figure 1: Top: Subdomain $[T_{p-1},T_{p}]$. Middle: Stage 1 discretization for
$[T_{p-1},T_{p}]$. Bottom: Stage 2 discretization for $[T_{p-1},T_{p}]$.

We now use the implicit Euler method on the coarse time step
$\widehat{\mathcal{I}}^{p}_{\hat{n}}$: given the approximate solution
$\widehat{U}_{n-1}$ at $\widehat{t}^{p}_{\hat{n}-1}$, compute approximation
$\widehat{U}_{n}\in{V}_{s}^{\widehat{q}_{s}}$ at
$\widehat{t}^{p}_{\hat{n}}$by,

$(\widehat{U}_{n},v)=(\widehat{U}_{n-1},v)+\Delta\hat{t}^{p}_{n}\left[-a(\widehat{U}_{n},v)+l(v)\right],$
(3)

for all $v\in{V}_{s}^{\widehat{q}_{s}}$,where
$\Delta\hat{t}^{p}_{n}=(\widehat{t}^{p}_{n}-\widehat{t}^{p}_{n-1})$. Moreover,
in the context of Parareal algorithm, we have
$\widehat{U}_{n-1}\in{V}_{s}^{\widehat{q}_{s}}$ for all ${\hat{n}}$ except
possibly for ${\hat{n}}=1$, where $\widehat{U}_{0}$ (the initial condition for
the coarse scale solver) may belong to ${V}_{s}^{q}$.

Similarly, the implicit Euler method on the fine time step
$\mathcal{I}^{p}_{n}$ is: given the approximate solution $U_{n-1}$ at
$t^{p}_{n-1}$, compute approximation $U_{n}\in{V}_{s}^{q}$ at $t^{p}_{n}$ by,

$(U_{n},v)=(U_{n-1},v)+\Delta t^{p}_{n}\left[-a(U_{n},v)+l(v)\right],$ (4)

for all $v\in{V}_{s}^{q}$, where $\Delta t^{p}_{n}=(t^{p}_{n}-t^{p}_{n-1})$.
In the context of Parareal algorithm, we have $U_{n-1}\in{V}_{s}^{q}$ for all
$n$ except for $n=1$ for which $U_{0}$ (the initial condition for the fine
scale solver) belongs to ${V}_{s}^{\widehat{q}_{s}}$.

In the TPA, we assume that equations 3 and 4 are solved exactly (up to
numerical precision).

###  2.2 Space-Time Parallel Algorithm (STPA)

The STPA introduces spatial parallelization (in addition to temporal
parallelization) by using a Schwarz domain decomposition iteration at the fine
scale. It is assumed in the spirit of the Parareal algorithm that the coarse
scale is solved implicitly without needing any spatial domain decomposition.

####  2.2.1 Schwarz Domain Decomposition at a single time step

The additive Schwarz Domain decomposition iteration is employed to solve the
spatial problem in (4). Similar to §2.1.1, we omit the subscript $p$ where
this is clear. We rewrite the problem (4) on subinterval $\mathcal{I}^{p}_{n}$
as: find $U_{n}\in{V}_{s}^{q_{s}}$ such that,

$B^{n}(U_{n},v)=l^{n}(v)$ (5)

for all $v\in{V}_{s}^{q}$, where the bilinear form $B^{n}(\cdot,\cdot)$ and
the linear form $l^{n}(\cdot)$ are defined by,

$\displaystyle B^{n}(u,v)$ $\displaystyle=(u,v)+\Delta t^{p}_{n}a(u,v),$ (6)
$\displaystyle l^{n}(v)$ $\displaystyle=(U_{n-1},v)+\Delta t^{p}_{n}l(v).$

Assume that we have $P_{s}$ overlapping subdomains
$\Omega_{1},\cdots,\Omega_{P_{s}}$ on $\Omega$, such that for any
$\Omega_{i}$, there exists a $\Omega_{j}$, $i\neq j$, for which
$\Omega_{i}\cap\Omega_{j}\neq\emptyset$ and $\cup_{i}\Omega_{i}=\Omega$. Let
$\mathcal{T}_{h,i}\equiv\mathcal{T}_{h}\cap\Omega_{i}$. We further assume that
the triangulation $\mathcal{T}_{h}$ is consistent with the domain
decomposition in the sense that for any $T_{s}\in\mathcal{T}_{h}$, if
$T_{s}\cap\Omega_{j}\neq\emptyset$ then $T_{s}\subset\Omega_{j}$. We let
$(\cdot,\cdot)_{ij}$ represent the $L_{2}(\Omega_{i}\cap\Omega_{j})$ inner
product. We denote by $B^{n}_{i}(\cdot,\cdot)$ the restriction of
$B^{n}(\cdot,\cdot)$ to $\Omega_{i}$ and $B^{n}_{ij}(\cdot,\cdot)$ the
restriction of $B^{n}(\cdot,\cdot)$ to $\Omega_{i}\cap\Omega_{j}$. Similarly,
we let $l^{n}_{i}(\cdot)$ be the restriction of $l^{n}(\cdot)$ to
$\Omega_{i}$. The additive Schwarz method for the solution of (5) is given in
Algorithm 2. The domain decomposition solution for $n=1,\ldots,N^{p}_{t}$ are
$U_{n}^{\\{K_{s}\\}}=$ ASDD($B^{n},l^{n},K_{s},P_{s},U_{n-1}^{\\{K_{s}\\}}$).
Note that the initial guess to the algorithm is the solution from the previous
time step $U_{n-1}^{\\{K_{s}\\}}$. In the algorithm $V_{{i}}^{q_{s}}$
represents the space of continuous piecewise polynomials of degree $q_{s}$ on
$\mathcal{T}_{h,i}$ such that if $U\in V_{{i}}^{q_{s}}$ then $U=0$ on
$\partial\Omega_{i}$, and $V_{{i},{k_{s}}}^{q_{s}}$ represents the space of
continuous piecewise polynomials of degree $q_{s}$ on $\mathcal{T}_{h,i}$ such
that if $U\in V_{{i},{k_{s}}}^{q_{s}}$ then $U=0$ on $\partial\Omega$ and
$U=U_{n}^{\\{k_{s}\\}}$ on $\partial\Omega_{i}\setminus\partial\Omega$. The
Richardson parameter $\tau$, is needed to ensure that the iteration converges
[37].

Algorithm 2 Overlapping additive Schwarz domain decomposition

procedure ASDD($B^{n},l^{n},K_{s},P_{s},U_{n}^{\\{0\\}}$)

for $k=0,1,2,\dots,K_{s}-1$ do

for $i=1,2,\dots,P_{s}$ do

Find $\widetilde{U}_{n,i}^{\\{k_{s}+1\\}}\in V_{{i},{k_{s}}}^{q_{s}}$ such
that

$B^{n}_{i}\big{(}\widetilde{U}_{n,i}^{\\{k_{s}+1\\}},v\big{)}=l^{n}_{i}(v),\quad\forall
v\in V_{{i}}^{q_{s}}.$ (7)

Let  $U_{n}^{\\{k_{s}+1\\}}=(1-\tau
P_{s})U_{n}^{\\{k_{s}\\}}+\tau\left(\sum_{i=1}^{P_{s}}{\Pi_{i}}\widetilde{U}_{n,i}^{\\{k_{s}+1\\}}\right)\;\hbox{where}\quad\Pi_{i}\widetilde{U}_{n,i}^{\\{k_{s}+1\\}}=\left\\{\begin{gathered}\begin{aligned}
&\widetilde{U}_{n,i}^{\\{k_{s}+1\\}},\;&&\hbox{on }\overline{\Omega}_{i},\\\
&U_{n}^{\\{k_{s}\\}},&&\hbox{on
}\Omega\backslash\overline{\Omega}_{i}.\end{aligned}\end{gathered}\right.$ (8)

end for

end for

return $U_{n}^{\\{K_{s}\\}}$

end procedure

####  2.2.2 The global STPA solution

In the context of the Parareal algorithm (Algorithm 1), we denote the solution
by $U^{p,\\{k_{t},K_{s}\\}}$, where we have now added $K_{s}$ to the
superscript to indicate the dependence of the solution on domain decomposition
iterations at the $k_{t}$th Parareal iteration. The global STPA solution on
$[0,T]$ after $k_{t}$ Parareal iterations and $K_{s}$ domain decomposition
iterations is defined as $U^{\\{k_{t},K_{s}\\}}=U^{p,\\{k_{t},K_{s}\\}}$ for
$t\in(T_{p-1},T_{p}]$.

##  3 Variational methods and preliminary _a posteriori_ error analysis

In this section, we introduce basic elements needed for the error analysis of
the numerical value of the QoI (see (2)) computed using the fine scale
solution from either the the Time Parallel Algorithm or the Space-Time
Parallel Algorithm. The main tools used in the analysis are adjoint problems
coupled with an appropriate definition of residuals. This kind of analysis
relies on the solution being in a variational format. Hence, we present the
implicit Euler method as a variational method in §3.3, and point out a
fundamental property of variational methods which is used subsequently in the
analysis of TPA and STPA.

###  3.1 Variational methods

We define two variational methods; the continuous Galerkin method and the
discontinuous Galerkin method for the fine scale; the extension to the coarse
scale is analogous. Note that the usage of “continuous” or “discontinuous”
refers to the property of the solution in time only; in space we always use
continuous finite elements.

On the space-time slab $S^{p}_{n}=\Omega\times\mathcal{I}^{p}_{n}$ we define
the space

${W}_{p,n}^{q_{t},q_{s}}:=\\{w(x,t):w(x,t)=\sum_{j=0}^{\widehat{q}_{t}}t^{j}v_{j}(x),v_{j}\in{V}_{s}^{q_{s}},(x,t)\in
S^{p}_{n}\\},$ (9)

The continuous Galerkin method, cG($q_{t}$) is: find $U_{CG}$ such that it is
continuous on $S^{p}_{n}$, its restriction to $S^{p}_{n}$ belongs to
${W}_{p,n}^{q_{t},q_{s}}$ and also satisfies,

$\langle(\dot{U}_{CG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}}=\langle-a(U_{CG},v)+l(v)\rangle_{{\mathcal{I}}^{{p}}_{n}}$
(10)

for all $v\in{W}_{p,n}^{q_{t}-1,q_{s}}$.

We denote the jump across $t^{p}_{n}$ as $[w]_{n,p}=w_{n,p}^{+}-w_{n,p}^{-}$,
where $w_{n,p}^{\pm}=\lim_{s\rightarrow(t^{p}_{n})^{\pm}}w(s)$. The
discontinuous Galerkin method, dG($q_{t}$) is: find $U_{CG}$ its restriction
to $S^{p}_{n}$ belongs to ${W}_{p,n}^{q_{t},q_{s}}$ and also satisfies,

$\langle(\dot{U}_{DG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}}+([U_{DG}]_{n-1,p},v_{n-1,p}^{+})=\langle-a(U_{DG},v)+l(v)\rangle_{{\mathcal{I}}^{{p}}_{n}}$
(11)

for all $v\in{W}_{p,n}^{q_{t},q_{s}}$. Note that the solution may be
discontinuous in at the time nodes $t^{p}_{n}$.

###  3.2  A posteriori error analysis for variational methods on a single time
subdomain

Variational methods satisfy an important property that underpins the _a
posteriori_ error analysis presented here. Let $\phi\in
L^{2}(T_{p-1},T_{p};H^{1}_{0}(\Omega))$ satisfy,

$(-\dot{\phi},v)=-a(v,\phi)\quad\forall v\in H^{1}_{0}(\Omega)\text{ and
}t\in(T_{p-1},T_{p}]$ (12)

Equation (12) is referred to as an adjoint equation. Next we quantify the
accumulation of error contributions.

######  Lemma 1.

Let $e_{\mathrm{Var}}=u-U_{\mathrm{Var}}$, where $u$ is the solution to (1),
$\mathrm{Var}\in\\{CG,DG\\}$ and $U_{\mathrm{Var}}$ is either the cG($q_{t}$)
or the dG($q_{t}$) solution on a time subdomain $p$. Then,

$(\phi,e_{\mathrm{Var}})(T_{p})=(\phi,e_{\mathrm{Var}})(T_{p-1})+\mathcal{R}_{\mathrm{Var}}^{p}(U_{\mathrm{Var}},\phi)$
(13)

where $\phi$ is defined in (12) and
$\mathcal{R}_{\mathrm{Var}}^{p}(U_{\mathrm{Var}},\phi)$ is the adjoint-
weighted space-time residual of the solution $U_{\mathrm{Var}}$ on time
subdomain $p$ defined as

$\mathcal{R}_{\mathrm{Var}}^{p}(U_{\mathrm{Var}},\phi)=\begin{cases}\sum_{n=1}^{{N}_{p}}\left[\langle
l(v)-a(U_{CG},v\rangle_{{\mathcal{I}}^{{p}}_{n}}-\langle(\dot{U}_{CG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}}\right],\qquad&\text{
if }\mathrm{Var}=CG,\\\ \sum_{n=1}^{{N}_{p}}\left[\langle
l(v)-a(U_{DG},v\rangle_{{\mathcal{I}}^{{p}}_{n}}-\langle(\dot{U}_{DG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}}-([U_{DG}]_{n-1,p},\phi_{n-1,p}^{+})\right],\;&\text{
if }\mathrm{Var}=DG.\end{cases}$ (14)

###### Proof.

We prove property (13) for the dG method, the proof for the cG method is
similar. In the notation of the dG method, we have

$e_{\mathrm{Var}}(T_{p-1})=e_{\mathrm{Var}}^{+}(t^{p}_{0})\qquad\text{and}\qquad
e_{\mathrm{Var}}(T_{p})=e_{\mathrm{Var}}^{-}(t^{p}_{{N}_{p}}).$ (15)

Consider

$\langle(-\dot{\phi},e_{DG})\rangle_{{\mathcal{I}}^{{p}}_{n}}=-\langle
a(e_{DG},\phi)\rangle_{{\mathcal{I}}^{{p}}_{n}}.$ (16)

The left-hand side is

$\langle(-\dot{\phi},e_{DG})\rangle_{{\mathcal{I}}^{{p}}_{n}}=-(\phi,e_{DG}^{-})(t^{p}_{n})+(\phi,e_{DG}^{+})(t^{p}_{n-1})+\langle(\phi,\dot{u}-\dot{U}_{DG})\rangle_{{\mathcal{I}}^{{p}}_{n}}.$
(17)

Combining (16) and (17) with (1) leads to

$\displaystyle(\phi,e_{DG}^{-})(t^{p}_{n})$
$\displaystyle=(\phi,e_{DG}^{+})(t^{p}_{n-1})+\langle
l(v)-a(U_{DG}^{p},v\rangle_{{\mathcal{I}}^{{p}}_{n}}-\langle(\dot{U}_{DG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}},$
(18)
$\displaystyle=(\phi,e_{DG}^{-})(t^{p}_{n-1})-(\phi_{n-1,p}^{+},[U_{DG}]_{n,p})+\langle
l(v)-a(U_{DG}^{p},v\rangle_{{\mathcal{I}}^{{p}}_{n}}-\langle(\dot{U}_{DG},v)\rangle_{{\mathcal{I}}^{{p}}_{n}}$

Summing (18) over $n=1$ to $n=N^{p}_{t}$ and using (15),

$\displaystyle(\phi,e_{DG})(T_{p})$
$\displaystyle=(\phi,e_{DG})(T_{p-1})+\mathcal{R}^{p}(U_{DG},\phi),$ (19)

which is the desired result. ∎

Note that a similar result, with similar derivation, also holds for the coarse
scale.

###  3.3 _A posteriori_ error analysis for implicit Euler on a single time
subdomain

It is well known that the implicit Euler method with a certain choice of
quadrature is nodally equivalent to the $dG(0)$ method and hence may be
considered a variational method [16]. For completeness, we show this
equivalence for the $dG(0)$ in (11) with the implicit Euler method in (4). For
the $dG(0)$ method, we have $\dot{U}_{DG}=0$, since the solution is constant
in time over $\mathcal{I}^{p}_{n}$ for a fixed spatial point. Moreover,
${W}_{p,n}^{0,q_{s}}$ is simply ${V}_{s}^{q_{s}}$. Now, using these facts,
along with the right-hand rule of evaluating quadrature in (11),

$((U_{DG})_{n,p}^{-}-(U_{DG})_{n,p}^{+},v_{n-1,p}^{+})=\Delta
t^{p}_{n}\left[a((U_{DG})_{n,p}^{-},v_{n,p}^{-})+l(v)\right]$ (20)

Now, since $v\in{W}_{p,n}^{0,q_{s}}$ is constant in time, we have
$v_{n-1,p}^{+}=v_{n,p}^{-}=v$. Finally, identifying $(U_{DG})_{n,p}^{-}$ with
$U_{n}$ and $(U_{DG})_{n,p}^{+}$ with $U_{n-1}$ completes the equivalence of
the two equations.

Now we define coarse and fine scale analogues to (14) for the implicit Euler
method. The fine scale residuals are,

$\displaystyle\mathcal{R}^{p}_{n}(w,v)$ $\displaystyle=\langle
l(v)-a(w,v\rangle_{{\mathcal{I}}^{{p}}_{n}}-([w]_{n-1,p},v^{+}),$ (21)
$\displaystyle\mathcal{R}^{p}(w,v)$
$\displaystyle=\sum_{n=1}^{{N}_{p}}\mathcal{R}^{p}_{n}(w,v).$ (22)

The coarse scale residuals are,

$\displaystyle\widehat{\mathcal{R}}^{p}_{\hat{n}}(w,v)$ $\displaystyle=\langle
l(v)-a(w,v\rangle_{\widehat{\mathcal{I}}^{p}_{\hat{n}}}-([w]_{\hat{n}-1,p},v^{+}),$
(23) $\displaystyle\widehat{\mathcal{R}}^{p}(w,v)$
$\displaystyle=\sum_{\hat{n}=1}^{\hat{N}_{p}}\widehat{\mathcal{R}}^{p}_{\hat{n}}(w,v).$
(24)

Notice that there is no time derivative in the residuals because, as discussed
earlier, those terms are 0 for the $dG(0)$ method.

The analogues of property (13) for the coarse and fine scale solutions in
Algorithm 1 follow immediately. Let $U^{p}=F^{p}[U^{p}(T_{p-1})]$ be the fine-
scale solution on time domain $p$, and $e^{p}=u-U^{p}$. Then,

$(\phi,e^{p})(T_{p})=(\phi,e^{p})(T_{p-1})+\mathcal{R}^{p}(U^{p},\phi).$ (25)

Let $\widehat{U}^{p}=\widehat{G}[\widehat{U}^{p}(T_{p-1})]$ be the coarse
scale solution on time domain $p$. Let $\widehat{e}^{p}=u-\widehat{U}^{p}$.
Then,

$(\phi,\widehat{e}^{p})(T_{p})=(\phi,\widehat{e}^{p})(T_{p-1})+\widehat{\mathcal{R}}^{p}(\widehat{U}^{p},\phi).$
(26)

##  4 _A posteriori_ error analysis for the Time Parallel Algorithm (TPA)

In this section we derive estimates for the error in the QoI computed from the
fine scale solution of the TPA presented in §2.1. That is, we seek to estimate
$\mathcal{Q}(u-U^{\\{K_{t}\\}})=(\psi,u(T)-U^{\\{K_{t}\\}}(T))$. The analysis
for TPA is quite similar to the analysis for ODEs appearing in [7], except
that the analysis for PDEs includes the error in the initial conditions, which
is assumed to be zero for ODEs.

###  4.1 Coarse and fine scale adjoint problems

The coarse scale adjoint problem is: find $\widehat{\phi}\in
L^{2}(0,T;H^{1}_{0}(\Omega))$ such that

$\left.\begin{gathered}\begin{aligned}
(-\dot{\widehat{\phi}},v)&=-a(v,\widehat{\phi})\quad\forall v\in
H^{1}_{0}(\Omega)\text{ and }t\in(0,T],\\\
\widehat{\phi}(T)&=\psi.\end{aligned}\end{gathered}\qquad\right\\}$ (27)

Here, the term “coarse” indicates that in numerical experiments, this adjoint
is approximated using a coarse (spatial and temporal) scale space.

For $p=1,2,\ldots,P_{t}$, the fine scale adjoint problem for the $p$th time
domain $[T_{p-1},T_{p}]$ is: find $\phi^{p}\in
L^{2}(T_{p-1},T_{p};H^{1}_{0}(\Omega))$ such that

$\left.\begin{gathered}\begin{aligned}
(-\dot{\phi}^{p},v)&=-a(v,\phi^{p})\quad\forall v\in H^{1}_{0}(\Omega)\text{
and }t\in(T_{p-1},T_{p}],\\\
\phi^{p}(T_{p})&=\widehat{\phi}(T_{p}).\end{aligned}\end{gathered}\qquad\right\\}$
(28)

Notice that the fine scale adjoint problem on the $p$th time subdomain
receives initial conditions from the coarse scale adjoint problem and
therefore the adjoint on the fine temporal scale can be solved independently
on each of the $p$ time domains. However, while reducing solution cost, the
resulting discontinuities in the fine temporal scale adjoint solution
occurring at the boundaries of temporal subdomains give rise to an additional
contribution to the error. In numerical experiments experiments, the fine
scale adjoints are approximated using a fine (relative to the space used for
the approximation of the coarse scale adjoints) scale space.

We seek the error in the fine scale solution below. However, since coarse
scale solutions are used as initial conditions for fine scale integration on
each subdomain, the analysis indicates the expected interaction between errors
at the two different scales.

###  4.2 Error representations

Let $e^{\\{k_{t}\\}}=u-U^{\\{k_{t}\\}}$ and
$\widehat{e}^{\\{k_{t}\\}}=u-\widehat{U}^{\\{k_{t}\\}}$.

######  Lemma 2.

$\displaystyle(\psi,e^{\\{k_{t}\\}}(T))$
$\displaystyle=\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},\widehat{e}^{1,\\{k_{t}\\}})(0)$
(29)
$\displaystyle+\sum_{p=2}^{P_{t}}\bigg{[}(\phi^{p}-\widehat{\phi},\widehat{e}^{p-1,\\{k_{t}\\}})+(\phi^{p}-\widehat{\phi},\widehat{U}^{p-1,\\{k_{t}\\}}-\widehat{U}^{p,\\{k_{t}\\}})+(\widehat{\phi},U^{p-1,\\{k_{t}\\}}-U^{p,\\{k_{t}\\}})\bigg{]}(T_{p-1}).$

###### Proof.

By definition of $U^{\\{k_{t}\\}}$, we have
$U^{\\{k_{t}\\}}(T)=U^{P_{t},\\{k_{t}\\}}(T)=U^{P_{t},\\{k_{t}\\}}(T_{P_{t}})$.
From (25) and $\phi^{p}(T_{p})=\widehat{\phi}(T_{p})$ (see (28)) we have

$(\widehat{\phi},e^{p,\\{k_{t}\\}})(T_{p})=(\phi^{p},e^{p,\\{k_{t}\\}})(T_{p-1})+\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p}).$
(30)

Summing over all time subdomains from $p=1$ to $P_{t}$, and isolating
$(\psi,e^{P_{t},\\{k_{t}\\}}(T))=(\psi,e^{\\{k_{t}\\}}(T))$ leads to

$\displaystyle(\psi,e^{\\{k_{t}\\}}(T))$
$\displaystyle=-\sum_{p=1}^{P_{t}-1}(\widehat{\phi},e^{p,\\{k_{t}\\}})(T_{p})+\sum_{p=2}^{P_{t}}(\phi^{p},e^{p,\\{k_{t}\\}})(T_{p-1})+\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},e^{1,\\{k_{t}\\}})(0).$
(31)

Using $U^{p,\\{k_{t}\\}}(T_{p-1})=\widehat{U}^{p,\\{k_{t}\\}}(T_{p-1})$ (see
Algorithm 1),

$(\psi,e^{\\{k_{t}\\}}(T))=-\sum_{p=1}^{P_{t}-1}(\widehat{\phi},e^{p,\\{k_{t}\\}})(T_{p})+\sum_{p=2}^{P_{t}}(\phi^{p},\widehat{e}^{p,\\{k_{t}\\}})(T_{p-1})+\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},\widehat{e}^{1,\\{k_{t}\\}})(0).$
(32)

Rearranging, adding and subtracting terms as necessary, we arrive at

$\displaystyle(\psi,e^{\\{k_{t}\\}}(T))$
$\displaystyle=\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},\widehat{e}^{1})(0)$
(33)
$\displaystyle+\sum_{p=2}^{P_{t}}\bigg{[}(\phi^{p}-\widehat{\phi},\widehat{e}^{p-1,\\{k_{t}\\}})+(\phi^{p}-\widehat{\phi},\widehat{U}^{p-1,\\{k_{t}\\}}-\widehat{U}^{p,\\{k_{t}\\}})+(\widehat{\phi},U^{p-1,\\{k_{t}\\}}-U^{p,\\{k_{t}\\}})\bigg{]}(T_{p-1}).$

∎

The terms $(\phi^{p}-\widehat{\phi},\widehat{e}^{p-1,\\{k_{t}\\}})(T_{p-1})$
in the error representation (29) are not directly computable as they contain
the unknown error $\widehat{e}^{p-1}$. We solve auxiliary adjoint problems to
account for these error terms. As noted earlier, this is the tradeoff between
solving the fine-scale adjoint problem from $T$ to $0$ and solving $P$ adjoint
problems on each subdomain independently using the coarse scale adjoint
problem as “initial” conditions. Consider $(P_{t}-1)$ auxiliary (coarse scale)
QoIs as

$\displaystyle\Psi(u)=(\phi^{p}-\widehat{\phi},u),$ (34)

and $(P_{t}-1)$ auxiliary adjoint problems: find $\dot{\widehat{\phi}}_{\rm
aux}^{p}\in L^{2}(0,T_{p-1};H^{1}_{0}(\Omega))$ such that

$\left.\begin{gathered}\begin{aligned} -\dot{\widehat{\phi}}_{\rm
aux}^{p}&=-a(v,\widehat{\phi}_{\rm aux}^{p}),\qquad\quad\forall v\in
H^{1}_{0}(\Omega)\text{ and }t\in(0,T_{p-1}]\\\ \widehat{\phi}_{\rm
aux}^{p}(x,T_{p-1})&=\phi^{p}(x,T_{p-1})-\widehat{\phi}(x,T_{p-1}).\end{aligned}\end{gathered}\qquad\right\\}$
(35)

Replacing $\phi$ by $\widehat{\phi}_{\rm aux}^{p}$ in (26) on $k$th time
subdomain we have

$(\widehat{\phi}_{\rm
aux}^{p},\widehat{e}^{k,\\{k_{t}\\}})(T_{k})=(\widehat{\phi}_{\rm
aux}^{p},\widehat{e}^{k,\\{k_{t}\\}})(T_{k-1})+\widehat{\mathcal{R}}^{k}(\widehat{U}^{k,\\{k_{t}\\}},\widehat{\phi}_{\rm
aux}^{p}).$ (36)

Summing from $k=1$ to $k=p-1$ and using (35),

$(\phi^{p}-\widehat{\phi},\widehat{e}^{p-1,\\{k_{t}\\}})(T_{p-1})=\sum_{k=2}^{p-1}(\widehat{\phi}_{\rm
aux}^{p},\widehat{U}^{k-1,\\{k_{t}\\}}-\widehat{U}^{k,\\{k_{t}\\}})(T_{k-1})+\sum_{k=1}^{p-1}\widehat{\mathcal{R}}^{k}(\widehat{U}^{k,\\{k_{t}\\}},\widehat{\phi}_{\rm
aux}^{p})+(\widehat{\phi}_{\rm aux}^{p},\widehat{e}^{1,\\{k_{t}\\}})(0).$ (37)

Combining (29) and (37) yields

$\displaystyle(\psi,e^{\\{k_{t}\\}}(T))$
$\displaystyle=\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},\widehat{e}^{1,\\{k_{t}\\}})(0)$
(38)
$\displaystyle+\sum_{p=2}^{P_{t}}\bigg{[}\sum_{k=2}^{p-1}(\widehat{\phi}_{\rm
aux}^{p},\widehat{U}^{k-1,\\{k_{t}\\}}-\widehat{U}^{k,\\{k_{t}\\}})(T_{k-1})+\sum_{k=1}^{p-1}\widehat{\mathcal{R}}^{k}(\widehat{U}^{k,\\{k_{t}\\}},\widehat{\phi}_{\rm
aux}^{p})+(\widehat{\phi}_{\rm
aux}^{p},\widehat{e}^{1,\\{k_{t}\\}})(0)\bigg{]}$
$\displaystyle+\sum_{p=2}^{P_{t}}\bigg{[}(\phi^{p}-\widehat{\phi},\widehat{U}^{p-1,\\{k_{t}\\}}-U^{p,\\{k_{t}\\}})+(\widehat{\phi},U^{p-1,\\{k_{t}\\}}-U^{p,\\{k_{t}\\}})\bigg{]}(T_{p-1}).$

We summarize as Theorem 1 below.

######  Theorem 1.

$(\psi,e^{\\{k_{t}\\}}(T))=\mathcal{D}+\mathcal{A}+\mathcal{C}+\mathcal{K},$
(39)

where

$\displaystyle\mathcal{D}$ $\displaystyle=$
$\displaystyle\sum_{p=1}^{P_{t}}\mathcal{R}^{p}(U^{p,\\{k_{t}\\}},\phi^{p})+(\phi^{1},\widehat{e}^{1,\\{k_{t}\\}}))(0),$
$\displaystyle\mathcal{A}$ $\displaystyle=$
$\displaystyle\sum_{p=2}^{P_{t}}\left[\sum_{k=2}^{p-1}(\widehat{\phi}_{\rm
aux}^{p},\widehat{U}^{k-1,\\{k_{t}\\}}-\widehat{U}^{k,\\{k_{t}\\}})(T_{k-1})+\sum_{k=1}^{p-1}\widehat{\mathcal{R}}^{k}(\widehat{U}^{k,\\{k_{t}\\}},\widehat{\phi}_{\rm
aux}^{p})+(\widehat{\phi}_{\rm
aux}^{p},\widehat{e}^{1,\\{k_{t}\\}})(0)\right],$ $\displaystyle\mathcal{C}$
$\displaystyle=$
$\displaystyle\sum_{p=2}^{P_{t}}(\phi^{p}-\widehat{\phi},\widehat{U}^{p-1,\\{k_{t}\\}}-\widehat{U}^{p,\\{k_{t}\\}})(T_{p-1}),$
$\displaystyle\mathcal{K}$ $\displaystyle=$
$\displaystyle\sum_{p=2}^{P_{t}}(\widehat{\phi},U^{p-1,\\{k_{t}\\}}-U^{p,\\{k_{t}\\}})(T_{p-1}).$

We identify the following components that comprise the total error.

  1. 1.

$\mathcal{D}$ is the “standard” fine scale discretization error contribution.

  2. 2.

$\mathcal{A}$ is the temporal auxiliary error contribution arising from the
discontinuities in the fine scale adjoint solution at the synchronization
points,

  3. 3.

$\mathcal{C}$ is the temporal coarse scale error contribution from the
discontinuities in the coarse scale solution and the fine scale adjoint
solution at synchronization points,

  4. 4.

$\mathcal{K}$ is the temporal iteration error contribution arising from the
difference between the (fine scale) solutions at the synchronization points.

Numerical examples demonstrating the accuracy of this error estimator are
provided in section 6.1. While the analysis above concerns the fine scale
error, we note that the coarse scale error may also be estimated as described
in C.

###  4.3 Extension to other time integrators

The analysis in the previous section assumed the implicit Euler method as the
time integration procedure for both the coarse and fine scale solutions.
However, the analysis remains valid provided the time integration method used
at the fine and coarse scales satisfies (25) and (26). These requirements are
satisfied for variational methods as demonstrated in §3. Moreover, these
requirements are also satisfied for a number of numerical schemes which may be
interpreted as a variational methods, e.g., Crank Nicolson, BDF, IMEX, etc
[12, 11, 8, 9, 5].

##  5 _A posteriori_ analysis of the Space-Time Parallel Algorithm (STPA)

In this section we derive error estimates for the error in the QoI computed
from the fine scale solution of the STPA presented in §2.2. That is, we seek
to estimate
$\mathcal{Q}(u-U^{p,\\{k_{t},K_{s}\\}})=(\psi,u(T)-U^{p,\\{k_{t},K_{s}\\}}(T))$.
Theorem 1 also gives the error in the STPA by replacing by replacing
$U^{p,\\{k_{t}\\}}$ with $U^{p,\\{k_{t},K_{s}\\}}$. However, this result does
not account for the error due to the spatial Schwarz domain decomposition
iteration. In this section we extend the analysis to account for this source
of error.

Let $U_{n}^{p,\\{k_{t},K_{s}\\}}=U^{p,\\{k_{t},K_{s}\\}}(t^{p}_{n})$, where
$U_{n}^{p,\\{k_{t},K_{s}\\}}$ is the domain decomposition solution to (5) on
time domain $p$ for Parareal iteration $k_{t}$. More precisely,
$U_{n}^{p,\\{k_{t},K_{s}\\}}$ is the domain decomposition solution to the
problem: find $U_{n}^{p,\\{k_{t}\\}}\in{V}_{s}^{q_{s}}$ such that

$B^{n}(U_{n}^{p,\\{k_{t}\\}},v)=l^{n}(v),$ (40)

for all $v\in{V}_{s}^{q}$, where the bilinear form $B^{n}(\cdot,\cdot)$ and
the linear form $l^{n}(\cdot)$ are given in (6), and
$U_{0}^{p,\\{k_{t}\\}}=U_{0}^{p,\\{k_{t},K_{s}\\}}=\widehat{U}^{p,\\{k_{t}\\}}(x,T_{p-1})$.

For analysis purposes we also introduce an analytical solution: find
$u_{n}^{p,\\{k_{t}\\}}\in H_{0}^{1}(\Omega)$ such that,

$B^{n}(u_{n}^{p,\\{k_{t}\\}},v)=l^{n}(v)\qquad\forall v\in H^{1}_{0}(\Omega).$
(41)

Based on the discussion in section 3.3, we may also consider
$u_{n}^{p,\\{k_{t}\\}}$ in variational form.

###  5.1 Overview of the strategy

Trivially,

$\displaystyle\left(u(t^{p}_{n})-U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p}(t_{n})\right)$
$\displaystyle=\underbrace{\left(u(t^{p}_{n})-u_{n}^{p,\\{k_{t}\\}},\phi^{p}(t_{n})\right)}_{\mathcal{E}_{I}}+\underbrace{\left(u_{n}^{p,\\{k_{t}\\}}-U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p}(t_{n})\right)}_{\mathcal{E}_{II}}$
(42)

  1. 1.

We estimate $\mathcal{E}_{I}$ by interpreting $u_{n}^{p,\\{k_{t}\\}}$ as a
variational solution and using _a posteriori_ techniques developed in section
3.

  2. 2.

We estimate $\mathcal{E}_{II}$ using the _a posteriori_ error analysis for
overlapping domain decomposition presented in [4]. This analysis allows us to
further split $\mathcal{E}_{II}$ into iterative and discretization components.

###  5.2 Estimating $\mathcal{E}_{I}$

In a derivation akin to that of (18), interpreting $u_{n}^{p,\\{k_{t}\\}}$ as
a space-time solution to the dG(0) method in time leads to,

$\mathcal{E}_{I}=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n})-u_{n}^{p,\\{k_{t}\\}})=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n-1})-U_{n-1}^{p,\\{k_{t},K_{s}\\}})+\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})$
(43)

where we also used the fact that
$u_{n}^{p,\\{k_{t}\\}}(t^{p}_{n-1})=U_{n-1}^{p,\\{k_{t},K_{s}\\}}$.

###  5.3 Estimating $\mathcal{E}_{II}$

The analysis of [4] is used to estimate the error term $\mathcal{E}_{II}$. We
define the global (spatial) adjoint problem: find $\phi_{{n}}^{p}\in
H^{1}_{0}(\Omega)$ such that,

$B^{n}(v,\phi_{{n}}^{p})=(\phi^{p},v),\quad\forall v\in H^{1}_{0}(\Omega),$
(44)

and adjoint problems on the spatial subdomains: find
$\phi_{{n},i}^{p,[k_{s}]}\in H^{1}_{0}(\Omega_{i})$ such that,

$B_{i}^{n}\left(v,\phi_{{n},i}^{p,[k_{s}]}\right)=\tau\sum_{j=1}^{P_{s}}\left\\{(\phi^{p}(t_{n})|_{\Omega_{j}},v)_{ij}-B_{ij}\left(v,\sum_{l=k_{s}+1}^{K_{s}}\phi_{{n},i}^{p,[l]}\right)\right\\},\quad\forall
v\in H^{1}_{0}(\Omega_{i}),$ (45)

where $\phi^{p}(t_{n})|_{\Omega_{i}}$ is the restriction of $\phi^{p}(t_{n})$
to $\Omega_{i}$, and $B_{i}^{n}$ and $B_{ij}^{n}$ are the restrictions of
$B^{n}$ to $\Omega_{i}$ and $\Omega_{ij}$ respectively. For a fixed $k_{s}$,
the adjoint problems (45) are independent for each $i$, so
$\phi_{{n},i}^{p,[k_{s}]}$ may be computed backwards from $K_{s},\ K_{s}-1,\
K_{s}-2,\cdots,1$ in parallel analogous to the solution strategy in the
additive Schwarz Algorithm 2.

######  Theorem 2.

$\mathcal{E}_{II}=\left(u_{n}^{p,\\{k_{t}\\}}-U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p}(t_{n})\right)=\mathcal{D}_{k,n}^{p}+\mathcal{D}_{s,n}^{p}$
(46)

where

$\displaystyle\mathcal{D}_{s,n}^{p}$
$\displaystyle=\sum_{i=1}^{P_{s}}\sum_{k_{s}=1}^{K_{s}}l^{n}(\phi_{{n},i}^{p,[k_{s}]})-B_{i}^{n}(\widetilde{U}_{n,i}^{p,\\{k_{s}\\}},\phi_{{n},i}^{p,[k_{s}]})$
(47) $\displaystyle\mathcal{D}_{k,n}^{p}$
$\displaystyle=l^{n}(\phi_{{n}}^{p})-B^{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi_{{n}}^{p})-\mathcal{D}_{s,n}^{p}$
(48)

where $\widetilde{U}_{n,i}^{p,\\{k_{s}\\}}$ are the spatial subdomain
solutions on spatial subdomain $i$ defined by (7).

The proof of Theorem 2 is provided in [4]. The term $\mathcal{D}_{s,n}^{p}$
quantifies the discretization error contribution (that is, due to using the
spaces $V_{{i},{k_{s}}}^{q_{s}}$) while the term $\mathcal{D}_{k,n}^{p}$
quantifies the domain decomposition iteration error contribution (that is, due
to using a finite number $K_{s}$ iterations).

###  5.4 Estimating $\mathcal{E}_{I}+\mathcal{E}_{II}$

Combining (43) and (46),

$\displaystyle\mathcal{E}_{I}+\mathcal{E}_{II}$
$\displaystyle=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n-1})-U_{n-1}^{p,\\{k_{t},K_{s}\\}})+\mathcal{R}^{p}_{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p})$
(49)
$\displaystyle=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n-1})-U_{n-1}^{p,\\{k_{t},K_{s}\\}})+\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})+\mathcal{D}_{k,n}^{p}+\mathcal{D}_{s,n}^{p},$

To provide a comparison with the two part construction above, following the
direct approach prescribed in (18) we have,

$\mathcal{E}_{I}+\mathcal{E}_{II}=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n})-U_{n}^{p,\\{k_{t},K_{s}\\}})=(\phi^{p}(t^{p}_{n}),u(t^{p}_{n-1})-U_{n-1}^{p,\\{k_{t},K_{s}\\}})+\mathcal{R}^{p}_{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p}).$
(50)

Comparing (50) and (49),

$\mathcal{R}^{p}_{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p})=\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})+\mathcal{D}_{k,n}^{p}+\mathcal{D}_{s,n}^{p}.$
(51)

Rearranging,

$\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})=\mathcal{R}^{p}_{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p})-\mathcal{D}_{k,n}^{p}-\mathcal{D}_{s,n}^{p}.$
(52)

Since all terms on the RHS of (52) are computable, so is
$\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})$.

###  5.5 Summing contributions over time subdomain $p$

Using (22), the residual over the time domain $[T_{p-1},T_{p}]$ is

$\mathcal{R}^{p}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p})=\sum_{n=1}^{{N}_{p}}\mathcal{R}^{p}_{n}(U_{n}^{p,\\{k_{t},K_{s}\\}},\phi^{p})=\sum_{n=1}^{{N}_{p}}\left(\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p})+\mathcal{D}_{k,n}^{p}+\mathcal{D}_{s,n}^{p}\right)=\mathcal{D}_{t}^{p}+\mathcal{D}_{s}^{p}+\mathcal{D}_{k}^{p}$
(53)

where

$\mathcal{D}_{t}^{p}=\sum_{n=1}^{{N}_{p}}\mathcal{R}^{p}_{n}(u_{n}^{p,\\{k_{t}\\}},\phi^{p}),\qquad\mathcal{D}_{k}^{p}=\sum_{n=1}^{{N}_{p}}\mathcal{D}_{k,n}^{p},\qquad\mathcal{D}_{s}^{p}=\sum_{n=1}^{{N}_{p}}\mathcal{D}_{s,n}^{p}.$
(54)

###  5.6 Error Representation for Space-Time Parallel Domain-Decomposition-
Parareal Algorithm

Let $e^{\\{k_{t},K_{s}\\}}=u(T)-U^{p,\\{k_{t},K_{s}\\}}(T)$.

######  Theorem 3.

$(\psi,e^{\\{k_{t},K_{s}\\}}(T))=\mathcal{D}_{t}+\mathcal{D}_{s}+\mathcal{D}_{k}+\mathcal{A}+\mathcal{C}+\mathcal{K},$
(55)

where $\mathcal{A},\mathcal{C}$ and $\mathcal{K}$ are defined in Theorem 1 by
(1) and,

$\mathcal{D}_{t}=\sum_{p=1}^{P_{t}}\mathcal{D}_{t}^{p}+(\phi^{1},\widehat{e}^{1,\\{k_{t}\\}})(0),\qquad\mathcal{D}_{s}=\sum_{p=1}^{P_{t}}\mathcal{D}_{s}^{p},\qquad\mathcal{D}_{k}=\sum_{p=1}^{P_{t}}\mathcal{D}_{k}^{p}.$
(56)

###### Proof.

The proof follows directly from Theorem 1 by replacing $U^{p,\\{k_{t}\\}}$
with $U^{p,\\{k_{t},K_{s}\\}}$ and then using (53). ∎

We identify the following additional error components.

  1. 1.

$\mathcal{D}_{t}$ is the error contribution due to the use of implicit Euler
time integration.

  2. 2.

$\mathcal{D}_{s}$ is the error contribution due to using a finite dimensional
space for the solution of domain decomposition.

  3. 3.

$\mathcal{D}_{k}$ is the error contribution due to using a finite number,
$K_{s}$, iterations of domain decomposition.

##  6 Numerical results

We present numerical results to support the analyses of both section 4 and
section 5, highlighting the accuracy of the error estimates developed.
Accordingly, we first consider the error analysis for the Time Parallel
Algorithm of section 6.1 and then consider the effect of a further spatial
domain decomposition iteration in the Space-Time Parallel Algorithm of section
6.2. The implicit Euler method is used in the fine and coarse scale solvers
for these two sections. In section 6.3 we demonstrate the accuracy of the
error estimates when a different time integration, a cG method, is employed in
the coarse and fine scale solvers.

The accuracy of the _a posteriori_ error estimates is measured by the
effectivity ratio $\gamma$, where

$\gamma=\frac{\hbox{Estimated error}}{\hbox{True error}}.$ (57)

The bilinear and linear forms considered in the numerical examples are

$a(u,v)=-(\nabla u,\nabla v)\quad\hbox{and}\quad l(v)=(f(x),v)$ (58)

where

$f(x)=\sin(\mu\pi x)(\mu^{2}\pi^{2}\cos(\nu\pi t)-\nu\pi\sin(\nu\pi t)).$ (59)

This choice of $f(x)$ corresponds to the true solution $u(x,t)=\cos(\nu\pi
t)\sin(\mu\pi x)$. In all examples, the final time $T=2$ and the quantity of
interest at the final time $T$ is defined by the “adjoint data”

$\psi(x)=\left\\{\qquad\begin{aligned}
10000[(x-0.2)^{2}(x-0.6)^{2}],\qquad&0.2<x<0.6,\\\ 0,\qquad\qquad\qquad&\text{
otherwise}.\end{aligned}\qquad\right.$ (60)

This choice of $\psi(x)$ corresponds to a local weighted average of the
solution around $x=0.4$ at the final time. The ratio of number of fine and
coarse time steps is denoted by $r=N_{t}/\widehat{N}_{t}$. For the spatial
discretization, the number of spatial finite elements at the two scales
remained fixed, i.e., $N_{s}=\widehat{N}_{s}$, but different degrees of
interpolation were employed according to whether the corresponding temporal
integration was performed on the coarse or fine scales. A lower degree of
spatial interpolation, $\widehat{q}_{s}$, was employed when constructing the
solution on the coarse scale, and a higher degree of spatial interpolation,
$q_{s}$, was employed on the same spatial mesh when constructing the solution
at the fine scale, i.e., $\hat{q}_{s}<q_{s}$. The notation used in this
section is also summarized in A.

Two different time steps were employed for the temporal integration, a coarse
time step and a smaller, fine time step. The implicit Euler method was chosen
as the time integrator for the coarse and fine scale solvers $\widehat{G}$ and
$F^{p}$ in §6.1 and §6.2. A continuous Galerkin method was used to obtain the
results in §6.3. The adjoint solutions in (27), (28), (35), (44) and (45) need
to be approximated. In all cases, the same temporal and spatial meshes were
used as those used to compute the numerical solutions, however, higher degree
approximations were employed to obtain accurate estimates. The adjoint
solutions corresponding to (27), (28)and (35) were approximated on the space-
time slab $S^{p}_{n}$ using the cG(3) method in the space ${W}_{p,n}^{3,3}$
(see (9) for the definition of this space). The adjoint solutions (44) and
(45), needed for the analysis of STPA, were approximated by using a spatial
finite element employing continuous piecewise cubic polynomial functions.

###  6.1 Paralellism in time only: Time Parallel Algorithm

We first consider the _a posteriori_ error estimate given by Theorem 1, Eq.
(1) in section 4. The Parareal algorithm presents numerous discretization
choices. We investigate the effect of the number of Parareal iterations in
section 6.1.1, and the effect of the number of temporal subdomains in section
6.1.2. We then consider the effect of the fine time scale and the coarse time
scale in sections 6.1.3 and 6.1.4 respectively. In sections 6.1.1 to 6.1.4, we
ensure the temporal errors dominate spatial errors by choosing a modestly
large number of spatial elements, so that the effects of changes to the
spatial discretization parameters can be observed. Finally we consider the
effect of the spatial discretization in section 6.1.5. In all cases the
effectivity ratio of the error estimator is 1.00. For these examples
$\nu=4,\mu=1$ in the function $f(x)$ in equation (59).

The examples not only provide convincing evidence of the accuracy of the error
estimator, but illustrate the importance of identifying and estimating
distinct error contributions. We observe fortuitous cancellation in the first
example so that further refinement actually increases the error. The second
example demonstrates a situation in which the refinement strategy has no
effect on the overall error since it does not address the dominant source of
error. Later examples show that when the dominant error term is targeted and
it is at least an order of magnitude larger than all the other error
contributions, refinements strategies have the anticipated effect on the
overall error.

####  6.1.1 Effect of the number of Parareal iterations

The total error in Table 1 initially decreases with the number of Parareal
iterations ($K_{t}$), but then increases somewhat. The initial decrease is
expected, since for $K_{t}=1$, the iteration error $\mathcal{K}$ is the
dominant source of error, and hence increasing the number of iterations leads
to a decrease in this iteration error, and hence the overall error as well.
However, after a single Parareal iteration, the iteration error $\mathcal{K}$
is no longer the dominant error component. Rather the discretization error
$\mathcal{D}$, which remains essentially constant during the iterative
process, becomes the dominant error. Moreover, the discretization and
iteration errors have opposite signs, and there is a fortuitous cancellation
of error for $K_{t}=2$ between these two terms. When the number of iterations
is increased to 3, $\mathcal{K}$ decreases as expected, but so does the
cancellation of error between this term and $\mathcal{D}$, and hence the total
error shows a modest increase.

$K_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
1 | -1.02e-01 | 1.00e+00 | 5.10e-02 | -1.53e-01 | 0.00e+00 | -2.30e-06
2 | 4.39e-02 | 1.00e+00 | 5.82e-02 | -1.43e-02 | 2.86e-06 | -2.91e-06
3 | 5.73e-02 | 1.00e+00 | 5.89e-02 | -1.59e-03 | 3.11e-06 | -2.96e-06
Table 1: $\widehat{N}_{t}$ = 20,$r$ = 16, $P_{t}$ = 10, $\widehat{N}_{s}$ =
20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.1.2 Effect of the number of temporal subdomains

Increasing the number of temporal subdomains ($P_{t}$) is not expected to
affect the discretization error $\mathcal{D}$ which is determined by the
spatial and temporal element scales, but is expected to increase the iteration
error $\mathcal{K}$. Both are supported by the results in Table 2 below. The
discretization error remains the largest error for all choices of the number
of temporal subdomains.

$P_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
2 | 1.16e-01 | 1.00e+00 | 1.16e-01 | 3.53e-07 | -5.81e-09 | -3.58e-08
5 | 1.16e-01 | 1.00e+00 | 1.16e-01 | 3.68e-05 | 6.18e-09 | 1.02e-07
10 | 1.12e-01 | 1.00e+00 | 1.15e-01 | -3.26e-03 | 1.52e-08 | 9.56e-08
Table 2: $\widehat{N}_{t}$ = 40,$r$ = 4, $K_{t}$ = 2, $\widehat{N}_{s}$ = 20,
$\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.1.3 Effect of the fine time scale

As is evident in Table 3, increasing the ratio $r$ between the temporal fine
and coarse scales reduces the discretization error $\mathcal{D}$. Since this
is the dominant error, this also leads to a decrease in the total error.

r | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
2 | 7.21e-01 | 1.00e+00 | 7.31e-01 | -9.49e-03 | 1.78e-04 | -3.77e-04
4 | 3.82e-01 | 1.00e+00 | 4.06e-01 | -2.31e-02 | 2.75e-04 | -4.08e-04
Table 3: $\widehat{N}_{t}$ = 10,$P_{t}$ = 10, $K_{t}$ = 2, $\widehat{q}_{t}$ =
1, $q_{t}$ = 1, $\widehat{N}_{s}$ = 20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.1.4 Effect of the coarse time scale

The results in Table 4 demonstrate that improving the accuracy of the coarse
temporal solution reduces _all_ components of the error, provided that the
temporal errors dominate the spatial errors.

$\widehat{N}_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
10 | 7.21e-01 | 1.00e+00 | 7.31e-01 | -9.49e-03 | 1.78e-04 | -3.77e-04
20 | 4.09e-01 | 1.00e+00 | 4.13e-01 | -3.91e-03 | 1.42e-06 | -2.94e-06
Table 4: $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2, $\widehat{N}_{s}$ = 20,
$\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.1.5 Effect of the spatial scale

For the numerical results presented in Table 5 we have increased the number of
coarse temporal elements to $100$, $r$ to $8$, and the number of Parareal
iterations to $6$, in order to ensure the temporal errors are small compared
with the spatial errors. While decreasing the discretization error as
anticipated, improving the spatial accuracy is also seen to decrease the
coarse temporal and auxiliary errors (though not monotonically) since these
both have a spatial component.

$\widehat{N}_{s}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
5 | 6.61e-02 | 1.00e+00 | 6.61e-02 | 1.00e-10 | 1.16e-07 | 9.08e-07
10 | 3.40e-02 | 1.00e+00 | 3.40e-02 | 9.78e-11 | 2.20e-07 | 1.92e-06
20 | 2.64e-02 | 1.00e+00 | 2.64e-02 | 9.72e-11 | 2.51e-09 | 6.32e-08
Table 5: $\widehat{N}_{t}$ = 100, $r$ = 8, $P_{t}$ = 10, $K_{t}$ = 6,
$\widehat{q}_{s}$ = 1, $q_{s}$ = 1

###  6.2 Space-Time Parallel Algorithm

The following results were obtained through a combination of the Parareal
integration in time and additive Schwarz domain decomposition in space. We
decompose the error $\mathcal{D}$ in to its various components as presented in
Theorem 3. The effects of varying the fine and coarse time scales are
considered in sections 6.2.1 and 6.2.2 respectively. The number of domain
decomposition iterations is varied in section 6.2.3, the number of spatial
subdomains in section 6.2.4, and the degree of overlap of the spatial
subdomains in 6.2.5. In all examples the effectivity ratio of the error
estimator is 1.00. For these examples, we set $\nu=4,\mu=2$ in the function
$f(x)$ defined by equation (59). The Richardson parameter used in spatial
domain decomposition iterations was always set to $\tau=0.4$.

Once again the examples not only provide convincing evidence of the accuracy
of the error estimator, but illustrate the complex interplay that can occur
between different contributions to the overall error. We frequently observe
second order effects of a refinement strategy on error components other than
those the strategy is directly targeting. These second order effects are
largely inconsequential if the contributions to the overall error have widely
different magnitude, but may become important when the error components are
similar in scale and particularly when they are opposite in sign.

####  6.2.1 Effect of the fine time scale

Consistent with the results in §6.1.3, decreasing the temporal time step
decreases the temporal component $\mathcal{D}_{t}$ of the discretization
error. Notice that for this example the number of spatial elements has been
increased so that the temporal discretization error is dominant. All other
error contributions in Table 6 are largely unaffected.

r | Est. Err. | $\gamma$ | $\mathcal{D}_{t}$ | $\mathcal{D}_{s}$ | $\mathcal{D}_{k}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---|---|---
2 | 2.42e-01 | 1.00e+00 | 2.45e-01 | 6.94e-08 | 1.11e-02 | 7.98e-04 | 1.23e-02 | -2.58e-02
4 | 1.55e-01 | 1.00e+00 | 1.48e-01 | 6.93e-08 | 1.37e-02 | 1.70e-03 | 1.95e-02 | -2.65e-02
8 | 1.05e-01 | 1.00e+00 | 8.34e-02 | 7.32e-08 | 2.51e-02 | 2.12e-03 | 2.30e-02 | -2.68e-02
Table 6: $\widehat{N}_{t}$ = 10,$P_{t}$ = 10, $K_{t}$ = 2, $\widehat{N}_{s}$ =
80, $P_{s}$ = 2, $K_{s}$=8, $\beta$ = 0.2, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.2.2 Effect of the coarse time scale

Again, consistent with the results in section 6.1.4, all temporal components
of the error decrease as the coarse time scale is decreased. The spatial error
components in Table 7 are seen to be largely unaffected since the number of
spatial elements has been chosen so that temporal errors dominate.

$\widehat{N}_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}_{t}$ | $\mathcal{D}_{s}$ | $\mathcal{D}_{k}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---|---|---
10 | 2.42e-01 | 1.00e+00 | 2.45e-01 | 6.94e-08 | 1.11e-02 | 7.98e-04 | 1.23e-02 | -2.58e-02
20 | 1.55e-01 | 1.00e+00 | 1.41e-01 | 6.55e-08 | 1.46e-02 | -9.88e-07 | -5.09e-04 | -4.81e-05
40 | 1.05e-01 | 1.00e+00 | 7.95e-02 | 6.89e-08 | 2.54e-02 | -4.12e-07 | -8.73e-06 | -3.05e-07
Table 7: $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2, $\widehat{N}_{s}$ = 80, $P_{s}$ =
2, $K_{s}$=8, $\beta$ = 0.2, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2

####  6.2.3 Effect of the number of domain decomposition iterations

The spatial iteration error $\mathcal{D}_{k}$ decreases with number of domain
decomposition iterations, as shown in Table 8, while the spatial and temporal
discretization errors remain essentially constant. There is a second order
effect of decreasing the temporal iterative and coarse time scale error
contributions since these error contributions also have a spatial component.

$K_{s}$ | Est. Err. | $\gamma$ | $\mathcal{D}_{t}$ | $\mathcal{D}_{s}$ | $\mathcal{D}_{k}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---|---|---
2 | 6.61e-01 | 1.00e+00 | 2.16e-01 | 1.95e-05 | 4.49e-01 | -1.10e-04 | -4.94e-03 | -4.52e-05
6 | 1.88e-01 | 1.00e+00 | 1.45e-01 | 1.71e-05 | 4.40e-02 | -1.84e-05 | -1.31e-03 | -4.87e-05
Table 8: $\widehat{N}_{t}$ = 20, $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2,
$\widehat{N}_{s}$ = 20, $P_{s}$ = 2, $\beta$ = 0.2, $\widehat{q}_{s}$ = 1,
$q_{s}$ = 2

####  6.2.4 Effect of the number of spatial subdomains

The spatial iteration error $\mathcal{D}_{k}$ increases with number of spatial
subdomains, while the spatial and temporal discretization errors remain
essentially constant. A second order effect is again apparent in Table 9,
where the temporal iterative error is see to also increase due to its spatial
component.

$P_{s}$ | Est. Err. | $\gamma$ | $\mathcal{D}_{t}$ | $\mathcal{D}_{s}$ | $\mathcal{D}_{k}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---|---|---
2 | 7.17e-01 | 1.00e+00 | 2.25e-01 | 1.02e-06 | 4.96e-01 | -8.83e-06 | -4.51e-03 | -4.46e-05
4 | 1.23e+00 | 1.00e+00 | 3.13e-01 | 1.26e-06 | 9.18e-01 | 3.84e-05 | -5.76e-03 | -4.43e-05
Table 9: $\widehat{N}_{t}$ = 20, $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2,
$\widehat{N}_{s}$ = 40, $K_{s}$=2, $\beta$ = 0.1, $\widehat{q}_{s}$ = 1,
$q_{s}$ = 2

####  6.2.5 Effect of spatial domain overlap

Increasing the degree of overlap between the spatial domains is expected to
decrease the spatial iterative error $\mathcal{D}_{k}$ while leaving the
spatial and temporal discretization errors largely unchanged. Slightly
different behavior is observed in Table 10 for this example, perhaps because
the temporal discretization error is orders of magnitude larger than the
spatial discretization error. Never-the-less, the error estimator is accurate
and the effectivity ratio is 1.00.

$\beta$ | Est. Err. | $\gamma$ | $\mathcal{D}_{t}$ | $\mathcal{D}_{s}$ | $\mathcal{D}_{k}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---|---|---
0.1 | 7.17e-01 | 1.00e+00 | 2.25e-01 | 1.02e-06 | 4.96e-01 | -8.83e-06 | -4.51e-03 | -4.46e-05
0.2 | 6.61e-01 | 1.00e+00 | 2.16e-01 | 1.95e-05 | 4.49e-01 | -1.10e-04 | -4.94e-03 | -4.52e-05
Table 10: $\widehat{N}_{t}$ = 20, $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2,
$\widehat{N}_{s}$ = 20, $P_{s}$ = 2, $K_{s}$=2, $\widehat{q}_{s}$ = 1, $q_{s}$
= 2

###  6.3 Results for a different time integrator for the Time parallel
Algorithm

We briefly demonstrate the accuracy of the _a posteriori_ error estimates when
the continuous Galerkin method, cG($q_{t}$) (see section 3), is employed as
the time integration scheme in the fine and coarse scale solvers for the TPA.
The approximation space for the coarse and fine scales on the the space-time
slab $S^{p}_{n}$ are ${W}_{p,n}^{\widehat{q}_{t},\widehat{q}_{s}}$ and
${W}_{p,n}^{q_{t},q_{s}}$ respectively. The results of Theorems 1 and 3 remain
valid; however, the definition of the residuals involved need to be modified
to reflect the cG method. The residual for the cG method is given in (14). The
results are qualitative similar to those in section 6.1, and hence we present
them without any comment. The aim is to show that the results are still
accurate, and that the analysis is applicable to a broad class of numerical
methods. The results are given in Tables 11, 12, 13, 14 and 15.

$K_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
1 | 1.02e-01 | 1.00e+00 | -4.46e-02 | 1.47e-01 | 0.00e+00 | 1.79e-04
2 | -6.34e-02 | 1.00e+00 | -3.86e-02 | -2.48e-02 | -1.99e-04 | 1.80e-04
3 | -3.81e-02 | 1.00e+00 | -3.96e-02 | 1.50e-03 | -1.67e-04 | 1.80e-04
Table 11: $\widehat{N}_{t}$ = 10, $r$ = 4, $P_{t}$ = 10, $\widehat{q}_{t}$ = 1, $q_{t}$ = 1, $\widehat{N}_{s}$ = 20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2  $P_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
2 | -3.95e-02 | 1.00e+00 | -3.95e-02 | 3.86e-07 | 1.28e-07 | 1.20e-08
5 | -3.96e-02 | 1.00e+00 | -3.95e-02 | -2.55e-05 | 6.94e-05 | -7.27e-05
10 | -6.34e-02 | 1.00e+00 | -3.86e-02 | -2.48e-02 | -1.99e-04 | 1.80e-04
Table 12: $\widehat{N}_{t}$ = 10, $r$ = 4, $K_{t}$ = 2, $\widehat{q}_{t}$ = 1, $q_{t}$ = 1, $\widehat{N}_{s}$ = 20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2  r | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
2 | -1.71e-01 | 1.00e+00 | -1.53e-01 | -1.82e-02 | -1.57e-04 | 1.78e-04
4 | -6.34e-02 | 1.00e+00 | -3.86e-02 | -2.48e-02 | -1.99e-04 | 1.80e-04
Table 13: $\widehat{N}_{t}$ = 10, $P_{t}$ = 10, $K_{t}$ = 2, $\widehat{q}_{t}$ = 1, $q_{t}$ = 1, $\widehat{N}_{s}$ = 20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2  $\widehat{N}_{t}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
10 | -1.71e-01 | 1.00e+00 | -1.53e-01 | -1.82e-02 | -1.57e-04 | 1.78e-04
20 | -4.10e-02 | 1.00e+00 | -3.95e-02 | -1.55e-03 | -5.17e-07 | 6.50e-07
Table 14: $r$ = 2, $P_{t}$ = 10, $K_{t}$ = 2, $\widehat{q}_{t}$ = 1, $q_{t}$ = 1, $\widehat{N}_{s}$ = 20, $\widehat{q}_{s}$ = 1, $q_{s}$ = 2  $\widehat{N}_{s}$ | Est. Err. | $\gamma$ | $\mathcal{D}$ | $\mathcal{K}$ | $\mathcal{C}$ | $\mathcal{A}$
---|---|---|---|---|---|---
5 | 3.73e-02 | 1.00e+00 | 3.69e-02 | 1.68e-10 | -1.95e-05 | 8.85e-04
10 | 5.55e-03 | 1.00e+00 | 5.54e-03 | 1.40e-10 | -7.21e-07 | 4.15e-05
20 | -1.93e-03 | 1.00e+00 | -1.93e-03 | 1.33e-10 | -6.45e-07 | 2.55e-07
Table 15: $\widehat{N}_{t}$ = 20, $r$ = 6, $P_{t}$ = 10, $K_{t}$ = 6,
$\widehat{q}_{t}$ = 1, $q_{t}$ = 1, $\widehat{q}_{s}$ = 1, $q_{s}$ = 1

##  7 Conclusions and future work

We first developed an accurate adjoint-based _a posteriori error_ analysis for
the Time Parallel Algorithm, which applies the Parareal method in time for the
solution of parabolic partial differential equations. This analysis does not
seek to separate spatial and temporal sources of error, but combines the two
as “discretization” error. Additional error contributions arise due to
incomplete iteration, discontinuities in the coarse forward solution, and the
fine adjoint solution when it is solved in parallel. We then extended this
analysis to the Space-Time Parallel Algorithm by assuming that the spatial
solution is determined through a second iterative method, in this case domain
decomposition. The combined analysis is able to disaggregate the spatial and
temporal discretization errors, as well as identify additional iterative
errors resulting from the domain decomposition iteration in space. Thus the
analysis presented here provides a basis for separating discretization and
iteration errors and for estimating the effects of incomplete iteration in
both space and time. Accurate error estimates provide a foundation for
adaptivity and are essential for accurate uncertainty quantification where it
is necessary to distinguish variation due to parameter changes from variation
due to numerical errors which can also vary across the parameter domain.

We have limited the analysis to linear problems and intend to extend these
results to nonlinear problems using linearization techniques that have
previously proven effective [7]. We also investigate more sophisticated
temporal solvers than backward Euler and a simple cG method. Parallel
iterative methods for solving PDEs require a large number of discretization
choices. The error analysis developed here, which accurately distinguishes
multiple sources of error, provides a sound foundation on which to make many
of these choices. Finally, we can investigate adaptive strategies, noting the
complex interaction between error components.

Waveform relaxation methods [23, 24, 25] make a fundamentally different choice
when combining the two iterative methods of domain decomposition and Parareal
iteration. Assume we wish to solve a parabolic partial differential equation
on $\Omega\times(0,T]$ and let $\Omega_{i}\subset\Omega,i=1,\dots p$ be a set
of overlapping (spatial) subdomains. Waveform relaxation methods consider
domain decomposition as the outer iteration and employs Parareal iteration (or
some other time integration technique) to solve subproblems on each spatio-
temporal subdomain $\Omega_{i}\times(0,T],i=1,\dots p$ independently, and then
perform a domain decomposition like iteration on the $p$ spatial-temporal
blocks. Eficient implementations of waveform relaxation require additional
computation to determine Robin conditions for the boundaries of subdomains. An
analysis of this approach is saved for future work, starting with an _a
posteriori_ error analysis for waveform relaxation implementing a simple,
discontinuous Galerkin method for time integration and then extending to
Parareal integration in time.

## Acknowledgments

J. Chaudhry’s work is supported by the NSF-DMS 1720402. S. Tavener’s work is
supported by NSF-DMS 1720473. D. Estep’s work is supported by NSF-DMS 1720473
and NSERC grants.

## References

  * [1] M. Ainsworth and T. Oden.  _A posteriori_ error estimation in finite element analysis.  John Wiley-Teubner, 2000.
  * [2] R. Becker and R. Rannacher.  An optimal control approach to _a posteriori_ error estimation in finite element methods.  Acta Numerica, pages 1–102, 2001.
  * [3] Tony F Chan and Tarek P Mathew.  Domain decomposition algorithms.  Acta numerica, 3:61–143, 1994.
  * [4] Jehanzeb Chaudhry, Don Estep, and Simon Tavener.  A posteriori error analysis for schwarz overlapping domain decomposition methods.  arXiv preprint arXiv:1907.01139, 2019.
  * [5] Jehanzeb H Chaudhry and JB Collins.  A posteriori error estimation for the spectral deferred correction method.  arXiv preprint arXiv:2003.03399, 2020.
  * [6] Jehanzeb H. Chaudhry, Donald Estep, , Victor Ginting, and Simon Tavener.  A posteriori analysis of an iterative multi-discretization method for reaction-diffusion systems.  Computer Methods in Applied Mechanics and Engineering, 267(0):1 – 22, 2013.
  * [7] Jehanzeb Hameed Chaudhry, Don Estep, Simon Tavener, Varis Carey, and Jeff Sandelin.  A posteriori error analysis of two-stage computation methods with application to efficient discretization and the parareal algorithm.  SIAM Journal on Numerical Analysis, 54(5):2974–3002, 2016.
  * [8] J.H. Chaudhry, D. Estep, V. Ginting, J.N. Shadid, and S.J. Tavener.  A posteriori error analysis of IMEX multi-step time integration methods for advection-diffusion-reaction equations.  Computer Methods in Applied Mechanics and Engineering, 285:730–751, Mar 2015.
  * [9] J.H. Chaudhry, J.N. Shadid, and T. Wildey.  A posteriori analysis of an IMEX entropy-viscosity formulation for hyperbolic conservation laws with dissipation.  Applied Numerical Mathematics, 135, 2019.
  * [10] J. H. Chaudry, D. Estep, V. Ginting, and S. Tavener.  A posteriori analysis for iterative solvers for non-autonomous evolution problems.  SIAM Journal on Uncertainty Quantification, 3, 2015.
  * [11] J. Collins, D. Estep, and S.J. Tavener.  A posteriori error estimates for explicit time integration methods.  BIT Numerical Mathematics, 2014.
  * [12] J.B. Collins, D. Estep, and S.J. Tavener.  A posteriori error estimation for the Lax–Wendroff finite difference scheme.  Journal of Computational and Applied Mathematics, 263:299–311, 2014\.
  * [13] Victorita Dolean, Pierre Jolivet, and Frédéric Nataf.  An introduction to domain decomposition methods: algorithms, theory, and parallel implementation.  SIAM, 2015.
  * [14] Matthew Emmett and Michael Minion.  Toward an efficient parallel in time method for partial differential equations.  Communications in Applied Mathematics and Computational Science, 7(1):105–132, 2012.
  * [15] Matthew Emmett and Michael L Minion.  Efficient implementation of a multi-level parallel in time algorithm.  In Domain Decomposition Methods in Science and Engineering XXI, pages 359–366. Springer, 2014.
  * [16] K. Eriksson, D. Estep, P. Hansbo, and C. Johnson.  Computational Differential Equations.  Cambridge University Press, Cambridge, 1996.
  * [17] Kenneth Eriksson, Don Estep, Peter Hansbo, and Claes Johnson.  Introduction to adaptive methods for differential equations.  Acta numerica, 4:105–158, 1995.
  * [18] D Estep, V Ginting, and S Tavener.  A posteriori analysis of a multirate numerical method for ordinary differential equations.  Computer Methods in Applied Mechanics and Engineering, 223:10–27, 2012.
  * [19] Robert D Falgout, Stephanie Friedhoff, Tz V Kolev, Scott P MacLachlan, and Jacob B Schroder.  Parallel time integration with multigrid.  SIAM Journal on Scientific Computing, 36(6):C635–C661, 2014.
  * [20] Robert D Falgout, Stephanie Friedhoff, Tz V Kolev, Scott P MacLachlan, Jacob B Schroder, and Stefan Vandewalle.  Multigrid methods with space–time concurrency.  Computing and Visualization in Science, 18(4):123–143, 2017.
  * [21] M. Gander and S. Vandewalle.  Analysis of the parareal time-parallel time-integration method.  SIAM Journal on Scientific Computing, 29(2):556–578, 2007.
  * [22] Martin J Gander.  50 years of time parallel time integration.  In Multiple shooting and time domain decomposition methods, pages 69–113. Springer, 2015.
  * [23] Martin J Gander and Laurence Halpern.  Optimized schwarz waveform relaxation methods for advection reaction diffusion problems.  SIAM Journal on Numerical Analysis, 45(2):666–697, 2007.
  * [24] Martin J Gander, Yao-Lin Jiang, and Rong-Jian Li.  Parareal schwarz waveform relaxation methods.  In Domain Decomposition Methods in Science and Engineering XX, pages 451–458. Springer, 2013.
  * [25] Martin J Gander and Andrew M Stuart.  Space-time continuous analysis of waveform relaxation for the heat equation.  SIAM Journal on Scientific Computing, 19(6):2014–2031, 1998.
  * [26] MartinJ. Gander and Ernst Hairer.  Nonlinear convergence analysis for the parareal algorithm.  In Domain Decomposition Methods in Science and Engineering XVII, volume 60 of Lecture Notes in Computational Science and Engineering, pages 45–56. Springer Berlin Heidelberg, 2008.
  * [27] M.B. Giles and E. Süli.  Adjoint methods for pdes: a posteriori error analysis and postprocessing by duality.  Acta Numerica, 11(1):145–236, 2002.
  * [28] SEBASTIAN Gotschel and Michael L Minion.  An efficient parallel-in-time method for optimization with parabolic pdes.  SIAM Journal on Scientific Computing, 41(6):C603–C626, 2019.
  * [29] Alexander J Howse, Hans De Sterck, Robert D Falgout, Scott MacLachlan, and Jacob Schroder.  Parallel-in-time multigrid with adaptive spatial coarsening for the linear advection and inviscid burgers equations.  SIAM Journal on Scientific Computing, 41(1):A538–A565, 2019.
  * [30] David E. Keyes, Youcef Saad, and Donald G. Truhlar, editors.  Domain-Based Parallelism and Problem Decomposition Methods in Computational Sciences and Engineering.  SIAM, 1995.
  * [31] C. Evans Lawrence.  Partial Differential Equations, volume 19 of Graduate Studies in Mathematics.  AMS, 2010.
  * [32] Pierre-Louis Lions.  On the schwarz alternating method. iii: a variant for nonoverlapping subdomains.  In Third international symposium on domain decomposition methods for partial differential equations, volume 6, pages 202–223. SIAM Philadelphia, PA, 1990.
  * [33] Pierre-Louis Lions et al.  On the schwarz alternating method. i.  In First international symposium on domain decomposition methods for partial differential equations, volume 1, page 42. Paris, France, 1988.
  * [34] Y. Maday.  The parareal in time algorithm.  2008\.
  * [35] Yvon Maday and Gabriel Turinici.  A parareal in time procedure for the control of partial differential equations.  Comptes Rendus Mathematique, 335(4):387–392, 2002.
  * [36] Yvon Maday and Gabriel Turinici.  A parareal in time procedure for the control of partial differential equations.  Comptes Rendus Mathematique, 335(4):387–392, 2002.
  * [37] Tarek Poonithara Abraham Mathew.  Domain Decomposition Methods for the Numerical Solution of Partial Differential Equations, volume 61 of Lecture Notes in Computational Science and Engineering.  Springer Berlin Heidelberg, Berlin, Heidelberg, 2008.
  * [38] Benjamin W Ong and Jacob B Schroder.  Applications of time parallelization.  Computing and Visualization in Science, 23(1):1–15, 2020.
  * [39] Barry F. Smith, Petter E. Bjørstad, and William Gropp.  Domain Decomposition: Parallel Multilevel Methods for Elliptic Partial Differential Equations.  Cambridge University Press, 1996.
  * [40] Andrea Toselli and Olof Widlund.  Domain Decomposition Methods - Algorithms and Theory, volume 34 of Springer Series in Computational Mathematics.  Springer, 2004.
  * [41] Jinchao Xu.  Iterative methods by space decomposition and subspace correction.  SIAM review, 34(4):581–613, 1992.

##  Appendix A Notation

Variables, errors and residuals

| Solutions
---|---
$u$ | Analytic solution
$U^{p}$ | Fine scale discrete solution on the $p$th temporal subdomain
$\widehat{U}^{p}$ | Coarse scale discrete solution on the $p$th temporal subdomain
$U^{p,\\{k_{t}\\}}$ | Fine scale discrete solution on the $p$th temporal subdomain after $k_{t}$ parareal iterations
$U^{p,\\{k_{t},k_{s}\\}}$ | Fine scale discrete solution on the $p$th temporal subdomain after $k_{t}$ parareal iterations
| and $k_{s}$ domain decomposition iterations
| Errors and residuals
$e^{p}$ | Fine scale error on the $p$th temporal subdomain
$\widehat{e}^{p}$ | Coarse scale error on the $p$th temporal subdomain
$\mathcal{R}^{p}$ | Fine scale residual on the $p$th temporal subdomain
$\widehat{\mathcal{R}}^{p}$ | Coarse scale residual on the $p$th temporal subdomain
$\widehat{\mathcal{R}}_{\rm aux}^{p}$ | Auxilary adjoint problem residual on the $p$th temporal subdomain
| Adjoints
$\phi$ | Fine scale adjoint solution
$\phi^{p}$ | Fine scale adjoint solution on the $p$th temporal subdomain
$\widehat{\phi}$ | Coarse scale adjoint solution
$\widehat{\phi}_{\rm aux}^{p}$ | Auxiliary adjoint solution on the $p$th temporal subdomain
Table 16: Notation: Variables, errors and residuals

Discretization choices

| Time
---|---
$N_{t}$ | Number of fine scale temporal finite elements
$\widehat{N}_{t}$ | Number of coarse scale temporal finite elements
$r$ |  $N_{t}/\widehat{N}_{t}$. Refinement factor for time
$q_{t}$ | Interpolation degree of fine scale temporal finite elements
$\widehat{q}_{t}$ | Interpolation degree of coarse scale temporal finite elements
$P_{t}$ | Number of parareal subdomains (synchronization intervals)
$K_{t}$ | Number of parareal iterations
| Space
$N_{s}$ | Number of fine scale spatial finite elements
$\widehat{N}_{s}$ | Number of coarse scale spatial finite elements
$q_{s}$ | Interpolation degree of fine scale spatial finite elements
$\widehat{q}_{s}$ | Interpolation degree of coarse scale spatial finite elements
$P_{s}$ | Number of spatial subdomains
$K_{s}$ | Number of domain decomposition iterations
$\beta$ | Amount of spatial overlap
Table 17: Notation: Discretization choices

Time discretization

$T_{p-1}$ | Left hand end of pth time subdomain
---|---
$T_{p}$ | Right hand end of $p$th time subdomain
| Course scale solution
$\widehat{N}^{p}_{t}$ | Number of coarse scale temporal finite elements in $p$th time subdomain
$\widehat{t}^{p}_{j-1}$ | Left hand end of $j$th coarse time interval in $p$th time subdomain
$\widehat{t}^{p}_{j}$ | Right hand end of $j$th coarse time interval in $p$th time subdomain
$\widehat{\mathcal{I}}^{p}_{j}$ | $[\widehat{t}^{p}_{j-1},\widehat{t}^{p}_{j}]$
$\widehat{\mathcal{T}}^{p}$ | $\left\\{\widehat{\mathcal{I}}^{p}_{1},\dots,\widehat{\mathcal{I}}^{p}_{j},\dots,\widehat{\mathcal{I}}^{p}_{\widehat{N}^{p}_{t}}\right\\}$
| Fine scale solution
$N^{p}_{t}$ | Number of fine scale temporal finite elements in $p$th time subdomain
$t^{p}_{j-1}$ | Left hand end of $j$th fine time interval in $p$th time subdomain
$t^{p}_{j}$ | Right hand end of $j$th fine time interval in $p$th time subdomain
$\mathcal{I}^{p}_{j}$ | $[t^{p}_{j-1},t^{p}_{j}]$
$\mathcal{T}^{p}$ | $\left\\{\mathcal{I}^{p}_{1},\dots,\mathcal{I}^{p}_{j},\dots,\mathcal{I}^{p}_{N^{p}_{t}}\right\\}$
Table 18: Notation: Time discretization

##  Appendix B Standard Parareal algorithm and equivalence

The standard Parareal algorithm only defines the solutions at the times
$T_{p}$ [34]. This algorithm is given in Algorithm 3 Here
$\widetilde{U}_{p}^{(k_{t})}\in{V}_{s}^{\widehat{q}_{s}}$ and
$\overline{U}_{p}^{(k_{t})}\in{V}_{s}^{q_{s}}$ are the coarse and fine scale
solutions at $T_{p}$ at iteration $k_{t}$.

Algorithm 3 Standard form of the Parareal algorithm

procedure PAR($P_{t},K_{t},\widehat{U}_{0}$)

$C^{p,\\{0\\}}=0,p=0,\ldots,P_{t}$ $\triangleright$ Initialize corrections

for $k_{t}=1,\ldots,K_{t}$ do

$\widetilde{U}_{0}^{(k_{t})}:=\widehat{U}_{0}$ $\triangleright$ Set initial
conditions

for $p=1,\ldots,P_{t}$ do

$\widetilde{U}_{p}^{(k_{t})}=\widehat{G}^{p}\left[\widetilde{U}_{p-1}^{(k_{t})}\right](T_{p})+C_{p}^{k_{t}-1}$
$\triangleright$ Serial computation

$\overline{U}_{p}^{(k_{t})}=F^{p}\left[\widetilde{U}_{p-1}^{(k_{t})}\right](T_{p})$
$\triangleright$ Parallel computation

$C_{p}^{k_{t}}=\overline{U}_{p}^{(k_{t})}-\widehat{G}^{p}\left[\widetilde{U}_{p-1}^{(k_{t})}\right](T_{p})$
$\triangleright$ Update corrections

end for

end for

end procedure

######  Theorem 4.

The standard Parareal algorithm in Algorithm 3 and the variational version in
Algorithm 1 are equivalent in the sense that,

$\displaystyle\widetilde{U}_{p}^{(K_{t})}=$
$\displaystyle\widehat{U}^{p,\\{K_{t}\\}}(x,T_{p})+C^{p,\\{K_{t}-1\\}}$ (61)
$\displaystyle\overline{U}_{p}^{(K_{t})}=$ $\displaystyle
U^{p,\\{K_{t}\\}}(x,T_{p})$ (62) $\displaystyle C^{p,\\{K_{t}\\}}=$
$\displaystyle C^{p,\\{K_{t}\\}}$ (63)

For proof, please see [7].

##  Appendix C Coarse scale error

While the focus of the article is on the fine scale error, we briefly outline
the result for the coarse scale error for completeness. Let
$\widehat{e}^{p,\\{k_{t}\\}}=u-\widehat{U}^{p,\\{k_{t}\\}}$. Then,

######  Theorem 5.

$\displaystyle(\psi_{T},\widehat{e}^{P_{t},\\{k_{t}\\}}(T))=\sum_{p=1}^{P_{t}}\widehat{\mathcal{R}}^{p}(\widehat{U}^{p,\\{k_{t}\\}},\widehat{\phi})-\sum_{p=1}^{P_{t}-1}(\widehat{\phi}(T_{p}),C^{p,\\{k_{t}-1\\}})+(\widehat{\phi}(0),\widehat{e}^{1,\\{k_{t}\\}}(0)).$
(64)

The proof is similar to the proof of Theorem 1.