www.mooseframework.org
EigenKernel.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 "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 
20 #include "libmesh/quadrature.h"
21 
22 template <>
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(dynamic_cast<MooseEigenSystem *>(&_fe_problem.getNonlinearSystemBase())),
39  _eigenvalue(NULL)
40 {
41  // The name to the postprocessor storing the eigenvalue
42  std::string eigen_pp_name;
43 
44  // If the "eigen_postprocessor" is given, use it. The isParamValid does not work here because of
45  // the default value, which
46  // you don't want to use if an EigenExecutioner exists.
47  if (hasPostprocessor("eigen_postprocessor"))
48  eigen_pp_name = getParam<PostprocessorName>("eigen_postprocessor");
49 
50  // Attempt to extract the eigenvalue postprocessor from the Executioner
51  else
52  {
54  if (exec)
55  eigen_pp_name = exec->getParam<PostprocessorName>("bx_norm");
56  }
57 
58  // If the postprocessor name was not provided and an EigenExecutionerBase is not being used,
59  // use the default value from the "eigen_postprocessor" parameter
60  if (eigen_pp_name.empty())
61  _eigenvalue = &getDefaultPostprocessorValue("eigen_postprocessor");
62 
63  // If the name does exist, then use the postprocessor value
64  else
65  {
66  if (_is_implicit)
67  _eigenvalue = &getPostprocessorValueByName(eigen_pp_name);
68  else
69  {
71  if (exec)
72  _eigenvalue = &exec->eigenvalueOld();
73  else
75  }
76  }
77 }
78 
79 void
81 {
82  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
83  _local_re.resize(re.size());
84  _local_re.zero();
85 
86  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
87  Real one_over_eigen = 1.0 / *_eigenvalue;
88  for (_i = 0; _i < _test.size(); _i++)
89  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
90  _local_re(_i) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpResidual();
91 
92  re += _local_re;
93 
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 
108  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
109  _local_ke.resize(ke.m(), ke.n());
110  _local_ke.zero();
111 
112  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
113  Real one_over_eigen = 1.0 / *_eigenvalue;
114  for (_i = 0; _i < _test.size(); _i++)
115  for (_j = 0; _j < _phi.size(); _j++)
116  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
117  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpJacobian();
118 
119  ke += _local_ke;
120 
121  if (_has_diag_save_in)
122  {
123  unsigned int rows = ke.m();
124  DenseVector<Number> diag(rows);
125  for (unsigned int i = 0; i < rows; i++)
126  diag(i) = _local_ke(i, i);
127 
128  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
129  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
130  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
131  }
132 }
133 
134 void
136 {
137  size_t jvar_num = jvar.number();
138  if (!_is_implicit)
139  return;
140 
141  if (jvar_num == _var.number())
142  computeJacobian();
143  else
144  {
145  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_num);
146  _local_ke.resize(ke.m(), ke.n());
147  _local_ke.zero();
148 
149  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
150  Real one_over_eigen = 1.0 / *_eigenvalue;
151  for (_i = 0; _i < _test.size(); _i++)
152  for (_j = 0; _j < jvar.phiSize(); _j++)
153  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
154  _local_ke(_i, _j) +=
155  _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpOffDiagJacobian(jvar_num);
156 
157  ke += _local_ke;
158  }
159 }
160 
161 bool
163 {
164  bool flag = MooseObject::enabled();
165  if (_eigen)
166  {
167  if (_is_implicit)
168  return flag && (!_eigen_sys->activeOnOld());
169  else
170  return flag && _eigen_sys->activeOnOld();
171  }
172  else
173  return flag;
174 }
MooseEigenSystem * _eigen_sys
EigenKernel always lives in EigenSystem.
Definition: EigenKernel.h:46
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: KernelBase.h:111
Executioner * getExecutioner() const
Retrieve the Executioner for this App.
Definition: MooseApp.h:255
bool hasPostprocessor(const std::string &name) const
Determine if the Postprocessor exists.
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:54
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:172
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:159
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:52
unsigned int number() const
Get variable number coming from libMesh.
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2019
const PostprocessorValue & getPostprocessorValueByName(const PostprocessorName &name)
Retrieve the value of the Postprocessor.
bool activeOnOld()
Return if eigen kernels should be on old solution.
const Real & eigenvalueOld()
The old eigenvalue used by inverse power iterations.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:176
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
Definition: MooseObject.h:191
virtual size_t phiSize() const =0
Return phi size.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:139
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: EigenKernel.C:80
const VariableTestValue & _test
the current test function
Definition: Kernel.h:57
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:177
const PostprocessorValue & getPostprocessorValueOldByName(const PostprocessorName &name)
This class provides reusable routines for eigenvalue executioners.
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:156
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:171
virtual void computeOffDiagJacobian(MooseVariableFEBase &) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: EigenKernel.C:135
virtual bool enabled() const override
Return the enabled status of the object.
Definition: EigenKernel.C:162
unsigned int _i
current index for the test function
Definition: KernelBase.h:165
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:162
unsigned int _j
current index for the shape function
Definition: KernelBase.h:168
InputParameters validParams< Kernel >()
Definition: Kernel.C:24
const PostprocessorValue & getDefaultPostprocessorValue(const std::string &name)
Return the default postprocessor value.
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:177
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
Definition: Kernel.h:20
InputParameters validParams< EigenKernel >()
Definition: EigenKernel.C:24
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 void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: EigenKernel.C:103
virtual Real computeQpOffDiagJacobian(unsigned int)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
Definition: KernelBase.h:116
DenseVector< Number > & residualBlock(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:720
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:63
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:99
bool _eigen
flag for as an eigen kernel or a normal kernel
Definition: EigenKernel.h:43
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:150