CUDA by example Chapter5 线程协作

并行线程块的分解

使用线程实现GPU上的矢量求和

在上一章代码下需要修改两个地方:

1. add<<<N,1>>>(dev_a,dev_b,dev_c) 修改为 add<<<1,N>>>(dev_a,dev_b,dev_c)
2. int tid = blockIdx.x 修改为 int tid = threadIdx.x

完整代码:

//代码5.2.1.1使用thread实现GPU上的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10

__global__ void add(int *a, int *b, int *c)
{
	int tid = threadIdx.x;//计算位于这个索引处的数据
	if (tid < N)
		c[tid] = a[tid] + b[tid];
}


int main()
{
	int a[N], b[N], c[N];
	int *dev_a, *dev_b, *dev_c;

	//在GPU上分配内存
	cudaMalloc((void **)&dev_a, N*sizeof(int));
	cudaMalloc((void **)&dev_b, N*sizeof(int));
	cudaMalloc((void **)&dev_c, N*sizeof(int));

	//在CPU上为数组'a'和‘b’赋值
	for (int i = 0; i < N; i++)
	{
		a[i] = -i;
		b[i] = i*i;
	}

	//将数组'a'和'b'复制到GPU
	cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
	cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

	add << <1, N >> >(dev_a, dev_b, dev_c);

	//将数组'c'从GPU复制到CPU
	cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

	//显示结果
	for (int i = 0; i < N; i++)
	{
		printf("%d + %d = %d\n", a[i], b[i], c[i]);
	}
	//释放在GPU上分配的内存
	cudaFree(dev_a);
	cudaFree(dev_b);
	cudaFree(dev_c);

	system("pause");
	return 0;
}

运行结果:
在这里插入图片描述

在GPU上对更长的矢量求和

原因:

  1. 在GPU中,每个Block可以拥有的thread数量存在上限,由设备属性结构中maxThreadsPerBlock指定
  2. 硬件将Block的数量限制为不超过65535个

因此,需要将Block与thread结合起来,以对更长的矢量进行求和。

在上一章的代码下需要修改以下3个地方:

1. int tid = blockIdx.x 修改为 int tid = threadIdx.x + blockIdx.x*blockDim.x
2. 假设将每个Block中线程数量固定为128,那么
	add<<<N,1>>>(dev_a,dev_b,dev_c) 修改为 add<<<(N+127)/128,128>>>(dev_a,dev_b,dev_c)
3. 由于(N+127)/128是将Block数量向上取整,所以申请的总线程数会大于等于需要的线程数,
	因此在访问输入数组和输出数组之前,必须检查线程的偏移是否位于0到N之间:
	if(tid < N)
		c[tid] = a[tid] + b[tid]

完整代码:

// 代码5.2.1.2使用block和thread实现GPU上的更大规模的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10000

__global__ void add(int *a, int *b, int *c)
{
	int tid = threadIdx.x + blockIdx.x*blockDim.x;//计算位于这个索引处的数据
	if (tid < N)
		c[tid] = a[tid] + b[tid];
}


int main()
{
	int a[N], b[N], c[N];
	int *dev_a, *dev_b, *dev_c;

	//在GPU上分配内存
	cudaMalloc((void **)&dev_a, N*sizeof(int));
	cudaMalloc((void **)&dev_b, N*sizeof(int));
	cudaMalloc((void **)&dev_c, N*sizeof(int));

	//在CPU上为数组'a'和‘b’赋值
	for (int i = 0; i < N; i++)
	{
		a[i] = -i;
		b[i] = i*i;
	}

	//将数组'a'和'b'复制到GPU
	cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
	cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

	add << <(N+127)/128, 128 >> >(dev_a, dev_b, dev_c);

	//将数组'c'从GPU复制到CPU
	cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

	//显示结果
	for (int i = 0; i < N; i++)
	{
		printf("%d + %d = %d\n", a[i], b[i], c[i]);
	}
	//释放在GPU上分配的内存
	cudaFree(dev_a);
	cudaFree(dev_b);
	cudaFree(dev_c);

	system("pause");
	return 0;
}

运行结果:
在这里插入图片描述

在GPU上对任意长度的矢量求和

原因:
GPU的Block数量由硬件限制,最大为65535;每个Block中thread的数量由设备属性结构中maxThreadsPerBlock指定,也存在一个上限,因此每个GPU可以提供的线程总数量具有上限
对任意大小的矢量进行并行计算,常用的trick就是在核函数中用一个while()循环对数据进行迭代计算,每次递增的数量 = blockDim.x * gridDim.x

在上一章的代码下需要修改以下几个地方:

1. 核函数
__global__ void add(int *a, int *b, int *c)
{
	int tid = threadIdx.x + blockIdx.x*blockDim.x;	//计算位于这个索引处的数据
	while (tid < N)
	{
		c[tid] = a[tid] + b[tid];
		tid += blockDim.x*gridDim.x;				//每次递增
	}
}

2. 以合适的固定大小的grid_size和block_size启动核函数
	add<<<128,128>>>(dev_a,dev_b,dev_c) 

完整代码:

// 代码5.2.1.3使用block和thread实现GPU上的任意规模的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10000

__global__ void add(int *a, int *b, int *c)
{
	int tid = threadIdx.x + blockIdx.x*blockDim.x;//计算位于这个索引处的数据
	while (tid < N)
	{
		c[tid] = a[tid] + b[tid];
		tid += blockDim.x*gridDim.x;
	}
		
}


int main()
{
	int a[N], b[N], c[N];
	int *dev_a, *dev_b, *dev_c;

	//在GPU上分配内存
	cudaMalloc((void **)&dev_a, N*sizeof(int));
	cudaMalloc((void **)&dev_b, N*sizeof(int));
	cudaMalloc((void **)&dev_c, N*sizeof(int));

	//在CPU上为数组'a'和‘b’赋值
	for (int i = 0; i < N; i++)
	{
		a[i] = -i;
		b[i] = i*i;
	}

	//将数组'a'和'b'复制到GPU
	cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
	cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

	add << <128, 128 >> >(dev_a, dev_b, dev_c);

	//将数组'c'从GPU复制到CPU
	cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

	//显示结果
	for (int i = 0; i < N; i++)
	{
		printf("%d + %d = %d\n", a[i], b[i], c[i]);
	}

	//验证GPU确实完成了我们要求的工作
	bool success = true;
	for (int i = 0; i < N; i++)
	{
		if ((a[i] + b[i]) != c[i])
		{
			printf("Error:%d + %d != %d\n", a[i], b[i], c[i]);
			success = false;
		}
	}
	if (success) printf("we did it!");

	//释放在GPU上分配的内存
	cudaFree(dev_a);
	cudaFree(dev_b);
	cudaFree(dev_c);

	system("pause");
	return 0;
}

运行结果:
在这里插入图片描述

在GPU上使用线程实现波纹效果

解释说明:

  1. 构建数据块结构体,用于组织CPU以及GPU端的数据
struct DataBlock
{
	unsigned char *dev_bitmap;
	CPUAnimBitmap *bitmap;//CPU中用于动画的一帧位图
};
  1. 像素映射
//将blockIdx/threadIdx映射到像素位置(二维)
int x = threadIdx.x + blockIdx.x*blockDim.x;
int y = threadIdx.y + blockIdx.y*blockDim.y;
//计算该像素(二维)的索引值
int offset = x + y*DIM;//blockDim.x*gridDim.x;

完整代码:

// 代码5.2.2使用block和thread实现波纹效果
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
#include "cpu_bitmap.h"
#include "cpu_anim.h"

#define DIM 1050

struct DataBlock
{
	unsigned char *dev_bitmap;
	CPUAnimBitmap *bitmap;//CPU中用于动画的一帧位图
};

//释放在GPU上分配的内存
void cleanup(DataBlock *d)
{
	cudaFree(d->dev_bitmap);
}


//生成图像的核函数是关键代码
//核函数需要知道当前的时间以便生成正确的帧
__global__ void kernal(unsigned char *ptr, int ticks)
{
	//将blockIdx/threadIdx映射到像素位置
	int x = threadIdx.x + blockIdx.x*blockDim.x;
	int y = threadIdx.y + blockIdx.y*blockDim.y;
	if (x < DIM && y < DIM)
	{
		int offset = x + y*DIM;//blockDim.x*gridDim.x;

		//现在计算这个位置上的值
		float fx = x - DIM / 2;
		float fy = y - DIM / 2;
		float d = sqrtf(fx*fx + fy*fy);
		unsigned char grey = (unsigned char)(128.0f + 127.0f*cos(d / 10.0f - ticks / 7.0f) / (d / 10.0f + 1.0f));

		ptr[offset * 4 + 0] = grey;//R
		ptr[offset * 4 + 1] = grey;//G
		ptr[offset * 4 + 2] = grey;//B
		ptr[offset * 4 + 3] = 255;//A
	}


}

//这里关心的是怎么生成一帧图像
void generate_frame(DataBlock *d, int ticks)
{
	dim3 blocks((DIM+15) / 16,(DIM+15) / 16);
	dim3 threads(16, 16);
	kernal << <blocks, threads >> >(d->dev_bitmap, ticks);
	cudaMemcpy(d->bitmap->get_ptr(), d->dev_bitmap, d->bitmap->image_size(), cudaMemcpyDeviceToHost);
}

int main()
{
	//data对于主机来说是全局的,如果要在kernal中使用,必须传递进去
	DataBlock data;
	CPUAnimBitmap bitmap(DIM, DIM, &data);
	data.bitmap = &bitmap;
	cudaMalloc((void **)&data.dev_bitmap, bitmap.image_size());
	bitmap.anim_and_exit((void(*)(void*, int))generate_frame, (void(*)(void*))cleanup);//动画怎么实现,这里并不关心
	//将generate_frame()传递进来,每当要生成一帧新的动画时,就要调用它一次
	//把cleanup()也传递进来,当动画程序外部中断结束时用于释放GPU内存

}

运行结果:
在这里插入图片描述

共享内存与同步

共享内存

  1. 用__shared__进行指定
__shared__ float cache[threadsPerBlock];
  1. 对于在GPU上启动的每一个线程块,CUDA C编译器都将创建共享内存中变量的一个副本,某个线程块中的每个线程都可以共享相应的共享内存变量副本,不同线程块之间的线程之间无法对共享内存变量进行共享
    在这里插入图片描述
  2. 共享内存大小一般设置成threadsPerBlock,这样Block中的每个thread都可以不受竞争的并行运行

线程同步

  1. 线程需要同步的原因
    由于引入共享内存,每个Block都拥有共享内存变量的一份copy,它们独立并行高效地运行。但是对所有Block中运行结果进行汇总的时候,必须要确定所有的Block都运行结束,这称之为线程同步否则将会产生竞态条件,使运行结果具有不确定性。
  2. 用__syncthreads()进行线程同步
    这个函数的调用将确保线程块中的每个线程都执行完__syncthreads()前面的语句后,才会执行下一条语句。当一个线程开始执行__syncthreads()之后的第一条语句时,线程块中的其他线程肯定都已经执行完了__syncthreads()

点积运算

1. 矢量点积运算举例
(x1,x2,x3,x4)·(y1,y2,y3,y4) = x1y1+x2y2+x3y3+x4y4
2. 规约运算
规约运算:对一个输入数组执行某种计算,然后产生一个更小的结果数组
并行规约运算:每个线程将cache[]中的两个值相加,然后将结果保存回cahce[].执行这种操作log2(threadsPerBlock)个步骤后,就能得到cache[]中所有值的总和。
3. 将block中所有线程计算的加和值保存到全局内存
因为一个block最终只产生一个结果,保存在cache[0]中.因此,只需要一个线程执行保存工作。这里选用了0号线程来执行它,其实每个线程都可以执行这个写入操作,书中讲到选用其他的线程来执行写入操作将会带来不必要的内存通信量
解释如下:
我没有去查英文原文,但我武断并且自信地认为,书中的意思是指如果不加if(cacheIndex==0),那么一个block中的每个线程都会将结果写入到全局内存中,显然,这么做在写入单个值时效率非常差。在本书示例中,选用任意线程来将一个block的最终结果(cache[0])写入到全局内存中所花费的代价是相同的,并不会产生额外的内存通信量。

4. 完整代码

// 代码5.3.1基于共享内存和同步的点积运算
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
//#include <cuda.h>
//#include <device_functions.h>

#define imin(a,b) (a<b?a:b)
#define sum_squares(x) (x*(x+1)*(2*x+1)/6)

const int N = 33 * 1024;
const int threadsPerBlock = 256;
const int blocksPerGrid = imin(32, (N + threadsPerBlock - 1) / threadsPerBlock);

__global__ void dot(float *a, float *b, float *c)
{
	__shared__ float cache[threadsPerBlock];
	int tid = threadIdx.x + blockIdx.x*blockDim.x;
	int cacheIndex = threadIdx.x;
	float temp = 0;
	while (tid < N)
	{
		temp += a[tid] * b[tid];
		tid += blockDim.x * gridDim.x;
	}

	//设置cache中相应位置上的值
	cache[cacheIndex] = temp;

	//对线程块中的线程进行同步
	__syncthreads();

	//规约
	//对于规约运算来说,以下代码要求threadPerBlock必须是2的整数
	int i = blockDim.x / 2;
	while (i != 0)
	{
		if (cacheIndex < i)
			cache[cacheIndex] += cache[cacheIndex + i];

		__syncthreads();//等待所有线程执行完上面的代码
		i /= 2;
	}
	//由第一个线程负责将这个block中计算的加和值保存到c中
	if (cacheIndex == 0)
	{
		c[blockIdx.x] = cache[0];
	}

}

int main()
{
	float *a, *b, c, *partial_c;
	float *dev_a, *dev_b, *dev_partial_c;
	//CPU上分配内存
	a = (float *)malloc(N*sizeof(float));
	b = (float *)malloc(N*sizeof(float));
	partial_c = (float *)malloc(blocksPerGrid*sizeof(float));

	//在GPU上分配内存
	cudaMalloc((void **)&dev_a, N*sizeof(float));
	cudaMalloc((void **)&dev_b, N*sizeof(float));
	cudaMalloc((void **)&dev_partial_c, blocksPerGrid*sizeof(float));

	//填充主机内存
	for (int i = 0; i < N; i++)
	{
		a[i] = i;
		b[i] = i * 2;
	}

	//将数组'a'和'b'复制到GPU上
	cudaMemcpy(dev_a, a, N*sizeof(float), cudaMemcpyHostToDevice);
	cudaMemcpy(dev_b, b, N*sizeof(float), cudaMemcpyHostToDevice);

	dot << <blocksPerGrid, threadsPerBlock >> >(dev_a, dev_b, dev_partial_c);

	//将数组'c'从GPU复制回CPU
	cudaMemcpy(partial_c, dev_partial_c, blocksPerGrid*sizeof(float), cudaMemcpyDeviceToHost);

	//在CPU上完成最终的求和运算
	c = 0;
	for (int i = 0; i < blocksPerGrid; i++)
	{
		c += partial_c[i];
	}

	//验证结果的正确性
	printf("Does GPU value %.6g = %.6g?\n", c, 2 * sum_squares((float)(N - 1)));

	//释放GPU上的内存
	cudaFree(dev_a);
	cudaFree(dev_b);
	cudaFree(dev_partial_c);

	//释放CPU上的内存
	delete[] a;
	delete[] b;
	delete[] partial_c;

	system("pause");
}

5. 运行结果
在这里插入图片描述
显然,计算结果通过了验证

不正确的点积运算优化

__syncthreads()函数将会等待一个block中所有的线程结束之后,才会运行它之后的内容!!!
对于规约运算:

	//规约
	//对于规约运算来说,以下代码要求threadPerBlock必须是2的整数
	int i = blockDim.x / 2;
	while (i != 0)
	{
		if (cacheIndex < i)
			cache[cacheIndex] += cache[cacheIndex + i];

		__syncthreads();//等待所有线程执行完上面的代码
		i /= 2;
	}

不能改为:

	//规约
	//对于规约运算来说,以下代码要求threadPerBlock必须是2的整数
	int i = blockDim.x / 2;
	while (i != 0)
	{
		if (cacheIndex < i)
		{
			cache[cacheIndex] += cache[cacheIndex + i];
			__syncthreads();//等待所有线程执行完上面的代码
		}
		i /= 2;
	}

这将导致线程发散,即某些线程需要执行__syncthreads()函数,某些线程不需要。显然这种线程发散,对于GPU来说是致命的,此时GPU会挂起,等待那些永远也不能有机会执行__syncthreads()函数的线程执行__syncthreads(),这是十分糟糕的。

基于共享内存的位图

// 代码5.3.3基于共享内存的位图
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
#include "cpu_bitmap.h"

#define DIM 1024
#define PI 3.1415926535897932f

__global__ void kernal(unsigned char *ptr)
{
	//将threadIdx/blockIdx映射到像素位置
	int x = threadIdx.x + blockIdx.x*blockDim.x;
	int y = threadIdx.y + blockIdx.y*blockDim.y;
	int offset = x + y*DIM;
	__shared__ float shared[16][16];

	//计算这个像素位置上的值
	const float period = 128.0f;
	shared[threadIdx.x][threadIdx.y] = 255 * (sinf(x*2.0f*PI / period) + 1.0f)*(sinf(y*2.0f*PI / period) + 1.0f) / 4.0f;
	
	__syncthreads();//这里的同步十分重要,如果没有同步,下面的写入操作 可能会在 上面的计算值操作之前完成,而这是不正确的
	
	//把值保存回像素,保留x和y的次序
	ptr[offset * 4 + 0] = 0;
	ptr[offset * 4 + 1] = shared[15 - threadIdx.x][15 - threadIdx.y];
	ptr[offset * 4 + 2] = 0;
	ptr[offset * 4 + 3] = 255;

}

int main()
{
	CPUBitmap bitmap(DIM, DIM);
	unsigned char *dev_bitmap;

	cudaMalloc((void **)&dev_bitmap, bitmap.image_size());

	dim3 grids(DIM / 16, DIM / 16);
	dim3 blocks(16, 16);
	kernal << <grids, blocks >> >(dev_bitmap);
	cudaMemcpy(bitmap.get_ptr(), dev_bitmap, bitmap.image_size(), cudaMemcpyDeviceToHost);

	bitmap.display_and_exit();
	cudaFree(dev_bitmap);
}

(添加了同步的)运行结果:
在这里插入图片描述
(未添加同步的)运行结果:
在这里插入图片描述
可以看到,那些随机出现的亮线和暗线就是写入操作在计算操作之前完成引起的,但是,共享内存的初始值不应该全部为0吗?有暗线是正常的,为什么会有亮线出现呢?出于好奇,我检查了共享内存shared[0][0]的值,没错啊,确实是0!百思不得其解,望着图书馆外校长的雕像发呆好久,然后将shared的所有初始值都打印出来:
在这里插入图片描述
看来,在GPU上申请共享内存时,它的初始值是不确定的!

基于共享内存位图的继续

为了更好的直观演示那些因为没有同步而出现的现象,我先将共享内存的初始值设置为255,然后再进行计算和赋值操作得到的结果如下:
在这里插入图片描述
这符合我们的预期。
完整代码:

/*
//代码5.2.1.1使用thread实现GPU上的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10

__global__ void add(int *a, int *b, int *c)
{
   int tid = threadIdx.x;//计算位于这个索引处的数据
   if (tid < N)
   	c[tid] = a[tid] + b[tid];
}


int main()
{
   int a[N], b[N], c[N];
   int *dev_a, *dev_b, *dev_c;

   //在GPU上分配内存
   cudaMalloc((void **)&dev_a, N*sizeof(int));
   cudaMalloc((void **)&dev_b, N*sizeof(int));
   cudaMalloc((void **)&dev_c, N*sizeof(int));

   //在CPU上为数组'a'和‘b’赋值
   for (int i = 0; i < N; i++)
   {
   	a[i] = -i;
   	b[i] = i*i;
   }

   //将数组'a'和'b'复制到GPU
   cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
   cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

   add << <1, N >> >(dev_a, dev_b, dev_c);

   //将数组'c'从GPU复制到CPU
   cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

   //显示结果
   for (int i = 0; i < N; i++)
   {
   	printf("%d + %d = %d\n", a[i], b[i], c[i]);
   }
   //释放在GPU上分配的内存
   cudaFree(dev_a);
   cudaFree(dev_b);
   cudaFree(dev_c);

   system("pause");
   return 0;
}
*/


/*
// 代码5.2.1.2使用block和thread实现GPU上的更大规模的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10000

__global__ void add(int *a, int *b, int *c)
{
   int tid = threadIdx.x + blockIdx.x*blockDim.x;//计算位于这个索引处的数据
   if (tid < N)
   	c[tid] = a[tid] + b[tid];
}


int main()
{
   int a[N], b[N], c[N];
   int *dev_a, *dev_b, *dev_c;

   //在GPU上分配内存
   cudaMalloc((void **)&dev_a, N*sizeof(int));
   cudaMalloc((void **)&dev_b, N*sizeof(int));
   cudaMalloc((void **)&dev_c, N*sizeof(int));

   //在CPU上为数组'a'和‘b’赋值
   for (int i = 0; i < N; i++)
   {
   	a[i] = -i;
   	b[i] = i*i;
   }

   //将数组'a'和'b'复制到GPU
   cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
   cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

   add << <(N+127)/128, 128 >> >(dev_a, dev_b, dev_c);

   //将数组'c'从GPU复制到CPU
   cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

   //显示结果
   for (int i = 0; i < N; i++)
   {
   	printf("%d + %d = %d\n", a[i], b[i], c[i]);
   }
   //释放在GPU上分配的内存
   cudaFree(dev_a);
   cudaFree(dev_b);
   cudaFree(dev_c);

   system("pause");
   return 0;
}
*/



/*
// 代码5.2.1.3使用block和thread实现GPU上的任意规模的矢量求和
//时间:2019.07.24
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>


#define N 10000

__global__ void add(int *a, int *b, int *c)
{
   int tid = threadIdx.x + blockIdx.x*blockDim.x;//计算位于这个索引处的数据
   while (tid < N)
   {
   	c[tid] = a[tid] + b[tid];
   	tid += blockDim.x*gridDim.x;
   }
   	
}


int main()
{
   int a[N], b[N], c[N];
   int *dev_a, *dev_b, *dev_c;

   //在GPU上分配内存
   cudaMalloc((void **)&dev_a, N*sizeof(int));
   cudaMalloc((void **)&dev_b, N*sizeof(int));
   cudaMalloc((void **)&dev_c, N*sizeof(int));

   //在CPU上为数组'a'和‘b’赋值
   for (int i = 0; i < N; i++)
   {
   	a[i] = -i;
   	b[i] = i*i;
   }

   //将数组'a'和'b'复制到GPU
   cudaMemcpy(dev_a, a, N*sizeof(int), cudaMemcpyHostToDevice);
   cudaMemcpy(dev_b, b, N*sizeof(int), cudaMemcpyHostToDevice);

   add << <128, 128 >> >(dev_a, dev_b, dev_c);

   //将数组'c'从GPU复制到CPU
   cudaMemcpy(c, dev_c, N*sizeof(int), cudaMemcpyDeviceToHost);

   //显示结果
   for (int i = 0; i < N; i++)
   {
   	printf("%d + %d = %d\n", a[i], b[i], c[i]);
   }

   //验证GPU确实完成了我们要求的工作
   bool success = true;
   for (int i = 0; i < N; i++)
   {
   	if ((a[i] + b[i]) != c[i])
   	{
   		printf("Error:%d + %d != %d\n", a[i], b[i], c[i]);
   		success = false;
   	}
   }
   if (success) printf("we did it!");

   //释放在GPU上分配的内存
   cudaFree(dev_a);
   cudaFree(dev_b);
   cudaFree(dev_c);

   system("pause");
   return 0;
}
*/



/*
// 代码5.2.2使用block和thread实现波纹效果
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
#include "cpu_bitmap.h"
#include "cpu_anim.h"

#define DIM 1050

struct DataBlock
{
   unsigned char *dev_bitmap;
   CPUAnimBitmap *bitmap;//CPU中用于动画的一帧位图
};

//释放在GPU上分配的内存
void cleanup(DataBlock *d)
{
   cudaFree(d->dev_bitmap);
}


//生成图像的核函数是关键代码
//核函数需要知道当前的时间以便生成正确的帧
__global__ void kernal(unsigned char *ptr, int ticks)
{
   //将blockIdx/threadIdx映射到像素位置
   int x = threadIdx.x + blockIdx.x*blockDim.x;
   int y = threadIdx.y + blockIdx.y*blockDim.y;
   if (x < DIM && y < DIM)
   {
   	int offset = x + y*DIM;//blockDim.x*gridDim.x;

   	//现在计算这个位置上的值
   	float fx = x - DIM / 2;
   	float fy = y - DIM / 2;
   	float d = sqrtf(fx*fx + fy*fy);
   	unsigned char grey = (unsigned char)(128.0f + 127.0f*cos(d / 10.0f - ticks / 7.0f) / (d / 10.0f + 1.0f));

   	ptr[offset * 4 + 0] = grey;//R
   	ptr[offset * 4 + 1] = grey;//G
   	ptr[offset * 4 + 2] = grey;//B
   	ptr[offset * 4 + 3] = 255;//A
   }


}

//这里关心的是怎么生成一帧图像
void generate_frame(DataBlock *d, int ticks)
{
   dim3 blocks((DIM+15) / 16,(DIM+15) / 16);
   dim3 threads(16, 16);
   kernal << <blocks, threads >> >(d->dev_bitmap, ticks);
   cudaMemcpy(d->bitmap->get_ptr(), d->dev_bitmap, d->bitmap->image_size(), cudaMemcpyDeviceToHost);
}

int main()
{
   //data对于主机来说是全局的,如果要在kernal中使用,必须传递进去
   DataBlock data;
   CPUAnimBitmap bitmap(DIM, DIM, &data);
   data.bitmap = &bitmap;
   cudaMalloc((void **)&data.dev_bitmap, bitmap.image_size());
   bitmap.anim_and_exit((void(*)(void*, int))generate_frame, (void(*)(void*))cleanup);//动画怎么实现,这里并不关心
   //将generate_frame()传递进来,每当要生成一帧新的动画时,就要调用它一次
   //把cleanup()也传递进来,当动画程序外部中断结束时用于释放GPU内存

}

*/


/*
// 代码5.3.1基于共享内存和同步的点积运算
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
//#include <cuda.h>
//#include <device_functions.h>

#define imin(a,b) (a<b?a:b)
#define sum_squares(x) (x*(x+1)*(2*x+1)/6)

const int N = 33 * 1024;
const int threadsPerBlock = 256;
const int blocksPerGrid = imin(32, (N + threadsPerBlock - 1) / threadsPerBlock);

__global__ void dot(float *a, float *b, float *c)
{
   __shared__ float cache[threadsPerBlock];
   int tid = threadIdx.x + blockIdx.x*blockDim.x;
   int cacheIndex = threadIdx.x;
   float temp = 0;
   while (tid < N)
   {
   	temp += a[tid] * b[tid];
   	tid += blockDim.x * gridDim.x;
   }

   //设置cache中相应位置上的值
   cache[cacheIndex] = temp;

   //对线程块中的线程进行同步
   __syncthreads();

   //规约
   //对于规约运算来说,以下代码要求threadPerBlock必须是2的整数
   int i = blockDim.x / 2;
   while (i != 0)
   {
   	if (cacheIndex < i)
   		cache[cacheIndex] += cache[cacheIndex + i];

   	__syncthreads();//等待所有线程执行完上面的代码
   	i /= 2;
   }
   //由第一个线程负责将这个block中计算的加和值保存到c中
   if (cacheIndex == 0)
   {
   	c[blockIdx.x] = cache[0];
   }

}

int main()
{
   float *a, *b, c, *partial_c;
   float *dev_a, *dev_b, *dev_partial_c;
   //CPU上分配内存
   a = (float *)malloc(N*sizeof(float));
   b = (float *)malloc(N*sizeof(float));
   partial_c = (float *)malloc(blocksPerGrid*sizeof(float));

   //在GPU上分配内存
   cudaMalloc((void **)&dev_a, N*sizeof(float));
   cudaMalloc((void **)&dev_b, N*sizeof(float));
   cudaMalloc((void **)&dev_partial_c, blocksPerGrid*sizeof(float));

   //填充主机内存
   for (int i = 0; i < N; i++)
   {
   	a[i] = i;
   	b[i] = i * 2;
   }

   //将数组'a'和'b'复制到GPU上
   cudaMemcpy(dev_a, a, N*sizeof(float), cudaMemcpyHostToDevice);
   cudaMemcpy(dev_b, b, N*sizeof(float), cudaMemcpyHostToDevice);

   dot << <blocksPerGrid, threadsPerBlock >> >(dev_a, dev_b, dev_partial_c);

   //将数组'c'从GPU复制回CPU
   cudaMemcpy(partial_c, dev_partial_c, blocksPerGrid*sizeof(float), cudaMemcpyDeviceToHost);

   //在CPU上完成最终的求和运算
   c = 0;
   for (int i = 0; i < blocksPerGrid; i++)
   {
   	c += partial_c[i];
   }

   //验证结果的正确性
   printf("Does GPU value %.6g = %.6g?\n", c, 2 * sum_squares((float)(N - 1)));

   //释放GPU上的内存
   cudaFree(dev_a);
   cudaFree(dev_b);
   cudaFree(dev_partial_c);

   //释放CPU上的内存
   delete[] a;
   delete[] b;
   delete[] partial_c;

   system("pause");
}
*/



/*
// 代码5.3.3基于共享内存的位图
//时间:2019.07.25
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
#include "cpu_bitmap.h"

#define DIM 1024
#define PI 3.1415926535897932f

__global__ void kernal(unsigned char *ptr)
{
//将threadIdx/blockIdx映射到像素位置
int x = threadIdx.x + blockIdx.x*blockDim.x;
int y = threadIdx.y + blockIdx.y*blockDim.y;
int offset = x + y*DIM;
__shared__ float shared[16][16];

//计算这个像素位置上的值
const float period = 128.0f;
shared[threadIdx.x][threadIdx.y] = 255 * (sinf(x*2.0f*PI / period) + 1.0f)*(sinf(y*2.0f*PI / period) + 1.0f) / 4.0f;

__syncthreads();//这里的同步十分重要,如果没有同步,下面的写入操作 可能会在 上面的计算值操作之前完成,而这是不正确的

//把值保存回像素,保留x和y的次序
ptr[offset * 4 + 0] = 0;
ptr[offset * 4 + 1] = shared[15 - threadIdx.x][15 - threadIdx.y];
ptr[offset * 4 + 2] = 0;
ptr[offset * 4 + 3] = 255;

}

int main()
{
CPUBitmap bitmap(DIM, DIM);
unsigned char *dev_bitmap;

cudaMalloc((void **)&dev_bitmap, bitmap.image_size());

dim3 grids(DIM / 16, DIM / 16);
dim3 blocks(16, 16);
kernal << <grids, blocks >> >(dev_bitmap);
cudaMemcpy(bitmap.get_ptr(), dev_bitmap, bitmap.image_size(), cudaMemcpyDeviceToHost);

bitmap.display_and_exit();
cudaFree(dev_bitmap);
}
*/






// 代码5.3.3基于共享内存的位图(添加查看共享内存初始值的代码,发现共享内存初始值为0)
//时间:2019.09.01
#include "cuda_runtime.h"
#include "device_launch_parameters.h"
#include <iostream>
#include "cpu_bitmap.h"

#define DIM 1024
#define PI 3.1415926535897932f

__global__ void kernal(unsigned char *ptr,float *dev_test)
{
   //将threadIdx/blockIdx映射到像素位置
   int x = threadIdx.x + blockIdx.x*blockDim.x;
   int y = threadIdx.y + blockIdx.y*blockDim.y;
   int offset = x + y*DIM;
   __shared__ float shared[16][16];


   dev_test[threadIdx.x * 16 + threadIdx.y] = shared[threadIdx.x][threadIdx.y];
   shared[threadIdx.x][threadIdx.y] = 255;
   __syncthreads();//当所有线程将shared初始值完全赋值为255时,才能继续执行下面的任务


   //计算这个像素位置上的值
   const float period = 128.0f;
   shared[threadIdx.x][threadIdx.y] = 255 * (sinf(x*2.0f*PI / period) + 1.0f)*(sinf(y*2.0f*PI / period) + 1.0f) / 4.0f;
   
   //__syncthreads();//这里的同步十分重要,如果没有同步,下面的写入操作 可能会在 上面的计算值操作之前完成,而这是不正确的
   
   //把值保存回像素,保留x和y的次序
   ptr[offset * 4 + 0] = 0;
   ptr[offset * 4 + 1] = shared[15 - threadIdx.x][15 - threadIdx.y];
   ptr[offset * 4 + 2] = 0;
   ptr[offset * 4 + 3] = 255;

}

int main()
{
   CPUBitmap bitmap(DIM, DIM);
   unsigned char *dev_bitmap;
   float *dev_test;//GPU共享内存初始值变量
   float test[16*16];//CPU共享内存初始值变量

   cudaMalloc((void **)&dev_bitmap, bitmap.image_size());
   cudaMalloc((void **)&dev_test, sizeof(float)*16*16);//给GPU共享内存初始值变量分配内存

   dim3 grids(DIM / 16, DIM / 16);
   dim3 blocks(16, 16);
   kernal << <grids, blocks >> >(dev_bitmap,dev_test);//核函数也做了相应改变
   cudaMemcpy(bitmap.get_ptr(), dev_bitmap, bitmap.image_size(), cudaMemcpyDeviceToHost);
   cudaMemcpy(test, dev_test, sizeof(float)*16*16, cudaMemcpyDeviceToHost);//将GPU共享内存初始值变量的值复制到CPU
   std::cout << "共享内存初始值为:" << std::endl;//输出共享内存初始值
   for (int i = 0; i < 16; i++)
   {
   	for (int j = 0; j < 16; j++)
   	{
   		std::cout << "  " << test[j+i*16];//输出共享内存初始值
   	}
   	std::cout << " " << std::endl;
   }

   cudaFree(dev_bitmap);
   cudaFree(dev_test);

   bitmap.display_and_exit();

}

注意,上述代码中将共享内存初始值赋值给全局内存、将共享内存初始化为255都是多线程的。

	dev_test[threadIdx.x * 16 + threadIdx.y] = shared[threadIdx.x][threadIdx.y];
	shared[threadIdx.x][threadIdx.y] = 255;
	__syncthreads();//当所有线程将shared初始值完全赋值为255时,才能继续执行下面的任务

如果在一个线程中,使用for循环来进行上述操作,即:

	/*
	//线程0将共享内存的初始值赋值给*dev_test
	//线程0将共享内存的初始值设置为255
	if (offset == 0)
	{
		for (int i = 0; i < 16; i++)
		{
			for (int j = 0; j < 16; j++)
			{
				dev_test[j + i * 16] = shared[i][j];
			}
		}
		for (int i = 0; i < 16; i++)
		{
			for (int j = 0; j < 16; j++)
			{
				shared[i][j] = 255;
			}
		}
	}
	__syncthreads();//所有线程必须等待线程0完成上述工作才能继续执行下面的任务
	*/

似乎看起来应该得到同样的效果,但是,亲身测试,结果并不符合预期,因此在CUDA编程中,不推荐在单个线程中使用for循环,具体原因目前还不清楚,后续会在补充说明。

评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值