Line data Source code
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 "PorousFlowPropertyAux.h"
11 :
12 : registerMooseObject("PorousFlowApp", PorousFlowPropertyAux);
13 : registerMooseObject("PorousFlowApp", ADPorousFlowPropertyAux);
14 :
15 : template <bool is_ad>
16 : InputParameters
17 30805 : PorousFlowPropertyAuxTempl<is_ad>::validParams()
18 : {
19 30805 : InputParameters params = AuxKernel::validParams();
20 61610 : params.addRequiredParam<UserObjectName>(
21 : "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names");
22 61610 : MooseEnum property_enum("pressure saturation temperature density viscosity mass_fraction relperm "
23 : "capillary_pressure enthalpy internal_energy secondary_concentration "
24 : "mineral_concentration mineral_reaction_rate porosity permeability "
25 : "hysteresis_order hysteresis_saturation_turning_point hysteretic_info");
26 61610 : params.addRequiredParam<MooseEnum>(
27 : "property", property_enum, "The fluid property that this auxillary kernel is to calculate");
28 61610 : params.addParam<unsigned int>("phase", 0, "The index of the phase this auxillary kernel acts on");
29 61610 : params.addParam<unsigned int>(
30 61610 : "liquid_phase", 0, "The index of the liquid phase (used for capillary pressure)");
31 61610 : params.addParam<unsigned int>(
32 61610 : "gas_phase", 1, "The index of the gas phase (used for capillary pressure)");
33 61610 : params.addParam<unsigned int>(
34 61610 : "fluid_component", 0, "The index of the fluid component this auxillary kernel acts on");
35 61610 : params.addParam<unsigned int>("secondary_species", 0, "The secondary chemical species number");
36 61610 : params.addParam<unsigned int>("mineral_species", 0, "The mineral chemical species number");
37 61610 : params.addParam<unsigned int>(
38 61610 : "hysteresis_turning_point", 0, "The hysteresis turning point number");
39 92415 : params.addRangeCheckedParam<unsigned int>(
40 61610 : "row", 0, "row>=0 & row<=2", "Row of permeability tensor to output");
41 92415 : params.addRangeCheckedParam<unsigned int>(
42 61610 : "column", 0, "column>=0 & column<=2", "Column of permeability tensor to output");
43 30805 : params.addClassDescription("AuxKernel to provide access to properties evaluated at quadpoints. "
44 : "Note that elemental AuxVariables must be used, so that these "
45 : "properties are integrated over each element.");
46 30805 : return params;
47 30805 : }
48 :
49 : template <bool is_ad>
50 16515 : PorousFlowPropertyAuxTempl<is_ad>::PorousFlowPropertyAuxTempl(const InputParameters & parameters)
51 : : AuxKernel(parameters),
52 16515 : _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
53 33030 : _property_enum(getParam<MooseEnum>("property").template getEnum<PropertyEnum>()),
54 33030 : _phase(getParam<unsigned int>("phase")),
55 33030 : _liquid_phase(getParam<unsigned int>("liquid_phase")),
56 33030 : _gas_phase(getParam<unsigned int>("gas_phase")),
57 33030 : _fluid_component(getParam<unsigned int>("fluid_component")),
58 33030 : _secondary_species(getParam<unsigned int>("secondary_species")),
59 33030 : _mineral_species(getParam<unsigned int>("mineral_species")),
60 33030 : _hysteresis_turning_point(getParam<unsigned int>("hysteresis_turning_point")),
61 33030 : _k_row(getParam<unsigned int>("row")),
62 49545 : _k_col(getParam<unsigned int>("column"))
63 : {
64 : // Check that the phase and fluid_component are valid
65 16515 : if (_phase >= _dictator.numPhases())
66 0 : paramError("phase",
67 : "Phase number entered is greater than the number of phases specified in the "
68 : "Dictator. Remember that indexing starts at 0");
69 :
70 16515 : if (_fluid_component >= _dictator.numComponents())
71 0 : paramError("fluid_component",
72 : "Fluid component number entered is greater than the number of fluid components "
73 : "specified in the Dictator. Remember that indexing starts at 0");
74 :
75 : // Check the parameters used to calculate capillary pressure
76 16515 : if (_property_enum == PropertyEnum::CAPILLARY_PRESSURE)
77 : {
78 44 : if (_liquid_phase >= _dictator.numPhases())
79 0 : paramError(
80 : "liquid_phase",
81 : "Liquid phase number entered is greater than the number of phases specified in the "
82 : "Dictator. Remember that indexing starts at 0");
83 :
84 44 : if (_gas_phase >= _dictator.numPhases())
85 0 : paramError("gas_phase",
86 : "Gas phase number entered is greater than the number of phases specified in the "
87 : "Dictator. Remember that indexing starts at 0");
88 :
89 44 : if (_liquid_phase == _gas_phase)
90 0 : paramError("liquid_phase", "Liquid phase number entered cannot be equal to gas_phase");
91 : }
92 :
93 16647 : if (_property_enum == PropertyEnum::SECONDARY_CONCENTRATION &&
94 132 : (_secondary_species >= _dictator.numAqueousEquilibrium()))
95 0 : paramError("secondary_species",
96 : "Secondary species number entered is greater than the number of aqueous equilibrium "
97 : "chemical reactions specified in the Dictator. Remember that indexing starts at 0");
98 :
99 16515 : if ((_property_enum == PropertyEnum::MINERAL_CONCENTRATION ||
100 16735 : _property_enum == PropertyEnum::MINERAL_REACTION_RATE) &&
101 220 : (_mineral_species >= _dictator.numAqueousKinetic()))
102 0 : paramError("mineral_species",
103 : "Mineral species number entered is greater than the number of aqueous "
104 : "precipitation-dissolution chemical reactions specified in the Dictator. Remember "
105 : "that indexing starts at 0");
106 :
107 16515 : if (_hysteresis_turning_point >= PorousFlowConstants::MAX_HYSTERESIS_ORDER)
108 2 : paramError("hysteresis_turning_point",
109 : "The maximum number of hysteresis turning points is ",
110 : PorousFlowConstants::MAX_HYSTERESIS_ORDER);
111 :
112 : // Only get material properties required by this instance of the AuxKernel
113 16513 : switch (_property_enum)
114 : {
115 1217 : case PropertyEnum::PRESSURE:
116 1217 : _pressure =
117 1217 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_porepressure_qp");
118 1217 : break;
119 :
120 2215 : case PropertyEnum::SATURATION:
121 2215 : _saturation =
122 2215 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_saturation_qp");
123 2215 : break;
124 :
125 154 : case PropertyEnum::TEMPERATURE:
126 154 : _temperature = &getGenericMaterialProperty<Real, is_ad>("PorousFlow_temperature_qp");
127 154 : break;
128 :
129 1008 : case PropertyEnum::DENSITY:
130 1008 : _fluid_density = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
131 : "PorousFlow_fluid_phase_density_qp");
132 1008 : break;
133 :
134 814 : case PropertyEnum::VISCOSITY:
135 814 : _fluid_viscosity =
136 814 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_viscosity_qp");
137 814 : break;
138 :
139 1379 : case PropertyEnum::MASS_FRACTION:
140 1379 : _mass_fractions = &getGenericMaterialProperty<std::vector<std::vector<Real>>, is_ad>(
141 : "PorousFlow_mass_frac_qp");
142 1379 : break;
143 :
144 858 : case PropertyEnum::RELPERM:
145 858 : _relative_permeability = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
146 : "PorousFlow_relative_permeability_qp");
147 858 : break;
148 :
149 44 : case PropertyEnum::CAPILLARY_PRESSURE:
150 44 : _pressure =
151 44 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_porepressure_qp");
152 44 : break;
153 :
154 528 : case PropertyEnum::ENTHALPY:
155 528 : _enthalpy = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
156 : "PorousFlow_fluid_phase_enthalpy_qp");
157 528 : break;
158 :
159 352 : case PropertyEnum::INTERNAL_ENERGY:
160 352 : _internal_energy = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
161 : "PorousFlow_fluid_phase_internal_energy_qp");
162 352 : break;
163 :
164 132 : case PropertyEnum::SECONDARY_CONCENTRATION:
165 132 : _sec_conc = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
166 : "PorousFlow_secondary_concentration_qp");
167 132 : break;
168 :
169 220 : case PropertyEnum::MINERAL_CONCENTRATION:
170 220 : _mineral_conc = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
171 : "PorousFlow_mineral_concentration_qp");
172 220 : break;
173 :
174 0 : case PropertyEnum::MINERAL_REACTION_RATE:
175 0 : _mineral_reaction_rate = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
176 : "PorousFlow_mineral_reaction_rate_qp");
177 0 : break;
178 :
179 647 : case PropertyEnum::POROSITY:
180 647 : _porosity = &getGenericMaterialProperty<Real, is_ad>("PorousFlow_porosity_qp");
181 647 : break;
182 :
183 1423 : case PropertyEnum::PERMEABILITY:
184 1423 : _permeability =
185 1423 : &getGenericMaterialProperty<RealTensorValue, is_ad>("PorousFlow_permeability_qp");
186 1423 : break;
187 :
188 2750 : case PropertyEnum::HYSTERESIS_ORDER:
189 2750 : _hys_order = &getMaterialProperty<unsigned int>("PorousFlow_hysteresis_order_qp");
190 2750 : break;
191 :
192 572 : case PropertyEnum::HYSTERESIS_SATURATION_TURNING_POINT:
193 572 : _hys_sat_tps =
194 572 : &getMaterialProperty<std::array<Real, PorousFlowConstants::MAX_HYSTERESIS_ORDER>>(
195 : "PorousFlow_hysteresis_saturation_tps_qp");
196 572 : break;
197 :
198 2200 : case PropertyEnum::HYSTERETIC_INFO:
199 2200 : _hys_info = &getMaterialProperty<Real>("PorousFlow_hysteretic_info_qp");
200 2200 : break;
201 : }
202 16513 : }
203 :
204 : template <bool is_ad>
205 : Real
206 4222917 : PorousFlowPropertyAuxTempl<is_ad>::computeValue()
207 : {
208 : Real property = 0.0;
209 :
210 4222917 : switch (_property_enum)
211 : {
212 452054 : case PropertyEnum::PRESSURE:
213 452054 : property = MetaPhysicL::raw_value((*_pressure)[_qp][_phase]);
214 452054 : break;
215 :
216 818332 : case PropertyEnum::SATURATION:
217 818332 : property = MetaPhysicL::raw_value((*_saturation)[_qp][_phase]);
218 818332 : break;
219 :
220 8684 : case PropertyEnum::TEMPERATURE:
221 8684 : property = MetaPhysicL::raw_value((*_temperature)[_qp]);
222 8684 : break;
223 :
224 327860 : case PropertyEnum::DENSITY:
225 327860 : property = MetaPhysicL::raw_value((*_fluid_density)[_qp][_phase]);
226 327860 : break;
227 :
228 248500 : case PropertyEnum::VISCOSITY:
229 248500 : property = MetaPhysicL::raw_value((*_fluid_viscosity)[_qp][_phase]);
230 248500 : break;
231 :
232 171364 : case PropertyEnum::MASS_FRACTION:
233 171364 : property = MetaPhysicL::raw_value((*_mass_fractions)[_qp][_phase][_fluid_component]);
234 171364 : break;
235 :
236 301694 : case PropertyEnum::RELPERM:
237 301694 : property = MetaPhysicL::raw_value((*_relative_permeability)[_qp][_phase]);
238 301694 : break;
239 :
240 16000 : case PropertyEnum::CAPILLARY_PRESSURE:
241 : property =
242 16000 : MetaPhysicL::raw_value((*_pressure)[_qp][_gas_phase] - (*_pressure)[_qp][_liquid_phase]);
243 16000 : break;
244 :
245 2268 : case PropertyEnum::ENTHALPY:
246 2268 : property = MetaPhysicL::raw_value((*_enthalpy)[_qp][_phase]);
247 2268 : break;
248 :
249 1668 : case PropertyEnum::INTERNAL_ENERGY:
250 1668 : property = MetaPhysicL::raw_value((*_internal_energy)[_qp][_phase]);
251 1668 : break;
252 :
253 200160 : case PropertyEnum::SECONDARY_CONCENTRATION:
254 200160 : property = MetaPhysicL::raw_value((*_sec_conc)[_qp][_secondary_species]);
255 200160 : break;
256 :
257 214276 : case PropertyEnum::MINERAL_CONCENTRATION:
258 214276 : property = MetaPhysicL::raw_value((*_mineral_conc)[_qp][_mineral_species]);
259 214276 : break;
260 :
261 0 : case PropertyEnum::MINERAL_REACTION_RATE:
262 0 : property = MetaPhysicL::raw_value((*_mineral_reaction_rate)[_qp][_mineral_species]);
263 0 : break;
264 :
265 291321 : case PropertyEnum::POROSITY:
266 291321 : property = MetaPhysicL::raw_value((*_porosity)[_qp]);
267 291321 : break;
268 :
269 432570 : case PropertyEnum::PERMEABILITY:
270 432570 : property = MetaPhysicL::raw_value((*_permeability)[_qp](_k_row, _k_col));
271 432570 : break;
272 :
273 384558 : case PropertyEnum::HYSTERESIS_ORDER:
274 384558 : property = (*_hys_order)[_qp];
275 384558 : break;
276 :
277 32396 : case PropertyEnum::HYSTERESIS_SATURATION_TURNING_POINT:
278 32396 : property = (*_hys_sat_tps)[_qp].at(_hysteresis_turning_point);
279 32396 : break;
280 :
281 319212 : case PropertyEnum::HYSTERETIC_INFO:
282 319212 : property = (*_hys_info)[_qp];
283 319212 : break;
284 : }
285 :
286 4222917 : return property;
287 : }
288 :
289 : template class PorousFlowPropertyAuxTempl<false>;
290 : template class PorousFlowPropertyAuxTempl<true>;
|
