24#ifndef DUMUX_POROUSMEDIUM_FLUXVARIABLESCACHE_FILLER_HH
25#define DUMUX_POROUSMEDIUM_FLUXVARIABLESCACHE_FILLER_HH
38template<
class TypeTag, DiscretizationMethod discMethod>
47template<
class TypeTag>
51template<
class TypeTag>
58 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
59 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
62 using Element =
typename GridView::template Codim<0>::Entity;
64 static constexpr bool advectionEnabled = ModelTraits::enableAdvection();
65 static constexpr bool diffusionEnabled = ModelTraits::enableMolecularDiffusion();
66 static constexpr bool heatConductionEnabled = ModelTraits::enableEnergyBalance();
68 static constexpr bool advectionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentAdvection>();
69 static constexpr bool diffusionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentMolecularDiffusion>();
70 static constexpr bool heatConductionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentHeatConduction>();
74 static constexpr bool isSolDependent = (advectionEnabled && advectionIsSolDependent) ||
75 (diffusionEnabled && diffusionIsSolDependent) ||
76 (heatConductionEnabled && heatConductionIsSolDependent);
80 : problemPtr_(&problem) {}
93 template<
class FluxVariablesCacheContainer,
class FluxVariablesCache>
94 void fill(FluxVariablesCacheContainer& fluxVarsCacheContainer,
95 FluxVariablesCache& scvfFluxVarsCache,
96 const Element& element,
97 const FVElementGeometry& fvGeometry,
98 const ElementVolumeVariables& elemVolVars,
99 const SubControlVolumeFace& scvf,
100 bool forceUpdateAll =
false)
105 if constexpr (advectionEnabled)
106 fillAdvection_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
107 if constexpr (diffusionEnabled)
108 fillDiffusion_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
109 if constexpr (heatConductionEnabled)
110 fillHeatConduction_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
114 if constexpr (advectionEnabled && advectionIsSolDependent)
115 fillAdvection_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
116 if constexpr (diffusionEnabled && diffusionIsSolDependent)
117 fillDiffusion_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
118 if constexpr (heatConductionEnabled && heatConductionIsSolDependent)
119 fillHeatConduction_(scvfFluxVarsCache, element, fvGeometry, elemVolVars, scvf);
125 const Problem& problem()
const
126 {
return *problemPtr_; }
129 template<
class FluxVariablesCache>
130 void fillAdvection_(FluxVariablesCache& scvfFluxVarsCache,
131 const Element& element,
132 const FVElementGeometry& fvGeometry,
133 const ElementVolumeVariables& elemVolVars,
134 const SubControlVolumeFace& scvf)
136 using AdvectionType = GetPropType<TypeTag, Properties::AdvectionType>;
137 using AdvectionFiller =
typename AdvectionType::Cache::Filler;
140 AdvectionFiller::fill(scvfFluxVarsCache, problem(), element, fvGeometry, elemVolVars, scvf, *
this);
144 template<
class FluxVariablesCache>
145 void fillDiffusion_(FluxVariablesCache& scvfFluxVarsCache,
146 const Element& element,
147 const FVElementGeometry& fvGeometry,
148 const ElementVolumeVariables& elemVolVars,
149 const SubControlVolumeFace& scvf)
151 using DiffusionType = GetPropType<TypeTag, Properties::MolecularDiffusionType>;
152 using DiffusionFiller =
typename DiffusionType::Cache::Filler;
153 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
155 static constexpr int numPhases = ModelTraits::numFluidPhases();
156 static constexpr int numComponents = ModelTraits::numFluidComponents();
159 if constexpr (FluidSystem::isTracerFluidSystem())
160 for (
unsigned int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
161 for (
unsigned int compIdx = 0; compIdx < numComponents; ++compIdx)
162 DiffusionFiller::fill(scvfFluxVarsCache, phaseIdx, compIdx, problem(), element, fvGeometry, elemVolVars, scvf, *
this);
164 for (
unsigned int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
165 for (
unsigned int compIdx = 0; compIdx < numComponents; ++compIdx)
166 if (compIdx != FluidSystem::getMainComponent(phaseIdx))
167 DiffusionFiller::fill(scvfFluxVarsCache, phaseIdx, compIdx, problem(), element, fvGeometry, elemVolVars, scvf, *
this);
171 template<
class FluxVariablesCache>
172 void fillHeatConduction_(FluxVariablesCache& scvfFluxVarsCache,
173 const Element& element,
174 const FVElementGeometry& fvGeometry,
175 const ElementVolumeVariables& elemVolVars,
176 const SubControlVolumeFace& scvf)
178 using HeatConductionType = GetPropType<TypeTag, Properties::HeatConductionType>;
179 using HeatConductionFiller =
typename HeatConductionType::Cache::Filler;
182 HeatConductionFiller::fill(scvfFluxVarsCache, problem(), element, fvGeometry, elemVolVars, scvf, *
this);
185 const Problem* problemPtr_;
189template<
class TypeTag>
195 using Element =
typename GridView::template Codim<0>::Entity;
198 using FVElementGeometry =
typename GridGeometry::LocalView;
199 using MpfaHelper =
typename GridGeometry::MpfaHelper;
200 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
205 using PrimaryDataHandle =
typename ElementFluxVariablesCache::PrimaryIvDataHandle;
206 using PrimaryLocalFaceData =
typename PrimaryInteractionVolume::Traits::LocalFaceData;
208 using SecondaryDataHandle =
typename ElementFluxVariablesCache::SecondaryIvDataHandle;
209 using SecondaryLocalFaceData =
typename SecondaryInteractionVolume::Traits::LocalFaceData;
211 static constexpr int dim = GridView::dimension;
212 static constexpr int dimWorld = GridView::dimensionworld;
214 static constexpr bool advectionEnabled = ModelTraits::enableAdvection();
215 static constexpr bool diffusionEnabled = ModelTraits::enableMolecularDiffusion();
216 static constexpr bool heatConductionEnabled = ModelTraits::enableEnergyBalance();
218 static constexpr bool advectionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentAdvection>();
219 static constexpr bool diffusionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentMolecularDiffusion>();
220 static constexpr bool heatConductionIsSolDependent = getPropValue<TypeTag, Properties::SolutionDependentHeatConduction>();
226 static constexpr bool isSolDependent =
true;
230 : problemPtr_(&problem) {}
244 template<
class FluxVarsCacheStorage,
class FluxVariablesCache,
class IVDataStorage>
245 void fill(FluxVarsCacheStorage& fluxVarsCacheStorage,
246 FluxVariablesCache& scvfFluxVarsCache,
247 IVDataStorage& ivDataStorage,
248 const Element& element,
249 const FVElementGeometry& fvGeometry,
250 const ElementVolumeVariables& elemVolVars,
251 const SubControlVolumeFace& scvf,
252 bool forceUpdateAll =
false)
255 elementPtr_ = &element;
256 fvGeometryPtr_ = &fvGeometry;
257 elemVolVarsPtr_ = &elemVolVars;
258 const auto& gridGeometry = fvGeometry.gridGeometry();
263 if (gridGeometry.vertexUsesSecondaryInteractionVolume(scvf.vertexIndex()))
268 const auto ivIndexInContainer = ivDataStorage.secondaryInteractionVolumes.size();
269 const auto& indexSet = gridGeometry.gridInteractionVolumeIndexSets().secondaryIndexSet(scvf);
270 ivDataStorage.secondaryInteractionVolumes.emplace_back();
271 secondaryIv_ = &ivDataStorage.secondaryInteractionVolumes.back();
272 secondaryIv_->bind(indexSet, problem(), fvGeometry);
275 ivDataStorage.secondaryDataHandles.emplace_back();
276 secondaryIvDataHandle_ = &ivDataStorage.secondaryDataHandles.back();
277 prepareDataHandle_(*secondaryIv_, *secondaryIvDataHandle_, forceUpdateAll);
280 fillCachesInInteractionVolume_<FluxVariablesCache>(fluxVarsCacheStorage, *secondaryIv_, ivIndexInContainer);
285 const auto ivIndexInContainer = scvfFluxVarsCache.ivIndexInContainer();
286 secondaryIv_ = &ivDataStorage.secondaryInteractionVolumes[ivIndexInContainer];
287 secondaryIvDataHandle_ = &ivDataStorage.secondaryDataHandles[ivIndexInContainer];
288 prepareDataHandle_(*secondaryIv_, *secondaryIvDataHandle_, forceUpdateAll);
291 fillCachesInInteractionVolume_<FluxVariablesCache>(fluxVarsCacheStorage, *secondaryIv_, ivIndexInContainer);
301 const auto ivIndexInContainer = ivDataStorage.primaryInteractionVolumes.size();
302 const auto& indexSet = gridGeometry.gridInteractionVolumeIndexSets().primaryIndexSet(scvf);
303 ivDataStorage.primaryInteractionVolumes.emplace_back();
304 primaryIv_ = &ivDataStorage.primaryInteractionVolumes.back();
305 primaryIv_->bind(indexSet, problem(), fvGeometry);
308 ivDataStorage.primaryDataHandles.emplace_back();
309 primaryIvDataHandle_ = &ivDataStorage.primaryDataHandles.back();
310 prepareDataHandle_(*primaryIv_, *primaryIvDataHandle_, forceUpdateAll);
313 fillCachesInInteractionVolume_<FluxVariablesCache>(fluxVarsCacheStorage, *primaryIv_, ivIndexInContainer);
318 const auto ivIndexInContainer = scvfFluxVarsCache.ivIndexInContainer();
319 primaryIv_ = &ivDataStorage.primaryInteractionVolumes[ivIndexInContainer];
320 primaryIvDataHandle_ = &ivDataStorage.primaryDataHandles[ivIndexInContainer];
321 prepareDataHandle_(*primaryIv_, *primaryIvDataHandle_, forceUpdateAll);
324 fillCachesInInteractionVolume_<FluxVariablesCache>(fluxVarsCacheStorage, *primaryIv_, ivIndexInContainer);
331 {
return *primaryIv_; }
335 {
return *secondaryIv_; }
339 {
return *primaryIvDataHandle_; }
343 {
return *secondaryIvDataHandle_; }
347 {
return *primaryLocalFaceData_; }
351 {
return *secondaryLocalFaceData_; }
355 const Problem& problem()
const {
return *problemPtr_; }
356 const Element& element()
const {
return *elementPtr_; }
357 const FVElementGeometry& fvGeometry()
const {
return *fvGeometryPtr_; }
358 const ElementVolumeVariables& elemVolVars()
const {
return *elemVolVarsPtr_; }
361 template<
class FluxVariablesCache,
class FluxVarsCacheStorage,
class InteractionVolume>
362 void fillCachesInInteractionVolume_(FluxVarsCacheStorage& fluxVarsCacheStorage,
363 InteractionVolume& iv,
364 unsigned int ivIndexInContainer)
367 static constexpr bool isSecondary = MpfaHelper::considerSecondaryIVs()
368 && std::is_same_v<InteractionVolume, SecondaryInteractionVolume>;
373 const auto numGlobalScvfs = iv.localFaceData().size();
374 std::vector<const SubControlVolumeFace*> ivScvfs(numGlobalScvfs);
375 std::vector<FluxVariablesCache*> ivFluxVarCaches(numGlobalScvfs);
378 for (
const auto& d : iv.localFaceData())
381 const auto& scvfJ = fvGeometry().scvf(d.gridScvfIndex());
383 ivFluxVarCaches[i] = &fluxVarsCacheStorage[scvfJ];
384 ivFluxVarCaches[i]->setIvIndexInContainer(ivIndexInContainer);
385 ivFluxVarCaches[i]->setUpdateStatus(
true);
386 ivFluxVarCaches[i]->setSecondaryIvUsage(isSecondary);
387 ivFluxVarCaches[i]->setIvLocalFaceIndex(d.ivLocalScvfIndex());
389 if (d.isOutsideFace())
390 ivFluxVarCaches[i]->setIndexInOutsideFaces(d.scvfLocalOutsideScvfIndex());
394 if constexpr (advectionEnabled)
395 fillAdvection_(iv, ivScvfs, ivFluxVarCaches);
396 if constexpr (diffusionEnabled)
397 fillDiffusion_(iv, ivScvfs, ivFluxVarCaches);
398 if constexpr (heatConductionEnabled)
399 fillHeatConduction_(iv, ivScvfs, ivFluxVarCaches);
403 template<
class InteractionVolume,
class FluxVariablesCache>
404 void fillAdvection_(InteractionVolume& iv,
405 const std::vector<const SubControlVolumeFace*>& ivScvfs,
406 const std::vector<FluxVariablesCache*>& ivFluxVarCaches)
408 using AdvectionType = GetPropType<TypeTag, Properties::AdvectionType>;
409 using AdvectionFiller =
typename AdvectionType::Cache::Filler;
412 for (
unsigned int i = 0; i < iv.localFaceData().size(); ++i)
416 if (std::is_same_v<PrimaryInteractionVolume, SecondaryInteractionVolume>)
419 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
420 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
422 else if (std::is_same_v<InteractionVolume, PrimaryInteractionVolume>)
423 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
425 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
428 AdvectionFiller::fill(*ivFluxVarCaches[i],
430 iv.element(iv.localFaceData()[i].ivLocalInsideScvIndex()),
439 template<
class InteractionVolume,
class FluxVariablesCache>
440 void fillDiffusion_(InteractionVolume& iv,
441 const std::vector<const SubControlVolumeFace*>& ivScvfs,
442 const std::vector<FluxVariablesCache*>& ivFluxVarCaches)
444 using DiffusionType = GetPropType<TypeTag, Properties::MolecularDiffusionType>;
445 using DiffusionFiller =
typename DiffusionType::Cache::Filler;
447 static constexpr int numPhases = ModelTraits::numFluidPhases();
448 static constexpr int numComponents = ModelTraits::numFluidComponents();
450 for (
unsigned int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
452 for (
unsigned int compIdx = 0; compIdx < numComponents; ++compIdx)
454 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
455 if constexpr (!FluidSystem::isTracerFluidSystem())
456 if (compIdx == FluidSystem::getMainComponent(phaseIdx))
460 for (
unsigned int i = 0; i < iv.localFaceData().size(); ++i)
464 if constexpr (std::is_same_v<PrimaryInteractionVolume, SecondaryInteractionVolume>)
467 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
468 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
470 else if constexpr (std::is_same_v<InteractionVolume, PrimaryInteractionVolume>)
471 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
473 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
476 DiffusionFiller::fill(*ivFluxVarCaches[i],
480 iv.element(iv.localFaceData()[i].ivLocalInsideScvIndex()),
491 template<
class InteractionVolume,
class FluxVariablesCache>
492 void fillHeatConduction_(InteractionVolume& iv,
493 const std::vector<const SubControlVolumeFace*>& ivScvfs,
494 const std::vector<FluxVariablesCache*>& ivFluxVarCaches)
496 using HeatConductionType = GetPropType<TypeTag, Properties::HeatConductionType>;
497 using HeatConductionFiller =
typename HeatConductionType::Cache::Filler;
500 for (
unsigned int i = 0; i < iv.localFaceData().size(); ++i)
504 if constexpr (std::is_same_v<PrimaryInteractionVolume, SecondaryInteractionVolume>)
507 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
508 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
510 else if constexpr (std::is_same_v<InteractionVolume, PrimaryInteractionVolume>)
511 primaryLocalFaceData_ = &(iv.localFaceData()[i]);
513 secondaryLocalFaceData_ = &(iv.localFaceData()[i]);
516 HeatConductionFiller::fill(*ivFluxVarCaches[i],
518 iv.element(iv.localFaceData()[i].ivLocalInsideScvIndex()),
527 template<
class InteractionVolume,
class DataHandle>
528 void prepareDataHandle_([[maybe_unused]] InteractionVolume& iv, [[maybe_unused]] DataHandle& handle, [[maybe_unused]]
bool forceUpdate)
531 if constexpr (advectionEnabled)
533 using AdvectionType = GetPropType<TypeTag, Properties::AdvectionType>;
535 prepareAdvectionHandle_(iv, handle, forceUpdate);
539 if constexpr (diffusionEnabled)
541 using DiffusionType = GetPropType<TypeTag, Properties::MolecularDiffusionType>;
543 prepareDiffusionHandles_(iv, handle, forceUpdate);
547 if constexpr (heatConductionEnabled)
549 using HeatConductionType = GetPropType<TypeTag, Properties::HeatConductionType>;
551 prepareHeatConductionHandle_(iv, handle, forceUpdate);
556 template<
class InteractionVolume,
class DataHandle>
557 void prepareAdvectionHandle_(InteractionVolume& iv, DataHandle& handle,
bool forceUpdateAll)
560 using Traits =
typename InteractionVolume::Traits;
561 using IvLocalAssembler =
typename Traits::template LocalAssembler<Problem, FVElementGeometry, ElementVolumeVariables>;
562 IvLocalAssembler localAssembler(problem(), fvGeometry(), elemVolVars());
565 auto getK = [] (
const auto& volVars) {
return volVars.permeability(); };
568 if (forceUpdateAll || advectionIsSolDependent)
569 localAssembler.assembleMatrices(handle.advectionHandle(), iv, getK);
572 for (
unsigned int pIdx = 0; pIdx < ModelTraits::numFluidPhases(); ++pIdx)
575 handle.advectionHandle().setPhaseIndex(pIdx);
578 auto getRho = [pIdx] (
const auto& volVars) {
return volVars.density(pIdx); };
579 static const bool enableGravity = getParamFromGroup<bool>(problem().paramGroup(),
"Problem.EnableGravity");
581 localAssembler.assembleGravity(handle.advectionHandle(), iv, getRho);
584 auto getPressure = [pIdx] (
const auto& volVars) {
return volVars.pressure(pIdx); };
585 localAssembler.assembleU(handle.advectionHandle(), iv, getPressure);
590 template<
class InteractionVolume,
class DataHandle>
591 void prepareDiffusionHandles_(InteractionVolume& iv,
595 for (
unsigned int phaseIdx = 0; phaseIdx < ModelTraits::numFluidPhases(); ++phaseIdx)
597 for (
unsigned int compIdx = 0; compIdx < ModelTraits::numFluidComponents(); ++compIdx)
600 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
601 if constexpr (!FluidSystem::isTracerFluidSystem())
602 if (compIdx == FluidSystem::getMainComponent(phaseIdx))
606 handle.diffusionHandle().setPhaseIndex(phaseIdx);
607 handle.diffusionHandle().setComponentIndex(compIdx);
609 using DiffusionType = GetPropType<TypeTag, Properties::MolecularDiffusionType>;
610 using EffDiffModel = GetPropType<TypeTag, Properties::EffectiveDiffusivityModel>;
613 using Traits =
typename InteractionVolume::Traits;
614 using IvLocalAssembler =
typename Traits::template LocalAssembler<Problem, FVElementGeometry, ElementVolumeVariables>;
615 IvLocalAssembler localAssembler(problem(), fvGeometry(), elemVolVars());
618 if (forceUpdateAll || diffusionIsSolDependent)
621 const auto getD = [&](
const auto& volVars)
623 if constexpr (FluidSystem::isTracerFluidSystem())
624 return Deprecated::template effectiveDiffusionCoefficient<EffDiffModel>(volVars, 0, 0, compIdx);
626 return Deprecated::template effectiveDiffusionCoefficient<EffDiffModel>(volVars, phaseIdx, FluidSystem::getMainComponent(phaseIdx), compIdx);
631 const auto& scv = *scvs(fvGeometry()).begin();
632 const auto& scvf = *scvfs(fvGeometry()).begin();
633 const auto& vv = elemVolVars()[scv];
634 const auto D = [&] ()
638 if constexpr (!FluidSystem::isTracerFluidSystem())
639 return max(1e-20, vv.diffusionCoefficient(phaseIdx, FluidSystem::getMainComponent(phaseIdx), compIdx));
641 return max(1e-20, vv.diffusionCoefficient(0, 0, compIdx));
645 localAssembler.assembleMatrices(handle.diffusionHandle(), iv, getD, eps);
649 auto getMassOrMoleFraction = [phaseIdx, compIdx] (
const auto& volVars)
652 volVars.moleFraction(phaseIdx, compIdx);
655 localAssembler.assembleU(handle.diffusionHandle(), iv, getMassOrMoleFraction);
661 template<
class InteractionVolume,
class DataHandle>
662 void prepareHeatConductionHandle_(InteractionVolume& iv, DataHandle& handle,
bool forceUpdateAll)
665 using Traits =
typename InteractionVolume::Traits;
666 using IvLocalAssembler =
typename Traits::template LocalAssembler<Problem, FVElementGeometry, ElementVolumeVariables>;
667 IvLocalAssembler localAssembler(problem(), fvGeometry(), elemVolVars());
670 auto getLambda = [] (
const auto& volVars) {
return volVars.effectiveThermalConductivity(); };
673 if (forceUpdateAll || heatConductionIsSolDependent)
674 localAssembler.assembleMatrices(handle.heatConductionHandle(), iv, getLambda);
677 auto getTemperature = [] (
const auto& volVars) {
return volVars.temperature(); };
678 localAssembler.assembleU(handle.heatConductionHandle(), iv, getTemperature);
681 const Problem* problemPtr_;
682 const Element* elementPtr_;
683 const FVElementGeometry* fvGeometryPtr_;
684 const ElementVolumeVariables* elemVolVarsPtr_;
689 PrimaryInteractionVolume* primaryIv_;
690 SecondaryInteractionVolume* secondaryIv_;
693 PrimaryDataHandle* primaryIvDataHandle_;
694 SecondaryDataHandle* secondaryIvDataHandle_;
700 const PrimaryLocalFaceData* primaryLocalFaceData_;
701 const SecondaryLocalFaceData* secondaryLocalFaceData_;
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
The available discretization methods in Dumux.
The reference frameworks and formulations available for splitting total fluxes into a advective and d...
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
Tensor::field_type computeTpfaTransmissibility(const SubControlVolumeFace &scvf, const SubControlVolume &scv, const Tensor &T, typename SubControlVolume::Traits::Scalar extrusionFactor)
Free function to evaluate the Tpfa transmissibility associated with the flux (in the form of flux = T...
Definition: tpfa/computetransmissibility.hh:47
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:149
Definition: fluxvariablescachefiller.hh:39
PorousMediumFluxVariablesCacheFillerImplementation(const Problem &problem)
The constructor. Sets the problem pointer.
Definition: fluxvariablescachefiller.hh:79
void fill(FluxVariablesCacheContainer &fluxVarsCacheContainer, FluxVariablesCache &scvfFluxVarsCache, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, bool forceUpdateAll=false)
function to fill the flux variables caches
Definition: fluxvariablescachefiller.hh:94
const PrimaryLocalFaceData & primaryIvLocalFaceData() const
returns the currently stored iv-local face data object
Definition: fluxvariablescachefiller.hh:346
void fill(FluxVarsCacheStorage &fluxVarsCacheStorage, FluxVariablesCache &scvfFluxVarsCache, IVDataStorage &ivDataStorage, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, bool forceUpdateAll=false)
function to fill the flux variables caches
Definition: fluxvariablescachefiller.hh:245
const PrimaryInteractionVolume & primaryInteractionVolume() const
returns the stored interaction volume pointer
Definition: fluxvariablescachefiller.hh:330
PorousMediumFluxVariablesCacheFillerImplementation(const Problem &problem)
The constructor. Sets problem pointer.
Definition: fluxvariablescachefiller.hh:229
const PrimaryDataHandle & primaryIvDataHandle() const
returns the stored data handle pointer
Definition: fluxvariablescachefiller.hh:338
const SecondaryDataHandle & secondaryIvDataHandle() const
returns the stored data handle pointer
Definition: fluxvariablescachefiller.hh:342
const SecondaryLocalFaceData & secondaryIvLocalFaceData() const
returns the currently stored iv-local face data object
Definition: fluxvariablescachefiller.hh:350
const SecondaryInteractionVolume & secondaryInteractionVolume() const
returns the stored interaction volume pointer
Definition: fluxvariablescachefiller.hh:334
Declares all properties used in Dumux.
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...