MPI 矩阵向量乘法 returns 有时正确有时奇怪的值
MPI matrix-vector-multiplication returns sometimes correct sometimes weird values
我有以下代码:
//Start MPI...
MPI_Init(&argc, &argv);
int size = atoi(argv[1]);
int delta = 10;
int rnk;
int p;
int root = 0;
MPI_Status mystatus;
MPI_Comm_rank(MPI_COMM_WORLD, &rnk);
MPI_Comm_size(MPI_COMM_WORLD, &p);
//Checking compatibility of size and number of processors
assert(size % p == 0);
//Initialize vector...
double *vector = NULL;
vector = malloc(size*sizeof(double));
double *matrix = NULL;
//Rank 0 -----------------------------------
if (rnk == 0) {
//Initialize vector...
srand(1);
for (int i = 0; i < size; i++) {
vector[i] = rand() % delta + 1;
}
printf("Initial vector:");
print_vector(vector, size);
//Initialize matrix...
matrix = malloc(size*size*sizeof(double));
srand(2);
for (int i = 0; i < (size*size); i++) {
matrix[i] = rand() % delta + 1;
}
//Print matrix...
printf("Initial matrix:");
print_flat_matrix(matrix, size);
}
//Calculating chunk_size...
int chunk_size = size/p;
//Initialize submatrix..
double *submatrix = malloc(size*chunk_size*sizeof(double));
//Initialize result vector...
double *result = malloc(chunk_size*sizeof(double));
//Broadcasting vector...
MPI_Bcast(vector, size, MPI_DOUBLE, root, MPI_COMM_WORLD);
MPI_Barrier(MPI_COMM_WORLD);
//Scattering matrix...
MPI_Scatter(matrix, (size*chunk_size), MPI_DOUBLE, submatrix, (size*chunk_size), MPI_DOUBLE, root, MPI_COMM_WORLD);
MPI_Barrier(MPI_COMM_WORLD);
printf("I am rank %d and first element of my vector is: %f and of my matrix1: %f/matrix2: %f/matrix3: %f/matrix4: %f\n", rnk, vector[0], submatrix[0], submatrix[1], submatrix[2], submatrix[3]);
//Calculating...
for (int i = 0; i < chunk_size; i++) {
for (int j = 0; j < size; j++) {
result[i] += (submatrix[(i*size)+j] * vector[j]);
printf("Rank %d; current result: %f, ", rnk, result[i]);
}
printf("\n");
printf("Rank %d; result: %f...\n", rnk, result[i]);
}
printf("Rank: %d; first result: %f\n", rnk, result[0]);
double *final_result = NULL;
//Rank 0 -----------------------------------
if (rnk == 0) {
final_result = malloc(size*sizeof(double));
}
//Gather...
MPI_Gather(result, chunk_size, MPI_DOUBLE, final_result, chunk_size, MPI_DOUBLE, root, MPI_COMM_WORLD);
//Rank 0 -----------------------------------
if (rnk == 0) {
printf("Final result:\n");
print_vector(final_result, size);
free(matrix);
free(final_result);
}
free(submatrix);
free(result);
free(vector);
MPI_Finalize();
当我 运行 程序时,它 运行 无误地完成了,但我最后打印的值并不总是正确的。有时我收到输出正确的矢量,有时部分正确,有时完全错误。错误的值要么恰好是 2 的错误值,要么是一些非常长的无用数字序列(在我看来,这必须是错误的内存访问,但我找不到任何东西而且也很奇怪,因为它有时有效)。
我也总是选择我的大小,以便它适合 mpi 创建的进程数。 mpi 在我的机器上创建了 4 个进程(测试和检查值),因此为了测试我的算法,我总是选择 4 作为大小值。更大的尺寸也会出现同样的问题。
期待您的帮助和意见,在此先感谢!
PS: 我在C
你熟悉 valgrind 吗?它会立即将您的注意力吸引到有问题的行上。
你的问题似乎是这一行:
result[i] += (submatrix[(i*size)+j] * vector[j]);
最初的结果[]是什么?它被拉下了堆。有时,如果幸运的话,它会是零。不要指望C.
的运气
初始化数组的方法有很多种。这里有一些方法,按最有可能被优化的顺序列出:
使用 calloc:
分配结果[]
double *result = calloc(chunk_size , sizeof(double));
或者,用 memset 初始化数组:
double *result = malloc(chunk_size *sizeof(double));
memset(result, 0, chunk_size *sizeof(double));
或者,可以遍历数组
for (i=0; i < chunk_size; i++)
result[i] = 0.0
我有以下代码:
//Start MPI...
MPI_Init(&argc, &argv);
int size = atoi(argv[1]);
int delta = 10;
int rnk;
int p;
int root = 0;
MPI_Status mystatus;
MPI_Comm_rank(MPI_COMM_WORLD, &rnk);
MPI_Comm_size(MPI_COMM_WORLD, &p);
//Checking compatibility of size and number of processors
assert(size % p == 0);
//Initialize vector...
double *vector = NULL;
vector = malloc(size*sizeof(double));
double *matrix = NULL;
//Rank 0 -----------------------------------
if (rnk == 0) {
//Initialize vector...
srand(1);
for (int i = 0; i < size; i++) {
vector[i] = rand() % delta + 1;
}
printf("Initial vector:");
print_vector(vector, size);
//Initialize matrix...
matrix = malloc(size*size*sizeof(double));
srand(2);
for (int i = 0; i < (size*size); i++) {
matrix[i] = rand() % delta + 1;
}
//Print matrix...
printf("Initial matrix:");
print_flat_matrix(matrix, size);
}
//Calculating chunk_size...
int chunk_size = size/p;
//Initialize submatrix..
double *submatrix = malloc(size*chunk_size*sizeof(double));
//Initialize result vector...
double *result = malloc(chunk_size*sizeof(double));
//Broadcasting vector...
MPI_Bcast(vector, size, MPI_DOUBLE, root, MPI_COMM_WORLD);
MPI_Barrier(MPI_COMM_WORLD);
//Scattering matrix...
MPI_Scatter(matrix, (size*chunk_size), MPI_DOUBLE, submatrix, (size*chunk_size), MPI_DOUBLE, root, MPI_COMM_WORLD);
MPI_Barrier(MPI_COMM_WORLD);
printf("I am rank %d and first element of my vector is: %f and of my matrix1: %f/matrix2: %f/matrix3: %f/matrix4: %f\n", rnk, vector[0], submatrix[0], submatrix[1], submatrix[2], submatrix[3]);
//Calculating...
for (int i = 0; i < chunk_size; i++) {
for (int j = 0; j < size; j++) {
result[i] += (submatrix[(i*size)+j] * vector[j]);
printf("Rank %d; current result: %f, ", rnk, result[i]);
}
printf("\n");
printf("Rank %d; result: %f...\n", rnk, result[i]);
}
printf("Rank: %d; first result: %f\n", rnk, result[0]);
double *final_result = NULL;
//Rank 0 -----------------------------------
if (rnk == 0) {
final_result = malloc(size*sizeof(double));
}
//Gather...
MPI_Gather(result, chunk_size, MPI_DOUBLE, final_result, chunk_size, MPI_DOUBLE, root, MPI_COMM_WORLD);
//Rank 0 -----------------------------------
if (rnk == 0) {
printf("Final result:\n");
print_vector(final_result, size);
free(matrix);
free(final_result);
}
free(submatrix);
free(result);
free(vector);
MPI_Finalize();
当我 运行 程序时,它 运行 无误地完成了,但我最后打印的值并不总是正确的。有时我收到输出正确的矢量,有时部分正确,有时完全错误。错误的值要么恰好是 2 的错误值,要么是一些非常长的无用数字序列(在我看来,这必须是错误的内存访问,但我找不到任何东西而且也很奇怪,因为它有时有效)。
我也总是选择我的大小,以便它适合 mpi 创建的进程数。 mpi 在我的机器上创建了 4 个进程(测试和检查值),因此为了测试我的算法,我总是选择 4 作为大小值。更大的尺寸也会出现同样的问题。
期待您的帮助和意见,在此先感谢!
PS: 我在C
你熟悉 valgrind 吗?它会立即将您的注意力吸引到有问题的行上。
你的问题似乎是这一行:
result[i] += (submatrix[(i*size)+j] * vector[j]);
最初的结果[]是什么?它被拉下了堆。有时,如果幸运的话,它会是零。不要指望C.
的运气初始化数组的方法有很多种。这里有一些方法,按最有可能被优化的顺序列出:
使用 calloc:
分配结果[]double *result = calloc(chunk_size , sizeof(double));
或者,用 memset 初始化数组:
double *result = malloc(chunk_size *sizeof(double));
memset(result, 0, chunk_size *sizeof(double));
或者,可以遍历数组
for (i=0; i < chunk_size; i++)
result[i] = 0.0