https://mooseframework.inl.gov
KernelGrad.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 "KernelGrad.h"
11 #include "Assembly.h"
12 #include "SystemBase.h"
13 #include "libmesh/quadrature.h"
14 
17 {
18  return Kernel::validParams();
19 }
20 
21 KernelGrad::KernelGrad(const InputParameters & parameters) : Kernel(parameters) {}
22 
23 void
25 {
27  const unsigned int n_test = _test.size();
28  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
29  {
31  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
33  }
35 
36  if (_has_save_in)
37  for (const auto & var : _save_in)
38  var->sys().solution().add_vector(_local_re, var->dofIndices());
39 }
40 
41 void
43 {
45 
46  const unsigned int n_test = _test.size();
47  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
48  for (_j = 0; _j < _phi.size(); _j++)
49  {
51  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
52  _local_ke(_i, _j) += value * _grad_test[_i][_qp];
53  }
54 
56 
58  {
59  const unsigned int rows = _local_ke.m();
60  DenseVector<Number> diag(rows);
61  for (unsigned int i = 0; i < rows; i++) // target for auto vectorization
62  diag(i) = _local_ke(i, i);
63 
64  for (const auto & var : _diag_save_in)
65  var->sys().solution().add_vector(diag, var->dofIndices());
66  }
67 }
68 
69 void
70 KernelGrad::computeOffDiagJacobian(const unsigned int jvar_num)
71 {
72  const auto & jvar = getVariable(jvar_num);
73 
74  if (jvar_num == _var.number())
76  else
77  {
78  prepareMatrixTag(_assembly, _var.number(), jvar_num);
79 
80  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
81  // on the displaced mesh
82  auto phi_size = jvar.dofIndices().size();
83 
84  for (_j = 0; _j < phi_size; _j++)
85  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
86  for (_i = 0; _i < _test.size(); _i++)
88 
90  }
91 }
92 
93 Real
95 {
96  return 0.0;
97 }
98 
101 {
102  return RealGradient(0.0);
103 }
static InputParameters validParams()
Definition: Kernel.C:24
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:72
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:65
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:52
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int)
For coupling standard variables.
Definition: Kernel.h:56
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: KernelGrad.C:24
virtual RealGradient precomputeQpResidual()=0
Called before forming the residual for an element.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int m() const
virtual Real computeQpResidual() override
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
Definition: KernelGrad.C:94
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.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
const VariableTestValue & _test
the current test function
Definition: Kernel.h:75
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: KernelGrad.C:42
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:70
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
virtual RealGradient precomputeQpJacobian()
Called before forming the jacobian for an element.
Definition: KernelGrad.C:100
unsigned int _i
current index for the test function
Definition: KernelBase.h:58
virtual const OutputTools< Real >::VariableValue & value()
The value of the variable this object is operating on.
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.
static InputParameters validParams()
Definition: KernelGrad.C:16
unsigned int _j
current index for the shape function
Definition: KernelBase.h:61
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const VariableTestGradient & _grad_test
gradient of the test function
Definition: Kernel.h:78
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
Definition: Kernel.h:15
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:81
KernelGrad(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: KernelGrad.C:21
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:43
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: KernelGrad.C:70