12#ifndef DUMUX_DISCRETIZATION_UPWINDSCHEME_HH 
   13#define DUMUX_DISCRETIZATION_UPWINDSCHEME_HH 
   21template<
class Gr
idGeometry, 
class DiscretizationMethod>
 
   29template<
class Gr
idGeometry>
 
   35template<
class ElemVolVars, 
class SubControlVolumeFace, 
class UpwindTermFunction, 
class Scalar>
 
   37                              const SubControlVolumeFace& scvf,
 
   38                              const UpwindTermFunction& upwindTerm,
 
   39                              Scalar flux, 
int phaseIdx)
 
   42    static const Scalar upwindWeight = 
getParamFromGroup<Scalar>(elemVolVars.gridVolVars().problem().paramGroup(), 
"Flux.UpwindWeight");
 
   44    const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
   45    const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
 
   49        return (upwindWeight*upwindTerm(outsideVolVars)
 
   50                + (1.0 - upwindWeight)*upwindTerm(insideVolVars));
 
   52        return (upwindWeight*upwindTerm(insideVolVars)
 
   53                + (1.0 - upwindWeight)*upwindTerm(outsideVolVars));
 
 
   59template<
class Gr
idGeometry, 
class DM>
 
   64    template<
class FluxVariables, 
class UpwindTermFunction, 
class Scalar>
 
   65    static Scalar 
apply(
const FluxVariables& fluxVars,
 
   66                        const UpwindTermFunction& upwindTerm,
 
   67                        Scalar flux, 
int phaseIdx)
 
   69        return apply(fluxVars.elemVolVars(), fluxVars.scvFace(), upwindTerm, flux, phaseIdx);
 
 
   73    template<
class ElemVolVars, 
class SubControlVolumeFace, 
class UpwindTermFunction, 
class Scalar>
 
   74    static Scalar 
apply(
const ElemVolVars& elemVolVars,
 
   75                        const SubControlVolumeFace& scvf,
 
   76                        const UpwindTermFunction& upwindTerm,
 
   77                        Scalar flux, 
int phaseIdx)
 
   79        return flux*
multiplier(elemVolVars, scvf, upwindTerm, flux, phaseIdx);
 
 
   83    template<
class ElemVolVars, 
class SubControlVolumeFace, 
class UpwindTermFunction, 
class Scalar>
 
   85                             const SubControlVolumeFace& scvf,
 
   86                             const UpwindTermFunction& upwindTerm,
 
   87                             Scalar flux, 
int phaseIdx)
 
 
 
   94template<
class Gr
idGeometry>
 
   97    using GridView = 
typename GridGeometry::GridView;
 
   98    static constexpr int dim = GridView::dimension;
 
   99    static constexpr int dimWorld = GridView::dimensionworld;
 
  103    template<
class ElemVolVars, 
class SubControlVolumeFace, 
class UpwindTermFunction, 
class Scalar>
 
  105                             const SubControlVolumeFace& scvf,
 
  106                             const UpwindTermFunction& upwindTerm,
 
  107                             Scalar flux, 
int phaseIdx)
 
  109        static_assert(dim == dimWorld, 
"Multiplier cannot be computed on surface grids using cell-centered scheme!");
 
 
  114    template<
class ElemVolVars, 
class SubControlVolumeFace, 
class UpwindTermFunction, 
class Scalar>
 
  115    static Scalar 
apply(
const ElemVolVars& elemVolVars,
 
  116                        const SubControlVolumeFace& scvf,
 
  117                        const UpwindTermFunction& upwindTerm,
 
  118                        Scalar flux, 
int phaseIdx)
 
  120        static_assert(dim == dimWorld, 
"This upwind scheme cannot be used on surface grids using cell-centered schemes!");
 
  121        return flux*
multiplier(elemVolVars, scvf, upwindTerm, flux, phaseIdx);
 
 
  125    template<
class FluxVariables, 
class UpwindTermFunction, 
class Scalar>
 
  126    static Scalar 
apply(
const FluxVariables& fluxVars,
 
  127                        const UpwindTermFunction& upwindTerm,
 
  128                        Scalar flux, 
int phaseIdx)
 
  130        const auto& scvf = fluxVars.scvFace();
 
  131        const auto& elemVolVars = fluxVars.elemVolVars();
 
  134        if constexpr (dim == dimWorld)
 
  135            return apply(elemVolVars, scvf, upwindTerm, flux, phaseIdx);
 
  138        if (scvf.numOutsideScvs() == 1)
 
  143        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  145        Scalar branchingPointUpwindTerm = 0.0;
 
  146        Scalar sumUpwindFluxes = 0.0;
 
  151            return upwindTerm(insideVolVars)*flux;
 
  153            sumUpwindFluxes += flux;
 
  155        for (
unsigned int i = 0; i < scvf.numOutsideScvs(); ++i)
 
  158            const auto& fvGeometry = fluxVars.fvGeometry();
 
  159            const auto outsideScvIdx = scvf.outsideScvIdx(i);
 
  160            const auto outsideElement = fvGeometry.gridGeometry().element(outsideScvIdx);
 
  161            const auto& flippedScvf = fvGeometry.flipScvf(scvf.index(), i);
 
  163            using AdvectionType = 
typename FluxVariables::AdvectionType;
 
  164            const auto outsideFlux = AdvectionType::flux(fluxVars.problem(),
 
  170                                                         fluxVars.elemFluxVarsCache());
 
  172            if (!signbit(outsideFlux))
 
  173                branchingPointUpwindTerm += upwindTerm(elemVolVars[outsideScvIdx])*outsideFlux;
 
  175                sumUpwindFluxes += outsideFlux;
 
  179        if (sumUpwindFluxes != 0.0)
 
  180            branchingPointUpwindTerm /= -sumUpwindFluxes;
 
  182            branchingPointUpwindTerm = 0.0;
 
  188            return flux*branchingPointUpwindTerm;
 
  190            return flux*upwindTerm(insideVolVars);
 
 
 
  195template<
class Gr
idGeometry>
 
static Scalar apply(const FluxVariables &fluxVars, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
Definition flux/upwindscheme.hh:126
static Scalar apply(const ElemVolVars &elemVolVars, const SubControlVolumeFace &scvf, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
applies a simple upwind scheme to the precalculated advective flux across the given scvf
Definition flux/upwindscheme.hh:115
static Scalar multiplier(const ElemVolVars &elemVolVars, const SubControlVolumeFace &scvf, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
returns the upwind factor which is multiplied to the advective flux across the given scvf
Definition flux/upwindscheme.hh:104
static Scalar multiplier(const ElemVolVars &elemVolVars, const SubControlVolumeFace &scvf, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
returns the upwind factor which is multiplied to the advective flux across the given scvf
Definition flux/upwindscheme.hh:84
static Scalar apply(const ElemVolVars &elemVolVars, const SubControlVolumeFace &scvf, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
applies a simple upwind scheme to the precalculated advective flux across the given scvf
Definition flux/upwindscheme.hh:74
static Scalar apply(const FluxVariables &fluxVars, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
applies a simple upwind scheme to the precalculated advective flux
Definition flux/upwindscheme.hh:65
Forward declaration of the upwind scheme implementation.
Definition flux/upwindscheme.hh:22
UpwindSchemeImpl< GridGeometry, typename GridGeometry::DiscretizationMethod > UpwindScheme
The upwind scheme used for the advective fluxes. This depends on the chosen discretization method.
Definition flux/upwindscheme.hh:30
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
The available discretization methods in Dumux.
Definition cellcentered/mpfa/elementvolumevariables.hh:24
Distance implementation details.
Definition cvfelocalresidual.hh:25
Scalar upwindSchemeMultiplier(const ElemVolVars &elemVolVars, const SubControlVolumeFace &scvf, const UpwindTermFunction &upwindTerm, Scalar flux, int phaseIdx)
returns the upwind factor which is multiplied to the advective flux across the given scvf
Definition flux/upwindscheme.hh:36
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.