blackoillocalresidual.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH
29#define EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH
30
31#include <opm/material/common/MathToolbox.hpp>
32#include <opm/material/fluidstates/BlackOilFluidState.hpp>
33
44
45#include <cassert>
46
47namespace Opm {
48
54template <class TypeTag>
55class BlackOilLocalResidual : public GetPropType<TypeTag, Properties::DiscLocalResidual>
56{
66
67 enum { conti0EqIdx = Indices::conti0EqIdx };
68 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
69 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
70 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
71
72 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
73 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
74 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
75
76 enum { gasCompIdx = FluidSystem::gasCompIdx };
77 enum { oilCompIdx = FluidSystem::oilCompIdx };
78 enum { waterCompIdx = FluidSystem::waterCompIdx };
79 enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
80
81 static constexpr bool waterEnabled = Indices::waterEnabled;
82 static constexpr bool gasEnabled = Indices::gasEnabled;
83 static constexpr bool oilEnabled = Indices::oilEnabled;
84 static constexpr bool compositionSwitchEnabled = (compositionSwitchIdx >= 0);
85
86 static constexpr bool blackoilConserveSurfaceVolume = getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>();
87 static constexpr bool enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>();
88 static constexpr bool enableDiffusion = getPropValue<TypeTag, Properties::EnableDiffusion>();
89 static constexpr bool enableConvectiveMixing = getPropValue<TypeTag, Properties::EnableConvectiveMixing>();
90
91 using Toolbox = MathToolbox<Evaluation>;
101
102public:
106 template <class LhsEval>
107 void computeStorage(Dune::FieldVector<LhsEval, numEq>& storage,
108 const ElementContext& elemCtx,
109 unsigned dofIdx,
110 unsigned timeIdx) const
111 {
112 // retrieve the intensive quantities for the SCV at the specified point in time
113 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
114 const auto& fs = intQuants.fluidState();
115
116 storage = 0.0;
117
118 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
119 if (!FluidSystem::phaseIsActive(phaseIdx)) {
120 if (Indices::numPhases == 3) { // add trivial equation for the pseudo phase
121 const unsigned activeCompIdx =
122 Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
123 if (timeIdx == 0) {
124 storage[conti0EqIdx + activeCompIdx] = variable<LhsEval>(0.0, conti0EqIdx + activeCompIdx);
125 }
126 else {
127 storage[conti0EqIdx + activeCompIdx] = 0.0;
128 }
129 }
130 continue;
131 }
132
133 const unsigned activeCompIdx =
134 Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
135 const LhsEval surfaceVolume =
136 Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx)) *
137 Toolbox::template decay<LhsEval>(fs.invB(phaseIdx)) *
138 Toolbox::template decay<LhsEval>(intQuants.porosity());
139
140 storage[conti0EqIdx + activeCompIdx] += surfaceVolume;
141
142 // account for dissolved gas
143 if (phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
144 const unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
145 storage[conti0EqIdx + activeGasCompIdx] +=
146 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rs()) *
147 surfaceVolume;
148 }
149
150 // account for dissolved gas in water phase
151 if (phaseIdx == waterPhaseIdx && FluidSystem::enableDissolvedGasInWater()) {
152 const unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
153 storage[conti0EqIdx + activeGasCompIdx] +=
154 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rsw()) *
155 surfaceVolume;
156 }
157
158 // account for vaporized oil
159 if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
160 const unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
161 storage[conti0EqIdx + activeOilCompIdx] +=
162 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rv()) *
163 surfaceVolume;
164 }
165
166 // account for vaporized water
167 if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedWater()) {
168 const unsigned activeWaterCompIdx = Indices::canonicalToActiveComponentIndex(waterCompIdx);
169 storage[conti0EqIdx + activeWaterCompIdx] +=
170 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rvw()) *
171 surfaceVolume;
172 }
173 }
174
175 adaptMassConservationQuantities_(storage, intQuants.pvtRegionIndex());
176
177 // deal with solvents (if present)
178 SolventModule::addStorage(storage, intQuants);
179
180 // deal with zFracton (if present)
181 ExtboModule::addStorage(storage, intQuants);
182
183 // deal with polymer (if present)
184 PolymerModule::addStorage(storage, intQuants);
185
186 // deal with energy (if present)
187 EnergyModule::addStorage(storage, intQuants);
188
189 // deal with foam (if present)
190 FoamModule::addStorage(storage, intQuants);
191
192 // deal with salt (if present)
193 BrineModule::addStorage(storage, intQuants);
194
195 // deal with micp (if present)
196 MICPModule::addStorage(storage, intQuants);
197 }
198
202 void computeFlux(RateVector& flux,
203 const ElementContext& elemCtx,
204 unsigned scvfIdx,
205 unsigned timeIdx) const
206 {
207 assert(timeIdx == 0);
208
209 flux = 0.0;
210
211 const ExtensiveQuantities& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
212 unsigned focusDofIdx = elemCtx.focusDofIndex();
213 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
214 if (!FluidSystem::phaseIsActive(phaseIdx)) {
215 continue;
216 }
217
218 const unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
219 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
220 const unsigned pvtRegionIdx = up.pvtRegionIndex();
221 if (upIdx == focusDofIdx) {
222 evalPhaseFluxes_<Evaluation>(flux, phaseIdx, pvtRegionIdx, extQuants, up.fluidState());
223 }
224 else {
225 evalPhaseFluxes_<Scalar>(flux, phaseIdx, pvtRegionIdx, extQuants, up.fluidState());
226 }
227 }
228 // deal with solvents (if present)
229 SolventModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
230
231 // deal with zFracton (if present)
232 ExtboModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
233
234 // deal with polymer (if present)
235 PolymerModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
236
237 // deal with energy (if present)
238 EnergyModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
239
240 // deal with foam (if present)
241 FoamModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
242
243 // deal with salt (if present)
244 BrineModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
245
246 // deal with micp (if present)
247 MICPModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
248
249 DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
250
251 // deal with convective mixing (if enabled)
252 ConvectiveMixingModule::addConvectiveMixingFlux(flux,elemCtx, scvfIdx, timeIdx);
253 }
254
258 void computeSource(RateVector& source,
259 const ElementContext& elemCtx,
260 unsigned dofIdx,
261 unsigned timeIdx) const
262 {
263 // retrieve the source term intrinsic to the problem
264 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
265
266 // deal with MICP (if present)
267 MICPModule::addSource(source, elemCtx, dofIdx, timeIdx);
268
269 // scale the source term of the energy equation
270 if constexpr (enableEnergy) {
271 source[Indices::contiEnergyEqIdx] *=
272 getPropValue<TypeTag, Properties::BlackOilEnergyScalingFactor>();
273 }
274 }
275
280 template <class UpEval, class FluidState>
281 static void evalPhaseFluxes_(RateVector& flux,
282 unsigned phaseIdx,
283 unsigned pvtRegionIdx,
284 const ExtensiveQuantities& extQuants,
285 const FluidState& upFs)
286 {
287 const auto& invB = getInvB_<FluidSystem, FluidState, UpEval>(upFs, phaseIdx, pvtRegionIdx);
288 const auto& surfaceVolumeFlux = invB*extQuants.volumeFlux(phaseIdx);
289 const unsigned activeCompIdx =
290 Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
291
292 if constexpr (blackoilConserveSurfaceVolume) {
293 flux[conti0EqIdx + activeCompIdx] += surfaceVolumeFlux;
294 }
295 else {
296 flux[conti0EqIdx + activeCompIdx] += surfaceVolumeFlux *
297 FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
298 }
299
300 if (phaseIdx == oilPhaseIdx) {
301 // dissolved gas (in the oil phase).
302 if (FluidSystem::enableDissolvedGas()) {
303 const auto& Rs = BlackOil::getRs_<FluidSystem, FluidState, UpEval>(upFs, pvtRegionIdx);
304
305 const unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
306 if constexpr (blackoilConserveSurfaceVolume) {
307 flux[conti0EqIdx + activeGasCompIdx] += Rs * surfaceVolumeFlux;
308 }
309 else {
310 flux[conti0EqIdx + activeGasCompIdx] +=
311 Rs * surfaceVolumeFlux *
312 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
313 }
314 }
315 } else if (phaseIdx == waterPhaseIdx) {
316 // dissolved gas (in the water phase).
317 if (FluidSystem::enableDissolvedGasInWater()) {
318 const auto& Rsw = BlackOil::getRsw_<FluidSystem, FluidState, UpEval>(upFs, pvtRegionIdx);
319
320 const unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
321 if constexpr (blackoilConserveSurfaceVolume) {
322 flux[conti0EqIdx + activeGasCompIdx] += Rsw * surfaceVolumeFlux;
323 }
324 else {
325 flux[conti0EqIdx + activeGasCompIdx] +=
326 Rsw * surfaceVolumeFlux *
327 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
328 }
329 }
330 }
331 else if (phaseIdx == gasPhaseIdx) {
332 // vaporized oil (in the gas phase).
333 if (FluidSystem::enableVaporizedOil()) {
334 const auto& Rv = BlackOil::getRv_<FluidSystem, FluidState, UpEval>(upFs, pvtRegionIdx);
335
336 const unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
337 if constexpr (blackoilConserveSurfaceVolume) {
338 flux[conti0EqIdx + activeOilCompIdx] += Rv * surfaceVolumeFlux;
339 }
340 else {
341 flux[conti0EqIdx + activeOilCompIdx] +=
342 Rv * surfaceVolumeFlux *
343 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
344 }
345 }
346 // vaporized water (in the gas phase).
347 if (FluidSystem::enableVaporizedWater()) {
348 const auto& Rvw = BlackOil::getRvw_<FluidSystem, FluidState, UpEval>(upFs, pvtRegionIdx);
349
350 const unsigned activeWaterCompIdx = Indices::canonicalToActiveComponentIndex(waterCompIdx);
351 if constexpr (blackoilConserveSurfaceVolume) {
352 flux[conti0EqIdx + activeWaterCompIdx] += Rvw * surfaceVolumeFlux;
353 }
354 else {
355 flux[conti0EqIdx + activeWaterCompIdx] +=
356 Rvw * surfaceVolumeFlux *
357 FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
358 }
359 }
360 }
361 }
362
374 template <class Scalar>
375 static void adaptMassConservationQuantities_(Dune::FieldVector<Scalar, numEq>& container,
376 unsigned pvtRegionIdx)
377 {
378 if constexpr (!blackoilConserveSurfaceVolume) {
379 // convert "surface volume" to mass. this is complicated a bit by the fact that
380 // not all phases are necessarily enabled. (we here assume that if a fluid phase
381 // is disabled, its respective "main" component is not considered as well.)
382
383 if constexpr (waterEnabled) {
384 const unsigned activeWaterCompIdx = Indices::canonicalToActiveComponentIndex(waterCompIdx);
385 container[conti0EqIdx + activeWaterCompIdx] *=
386 FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
387 }
388
389 if constexpr (gasEnabled) {
390 const unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
391 container[conti0EqIdx + activeGasCompIdx] *=
392 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
393 }
394
395 if constexpr (oilEnabled) {
396 const unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
397 container[conti0EqIdx + activeOilCompIdx] *=
398 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
399 }
400 }
401 }
402};
403
404} // namespace Opm
405
406#endif
Contains the classes required to extend the black-oil model by brine.
Classes required for dynamic convective mixing.
Classes required for molecular diffusion.
Contains the classes required to extend the black-oil model by energy.
Contains the classes required to extend the black-oil model by solvent component. For details,...
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by MICP.
Contains the classes required to extend the black-oil model by polymer.
Declares the properties required by the black oil model.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:55
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilbrinemodules.hh:195
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilbrinemodules.hh:162
Definition: blackoilconvectivemixingmodule.hh:64
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: blackoildiffusionmodule.hh:50
Contains the high level supplements required to extend the black oil model by energy.
Definition: blackoilenergymodules.hh:58
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilenergymodules.hh:181
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilenergymodules.hh:153
Contains the high level supplements required to extend the black oil model.
Definition: blackoilextbomodules.hh:62
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilextbomodules.hh:162
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilextbomodules.hh:197
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition: blackoilfoammodules.hh:58
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilfoammodules.hh:166
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilfoammodules.hh:206
Calculates the local residual of the black oil model.
Definition: blackoillocalresidual.hh:56
void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume.
Definition: blackoillocalresidual.hh:202
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: blackoillocalresidual.hh:258
static void evalPhaseFluxes_(RateVector &flux, unsigned phaseIdx, unsigned pvtRegionIdx, const ExtensiveQuantities &extQuants, const FluidState &upFs)
Helper function to calculate the flux of mass in terms of conservation quantities via specific fluid ...
Definition: blackoillocalresidual.hh:281
static void adaptMassConservationQuantities_(Dune::FieldVector< Scalar, numEq > &container, unsigned pvtRegionIdx)
Helper function to convert the mass-related parts of a Dune::FieldVector that stores conservation qua...
Definition: blackoillocalresidual.hh:375
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g. phase mass) within a finite sub-control volume.
Definition: blackoillocalresidual.hh:107
Contains the high level supplements required to extend the black oil model by MICP.
Definition: blackoilmicpmodules.hh:54
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilmicpmodules.hh:207
static void addSource(RateVector &source, const Problem &problem, const IntensiveQuantities &intQuants, unsigned globalSpaceIdex)
Definition: blackoilmicpmodules.hh:243
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilmicpmodules.hh:144
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilpolymermodules.hh:238
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:277
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:66
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilsolventmodules.hh:183
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilsolventmodules.hh:216
Definition: blackoilboundaryratevector.hh:39
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:233