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 "PorousFlowSinglePhaseBase.h"
11 :
12 : #include "FEProblem.h"
13 : #include "Conversion.h"
14 : #include "libmesh/string_to_enum.h"
15 :
16 : InputParameters
17 6401 : PorousFlowSinglePhaseBase::validParams()
18 : {
19 6401 : InputParameters params = PorousFlowActionBase::validParams();
20 12802 : params.addParam<bool>("add_darcy_aux", true, "Add AuxVariables that record Darcy velocity");
21 12802 : params.addParam<bool>("add_stress_aux", true, "Add AuxVariables that record effective stress");
22 19203 : params.addDeprecatedParam<bool>("use_brine",
23 12802 : false,
24 : "Whether to use a PorousFlowBrine Material",
25 : "This parameter should no longer be used. Instead use "
26 : "fluid_properties_type = PorousFlowBrine");
27 12802 : params.addRequiredParam<VariableName>("porepressure", "The name of the porepressure variable");
28 12802 : MooseEnum coupling_type("Hydro ThermoHydro HydroMechanical ThermoHydroMechanical", "Hydro");
29 12802 : params.addParam<MooseEnum>("coupling_type",
30 : coupling_type,
31 : "The type of simulation. For simulations involving Mechanical "
32 : "deformations, you will need to supply the correct Biot coefficient. "
33 : "For simulations involving Thermal flows, you will need an associated "
34 : "ConstantThermalExpansionCoefficient Material");
35 : MooseEnum fluid_properties_type("PorousFlowSingleComponentFluid PorousFlowBrine Custom",
36 12802 : "PorousFlowSingleComponentFluid");
37 12802 : params.addParam<MooseEnum>(
38 : "fluid_properties_type",
39 : fluid_properties_type,
40 : "Type of fluid properties to use. For 'PorousFlowSingleComponentFluid' you must provide a "
41 : "fp UserObject. For 'PorousFlowBrine' you must supply a nacl_name. For "
42 : "'Custom' your input file must include a Material that provides fluid properties such as "
43 : "density, viscosity, enthalpy and internal energy");
44 12802 : MooseEnum simulation_type_choice("steady transient", "transient");
45 12802 : params.addDeprecatedParam<MooseEnum>(
46 : "simulation_type",
47 : simulation_type_choice,
48 : "Whether a transient or steady-state simulation is being performed",
49 : "The execution type is now determined automatically. This parameter should no longer be "
50 : "used");
51 12802 : params.addParam<UserObjectName>(
52 : "fp",
53 : "The name of the user object for fluid "
54 : "properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid");
55 12802 : params.addCoupledVar("mass_fraction_vars",
56 : {},
57 : "List of variables that represent the mass fractions. With only one fluid "
58 : "component, this may be left empty. With N fluid components, the format is "
59 : "'f_0 f_1 f_2 ... f_(N-1)'. That is, the N^th component need not be "
60 : "specified because f_N = 1 - (f_0 + f_1 + ... + f_(N-1)). It is best "
61 : "numerically to choose the N-1 mass fraction variables so that they "
62 : "represent the fluid components with small concentrations. This Action "
63 : "will associated the i^th mass fraction variable to the equation for the "
64 : "i^th fluid component, and the pressure variable to the N^th fluid "
65 : "component.");
66 19203 : params.addDeprecatedParam<unsigned>(
67 : "nacl_index",
68 12802 : 0,
69 : "Index of NaCl variable in mass_fraction_vars, for "
70 : "calculating brine properties. Only required if use_brine is true.",
71 : "This parameter should no longer be used. Instead use nacl_name = the_nacl_variable_name");
72 12802 : params.addParam<VariableName>(
73 : "nacl_name",
74 : "Name of the NaCl variable. Only required if fluid_properties_type = PorousFlowBrine");
75 12802 : params.addParam<Real>(
76 : "biot_coefficient",
77 12802 : 1.0,
78 : "The Biot coefficient (relevant only for mechanically-coupled simulations)");
79 12802 : params.addParam<std::vector<AuxVariableName>>(
80 : "save_component_rate_in",
81 : {},
82 : "List of AuxVariables into which the rate-of-change of each fluid component at each node "
83 : "will be saved. There must be exactly N of these to match the N fluid components. The "
84 : "result will be measured in kg/s, where the kg is the mass of the fluid component at the "
85 : "node (or m^3/s if multiply_by_density=false). Note that this saves the result from the "
86 : "MassTimeDerivative Kernels, but NOT from the MassVolumetricExpansion Kernels.");
87 12802 : MooseEnum temp_unit_choice("Kelvin Celsius", "Kelvin");
88 12802 : params.addParam<MooseEnum>(
89 : "temperature_unit", temp_unit_choice, "The unit of the temperature variable");
90 12802 : MooseEnum p_unit_choice("Pa MPa", "Pa");
91 12802 : params.addParam<MooseEnum>(
92 : "pressure_unit",
93 : p_unit_choice,
94 : "The unit of the pressure variable used everywhere in the input file "
95 : "except for in the FluidProperties-module objects. This can be set to the non-default value "
96 : "only for fluid_properties_type = PorousFlowSingleComponentFluid");
97 12802 : MooseEnum time_unit_choice("seconds hours days years", "seconds");
98 12802 : params.addParam<MooseEnum>(
99 : "time_unit",
100 : time_unit_choice,
101 : "The unit of time used everywhere in the input file except for in the "
102 : "FluidProperties-module objects. This can be set to the non-default value only for "
103 : "fluid_properties_type = PorousFlowSingleComponentFluid");
104 12802 : params.addParam<std::string>("base_name",
105 : "The base_name used in the TensorMechanics strain calculator. This "
106 : "is only relevant for mechanically-coupled models.");
107 6401 : params.addClassDescription("Base class for single-phase simulations");
108 6401 : return params;
109 6401 : }
110 :
111 6401 : PorousFlowSinglePhaseBase::PorousFlowSinglePhaseBase(const InputParameters & params)
112 : : PorousFlowActionBase(params),
113 12802 : _pp_var(getParam<VariableName>("porepressure")),
114 12802 : _coupling_type(getParam<MooseEnum>("coupling_type").getEnum<CouplingTypeEnum>()),
115 6401 : _thermal(_coupling_type == CouplingTypeEnum::ThermoHydro ||
116 : _coupling_type == CouplingTypeEnum::ThermoHydroMechanical),
117 6401 : _mechanical(_coupling_type == CouplingTypeEnum::HydroMechanical ||
118 : _coupling_type == CouplingTypeEnum::ThermoHydroMechanical),
119 6401 : _fluid_properties_type(
120 12802 : getParam<MooseEnum>("fluid_properties_type").getEnum<FluidPropertiesTypeEnum>()),
121 12802 : _biot_coefficient(getParam<Real>("biot_coefficient")),
122 12802 : _add_darcy_aux(getParam<bool>("add_darcy_aux")),
123 12802 : _add_stress_aux(getParam<bool>("add_stress_aux")),
124 12802 : _save_component_rate_in(getParam<std::vector<AuxVariableName>>("save_component_rate_in")),
125 12802 : _temperature_unit(getParam<MooseEnum>("temperature_unit")),
126 12802 : _pressure_unit(getParam<MooseEnum>("pressure_unit")),
127 12802 : _time_unit(getParam<MooseEnum>("time_unit")),
128 19203 : _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : "")
129 : {
130 6401 : if (_thermal && _temperature_var.size() != 1)
131 0 : mooseError("PorousFlowSinglePhaseBase: You need to specify a temperature variable to perform "
132 : "non-isothermal simulations");
133 :
134 6401 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
135 : {
136 12596 : if (params.isParamValid("nacl_name"))
137 2 : paramError("nacl_name",
138 : "PorousFlowSinglePhaseBase: You should not specify a nacl_name when "
139 : "fluid_properties_type = PorousFlowSingleComponentFluid");
140 12592 : if (!params.isParamValid("fp"))
141 2 : paramError("fp",
142 : "PorousFlowSinglePhaseBase: You must specify fp when fluid_properties_type = "
143 : "PorousFlowSingleComponentFluid");
144 12588 : _fp = getParam<UserObjectName>("fp");
145 : }
146 :
147 6397 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
148 : {
149 206 : if (!params.isParamValid("nacl_name"))
150 2 : paramError("nacl_name",
151 : "PorousFlowSinglePhaseBase: You must specify nacl_name when "
152 : "fluid_properties_type = PorousFlowBrine");
153 202 : if (params.isParamValid("fp"))
154 2 : paramError("fp",
155 : "PorousFlowSinglePhaseBase: You should not specify fp when "
156 : "fluid_properties_type = PorousFlowBrine");
157 99 : if (_pressure_unit != "Pa")
158 2 : paramError("pressure_unit",
159 : "Must use pressure_unit = Pa for fluid_properties_type = PorousFlowBrine");
160 97 : if (_time_unit != "seconds")
161 2 : paramError("time_unit",
162 : "Must use time_unit = seconds for fluid_properties_type = PorousFlowBrine");
163 190 : _nacl_name = getParam<VariableName>("nacl_name");
164 : }
165 :
166 : auto save_component_rate_in_size = _save_component_rate_in.size();
167 6389 : if (save_component_rate_in_size && save_component_rate_in_size != _num_mass_fraction_vars + 1)
168 2 : paramError("save_component_rate_in",
169 : "The number of save_component_rate_in variables must be the number of fluid "
170 : "components + 1");
171 6387 : }
172 :
173 : void
174 6321 : PorousFlowSinglePhaseBase::addMaterialDependencies()
175 : {
176 6321 : PorousFlowActionBase::addMaterialDependencies();
177 :
178 : // Add necessary objects to list of PorousFlow objects added by this action
179 6321 : if (_mechanical)
180 : {
181 2320 : _included_objects.push_back("StressDivergenceTensors");
182 2320 : _included_objects.push_back("Gravity");
183 4640 : _included_objects.push_back("PorousFlowEffectiveStressCoupling");
184 : }
185 :
186 6321 : if (_thermal)
187 : {
188 2015 : _included_objects.push_back("PorousFlowHeatConduction");
189 2015 : if (_transient)
190 3460 : _included_objects.push_back("PorousFlowEnergyTimeDerivative");
191 : }
192 :
193 6321 : if (_thermal && _mechanical && _transient)
194 1220 : _included_objects.push_back("PorousFlowHeatVolumetricExpansion");
195 :
196 6321 : if (_add_darcy_aux)
197 10872 : _included_objects.push_back("PorousFlowDarcyVelocityComponent");
198 :
199 6321 : if (_add_stress_aux && _mechanical)
200 4070 : _included_objects.push_back("StressAux");
201 6321 : }
202 :
203 : void
204 1249 : PorousFlowSinglePhaseBase::addKernels()
205 : {
206 1249 : PorousFlowActionBase::addKernels();
207 :
208 1249 : if (_mechanical)
209 : {
210 1761 : for (unsigned i = 0; i < _ndisp; ++i)
211 : {
212 1297 : std::string kernel_name = "PorousFlowUnsaturated_grad_stress" + Moose::stringify(i);
213 1297 : std::string kernel_type = "StressDivergenceTensors";
214 1297 : if (_coord_system == Moose::COORD_RZ)
215 : kernel_type = "StressDivergenceRZTensors";
216 1297 : InputParameters params = _factory.getValidParams(kernel_type);
217 1297 : if (_subdomain_names_set)
218 1634 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
219 2594 : params.set<NonlinearVariableName>("variable") = _displacements[i];
220 1297 : params.set<std::vector<VariableName>>("displacements") = _coupled_displacements;
221 1297 : if (_thermal)
222 : {
223 922 : params.set<std::vector<VariableName>>("temperature") = _temperature_var;
224 922 : if (parameters().isParamValid("eigenstrain_names"))
225 : {
226 922 : params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") =
227 1383 : getParam<std::vector<MaterialPropertyName>>("eigenstrain_names");
228 : }
229 : }
230 1297 : params.set<unsigned>("component") = i;
231 2594 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
232 1297 : _problem->addKernel(kernel_type, kernel_name, params);
233 :
234 1297 : if (_gravity(i) != 0)
235 : {
236 130 : kernel_name = "PorousFlowUnsaturated_gravity" + Moose::stringify(i);
237 : kernel_type = "Gravity";
238 65 : params = _factory.getValidParams(kernel_type);
239 65 : if (_subdomain_names_set)
240 76 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
241 130 : params.set<NonlinearVariableName>("variable") = _displacements[i];
242 65 : params.set<Real>("value") = _gravity(i);
243 130 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
244 65 : _problem->addKernel(kernel_type, kernel_name, params);
245 : }
246 :
247 2594 : kernel_name = "PorousFlowUnsaturated_EffStressCoupling" + Moose::stringify(i);
248 : kernel_type = "PorousFlowEffectiveStressCoupling";
249 1297 : params = _factory.getValidParams(kernel_type);
250 1297 : if (_subdomain_names_set)
251 1634 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
252 2594 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
253 2594 : params.set<NonlinearVariableName>("variable") = _displacements[i];
254 1297 : params.set<Real>("biot_coefficient") = _biot_coefficient;
255 1297 : params.set<unsigned>("component") = i;
256 2594 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
257 1297 : _problem->addKernel(kernel_type, kernel_name, params);
258 1297 : }
259 : }
260 :
261 1249 : if (_thermal)
262 : {
263 403 : std::string kernel_name = "PorousFlowUnsaturated_HeatConduction";
264 403 : std::string kernel_type = "PorousFlowHeatConduction";
265 403 : InputParameters params = _factory.getValidParams(kernel_type);
266 806 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
267 806 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
268 403 : _problem->addKernel(kernel_type, kernel_name, params);
269 :
270 403 : if (_transient)
271 : {
272 : kernel_name = "PorousFlowUnsaturated_EnergyTimeDerivative";
273 : kernel_type = "PorousFlowEnergyTimeDerivative";
274 346 : params = _factory.getValidParams(kernel_type);
275 692 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
276 692 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
277 346 : params.set<bool>("strain_at_nearest_qp") = _strain_at_nearest_qp;
278 346 : if (!_base_name.empty())
279 0 : params.set<std::string>("base_name") = _base_name;
280 346 : _problem->addKernel(kernel_type, kernel_name, params);
281 : }
282 403 : }
283 :
284 1249 : if (_thermal && _mechanical && _transient)
285 : {
286 122 : std::string kernel_name = "PorousFlowUnsaturated_HeatVolumetricExpansion";
287 122 : std::string kernel_type = "PorousFlowHeatVolumetricExpansion";
288 122 : InputParameters params = _factory.getValidParams(kernel_type);
289 244 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
290 244 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
291 122 : params.set<bool>("strain_at_nearest_qp") = _strain_at_nearest_qp;
292 122 : _problem->addKernel(kernel_type, kernel_name, params);
293 122 : }
294 1249 : }
295 :
296 : void
297 2528 : PorousFlowSinglePhaseBase::addAuxObjects()
298 : {
299 2528 : PorousFlowActionBase::addAuxObjects();
300 :
301 2528 : if (_add_darcy_aux)
302 2174 : addDarcyAux(_gravity);
303 :
304 2528 : if (_add_stress_aux && _mechanical)
305 814 : addStressAux();
306 2528 : }
307 :
308 : void
309 1267 : PorousFlowSinglePhaseBase::addMaterials()
310 : {
311 1267 : PorousFlowActionBase::addMaterials();
312 :
313 : // add Materials
314 2534 : if (_deps.dependsOn(_included_objects, "temperature_qp"))
315 1267 : addTemperatureMaterial(false);
316 :
317 2534 : if (_deps.dependsOn(_included_objects, "temperature_nodal"))
318 1029 : addTemperatureMaterial(true);
319 :
320 2534 : if (_deps.dependsOn(_included_objects, "mass_fraction_qp"))
321 439 : addMassFractionMaterial(false);
322 :
323 2534 : if (_deps.dependsOn(_included_objects, "mass_fraction_nodal"))
324 1018 : addMassFractionMaterial(true);
325 :
326 1267 : const bool compute_rho_mu_qp = _deps.dependsOn(_included_objects, "density_qp") ||
327 1267 : _deps.dependsOn(_included_objects, "viscosity_qp");
328 1267 : const bool compute_e_qp = _deps.dependsOn(_included_objects, "internal_energy_qp");
329 1267 : const bool compute_h_qp = _deps.dependsOn(_included_objects, "enthalpy_qp");
330 :
331 1267 : if (compute_rho_mu_qp || compute_e_qp || compute_h_qp)
332 : {
333 1267 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
334 38 : addBrineMaterial(
335 19 : _nacl_name, false, 0, compute_rho_mu_qp, compute_e_qp, compute_h_qp, _temperature_unit);
336 1248 : else if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
337 1248 : addSingleComponentFluidMaterial(false,
338 : 0,
339 : compute_rho_mu_qp,
340 : compute_e_qp,
341 : compute_h_qp,
342 1248 : _fp,
343 1248 : _temperature_unit,
344 1248 : _pressure_unit,
345 1248 : _time_unit);
346 : }
347 :
348 1267 : const bool compute_rho_mu_nodal = _deps.dependsOn(_included_objects, "density_nodal") ||
349 1562 : _deps.dependsOn(_included_objects, "viscosity_nodal");
350 1267 : const bool compute_e_nodal = _deps.dependsOn(_included_objects, "internal_energy_nodal");
351 1267 : const bool compute_h_nodal = _deps.dependsOn(_included_objects, "enthalpy_nodal");
352 :
353 1267 : if (compute_rho_mu_nodal || compute_e_nodal || compute_h_nodal)
354 : {
355 1029 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
356 38 : addBrineMaterial(_nacl_name,
357 : true,
358 : 0,
359 : compute_rho_mu_nodal,
360 : compute_e_nodal,
361 : compute_h_nodal,
362 19 : _temperature_unit);
363 1010 : else if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
364 1010 : addSingleComponentFluidMaterial(true,
365 : 0,
366 : compute_rho_mu_nodal,
367 : compute_e_nodal,
368 : compute_h_nodal,
369 1010 : _fp,
370 1010 : _temperature_unit,
371 1010 : _pressure_unit,
372 1010 : _time_unit);
373 : }
374 :
375 2534 : if (_deps.dependsOn(_included_objects, "effective_pressure_qp"))
376 1267 : addEffectiveFluidPressureMaterial(false);
377 :
378 2534 : if (_deps.dependsOn(_included_objects, "effective_pressure_nodal"))
379 872 : addEffectiveFluidPressureMaterial(true);
380 1267 : }
381 :
382 : void
383 1277 : PorousFlowSinglePhaseBase::addDictator()
384 : {
385 1277 : const std::string uo_name = _dictator_name;
386 1277 : const std::string uo_type = "PorousFlowDictator";
387 1277 : InputParameters params = _factory.getValidParams(uo_type);
388 1277 : std::vector<VariableName> pf_vars = _mass_fraction_vars;
389 1277 : pf_vars.push_back(_pp_var);
390 : // Only couple if in the same nonlinear system
391 2554 : params.set<SolverSystemName>("solver_sys") =
392 1277 : _problem->getSystemBase(_problem->getSystem(_pp_var).name()).name();
393 1680 : if (_thermal &&
394 403 : (_problem->getSystem(_temperature_var[0]).number() == _problem->getSystem(_pp_var).number()))
395 403 : pf_vars.push_back(_temperature_var[0]);
396 1741 : if (_mechanical && _coupled_displacements.size() &&
397 464 : (_problem->getSystem(_coupled_displacements[0]).number() ==
398 464 : _problem->getSystem(_pp_var).number()))
399 445 : pf_vars.insert(pf_vars.end(), _coupled_displacements.begin(), _coupled_displacements.end());
400 1277 : params.set<std::vector<VariableName>>("porous_flow_vars") = pf_vars;
401 1277 : params.set<unsigned int>("number_fluid_phases") = 1;
402 1277 : params.set<unsigned int>("number_fluid_components") = _num_mass_fraction_vars + 1;
403 1277 : params.set<unsigned int>("number_aqueous_equilibrium") = _num_aqueous_equilibrium;
404 1277 : params.set<unsigned int>("number_aqueous_kinetic") = _num_aqueous_kinetic;
405 1277 : _problem->addUserObject(uo_type, uo_name, params);
406 2554 : }
|
