energymodule.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*/
29#ifndef EWOMS_ENERGY_MODULE_HH
30#define EWOMS_ENERGY_MODULE_HH
31
32#include <dune/common/fvector.hh>
33
34#include <opm/material/common/Valgrind.hpp>
35
38
40
41#include <cassert>
42#include <cmath>
43#include <stdexcept>
44#include <string>
45#include <type_traits>
46
47namespace Opm {
53template <class TypeTag, bool enableEnergy>
55
59template <class TypeTag>
60class EnergyModule<TypeTag, /*enableEnergy=*/false>
61{
70
71 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
72
73 using EvalEqVector = Dune::FieldVector<Evaluation, numEq>;
74
75public:
79 static void registerParameters()
80 {}
81
87 static std::string primaryVarName(unsigned)
88 { return ""; }
89
95 static std::string eqName(unsigned)
96 { return ""; }
97
102 static Scalar primaryVarWeight(const Model&,
103 unsigned,
104 unsigned)
105 { return -1; }
106
110 static Scalar eqWeight(const Model&,
111 unsigned,
112 unsigned)
113 { return -1; }
114
118 template <class FluidState>
119 static void setPriVarTemperatures(PrimaryVariables&,
120 const FluidState&)
121 {}
122
127 template <class RateVector, class FluidState>
128 static void setEnthalpyRate(RateVector&,
129 const FluidState&,
130 unsigned,
131 const Evaluation&)
132 {}
133
137 static void setEnthalpyRate(RateVector&,
138 const Evaluation&)
139 {}
140
144 static void addToEnthalpyRate(RateVector&,
145 const Evaluation&)
146 {}
147
151 static Scalar thermalConductionRate(const ExtensiveQuantities&)
152 { return 0.0; }
153
158 template <class LhsEval>
159 static void addPhaseStorage(Dune::FieldVector<LhsEval, numEq>&,
160 const IntensiveQuantities&,
161 unsigned)
162 {}
163
168 template <class LhsEval, class Scv>
169 static void addFracturePhaseStorage(Dune::FieldVector<LhsEval, numEq>&,
170 const IntensiveQuantities&,
171 const Scv&,
172 unsigned)
173 {}
174
179 template <class LhsEval>
180 static void addSolidEnergyStorage(Dune::FieldVector<LhsEval, numEq>&,
181 const IntensiveQuantities&)
182 {}
183
190 template <class Context>
191 static void addAdvectiveFlux(RateVector&,
192 const Context&,
193 unsigned,
194 unsigned)
195 {}
196
202 template <class Context>
203 static void handleFractureFlux(RateVector&,
204 const Context&,
205 unsigned,
206 unsigned)
207 {}
208
215 template <class Context>
216 static void addDiffusiveFlux(RateVector&,
217 const Context&,
218 unsigned,
219 unsigned)
220 {}
221};
222
226template <class TypeTag>
227class EnergyModule<TypeTag, /*enableEnergy=*/true>
228{
239
240 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
241 enum { numPhases = FluidSystem::numPhases };
242 enum { energyEqIdx = Indices::energyEqIdx };
243 enum { temperatureIdx = Indices::temperatureIdx };
244
245 using EvalEqVector = Dune::FieldVector<Evaluation, numEq>;
246 using Toolbox = Opm::MathToolbox<Evaluation>;
247
248public:
252 static void registerParameters()
253 {}
254
260 static std::string primaryVarName(unsigned pvIdx)
261 {
262 if (pvIdx == temperatureIdx) {
263 return "temperature";
264 }
265 return "";
266 }
267
273 static std::string eqName(unsigned eqIdx)
274 {
275 if (eqIdx == energyEqIdx) {
276 return "energy";
277 }
278 return "";
279 }
280
285 static Scalar primaryVarWeight(const Model& model, unsigned globalDofIdx, unsigned pvIdx)
286 {
287 if (pvIdx != temperatureIdx) {
288 return -1;
289 }
290
291 // make the weight of the temperature primary variable inversly proportional to its value
292 return std::max(1.0 / 1000,
293 1.0 / model.solution(/*timeIdx=*/0)[globalDofIdx][temperatureIdx]);
294 }
295
299 static Scalar eqWeight(const Model&,
300 unsigned,
301 unsigned eqIdx)
302 {
303 if (eqIdx != energyEqIdx) {
304 return -1;
305 }
306
307 // approximate change of internal energy of 1kg of liquid water for a temperature
308 // change of 30K
309 return 1.0 / (4.184e3 * 30.0);
310 }
311
315 static void setEnthalpyRate(RateVector& rateVec, const Evaluation& rate)
316 { rateVec[energyEqIdx] = rate; }
317
321 static void addToEnthalpyRate(RateVector& rateVec, const Evaluation& rate)
322 { rateVec[energyEqIdx] += rate; }
323
327 static Evaluation thermalConductionRate(const ExtensiveQuantities& extQuants)
328 { return -extQuants.temperatureGradNormal() * extQuants.thermalConductivity(); }
329
334 template <class RateVector, class FluidState>
335 static void setEnthalpyRate(RateVector& rateVec,
336 const FluidState& fluidState,
337 unsigned phaseIdx,
338 const Evaluation& volume)
339 {
340 rateVec[energyEqIdx] =
341 volume *
342 fluidState.density(phaseIdx) *
343 fluidState.enthalpy(phaseIdx);
344 }
345
349 template <class FluidState>
350 static void setPriVarTemperatures(PrimaryVariables& priVars,
351 const FluidState& fs)
352 {
353 priVars[temperatureIdx] = Toolbox::value(fs.temperature(/*phaseIdx=*/0));
354#ifndef NDEBUG
355 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
356 assert(std::abs(Toolbox::value(fs.temperature(/*phaseIdx=*/0))
357 - Toolbox::value(fs.temperature(phaseIdx))) < 1e-30);
358 }
359#endif
360 }
361
366 template <class LhsEval>
367 static void addPhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage,
368 const IntensiveQuantities& intQuants,
369 unsigned phaseIdx)
370 {
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());
377 }
378
383 template <class Scv, class LhsEval>
384 static void addFracturePhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage,
385 const IntensiveQuantities& intQuants,
386 const Scv& scv,
387 unsigned phaseIdx)
388 {
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();
396 }
397
402 template <class LhsEval>
403 static void addSolidEnergyStorage(Dune::FieldVector<LhsEval, numEq>& storage,
404 const IntensiveQuantities& intQuants)
405 { storage[energyEqIdx] += Opm::decay<LhsEval>(intQuants.solidInternalEnergy()); }
406
413 template <class Context>
414 static void addAdvectiveFlux(RateVector& flux,
415 const Context& context,
416 unsigned spaceIdx,
417 unsigned timeIdx)
418 {
419 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
420
421 // advective energy flux in all phases
422 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
423 if (!context.model().phaseIsConsidered(phaseIdx)) {
424 continue;
425 }
426
427 // intensive quantities of the upstream and the downstream DOFs
428 const unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
429 const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
430
431 flux[energyEqIdx] +=
432 extQuants.volumeFlux(phaseIdx) *
433 up.fluidState().enthalpy(phaseIdx) *
434 up.fluidState().density(phaseIdx);
435 }
436 }
437
443 template <class Context>
444 static void handleFractureFlux(RateVector& flux,
445 const Context& context,
446 unsigned spaceIdx,
447 unsigned timeIdx)
448 {
449 const auto& scvf = context.stencil(timeIdx).interiorFace(spaceIdx);
450 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
451
452 // reduce the energy flux in the matrix by the half the width occupied by the
453 // fracture
454 flux[energyEqIdx] *= 1 - extQuants.fractureWidth() / (2 * scvf.area());
455
456 // advective energy flux in all phases
457 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
458 if (!context.model().phaseIsConsidered(phaseIdx)) {
459 continue;
460 }
461
462 // intensive quantities of the upstream and the downstream DOFs
463 const unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
464 const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
465
466 flux[energyEqIdx] +=
467 extQuants.fractureVolumeFlux(phaseIdx) *
468 up.fluidState().enthalpy(phaseIdx) *
469 up.fluidState().density(phaseIdx);
470 }
471 }
472
479 template <class Context>
480 static void addDiffusiveFlux(RateVector& flux,
481 const Context& context,
482 unsigned spaceIdx,
483 unsigned timeIdx)
484 {
485 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
486
487 // diffusive energy flux
488 flux[energyEqIdx] += -extQuants.temperatureGradNormal() * extQuants.thermalConductivity();
489 }
490};
491
498template <unsigned PVOffset, bool enableEnergy>
500
504template <unsigned PVOffset>
505struct EnergyIndices<PVOffset, /*enableEnergy=*/false>
506{
508 enum { temperatureIdx = -1000 };
509
511 enum { energyEqIdx = -1000 };
512
513protected:
514 enum { numEq_ = 0 };
515};
516
520template <unsigned PVOffset>
521struct EnergyIndices<PVOffset, /*enableEnergy=*/true>
522{
524 enum { temperatureIdx = PVOffset };
525
527 enum { energyEqIdx = PVOffset };
528
529protected:
530 enum { numEq_ = 1 };
531};
532
539template <class TypeTag, bool enableEnergy>
541
545template <class TypeTag>
546class EnergyIntensiveQuantities<TypeTag, /*enableEnergy=*/false>
547{
552
553 using Toolbox = Opm::MathToolbox<Evaluation>;
554
555public:
560 Evaluation solidInternalEnergy() const
561 {
562 throw std::logic_error("solidInternalEnergy() does not make sense for isothermal models");
563 }
564
569 Evaluation thermalConductivity() const
570 {
571 throw std::logic_error("thermalConductivity() does not make sense for isothermal models");
572 }
573
574protected:
578 template <class FluidState, class Context>
579 static void updateTemperatures_(FluidState& fluidState,
580 const Context& context,
581 unsigned spaceIdx,
582 unsigned timeIdx)
583 {
584 const Scalar T = context.problem().temperature(context, spaceIdx, timeIdx);
585 fluidState.setTemperature(Toolbox::createConstant(T));
586 }
587
592 template <class FluidState>
593 void update_(FluidState&,
594 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>&,
595 const ElementContext&,
596 unsigned,
597 unsigned)
598 {}
599};
600
604template <class TypeTag>
605class EnergyIntensiveQuantities<TypeTag, /*enableEnergy=*/true>
606{
614
615 enum { numPhases = FluidSystem::numPhases };
616 enum { energyEqIdx = Indices::energyEqIdx };
617 enum { temperatureIdx = Indices::temperatureIdx };
618
619 using Toolbox = Opm::MathToolbox<Evaluation>;
620
621protected:
625 template <class FluidState, class Context>
626 static void updateTemperatures_(FluidState& fluidState,
627 const Context& context,
628 unsigned spaceIdx,
629 unsigned timeIdx)
630 {
631 const auto& priVars = context.primaryVars(spaceIdx, timeIdx);
632 if constexpr (std::is_same_v<Evaluation, Scalar>) { // finite differences
633 fluidState.setTemperature(Toolbox::createConstant(priVars[temperatureIdx]));
634 }
635 else {
636 // automatic differentiation
637 if (timeIdx == 0) {
638 fluidState.setTemperature(Toolbox::createVariable(priVars[temperatureIdx], temperatureIdx));
639 }
640 else {
641 fluidState.setTemperature(Toolbox::createConstant(priVars[temperatureIdx]));
642 }
643 }
644 }
645
650 template <class FluidState>
651 void update_(FluidState& fs,
652 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
653 const ElementContext& elemCtx,
654 unsigned dofIdx,
655 unsigned timeIdx)
656 {
657 // set the specific enthalpies of the fluids
658 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
659 if (!elemCtx.model().phaseIsConsidered(phaseIdx)) {
660 continue;
661 }
662
663 fs.setEnthalpy(phaseIdx,
664 FluidSystem::enthalpy(fs, paramCache, phaseIdx));
665 }
666
667 // compute and set the volumetric internal energy of the solid phase
668 const auto& problem = elemCtx.problem();
669 const auto& solidEnergyParams = problem.solidEnergyLawParams(elemCtx, dofIdx, timeIdx);
670 const auto& thermalCondParams = problem.thermalConductionLawParams(elemCtx, dofIdx, timeIdx);
671
672 solidInternalEnergy_ = SolidEnergyLaw::solidInternalEnergy(solidEnergyParams, fs);
673 thermalConductivity_ = ThermalConductionLaw::thermalConductivity(thermalCondParams, fs);
674
675 Opm::Valgrind::CheckDefined(solidInternalEnergy_);
676 Opm::Valgrind::CheckDefined(thermalConductivity_);
677 }
678
679public:
684 const Evaluation& solidInternalEnergy() const
685 { return solidInternalEnergy_; }
686
691 const Evaluation& thermalConductivity() const
692 { return thermalConductivity_; }
693
694private:
695 Evaluation solidInternalEnergy_;
696 Evaluation thermalConductivity_;
697};
698
705template <class TypeTag, bool enableEnergy>
707
711template <class TypeTag>
712class EnergyExtensiveQuantities<TypeTag, /*enableEnergy=*/false>
713{
716
717protected:
722 void update_(const ElementContext&,
723 unsigned,
724 unsigned)
725 {}
726
727 template <class Context, class FluidState>
728 void updateBoundary_(const Context&,
729 unsigned,
730 unsigned,
731 const FluidState&)
732 {}
733
734public:
739 {
740 throw std::logic_error("Calling temperatureGradNormal() does not make sense "
741 "for isothermal models");
742 }
743
747 Scalar thermalConductivity() const
748 {
749 throw std::logic_error("Calling thermalConductivity() does not make sense for "
750 "isothermal models");
751 }
752};
753
757template <class TypeTag>
758class EnergyExtensiveQuantities<TypeTag, /*enableEnergy=*/true>
759{
764
765 enum { dimWorld = GridView::dimensionworld };
766 using EvalDimVector = Dune::FieldVector<Evaluation, dimWorld>;
767 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
768
769protected:
774 void update_(const ElementContext& elemCtx, unsigned faceIdx, unsigned timeIdx)
775 {
776 const auto& gradCalc = elemCtx.gradientCalculator();
777 TemperatureCallback<TypeTag> temperatureCallback(elemCtx);
778
779 EvalDimVector temperatureGrad;
780 gradCalc.calculateGradient(temperatureGrad,
781 elemCtx,
782 faceIdx,
783 temperatureCallback);
784
785 // scalar product of temperature gradient and scvf normal
786 const auto& face = elemCtx.stencil(/*timeIdx=*/0).interiorFace(faceIdx);
787
788 temperatureGradNormal_ = 0;
789 for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
790 temperatureGradNormal_ += (face.normal()[dimIdx]*temperatureGrad[dimIdx]);
791 }
792
793 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, timeIdx);
794 const auto& intQuantsInside = elemCtx.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
795 const auto& intQuantsOutside = elemCtx.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
796
797 // arithmetic mean
798 thermalConductivity_ = 0.5 * (intQuantsInside.thermalConductivity() +
799 intQuantsOutside.thermalConductivity());
800 Opm::Valgrind::CheckDefined(thermalConductivity_);
801 }
802
803 template <class Context, class FluidState>
804 void updateBoundary_(const Context& context, unsigned bfIdx, unsigned timeIdx, const FluidState& fs)
805 {
806 const auto& stencil = context.stencil(timeIdx);
807 const auto& face = stencil.boundaryFace(bfIdx);
808
809 const auto& elemCtx = context.elementContext();
810 const unsigned insideScvIdx = face.interiorIndex();
811 const auto& insideScv = elemCtx.stencil(timeIdx).subControlVolume(insideScvIdx);
812
813 const auto& intQuantsInside = elemCtx.intensiveQuantities(insideScvIdx, timeIdx);
814 const auto& fsInside = intQuantsInside.fluidState();
815
816 // distance between the center of the SCV and center of the boundary face
817 DimVector distVec = face.integrationPos();
818 distVec -= insideScv.geometry().center();
819
820 Scalar dist = 0;
821 for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
822 dist += distVec[dimIdx] * face.normal()[dimIdx];
823 }
824
825 // if the following assertation triggers, the center of the
826 // center of the interior SCV was not inside the element!
827 assert(dist > 0);
828
829 // calculate the temperature gradient using two-point gradient
830 // appoximation
831 temperatureGradNormal_ =
832 (fs.temperature(/*phaseIdx=*/0) - fsInside.temperature(/*phaseIdx=*/0)) / dist;
833
834 // take the value for thermal conductivity from the interior finite volume
835 thermalConductivity_ = intQuantsInside.thermalConductivity();
836 }
837
838public:
842 const Evaluation& temperatureGradNormal() const
843 { return temperatureGradNormal_; }
844
849 const Evaluation& thermalConductivity() const
850 { return thermalConductivity_; }
851
852private:
853 Evaluation temperatureGradNormal_;
854 Evaluation thermalConductivity_;
855};
856
857} // namespace Opm
858
859#endif
void updateBoundary_(const Context &, unsigned, unsigned, const FluidState &)
Definition: energymodule.hh:728
Scalar temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition: energymodule.hh:738
Scalar thermalConductivity() const
The total thermal conductivity at the face .
Definition: energymodule.hh:747
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:722
const Evaluation & temperatureGradNormal() const
The temperature gradient times the face normal [K m^2 / m].
Definition: energymodule.hh:842
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:774
const Evaluation & thermalConductivity() const
The total thermal conductivity at the face .
Definition: energymodule.hh:849
void updateBoundary_(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fs)
Definition: energymodule.hh:804
Provides the quantities required to calculate energy fluxes.
Definition: energymodule.hh:706
Evaluation thermalConductivity() const
Returns the total thermal conductivity of the solid matrix in the sub-control volume.
Definition: energymodule.hh:569
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:579
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:593
Evaluation solidInternalEnergy() const
Returns the volumetric internal energy of the solid matrix in the sub-control volume.
Definition: energymodule.hh:560
const Evaluation & solidInternalEnergy() const
Returns the volumetric internal energy of the solid matrix in the sub-control volume.
Definition: energymodule.hh:684
const Evaluation & thermalConductivity() const
Returns the total conductivity capacity of the solid matrix in the sub-control volume.
Definition: energymodule.hh:691
void update_(FluidState &fs, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > &paramCache, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate energy fluxes.
Definition: energymodule.hh:651
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:626
Provides the volumetric quantities required for the energy equation.
Definition: energymodule.hh:540
static Scalar eqWeight(const Model &, unsigned, unsigned)
Returns the relative weight of a equation of the residual.
Definition: energymodule.hh:110
static void setPriVarTemperatures(PrimaryVariables &, const FluidState &)
Given a fluid state, set the temperature in the primary variables.
Definition: energymodule.hh:119
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:203
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:159
static void addToEnthalpyRate(RateVector &, const Evaluation &)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:144
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:216
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:180
static void setEnthalpyRate(RateVector &, const Evaluation &)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:137
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:191
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:169
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:87
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:95
static Scalar thermalConductionRate(const ExtensiveQuantities &)
Add the rate of the conductive energy flux to a rate vector.
Definition: energymodule.hh:151
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:128
static Scalar primaryVarWeight(const Model &, unsigned, unsigned)
Returns the relative weight of a primary variable for calculating relative errors.
Definition: energymodule.hh:102
static void registerParameters()
Register all run-time parameters for the energy module.
Definition: energymodule.hh:79
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:444
static void addToEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Add the rate of the enthalpy flux to a rate vector.
Definition: energymodule.hh:321
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
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
static void setEnthalpyRate(RateVector &rateVec, const Evaluation &rate)
Set the rate of energy flux of a rate vector.
Definition: energymodule.hh:315
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 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:480
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 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:273
static void registerParameters()
Register all run-time parameters for the energy module.
Definition: energymodule.hh:252
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:403
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:414
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:260
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:285
static Scalar eqWeight(const Model &, unsigned, unsigned eqIdx)
Returns the relative weight of a equation.
Definition: energymodule.hh:299
static Evaluation thermalConductionRate(const ExtensiveQuantities &extQuants)
Returns the conductive energy flux for a given flux integration point.
Definition: energymodule.hh:327
Provides the auxiliary methods required for consideration of the energy equation.
Definition: energymodule.hh:54
Callback class for temperature.
Definition: quantitycallbacks.hh:49
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.
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
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:499