Cuda编程Thrust库
傅奕
1.Thrust介绍
Thrust是一个C++模板库，为GPU计算提供高效的算法和数据结构。它是由NVIDIA开发的，旨在提供一个简单易用的接口来开发高性能GPU程序。
Thrust库基于CUDA，并且是在CUDA C++语言扩展的基础上构建的。Thrust为GPU编程提供了一种高级编程范式，使得开发人员可以使用类似于STL的算法和数据结构来加速应用程序。
Thrust库提供了一个广泛的算法库，包括排序、归约、扫描、拷贝、填充等等。这些算法都是为了在GPU上获得最佳的性能而进行了优化。此外，Thrust还提供了多种数据结构，包括向量、列表、映射和集合等等。与其他GPU编程框架相比，Thrust的优势在于它采用了C++语言的通用性，并且提供了一个高度模板化的API，使得开发人员可以轻松地将其应用于现有的代码中。同时，Thrust还提供了高度抽象化的接口，使得开发人员可以专注于算法本身，而不需要过多地关注底层的GPU架构和细节。（尤其在深度学习部署时，可以让前处理与后处理部分也运行在GPU上，对算法最终速度得到提升）
2.环境
Thrust的配置环境与CUDA编程的配置环境相同，可以参考我之前写的文章：CUDA编程（一）安装。
3.vector
3.1host_vector 和device_vector介绍
thrust提供两个矢量容器host_vector 和device_vector。
host_vector 存储在主机内存中。
device_vector在GPU设备内存中，必须运行在.cu文件，否则将无法运行。
#include <iostream>
using namespace std;
#include <thrust/host_vector.h>
#include <thrust/device_vector.h>
int main()
{
// H has storage for 4 integers
thrust::host_vector<int> H(4);
// initialize individual elements
H[0] = 14;
H[1] = 20;
H[2] = 38;
H[3] = 46;
// H.size() 显示尺寸
std::cout << "H has size " << H.size() << std::endl;
// 打印数据
for (int i = 0; i < H.size(); i++)
std::cout << "H[" << i << "] = " << H[i] << std::endl;
// 重设H尺寸
H.resize(2);
std::cout << "H now has size " << H.size() << std::endl;
// 拷贝 host_vector H 到 device_vector D，CUDA 与cpu可以互通
thrust::device_vector<int> D = H;
// 修改D数据
D[0] = 99;
D[1] = 88;
// 打印D
for (int i = 0; i < D.size(); i++)
std::cout << "D[" << i << "] = " << D[i] << std::endl;
// 函数结束时H,D会自动释放
return 0;
}这个例子说了 “=”可以用于显存上的数组与内存上的数组内容互相拷贝。要注意一个单独的元素device_vector可以使用标准括号表示来访问。但是，因为每个访问都需要呼叫cudaMemcpy，他们应该谨慎使用。稍后会给出一些更高效的方式。
3.2初始化与复制的方法
#include <thrust/host_vector.h>
#include <thrust/device_vector.h>
#include <thrust/copy.h>
#include <thrust/fill.h>
#include <thrust/sequence.h>
#include <iostream>
using namespace std;
int main(void)
{
//初始化10个1
thrust::device_vector<int> D(10, 1);
// 将索引0--6的的值变为9
thrust::fill(D.begin(), D.begin() + 7, 9);
//用D的前5个数据对H进行初始化
thrust::host_vector<int> H(D.begin(), D.begin() + 5);
// 将H变为 0, 1, 2, 3, ...
thrust::sequence(H.begin(), H.end());
// H拷贝到D的开始位置
thrust::copy(H.begin(), H.end(), D.begin());
// print D
for (int i = 0; i < D.size(); i++)
std::cout << "D[" << i << "] = " << D[i] << std::endl;
cout << "柯西的笔" << endl;
return 0;
}3.3原始指针与device_ptr之间的相互转化
#include <iostream>
using namespace std;
#include <thrust/host_vector.h>
#include <thrust/device_vector.h>
#include <thrust/device_malloc.h>
int main()
{
size_t N = 10;
// 申请gpu内存
int* raw_ptr;
cudaMalloc((void**)&raw_ptr, N * sizeof(int));
//用device_ptr封装该指针
thrust::device_ptr<int> dev_ptr(raw_ptr);
// 全部初始化为0
thrust::fill(dev_ptr, dev_ptr + N, (int)0);
//打印
for (int i = 0; i < N; i++)
std::cout << "dev_ptr[" << i << "] = " << dev_ptr[i] << std::endl;
// 将device_ptr变回普通指针
int* raw_ptr2 = thrust::raw_pointer_cast(dev_ptr);
cout << "柯西的笔" << endl;
return 0;
}4.算法
thrust提供了很多公用的并行算法，所有算法都有device和host两种实现。
4.1Transformations
下面将介绍一些算法的使用：thrust::negate<int>()取反操作，thrust::modulus<int>()相除余数。thrust还提供了许多其他算法有需要的可以读官方文档。
#include <thrust/device_vector.h>
#include <thrust/transform.h>
#include <thrust/sequence.h>
#include <thrust/copy.h>
#include <thrust/fill.h>
#include <thrust/replace.h>
#include <thrust/functional.h>
#include <iostream>
int main(void)
{
// 定义3个device_vectors
thrust::device_vector<int> X(10);
thrust::device_vector<int> Y(10);
thrust::device_vector<int> Z(10);
// X初始化0,1,2,3, ....
thrust::sequence(X.begin(), X.end());
// compute Y = -X thrust::negate<int>()取反操作
thrust::transform(X.begin(), X.end(), Y.begin(), thrust::negate<int>());
for (int i = 0; i < Y.size(); i++)
std::cout << "Y[" << i << "] = " << Y[i] << std::endl;
// Z全部用2填充
thrust::fill(Z.begin(), Z.end(), 2);
// compute Y = X mod Z thrust::modulus<int>()相除余数
thrust::transform(X.begin(), X.end(), Z.begin(), Y.begin(), thrust::modulus<int>());
for (int i = 0; i < Y.size(); i++)
std::cout << "Y[" << i << "] = " << Y[i] << std::endl;
// 将Y中所有的1替换成10
thrust::replace(Y.begin(), Y.end(), 1, 10);
for (int i = 0; i < Y.size(); i++)
std::cout << "Y[" << i << "] = " << Y[i] << std::endl;
return 0;
}实现Y = A*X+Y的两种方法
#include <thrust/device_vector.h>
#include <thrust/transform.h>
#include <thrust/fill.h>
#include <thrust/functional.h>
#include <thrust/sequence.h>
#include <iostream>
using namespace std;
//实现方式1
void saxpy_slow(float A, thrust::device_vector<float>& X, thrust::device_vector<float>& Y)
{
thrust::device_vector<float> temp(X.size());
// 用数字A填充temp
thrust::fill(temp.begin(), temp.end(), A);
// temp = A * X thrust::multiplies<float>()相乘
thrust::transform(X.begin(), X.end(), temp.begin(), temp.begin(), thrust::multiplies<float>());
// Y = A * X + Y thrust::plus<float>()相加
thrust::transform(temp.begin(), temp.end(), Y.begin(), Y.begin(), thrust::plus<float>());
}
//实现方式2
struct saxpy_functor
{
const float a;
saxpy_functor(float _a) : a(_a) {}
__host__ __device__
float operator()(const float& x, const float& y) const {
return a * x + y;
}
};
void saxpy_fast(float A, thrust::device_vector<float>& X, thrust::device_vector<float>& Y)
{
// Y = A * X + Y
thrust::transform(X.begin(), X.end(), Y.begin(), Y.begin(), saxpy_functor(A));
}
int main()
{
//定义初始数据
float A = 10;
thrust::device_vector<float> X(10);
thrust::device_vector<float> Y(10);
thrust::sequence(X.begin(), X.end());
thrust::fill(Y.begin(), Y.end(), 5);
/*saxpy_slow(A, X, Y);
for (int i = 0; i < Y.size(); i++)
std::cout << "Y[" << i << "] = " << Y[i] << std::endl;*/
saxpy_fast(A, X, Y);
for (int i = 0; i < Y.size(); i++)
std::cout << "Y[" << i << "] = " << Y[i] << std::endl;
cout << "柯西的笔" << endl;
return 0;
} saxpy_fast共有2N次读,1N写的开销
saxpy_slow有4N次读,3N次写的开销
在像这样的算法中，通常值得应用kernel fusion（将多个操作组合到单个kernel中）以最小化内存读写的开销。
4.2Reductions(归约)
归约算法使用二元操作符将一个序列的值转换为单个值，可以求和，求最大值等
#include <thrust/device_vector.h>
#include <thrust/functional.h>
#include <thrust/sequence.h>
#include <iostream>
using namespace std;
int main()
{
thrust::device_vector<int> D(10);
thrust::sequence(D.begin(), D.end());
//这里前两个参数提供归约序列,后两个提供初始值和归约运算符
int sum = thrust::reduce(D.begin(), D.end(), (int)0, thrust::plus<int>());
/*int sum = thrust::reduce(D.begin(), D.end(), (int)0);
int sum = thrust::reduce(D.begin(), D.end());*/
cout << sum << endl;
return 0;
}4.3其他算法
最大值、最大值索引、最小值、最小值索引、相同数据计数、累积和、排序、关键词排序等等。
#include <thrust/device_vector.h>
#include <thrust/transform.h>
#include <thrust/functional.h>
#include <thrust/sequence.h>
#include <thrust/count.h>
#include <thrust/extrema.h>
#include <thrust/sort.h>
#include <thrust/remove.h>
#include <thrust/unique.h>
#include <thrust/scan.h>
#include <iostream>
using namespace std;
int main()
{
thrust::device_vector<int> D(10);
thrust::sequence(D.begin(), D.end());
/**************************************/
//计算对应数据数量
D[2] = 5;
D[7] = 5;
D[9] = 5;
int result = thrust::count(D.begin(), D.end(), 5);
cout << result << endl; //有4个5
//输出最小值和最小位置
thrust::device_vector<int>::iterator minindex = thrust::min_element(D.begin(), D.end());
std::cout << "min value : " << *minindex << " at position " << (minindex - D.begin()) << std::endl;
//输出最大值和最大位置
thrust::device_vector<int>::iterator maxindex = thrust::max_element(D.begin(), D.end());
std::cout << "max value : " << *maxindex << " at position " << (maxindex - D.begin()) << std::endl;
/*********************************************/
//去除对应元素
thrust::device_vector<int> F(5, 6);
F[2] = 4;
auto NewF_End = thrust::remove(F.begin(), F.end(), 4);
for (auto it = F.begin(); it != NewF_End; it++)
std::cout << "NewF = " << *it << std::endl;
//去除相邻的重复元素
thrust::device_vector<int> C(10);
thrust::sequence(C.begin(), C.end());
C[0] = 1;
C[6] = 1;
auto NewC_End = thrust::unique(C.begin(), C.end());
for (auto it = C.begin(); it != NewC_End; it++)
std::cout << "NewC = " << *it << std::endl;
//求累积和
thrust::device_vector<int> Input(5);
thrust::device_vector<int> Output(5);
thrust::sequence(Input.begin(), Input.end());
thrust::inclusive_scan(Input.begin(), Input.end(), Output.begin());
for (int i =0;i<Output.size();i++)
std::cout << "Output[" << i << "] = " << Output[i] << std::endl;
/**********************************************/
//排序
/*thrust::sort(D.begin(), D.end());
for (int i = 0; i < D.size(); i++)
std::cout << "D[" << i << "] = " << D[i] << std::endl;*/
//稳定排序算法(从小到大
thrust::stable_sort(D.begin(), D.end());
for (int i = 0; i < D.size(); i++)
std::cout << "D[" << i << "] = " << D[i] << std::endl;
//稳定排序算法(从大到小
thrust::stable_sort(D.begin(), D.end(), thrust::greater<int>());
for (int i = 0; i < D.size(); i++)
std::cout << "D[" << i << "] = " << D[i] << std::endl;
//通过key排序
vector<int> keys = { 0,2,4,1,3 };
vector<string> values = { "0","2","4","1","3" };
/*thrust::sort_by_key(keys.begin(), keys.end(), values.begin());
for (int i = 0; i < keys.size(); i++)
std::cout << "keys[" << i << "] = " << keys[i] <<" " << "values[" << i << "] ="<<values[i]<<std::endl;*/
//稳定key排序
thrust::stable_sort_by_key(keys.begin(), keys.end(), values.begin());
for (int i = 0; i < keys.size(); i++)
std::cout << "keys[" << i << "] = " << keys[i] << " " << "values[" << i << "] =" << values[i] << std::endl;
cout << "柯西的笔" << endl;
return 0;
}5.小结
Thrust是一个易用、高性能、可扩展和可移植的GPU编程库，可以帮助程序员更轻松地利用GPU的并行性能来加速通用计算。
参考文献
https://docs.nvidia.com/cuda/thrust/
类似资料：