From 68cb2878c1815370e884d29686f058a6da83cd29 Mon Sep 17 00:00:00 2001
From: Timo Koch <timokoch@uio.no>
Date: Sat, 9 Dec 2023 18:11:39 +0100
Subject: [PATCH] [ff][brinkman] Add viscosity scaling
---
dumux/freeflow/navierstokes/momentum/brinkman.hh | 9 +++++----
dumux/freeflow/navierstokes/momentum/problem.hh | 14 ++++++++++++++
.../freeflow/couplingmanager_staggered.hh | 14 +++++++++++++-
3 files changed, 32 insertions(+), 5 deletions(-)
diff --git a/dumux/freeflow/navierstokes/momentum/brinkman.hh b/dumux/freeflow/navierstokes/momentum/brinkman.hh
index 8ec503514c..7cdb3fc020 100644
--- a/dumux/freeflow/navierstokes/momentum/brinkman.hh
+++ b/dumux/freeflow/navierstokes/momentum/brinkman.hh
@@ -33,7 +33,7 @@ namespace Dumux {
* The Navier-Stokes model can be extended to a Darcy-Brinkman model by adding
* the term:
* \f[
- + \epsilon_B \mathbf{K}^{-1} \mathbf{v}
+ + \epsilon_B \mu \mathbf{K}^{-1} \mathbf{v}
* \f]
* to the momentum balance. This can be achieved with the helper function \ref addBrinkmanTerm.
* The function relies on the spatial parameters class being based on <code>BrinkmanSpatialParams</code>
@@ -55,26 +55,27 @@ void addBrinkmanTerm(
{
const auto brinkmanEpsilon = problem.spatialParams().brinkmanEpsilon(element, fvGeometry, scv);
const auto invK = problem.spatialParams().inversePermeability(element, fvGeometry, scv);
+ const auto mu = problem.effectiveViscosity(element, fvGeometry, scv);
using DM = typename FVElementGeometry::GridGeometry::DiscretizationMethod;
if constexpr (DiscretizationMethods::isCVFE<DM>)
{
const auto velocity = elemVolVars[scv].velocity();
- source -= brinkmanEpsilon * mv(invK, velocity); // eps K^-1 velocity
+ source -= mu * brinkmanEpsilon * mv(invK, velocity);
}
else if constexpr (DM{} == DiscretizationMethods::fcstaggered)
{
if constexpr (Dune::IsNumber<std::decay_t<decltype(invK)>>::value)
{
const auto velocity = elemVolVars[scv].velocity();
- source -= brinkmanEpsilon * invK * velocity;
+ source -= mu * brinkmanEpsilon * invK * velocity;
}
else
{
// permeability is tensor-valued, use velocity vector reconstruction
const auto getVelocitySCV = [&](const auto& scv){ return elemVolVars[scv].velocity(); };
const auto velocity = StaggeredVelocityReconstruction::faceVelocity(scv, fvGeometry, getVelocitySCV);
- source -= brinkmanEpsilon * mv(invK, velocity)[scv.dofAxis()]; // eps K^-1 velocity
+ source -= mu * brinkmanEpsilon * mv(invK, velocity)[scv.dofAxis()];
}
}
else
diff --git a/dumux/freeflow/navierstokes/momentum/problem.hh b/dumux/freeflow/navierstokes/momentum/problem.hh
index 772874be2f..6a39312e6f 100644
--- a/dumux/freeflow/navierstokes/momentum/problem.hh
+++ b/dumux/freeflow/navierstokes/momentum/problem.hh
@@ -317,6 +317,20 @@ public:
return asImp_().effectiveViscosityAtPos(scvf.ipGlobal());
}
+ /*!
+ * \brief Returns the effective dynamic viscosity at a given sub control volume.
+ * \note Overload this if a fixed viscosity shall be prescribed.
+ */
+ Scalar effectiveViscosity(const Element& element,
+ const FVElementGeometry& fvGeometry,
+ const SubControlVolume& scv) const
+ {
+ if constexpr (isCoupled_)
+ return couplingManager_->effectiveViscosity(element, fvGeometry, scv);
+ else
+ return asImp_().effectiveViscosityAtPos(scv.dofPosition());
+ }
+
/*!
* \brief Returns the effective dynamic viscosity at a given position.
*/
diff --git a/dumux/multidomain/freeflow/couplingmanager_staggered.hh b/dumux/multidomain/freeflow/couplingmanager_staggered.hh
index bb0997085f..62a73a382c 100644
--- a/dumux/multidomain/freeflow/couplingmanager_staggered.hh
+++ b/dumux/multidomain/freeflow/couplingmanager_staggered.hh
@@ -324,7 +324,19 @@ public:
return mu(momentumCouplingContext_()[0].curElemVolVars);
}
- /*!
+ /*!
+ * \brief Returns the pressure at a given sub control volume.
+ */
+ Scalar effectiveViscosity(const Element<freeFlowMomentumIndex>& element,
+ const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
+ const SubControlVolume<freeFlowMomentumIndex>& scv) const
+ {
+ bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex>(), element, fvGeometry.elementIndex());
+ const auto& insideMassScv = momentumCouplingContext_()[0].fvGeometry.scv(scv.elementIndex());
+ return momentumCouplingContext_()[0].curElemVolVars[insideMassScv].viscosity();
+ }
+
+ /*!
* \brief Returns the velocity at a given sub control volume face.
*/
VelocityVector faceVelocity(const Element<freeFlowMassIndex>& element,
--
GitLab