What can I do to optimize this cuda kernel (generate random numbers)

I am trying to generate array of random numbers using cuda cuRand.

My questions:

1. How can I calculate how fast my kernel is? (with data below) (in terms of GB/s)
   (Right now i think its about 3Gb/s)
2. What optimizations can I use here to reach better speed?

My specifications:

1. OS: Windows 10 x64 Enterprise.
2. CPU: Intel Core i7-4720HQ @ 260 GHz
3. GPU: NVIDIA GeForce GTX 960M
4. IDE: Visual Studio 2015 update 3 + cuda toolkit.
5. CUDA toolkit version: 8.0.61
6. GPU Bus Width: 128 Bit (80.2 GB/s)
7. GPU Compute Capability: 5.0
8. Test array size: 87654321
   (I am calling void cudaGenRandNum(const long long int & size, T * const od_values) with size == 87654321 and od_values is empty array of doubles on devcie)
9. Test array size in Mb: 668.75 Mb

Code:

__global__
void init_kernel_1D(const long long int size, const long long int seed, curandState_t * const od_state) {
    const long long int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < size) {
        curand_init(seed, idx, 0, od_state + idx);
    }
}

template <typename T> __global__
void rand_kernel_1D(const long long int size, const long long int lim, curandState_t * const state, T * const od_value) {
    const long long int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if ((idx < size) && (idx < lim)) {
        od_value[idx] = curand_normal_double(state + idx);
    }
}

template <typename T> inline __forceinline__
void cudaGenRandNum(const long long int & size, T * const od_values) {

    // Size of curandState array:
    const long long int step = 100000;
    dim3 block(256, 1, 1);
    dim3 grid(((step + block.x - 1) / block.x), 1, 1);

    // Device memory: 
    curandState * d_states;
    ERRCHK(cudaMalloc(&d_states, step * sizeof(curandState)); // sizeof( curandState ) * step = about 3.89 Mb

    // Initialize states:
    init_kernel_1D <<<grid, block>>> (step, time(NULL), d_states);
    ERRCHK(cudaGetLastError());
    ERRCHK(cudaDeviceSynchronize());

    // Gen. rand. nubers:
    long long int remained = size;
    long long int offest = 0;
    while (remained > 0) {
        rand_kernel_1D <T> <<<grid, block>>> (step, size - offest, d_states, od_values + offest);
        ERRCHK(cudaGetLastError());

        remained -= step;
        offest += step;
    }
    ERRCHK(cudaDeviceSynchronize());
    ERRCHK(cudaFree(d_states));
}

Profiling result:

Time(%)      Time     Calls       Avg       Min       Max  Name
 85.55%  4.27294s     20171  211.84us  127.87us  332.64us  void rand_kernel_1D<double>(__int64, __int64, curandStateXORWOW*, double*)
 13.31%  664.75ms        23  28.902ms  28.721ms  29.954ms  init_kernel_1D(__int64, __int64, curandStateXORWOW*)
  1.14%  56.878ms         1  56.878ms  56.878ms  56.878ms  [CUDA memcpy DtoH]

Result from custom benchmark:

My custom benchmark function is repeatedly calling cudaGenRandNum<double>(size, d_values); for about 10 seconds (inside while loop) and prints average time.

Time measurements:
 0) Nanoseconds:   219925721.174
 1) Microseconds:  219925.721
 2) Milliseconds:  219.926
 3) Seconds:       0.220
 4) Minutes:       0.004
 5) Hours:         0.000