15#ifndef DUMUX_MP_LINEAR_MATERIAL_HH 
   16#define DUMUX_MP_LINEAR_MATERIAL_HH 
   19#include <dune/common/fvector.hh> 
   32template<
class S, 
int numFlu
idPhases>
 
   34: 
public Adapter<MPLinearMaterial<S, numFluidPhases>, MultiPhasePcKrSw>
 
   38        Params(
const std::array<S, numFluidPhases>& pcMaxSat,
 
   39               const std::array<S, numFluidPhases>& pcMinSat)
 
   40        : pcMaxSat_(pcMaxSat), pcMinSat_(pcMinSat) {}
 
   42        S pcMaxSat(
int phaseIdx)
 const { 
return pcMaxSat_[phaseIdx]; }
 
   43        S pcMinSat(
int phaseIdx)
 const { 
return pcMinSat_[phaseIdx]; }
 
   45        std::array<S, numFluidPhases> pcMaxSat_;
 
   46        std::array<S, numFluidPhases> pcMinSat_;
 
   54    : basicParams_(basicParams)
 
 
   68    template <
class Flu
idState>
 
   71        static_assert(FluidState::numPhases == numFluidPhases, 
"FluidState doesn't match the number of fluid phases!");
 
   72        Dune::FieldVector<typename FluidState::Scalar, numFluidPhases> values;
 
   73        for (
int phaseIdx = 0; phaseIdx < numFluidPhases; ++phaseIdx)
 
   75            const Scalar saturation = state.saturation(phaseIdx);
 
   77                saturation*basicParams_.pcMaxSat(phaseIdx) +
 
   78                (1 - saturation)*basicParams_.pcMinSat(phaseIdx);
 
 
   88    template <
class Flu
idState>
 
   91        static_assert(FluidState::numPhases == numFluidPhases, 
"FluidState doesn't match the number of fluid phases!");
 
   93        Dune::FieldVector<typename FluidState::Scalar, FluidState::numPhases> values;
 
   94        for (
int phaseIdx = 0; phaseIdx < FluidState::numPhases; ++phaseIdx)
 
   95            values[phaseIdx] = clamp(state.saturation(phaseIdx), 0.0, 1.0);
 
 
 
auto relativePermeabilities(const FluidState &state, int wPhaseIdx=0) const
The relative permeability of all phases.
Definition mplinearmaterial.hh:89
S Scalar
Definition mplinearmaterial.hh:51
Params BasicParams
Definition mplinearmaterial.hh:50
auto capillaryPressures(const FluidState &state, int wPhaseIdx=0) const
The linear capillary pressure-saturation curve.
Definition mplinearmaterial.hh:69
MPLinearMaterial(const BasicParams &basicParams)
Definition mplinearmaterial.hh:53
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
Definition brookscorey.hh:23
Adapter to inherit from, allowing the inheriting class to be wrapped by the makeFluidMatrixInteractio...
Definition fluidmatrixinteraction.hh:56