 
         // main loop
         for (int t = 0; t < nt; t++) {
+
+            printf("t = %i\n", t);
             // offset by kernels in previous
             int kx_offset = 0;
             int ky_offset = kernels[0].y * rank * blocks[0].y * threads[0].y;

    int y = ky_offset + blockIdx.y * blockDim.y + threadIdx.y;
    int l = threadIdx.z;

    if ((x > 0) && (x < (nx-1)) && (y > 0) && (y < (ny-1)) && (l < nlayers)) {

        float a = dt * alpha *
            U_half[((y * nx + x) * nlayers + l)*4] * (0.5 / dx) * (sum_phs[(y * nx + x+1) * nlayers + l] -

    int count;
    cudaGetDeviceCount(&count);

    if (rank == 0) {
        cudaError_t err = cudaGetLastError();
        // check that we actually have some GPUS
        if (err != cudaSuccess) {
            printf("Error: %s\n", cudaGetErrorString(err));
            printf("Aborting program.\n");
            return;
        }

        printf("Found %i CUDA devices\n", count);
    }

    int maxThreads = 256;
    int maxBlocks = 256; //64;
