diff --git a/dumux/freeflow/navierstokes/momentum/brinkman.hh b/dumux/freeflow/navierstokes/momentum/brinkman.hh
index 8ec503514c74063bc80c1b280999af9ebef7c190..7cdb3fc020a5b5aa3ec9dda5d7be3c5318658330 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 772874be2ff5548b04e1044407dd01b8b56a1fc3..6a39312e6f6350afa238984990f2d9d565420698 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 bb0997085fcc1d884cf5896a7de56644d8c3c6b0..62a73a382c40d847c19d40dfe32648843698b8ee 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,