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 "PorousFlowHystereticRelativePermeabilityLiquid.h"
11 : #include "PorousFlowVanGenuchten.h"
12 :
13 : registerMooseObject("PorousFlowApp", PorousFlowHystereticRelativePermeabilityLiquid);
14 :
15 : InputParameters
16 534 : PorousFlowHystereticRelativePermeabilityLiquid::validParams()
17 : {
18 534 : InputParameters params = PorousFlowHystereticRelativePermeabilityBase::validParams();
19 1602 : params.addRangeCheckedParam<Real>(
20 : "liquid_modification_range",
21 1068 : 0.9,
22 : "liquid_modification_range > 0 & liquid_modification_range <= 1",
23 : "The wetting liquid relative permeability will be a cubic between S_l = "
24 : "liquid_modification_range * (1 - S_gr_Del) and S_l = 1.0 - 0.5 * S_gr_Del");
25 534 : params.addClassDescription(
26 : "PorousFlow material that computes relative permeability of the liquid phase in 1-phase or "
27 : "2-phase models that include hysteresis. You should ensure that the 'phase' for this "
28 : "Material does indeed represent the liquid phase");
29 534 : return params;
30 0 : }
31 :
32 416 : PorousFlowHystereticRelativePermeabilityLiquid::PorousFlowHystereticRelativePermeabilityLiquid(
33 416 : const InputParameters & parameters)
34 : : PorousFlowHystereticRelativePermeabilityBase(parameters),
35 414 : _liquid_phase(_phase_num),
36 414 : _liquid_modification_range(getParam<Real>("liquid_modification_range")),
37 414 : _kl_begin(_nodal_material ? declareProperty<Real>("PorousFlow_hysteresis_kr_l_begin_nodal")
38 810 : : declareProperty<Real>("PorousFlow_hysteresis_kr_l_begin_qp")),
39 414 : _klp_begin(_nodal_material ? declareProperty<Real>("PorousFlow_hysteresis_krp_l_begin_nodal")
40 810 : : declareProperty<Real>("PorousFlow_hysteresis_krp_l_begin_qp")),
41 414 : _kl_end(_nodal_material ? declareProperty<Real>("PorousFlow_hysteresis_kr_l_end_nodal")
42 810 : : declareProperty<Real>("PorousFlow_hysteresis_kr_l_end_qp")),
43 414 : _klp_end(_nodal_material ? declareProperty<Real>("PorousFlow_hysteresis_krp_l_end_nodal")
44 1226 : : declareProperty<Real>("PorousFlow_hysteresis_krp_l_end_qp"))
45 : {
46 414 : }
47 :
48 : void
49 38728 : PorousFlowHystereticRelativePermeabilityLiquid::computeRelPermQp()
50 : {
51 38728 : const Real sl = _saturation[_qp][_liquid_phase];
52 :
53 38728 : if (_hys_order[_qp] == 0)
54 : {
55 20342 : _relative_permeability[_qp] = PorousFlowVanGenuchten::relativePermeabilityHys(
56 20342 : sl, _s_lr, 0.0, _s_gr_max, 1.0, _m, _liquid_modification_range, 0.0, 0.0, 0.0, 0.0);
57 20342 : _drelative_permeability_ds[_qp] = PorousFlowVanGenuchten::drelativePermeabilityHys(
58 20342 : sl, _s_lr, 0.0, _s_gr_max, 1.0, _m, _liquid_modification_range, 0.0, 0.0, 0.0, 0.0);
59 : }
60 : else
61 : {
62 : // following ternary deals with the case where the turning-point saturation occurs in the
63 : // low-saturation region (tp_sat < _s_lr). There is "no hysteresis along the extension"
64 : // according to Doughty2008, so assume that the wetting curve is the same as would occur if
65 : // the turning-point saturation occured at _s_lr
66 : const Real effective_liquid_tp =
67 18386 : (_hys_sat_tps[_qp].at(0) < _s_lr) ? _s_lr : _hys_sat_tps[_qp].at(0);
68 18386 : const Real s_gas_max = (_hys_sat_tps[_qp].at(0) < _s_lr) ? _s_gr_max : _s_gr_tp0[_qp];
69 18386 : _relative_permeability[_qp] =
70 18386 : PorousFlowVanGenuchten::relativePermeabilityHys(sl,
71 : _s_lr,
72 : s_gas_max,
73 18386 : _s_gr_max,
74 : effective_liquid_tp,
75 18386 : _m,
76 18386 : _liquid_modification_range,
77 18386 : _kl_begin[_qp],
78 18386 : _klp_begin[_qp],
79 18386 : _kl_end[_qp],
80 18386 : _klp_end[_qp]);
81 18386 : _drelative_permeability_ds[_qp] =
82 18386 : PorousFlowVanGenuchten::drelativePermeabilityHys(sl,
83 18386 : _s_lr,
84 : s_gas_max,
85 18386 : _s_gr_max,
86 : effective_liquid_tp,
87 18386 : _m,
88 18386 : _liquid_modification_range,
89 18386 : _kl_begin[_qp],
90 18386 : _klp_begin[_qp],
91 18386 : _kl_end[_qp],
92 18386 : _klp_end[_qp]);
93 : }
94 38728 : }
95 :
96 : void
97 2006 : PorousFlowHystereticRelativePermeabilityLiquid::computeTurningPoint0Info(Real tp_sat)
98 : {
99 2006 : PorousFlowHystereticRelativePermeabilityBase::computeTurningPoint0Info(tp_sat);
100 :
101 : // following ternaries deal with the case where the turning-point saturation occurs in the
102 : // low-saturation region (tp_sat < _s_lr). There is "no hysteresis along the extension"
103 : // according to Doughty2008, so assume that the wetting curve is the same as would occur if
104 : // the turning-point saturation occured at _s_lr
105 2006 : const Real effective_liquid_tp = (tp_sat < _s_lr) ? _s_lr : tp_sat;
106 2006 : const Real s_gas_max = (tp_sat < _s_lr) ? _s_gr_max : _s_gr_tp0[_qp];
107 :
108 : // compute and record the wetting-curve value and derivative at s_begin
109 2006 : const Real s_begin = _liquid_modification_range * (1.0 - s_gas_max);
110 2006 : _kl_begin[_qp] = PorousFlowVanGenuchten::relativePermeabilityHys(s_begin,
111 : _s_lr,
112 : s_gas_max,
113 2006 : _s_gr_max,
114 : effective_liquid_tp,
115 2006 : _m,
116 : _liquid_modification_range,
117 : 0.0,
118 : 0.0,
119 : 0.0,
120 : 0.0);
121 2006 : _klp_begin[_qp] = PorousFlowVanGenuchten::drelativePermeabilityHys(s_begin,
122 2006 : _s_lr,
123 : s_gas_max,
124 2006 : _s_gr_max,
125 : effective_liquid_tp,
126 2006 : _m,
127 2006 : _liquid_modification_range,
128 : 0.0,
129 : 0.0,
130 : 0.0,
131 : 0.0);
132 : // compute and record the drying curve information at s_end (the drying curve is used because
133 : // s_end = 1.0 - 0.5 * s_gas_max in PorousFlowVanGenuchten::relativePermeabilityHys)
134 2006 : const Real s_end = 1.0 - 0.5 * s_gas_max;
135 2006 : _kl_end[_qp] = PorousFlowVanGenuchten::relativePermeabilityHys(s_end,
136 2006 : _s_lr,
137 : s_gas_max,
138 2006 : _s_gr_max,
139 : effective_liquid_tp,
140 2006 : _m,
141 2006 : _liquid_modification_range,
142 : 0.0,
143 : 0.0,
144 : 0.0,
145 : 0.0);
146 4012 : _klp_end[_qp] = PorousFlowVanGenuchten::drelativePermeabilityHys(s_end,
147 2006 : _s_lr,
148 : s_gas_max,
149 2006 : _s_gr_max,
150 : effective_liquid_tp,
151 2006 : _m,
152 2006 : _liquid_modification_range,
153 : 0.0,
154 : 0.0,
155 : 0.0,
156 : 0.0);
157 2006 : }
|
