www.mooseframework.org
KernelValue.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 "KernelValue.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 #include "SubProblem.h"
21 #include "SystemBase.h"
22 
23 #include "libmesh/quadrature.h"
24 
25 template <>
28 {
30  return params;
31 }
32 
33 KernelValue::KernelValue(const InputParameters & parameters) : Kernel(parameters) {}
34 
35 void
37 {
38  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
39  _local_re.resize(re.size());
40  _local_re.zero();
41 
42  const unsigned int n_test = _test.size();
43  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
44  {
46  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
47  _local_re(_i) += value * _test[_i][_qp];
48  }
49 
50  re += _local_re;
51 
52  if (_has_save_in)
53  {
54  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
55  for (const auto & var : _save_in)
56  var->sys().solution().add_vector(_local_re, var->dofIndices());
57  }
58 }
59 
60 void
62 {
63  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
64  _local_ke.resize(ke.m(), ke.n());
65  _local_ke.zero();
66 
67  const unsigned int n_test = _test.size();
68  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
69  for (_j = 0; _j < _phi.size(); _j++)
70  {
72  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
73  _local_ke(_i, _j) += value * _test[_i][_qp];
74  }
75 
76  ke += _local_ke;
77 
79  {
80  unsigned int rows = ke.m();
81  DenseVector<Number> diag(rows);
82  for (unsigned int i = 0; i < rows; i++) // target for auto vectorization
83  diag(i) = _local_ke(i, i);
84 
85  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
86  for (const auto & var : _diag_save_in)
87  var->sys().solution().add_vector(diag, var->dofIndices());
88  }
89 }
90 
91 void
93 {
94  if (jvar == _var.number())
96  else
97  {
98  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
99 
100  for (_j = 0; _j < _phi.size(); _j++)
101  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
102  for (_i = 0; _i < _test.size(); _i++)
103  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
104  }
105 }
106 
107 Real
109 {
110  return 0.0;
111 }
112 
113 Real
115 {
116  return 0.0;
117 }
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
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...
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: KernelValue.C:92
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
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:117
KernelValue(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: KernelValue.C:33
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
virtual Real computeQpResidual() override
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
Definition: KernelValue.C:108
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: KernelValue.C:61
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
virtual Real precomputeQpJacobian()
Called before forming the jacobian for an element.
Definition: KernelValue.C:114
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
InputParameters validParams< Kernel >()
Definition: Kernel.C:30
unsigned int number() const
Get variable number coming from libMesh.
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: KernelValue.C:36
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
virtual const VariableValue & value()
The value of the variable this object is operating on.
Definition: Kernel.h:25
InputParameters validParams< KernelValue >()
Definition: KernelValue.C:27
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
virtual Real precomputeQpResidual()=0
Called before forming the residual for an element.
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131