diff --git a/PoissonSolver/2D-omp/src/solver.c b/PoissonSolver/2D-omp/src/solver.c index d4b52f9..b306265 100644 --- a/PoissonSolver/2D-omp/src/solver.c +++ b/PoissonSolver/2D-omp/src/solver.c @@ -62,30 +62,53 @@ void initSolver(Solver* solver, Parameter* params, int problem) solver->p = allocate(64, bytesize); solver->rhs = allocate(64, bytesize); - double dx = solver->dx; - double dy = solver->dy; - double* p = solver->p; - double* rhs = solver->rhs; -#pragma omp parallel for collapse(2) - for (int i = 0; i < imax + 2; i++) { - for (int j = 0; j < jmax + 2; j++) { - // P(i, j) = sin(2.0 * PI * i * dx * 2.0) + sin(2.0 * PI * j * dy * 2.0); - P(i, j) = 0; - } + double dx = solver->dx; + double dy = solver->dy; + double* p = solver->p; + double* rhs = solver->rhs; + int dim[2] = { 0 }; + int num_threads = 1; +#pragma omp parallel + { +#pragma omp critical + num_threads = omp_get_num_threads(); } + omp_create_dim(num_threads, dim); + printf("%d: { %d, %d}\n", num_threads, dim[0], dim[1]); +#pragma omp parallel + { + int jsw, isw; + double local_res = 0.0; + int li_start = get_dim_start(get_x_choord(omp_get_thread_num(), dim), + dim[0], + solver->imax); + int lj_start = get_dim_start(get_y_choord(omp_get_thread_num(), dim), + dim[1], + solver->jmax); + int limax = li_start + distribute_dim(get_x_choord(omp_get_thread_num(), dim), + dim[0], + solver->imax); + int ljmax = lj_start + distribute_dim(get_y_choord(omp_get_thread_num(), dim), + dim[1], + solver->jmax); - if (problem == 2) { -#pragma omp parallel for collapse(2) - for (int i = 0; i < imax + 2; i++) { - for (int j = 0; j < jmax + 2; j++) { - RHS(i, j) = sin(2.0 * PI * i * dx); + for (int i = li_start; i < limax + 2; i++) { + for (int j = lj_start; j < ljmax + 2; j++) { + P(i, j) = sin(2.0 * PI * i * dx * 2.0) + sin(2.0 * PI * j * dy * 2.0); } } - } else { -#pragma omp parallel for collapse(2) - for (int i = 0; i < imax + 2; i++) { - for (int j = 0; j < jmax + 2; j++) { - RHS(i, j) = 0.0; + + if (problem == 2) { + for (int i = li_start; i < limax + 2; i++) { + for (int j = lj_start; j < ljmax + 2; j++) { + RHS(i, j) = sin(2.0 * PI * i * dx); + } + } + } else { + for (int i = li_start; i < limax + 2; i++) { + for (int j = lj_start; j < ljmax + 2; j++) { + RHS(i, j) = 0.0; + } } } } @@ -117,9 +140,9 @@ void solveRB(Solver* solver) omp_create_dim(num_threads, dim); printf("%d: { %d, %d}\n", num_threads, dim[0], dim[1]); while ((res >= epssq) && (it < itermax)) { - res = 0.0; #pragma omp parallel { + res = 0.0; int jsw, isw; double local_res = 0.0; int li_start = get_dim_start(get_x_choord(omp_get_thread_num(), dim),