                      cusp::array1d<ValueType,cusp::device_memory>& x) const
     {
         // TODO check sizes
-        cusp::coo_matrix<IndexType,ValueType,MemorySpace> & A = levels[0].A;
+        const cusp::coo_matrix<IndexType,ValueType,MemorySpace> & A = levels[0].A;
 
         cusp::array1d<ValueType,MemorySpace> residual(A.num_rows);  // TODO eliminate temporaries
             

        ASSERT_EQUAL(P.row_indices[4], 2); ASSERT_EQUAL(P.column_indices[4], 1);
        ASSERT_EQUAL(P.row_indices[5], 3); ASSERT_EQUAL(P.column_indices[5], 1);

        ASSERT_ALMOST_EQUAL(P.values[0],  0.68426213);
        ASSERT_ALMOST_EQUAL(P.values[1],  0.68426213);
        ASSERT_ALMOST_EQUAL(P.values[2], -0.18426213);
        ASSERT_ALMOST_EQUAL(P.values[3], -0.18426213);
        ASSERT_ALMOST_EQUAL(P.values[4],  0.68426213);
        ASSERT_ALMOST_EQUAL(P.values[5],  0.68426213);
    }
}
DECLARE_HOST_DEVICE_UNITTEST(TestSmoothProlongator);

            cusp::blas::axpby(b, residual, residual, 1.0f, -1.0f);
            ValueType norm = cusp::blas::nrm2(residual);

            //printf("%10.8f  %6.4f\n", norm, norm/last_norm);

            last_norm = norm;
        }
