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 "PorousFlowPorosityLinear.h"
11 :
12 : registerMooseObject("PorousFlowApp", PorousFlowPorosityLinear);
13 :
14 : InputParameters
15 455 : PorousFlowPorosityLinear::validParams()
16 : {
17 455 : InputParameters params = PorousFlowPorosityBase::validParams();
18 910 : params.addParam<bool>("strain_at_nearest_qp",
19 910 : false,
20 : "When calculating nodal porosity that depends on strain, use the strain at "
21 : "the nearest quadpoint. This adds a small extra computational burden, and "
22 : "is not necessary for simulations involving only linear lagrange elements. "
23 : " If you set this to true, you will also want to set the same parameter to "
24 : "true for related Kernels and Materials");
25 910 : params.addRequiredCoupledVar(
26 : "porosity_ref",
27 : "The porosity at reference effective porepressure, temperature and volumetric strain. This "
28 : "must be a real number or a constant monomial variable (not a linear lagrange or other type "
29 : "of variable)");
30 910 : params.addCoupledVar(
31 : "P_ref",
32 : 0.0,
33 : "The reference effective porepressure. This should usually be a linear-lagrange variable");
34 910 : params.addCoupledVar(
35 : "T_ref",
36 : 0.0,
37 : "The reference temperature. This should usually be a linear-lagrange variable");
38 910 : params.addCoupledVar("epv_ref",
39 : 0.0,
40 : "The reference volumetric strain. This must be a real number or a constant "
41 : "monomial variable (not a linear lagrange or other type of variable)");
42 910 : params.addParam<Real>("P_coeff", 0.0, "Effective porepressure coefficient");
43 910 : params.addParam<Real>("T_coeff", 0.0, "Temperature coefficient");
44 910 : params.addParam<Real>("epv_coeff", 0.0, "Volumetric-strain coefficient");
45 1365 : params.addRangeCheckedParam<Real>(
46 : "porosity_min",
47 910 : 0.0,
48 : "porosity_min >= 0",
49 : "Minimum allowed value of the porosity: if the linear relationship gives "
50 : "values less than this value, then porosity is set to this value instead");
51 910 : params.addParam<Real>(
52 : "zero_modifier",
53 910 : 1E-3,
54 : "If the linear relationship produces porosity < porosity_min, then porosity is set "
55 : "porosity_min. This means the derivatives of it will be zero. However, these zero "
56 : "derivatives often result in poor NR convergence, so the derivatives are set to "
57 : "_zero_modifier * (values that are relevant for min porosity) to hint to the NR that "
58 : "porosity is not always constant");
59 910 : params.addParamNamesToGroup("zero_modifier", "Advanced");
60 455 : params.addClassDescription(
61 : "This Material calculates the porosity in PorousFlow simulations using the relationship "
62 : "porosity_ref + P_coeff * (P - P_ref) + T_coeff * (T - T_ref) + epv_coeff * (epv - epv_ref), "
63 : "where P is the effective porepressure, T is the temperature and epv is the volumetric "
64 : "strain");
65 455 : return params;
66 0 : }
67 :
68 357 : PorousFlowPorosityLinear::PorousFlowPorosityLinear(const InputParameters & parameters)
69 : : PorousFlowPorosityBase(parameters),
70 357 : _strain_at_nearest_qp(getParam<bool>("strain_at_nearest_qp")),
71 714 : _porosity_min(getParam<Real>("porosity_min")),
72 357 : _phi_ref(coupledValue("porosity_ref")),
73 615 : _P_ref(_nodal_material ? coupledDofValues("P_ref") : coupledValue("P_ref")),
74 615 : _T_ref(_nodal_material ? coupledDofValues("T_ref") : coupledValue("T_ref")),
75 357 : _epv_ref(coupledValue("epv_ref")),
76 714 : _P_coeff(getParam<Real>("P_coeff")),
77 714 : _T_coeff(getParam<Real>("T_coeff")),
78 714 : _epv_coeff(getParam<Real>("epv_coeff")),
79 :
80 507 : _uses_volstrain(parameters.isParamSetByUser("epv_coeff") && _epv_coeff != 0.0),
81 714 : _vol_strain_qp(getOptionalMaterialProperty<Real>("PorousFlow_total_volumetric_strain_qp")),
82 714 : _dvol_strain_qp_dvar(getOptionalMaterialProperty<std::vector<RealGradient>>(
83 : "dPorousFlow_total_volumetric_strain_qp_dvar")),
84 :
85 375 : _uses_pf(parameters.isParamSetByUser("P_coeff") && _P_coeff != 0.0),
86 714 : _pf_nodal(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_nodal")),
87 714 : _pf_qp(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_qp")),
88 714 : _dpf_dvar_nodal(getOptionalMaterialProperty<std::vector<Real>>(
89 : "dPorousFlow_effective_fluid_pressure_nodal_dvar")),
90 714 : _dpf_dvar_qp(getOptionalMaterialProperty<std::vector<Real>>(
91 : "dPorousFlow_effective_fluid_pressure_qp_dvar")),
92 :
93 507 : _uses_T(parameters.isParamSetByUser("T_coeff") && _T_coeff != 0.0),
94 714 : _temperature_nodal(getOptionalMaterialProperty<Real>("PorousFlow_temperature_nodal")),
95 714 : _temperature_qp(getOptionalMaterialProperty<Real>("PorousFlow_temperature_qp")),
96 357 : _dtemperature_dvar_nodal(
97 357 : getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_nodal_dvar")),
98 357 : _dtemperature_dvar_qp(
99 357 : getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_qp_dvar")),
100 1071 : _zero_modifier(getParam<Real>("zero_modifier"))
101 : {
102 357 : }
103 :
104 : void
105 336 : PorousFlowPorosityLinear::initialSetup()
106 : {
107 336 : _pf = _nodal_material ? _pf_nodal.get() : _pf_qp.get();
108 336 : _dpf_dvar = _nodal_material ? _dpf_dvar_nodal.get() : _dpf_dvar_qp.get();
109 336 : if (!(_pf && _dpf_dvar) && _uses_pf)
110 : {
111 2 : mooseError("PorousFlowPorosityLinear: When P_coeff is given you must have an effective fluid "
112 : "pressure Material");
113 : }
114 :
115 334 : _temperature = _nodal_material ? _temperature_nodal.get() : _temperature_qp.get();
116 334 : _dtemperature_dvar =
117 334 : _nodal_material ? _dtemperature_dvar_nodal.get() : _dtemperature_dvar_qp.get();
118 334 : if (!_temperature && _uses_T)
119 2 : mooseError(
120 : "PorousFlowPorosityLinear: When T_coeff is given you must have a temperature Material");
121 :
122 332 : if (!(_vol_strain_qp && _dvol_strain_qp_dvar) && _uses_volstrain)
123 2 : mooseError("PorousFlowPorosityLinear: When epv_coeff is given you must have a "
124 : "volumetric-strain Material");
125 :
126 330 : PorousFlowPorosityBase::initialSetup();
127 330 : }
128 :
129 : void
130 179731 : PorousFlowPorosityLinear::initQpStatefulProperties()
131 : {
132 179731 : _porosity[_qp] = _phi_ref[0];
133 179731 : (*_dporosity_dvar)[_qp].assign(_num_var, 0.0);
134 179731 : (*_dporosity_dgradvar)[_qp].assign(_num_var, 0.0);
135 179731 : }
136 :
137 : void
138 1028640 : PorousFlowPorosityLinear::computeQpProperties()
139 : {
140 : // note the [0] below: _phi_ref is a constant monomial and we use [0] regardless of
141 : // _nodal_material
142 1028640 : Real porosity = _phi_ref[0];
143 1028640 : if (_uses_pf)
144 1028640 : porosity += _P_coeff * ((*_pf)[_qp] - _P_ref[_qp]);
145 1028640 : if (_uses_T)
146 2440 : porosity += _T_coeff * ((*_temperature)[_qp] - _T_ref[_qp]);
147 1028640 : if (_uses_volstrain)
148 : {
149 : // Note that in the following _strain[_qp] is evaluated at q quadpoint
150 : // So _porosity_nodal[_qp], which should be the nodal value of porosity
151 : // actually uses the strain at a quadpoint. This
152 : // is OK for LINEAR elements, as strain is constant over the element anyway.
153 : const unsigned qp_to_use =
154 2440 : (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
155 2440 : porosity += _epv_coeff * (_vol_strain_qp[qp_to_use] - _epv_ref[0]);
156 : }
157 :
158 1028640 : if (porosity < _porosity_min)
159 205 : _porosity[_qp] = _porosity_min;
160 : else
161 1028435 : _porosity[_qp] = porosity;
162 :
163 1028640 : (*_dporosity_dvar)[_qp].resize(_num_var);
164 1028640 : (*_dporosity_dgradvar)[_qp].resize(_num_var);
165 3092760 : for (unsigned int pvar = 0; pvar < _num_var; ++pvar)
166 : {
167 : Real deriv = 0.0;
168 2064120 : if (_uses_pf)
169 2064120 : deriv += _P_coeff * (*_dpf_dvar)[_qp][pvar];
170 2064120 : if (_uses_T)
171 11720 : deriv += _T_coeff * (*_dtemperature_dvar)[_qp][pvar];
172 2064120 : if (porosity < _porosity_min)
173 590 : (*_dporosity_dvar)[_qp][pvar] = _zero_modifier * deriv;
174 : else
175 2063530 : (*_dporosity_dvar)[_qp][pvar] = deriv;
176 :
177 : RealGradient deriv_grad(0.0, 0.0, 0.0);
178 2064120 : if (_uses_volstrain)
179 : {
180 : const unsigned qp_to_use =
181 11720 : (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
182 11720 : deriv_grad += _epv_coeff * _dvol_strain_qp_dvar[qp_to_use][pvar];
183 : }
184 2064120 : if (porosity < _porosity_min)
185 590 : (*_dporosity_dgradvar)[_qp][pvar] = _zero_modifier * deriv_grad;
186 : else
187 2063530 : (*_dporosity_dgradvar)[_qp][pvar] = deriv_grad;
188 : }
189 1028640 : }
|