fixed bug in adaptive 2p model after change of initialization order

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@8885 2fb0f335-1f38-0410-981e-8018bf24f1b0

dumux/decoupled/2p/diffusion/fv/fvpressure2padaptive.hh

@@ -106,6 +106,29 @@ public:
     // Function which calculates the flux entry
     void getFlux(EntryType&, const Intersection&, const CellData&, const bool);
 
+    //!Initialize
+    void initialize()
+    {
+        ParentType::initialize();
+
+        if (!compressibility_)
+        {
+            ElementIterator element = problem_.gridView().template begin<0>();
+            FluidState fluidState;
+            fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
+            fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
+            fluidState.setTemperature(problem_.temperature(*element));
+            fluidState.setSaturation(wPhaseIdx, 1.);
+            fluidState.setSaturation(nPhaseIdx, 0.);
+            density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
+            density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
+            viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
+            viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
+        }
+
+        velocity_.initialize();
+    }
+
     /*! \brief Pressure update
      *
      *  \copydetails FVPressure2P::update()
@@ -166,20 +189,10 @@ public:
             DUNE_THROW(Dune::NotImplemented, "Saturation type not supported!");
         }
 
-        if (!compressibility_)
-        {
-            ElementIterator element = problem_.gridView().template begin<0>();
-            FluidState fluidState;
-            fluidState.setPressure(wPhaseIdx, problem_.referencePressure(*element));
-            fluidState.setPressure(nPhaseIdx, problem_.referencePressure(*element));
-            fluidState.setTemperature(problem_.temperature(*element));
-            fluidState.setSaturation(wPhaseIdx, 1.);
-            fluidState.setSaturation(nPhaseIdx, 0.);
-            density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
-            density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
-            viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
-            viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
-        }
+        density_[wPhaseIdx] = 0.;
+        density_[nPhaseIdx] = 0.;
+        viscosity_[wPhaseIdx] = 0.;
+        viscosity_[nPhaseIdx] = 0.;
     }
 
 private:

dumux/decoupled/2p/diffusion/fv/fvvelocity2padaptive.hh

@@ -116,6 +116,17 @@ public:
             DUNE_THROW(Dune::NotImplemented, "Velocity type not supported!");
         }
 
+        density_[wPhaseIdx] = 0.;
+        density_[nPhaseIdx] = 0.;
+        viscosity_[wPhaseIdx] = 0.;
+        viscosity_[nPhaseIdx] = 0.;
+    }
+
+    //!For initialization
+    void initialize()
+    {
+        ParentType::initialize();
+
         if (!compressibility_)
         {
             ElementIterator element = problem_.gridView().template begin<0> ();
@@ -197,6 +208,9 @@ void FVVelocity2PAdaptive<TypeTag>::calculateVelocity(const Intersection& inters
         //get face index
         int isIndexI = intersection.indexInInside();
 
+        Scalar faceArea = intersection.geometry().volume();
+        Scalar faceAreaSum = faceArea;
+
         //get face normal
         const Dune::FieldVector<Scalar, dim>& unitOuterNormal = intersection.centerUnitOuterNormal();
 
@@ -225,6 +239,8 @@ void FVVelocity2PAdaptive<TypeTag>::calculateVelocity(const Intersection& inters
                 {
                     globalIdxK = problem_.variables().index(*neighborPointer2);
                     elementK = neighborPointer2;
+                    faceAreaSum += isItI->geometry().volume();
+
                     break;
                 }
             }
@@ -440,8 +456,8 @@ void FVVelocity2PAdaptive<TypeTag>::calculateVelocity(const Intersection& inters
         cellDataJ.fluxData().setVelocityMarker(isIndexJ);
 
         //times 0.5 because cell face with hanging node is called twice! Do not set marker because it should be called twice!
-        velocityW *= 0.5;
-        velocityNW *= 0.5;
+        velocityW *= faceArea/faceAreaSum;
+        velocityNW *= faceArea/faceAreaSum;
         cellData.fluxData().addVelocity(wPhaseIdx, isIndexI, velocityW);
         cellData.fluxData().addVelocity(nPhaseIdx, isIndexI, velocityNW);
     }

dumux/decoupled/common/fv/fvvelocity.hh

@@ -69,6 +69,7 @@ template<class TypeTag, class Velocity> class FVVelocity
 
 public:
 
+    //!Initialize velocity implementation
     void initialize()
     {
         velocity_.initialize();