             // read the current neighbor index (MEM TRANSFER: 4 bytes)
             cur_j = next_j;
             if (neigh_idx+tpp < n_neigh)
-                next_j = d_nlist[my_head + neigh_idx+tpp];
+                {
+                #if __CUDA_ARCH__ >= 350
+                next_j = __ldg(d_nlist + my_head + neigh_idx + tpp);
+                #else
+                next_j = tex1Dfetch(pdata_nlist_tex, my_head + neigh_idx + tpp);
+                #endif
+                }
+//                 next_j = d_nlist[my_head + neigh_idx+tpp];
 
             // get the neighbor's position (MEM TRANSFER: 16 bytes)
             Scalar4 postypej = texFetchScalar4(d_pos, pdata_pos_tex, cur_j);

    pdata_charge_tex.filterMode = cudaFilterModePoint;
    cudaBindTexture(0, pdata_charge_tex, pair_args.d_charge, sizeof(Scalar) * pair_args.n_max);
    
    pdata_nlist_tex.normalized = false;
    pdata_nlist_tex.filterMode = cudaFilterModePoint;
    cudaBindTexture(0, pdata_nlist_tex, pair_args.d_nlist, sizeof(unsigned int) * pair_args.size_neigh_list);

    Scalar virialyz = Scalar(0.0);
    Scalar virialzz = Scalar(0.0);

//     #if __CUDA_ARCH__ >= 350
//     unsigned int my_head = __ldg(d_head_list + idx);
//     #else
//     unsigned int my_head = tex1Dfetch(pdata_head_list_tex, idx);
//     #endif
    unsigned int my_head = d_head_list[idx];
    unsigned int cur_j = 0;
    
    // first read only neighborlist fetch
    #if __CUDA_ARCH__ >= 350
    unsigned int next_j = threadIdx.x%tpp < n_neigh ?
        __ldg(d_nlist + my_head + threadIdx.x%tpp) : 0;
    #else
    unsigned int next_j = threadIdx.x%tpp < n_neigh ?
        tex1Dfetch(pdata_nlist_tex, my_head + threadIdx.x%tpp) : 0;
    #endif
    
//     unsigned int next_j = threadIdx.x%tpp < n_neigh ?
//         d_nlist[my_head + threadIdx.x%tpp] : 0;
    // loop over neighbors
    for (int neigh_idx = threadIdx.x%tpp; neigh_idx < n_neigh; neigh_idx+=tpp)
        {
