137 enum { conti0EqIdx = Indices::conti0EqIdx };
138 enum { dimWorld = GridView::dimensionworld };
139 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
140 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
142 static constexpr unsigned contiEnergyEqIdx = Indices::contiEnergyEqIdx;
146 struct ConvectiveMixingModuleParam
148 std::vector<bool> active_;
149 std::vector<Scalar> Xhi_;
150 std::vector<Scalar> Psi_;
154 static void beginEpisode(
const EclipseState& eclState,
157 ConvectiveMixingModuleParam& info)
160 std::size_t
numRegions = eclState.runspec().tabdims().getNumPVTTables();
162 if (info.active_.empty()) {
168 info.active_[i] =
control.drsdtConvective(i);
169 if (
control.drsdtConvective(i)) {
170 info.Xhi_[i] =
control.getMaxDRSDT(i);
171 info.Psi_[i] =
control.getPsi(i);
177 static void modifyAvgDensity(Evaluation&
rhoAvg,
181 const ConvectiveMixingModuleParam& info) {
183 if (info.active_.empty()) {
190 if (
phaseIdx == FluidSystem::gasPhaseIdx) {
194 const auto& liquidPhaseIdx = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
195 FluidSystem::waterPhaseIdx :
196 FluidSystem::oilPhaseIdx;
199 const auto&
t_in =
intQuantsIn.fluidState().temperature(liquidPhaseIdx);
203 const auto&
bLiquidIn = (FluidSystem::phaseIsActive(waterPhaseIdx)) ?
208 FluidSystem::waterPvt().waterReferenceDensity(
intQuantsIn.pvtRegionIndex()) :
213 const auto t_ex = Toolbox::value(
intQuantsEx.fluidState().temperature(liquidPhaseIdx));
214 const auto p_ex = Toolbox::value(
intQuantsEx.fluidState().pressure(liquidPhaseIdx));
217 const auto bLiquidEx = (FluidSystem::phaseIsActive(waterPhaseIdx)) ?
218 FluidSystem::waterPvt().inverseFormationVolumeFactor(
intQuantsEx.pvtRegionIndex(),
224 FluidSystem::waterPvt().waterReferenceDensity(
intQuantsEx.pvtRegionIndex()) :
232 template <
class Context>
233 static void addConvectiveMixingFlux(RateVector&
flux,
234 const Context& elemCtx,
238 const auto& problem = elemCtx.problem();
239 const auto& stencil = elemCtx.stencil(
timeIdx);
249 Scalar faceArea =
scvf.area();
250 const Scalar
g = problem.gravity()[dimWorld - 1];
256 addConvectiveMixingFlux(
flux,
264 problem.moduleParams().convectiveMixingModuleParam);
278 const Scalar faceArea,
279 const ConvectiveMixingModuleParam& info)
281 if (info.active_.empty()) {
289 const auto& liquidPhaseIdx = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
290 FluidSystem::waterPhaseIdx :
291 FluidSystem::oilPhaseIdx;
292 const Evaluation
SoMax = 0.0;
295 const auto&
t_in =
intQuantsIn.fluidState().temperature(liquidPhaseIdx);
303 const auto bLiquidSatIn = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
307 const auto&
densityLiquidIn = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
308 FluidSystem::waterPvt().waterReferenceDensity(
intQuantsIn.pvtRegionIndex()) :
309 FluidSystem::oilPvt().oilReferenceDensity(
intQuantsIn.pvtRegionIndex());
323 const auto bLiquidSatEx = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
328 const auto&
densityLiquidEx = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
329 FluidSystem::waterPvt().waterReferenceDensity(
intQuantsEx.pvtRegionIndex()) :
330 FluidSystem::oilPvt().oilReferenceDensity(
intQuantsEx.pvtRegionIndex());
354 const auto&
Rsup = (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) ?
355 up.fluidState().Rsw() :
356 up.fluidState().Rs();
358 const Evaluation&
transMult =
up.rockCompTransMultiplier();
359 const auto&
invB =
up.fluidState().invB(liquidPhaseIdx);
360 const auto&
visc =
up.fluidState().viscosity(liquidPhaseIdx);
366 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
372 if constexpr (enableEnergy) {
373 const auto&
h =
up.fluidState().enthalpy(liquidPhaseIdx) * FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx,
up.pvtRegionIndex());
static void addConvectiveMixingFlux(RateVector &, const IntensiveQuantities &, const IntensiveQuantities &, const unsigned, const unsigned, const Scalar, const Scalar, const Scalar, const ConvectiveMixingModuleParam &)
Adds the convective mixing mass flux flux to the flux vector over a flux integration point.
Definition blackoilconvectivemixingmodule.hh:113
static void addConvectiveMixingFlux(RateVector &flux, const IntensiveQuantities &intQuantsIn, const IntensiveQuantities &intQuantsEx, const unsigned globalIndexIn, const unsigned globalIndexEx, const Scalar distZg, const Scalar trans, const Scalar faceArea, const ConvectiveMixingModuleParam &info)
Adds the convective mixing mass flux flux to the flux vector over a flux integration point.
Definition blackoilconvectivemixingmodule.hh:271
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235