www.mooseframework.org
VectorKernel.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 "VectorKernel.h"
11 #include "Assembly.h"
12 #include "MooseVariableFE.h"
13 #include "MooseVariableScalar.h"
14 #include "SubProblem.h"
15 #include "NonlinearSystem.h"
16 
17 #include "libmesh/threads.h"
18 #include "libmesh/quadrature.h"
19 
20 template <>
23 {
25  params.registerBase("VectorKernel");
26  return params;
27 }
28 
30  : KernelBase(parameters),
32  false,
33  "variable",
36  _var(*mooseVariable()),
37  _test(_var.phi()),
38  _grad_test(_var.gradPhi()),
39  _phi(_assembly.phi(_var)),
40  _grad_phi(_assembly.gradPhi(_var)),
41  _u(_is_implicit ? _var.sln() : _var.slnOld()),
42  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
43 {
45 }
46 
47 void
49 {
50  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
51  _local_re.resize(re.size());
52  _local_re.zero();
53 
55  for (_i = 0; _i < _test.size(); _i++)
56  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
57  {
58  Real residual = _JxW[_qp] * _coord[_qp] * computeQpResidual();
59  _local_re(_i) += residual;
60  }
61 
62  re += _local_re;
63 
64  if (_has_save_in)
65  {
66  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
67  for (const auto & var : _save_in)
68  var->sys().solution().add_vector(_local_re, var->dofIndices());
69  }
70 }
71 
72 void
74 {
75  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
76  _local_ke.resize(ke.m(), ke.n());
77  _local_ke.zero();
78 
80  for (_i = 0; _i < _test.size(); _i++)
81  for (_j = 0; _j < _phi.size(); _j++)
82  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
84 
85  ke += _local_ke;
86 
88  {
89  unsigned int rows = ke.m();
90  DenseVector<Number> diag(rows);
91  for (unsigned int i = 0; i < rows; i++)
92  diag(i) = _local_ke(i, i);
93 
94  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
95  for (const auto & var : _diag_save_in)
96  var->sys().solution().add_vector(diag, var->dofIndices());
97  }
98 }
99 
100 void
102 {
103  size_t jvar_num = jvar.number();
104  if (jvar_num == _var.number())
105  computeJacobian();
106  else
107  {
108  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_num);
109 
110  precalculateOffDiagJacobian(jvar_num);
111  for (_i = 0; _i < _test.size(); _i++)
112  for (_j = 0; _j < jvar.phiSize(); _j++)
113  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
114  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
115  }
116 }
117 
118 void
120 {
121  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
123 
124  for (_i = 0; _i < _test.size(); _i++)
125  for (_j = 0; _j < jv.order(); _j++)
126  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
127  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
128 }
VarFieldType
Definition: MooseTypes.h:488
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: KernelBase.h:111
const VectorVariableTestValue & _test
the current test function
Definition: VectorKernel.h:52
virtual void computeResidual() override
Compute this VectorKernel&#39;s contribution to the residual.
Definition: VectorKernel.C:48
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:172
virtual void precalculateOffDiagJacobian(unsigned int)
Definition: KernelBase.h:123
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:159
VectorValue< Real > RealVectorValue
Definition: Assembly.h:31
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
SystemBase & _sys
Reference to the EquationSystem object.
Definition: KernelBase.h:133
THREAD_ID _tid
The thread ID for this kernel.
Definition: KernelBase.h:136
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
VectorKernel(const InputParameters &parameters)
Definition: VectorKernel.C:29
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.
virtual void precalculateJacobian()
Definition: KernelBase.h:122
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
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
MooseVariableFE< RealVectorValue > * mooseVariable() const
Get the variable that this object is using.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
This is the common base class for the two main kernel types implemented in MOOSE, EigenKernel and Ker...
Definition: KernelBase.h:44
virtual void computeJacobian() override
Compute this VectorKernel&#39;s contribution to the diagonal Jacobian entries.
Definition: VectorKernel.C:73
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:177
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
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
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
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
unsigned int _j
current index for the shape function
Definition: KernelBase.h:168
virtual void precalculateResidual()
Following methods are used for Kernels that need to perform a per-element calculation.
Definition: KernelBase.h:121
VectorMooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: VectorKernel.h:49
const VectorVariablePhiValue & _phi
the current shape functions
Definition: VectorKernel.h:58
InputParameters validParams< VectorKernel >()
Definition: VectorKernel.C:22
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
Interface for objects that need to get values of MooseVariables.
Class for scalar variables (they are different).
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
Definition: Moose.h:112
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
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:149
InputParameters validParams< KernelBase >()
Definition: KernelBase.C:22
virtual void computeOffDiagJacobian(MooseVariableFEBase &jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: VectorKernel.C:101
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: VectorKernel.C:119
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:150