Isogeometric residual minimization (iGRM) for non-stationary Stokes and Navier–Stokes problems

https://doi.org/10.1016/j.camwa.2020.11.013Get rights and content

Abstract

We show that it is possible to obtain a linear computational cost FEM-based solver for non-stationary Stokes and Navier–Stokes equations. Our method employs a technique developed by Guermond and Minev (2011), which consists of singular perturbation plus a splitting scheme. While the time-integration schemes are implicit, we use finite elements to discretize the spatial counterparts. At each time-step, we solve a PDE having weak-derivatives in one direction only (which allows for the linear computational cost), at the expense of handling strong second-order derivatives of the previous time step solution, on the right-hand side of these PDEs. This motivates the use of smooth functions such as B-splines. For high Reynolds numbers, some of these PDEs become unstable. To deal robustly with these instabilities, we propose to use a residual minimization technique. We test our method on problems having manufactured solutions, as well as on the cavity flow problem.

Introduction

In this paper, we focus on non-stationary Stokes and Navier–Stokes problems, which requires special stabilization effort, for large Reynolds numbers. There are multiple stabilization techniques developed for standard finite element methods [1], [2], [3], [4], [5]. Here, we apply a residual minimization approach for stabilization (see, e.g., [6], [7]), together with B-spline basis functions from isogeometric analysis (IGA) [8] for the discretization in space.

Minimum residual methods aim to find uhUh such that uh=argminwhUhb(wh,)()V,where U and V are Hilbert spaces, b:U×VR is a continuous bilinear (weak) form, UhU is a discrete trial space, and V (the dual space of V) is a given right-hand side. Several discretization techniques are particular incarnations of this wide-class of residual minimization methods. These include: the least-squares finite element method [9], the discontinuous Petrov–Galerkin method (DPG) with optimal test functions [10], [11], or the variational stabilization method [7]. We approach the residual minimization method using its saddle point (mixed) formulation, e.g., as described in [6].

We employ the splitting scheme from Jean-Luc Guermond and Petar Minev [12] to express the non-stationary Stokes or Navier–Stokes problem as a sequence of pressure predictor, velocity update, and pressure corrector. Each of these steps can be solved in a linear computational cost O(N) using the Kronecker product structure of the underlying matrices. Instead of using finite differences, we use finite elements for space discretizations.

The method presented in this paper is an extension of the linear computational cost Kronecker product solver that we proposed for the advection–diffusion problem [13]. The algorithm we propose uses the method of lines (discretization in space with time iterations) and delivers a linear computational cost O(N) solver for each time step. In this sense, it is an alternative to the space–time formulation [14], where an iterative solver is employed. In that configuration, the size of the mesh is equal to M=N×k (for the uniform case), where k is the number of time steps. Hence, if c is the number of iterations of the solver, its total cost becomes M×c=N×k×c, while our cost is just N×k. However, the space–time formulation allows for adaptation in both space and time, where some parts of the space–time mesh can use “larger” time steps than the others. Thus, the case of space–time adaptivity can be indeed competitive to our method. We also differ from [14] in the sense that we deliver automatic stabilization, with a linear computational cost.

Our paper is also different from [12]. Even though they use finite differences discretization with a linear computational cost direction splitting algorithm, they do not incorporate residual minimization stabilization. On the other hand, we discretize the space with B-splines over the patch of elements, and we use a weak form of the splitting scheme. We add the residual minimization procedure in a way that preserves the linear computational cost of the solver. Our method delivers high continuity approximations and automatic stabilization in every time step.

The residual minimization method used here is very similar to the DPG method developed by [15]. Our spaces, trial, and test are continuous, while they are broken in the DPG method. The motivation behind breaking the spaces is to obtain a block-diagonal structure of the matrix. Thus, static condensation practically eliminates the inner product matrix, leaving only fluxes or traces over the edges between elements. However, breaking the spaces makes the linear cost factorization impossible since it destroys the Kronecker product structure of the matrix. So does the mesh adaptation. For the DPG method, there are some modern multi-grid solvers developed, allowing for fast factorizations of the system [16].

We show that the Navier–Stokes simulation using the alternating directions implicit solver without the residual minimization stabilization generates some unexpected oscillations for high Reynolds number. In our paper, we use the splitting scheme from Jean-Luc Guermond and Petar Minev [12], which translates the Navier–Stokes problem into the velocity updates involving reaction-dominated diffusion kind of equation, followed by the explicit pressure updates. Instead of using finite differences, we use standard and higher-order finite elements for space discretizations. This kind of setup is unstable for high Reynolds number in a similar manner as to how the reaction-dominated diffusion equations are unstable for a very small diffusion coefficient. We show how the stabilization for high Reynolds number of the Navier–Stokes formulation is performed by the residual minimization method.

The novel contributions of our paper with respect to [17], [18] are the following. In [17], the authors propose and study Stokes problem formulation with IGA discretization. They use the Taylor–Hood isogeometric element, and the NURBS basis functions. They augment the paper with 2D and 3D numerical results for the Stokes flow problem. In [18], the authors introduce the Raviart–Thomas, Nedelec, and Taylor–Hood elements within the IGA formulation. They test the numerical properties of these formulations on the 2D Stokes problem using simple NURBS geometries.

In summary, papers [17], [18] deal with the stationary Stokes problem, while the main focus of our paper is on the non-stationary Stokes a Navier–Stokes problem. Our paper shows that for the implicit time integration scheme with non-stationary Stokes or Navier–Stokes equations, it is possible to perform direction splitting, including the residual minimization method. The computational cost of the resulting direct solver is linear. Our paper’s main scientific contribution is the linear computational cost solver for non-stationary Stokes and Navier–Stokes problems, stabilized with the residual minimization method. We perform numerical experiments for several trial and test spaces.

The common understanding is that the direction splitting solver is limited to tensor product grids. It has been shown in [19], [20] that this solver can be extended to complicated geometries for finite difference simulations. The investigation on the possibility of applying this method for complicated geometries in the context of the finite element method (classical and isogeometric) will be the subject of our future work, either by building iterative solvers or using similar tricks as presented in [19], [20].

In our paper, we advocate linear computational cost alternating direction, implicit solver, on tensor product grids. We show that higher continuity spaces for this kind of solvers reduce the computational cost of the simulations, maintaining high numerical accuracy. However, the standard finite element spaces can be augmented with hp-adaptivity, as presented in DPG framework [10], resulting in exponential convergence of the numerical error with respect to the mesh size. In this paper, we advocate an alternative approach, and we are aware that higher continuity spaces are not suitable for hp adaptivity.

The structure of this paper is as follows: We describe in Section 2 the residual minimization method in an abstract way, and we introduce the functional settings that we use in the paper. Section 3 introduces the non-stationary Stokes and Navier–Stokes problems, together with alternating direction splitting and residual minimization methods. Finally, Section 4 presents numerical results using manufactured solutions and the cavity flow problem [21], [22]. We conclude the paper in Section 5.

Section snippets

The residual minimization method

Let U and V be Hilbert spaces and let b:U×VR be a continuous bilinear form. Given V,1 an abstraction of the variational problems that we will find throughout this paper is: Find uU such thatb(u,v)=(v),vV.We assume that problem (1) is well-posed in the sense that for any V, there exists a unique solution uU, and uUV=b(u,)V with a (stability) constant that only depends on b(,). Let

Time-dependent Stokes and Navier–Stokes problems

Let Ω=Ωx×Ωy be a square domain and I=(0,T]R. We consider the following two-dimensional Navier–Stokes equation: tv+(v)v1ReΔv+p=finΩ×I,v=0inΩ×I,v=0inΓ×I,v(0)=v0inΩ,where v=(v1,v2) is the unknown velocity field and p is the unknown pressure. Here, Γ=ΓxΓy denotes the boundary of the spatial domain Ω, f is a given source, v0 is a given initial condition, and Re stands for the dimensionless Reynolds number. When we drop the advection term (v)v, we obtain the non-stationary Stokes problem,

Numerical results

In this section, we provide several numerical results for non-stationary Stokes and Navier–Stokes problems. First, in Section 4.1, we show that the Navier–Stokes simulation using the alternating directions implicit solver without the RM stabilization produces some oscillations for high Reynolds number. In our paper, we use the splitting scheme from Jean-Luc Guermond and Petar Minev [12], which translates the Navier–Stokes problem into the velocity updates involving reaction-dominated diffusion

Conclusions

We applied residual minimization procedures to non-stationary Stokes and Navier–Stokes problems in simple, functional settings. We experiment with different polynomial order and continuity of both trial and test spaces. For non-stationary Stokes and Navier–Stokes problems, we employed the splitting scheme originally proposed in [12] for finite-difference simulations. We applied the residual minimization (RM) in every time step for stabilization. We obtained a linear computational cost solver

CRediT authorship contribution statement

M. Łoś: Methodology, Software, Visualization. I. Muga: Formal analysis, Methodology, Funding acquisition, Project administration, Supervision, Writing - review & editing. J. Muñoz-Matute: Formal analysis, Methodology, Writing - review & editing. M. Paszyński: Formal analysis, Methodology, Funding acquisition, Project administration, Supervision, Writing - original draft, Writing - review & editing.

Acknowledgements

The work of Maciej Paszyński and Marcin Łoś, and the visit of Judit Muñoz-Matute and Ignacio Muga in Kraków is supported by National Science Centre, Poland Grant No. 2017/26/M/ST1/00281. Ignacio Muga has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie Grant Agreement No. 777778 (MATHROCKS). The authors also want to thank the Chilean Conicyt project PAI80160025, and the unknown reviewers for their insightful suggestions,

References (37)

  • WoźniakM. et al.

    Computational cost estimates for parallel shared memory isogeometric multi-frontal solvers

    Comput. Math. Appl.

    (2014)
  • GurgulG. et al.

    Object-oriented implementation of the alternating directions implicit solver for isogeometric analysis

    Adv. Eng. Softw.

    (2019)
  • LipingG.

    Stability and super convergence analysis of ADI-FDTD for the 2D Maxwell equations in a lossy medium

    Acta Math. Sci.

    (2012)
  • PaszyńskiM. et al.

    Verification of goal-oriented hp-adaptivity

    Comput. Math. Appl.

    (2005)
  • JansenK.E. et al.

    A better consistency for low-order stabilized finite element methods

    Comput. Methods Appl. Mech. Engrg.

    (1997)
  • MatuszykP. et al.
  • ChanJ. et al.

    A Minimum-Residual Finite Element Method for the Convection-Diffusion EquationsICES-Report 13–12

    (2013)
  • CohenA. et al.

    Adaptivity and variational stabilization for convection-diffusion equations

    Math. Modelling Numer. Anal.

    (2012)
  • Cited by (0)

    View full text