diff --git a/install b/install
index 217830e0b630bc3bb0a67552eec438508640725c..0c15ff9b1e3e6f64833becbfd12d2b004e53e2c2 100755
--- a/install
+++ b/install
@@ -2,7 +2,7 @@
 
 TARGET=TNL
 INSTALL_PREFIX=${HOME}/local
-WITH_CUDA=yes
+WITH_CUDA=no
 TEMPLATE_EXPLICIT_INSTANTIATION=yes
 #VERBOSE="VERBOSE=1"
 
diff --git a/src/functions/tnlConstantFunction.h b/src/functions/tnlConstantFunction.h
index e7fcfd2aa474563627ef9d1e3c71f46c4d844117..5188269769979fcf0ecf7d2a8cf462e2424991e9 100644
--- a/src/functions/tnlConstantFunction.h
+++ b/src/functions/tnlConstantFunction.h
@@ -48,14 +48,14 @@ class tnlConstantFunction
              int ZDiffOrder,
              typename Vertex >
 #else
-   template< int XDiffOrder = 0,
+   template< int XDiffOrder,
              int YDiffOrder = 0,
              int ZDiffOrder = 0,
              typename Vertex = VertexType >
 #endif
    RealType getValue( const Vertex& v,
                       const Real& time = 0.0 ) const;
-#ifdef HAVE_NOT_CXX11
+
    template< typename Vertex >
    RealType getValue( const Vertex& v,
                       const Real& time = 0.0 ) const
@@ -63,8 +63,6 @@ class tnlConstantFunction
       return getValue< 0, 0, 0, Vertex >( v, time );
    }
 
-#endif
-
    protected:
 
    RealType value;
diff --git a/src/implementation/core/arrays/tnlMultiArray1D_impl.h b/src/implementation/core/arrays/tnlMultiArray1D_impl.h
index fe408d60e711c5493f36f318edc4554c28300217..599aaaff8ebada962042c638afd2a95f647f1c61 100644
--- a/src/implementation/core/arrays/tnlMultiArray1D_impl.h
+++ b/src/implementation/core/arrays/tnlMultiArray1D_impl.h
@@ -103,7 +103,7 @@ template< typename Element, typename Device, typename Index >
 Index tnlMultiArray< 1, Element, Device, Index > :: getElementIndex( const Index i ) const
 {
    tnlAssert( i >= 0 && i < this -> dimensions[ 0 ],
-              printf( "i = %d this -> dimensions[ 0 ] = %d \n", ( int ) i, ( int ) this -> dimensions[ 0 ] ) );
+              cerr << "i = " << i << " this -> dimensions[ 0 ] = " <<  this -> dimensions[ 0 ] );
    return i;
 }
 
diff --git a/src/implementation/core/arrays/tnlMultiArray2D_impl.h b/src/implementation/core/arrays/tnlMultiArray2D_impl.h
index 62781a9fa9bc68f7b4074b764d7570066c59fd41..da5eb0e1ee878313f9df3ecfc1182cea8cf8832f 100644
--- a/src/implementation/core/arrays/tnlMultiArray2D_impl.h
+++ b/src/implementation/core/arrays/tnlMultiArray2D_impl.h
@@ -115,8 +115,8 @@ template< typename Element, typename Device, typename Index >
 Index tnlMultiArray< 2, Element, Device, Index > :: getElementIndex( const Index j, const Index i ) const
 {
    tnlAssert( i >= 0 && i < this -> dimensions[ 0 ] && j >= 0 && j < this -> dimensions[ 1 ],
-              printf( "i = %d j = %d this -> dimensions[ 0 ] = %d this -> dimensions[ 1 ] = %d \n",
-               i, j, this -> dimensions[ 0 ], this -> dimensions[ 1 ] ) );
+              cerr << "i = " << i << " j = " << j << " this -> dimensions[ 0 ] = " <<  this -> dimensions[ 0 ]
+                   << " this -> dimensions[ 1 ] = " << this -> dimensions[ 1 ] );
    return j * this -> dimensions[ 0 ] + i;
 }
 
diff --git a/src/implementation/functions/tnlConstantFunction_impl.h b/src/implementation/functions/tnlConstantFunction_impl.h
index 9c83b654c3021f16d9f4cca93c244babe191b2bd..f4c2c6b7de691227caabd0a0f076887daa1f0564 100644
--- a/src/implementation/functions/tnlConstantFunction_impl.h
+++ b/src/implementation/functions/tnlConstantFunction_impl.h
@@ -51,7 +51,7 @@ tnlConstantFunction< Dimensions, Real >::
 setup( const tnlParameterContainer& parameters,
       const tnlString& prefix )
 {
-   this->setValue( parameters.GetParameter< double >( prefix + "-value") );
+   this->setValue( parameters.GetParameter< double >( prefix + "value") );
    return true;
 }
 
diff --git a/src/implementation/functions/tnlTestFunction_impl.h b/src/implementation/functions/tnlTestFunction_impl.h
index 901415b6da2feb65b6cf0db86fe72f77456551a6..fe15d1f7e1460ce5de9ad910f14d84edb38f4809 100644
--- a/src/implementation/functions/tnlTestFunction_impl.h
+++ b/src/implementation/functions/tnlTestFunction_impl.h
@@ -44,9 +44,10 @@ configSetup( tnlConfigDescription& config,
              const tnlString& prefix )
 {
    config.addRequiredEntry< tnlString >( prefix + "test-function", "Testing function." );
+      config.addEntryEnum( "constant" );
+      config.addEntryEnum( "exp-bump" );
       config.addEntryEnum( "sin-wave" );
       config.addEntryEnum( "sin-bumps" );
-      config.addEntryEnum( "exp-bump" );
    config.addEntry     < double >( prefix + "value", "Value of the constant function.", 0.0 );
    config.addEntry     < double >( prefix + "wave-length", "Wave length of the sine based test functions.", 1.0 );
    config.addEntry     < double >( prefix + "wave-length-x", "Wave length of the sine based test functions.", 1.0 );
diff --git a/src/implementation/matrices/tnlCSRMatrix_impl.h b/src/implementation/matrices/tnlCSRMatrix_impl.h
index 8ae0f660fc1bd137a947dfd8c7cab2b1397c9794..400e3a3ae017abd25007b7b8735ea8d2809733bc 100644
--- a/src/implementation/matrices/tnlCSRMatrix_impl.h
+++ b/src/implementation/matrices/tnlCSRMatrix_impl.h
@@ -426,8 +426,7 @@ typename Vector::RealType tnlCSRMatrix< Real, Device, Index >::rowVectorProduct(
    const IndexType rowEnd = this->rowPointers[ row + 1 ];
    IndexType column;
    while( elementPtr < rowEnd &&
-          ( column = this->columnIndexes[ elementPtr ] ) < this->columns &&
-          column != this->getPaddingIndex() )
+          ( column = this->columnIndexes[ elementPtr ] ) != this->getPaddingIndex() )
       result += this->values[ elementPtr++ ] * vector[ column ];
    return result;
 }
@@ -472,9 +471,9 @@ template< typename Real,
           typename Index >
    template< typename Vector >
 bool tnlCSRMatrix< Real, Device, Index >::performSORIteration( const Vector& b,
-                                                                                    const IndexType row,
-                                                                                    Vector& x,
-                                                                                    const RealType& omega ) const
+                                                               const IndexType row,
+                                                               Vector& x,
+                                                               const RealType& omega ) const
 {
    tnlAssert( row >=0 && row < this->getRows(),
               cerr << "row = " << row
@@ -487,12 +486,12 @@ bool tnlCSRMatrix< Real, Device, Index >::performSORIteration( const Vector& b,
    IndexType elementPtr = this->rowPointers[ row ];
    const IndexType rowEnd = this->rowPointers[ row + 1 ];
    IndexType column;
-   while( elementPtr < rowEnd && ( column = this->columnIndexes[ elementPtr ] ) < this->columns )
+   while( elementPtr < rowEnd && ( column = this->columnIndexes[ elementPtr ] ) != this->getPaddingIndex() )
    {
       if( column == row )
-         diagonalValue = this->values.getElement( elementPtr );
+         diagonalValue = this->values[ elementPtr ];
       else
-         sum += this->values.getElement( row * this->diagonalsShift.getSize() + elementPtr ) * x. getElement( column );
+         sum += this->values[ elementPtr ] * x[ column ];
       elementPtr++;
    }
    if( diagonalValue == ( Real ) 0.0 )
@@ -500,7 +499,7 @@ bool tnlCSRMatrix< Real, Device, Index >::performSORIteration( const Vector& b,
       cerr << "There is zero on the diagonal in " << row << "-th row of the matrix " << this->getName() << ". I cannot perform SOR iteration." << endl;
       return false;
    }
-   x. setElement( row, x[ row ] + omega / diagonalValue * ( b[ row ] - sum ) );
+   x[ row ] = ( 1.0 - omega ) * x[ row ] + omega / diagonalValue * ( b[ row ] - sum );
    return true;
 }
 
diff --git a/src/implementation/matrices/tnlDenseMatrix_impl.h b/src/implementation/matrices/tnlDenseMatrix_impl.h
index 4d41844018405abcd8937b011ff4960d20018721..551804c99fbe1932f0433e488eb107562f57e683 100644
--- a/src/implementation/matrices/tnlDenseMatrix_impl.h
+++ b/src/implementation/matrices/tnlDenseMatrix_impl.h
@@ -142,7 +142,8 @@ bool tnlDenseMatrix< Real, Device, Index >::setElementFast( const IndexType row,
 {
    tnlAssert( row >= 0 && row < this->getRows() &&
               column >= 0 && column < this->getColumns(),
-              printf( " row = %d, column = %d, this->getRows = %d, this->getColumns() = %d \n", row, column, this->getRows(), this->getColumns() ) );
+              cerr << " row = " << row << " column = " << column << " this->getRows() = " << this->getRows()
+                   << " this->getColumns() = " << this->getColumns() );
    this->values.operator[]( this->getElementIndex( row, column ) ) = value;
    return true;
 }
diff --git a/src/implementation/solvers/linear/stationary/tnlSORSolver_impl.h b/src/implementation/solvers/linear/stationary/tnlSORSolver_impl.h
index 54c121d9eb3dbf12e341dd39abe32f1d13c913db..73e9a0865b901cda34278be2e740a68476eefa70 100644
--- a/src/implementation/solvers/linear/stationary/tnlSORSolver_impl.h
+++ b/src/implementation/solvers/linear/stationary/tnlSORSolver_impl.h
@@ -51,7 +51,7 @@ setup( const tnlParameterContainer& parameters,
        const tnlString& prefix )
 {
    tnlIterativeSolver< RealType, IndexType >::setup( parameters, prefix );
-   this->setOmega( parameters.GetParameter< int >( "sor-omega" ) );
+   this->setOmega( parameters.GetParameter< double >( prefix + "sor-omega" ) );
 }
 
 
@@ -94,23 +94,25 @@ bool tnlSORSolver< Matrix, Preconditioner > :: solve( const Vector& b,
 
    RealType bNorm = b. lpNorm( ( RealType ) 2.0 );
 
-   while( this -> getIterations() < this -> getMaxIterations() &&
-          this -> getResidue() > this -> getConvergenceResidue() )
+   while( this->nextIteration() )
    {
-      /*matrix -> performSORIteration( this -> getOmega(),
-                                     b,
-                                     x,
-                                     0,
-                                     size );*/
-      if( this -> getIterations() % 10 == 0 )
-         this -> setResidue( ResidueGetter :: getResidue( *matrix, b, x, bNorm ) );
-      if( ! this -> nextIteration() )
-         return false;
+      for( IndexType row = 0; row < size; row ++ )
+         matrix->performSORIteration( b,
+                                      row,
+                                      x,
+                                      this->getOmega() );
+      //if( this -> getIterations() % 10 == 0 )
+         this -> setResidue( ResidueGetter :: getResidue( *matrix, x, b, bNorm ) );
       this -> refreshSolverMonitor();
    }
-   this -> setResidue( ResidueGetter :: getResidue( *matrix, b, x, bNorm ) );
+   this -> setResidue( ResidueGetter :: getResidue( *matrix, x, b, bNorm ) );
    this -> refreshSolverMonitor();
-      if( this -> getResidue() > this -> getConvergenceResidue() ) return false;
+   if( this -> getResidue() > this -> getConvergenceResidue() ||
+       std::isnan( this->getResidue() ) )
+   {
+      cerr << "The residue ( " << this->getResidue() << " ) is over the convergence residue ( " << this -> getConvergenceResidue() << ")." << endl;
+      return false;
+   }
    return true;
 };
 
diff --git a/src/implementation/solvers/pde/tnlSemiImplicitTimeStepper_impl.h b/src/implementation/solvers/pde/tnlSemiImplicitTimeStepper_impl.h
index 785866a56b0d8e870ff1047bfd51a6f462e459fc..becc60c4a2f7c85249c5b3f873f0012d2f3a2cf5 100644
--- a/src/implementation/solvers/pde/tnlSemiImplicitTimeStepper_impl.h
+++ b/src/implementation/solvers/pde/tnlSemiImplicitTimeStepper_impl.h
@@ -159,9 +159,11 @@ solve( const RealType& time,
          cout << "                                                                  Solving the linear system for time " << t << "             \r" << flush;
       if( ! this->linearSystemSolver->template solve< DofVectorType, tnlLinearResidueGetter< MatrixType, DofVectorType > >( this->rightHandSide, dofVector ) )
       {
-         cerr << "The linear system solver did not converge." << endl;
+         cerr << endl << "The linear system solver did not converge." << endl;
          return false;
       }
+      //if( verbose )
+      //   cout << endl;
       if( ! this->problem->postIterate( t,
                                         currentTau,
                                         mesh,
diff --git a/src/implementation/solvers/tnlIterativeSolver_impl.h b/src/implementation/solvers/tnlIterativeSolver_impl.h
index e423b30e5af2a2bb5e5da47aa411f248acb2f966..ff700563d32f72768884437b13bd3f79070cc7e8 100644
--- a/src/implementation/solvers/tnlIterativeSolver_impl.h
+++ b/src/implementation/solvers/tnlIterativeSolver_impl.h
@@ -97,10 +97,23 @@ bool tnlIterativeSolver< Real, Index> :: nextIteration()
        this->currentIteration % this->refreshRate == 0 )
       solverMonitor->refresh();
 
-   if( std::isnan( this->getResidue() ) ||
-       ( this->getResidue() > this->getDivergenceResidue() &&
-         this->getIterations() > this->minIterations ) ||
-       this->getIterations() > this->getMaxIterations() )
+   if( std::isnan( this->getResidue() ) )
+   {
+      cerr << "The residue is NaN." << endl;
+      return false;
+   }
+   if(( this->getResidue() > this->getDivergenceResidue() &&
+         this->getIterations() > this->minIterations ) )
+   {
+      cerr << "The residue has exceeded allowed tolerance " << this->getDivergenceResidue() << "." << endl;
+      return false;
+   }
+   if( this->getIterations() > this->getMaxIterations() )
+   {
+      cerr << "The solver has exceeded maximal allowed number of iterations " << this->getMaxIterations() << "." << endl;
+      return false;
+   }
+   if( this->getResidue() < this->getConvergenceResidue() )
       return false;
    return true;
 }
diff --git a/tools/src/tnl-init.h b/tools/src/tnl-init.h
index 088f16fbecb7ac5c93ad2a162f9d249eabe65fc9..90512ef3b076c2742c7bbcb079a7afdc7802bc8d 100644
--- a/tools/src/tnl-init.h
+++ b/tools/src/tnl-init.h
@@ -55,7 +55,7 @@ bool renderFunction( const tnlParameterContainer& parameters )
    double tau = parameters.GetParameter< double >( "snapshot-period" );
    bool numericalDifferentiation = parameters.GetParameter< bool >( "numerical-differentiation" );
    int step( 0 );
-   const int steps = ceil( finalTime / tau );
+   const int steps = tau > 0 ? ceil( finalTime / tau ): 0;
 
    while( step <= steps )
    {