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 243 : PorousFlowPorosityLinear::validParams()
16 : {
17 243 : InputParameters params = PorousFlowPorosityBase::validParams();
18 486 : params.addParam<bool>("strain_at_nearest_qp",
19 486 : 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 486 : 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 486 : params.addCoupledVar(
31 : "P_ref",
32 : 0.0,
33 : "The reference effective porepressure. This should usually be a linear-lagrange variable");
34 486 : params.addCoupledVar(
35 : "T_ref",
36 : 0.0,
37 : "The reference temperature. This should usually be a linear-lagrange variable");
38 486 : 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 486 : params.addParam<Real>("P_coeff", 0.0, "Effective porepressure coefficient");
43 486 : params.addParam<Real>("T_coeff", 0.0, "Temperature coefficient");
44 486 : params.addParam<Real>("epv_coeff", 0.0, "Volumetric-strain coefficient");
45 729 : params.addRangeCheckedParam<Real>(
46 : "porosity_min",
47 486 : 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 486 : params.addParam<Real>(
52 : "zero_modifier",
53 486 : 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 486 : params.addParamNamesToGroup("zero_modifier", "Advanced");
60 243 : 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 243 : return params;
66 0 : }
67 :
68 189 : PorousFlowPorosityLinear::PorousFlowPorosityLinear(const InputParameters & parameters)
69 : : PorousFlowPorosityBase(parameters),
70 189 : _strain_at_nearest_qp(getParam<bool>("strain_at_nearest_qp")),
71 378 : _porosity_min(getParam<Real>("porosity_min")),
72 189 : _phi_ref(coupledValue("porosity_ref")),
73 327 : _P_ref(_nodal_material ? coupledDofValues("P_ref") : coupledValue("P_ref")),
74 327 : _T_ref(_nodal_material ? coupledDofValues("T_ref") : coupledValue("T_ref")),
75 189 : _epv_ref(coupledValue("epv_ref")),
76 378 : _P_coeff(getParam<Real>("P_coeff")),
77 378 : _T_coeff(getParam<Real>("T_coeff")),
78 378 : _epv_coeff(getParam<Real>("epv_coeff")),
79 :
80 267 : _uses_volstrain(parameters.isParamSetByUser("epv_coeff") && _epv_coeff != 0.0),
81 378 : _vol_strain_qp(getOptionalMaterialProperty<Real>("PorousFlow_total_volumetric_strain_qp")),
82 378 : _dvol_strain_qp_dvar(getOptionalMaterialProperty<std::vector<RealGradient>>(
83 : "dPorousFlow_total_volumetric_strain_qp_dvar")),
84 :
85 207 : _uses_pf(parameters.isParamSetByUser("P_coeff") && _P_coeff != 0.0),
86 378 : _pf_nodal(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_nodal")),
87 378 : _pf_qp(getOptionalMaterialProperty<Real>("PorousFlow_effective_fluid_pressure_qp")),
88 378 : _dpf_dvar_nodal(getOptionalMaterialProperty<std::vector<Real>>(
89 : "dPorousFlow_effective_fluid_pressure_nodal_dvar")),
90 378 : _dpf_dvar_qp(getOptionalMaterialProperty<std::vector<Real>>(
91 : "dPorousFlow_effective_fluid_pressure_qp_dvar")),
92 :
93 267 : _uses_T(parameters.isParamSetByUser("T_coeff") && _T_coeff != 0.0),
94 378 : _temperature_nodal(getOptionalMaterialProperty<Real>("PorousFlow_temperature_nodal")),
95 378 : _temperature_qp(getOptionalMaterialProperty<Real>("PorousFlow_temperature_qp")),
96 189 : _dtemperature_dvar_nodal(
97 189 : getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_nodal_dvar")),
98 189 : _dtemperature_dvar_qp(
99 189 : getOptionalMaterialProperty<std::vector<Real>>("dPorousFlow_temperature_qp_dvar")),
100 567 : _zero_modifier(getParam<Real>("zero_modifier"))
101 : {
102 189 : }
103 :
104 : void
105 168 : PorousFlowPorosityLinear::initialSetup()
106 : {
107 168 : _pf = _nodal_material ? _pf_nodal.get() : _pf_qp.get();
108 168 : _dpf_dvar = _nodal_material ? _dpf_dvar_nodal.get() : _dpf_dvar_qp.get();
109 168 : 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 166 : _temperature = _nodal_material ? _temperature_nodal.get() : _temperature_qp.get();
116 166 : _dtemperature_dvar =
117 166 : _nodal_material ? _dtemperature_dvar_nodal.get() : _dtemperature_dvar_qp.get();
118 166 : if (!_temperature && _uses_T)
119 2 : mooseError(
120 : "PorousFlowPorosityLinear: When T_coeff is given you must have a temperature Material");
121 :
122 164 : 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 162 : PorousFlowPorosityBase::initialSetup();
127 162 : }
128 :
129 : void
130 111291 : PorousFlowPorosityLinear::initQpStatefulProperties()
131 : {
132 111291 : _porosity[_qp] = _phi_ref[0];
133 111291 : (*_dporosity_dvar)[_qp].assign(_num_var, 0.0);
134 111291 : (*_dporosity_dgradvar)[_qp].assign(_num_var, 0.0);
135 111291 : }
136 :
137 : void
138 720280 : PorousFlowPorosityLinear::computeQpProperties()
139 : {
140 : // note the [0] below: _phi_ref is a constant monomial and we use [0] regardless of
141 : // _nodal_material
142 720280 : Real porosity = _phi_ref[0];
143 720280 : if (_uses_pf)
144 720280 : porosity += _P_coeff * ((*_pf)[_qp] - _P_ref[_qp]);
145 720280 : if (_uses_T)
146 1940 : porosity += _T_coeff * ((*_temperature)[_qp] - _T_ref[_qp]);
147 720280 : 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 1940 : (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
155 1940 : porosity += _epv_coeff * (_vol_strain_qp[qp_to_use] - _epv_ref[0]);
156 : }
157 :
158 720280 : if (porosity < _porosity_min)
159 143 : _porosity[_qp] = _porosity_min;
160 : else
161 720137 : _porosity[_qp] = porosity;
162 :
163 720280 : (*_dporosity_dvar)[_qp].resize(_num_var);
164 720280 : (*_dporosity_dgradvar)[_qp].resize(_num_var);
165 2166324 : for (unsigned int pvar = 0; pvar < _num_var; ++pvar)
166 : {
167 : Real deriv = 0.0;
168 1446044 : if (_uses_pf)
169 1446044 : deriv += _P_coeff * (*_dpf_dvar)[_qp][pvar];
170 1446044 : if (_uses_T)
171 9364 : deriv += _T_coeff * (*_dtemperature_dvar)[_qp][pvar];
172 1446044 : if (porosity < _porosity_min)
173 436 : (*_dporosity_dvar)[_qp][pvar] = _zero_modifier * deriv;
174 : else
175 1445608 : (*_dporosity_dvar)[_qp][pvar] = deriv;
176 :
177 : RealGradient deriv_grad(0.0, 0.0, 0.0);
178 1446044 : if (_uses_volstrain)
179 : {
180 : const unsigned qp_to_use =
181 9364 : (_nodal_material && (_bnd || _strain_at_nearest_qp) ? nearestQP(_qp) : _qp);
182 9364 : deriv_grad += _epv_coeff * _dvol_strain_qp_dvar[qp_to_use][pvar];
183 : }
184 1446044 : if (porosity < _porosity_min)
185 436 : (*_dporosity_dgradvar)[_qp][pvar] = _zero_modifier * deriv_grad;
186 : else
187 1445608 : (*_dporosity_dgradvar)[_qp][pvar] = deriv_grad;
188 : }
189 720280 : }
|
