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