blackoilenergymodules.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*/
28#ifndef EWOMS_BLACK_OIL_ENERGY_MODULE_HH
29#define EWOMS_BLACK_OIL_ENERGY_MODULE_HH
30
31#include "blackoilproperties.hh"
35
36#include <opm/material/common/Tabulated1DFunction.hpp>
37
38#include <opm/material/common/Valgrind.hpp>
39
40#include <dune/common/fvector.hh>
41
42#include <string>
43
44namespace Opm {
50template <class TypeTag, bool enableEnergyV = getPropValue<TypeTag, Properties::EnableEnergy>()>
52{
64
65 static constexpr unsigned temperatureIdx = Indices::temperatureIdx;
66 static constexpr unsigned contiEnergyEqIdx = Indices::contiEnergyEqIdx;
67
68 static constexpr unsigned enableEnergy = enableEnergyV;
69 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
70 static constexpr unsigned numPhases = FluidSystem::numPhases;
71
72public:
77 static void registerParameters()
78 {
79 if constexpr (enableEnergy)
81 }
82
86 static void registerOutputModules(Model& model,
87 Simulator& simulator)
88 {
89 if constexpr (enableEnergy)
90 model.addOutputModule(new VtkBlackOilEnergyModule<TypeTag>(simulator));
91 }
92
93 static bool primaryVarApplies(unsigned pvIdx)
94 {
95 if constexpr (enableEnergy)
96 return pvIdx == temperatureIdx;
97 else
98 return false;
99
100 }
101
102 static std::string primaryVarName([[maybe_unused]] unsigned pvIdx)
103 {
104 assert(primaryVarApplies(pvIdx));
105
106 return "temperature";
107 }
108
109 static Scalar primaryVarWeight([[maybe_unused]] unsigned pvIdx)
110 {
111 assert(primaryVarApplies(pvIdx));
112
113 // TODO: it may be beneficial to chose this differently.
114 return static_cast<Scalar>(1.0);
115 }
116
117 static bool eqApplies(unsigned eqIdx)
118 {
119 if constexpr (enableEnergy)
120 return eqIdx == contiEnergyEqIdx;
121 else
122 return false;
123 }
124
125 static std::string eqName([[maybe_unused]] unsigned eqIdx)
126 {
127 assert(eqApplies(eqIdx));
128
129 return "conti^energy";
130 }
131
132 static Scalar eqWeight([[maybe_unused]] unsigned eqIdx)
133 {
134 assert(eqApplies(eqIdx));
135
136 return 1.0;
137 }
138
139 // must be called after water storage is computed
140 template <class LhsEval>
141 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
142 const IntensiveQuantities& intQuants)
143 {
144 if constexpr (enableEnergy) {
145 const auto& poro = decay<LhsEval>(intQuants.porosity());
146
147 // accumulate the internal energy of the fluids
148 const auto& fs = intQuants.fluidState();
149 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
150 if (!FluidSystem::phaseIsActive(phaseIdx))
151 continue;
152
153 const auto& u = decay<LhsEval>(fs.internalEnergy(phaseIdx));
154 const auto& S = decay<LhsEval>(fs.saturation(phaseIdx));
155 const auto& rho = decay<LhsEval>(fs.density(phaseIdx));
156
157 storage[contiEnergyEqIdx] += poro*S*u*rho;
158 }
159
160 // add the internal energy of the rock
161 Scalar rockFraction = intQuants.rockFraction();
162 const auto& uRock = decay<LhsEval>(intQuants.rockInternalEnergy());
163 storage[contiEnergyEqIdx] += rockFraction*uRock;
164 storage[contiEnergyEqIdx] *= getPropValue<TypeTag, Properties::BlackOilEnergyScalingFactor>();
165 }
166 }
167
168 static void computeFlux([[maybe_unused]] RateVector& flux,
169 [[maybe_unused]] const ElementContext& elemCtx,
170 [[maybe_unused]] unsigned scvfIdx,
171 [[maybe_unused]] unsigned timeIdx)
172 {
173 if constexpr (enableEnergy) {
174 flux[contiEnergyEqIdx] = 0.0;
175
176 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
177 unsigned focusIdx = elemCtx.focusDofIndex();
178 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
179 if (!FluidSystem::phaseIsActive(phaseIdx))
180 continue;
181
182 unsigned upIdx = extQuants.upstreamIndex(phaseIdx);
183 if (upIdx == focusIdx)
184 addPhaseEnthalpyFlux_<Evaluation>(flux, phaseIdx, elemCtx, scvfIdx, timeIdx);
185 else
186 addPhaseEnthalpyFlux_<Scalar>(flux, phaseIdx, elemCtx, scvfIdx, timeIdx);
187 }
188
189 // diffusive energy flux
190 flux[contiEnergyEqIdx] += extQuants.energyFlux();
191 flux[contiEnergyEqIdx] *= getPropValue<TypeTag, Properties::BlackOilEnergyScalingFactor>();
192 }
193 }
194
195 static void addHeatFlux(RateVector& flux,
196 const Evaluation& heatFlux)
197 {
198 if constexpr (enableEnergy) {
199 // diffusive energy flux
200 flux[contiEnergyEqIdx] += heatFlux;
201 flux[contiEnergyEqIdx] *= getPropValue<TypeTag, Properties::BlackOilEnergyScalingFactor>();
202 }
203 }
204
205
206
207 template <class UpEval, class Eval, class FluidState>
208 static void addPhaseEnthalpyFluxes_(RateVector& flux,
209 unsigned phaseIdx,
210 const Eval& volumeFlux,
211 const FluidState& upFs)
212 {
213 flux[contiEnergyEqIdx] +=
214 decay<UpEval>(upFs.enthalpy(phaseIdx))
215 * decay<UpEval>(upFs.density(phaseIdx))
216 * volumeFlux;
217 }
218
219 template <class UpstreamEval>
220 static void addPhaseEnthalpyFlux_(RateVector& flux,
221 unsigned phaseIdx,
222 const ElementContext& elemCtx,
223 unsigned scvfIdx,
224 unsigned timeIdx)
225 {
226 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
227 unsigned upIdx = extQuants.upstreamIndex(phaseIdx);
228 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
229 const auto& fs = up.fluidState();
230 const auto& volFlux = extQuants.volumeFlux(phaseIdx);
231 addPhaseEnthalpyFluxes_<UpstreamEval>(flux,
232 phaseIdx,
233 volFlux,
234 fs);
235 }
236
237 static void addToEnthalpyRate(RateVector& flux,
238 const Evaluation& hRate)
239 {
240 if constexpr (enableEnergy)
241 flux[contiEnergyEqIdx] += hRate;
242 }
243
247 static void assignPrimaryVars(PrimaryVariables& priVars,
248 Scalar)
249 {
250 if constexpr (enableEnergy)
251 priVars[temperatureIdx] = temperatureIdx;
252 }
253
257 template <class FluidState>
258 static void assignPrimaryVars(PrimaryVariables& priVars,
259 const FluidState& fluidState)
260 {
261 if constexpr (enableEnergy)
262 priVars[temperatureIdx] = fluidState.temperature(/*phaseIdx=*/0);
263 }
264
268 static void updatePrimaryVars(PrimaryVariables& newPv,
269 const PrimaryVariables& oldPv,
270 const EqVector& delta)
271 {
272 if constexpr (enableEnergy)
273 // do a plain unchopped Newton update
274 newPv[temperatureIdx] = oldPv[temperatureIdx] - delta[temperatureIdx];
275 }
276
280 static Scalar computeUpdateError(const PrimaryVariables&,
281 const EqVector&)
282 {
283 // do not consider consider the cange of energy primary variables for
284 // convergence
285 // TODO: maybe this should be changed
286 return static_cast<Scalar>(0.0);
287 }
288
292 static Scalar computeResidualError(const EqVector& resid)
293 {
294 // do not weight the residual of energy when it comes to convergence
295 return std::abs(scalarValue(resid[contiEnergyEqIdx]));
296 }
297
298 template <class DofEntity>
299 static void serializeEntity(const Model& model, std::ostream& outstream, const DofEntity& dof)
300 {
301 if constexpr (enableEnergy) {
302 unsigned dofIdx = model.dofMapper().index(dof);
303 const PrimaryVariables& priVars = model.solution(/*timeIdx=*/0)[dofIdx];
304 outstream << priVars[temperatureIdx];
305 }
306 }
307
308 template <class DofEntity>
309 static void deserializeEntity(Model& model, std::istream& instream, const DofEntity& dof)
310 {
311 if constexpr (enableEnergy) {
312 unsigned dofIdx = model.dofMapper().index(dof);
313 PrimaryVariables& priVars0 = model.solution(/*timeIdx=*/0)[dofIdx];
314 PrimaryVariables& priVars1 = model.solution(/*timeIdx=*/1)[dofIdx];
315
316 instream >> priVars0[temperatureIdx];
317
318 // set the primary variables for the beginning of the current time step.
319 priVars1 = priVars0[temperatureIdx];
320 }
321 }
322};
323
331template <class TypeTag, bool enableEnergyV = getPropValue<TypeTag, Properties::EnableEnergy>()>
333{
335
345
347
348 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
349 static constexpr int temperatureIdx = Indices::temperatureIdx;
350 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
351
352
353public:
358 void updateTemperature_(const ElementContext& elemCtx,
359 unsigned dofIdx,
360 unsigned timeIdx)
361 {
362 auto& fs = asImp_().fluidState_;
363 const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
364
365 // set temperature
366 fs.setTemperature(priVars.makeEvaluation(temperatureIdx, timeIdx, elemCtx.linearizationType()));
367 }
368
373 void updateTemperature_([[maybe_unused]] const Problem& problem,
374 const PrimaryVariables& priVars,
375 [[maybe_unused]] unsigned globalDofIdx,
376 const unsigned timeIdx,
377 const LinearizationType& lintype)
378 {
379 auto& fs = asImp_().fluidState_;
380 fs.setTemperature(priVars.makeEvaluation(temperatureIdx, timeIdx, lintype));
381 }
382
387 void updateEnergyQuantities_(const ElementContext& elemCtx,
388 unsigned dofIdx,
389 unsigned timeIdx,
390 const typename FluidSystem::template ParameterCache<Evaluation>& paramCache)
391 {
392 auto& fs = asImp_().fluidState_;
393
394 // compute the specific enthalpy of the fluids, the specific enthalpy of the rock
395 // and the thermal condictivity coefficients
396 for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
397 if (!FluidSystem::phaseIsActive(phaseIdx)) {
398 continue;
399 }
400
401 const auto& h = FluidSystem::enthalpy(fs, paramCache, phaseIdx);
402 fs.setEnthalpy(phaseIdx, h);
403 }
404
405 const auto& solidEnergyLawParams = elemCtx.problem().solidEnergyLawParams(elemCtx, dofIdx, timeIdx);
406 rockInternalEnergy_ = SolidEnergyLaw::solidInternalEnergy(solidEnergyLawParams, fs);
407
408 const auto& thermalConductionLawParams = elemCtx.problem().thermalConductionLawParams(elemCtx, dofIdx, timeIdx);
409 totalThermalConductivity_ = ThermalConductionLaw::thermalConductivity(thermalConductionLawParams, fs);
410
411 // Retrieve the rock fraction from the problem
412 // Usually 1 - porosity, but if pvmult is used to modify porosity
413 // we will apply the same multiplier to the rock fraction
414 // i.e. pvmult*(1 - porosity) and thus interpret multpv as a volume
415 // multiplier. This is to avoid negative rock volume for pvmult*porosity > 1
416 const unsigned cell_idx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
417 rockFraction_ = elemCtx.problem().rockFraction(cell_idx, timeIdx);
418 }
419
420 const Evaluation& rockInternalEnergy() const
421 { return rockInternalEnergy_; }
422
423 const Evaluation& totalThermalConductivity() const
424 { return totalThermalConductivity_; }
425
426 const Scalar& rockFraction() const
427 { return rockFraction_; }
428
429protected:
430 Implementation& asImp_()
431 { return *static_cast<Implementation*>(this); }
432
436};
437
438template <class TypeTag>
440{
447
449 static constexpr bool enableTemperature = getPropValue<TypeTag, Properties::EnableTemperature>();
450
451public:
452 void updateTemperature_([[maybe_unused]] const ElementContext& elemCtx,
453 [[maybe_unused]] unsigned dofIdx,
454 [[maybe_unused]] unsigned timeIdx)
455 {
456 if constexpr (enableTemperature) {
457 // even if energy is conserved, the temperature can vary over the spatial
458 // domain if the EnableTemperature property is set to true
459 auto& fs = asImp_().fluidState_;
460 const Scalar T = elemCtx.problem().temperature(elemCtx, dofIdx, timeIdx);
461 fs.setTemperature(T);
462 }
463 }
464
465 template<class Problem>
466 void updateTemperature_([[maybe_unused]] const Problem& problem,
467 [[maybe_unused]] const PrimaryVariables& priVars,
468 [[maybe_unused]] unsigned globalDofIdx,
469 [[maybe_unused]] unsigned timeIdx,
470 [[maybe_unused]] const LinearizationType& lintype
471 )
472 {
473 if constexpr (enableTemperature) {
474 auto& fs = asImp_().fluidState_;
475 // even if energy is conserved, the temperature can vary over the spatial
476 // domain if the EnableTemperature property is set to true
477 const Scalar T = problem.temperature(globalDofIdx, timeIdx);
478 fs.setTemperature(T);
479 }
480 }
481
482 void updateEnergyQuantities_(const ElementContext&,
483 unsigned,
484 unsigned,
485 const typename FluidSystem::template ParameterCache<Evaluation>&)
486 { }
487
488 const Evaluation& rockInternalEnergy() const
489 { throw std::logic_error("Requested the rock internal energy, which is "
490 "unavailable because energy is not conserved"); }
491
492 const Evaluation& totalThermalConductivity() const
493 { throw std::logic_error("Requested the total thermal conductivity, which is "
494 "unavailable because energy is not conserved"); }
495
496protected:
497 Implementation& asImp_()
498 { return *static_cast<Implementation*>(this); }
499};
500
501
509template <class TypeTag, bool enableEnergyV = getPropValue<TypeTag, Properties::EnableEnergy>()>
511{
513
521
522 using Toolbox = MathToolbox<Evaluation>;
523
525
526 static const int dimWorld = GridView::dimensionworld;
527 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
528 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
529public:
530 template<class FluidState>
531 static void updateEnergy(Evaluation& energyFlux,
532 const unsigned& focusDofIndex,
533 const unsigned& inIdx,
534 const unsigned& exIdx,
535 const IntensiveQuantities& inIq,
536 const IntensiveQuantities& exIq,
537 const FluidState& inFs,
538 const FluidState& exFs,
539 const Scalar& inAlpha,
540 const Scalar& outAlpha,
541 const Scalar& faceArea)
542 {
543 Evaluation deltaT;
544 if (focusDofIndex == inIdx)
545 deltaT =
546 decay<Scalar>(exFs.temperature(/*phaseIdx=*/0))
547 - inFs.temperature(/*phaseIdx=*/0);
548 else if (focusDofIndex == exIdx)
549 deltaT =
550 exFs.temperature(/*phaseIdx=*/0)
551 - decay<Scalar>(inFs.temperature(/*phaseIdx=*/0));
552 else
553 deltaT =
554 decay<Scalar>(exFs.temperature(/*phaseIdx=*/0))
555 - decay<Scalar>(inFs.temperature(/*phaseIdx=*/0));
556
557 Evaluation inLambda;
558 if (focusDofIndex == inIdx)
559 inLambda = inIq.totalThermalConductivity();
560 else
561 inLambda = decay<Scalar>(inIq.totalThermalConductivity());
562
563 Evaluation exLambda;
564 if (focusDofIndex == exIdx)
565 exLambda = exIq.totalThermalConductivity();
566 else
567 exLambda = decay<Scalar>(exIq.totalThermalConductivity());
568
569 Evaluation H;
570 const Evaluation& inH = inLambda*inAlpha;
571 const Evaluation& exH = exLambda*outAlpha;
572 if (inH > 0 && exH > 0) {
573 // compute the "thermal transmissibility". In contrast to the normal
574 // transmissibility this cannot be done as a preprocessing step because the
575 // average thermal conductivity is analogous to the permeability but
576 // depends on the solution.
577 H = 1.0/(1.0/inH + 1.0/exH);
578 }
579 else
580 H = 0.0;
581
582 energyFlux = deltaT * (-H/faceArea);
583 }
584
585 void updateEnergy(const ElementContext& elemCtx,
586 unsigned scvfIdx,
587 unsigned timeIdx)
588 {
589 const auto& stencil = elemCtx.stencil(timeIdx);
590 const auto& scvf = stencil.interiorFace(scvfIdx);
591
592 const Scalar faceArea = scvf.area();
593 const unsigned inIdx = scvf.interiorIndex();
594 const unsigned exIdx = scvf.exteriorIndex();
595 const auto& inIq = elemCtx.intensiveQuantities(inIdx, timeIdx);
596 const auto& exIq = elemCtx.intensiveQuantities(exIdx, timeIdx);
597 const auto& inFs = inIq.fluidState();
598 const auto& exFs = exIq.fluidState();
599 const Scalar inAlpha = elemCtx.problem().thermalHalfTransmissibilityIn(elemCtx, scvfIdx, timeIdx);
600 const Scalar outAlpha = elemCtx.problem().thermalHalfTransmissibilityOut(elemCtx, scvfIdx, timeIdx);
601 updateEnergy(energyFlux_,
602 elemCtx.focusDofIndex(),
603 inIdx,
604 exIdx,
605 inIq,
606 exIq,
607 inFs,
608 exFs,
609 inAlpha,
610 outAlpha,
611 faceArea);
612 }
613
614 template <class Context, class BoundaryFluidState>
615 void updateEnergyBoundary(const Context& ctx,
616 unsigned scvfIdx,
617 unsigned timeIdx,
618 const BoundaryFluidState& boundaryFs)
619 {
620 const auto& stencil = ctx.stencil(timeIdx);
621 const auto& scvf = stencil.boundaryFace(scvfIdx);
622
623 unsigned inIdx = scvf.interiorIndex();
624 const auto& inIq = ctx.intensiveQuantities(inIdx, timeIdx);
625 const auto& focusDofIdx = ctx.focusDofIndex();
626 const Scalar alpha = ctx.problem().thermalHalfTransmissibilityBoundary(ctx, scvfIdx);
627 updateEnergyBoundary(energyFlux_, inIq, focusDofIdx, inIdx, alpha, boundaryFs);
628 }
629
630 template <class BoundaryFluidState>
631 static void updateEnergyBoundary(Evaluation& energyFlux,
632 const IntensiveQuantities& inIq,
633 unsigned focusDofIndex,
634 unsigned inIdx,
635 Scalar alpha,
636 const BoundaryFluidState& boundaryFs)
637 {
638 const auto& inFs = inIq.fluidState();
639 Evaluation deltaT;
640 if (focusDofIndex == inIdx)
641 deltaT =
642 boundaryFs.temperature(/*phaseIdx=*/0)
643 - inFs.temperature(/*phaseIdx=*/0);
644 else
645 deltaT =
646 decay<Scalar>(boundaryFs.temperature(/*phaseIdx=*/0))
647 - decay<Scalar>(inFs.temperature(/*phaseIdx=*/0));
648
649 Evaluation lambda;
650 if (focusDofIndex == inIdx)
651 lambda = inIq.totalThermalConductivity();
652 else
653 lambda = decay<Scalar>(inIq.totalThermalConductivity());
654
655
656 if (lambda > 0.0) {
657 // compute the "thermal transmissibility". In contrast to the normal
658 // transmissibility this cannot be done as a preprocessing step because the
659 // average thermal conductivity is analogous to the permeability but depends
660 // on the solution.
661 energyFlux = deltaT*lambda*(-alpha);
662 }
663 else
664 energyFlux = 0.0;
665 }
666
667 const Evaluation& energyFlux() const
668 { return energyFlux_; }
669
670private:
671 Implementation& asImp_()
672 { return *static_cast<Implementation*>(this); }
673
674 Evaluation energyFlux_;
675};
676
677template <class TypeTag>
679{
684public:
685 template<class FluidState>
686 static void updateEnergy(Evaluation& /*energyFlux*/,
687 const unsigned& /*focusDofIndex*/,
688 const unsigned& /*inIdx*/,
689 const unsigned& /*exIdx*/,
690 const IntensiveQuantities& /*inIq*/,
691 const IntensiveQuantities& /*exIq*/,
692 const FluidState& /*inFs*/,
693 const FluidState& /*exFs*/,
694 const Scalar& /*inAlpha*/,
695 const Scalar& /*outAlpha*/,
696 const Scalar& /*faceArea*/)
697 {}
698
699 void updateEnergy(const ElementContext&,
700 unsigned,
701 unsigned)
702 {}
703
704 template <class Context, class BoundaryFluidState>
705 void updateEnergyBoundary(const Context&,
706 unsigned,
707 unsigned,
708 const BoundaryFluidState&)
709 {}
710
711 template <class BoundaryFluidState>
712 static void updateEnergyBoundary(Evaluation& /*heatFlux*/,
713 const IntensiveQuantities& /*inIq*/,
714 unsigned /*focusDofIndex*/,
715 unsigned /*inIdx*/,
716 unsigned /*timeIdx*/,
717 Scalar /*alpha*/,
718 const BoundaryFluidState& /*boundaryFs*/)
719 {}
720 const Evaluation& energyFlux() const
721 { throw std::logic_error("Requested the energy flux, but energy is not conserved"); }
722};
723
724
725} // namespace Opm
726
727#endif
Declares the properties required by the black oil model.
void updateEnergyBoundary(const Context &, unsigned, unsigned, const BoundaryFluidState &)
Definition: blackoilenergymodules.hh:705
static void updateEnergy(Evaluation &, const unsigned &, const unsigned &, const unsigned &, const IntensiveQuantities &, const IntensiveQuantities &, const FluidState &, const FluidState &, const Scalar &, const Scalar &, const Scalar &)
Definition: blackoilenergymodules.hh:686
static void updateEnergyBoundary(Evaluation &, const IntensiveQuantities &, unsigned, unsigned, unsigned, Scalar, const BoundaryFluidState &)
Definition: blackoilenergymodules.hh:712
void updateEnergy(const ElementContext &, unsigned, unsigned)
Definition: blackoilenergymodules.hh:699
const Evaluation & energyFlux() const
Definition: blackoilenergymodules.hh:720
Provides the energy specific extensive quantities to the generic black-oil module's extensive quantit...
Definition: blackoilenergymodules.hh:511
const Evaluation & energyFlux() const
Definition: blackoilenergymodules.hh:667
static void updateEnergy(Evaluation &energyFlux, const unsigned &focusDofIndex, const unsigned &inIdx, const unsigned &exIdx, const IntensiveQuantities &inIq, const IntensiveQuantities &exIq, const FluidState &inFs, const FluidState &exFs, const Scalar &inAlpha, const Scalar &outAlpha, const Scalar &faceArea)
Definition: blackoilenergymodules.hh:531
static void updateEnergyBoundary(Evaluation &energyFlux, const IntensiveQuantities &inIq, unsigned focusDofIndex, unsigned inIdx, Scalar alpha, const BoundaryFluidState &boundaryFs)
Definition: blackoilenergymodules.hh:631
void updateEnergyBoundary(const Context &ctx, unsigned scvfIdx, unsigned timeIdx, const BoundaryFluidState &boundaryFs)
Definition: blackoilenergymodules.hh:615
void updateEnergy(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilenergymodules.hh:585
void updateEnergyQuantities_(const ElementContext &, unsigned, unsigned, const typename FluidSystem::template ParameterCache< Evaluation > &)
Definition: blackoilenergymodules.hh:482
const Evaluation & rockInternalEnergy() const
Definition: blackoilenergymodules.hh:488
const Evaluation & totalThermalConductivity() const
Definition: blackoilenergymodules.hh:492
void updateTemperature_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: blackoilenergymodules.hh:452
void updateTemperature_(const Problem &problem, const PrimaryVariables &priVars, unsigned globalDofIdx, unsigned timeIdx, const LinearizationType &lintype)
Definition: blackoilenergymodules.hh:466
Implementation & asImp_()
Definition: blackoilenergymodules.hh:497
Provides the volumetric quantities required for the equations needed by the energys extension of the ...
Definition: blackoilenergymodules.hh:333
Scalar rockFraction_
Definition: blackoilenergymodules.hh:435
Implementation & asImp_()
Definition: blackoilenergymodules.hh:430
Evaluation totalThermalConductivity_
Definition: blackoilenergymodules.hh:434
const Evaluation & rockInternalEnergy() const
Definition: blackoilenergymodules.hh:420
Evaluation rockInternalEnergy_
Definition: blackoilenergymodules.hh:433
void updateTemperature_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the temperature of the intensive quantity's fluid state.
Definition: blackoilenergymodules.hh:358
void updateEnergyQuantities_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, const typename FluidSystem::template ParameterCache< Evaluation > &paramCache)
Compute the intensive quantities needed to handle energy conservation.
Definition: blackoilenergymodules.hh:387
void updateTemperature_(const Problem &problem, const PrimaryVariables &priVars, unsigned globalDofIdx, const unsigned timeIdx, const LinearizationType &lintype)
Update the temperature of the intensive quantity's fluid state.
Definition: blackoilenergymodules.hh:373
const Scalar & rockFraction() const
Definition: blackoilenergymodules.hh:426
const Evaluation & totalThermalConductivity() const
Definition: blackoilenergymodules.hh:423
Contains the high level supplements required to extend the black oil model by energy.
Definition: blackoilenergymodules.hh:52
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilenergymodules.hh:168
static void registerParameters()
Register all run-time parameters for the black-oil energy module.
Definition: blackoilenergymodules.hh:77
GetPropType< TypeTag, Properties::ExtensiveQuantities > ExtensiveQuantities
Definition: blackoilenergymodules.hh:73
static Scalar computeResidualError(const EqVector &resid)
Return how much a residual is considered an error.
Definition: blackoilenergymodules.hh:292
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilenergymodules.hh:309
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilenergymodules.hh:141
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilenergymodules.hh:93
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilenergymodules.hh:102
static std::string eqName(unsigned eqIdx)
Definition: blackoilenergymodules.hh:125
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilenergymodules.hh:299
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilenergymodules.hh:109
static void updatePrimaryVars(PrimaryVariables &newPv, const PrimaryVariables &oldPv, const EqVector &delta)
Do a Newton-Raphson update the primary variables of the energys.
Definition: blackoilenergymodules.hh:268
static void addToEnthalpyRate(RateVector &flux, const Evaluation &hRate)
Definition: blackoilenergymodules.hh:237
static void addHeatFlux(RateVector &flux, const Evaluation &heatFlux)
Definition: blackoilenergymodules.hh:195
static void registerOutputModules(Model &model, Simulator &simulator)
Register all energy specific VTK and ECL output modules.
Definition: blackoilenergymodules.hh:86
static void addPhaseEnthalpyFlux_(RateVector &flux, unsigned phaseIdx, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilenergymodules.hh:220
static void addPhaseEnthalpyFluxes_(RateVector &flux, unsigned phaseIdx, const Eval &volumeFlux, const FluidState &upFs)
Definition: blackoilenergymodules.hh:208
static bool eqApplies(unsigned eqIdx)
Definition: blackoilenergymodules.hh:117
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar)
Assign the energy specific primary variables to a PrimaryVariables object.
Definition: blackoilenergymodules.hh:247
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilenergymodules.hh:280
static void assignPrimaryVars(PrimaryVariables &priVars, const FluidState &fluidState)
Assign the energy specific primary variables to a PrimaryVariables object.
Definition: blackoilenergymodules.hh:258
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilenergymodules.hh:132
VTK output module for the black oil model's energy related quantities.
Definition: vtkblackoilenergymodule.hh:62
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
This method contains all callback classes for quantities that are required by some extensive quantiti...
Definition: linearizationtype.hh:35