Программирование на видеокартах GPGPU (1184391), страница 5
Текст из файла (страница 5)
Осуществляется это путём поочерёдного сравненияи перестановки сначала соседних пар (с чередованием порядка следования), затем четвёрок(сначала — удалённых пар, потом — соседних), затем аналогичным образом восьмёрок и такдалее.Пары, используемые для сравнения и перестановки на каждом шаге, показаны на рисунке:Иллюстрация действий в Bitonic Sort (http://en.wikipedia.org/wiki/Bitonic_sorter)Оказывается, все эти действия укладываются в однотипную последовательность шагов,описываемую двумя вложенными циклами (см.
код далее).Параметры для Bitonic Sort с помощью CUDA:/* Каждая нить "работает" с одной величиной */#define THREADS 512 // 2^9#define BLOCKS 32768 // 2^15#define NUM_VALS THREADS*BLOCKS17Функция-ядро:__global__ void bitonic_sort_step(float *dev_values, int j, int k){unsigned int i, ixj; /* Сортируемые "партнёры": i and ixj */i = threadIdx.x + blockDim.x * blockIdx.x;ixj = i^j;if ((ixj)>i) {if ((i&k)==0) {/* Сортировка по возрастанию */if (dev_values[i]>dev_values[ixj]) {/* Обмен i и ixj */float temp = dev_values[i];dev_values[i] = dev_values[ixj];dev_values[ixj] = temp;}}if ((i&k)!=0) {/* Сортировка по убыванию */if (dev_values[i]<dev_values[ixj]) {/* Обмен i и ixj */float temp = dev_values[i];dev_values[i] = dev_values[ixj];dev_values[ixj] = temp;}}}}А вот сама процедура сортировки, включающая формирование битонической последовательностииз произвольной заданной:/*** Bitonic Sort - сортировка происходит "на месте"!*/void bitonic_sort(float *values){float *dev_values;size_t size = NUM_VALS * sizeof(float);cudaMalloc((void**) &dev_values, size);cudaMemcpy(dev_values, values, size, cudaMemcpyHostToDevice);dim3 blocks(BLOCKS,1);dim3 threads(THREADS,1);/* Number of blocks/* Number of threads*/*/int j, k;/* k - "длина" первого шага сравнения-обмена; он удваивается */for (k = 2; k <= NUM_VALS; k <<= 1) {18/* Перестановки должны продолжаться пока шаг не уменьшился */for (j=k>>1; j>0; j=j>>1) {bitonic_sort_step<<<blocks, threads>>>(dev_values, j, k);}}cudaMemcpy(values, dev_values, size, cudaMemcpyDeviceToHost);cudaFree(dev_values);}Выделяется память для сортируемых величин на видеокарте (GPU) и они копируются туда.
Затемнесколько раз вызывается функция-ядро для различных параметров j, k. И всё оказываетсяотсортированным, причём прямо "на месте"! Результат копируется с видеокарты на компьютер,выделенная ранее память на карте — освобождается."Косвенное" использование возможностей GPUНе всегда у программирующего есть возможность (и желание) изучать "прямое"программирование графической карты. Это может быть связано с изменчивостью архитектуры ичастой сменой поколений графических карт, завязанностью такого выбора на одного конкретногопроизводителя карт, сложностью самой архитектуры.
К счастью, довольно часто можно неуглубляться в изучение "железок", а воспользоваться тем, что наиболее употребительныеалгоритмы уже реализованы в рамках либо поставляемых с CUDA, либо каких-то стороннихбиблиотек. В состав CUDA входят такие специализированные библиотеки, как Thrust, cuBLAS,cuFFT, cuRAND и т.п. Библиотека Thrust — это средства работы с контейнером vector награфической карте (GPU) в стиле библиотеки STL. Здесь можно не отвлекаться навыделение/освобождение памяти на компьютере (CPU) или на графической карте (GPU), а простообъявлять необходимые вектора, причём освобождаться они будут автоматически.Выделение памяти для вектора в памяти компьютера:thrust::host_vector<int> h_vec(16*1024*1024);Генерация случайных значений (в стиле библиотеки STL):thrust::generate(h_vec.begin(), h_vec.end(), rand);Создание места для вектора на графической карте и его копирование:thrust::device_vector<int> d_vec = h_vec;Сортировка (на видеокарте; прямо в месте расположения вектора):thrust::sort(d_vec.begin(), d_vec.end());Обратное копирование результата с видеокарты в память компьютера:thrust::copy(d_vec.begin(), d_vec.end(), h_vec.begin());Видно, что здесь используется явное указание пространства имён, поскольку названия совпадаютс аналогичными названиями из библиотеки STL.19Ещё один пример: умножение вектора на число и суммирование с другим вектором.В стиле CUDA C:__global__void saxpy_kernel(int n, float a, float *x, float *y){const int i = blockDim.x*blockIdx.x+threadIdx.x;if (i < n)y[i] = a*x[i]+y[i];}void saxpy(int n, float a, float *x, float *y){// параметры запуска ядраint block_size = 256;int grid_size = (n + block_size-1)/block_size;// запуск ядра saxpy_kernelsaxpy_kernel<<<grid_size,block_size>>>(n, a, x, y);}В стиле Thrust:struct saxpy_functor{const float a;saxpy_functor(float _a) : a(_a) {}__host__ __device__float operator()(float x, float y) { return a*x+y; }};void saxpy(float a, device_vector<float>& x, device_vector<float>& y){// определение функтораsaxpy_functor func(a);// вызов преобразованияtransform(x.begin(), x.end(), y.begin(), y.begin(), func);}Напомним, что функтор — это переменная, к которой можно приписать круглые скобки сосписком параметров в них, что приведёт к вызову оператора, обозначаемого круглыми скобками(operator()).
Может быть также именем класса (или структуры). Если функтор — шаблонный,то имя шаблона с типом тоже будет функтором.Обратите внимание, что оператор "круглые скобки" помечен спецификаторами и __host__, и__device__, что означает, что он может применяться как в host-векторах, так и в deviceвекторах.Вывод результирующего вектора V (на компьютере) легко делается с помощью STL:std::copy(V.begin(), V.end(), std::ostream_iterator<float>(std::cout," "));std::cout<<std::endl;205Библиотеки, использующие CUDA GPUПопытаемся вкратце охарактеризовать вычислительные библиотеки, использующие CUDA-карты.Некоторые из таких библиотек в настоящее время уже являются частью CUDA, поэтому имеетсмысл сначала поговорить о них. Разберёмся сначала с Thrust, cuBLAS и cuSPARSE.Thrust: состав и возможностиБиблиотека Thrust является библиотекой шаблонов C++ для CUDA, аналогичной STL, но покаболее простой.
Она содержит только один вид контейнеров — векторы, однако два типа их:thrust::host_vector<T> и thrust::device_vector<T>. Первый (как следует изназвания) предназначен для работы в памяти компьютера, второй — на графической карте.Подобно вектору в STL — это общие контейнеры, хранящие любые типы данных и избавляющиепрограммиста от забот с выделением/освобождением памяти, а также упрощающие обменданными между CPU и GPU.С помощью Thrust описываются собственно вычисления, а не то, как будут храниться данные илипроизводиться эти вычисления.
В библиотеке обеспечивается абстрактный интерфейс кфундаментальным параллельным алгоритмам, таким как сортировка (thrust::sort()),редуцирование (thrust::reduce()) и др. Широко используется принятое в STL указаниедиапазона как пары итераторов. Сами итераторы могут использоваться и как указатели, в томчисле на память в GPU (например, для передачи в функцию-ядро):// выделяем память на графическом устройствеthrust::device_vector<int> d_vec(M);// получаем указатель на этот вектор в памяти GPUint * ptr = thrust::raw_pointer_cast(&d_vec[0]);// используем указатель в функции-ядреmy_kernel<<<N/256, 256>>>(N, ptr);// Операция разыменования в памяти CPU не имеет смысла!Видно, что библиотека использует собственное пространство имён (thrust), что позволяетизбежать коллизий с аналогичными именами STL (например, thrust::sort() иstd::sort()).Существуют в библиотеке также так называемые fancy iterators ("воображаемые" итераторы),ведущие себя в алгоритмах как "настоящие" итераторы:constant_iterator(подобен итератору в бесконечном массиве, заполненном одной величиной)counting_iterator(подобен итератору в бесконечном массиве, заполненном последовательными величинами)transform_iterator(производит итерирование по последовательности, преобразованной с помощью заданной функции)zip_iterator(производит итерирование как бы по массиву структур, хотя принимает массивы отдельных величин,составляющих структуру)permutation_iterator21Для освоения этих возможностей полезно поэкспериментировать с некоторыми примерами,найденными в Сети.В интересной статье "odeint v2 — Solving ordinary differential equations in C++"(http://www.codeproject.com/Articles/268589/) есть пример решения дифференциальныхуравненийвподобномстиле(сиспользованиемthrust::device_vector,thrust::for_each, thrust::get<i>, thrust::make_zip_iterator).
Правда, полноготекста программы там вроде бы нет.Пример monte_carlo.cu — это известный способ оценивания константы π методом МонтеКарло. В нём используются формирование равномерно распределённой случайной величины(thrust::uniform_real_distribution<float>),редуцирующеепреобразованиеthrust::transform_reduce(), а также упомянутый выше "считающий" итераторthrust::counting_iterator<int>().Интересный пример argmin_row_space.cu использования изначально векторной природыопераций в Thrust для работы с матрицей (отыскание положения минимальных элементов в еёстроках) приведён в заметке:http://peterwittek.com/2013/04/argmin-on-the-rows-of-a-matrix-with-thrust/.Тамиспользуютсяthrust::device_vector,thrust::make_zip_iterator(),thrust::counting_iterator(), thrust::make_transform_iterator().Стоит попробовать приводимый в статье "A Brief Test on the Code Efficiency of CUDA andThrust"(http://www.codeproject.com/Articles/83757/)полезныйпримерthrustExample.cu (он располагается в прилагаемом к статье архиве с исходным кодом:http://www.codeproject.com/KB/Parallel_Programming/test-on-thrustefficiency/thrustExample.zip).Там сравниваются скорости вычисления суммы квадратов вектора просто в CUDA, с Thrust, а такжена CPU.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
