diff --git a/src/navier-stokes.c b/src/navier-stokes.c
index 705a302..8cfaa7a 100644
--- a/src/navier-stokes.c
+++ b/src/navier-stokes.c
@@ -615,7 +615,7 @@ double ns_step_rkf45(
   bool irreversible
 ){
   int kx,ky;
-  double err;
+  double err,relative;
 
   // k1: u(t)
   ns_rhs(k1, u, K1, K2, N1, N2, nu, L, g, fft1, fft2, ifft, irreversible);
@@ -662,21 +662,25 @@ double ns_step_rkf45(
   
   // compute error
   err=0;
+  relative=0;
   for(kx=0;kx<=K1;kx++){
     for(ky=(kx>0 ? -K2 : 1);ky<=K2;ky++){
       // difference between 5th order and 4th order
       err+=cabs(delta*(1./360*k1[klookup_sym(kx,ky,K2)]-128./4275*k3[klookup_sym(kx,ky,K2)]-2197./75240*k4[klookup_sym(kx,ky,K2)]+1./50*k5[klookup_sym(kx,ky,K2)]+2./55*k6[klookup_sym(kx,ky,K2)]));
+      // next step
+      tmp[klookup_sym(kx,ky,K2)]=delta*(16./135*k1[klookup_sym(kx,ky,K2)]+6656./12825*k3[klookup_sym(kx,ky,K2)]+28561./56430*k4[klookup_sym(kx,ky,K2)]-9./50*k5[klookup_sym(kx,ky,K2)]+2./55*k6[klookup_sym(kx,ky,K2)]);
+      relative+=cabs(tmp[klookup_sym(kx,ky,K2)]);
     }
   }
 
   // new delta
-  delta=factor*delta*pow(tolerance/err,0.2);
+  delta=factor*delta*pow(relative*tolerance/err,0.2);
   // compare relative error with tolerance
-  if(err<tolerance){
+  if(err<relative*tolerance){
     // add to output
     for(kx=0;kx<=K1;kx++){
       for(ky=(kx>0 ? -K2 : 1);ky<=K2;ky++){
-	u[klookup_sym(kx,ky,K2)]+=delta*(16./135*k1[klookup_sym(kx,ky,K2)]+6656./12825*k3[klookup_sym(kx,ky,K2)]+28561./56430*k4[klookup_sym(kx,ky,K2)]-9./50*k5[klookup_sym(kx,ky,K2)]+2./55*k6[klookup_sym(kx,ky,K2)]);
+	u[klookup_sym(kx,ky,K2)]+=tmp[klookup_sym(kx,ky,K2)];
       }
     }
     return delta;