12#ifndef DUMUX_NAVIERSTOKES_SCALAR_BOUNDARY_FLUXHELPER_HH 
   13#define DUMUX_NAVIERSTOKES_SCALAR_BOUNDARY_FLUXHELPER_HH 
   15#include <dune/common/float_cmp.hh> 
   16#include <dune/common/std/type_traits.hh> 
   26template <
class Indices>
 
   27using NonisothermalDetector = 
decltype(std::declval<Indices>().energyEqIdx);
 
   29template<
class Indices>
 
   30static constexpr bool isNonIsothermal()
 
   31{ 
return Dune::Std::is_detected<NonisothermalDetector, Indices>::value; }
 
   40template<
class AdvectiveFlux>
 
   46    template<
class VolumeVariables, 
class SubControlVolumeFace, 
class Scalar, 
class UpwindTerm>
 
   48                                            const VolumeVariables& outsideVolVars,
 
   49                                            const SubControlVolumeFace& scvf,
 
   50                                            const Scalar volumeFlux,
 
   51                                            const Scalar upwindWeight,
 
   52                                            UpwindTerm upwindTerm)
 
   55        const bool insideIsUpstream = !signbit(volumeFlux);
 
   57        const auto& upstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
 
   58        const auto& downstreamVolVars = insideIsUpstream ? outsideVolVars : insideVolVars;
 
   60        return (upwindWeight * upwindTerm(upstreamVolVars) +
 
   61               (1.0 - upwindWeight) * upwindTerm(downstreamVolVars))
 
 
   65    template<
class Indices, 
class NumEqVector, 
class UpwindFunction>
 
   67                                              const UpwindFunction& upwind)
 
   70        AdvectiveFlux::addAdvectiveFlux(flux, upwind);
 
   73        if constexpr (Detail::isNonIsothermal<Indices>())
 
   75            auto upwindTerm = [](
const auto& volVars) { 
return volVars.density()*volVars.enthalpy(); };
 
   76            flux[Indices::energyEqIdx] = upwind(upwindTerm);
 
 
   84    template<
class Problem, 
class Element, 
class FVElementGeometry, 
class ElementVolumeVariables>
 
   86                                  const Element& element,
 
   87                                  const FVElementGeometry& fvGeometry,
 
   88                                  const typename FVElementGeometry::SubControlVolumeFace& scvf,
 
   89                                  const ElementVolumeVariables& elemVolVars,
 
   90                                  typename ElementVolumeVariables::VolumeVariables::PrimaryVariables&& outsideBoundaryPriVars,
 
   91                                  const typename ElementVolumeVariables::VolumeVariables::PrimaryVariables::value_type upwindWeight = 1.0)
 
   93        using VolumeVariables = 
typename ElementVolumeVariables::VolumeVariables;
 
   94        using PrimaryVariables = 
typename VolumeVariables::PrimaryVariables;
 
   97        const auto velocity = problem.faceVelocity(element,fvGeometry, scvf);
 
   98        const auto volumeFlux = velocity * scvf.unitOuterNormal();
 
  100        const bool insideIsUpstream = !signbit(volumeFlux);
 
  101        const VolumeVariables& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  102        const VolumeVariables& outsideVolVars = [&]()
 
  105            if (insideIsUpstream && Dune::FloatCmp::eq(upwindWeight, 1.0, 1e-6))
 
  106                return insideVolVars;
 
  110                VolumeVariables boundaryVolVars;
 
  114                                       fvGeometry.scv(scvf.insideScvIdx()));
 
  115                return boundaryVolVars;
 
  119        auto upwindFuntion = [&](
const auto& upwindTerm)
 
 
  134    template<
class Problem, 
class Element, 
class FVElementGeometry, 
class ElementVolumeVariables>
 
  136                                  const Element& element,
 
  137                                  const FVElementGeometry& fvGeometry,
 
  138                                  const typename FVElementGeometry::SubControlVolumeFace& scvf,
 
  139                                  const ElementVolumeVariables& elemVolVars)
 
  141        using VolumeVariables = 
typename ElementVolumeVariables::VolumeVariables;
 
  142        using NumEqVector = 
typename VolumeVariables::PrimaryVariables;
 
  144        const auto velocity = problem.faceVelocity(element,fvGeometry, scvf);
 
  145        const auto volumeFlux = velocity * scvf.unitOuterNormal();
 
  147        const bool insideIsUpstream = !signbit(volumeFlux);
 
  148        const VolumeVariables& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  150        if constexpr (VolumeVariables::FluidSystem::isCompressible(0)  || NumEqVector::size() > 1)
 
  152            static const bool verbose = 
getParamFromGroup<bool>(problem.paramGroup(), 
"Flux.EnableOutflowReversalWarning", 
true);
 
  154            if (verbose && !insideIsUpstream && abs(volumeFlux) > 1e-10)
 
  156                std::cout << 
"velo " << velocity << 
", flux " << volumeFlux << std::endl;
 
  157                std::cout << 
"\n ********** WARNING ********** \n\n" 
  158                "Outflow condition set at " << scvf.center() << 
" might be invalid due to flow reversal. " 
  160                "scalarOutflowFlux(problem, element, fvGeometry, scvf, elemVolVars, outsideBoundaryPriVars, upwindWeight) \n" 
  161                "instead where you can specify primary variables for inflow situations.\n" 
  162                "\n ***************************** \n" << std::endl;
 
  166        auto upwindFuntion = [&](
const auto& upwindTerm)
 
 
 
Element solution classes and factory functions.
typename NumEqVectorTraits< PrimaryVariables >::type NumEqVector
A vector with the same size as numbers of equations This is the default implementation and has to be ...
Definition numeqvector.hh:34
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethods::cctpfa||GridGeometry::discMethod==DiscretizationMethods::ccmpfa, CCElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for cell-centered schemes.
Definition cellcentered/elementsolution.hh:101
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
Define some often used mathematical functions.
Distance implementation details.
Definition cvfelocalresidual.hh:25
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Navier Stokes scalar boundary flux helper.
Definition scalarfluxhelper.hh:42
static Scalar advectiveScalarUpwindFlux(const VolumeVariables &insideVolVars, const VolumeVariables &outsideVolVars, const SubControlVolumeFace &scvf, const Scalar volumeFlux, const Scalar upwindWeight, UpwindTerm upwindTerm)
Return the area-specific, weighted advective flux of a scalar quantity.
Definition scalarfluxhelper.hh:47
static auto scalarOutflowFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars)
Return the area-specific outflow fluxes for all scalar balance equations. This should only be used if...
Definition scalarfluxhelper.hh:135
static auto scalarOutflowFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, typename ElementVolumeVariables::VolumeVariables::PrimaryVariables &&outsideBoundaryPriVars, const typename ElementVolumeVariables::VolumeVariables::PrimaryVariables::value_type upwindWeight=1.0)
Return the area-specific outflow fluxes for all scalar balance equations. The values specified in out...
Definition scalarfluxhelper.hh:85
static void addModelSpecificAdvectiveFlux(NumEqVector &flux, const UpwindFunction &upwind)
Definition scalarfluxhelper.hh:66