Time-dependent PDEs arise very often in many scientific areas, such as mechanics, biology, economics, or chemistry, just to name a few. The lack of a closed form solution for general time-dependent PDEs requires one to employ numerical methods in order to find an approximation of it. These methods consider suitable discretizations of the quantity involved; in particular, they employ suitable time-stepping schemes as discretization of the time derivative. The majority of the solvers for time-dependent PDEs is based on classical linear multistep methods. Within this framework, the approximation of the solution at time \(t_n\) is evaluated as a linear combination of the \(s\) previous steps. The wide use of multistep methods is due to their simplicity; in fact, the structure of the discretized problem allows one to employ solvers for the corresponding stationary PDE as a solver for the time-dependent one. Despite this favourable quality, multistep methods have a drowback: they are (in general) not A-stable, a property that allows one to choose an arbitrary time-step for the integration. In fact, as stated by the second Dahlquist barrier, see for example [4, Theorem 6.6], an A-stable linear multistep method cannot have order of convergence greater than two. By contrast, one can devise an implicit Runge–Kutta method so that not only it is A-stable, but also such that the method possesses more desirable stability properties (e.g., L- or B-stability, see for instance [2, 3, 4]). However, the better stability properties of Runge–Kutta methods come to a price: the discretization results in a non-linear block system, to be solved for the so called stages of the method at each time step. For this reason, in recent years researchers have devoted their effort in devising efficient and robust linear solvers for the solution of block systems arising from the Runge–Kutta discretization of a time-dependent PDE, see for example [1, 5, 6, 7].

Consider the integration of an ODE of the form \(v'(t)=f(v(t),t)\) between \(0\) and a final time \(t_f>0\), given the initial condition \(v(0)=v_0\). By employing a constant time-step \(\tau \), an \(s\)-stage Runge–Kutta time-stepping scheme applied to \(v'(t)=f(v(t),t)\) reads as follows:

\[v_{n+1}= v_{n} + \tau \sum _{i=1}^{s} b_i k_{i,n}, \quad n = 0, \ldots , n_t-1,\]

where the stages \(k_{i,n}\) are given by

\begin{equation} \label {stages} k_{i,n} = f\left (v_n + \tau \sum _{j=1}^{s} a_{i,j} k_{j,n}, t_n + c_i \tau \right ), \qquad i=1,\ldots ,s, \end{equation}

with \(t_n=n \tau \). The Runge–Kutta method is uniquely defined by the coefficients \(a_{i,j}\), the weights \(b_{i}\), and the nodes \(c_{i}\), for \(i,j=1,\ldots ,s\). For this reason, an \(s\)-stage Runge–Kutta method is represented by the following Butcher tableau:

\[\def \arraystretch {1.2} \begin {array}{c|c} \mathbf {c}_{\mathrm {RK}} & A_{\mathrm {RK}}\\ \hline & \mathbf {b}_{\mathrm {RK}}^\top \end {array}\]

where \(A_{\mathrm {RK}} =\{a_{i,j} \}_{i,j=1}^s\), \(\mathbf {b}_{\mathrm {RK}} =[b_1, \ldots , b_s ]^\top \), and \(\mathbf {c}_{\mathrm {RK}} =[c_1, \ldots , c_s ]^\top \).

In this talk, we consider the following instationary incompressible Navier–Stokes equations

\[ \left \{ \begin {array}{rl} \frac {\partial \vec {v}}{\partial t} - \nu \nabla ^2 \vec {v} + \vec {v} \cdot \nabla \vec {v} + \nabla {p} = \vec {f}(\mathbf {x},t) & \quad \mathrm {in} \; \Omega \times (0,t_{f}), \\ - \nabla \cdot \vec {v} = 0 & \quad \mathrm {in} \; \Omega \times (0,t_{f}), \\ \vec {v}(\mathbf {x},t) = \vec {g}(\mathbf {x},t) & \quad \mathrm {on} \; \partial \Omega \times (0,t_{f}),\\ \vec {v}(\mathbf {x},0) = \vec {v}_0(\mathbf {x}) & \quad \mathrm {in} \; \Omega . \end {array} \right . \]

The functions \(\vec {f}\) and \(\vec {g}\) as well as the initial condition \(\vec {v}_0(\mathbf {x})\) are known. Further, \(\nu \) is the viscosity of the fluid. We integrate the problem with a Runge–Kutta scheme in time. The time discretization results in a non-linear system to be solved for the stages of the method at each time step. Specifically, introducing the variables

\[\begin {array}{l} \vec {w}^{\,v}_{n,i}= \vec {v}_n + \tau \sum _{j=1}^s a_{i,j} \vec {k}^{\, v}_{j,n}, \quad i = 1, \ldots , s,\\ w^p_{n,i}= p_n + \tau \sum _{j=1}^s a_{i,j} k^p_{j,n}, \quad i = 1, \ldots , s, \end {array}\]

at each time step the non-linear system (1) that characterizes the stages of the Runge–Kutta method reads as follows:

\begin{equation} \label {RK_NavierStokes} \left \{ \begin{array}{rl} \vec {w}^{\,v}_{n,i} - \nu \nabla ^2 \vec {w}^{\,v}_{n,i} + \vec {w}^{\,v}_{n,i} \cdot \nabla \vec {w}^{\,v}_{n,i} + \nabla w^p_{n,i} = \vec {f}(\mathbf {x},t) & \quad i = 1, \ldots , s, \\ - \nabla \cdot \vec {w}^{\,v}_{n,i} = 0 & \quad i = 1, \ldots , s. \end {array} \right . \end{equation}

Then, the solutions at time \(t_n + \tau \) are given by

\[\begin {array}{l} \vec {v}_{n+1} = \vec {v}_n + \tau \sum _{i=1}^{s} b_i \vec {k}^{\,v}_{i,n},\\ p_{n+1} = p_n + \tau \sum _{i=1}^{s} b_i k^p_{i,n}. \end {array}\]

In order to find a numerical solution, we consider a Newton linearization of the non-linear problem in (2), which is then discretized with suitable finite elements. The resulting linear system presents a saddle-point block structure, and can be very large and sparse in real-life applications. For this reason, in order to find a solution one requires the use of preconditioned iterative methods. We adopt an augmented Lagrangian-based preconditioner, and employ saddle-point theory for deriving approximations of the \((1,1)\)-block and the Schur complement. Numerical experiments show the effectiveness and robustness of our approach, for a range of problem parameters.

References

• [1]  O. Axelsson, I. Dravins, and M. Neytcheva, Stage-parallel preconditioners for implicit Runge–Kutta methods of arbitrarily high order. Linear problems, Numer. Linear Algebra Appl., 31 (2024), p. e2532.

• [2]  J. C. Butcher, Numerical methods for ordinary differential equations, John Wiley & Sons, Ltd, 3rd ed., 2016.

• [3]  E. Hairer and G. Wanner, Solving ordinary differential equations II: stiff and differential-algebraic problems, Springer Berlin, Heidelberg, 2nd ed., 1996.

• [4]  J. D. Lambert, Numerical method for ordinary differential systems: the initial value problem, John Wiley & Sons, Ltd, New York, 1991.

• [5]  K.-A. Mardal, T. K. Nilssen, and G. A. Staff, Order-optimal preconditioners for implicit Runge–Kutta schemes applied to parabolic PDEs, SIAM J. Sci. Comput., 29 (2007), pp. 361–375.

• [6]  M. M. Rana, V. E. Howle, K. Long, A. Meek, and W. Milestone, A new block preconditioner for implicit Runge–Kutta methods for parabolic PDE problems, SIAM J. Sci. Comput., 43 (2021), pp. S475–S495.

• [7]  J. Van Lent and S. Vandewalle, Multigrid methods for implicit Runge–Kutta and boundary value method discretizations of parabolic PDEs, SIAM J. Sci. Comput., 27 (2005), pp. 67–92.