www.mooseframework.org
FieldSplitPreconditioner.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "libmesh/petsc_macro.h"
16 #if defined(LIBMESH_HAVE_PETSC) && !PETSC_VERSION_LESS_THAN(3, 3, 0)
18 
19 // MOOSE includes
20 #include "FEProblem.h"
21 #include "MooseEnum.h"
22 #include "MooseVariable.h"
23 #include "NonlinearSystem.h"
24 #include "PetscSupport.h"
25 
26 #include "libmesh/libmesh_common.h"
27 #include "libmesh/petsc_nonlinear_solver.h"
28 #include "libmesh/coupling_matrix.h"
29 
30 template <>
33 {
35 
36  params.addParam<std::vector<std::string>>(
37  "off_diag_row",
38  "The off diagonal row you want to add into the matrix, it will be associated "
39  "with an off diagonal column from the same position in off_diag_colum.");
40  params.addParam<std::vector<std::string>>("off_diag_column",
41  "The off diagonal column you want to add into the "
42  "matrix, it will be associated with an off diagonal "
43  "row from the same position in off_diag_row.");
44  // We should use full coupling Jacobian matrix by default
45  params.addParam<bool>("full",
46  true,
47  "Set to true if you want the full set of couplings. Simply "
48  "for convenience so you don't have to set every off_diag_row "
49  "and off_diag_column combination.");
50  params.addRequiredParam<std::vector<std::string>>(
51  "topsplit", "entrance to splits, the top split will specify how splits will go.");
52  return params;
53 }
54 
56  : MoosePreconditioner(parameters),
57  _top_split(getParam<std::vector<std::string>>("topsplit")),
58  _nl(_fe_problem.getNonlinearSystemBase())
59 {
60  // number of variables
61  unsigned int n_vars = _nl.nVariables();
62  // if we want to construct a full Jacobian?
63  // it is recommended to have a full Jacobian for using
64  // the fieldSplit preconditioner
65  bool full = getParam<bool>("full");
66  // how variables couple
67  std::unique_ptr<CouplingMatrix> cm = libmesh_make_unique<CouplingMatrix>(n_vars);
68  if (!full)
69  {
70  // put 1s on diagonal
71  for (unsigned int i = 0; i < n_vars; i++)
72  (*cm)(i, i) = 1;
73 
74  // off-diagonal entries
75  std::vector<std::vector<unsigned int>> off_diag(n_vars);
76  for (unsigned int i = 0; i < getParam<std::vector<std::string>>("off_diag_row").size(); i++)
77  {
78  unsigned int row =
79  _nl.getVariable(0, getParam<std::vector<std::string>>("off_diag_row")[i]).number();
80  unsigned int column =
81  _nl.getVariable(0, getParam<std::vector<std::string>>("off_diag_column")[i]).number();
82  (*cm)(row, column) = 1;
83  }
84  }
85  else
86  {
87  for (unsigned int i = 0; i < n_vars; i++)
88  for (unsigned int j = 0; j < n_vars; j++)
89  (*cm)(i, j) = 1; // full coupling
90  }
91  _fe_problem.setCouplingMatrix(std::move(cm));
92 
93  // turn on a flag
95 
96  // set a top splitting
98 
99  // apply prefix and store PETSc options
101 }
102 
103 #endif
NonlinearSystemBase & getNonlinearSystemBase()
FieldSplitPreconditioner(const InputParameters &parameters)
Constructor.
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 setCouplingMatrix(std::unique_ptr< CouplingMatrix > cm)
Set custom coupling matrix.
void useFieldSplitPreconditioner(bool use=true)
If called with true this system will use a field split preconditioner matrix.
FEProblemBase & _fe_problem
Subproblem this preconditioner is part of.
InputParameters validParams< MoosePreconditioner >()
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)
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
InputParameters validParams< FieldSplitPreconditioner >()
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
Definition: MooseObject.h:122
std::vector< std::string > _top_split
top split
unsigned int number() const
Get variable number coming from libMesh.
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...
virtual unsigned int nVariables()
Get the number of variables in this system.
Definition: SystemBase.C:580