 #endif//DG_BENCHMARK
         for( unsigned j=0; j<p.itstp; j++)
         {
-            karniadakis.step( exp, imp, time, y1);
+            //karniadakis.step( exp, imp, time, y0);
+            stepper.step( exp, time, y0, time, y0, dt, dg::pid_control, dg::l2norm, 1e-5, 1e-10);
+            if ( stepper.failed() ) failed_counter ++;
+            DG_RANK0 std::cout << "Time "<<time<<" dt "<<dt<<" failed counter "<<failed_counter<<"\n";
             //store accuracy details
             {
                 DG_RANK0 std::cout << "(m_tot-m_0)/m_0: "<< (exp.mass()-mass0)/mass_blob0<<"\t";

    ///////////////////////////////////////Timeloop/////////////////////////////////
    const double mass0 = exp.mass(), mass_blob0 = mass0 - grid.lx()*grid.ly();
    double E0 = exp.energy(), E1 = 0, diff = 0;
    dg::Timer t;
    t.tic();
    try

    std::vector<dg::x::DVec> transferD(4, dg::evaluate(dg::zero, grid_out));
    dg::x::HVec transferH(dg::evaluate(dg::zero, grid_out));
    dg::x::IDMatrix interpolate = dg::create::interpolation( grid_out, grid);
    dg::blas2::symv( interpolate, y0[0], transferD[0]);
    dg::blas2::symv( interpolate, y0[1], transferD[1]);
    dg::blas2::symv( interpolate, exp.potential()[0], transferD[2]);
    dg::blas2::symv( imp.laplacianM(), exp.potential()[0], transfer);
    dg::blas2::symv( interpolate, transfer, transferD[3]);
