www.mooseframework.org
FieldSplitPreconditioner.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "libmesh/petsc_macro.h"
12 
13 // MOOSE includes
14 #include "FEProblem.h"
15 #include "MooseEnum.h"
16 #include "MooseVariableFE.h"
17 #include "NonlinearSystem.h"
18 #include "PetscSupport.h"
19 
20 #include "libmesh/libmesh_common.h"
21 #include "libmesh/petsc_nonlinear_solver.h"
22 #include "libmesh/coupling_matrix.h"
23 
25 
28 {
30  params.addClassDescription("Preconditioner designed to map onto PETSc's PCFieldSplit.");
31 
32  params.addRequiredParam<std::vector<std::string>>(
33  "topsplit", "Entrance to splits, the top split will specify how splits will go.");
34  // We should use full coupling Jacobian matrix by default
35  params.addParam<bool>("full",
36  true,
37  "Set to true if you want the full set of couplings between variables "
38  "simply for convenience so you don't have to set every off_diag_row "
39  "and off_diag_column combination.");
40  return params;
41 }
42 
44  : MoosePreconditioner(parameters),
45  _top_split(getParam<std::vector<std::string>>("topsplit")),
46  _nl(_fe_problem.getNonlinearSystemBase(_nl_sys_num))
47 {
48  // number of variables
49  unsigned int n_vars = _nl.nVariables();
50  // if we want to construct a full Jacobian?
51  // it is recommended to have a full Jacobian for using
52  // the fieldSplit preconditioner
53  bool full = getParam<bool>("full");
54 
55  // how variables couple
56  std::unique_ptr<CouplingMatrix> cm = std::make_unique<CouplingMatrix>(n_vars);
57  if (!full)
58  {
59  if (isParamValid("off_diag_row") && isParamValid("off_diag_column"))
60  {
61 
62  const auto off_diag_rows = getParam<std::vector<NonlinearVariableName>>("off_diag_row");
63  const auto off_diag_columns = getParam<std::vector<NonlinearVariableName>>("off_diag_column");
64 
65  // put 1s on diagonal
66  for (unsigned int i = 0; i < n_vars; i++)
67  (*cm)(i, i) = 1;
68 
69  // off-diagonal entries
70  std::vector<std::vector<unsigned int>> off_diag(n_vars);
71  if (off_diag_rows.size() * off_diag_columns.size() != 0 &&
72  off_diag_rows.size() == off_diag_columns.size())
73  for (const auto i : index_range(off_diag_rows))
74  {
75  unsigned int row = _nl.getVariable(0, off_diag_rows[i]).number();
76  unsigned int column = _nl.getVariable(0, off_diag_columns[i]).number();
77  (*cm)(row, column) = 1;
78  }
79  }
80  }
81  else
82  {
83  for (unsigned int i = 0; i < n_vars; i++)
84  for (unsigned int j = 0; j < n_vars; j++)
85  (*cm)(i, j) = 1; // full coupling
86  }
87  setCouplingMatrix(std::move(cm));
88 
89  // turn on a flag
91 
92  // set a top splitting
94 
95  // apply prefix and store PETSc options
97 }
Implements a preconditioner designed to map onto PETSc&#39;s PCFieldSplit.
unsigned int number() const
Get variable number coming from libMesh.
FieldSplitPreconditioner(const InputParameters &parameters)
void setDecomposition(const std::vector< std::string > &decomposition)
If called with a single string, it is used as the name of a the top-level decomposition split...
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void useFieldSplitPreconditioner(bool use=true)
If called with true this system will use a field split preconditioner matrix.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
Base class for MOOSE preconditioners.
NonlinearSystemBase & _nl
The nonlinear system this FSP is associated with (convenience reference)
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
virtual unsigned int nVariables() const
Get the number of variables in this system.
Definition: SystemBase.C:840
unsigned int n_vars
static InputParameters validParams()
std::vector< std::string > _top_split
top split
static InputParameters validParams()
Constructor.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
registerMooseObjectAliased("MooseApp", FieldSplitPreconditioner, "FSP")
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:79
void setCouplingMatrix(std::unique_ptr< CouplingMatrix > cm)
Setup the coupling matrix on the finite element problem.
auto index_range(const T &sizable)