我是CUDA的新手。我用的是几个1D的线块。只是线的一个维度,和块的一个维度。问题是,我有一个2维的数组,接下来的程序部分一定是这样的。所以,我做了cudaMemCpy将这个数组传递给一个一维数组,也就是我将在内核中使用的数组。在内核中,我使用了高度块,和宽度线程(这些都是数组的量度。depth[height][width]
).
所以现在在内核中,我已经得到了int i,这将是块的id,我想用它来获得高度。然后,j将是区块内的线程id,并以 int j = threadIdx.x + blockIdx.x * width;
我想我会得到宽度。块ID由块中的线程数加上对应的线程。像模拟一个 对于. 然后,在1D数组中,我会用j来保存我需要的结果 depth[j]
我不知道我是否误解了如何使用CUDA的块和线程,甚至在我做cudaMemCpy回主机数组的部分我做错了,但我没有得到规定的结果,我不知道如何解决这个问题。
我必须用OpenMP的模板进行paralelize。主要的问题,是解决图像的一些点,看看他们是否是曼德尔布罗特集的一部分,它应该是打印的图像显示在白色的部分(depth[j] = -1
)或其他颜色(depth[j] = z
),这只是为了解释程序的其他部分,如果有帮助的话,但我在这里不需要做任何事情,只是在CUDA部分要paralelized。
这里是库达的零件。
设备。
__global__ void BAdd(int width, int* depth, float ymax, float xmin, double k, int max_iter){
int i = blockIdx.x;
int j = threadIdx.x + blockIdx.x * width;
double y = ymax-i/k;
double x = xmin+j/k;
//mandelbrot function
int z = 1;
double temp_real, cr, p_real = x;
double ci, p_imag = y;
cr = x*x; ci = y*y; // cr+ci is the square of the module of x + i·y
while (((cr+ci)<=4) && (z < max_iter)){
temp_real = cr - ci + x;
p_imag = 2 * p_real * p_imag + y;
p_real = temp_real;
z++; // z(z) = p_real + i · p_imag
cr = p_real * p_real; ci = p_imag * p_imag;
}
if((cr+ci) <= 4){
depth[j] = -1;
}
else{
depth[j] = z;
}
}
主机。
int *device_depth;
int size = height*width*sizeof(int);;
cudaMalloc((void**)&device_depth, size);
cudaMemcpy(device_depth, depth, size, cudaMemcpyHostToDevice);
int blocks = height, threads = width;
BAdd<<<blocks, threads>>>(width, device_depth, ymax, xmin, k, MAX_ITER);
cudaMemcpy(depth, device_depth, size, cudaMemcpyDeviceToHost);
cudaFree(device_depth);
cudaDeviceSynchronize();
这就是剩下的程序,因为如果能帮到你。
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>
#include <math.h>
#include <string.h>
#define MAX_ITER 1000 // Máximo número de términos calculados en cada sucesión
// Cuanto mayor sea MAX_ITER más calidad de imagen y mayor tiempo de ejecución
#define RK 250 // RK*RK número aproximado de puntos de la imagen
这里的这部分,mandel_val函数,应该删掉,在内核中引入。
/* La función mandel_val determina si el número x+i·y
* pertenece (probablemente) o no (con seguridad) al conjunto de Mandelbrot */
int mandel_val(double x, double y, int max_iter){
int j = 1;
double temp_real, cr, p_real = x;
double ci, p_imag = y;
cr = x*x; ci = y*y; // cr+ci es el cuadrado del módulo de x + i·y
while (((cr+ci)<=4) && (j < max_iter)){
temp_real = cr - ci + x;
p_imag = 2 * p_real * p_imag + y;
p_real = temp_real;
j++; // z(j) = p_real + i · p_imag
cr = p_real * p_real; ci = p_imag * p_imag;
}
if((cr+ci) <= 4) return -1; else return j;
};
int main(int argc, char *argv[]) {
int n_imag, imag, t_imag, i, j, value, valuer, valueg, valueb;
float xmin, xmax, ymin, ymax;
double x,y,l,t;
char archivoppm[20];
FILE* filein;
// A. Abrir archivo de entrada y leer número de imágenes
if (argc == 2) {
filein = fopen(argv[1], "r");
if (filein == NULL) {
fprintf(stderr, "No puedo abrir %s\n", argv[1]);
fprintf(stderr, "uso: %s <archivo de entrada> \n",argv[0]);
exit(0);
}
fscanf(filein,"%d",&n_imag);
} else {
fprintf(stderr, "uso: %s <archivo de entrada> \n",argv[0]);
exit(0);
}
// B. Obtener una a una las imágenes y los tiempos respectivos
for (imag = 0; imag <n_imag; imag ++) {
// B1. Leer tipo de imagen (1 = rectangular, 2 = cuadrada)
fscanf(filein,"%d",&t_imag);
if (t_imag == 1) {
fscanf(filein,"%f %f %f %f %s", &xmin, &xmax, &ymin, &ymax, archivoppm);
} else if (t_imag == 2) {
fscanf(filein,"%lf %lf %lf %s", &x, &y, &l, archivoppm);
l = l/2; xmin = x-l; xmax = x+l; ymin = y-l; ymax = y+l;
} else {fprintf(stderr, "\n Error en formato del archivo de entrada\n");exit(0);}
// B2. Calcular, en pixels, la anchura, width, y la altura, height, de la imagen
double k = RK/sqrt((xmax-xmin)*(ymax-ymin));
int width = k*(xmax-xmin); // Hay conversión de double a int
int height = k*(ymax-ymin); // Íd.
// width x height = RK² aproximadamente
// k = width/(xmax-xmin) = height/(ymax-ymin);
// B3. Se abre y configura el archivo de imagen ppm
strcat(archivoppm,".ppm"); //archivoppm = archivoppm + ".ppm";
FILE *f_imag = fopen(archivoppm, "wb");
fprintf(f_imag, "P3\n%d %d\n255\n", width, height);
/* Formato P3: cada pixel se representa por sus valores de rojo, verde y azul (entre 0 y 255);
* width x height: dimensión imagen en pixels */
int depth [height][width];
// Necesitamos esta matriz para colocar secuencialmente los pixels en el archivo ppm
t = omp_get_wtime();
而这一点,我应该改成cuda。
// B4. Se calcula el punto (x,y) y se determina si incluir (x + i·y) en el conjunto de Mandelbrot
#pragma omp parallel for private(j,x,y) schedule(dynamic,height/omp_get_num_threads()/16)
for(i = 0; i < height; i++){
y = ymax - i/k;
for(j = 0; j < width; j++) {
x = xmin + j/k;
// y = ymax - i/k; se ha sacado del for j
depth[i][j] = mandel_val(x, y, MAX_ITER);
}
}
for(i = 0; i < height; i++){
for(j = 0; j < width; j++) {
printf("In i:%d, j:%d we've got: %d \n",i, j, depth[i][j]);
}
}
而在这里它完成了库达的部分
// B5. Se determina el color de cada pixel y se imprime secuencialmente en el archivo
for(i = 0; i < height; i++){
for(j = 0; j < width; j++) {
value = depth[i][j];
if(value == -1) // se supone que (x+i·y) pertenece a M --> pixel blanco
fprintf(f_imag, " 255 255 255 ");
else
// se sabe que (x+i·y) no pertenece a M --> color del pixel según número de iteraciones
{
value = ((float)value/MAX_ITER) * 16777215;
// Hacemos que value quede entre 0 y 2²⁴ y extraemos sus componentes rgb
valueg = value >>8; valueb = value % 256;
valuer = valueg >> 8; valueg = valueg % 256;
// valuer = valuer + (255-valuer)*0.4; // potenciamos los rojos (?)
fprintf(f_imag," %d %d %d ",valuer,valueg,valueb);
}
}
}
fprintf(f_imag, "\n"); fclose(f_imag);
// B6. Informamos al usuario
printf("\n En el archivo %s se ha creado una imagen de la zona explorada", archivoppm);
printf("\n\n Tiempo empleado %0.2f segundos\n\n",omp_get_wtime()-t);
} // fin de: for (imag = 0; imag <n; imag ++)
fclose(filein);
}