16 #include "libmesh/implicit_system.h" 17 #include "libmesh/static_condensation_dof_map.h" 26 params.
set<
bool>(
"full") =
true;
28 params.
addParam<std::vector<NonlinearVariableName>>(
31 "A list of variables for whom to not statically condense their degrees of freedom out of the " 32 "system. By default all degrees of freedom on element interiors are condensed out.");
43 return _nl.
name() +
"_condensed_";
50 auto check_param = [
this](
const auto & param_name)
54 "This class prefixes every PETSc option so that it applies to the condensed " 55 "system. Given that, there are multiple issues with setting '",
57 "': 1) it applies to the nonlinear solver algorithm whereas the prefix this class " 58 "applies makes PETSc options conceptually applicable only to the linear solve of " 59 "the statically condensed system 2) these are singleton MOOSE-wrapped PETSc " 60 "options. Consequently even if having multiple prefixes for a system's nonlinear " 61 "solver options made sense, we don't support it. E.g. if you specify '",
63 "' in both this object's block and the Executioner block, then there will be a " 66 check_param(
"mffd_type");
67 check_param(
"solve_type");
72 "Static condensation preconditioning should use a PJFNK solve type. This is because it " 73 "provides a means to compute the action of the Jacobian on vectors, which we otherwise " 74 "would not have because when using static condensation, the Jacobian is never formed. Note " 75 "that actions of the Jacobian on a vector are necessary for things like: GMRES, printing " 76 "of linear residuals, or application of line searches like cubic backtracking. These " 77 "particular operations can be avoided by using '-ksp_type preonly', disabling printing of " 78 "linear residuals, and using the 'cp' or 'basic' line search.");
82 mooseError(
"Static condensation can only be used with implicit systems");
83 implicit_sys->create_static_condensation();
84 _sc_dof_map = &implicit_sys->get_dof_map().get_static_condensation();
86 std::unordered_set<unsigned int> uncondensed_vars;
87 for (
auto & nl_var_name :
getParam<std::vector<NonlinearVariableName>>(
"dont_condense_vars"))
NonlinearSystemBase & _nl
The nonlinear system whose linearization this preconditioner should be applied to.
virtual void initialSetup() override
Perform some setup tasks such as storing the PETSc options.
void paramError(const std::string ¶m, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
unsigned int number() const
Get variable number coming from libMesh.
static InputParameters validParams()
FEProblemBase & _fe_problem
Subproblem this preconditioner is part of.
static InputParameters validParams()
libMesh::StaticCondensationDofMap * _sc_dof_map
Pointer to the libMesh static condensation dof map object.
virtual const std::string & name() const
void dont_condense_vars(const std::unordered_set< unsigned int > &vars)
std::string prefix() const
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
unsigned int number() const
Gets the number of this system.
libMesh::StaticCondensation * _sc_system_matrix
Pointer to the libMesh static condensation system matrix object.
Static condensation preconditioner.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type and optionally a file path to the top-level block p...
SolverParams & solverParams(unsigned int solver_sys_num=0)
Get the solver parameters.
Single matrix preconditioner.
Preconditioned Jacobian-Free Newton Krylov.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
registerMooseObjectAliased("MooseApp", MooseStaticCondensationPreconditioner, "StaticCondensation")
MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a variable of with specified name.
void storePetscOptions(FEProblemBase &fe_problem, const std::string &prefix, const ParallelParamObject ¶m_object)
Stores the PETSc options supplied from the parameter object on the problem.
MooseStaticCondensationPreconditioner(const InputParameters ¶ms)
const std::string & _type
The type of this class.
virtual libMesh::System & system() override
Get the reference to the libMesh system.