12#ifndef DUMUX_FREEFLOW_NAVIERSTOKES_SLIPCONDITION_HH 
   13#define DUMUX_FREEFLOW_NAVIERSTOKES_SLIPCONDITION_HH 
   15#include <dune/common/fvector.hh> 
   27    static std::string 
name() { 
return "Beavers-Joseph"; }
 
 
   35    static std::string 
name() { 
return "Beavers-Joseph-Saffman"; }
 
 
 
   58template<
class DiscretizationMethod, 
class SlipCondition>
 
   69template<
class SlipCondition>
 
   84    template<
class Problem, 
class FVElementGeometry, 
class ElementVolumeVariables, 
class Scalar>
 
   86                           const FVElementGeometry& fvGeometry,
 
   87                           const typename FVElementGeometry::SubControlVolumeFace& scvf,
 
   88                           const ElementVolumeVariables& elemVolVars,
 
   89                           Scalar tangentialVelocityDeriv)
 
   91        assert(scvf.isLateral());
 
   92        assert(scvf.boundary());
 
   94        static constexpr int dimWorld = FVElementGeometry::GridGeometry::GridView::dimensionworld;
 
   95        using Vector = Dune::FieldVector<Scalar, dimWorld>;
 
   97        Vector porousMediumVelocity(0.0);
 
  100            porousMediumVelocity = problem.porousMediumVelocity(fvGeometry, scvf);
 
  102            DUNE_THROW(Dune::NotImplemented, 
"Fcstaggered currently only implements BJ or BJS slip conditions");
 
  104        const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
  108        tangent[scv.dofAxis()] = 1.0;
 
  112        const Scalar betaBJ = problem.betaBJ(fvGeometry, scvf, tangent);
 
  113        const Scalar distanceNormalToBoundary = (scvf.ipGlobal() - scv.dofPosition()).two_norm();
 
  115        static const bool onlyNormalGradient = 
getParamFromGroup<bool>(problem.paramGroup(), 
"FreeFlow.EnableUnsymmetrizedVelocityGradientForBeaversJoseph", 
false);
 
  116        if (onlyNormalGradient)
 
  117            tangentialVelocityDeriv = 0.0;
 
  119        const Scalar scalarSlipVelocity = (tangentialVelocityDeriv*distanceNormalToBoundary
 
  120            + porousMediumVelocity * tangent * betaBJ * distanceNormalToBoundary
 
  121            + elemVolVars[scv].velocity()) / (betaBJ*distanceNormalToBoundary + 1.0);
 
  123        return scalarSlipVelocity*tangent;
 
 
 
static constexpr SlipCondition slipCondition
Definition slipcondition.hh:73
static auto velocity(const Problem &problem, const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, Scalar tangentialVelocityDeriv)
Returns the slip velocity at a porous boundary based on the Beavers-Joseph(-Saffman) condition.
Definition slipcondition.hh:85
Navier Stokes slip velocity policy.
Definition slipcondition.hh:59
constexpr BJS bjs
Tag for the Beavers-Joseph-Saffman slip condition.
Definition slipcondition.hh:48
constexpr BJ bj
Tag for the Beavers-Joseph slip condition.
Definition slipcondition.hh:42
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 slipcondition.hh:20
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Tag for the Beavers-Joseph slip condition.
Definition slipcondition.hh:26
static std::string name()
Definition slipcondition.hh:27
Tag for the Beavers-Joseph-Saffman slip condition.
Definition slipcondition.hh:34
static std::string name()
Definition slipcondition.hh:35
Helper class to create (named and comparable) tagged types Tags any given type. The tagged type is eq...
Definition tag.hh:30
Helper class to create (named and comparable) tagged types.