www.mooseframework.org
RichardsMaterial.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 #ifndef RICHARDSMATERIAL_H
9 #define RICHARDSMATERIAL_H
10 
11 #include "Material.h"
12 
13 #include "RichardsVarNames.h"
14 #include "RichardsDensity.h"
15 #include "RichardsRelPerm.h"
16 #include "RichardsSeff.h"
17 #include "RichardsSat.h"
18 #include "RichardsSUPG.h"
19 
20 // Forward Declarations
21 class RichardsMaterial;
22 
23 template <>
24 InputParameters validParams<RichardsMaterial>();
25 
26 class RichardsMaterial : public Material
27 {
28 public:
29  RichardsMaterial(const InputParameters & parameters);
30 
31 protected:
34 
36  const VariableValue & _por_change;
37  const VariableValue & _por_change_old;
38 
40  RealTensorValue _material_perm;
41 
43  RealVectorValue _material_gravity;
44 
46  MaterialProperty<Real> & _porosity_old;
47  MaterialProperty<Real> & _porosity;
48  MaterialProperty<RealTensorValue> & _permeability;
49  MaterialProperty<RealVectorValue> & _gravity;
50 
53  unsigned int _num_p;
54 
55  std::vector<const RichardsRelPerm *> _material_relperm_UO;
56  std::vector<const RichardsSeff *> _material_seff_UO;
57  std::vector<const RichardsSat *> _material_sat_UO;
58  std::vector<const RichardsDensity *> _material_density_UO;
59  std::vector<const RichardsSUPG *> _material_SUPG_UO;
60 
61  std::vector<const VariableValue *> _perm_change;
62 
63  virtual void computeProperties();
64 
72  void computePandSeff();
73 
80  void computeDerivedQuantities(unsigned int qp);
81 
87  void compute2ndDerivedQuantities(unsigned int qp);
88 
94  void zeroSUPG(unsigned int qp);
95 
97  void computeSUPG();
98 
99 private:
102 
103  std::vector<Real> _material_viscosity;
104 
106  MaterialProperty<std::vector<Real>> & _pp_old;
107 
109  MaterialProperty<std::vector<Real>> & _pp;
110 
112  MaterialProperty<std::vector<std::vector<Real>>> & _dpp_dv;
113 
115  MaterialProperty<std::vector<std::vector<std::vector<Real>>>> & _d2pp_dv;
116 
118  MaterialProperty<std::vector<Real>> & _viscosity;
119 
121  MaterialProperty<std::vector<Real>> & _density_old;
122 
124  MaterialProperty<std::vector<Real>> & _density;
125 
127  MaterialProperty<std::vector<std::vector<Real>>> & _ddensity_dv;
128 
130  MaterialProperty<std::vector<Real>> & _seff_old;
131 
133  MaterialProperty<std::vector<Real>> & _seff; // effective saturation
134 
136  MaterialProperty<std::vector<std::vector<Real>>> & _dseff_dv; // d(seff)/dp
137 
139  MaterialProperty<std::vector<std::vector<std::vector<Real>>>> & _d2seff_dv;
140 
142  MaterialProperty<std::vector<Real>> & _sat_old;
143 
145  MaterialProperty<std::vector<Real>> & _sat;
146 
148  MaterialProperty<std::vector<std::vector<Real>>> & _dsat_dv;
149 
151  MaterialProperty<std::vector<Real>> & _rel_perm;
152 
154  MaterialProperty<std::vector<std::vector<Real>>> & _drel_perm_dv;
155 
157  MaterialProperty<std::vector<Real>> & _mass_old;
158 
160  MaterialProperty<std::vector<Real>> & _mass;
161 
163  MaterialProperty<std::vector<std::vector<Real>>> & _dmass;
164 
166  MaterialProperty<std::vector<RealVectorValue>> & _flux_no_mob;
167 
169  MaterialProperty<std::vector<std::vector<RealVectorValue>>> & _dflux_no_mob_dv;
170 
172  MaterialProperty<std::vector<std::vector<RealTensorValue>>> & _dflux_no_mob_dgradv;
173 
175  MaterialProperty<std::vector<RealVectorValue>> & _flux;
176 
178  MaterialProperty<std::vector<std::vector<RealVectorValue>>> & _dflux_dv;
179 
181  MaterialProperty<std::vector<std::vector<RealTensorValue>>> & _dflux_dgradv;
182 
184  MaterialProperty<std::vector<std::vector<std::vector<RealVectorValue>>>> & _d2flux_dvdv;
185 
187  MaterialProperty<std::vector<std::vector<std::vector<RealTensorValue>>>> & _d2flux_dgradvdv;
188 
190  MaterialProperty<std::vector<std::vector<std::vector<RealTensorValue>>>> & _d2flux_dvdgradv;
191 
192  MaterialProperty<std::vector<RealVectorValue>> & _tauvel_SUPG; // tauSUPG * velSUPG
193  MaterialProperty<std::vector<std::vector<RealTensorValue>>> &
194  _dtauvel_SUPG_dgradp; // d (_tauvel_SUPG_i)/d(_grad_variable_j)
195  MaterialProperty<std::vector<std::vector<RealVectorValue>>> &
196  _dtauvel_SUPG_dp; // d (_tauvel_SUPG_i)/d(variable_j)
197 
199  std::vector<std::vector<std::vector<Real>>> _d2density;
200 
202  std::vector<std::vector<std::vector<Real>>> _d2rel_perm_dv;
203 
204  std::vector<const VariableValue *> _pressure_vals;
205  std::vector<const VariableValue *> _pressure_old_vals;
206  std::vector<const VariableGradient *> _grad_p;
207 
215  void zero2ndDerivedQuantities(unsigned int qp);
216 };
217 
218 #endif // RICHARDSMATERIAL_H
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _dflux_no_mob_dgradv
d(_flux_no_mob_i)/d(grad(variable_j))
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > & _d2pp_dv
d^2(porepressure_i)/d(variable_j)/d(variable_k)
MaterialProperty< std::vector< Real > > & _pp
porepressure(s)
std::vector< const VariableValue * > _perm_change
MaterialProperty< RealVectorValue > & _gravity
MaterialProperty< std::vector< std::vector< std::vector< RealTensorValue > > > > & _d2flux_dvdgradv
d^2(Richards flux_i)/d(variable_j)/d(grad(variable_k)), here flux_i is the i_th flux, which is itself a RealVectorValue. We should have _d2flux_dvdgradv[i][j][k] = _d2flux_dgradvdv[i][k][j], but i think it is more clear having both, and hopefully not a blowout on memory/CPU.
std::vector< std::vector< std::vector< Real > > > _d2rel_perm_dv
d^2(relperm_i)/dP_j/dP_k - used in various derivative calculations
RichardsMaterial(const InputParameters &parameters)
MaterialProperty< std::vector< std::vector< std::vector< RealTensorValue > > > > & _d2flux_dgradvdv
d^2(Richards flux_i)/d(grad(variable_j))/d(variable_k), here flux_i is the i_th flux, which is itself a RealVectorValue
std::vector< const RichardsRelPerm * > _material_relperm_UO
MaterialProperty< std::vector< Real > > & _pp_old
old values of porepressure(s)
const RichardsVarNames & _richards_name_UO
The variable names userobject for the Richards variables.
MaterialProperty< std::vector< RealVectorValue > > & _flux
fluid flux (a vector of fluxes for multicomponent)
std::vector< Real > _material_viscosity
unsigned int _num_p
This holds maps between pressure_var or pressure_var, sat_var used in RichardsMaterial and kernels...
MaterialProperty< std::vector< Real > > & _sat
saturation (vector of saturations in case of multiphase)
std::vector< const RichardsSUPG * > _material_SUPG_UO
MaterialProperty< std::vector< std::vector< Real > > > & _drel_perm_dv
d(relperm_i)/d(variable_j)
Real _material_por
porosity as entered by the user
const VariableValue & _por_change
porosity changes. if not entered they default to zero
void zero2ndDerivedQuantities(unsigned int qp)
Zeroes 2nd derivatives of the flux.
MaterialProperty< std::vector< Real > > & _seff_old
old effective saturation
MaterialProperty< std::vector< std::vector< Real > > > & _dpp_dv
d(porepressure_i)/d(variable_j)
void computePandSeff()
computes the quadpoint values of the porepressure(s) and effective saturation(s), and their derivativ...
MaterialProperty< std::vector< std::vector< std::vector< RealVectorValue > > > > & _d2flux_dvdv
d^2(Richards flux_i)/d(variable_j)/d(variable_k), here flux_i is the i_th flux, which is itself a Rea...
MaterialProperty< RealTensorValue > & _permeability
MaterialProperty< std::vector< std::vector< RealVectorValue > > > & _dflux_no_mob_dv
d(_flux_no_mob_i)/d(variable_j)
InputParameters validParams< RichardsMaterial >()
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _dtauvel_SUPG_dgradp
std::vector< const VariableGradient * > _grad_p
MaterialProperty< std::vector< std::vector< RealVectorValue > > > & _dflux_dv
d(Richards flux_i)/d(variable_j), here flux_i is the i_th flux, which is itself a RealVectorValue ...
std::vector< const VariableValue * > _pressure_vals
MaterialProperty< std::vector< Real > > & _rel_perm
relative permeability (vector of relative permeabilities in case of multiphase)
MaterialProperty< std::vector< Real > > & _mass
fluid mass (a vector of masses for multicomponent)
RealVectorValue _material_gravity
gravity as entered by user
MaterialProperty< std::vector< RealVectorValue > > & _flux_no_mob
permeability*(grad(P) - density*gravity) (a vector of these for multicomponent)
std::vector< const RichardsDensity * > _material_density_UO
void computeDerivedQuantities(unsigned int qp)
Computes the "derived" quantities — those that depend on porepressure(s) and effective saturation(s)...
std::vector< const VariableValue * > _pressure_old_vals
MaterialProperty< std::vector< Real > > & _density
fluid density (or densities for multiphase problems)
MaterialProperty< std::vector< Real > > & _sat_old
old saturation
std::vector< const RichardsSeff * > _material_seff_UO
MaterialProperty< std::vector< Real > > & _mass_old
old value of fluid mass (a vector of masses for multicomponent)
MaterialProperty< Real > & _porosity
void computeSUPG()
Computes the tauvel_SUPG and its derivatives.
MaterialProperty< std::vector< RealVectorValue > > & _tauvel_SUPG
MaterialProperty< std::vector< std::vector< std::vector< Real > > > > & _d2seff_dv
d^2(Seff_i)/d(variable_j)/d(variable_k)
void compute2ndDerivedQuantities(unsigned int qp)
Computes 2nd derivatives of the flux.
std::vector< const RichardsSat * > _material_sat_UO
void zeroSUPG(unsigned int qp)
Assigns and zeroes the MaterialProperties associated with SUPG.
MaterialProperty< std::vector< Real > > & _density_old
old fluid density (or densities for multiphase problems)
MaterialProperty< Real > & _porosity_old
material properties
const VariableValue & _por_change_old
MaterialProperty< std::vector< Real > > & _seff
effective saturation (vector of effective saturations in case of multiphase)
RealTensorValue _material_perm
permeability as entered by the user
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _dflux_dgradv
d(Richards flux_i)/d(grad(variable_j)), here flux_i is the i_th flux, which is itself a RealVectorVal...
Real _trace_perm
trace of permeability tensor
MaterialProperty< std::vector< std::vector< Real > > > & _ddensity_dv
d(density_i)/d(variable_j)
MaterialProperty< std::vector< std::vector< Real > > > & _dsat_dv
d(saturation_i)/d(variable_j)
MaterialProperty< std::vector< std::vector< Real > > > & _dseff_dv
d(Seff_i)/d(variable_j)
MaterialProperty< std::vector< std::vector< Real > > > & _dmass
d(fluid mass_i)/dP_j (a vector of masses for multicomponent)
MaterialProperty< std::vector< std::vector< RealVectorValue > > > & _dtauvel_SUPG_dp
std::vector< std::vector< std::vector< Real > > > _d2density
d^2(density)/dp_j/dP_k - used in various derivative calculations
MaterialProperty< std::vector< Real > > & _viscosity
fluid viscosity (or viscosities in the multiphase case)
virtual void computeProperties()