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 
34 
35 #include <opm/common/Valgrind.hpp>
36 #include <opm/common/Unused.hpp>
37 #include <opm/common/ErrorMacros.hpp>
38 #include <opm/common/Exceptions.hpp>
39 
40 #include <dune/common/fvector.hh>
41 
42 #include <string>
43 
44 namespace Ewoms {
45 namespace Properties {
46 NEW_PROP_TAG(Indices);
47 NEW_PROP_TAG(EnableEnergy);
48 NEW_PROP_TAG(HeatConductionLaw);
49 NEW_PROP_TAG(HeatConductionLawParams);
50 }}
51 
52 namespace Ewoms {
58 template <class TypeTag, bool enableEnergy>
60 
64 template <class TypeTag>
65 class EnergyModule<TypeTag, /*enableEnergy=*/false>
66 {
67  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
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;
74  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
75 
76  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
77 
78  typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
79 
80 public:
84  static void registerParameters()
85  {}
86 
92  static std::string primaryVarName(unsigned pvIdx OPM_UNUSED)
93  { return ""; }
94 
100  static std::string eqName(unsigned eqIdx OPM_UNUSED)
101  { return ""; }
102 
107  static Scalar primaryVarWeight(const Model& model OPM_UNUSED,
108  unsigned globalDofIdx OPM_UNUSED,
109  unsigned pvIdx OPM_UNUSED)
110  { return -1; }
111 
115  static Scalar eqWeight(const Model& model OPM_UNUSED,
116  unsigned globalDofIdx OPM_UNUSED,
117  unsigned eqIdx OPM_UNUSED)
118  { return -1; }
119 
123  template <class FluidState>
124  static void setPriVarTemperatures(PrimaryVariables& priVars OPM_UNUSED,
125  const FluidState& fs OPM_UNUSED)
126  {}
127 
132  template <class RateVector, class FluidState>
133  static void setEnthalpyRate(RateVector& rateVec OPM_UNUSED,
134  const FluidState& fluidState OPM_UNUSED,
135  unsigned phaseIdx OPM_UNUSED,
136  const Evaluation& volume OPM_UNUSED)
137  {}
138 
142  static void setEnthalpyRate(RateVector& rateVec OPM_UNUSED,
143  const Evaluation& rate OPM_UNUSED)
144  {}
145 
149  static void addToEnthalpyRate(RateVector& rateVec OPM_UNUSED,
150  const Evaluation& rate OPM_UNUSED)
151  {}
152 
156  static Scalar heatConductionRate(const ExtensiveQuantities& extQuants OPM_UNUSED)
157  { return 0.0; }
158 
163  template <class LhsEval>
164  static void addPhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage OPM_UNUSED,
165  const IntensiveQuantities& intQuants OPM_UNUSED,
166  unsigned phaseIdx OPM_UNUSED)
167  {}
168 
173  template <class LhsEval, class Scv>
174  static void addFracturePhaseStorage(Dune::FieldVector<LhsEval, numEq>& storage OPM_UNUSED,
175  const IntensiveQuantities& intQuants OPM_UNUSED,
176  const Scv& scv OPM_UNUSED,
177  unsigned phaseIdx OPM_UNUSED)
178  {}
179 
184  template <class LhsEval>
185  static void addSolidHeatStorage(Dune::FieldVector<LhsEval, numEq>& storage OPM_UNUSED,
186  const IntensiveQuantities& intQuants OPM_UNUSED)
187  {}
188 
195  template <class Context>
196  static void addAdvectiveFlux(RateVector& flux OPM_UNUSED,
197  const Context& context OPM_UNUSED,
198  unsigned spaceIdx OPM_UNUSED,
199  unsigned timeIdx OPM_UNUSED)
200  {}
201 
207  template <class Context>
208  static void handleFractureFlux(RateVector& flux OPM_UNUSED,
209  const Context& context OPM_UNUSED,
210  unsigned spaceIdx OPM_UNUSED,
211  unsigned timeIdx OPM_UNUSED)
212  {}
213 
220  template <class Context>
221  static void addDiffusiveFlux(RateVector& flux OPM_UNUSED,
222  const Context& context OPM_UNUSED,
223  unsigned spaceIdx OPM_UNUSED,
224  unsigned timeIdx OPM_UNUSED)
225  {}
226 };
227 
231 template <class TypeTag>
232 class EnergyModule<TypeTag, /*enableEnergy=*/true>
233 {
234  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
235  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
236  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
237  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
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;
242  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
243  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
244 
245  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
246  enum { numPhases = FluidSystem::numPhases };
247  enum { energyEqIdx = Indices::energyEqIdx };
248  enum { temperatureIdx = Indices::temperatureIdx };
249 
250  typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
251  typedef Opm::MathToolbox<Evaluation> Toolbox;
252 
253 public:
257  static void registerParameters()
258  {}
259 
265  static std::string primaryVarName(unsigned pvIdx)
266  {
267  if (pvIdx == temperatureIdx)
268  return "temperature";
269  return "";
270  }
271 
277  static std::string eqName(unsigned eqIdx)
278  {
279  if (eqIdx == energyEqIdx)
280  return "energy";
281  return "";
282  }
283 
288  static Scalar primaryVarWeight(const Model& model, unsigned globalDofIdx, unsigned pvIdx)
289  {
290  if (pvIdx != temperatureIdx)
291  return -1;
292 
293  // make the weight of the temperature primary variable inversly proportional to its value
294  return std::max(1.0/1000, 1.0/model.solution(/*timeIdx=*/0)[globalDofIdx][temperatureIdx]);
295  }
296 
300  static Scalar eqWeight(const Model& model OPM_UNUSED,
301  unsigned globalDofIdx OPM_UNUSED,
302  unsigned eqIdx)
303  {
304  if (eqIdx != energyEqIdx)
305  return -1;
306 
307  // approximate heat capacity of 1kg of air
308  return 1.0 / 1.0035e3;
309  }
310 
314  static void setEnthalpyRate(RateVector& rateVec, const Evaluation& rate)
315  { rateVec[energyEqIdx] = rate; }
316 
320  static void addToEnthalpyRate(RateVector& rateVec, const Evaluation& rate)
321  { rateVec[energyEqIdx] += rate; }
322 
327  static Evaluation heatConductionRate(const ExtensiveQuantities& extQuants)
328  { return -extQuants.temperatureGradNormal() * extQuants.heatConductivity(); }
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 addSolidHeatStorage(Dune::FieldVector<LhsEval, numEq>& storage,
404  const IntensiveQuantities& intQuants)
405  {
406  storage[energyEqIdx] +=
407  Toolbox::template decay<LhsEval>(intQuants.heatCapacitySolid())
408  * Toolbox::template decay<LhsEval>(intQuants.fluidState().temperature(/*phaseIdx=*/0));
409  }
410 
417  template <class Context>
418  static void addAdvectiveFlux(RateVector& flux,
419  const Context& context,
420  unsigned spaceIdx,
421  unsigned timeIdx)
422  {
423  const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
424 
425  // advective heat flux in all phases
426  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
427  if (!context.model().phaseIsConsidered(phaseIdx))
428  continue;
429 
430  // intensive quantities of the upstream and the downstream DOFs
431  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
432  const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
433 
434  flux[energyEqIdx] +=
435  extQuants.volumeFlux(phaseIdx)
436  * up.fluidState().enthalpy(phaseIdx)
437  * up.fluidState().density(phaseIdx);
438  }
439  }
440 
446  template <class Context>
447  static void handleFractureFlux(RateVector& flux,
448  const Context& context,
449  unsigned spaceIdx,
450  unsigned timeIdx)
451  {
452  const auto& scvf = context.stencil(timeIdx).interiorFace(spaceIdx);
453  const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
454 
455  // reduce the heat flux in the matrix by the half the width
456  // occupied by the fracture
457  flux[energyEqIdx] *=
458  1 - extQuants.fractureWidth()/(2*scvf.area());
459 
460  // advective heat flux in all phases
461  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
462  if (!context.model().phaseIsConsidered(phaseIdx))
463  continue;
464 
465  // intensive quantities of the upstream and the downstream DOFs
466  unsigned upIdx = static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
467  const IntensiveQuantities& up = context.intensiveQuantities(upIdx, timeIdx);
468 
469  flux[energyEqIdx] +=
470  extQuants.fractureVolumeFlux(phaseIdx)
471  * up.fluidState().enthalpy(phaseIdx)
472  * up.fluidState().density(phaseIdx);
473  }
474  }
475 
482  template <class Context>
483  static void addDiffusiveFlux(RateVector& flux,
484  const Context& context,
485  unsigned spaceIdx,
486  unsigned timeIdx)
487  {
488  const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
489 
490  // diffusive heat flux
491  flux[energyEqIdx] +=
492  - extQuants.temperatureGradNormal()
493  * extQuants.heatConductivity();
494  }
495 };
496 
503 template <unsigned PVOffset, bool enableEnergy>
505 
509 template <unsigned PVOffset>
510 struct EnergyIndices<PVOffset, /*enableEnergy=*/false>
511 {
512 protected:
513  enum { numEq_ = 0 };
514 };
515 
519 template <unsigned PVOffset>
520 struct EnergyIndices<PVOffset, /*enableEnergy=*/true>
521 {
523  enum { temperatureIdx = PVOffset };
524 
526  enum { energyEqIdx = PVOffset };
527 
528 protected:
529  enum { numEq_ = 1 };
530 };
531 
538 template <class TypeTag, bool enableEnergy>
540 
544 template <class TypeTag>
545 class EnergyIntensiveQuantities<TypeTag, /*enableEnergy=*/false>
546 {
547  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
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;
551 
552  typedef Opm::MathToolbox<Evaluation> Toolbox;
553 
554 public:
560  Evaluation heatCapacitySolid() const
561  {
562  OPM_THROW(std::logic_error, "Method heatCapacitySolid() does not make "
563  "sense for isothermal models");
564  }
565 
571  Evaluation heatConductivity() const
572  {
573  OPM_THROW(std::logic_error, "Method heatConductivity() does not make "
574  "sense for isothermal models");
575  }
576 
577 protected:
581  template <class FluidState, class Context>
582  static void updateTemperatures_(FluidState& fluidState,
583  const Context& context,
584  unsigned spaceIdx,
585  unsigned timeIdx)
586  {
587  Scalar T = context.problem().temperature(context, spaceIdx, timeIdx);
588  fluidState.setTemperature(Toolbox::createConstant(T));
589  }
590 
595  template <class FluidState>
596  void update_(FluidState& fs OPM_UNUSED,
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)
601  { }
602 };
603 
607 template <class TypeTag>
608 class EnergyIntensiveQuantities<TypeTag, /*enableEnergy=*/true>
609 {
610  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
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;
615  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
616 
617  enum { numPhases = FluidSystem::numPhases };
618  enum { energyEqIdx = Indices::energyEqIdx };
619  enum { temperatureIdx = Indices::temperatureIdx };
620 
621  typedef Opm::MathToolbox<Evaluation> Toolbox;
622 
623 protected:
627  template <class FluidState, class Context>
628  static void updateTemperatures_(FluidState& fluidState,
629  const Context& context,
630  unsigned spaceIdx,
631  unsigned timeIdx)
632  {
633  const auto& priVars = context.primaryVars(spaceIdx, timeIdx);
634  Evaluation val;
635  if (timeIdx == 0)
636  val = Toolbox::createVariable(priVars[temperatureIdx], temperatureIdx);
637  else
638  val = Toolbox::createConstant(priVars[temperatureIdx]);
639 
640  fluidState.setTemperature(val);
641  }
642 
647  template <class FluidState>
648  void update_(FluidState& fs,
649  typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
650  const ElementContext& elemCtx,
651  unsigned dofIdx,
652  unsigned timeIdx)
653  {
654  // set the specific enthalpies of the fluids
655  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
656  if (!elemCtx.model().phaseIsConsidered(phaseIdx))
657  continue;
658 
659  fs.setEnthalpy(phaseIdx,
660  FluidSystem::enthalpy(fs, paramCache, phaseIdx));
661  }
662 
663  // compute and set the heat capacity of the solid phase
664  const auto& problem = elemCtx.problem();
665  const auto& heatCondParams = problem.heatConductionParams(elemCtx, dofIdx, timeIdx);
666 
667  heatCapacitySolid_ = problem.heatCapacitySolid(elemCtx, dofIdx, timeIdx);
668  heatConductivity_ =
669  HeatConductionLaw::template heatConductivity<FluidState, Evaluation>(heatCondParams, fs);
670 
671  Opm::Valgrind::CheckDefined(heatCapacitySolid_);
672  Opm::Valgrind::CheckDefined(heatConductivity_);
673  }
674 
675 public:
681  const Evaluation& heatCapacitySolid() const
682  { return heatCapacitySolid_; }
683 
689  const Evaluation& heatConductivity() const
690  { return heatConductivity_; }
691 
692 private:
693  Evaluation heatCapacitySolid_;
694  Evaluation heatConductivity_;
695 };
696 
703 template <class TypeTag, bool enableEnergy>
705 
709 template <class TypeTag>
710 class EnergyExtensiveQuantities<TypeTag, /*enableEnergy=*/false>
711 {
712  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
713  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
714 
715 protected:
720  void update_(const ElementContext& elemCtx OPM_UNUSED,
721  unsigned faceIdx OPM_UNUSED,
722  unsigned timeIdx OPM_UNUSED)
723  {}
724 
725  template <class Context, class FluidState>
726  void updateBoundary_(const Context& context OPM_UNUSED,
727  unsigned bfIdx OPM_UNUSED,
728  unsigned timeIdx OPM_UNUSED,
729  const FluidState& fs OPM_UNUSED)
730  {}
731 
732 public:
736  Scalar temperatureGradNormal() const
737  {
738  OPM_THROW(std::logic_error, "Method temperatureGradNormal() does not "
739  "make sense for isothermal models");
740  }
741 
746  Scalar heatConductivity() const
747  {
748  OPM_THROW(std::logic_error, "Method heatConductivity() does not make "
749  "sense for isothermal models");
750  }
751 };
752 
756 template <class TypeTag>
757 class EnergyExtensiveQuantities<TypeTag, /*enableEnergy=*/true>
758 {
759  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
760  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
761  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
762  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
763 
764  enum { dimWorld = GridView::dimensionworld };
765  typedef Dune::FieldVector<Evaluation, dimWorld> EvalDimVector;
766  typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
767 
768 protected:
773  void update_(const ElementContext& elemCtx, unsigned faceIdx, unsigned timeIdx)
774  {
775  const auto& gradCalc = elemCtx.gradientCalculator();
776  Ewoms::TemperatureCallback<TypeTag> temperatureCallback(elemCtx);
777 
778  EvalDimVector temperatureGrad;
779  gradCalc.calculateGradient(temperatureGrad,
780  elemCtx,
781  faceIdx,
782  temperatureCallback);
783 
784  // scalar product of temperature gradient and scvf normal
785  const auto& face = elemCtx.stencil(/*timeIdx=*/0).interiorFace(faceIdx);
786 
787  temperatureGradNormal_ = 0;
788  for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
789  temperatureGradNormal_ += (face.normal()[dimIdx]*temperatureGrad[dimIdx]);
790 
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);
794 
795  // arithmetic mean
796  heatConductivity_ =
797  0.5 * (intQuantsInside.heatConductivity() + intQuantsOutside.heatConductivity());
798  Opm::Valgrind::CheckDefined(heatConductivity_);
799  }
800 
801  template <class Context, class FluidState>
802  void updateBoundary_(const Context& context, unsigned bfIdx, unsigned timeIdx, const FluidState& fs)
803  {
804  const auto& stencil = context.stencil(timeIdx);
805  const auto& face = stencil.boundaryFace(bfIdx);
806 
807  const auto& elemCtx = context.elementContext();
808  unsigned insideScvIdx = face.interiorIndex();
809  const auto& insideScv = elemCtx.stencil(timeIdx).subControlVolume(insideScvIdx);
810 
811  const auto& intQuantsInside = elemCtx.intensiveQuantities(insideScvIdx, timeIdx);
812  const auto& fsInside = intQuantsInside.fluidState();
813 
814  // distance between the center of the SCV and center of the boundary face
815  DimVector distVec = face.integrationPos();
816  distVec -= insideScv.geometry().center();
817 
818  Scalar tmp = 0;
819  for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
820  tmp += distVec[dimIdx] * face.normal()[dimIdx];
821  Scalar dist = tmp;
822 
823  // if the following assertation triggers, the center of the
824  // center of the interior SCV was not inside the element!
825  assert(dist > 0);
826 
827  // calculate the temperature gradient using two-point gradient
828  // appoximation
829  temperatureGradNormal_ =
830  (fs.temperature(/*phaseIdx=*/0) - fsInside.temperature(/*phaseIdx=*/0)) / dist;
831 
832  // take the value for heat conductivity from the interior finite volume
833  heatConductivity_ = intQuantsInside.heatConductivity();
834  }
835 
836 public:
840  const Evaluation& temperatureGradNormal() const
841  { return temperatureGradNormal_; }
842 
847  const Evaluation& heatConductivity() const
848  { return heatConductivity_; }
849 
850 private:
851  Evaluation temperatureGradNormal_;
852  Evaluation heatConductivity_;
853 };
854 
855 } // namespace Ewoms
856 
857 #endif
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 > &paramCache 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 > &paramCache, 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