https://mooseframework.inl.gov
EigenKernel.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 "EigenKernel.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "EigenExecutionerBase.h"
15 #include "Executioner.h"
16 #include "MooseApp.h"
17 #include "MooseEigenSystem.h"
18 #include "MooseVariableFE.h"
19 #include "FEProblemBase.h"
20 
21 #include "libmesh/quadrature.h"
22 
25 {
27  params.addParam<bool>(
28  "eigen", true, "Use for eigenvalue problem (true) or source problem (false)");
29  params.addParam<PostprocessorName>(
30  "eigen_postprocessor", 1.0, "The name of the postprocessor that provides the eigenvalue.");
31  params.registerBase("EigenKernel");
32  return params;
33 }
34 
36  : Kernel(parameters),
37  _eigen(getParam<bool>("eigen")),
38  _eigen_sys(
39  dynamic_cast<MooseEigenSystem *>(&_fe_problem.getNonlinearSystemBase(_sys.number()))),
40  _eigenvalue(NULL)
41 {
42  // The name to the postprocessor storing the eigenvalue
43  std::string eigen_pp_name;
44 
45  // If the "eigen_postprocessor" is given, use it. The isParamValid does not work here because of
46  // the default value, which
47  // you don't want to use if an EigenExecutioner exists.
48  if (!isDefaultPostprocessorValue("eigen_postprocessor"))
49  eigen_pp_name = getPostprocessorName("eigen_postprocessor");
50 
51  // Attempt to extract the eigenvalue postprocessor from the Executioner
52  else
53  {
55  if (exec)
56  eigen_pp_name = exec->getParam<PostprocessorName>("bx_norm");
57  }
58 
59  // If the postprocessor name was not provided and an EigenExecutionerBase is not being used,
60  // use the default value from the "eigen_postprocessor" parameter
61  if (eigen_pp_name.empty())
62  _eigenvalue = &getPostprocessorValue("eigen_postprocessor");
63 
64  // If the name does exist, then use the postprocessor value
65  else
66  {
67  if (_is_implicit)
68  _eigenvalue = &getPostprocessorValueByName(eigen_pp_name);
69  else
70  {
72  if (exec)
73  _eigenvalue = &exec->eigenvalueOld();
74  else
76  }
77  }
78 }
79 
80 void
82 {
84 
86 
87  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
88  Real one_over_eigen = 1.0 / *_eigenvalue;
89  for (_i = 0; _i < _test.size(); _i++)
90  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
91  _local_re(_i) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpResidual();
92 
94  if (_has_save_in)
95  {
96  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
97  for (const auto & var : _save_in)
98  var->sys().solution().add_vector(_local_re, var->dofIndices());
99  }
100 }
101 
102 void
104 {
105  if (!_is_implicit)
106  return;
107 
109 
111  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
112  Real one_over_eigen = 1.0 / *_eigenvalue;
113  for (_i = 0; _i < _test.size(); _i++)
114  for (_j = 0; _j < _phi.size(); _j++)
115  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
116  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpJacobian();
117 
119 
121  {
123  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
124  for (const auto & var : _diag_save_in)
125  var->sys().solution().add_vector(diag, var->dofIndices());
126  }
127 }
128 
129 void
130 EigenKernel::computeOffDiagJacobian(const unsigned int jvar_num)
131 {
132  if (!_is_implicit)
133  return;
134 
135  if (jvar_num == _var.number())
136  computeJacobian();
137  else
138  {
139  const auto & jvar = getVariable(jvar_num);
140 
141  prepareMatrixTag(_assembly, _var.number(), jvar_num);
142 
143  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
144  // on the displaced mesh
145  auto phi_size = jvar.dofIndices().size();
146  mooseAssert(
147  phi_size * jvar.count() == _local_ke.n(),
148  "The size of the phi container does not match the number of local Jacobian columns");
149 
150  if (_local_ke.m() != _test.size())
151  return;
152 
153  precalculateOffDiagJacobian(jvar_num);
154  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
155  Real one_over_eigen = 1.0 / *_eigenvalue;
156  if (jvar.count() == 1)
157  {
158  for (_i = 0; _i < _test.size(); _i++)
159  for (_j = 0; _j < phi_size; _j++)
160  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
161  _local_ke(_i, _j) +=
162  _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpOffDiagJacobian(jvar.number());
163  }
164  else
165  {
166  unsigned int n = phi_size;
167  for (_i = 0; _i < _test.size(); _i++)
168  for (_j = 0; _j < n; _j++)
169  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
170  {
171  RealEigenVector v = _JxW[_qp] * _coord[_qp] * one_over_eigen *
172  computeQpOffDiagJacobianArray(static_cast<ArrayMooseVariable &>(
173  const_cast<MooseVariableFieldBase &>(jvar)));
174  for (unsigned int k = 0; k < v.size(); ++k)
175  _local_ke(_i, _j + k * n) += v(k);
176  }
177  }
178 
180  }
181 }
182 
183 bool
185 {
186  bool flag = MooseObject::enabled();
187  if (_eigen)
188  {
189  if (!_eigen_sys)
190  mooseError("Eigen kernel ",
191  name(),
192  " requires a MooseEigenSystem and was designed to work with old eigenvalue",
193  " executioners such as 'NonlinearEigen'. It is suggested to use the new",
194  " eigenvalue executioner 'Eigenvalue' along with kernel tagging");
195 
196  if (_is_implicit)
197  return flag && (!_eigen_sys->activeOnOld());
198  else
199  return flag && _eigen_sys->activeOnOld();
200  }
201  else
202  return flag;
203 }
MooseEigenSystem * _eigen_sys
EigenKernel always lives in EigenSystem.
Definition: EigenKernel.h:41
static InputParameters validParams()
Definition: Kernel.C:24
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
virtual RealEigenVector computeQpOffDiagJacobianArray(const ArrayMooseVariable &jvar)
For coupling array variables.
Definition: Kernel.h:66
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:72
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:65
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:52
virtual void precalculateOffDiagJacobian(unsigned int)
const Real * _eigenvalue
A pointer to the eigenvalue that is stored in a postprocessor This is a pointer so that the method fo...
Definition: EigenKernel.h:47
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int)
For coupling standard variables.
Definition: Kernel.h:56
bool activeOnOld()
Return if eigen kernels should be on old solution.
const Real & eigenvalueOld()
The old eigenvalue used by inverse power iterations.
virtual void precalculateResidual()
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool computingScalingJacobian() const
Whether we are computing an initial Jacobian for automatic variable scaling.
Definition: SystemBase.C:1519
unsigned int m() const
const PostprocessorValue & getPostprocessorValueOldByName(const PostprocessorName &name) const
const PostprocessorValue & getPostprocessorValue(const std::string &param_name, const unsigned int index=0) const
doco-normal-methods-begin Retrieve the value of a Postprocessor or one of it&#39;s old or older values ...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:69
DenseVector< Real > getJacobianDiagonal(DenseMatrix< Number > &ke)
Definition: Assembly.h:1857
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:57
static InputParameters validParams()
Definition: EigenKernel.C:24
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
SystemBase & _sys
Reference to the EquationSystem object.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: EigenKernel.C:81
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
const VariableTestValue & _test
the current test function
Definition: Kernel.h:75
virtual void precalculateJacobian()
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:70
This class provides reusable routines for eigenvalue executioners.
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:49
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:64
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: Kernel.h:51
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
virtual bool enabled() const override
Return the enabled status of the object.
Definition: EigenKernel.C:184
unsigned int _i
current index for the test function
Definition: KernelBase.h:58
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: EigenKernel.C:130
MooseApp & _app
The MOOSE application this is associated with.
Definition: MooseBase.h:84
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
bool isDefaultPostprocessorValue(const std::string &param_name, const unsigned int index=0) const
Determine whether or not the Postprocessor is a default value.
Executioner * getExecutioner() const
Retrieve the Executioner for this App.
Definition: MooseApp.C:2085
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:55
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual const PostprocessorValue & getPostprocessorValueByName(const PostprocessorName &name) const
Retrieve the value of the Postprocessor.
unsigned int _j
current index for the shape function
Definition: KernelBase.h:61
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
Definition: Kernel.h:15
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: EigenKernel.C:103
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
unsigned int n() const
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:81
EigenKernel(const InputParameters &parameters)
Definition: EigenKernel.C:35
bool _is_implicit
If the object is using implicit or explicit form.
virtual bool enabled() const
Return the enabled status of the object.
Definition: MooseObject.h:40
const PostprocessorName & getPostprocessorName(const std::string &param_name, const unsigned int index=0) const
Get the name of a postprocessor.
bool _eigen
flag for as an eigen kernel or a normal kernel
Definition: EigenKernel.h:38
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:43