From c719d074d4aceb6d136dbf995333dad561ad262d Mon Sep 17 00:00:00 2001
From: Sina Ackermann <sina.ackermann@iws.uni-stuttgart.de>
Date: Fri, 21 Apr 2017 14:02:50 +0200
Subject: [PATCH] [staggeredGrid] Revise Fourier's law
* correct distance computation for gradient
---
.../staggered/freeflow/fourierslaw.hh | 32 ++++++++-----------
1 file changed, 14 insertions(+), 18 deletions(-)
diff --git a/dumux/discretization/staggered/freeflow/fourierslaw.hh b/dumux/discretization/staggered/freeflow/fourierslaw.hh
index 06214e09f3..6a7af1a3df 100644
--- a/dumux/discretization/staggered/freeflow/fourierslaw.hh
+++ b/dumux/discretization/staggered/freeflow/fourierslaw.hh
@@ -51,6 +51,7 @@ class FouriersLawImplementation<TypeTag, DiscretizationMethods::Staggered >
{
using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+ using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
@@ -64,7 +65,9 @@ class FouriersLawImplementation<TypeTag, DiscretizationMethods::Staggered >
static const int dimWorld = GridView::dimensionworld;
enum {
- energyBalanceIdx = Indices::energyBalanceIdx
+ massBalanceIdx = Indices::massBalanceIdx,
+ energyBalanceIdx = Indices::energyBalanceIdx,
+ phaseIdx = Indices::phaseIdx
};
public:
@@ -98,35 +101,28 @@ public:
outsideLambda *= outsideVolVars.extrusionFactor();
// the resulting averaged conductivity tensor
-// const auto lambda = problem.spatialParams().harmonicMean(insideLambda, outsideLambda, scvf.unitOuterNormal());
- const auto lambda = harmonicMean(insideLambda, outsideLambda); // TODO unitOuterNormal?!
+ const auto lambda = harmonicMean(insideLambda, outsideLambda);
const Scalar insideTemp = insideVolVars.temperature();
- Scalar distance(0.0), outsideTemp(insideTemp);
+ const Scalar outsideTemp = scvf.boundary() ? problem.dirichletAtPos(scvf.center())[energyBalanceIdx]
+ : outsideVolVars.temperature();
+ const Scalar distance = scvf.boundary() ? (insideScv.dofPosition() - scvf.ipGlobal()).two_norm()
+ : (outsideScv.dofPosition() - scvf.ipGlobal()).two_norm();
if(scvf.boundary())
{
const auto bcTypes = problem.boundaryTypesAtPos(scvf.center());
- if(bcTypes.isOutflow(energyBalanceIdx))
- return flux; // TODO flux = 0??
- else if(bcTypes.isNeumann(energyBalanceIdx))
- return flux; // TODO: implement neumann
- else
- {
- distance = (insideScv.dofPosition() - scvf.ipGlobal()).two_norm();
- outsideTemp = problem.dirichletAtPos(scvf.center())[energyBalanceIdx];
- }
- }
- else
- {
- distance = (insideScv.dofPosition() - outsideScv.dofPosition()).two_norm();
- outsideTemp = outsideVolVars.temperature();
+ if(bcTypes.isOutflow(energyBalanceIdx) || bcTypes.isNeumann(energyBalanceIdx))
+ return flux;
}
flux[energyBalanceIdx] = -1.0 * (insideTemp - outsideTemp);
flux[energyBalanceIdx] *= lambda / distance;
+ flux *= scvf.area() * sign(scvf.outerNormalScalar());
return flux;
}
+
+
};
} // end namespace
--
GitLab