thrust :: max_element 비교가 느리다 cublasIsamax -보다 효율적인 구현입니까?

Question

CUDA의 배열에서 최대 값의 인덱스를 찾는 데는 빠르고 효율적인 구현이 필요합니다. 이 작업은 여러 번 수행해야합니다. 원래 cublasIsamax를 사용했지만 슬프게도 최대 절대 값의 인덱스를 반환합니다. 이는 절대 원하는 값이 아닙니다. 대신 thrust :: max_element를 사용하고 있지만 cublasIsamax와 비교할 때 속도가 느립니다. 다음과 같은 방식으로 사용합니다 :

//d_vector is a pointer on the device pointing to the beginning of the vector, containing nrElements floats. 
thrust::device_ptr<float> d_ptr = thrust::device_pointer_cast(d_vector); 
thrust::device_vector<float>::iterator d_it = thrust::max_element(d_ptr, d_ptr + nrElements); 
max_index = d_it - (thrust::device_vector<float>::iterator)d_ptr; 

벡터 범위의 요소 수는 10'000에서 20'000 사이입니다. thrust :: max_element와 cublasIsamax 사이의 속도 차이는 다소 크다. 아마도 나는 모른 채 여러 가지 메모리 트랜잭션을 수행하고있을 것입니다.

Answer 1

보다 효율적인 구현은 CUDA에 고유 한 최대 색인 축소 코드를 작성하는 것입니다. cublasIsamax이 이런 식으로 뭔가를 사용하고있을 가능성이 큽니다.

우리는 3 개 접근법을 비교할 수 있습니다

1. thrust::max_element
2. cublasIsamax
3. 사용자 정의 CUDA 커널

을 여기에 완벽하게 일을 예입니다 :

$ cat t665.cu 
#include <cublas_v2.h> 
#include <thrust/extrema.h> 
#include <thrust/device_ptr.h> 
#include <thrust/device_vector.h> 
#include <iostream> 
#include <stdlib.h> 

#define DSIZE 10000 
// nTPB should be a power-of-2 
#define nTPB 256 
#define MAX_KERNEL_BLOCKS 30 
#define MAX_BLOCKS ((DSIZE/nTPB)+1) 
#define MIN(a,b) ((a>b)?b:a) 
#define FLOAT_MIN -1.0f 

#include <time.h> 
#include <sys/time.h> 

unsigned long long dtime_usec(unsigned long long prev){ 
#define USECPSEC 1000000ULL 
    timeval tv1; 
    gettimeofday(&tv1,0); 
    return ((tv1.tv_sec * USECPSEC)+tv1.tv_usec) - prev; 
} 

__device__ volatile float blk_vals[MAX_BLOCKS]; 
__device__ volatile int blk_idxs[MAX_BLOCKS]; 
__device__ int blk_num = 0; 

template <typename T> 
__global__ void max_idx_kernel(const T *data, const int dsize, int *result){ 

    __shared__ volatile T vals[nTPB]; 
    __shared__ volatile int idxs[nTPB]; 
    __shared__ volatile int last_block; 
    int idx = threadIdx.x+blockDim.x*blockIdx.x; 
    last_block = 0; 
    T my_val = FLOAT_MIN; 
    int my_idx = -1; 
    // sweep from global memory 
    while (idx < dsize){ 
    if (data[idx] > my_val) {my_val = data[idx]; my_idx = idx;} 
    idx += blockDim.x*gridDim.x;} 
    // populate shared memory 
    vals[threadIdx.x] = my_val; 
    idxs[threadIdx.x] = my_idx; 
    __syncthreads(); 
    // sweep in shared memory 
    for (int i = (nTPB>>1); i > 0; i>>=1){ 
    if (threadIdx.x < i) 
     if (vals[threadIdx.x] < vals[threadIdx.x + i]) {vals[threadIdx.x] = vals[threadIdx.x+i]; idxs[threadIdx.x] = idxs[threadIdx.x+i]; } 
    __syncthreads();} 
    // perform block-level reduction 
    if (!threadIdx.x){ 
    blk_vals[blockIdx.x] = vals[0]; 
    blk_idxs[blockIdx.x] = idxs[0]; 
    if (atomicAdd(&blk_num, 1) == gridDim.x - 1) // then I am the last block 
     last_block = 1;} 
    __syncthreads(); 
    if (last_block){ 
    idx = threadIdx.x; 
    my_val = FLOAT_MIN; 
    my_idx = -1; 
    while (idx < gridDim.x){ 
     if (blk_vals[idx] > my_val) {my_val = blk_vals[idx]; my_idx = blk_idxs[idx]; } 
     idx += blockDim.x;} 
    // populate shared memory 
    vals[threadIdx.x] = my_val; 
    idxs[threadIdx.x] = my_idx; 
    __syncthreads(); 
    // sweep in shared memory 
    for (int i = (nTPB>>1); i > 0; i>>=1){ 
     if (threadIdx.x < i) 
     if (vals[threadIdx.x] < vals[threadIdx.x + i]) {vals[threadIdx.x] = vals[threadIdx.x+i]; idxs[threadIdx.x] = idxs[threadIdx.x+i]; } 
     __syncthreads();} 
    if (!threadIdx.x) 
     *result = idxs[0]; 
    } 
} 

int main(){ 

    int nrElements = DSIZE; 
    float *d_vector, *h_vector; 
    h_vector = new float[DSIZE]; 
    for (int i = 0; i < DSIZE; i++) h_vector[i] = rand()/(float)RAND_MAX; 
    h_vector[10] = 10; // create definite max element 
    cublasHandle_t my_handle; 
    cublasStatus_t my_status = cublasCreate(&my_handle); 
    cudaMalloc(&d_vector, DSIZE*sizeof(float)); 
    cudaMemcpy(d_vector, h_vector, DSIZE*sizeof(float), cudaMemcpyHostToDevice); 
    int max_index = 0; 
    unsigned long long dtime = dtime_usec(0); 
    //d_vector is a pointer on the device pointing to the beginning of the vector, containing nrElements floats. 
    thrust::device_ptr<float> d_ptr = thrust::device_pointer_cast(d_vector); 
    thrust::device_vector<float>::iterator d_it = thrust::max_element(d_ptr, d_ptr + nrElements); 
    max_index = d_it - (thrust::device_vector<float>::iterator)d_ptr; 
    cudaDeviceSynchronize(); 
    dtime = dtime_usec(dtime); 
    std::cout << "thrust time: " << dtime/(float)USECPSEC << " max index: " << max_index << std::endl; 
    max_index = 0; 
    dtime = dtime_usec(0); 
    my_status = cublasIsamax(my_handle, DSIZE, d_vector, 1, &max_index); 
    cudaDeviceSynchronize(); 
    dtime = dtime_usec(dtime); 
    std::cout << "cublas time: " << dtime/(float)USECPSEC << " max index: " << max_index << std::endl; 
    max_index = 0; 
    int *d_max_index; 
    cudaMalloc(&d_max_index, sizeof(int)); 
    dtime = dtime_usec(0); 
    max_idx_kernel<<<MIN(MAX_KERNEL_BLOCKS, ((DSIZE+nTPB-1)/nTPB)), nTPB>>>(d_vector, DSIZE, d_max_index); 
    cudaMemcpy(&max_index, d_max_index, sizeof(int), cudaMemcpyDeviceToHost); 
    dtime = dtime_usec(dtime); 
    std::cout << "kernel time: " << dtime/(float)USECPSEC << " max index: " << max_index << std::endl; 


    return 0; 
} 
$ nvcc -O3 -arch=sm_20 -o t665 t665.cu -lcublas 
$ ./t665 
thrust time: 0.00075 max index: 10 
cublas time: 6.3e-05 max index: 11 
kernel time: 2.5e-05 max index: 10 
$ 

참고 :

1. CUBLAS uses 1-based indexing이므로 CUBLAS는 다른 인덱스보다 1 높은 인덱스를 반환합니다.
2. CUBLAS_POINTER_MODE_DEVICE을 사용하면 CUBLAS might be quicker이지만 유효성 검사를 위해 결과를 여전히 호스트에 복사해야합니다.
3. CUBLAS : CUBLAS_POINTER_MODE_DEVICE은 비동기식이어야하므로 여기에 표시된 호스트 기반 타이밍에는 cudaDeviceSynchronize()이 바람직합니다. 경우에 따라 추력도 비동기식이 될 수 있습니다.
4. CUBLAS와 다른 방법을 비교하고 편리하게하기 위해 데이터에 음수가 아닌 값을 사용하고 있습니다. 음수 값을 사용하는 경우 FLOAT_MIN 값을 조정할 수도 있습니다.
5. 성능에 이상이있을 경우 nTPB 및 MAX_KERNEL_BLOCKS 매개 변수를 조정하여 특정 GPU에서 성능을 최대화 할 수 있는지 확인할 수 있습니다. 또한 커널 코드는 (2 개의) 스레드 차단 감소의 마지막 단계에 대해 워프 - 동기 모드로 신중하게 전환하지 않아서 테이블에 성능을 남겨 놓을 수 있습니다.
6. 스레드 블로킹 감소 커널은 블록 제거/마지막 블록 전략을 사용하여 추가 커널 시작의 오버 헤드를 줄여 최종 축소를 수행합니다.