     std::cout << "eps \t # iterations \t error \t hx_max\t hy_max \t time/iteration \n";
     std::cout << eps<<"\t";
     t.tic();
-    //pol.compute_rhs( b, bmod);
-    //unsigned number = invert(pol, x,bmod);// vol3d, v3d );
-    unsigned number = invert(pol, x_fine ,b_fine, vol3d, v3d );
+    pol_refined.compute_rhs( b, bmod);
+    //unsigned number = invert(pol_refined, x,bmod);// vol3d, v3d );
+    unsigned number = invert(pol, x_fine ,b_fine, vol3dFINE, v3dFINE );
     dg::blas2::gemv( P, x_fine, x);
     std::cout <<number<<"\t";
     t.toc();

    dg::blas2::gemv( Q, b, b_fine);
    //pol.set_chi( chi);
    pol.set_chi( chi_fine);
    //compute error
    dg::DVec error( solution);
    const double eps = 1e-10;

    orthogonal::refined::GridX3d<dg::DVec> g3d(add_x, add_y, howmanyX, howmanyY, gp, psi_0, 0.25, 1./22., n_ref, n, Nx, Ny,Nz, dg::DIR, dg::NEU);
    orthogonal::refined::GridX2d<dg::DVec> g2d = g3d.perp_grid();
    dg::Elliptic<orthogonal::refined::GridX3d<dg::DVec>, dg::Composite<dg::DMatrix>, dg::DVec> pol( g3d, dg::not_normed, dg::centered);
    dg::RefinedElliptic<orthogonal::refined::GridX3d<dg::DVec>, dg::IDMatrix, dg::Composite<dg::DMatrix>, dg::DVec> pol_refined( g3d, dg::not_normed, dg::centered);
    psi_1 = g3d.psi1();
    std::cout << "psi 1 is          "<<psi_1<<"\n";

