Unconditionally Strong Stability Preserving Extensions of the TR-BDF2 Method

Abstract

We analyze monotonicity, strong stability and positivity of the TR-BDF2 method, interpreting these properties in the framework of absolute monotonicity. The radius of absolute monotonicity is computed and it is shown that the parameter value which makes the method L-stable is also the value which maximizes the radius of monotonicity. In order to achieve unconditional monotonicity, hybrid variants of TR-BDF2 are proposed, that reduce the formal order of accuracy, while keeping the native L-stability property, which is useful for the application to stiff problems. Numerical experiments compare these different hybridization strategies to other methods used in stiff and mildly stiff problems. The results show that the proposed strategies provide a good compromise between accuracy and robustness at high CFL numbers, without suffering from the limitations of alternative approaches already available in literature.

Appendices

Appendix 1: Additional Results on Monotonicity

Discrete maximum principle. The discrete maximum principle is closely related to the property of monotonicity.

Definition 7

Discrete Maximum Principle. The RK method (2) respects the discrete maximum principle if $\forall u \,{\in }\, \mathbb {R}^m$ it guarantees that

$$\begin{aligned} \min _{1 \le j \le m} u_{j}^{0} \le u_{i}^{n+1} \le \max _{1 \le j \le m} u_{j}^{0} \quad (1 \le i \le m) \end{aligned}$$

under the assumption that for $0 \,{<}\, h \,{\le }\, \tau _{0}$ and $\forall u \,{\in }\, \mathbb {R}^m$ with components $u_p$

$$\begin{aligned} \min _{1 \le q \le m} u_q \le u_p + \tau _0 f_p(t,u(t)) \le \max _{1 \le q \le m} u_q, \quad (1 \le p \le m). \end{aligned}$$

(47)

Similarly to the other properties, range boundedness may be verified under a step size restriction analogous to (5). Again following [29] and [46], we introduce two relevant sublinear functionals, denoted as max and min functional

$$\begin{aligned} \Vert u \Vert _>&= \max _j u_j \end{aligned}$$

(48a)

$$\begin{aligned} \Vert u \Vert _<&= -\min _j u_j \end{aligned}$$

(48b)

which allow us to write the assumption (47) in the form (3). By assuming the monotonicity property 1 under the functionals (48) the discrete maximum principle 7 directly follows.

Contractivity. Contractivity of numerical approximations has been extensively studied. Relevant conclusions on step size conditions for contractivity have been given in [45], while contractivity of RK for nonlinear problems was thoroughly examined in [37].

Definition 8

Contractivity. The RK method (2) is contractive if $\Vert \tilde{u}^{n} - u^{n}\Vert \,{\le }\, \Vert \tilde{u}^{n-1} - u^{n-1}\Vert$ under the assumption that

$$\begin{aligned} \Vert \tilde{u} - u + h (f(t,\tilde{u}) - f(t,u))\Vert \le \Vert \tilde{u} - u\Vert \quad \text {for}\quad 0 < h \le \tau _{0}. \end{aligned}$$

(49)

Usually, Definition 8 is verified under a step size restriction in the form of (5). For conditional contractivity, the circle condition that was originally assumed in [37] is

$$\begin{aligned} \Vert f(t,\tilde{u}) - f(t,u) + \rho (\tilde{u} - u) \Vert \le \rho \Vert \tilde{u} - u\Vert \end{aligned}$$

(50)

It was shown later in [24] that this condition can be considered as a special form of (3). By introducing the auxiliary space $\mathbb {V} \,{=}\, \mathbb {R}^m\times \mathbb {R}^m$ and considering $\Vert G_i\Vert \,{=}\,\Vert \tilde{g}^i - g^i\Vert$ and $\Vert U_i\Vert \,{=}\,\Vert \tilde{u}_i^{n} - u_i^{n}\Vert$ the circle condition (49) can be reformulated as (3) in the space of perturbations. See also [46] for additional considerations on these issues.

Appendix 2: Tables of Errors and Workload in Numerical Experiments

Here we report the Tables 8, 9, 10, 11, and 12 of the errors and workload as measured during the numerical experiments to show the typical error magnitude and workload to be expected from the different methods. We dot not claim that the error-workload tables reported here are immediately relevant for the selection of step sizes or numerical methods, but we consider them as representative of the relative workload to be expected from the methods assessed.

Table 8 Error and workload in CPU time (s) for the Brusselator problem

Table 9 Error and workload in CPU time (s) for the advection problem with non smooth initial condition

Table 10 Error and workload in CPU time (s) for the advection diffusion reaction problem with non smooth initial condition

Table 11 Error and workload in CPU time (s) for the Burgers equation with van Leer limiter and non smooth initial condition

Table 12 Error and workload in CPU time (s) for the Buckley–Leverett equation with Koren limiter and non smooth initial condition