29 #ifndef EWOMS_ENERGY_MODULE_HH 30 #define EWOMS_ENERGY_MODULE_HH 35 #include <opm/common/Valgrind.hpp> 36 #include <opm/common/Unused.hpp> 37 #include <opm/common/ErrorMacros.hpp> 38 #include <opm/common/Exceptions.hpp> 40 #include <dune/common/fvector.hh> 45 namespace Properties {
58 template <
class TypeTag,
bool enableEnergy>
64 template <
class TypeTag>
68 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
69 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
70 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
71 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
72 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
73 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
78 typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
100 static std::string
eqName(
unsigned eqIdx OPM_UNUSED)
108 unsigned globalDofIdx OPM_UNUSED,
109 unsigned pvIdx OPM_UNUSED)
115 static Scalar
eqWeight(
const Model& model OPM_UNUSED,
116 unsigned globalDofIdx OPM_UNUSED,
117 unsigned eqIdx OPM_UNUSED)
123 template <
class Flu
idState>
125 const FluidState& fs OPM_UNUSED)
132 template <
class RateVector,
class Flu
idState>
134 const FluidState& fluidState OPM_UNUSED,
135 unsigned phaseIdx OPM_UNUSED,
136 const Evaluation& volume OPM_UNUSED)
143 const Evaluation& rate OPM_UNUSED)
150 const Evaluation& rate OPM_UNUSED)
163 template <
class LhsEval>
165 const IntensiveQuantities& intQuants OPM_UNUSED,
166 unsigned phaseIdx OPM_UNUSED)
173 template <
class LhsEval,
class Scv>
175 const IntensiveQuantities& intQuants OPM_UNUSED,
176 const Scv& scv OPM_UNUSED,
177 unsigned phaseIdx OPM_UNUSED)
184 template <
class LhsEval>
186 const IntensiveQuantities& intQuants OPM_UNUSED)
195 template <
class Context>
197 const Context& context OPM_UNUSED,
198 unsigned spaceIdx OPM_UNUSED,
199 unsigned timeIdx OPM_UNUSED)
207 template <
class Context>
209 const Context& context OPM_UNUSED,
210 unsigned spaceIdx OPM_UNUSED,
211 unsigned timeIdx OPM_UNUSED)
220 template <
class Context>
222 const Context& context OPM_UNUSED,
223 unsigned spaceIdx OPM_UNUSED,
224 unsigned timeIdx OPM_UNUSED)
231 template <
class TypeTag>
235 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
236 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
238 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
239 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
240 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
241 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
246 enum { numPhases = FluidSystem::numPhases };
247 enum { energyEqIdx = Indices::energyEqIdx };
248 enum { temperatureIdx = Indices::temperatureIdx };
250 typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
251 typedef Opm::MathToolbox<Evaluation> Toolbox;
267 if (pvIdx == temperatureIdx)
268 return "temperature";
277 static std::string
eqName(
unsigned eqIdx)
279 if (eqIdx == energyEqIdx)
290 if (pvIdx != temperatureIdx)
294 return std::max(1.0/1000, 1.0/model.solution(0)[globalDofIdx][temperatureIdx]);
300 static Scalar
eqWeight(
const Model& model OPM_UNUSED,
301 unsigned globalDofIdx OPM_UNUSED,
304 if (eqIdx != energyEqIdx)
308 return 1.0 / 1.0035e3;
315 { rateVec[energyEqIdx] = rate; }
321 { rateVec[energyEqIdx] += rate; }
328 {
return -extQuants.temperatureGradNormal() * extQuants.heatConductivity(); }
334 template <
class RateVector,
class Flu
idState>
336 const FluidState& fluidState,
338 const Evaluation& volume)
340 rateVec[energyEqIdx] =
342 * fluidState.density(phaseIdx)
343 * fluidState.enthalpy(phaseIdx);
349 template <
class Flu
idState>
351 const FluidState& fs)
353 priVars[temperatureIdx] = Toolbox::value(fs.temperature(0));
355 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
356 assert(std::abs(Toolbox::value(fs.temperature(0))
357 - Toolbox::value(fs.temperature(phaseIdx))) < 1e-30);
366 template <
class LhsEval>
368 const IntensiveQuantities& intQuants,
371 const auto& fs = intQuants.fluidState();
372 storage[energyEqIdx] +=
373 Toolbox::template decay<LhsEval>(fs.density(phaseIdx))
374 * Toolbox::template decay<LhsEval>(fs.internalEnergy(phaseIdx))
375 * Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
376 * Toolbox::template decay<LhsEval>(intQuants.porosity());
383 template <
class Scv,
class LhsEval>
385 const IntensiveQuantities& intQuants,
389 const auto& fs = intQuants.fractureFluidState();
390 storage[energyEqIdx] +=
391 Toolbox::template decay<LhsEval>(fs.density(phaseIdx))
392 * Toolbox::template decay<LhsEval>(fs.internalEnergy(phaseIdx))
393 * Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
394 * Toolbox::template decay<LhsEval>(intQuants.fracturePorosity())
395 * Toolbox::template decay<LhsEval>(intQuants.fractureVolume())/scv.volume();
402 template <
class LhsEval>
404 const IntensiveQuantities& intQuants)
406 storage[energyEqIdx] +=
407 Toolbox::template decay<LhsEval>(intQuants.heatCapacitySolid())
408 * Toolbox::template decay<LhsEval>(intQuants.fluidState().temperature(0));
417 template <
class Context>
419 const Context& context,
423 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
426 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
427 if (!context.model().phaseIsConsidered(phaseIdx))
431 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
432 const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
435 extQuants.volumeFlux(phaseIdx)
436 * up.fluidState().enthalpy(phaseIdx)
437 * up.fluidState().density(phaseIdx);
446 template <
class Context>
448 const Context& context,
452 const auto& scvf = context.stencil(timeIdx).interiorFace(spaceIdx);
453 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
458 1 - extQuants.fractureWidth()/(2*scvf.area());
461 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
462 if (!context.model().phaseIsConsidered(phaseIdx))
466 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
467 const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
470 extQuants.fractureVolumeFlux(phaseIdx)
471 * up.fluidState().enthalpy(phaseIdx)
472 * up.fluidState().density(phaseIdx);
482 template <
class Context>
484 const Context& context,
488 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
492 - extQuants.temperatureGradNormal()
493 * extQuants.heatConductivity();
503 template <
unsigned PVOffset,
bool enableEnergy>
509 template <
unsigned PVOffset>
519 template <
unsigned PVOffset>
523 enum { temperatureIdx = PVOffset };
526 enum { energyEqIdx = PVOffset };
538 template <
class TypeTag,
bool enableEnergy>
544 template <
class TypeTag>
548 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
549 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
550 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
552 typedef Opm::MathToolbox<Evaluation> Toolbox;
562 OPM_THROW(std::logic_error,
"Method heatCapacitySolid() does not make " 563 "sense for isothermal models");
573 OPM_THROW(std::logic_error,
"Method heatConductivity() does not make " 574 "sense for isothermal models");
581 template <
class Flu
idState,
class Context>
583 const Context& context,
587 Scalar T = context.problem().temperature(context, spaceIdx, timeIdx);
588 fluidState.setTemperature(Toolbox::createConstant(T));
595 template <
class Flu
idState>
597 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache OPM_UNUSED,
598 const ElementContext& elemCtx OPM_UNUSED,
599 unsigned dofIdx OPM_UNUSED,
600 unsigned timeIdx OPM_UNUSED)
607 template <
class TypeTag>
611 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
612 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
613 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
614 typedef typename GET_PROP_TYPE(TypeTag, HeatConductionLaw) HeatConductionLaw;
617 enum { numPhases = FluidSystem::numPhases };
618 enum { energyEqIdx = Indices::energyEqIdx };
619 enum { temperatureIdx = Indices::temperatureIdx };
621 typedef Opm::MathToolbox<Evaluation> Toolbox;
627 template <
class Flu
idState,
class Context>
629 const Context& context,
633 const auto& priVars = context.primaryVars(spaceIdx, timeIdx);
636 val = Toolbox::createVariable(priVars[temperatureIdx], temperatureIdx);
638 val = Toolbox::createConstant(priVars[temperatureIdx]);
640 fluidState.setTemperature(val);
647 template <
class Flu
idState>
649 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
650 const ElementContext& elemCtx,
655 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
656 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
659 fs.setEnthalpy(phaseIdx,
660 FluidSystem::enthalpy(fs, paramCache, phaseIdx));
664 const auto& problem = elemCtx.problem();
665 const auto& heatCondParams = problem.heatConductionParams(elemCtx, dofIdx, timeIdx);
667 heatCapacitySolid_ = problem.heatCapacitySolid(elemCtx, dofIdx, timeIdx);
669 HeatConductionLaw::template heatConductivity<FluidState, Evaluation>(heatCondParams, fs);
671 Opm::Valgrind::CheckDefined(heatCapacitySolid_);
672 Opm::Valgrind::CheckDefined(heatConductivity_);
682 {
return heatCapacitySolid_; }
690 {
return heatConductivity_; }
693 Evaluation heatCapacitySolid_;
694 Evaluation heatConductivity_;
703 template <
class TypeTag,
bool enableEnergy>
709 template <
class TypeTag>
713 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
720 void update_(
const ElementContext& elemCtx OPM_UNUSED,
721 unsigned faceIdx OPM_UNUSED,
722 unsigned timeIdx OPM_UNUSED)
725 template <
class Context,
class Flu
idState>
726 void updateBoundary_(
const Context& context OPM_UNUSED,
727 unsigned bfIdx OPM_UNUSED,
728 unsigned timeIdx OPM_UNUSED,
729 const FluidState& fs OPM_UNUSED)
738 OPM_THROW(std::logic_error,
"Method temperatureGradNormal() does not " 739 "make sense for isothermal models");
748 OPM_THROW(std::logic_error,
"Method heatConductivity() does not make " 749 "sense for isothermal models");
756 template <
class TypeTag>
759 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
761 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
764 enum { dimWorld = GridView::dimensionworld };
765 typedef Dune::FieldVector<Evaluation, dimWorld> EvalDimVector;
766 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
773 void update_(
const ElementContext& elemCtx,
unsigned faceIdx,
unsigned timeIdx)
775 const auto& gradCalc = elemCtx.gradientCalculator();
778 EvalDimVector temperatureGrad;
779 gradCalc.calculateGradient(temperatureGrad,
782 temperatureCallback);
785 const auto& face = elemCtx.stencil(0).interiorFace(faceIdx);
787 temperatureGradNormal_ = 0;
788 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
789 temperatureGradNormal_ += (face.normal()[dimIdx]*temperatureGrad[dimIdx]);
791 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, timeIdx);
792 const auto& intQuantsInside = elemCtx.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
793 const auto& intQuantsOutside = elemCtx.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
797 0.5 * (intQuantsInside.heatConductivity() + intQuantsOutside.heatConductivity());
798 Opm::Valgrind::CheckDefined(heatConductivity_);
801 template <
class Context,
class Flu
idState>
802 void updateBoundary_(
const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fs)
804 const auto& stencil = context.stencil(timeIdx);
805 const auto& face = stencil.boundaryFace(bfIdx);
807 const auto& elemCtx = context.elementContext();
808 unsigned insideScvIdx = face.interiorIndex();
809 const auto& insideScv = elemCtx.stencil(timeIdx).subControlVolume(insideScvIdx);
811 const auto& intQuantsInside = elemCtx.intensiveQuantities(insideScvIdx, timeIdx);
812 const auto& fsInside = intQuantsInside.fluidState();
815 DimVector distVec = face.integrationPos();
816 distVec -= insideScv.geometry().center();
819 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
820 tmp += distVec[dimIdx] * face.normal()[dimIdx];
829 temperatureGradNormal_ =
830 (fs.temperature(0) - fsInside.temperature(0)) / dist;
833 heatConductivity_ = intQuantsInside.heatConductivity();
841 {
return temperatureGradNormal_; }
848 {
return heatConductivity_; }
851 Evaluation temperatureGradNormal_;
852 Evaluation heatConductivity_;
Scalar heatConductivity() const
The total heat conductivity at the face .
Definition: energymodule.hh:746
static void addAdvectiveFlux(RateVector &flux OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED)
Evaluates the advective energy fluxes over a face of a subcontrol volume and adds the result in the f...
Definition: energymodule.hh:196
static void handleFractureFlux(RateVector &flux OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED)
Evaluates the advective energy fluxes over a fracture which should be attributed to a face of a subco...
Definition: energymodule.hh:208
static void addSolidHeatStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Add the energy storage term for a fluid phase to an equation vector.
Definition: energymodule.hh:403
static Scalar heatConductionRate(const ExtensiveQuantities &extQuants OPM_UNUSED)
Add the rate of the conductive heat flux to a rate vector.
Definition: energymodule.hh:156
Provides the quantities required to calculate energy fluxes.
Definition: energymodule.hh:704
static void setPriVarTemperatures(PrimaryVariables &priVars, const FluidState &fs)
Given a fluid state, set the temperature in the primary variables.
Definition: energymodule.hh:350
static void setEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:314
static std::string eqName(unsigned eqIdx)
Returns the name of an equation or an empty string if the specified equation index does not belong to...
Definition: energymodule.hh:277
Definition: baseauxiliarymodule.hh:37
static std::string primaryVarName(unsigned pvIdx)
Returns the name of a primary variable or an empty string if the specified primary variable index doe...
Definition: energymodule.hh:265
static void addFracturePhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage OPM_UNUSED, const IntensiveQuantities &intQuants OPM_UNUSED, const Scv &scv OPM_UNUSED, unsigned phaseIdx OPM_UNUSED)
Add the energy storage term for a fluid phase to an equation vector.
Definition: energymodule.hh:174
static void addFracturePhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants, const Scv &scv, unsigned phaseIdx)
Add the energy storage term for a fluid phase to an equation vector.
Definition: energymodule.hh:384
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:59
void update_(FluidState &fs OPM_UNUSED, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache OPM_UNUSED, const ElementContext &elemCtx OPM_UNUSED, unsigned dofIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:596
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
static void addDiffusiveFlux(RateVector &flux OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED)
Adds the diffusive heat flux to the flux vector over the face of a sub-control volume.
Definition: energymodule.hh:221
static Scalar eqWeight(const Model &model OPM_UNUSED, unsigned globalDofIdx OPM_UNUSED, unsigned eqIdx OPM_UNUSED)
Returns the relative weight of a equation of the residual.
Definition: energymodule.hh:115
const Evaluation & heatConductivity() const
The total heat conductivity at the face .
Definition: energymodule.hh:847
static void registerParameters()
Register all run-time parameters for the energy module.
Definition: energymodule.hh:84
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Provides the indices required for the energy equation.
Definition: energymodule.hh:504
Evaluation heatCapacitySolid() const
Returns the total heat capacity of the rock matrix in the sub-control volume.
Definition: energymodule.hh:560
static Evaluation heatConductionRate(const ExtensiveQuantities &extQuants)
Returns the rate of the conductive heat flux for a given flux integration point.
Definition: energymodule.hh:327
static void registerParameters()
Register all run-time parameters for the energy module.
Definition: energymodule.hh:257
void update_(FluidState &fs, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:648
static void addToEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:320
Declare the properties used by the infrastructure code of the finite volume discretizations.
static void setEnthalpyRate(RateVector &rateVec, const FluidState &fluidState, unsigned phaseIdx, const Evaluation &volume)
Given a fluid state, set the enthalpy rate which emerges from a volumetric rate.
Definition: energymodule.hh:335
static void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage OPM_UNUSED, const IntensiveQuantities &intQuants OPM_UNUSED, unsigned phaseIdx OPM_UNUSED)
Add the energy storage term for a fluid phase to an equation vector.
Definition: energymodule.hh:164
static std::string primaryVarName(unsigned pvIdx OPM_UNUSED)
Returns the name of a primary variable or an empty string if the specified primary variable index doe...
Definition: energymodule.hh:92
Provides the volumetric quantities required for the energy equation.
Definition: energymodule.hh:539
Scalar temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition: energymodule.hh:736
This method contains all callback classes for quantities that are required by some extensive quantiti...
static Scalar primaryVarWeight(const Model &model, unsigned globalDofIdx, unsigned pvIdx)
Returns the relative weight of a primary variable for calculating relative errors.
Definition: energymodule.hh:288
static void setPriVarTemperatures(PrimaryVariables &priVars OPM_UNUSED, const FluidState &fs OPM_UNUSED)
Given a fluid state, set the temperature in the primary variables.
Definition: energymodule.hh:124
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:773
static Scalar primaryVarWeight(const Model &model OPM_UNUSED, unsigned globalDofIdx OPM_UNUSED, unsigned pvIdx OPM_UNUSED)
Returns the relative weight of a primary variable for calculating relative errors.
Definition: energymodule.hh:107
static void addSolidHeatStorage(Dune::FieldVector< LhsEval, numEq > &storage OPM_UNUSED, const IntensiveQuantities &intQuants OPM_UNUSED)
Add the energy storage term for the fracture part a fluid phase to an equation vector.
Definition: energymodule.hh:185
static void handleFractureFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Evaluates the advective energy fluxes over a fracture which should be attributed to a face of a subco...
Definition: energymodule.hh:447
static void updateTemperatures_(FluidState &fluidState, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Update the temperatures of the fluids of a fluid state.
Definition: energymodule.hh:628
static void addAdvectiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Evaluates the advective energy fluxes for a flux integration point and adds the result in the flux ve...
Definition: energymodule.hh:418
static void setEnthalpyRate(RateVector &rateVec OPM_UNUSED, const Evaluation &rate OPM_UNUSED)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:142
static void setEnthalpyRate(RateVector &rateVec OPM_UNUSED, const FluidState &fluidState OPM_UNUSED, unsigned phaseIdx OPM_UNUSED, const Evaluation &volume OPM_UNUSED)
Given a fluid state, set the enthalpy rate which emerges from a volumetric rate.
Definition: energymodule.hh:133
static Scalar eqWeight(const Model &model OPM_UNUSED, unsigned globalDofIdx OPM_UNUSED, unsigned eqIdx)
Returns the relative weight of a equation.
Definition: energymodule.hh:300
Callback class for temperature.
Definition: quantitycallbacks.hh:47
static std::string eqName(unsigned eqIdx OPM_UNUSED)
Returns the name of an equation or an empty string if the specified equation index does not belong to...
Definition: energymodule.hh:100
const Evaluation & temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition: energymodule.hh:840
static void addDiffusiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Adds the diffusive heat flux to the flux vector over the face of a sub-control volume.
Definition: energymodule.hh:483
const Evaluation & heatCapacitySolid() const
Returns the total heat capacity of the rock matrix in the sub-control volume.
Definition: energymodule.hh:681
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
static void updateTemperatures_(FluidState &fluidState, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Update the temperatures of the fluids of a fluid state.
Definition: energymodule.hh:582
const Evaluation & heatConductivity() const
Returns the total conductivity capacity of the rock matrix in the sub-control volume.
Definition: energymodule.hh:689
Evaluation heatConductivity() const
Returns the total conductivity capacity of the rock matrix in the sub-control volume.
Definition: energymodule.hh:571
static void addToEnthalpyRate(RateVector &rateVec OPM_UNUSED, const Evaluation &rate OPM_UNUSED)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:149
static void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants, unsigned phaseIdx)
Add the energy storage term for a fluid phase to an equation vector.
Definition: energymodule.hh:367
void update_(const ElementContext &elemCtx OPM_UNUSED, unsigned faceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:720