29#ifndef EWOMS_ENERGY_MODULE_HH
30#define EWOMS_ENERGY_MODULE_HH
32#include <dune/common/fvector.hh>
34#include <opm/material/common/Valgrind.hpp>
49template <
class TypeTag,
bool enableEnergy>
55template <
class TypeTag>
69 using EvalEqVector = Dune::FieldVector<Evaluation, numEq>;
114 template <
class Flu
idState>
123 template <
class RateVector,
class Flu
idState>
154 template <
class LhsEval>
156 const IntensiveQuantities&,
164 template <
class LhsEval,
class Scv>
166 const IntensiveQuantities&,
175 template <
class LhsEval>
177 const IntensiveQuantities&)
186 template <
class Context>
198 template <
class Context>
211 template <
class Context>
222template <
class TypeTag>
237 enum { numPhases = FluidSystem::numPhases };
238 enum { energyEqIdx = Indices::energyEqIdx };
239 enum { temperatureIdx = Indices::temperatureIdx };
241 using EvalEqVector = Dune::FieldVector<Evaluation, numEq>;
258 if (
pvIdx == temperatureIdx)
259 return "temperature";
270 if (
eqIdx == energyEqIdx)
281 if (
pvIdx != temperatureIdx)
285 return std::max(1.0/1000, 1.0/model.solution(0)[
globalDofIdx][temperatureIdx]);
295 if (
eqIdx != energyEqIdx)
300 return 1.0 / (4.184e3 * 30.0);
325 template <
class RateVector,
class Flu
idState>
327 const FluidState& fluidState,
329 const Evaluation& volume)
340 template <
class Flu
idState>
342 const FluidState&
fs)
344 priVars[temperatureIdx] = Toolbox::value(
fs.temperature(0));
347 assert(std::abs(Toolbox::value(
fs.temperature(0))
348 - Toolbox::value(
fs.temperature(
phaseIdx))) < 1
e-30);
357 template <
class LhsEval>
374 template <
class Scv,
class LhsEval>
393 template <
class LhsEval>
404 template <
class Context>
414 if (!context.model().phaseIsConsidered(
phaseIdx))
419 const IntensiveQuantities&
up = context.intensiveQuantities(
upIdx,
timeIdx);
433 template <
class Context>
449 if (!context.model().phaseIsConsidered(
phaseIdx))
454 const IntensiveQuantities&
up = context.intensiveQuantities(
upIdx,
timeIdx);
469 template <
class Context>
490template <
unsigned PVOffset,
bool enableEnergy>
496template <
unsigned PVOffset>
500 enum { temperatureIdx = -1000 };
503 enum { energyEqIdx = -1000 };
512template <
unsigned PVOffset>
531template <
class TypeTag,
bool enableEnergy>
537template <
class TypeTag>
554 throw std::logic_error(
"solidInternalEnergy() does not make sense for isothermal models");
563 throw std::logic_error(
"thermalConductivity() does not make sense for isothermal models");
570 template <
class Flu
idState,
class Context>
577 fluidState.setTemperature(Toolbox::createConstant(T));
584 template <
class Flu
idState>
587 const ElementContext&,
596template <
class TypeTag>
607 enum { numPhases = FluidSystem::numPhases };
608 enum { energyEqIdx = Indices::energyEqIdx };
609 enum { temperatureIdx = Indices::temperatureIdx };
617 template <
class Flu
idState,
class Context>
625 if (std::is_same<Evaluation, Scalar>::value)
626 val = Toolbox::createConstant(priVars[temperatureIdx]);
630 val = Toolbox::createVariable(priVars[temperatureIdx], temperatureIdx);
632 val = Toolbox::createConstant(priVars[temperatureIdx]);
634 fluidState.setTemperature(
val);
641 template <
class Flu
idState>
658 const auto& problem =
elemCtx.problem();
662 solidInternalEnergy_ = SolidEnergyLaw::solidInternalEnergy(solidEnergyParams,
fs);
665 Opm::Valgrind::CheckDefined(solidInternalEnergy_);
666 Opm::Valgrind::CheckDefined(thermalConductivity_);
675 {
return solidInternalEnergy_; }
682 {
return thermalConductivity_; }
685 Evaluation solidInternalEnergy_;
686 Evaluation thermalConductivity_;
695template <
class TypeTag,
bool enableEnergy>
701template <
class TypeTag>
717 template <
class Context,
class Flu
idState>
718 void updateBoundary_(
const Context&,
730 throw std::logic_error(
"Calling temperatureGradNormal() does not make sense "
731 "for isothermal models");
739 throw std::logic_error(
"Calling thermalConductivity() does not make sense for "
740 "isothermal models");
747template <
class TypeTag>
755 enum { dimWorld = GridView::dimensionworld };
756 using EvalDimVector = Dune::FieldVector<Evaluation, dimWorld>;
757 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
778 temperatureGradNormal_ = 0;
787 thermalConductivity_ =
789 Opm::Valgrind::CheckDefined(thermalConductivity_);
792 template <
class Context,
class Flu
idState>
795 const auto& stencil = context.stencil(
timeIdx);
796 const auto&
face = stencil.boundaryFace(
bfIdx);
798 const auto&
elemCtx = context.elementContext();
820 temperatureGradNormal_ =
832 {
return temperatureGradNormal_; }
839 {
return thermalConductivity_; }
842 Evaluation temperatureGradNormal_;
843 Evaluation thermalConductivity_;
Scalar temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition energymodule.hh:728
Scalar thermalConductivity() const
The total thermal conductivity at the face .
Definition energymodule.hh:737
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate energy fluxes.
Definition energymodule.hh:712
const Evaluation & temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition energymodule.hh:831
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate energy fluxes.
Definition energymodule.hh:764
const Evaluation & thermalConductivity() const
The total thermal conductivity at the face .
Definition energymodule.hh:838
Provides the quantities required to calculate energy fluxes.
Definition energymodule.hh:696
Evaluation thermalConductivity() const
Returns the total thermal conductivity of the solid matrix in the sub-control volume.
Definition energymodule.hh:561
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:571
void update_(FluidState &, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > &, const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate energy fluxes.
Definition energymodule.hh:585
Evaluation solidInternalEnergy() const
Returns the volumetric internal energy of the solid matrix in the sub-control volume.
Definition energymodule.hh:552
const Evaluation & solidInternalEnergy() const
Returns the volumetric internal energy of the solid matrix in the sub-control volume.
Definition energymodule.hh:674
const Evaluation & thermalConductivity() const
Returns the total conductivity capacity of the solid matrix in the sub-control volume.
Definition energymodule.hh:681
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:642
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:618
Provides the volumetric quantities required for the energy equation.
Definition energymodule.hh:532
static Scalar eqWeight(const Model &, unsigned, unsigned)
Returns the relative weight of a equation of the residual.
Definition energymodule.hh:106
static void setPriVarTemperatures(PrimaryVariables &, const FluidState &)
Given a fluid state, set the temperature in the primary variables.
Definition energymodule.hh:115
static void handleFractureFlux(RateVector &, const Context &, unsigned, unsigned)
Evaluates the advective energy fluxes over a fracture which should be attributed to a face of a subco...
Definition energymodule.hh:199
static void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &, const IntensiveQuantities &, unsigned)
Add the energy storage term for a fluid phase to an equation vector.
Definition energymodule.hh:155
static void addToEnthalpyRate(RateVector &, const Evaluation &)
Add the rate of the enthalpy flux to a rate vector.
Definition energymodule.hh:140
static void addDiffusiveFlux(RateVector &, const Context &, unsigned, unsigned)
Adds the diffusive energy flux to the flux vector over the face of a sub-control volume.
Definition energymodule.hh:212
static void addSolidEnergyStorage(Dune::FieldVector< LhsEval, numEq > &, const IntensiveQuantities &)
Add the energy storage term for the fracture part a fluid phase to an equation vector.
Definition energymodule.hh:176
static void setEnthalpyRate(RateVector &, const Evaluation &)
Add the rate of the enthalpy flux to a rate vector.
Definition energymodule.hh:133
static void addAdvectiveFlux(RateVector &, const Context &, unsigned, unsigned)
Evaluates the advective energy fluxes over a face of a subcontrol volume and adds the result in the f...
Definition energymodule.hh:187
static void addFracturePhaseStorage(Dune::FieldVector< LhsEval, numEq > &, const IntensiveQuantities &, const Scv &, unsigned)
Add the energy storage term for a fluid phase to an equation vector.
Definition energymodule.hh:165
static std::string primaryVarName(unsigned)
Returns the name of a primary variable or an empty string if the specified primary variable index doe...
Definition energymodule.hh:83
static std::string eqName(unsigned)
Returns the name of an equation or an empty string if the specified equation index does not belong to...
Definition energymodule.hh:91
static Scalar thermalConductionRate(const ExtensiveQuantities &)
Add the rate of the conductive energy flux to a rate vector.
Definition energymodule.hh:147
static void setEnthalpyRate(RateVector &, const FluidState &, unsigned, const Evaluation &)
Given a fluid state, set the enthalpy rate which emerges from a volumetric rate.
Definition energymodule.hh:124
static Scalar primaryVarWeight(const Model &, unsigned, unsigned)
Returns the relative weight of a primary variable for calculating relative errors.
Definition energymodule.hh:98
static void registerParameters()
Register all run-time parameters for the energy module.
Definition energymodule.hh:75
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:434
static void addToEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Add the rate of the enthalpy flux to a rate vector.
Definition energymodule.hh:312
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:375
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:358
static void setEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Set the rate of energy flux of a rate vector.
Definition energymodule.hh:306
static void setPriVarTemperatures(PrimaryVariables &priVars, const FluidState &fs)
Given a fluid state, set the temperature in the primary variables.
Definition energymodule.hh:341
static void addDiffusiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Adds the diffusive energy flux to the flux vector over the face of a sub-control volume.
Definition energymodule.hh:470
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:326
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:268
static void registerParameters()
Register all run-time parameters for the energy module.
Definition energymodule.hh:248
static void addSolidEnergyStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Add the energy storage term for a fluid phase to an equation vector.
Definition energymodule.hh:394
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:405
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:256
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:279
static Scalar eqWeight(const Model &, unsigned, unsigned eqIdx)
Returns the relative weight of a equation.
Definition energymodule.hh:291
static Evaluation thermalConductionRate(const ExtensiveQuantities &extQuants)
Returns the conductive energy flux for a given flux integration point.
Definition energymodule.hh:318
Provides the auxiliary methods required for consideration of the energy equation.
Definition energymodule.hh:50
Callback class for temperature.
Definition quantitycallbacks.hh:48
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
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...
Provides the indices required for the energy equation.
Definition energymodule.hh:491