https://mooseframework.inl.gov
PorousFlowFluidStateBaseMaterial.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 
12 
13 template <bool is_ad>
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 }
28 
29 template <bool is_ad>
31  const InputParameters & parameters)
32  : PorousFlowVariableBaseTempl<is_ad>(parameters),
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)),
39  _grad_mass_frac_qp(
40  _nodal_material
41  ? nullptr
42  : &this->template declareGenericProperty<std::vector<std::vector<RealGradient>>, is_ad>(
43  "PorousFlow_grad_mass_frac" + _sfx)),
44  _dmass_frac_dvar(
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)),
65  _dfluid_internal_energy_dvar(
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 }
73 
74 template <>
77 {
79 }
80 
81 template <>
84 {
85  return value;
86 }
87 
88 template <bool is_ad>
89 void
91 {
92  _is_initqp = true;
93  // Set the size of pressure and saturation vectors
95 
96  // Set the size of all other vectors
97  setMaterialVectorSize();
98 
99  thermophysicalProperties();
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);
112  _fluid_internal_energy[_qp][ph] = genericValue(_fsp[ph].internal_energy);
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 }
118 
119 template <bool is_ad>
120 void
122 {
123  _is_initqp = false;
124  // Size pressure and saturation and prepare the derivative vectors
126 
127  // Set the size of all other vectors
128  setMaterialVectorSize();
129 
130  // Calculate all required thermophysical properties
131  thermophysicalProperties();
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);
140  _fluid_internal_energy[_qp][ph] = genericValue(_fsp[ph].internal_energy);
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 }
146 
147 template <bool is_ad>
148 void
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 }
195 
virtual void computeQpProperties() override
Moose::GenericType< Real, is_ad > GenericReal
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 const std::string density
Definition: NS.h:33
auto raw_value(const Eigen::Map< T > &in)
static InputParameters validParams()
GenericReal< is_ad > genericValue(const ADReal &value)
Templated method that takes an ADReal value, and returns the raw value when is_ad = false...
void addRequiredParam(const std::string &name, const std::string &doc_string)
Base class for capillary pressure for multiphase flow in porous media.
Base class for thermophysical variable materials, which assemble materials for primary variables such...
virtual void setMaterialVectorSize() const
Size material property vectors and initialise with zeros.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
AD data structure to pass calculated thermophysical properties.
const unsigned int _num_components
Number of components.
virtual void initQpStatefulProperties() override
Fluid state base class using a persistent set of primary variables for multiphase, single and multicomponent cases.
const unsigned int _num_phases
Number of phases.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string pressure
Definition: NS.h:56
std::vector< FluidStateProperties > _fsp
FluidStateProperties data for each fluid component in each fluid phase.
void addClassDescription(const std::string &doc_string)
static const std::string internal_energy
Definition: NS.h:61
PorousFlowFluidStateBaseMaterialTempl(const InputParameters &parameters)