   // Always use device 0
   cudaSetDevice(0);
 
+  CHECK_CALL(cudaDeviceSetCacheConfig(cudaFuncCachePreferShared));
+
   // Parse options
   int opt;
   while((opt = getopt_long(argc, argv, "cen:pr:t:vx:y:", commandLineOptions, NULL)) != -1) {

    // Sync before overwriting shared
    __syncthreads();

    buffer[threadIdx.x] = ((baseOffset >= 0) && (baseOffset <= (problemSize-1))) ? average : buffer[threadIdx.x];

    // Sync before re-reading shared
    __syncthreads();
  }

  if(threadIdx.x >= (TIME_TILE_SIZE-1) &&
     threadIdx.x <= (blockDim.x-1-(TIME_TILE_SIZE-1))) {
    output[baseOffset] = buffer[threadIdx.x];
#ifdef DEBUG
    printf("[%d, %d]: Write output at %d\n", blockIdx.x, threadIdx.x, baseOffset);
#endif

__global__
void jacobi_2d_kernel_overlapped(FloatType* input, FloatType* output,
                                 int32_t problemSize) {
  // Determine our start position
  int offsetI = blockIdx.x * (blockDim.x-2*(TIME_TILE_SIZE-1)) + threadIdx.x;
  offsetI -= TIME_TILE_SIZE-1;
  int offsetJ = blockIdx.y * (blockDim.y-2*(TIME_TILE_SIZE-1)) + threadIdx.y;
  offsetJ -= TIME_TILE_SIZE-1;

  // Load data into shared
  //FloatType c = ((baseOffset >= 0) && (baseOffset <= (problemSize-1))) ? input[baseOffset] : ZERO;
  //FloatType l = (baseOffset > 0) ? input[baseOffset-1] : ZERO;
  //FloatType r = (baseOffset < (problemSize-1)) ? input[baseOffset+1] : ZERO;

  FloatType c = ((offsetI >= 0) && (offsetI <= (problemSize-1)) &&
                 (offsetJ >= 0) && (offsetJ <= (problemSize-1)))
                 ? input[offsetI*problemSize+offsetJ] : ZERO;
  FloatType l = ((offsetI-1 >= 0) && (offsetI-1 <= (problemSize-1)) &&
                 (offsetJ >= 0) && (offsetJ <= (problemSize-1)))
                 ? input[(offsetI-1)*problemSize+offsetJ] : ZERO;
  FloatType r = ((offsetI+1 >= 0) && (offsetI+1 <= (problemSize-1)) &&
                 (offsetJ >= 0) && (offsetJ <= (problemSize-1)))
                 ? input[(offsetI+1)*problemSize+offsetJ] : ZERO;
  FloatType t = ((offsetI >= 0) && (offsetI <= (problemSize-1)) &&
                 (offsetJ-1 >= 0) && (offsetJ-1 <= (problemSize-1)))
                 ? input[offsetI*problemSize+offsetJ-1] : ZERO;
  FloatType b = ((offsetI >= 0) && (offsetI <= (problemSize-1)) &&
                 (offsetJ+1 >= 0) && (offsetJ+1 <= (problemSize-1)))
                 ? input[offsetI*problemSize+offsetJ+1] : ZERO;

  FloatType average = (l + c + r + t + b) / FIVE;
  buffer[threadIdx.x*blockDim.x+threadIdx.y] = 
    ((offsetI >= 0) && (offsetI <= (problemSize-1)) &&
    (offsetJ >= 0) && (offsetJ <= (problemSize-1))) ? average : buffer[threadIdx.x*blockDim.x+threadIdx.y];

  __syncthreads();

