www.mooseframework.org
FiniteDifferencePreconditioner.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 
11 
12 // MOOSE includes
13 #include "FEProblem.h"
14 #include "MooseVariableFE.h"
15 #include "NonlinearSystem.h"
16 
17 #include "libmesh/coupling_matrix.h"
18 
20 
23 {
25 
26  params.addClassDescription("Finite difference preconditioner (FDP) builds a numerical Jacobian "
27  "for preconditioning, only use for testing and verification.");
28 
29  params.addParam<bool>("implicit_geometric_coupling",
30  false,
31  "Set to true if you want to add entries into the "
32  "matrix for degrees of freedom that might be coupled "
33  "by inspection of the geometric search objects.");
34 
35  MooseEnum finite_difference_type("standard coloring", "coloring");
36  params.addParam<MooseEnum>("finite_difference_type",
37  finite_difference_type,
38  "standard: standard finite difference"
39  "coloring: finite difference based on coloring");
40 
41  return params;
42 }
43 
45  : MoosePreconditioner(params),
46  _finite_difference_type(getParam<MooseEnum>("finite_difference_type"))
47 {
48  if (n_processors() > 1)
49  mooseWarning("Finite differencing to assemble the Jacobian is MUCH MUCH slower than forming "
50  "the Jacobian by hand, so don't complain about performance if you use it!");
51 
53  unsigned int n_vars = nl.nVariables();
54 
55  std::unique_ptr<CouplingMatrix> cm = std::make_unique<CouplingMatrix>(n_vars);
56 
57  bool full = getParam<bool>("full");
58 
59  // standard finite difference method will add off-diagonal entries
60  if (_finite_difference_type == "standard")
61  full = true;
62 
63  if (!full)
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 (isParamValid("off_diag_row") && isParamValid("off_diag_column"))
72 
73  for (const auto i : index_range(getParam<std::vector<NonlinearVariableName>>("off_diag_row")))
74  {
75  unsigned int row =
76  nl.getVariable(0, getParam<std::vector<NonlinearVariableName>>("off_diag_row")[i])
77  .number();
78  unsigned int column =
79  nl.getVariable(0, getParam<std::vector<NonlinearVariableName>>("off_diag_column")[i])
80  .number();
81  (*cm)(row, column) = 1;
82  }
83 
84  // TODO: handle coupling entries between NL-vars and SCALAR-vars
85  }
86  else
87  {
88  for (unsigned int i = 0; i < n_vars; i++)
89  for (unsigned int j = 0; j < n_vars; j++)
90  (*cm)(i, j) = 1;
91  }
92 
93  setCouplingMatrix(std::move(cm));
94 
95  bool implicit_geometric_coupling = getParam<bool>("implicit_geometric_coupling");
96 
97  nl.addImplicitGeometricCouplingEntriesToJacobian(implicit_geometric_coupling);
98 
99  // Set the jacobian to null so that libMesh won't override our finite differenced jacobian
101 }
void addImplicitGeometricCouplingEntriesToJacobian(bool add=true)
If called with true this will add entries into the jacobian to link together degrees of freedom that ...
unsigned int number() const
Get variable number coming from libMesh.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const unsigned int _nl_sys_num
The nonlinear system number whose linearization this preconditioner should be applied to...
FEProblemBase & _fe_problem
Subproblem this preconditioner is part of.
FiniteDifferencePreconditioner(const InputParameters &params)
void mooseWarning(Args &&... args) const
Emits a warning prefixed with object name and type.
Base class for MOOSE preconditioners.
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
Nonlinear system to be solved.
processor_id_type n_processors() const
void useFiniteDifferencedPreconditioner(bool use=true)
If called with true this system will use a finite differenced form of the Jacobian as the preconditio...
unsigned int n_vars
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
static InputParameters validParams()
NonlinearSystemBase & getNonlinearSystemBase(const unsigned int sys_num)
Finite difference preconditioner.
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...
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...
registerMooseObjectAliased("MooseApp", FiniteDifferencePreconditioner, "FDP")
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)