NavierStokes base problems changed

extensions/base_problems/NavierStokesPhase_TaylorHood.cc

@@ -19,7 +19,8 @@ NavierStokesPhase_TaylorHood::NavierStokesPhase_TaylorHood(const std::string &na
   Parameters::get(name + "->epsilon", epsilon);
   Parameters::get(name + "->alpha", alpha);
   
-};
+}
+
 
 NavierStokesPhase_TaylorHood::~NavierStokesPhase_TaylorHood()
 {
@@ -27,6 +28,7 @@ NavierStokesPhase_TaylorHood::~NavierStokesPhase_TaylorHood()
     delete fileWriterPhase;
 }
 
+
 void NavierStokesPhase_TaylorHood::initData()
 {
   super::initData();
@@ -52,8 +54,15 @@ void NavierStokesPhase_TaylorHood::transferInitialSolution(AdaptInfo *adaptInfo)
 void NavierStokesPhase_TaylorHood::closeTimestep(AdaptInfo *adaptInfo)
 {
   super::closeTimestep(adaptInfo);
-  
   phaseOld->interpol(phase);
+  writeFiles(adaptInfo, false);
+}
+
+
+void NavierStokesPhase_TaylorHood::writeFiles(AdaptInfo *adaptInfo, bool force)
+{ FUNCNAME("NavierStokesPhase_TaylorHood::closeTimestep()");
+
+  super::writeFiles(adaptInfo, force);
   fileWriterPhase->writeFiles(adaptInfo, false);
 }
 
@@ -76,7 +85,7 @@ void NavierStokesPhase_TaylorHood::fillOperators()
     /// < (1/tau)*u_i^old , psi >
     Operator *opTimeOld = new Operator(prob->getFeSpace(i), prob->getFeSpace(i));
     opTimeOld->addTerm(new Phase_ZOT(phaseOld, density));
-    opTimeOld->setUhOld(prob->getSolution()->getDOFVector(i));
+    opTimeOld->setUhOld(getOldSolution(i));
     prob->addVectorOperator(*opTimeOld, i, getInvTau());
  
     /// < u^old*grad(u_i^old) , psi >
@@ -148,7 +157,7 @@ void NavierStokesPhase_TaylorHood::addLaplaceTerm(int i)
     if (laplaceType == 1) {
       for (unsigned j = 0; j < dow; ++j) {
         Operator *opLaplaceUi1 = new Operator(prob->getFeSpace(i), prob->getFeSpace(j));
-        opLaplaceUi1->addTerm(new MatrixIJPhase_SOT(phase, 1-i, 1-j, viscosity));
+        opLaplaceUi1->addTerm(new MatrixIJPhase_SOT(phase, j, i, viscosity));
         prob->addMatrixOperator(*opLaplaceUi1, i, j, &theta);
 
         opLaplaceUi1->setUhOld(prob->getSolution()->getDOFVector(j));

extensions/base_problems/NavierStokesPhase_TaylorHood.h

@@ -29,6 +29,8 @@ public: // methods
   virtual void transferInitialSolution(AdaptInfo *adaptInfo);  
   virtual void closeTimestep(AdaptInfo *adaptInfo);
   
+  virtual void writeFiles(AdaptInfo *adaptInfo, bool force = false);
+  
   // === getting/setting methods ===
   
   DOFVector<double> *getPhase()
@@ -51,6 +53,11 @@ public: // methods
     return epsilon;
   }
   
+  DOFVector<WorldVector<double> >* getVelocity()
+  {
+    return velocity;
+  }
+  
   void setPhase(DOFVector<double>* phase_)
   {
     phase = phase_;

extensions/base_problems/NavierStokes_TaylorHood.cc

@@ -16,7 +16,10 @@ NavierStokes_TaylorHood::NavierStokes_TaylorHood(const std::string &name_) :
   theta(0.5),
   theta1(0.5),
   minusTheta1(-0.5),
-  velocity(NULL)
+  velocity(NULL),
+  fileWriter(NULL),
+  initialX(NULL),
+  initialY(NULL)
 {  
   // force the homogeniouse dirichlet BC if combination of dirichlet and neumann BC
   // are set and AMDiS can not handle this combination automatically
@@ -50,6 +53,10 @@ NavierStokes_TaylorHood::NavierStokes_TaylorHood(const std::string &name_) :
     Initfile::get("mesh->scale",dimension);
   } else
     dimension.set(1.0);
+  
+  
+  for (size_t i = 0; i < dow; i++)
+    oldSolution[i] = NULL;
 }
 
 
@@ -60,8 +67,18 @@ NavierStokes_TaylorHood::~NavierStokes_TaylorHood()
     delete initialX;
     delete initialY;
   }
-  if (velocity != NULL)
+  
+  if (velocity != NULL) {
     delete velocity;
+    velocity = NULL;
+  }
+  
+  for (size_t i = 0; i < dow; i++) {
+    if (oldSolution[i] != NULL)
+      delete oldSolution[i];
+    oldSolution[i] = NULL;
+  }
+  
   delete fileWriter;
 }
 
@@ -71,6 +88,10 @@ void NavierStokes_TaylorHood::initData()
 
   if (velocity == NULL)
     velocity = new DOFVector<WorldVector<double> >(getFeSpace(0), "velocity");
+  
+  for (size_t i = 0; i < dow; i++)
+    oldSolution[i] = new DOFVector<double>(getFeSpace(i), "old(v_"+Helpers::toString(i)+")");
+  
   fileWriter = new FileVectorWriter(name + "->velocity->output", getFeSpace()->getMesh(), velocity);
 
   super::initData();
@@ -111,7 +132,6 @@ void NavierStokes_TaylorHood::solveInitialProblem(AdaptInfo *adaptInfo)
     prob->getSolution()->getDOFVector(i)->interpol(initialFcts[i]);
     prob->setExactSolutionFct(initialFcts[i], i);
   }
-
 }
 
 
@@ -123,6 +143,9 @@ void NavierStokes_TaylorHood::transferInitialSolution(AdaptInfo *adaptInfo)
   for (int i = 0; i < dow; i++)
     prob->setExactSolution(prob->getSolution()->getDOFVector(i), i);
 
+  for (size_t i = 0; i < dow; i++)
+    oldSolution[i]->copy(*prob->getSolution()->getDOFVector(i));
+  
   fileWriter->writeFiles(adaptInfo, false);
   writeFiles(adaptInfo, false);
 
@@ -211,7 +234,7 @@ void NavierStokes_TaylorHood::addLaplaceTerm(int i)
     if (laplaceType == 1) {
       for (unsigned j = 0; j < dow; ++j) {
         Operator *opLaplaceUi1 = new Operator(getFeSpace(i), getFeSpace(j));
-        opLaplaceUi1->addTerm(new MatrixIJ_SOT(1-i, 1-j, viscosity));
+        opLaplaceUi1->addTerm(new MatrixIJ_SOT(j, i, viscosity));
         prob->addMatrixOperator(*opLaplaceUi1, i, j, &theta, &theta);
 
 	if (abs(minusTheta1) > DBL_TOL) {
@@ -244,6 +267,14 @@ void NavierStokes_TaylorHood::closeTimestep(AdaptInfo *adaptInfo)
 { FUNCNAME("NavierStokes_TaylorHood::closeTimestep()");
 
   calcVelocity();
+  for (size_t i = 0; i < dow; i++)
+    oldSolution[i]->copy(*prob->getSolution()->getDOFVector(i));
+}
+
+
+void NavierStokes_TaylorHood::writeFiles(AdaptInfo *adaptInfo, bool force)
+{ FUNCNAME("NavierStokesPhase_TaylorHood::closeTimestep()");
+
+  super::writeFiles(adaptInfo, force);
   fileWriter->writeFiles(adaptInfo, false);
-  writeFiles(adaptInfo, false);
 }

extensions/base_problems/NavierStokes_TaylorHood.h

@@ -43,6 +43,8 @@ public: // methods
   virtual void transferInitialSolution(AdaptInfo *adaptInfo);  
 
   virtual void closeTimestep(AdaptInfo *adaptInfo);
+  
+  virtual void writeFiles(AdaptInfo *adaptInfo, bool force = false);
 
   // === getting/setting methods ===
 
@@ -51,10 +53,15 @@ public: // methods
 
   DOFVector<WorldVector<double> >* getVelocity()
   {
-    if (velocity == NULL)
-      velocity = new DOFVector<WorldVector<double> >(getFeSpace(0), "velocity");
+    TEST_EXIT_DBG(velocity != NULL)("velocity is NULL!");
     return velocity;
   }
+  
+  DOFVector<double> *getOldSolution(int i)
+  {
+    TEST_EXIT_DBG(oldSolution[i] != NULL)("Für die angegebene Komponente ist keine OldSolution erstellt!");
+    return oldSolution[i];
+  }
 
 // protected: // methods
 
@@ -99,5 +106,6 @@ protected: // variables
   WorldVector<double> dimension;
 
   FileVectorWriter *fileWriter;
+  std::vector<DOFVector<double>*> oldSolution;
 };
 #endif // NAVIER_STOKES_TAYLOR_HOOD_H