From 3fe96c43f3216012c8228d058ad8b72015067e55 Mon Sep 17 00:00:00 2001
From: Kilian <kilian.weishaupt@iws.uni-stuttgart.de>
Date: Tue, 2 Jul 2019 13:39:30 +0200
Subject: [PATCH] [staggered][upwindfluxvars] Use new helper class and new BJ
function
* not done for getParallelVelocityFromOtherBoundary -> will raise
deprecation warning
---
.../navierstokes/staggered/fluxvariables.hh | 8 ++--
.../staggered/staggeredupwindfluxvariables.hh | 43 +++++++++++--------
2 files changed, 29 insertions(+), 22 deletions(-)
diff --git a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
index ba330448da..09a4a30611 100644
--- a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
+++ b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
@@ -387,8 +387,10 @@ public:
// If none of the above boundary conditions apply for the given sub face, proceed to calculate the tangential momentum flux.
if (problem.enableInertiaTerms())
lateralFlux += computeAdvectivePartOfLateralMomentumFlux_(problem, fvGeometry, element,
- scvf, elemVolVars, faceVars,
- gridFluxVarsCache, lateralFaceBoundaryTypes, localSubFaceIdx);
+ scvf, elemVolVars, faceVars,
+ gridFluxVarsCache,
+ currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
+ localSubFaceIdx);
lateralFlux += computeDiffusivePartOfLateralMomentumFlux_(problem, fvGeometry, element,
scvf, elemVolVars, faceVars,
@@ -487,7 +489,7 @@ private:
const Scalar transportingVelocity = faceVars.velocityLateralInside(localSubFaceIdx);
return StaggeredUpwindFluxVariables<TypeTag, upwindSchemeOrder>::computeUpwindedLateralMomentum(problem, fvGeometry, element, scvf, elemVolVars, faceVars,
- gridFluxVarsCache, localSubFaceIdx, lateralFaceBoundaryTypes)
+ gridFluxVarsCache, localSubFaceIdx, currentScvfBoundaryTypes, lateralFaceBoundaryTypes)
* transportingVelocity * lateralFace.directionSign() * lateralFace.area() * 0.5 * extrusionFactor_(elemVolVars, lateralFace);
}
diff --git a/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh b/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh
index 29baf9ae1b..d8485ec10c 100644
--- a/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh
+++ b/dumux/freeflow/navierstokes/staggered/staggeredupwindfluxvariables.hh
@@ -33,6 +33,7 @@
#include <dumux/discretization/method.hh>
#include <dumux/freeflow/staggeredupwindmethods.hh>
+#include "velocitygradients.hh"
namespace Dumux {
@@ -69,6 +70,7 @@ class StaggeredUpwindFluxVariables
using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
using FacePrimaryVariables = GetPropType<TypeTag, Properties::FacePrimaryVariables>;
using BoundaryTypes = GetPropType<TypeTag, Properties::BoundaryTypes>;
+ using VelocityGradients = StaggeredVelocityGradients<Scalar, FVGridGeometry, BoundaryTypes, Indices>;
using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
static constexpr bool useHigherOrder = upwindSchemeOrder > 1;
@@ -121,6 +123,7 @@ public:
const FaceVariables& faceVars,
const GridFluxVariablesCache& gridFluxVarsCache,
const int localSubFaceIdx,
+ const Dune::Std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
const Dune::Std::optional<BoundaryTypes>& lateralFaceBoundaryTypes)
{
const auto eIdx = scvf.insideScvIdx();
@@ -138,15 +141,15 @@ public:
if (canHigherOrder)
{
const auto parallelUpwindingMomenta = getLateralUpwindingMomenta_(problem, fvGeometry, element, scvf, elemVolVars, faceVars,
- transportingVelocity, localSubFaceIdx, lateralFaceBoundaryTypes,
- std::integral_constant<bool, useHigherOrder>{});
+ transportingVelocity, localSubFaceIdx, currentScvfBoundaryTypes,
+ lateralFaceBoundaryTypes, std::integral_constant<bool, useHigherOrder>{});
return doLateralMomentumUpwinding_(fvGeometry, scvf, parallelUpwindingMomenta, transportingVelocity, localSubFaceIdx, gridFluxVarsCache);
}
else
{
const auto parallelUpwindingMomenta = getLateralUpwindingMomenta_(problem, fvGeometry, element, scvf, elemVolVars, faceVars,
- transportingVelocity, localSubFaceIdx, lateralFaceBoundaryTypes,
- std::integral_constant<bool, false>{});
+ transportingVelocity, localSubFaceIdx, currentScvfBoundaryTypes,
+ lateralFaceBoundaryTypes, std::integral_constant<bool, false>{});
return doLateralMomentumUpwinding_(fvGeometry, scvf, parallelUpwindingMomenta, transportingVelocity, localSubFaceIdx, gridFluxVarsCache);
}
}
@@ -186,11 +189,7 @@ private:
{
// Depending on selfIsUpstream I have to check if I have a forward or a backward neighbor to retrieve
const bool selfIsUpstream = ownScvf.directionSign() != sign(transportingVelocity);
-
- if (selfIsUpstream)
- return ownScvf.hasForwardNeighbor(0);
- else
- return ownScvf.hasBackwardNeighbor(0);
+ return selfIsUpstream ? ownScvf.hasForwardNeighbor(0) : ownScvf.hasBackwardNeighbor(0);
}
/*!
@@ -369,6 +368,7 @@ private:
const FaceVariables& faceVars,
const Scalar transportingVelocity,
const int localSubFaceIdx,
+ const Dune::Std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
const Dune::Std::optional<BoundaryTypes>& lateralFaceBoundaryTypes,
std::true_type)
{
@@ -383,8 +383,8 @@ private:
if (!ownScvf.hasParallelNeighbor(localSubFaceIdx, 0))
{
const Scalar boundaryMomentum = getParallelVelocityFromBoundary_(problem, element, fvGeometry, ownScvf,
- faceVars.velocitySelf(), localSubFaceIdx,
- lateralFaceBoundaryTypes) * insideVolVars.density();
+ faceVars, currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
+ localSubFaceIdx) * insideVolVars.density();
return std::array<Scalar, 3>{boundaryMomentum, boundaryMomentum, boundaryMomentum};
}
@@ -402,8 +402,8 @@ private:
momenta[0] = faceVars.velocityParallel(localSubFaceIdx, 0) * insideVolVars.density();
else
momenta[0] = getParallelVelocityFromBoundary_(problem, element, fvGeometry, ownScvf,
- faceVars.velocitySelf(), localSubFaceIdx,
- lateralFaceBoundaryTypes) * insideVolVars.density();
+ faceVars, currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
+ localSubFaceIdx) * insideVolVars.density();
// The local index of the faces that is opposite to localSubFaceIdx
const int oppositeSubFaceIdx = localSubFaceIdx % 2 ? localSubFaceIdx - 1 : localSubFaceIdx + 1;
@@ -450,6 +450,7 @@ private:
const FaceVariables& faceVars,
const Scalar transportingVelocity,
const int localSubFaceIdx,
+ const Dune::Std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
const Dune::Std::optional<BoundaryTypes>& lateralFaceBoundaryTypes,
std::false_type)
{
@@ -463,8 +464,8 @@ private:
const Scalar momentumParallel = ownScvf.hasParallelNeighbor(localSubFaceIdx, 0)
? faceVars.velocityParallel(localSubFaceIdx, 0) * outsideVolVars.density()
: (getParallelVelocityFromBoundary_(problem, element, fvGeometry, ownScvf,
- faceVars.velocitySelf(), localSubFaceIdx,
- lateralFaceBoundaryTypes)
+ faceVars, currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
+ localSubFaceIdx)
* insideVolVars.density());
// If the lateral face lies on a boundary, we assume that the parallel velocity on the boundary is actually known,
@@ -572,13 +573,15 @@ private:
const Element& element,
const FVElementGeometry& fvGeometry,
const SubControlVolumeFace& scvf,
- const Scalar velocitySelf,
- const int localSubFaceIdx,
- const Dune::Std::optional<BoundaryTypes>& lateralFaceBoundaryTypes)
+ const FaceVariables& faceVars,
+ const Dune::Std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
+ const Dune::Std::optional<BoundaryTypes>& lateralFaceBoundaryTypes,
+ const int localSubFaceIdx)
{
// Find out what boundary type is set on the lateral face
const bool lateralFaceHasDirichletPressure = lateralFaceBoundaryTypes && lateralFaceBoundaryTypes->isDirichlet(Indices::pressureIdx);
const bool lateralFaceHasBJS = lateralFaceBoundaryTypes && lateralFaceBoundaryTypes->isBJS(Indices::velocity(scvf.directionIndex()));
+ const Scalar velocitySelf = faceVars.velocitySelf();
// If there is a Dirichlet condition for the pressure we assume zero gradient for the velocity,
// so the velocity at the boundary equal to that on the scvf.
@@ -589,9 +592,11 @@ private:
{
const auto eIdx = scvf.insideScvIdx();
const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx);
+ const Scalar velocityGrad_ji = VelocityGradients::velocityGradJI(problem, element, fvGeometry, scvf, faceVars,
+ currentScvfBoundaryTypes, lateralFaceBoundaryTypes, localSubFaceIdx);
const SubControlVolume& scv = fvGeometry.scv(scvf.insideScvIdx());
- return problem.bjsVelocity(element, scv, lateralFace, velocitySelf);
+ return problem.beaversJosephVelocity(element, scv, lateralFace, velocitySelf, velocityGrad_ji);
}
else
{
--
GitLab