diff --git a/py_template/test_step.cpp b/py_template/test_step.cpp
index f76cda73ecd6d76c6940e132dfbe13ee5a8e49c7..deca5c8ccf21fedf0c8a7bbaabbb9e01490205e2 100644
--- a/py_template/test_step.cpp
+++ b/py_template/test_step.cpp
@@ -6,11 +6,9 @@ fs->cd->read(
         "Kdata0",
         (void*)fs->cvorticity);
 DEBUG_MSG("field read\n");
+DEBUG_MSG("######### %d %g\n", fs->iteration, fs->correl_vec(fs->cvorticity, fs->cvorticity));
 fs->step(0.0001);
-DEBUG_MSG("after time step\n");
-fs->cd->write(
-        "Kdata1",
-        (void*)fs->cvorticity);
+DEBUG_MSG("######### %d %g\n", fs->iteration, fs->correl_vec(fs->cvorticity, fs->cvorticity));
 
 delete fs;
 DEBUG_MSG("fluid_solver object deleted\n");
diff --git a/src/fluid_solver.cpp b/src/fluid_solver.cpp
index 88c522affbe82b3bcbb5d0c3276b2b5eb4cb823e..fbfe6ce5e245b0a211b1e2aabe41d1b75509085b 100644
--- a/src/fluid_solver.cpp
+++ b/src/fluid_solver.cpp
@@ -250,6 +250,7 @@ void fluid_solver<R>::step(double dt) \
             this->cv[3][cindex*3+2][1] = (this->cv[0][cindex*3+2][1]*factor0 + 2*(this->cv[1][cindex*3+2][1] + dt*this->cu[cindex*3+2][1])*factor1)/3; \
             ); \
  \
+    this->iteration++; \
 }
 /*****************************************************************************/
 
diff --git a/src/fluid_solver_base.cpp b/src/fluid_solver_base.cpp
index fadf5d79bf5d0dcd19ca8cb605bd510773bad84a..645aa67591e25a1936d19d434804f193b66b05d9 100644
--- a/src/fluid_solver_base.cpp
+++ b/src/fluid_solver_base.cpp
@@ -40,6 +40,8 @@ fluid_solver_base<R>::fluid_solver_base( \
         double DKY, \
         double DKZ) \
 { \
+    this->iteration = 0; \
+ \
     int ntmp[4]; \
     ntmp[0] = nz; \
     ntmp[1] = ny; \
@@ -119,6 +121,27 @@ fluid_solver_base<R>::~fluid_solver_base() \
  \
     delete this->cd; \
     delete this->rd; \
+} \
+ \
+template<> \
+R fluid_solver_base<R>::correl_vec(C *a, C *b) \
+{ \
+    double val_process = 0.0, val; \
+    double k2; \
+    int factor = 1;\
+    CLOOP( \
+            k2 = (this->kx[xindex]*this->kx[xindex] + \
+                  this->ky[yindex]*this->ky[yindex] + \
+                  this->kz[zindex]*this->kz[zindex]); \
+            factor = (xindex == 0) ? 1 : 2; \
+            val_process += factor * ((*(a + cindex))[0] * (*(b + cindex))[0] + \
+                                     (*(a + cindex))[1] * (*(b + cindex))[1]); \
+            );\
+    MPI_Allreduce( \
+            (void*)(&val_process), \
+            (void*)(&val), \
+            1, MPI_DOUBLE_PRECISION, MPI_SUM, MPI_COMM_WORLD); \
+    return R(val); \
 }
 
 /*****************************************************************************/
diff --git a/src/fluid_solver_base.hpp b/src/fluid_solver_base.hpp
index dd284139f5b76c5a1fe0db9f9b17eda4b7f8759b..f0001e5323fd06731e262f4074c222346f5de6a4 100644
--- a/src/fluid_solver_base.hpp
+++ b/src/fluid_solver_base.hpp
@@ -68,6 +68,8 @@ class fluid_solver_base
                 double DKY = 1.0,
                 double DKZ = 1.0);
         ~fluid_solver_base();
+
+        rnumber correl_vec(cnumber *a, cnumber *b);
 };
 
 
diff --git a/test.py b/test.py
index b2bdd40f29b451ffb6b72a8dbb24e3e6f91c8056..48432b616142e32f291922b1a56253d1a3de4ab9 100755
--- a/test.py
+++ b/test.py
@@ -60,19 +60,5 @@ Kdata0[..., 0] = Kdata00
 Kdata0[..., 1] = Kdata01
 Kdata0[..., 2] = Kdata02
 Kdata0.tofile("Kdata0")
-run_test('test_step')
-Kdata1 = np.fromfile('Kdata1', dtype = np.complex64).reshape(Kdata0.shape)
-
-#print np.max(np.abs(Kdata0 - Kdata1))
-print np.max(np.abs(Kdata0))
-print np.max(np.abs(Kdata1))
-#
-#fig = plt.figure(figsize=(12, 6))
-#a = fig.add_subplot(121)
-#a.imshow(abs(Kdata0[4, :, :, 2]), interpolation = 'none')
-#a = fig.add_subplot(122)
-#a.imshow(abs(Kdata1[4, :, :, 2]), interpolation = 'none')
-#fig.savefig('tmp.pdf', format = 'pdf')
-
-
+run_test('test_step', ncpu = 1)