 
   int *directS = shmem;                              //  0*DIM,  1*DIM,  1*DIM
   int *nodesS  = directS + DIM;                      //  1*DIM, 10*DIM,  9*DIM
-  int *prefix  = nodesS  + DIM*9;                    // 10*DIM, 11*DIM,  1*DIM
+  int *prefix  = nodesS  + DIM*8;                    //  9*DIM, 10*DIM,  1*DIM
 
   const int NJMAX = DIM*3;
   int    *body_list = (int*   )&nodesS   [  DIM]; //  2*DIM,   5*DIM,  2*DIM

  __syncthreads();


  //Sum the interaction counters
  float  *sh_ds2 = (float*)&sh_acc[DIM];
  int    *sh_ngb = (int*  )&sh_ds2[DIM];
  sh_ds2[tid] = direCount;
  sh_ngb[tid] = apprCount;

  __syncthreads();


  if (ty == 0) {
#pragma unroll
    for (int i = 1; i < DIMy; i++){
      int idx = (i << DIM2x) + tx;
      direCount  += sh_ds2[idx];
      apprCount  += sh_ngb[idx];
    }
  }
  __syncthreads();

  return acc_i;
}



    const int blockDim2 = NTHREAD2;
    __shared__ int shmem[10*(1 << blockDim2)];
    //    __shared__ int shmem[24*(1 << blockDim2)]; is possible on FERMI
    //    int             lmem[LMEM_STACK_SIZE];

