 
         __syncthreads();
 
-#ifdef _ORIG_SHMEM_
-        lmem = &MEM_BUF[gridDim.x*LMEM_STACK_SIZE*blockDim.x];    //Use the extra large buffer
-#else
         lmem = &MEM_BUF[gridDim.x*(LMEM_STACK_SIZE*blockDim.x + LMEM_EXTRA_SIZE)];    //Use the extra large buffer
-#endif
         apprCount = direCount = 0;
         acc_i = approximate_gravity<blockDim2, 8>( DIM2x, DIM2y, tid, tx, ty,
             body_i, pos_i, group_pos,

            multipole_data, body_pos,
            shmem, lmem, ngb_i, apprCount, direCount, boxSizeInfo, curGroupSize, boxCenterInfo,
            group_eps, body_vel);

#ifdef _ORIG_SHMEM_
        lmem = &MEM_BUF[blockIdx.x* LMEM_STACK_SIZE*blockDim.x]; //Back to normal location
#else
        lmem = &MEM_BUF[blockIdx.x*(LMEM_STACK_SIZE*blockDim.x + LMEM_EXTRA_SIZE)];
#endif

        if(threadIdx.x == 0)
        {

  return y;
} 

/************************************/
/********* SEGMENTED SCAN ***********/
/************************************/

 
#ifdef _KEPLERCODE_
  template<const int SIZE2>
__device__ __forceinline__ int inclusive_segscan_warp(int value, const int distance)
{
#if 0
  const unsigned int laneId = threadIdx.x & (WARP_SIZE - 1);
#endif
  
  const int SIZE = 1 << SIZE2; 

#if 0
  for (int i = 1; i <= SIZE; i <<= 1) 
    value += __shfl_up(value, i, SIZE) & BTEST(laneId >= i && i <= distance);
#else
  for (int i = 0; i < SIZE2; i++)
    value += __shfl_up(value, 1 << i, SIZE) & BTEST(laneId >= (1<<i)) & BTEST((1<<i) <= distance);
#endif

  return value;
}
#else /* !_KEPLERCODE_ */
#error "Only _KEPLERCODE_ is supported at this time"
#endif


__device__ __forceinline__ int bfi(const int x, const int y, const int bit, const int numBits) 
{
	int ret;
	asm("bfi.b32 %0, %1, %2, %3, %4;" : 
		"=r"(ret) : "r"(y), "r"(x), "r"(bit), "r"(numBits));
	return ret;
}

template<const int DIM2>
__device__ __forceinline__ int inclusive_segscan_block64(
    int *shmem, const int tid, const int packed_value, int &dist_block, int &nseg)
{
#if 0
  const int laneId = tid & (WARP_SIZE - 1);
  const int warpId = tid >> WARP_SIZE2;
#endif

  const int  flag = packed_value < 0;
  const int  mask = BTEST(flag);
  const int value = (mask & (-1-packed_value)) + (~mask & 1);
  
  shmem[(warpId << WARP_SIZE2) + WARP_SIZE - 1 - laneId] = flag;
  const int flags1 = __ballot(shmem[tid]);
  shmem[tid] = __popc(flags1);
  __syncthreads();
  nseg += shmem[0] + shmem[WARP_SIZE];
  __syncthreads();
  
  shmem[tid] = __clz (flags1);
  __syncthreads();

  const int flags = __ballot(flag) & bfi(0, 0xffffffff, 0, laneId + 1);
	const int distance = __clz(flags) + laneId - 31;

  int dist0  = shmem[WARP_SIZE];
  dist_block = shmem[0] + (BTEST(shmem[0] == WARP_SIZE) & dist0);
  
  __syncthreads();

  shmem[tid] = inclusive_segscan_warp<WARP_SIZE2>(value, distance);
  __syncthreads();
 
  shmem[tid] += shmem[WARP_SIZE - 1] & BTEST(tid >= WARP_SIZE) & BTEST(tid < WARP_SIZE + dist0);
  __syncthreads();

  const int val = shmem[(1 << DIM2) - 1];
  __syncthreads();

  return val;
}

/* does not work if segment size > 64 (= thead block size) */
template<const int DIM2>
__device__ __forceinline__ int inclusive_segscan_array(int *shmem_in, const int N, const int tid)
{
  const int DIM = 1 << DIM2;

  int dist, nseg = 0;
  int y  = inclusive_segscan_block64<DIM2>(shmem_in, tid, shmem_in[tid], dist, nseg);
  if (N <= DIM) return nseg;

  for (int p = DIM; p < N; p += DIM)
  {
    int *shmem = shmem_in + p;
    int y1  = inclusive_segscan_block64<DIM2>(shmem, tid, shmem[tid], dist, nseg);
    shmem[tid] += y & BTEST(tid < dist);
    y = y1;
  }

  return nseg;
}


/*************** Tree walk ************/


  template<int SHIFT>
