diff --git a/dumux/discretization/cellcentered/tpfa/computetransmissibility.hh b/dumux/discretization/cellcentered/tpfa/computetransmissibility.hh
new file mode 100644
index 0000000000000000000000000000000000000000..9c73844cf03eba1f8902e53227d18a46921e75d1
--- /dev/null
+++ b/dumux/discretization/cellcentered/tpfa/computetransmissibility.hh
@@ -0,0 +1,86 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ * See the file COPYING for full copying permissions. *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 2 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ *****************************************************************************/
+/*!
+ * \file
+ * \brief This file contains free functions to evaluate the transmissibilities
+ * associated with flux evaluations across sub-control volume faces
+ * in the context of the cell-centered TPFA scheme.
+ */
+#ifndef DUMUX_DISCRETIZATION_CC_TPFA_COMPUTE_TRANSMISSIBILITY_HH
+#define DUMUX_DISCRETIZATION_CC_TPFA_COMPUTE_TRANSMISSIBILITY_HH
+
+namespace Dumux
+{
+
+/*!
+ * \ingroup Tpfa
+ *
+ * \brief Free function to evaluate the Tpfa transmissibility associated
+ * with the flux (in the form of flux = T*gradU) across a
+ * sub-control volume face stemming from a given sub-control
+ * volume with corresponding tensor T.
+ *
+ * \param scvf The sub-control volume face
+ * \param scv The neighboring sub-control volume
+ * \param K The tensor living in the neighboring scv
+ * \param extrusionFactor The extrusion factor of the scv
+ */
+template< class SubControlVolumeFace, class SubControlVolume, class FieldScalar, int dimWorld >
+FieldScalar computeTpfaTransmissibility(const SubControlVolumeFace& scvf,
+ const SubControlVolume& scv,
+ const Dune::FieldMatrix<FieldScalar, dimWorld, dimWorld>& T,
+ typename SubControlVolume::Traits::Scalar extrusionFactor)
+{
+ using GlobalPosition = typename SubControlVolumeFace::Traits::GlobalPosition;
+ GlobalPosition Knormal;
+ T.mv(scvf.unitOuterNormal(), Knormal);
+
+ auto distanceVector = scvf.ipGlobal();
+ distanceVector -= scv.center();
+ distanceVector /= distanceVector.two_norm2();
+
+ return (Knormal*distanceVector) * extrusionFactor;
+}
+
+/*!
+ * \ingroup Tpfa
+ *
+ * \brief Specialization of the above function for scalar T.
+ *
+ * \param scvf The sub-control volume face
+ * \param scv The neighboring sub-control volume
+ * \param t The scalar quantity living in the neighboring scv
+ * \param extrusionFactor The extrusion factor of the scv
+ */
+template< class SubControlVolumeFace, class SubControlVolume, class FieldScalar >
+FieldScalar computeTpfaTransmissibility(const SubControlVolumeFace& scvf,
+ const SubControlVolume &scv,
+ FieldScalar t,
+ typename SubControlVolumeFace::Traits::Scalar extrusionFactor)
+{
+ auto distanceVector = scvf.ipGlobal();
+ distanceVector -= scv.center();
+ distanceVector /= distanceVector.two_norm2();
+
+ return t * extrusionFactor * (distanceVector * scvf.unitOuterNormal());
+}
+
+} // end namespace Dumux
+
+#endif
diff --git a/dumux/discretization/cellcentered/tpfa/darcyslaw.hh b/dumux/discretization/cellcentered/tpfa/darcyslaw.hh
index 3a00ac3028abdf840c8b3122bada60e09028ccbc..197ca72ce125257534b25e79d4e705a4af0b1bd9 100644
--- a/dumux/discretization/cellcentered/tpfa/darcyslaw.hh
+++ b/dumux/discretization/cellcentered/tpfa/darcyslaw.hh
@@ -25,14 +25,10 @@
#ifndef DUMUX_DISCRETIZATION_CC_TPFA_DARCYS_LAW_HH
#define DUMUX_DISCRETIZATION_CC_TPFA_DARCYS_LAW_HH
-#include <memory>
-
-#include <dune/common/float_cmp.hh>
-
-#include <dumux/common/math.hh>
#include <dumux/common/parameters.hh>
#include <dumux/discretization/methods.hh>
+#include <dumux/discretization/cellcentered/tpfa/computetransmissibility.hh>
namespace Dumux
{
@@ -261,8 +257,8 @@ public:
return volVars.permeability();
};
- const Scalar ti = calculateOmega_(scvf, getPermeability(insideVolVars, scvf.ipGlobal()),
- insideScv, insideVolVars.extrusionFactor());
+ const Scalar ti = computeTpfaTransmissibility(scvf, insideScv, getPermeability(insideVolVars, scvf.ipGlobal()),
+ insideVolVars.extrusionFactor());
// for the boundary (dirichlet) or at branching points we only need ti
if (scvf.boundary() || scvf.numOutsideScvs() > 1)
@@ -283,14 +279,14 @@ public:
{
// normal grids
if (dim == dimWorld)
- return -1.0*calculateOmega_(scvf, getPermeability(outsideVolVars, scvf.ipGlobal()),
- outsideScv, outsideVolVars.extrusionFactor());
+ return -1.0*computeTpfaTransmissibility(scvf, outsideScv, getPermeability(outsideVolVars, scvf.ipGlobal()),
+ outsideVolVars.extrusionFactor());
// embedded surface and network grids
//(the outside normal vector might differ from the inside normal vector)
else
- return calculateOmega_(fvGeometry.flipScvf(scvf.index()), getPermeability(outsideVolVars, scvf.ipGlobal()),
- outsideScv, outsideVolVars.extrusionFactor());
+ return computeTpfaTransmissibility(fvGeometry.flipScvf(scvf.index()), outsideScv, getPermeability(outsideVolVars, scvf.ipGlobal()),
+ outsideVolVars.extrusionFactor());
}();
@@ -303,42 +299,6 @@ public:
return tij;
}
-
-private:
- //! compute the transmissibility ti, overload for tensor permeabilites
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- const DimWorldMatrix &K,
- const SubControlVolume &scv,
- Scalar extrusionFactor)
- {
- GlobalPosition Knormal;
- K.mv(scvf.unitOuterNormal(), Knormal);
-
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = Knormal * distanceVector;
- omega *= extrusionFactor;
-
- return omega;
- }
-
- //! compute the transmissibility ti, overload for scalar permeabilites
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- const Scalar K,
- const SubControlVolume &scv,
- Scalar extrusionFactor)
- {
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = K * (distanceVector * scvf.unitOuterNormal());
- omega *= extrusionFactor;
-
- return omega;
- }
};
} // end namespace Dumux
diff --git a/dumux/discretization/cellcentered/tpfa/fickslaw.hh b/dumux/discretization/cellcentered/tpfa/fickslaw.hh
index 8f8696714773f4e5b18d7f0ec81e733812c9e229..36667762470ce5b9685cff50cbd156f20e51b945 100644
--- a/dumux/discretization/cellcentered/tpfa/fickslaw.hh
+++ b/dumux/discretization/cellcentered/tpfa/fickslaw.hh
@@ -24,13 +24,11 @@
#ifndef DUMUX_DISCRETIZATION_CC_TPFA_FICKS_LAW_HH
#define DUMUX_DISCRETIZATION_CC_TPFA_FICKS_LAW_HH
-#include <dune/common/float_cmp.hh>
-#include <dumux/common/math.hh>
#include <dumux/common/parameters.hh>
#include <dumux/discretization/methods.hh>
-#include <dumux/discretization/fluxvariablescaching.hh>
+#include <dumux/discretization/cellcentered/tpfa/computetransmissibility.hh>
namespace Dumux
{
@@ -180,10 +178,7 @@ public:
using EffDiffModel = typename GET_PROP_TYPE(TypeTag, EffectiveDiffusivityModel);
auto insideD = insideVolVars.diffusionCoefficient(phaseIdx, compIdx);
insideD = EffDiffModel::effectiveDiffusivity(insideVolVars.porosity(), insideVolVars.saturation(phaseIdx), insideD);
- const Scalar ti = calculateOmega_(scvf,
- insideD,
- insideScv,
- insideVolVars.extrusionFactor());
+ const Scalar ti = computeTpfaTransmissibility(scvf, insideScv, insideD, insideVolVars.extrusionFactor());
// for the boundary (dirichlet) or at branching points we only need ti
if (scvf.boundary() || scvf.numOutsideScvs() > 1)
@@ -203,15 +198,12 @@ public:
Scalar tj;
if (dim == dimWorld)
// assume the normal vector from outside is anti parallel so we save flipping a vector
- tj = -1.0*calculateOmega_(scvf,
- outsideD,
- outsideScv,
- outsideVolVars.extrusionFactor());
+ tj = -1.0*computeTpfaTransmissibility(scvf, outsideD, outsideScv, outsideVolVars.extrusionFactor());
else
- tj = calculateOmega_(fvGeometry.flipScvf(scvf.index()),
- outsideD,
- outsideScv,
- outsideVolVars.extrusionFactor());
+ tj = computeTpfaTransmissibility(fvGeometry.flipScvf(scvf.index()),
+ outsideD,
+ outsideScv,
+ outsideVolVars.extrusionFactor());
// check if we are dividing by zero!
if (ti*tj <= 0.0)
@@ -264,41 +256,6 @@ public:
}
return rho/(scvf.numOutsideScvs()+1);
}
-
-private:
-
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- const DimWorldMatrix &D,
- const SubControlVolume &scv,
- const Scalar extrusionFactor)
- {
- GlobalPosition Dnormal;
- D.mv(scvf.unitOuterNormal(), Dnormal);
-
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = Dnormal * distanceVector;
- omega *= extrusionFactor;
-
- return omega;
- }
-
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- const Scalar D,
- const SubControlVolume &scv,
- const Scalar extrusionFactor)
- {
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = D * (distanceVector * scvf.unitOuterNormal());
- omega *= extrusionFactor;
-
- return omega;
- }
};
} // end namespace
diff --git a/dumux/discretization/cellcentered/tpfa/fourierslaw.hh b/dumux/discretization/cellcentered/tpfa/fourierslaw.hh
index a7a0e3f8af26fd3315a1ad3dc0176370ffb4dee8..62833620c2bb7c7f6ca41951d6c1de7339da5798 100644
--- a/dumux/discretization/cellcentered/tpfa/fourierslaw.hh
+++ b/dumux/discretization/cellcentered/tpfa/fourierslaw.hh
@@ -24,12 +24,10 @@
#ifndef DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
#define DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
-#include <dune/common/float_cmp.hh>
-#include <dumux/common/math.hh>
#include <dumux/common/parameters.hh>
#include <dumux/discretization/methods.hh>
-#include <dumux/discretization/fluxvariablescaching.hh>
+#include <dumux/discretization/cellcentered/tpfa/computetransmissibility.hh>
namespace Dumux
{
@@ -146,7 +144,7 @@ public:
const auto& insideVolVars = elemVolVars[insideScvIdx];
auto insideLambda = ThermalConductivityModel::effectiveThermalConductivity(insideVolVars, problem.spatialParams(), element, fvGeometry, insideScv);
- Scalar ti = calculateOmega_(scvf, insideLambda, insideScv, insideVolVars.extrusionFactor());
+ Scalar ti = computeTpfaTransmissibility(scvf, insideScv, insideLambda, insideVolVars.extrusionFactor());
// for the boundary (dirichlet) or at branching points we only need ti
if (scvf.boundary() || scvf.numOutsideScvs() > 1)
@@ -169,9 +167,9 @@ public:
Scalar tj;
if (dim == dimWorld)
// assume the normal vector from outside is anti parallel so we save flipping a vector
- tj = -1.0*calculateOmega_(scvf, outsideLambda, outsideScv, outsideVolVars.extrusionFactor());
+ tj = -1.0*computeTpfaTransmissibility(scvf, outsideScv, outsideLambda, outsideVolVars.extrusionFactor());
else
- tj = calculateOmega_(fvGeometry.flipScvf(scvf.index()), outsideLambda, outsideScv, outsideVolVars.extrusionFactor());
+ tj = computeTpfaTransmissibility(fvGeometry.flipScvf(scvf.index()), outsideScv, outsideLambda, outsideVolVars.extrusionFactor());
// check for division by zero!
if (ti*tj <= 0.0)
@@ -211,35 +209,6 @@ private:
}
return sumTempTi/sumTi;
}
-
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- const DimWorldMatrix& lambda,
- const SubControlVolume& scv,
- Scalar extrusionFactor)
- {
- GlobalPosition lambdaNormal;
- lambda.mv(scvf.unitOuterNormal(), lambdaNormal);
-
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = lambdaNormal * distanceVector;
- return omega*extrusionFactor;
- }
-
- static Scalar calculateOmega_(const SubControlVolumeFace& scvf,
- Scalar lambda,
- const SubControlVolume &scv,
- Scalar extrusionFactor)
- {
- auto distanceVector = scvf.ipGlobal();
- distanceVector -= scv.center();
- distanceVector /= distanceVector.two_norm2();
-
- Scalar omega = lambda * (distanceVector * scvf.unitOuterNormal());
- return omega*extrusionFactor;
- }
};
} // end namespace Dumux