From 9951e9b28610b59675a3af3f4140e7c2b38ef0d0 Mon Sep 17 00:00:00 2001
From: Markus Wolff <markus.wolff@twt-gmbh.de>
Date: Wed, 30 May 2012 10:01:43 +0000
Subject: [PATCH] some new functionality of gridadapt class, changed call of
 updateMaterialLaw function in base problem

   - GridAdapt: new function wasAdapted() -> returns true or false after
     call of adaptGrid(), new function init() -> calls init() function
     of grid adaption indicator (all indicators need an init() function!), new parameter adaptationInterval -> allows to adapt only at certain time steps,
   - ImpetProblem: updateMaterialLaw() is now called in the
     preTimeStepFunction() per default

   reviewed by Benjamin



git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@8403 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
 dumux/decoupled/common/gridadapt.hh           | 42 ++++++++++++++++---
 dumux/decoupled/common/gridadaptproperties.hh |  4 ++
 dumux/decoupled/common/impetproblem.hh        |  9 ++--
 3 files changed, 45 insertions(+), 10 deletions(-)

diff --git a/dumux/decoupled/common/gridadapt.hh b/dumux/decoupled/common/gridadapt.hh
index 7f16365f77..12a83ed3dc 100644
--- a/dumux/decoupled/common/gridadapt.hh
+++ b/dumux/decoupled/common/gridadapt.hh
@@ -66,16 +66,22 @@ public:
      * @param problem The problem
      */
     GridAdapt (Problem& problem)
-    : problem_(problem), adaptionIndicator_(problem)
+    : problem_(problem), adaptionIndicator_(problem), marked_(0), coarsened_(0)
     {
         levelMin_ = GET_PARAM(TypeTag, int, MinLevel);
         levelMax_ = GET_PARAM(TypeTag, int, MaxLevel);
+        adaptationInterval_ = GET_PARAM(TypeTag, int, AdaptionInterval);
 
         if (levelMin_ < 0)
             Dune::dgrave <<  __FILE__<< ":" <<__LINE__
             << " :  Dune cannot coarsen to gridlevels smaller 0! "<< std::endl;
     }
 
+    void init()
+    {
+        adaptionIndicator_.init();
+    }
+
     /*!
      * @brief Standard method to adapt the grid
      *
@@ -94,6 +100,10 @@ public:
         // reset internal counter for marked elements
         marked_ = coarsened_ = 0;
 
+        // check for adaption interval: Adapt only at certain time step indices
+        if (problem_.timeManager().timeStepIndex() % adaptationInterval_ != 0)
+            return;
+
         /**** 1) determine refining parameter if standard is used ***/
         // if not, the indicatorVector and refinement Bounds have to
         // specified by the problem through setIndicator()
@@ -133,9 +143,6 @@ public:
         // delete markers in grid
         problem_.grid().postAdapt();
 
-        // adapt secondary variables
-        problem_.pressureModel().updateMaterialLaws();
-
         // write out new grid
 //        Dune::VTKWriter<LeafGridView> vtkwriter(leafView);
 //        vtkwriter.write("latestgrid",Dune::VTKOptions::binaryappended);
@@ -211,6 +218,11 @@ public:
         return marked_;
     }
 
+    bool wasAdapted()
+    {
+        return (marked_ != 0 || coarsened_ != 0);
+    }
+
     /*!
      * Sets minimum and maximum refinement levels
      *
@@ -245,6 +257,16 @@ public:
         return levelMin_;
     }
 
+    AdaptionIndicator& adaptionIndicator()
+    {
+        return adaptionIndicator_;
+    }
+
+    AdaptionIndicator& adaptionIndicator() const
+    {
+        return adaptionIndicator_;
+    }
+
 private:
     /*!
      * @brief Method ensuring the refinement ratio of 2:1
@@ -336,11 +358,13 @@ private:
     Problem& problem_;
     AdaptionIndicator adaptionIndicator_;
 
+    int marked_;
+    int coarsened_;
+
     int levelMin_;
     int levelMax_;
 
-    int marked_;
-    int coarsened_;
+    int adaptationInterval_;
 };
 
 /*!
@@ -359,8 +383,14 @@ class GridAdapt<TypeTag, false>
     typedef typename SolutionTypes::ScalarSolution ScalarSolutionType;
 
 public:
+    void init()
+    {};
     void adaptGrid()
     {};
+    bool wasAdapted()
+    {
+        return false;
+    }
     void setLevels(int, int)
     {};
     void setTolerance(int, int)
diff --git a/dumux/decoupled/common/gridadaptproperties.hh b/dumux/decoupled/common/gridadaptproperties.hh
index 08af9109cd..17def39553 100644
--- a/dumux/decoupled/common/gridadaptproperties.hh
+++ b/dumux/decoupled/common/gridadaptproperties.hh
@@ -63,6 +63,9 @@ NEW_PROP_TAG(RefineThreshold);
 //! Tolerance for coarsening
 NEW_PROP_TAG(CoarsenThreshold);
 
+//! Time step interval for adaption
+NEW_PROP_TAG(AdaptionInterval);
+
 //no adaptive grid
 SET_BOOL_PROP(GridAdaptTypeTag, AdaptiveGrid, false);
 
@@ -73,6 +76,7 @@ SET_SCALAR_PROP(GridAdaptTypeTag, RefineTolerance, 0.05);
 SET_SCALAR_PROP(GridAdaptTypeTag, CoarsenTolerance, 0.001);
 SET_SCALAR_PROP(GridAdaptTypeTag, RefineThreshold, 0.0);
 SET_SCALAR_PROP(GridAdaptTypeTag, CoarsenThreshold, 0.0);
+SET_INT_PROP(GridAdaptTypeTag, AdaptionInterval, 1);
 
 } // namespace Properties
 } // namespace Dumux
diff --git a/dumux/decoupled/common/impetproblem.hh b/dumux/decoupled/common/impetproblem.hh
index e6999b5bed..5b4383d02c 100644
--- a/dumux/decoupled/common/impetproblem.hh
+++ b/dumux/decoupled/common/impetproblem.hh
@@ -340,6 +340,8 @@ public:
         // set the initial condition of the model
         variables_.initialize();
         model().initialize();
+        if (adaptiveGrid)
+            gridAdapt().init();
     }
 
     /*!
@@ -352,6 +354,7 @@ public:
         // if it is not used, this method does nothing.
         if (adaptiveGrid)
             this->gridAdapt().adaptGrid();
+        asImp_().pressureModel().updateMaterialLaws();
     }
 
     /*!
@@ -417,9 +420,7 @@ public:
      * current solution to disk.
      */
     void postTimeStep()
-    {
-        asImp_().pressureModel().updateMaterialLaws();
-    };
+    {};
 
     /*!
      * \brief Called by the time manager after everything which can be
@@ -479,7 +480,7 @@ public:
      */
     void setOutputTimeInterval(const Scalar timeInterval)
     {
-        outputTimeInterval_ = timeInterval;
+        outputTimeInterval_ = (timeInterval > 0.0) ? timeInterval : 1e100;
         timeManager().startNextEpisode(outputTimeInterval_);
     }
 
-- 
GitLab