www.mooseframework.org
EigenKernel.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 "EigenKernel.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "EigenExecutionerBase.h"
20 #include "Executioner.h"
21 #include "MooseApp.h"
22 #include "MooseEigenSystem.h"
23 #include "MooseVariable.h"
24 
25 #include "libmesh/quadrature.h"
26 
27 template <>
30 {
32  params.addParam<bool>(
33  "eigen", true, "Use for eigenvalue problem (true) or source problem (false)");
34  params.addParam<PostprocessorName>(
35  "eigen_postprocessor", 1.0, "The name of the postprocessor that provides the eigenvalue.");
36  params.registerBase("EigenKernel");
37  return params;
38 }
39 
41  : KernelBase(parameters),
42  _u(_is_implicit ? _var.sln() : _var.slnOld()),
43  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
44  _eigen(getParam<bool>("eigen")),
45  _eigen_sys(dynamic_cast<MooseEigenSystem *>(&_fe_problem.getNonlinearSystemBase())),
46  _eigenvalue(NULL)
47 {
48  // The name to the postprocessor storing the eigenvalue
49  std::string eigen_pp_name;
50 
51  // If the "eigen_postprocessor" is given, use it. The isParamValid does not work here because of
52  // the default value, which
53  // you don't want to use if an EigenExecutioner exists.
54  if (hasPostprocessor("eigen_postprocessor"))
55  eigen_pp_name = getParam<PostprocessorName>("eigen_postprocessor");
56 
57  // Attempt to extract the eigenvalue postprocessor from the Executioner
58  else
59  {
61  if (exec)
62  eigen_pp_name = exec->getParam<PostprocessorName>("bx_norm");
63  }
64 
65  // If the postprocessor name was not provided and an EigenExecutionerBase is not being used,
66  // use the default value from the "eigen_postprocessor" parameter
67  if (eigen_pp_name.empty())
68  _eigenvalue = &getDefaultPostprocessorValue("eigen_postprocessor");
69 
70  // If the name does exist, then use the postprocessor value
71  else
72  {
73  if (_is_implicit)
74  _eigenvalue = &getPostprocessorValueByName(eigen_pp_name);
75  else
76  {
78  if (exec)
79  _eigenvalue = &exec->eigenvalueOld();
80  else
82  }
83  }
84 }
85 
86 void
88 {
89  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
90  _local_re.resize(re.size());
91  _local_re.zero();
92 
93  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
94  Real one_over_eigen = 1.0 / *_eigenvalue;
95  for (_i = 0; _i < _test.size(); _i++)
96  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
97  _local_re(_i) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpResidual();
98 
99  re += _local_re;
100 
101  if (_has_save_in)
102  {
103  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
104  for (const auto & var : _save_in)
105  var->sys().solution().add_vector(_local_re, var->dofIndices());
106  }
107 }
108 
109 void
111 {
112  if (!_is_implicit)
113  return;
114 
115  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
116  _local_ke.resize(ke.m(), ke.n());
117  _local_ke.zero();
118 
119  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
120  Real one_over_eigen = 1.0 / *_eigenvalue;
121  for (_i = 0; _i < _test.size(); _i++)
122  for (_j = 0; _j < _phi.size(); _j++)
123  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
124  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpJacobian();
125 
126  ke += _local_ke;
127 
128  if (_has_diag_save_in)
129  {
130  unsigned int rows = ke.m();
131  DenseVector<Number> diag(rows);
132  for (unsigned int i = 0; i < rows; i++)
133  diag(i) = _local_ke(i, i);
134 
135  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
136  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
137  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
138  }
139 }
140 
141 void
143 {
144  if (!_is_implicit)
145  return;
146 
147  if (jvar == _var.number())
148  computeJacobian();
149  else
150  {
151  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
152  _local_ke.resize(ke.m(), ke.n());
153  _local_ke.zero();
154 
155  mooseAssert(*_eigenvalue != 0.0, "Can't divide by zero eigenvalue in EigenKernel!");
156  Real one_over_eigen = 1.0 / *_eigenvalue;
157  for (_i = 0; _i < _test.size(); _i++)
158  for (_j = 0; _j < _phi.size(); _j++)
159  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
160  _local_ke(_i, _j) +=
161  _JxW[_qp] * _coord[_qp] * one_over_eigen * computeQpOffDiagJacobian(jvar);
162 
163  ke += _local_ke;
164  }
165 }
166 
167 bool
169 {
170  bool flag = MooseObject::enabled();
171  if (_eigen)
172  {
173  if (_is_implicit)
174  return flag && (!_eigen_sys->activeOnOld());
175  else
176  return flag && _eigen_sys->activeOnOld();
177  }
178  else
179  return flag;
180 }
MooseEigenSystem * _eigen_sys
EigenKernel always lives in EigenSystem.
Definition: EigenKernel.h:60
std::vector< MooseVariable * > _diag_save_in
Definition: KernelBase.h:176
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:140
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
Definition: KernelBase.h:167
virtual void computeOffDiagJacobian(unsigned int) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: EigenKernel.C:142
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:66
const VariableTestValue & _test
the current test function
Definition: KernelBase.h:152
bool hasPostprocessor(const std::string &name) const
Determine if the Postprocessor exists.
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.
QBase *& _qrule
active quadrature rule
Definition: KernelBase.h:137
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
Definition: KernelBase.h:164
virtual bool enabled()
Return the enabled status of the object.
Definition: MooseObject.h:77
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:175
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
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:117
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: EigenKernel.C:87
This is the common base class for the two main kernel types implemented in MOOSE, EigenKernel and Ker...
Definition: KernelBase.h:47
const PostprocessorValue & getPostprocessorValueOldByName(const PostprocessorName &name)
This class provides reusable routines for eigenvalue executioners.
virtual bool enabled() override
Return the enabled status of the object.
Definition: EigenKernel.C:168
std::vector< MooseVariable * > _save_in
Definition: KernelBase.h:171
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:170
unsigned int _i
current index for the test function
Definition: KernelBase.h:146
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
Definition: MooseObject.h:122
virtual Real computeQpOffDiagJacobian(unsigned int)
Definition: EigenKernel.h:48
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:143
const VariablePhiValue & _phi
the current shape functions
Definition: KernelBase.h:158
unsigned int _j
current index for the shape function
Definition: KernelBase.h:149
unsigned int number() const
Get variable number coming from libMesh.
const PostprocessorValue & getDefaultPostprocessorValue(const std::string &name)
Return the default postprocessor value.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:108
InputParameters validParams< EigenKernel >()
Definition: EigenKernel.C:29
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...
Executioner * getExecutioner()
Retrieve the Executioner for this App.
Definition: MooseApp.h:233
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: EigenKernel.C:110
virtual Real computeQpResidual()=0
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
EigenKernel(const InputParameters &parameters)
Definition: EigenKernel.C:40
bool _is_implicit
If the object is using implicit or explicit form.
InputParameters validParams< KernelBase >()
Definition: KernelBase.C:27
bool _eigen
flag for as an eigen kernel or a normal kernel
Definition: EigenKernel.h:57
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131
virtual Real computeQpJacobian()
Definition: EigenKernel.h:47