From 440b5db1064c2d9764f509bb05101102d2a16ba6 Mon Sep 17 00:00:00 2001
From: Markus Wolff <markus.wolff@twt-gmbh.de>
Date: Thu, 8 Dec 2011 10:45:41 +0000
Subject: [PATCH] corrected sources and neumann boundaries

   - added division by phase densities



git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@6966 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
 .../2p/diffusion/mimetic/mimeticgroundwater.hh         | 10 ++++++++--
 1 file changed, 8 insertions(+), 2 deletions(-)

diff --git a/dumux/decoupled/2p/diffusion/mimetic/mimeticgroundwater.hh b/dumux/decoupled/2p/diffusion/mimetic/mimeticgroundwater.hh
index 0757ff7aba..e089be5cde 100644
--- a/dumux/decoupled/2p/diffusion/mimetic/mimeticgroundwater.hh
+++ b/dumux/decoupled/2p/diffusion/mimetic/mimeticgroundwater.hh
@@ -332,9 +332,11 @@ public:
         Pi.umv(Wc, F);
         //      std::cout << "Pi = \dim" << Pi << "c = " << c << ", F = " << F << std::endl;
 
+        Scalar densityW = problem_.variables().densityWetting(globalIdx);
+        Scalar densityNW = problem_.variables().densityNonwetting(globalIdx);
         PrimaryVariables sourceVec(0.0);
         problem_.source(sourceVec, element);
-        qmean = volume*(sourceVec[wetting] + sourceVec[nonwetting]);
+        qmean = volume*(sourceVec[wetting]/densityW + sourceVec[nonwetting]/densityNW);
     }
 
 private:
@@ -402,8 +404,12 @@ private:
 
                 if (bcType.isNeumann(pressEqIdx))
                 {
+                    int globalIdx = problem_.variables().index(element);
+                    Scalar densityW = problem_.variables().densityWetting(globalIdx);
+                    Scalar densityNW = problem_.variables().densityNonwetting(globalIdx);
+
                     problem_.neumann(boundValues, *it);
-                    Scalar J = (boundValues[wetting]+boundValues[nonwetting]);
+                    Scalar J = (boundValues[wetting]/densityW + boundValues[nonwetting]/densityNW);
                     this->b[faceIndex] -= J*it->geometry().volume();
                 }
                 else if (bcType.isDirichlet(pressEqIdx))
-- 
GitLab