From 6703d0ee10aa820cb0bd2f77c18557596033d282 Mon Sep 17 00:00:00 2001 From: Peter Hill Date: Mon, 1 Dec 2025 14:21:42 +0000 Subject: [PATCH] Remove explicit use of `MPI_COMM_WORLD` Fixes #3215 --- src/invert/laplace/impls/cyclic/cyclic_laplace.cxx | 14 ++++++++------ src/invert/laplace/impls/pcr/pcr.cxx | 8 ++++---- src/invert/laplace/impls/pcr_thomas/pcr_thomas.cxx | 8 ++++---- src/mesh/interpolation/hermite_spline_xz.cxx | 10 ++-------- 4 files changed, 18 insertions(+), 22 deletions(-) diff --git a/src/invert/laplace/impls/cyclic/cyclic_laplace.cxx b/src/invert/laplace/impls/cyclic/cyclic_laplace.cxx index 5ce4e540b7..837ac3755f 100644 --- a/src/invert/laplace/impls/cyclic/cyclic_laplace.cxx +++ b/src/invert/laplace/impls/cyclic/cyclic_laplace.cxx @@ -38,8 +38,10 @@ #if not BOUT_USE_METRIC_3D #include "cyclic_laplace.hxx" -#include "bout/assert.hxx" -#include "bout/bout_types.hxx" + +#include +#include +#include #include #include #include @@ -591,21 +593,21 @@ void LaplaceCyclic ::verify_solution(const Matrix& a_ver, } if (xproc > 0) { - MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, MPI_COMM_WORLD, + MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, BoutComm::get(), &request[1]); for (int kz = 0; kz < nsys; kz++) { sbufdown[kz] = x_ver(kz, 1); } - MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, MPI_COMM_WORLD, + MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, BoutComm::get(), &request[0]); } if (xproc < nprocs - 1) { - MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, MPI_COMM_WORLD, + MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, BoutComm::get(), &request[3]); for (int kz = 0; kz < nsys; kz++) { sbufup[kz] = x_ver(kz, nx); } - MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, MPI_COMM_WORLD, + MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, BoutComm::get(), &request[2]); } diff --git a/src/invert/laplace/impls/pcr/pcr.cxx b/src/invert/laplace/impls/pcr/pcr.cxx index 48bbdbac4b..b619d80100 100644 --- a/src/invert/laplace/impls/pcr/pcr.cxx +++ b/src/invert/laplace/impls/pcr/pcr.cxx @@ -1041,21 +1041,21 @@ void LaplacePCR ::verify_solution(const Matrix& a_ver, } if (xproc > 0) { - MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, MPI_COMM_WORLD, + MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, BoutComm::get(), &request[1]); for (int kz = 0; kz < nsys; kz++) { sbufdown[kz] = x_ver(kz, 1); } - MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, MPI_COMM_WORLD, + MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, BoutComm::get(), &request[0]); } if (xproc < nprocs - 1) { - MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, MPI_COMM_WORLD, + MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, BoutComm::get(), &request[3]); for (int kz = 0; kz < nsys; kz++) { sbufup[kz] = x_ver(kz, nx); } - MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, MPI_COMM_WORLD, + MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, BoutComm::get(), &request[2]); } diff --git a/src/invert/laplace/impls/pcr_thomas/pcr_thomas.cxx b/src/invert/laplace/impls/pcr_thomas/pcr_thomas.cxx index 61c8f58694..b48b20f246 100644 --- a/src/invert/laplace/impls/pcr_thomas/pcr_thomas.cxx +++ b/src/invert/laplace/impls/pcr_thomas/pcr_thomas.cxx @@ -1140,21 +1140,21 @@ void LaplacePCR_THOMAS ::verify_solution(const Matrix& a_ver, } if (xproc > 0) { - MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, MPI_COMM_WORLD, + MPI_Irecv(&rbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 901, BoutComm::get(), &request[1]); for (int kz = 0; kz < nsys; kz++) { sbufdown[kz] = x_ver(kz, 1); } - MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, MPI_COMM_WORLD, + MPI_Isend(&sbufdown[0], nsys, MPI_DOUBLE_COMPLEX, myrank - 1, 900, BoutComm::get(), &request[0]); } if (xproc < nprocs - 1) { - MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, MPI_COMM_WORLD, + MPI_Irecv(&rbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 900, BoutComm::get(), &request[3]); for (int kz = 0; kz < nsys; kz++) { sbufup[kz] = x_ver(kz, nx); } - MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, MPI_COMM_WORLD, + MPI_Isend(&sbufup[0], nsys, MPI_DOUBLE_COMPLEX, myrank + 1, 901, BoutComm::get(), &request[2]); } diff --git a/src/mesh/interpolation/hermite_spline_xz.cxx b/src/mesh/interpolation/hermite_spline_xz.cxx index c0040d096e..5064c7d4a1 100644 --- a/src/mesh/interpolation/hermite_spline_xz.cxx +++ b/src/mesh/interpolation/hermite_spline_xz.cxx @@ -21,6 +21,7 @@ **************************************************************************/ #include "../impls/bout/boutmesh.hxx" +#include "bout/bout.hxx" #include "bout/globals.hxx" #include "bout/index_derivs_interface.hxx" #include "bout/interpolation_xz.hxx" @@ -133,16 +134,9 @@ XZHermiteSpline::XZHermiteSpline(int y_offset, Mesh* mesh) #ifdef HS_USE_PETSC petsclib = new PetscLib( &Options::root()["mesh:paralleltransform:xzinterpolation:hermitespline"]); - // MatCreate(MPI_Comm comm,Mat *A) - // MatCreate(MPI_COMM_WORLD, &petscWeights); - // MatSetSizes(petscWeights, m, m, M, M); - // PetscErrorCode MatCreateAIJ(MPI_Comm comm, PetscInt m, PetscInt n, PetscInt M, - // PetscInt N, PetscInt d_nz, const PetscInt d_nnz[], - // PetscInt o_nz, const PetscInt o_nnz[], Mat *A) - // MatSetSizes(Mat A,PetscInt m,PetscInt n,PetscInt M,PetscInt N) const int m = localmesh->LocalNx * localmesh->LocalNy * localmesh->LocalNz; const int M = m * localmesh->getNXPE() * localmesh->getNYPE(); - MatCreateAIJ(MPI_COMM_WORLD, m, m, M, M, 16, nullptr, 16, nullptr, &petscWeights); + MatCreateAIJ(BoutComm::get(), m, m, M, M, 16, nullptr, 16, nullptr, &petscWeights); #endif #endif #ifndef HS_USE_PETSC