index.rst

FENaPack

Mathematical background

Navier-Stokes equations

$$\begin{aligned} \left[\begin{array}{cc} -\nu\Delta + \mathbf{v}\cdot\nabla & \nabla \\\ -\operatorname{div} & 0 \end{array}\right] \left[\begin{array}{c} \mathbf{u} \\\ p \end{array}\right] = \left[\begin{array}{c} \mathbf{f} \\\ 0 \end{array}\right] \end{aligned}$$

solved by GMRES preconditioned from right by

$$\begin{aligned} \mathbb{P} := \left[\begin{array}{cc} -\nu\Delta + \mathbf{v}\cdot\nabla & \nabla \\\ & -\mathbb{S} \end{array}\right] \end{aligned}$$

with Schur complement 𝕊 =  − div (−νΔ+v⋅∇)^− 1∇ would converge in two iterations. Unfortunately 𝕊 is dense. Possible trick is to approximate 𝕊 by swapping the order of the operators

𝕊 ≈ 𝕏_BRM1 := ( − Δ)(−νΔ+v⋅∇)^− 1

or

𝕊 ≈ 𝕏_BRM2 := (−νΔ+v⋅∇)^− 1( − Δ).

This gives rise to the action of 11-block of preconditioner ℙ^− 1 given by

𝕏_BRM1^− 1 := (−νΔ+v⋅∇)( − Δ)^− 1.

or

𝕏_BRM2^− 1 := ( − Δ)^− 1(−νΔ+v⋅∇).

Obviously additional artificial boundary condition for Laplacian solve  − Δ^− 1 is needed in the action of preconditioner. Modifying the approach from we implement 𝕏_BRM1^− 1 as

𝕏_BRM1^− 1 := 𝕄_p^− 1(𝕀 + 𝕂_p𝔸_p^− 1)

where 𝕄_p is ν^− 1-multiple of mass matrix on pressure, 𝕂_p ≈ ν^− 1v ⋅ ∇ is a pressure convection matrix, and 𝔸_p^− 1 ≈ ( − Δ)^− 1 is a pressure Laplacian solve with zero boundary condition on inlet. This is implemented by :pyfenapack.preconditioners.PCDPC_BRM1 and demos/navier-stokes-pcd.

Analogically we prefer to express BRM2 approach as

𝕏_BRM2^− 1 := (𝕀 + 𝔸_p^− 1𝕂_p)𝕄_p^− 1

now with zero boundary condition on outlet for Laplacian solve and additional Robin term in convection matrix 𝕂_p roughly as stated in, section 9.2.2. See also demos/navier-stokes-pcd and :pyfenapack.preconditioners.PCDPC_BRM2.

Extension to time-dependent problems (PCDR preconditioners)

Time disretization applied in unsteady problems typically leads to the need to incorporate a reaction term into the preconditioner. Typically, we end up with

$$\mathbb{X}_\mathrm{BRM1}^{-1} := \left(\frac{1}{\tau}-\nu\Delta+\mathbf{v}\cdot\nabla\right)(-\Delta)^{-1},$$

or

$$\mathbb{X}_\mathrm{BRM2}^{-1} := (-\Delta)^{-1}\left(\frac{1}{\tau}-\nu\Delta+\mathbf{v}\cdot\nabla\right),$$

where τ denotes a fixed time step and the original PCD preconditioner thus becomes PCDR (pressure-convection-diffusion-reaction) preconditioner. A straightforward way of how to implement the above actions is to update the pressure convection matrix 𝕂_p by a contribution corresponding to the scaled pressure mass matrix, namely

𝕏_BRM1^− 1 := 𝕄_p^− 1(𝕀+(𝕂_p+τ^− 1𝕄_p)𝔸_p^− 1),

or

𝕏_BRM2^− 1 := (𝕀+𝔸_p^− 1(𝕂_p+τ^− 1𝕄_p))𝕄_p^− 1.

However, for unsteady problems we prefer to use the following elaborated implementation of PCDR preconditioners, namely

𝕏_BRM1^− 1 := ℝ_p^− 1 + 𝕄_p^− 1(𝕀 + 𝕂_p𝔸_p^− 1),

or

𝕏_BRM2^− 1 := ℝ_p^− 1 + (𝕀 + 𝔸_p^− 1𝕂_p)𝕄_p^− 1,

where $\mathbb{R}_p^{-1} \approx \frac{1}{\tau} (-\Delta)^{-1}$, while ℝ_p itself is approximated and implemented as

ℝ_p := 𝔹(τ^− 1𝔻_M)^− 1𝔹^T,

Here, 𝔻_M is the diagonal of the velocity mass matrix, 𝔻_M = diag (𝕄_u), and 𝔹^T corresponds to the discrete pressure gradient which is obtained as the 01-block of the original system matrix. Let us emphasize that this submatrix is extracted from the system matrix with velocity Dirichlet boundary conditions being applied on it.

The choice of ℝ_p as above can be justified especially in the case of τ → 0₊, for which

$$\mathbb{S}^{-1} := \left(-\operatorname{div}\left(\frac{1}{\tau} - \nu\Delta+\mathbf{v}\cdot\nabla\right)^{-1}\nabla\right)^{-1} \approx \frac{1}{\tau}\left(\mathbb{B} \mathbb{M}_{\mathbf{u}}^{-1} \mathbb{B}^T\right)^{-1},$$

and simultaneously $\mathbb{X}^{-1} \approx \mathbb{R}_p^{-1} = \frac{1}{\tau} \left(\mathbb{B} \mathbb{D}_\mathrm{M}^{-1} \mathbb{B}^T\right)^{-1}$. The same approximation of the minus Laplacian operator was previously used also in, see Remark 9.6 therein.

demos/navier-stokes-pcd demos/unsteady-navier-stokes-pcd

circle

Manual and API Reference

API Reference <api-doc/fenapack>

• modindex
• genindex
• search