47    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   48    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   53    using Element = 
typename GridView::template Codim<0>::Entity;
 
   58    using ComponentFluxVector = Dune::FieldVector<Scalar, numComponents>;
 
   59    using ReducedComponentVector = Dune::FieldVector<Scalar, numComponents-1>;
 
   60    using ReducedComponentMatrix = Dune::FieldMatrix<Scalar, numComponents-1, numComponents-1>;
 
   69    { 
return referenceSystem; }
 
 
   76    static ComponentFluxVector 
flux(
const Problem& problem,
 
   77                                    const Element& element,
 
   78                                    const FVElementGeometry& fvGeometry,
 
   79                                    const ElementVolumeVariables& elemVolVars,
 
   80                                    const SubControlVolumeFace& scvf,
 
   82                                    const ElementFluxVariablesCache& elemFluxVarsCache)
 
   86        ComponentFluxVector componentFlux(0.0);
 
   87        ReducedComponentMatrix reducedDiffusionMatrix(0.0);
 
   88        ReducedComponentVector reducedFlux(0.0);
 
   89        ComponentFluxVector moleFrac(0.0);
 
   90        ReducedComponentVector normalX(0.0);
 
   93        const auto& fluxVarsCache = elemFluxVarsCache[scvf];
 
   94        const auto& shapeValues = fluxVarsCache.shapeValues();
 
   97        Scalar avgMolarMass(0.0);
 
   98        for (
auto&& scv : scvs(fvGeometry))
 
  100            const auto& volVars = elemVolVars[scv];
 
  103            rho +=  volVars.density(phaseIdx)*shapeValues[scv.indexInElement()][0];
 
  105            avgMolarMass += volVars.averageMolarMass(phaseIdx)*shapeValues[scv.indexInElement()][0];
 
  107            for (
int compIdx = 0; compIdx < numComponents; compIdx++)
 
  109              moleFrac[compIdx] += volVars.moleFraction(phaseIdx, compIdx)*shapeValues[scv.indexInElement()][0];
 
  113        reducedDiffusionMatrix = setupMSMatrix_(problem, element, fvGeometry, elemVolVars, scvf, phaseIdx, moleFrac, avgMolarMass);
 
  115        for (
int compIdx = 0; compIdx < numComponents-1; compIdx++)
 
  117            Dune::FieldVector<Scalar, GridView::dimensionworld> gradX(0.0);
 
  118            for (
auto&& scv : scvs(fvGeometry))
 
  120                const auto& volVars = elemVolVars[scv];
 
  123                gradX.axpy(volVars.moleFraction(phaseIdx, compIdx), fluxVarsCache.gradN(scv.indexInElement()));
 
  126           normalX[compIdx] = gradX *scvf.unitOuterNormal();
 
  128         reducedDiffusionMatrix.solve(reducedFlux,normalX);
 
  129         reducedFlux *= -1.0*rho*Extrusion::area(fvGeometry, scvf);
 
  131        for (
int compIdx = 0; compIdx < numComponents-1; compIdx++)
 
  133            componentFlux[compIdx] = reducedFlux[compIdx];
 
  134            componentFlux[numComponents-1] -= reducedFlux[compIdx];
 
  136        return componentFlux ;
 
 
  141    static ReducedComponentMatrix setupMSMatrix_(
const Problem& problem,
 
  142                                                 const Element& element,
 
  143                                                 const FVElementGeometry& fvGeometry,
 
  144                                                 const ElementVolumeVariables& elemVolVars,
 
  145                                                 const SubControlVolumeFace& scvf,
 
  147                                                 const ComponentFluxVector moleFrac,
 
  148                                                 const Scalar avgMolarMass)
 
  150        ReducedComponentMatrix reducedDiffusionMatrix(0.0);
 
  152        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  153        const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
 
  156        if(Dune::FloatCmp::eq<Scalar>(insideVolVars.saturation(phaseIdx), 0) || Dune::FloatCmp::eq<Scalar>(outsideVolVars.saturation(phaseIdx), 0))
 
  157            return reducedDiffusionMatrix;
 
  159        for (
int compIIdx = 0; compIIdx < numComponents-1; compIIdx++)
 
  162            const auto xi = moleFrac[compIIdx];
 
  163            const auto Mn = FluidSystem::molarMass(numComponents-1);
 
  165            auto tinInside = insideVolVars.effectiveDiffusionCoefficient(phaseIdx, compIIdx, numComponents-1);
 
  166            auto tinOutside = outsideVolVars.effectiveDiffusionCoefficient(phaseIdx, compIIdx, numComponents-1);
 
  169            tinInside *= insideVolVars.extrusionFactor();
 
  170            tinOutside *= outsideVolVars.extrusionFactor();
 
  173            const auto tin = 
faceTensorAverage(tinInside, tinOutside, scvf.unitOuterNormal());
 
  176            reducedDiffusionMatrix[compIIdx][compIIdx] += xi*avgMolarMass/(tin*Mn);
 
  179            for (
int compJIdx = 0; compJIdx < numComponents; compJIdx++)
 
  182                if (compIIdx == compJIdx)
 
  186                const auto xj = moleFrac[compJIdx];
 
  187                const auto Mi = FluidSystem::molarMass(compIIdx);
 
  188                const auto Mj = FluidSystem::molarMass(compJIdx);
 
  191                auto tijInside = insideVolVars.effectiveDiffusionCoefficient(phaseIdx, compIIdx, compJIdx);
 
  192                auto tijOutside = outsideVolVars.effectiveDiffusionCoefficient(phaseIdx, compIIdx, compJIdx);
 
  195                tijInside *= insideVolVars.extrusionFactor();
 
  196                tijOutside *= outsideVolVars.extrusionFactor();
 
  199                const auto tij = 
faceTensorAverage(tijInside, tijOutside, scvf.unitOuterNormal());
 
  201                reducedDiffusionMatrix[compIIdx][compIIdx] += xj*avgMolarMass/(tij*Mi);
 
  202                if (compJIdx < numComponents-1)
 
  203                    reducedDiffusionMatrix[compIIdx][compJIdx] +=xi*(avgMolarMass/(tin*Mn) - avgMolarMass/(tij*Mj));
 
  206        return reducedDiffusionMatrix;
 
 
static ComponentFluxVector flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const ElementFluxVariablesCache &elemFluxVarsCache)
Returns the diffusive fluxes of all components within a fluid phase across the given sub-control volu...
Definition box/maxwellstefanslaw.hh:76
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Scalar faceTensorAverage(const Scalar T1, const Scalar T2, const Dune::FieldVector< Scalar, dim > &normal)
Average of a discontinuous scalar field at discontinuity interface (for compatibility reasons with th...
Definition facetensoraverage.hh:29
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166