28#ifndef EWOMS_BLACK_OIL_ENERGY_MODULE_HH
29#define EWOMS_BLACK_OIL_ENERGY_MODULE_HH
36#include <opm/material/common/Tabulated1DFunction.hpp>
38#include <opm/material/common/Valgrind.hpp>
40#include <dune/common/fvector.hh>
65 static constexpr unsigned temperatureIdx = Indices::temperatureIdx;
66 static constexpr unsigned contiEnergyEqIdx = Indices::contiEnergyEqIdx;
70 static constexpr unsigned numPhases = FluidSystem::numPhases;
79 if constexpr (enableEnergy)
89 if constexpr (enableEnergy)
93 static bool primaryVarApplies(
unsigned pvIdx)
95 if constexpr (enableEnergy)
96 return pvIdx == temperatureIdx;
106 return "temperature";
114 return static_cast<Scalar
>(1.0);
117 static bool eqApplies(
unsigned eqIdx)
119 if constexpr (enableEnergy)
120 return eqIdx == contiEnergyEqIdx;
129 return "conti^energy";
140 template <
class LhsEval>
141 static void addStorage(Dune::FieldVector<LhsEval, numEq>&
storage,
144 if constexpr (enableEnergy) {
150 if (!FluidSystem::phaseIsActive(
phaseIdx))
157 storage[contiEnergyEqIdx] += poro*
S*
u*rho;
161 Scalar rockFraction =
intQuants.rockFraction();
173 if constexpr (enableEnergy) {
174 flux[contiEnergyEqIdx] = 0.0;
179 if (!FluidSystem::phaseIsActive(
phaseIdx))
195 static void addHeatFlux(RateVector&
flux,
198 if constexpr (enableEnergy) {
207 template <
class UpEval,
class Eval,
class Flu
idState>
208 static void addPhaseEnthalpyFluxes_(RateVector&
flux,
210 const Eval& volumeFlux,
211 const FluidState&
upFs)
213 flux[contiEnergyEqIdx] +=
219 template <
class UpstreamEval>
220 static void addPhaseEnthalpyFlux_(RateVector&
flux,
229 const auto&
fs =
up.fluidState();
237 static void addToEnthalpyRate(RateVector&
flux,
238 const Evaluation&
hRate)
240 if constexpr (enableEnergy)
250 if constexpr (enableEnergy)
251 priVars[temperatureIdx] = temperatureIdx;
257 template <
class Flu
idState>
259 const FluidState& fluidState)
261 if constexpr (enableEnergy)
262 priVars[temperatureIdx] = fluidState.temperature(0);
269 const PrimaryVariables&
oldPv,
270 const EqVector&
delta)
272 if constexpr (enableEnergy)
286 return static_cast<Scalar
>(0.0);
298 template <
class DofEntity>
301 if constexpr (enableEnergy) {
302 unsigned dofIdx = model.dofMapper().index(
dof);
303 const PrimaryVariables& priVars = model.solution(0)[
dofIdx];
308 template <
class DofEntity>
311 if constexpr (enableEnergy) {
312 unsigned dofIdx = model.dofMapper().index(
dof);
349 static constexpr int temperatureIdx = Indices::temperatureIdx;
350 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
362 auto&
fs = asImp_().fluidState_;
366 fs.setTemperature(priVars.makeEvaluation(temperatureIdx,
timeIdx,
elemCtx.linearizationType()));
374 const PrimaryVariables& priVars,
379 auto&
fs = asImp_().fluidState_;
380 fs.setTemperature(priVars.makeEvaluation(temperatureIdx,
timeIdx,
lintype));
392 auto&
fs = asImp_().fluidState_;
397 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
406 rockInternalEnergy_ = SolidEnergyLaw::solidInternalEnergy(solidEnergyLawParams,
fs);
409 totalThermalConductivity_ = ThermalConductionLaw::thermalConductivity(thermalConductionLawParams,
fs);
420 const Evaluation& rockInternalEnergy()
const
421 {
return rockInternalEnergy_; }
423 const Evaluation& totalThermalConductivity()
const
424 {
return totalThermalConductivity_; }
426 const Scalar& rockFraction()
const
427 {
return rockFraction_; }
430 Implementation& asImp_()
431 {
return *
static_cast<Implementation*
>(
this); }
433 Evaluation rockInternalEnergy_;
434 Evaluation totalThermalConductivity_;
435 Scalar rockFraction_;
438template <
class TypeTag>
456 if constexpr (enableTemperature) {
459 auto&
fs = asImp_().fluidState_;
461 fs.setTemperature(T);
465 template<
class Problem>
473 if constexpr (enableTemperature) {
474 auto&
fs = asImp_().fluidState_;
478 fs.setTemperature(T);
488 const Evaluation& rockInternalEnergy()
const
489 {
throw std::logic_error(
"Requested the rock internal energy, which is "
490 "unavailable because energy is not conserved"); }
492 const Evaluation& totalThermalConductivity()
const
493 {
throw std::logic_error(
"Requested the total thermal conductivity, which is "
494 "unavailable because energy is not conserved"); }
497 Implementation& asImp_()
498 {
return *
static_cast<Implementation*
>(
this); }
526 static const int dimWorld = GridView::dimensionworld;
527 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
528 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
530 template<
class Flu
idState>
531 static void updateEnergy(Evaluation& energyFlux,
532 const unsigned& focusDofIndex,
533 const unsigned&
inIdx,
534 const unsigned&
exIdx,
535 const IntensiveQuantities&
inIq,
536 const IntensiveQuantities&
exIq,
537 const FluidState&
inFs,
538 const FluidState&
exFs,
539 const Scalar& inAlpha,
540 const Scalar& outAlpha,
541 const Scalar& faceArea)
544 if (focusDofIndex ==
inIdx)
547 -
inFs.temperature(0);
548 else if (focusDofIndex ==
exIdx)
558 if (focusDofIndex ==
inIdx)
564 if (focusDofIndex ==
exIdx)
582 energyFlux =
deltaT * (-
H/faceArea);
585 void updateEnergy(
const ElementContext&
elemCtx,
592 const Scalar faceArea =
scvf.area();
593 const unsigned inIdx =
scvf.interiorIndex();
594 const unsigned exIdx =
scvf.exteriorIndex();
597 const auto&
inFs =
inIq.fluidState();
598 const auto&
exFs =
exIq.fluidState();
601 updateEnergy(energyFlux_,
614 template <
class Context,
class BoundaryFlu
idState>
626 const Scalar alpha =
ctx.problem().thermalHalfTransmissibilityBoundary(
ctx,
scvfIdx);
630 template <
class BoundaryFlu
idState>
631 static void updateEnergyBoundary(Evaluation& energyFlux,
632 const IntensiveQuantities&
inIq,
633 unsigned focusDofIndex,
638 const auto&
inFs =
inIq.fluidState();
640 if (focusDofIndex ==
inIdx)
643 -
inFs.temperature(0);
650 if (focusDofIndex ==
inIdx)
667 const Evaluation& energyFlux()
const
668 {
return energyFlux_; }
671 Implementation& asImp_()
672 {
return *
static_cast<Implementation*
>(
this); }
674 Evaluation energyFlux_;
677template <
class TypeTag>
685 template<
class Flu
idState>
686 static void updateEnergy(Evaluation& ,
690 const IntensiveQuantities& ,
691 const IntensiveQuantities& ,
699 void updateEnergy(
const ElementContext&,
704 template <
class Context,
class BoundaryFlu
idState>
705 void updateEnergyBoundary(
const Context&,
711 template <
class BoundaryFlu
idState>
712 static void updateEnergyBoundary(Evaluation& ,
713 const IntensiveQuantities& ,
720 const Evaluation& energyFlux()
const
721 {
throw std::logic_error(
"Requested the energy flux, but energy is not conserved"); }
Declares the properties required by the black oil model.
Provides the energy specific extensive quantities to the generic black-oil module's extensive quantit...
Definition blackoilenergymodules.hh:511
Provides the volumetric quantities required for the equations needed by the energys extension of the ...
Definition blackoilenergymodules.hh:333
void updateTemperature_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the temperature of the intensive quantity's fluid state.
Definition blackoilenergymodules.hh:358
void updateEnergyQuantities_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, const typename FluidSystem::template ParameterCache< Evaluation > ¶mCache)
Compute the intensive quantities needed to handle energy conservation.
Definition blackoilenergymodules.hh:387
void updateTemperature_(const Problem &problem, const PrimaryVariables &priVars, unsigned globalDofIdx, const unsigned timeIdx, const LinearizationType &lintype)
Update the temperature of the intensive quantity's fluid state.
Definition blackoilenergymodules.hh:373
Contains the high level supplements required to extend the black oil model by energy.
Definition blackoilenergymodules.hh:52
static void registerParameters()
Register all run-time parameters for the black-oil energy module.
Definition blackoilenergymodules.hh:77
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition blackoilenergymodules.hh:292
static void updatePrimaryVars(PrimaryVariables &newPv, const PrimaryVariables &oldPv, const EqVector &delta)
Do a Newton-Raphson update the primary variables of the energys.
Definition blackoilenergymodules.hh:268
static void registerOutputModules(Model &model, Simulator &simulator)
Register all energy specific VTK and ECL output modules.
Definition blackoilenergymodules.hh:86
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar)
Assign the energy specific primary variables to a PrimaryVariables object.
Definition blackoilenergymodules.hh:247
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition blackoilenergymodules.hh:280
static void assignPrimaryVars(PrimaryVariables &priVars, const FluidState &fluidState)
Assign the energy specific primary variables to a PrimaryVariables object.
Definition blackoilenergymodules.hh:258
VTK output module for the black oil model's energy related quantities.
Definition vtkblackoilenergymodule.hpp:54
The common code for the linearizers of non-linear systems of equations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
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
This method contains all callback classes for quantities that are required by some extensive quantiti...
Definition linearizationtype.hh:35
VTK output module for the black oil model's energy related quantities.