Loading [MathJax]/extensions/tex2jax.js
https://mooseframework.inl.gov
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
VectorKernel.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 
22 {
24  params.registerBase("VectorKernel");
25  return params;
26 }
27 
29  : KernelBase(parameters),
31  false,
32  "variable",
35  _var(*mooseVariable()),
36  _test(_var.phi()),
37  _grad_test(_var.gradPhi()),
38  _phi(_assembly.phi(_var)),
39  _grad_phi(_assembly.gradPhi(_var)),
40  _u(_is_implicit ? _var.sln() : _var.slnOld()),
41  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
42 {
44 }
45 
46 void
48 {
51  for (_i = 0; _i < _test.size(); _i++)
52  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
53  {
54  Real residual = _JxW[_qp] * _coord[_qp] * computeQpResidual();
55  _local_re(_i) += residual;
56  }
58 
59  if (_has_save_in)
60  for (const auto & var : _save_in)
61  var->sys().solution().add_vector(_local_re, var->dofIndices());
62 }
63 
64 void
66 {
69  for (_i = 0; _i < _test.size(); _i++)
70  for (_j = 0; _j < _phi.size(); _j++)
71  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
74 
76  {
77  unsigned int rows = _local_ke.m();
78  DenseVector<Number> diag(rows);
79  for (unsigned int i = 0; i < rows; i++)
80  diag(i) = _local_ke(i, i);
81 
82  for (const auto & var : _diag_save_in)
83  var->sys().solution().add_vector(diag, var->dofIndices());
84  }
85 }
86 
87 void
88 VectorKernel::computeOffDiagJacobian(const unsigned int jvar_num)
89 {
90  if (jvar_num == _var.number())
92  else
93  {
94  const auto & jvar = getVariable(jvar_num);
95 
96  prepareMatrixTag(_assembly, _var.number(), jvar_num);
97 
98  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
99  // on the displaced mesh
100  const auto phi_size = jvar.dofIndices().size();
101 
102  precalculateOffDiagJacobian(jvar_num);
103  for (_i = 0; _i < _test.size(); _i++)
104  for (_j = 0; _j < phi_size; _j++)
105  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
106  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
107 
109  }
110 }
111 
112 void
114 {
117 
118  for (_i = 0; _i < _test.size(); _i++)
119  for (_j = 0; _j < jv.order(); _j++)
120  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
123 }
VarFieldType
Definition: MooseTypes.h:721
const VectorVariableTestValue & _test
the current test function
Definition: VectorKernel.h:65
virtual void computeResidual() override
Compute this VectorKernel&#39;s contribution to the residual.
Definition: VectorKernel.C:47
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:65
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:52
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
virtual void precalculateResidual()
virtual Real computeQpOffDiagJacobian(unsigned int)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
Definition: VectorKernel.h:54
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool computingScalingJacobian() const
Whether we are computing an initial Jacobian for automatic variable scaling.
Definition: SystemBase.C:1519
unsigned int m() const
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: VectorKernel.h:48
VectorKernel(const InputParameters &parameters)
Definition: VectorKernel.C:28
THREAD_ID _tid
The thread ID for this kernel.
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:69
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
SystemBase & _sys
Reference to the EquationSystem object.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
MooseVariableFE< RealVectorValue > * mooseVariable() const
Return the MooseVariableFE object that this interface acts on.
virtual Real computeQpOffDiagJacobianScalar(unsigned int)
For coupling scalar variables.
Definition: VectorKernel.h:59
This is the common base class for the three main kernel types implemented in MOOSE, Kernel, VectorKernel and ArrayKernel.
Definition: KernelBase.h:23
virtual void computeJacobian() override
Compute this VectorKernel&#39;s contribution to the diagonal Jacobian entries.
Definition: VectorKernel.C:65
virtual void precalculateJacobian()
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:70
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:714
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:49
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:64
unsigned int _i
current index for the test function
Definition: KernelBase.h:58
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:55
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
unsigned int _j
current index for the shape function
Definition: KernelBase.h:61
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:144
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
const VectorVariablePhiValue & _phi
the current shape functions
Definition: VectorKernel.h:71
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static InputParameters validParams()
Definition: VectorKernel.C:21
DenseVector< Number > _local_re
Holds local 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 void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: VectorKernel.C:88
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
static InputParameters validParams()
Definition: KernelBase.C:21
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: VectorKernel.C:113
const VectorMooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: VectorKernel.h:62
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:43