Abstract
We construct and analyze a second-order implicit–explicit (IMEX) scheme for the time integration of semilinear second-order wave equations. The scheme treats the stiff linear part of the problem implicitly and the nonlinear part explicitly. This makes the scheme unconditionally stable and at the same time very efficient, since it only requires the solution of one linear system of equations per time step. For the combination of the IMEX scheme with a general, abstract, nonconforming space discretization we prove a full discretization error bound. We then apply the method to a nonconforming finite element discretization of an acoustic wave equation with a kinetic boundary condition. This yields a fully discrete scheme and a corresponding a-priori error estimate.
Similar content being viewed by others
1 Introduction
In this paper we construct and analyze an implicit–explicit (IMEX) time integration scheme for second-order semilinear wave equations of the form
in a suitable Hilbert space. Here, A and B are unbounded operators and f is a locally Lipschitz continuous nonlinearity. The IMEX scheme is constructed as a combination of the explicit leapfrog method and the implicit Crank–Nicolson scheme. It treats the unbounded linear part of the differential equation implicitly and the nonlinear part explicitly. We show that the scheme is unconditionally stable in the sense that the time-step size is only restricted by the Lipschitz constant of f but not by the linear operators A and B.
We combine this IMEX scheme with an abstract, nonconforming space discretization within the framework of [11,12,13]. These papers provide a unified error analysis (UEA) which allows one to analyze nonconforming space discretizations of wave equations in a systematic way. Our main result is an error bound that is second order in time and contains abstract space-discretization errors. The error result can then be used to prove convergence rates for specific problems and discretizations by plugging in geometric and interpolation error results. The fully discrete scheme is very efficient. In fact, we show that one time step only requires the solution of one linear system and one application of discretizations of A, B, and f, respectively. Since the construction of the scheme is based on two second-order methods and our analysis makes use of the specific form of the method, the generalization to higher order is not straightforward and out of the scope of this paper. Higher-order IMEX schemes for second-order equations will be part of future work.
There is a rich literature on IMEX schemes for first-order equations, in particular, there is a well-developed theory for IMEX Runge–Kutta schemes [3, 6] or IMEX multistep schemes [1, 4, 8, 15], for instance. In [6, 15] an error analysis for ODEs is presented, while [1] contains discretization errors for IMEX schemes applied to conformal space discretizations of quasilinear parabolic evolution equations. IMEX schemes are used in applications, e.g., in structural dynamics and fluid-structure interaction [22], hydrodynamics [16], sea-ice dynamics [20], or atmospheric dynamics, see, e.g., [9], to mention just a few examples. There exists also a so-called Crank–Nicolson-leapfrog IMEX scheme which is obtained from a combination of the Crank–Nicolson and the leapfrog scheme for first-order equations, cf. [17, 18], and references therein. However, this scheme is not equivalent to the scheme we construct here, since the leapfrog schemes for first- and second-order equations are not equivalent and indeed have completely different stability properties. More precisely, the Crank–Nicolson-leapfrog scheme is only stable, if the explicitly treated part is skew symmetric.
To solve a second-order equation one can either reformulate it equivalently into a first-order equation and apply an IMEX scheme to it or one can design a scheme for the original second-order form. An example for the first option is the Crank–Nicolson-leapfrog IMEX scheme presented in [17, 18]. In contrast, we decided to take the second option and present a new scheme which to the best of our knowledge was not considered in literature so far. In fact, we are not aware of any IMEX scheme exploiting the special structure of second-order equations. The advantage of this approach is the efficiency of the scheme, as will be discussed in detail in Sect. 2.2.
The main contribution of this paper is to provide a full discretization error analysis of semilinear wave equations in a quite general framework. So far, such an error analysis does not even exist for the Crank–Nicolson scheme, which is also covered as a byproduct of our analysis of the IMEX scheme. The challenge of such a rigorous analysis is that it applies to abstract, non-conforming space discretizations of semilinear wave-type equations.
As an application of our abstract theory, we consider an acoustic wave equation with kinetic boundary conditions that fits into the abstract setting, cf., [13]. Kinetic boundary conditions are a special case of dynamic boundary conditions that contain tangential derivatives and are intrinsically posed on domains with (piecewise) smooth and therefore possibly curved boundaries. Hence, the spatial discretization has to be done on an approximated domain rendering the discretization nonconforming.
The paper is organized as follows: in Sect. 2 we present the problem setting, introduce the IMEX scheme for second-order wave equations, and state a second-order error bound for the time discretization error. In Sect. 3 we briefly recall the UEA and present the fully discrete scheme as a combination of the IMEX scheme with a general space discretization. Afterwards we state and prove the main result, namely the abstract full discretization error bound. Finally, in Sect. 4, we consider a semilinear acoustic wave equation with a kinetic boundary conditions as an example fitting into the abstract setting. We present a finite element space discretization and the full discretization error bound. We finish the paper with numerical experiments underlining the theoretical error bounds and the efficiency of the IMEX scheme.
2 The implicit–explicit (IMEX) scheme
In this section we first introduce the problem setting and then present the IMEX scheme and its properties.
2.1 Continuous problem
Let V, H be Hilbert spaces with \(V \subset H\). We consider the following second-order variational equation as a prototype of second-order wave equations in weak formulation: find \(u{:}[0,T] \rightarrow V\) s.t.
with bilinear forms
and \(f {:}[0,T]\times V \rightarrow H\).
For the rest of the paper we require the following assumptions without further mentioning it everywhere:
Assumption 2.1
-
(a)
The bilinear form \(m\) is a scalar product on H. In the following, we equip H with \(m\) and denote the corresponding norm by \({\Vert {\cdot } \Vert _{m}}\).
-
(b)
The bilinear form \(a\) is symmetric and there exists a constant \(\alpha _\text {G}> 0\) s.t.
$$\begin{aligned} {\widetilde{a}}= a+ \alpha _\text {G}m\end{aligned}$$is a scalar product on V. In the following we equip V with \({\widetilde{a}}\) and denote the corresponding norm by \({\Vert {\cdot } \Vert _{{\widetilde{a}}}}\).
-
(c)
The space V is densely embedded in H, i.e., there exists an embedding constant \(C_{\text {emb}}\), s.t.
$$\begin{aligned} {\Vert {v} \Vert _{m}} \le C_{\text {emb}}{\Vert {v} \Vert _{{\widetilde{a}}}} \qquad \text {for all }v \in V. \end{aligned}$$(2) -
(d)
The bilinear form \(b\) is bounded and quasi-monotone, i.e., there exist constants \(C_{B}> 0\) and \(\beta _{\text {qm}}\) s.t.
$$\begin{aligned} b\big ( v, w \big )&\le C_{B}{\Vert {v} \Vert _{{\widetilde{a}}}}{\Vert {w} \Vert _{m}}&\text {for all } v \in V,w \in H, \end{aligned}$$(3a)$$\begin{aligned} b\big ( v, v \big ) + \beta _{\text {qm}}m\big ( v, v \big )&\ge 0&\text {for all } v \in V. \end{aligned}$$(3b) -
(e)
The nonlinearity \(f\) satisfies \(f\in C^1([0,T]\times V;H)\) and is locally Lipschitz-continuous on V with constant \(L_{T,\rho }\), i.e., for all \(t\le T\) and \(v,w \in V\) with \({\Vert {v} \Vert _{{\widetilde{a}}}},{\Vert {w} \Vert _{{\widetilde{a}}}}\le \rho \) we have
$$\begin{aligned} {\Vert {f(t,v) - f(t,w)} \Vert _{m}} \le L_{T,\rho } {\Vert {v-w} \Vert _{{\widetilde{a}}}}. \end{aligned}$$(4)
Example 2.2
We consider the semilinear damped wave equation
in a domain \(\varOmega \subset \mathbb {R}^d, d = 2,3\).
-
(a)
It is well known, that the weak formulation of (5) equipped with Dirichlet boundary condition \(u = 0\) on \(\partial \varOmega \) fits in the setting presented above with
$$\begin{aligned} V = H_{0}^1(\varOmega ),\qquad H = L^2(\varOmega ). \end{aligned}$$ -
(b)
If \(\varOmega \) has a \(C^2\) boundary, we can equip (5) with semilinear acoustic boundary conditions that have the form
$$\begin{aligned} u_{tt} + \partial _\mathbf{n}u - \Delta _{\partial \varOmega } u = |u|^2u\qquad \text {in } (0,T) \times \partial \varOmega . \end{aligned}$$We discuss a more general form of this example in Sect. 4, and show that it also fits in the setting of Sect. 2.1.
In the following we consider (1) as evolution equation on H:
with operators \(A {:}D(A)\rightarrow H\) and \(B{:}V\rightarrow H\) induced by \(a\) and \(b\), i.e.,
with
The wellposedness of (6) can be obtained by classical semigroup theory applied to the first-order formulation (7) of the equation, as detailed in [13].
Theorem 2.3
The problem (6) is locally wellposed, i.e., for all \(u^0\in D(A), v^0\in V\) there exists a time \(t^*= t^*(u^0,v^0)\) s.t. for all \(T<t^*\), (6) has a unique solution
2.2 Construction of the IMEX scheme
We modify the Crank–Nicolson scheme such that it treats the nonlinear part of (6) explicitly and thus avoids the solution of nonlinear equations.
To derive and analyze the scheme we state (6) in a first-order formulation. Set \(u' = v\) and
Then (6) can be written as
We consider this equation in the Hilbert space \((X,p) {:=}(V,{\widetilde{a}}) \times (H,m)\), where p is the natural inner product, and \(D(S) = D(A) \times V\).
Let \(\tau > 0\) be the time step and \(t_n{:=}\tau n\) for \(n\in \mathbb {N}\). By
we denote the numerical approximation of the exact solution of (7) at time \(t_n\) and we further use the short notations
The Crank–Nicolson scheme applied to (7) has the form
and can be written as
By [11, Lemma 2.14 and Lemma 4.2] we have the following properties of \(R_\pm \):
Lemma 2.4
Let \(c_{\text {qm}}= \frac{1}{2} \alpha _\text {G}C_{\text {emb}}+ \beta _{\text {qm}}\) with \(C_{\text {emb}}\) defined in (2) and \(\alpha _\text {G},\beta _{\text {qm}}\) from Assumption 2.1. Then, for \(\tau c_{\text {qm}}< 2\), the following assertions hold true:
-
(a)
\(R_+\) is invertible with \({\Vert {R_+^{-1}} \Vert _{X \leftarrow X}} \le 1\) and \(R_+^{-1}x \in D(S)\) for all \(x \in X\).
-
(b)
\(R{:=}R_+^{-1}R_-\) has a continuous extension on X satisfying \({\Vert {R} \Vert _{X \leftarrow X}} \le \mathrm e^{\tau c_{\text {qm}}}\).
By applying \(R_+^{-1}\) to (8b), the Crank–Nicolson scheme reads
Lemma 2.5
The Crank–Nicolson scheme (8a) can equivalently be rewritten in a half-full-half step formulation
Proof
Using
we write (8a) component wise:
The first equation gives (9b). By eliminating \(u^{n+1}\) in the second equation we obtain
which is equivalent to the two half steps (9a) and (9c). \(\square \)
By replacing in (2.5) the trapezoidal rule for the nonlinear part by the left/right rectangular rule, respectively, we obtain the following IMEX scheme:
It can be interpreted as a combination of the Crank–Nicolson scheme for the linear and the leapfrog scheme for the nonlinear part, respectively. Obviously, it is time reversible.
Remark 2.6
An equivalent representation of \(v^{n+1}\) is obtained by subtracting (11a) from (11c), namely
It is computationally more efficient because of the elimination of the operators A and B.
The implementation is comprised by solving the linear system in (11a), and then computing (11b), and (11d). Altogether, each time step requires the solution of one linear system, one application of A and one evaluation of the nonlinearity (note that \(f^{n+1}\) can be reused in the next time step).
2.3 Wellposedness of the IMEX scheme
By (11a) we have
with \(Q_\pm {:}D(A) \rightarrow H\) given by
Since these operators play an important role in the analysis of the method, we collect some of their properties.
Lemma 2.7
(Properties of \(Q_\pm \)) Let
Then, the operator \(Q_+\) is invertible and its inverse \({Q_+^{-1}{:}H \rightarrow D(A)}\) satisfies
Proof
(a) We first prove that \(Q_+\) is invertible. By (13) and Assumption 2.1, the bilinear form
is V-elliptic. Hence, by the Lax–Milgram lemma, for a given \(v \in H \subset V'\) there exists a unique \(z \in V\) such that
or equivalently
This yields \(z \in D(A)\) and \(Q_+z = v\), hence \(Q_+\) is invertible.
(b) Proof of the bounds (14): Let \(v \in H\) and set \(z = Q_+^{-1}v \in D(A)\). Then we have
due to the quasi-monotonicity of \(b\) (3b) and (13). This immediately yields (14a) and (14b). The bound (14c) can be shown similar to the bound for \(R_+^{-1}R_-\) in the proof of [11, Lemma 2.14]. \(\square \)
Corollary 2.8
The IMEX scheme is wellposed in \(D(A)\times H\), i.e., for \(u^0\in D(A)\) and \(v^0\in H\) the numerical approximations satisfy
Proof
We use induction. The statement holds for \(n=0\) by assumption, hence we assume that \(u^n\in D(A), v^n\in H\) for some \(n\ge 0\). By Lemma 2.7, \(Q_+\) is invertible and (12) implies \(v^{n+\frac{1}{2}}\in D(A)\). From (11b) and (11c) we then get \(u^{n+1}\in D(A), v^{n+1}\in H\). \(\square \)
2.4 Error bound for the IMEX scheme
We now state a second-order error bound for the IMEX scheme.
Theorem 2.9
Let the assumptions of Theorem 2.3 be satisfied, \(T<t^*(u^0,v^0)\), and let the exact solution u of (6) satisfy \(u \in C^4([0,T],H) \cap C^3([0,T],V)\) and \(f(u) \in C^2([0,T],H)\). Then for all \(\tau \) sufficiently small and all \(t_n< T\), the approximations \(u^n\) from the IMEX scheme (11) are bounded by
Moreover, \(u^n,v^n\) satisfy the error bound
with \(M= c_{\text {qm}}+ \frac{ \big (1 + ( 3/2 )^{1/2}\big )L_{T,\rho }}{1- \big (1 + ( 3/2 )^{1/2}\big )L_{T,\rho }\tau }\), \(c_{\text {qm}}= \frac{1}{2}\alpha _\text {G}C_{\text {emb}}+ \beta _{\text {qm}}\),
and a constant C that only depends on T but is independent of \(\tau \), L, and u.
Since the proof works with the same arguments as the the more complicated proof of Theorem 3.3 for the full discretization error of the IMEX scheme, we do not present it here, cf., also Remark 3.4 a).
2.5 The IMEX scheme in first-order formulation
For the error analysis we rewrite the IMEX scheme (11) as a perturbation of the one-step formulation of the CN scheme (8b). A similar idea was used in [14] for analyzing the leapfrog scheme and locally implicit schemes for Maxwell equations.
The formulation (8c) of the Crank–Nicolson scheme can be used to prove stability. For the IMEX scheme we now derive a similar formulation.
Lemma 2.10
-
(a)
We have the following representations of \(R_+^{-1}\) and R:
$$\begin{aligned} R_+^{-1}&= \begin{bmatrix} Q_+^{-1}\left( {\text {I}}+ \frac{\tau }{2}B\right) &{} \frac{\tau }{2}Q_+^{-1}\\ -\frac{2}{\tau } + \frac{2}{\tau } Q_+^{-1}({\text {I}}+ \frac{\tau }{2}B)\, &{} Q_+^{-1} \end{bmatrix}, \end{aligned}$$(15a)$$\begin{aligned} R&= \begin{bmatrix} -{\text {I}}+ Q_+^{-1}\left( 2{\text {I}}+ \tau B\right) &{} \tau Q_+^{-1}\\ -\tfrac{4}{\tau }{\text {I}}+ \tfrac{1}{\tau } Q_+^{-1}\left( 4{\text {I}}+ 2\tau B\right) \; &{} Q_-Q_+^{-1} \end{bmatrix}. \end{aligned}$$(15b) -
(b)
For all \(w \in V\) we have
$$\begin{aligned} R_+^{-1}\begin{bmatrix} w \\ -Bw \end{bmatrix} = \begin{bmatrix} Q_+^{-1}w \\ - \left( B + \frac{\tau }{2}A \right) Q_+^{-1}w \end{bmatrix}. \end{aligned}$$(15c) -
(c)
The IMEX scheme (11) is equivalent to
$$\begin{aligned} \begin{aligned} x^{n+1}= Rx^n&+ \frac{\tau }{2}R_+^{-1}\big (g^n+ g^{n+1}\big ) + \frac{\tau ^2}{4}\begin{bmatrix} Q_+^{-1}\big ( f^n- f^{n+1}\big ) \\ - \left( B + \frac{\tau }{2}A\right) Q_+^{-1}\big ( f^n- f^{n+1}\big ) \end{bmatrix}. \end{aligned} \end{aligned}$$(16)
Proof
First note that the right-hand side of (15a) is a well-defined mapping from X to D(S) by Lemma 2.7. The identities (15) can be verified by straightforward calculations.
(c) To make the following calculations well defined, we assume for the moment that \(v^n, v^{n+1}\in V\). We eliminate \(v^{n+\frac{1}{2}}\) from the scheme (11), by subtracting (11c) from (11a). This yields
which differs from the Crank–Nicolson scheme by the contributions of the nonlinearity \(f\), cf. (10). Note that we have \(f^n- f^{n+1}\in V\) since, by Corollary 2.8, \(v^n, v^{n+1}\in V\) and \(v^{n+\frac{1}{2}}\in D(A)\) . Inserting (17) into (11b) gives
On the other hand, by adding (11a) and (11c) and inserting (11b) and (17) we get
With the definition (8b) of \(R_\pm \) we can express (18) and (19) in first-order formulation as
Multiplying by \(R_+^{-1}\) and using (15c) shows that the IMEX scheme is equivalent to (16) under the additional assumption \(v^n, v^{n+1}\in V\). Since both formulations are also well defined for \(v^n, v^{n+1}\in H\) and since V is dense in H, we also get their equivalence for \(v^n, v^{n+1}\in H\). \(\square \)
3 Full discretization
In this section we combine the IMEX scheme with an abstract space discretization to obtain a fully discrete scheme. We use the framework introduced in [12] for linear equations and extended in [13] to the semilinear case. It is rather general and allows one to cover conforming as well es nonconforming space discretizations, the latter being relevant for the discretization of equations with dynamic boundary conditions.
3.1 Framework
Let \((V_h)_h\) be a family of finite dimensional vector spaces for the spatial approximation related to a discretization parameter h, e.g., the mesh width of a spatial grid.
For all \(V_h\in (V_h)_h\) we consider the following discretization of (1): find \(u_h{:}[0,T] \rightarrow V_h\) s.t.
Here, \(m_{h}, a_h, b_h, f_h, u_h^0\) and \(v_h^0\) are approximations of their corresponding continuous counterparts and satisfy similar properties as in Assumption 2.1.
Assumption 3.1
In the following statements, all constants are independent of h.
-
(a)
The bilinear form \(m_{h}\) is a scalar product on \(V_h\) and we denote \(V_h\) equipped with this scalar product by \(H_h\).
-
(b)
The bilinear form \(a_h\) is symmetric and there exists a constant \(\widehat{\alpha }_\text {G}\ge 0\) s.t.
$$\begin{aligned} {\widetilde{a}}_h= a_h+ \widehat{\alpha }_\text {G}m_{h}\end{aligned}$$is a scalar product on \(V_h\). In the following we equip \(V_h\) with \({\widetilde{a}}_h\).
-
(c)
There exists a constant \(\widehat{C}_{\text {emb}}>0\) s.t. \({\Vert {v_h} \Vert _{m_{h}}}\le \widehat{C}_{\text {emb}}{\Vert {v_h} \Vert _{{\widetilde{a}}_h}}\) for all \(v_h\in V_h\).
-
(d)
The bilinear form \(b_h{:}V_h\times H_h\rightarrow \mathbb {R}\) is continuous and there exists a \(\widehat{\beta }_{\text {qm}}\ge 0\) s.t.
$$\begin{aligned} b_h\big ( v, v \big ) + \widehat{\beta }_{\text {qm}}{\Vert {v} \Vert _{m_{h}}}^2 \ge 0 \quad \text {for all } v_h\in V_h. \end{aligned}$$ -
(e)
The discrete nonlinearity \(f_h{:}[0,T] \times V_h\rightarrow H_h\) is locally Lipschitz continuous on \(V_h\) with constant \(\widehat{L}_M\), analogously to (4).
To reformulate (20) as an evolution equation on \(V_h\), we define \(A_h, B_h\in \mathcal {L}(V_h;V_h)\) via
Then (20) is equivalent to
Analogously to the continuous case we can rewrite this in a first-order formulation. With the Hilbert space \(X_h= V_h\times H_h\) and
(21) is equivalent to
The IMEX scheme (11) applied to (21), reads
where we used the short notation \(f_h^n{:=}f_h(t_n,u_h^n)\). As in the continuous case, (23c) can be replaced by the more efficient update
3.2 Error analysis
We first introduce some notation that is required for the unified error analysis presented in [12, 13].
To relate the discrete and the continuous solution we assume that there exists a lift operator \({\mathcal {L}_h^V\in \mathcal {L}(V_h;V)}\) which satisfies
for all \(v_h\in V_h\) with constants \(C_H, C_V> 0\) which are independent of h.
The adjoints
of \(\mathcal {L}_h^V\) play an important role in the error analysis. They are defined via
The corresponding first-order operators \(\mathcal {L}_h{:}X_h\rightarrow X\) and \({\mathcal {L}}_h^*{:}X \rightarrow X_h\) are defined as
Let \(Z^V\overset{d}{\hookrightarrow }V\) be a dense subspace and \({I_h\in \mathcal {L}(Z^V;V_h)}\) be an interpolation operator satisfying
with \(\widehat{C}_{I}> 0\) independent of h. We define \(Z = V \times Z^V\overset{d}{\hookrightarrow }X\) and the first-order reference operator \({J_h{:}Z \rightarrow X_h}\) by
Furthermore, for \(v_h, w_h\in V_h\), the errors in the scalar products are defined via.
For \(z = \begin{bmatrix} v \\ w \end{bmatrix} \in Z\), the discretization errors in the linear operator and the nonlinearity are given by the following remainder terms:
Analogously to the continuous case we define the following operators:
Since the setting is the same as in the continuous case with constants independent of h, Lemmas 2.4 and 2.7 transfer directly to the discrete case with the continuous constants replaced by the discrete ones.
Our analysis relies on the following regularity assumptions.
Assumption 3.2
Let the assumptions of Theorem 2.3 be satisfied and let u be the classical solution of (6) that satisfies additionally
for a \(T < t^*(u^0, v^0)\).
We now state our main results. First we present an abstract error bound depending on the radius of a ball containing the exact and the numerical solution. Here, this radius depends on \(\tau \) and h. Afterwards, in Corollary 3.5, we show that under additional consistency assumptions for the space discretization and for sufficiently small \(\tau \) and h, the fully discrete approximations \(u_h^n\) are bounded in terms of the exact solution only.
Theorem 3.3
Let Assumptions 3.1 and 3.2 be satisfied. For T given in Assumption 3.2 we set
and define \(\widehat{c}_{\text {qm}}= \widehat{\alpha }_\text {G}\widehat{C}_{\text {emb}}/2 + \widehat{\beta }_{\text {qm}}\). If \(\tau \) satisfies the step-size restriction
then for all \(n > 0\) with \(t_n< T\), the fully discrete approximations \(u_h^n, v_h^n\) given by the scheme (23) satisfy the error bound
with
and a constant C that only depends on T, but which is independent of \(\tau \), h, \(\widehat{L}\), and u. The constants \(E_{h,i} = E_{h,i}(u)\) contain the abstract space discretization errors and are given by
and \(E= E(u)\) is given in Theorem 2.9.
Proof
All error terms arising from the space discretization can be expressed within the unified error analysis and were bounded against \(E_{h,i}, i=0,\ldots ,4\) in [12, Theorem 4.8], and [13, Theorem 3.9], respectively.
For the proof of the error bound (26), we use the first-order formulation of the IMEX scheme. We use the notation
The proof consists of four main steps.
(a) Splitting of the error. The error can be split via
Due to the continuity of the lift operator and [12, Theorem 4.8] we have
Hence, it remains to bound the discrete error \({\Vert {e_h^n} \Vert _{X_h}}\).
(b) Derivation of an error recursion for \(e_h^n\). Since the discrete operators share the properties of their continuous counterparts, we can rewrite the fully discrete scheme (23) analogously to Lemma 2.10 as
To derive a recursion for the discrete error, we insert \(J_hx\) into the fully discrete scheme (28) and obtain
with a defect \(\Delta _h^{n+1}\) which is yet to be determined. We can interpret \(\Delta _h^{n+1}\) as a perturbation of the defect \(\Delta _{h,{\text {CN}}}^{n+1}\) of the fully discrete Crank–Nicolson scheme given by
In fact we have
A simple calculation shows that \(J_hx\) inserted in the Crank–Nicolson scheme satisfies
where
contains the abstract space discretization errors and
is the defect of the Crank–Nicolson scheme applied to the continuous equation (7). By applying \(\widehat{R}_+^{-1}\) to (32) we obtain with (30) and (31)
Subtracting (29) from (28) yields the error recursion
(c) Stability. Solving (35) gives
Taking the norm, using the triangle inequality, and \({\Vert {\widehat{R}} \Vert }{_{X_h\leftarrow X_h}} \le e^{\tau \widehat{c}_{\text {qm}}}\) yields
We investigate the different terms in the first sum separately. For this we use the Lipschitz-continuity of the discrete nonlinearity and the bounds from Lemmas 2.7 and 2.4 for the discrete case. Note that by (38) we have \({\Vert {\mathcal {L}_h^{V* }u(t)} \Vert _{{\widetilde{a}}_h}}\le \rho \) for all \(t\le T\). With \(\tau \widehat{c}_{\text {qm}}< 2\) and \({\Vert {\widehat{R}_+^{-1}} \Vert }{_{X_h\leftarrow X_h}} \le 1\), we obtain:
Using
we have for \(\frac{\tau ^2}{2}\widehat{\alpha }_\text {G}+ \tau \widehat{\beta }_{\text {qm}}\le 1\)
With \(C_{3/2}= 1 + ( 3/2 )^{1/2}\), this yields
By applying a discrete Gronwall Lemma, multiplying by \(\mathrm e^{n\tau \widehat{c}_{\text {qm}}}\), and inserting (34), we obtain for
(d) Defects. The initial error \(e_h^0\) is given by the discretization errors of the initial values and bounded by
From (33) we obtain for the defect containing the space discretization errors
By [12, Theorem 4.8] and [13, Theorem 3.9] we have
The Crank–Nicolson defect was bounded in [11, Lemma 2.15] by
To bound the additional IMEX defect we split it into
with
where we used the additional notation \(\widetilde{f}^m= f(t_m,\widetilde{u}^m)\). The first term is bounded by \({\Vert {\widetilde{\delta }_{h,1}^{m}} \Vert _{X_h}} \le C \tau E_{h,1}\).
The terms \(\widetilde{\delta }_{h,2}^{m}\) and \(\widetilde{\delta }_{h,3}^{m}\) are only of order \(\tau ^2\), which is not sufficient to obtain a global error of order two. Moreover, a combination of both terms from two successive time steps allows to gain an additional factor of \(\tau \). With the explicit representation of \(\widehat{R}\) analogous to that of R in (15b), we obtain
Using this together with the bound bounds for \(\widehat{Q}_+^{-1}\) and \(\widehat{Q}_-\widehat{Q}_+^{-1}\) from Lemma 2.7 and the continuity of the adjoint lift operator, leads to the bound
and hence,
Inserting all bounds into (36) yields
Finally, the error bound (26) follows from
and (27). \(\square \)
Remark 3.4
-
(a)
Theorem 2.9 can be proven by replacing all discrete quantities by their continuous counterparts in the proof of Theorem 3.3. Since all assumptions in the discrete setting remain valid in the continuous case, the proof applies to the latter as well and all space discretization errors vanish.
-
(b)
The step-size restriction in Theorem 2.9 is not a CFL condition, since it only depends on constants that are independent of the mesh width h. Note that in the special case of a linear problem with \(\widehat{\alpha }_\text {G}= \widehat{\beta }_{\text {qm}}= 0\), the scheme is unconditionally stable.
-
(c)
Theorem 3.3 and Corollary 3.5 are also valid for the Crank–Nicolson scheme. In this case, the error recursion (35) simplifies to
$$\begin{aligned} e_h^{n+1}= \widehat{R}e_h^n+ \frac{\tau }{2}\widehat{R}_+^{-1}\left( g_h^{n+1}- \widetilde{g}_h^{n+1}+ g_h^n- \widetilde{g}_h^n\right) + \delta _h^{n+1} + \widehat{R}_+^{-1}{\mathcal {L}}_h^*\delta _{{\text {CN}}}^{n+1}\end{aligned}$$and the error bound (26) holds with \(E_\text {CN}(u)\! =\! {\Vert {u^{(4)}} \Vert _{L^\infty \left( [0,T];H \right) }} \!{+}\! {\Vert {u^{(3)}} \Vert _{L^\infty \left( [0,T];V \right) }}\) (instead of \(E\)) and \({1 + ( 3/2 )^{1/2}}\) replaced by 1 in the CFL condition and the error bound.
Under additional consistency assumptions for the space discretization, the following corollary states that for sufficiently small \(\tau \) and h, the fully discrete approximations are bounded in terms of the exact solution and converge.
Corollary 3.5
Let Assumptions 3.1 and 3.2 be satisfied.
- a):
-
For T given in Assumption 3.2 define
$$\begin{aligned} \rho = 2C_V{\Vert {u} \Vert _{L^\infty \left( [0,T];V \right) }}. \end{aligned}$$(38)If \(E_{h,i}\overset{h \rightarrow 0}{\longrightarrow }0\) for \(i = 0,\ldots , 4\), then there exist \(h^*>0\) and \(\tau ^*> 0\) with \(\tau ^*\) independent of h s.t. for all \(h< h^*, \tau < \tau ^*\) we have
$$\begin{aligned} \max _{t_n\le T} {\Vert {u_h^n} \Vert _{{\widetilde{a}}_h}} \le \rho , \end{aligned}$$(39)and the fully discrete solution converges, i.e.,
$$\begin{aligned} \max _{t_n\le T}\lbrace {\Vert {\mathcal {L}_h^Vu_h^n- u(t_n)} \Vert _{{\widetilde{a}}}} + {\Vert {\mathcal {L}_h^Vv_h^n- u'(t_n)} \Vert _{m}}\rbrace \rightarrow 0, \qquad \tau , h \rightarrow 0. \end{aligned}$$ - b):
-
If additionally \( E_{h,i}\le h^k\) for a \(k>0\) and \(i = 0,\ldots , 4\), we obtain the error bound
$$\begin{aligned} \max _{t_n\le T}\lbrace {\Vert {\mathcal {L}_h^Vu_h^n- u(t_n)} \Vert _{{\widetilde{a}}}} + {\Vert {\mathcal {L}_h^Vv_h^n- u'(t_n)} \Vert _{m}}\rbrace \le C\mathrm e^{\widehat{M}t_n}(\tau ^2 +h^k) \end{aligned}$$(40)with \(\widehat{M}\) defined in Theorem 3.3 and a constant C independent of \(\tau \) and h.
Proof
It remains to prove the bound (39) for \(\tau \) and h sufficiently small, since all other assertions then follow immediately from Theorem 3.3.
Let \(f_h^\rho {:}[0,T]\times V_h\rightarrow H_h\) be a function that is globally Lipschitz-continuous on \(V_h\) with Lipschitz constant \(\widehat{L}_{T,\rho }\) and satisfies
Further let \(u_h^{\rho ,n}\) be the iterates of the IMEX scheme (23) with \(f_h\) replaced by \(f_h^\rho \). Note that due to (38) we have
Hence, as in the proof of Theorem 3.3, we obtain similar to the bound of the first component in (37)
for all \(t_n\le T\). Since \(E_{h,i}\overset{h \rightarrow 0}{\longrightarrow }0\), and E is independent of h and \(\tau \), we can choose \(h^*,\tau ^*>0\) s.t. for all \(h<h^*, \tau <\tau ^*\) we have
and hence together with (38)
This implies, that for all \(t_n\le T\) the iterates \(u_h^{\rho ,n}\) coincide with the original iterates \(u_h^n\) and thereby \({\Vert {u_h^n} \Vert _{{\widetilde{a}}_h}} = {\Vert {u_h^{\rho ,n}} \Vert _{{\widetilde{a}}_h}} \le \rho \). \(\square \)
4 Application: semilinear wave equation with kinetic boundary conditions
In this section we consider the IMEX scheme applied to a finite element discretization of a semilinear acoustic wave equation with a kinetic boundary condition. Kinetic boundary conditions serve as an effective model for the interaction of waves with a boundary covered by a thin layer. A derivation can be found in, e.g., [10], and the wellposedness was proven in [21]. The space discretization we present in this section was analyzed in [11, 12, 19].
We show that this example fits into the abstract theory presented in the previous sections.
4.1 Formulation of the equations
Let \(\varOmega \subset \mathbb {R}^d, d = 2,3\), be a bounded domain with smooth boundary \(\Gamma = \partial \varOmega \). With \(\Delta _\Gamma \) we denote the Laplace-Beltrami operator on \(\Gamma \) and with \(\mathbf{n}\) the outer normal vector.
We consider the following semilinear acoustic wave equations with kinetic boundary conditions. Seek \(u {:}[0,T] \times \overline{\varOmega } \rightarrow \mathbb {R}\) satisfying
where the nonlinearities and the coefficients satisfy the following conditions.
Assumption 4.1
-
(a)
The nonlinearities satisfy
$$\begin{aligned} f_{\varOmega }\in C^1 ([0,T]\times \overline{\varOmega } \times \mathbb {R}; \mathbb {R}), \qquad f_{\Gamma }\in C^1 ( [0,T]\times \Gamma \times \mathbb {R}; \mathbb {R}). \end{aligned}$$(42)Moreover, they satisfy the following growth condition, that there exist
$$\begin{aligned} \zeta _\varOmega {\left\{ \begin{array}{ll}< \infty , &{} d=2, \\ \le \frac{d}{d-2} , &{} d\ge 3, \end{array}\right. } \qquad \text {and} \qquad \zeta _\Gamma {\left\{ \begin{array}{ll} < \infty , &{} d=2,3, \\ \le \frac{d-1}{d-3} , &{} d\ge 4 \end{array}\right. } \end{aligned}$$(43)such that for all \((t,\mathbf {x},u)\in [0,T]\times \varOmega \times \mathbb {R}\)
$$\begin{aligned} \begin{aligned} |f_{\varOmega }(t,\mathbf {x},u)|&\le C (1 + |u|^{\zeta _\varOmega }),\qquad |\nabla f_{\varOmega }(t,\mathbf {x},u)|\le C (1 + |u|^{\zeta _\varOmega - 1}), \end{aligned} \end{aligned}$$(44)and for all \((t,\mathbf {x},u)\in [0,T]\times \Gamma \times \mathbb {R}\)
$$\begin{aligned} \begin{aligned} |f_{\Gamma }(t,\mathbf {x},u)|&\le C (1 + |u|^{\zeta _\Gamma }),\qquad |\nabla f_{\Gamma }(t,\mathbf {x},u)|\le C (1 + |u|^{\zeta _\Gamma - 1}) \end{aligned} \end{aligned}$$hold true.
-
(b)
The coefficients \(\alpha _\varOmega \in C(\overline{\varOmega }), \beta _\varOmega \in C^1(\overline{\varOmega })^d, \alpha _\Gamma \in C(\Gamma )\) and \(\beta _\Gamma \in C^1(\Gamma )^d\) are non-negative and satisfy
$$\begin{aligned} \alpha _\varOmega - \frac{1}{2}{\text {div}}\beta _\varOmega \ge 0 \quad \text { in }\varOmega , \qquad \alpha _\Gamma + \frac{1}{2}\left( \beta _\varOmega \cdot n - {\text {div}}_\Gamma \beta _\Gamma \right) \ge 0 \quad \text { on }\Gamma . \end{aligned}$$
In [13] was shown that the weak formulation of (41) is of the form (2) with
and \(f{:}[0,T]\times V\rightarrow H\) defined via
Note that for the subset relation \(V\subset H\), each element \(v\in V\) is identified with \((v,\gamma (v))\in H\). Furthermore, Assumption 2.1 is satisfied and we have \(D(A) = H^2(\varOmega ;\Gamma )\). Thus, Theorem 2.3 yields the existence of a solution u of (41) .
4.2 Space discretization
As in [13], we use the bulk-surface finite element method presented in [7] to discretize (41) in space. This discretization was also considered in [11, 12] for linear problems.
We start by giving a short summary of this method and refer to [7, 11] for more details.
Bulk-surface finite element method
Let \((\mathcal T_h)_h\) be a quasi-uniform family of consistent meshes of order p isoparametric elements with maximal mesh width h. For each \(\mathcal T_h\in (\mathcal T_h)_h\) the discretized domain and its boundary are denoted by
We define the bulk and the surface finite element space of order \(p \ge 1\) via
Here \(\mathbb {P}_p(\widehat{K})\) denotes the space of polynomials of degree p on a reference triangle \(\widehat{K}\) and \(F_K\) is a transformation from \(\widehat{K}\) to \(K\). This discretization is nonconforming because \(\varOmega _h\ne \varOmega \). In [7], an elementwise smooth homeomorphism \(G_h{:}\varOmega _h\rightarrow \varOmega \) with
is constructed. This allows us to define lifted versions of \(v_h\in V_{h,p}^\varOmega \) and \(\vartheta _h \in V_{h,p}^\Gamma \) as
The mapping \(G_h\) is constructed in such a way, that \(G_h(a_i)=a_i\), \(i=1,\ldots ,N=\dim V_h\), where \(a_1,\ldots ,a_N\in \varOmega _h\) are the nodes corresponding to the finite element discretization. This implies \(v_h^\ell (a_i)=v_h(a_i)\) for \(i=1,\ldots ,N\) and for all \(v_h\in V_{h,p}^\varOmega \).
By \(I_{h,\varOmega }{:}C(\overline{\varOmega })\rightarrow V_{h,p}^\varOmega \) and \(I_{h,\Gamma }{:}C(\Gamma ) \rightarrow V_{h,p}^\Gamma \) we denote the nodal interpolation operator in \(\varOmega \) and on \(\Gamma \), respectively. By construction, the nodes on the surface coincide with the bulk nodes and therefore we have
The semidiscrete equation
We now present the space discretization in the framework of Sect. 3. As finite element space we choose \(V_h= V_{h,p}^\varOmega \). Furthermore we set \(I_h{:=}I_{h,\varOmega }\big |_{Z^V}{:}Z^V\rightarrow V_h\), where
and define the lift operator via
with \(v^\ell \) given in (46). The spatial discretization of (41) can then be written as (20) where the discretizations \({m_{h},b_h,a_h{:}V_h\times V_h\rightarrow \mathbb {R}}\) of \(m,b,\) and \(a\) are defined via
The discretized nonlinearity \(f_h{:}[0,T]\times V_h\rightarrow H_h\) is given by
for all \(\varphi _h \in V_h\).
In [11, 13] it was shown, that Assumption 3.1 is satisfied.
Remark 4.2
The nodal interpolation only requires function evaluations in the nodes
\(a_1,\ldots ,a_N\). Since they are invariant under the lift operator, the computation of \(v_h^\ell \) is not necessary. It is only needed for the definition of \(f_h\) since the interpolation operator acts on functions defined on \(\varOmega \).
4.3 Full discretization error bound
We now state an error bound for the full discretization of (41) with the bulk-surface finite element method and the IMEX scheme (23).
Corollary 4.3
Let \(1\le p \le k\) and \(\Gamma \in C^{k+1}\). Furthermore let Assumption 4.1 be satisfied and let u be a solution of (41) on [0, T] satisfying
Then there exist \(\tau ^*, h^*, \rho > 0\) s.t. for all \(0< h< h^*, 0< \tau < \tau ^*\), and \(t_n\le T\), the approximations \(u_h^n\) and \(v_h^n\) given by (23) with bulk-surface elements of order p, satisfy
with a constant C independent of \(\tau \) and h. The Lipschitz constant of the discretized nonlinearity is given by
where \(\sigma (\varOmega )\) and \(\sigma (\Gamma )\) denote the measure of \(\varOmega \) and \(\Gamma \), respectively, and \(\zeta _\Gamma \) and \(\zeta _\varOmega \) are given in Assumption 4.1.
Proof
In [13] it was shown, that Assumptions 2.1 and 3.1 are satisfied with \(\widehat{C}_{\text {emb}}= \widehat{\alpha }_\text {G}= 1\), \(\widehat{\beta }_{\text {qm}}=0\), and \(\widehat{L}_{T,\rho }\) given in (50). The regularity assumptions on u are such that \(u \in C^4\left( [0,T],H\right) \cap C^3\left( [0,T],V\right) \cap C^2\left( [0,T],Z^V\right) \), cf. (45) and (47). Since additionally \(Z^V= D(A) = H^2(\varOmega ;\Gamma )\), we have that \(Au\in C^2([0,T];H)\) und hence \(f(u) = u'' + Bu' + Au\in C^2([0,T];H)\).
Thus, also Assumption 3.2 is satisfied, and we can apply Corollary 3.5.
Under the above assumptions, in [12, 13] it was shown that the space discretization error terms are bounded by \(E_{h,i}\le Ch^p\). So the bound (49) follows then directly by (40). \(\square \)
4.4 Implementation
In the following numerical experiments we compare the IMEX scheme with the Crank–Nicolson and the explicit classical Runge–Kutta scheme of order 4 applied to the space discretized wave equation with kinetic boundary conditions. For the implementation we used the C++ finite element-library deal.II (version 9.2) [2, 5]. The codes to reproduce the experiments are available on https://doi.org/10.5445/IR/1000127003. We ran the experiments on a computer with an i5 processor (3.5 GHz) and 16 GB RAM.
To comment on the implementation we first introduce some additional notation. For a finite element function \(u_h\in V_h\) we denote by \(\mathbf {u}\in \mathbb {R}^{N}\), the corresponding coefficient vector in the finite element basis. Furthermore, \(\mathbf {M}\in \mathbb {R}^{N\times N}\) is the mass matrix, \(\mathbf {A}, \mathbf {B}\in \mathbb {R}^{N\times N}\) are the stiffness matrices related to \(A_h\) and \(B_h\), respectively, and \(\mathbf {f}^n\in \mathbb {R}^N\) denotes the load vector corresponding to \(f_h^n= f_h(t_n,u^n)\), \(n\in \mathbb {N}\).
IMEX scheme
The fully discrete IMEX scheme (23) reads
The linear system in (51a) has the form
We solve this system with the GMRES solver provided by deal.II and either a sparse incomplete LU or a geometric multigrid preconditioner. For the measure of the error in the GMRES iterations, the residual r with corresponding coefficient vector \(\mathbf{r}\) is used. A suitable stopping criteria would be
where \(\text {tol}\) is a given tolerance, since then in (51b) the error in \(\mathbf {u}^{n+1}\) caused by the solution of the linear system in the \({\Vert {\cdot } \Vert _{{\widetilde{a}}_h}}\) norm is of order \(\tau ^3\), which corresponds to the local error of the IMEX scheme. In practice, the computation of \({\Vert {r} \Vert _{{\widetilde{a}}_h}}\) is quite expensive and we thus used the stopping criterion
in a grid dependent scaled Euclidean norm \(\Vert \cdot \Vert _{h,2} = h^{d/2}\Vert \cdot \Vert _{2}\). This is much more efficient, since this norm is available within the GMRES code at no additional cost. The criterion worked well in our numerical experiments as we will show in Sect. 4.5. We always use \(\text {tol}= 0.01\) in our numerical examples, which is chosen s.t. the errors in solving the linear systems do not destroy the overall order of convergence.
Note that in the IMEX scheme (23), only \(\mathbf {M}\mathbf {v}^{n+1}\) is required so that we neither compute nor store \(\mathbf {v}^{n+1}\) itself.
Crank–Nicolson scheme
The fully discrete Crank–Nicolson scheme can be written in the form
We solve the nonlinear equation (52a) with a simplified Newton method where we use \(\mathbf {Q}_+\) as an approximation to the Jacobian. The linear equations are solved as in the IMEX scheme. We stop the Newton scheme when the update \(\Delta u\) satisfies \(\Vert \Delta u\Vert _{h,2} \le \tau ^3 \widetilde{\text {tol}}\) with a given tolerance \(\widetilde{\text {tol}}\). In the numerical examples we use \(\widetilde{\text {tol}}= 0.1\), which is again chosen s.t. the Newton errors do not destroy the overall order of convergence. All matrix vector products appearing in (52a) and (52b) are computed only ones and saved in temporary vectors, as well as all terms that can be reused in the next time step. As in the IMEX scheme we do not compute \(\mathbf {v}^{n+1}\) but only \(\mathbf {M}\mathbf {v}^{n+1}\).
Classical Runge–Kutta scheme
The classical Runge–Kutta scheme is an explicit scheme of order four that is suited for hyperbolic problems because its stability region contains an interval on the imaginary axis. This is in contrast to the second-order schemes by Heun and Runge, which intersect with the imaginary axis in the origin only. We implemented it using mass lumping to obtain a fully explicit scheme. Note that the space discretization with mass lumping also fits into the setting of Sect. 3, as it was shown in [12] for a linear acoustic wave equation.
4.5 Numerical examples
We consider the semilinear wave equation with kinetic boundary conditions (41) with \(\alpha _\varOmega = 1, \beta _\varOmega (\mathbf {x}) = \mathbf {x}\) and \(\alpha _\Gamma = \beta _\Gamma = 0\) on the unit disc \(\varOmega = B(0,1) \subset \mathbb {R}^2\). As nonlinearities, we choose
Error \(E_h(0.8)\) of the IMEX scheme (solved with GMRES and ILU preconditioner), the Crank–Nicolson scheme, and the classical Runge–Kutta method plotted against step size \(\tau \) for coarse space discretization (328,193 degrees of freedom, left) and fine space discretization (18,3118,745 degrees of freedom, right)
Error \(E_h(0.8)\) of the IMEX scheme, solved with GMRES and ILU/Multigrid(MG, F-cycle with 8 levels) preconditioner, the Crank–Nicolson scheme, and the classical Runge–Kutta method plotted against runtime for coarse space discretization (3288,193 degrees of freedom, top) and fine space discretization (18,3118,745 degrees of freedom, bottom)
with
and as initial values
The example is chosen such that the exact solution is given by
which allows us to compute reliable errors. For the space discretization we use isoparametric elements of order \(p=2\), and choose \(u_h^0= I_{h,\varOmega }u^0\) and \(v_h^0= I_{h,\varOmega }v^0\) as discrete initial values.
As we cannot compute the lift of a finite element function exactly, we consider the error
We evaluate the integrals with a quadrature rule of order 4 such that the quadrature error is negligible. The restriction of u to \(\varOmega _h\) is possible since for convex domains we have \(\varOmega _h\subset \varOmega \).
Error \(E_h(0.8)\) of the IMEX scheme plotted against the runtime when using the two different error estimates as stopping criteria for the GMRES scheme as discussed in Sect. 4.4 for a coarse space discretization (3288,193 degrees of freedom)
In Fig. 1 the errors of the IMEX, the Crank–Nicolson, and the classical Runge–Kutta scheme are plotted against the time-step size \(\tau \) for a coarse (\(h\approx 0.014\)) and a fine (\(h\approx 0.007\)) space discretization, respectively. As predicted by Corollary 1, the IMEX and the Crank–Nicolson scheme converge with order two until the error of the space discretization is reached. The explicit Runge–Kutta scheme is only stable under a strong CFL condition and then the error reaches immediately the space discretization error plateau.
Figure 2 shows the errors of the different schemes plotted against the runtime for the same coarse and fine space discretization as in Fig. 1. It can be observed, that the IMEX scheme is significantly faster than the Crank–Nicolson scheme. For errors of the magnitude of the space discretization error plateau, the classical Runge–Kutta scheme is more efficient than the IMEX scheme, but the IMEX scheme outperforms the Runge–Kutta scheme if less accuracy is sufficient. The Runge–Kutta method has the disadvantage that the stability limit in applications is not exactly known, and therefore there is a risk that it will not be stable if a too large time-step size is chosen, or the effort is unnecessarily high if the time-step size is too small. For the large system obtained by the fine space discretization and large time-step sizes, it can be observed, that the use of the multigrid preconditioner is quite efficient.
Finally, Fig. 3 shows a comparison of the runtimes of the IMEX scheme when using the different stopping criteria for the GMRES solver discussed in Sect. 4.4, namely using \({\Vert {r} \Vert _{{\widetilde{a}}_h}}\) or \(\Vert r\Vert _{h,2}\) as estimate for the error, respectively. It can be seen that the afford of computing the (better suited) \({\Vert {r} \Vert _{{\widetilde{a}}_h}}\) is too high and does not pay off.
Change history
22 March 2021
A Correction to this paper has been published: https://doi.org/10.1007/s00211-021-01190-y
References
Akrivis, G., Crouzeix, M., Makridakis, C.: Implicit–explicit multistep methods for quasilinear parabolic equations. Numer. Math. 82(4), 521–541 (1999). https://doi.org/10.1007/s002110050429
Arndt, D., Bangerth, W., Blais, B., Clevenger, T.C., Fehling, M., Grayver, A.V., Heister, T., Heltai, L., Kronbichler, M., Maier, M., Munch, P., Pelteret, J.P., Rastak, R., Thomas, I., Turcksin, B., Wang, Z., Wells, D.: The deal.II library, version 9.2. J. Numer. Math. 28(3), 131–146 (2020). https://doi.org/10.1515/jnma-2020-0043
Ascher, U.M., Ruuth, S.J., Spiteri, R.J.: Implicit–explicit Runge–Kutta methods for time-dependent partial differential equations. Appl. Numer. Math. 25(2–3), 151–167 (1997). https://doi.org/10.1016/S0168-9274(97)00056-1. ((Special issue on time integration (Amsterdam, 1996)))
Ascher, U.M., Ruuth, S.J., Wetton, B.T.R.: Implicit–explicit methods for time-dependent partial differential equations. SIAM J. Numer. Anal. 32(3), 797–823 (1995). https://doi.org/10.1137/0732037
Bangerth, W., Hartmann, R., Kanschat, G.: deal.II—a general-purpose object-oriented finite element library. ACM Trans. Math. Softw. (2007). https://doi.org/10.1145/1268776.1268779
Boscarino, S.: Error analysis of IMEX Runge–Kutta methods derived from differential-algebraic systems. SIAM J. Numer. Anal. 45(4), 1600–1621 (2007). https://doi.org/10.1137/060656929
Elliott, C.M., Ranner, T.: Finite element analysis for a coupled bulk-surface partial differential equation. IMA J. Numer. Anal. 33(2), 377–402 (2013). https://doi.org/10.1093/imanum/drs022
Frank, J., Hundsdorfer, W., Verwer, J.G.: On the stability of implicit–explicit linear multistep methods. Appl. Numer. Math. 25(2–3), 193–205 (1997). https://doi.org/10.1016/S0168-9274(97)00059-7. ((Special issue on time integration (Amsterdam, 1996)))
Gardner, D.J., Guerra, J.E., Hamon, F.P., Reynolds, D.R., Ullrich, P.A., Woodward, C.S.: Implicit–explicit (IMEX) Runge–Kutta methods for non-hydrostatic atmospheric models. Geosci. Model Dev. 11(4), 1497 (2018)
Goldstein, G.R.: Derivation and physical interpretation of general boundary conditions. Adv. Differ. Equ. 11(4), 457–480 (2006)
Hipp, D.: A unified error analysis for spatial discretizations of wave-type equations with applications to dynamic boundary conditions. Ph.D. thesis, Karlsruhe Institute of Technology (2017). https://publikationen.bibliothek.kit.edu/1000070952
Hipp, D., Hochbruck, M., Stohrer, C.: Unified error analysis for nonconforming space discretizations of wave-type equations. IMA J. Numer. Anal. 39(3), 1206–1245 (2019). https://doi.org/10.1093/imanum/dry036
Hochbruck, M., Leibold, J.: Finite element discretization of semilinear acoustic wave equations with kinetic boundary conditions. Electron. Trans. Numer. Anal. 53, 522–540 (2020). https://doi.org/10.1553/etna_vol53s522
Hochbruck, M., Sturm, A.: Error analysis of a second-order locally implicit method for linear Maxwell’s equations. SIAM J. Numer. Anal. 54(5), 3167–3191 (2016). https://doi.org/10.1137/15M1038037
Hundsdorfer, W., Ruuth, S.J.: IMEX extensions of linear multistep methods with general monotonicity and boundedness properties. J. Comput. Phys. 225(2), 2016–2042 (2007). https://doi.org/10.1016/j.jcp.2007.03.003
Kadioglu, S.Y., Knoll, D.A., Lowrie, R.B., Rauenzahn, R.M.: A second order self-consistent IMEX method for radiation hydrodynamics. J. Comput. Phys. 229(22), 8313–8332 (2010). https://doi.org/10.1016/j.jcp.2010.07.019
Layton, W., Li, Y., Trenchea, C.: Recent developments in IMEX methods with time filters for systems of evolution equations. J. Comput. Appl. Math. 299, 50–67 (2016). https://doi.org/10.1016/j.cam.2015.09.038
Layton, W., Trenchea, C.: Stability of two IMEX methods, CNLF and BDF2-AB2, for uncoupling systems of evolution equations. Appl. Numer. Math. 62(2), 112–120 (2012). https://doi.org/10.1016/j.apnum.2011.10.006
Leibold, J.: Semilineare Wellengleichungen mit dynamischen Randbedingungen. Master’s thesis, Karlsruhe Institute of Technology (2017). http://na.math.kit.edu/download/thesis/2017-Leibold.pdf
Lemieux, J.F., Knoll, D.A., Losch, M., Girard, C.: A second-order accurate in time IMplicit–EXplicit (IMEX) integration scheme for sea ice dynamics. J. Comput. Phys. 263, 375–392 (2014). https://doi.org/10.1016/j.jcp.2014.01.010
Vitillaro, E.: Strong solutions for the wave equation with a kinetic boundary condition. In: Recent Trends in Nonlinear Partial Differential Equations. I. Evolution Problems, Contemp. Math., vol. 594, pp. 295–307. Am. Math. Soc. Providence, RI (2013). https://doi.org/10.1090/conm/594/11793
van Zuijlen, A.H., Bijl, H.: Implicit and explicit higher order time integration schemes for structural dynamics and fluid–structure interaction computations. Comput. Struct. 83(2–3), 93–105 (2005)
Funding
Open Access funding enabled and organized by Projekt DEAL.
Author information
Authors and Affiliations
Corresponding author
Additional information
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—Project-ID 258734477—SFB 1173.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Hochbruck, M., Leibold, J. An implicit–explicit time discretization scheme for second-order semilinear wave equations with application to dynamic boundary conditions. Numer. Math. 147, 869–899 (2021). https://doi.org/10.1007/s00211-021-01184-w
Received:
Revised:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1007/s00211-021-01184-w
Keywords
- Implicit–explicit time integration
- IMEX
- Dynamic boundary conditions
- Semilinear wave equation
- Nonconforming space discretization
- Error analysis
- A-priori error bounds
- Semilinear evolution equations
- Operator semigroups