www.mooseframework.org
Kernel.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 "Kernel.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 #include "MooseVariableScalar.h"
21 #include "Problem.h"
22 #include "SubProblem.h"
23 #include "SystemBase.h"
24 
25 #include "libmesh/threads.h"
26 #include "libmesh/quadrature.h"
27 
28 template <>
31 {
33  params.registerBase("Kernel");
34  return params;
35 }
36 
37 Kernel::Kernel(const InputParameters & parameters)
38  : KernelBase(parameters),
39  _u(_is_implicit ? _var.sln() : _var.slnOld()),
40  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
41  _u_dot(_var.uDot()),
42  _du_dot_du(_var.duDotDu())
43 {
44 }
45 
46 void
48 {
49  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
50  _local_re.resize(re.size());
51  _local_re.zero();
52 
54  for (_i = 0; _i < _test.size(); _i++)
55  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
57 
58  re += _local_re;
59 
60  if (_has_save_in)
61  {
62  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
63  for (const auto & var : _save_in)
64  var->sys().solution().add_vector(_local_re, var->dofIndices());
65  }
66 }
67 
68 void
70 {
71  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
72  _local_ke.resize(ke.m(), ke.n());
73  _local_ke.zero();
74 
76  for (_i = 0; _i < _test.size(); _i++)
77  for (_j = 0; _j < _phi.size(); _j++)
78  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
80 
81  ke += _local_ke;
82 
84  {
85  unsigned int rows = ke.m();
86  DenseVector<Number> diag(rows);
87  for (unsigned int i = 0; i < rows; i++)
88  diag(i) = _local_ke(i, i);
89 
90  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
91  for (const auto & var : _diag_save_in)
92  var->sys().solution().add_vector(diag, var->dofIndices());
93  }
94 }
95 
96 void
98 {
99  if (jvar == _var.number())
100  computeJacobian();
101  else
102  {
103  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
104 
106  for (_i = 0; _i < _test.size(); _i++)
107  for (_j = 0; _j < _phi.size(); _j++)
108  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
109  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
110  }
111 }
112 
113 void
115 {
116  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
118 
119  for (_i = 0; _i < _test.size(); _i++)
120  for (_j = 0; _j < jv.order(); _j++)
121  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
122  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
123 }
124 
125 Real
127 {
128  return 0;
129 }
130 
131 Real
132 Kernel::computeQpOffDiagJacobian(unsigned int /*jvar*/)
133 {
134  return 0;
135 }
136 
137 void
139 {
140 }
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: Kernel.C:47
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: Kernel.C:126
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
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
const VariableTestValue & _test
the current test function
Definition: KernelBase.h:152
SystemBase & _sys
Reference to the EquationSystem object.
Definition: KernelBase.h:111
InputParameters validParams< Kernel >()
Definition: Kernel.C:30
THREAD_ID _tid
The thread ID for this kernel.
Definition: KernelBase.h:114
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
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
Kernel(const InputParameters &parameters)
Definition: Kernel.C:37
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
This is the common base class for the two main kernel types implemented in MOOSE, EigenKernel and Ker...
Definition: KernelBase.h:47
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
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
Definition: Kernel.C:132
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
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: Kernel.C:69
unsigned int _j
current index for the shape function
Definition: KernelBase.h:149
unsigned int number() const
Get variable number coming from libMesh.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
virtual void precalculateOffDiagJacobian(unsigned int)
Definition: Kernel.h:48
Class for scalar variables (they are different).
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
virtual void precalculateResidual()
Following methods are used for Kernels that need to perform a per-element calculation.
Definition: Kernel.C:138
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: Kernel.C:97
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:121
InputParameters validParams< KernelBase >()
Definition: KernelBase.C:27
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: Kernel.C:114
virtual void precalculateJacobian()
Definition: Kernel.h:47
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131