https://mooseframework.inl.gov
Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
PorousFlowFluidStateBaseMaterialTempl< is_ad > Class Template Referenceabstract

Fluid state base class using a persistent set of primary variables for multiphase, single and multicomponent cases. More...

#include <PorousFlowFluidStateBaseMaterial.h>

Inheritance diagram for PorousFlowFluidStateBaseMaterialTempl< is_ad >:
[legend]

Public Types

typedef DerivativeMaterialPropertyNameInterface::SymbolName SymbolName
 

Public Member Functions

 PorousFlowFluidStateBaseMaterialTempl (const InputParameters &parameters)
 
const GenericMaterialProperty< U, is_ad > & getDefaultMaterialProperty (const std::string &name)
 
const GenericMaterialProperty< U, is_ad > & getDefaultMaterialPropertyByName (const std::string &name)
 
void validateDerivativeMaterialPropertyBase (const std::string &base)
 
const MaterialPropertyName derivativePropertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName derivativePropertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName derivativePropertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName derivativePropertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
void validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)
 
void validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)
 
void validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
const MaterialPropertyName propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 
const MaterialPropertyName propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

virtual void initQpStatefulProperties () override
 
virtual void computeQpProperties () override
 
virtual void setMaterialVectorSize () const
 Size material property vectors and initialise with zeros. More...
 
virtual void thermophysicalProperties ()=0
 Calculates all required thermophysical properties and derivatives for each phase and fluid component. More...
 
GenericReal< is_ad > genericValue (const ADReal &value)
 Templated method that takes an ADReal value, and returns the raw value when is_ad = false, and the AD value when is_ad = true. More...
 
template<>
GenericReal< false > genericValue (const ADReal &value)
 
template<>
GenericReal< true > genericValue (const ADReal &value)
 

Protected Attributes

const Real _T_c2k
 Conversion from degrees Celsius to degrees Kelvin. More...
 
bool _is_initqp
 Flag to indicate whether stateful properties should be computed. More...
 
std::vector< FluidStateProperties_fsp
 FluidStateProperties data for each fluid component in each fluid phase. More...
 
const PorousFlowCapillaryPressure_pc
 Capillary pressure UserObject. More...
 
const std::string _sfx
 Suffix to append to material property names (either qp or nodal as required) More...
 
GenericMaterialProperty< std::vector< std::vector< Real > >, is_ad > & _mass_frac
 Mass fraction matrix (indexing is fluid component in fluid phase) More...
 
GenericMaterialProperty< std::vector< std::vector< RealGradient > >, is_ad > *const _grad_mass_frac_qp
 Gradient of the mass fraction matrix (only defined at the qps) More...
 
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > *const _dmass_frac_dvar
 Derivative of the mass fraction matrix with respect to the Porous Flow variables. More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_density
 Fluid density of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_density_dvar
 Derivative of the fluid density for each phase wrt PorousFlow variables. More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_viscosity
 Viscosity of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_viscosity_dvar
 Derivative of the fluid viscosity for each phase wrt PorousFlow variables. More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_enthalpy
 Enthalpy of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_enthalpy_dvar
 Derivative of the fluid enthalpy for each phase wrt PorousFlow variables. More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_internal_energy
 Internal energy of each phase. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_internal_energy_dvar
 Derivative of the fluid internal energy for each phase wrt PorousFlow variables. More...
 
 usingPorousFlowVariableBaseMembers
 
const unsigned int _num_phases
 Number of phases. More...
 
const unsigned int _num_components
 Number of components. More...
 
const unsigned int _num_pf_vars
 Number of PorousFlow variables. More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _porepressure
 Computed nodal or quadpoint values of porepressure of the phases. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dporepressure_dvar
 d(porepressure)/d(PorousFlow variable) More...
 
GenericMaterialProperty< std::vector< RealGradient >, is_ad > *const _gradp_qp
 Grad(p) at the quadpoints. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dgradp_qp_dgradv
 d(grad porepressure)/d(grad PorousFlow variable) at the quadpoints More...
 
MaterialProperty< std::vector< std::vector< RealGradient > > > *const _dgradp_qp_dv
 d(grad porepressure)/d(PorousFlow variable) at the quadpoints More...
 
GenericMaterialProperty< std::vector< Real >, is_ad > & _saturation
 Computed nodal or qp saturation of the phases. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dsaturation_dvar
 d(saturation)/d(PorousFlow variable) More...
 
GenericMaterialProperty< std::vector< RealGradient >, is_ad > *const _grads_qp
 Grad(s) at the quadpoints. More...
 
MaterialProperty< std::vector< std::vector< Real > > > *const _dgrads_qp_dgradv
 d(grad saturation)/d(grad PorousFlow variable) at the quadpoints More...
 
MaterialProperty< std::vector< std::vector< RealGradient > > > *const _dgrads_qp_dv
 d(grad saturation)/d(PorousFlow variable) at the quadpoints More...
 

Detailed Description

template<bool is_ad>
class PorousFlowFluidStateBaseMaterialTempl< is_ad >

Fluid state base class using a persistent set of primary variables for multiphase, single and multicomponent cases.

Definition at line 22 of file PorousFlowFluidStateBaseMaterial.h.

Constructor & Destructor Documentation

◆ PorousFlowFluidStateBaseMaterialTempl()

Definition at line 30 of file PorousFlowFluidStateBaseMaterial.C.

33  _T_c2k(this->template getParam<MooseEnum>("temperature_unit") == 0 ? 0.0 : 273.15),
34  _is_initqp(false),
35  _pc(this->template getUserObject<PorousFlowCapillaryPressure>("capillary_pressure")),
36  _sfx(_nodal_material ? "_nodal" : "_qp"),
37  _mass_frac(this->template declareGenericProperty<std::vector<std::vector<Real>>, is_ad>(
38  "PorousFlow_mass_frac" + _sfx)),
40  _nodal_material
41  ? nullptr
42  : &this->template declareGenericProperty<std::vector<std::vector<RealGradient>>, is_ad>(
43  "PorousFlow_grad_mass_frac" + _sfx)),
45  is_ad ? nullptr
46  : &this->template declareProperty<std::vector<std::vector<std::vector<Real>>>>(
47  "dPorousFlow_mass_frac" + _sfx + "_dvar")),
48  _fluid_density(this->template declareGenericProperty<std::vector<Real>, is_ad>(
49  "PorousFlow_fluid_phase_density" + _sfx)),
50  _dfluid_density_dvar(is_ad ? nullptr
51  : &this->template declareProperty<std::vector<std::vector<Real>>>(
52  "dPorousFlow_fluid_phase_density" + _sfx + "_dvar")),
53  _fluid_viscosity(this->template declareGenericProperty<std::vector<Real>, is_ad>(
54  "PorousFlow_viscosity" + _sfx)),
55  _dfluid_viscosity_dvar(is_ad ? nullptr
56  : &this->template declareProperty<std::vector<std::vector<Real>>>(
57  "dPorousFlow_viscosity" + _sfx + "_dvar")),
58  _fluid_enthalpy(this->template declareGenericProperty<std::vector<Real>, is_ad>(
59  "PorousFlow_fluid_phase_enthalpy" + _sfx)),
60  _dfluid_enthalpy_dvar(is_ad ? nullptr
61  : &this->template declareProperty<std::vector<std::vector<Real>>>(
62  "dPorousFlow_fluid_phase_enthalpy" + _sfx + "_dvar")),
63  _fluid_internal_energy(this->template declareGenericProperty<std::vector<Real>, is_ad>(
64  "PorousFlow_fluid_phase_internal_energy" + _sfx)),
66  is_ad ? nullptr
67  : &this->template declareProperty<std::vector<std::vector<Real>>>(
68  "dPorousFlow_fluid_phase_internal_energy" + _sfx + "_dvar"))
69 {
70  // Set the size of the FluidStateProperties vector
72 }
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_internal_energy_dvar
Derivative of the fluid internal energy for each phase wrt PorousFlow variables.
const PorousFlowCapillaryPressure & _pc
Capillary pressure UserObject.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_viscosity
Viscosity of each phase.
const Real _T_c2k
Conversion from degrees Celsius to degrees Kelvin.
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_viscosity_dvar
Derivative of the fluid viscosity for each phase wrt PorousFlow variables.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_internal_energy
Internal energy of each phase.
GenericMaterialProperty< std::vector< std::vector< Real > >, is_ad > & _mass_frac
Mass fraction matrix (indexing is fluid component in fluid phase)
GenericMaterialProperty< std::vector< std::vector< RealGradient > >, is_ad > *const _grad_mass_frac_qp
Gradient of the mass fraction matrix (only defined at the qps)
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > *const _dmass_frac_dvar
Derivative of the mass fraction matrix with respect to the Porous Flow variables. ...
Base class for thermophysical variable materials, which assemble materials for primary variables such...
const std::string _sfx
Suffix to append to material property names (either qp or nodal as required)
AD data structure to pass calculated thermophysical properties.
const unsigned int _num_components
Number of components.
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_enthalpy_dvar
Derivative of the fluid enthalpy for each phase wrt PorousFlow variables.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_density
Fluid density of each phase.
bool _is_initqp
Flag to indicate whether stateful properties should be computed.
const unsigned int _num_phases
Number of phases.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_enthalpy
Enthalpy of each phase.
std::vector< FluidStateProperties > _fsp
FluidStateProperties data for each fluid component in each fluid phase.
MaterialProperty< std::vector< std::vector< Real > > > *const _dfluid_density_dvar
Derivative of the fluid density for each phase wrt PorousFlow variables.

Member Function Documentation

◆ computeQpProperties()

template<bool is_ad>
void PorousFlowFluidStateBaseMaterialTempl< is_ad >::computeQpProperties ( )
overrideprotectedvirtual

Reimplemented from PorousFlowVariableBaseTempl< is_ad >.

Reimplemented in PorousFlowFluidStateTempl< is_ad >, and PorousFlowFluidStateSingleComponentTempl< is_ad >.

Definition at line 121 of file PorousFlowFluidStateBaseMaterial.C.

Referenced by PorousFlowFluidStateSingleComponentTempl< is_ad >::computeQpProperties(), and PorousFlowFluidStateTempl< is_ad >::computeQpProperties().

122 {
123  _is_initqp = false;
124  // Size pressure and saturation and prepare the derivative vectors
126 
127  // Set the size of all other vectors
129 
130  // Calculate all required thermophysical properties
132 
133  for (unsigned int ph = 0; ph < _num_phases; ++ph)
134  {
135  _saturation[_qp][ph] = genericValue(_fsp[ph].saturation);
136  _porepressure[_qp][ph] = genericValue(_fsp[ph].pressure);
137  _fluid_density[_qp][ph] = genericValue(_fsp[ph].density);
138  _fluid_viscosity[_qp][ph] = genericValue(_fsp[ph].viscosity);
139  _fluid_enthalpy[_qp][ph] = genericValue(_fsp[ph].enthalpy);
141 
142  for (unsigned int comp = 0; comp < _num_components; ++comp)
143  _mass_frac[_qp][ph][comp] = genericValue(_fsp[ph].mass_fraction[comp]);
144  }
145 }
virtual void computeQpProperties() override
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_viscosity
Viscosity of each phase.
GenericMaterialProperty< std::vector< Real >, is_ad > & _porepressure
Computed nodal or quadpoint values of porepressure of the phases.
static const std::string density
Definition: NS.h:33
virtual void thermophysicalProperties()=0
Calculates all required thermophysical properties and derivatives for each phase and fluid component...
GenericReal< is_ad > genericValue(const ADReal &value)
Templated method that takes an ADReal value, and returns the raw value when is_ad = false...
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_internal_energy
Internal energy of each phase.
GenericMaterialProperty< std::vector< std::vector< Real > >, is_ad > & _mass_frac
Mass fraction matrix (indexing is fluid component in fluid phase)
virtual void setMaterialVectorSize() const
Size material property vectors and initialise with zeros.
const unsigned int _num_components
Number of components.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_density
Fluid density of each phase.
bool _is_initqp
Flag to indicate whether stateful properties should be computed.
const unsigned int _num_phases
Number of phases.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_enthalpy
Enthalpy of each phase.
static const std::string pressure
Definition: NS.h:56
std::vector< FluidStateProperties > _fsp
FluidStateProperties data for each fluid component in each fluid phase.
static const std::string internal_energy
Definition: NS.h:61
GenericMaterialProperty< std::vector< Real >, is_ad > & _saturation
Computed nodal or qp saturation of the phases.

◆ genericValue() [1/3]

template<bool is_ad>
GenericReal<is_ad> PorousFlowFluidStateBaseMaterialTempl< is_ad >::genericValue ( const ADReal value)
protected

Templated method that takes an ADReal value, and returns the raw value when is_ad = false, and the AD value when is_ad = true.

Parameters
valueAD value
Returns
raw or AD value depending on is_ad

◆ genericValue() [2/3]

template<>
GenericReal< false > PorousFlowFluidStateBaseMaterialTempl< false >::genericValue ( const ADReal value)
protected

Definition at line 76 of file PorousFlowFluidStateBaseMaterial.C.

77 {
78  return MetaPhysicL::raw_value(value);
79 }
auto raw_value(const Eigen::Map< T > &in)

◆ genericValue() [3/3]

template<>
GenericReal< true > PorousFlowFluidStateBaseMaterialTempl< true >::genericValue ( const ADReal value)
protected

Definition at line 83 of file PorousFlowFluidStateBaseMaterial.C.

84 {
85  return value;
86 }
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)

◆ initQpStatefulProperties()

template<bool is_ad>
void PorousFlowFluidStateBaseMaterialTempl< is_ad >::initQpStatefulProperties ( )
overrideprotectedvirtual

Reimplemented from PorousFlowVariableBaseTempl< is_ad >.

Reimplemented in PorousFlowFluidStateTempl< is_ad >, and PorousFlowFluidStateSingleComponentTempl< is_ad >.

Definition at line 90 of file PorousFlowFluidStateBaseMaterial.C.

Referenced by PorousFlowFluidStateSingleComponentTempl< is_ad >::initQpStatefulProperties(), and PorousFlowFluidStateTempl< is_ad >::initQpStatefulProperties().

91 {
92  _is_initqp = true;
93  // Set the size of pressure and saturation vectors
95 
96  // Set the size of all other vectors
98 
100 
101  // Set the initial values of the properties at the nodes.
102  // Note: not required for qp materials as no old values at the qps are requested
103  // unless the material is AD (for the FV case there are no nodal materials)
104  if (_nodal_material || is_ad)
105  for (unsigned int ph = 0; ph < _num_phases; ++ph)
106  {
107  _saturation[_qp][ph] = genericValue(_fsp[ph].saturation);
108  _porepressure[_qp][ph] = genericValue(_fsp[ph].pressure);
109  _fluid_density[_qp][ph] = genericValue(_fsp[ph].density);
110  _fluid_viscosity[_qp][ph] = genericValue(_fsp[ph].viscosity);
111  _fluid_enthalpy[_qp][ph] = genericValue(_fsp[ph].enthalpy);
113 
114  for (unsigned int comp = 0; comp < _num_components; ++comp)
115  _mass_frac[_qp][ph][comp] = genericValue(_fsp[ph].mass_fraction[comp]);
116  }
117 }
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_viscosity
Viscosity of each phase.
GenericMaterialProperty< std::vector< Real >, is_ad > & _porepressure
Computed nodal or quadpoint values of porepressure of the phases.
static const std::string density
Definition: NS.h:33
virtual void thermophysicalProperties()=0
Calculates all required thermophysical properties and derivatives for each phase and fluid component...
GenericReal< is_ad > genericValue(const ADReal &value)
Templated method that takes an ADReal value, and returns the raw value when is_ad = false...
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_internal_energy
Internal energy of each phase.
GenericMaterialProperty< std::vector< std::vector< Real > >, is_ad > & _mass_frac
Mass fraction matrix (indexing is fluid component in fluid phase)
virtual void setMaterialVectorSize() const
Size material property vectors and initialise with zeros.
const unsigned int _num_components
Number of components.
virtual void initQpStatefulProperties() override
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_density
Fluid density of each phase.
bool _is_initqp
Flag to indicate whether stateful properties should be computed.
const unsigned int _num_phases
Number of phases.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_enthalpy
Enthalpy of each phase.
static const std::string pressure
Definition: NS.h:56
std::vector< FluidStateProperties > _fsp
FluidStateProperties data for each fluid component in each fluid phase.
static const std::string internal_energy
Definition: NS.h:61
GenericMaterialProperty< std::vector< Real >, is_ad > & _saturation
Computed nodal or qp saturation of the phases.

◆ setMaterialVectorSize()

template<bool is_ad>
void PorousFlowFluidStateBaseMaterialTempl< is_ad >::setMaterialVectorSize ( ) const
protectedvirtual

Size material property vectors and initialise with zeros.

Reimplemented in PorousFlowFluidStateSingleComponentTempl< is_ad >.

Definition at line 149 of file PorousFlowFluidStateBaseMaterial.C.

Referenced by PorousFlowFluidStateSingleComponentTempl< is_ad >::setMaterialVectorSize().

150 {
151  _fluid_density[_qp].resize(_num_phases, 0.0);
152  _fluid_viscosity[_qp].resize(_num_phases, 0.0);
153  _fluid_enthalpy[_qp].resize(_num_phases, 0.0);
154  _fluid_internal_energy[_qp].resize(_num_phases, 0.0);
155  _mass_frac[_qp].resize(_num_phases);
156 
157  for (unsigned int ph = 0; ph < _num_phases; ++ph)
158  _mass_frac[_qp][ph].resize(_num_components);
159 
160  // Derivatives and gradients are not required in initQpStatefulProperties
161  if (!_is_initqp)
162  {
163  // The gradient of the mass fractions are needed for qp materials and AD materials
164  if (!_nodal_material || is_ad)
165  {
166  (*_grad_mass_frac_qp)[_qp].resize(_num_phases);
167 
168  for (unsigned int ph = 0; ph < _num_phases; ++ph)
169  (*_grad_mass_frac_qp)[_qp][ph].assign(_num_components, RealGradient());
170  }
171 
172  // No derivatives are required for AD materials
173  if (!is_ad)
174  {
175  (*_dfluid_density_dvar)[_qp].resize(_num_phases);
176  (*_dfluid_viscosity_dvar)[_qp].resize(_num_phases);
177  (*_dfluid_enthalpy_dvar)[_qp].resize(_num_phases);
178  (*_dfluid_internal_energy_dvar)[_qp].resize(_num_phases);
179  (*_dmass_frac_dvar)[_qp].resize(_num_phases);
180 
181  for (unsigned int ph = 0; ph < _num_phases; ++ph)
182  {
183  (*_dfluid_density_dvar)[_qp][ph].assign(_num_pf_vars, 0.0);
184  (*_dfluid_viscosity_dvar)[_qp][ph].assign(_num_pf_vars, 0.0);
185  (*_dfluid_enthalpy_dvar)[_qp][ph].assign(_num_pf_vars, 0.0);
186  (*_dfluid_internal_energy_dvar)[_qp][ph].assign(_num_pf_vars, 0.0);
187  (*_dmass_frac_dvar)[_qp][ph].resize(_num_components);
188 
189  for (unsigned int comp = 0; comp < _num_components; ++comp)
190  (*_dmass_frac_dvar)[_qp][ph][comp].assign(_num_pf_vars, 0.0);
191  }
192  }
193  }
194 }
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_viscosity
Viscosity of each phase.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_internal_energy
Internal energy of each phase.
GenericMaterialProperty< std::vector< std::vector< Real > >, is_ad > & _mass_frac
Mass fraction matrix (indexing is fluid component in fluid phase)
GenericMaterialProperty< std::vector< std::vector< RealGradient > >, is_ad > *const _grad_mass_frac_qp
Gradient of the mass fraction matrix (only defined at the qps)
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > *const _dmass_frac_dvar
Derivative of the mass fraction matrix with respect to the Porous Flow variables. ...
const unsigned int _num_components
Number of components.
const unsigned int _num_pf_vars
Number of PorousFlow variables.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_density
Fluid density of each phase.
bool _is_initqp
Flag to indicate whether stateful properties should be computed.
const unsigned int _num_phases
Number of phases.
GenericMaterialProperty< std::vector< Real >, is_ad > & _fluid_enthalpy
Enthalpy of each phase.

◆ thermophysicalProperties()

template<bool is_ad>
virtual void PorousFlowFluidStateBaseMaterialTempl< is_ad >::thermophysicalProperties ( )
protectedpure virtual

Calculates all required thermophysical properties and derivatives for each phase and fluid component.

Implemented in PorousFlowFluidStateTempl< is_ad >, and PorousFlowFluidStateSingleComponentTempl< is_ad >.

◆ validParams()

template<bool is_ad>
InputParameters PorousFlowFluidStateBaseMaterialTempl< is_ad >::validParams ( )
static

Definition at line 15 of file PorousFlowFluidStateBaseMaterial.C.

Referenced by PorousFlowFluidStateSingleComponentTempl< is_ad >::validParams(), and PorousFlowFluidStateTempl< is_ad >::validParams().

16 {
18  MooseEnum unit_choice("Kelvin=0 Celsius=1", "Kelvin");
19  params.addParam<MooseEnum>(
20  "temperature_unit", unit_choice, "The unit of the temperature variable");
21  params.addRequiredParam<UserObjectName>("capillary_pressure",
22  "Name of the UserObject defining the capillary pressure");
23  params.addPrivateParam<std::string>("pf_material_type", "fluid_state");
24  params.addClassDescription("Base class for fluid state calculations using persistent primary "
25  "variables and a vapor-liquid flash");
26  return params;
27 }
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
void addPrivateParam(const std::string &name, const T &value)
static InputParameters validParams()
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _dfluid_density_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_dfluid_density_dvar
protected

Derivative of the fluid density for each phase wrt PorousFlow variables.

Definition at line 70 of file PorousFlowFluidStateBaseMaterial.h.

◆ _dfluid_enthalpy_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_dfluid_enthalpy_dvar
protected

Derivative of the fluid enthalpy for each phase wrt PorousFlow variables.

Definition at line 78 of file PorousFlowFluidStateBaseMaterial.h.

◆ _dfluid_internal_energy_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_dfluid_internal_energy_dvar
protected

Derivative of the fluid internal energy for each phase wrt PorousFlow variables.

Definition at line 82 of file PorousFlowFluidStateBaseMaterial.h.

◆ _dfluid_viscosity_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_dfluid_viscosity_dvar
protected

Derivative of the fluid viscosity for each phase wrt PorousFlow variables.

Definition at line 74 of file PorousFlowFluidStateBaseMaterial.h.

◆ _dgradp_qp_dgradv

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBaseTempl< is_ad >::_dgradp_qp_dgradv
protectedinherited

d(grad porepressure)/d(grad PorousFlow variable) at the quadpoints

Definition at line 51 of file PorousFlowVariableBase.h.

◆ _dgradp_qp_dv

template<bool is_ad>
MaterialProperty<std::vector<std::vector<RealGradient> > >* const PorousFlowVariableBaseTempl< is_ad >::_dgradp_qp_dv
protectedinherited

d(grad porepressure)/d(PorousFlow variable) at the quadpoints

Definition at line 54 of file PorousFlowVariableBase.h.

◆ _dgrads_qp_dgradv

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBaseTempl< is_ad >::_dgrads_qp_dgradv
protectedinherited

d(grad saturation)/d(grad PorousFlow variable) at the quadpoints

Definition at line 66 of file PorousFlowVariableBase.h.

◆ _dgrads_qp_dv

template<bool is_ad>
MaterialProperty<std::vector<std::vector<RealGradient> > >* const PorousFlowVariableBaseTempl< is_ad >::_dgrads_qp_dv
protectedinherited

d(grad saturation)/d(PorousFlow variable) at the quadpoints

Definition at line 69 of file PorousFlowVariableBase.h.

◆ _dmass_frac_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<std::vector<Real> > > >* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_dmass_frac_dvar
protected

Derivative of the mass fraction matrix with respect to the Porous Flow variables.

Definition at line 66 of file PorousFlowFluidStateBaseMaterial.h.

◆ _dporepressure_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBaseTempl< is_ad >::_dporepressure_dvar
protectedinherited

d(porepressure)/d(PorousFlow variable)

Definition at line 45 of file PorousFlowVariableBase.h.

◆ _dsaturation_dvar

template<bool is_ad>
MaterialProperty<std::vector<std::vector<Real> > >* const PorousFlowVariableBaseTempl< is_ad >::_dsaturation_dvar
protectedinherited

d(saturation)/d(PorousFlow variable)

Definition at line 60 of file PorousFlowVariableBase.h.

◆ _fluid_density

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_fluid_density
protected

Fluid density of each phase.

Definition at line 68 of file PorousFlowFluidStateBaseMaterial.h.

◆ _fluid_enthalpy

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_fluid_enthalpy
protected

Enthalpy of each phase.

Definition at line 76 of file PorousFlowFluidStateBaseMaterial.h.

◆ _fluid_internal_energy

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_fluid_internal_energy
protected

Internal energy of each phase.

Definition at line 80 of file PorousFlowFluidStateBaseMaterial.h.

◆ _fluid_viscosity

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_fluid_viscosity
protected

Viscosity of each phase.

Definition at line 72 of file PorousFlowFluidStateBaseMaterial.h.

◆ _fsp

template<bool is_ad>
std::vector<FluidStateProperties> PorousFlowFluidStateBaseMaterialTempl< is_ad >::_fsp
protected

FluidStateProperties data for each fluid component in each fluid phase.

Definition at line 55 of file PorousFlowFluidStateBaseMaterial.h.

Referenced by PorousFlowFluidStateBaseMaterialTempl< is_ad >::PorousFlowFluidStateBaseMaterialTempl().

◆ _grad_mass_frac_qp

template<bool is_ad>
GenericMaterialProperty<std::vector<std::vector<RealGradient> >, is_ad>* const PorousFlowFluidStateBaseMaterialTempl< is_ad >::_grad_mass_frac_qp
protected

Gradient of the mass fraction matrix (only defined at the qps)

Definition at line 64 of file PorousFlowFluidStateBaseMaterial.h.

◆ _gradp_qp

template<bool is_ad>
GenericMaterialProperty<std::vector<RealGradient>, is_ad>* const PorousFlowVariableBaseTempl< is_ad >::_gradp_qp
protectedinherited

Grad(p) at the quadpoints.

Definition at line 48 of file PorousFlowVariableBase.h.

◆ _grads_qp

template<bool is_ad>
GenericMaterialProperty<std::vector<RealGradient>, is_ad>* const PorousFlowVariableBaseTempl< is_ad >::_grads_qp
protectedinherited

Grad(s) at the quadpoints.

Definition at line 63 of file PorousFlowVariableBase.h.

◆ _is_initqp

template<bool is_ad>
bool PorousFlowFluidStateBaseMaterialTempl< is_ad >::_is_initqp
protected

Flag to indicate whether stateful properties should be computed.

Definition at line 53 of file PorousFlowFluidStateBaseMaterial.h.

◆ _mass_frac

template<bool is_ad>
GenericMaterialProperty<std::vector<std::vector<Real> >, is_ad>& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_mass_frac
protected

Mass fraction matrix (indexing is fluid component in fluid phase)

Definition at line 62 of file PorousFlowFluidStateBaseMaterial.h.

◆ _num_components

template<bool is_ad>
const unsigned int PorousFlowVariableBaseTempl< is_ad >::_num_components
protectedinherited

◆ _num_pf_vars

template<bool is_ad>
const unsigned int PorousFlowVariableBaseTempl< is_ad >::_num_pf_vars
protectedinherited

Number of PorousFlow variables.

Definition at line 39 of file PorousFlowVariableBase.h.

◆ _num_phases

template<bool is_ad>
const unsigned int PorousFlowVariableBaseTempl< is_ad >::_num_phases
protectedinherited

◆ _pc

template<bool is_ad>
const PorousFlowCapillaryPressure& PorousFlowFluidStateBaseMaterialTempl< is_ad >::_pc
protected

Capillary pressure UserObject.

Definition at line 57 of file PorousFlowFluidStateBaseMaterial.h.

◆ _porepressure

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowVariableBaseTempl< is_ad >::_porepressure
protectedinherited

◆ _saturation

template<bool is_ad>
GenericMaterialProperty<std::vector<Real>, is_ad>& PorousFlowVariableBaseTempl< is_ad >::_saturation
protectedinherited

◆ _sfx

template<bool is_ad>
const std::string PorousFlowFluidStateBaseMaterialTempl< is_ad >::_sfx
protected

Suffix to append to material property names (either qp or nodal as required)

Definition at line 59 of file PorousFlowFluidStateBaseMaterial.h.

◆ _T_c2k

template<bool is_ad>
const Real PorousFlowFluidStateBaseMaterialTempl< is_ad >::_T_c2k
protected

Conversion from degrees Celsius to degrees Kelvin.

Definition at line 51 of file PorousFlowFluidStateBaseMaterial.h.

◆ usingPorousFlowVariableBaseMembers

template<bool is_ad>
PorousFlowFluidStateBaseMaterialTempl< is_ad >::usingPorousFlowVariableBaseMembers
protected

Definition at line 84 of file PorousFlowFluidStateBaseMaterial.h.


The documentation for this class was generated from the following files: