www.mooseframework.org
KernelGrad.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 "KernelGrad.h"
11 #include "Assembly.h"
12 #include "libmesh/quadrature.h"
13 
14 template <>
17 {
19  return params;
20 }
21 
22 KernelGrad::KernelGrad(const InputParameters & parameters) : Kernel(parameters) {}
23 
24 void
26 {
27  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
28  _local_re.resize(re.size());
29  _local_re.zero();
30 
31  const unsigned int n_test = _test.size();
32  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
33  {
34  RealGradient value = precomputeQpResidual() * _JxW[_qp] * _coord[_qp];
35  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
37  }
38 
39  re += _local_re;
40 
41  if (_has_save_in)
42  {
43  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
44  for (const auto & var : _save_in)
45  var->sys().solution().add_vector(_local_re, var->dofIndices());
46  }
47 }
48 
49 void
51 {
52  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
53  _local_ke.resize(ke.m(), ke.n());
54  _local_ke.zero();
55 
56  const unsigned int n_test = _test.size();
57  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
58  for (_j = 0; _j < _phi.size(); _j++)
59  {
60  RealGradient value = precomputeQpJacobian() * _JxW[_qp] * _coord[_qp];
61  for (_i = 0; _i < n_test; _i++) // target for auto vectorization
62  _local_ke(_i, _j) += value * _grad_test[_i][_qp];
63  }
64 
65  ke += _local_ke;
66 
68  {
69  const unsigned int rows = ke.m();
70  DenseVector<Number> diag(rows);
71  for (unsigned int i = 0; i < rows; i++) // target for auto vectorization
72  diag(i) = _local_ke(i, i);
73 
74  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
75  for (const auto & var : _diag_save_in)
76  var->sys().solution().add_vector(diag, var->dofIndices());
77  }
78 }
79 
80 void
82 {
83  size_t jvar_num = jvar.number();
84  if (jvar_num == _var.number())
86  else
87  {
88  DenseMatrix<Number> & Ke = _assembly.jacobianBlock(_var.number(), jvar_num);
89 
90  for (_j = 0; _j < jvar.phiSize(); _j++)
91  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
92  for (_i = 0; _i < _test.size(); _i++)
93  Ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
94  }
95 }
96 
97 Real
99 {
100  return 0.0;
101 }
102 
103 RealGradient
105 {
106  return RealGradient(0.0);
107 }
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:54
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:172
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:159
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
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: KernelGrad.C:25
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...
virtual Real computeQpResidual() override
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
Definition: KernelGrad.C:98
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 size_t phiSize() const =0
Return phi size.
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:139
const VariableTestValue & _test
the current test function
Definition: Kernel.h:57
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: KernelGrad.C:50
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:177
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
virtual RealGradient precomputeQpJacobian()
Called before forming the jacobian for an element.
Definition: KernelGrad.C:104
unsigned int _i
current index for the test function
Definition: KernelBase.h:165
virtual const OutputTools< Real >::VariableValue & value()
The value of the variable this object is operating on.
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:162
unsigned int _j
current index for the shape function
Definition: KernelBase.h:168
InputParameters validParams< Kernel >()
Definition: Kernel.C:24
const VariableTestGradient & _grad_test
gradient of the test function
Definition: Kernel.h:60
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
Definition: Kernel.h:20
virtual void computeOffDiagJacobian(MooseVariableFEBase &jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: KernelGrad.C:81
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
InputParameters validParams< KernelGrad >()
Definition: KernelGrad.C:16
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:63
KernelGrad(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: KernelGrad.C:22
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:150