[freeflow][fluxVars] Add generic advectiveFlux function

* similiar to PM flow, takes lambda to get transported quantity

dumux/freeflow/navierstokes/staggered/fluxvariables.hh

@@ -91,6 +91,30 @@ class NavierStokesFluxVariablesImpl<TypeTag, DiscretizationMethods::Staggered>
 
 public:
 
+    template<class UpwindTerm>
+    static Scalar advectiveFluxForCellCenter(const ElementVolumeVariables& elemVolVars,
+                                             const ElementFaceVariables& elemFaceVars,
+                                             const SubControlVolumeFace &scvf,
+                                             UpwindTerm upwindTerm,
+                                             bool isOutflow = false)
+    {
+        const Scalar velocity = elemFaceVars[scvf].velocitySelf();
+        const bool insideIsUpstream = sign(scvf.outerNormalScalar()) == sign(velocity);
+        static const Scalar upWindWeight = getParamFromGroup<Scalar>(GET_PROP_VALUE(TypeTag, ModelParameterGroup), "Implicit.UpwindWeight");
+
+        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
+        const auto& outsideVolVars = isOutflow ? insideVolVars : elemVolVars[scvf.outsideScvIdx()];
+
+        const auto& upstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
+        const auto& downstreamVolVars = insideIsUpstream ? outsideVolVars : insideVolVars;
+
+        const Scalar flux = (upWindWeight * upwindTerm(upstreamVolVars) +
+                            (1.0 - upWindWeight) * upwindTerm(downstreamVolVars))
+                            * velocity * scvf.area() * sign(scvf.outerNormalScalar());
+
+        return flux;
+    }
+
     CellCenterPrimaryVariables computeFluxForCellCenter(const Problem& problem,
                                                         const Element &element,
                                                         const FVElementGeometry& fvGeometry,
@@ -99,27 +123,14 @@ public:
                                                         const SubControlVolumeFace &scvf,
                                                         const FluxVariablesCache& fluxVarsCache)
     {
-        const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
-        const auto& insideVolVars = elemVolVars[insideScv];
+        auto upwindTerm = [](const auto& volVars) { return volVars.density(); };
 
-        // if we are on an inflow/outflow boundary, use the volVars of the element itself
-        const auto& outsideVolVars = scvf.boundary() ?  insideVolVars : elemVolVars[scvf.outsideScvIdx()];
+        const Scalar flux = advectiveFluxForCellCenter(elemVolVars, elemFaceVars, scvf, upwindTerm, false);
 
-        CellCenterPrimaryVariables flux(0.0);
-        const Scalar velocity = elemFaceVars[scvf].velocitySelf();
+        CellCenterPrimaryVariables result(0.0);
+        result[massBalanceIdx] = flux;
 
-        const bool insideIsUpstream = sign(scvf.outerNormalScalar()) == sign(velocity) ? true : false;
-        const auto& upstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
-        const auto& downstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
-
-        static const Scalar upWindWeight = getParamFromGroup<Scalar>(GET_PROP_VALUE(TypeTag, ModelParameterGroup), "Implicit.UpwindWeight");
-
-        flux = (upWindWeight * upstreamVolVars.density() +
-               (1.0 - upWindWeight) * downstreamVolVars.density()) * velocity;
-
-        flux *= scvf.area() * sign(scvf.outerNormalScalar());
-
-        return flux;
+        return result;
     }
 
     void computeCellCenterToCellCenterStencil(Stencil& stencil,