diff --git a/AMDiS/src/AdaptInfo.h b/AMDiS/src/AdaptInfo.h
index d4be7b40ec64c6bfa3424c084985ce1a99a279cc..0bf7becadf0e1a81ad5a52b5b939be38b2af8ed2 100644
--- a/AMDiS/src/AdaptInfo.h
+++ b/AMDiS/src/AdaptInfo.h
@@ -213,6 +213,8 @@ namespace AMDiS {
       timestepNumber = 0;
       solverIterations = 0;
       solverResidual = 0.0;
+
+      GET_PARAMETER(0, name + "->timestep", "%f", &timestep);
     };
 
     /** \brief
diff --git a/AMDiS/src/Operator.cc b/AMDiS/src/Operator.cc
index ab73a1f8a3915ee10d48c933f47413d1d77b4a1d..af4800640c06a4dc70e0289130cecb33e9483406 100644
--- a/AMDiS/src/Operator.cc
+++ b/AMDiS/src/Operator.cc
@@ -407,8 +407,6 @@ namespace AMDiS {
     coordsAtQPs = subAssembler->getCoordsAtQPs(elInfo, quad);
   }
  
-  // bis hierhin erstmal
-
   void VecAndGradAtQP_ZOT::initElement(const ElInfo* elInfo, 
 				       SubAssembler* subAssembler,
 				       Quadrature *quad) 
@@ -2157,4 +2155,34 @@ namespace AMDiS {
     }
   }
 
+  void VecGrad_FOT::initElement(const ElInfo* elInfo, 
+				SubAssembler* subAssembler,
+				Quadrature *quad) 
+  {
+    vecAtQPs = subAssembler->getVectorAtQPs(vec1, elInfo, quad);
+    gradAtQPs = subAssembler->getGradientsAtQPs(vec2, elInfo, quad);
+  }
+
+  void VecGrad_FOT::getLb(const ElInfo *elInfo, int numPoints, VectorOfFixVecs<DimVec<double> >& Lb) const {
+    const DimVec<WorldVector<double> > &Lambda = elInfo->getGrdLambda();
+    for (int iq = 0; iq < numPoints; iq++) {
+      lb(Lambda, (*vecFct)(vecAtQPs[iq], gradAtQPs[iq]), Lb[iq], 1.0);
+    }
+  }
+  
+  void VecGrad_FOT::eval(int numPoints,
+			 const double *uhAtQP,
+			 const WorldVector<double> *grdUhAtQP,
+			 const WorldMatrix<double> *D2UhAtQP,
+			 double *result,
+			 double factor) const
+  {
+    if (grdUhAtQP) {
+      for (int iq = 0; iq < numPoints; iq++) {
+	WorldVector<double> b = (*vecFct)(vecAtQPs[iq], gradAtQPs[iq]);
+	result[iq] += b * grdUhAtQP[iq] * factor;
+      }
+    }
+  }
+
 }
diff --git a/AMDiS/src/Operator.h b/AMDiS/src/Operator.h
index 24c70d9cfff51b711e82b0ec32e7ffcdb888daf7..c99b5a0863a1233b685a6d0bdb2589775f02f7fa 100644
--- a/AMDiS/src/Operator.h
+++ b/AMDiS/src/Operator.h
@@ -1900,6 +1900,67 @@ namespace AMDiS {
     AbstractFunction<WorldVector<double>, WorldVector<double> > *g;
   };
 
+
+  /**
+   * \ingroup Assembler
+   *
+   * \brief
+   * First order term: \f$ b(v(\vec{x}), \nabla w(\vec{x})) \cdot \nabla u(\vec{x})\f$.
+   */
+  class VecGrad_FOT : public FirstOrderTerm
+  {
+  public:
+    /** \brief
+     * Constructor.
+     */
+    VecGrad_FOT(DOFVectorBase<double> *dv1, DOFVectorBase<double> *dv2,
+		BinaryAbstractFunction<WorldVector<double>, double,
+		WorldVector<double> > *vecFct_)
+    : FirstOrderTerm(vecFct_->getDegree()), vec1(dv1), vec2(dv2), vecFct(vecFct_) 
+    {};
+
+    /** \brief
+     * Implementation of \ref OperatorTerm::initElement().
+     */
+    void initElement(const ElInfo* elInfo, SubAssembler* subAssembler,
+		     Quadrature *quad = NULL);
+    
+    /** \brief
+     * Implements FirstOrderTerm::getLb().
+     */
+    void getLb(const ElInfo *elInfo, int qPoint, 
+	       VectorOfFixVecs<DimVec<double> >& Lb) const; 
+    
+    /** \brief
+     * Implements FirstOrderTerm::eval().
+     */
+    void eval(int numPoints,
+	      const double *uhAtQP,
+	      const WorldVector<double> *grdUhAtQP,
+	      const WorldMatrix<double> *D2UhAtQP,
+	      double *result,
+	      double factor) const;
+    
+  protected:
+    /** \brief
+     * DOFVector to be evaluated at quadrature points.
+     */
+    DOFVectorBase<double>* vec1;
+    DOFVectorBase<double>* vec2;
+    
+    /** \brief
+     * Vector v at quadrature points.
+     */
+    double *vecAtQPs;
+
+    WorldVector<double> *gradAtQPs;
+    
+    /** \brief
+     * Function for b.
+     */
+    BinaryAbstractFunction<WorldVector<double>, double, WorldVector<double> > *vecFct;
+  };
+
   // ============================================================================
 
   class General_FOT : public FirstOrderTerm
diff --git a/AMDiS/src/parareal/AdaptParaReal.cc b/AMDiS/src/parareal/AdaptParaReal.cc
index f02e575475b5b2cb1f809f2a764c9c832bebd882..9ddf1176bdb8cccd766ff2f0804a4349d6e4cad5 100644
--- a/AMDiS/src/parareal/AdaptParaReal.cc
+++ b/AMDiS/src/parareal/AdaptParaReal.cc
@@ -12,6 +12,7 @@ namespace AMDiS {
     adaptInfo_->setTimestep(coarseTimestep);
     pararealProb->setStoreSolutions(false, true);
     pararealProb->setStoreInitSolution(true);
+    pararealProb->setInitSolutionVec(NULL);
     AdaptInstationary::adapt();
 
     // ParaReal iterations
diff --git a/AMDiS/src/parareal/AdaptParaReal.h b/AMDiS/src/parareal/AdaptParaReal.h
index c6a66416aca885fec6a8d80716aacd3008ca85dd..943ddb5f2e0e70ca0862829b043f9ce9dce169d5 100644
--- a/AMDiS/src/parareal/AdaptParaReal.h
+++ b/AMDiS/src/parareal/AdaptParaReal.h
@@ -52,7 +52,7 @@ namespace AMDiS {
       TEST_EXIT(fineTimestep > 0.0)("fineTimestep must be greater than zero!\n");
       TEST_EXIT(coarseTimestep > 0.0)("coarseTimestep must be greater than zero!\n");
       TEST_EXIT(fineTimestep < coarseTimestep)("fineTimestep must be smaller than coarseTimestep!\n");
-      TEST_EXIT(pararealIter > 0)("ParaReal Iterations must be greater than zero!\n");
+      TEST_EXIT(pararealIter >= 0)("ParaReal Iterations must be >= 0!\n");
     }
 
     /** \brief
diff --git a/AMDiS/src/parareal/ProblemBase.h b/AMDiS/src/parareal/ProblemBase.h
index d129a5d05014d9bcd52a934595cb1fe17e5aa5a5..a47b2707c2a48f8fa172b7b7bb217c8ded04db01 100644
--- a/AMDiS/src/parareal/ProblemBase.h
+++ b/AMDiS/src/parareal/ProblemBase.h
@@ -165,7 +165,11 @@ namespace AMDiS {
     }
     
     virtual void solveInitialProblem(AdaptInfo *adaptInfo) {
-      ProblemInstatScal::solveInitialProblem(adaptInfo);
+      if (initSolutionVec == NULL) {
+	ProblemInstatScal::solveInitialProblem(adaptInfo);
+      } else {
+	*(problemStat->getSolution()) = *initSolutionVec;
+      }
       
       if (storeInitSolution) {
 	storeSolution(problemStat->getSolution());
@@ -194,7 +198,11 @@ namespace AMDiS {
     }
     
     virtual void solveInitialProblem(AdaptInfo *adaptInfo) {
-      ProblemInstatVec::solveInitialProblem(adaptInfo);
+      if (initSolutionVec == NULL) {
+	ProblemInstatVec::solveInitialProblem(adaptInfo);
+      } else {
+	*(problemStat->getSolution()) = *initSolutionVec;
+      }
       
       if (storeInitSolution) {
 	storeSolution(problemStat->getSolution());