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 3236 : PorousFlowSinglePhaseBase::validParams()
18 : {
19 3236 : InputParameters params = PorousFlowActionBase::validParams();
20 6472 : params.addParam<bool>("add_darcy_aux", true, "Add AuxVariables that record Darcy velocity");
21 6472 : params.addParam<bool>("add_stress_aux", true, "Add AuxVariables that record effective stress");
22 9708 : params.addDeprecatedParam<bool>("use_brine",
23 6472 : false,
24 : "Whether to use a PorousFlowBrine Material",
25 : "This parameter should no longer be used. Instead use "
26 : "fluid_properties_type = PorousFlowBrine");
27 6472 : params.addRequiredParam<VariableName>("porepressure", "The name of the porepressure variable");
28 6472 : MooseEnum coupling_type("Hydro ThermoHydro HydroMechanical ThermoHydroMechanical", "Hydro");
29 6472 : 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 6472 : "PorousFlowSingleComponentFluid");
37 6472 : 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 6472 : MooseEnum simulation_type_choice("steady transient", "transient");
45 6472 : 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 6472 : params.addParam<UserObjectName>(
52 : "fp",
53 : "The name of the user object for fluid "
54 : "properties. Only needed if fluid_properties_type = PorousFlowSingleComponentFluid");
55 6472 : 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 9708 : params.addDeprecatedParam<unsigned>(
67 : "nacl_index",
68 6472 : 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 6472 : params.addParam<VariableName>(
73 : "nacl_name",
74 : "Name of the NaCl variable. Only required if fluid_properties_type = PorousFlowBrine");
75 6472 : params.addParam<Real>(
76 : "biot_coefficient",
77 6472 : 1.0,
78 : "The Biot coefficient (relevant only for mechanically-coupled simulations)");
79 6472 : 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 6472 : MooseEnum temp_unit_choice("Kelvin Celsius", "Kelvin");
88 6472 : params.addParam<MooseEnum>(
89 : "temperature_unit", temp_unit_choice, "The unit of the temperature variable");
90 6472 : MooseEnum p_unit_choice("Pa MPa", "Pa");
91 6472 : 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 6472 : MooseEnum time_unit_choice("seconds hours days years", "seconds");
98 6472 : 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 6472 : 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 3236 : params.addClassDescription("Base class for single-phase simulations");
108 3236 : return params;
109 3236 : }
110 :
111 3236 : PorousFlowSinglePhaseBase::PorousFlowSinglePhaseBase(const InputParameters & params)
112 : : PorousFlowActionBase(params),
113 6472 : _pp_var(getParam<VariableName>("porepressure")),
114 6472 : _coupling_type(getParam<MooseEnum>("coupling_type").getEnum<CouplingTypeEnum>()),
115 3236 : _thermal(_coupling_type == CouplingTypeEnum::ThermoHydro ||
116 : _coupling_type == CouplingTypeEnum::ThermoHydroMechanical),
117 3236 : _mechanical(_coupling_type == CouplingTypeEnum::HydroMechanical ||
118 : _coupling_type == CouplingTypeEnum::ThermoHydroMechanical),
119 3236 : _fluid_properties_type(
120 6472 : getParam<MooseEnum>("fluid_properties_type").getEnum<FluidPropertiesTypeEnum>()),
121 6472 : _biot_coefficient(getParam<Real>("biot_coefficient")),
122 6472 : _add_darcy_aux(getParam<bool>("add_darcy_aux")),
123 6472 : _add_stress_aux(getParam<bool>("add_stress_aux")),
124 6472 : _save_component_rate_in(getParam<std::vector<AuxVariableName>>("save_component_rate_in")),
125 6472 : _temperature_unit(getParam<MooseEnum>("temperature_unit")),
126 6472 : _pressure_unit(getParam<MooseEnum>("pressure_unit")),
127 6472 : _time_unit(getParam<MooseEnum>("time_unit")),
128 9708 : _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : "")
129 : {
130 3236 : 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 3236 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
135 : {
136 6366 : 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 6362 : if (!params.isParamValid("fp"))
141 2 : paramError("fp",
142 : "PorousFlowSinglePhaseBase: You must specify fp when fluid_properties_type = "
143 : "PorousFlowSingleComponentFluid");
144 6358 : _fp = getParam<UserObjectName>("fp");
145 : }
146 :
147 3232 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
148 : {
149 106 : if (!params.isParamValid("nacl_name"))
150 2 : paramError("nacl_name",
151 : "PorousFlowSinglePhaseBase: You must specify nacl_name when "
152 : "fluid_properties_type = PorousFlowBrine");
153 102 : if (params.isParamValid("fp"))
154 2 : paramError("fp",
155 : "PorousFlowSinglePhaseBase: You should not specify fp when "
156 : "fluid_properties_type = PorousFlowBrine");
157 49 : if (_pressure_unit != "Pa")
158 2 : paramError("pressure_unit",
159 : "Must use pressure_unit = Pa for fluid_properties_type = PorousFlowBrine");
160 47 : if (_time_unit != "seconds")
161 2 : paramError("time_unit",
162 : "Must use time_unit = seconds for fluid_properties_type = PorousFlowBrine");
163 90 : _nacl_name = getParam<VariableName>("nacl_name");
164 : }
165 :
166 : auto save_component_rate_in_size = _save_component_rate_in.size();
167 3224 : 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 3222 : }
172 :
173 : void
174 3156 : PorousFlowSinglePhaseBase::addMaterialDependencies()
175 : {
176 3156 : PorousFlowActionBase::addMaterialDependencies();
177 :
178 : // Add necessary objects to list of PorousFlow objects added by this action
179 3156 : if (_mechanical)
180 : {
181 1110 : _included_objects.push_back("StressDivergenceTensors");
182 1110 : _included_objects.push_back("Gravity");
183 2220 : _included_objects.push_back("PorousFlowEffectiveStressCoupling");
184 : }
185 :
186 3156 : if (_thermal)
187 : {
188 1050 : _included_objects.push_back("PorousFlowHeatConduction");
189 1050 : if (_transient)
190 1830 : _included_objects.push_back("PorousFlowEnergyTimeDerivative");
191 : }
192 :
193 3156 : if (_thermal && _mechanical && _transient)
194 600 : _included_objects.push_back("PorousFlowHeatVolumetricExpansion");
195 :
196 3156 : if (_add_darcy_aux)
197 5352 : _included_objects.push_back("PorousFlowDarcyVelocityComponent");
198 :
199 3156 : if (_add_stress_aux && _mechanical)
200 1950 : _included_objects.push_back("StressAux");
201 3156 : }
202 :
203 : void
204 616 : PorousFlowSinglePhaseBase::addKernels()
205 : {
206 616 : PorousFlowActionBase::addKernels();
207 :
208 616 : if (_mechanical)
209 : {
210 843 : for (unsigned i = 0; i < _ndisp; ++i)
211 : {
212 621 : std::string kernel_name = "PorousFlowUnsaturated_grad_stress" + Moose::stringify(i);
213 621 : std::string kernel_type = "StressDivergenceTensors";
214 621 : if (_coord_system == Moose::COORD_RZ)
215 : kernel_type = "StressDivergenceRZTensors";
216 621 : InputParameters params = _factory.getValidParams(kernel_type);
217 621 : if (_subdomain_names_set)
218 774 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
219 1242 : params.set<NonlinearVariableName>("variable") = _displacements[i];
220 621 : params.set<std::vector<VariableName>>("displacements") = _coupled_displacements;
221 621 : if (_thermal)
222 : {
223 450 : params.set<std::vector<VariableName>>("temperature") = _temperature_var;
224 450 : if (parameters().isParamValid("eigenstrain_names"))
225 : {
226 450 : params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") =
227 675 : getParam<std::vector<MaterialPropertyName>>("eigenstrain_names");
228 : }
229 : }
230 621 : params.set<unsigned>("component") = i;
231 1242 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
232 621 : _problem->addKernel(kernel_type, kernel_name, params);
233 :
234 621 : if (_gravity(i) != 0)
235 : {
236 66 : kernel_name = "PorousFlowUnsaturated_gravity" + Moose::stringify(i);
237 : kernel_type = "Gravity";
238 33 : params = _factory.getValidParams(kernel_type);
239 33 : if (_subdomain_names_set)
240 36 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
241 66 : params.set<NonlinearVariableName>("variable") = _displacements[i];
242 33 : params.set<Real>("value") = _gravity(i);
243 66 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
244 33 : _problem->addKernel(kernel_type, kernel_name, params);
245 : }
246 :
247 1242 : kernel_name = "PorousFlowUnsaturated_EffStressCoupling" + Moose::stringify(i);
248 : kernel_type = "PorousFlowEffectiveStressCoupling";
249 621 : params = _factory.getValidParams(kernel_type);
250 621 : if (_subdomain_names_set)
251 774 : params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
252 1242 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
253 1242 : params.set<NonlinearVariableName>("variable") = _displacements[i];
254 621 : params.set<Real>("biot_coefficient") = _biot_coefficient;
255 621 : params.set<unsigned>("component") = i;
256 1242 : params.set<bool>("use_displaced_mesh") = getParam<bool>("use_displaced_mesh");
257 621 : _problem->addKernel(kernel_type, kernel_name, params);
258 621 : }
259 : }
260 :
261 616 : if (_thermal)
262 : {
263 210 : std::string kernel_name = "PorousFlowUnsaturated_HeatConduction";
264 210 : std::string kernel_type = "PorousFlowHeatConduction";
265 210 : InputParameters params = _factory.getValidParams(kernel_type);
266 420 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
267 420 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
268 210 : _problem->addKernel(kernel_type, kernel_name, params);
269 :
270 210 : if (_transient)
271 : {
272 : kernel_name = "PorousFlowUnsaturated_EnergyTimeDerivative";
273 : kernel_type = "PorousFlowEnergyTimeDerivative";
274 183 : params = _factory.getValidParams(kernel_type);
275 366 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
276 366 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
277 183 : params.set<bool>("strain_at_nearest_qp") = _strain_at_nearest_qp;
278 183 : if (!_base_name.empty())
279 0 : params.set<std::string>("base_name") = _base_name;
280 183 : _problem->addKernel(kernel_type, kernel_name, params);
281 : }
282 210 : }
283 :
284 616 : if (_thermal && _mechanical && _transient)
285 : {
286 60 : std::string kernel_name = "PorousFlowUnsaturated_HeatVolumetricExpansion";
287 60 : std::string kernel_type = "PorousFlowHeatVolumetricExpansion";
288 60 : InputParameters params = _factory.getValidParams(kernel_type);
289 120 : params.set<UserObjectName>("PorousFlowDictator") = _dictator_name;
290 120 : params.set<NonlinearVariableName>("variable") = _temperature_var[0];
291 60 : params.set<bool>("strain_at_nearest_qp") = _strain_at_nearest_qp;
292 60 : _problem->addKernel(kernel_type, kernel_name, params);
293 60 : }
294 616 : }
295 :
296 : void
297 1262 : PorousFlowSinglePhaseBase::addAuxObjects()
298 : {
299 1262 : PorousFlowActionBase::addAuxObjects();
300 :
301 1262 : if (_add_darcy_aux)
302 1070 : addDarcyAux(_gravity);
303 :
304 1262 : if (_add_stress_aux && _mechanical)
305 390 : addStressAux();
306 1262 : }
307 :
308 : void
309 634 : PorousFlowSinglePhaseBase::addMaterials()
310 : {
311 634 : PorousFlowActionBase::addMaterials();
312 :
313 : // add Materials
314 1268 : if (_deps.dependsOn(_included_objects, "temperature_qp"))
315 634 : addTemperatureMaterial(false);
316 :
317 1268 : if (_deps.dependsOn(_included_objects, "temperature_nodal"))
318 510 : addTemperatureMaterial(true);
319 :
320 1268 : if (_deps.dependsOn(_included_objects, "mass_fraction_qp"))
321 213 : addMassFractionMaterial(false);
322 :
323 1268 : if (_deps.dependsOn(_included_objects, "mass_fraction_nodal"))
324 503 : addMassFractionMaterial(true);
325 :
326 634 : const bool compute_rho_mu_qp = _deps.dependsOn(_included_objects, "density_qp") ||
327 634 : _deps.dependsOn(_included_objects, "viscosity_qp");
328 634 : const bool compute_e_qp = _deps.dependsOn(_included_objects, "internal_energy_qp");
329 634 : const bool compute_h_qp = _deps.dependsOn(_included_objects, "enthalpy_qp");
330 :
331 634 : if (compute_rho_mu_qp || compute_e_qp || compute_h_qp)
332 : {
333 634 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
334 18 : addBrineMaterial(
335 9 : _nacl_name, false, 0, compute_rho_mu_qp, compute_e_qp, compute_h_qp, _temperature_unit);
336 625 : else if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
337 625 : addSingleComponentFluidMaterial(false,
338 : 0,
339 : compute_rho_mu_qp,
340 : compute_e_qp,
341 : compute_h_qp,
342 625 : _fp,
343 625 : _temperature_unit,
344 625 : _pressure_unit,
345 625 : _time_unit);
346 : }
347 :
348 634 : const bool compute_rho_mu_nodal = _deps.dependsOn(_included_objects, "density_nodal") ||
349 785 : _deps.dependsOn(_included_objects, "viscosity_nodal");
350 634 : const bool compute_e_nodal = _deps.dependsOn(_included_objects, "internal_energy_nodal");
351 634 : const bool compute_h_nodal = _deps.dependsOn(_included_objects, "enthalpy_nodal");
352 :
353 634 : if (compute_rho_mu_nodal || compute_e_nodal || compute_h_nodal)
354 : {
355 510 : if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowBrine)
356 18 : addBrineMaterial(_nacl_name,
357 : true,
358 : 0,
359 : compute_rho_mu_nodal,
360 : compute_e_nodal,
361 : compute_h_nodal,
362 9 : _temperature_unit);
363 501 : else if (_fluid_properties_type == FluidPropertiesTypeEnum::PorousFlowSingleComponentFluid)
364 501 : addSingleComponentFluidMaterial(true,
365 : 0,
366 : compute_rho_mu_nodal,
367 : compute_e_nodal,
368 : compute_h_nodal,
369 501 : _fp,
370 501 : _temperature_unit,
371 501 : _pressure_unit,
372 501 : _time_unit);
373 : }
374 :
375 1268 : if (_deps.dependsOn(_included_objects, "effective_pressure_qp"))
376 634 : addEffectiveFluidPressureMaterial(false);
377 :
378 1268 : if (_deps.dependsOn(_included_objects, "effective_pressure_nodal"))
379 435 : addEffectiveFluidPressureMaterial(true);
380 634 : }
381 :
382 : void
383 644 : PorousFlowSinglePhaseBase::addDictator()
384 : {
385 644 : const std::string uo_name = _dictator_name;
386 644 : const std::string uo_type = "PorousFlowDictator";
387 644 : InputParameters params = _factory.getValidParams(uo_type);
388 644 : std::vector<VariableName> pf_vars = _mass_fraction_vars;
389 644 : pf_vars.push_back(_pp_var);
390 : // Only couple if in the same nonlinear system
391 1288 : params.set<SolverSystemName>("solver_sys") =
392 644 : _problem->getSystemBase(_problem->getSystem(_pp_var).name()).name();
393 854 : if (_thermal &&
394 210 : (_problem->getSystem(_temperature_var[0]).number() == _problem->getSystem(_pp_var).number()))
395 210 : pf_vars.push_back(_temperature_var[0]);
396 866 : if (_mechanical && _coupled_displacements.size() &&
397 222 : (_problem->getSystem(_coupled_displacements[0]).number() ==
398 222 : _problem->getSystem(_pp_var).number()))
399 213 : pf_vars.insert(pf_vars.end(), _coupled_displacements.begin(), _coupled_displacements.end());
400 644 : params.set<std::vector<VariableName>>("porous_flow_vars") = pf_vars;
401 644 : params.set<unsigned int>("number_fluid_phases") = 1;
402 644 : params.set<unsigned int>("number_fluid_components") = _num_mass_fraction_vars + 1;
403 644 : params.set<unsigned int>("number_aqueous_equilibrium") = _num_aqueous_equilibrium;
404 644 : params.set<unsigned int>("number_aqueous_kinetic") = _num_aqueous_kinetic;
405 644 : _problem->addUserObject(uo_type, uo_name, params);
406 1288 : }
|
