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 14593 : PorousFlowPropertyAuxTempl<is_ad>::validParams()
18 : {
19 14593 : InputParameters params = AuxKernel::validParams();
20 29186 : params.addRequiredParam<UserObjectName>(
21 : "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names");
22 29186 : 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 29186 : params.addRequiredParam<MooseEnum>(
27 : "property", property_enum, "The fluid property that this auxillary kernel is to calculate");
28 29186 : params.addParam<unsigned int>("phase", 0, "The index of the phase this auxillary kernel acts on");
29 29186 : params.addParam<unsigned int>(
30 29186 : "liquid_phase", 0, "The index of the liquid phase (used for capillary pressure)");
31 29186 : params.addParam<unsigned int>(
32 29186 : "gas_phase", 1, "The index of the gas phase (used for capillary pressure)");
33 29186 : params.addParam<unsigned int>(
34 29186 : "fluid_component", 0, "The index of the fluid component this auxillary kernel acts on");
35 29186 : params.addParam<unsigned int>("secondary_species", 0, "The secondary chemical species number");
36 29186 : params.addParam<unsigned int>("mineral_species", 0, "The mineral chemical species number");
37 29186 : params.addParam<unsigned int>(
38 29186 : "hysteresis_turning_point", 0, "The hysteresis turning point number");
39 43779 : params.addRangeCheckedParam<unsigned int>(
40 29186 : "row", 0, "row>=0 & row<=2", "Row of permeability tensor to output");
41 43779 : params.addRangeCheckedParam<unsigned int>(
42 29186 : "column", 0, "column>=0 & column<=2", "Column of permeability tensor to output");
43 14593 : 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 14593 : return params;
47 14593 : }
48 :
49 : template <bool is_ad>
50 7677 : PorousFlowPropertyAuxTempl<is_ad>::PorousFlowPropertyAuxTempl(const InputParameters & parameters)
51 : : AuxKernel(parameters),
52 7677 : _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
53 15354 : _property_enum(getParam<MooseEnum>("property").template getEnum<PropertyEnum>()),
54 15354 : _phase(getParam<unsigned int>("phase")),
55 15354 : _liquid_phase(getParam<unsigned int>("liquid_phase")),
56 15354 : _gas_phase(getParam<unsigned int>("gas_phase")),
57 15354 : _fluid_component(getParam<unsigned int>("fluid_component")),
58 15354 : _secondary_species(getParam<unsigned int>("secondary_species")),
59 15354 : _mineral_species(getParam<unsigned int>("mineral_species")),
60 15354 : _hysteresis_turning_point(getParam<unsigned int>("hysteresis_turning_point")),
61 15354 : _k_row(getParam<unsigned int>("row")),
62 23031 : _k_col(getParam<unsigned int>("column"))
63 : {
64 : // Check that the phase and fluid_component are valid
65 7677 : 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 7677 : 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 7677 : if (_property_enum == PropertyEnum::CAPILLARY_PRESSURE)
77 : {
78 20 : 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 20 : 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 20 : if (_liquid_phase == _gas_phase)
90 0 : paramError("liquid_phase", "Liquid phase number entered cannot be equal to gas_phase");
91 : }
92 :
93 7737 : if (_property_enum == PropertyEnum::SECONDARY_CONCENTRATION &&
94 60 : (_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 7677 : if ((_property_enum == PropertyEnum::MINERAL_CONCENTRATION ||
100 7777 : _property_enum == PropertyEnum::MINERAL_REACTION_RATE) &&
101 100 : (_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 7677 : 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 7675 : switch (_property_enum)
114 : {
115 557 : case PropertyEnum::PRESSURE:
116 557 : _pressure =
117 557 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_porepressure_qp");
118 557 : break;
119 :
120 1035 : case PropertyEnum::SATURATION:
121 1035 : _saturation =
122 1035 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_saturation_qp");
123 1035 : break;
124 :
125 70 : case PropertyEnum::TEMPERATURE:
126 70 : _temperature = &getGenericMaterialProperty<Real, is_ad>("PorousFlow_temperature_qp");
127 70 : break;
128 :
129 470 : case PropertyEnum::DENSITY:
130 470 : _fluid_density = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
131 : "PorousFlow_fluid_phase_density_qp");
132 470 : break;
133 :
134 370 : case PropertyEnum::VISCOSITY:
135 370 : _fluid_viscosity =
136 370 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_viscosity_qp");
137 370 : break;
138 :
139 631 : case PropertyEnum::MASS_FRACTION:
140 631 : _mass_fractions = &getGenericMaterialProperty<std::vector<std::vector<Real>>, is_ad>(
141 : "PorousFlow_mass_frac_qp");
142 631 : break;
143 :
144 390 : case PropertyEnum::RELPERM:
145 390 : _relative_permeability = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
146 : "PorousFlow_relative_permeability_qp");
147 390 : break;
148 :
149 20 : case PropertyEnum::CAPILLARY_PRESSURE:
150 20 : _pressure =
151 20 : &getGenericMaterialProperty<std::vector<Real>, is_ad>("PorousFlow_porepressure_qp");
152 20 : break;
153 :
154 240 : case PropertyEnum::ENTHALPY:
155 240 : _enthalpy = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
156 : "PorousFlow_fluid_phase_enthalpy_qp");
157 240 : break;
158 :
159 160 : case PropertyEnum::INTERNAL_ENERGY:
160 160 : _internal_energy = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
161 : "PorousFlow_fluid_phase_internal_energy_qp");
162 160 : break;
163 :
164 60 : case PropertyEnum::SECONDARY_CONCENTRATION:
165 60 : _sec_conc = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
166 : "PorousFlow_secondary_concentration_qp");
167 60 : break;
168 :
169 100 : case PropertyEnum::MINERAL_CONCENTRATION:
170 100 : _mineral_conc = &getGenericMaterialProperty<std::vector<Real>, is_ad>(
171 : "PorousFlow_mineral_concentration_qp");
172 100 : 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 311 : case PropertyEnum::POROSITY:
180 311 : _porosity = &getGenericMaterialProperty<Real, is_ad>("PorousFlow_porosity_qp");
181 311 : break;
182 :
183 751 : case PropertyEnum::PERMEABILITY:
184 751 : _permeability =
185 751 : &getGenericMaterialProperty<RealTensorValue, is_ad>("PorousFlow_permeability_qp");
186 751 : break;
187 :
188 1250 : case PropertyEnum::HYSTERESIS_ORDER:
189 1250 : _hys_order = &getMaterialProperty<unsigned int>("PorousFlow_hysteresis_order_qp");
190 1250 : break;
191 :
192 260 : case PropertyEnum::HYSTERESIS_SATURATION_TURNING_POINT:
193 260 : _hys_sat_tps =
194 260 : &getMaterialProperty<std::array<Real, PorousFlowConstants::MAX_HYSTERESIS_ORDER>>(
195 : "PorousFlow_hysteresis_saturation_tps_qp");
196 260 : break;
197 :
198 1000 : case PropertyEnum::HYSTERETIC_INFO:
199 1000 : _hys_info = &getMaterialProperty<Real>("PorousFlow_hysteretic_info_qp");
200 1000 : break;
201 : }
202 7675 : }
203 :
204 : template <bool is_ad>
205 : Real
206 3032082 : PorousFlowPropertyAuxTempl<is_ad>::computeValue()
207 : {
208 : Real property = 0.0;
209 :
210 3032082 : switch (_property_enum)
211 : {
212 303494 : case PropertyEnum::PRESSURE:
213 303494 : property = MetaPhysicL::raw_value((*_pressure)[_qp][_phase]);
214 303494 : break;
215 :
216 661522 : case PropertyEnum::SATURATION:
217 661522 : property = MetaPhysicL::raw_value((*_saturation)[_qp][_phase]);
218 661522 : break;
219 :
220 5846 : case PropertyEnum::TEMPERATURE:
221 5846 : property = MetaPhysicL::raw_value((*_temperature)[_qp]);
222 5846 : break;
223 :
224 227450 : case PropertyEnum::DENSITY:
225 227450 : property = MetaPhysicL::raw_value((*_fluid_density)[_qp][_phase]);
226 227450 : break;
227 :
228 173842 : case PropertyEnum::VISCOSITY:
229 173842 : property = MetaPhysicL::raw_value((*_fluid_viscosity)[_qp][_phase]);
230 173842 : break;
231 :
232 115576 : case PropertyEnum::MASS_FRACTION:
233 115576 : property = MetaPhysicL::raw_value((*_mass_fractions)[_qp][_phase][_fluid_component]);
234 115576 : break;
235 :
236 203750 : case PropertyEnum::RELPERM:
237 203750 : property = MetaPhysicL::raw_value((*_relative_permeability)[_qp][_phase]);
238 203750 : break;
239 :
240 10948 : case PropertyEnum::CAPILLARY_PRESSURE:
241 : property =
242 10948 : MetaPhysicL::raw_value((*_pressure)[_qp][_gas_phase] - (*_pressure)[_qp][_liquid_phase]);
243 10948 : break;
244 :
245 1542 : case PropertyEnum::ENTHALPY:
246 1542 : property = MetaPhysicL::raw_value((*_enthalpy)[_qp][_phase]);
247 1542 : break;
248 :
249 1134 : case PropertyEnum::INTERNAL_ENERGY:
250 1134 : property = MetaPhysicL::raw_value((*_internal_energy)[_qp][_phase]);
251 1134 : break;
252 :
253 134160 : case PropertyEnum::SECONDARY_CONCENTRATION:
254 134160 : property = MetaPhysicL::raw_value((*_sec_conc)[_qp][_secondary_species]);
255 134160 : break;
256 :
257 143716 : case PropertyEnum::MINERAL_CONCENTRATION:
258 143716 : property = MetaPhysicL::raw_value((*_mineral_conc)[_qp][_mineral_species]);
259 143716 : 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 206550 : case PropertyEnum::POROSITY:
266 206550 : property = MetaPhysicL::raw_value((*_porosity)[_qp]);
267 206550 : break;
268 :
269 330024 : case PropertyEnum::PERMEABILITY:
270 330024 : property = MetaPhysicL::raw_value((*_permeability)[_qp](_k_row, _k_col));
271 330024 : break;
272 :
273 267582 : case PropertyEnum::HYSTERESIS_ORDER:
274 267582 : property = (*_hys_order)[_qp];
275 267582 : break;
276 :
277 21638 : case PropertyEnum::HYSTERESIS_SATURATION_TURNING_POINT:
278 21638 : property = (*_hys_sat_tps)[_qp].at(_hysteresis_turning_point);
279 21638 : break;
280 :
281 223308 : case PropertyEnum::HYSTERETIC_INFO:
282 223308 : property = (*_hys_info)[_qp];
283 223308 : break;
284 : }
285 :
286 3032082 : return property;
287 : }
288 :
289 : template class PorousFlowPropertyAuxTempl<false>;
290 : template class PorousFlowPropertyAuxTempl<true>;
|
