INSADMomentumGradDiv

This kernel implements grad-div stabilization for the purpose of building scalable preconditioners at high Reynolds numbers for Navier-Stokes as introduced in (Benzi and Olshanskii, 2006). Its weak form is given by

γ(div uh,div vh)\gamma\left(\textrm{div}\ \mathbf{u}_h, \textrm{div}\ \mathbf{v}_h\right)

where γ\gamma is a penalty parameter, uh\mathbf{u}_h is the velocity vector approximate solution, and vh\mathbf{v}_h are its associated test functions. For γ\gamma not too small, the Schur complement inverse of a linearized (via Newton or Picard) stable (no pressure on-diagonal) Navier-Stokes discretization is well approximated by

S1=(ν+γ)Mp1S^{-1} = -\left(\nu + \gamma\right)M_p^{-1}

where ν\nu is the kinematic viscosity and MpM_p is the pressure mass matrix ((Farrell et al., 2019)). In (Benzi and Olshanskii, 2006) the authors found that γ\gamma on the order of the velocity magnitude provides ν\nu- and hh-independent (the latter denoting element size) convergence for a wide range of values for ν\nu and hh. For using the (scaled) pressure mass matrix as the preconditioner for the Schur complement, the NavierStokesProblem should be used in conjunction with FSP. A demonstration of this combination is given in the steady_vector_fsp_al.i input

[Problem]
  type = NavierStokesProblem
  mass_matrix = 'mass'
  extra_tag_matrices = 'mass'
  use_pressure_mass_matrix = true
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)

where we have indicated that we should use the pressure mass matrix as the preconditioner for the Schur complement using the "use_pressure_mass_matrix" parameter. The scaled pressure mass matrix is built using the MassMatrix kernel with block

[Kernels]
  [mass_kernel]
    type = MassMatrix
    variable = p
    matrix_tags = 'mass'
    density = '${fparse -1/(gamma + mu)}'
  []
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)

The field split options are shown in

[Preconditioning]
  [FSP]
    type = FSP
    topsplit = 'up'
    [up]
      splitting = 'u p'
      splitting_type = schur
      petsc_options_iname = '-pc_fieldsplit_schur_fact_type  -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
      petsc_options_value = 'full                            self                              300                fgmres    right        1e-4'
    []
    [u]
      vars = 'vel'
      petsc_options = '-ksp_converged_reason'
      petsc_options_iname = '-pc_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side -pc_factor_mat_solver_type'
      petsc_options_value = 'ilu      gmres     1e-2      300                right        strumpack'
    []
    [p]
      vars = 'p'
      petsc_options = '-ksp_converged_reason'
      petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side'
      petsc_options_value = 'gmres     300                1e-2      ilu      right'
    []
  []
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)

In this case we only wish to display the effectiveness of the pressure mass matrix as a preconditioner for the Schur complement, so we use LU solves for inverting the velocity block and MpM_p. Indeed for ν=103\nu=10^{-3} (equivalent to Re=1000\textrm{Re}=1000), the solution of the Schur complement system to a relative tolerance of (10210^{-2}) takes no more than 4 Krylov iterations when running on grids ranging from 8x8 to 32x32. Scalable solves for the velocity block in the context of augmented Lagrange formulations are the subject of Benzi and Faccio (2023).

Input Parameters

  • gammaThe grad-div stabilization coefficient

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The grad-div stabilization coefficient

  • variableThe name of the variable that this residual object operates on

    C++ Type:NonlinearVariableName

    Unit:(no unit assumed)

    Controllable:No

    Description:The name of the variable that this residual object operates on

Required Parameters

  • blockThe list of blocks (ids or names) that this object will be applied

    C++ Type:std::vector<SubdomainName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The list of blocks (ids or names) that this object will be applied

  • displacementsThe displacements

    C++ Type:std::vector<VariableName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The displacements

  • prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.

    C++ Type:MaterialPropertyName

    Unit:(no unit assumed)

    Controllable:No

    Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.

  • use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

Optional Parameters

  • absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution

    C++ Type:std::vector<TagName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution

  • extra_matrix_tagsThe extra tags for the matrices this Kernel should fill

    C++ Type:std::vector<TagName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The extra tags for the matrices this Kernel should fill

  • extra_vector_tagsThe extra tags for the vectors this Kernel should fill

    C++ Type:std::vector<TagName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The extra tags for the vectors this Kernel should fill

  • matrix_tagssystemThe tag for the matrices this Kernel should fill

    Default:system

    C++ Type:MultiMooseEnum

    Unit:(no unit assumed)

    Options:nontime, system

    Controllable:No

    Description:The tag for the matrices this Kernel should fill

  • vector_tagsnontimeThe tag for the vectors this Kernel should fill

    Default:nontime

    C++ Type:MultiMooseEnum

    Unit:(no unit assumed)

    Options:nontime, time

    Controllable:No

    Description:The tag for the vectors this Kernel should fill

Tagging Parameters

  • control_tagsAdds user-defined labels for accessing object parameters via control logic.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Adds user-defined labels for accessing object parameters via control logic.

  • diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)

    C++ Type:std::vector<AuxVariableName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)

  • enableTrueSet the enabled status of the MooseObject.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:Yes

    Description:Set the enabled status of the MooseObject.

  • implicitTrueDetermines whether this object is calculated using an implicit or explicit form

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Determines whether this object is calculated using an implicit or explicit form

  • save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)

    C++ Type:std::vector<AuxVariableName>

    Unit:(no unit assumed)

    Controllable:No

    Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)

  • seed0The seed for the master random number generator

    Default:0

    C++ Type:unsigned int

    Unit:(no unit assumed)

    Controllable:No

    Description:The seed for the master random number generator

  • use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

References

  1. Michele Benzi and Chiara Faccio. Solving linear systems of the form (a + uu^t)x = b by preconditioned iterative methods. SIAM Journal on Scientific Computing, pages S51–S70, 2023.[BibTeX]
  2. Michele Benzi and Maxim A Olshanskii. An augmented lagrangian-based approach to the oseen problem. SIAM Journal on Scientific Computing, 28(6):2095–2113, 2006.[BibTeX]
  3. Patrick E Farrell, Lawrence Mitchell, and Florian Wechsung. An augmented lagrangian preconditioner for the 3d stationary incompressible navier–stokes equations at high reynolds number. SIAM Journal on Scientific Computing, 41(5):A3073–A3096, 2019.[BibTeX]