         time += dt;
     }
     ////////////////////////////////////////////////////////////////////
-    cout << "Analytic formula enstrophy "<<lamb.enstrophy()<<endl;
-    cout << "Analytic formula energy    "<<lamb.energy()<<endl;
+    //cout << "Analytic formula enstrophy "<<lamb.enstrophy()<<endl;
+    //cout << "Analytic formula energy    "<<lamb.energy()<<endl;
     cout << "Total vorticity           is: "<<blas2::dot( stencil.data(), S2D<double, n>(hx, hy), y0) << "\n";
     cout << "Relative enstrophy error  is: "<<(0.5*blas2::dot( S2D<double, n>(hx, hy), y0) - enstrophy)/enstrophy<<"\n";
     test( y0, y1); //get the potential ready
     cout << "Relative energy error     is: "<<(0.5*blas2::dot( test.potential(), S2D<double, n>(hx, hy), y0) - energy)/energy<<"\n";
+
+    blas1::axpby( 1., sol, -1., y0);
+    cout << "Distance to solution "<<sqrt( blas2::dot( S2D<double,n >(hx,hy), y0)) << endl;
 
     //energy and enstrophy errrors are due to timestep only ( vorticity is exactly conserved)
     // k = 2 | p = 3
     // k = 3 | p = 4
     // k = 4 | p = 5
 
+    //solution to sin(x)sin(y) 
+    // n = 1 
+    // n = 2 | p = 2
+    // n = 3 | p = 2.6
+    // n = 4 | p = 4
+
     return 0;
 
 }

using namespace std;
using namespace dg;

const unsigned n = 3;
const unsigned Nx = 33;
const unsigned Ny = 33;
const double lx = 1.;
const double ly = 1.;

const unsigned k = 3;
const double D = 0.0;


typedef cusp::device_memory Memory;

double initial( double x, double y){ return 2.*sin(x)*sin(y);}
double solution( double x, double y){ return 2.*sin(x)*sin(y)*exp(-2.*D*T);}
using namespace std;

int main()
