www.mooseframework.org
MatDiffusionBase.h
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 #pragma once
11 
12 #include "Kernel.h"
13 #include "JvarMapInterface.h"
15 
25 template <typename T>
26 class MatDiffusionBase : public DerivativeMaterialInterface<JvarMapKernelInterface<Kernel>>
27 {
28 public:
30 
31  virtual void initialSetup();
32 
35 
36 protected:
37  virtual Real computeQpResidual();
38  virtual Real computeQpJacobian();
39  virtual Real computeQpOffDiagJacobian(unsigned int jvar);
40  virtual Real computeQpCJacobian();
41 
44 
47 
49  std::vector<const MaterialProperty<T> *> _dDdarg;
50 
52  const bool _is_coupled;
53 
55  unsigned int _conc_var;
56 
59 };
60 
61 template <typename T>
64 {
66  params.addParam<MaterialPropertyName>("D_name", "D", "The name of the diffusivity");
67  params.addCoupledVar("args", "Vector of arguments of the diffusivity");
68  params.addCoupledVar("conc",
69  "Coupled concentration variable for kernel to operate on; if this "
70  "is not specified, the kernel's nonlinear variable will be used as "
71  "usual");
72  return params;
73 }
74 
75 template <typename T>
78  _D(getMaterialProperty<T>("D_name")),
79  _dDdc(getMaterialPropertyDerivative<T>("D_name", _var.name())),
80  _dDdarg(_coupled_moose_vars.size()),
81  _is_coupled(isCoupled("conc")),
82  _conc_var(_is_coupled ? coupled("conc") : _var.number()),
83  _grad_conc(_is_coupled ? coupledGradient("conc") : _grad_u)
84 {
85  // fetch derivatives
86  for (unsigned int i = 0; i < _dDdarg.size(); ++i)
87  _dDdarg[i] = &getMaterialPropertyDerivative<T>("D_name", _coupled_moose_vars[i]->name());
88 }
89 
90 template <typename T>
91 void
93 {
94  validateNonlinearCoupling<Real>("D_name");
95 }
96 
97 template <typename T>
98 Real
100 {
101  return _D[_qp] * _grad_conc[_qp] * _grad_test[_i][_qp];
102 }
103 
104 template <typename T>
105 Real
107 {
108  Real sum = _phi[_j][_qp] * _dDdc[_qp] * _grad_conc[_qp] * _grad_test[_i][_qp];
109  if (!_is_coupled)
110  sum += computeQpCJacobian();
111 
112  return sum;
113 }
114 
115 template <typename T>
116 Real
118 {
119  // get the coupled variable jvar is referring to
120  const unsigned int cvar = mapJvarToCvar(jvar);
121 
122  Real sum = (*_dDdarg[cvar])[_qp] * _phi[_j][_qp] * _grad_conc[_qp] * _grad_test[_i][_qp];
123  if (_conc_var == jvar)
124  sum += computeQpCJacobian();
125 
126  return sum;
127 }
128 
129 template <typename T>
130 Real
132 {
133  return _D[_qp] * _grad_phi[_j][_qp] * _grad_test[_i][_qp];
134 }
OutputTools< Real >::VariableGradient VariableGradient
Definition: MooseTypes.h:198
const MaterialProperty< T > & _D
diffusion coefficient
const MaterialProperty< T > & _dDdc
diffusion coefficient derivative w.r.t. the kernel variable
const VariableGradient & _grad_conc
Gradient of the concentration.
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
virtual Real computeQpCJacobian()
unsigned int _conc_var
int label for the Concentration
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::vector< const MaterialProperty< T > * > _dDdarg
diffusion coefficient derivatives w.r.t. coupled variables
virtual Real computeQpResidual()
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
Interface class ("Veneer") for Kernel to provide a mapping from &#39;jvar&#39; in computeQpOffDiagJacobian in...
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
MatDiffusionBase(const InputParameters &parameters)
static InputParameters validParams()
in class templates this function has to be a static member
This class template implements a diffusion kernel with a mobility that can vary spatially and can dep...
void addCoupledVar(const std::string &name, const std::string &doc_string)
This method adds a coupled variable name pair.
InputParameters validParams< Kernel >()
Definition: Kernel.C:24
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Interface class ("Veneer") to provide generator methods for derivative material property names...
std::vector< MooseVariableFEBase * > _coupled_moose_vars
Vector of all coupled variables.
Definition: Coupleable.h:740
Concrete definition of a parameter value for a specified type.
Definition: Kernel.h:20
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:59
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
const bool _is_coupled
is the kernel used in a coupled form?
virtual void initialSetup()
Gets called at the beginning of the simulation before this object is asked to do its job...