После написания последовательной программы мне нужно распараллелить ее.Вот небольшая часть, которая почему-то не работает.Когда N > 64 и 4 потока, программа начинает выдавать ошибку сегментации.И с 2-мя потоками все отлично работает.Я пытался установить переменную окружения KMP_STACK_SIZE = 128m, но это мне не помогло.В чем может быть проблема?
#include <omp.h>
void setMatrix(double *matrix, int size) {
for (int i = 0; i < size; i++)
for (int j = 0; j < size; j++) {
if (i == j) {
matrix[i * size + j] = 2;
} else
matrix[i * size + j] = 1;
}
}
void setVector(double *vector, int size, int value) {
for (int i = 0; i < size; i++) {
vector[i] = value;
}
}
void clearVec(double *vector, int size) {
for (int i = 0; i < size; i++) {
vector[i] = 0;
}
}
void mulMatrAndVec(double *result, const double *matrix, const double *vector, int size) {
clearVec(result, size);
for (int i = 0; i < size; i++) {
for (int j = 0; j < size; j++) {
result[i] += matrix[i * size + j] * vector[j];
}
}
}
void subVectors(double *result, const double *vector1, const double *vector2, int size) {
for (int i = 0; i < size; i++) {
result[i] = vector1[i] - vector2[i];
}
}
void mulMatrAndVecMP(double *result, const double *matrix, const double *vector, int place, int blockSize, int size) {
for (int i = place; i < (place + blockSize); i++) {
for (int j = 0; j < size; j++) {
result[i] += matrix[i * size + j] * vector[j];
}
}
}
int main(int argc, char *argv[]) {
int N = 68;
double *A = (double *) malloc(N * N * sizeof(double));
double *x = (double *) malloc(N * sizeof(double));
double *b = (double *) malloc(N * sizeof(double));
double *u = (double *) malloc(N * sizeof(double));
double *r1 = (double *) malloc(N * sizeof(double));
double *r2 = (double *) malloc(N * sizeof(double));
double *z = (double *) malloc(N * sizeof(double));
double *vec1 = (double *) malloc(N*N * sizeof(double));
double *vec2 = (double *) malloc(N * sizeof(double));
double a = 0;
double bt = 0;
int threadNum, threadCount;
setMatrix(A, N);
setVector(x, N, 0);
setVector(b, N, N + 1);
mulMatrAndVec(vec1, A, x, N);
subVectors(r1, b, vec1, N);
clearVec(vec1, N);
memcpy(z, r1, N * sizeof(double));
memcpy(r2, r1, N * sizeof(double));
omp_set_num_threads(4);
threadCount = omp_get_num_threads();
#pragma omp parallel private(threadNum) shared(threadCount, vec1, A, z)
{
threadNum = omp_get_thread_num();
mulMatrAndVecMP(vec1, A, z, (threadNum * N) / threadCount, N / threadCount, N);
}
free(A);
free(x);
free(b);
free(r1);
free(r2);
free(z);
free(vec1);
free(vec2);
free(u);
return 0;
}