Mersenne Twister CUDA для произвольного количества потоков

Реализация CUDA генератора случайных чисел Mersenne Twister (MT) ограничена максимальным количеством потоков/блоков из 256 и 200 блоков/сетки, т. е. максимальное количество потоков равно 51200.

Поэтому невозможно запустить ядро, использующее МТ с

kernel<<<blocksPerGrid, threadsPerBlock>>>(devMTGPStates, ...)

куда

int blocksPerGrid = (n+threadsPerBlock-1)/threadsPerBlock;

а n — общее количество потоков.

Как лучше всего использовать MT для threads > 51200?

Мой подход, если использовать постоянные значения для blocksPerGrid и threadsPerBlock, например. <<<128,128>>> и используйте в коде ядра следующее:

__global__ void kernel(curandStateMtgp32 *state, int n, ...) { 

    int id = threadIdx.x+blockIdx.x*blockDim.x;

    while (id < n) {

        float x = curand_normal(&state[blockIdx.x]);
        /* some more calls to curand_normal() followed
           by the algorithm that works with the data */

        id += blockDim.x*gridDim.x; 
    }
}

Я не уверен, что это правильный путь или он может нежелательным образом повлиять на статус MT?

Спасибо.


person PhillipD    schedule 21.10.2013    source источник


Ответы (2)


Я предлагаю вам прочитать документацию по CURAND. тщательно и основательно.

API MT будет наиболее эффективным при использовании 256 потоков на блок и до 64 блоков для генерации чисел.

Если вам нужно больше, у вас есть множество вариантов:

  1. просто сгенерируйте больше чисел из существующего набора состояний (например, 64 блока, 256 потоков) и распределите эти числа среди потоков, которым они нужны.
  2. Используйте более одного состояния на блок (но это не позволяет вам превысить общий лимит в наборе состояний, это просто устраняет необходимость в одном блоке).
  3. Создайте несколько генераторов MT с независимыми начальными значениями (и, следовательно, независимыми наборами состояний).

В целом, я не вижу проблемы с ядром, которое вы описали, и оно примерно соответствует варианту 1 выше. Однако это не позволяет вам превышать 51200 потоков. (ваш пример имеет <<<128, 128>>>, поэтому 16384 потока)

person Robert Crovella    schedule 21.10.2013
comment
Спасибо за ваш ответ. Между тем я обнаружил, что, например. результат для потока 0 такой же, как для потока 16384 и так далее, поэтому моя идея не будет работать для многих потоков. Я буду размышлять над вашими предложениями и над тем, как применить их к моей программе. Кстати: согласно документации я делаю вывод, что не должно быть проблем с использованием XORWOW, скажем, для потоков 1e6 (т.е. статусов 1e6) одновременно. Или возможны какие-то ограничения? - person PhillipD; 22.10.2013
comment
Единственные ограничения этого типа, о которых я знаю, относятся к МП. MRG и XORWOW не должны иметь таких ограничений. - person Robert Crovella; 24.10.2013

Следуя ответу Роберта, ниже я привожу полностью проработанный пример использования Mersenne Twister cuRAND для произвольного количества потоков. Я использую первый вариант Роберта, чтобы генерировать больше чисел из существующего набора состояний и распределять эти числа среди потоков, которые в них нуждаются.

// --- Generate random numbers with cuRAND's Mersenne Twister

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

#include <cuda.h>
#include <curand_kernel.h>
/* include MTGP host helper functions */
#include <curand_mtgp32_host.h>

#define BLOCKSIZE   256
#define GRIDSIZE    64

/*******************/
/* GPU ERROR CHECK */
/*******************/
#define gpuErrchk(x) do { if((x) != cudaSuccess) { \
    printf("Error at %s:%d\n",__FILE__,__LINE__); \
    return EXIT_FAILURE;}} while(0)

#define CURAND_CALL(x) do { if((x) != CURAND_STATUS_SUCCESS) { \
    printf("Error at %s:%d\n",__FILE__,__LINE__); \
    return EXIT_FAILURE;}} while(0)

/*******************/
/* iDivUp FUNCTION */
/*******************/
__host__ __device__ int iDivUp(int a, int b) { return ((a % b) != 0) ? (a / b + 1) : (a / b); }

/*********************/
/* GENERATION KERNEL */
/*********************/
__global__ void generate_kernel(curandStateMtgp32 * __restrict__ state, float * __restrict__ result, const int N)
{
    int tid = threadIdx.x + blockIdx.x * blockDim.x;
    for (int k = tid; k < N; k += blockDim.x * gridDim.x)
        result[k] = curand_uniform(&state[blockIdx.x]);
}

/********/
/* MAIN */
/********/
int main()
{
    const int N = 217 * 123;

    // --- Allocate space for results on host
    float *hostResults = (float *)malloc(N * sizeof(float));

    // --- Allocate and initialize space for results on device 
    float *devResults; gpuErrchk(cudaMalloc(&devResults, N * sizeof(float)));
    gpuErrchk(cudaMemset(devResults, 0, N * sizeof(float)));

    // --- Setup the pseudorandom number generator
    curandStateMtgp32 *devMTGPStates; gpuErrchk(cudaMalloc(&devMTGPStates, GRIDSIZE * sizeof(curandStateMtgp32)));
    mtgp32_kernel_params *devKernelParams; gpuErrchk(cudaMalloc(&devKernelParams, sizeof(mtgp32_kernel_params)));
    CURAND_CALL(curandMakeMTGP32Constants(mtgp32dc_params_fast_11213, devKernelParams));
    //CURAND_CALL(curandMakeMTGP32KernelState(devMTGPStates, mtgp32dc_params_fast_11213, devKernelParams, GRIDSIZE, 1234));
    CURAND_CALL(curandMakeMTGP32KernelState(devMTGPStates, mtgp32dc_params_fast_11213, devKernelParams, GRIDSIZE, time(NULL)));

    // --- Generate pseudo-random sequence and copy to the host
    generate_kernel << <GRIDSIZE, BLOCKSIZE >> >(devMTGPStates, devResults, N);
    gpuErrchk(cudaPeekAtLastError());
    gpuErrchk(cudaDeviceSynchronize());
    gpuErrchk(cudaMemcpy(hostResults, devResults, N * sizeof(float), cudaMemcpyDeviceToHost));

    // --- Print results
    //for (int i = 0; i < N; i++) {
    for (int i = 0; i < 10; i++) {
        printf("%f\n", hostResults[i]);
    }

    // --- Cleanup
    gpuErrchk(cudaFree(devMTGPStates));
    gpuErrchk(cudaFree(devResults));
    free(hostResults);

    return 0;
}
person Vitality    schedule 30.03.2019