https://mooseframework.inl.gov
ArrayReaction.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 "ArrayReaction.h"
11 
13 
16 {
18  params.addRequiredParam<MaterialPropertyName>(
19  "reaction_coefficient",
20  "The name of the reactivity, can be scalar, vector, or matrix material property.");
21  params.addClassDescription("The array reaction operator with the weak "
22  "form of $(\\psi_i, u_h)$.");
23  return params;
24 }
25 
27  : ArrayKernel(parameters),
28  _r(hasMaterialProperty<Real>("reaction_coefficient")
29  ? &getMaterialProperty<Real>("reaction_coefficient")
30  : nullptr),
31  _r_array(hasMaterialProperty<RealEigenVector>("reaction_coefficient")
32  ? &getMaterialProperty<RealEigenVector>("reaction_coefficient")
33  : nullptr),
34  _r_2d_array(hasMaterialProperty<RealEigenMatrix>("reaction_coefficient")
35  ? &getMaterialProperty<RealEigenMatrix>("reaction_coefficient")
36  : nullptr)
37 {
38  if (!_r && !_r_array && !_r_2d_array)
39  {
40  MaterialPropertyName mat = getParam<MaterialPropertyName>("reaction_coefficient");
41  mooseError("Property " + mat + " is of unsupported type for ArrayReaction");
42  }
43 }
44 
45 void
47 {
48 
49  if (_r)
50  residual = (*_r)[_qp] * _u[_qp] * _test[_i][_qp];
51 
52  else if (_r_array)
53  {
54  mooseAssert((*_r_array)[_qp].size() == _var.count(),
55  "reaction_coefficient size is inconsistent with the number of components of array "
56  "variable");
57  // WARNING: use noalias() syntax with caution. See ArrayDiffusion.C for more details.
58  residual.noalias() = (*_r_array)[_qp].asDiagonal() * _u[_qp] * _test[_i][_qp];
59  }
60 
61  else
62  {
63  mooseAssert((*_r_2d_array)[_qp].cols() == _var.count(),
64  "reaction_coefficient size is inconsistent with the number of components of array "
65  "variable");
66  mooseAssert((*_r_2d_array)[_qp].rows() == _var.count(),
67  "reaction_coefficient size is inconsistent with the number of components of array "
68  "variable");
69  // WARNING: use noalias() syntax with caution. See ArrayDiffusion.C for more details.
70  residual.noalias() = (*_r_2d_array)[_qp] * _u[_qp] * _test[_i][_qp];
71  }
72 }
73 
76 {
77  if (_r)
78  return RealEigenVector::Constant(_var.count(), _phi[_j][_qp] * _test[_i][_qp] * (*_r)[_qp]);
79  else if (_r_array)
80  return _phi[_j][_qp] * _test[_i][_qp] * (*_r_array)[_qp];
81 
82  else
83  return _phi[_j][_qp] * _test[_i][_qp] * (*_r_2d_array)[_qp].diagonal();
84 }
85 
88 {
89  if (jvar.number() == _var.number() && _r_2d_array)
90  return _phi[_j][_qp] * _test[_i][_qp] * (*_r_2d_array)[_qp];
91  else
93 }
virtual RealEigenVector computeQpJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian at the current quadrature point...
Definition: ArrayReaction.C:75
ArrayReaction(const InputParameters &parameters)
Definition: ArrayReaction.C:26
unsigned int number() const
Get variable number coming from libMesh.
unsigned int count() const
Get the number of components Note: For standard and vector variables, the number is one...
const MaterialProperty< RealEigenMatrix > *const _r_2d_array
matrix diffusion coefficient
Definition: ArrayReaction.h:31
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
This class provides an interface for common operations on field variables of both FE and FV types wit...
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
static InputParameters validParams()
Definition: ArrayKernel.C:23
virtual RealEigenMatrix computeQpOffDiagJacobian(const MooseVariableFEBase &jvar) override
This is the virtual that derived classes should override for computing a full Jacobian component...
Definition: ArrayReaction.C:87
virtual void computeQpResidual(RealEigenVector &residual) override
Compute this Kernel&#39;s contribution to the residual at the current quadrature point, to be filled in residual.
Definition: ArrayReaction.C:46
const ArrayVariableTestValue & _test
the current test function
Definition: ArrayKernel.h:88
unsigned int _i
current index for the test function
Definition: KernelBase.h:58
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
Definition: MooseTypes.h:149
const MaterialProperty< RealEigenVector > *const _r_array
array diffusion coefficient
Definition: ArrayReaction.h:29
virtual RealEigenMatrix computeQpOffDiagJacobian(const MooseVariableFEBase &jvar)
This is the virtual that derived classes should override for computing a full Jacobian component...
Definition: ArrayKernel.C:215
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 MaterialProperty< Real > *const _r
scalar diffusion coefficient
Definition: ArrayReaction.h:27
const ArrayVariablePhiValue & _phi
the current shape functions
Definition: ArrayKernel.h:95
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
const ArrayVariableValue & _u
Holds the solution at current quadrature points.
Definition: ArrayKernel.h:101
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
ArrayMooseVariable & _var
This is an array kernel so we cast to a ArrayMooseVariable.
Definition: ArrayKernel.h:85
registerMooseObject("MooseApp", ArrayReaction)
static InputParameters validParams()
Definition: ArrayReaction.C:15
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:43