Как я могу пропустить четвертый элемент в float4 при использовании cublas sgemv?

Часть кода, над которым я работаю, требует как можно более быстрого выполнения векторного умножения матриц, то есть с использованием оптимизированной сторонней библиотеки, такой как cublas (хотя тот же принцип применим к любому процессору).

Проблема в том, что между элементами в векторе есть своего рода шаг:

Матрица хранится в виде трехмерного массива с плавающей запятой 3Nx3N.

Вектор сохраняется как массив N 1D float4s, но должны использоваться только первые три элемента каждого float4, четвертый должен игнорироваться.

N составляет порядка миллионов.

Если бы вектор хранился как float3 вместо float4, я мог бы просто привести указатель к float, как в этом рабочем примере:

//Compile with nvcc test.cu -O3 -lcublas -o test

/*
Multiply a 3Nx3N float matrix, M,  by a vector, X, of N float3 elements

The result, Y, is a 3N float vector
-----------------------

What if X is a vector of N float4?

How can I tell cublas to skip the forth element?

*/

#include<iostream>
#include<thrust/device_vector.h>
#include<cuda_runtime.h>
#include<cublas_v2.h>

using namespace std;

int main(){

int N = 3;

thrust::device_vector<float3> X(N);

thrust::device_vector<float> Y(3*N);

for(int i=0; i<N; i++)
X[i] = make_float3(1,1,1); //make_float4(1,1,1,0); //in the case of float4 i.e., The result should be the same

thrust::device_vector<float> M(3*N*3*N, 1);cublasHandle_t handle;

cublasCreate(&handle);

float beta = 0.0f;
float alpha = 1.0f;
cublasSgemv(handle, CUBLAS_OP_T,
3*N, 3*N,
&alpha,
thrust::raw_pointer_cast(&M[0]), 3*N,
(float*) thrust::raw_pointer_cast(&X[0]), 1,
&beta,
thrust::raw_pointer_cast(&Y[0]), 1);

cout<<"Performed Y = M·X\n\tX = ";
for(int i=0; i<N; i++){
float3 Xi = X[i];
cout<<Xi.x<<" "<<Xi.y<<" "<<Xi.z<<" ";
}
cout<<"\n\tY = ";
for(int i=0; i<3*N; i++){
cout<<Y[i]<<" ";
}
cout<<endl;

return 0;
}

Но как я могу выполнить эту операцию, если вектор X хранится как float4 s?

Учитывая, что float4 * можно интерпретировать как float * с в 4 раза большим количеством элементов, вопрос может быть более общим (хотя меня интересует только случай float4);
Если между 3 «полезными» элементами есть шаг. Я хочу сказать кублам, что массив не слипается в памяти. Но что-то вроде: есть 3 элемента в начале, следующие три — это элементы «шага» после этого и т. Д.
Подобно тому, что вы можете сделать в OpenGL с объектами массива вершин.

РЕДАКТИРОВАТЬ:

Ответы показали, что наиболее жизнеспособный метод — просто скопировать выделенный массив во временный преобразованный массив float3, понятный cublas..

На данный момент есть два варианта:

1. Use cudaMemcpy2D
2. Use a thrust transformation
3. Use a custom copy kernel

Я написал этот код для проверки трех случаев:

//Compile with Compile with: nvcc test.cu -O3 -lcublas -o test
#include<iostream>
#include<thrust/device_vector.h>
#include<cuda.h>
#include<cuda_runtime.h>
#include<cublas_v2.h>

using namespace std;struct Timer{
cudaEvent_t start, stop;
float time;

void tic(){
cudaEventCreate(&start);
cudaEventCreate(&stop);
cudaEventRecord(start, 0);
}
float toc(){
cudaEventRecord(stop, 0);
cudaEventSynchronize(stop);
cudaEventElapsedTime(&time, start, stop);

cudaEventDestroy(start);
cudaEventDestroy(stop);

return time;
}

};struct copy_functor{
copy_functor(){}
__device__ float3 operator() (const float4& X4){
return make_float3(X4.x, X4.y, X4.z);
}
};__global__ void copy_kernel(const float4* __restrict__ X4, float3* __restrict__ X3, int N){
int id = blockIdx.x*blockDim.x + threadIdx.x;
if(id < N){
float4 x4 = X4[id];
X3[id] = make_float3(x4.x, x4.y, x4.z);
}
}

int main(){

int N = 1000000;
int Ntest = 1000;

Timer t;

thrust::device_vector<float3> X3(N, make_float3(0,0,0));
thrust::device_vector<float4> X4(N, make_float4(1,1,1,10));/*************************CUDAMEMCPY2D*******************/
t.tic();

for(int i= 0; i<Ntest; i++){
cudaMemcpy2DAsync(thrust::raw_pointer_cast(&X3[0]),
3*sizeof(float),
thrust::raw_pointer_cast(&X4[0]),
4*sizeof(float),
3*sizeof(float),
N,
cudaMemcpyDeviceToDevice);
cudaDeviceSynchronize();
}
printf ("Time for cudaMemcpy2DAsync: %f ms\n", t.toc()/(float)Ntest);/************************THRUST***********************/
t.tic();

for(int i= 0; i<Ntest; i++){
transform(X4.begin(), X4.end(), X3.begin(), copy_functor());
cudaDeviceSynchronize();
}

printf ("Time for thrust transformation: %f ms\n", t.toc()/(float)Ntest);

/*********************COPY KERNEL*****************************/

t.tic();
for(int i= 0; i<Ntest; i++){
copy_kernel<<< N/128 + 1, 128 >>>(thrust::raw_pointer_cast(&X4[0]),
thrust::raw_pointer_cast(&X3[0]), N);
cudaDeviceSynchronize();
}
printf ("Time for copy kernel: %f ms\n", t.toc()/(float)Ntest);return 0;
}

Обратите внимание, что я выполняю в среднем 1000 копий.

Вывод этого кода в GTX 980 следующий:

Time for cudaMemcpy2DAsync: 1.465522 ms
Time for thrust transformation: 0.178745 ms
Time for copy kernel: 0.168507 ms

cudaMemcpy2D на порядок медленнее остальных.

тяги и копирование ядра очень похожи и самый быстрый способ

Такое поведение остается с любым количеством элементов.

EDIT2:

Другие ответы предполагают, что GEMM может быть использован для передачи шага. Без необходимости временного массива.

Интерпретация матричного вектора мул. как матрица матрица мул. будет сделано так:

 cublasSgemm(handle, CUBLAS_OP_T, CUBLAS_OP_T,
3*N, 1 /*m*/, 3*N,
&alpha,
thrust::raw_pointer_cast(&M[0]), 3*N,
(float*) thrust::raw_pointer_cast(&X3[0]), 1 /*ldb*/,
&beta,
thrust::raw_pointer_cast(&Y[0]), 3*N);

Однако, на данный момент, я не знаю, как пройти X4 вместо X3. Похоже, что решение находится в параметрах m и ldb.

1

Решение

Вы можете рассматривать ваш 1-D вектор float4 как матрицу Nx3 2-D с плавающей строкой 4 и использовать cudaMemcpy2DAsync изменить шаг с 4 до 3 с

cudaMemcpy2DAsync(dst,
3*sizeof(float),
src,
4*sizeof(float),
3*sizeof(float),
N,
cudaMemcpyDeviceToDevice);

Тогда dst может рассматриваться как трехмерный вектор с плавающей запятой и передаваться в gemv() непосредственно.

Учитывая масштаб вашего Nвремя копирования не заметно по сравнению с gemv(),


РЕДАКТИРОВАТЬ

Результат теста @Apo показывает, что лучше использовать ядро ​​копии вместо cudaMemcpy2DAsync, Я был слишком ожидаем cudaMemcpy2DAsync и думал, что это будет хорошо оптимизировано и будет иметь лучшую производительность для всех случаев.

1

Другие решения

Других решений пока нет …

По вопросам рекламы ammmcru@yandex.ru
Adblock
detector