diff --git a/src/fluid/boundary/axis.cpp b/src/fluid/boundary/axis.cpp index 47f0c295..a6762193 100644 --- a/src/fluid/boundary/axis.cpp +++ b/src/fluid/boundary/axis.cpp @@ -394,11 +394,14 @@ void Axis::ExchangeMPI(int side) { // Load the buffers with data int ibeg,iend,jbeg,jend,kbeg,kend,offset; int nx,ny,nz; - auto bufferSend = this->bufferSend; + Buffer bufferSend = this->bufferSend; IdefixArray1D map = this->mapVars; IdefixArray4D Vc = this->Vc; IdefixArray4D Vs = this->Vs; + //reset pointer + bufferSend.ResetPointer(); + // If MPI Persistent, start receiving even before the buffers are filled MPI_Status sendStatus; @@ -419,54 +422,56 @@ void Axis::ExchangeMPI(int side) { kbeg = data->beg[KDIR]; kend = data->end[KDIR]; nz = kend - kbeg; + + // Create the base box to be patched by the ops + BoundingBox baseBox; + baseBox[IDIR][0] = ibeg; + baseBox[IDIR][1] = iend; + baseBox[JDIR][0] = jbeg; + baseBox[JDIR][1] = jend; + baseBox[KDIR][0] = kbeg; + baseBox[KDIR][1] = kend; + if(side==left) { - idefix_for("LoadBufferX2Vc",0,mapNVars,kbeg,kend,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int n, int k, int j, int i) { - bufferSend(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + n*nx*ny*nz ) = - Vc(map(n),k,j+ny,i); - } - ); + //shift by NY + BoundingBox sendBoxVc = baseBox; + sendBoxVc[JDIR][0] += ny; + sendBoxVc[JDIR][1] += ny; + bufferSend.Pack(Vc, map, sendBoxVc); + if (haveMHD) { - int VsIndex = mapNVars*nx*ny*nz; - idefix_for("LoadBufferX2IDIR",kbeg,kend,jbeg,jend,ibeg,iend+1, - KOKKOS_LAMBDA (int k, int j, int i) { - bufferSend(i + (j-jbeg)*(nx+1) + (k-kbeg)*(nx+1)*ny + VsIndex ) = - Vs(IDIR,k,j+ny,i); - } - ); - VsIndex = mapNVars*nx*ny*nz + (nx+1)*ny*nz; - idefix_for("LoadBufferX2KDIR",kbeg,kend+1,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int k, int j, int i) { - bufferSend(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + VsIndex ) = - Vs(KDIR,k,j+ny,i); - } - ); + //extend by one the end on idir + BoundingBox sendBoxVsIdir = sendBoxVc; + sendBoxVsIdir[IDIR][1] += 1; + bufferSend.Pack(Vs, IDIR, sendBoxVsIdir); + + //extend by one the end on kdir + BoundingBox sendBoxVsKdir = sendBoxVc; + sendBoxVsKdir[KDIR][1] += 1; + bufferSend.Pack(Vs, KDIR, sendBoxVsKdir); } // MHD } else if(side==right) { - idefix_for("LoadBufferX2Vc",0,mapNVars,kbeg,kend,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int n, int k, int j, int i) { - bufferSend(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + n*nx*ny*nz ) = - Vc(map(n),k,j+offset-ny,i); - } - ); + //shift by offset - NY (to select the right part) + BoundingBox sendBoxVc = baseBox; + sendBoxVc[JDIR][0] += offset - ny; + sendBoxVc[JDIR][1] += offset - ny; + bufferSend.Pack(Vc, map, sendBoxVc); // Load face-centered field in the buffer if (haveMHD) { - int VsIndex = mapNVars*nx*ny*nz; - idefix_for("LoadBufferX2IDIR",kbeg,kend,jbeg,jend,ibeg,iend+1, - KOKKOS_LAMBDA (int k, int j, int i) { - bufferSend(i + (j-jbeg)*(nx+1) + (k-kbeg)*(nx+1)*ny + VsIndex ) = - Vs(IDIR,k,j+offset-ny,i); - } - ); - VsIndex = mapNVars*nx*ny*nz + (nx+1)*ny*nz; - - idefix_for("LoadBufferX2KDIR",kbeg,kend+1,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int k, int j, int i) { - bufferSend(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + VsIndex ) = - Vs(KDIR,k,j+offset-ny,i); - } - ); + //extend by one the end on idir + take the right part or the mesh on jdir + BoundingBox sendBoxVsIdir = baseBox; + sendBoxVsIdir[IDIR][1] += 1; + sendBoxVsIdir[JDIR][0] += offset - ny; + sendBoxVsIdir[JDIR][1] += offset - ny; + bufferSend.Pack(Vs, IDIR, sendBoxVsIdir); + + //extend by one the end on kdir + take the right part or the mesh on jdir + BoundingBox sendBoxVsKdir = baseBox; + sendBoxVsKdir[KDIR][1] += 1; + sendBoxVsKdir[JDIR][0] += offset - ny; + sendBoxVsKdir[JDIR][1] += offset - ny; + bufferSend.Pack(Vs, KDIR, sendBoxVsKdir); } // MHD } // side==right @@ -478,63 +483,55 @@ void Axis::ExchangeMPI(int side) { idfx::mpiCallsTimer += MPI_Wtime() - tStart; // Unpack - auto bufferRecv=this->bufferRecv; + Buffer bufferRecv=this->bufferRecv; + bufferRecv.ResetPointer(); auto sVc = this->symmetryVc; if(side==left) { - idefix_for("StoreBufferAxis",0,mapNVars,kbeg,kend,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int n, int k, int j, int i) { - Vc(map(n),k,jend-(j-jbeg)-1,i) = - sVc(map(n))*bufferRecv(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + n*nx*ny*nz ); - } - ); + //unpack Vc + BoundingBox recvBoxVc = baseBox; + bufferRecv.UnpackJDirSymetric(Vc, map, sVc, recvBoxVc); // Load face-centered field in the buffer if (haveMHD) { int VsIndex = mapNVars*nx*ny*nz; auto sVs = this->symmetryVs; - idefix_for("StoreBufferX2IDIR",kbeg,kend,jbeg,jend,ibeg,iend+1, - KOKKOS_LAMBDA (int k, int j, int i) { - Vs(IDIR,k,jend-(j-jbeg)-1,i) = - sVs(IDIR)*bufferRecv(i + (j-jbeg)*(nx+1) + (k-kbeg)*(nx+1)*ny + VsIndex ); - } - ); - VsIndex = mapNVars*nx*ny*nz + (nx+1)*ny*nz; + //unpack Vs face-centered + BoundingBox recvBoxVsIdir = baseBox; + recvBoxVsIdir[IDIR][1] += 1; + bufferRecv.UnpackJDirSymetric(Vs, IDIR, sVs(IDIR), recvBoxVsIdir); - idefix_for("StoreBufferX2KDIR",kbeg,kend+1,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA ( int k, int j, int i) { - Vs(KDIR,k,jend-(j-jbeg)-1,i) = - sVs(KDIR)*bufferRecv(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + VsIndex ); - } - ); + //unpack Vs face-centered + BoundingBox recvBoxVsKdir = baseBox; + recvBoxVsKdir[KDIR][1] += 1; + bufferRecv.UnpackJDirSymetric(Vs, KDIR, sVs(KDIR), recvBoxVsKdir); } } else if(side==right) { - idefix_for("StoreBufferAxis",0,mapNVars,kbeg,kend,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA (int n, int k, int j, int i) { - Vc(map(n),k,jend-(j-jbeg)-1+offset,i) = - sVc(map(n))*bufferRecv(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + n*nx*ny*nz ); - } - ); + //unpack Vc on right part + BoundingBox recvBoxVc = baseBox; + recvBoxVc[JDIR][0] += offset; + recvBoxVc[JDIR][1] += offset; + bufferRecv.UnpackJDirSymetric(Vc, map, sVc, recvBoxVc); // Load face-centered field in the buffer if (haveMHD) { int VsIndex = mapNVars*nx*ny*nz; auto sVs = this->symmetryVs; - idefix_for("StoreBufferX2IDIR",kbeg,kend,jbeg,jend,ibeg,iend+1, - KOKKOS_LAMBDA (int k, int j, int i) { - Vs(IDIR,k,jend-(j-jbeg)-1+offset,i) = - sVs(IDIR)*bufferRecv(i + (j-jbeg)*(nx+1) + (k-kbeg)*(nx+1)*ny + VsIndex ); - } - ); - VsIndex = mapNVars*nx*ny*nz + (nx+1)*ny*nz; - idefix_for("StoreBufferX2KDIR",kbeg,kend+1,jbeg,jend,ibeg,iend, - KOKKOS_LAMBDA ( int k, int j, int i) { - Vs(KDIR,k,jend-(j-jbeg)-1+offset,i) = - sVs(KDIR)*bufferRecv(i + (j-jbeg)*nx + (k-kbeg)*nx*ny + VsIndex ); - } - ); + //unpack Vs face-centered on right part + BoundingBox recvBoxVsIdir = baseBox; + recvBoxVsIdir[IDIR][1] += 1; + recvBoxVsIdir[JDIR][0] += offset; + recvBoxVsIdir[JDIR][1] += offset; + bufferRecv.UnpackJDirSymetric(Vs, IDIR, sVs(IDIR), recvBoxVsIdir); + + //unpack Vs face-centered on right part + BoundingBox recvBoxVsKdir = baseBox; + recvBoxVsKdir[KDIR][1] += 1; + recvBoxVsKdir[JDIR][0] += offset; + recvBoxVsKdir[JDIR][1] += offset; + bufferRecv.UnpackJDirSymetric(Vs, KDIR, sVs(KDIR), recvBoxVsKdir); } // MHD } @@ -594,8 +591,9 @@ void Axis::InitMPI() { #endif // DIMENSIONS } - this->bufferRecv = IdefixArray1D("bufferRecvAxis", bufferSize); - this->bufferSend = IdefixArray1D("bufferSendAxis", bufferSize); + //build buffers + this->bufferRecv = Buffer(bufferSize); + this->bufferSend = Buffer(bufferSize); // init persistent communications // We receive from procRecv, and we send to procSend @@ -605,10 +603,10 @@ void Axis::InitMPI() { MPI_SAFE_CALL(MPI_Cart_shift(data->mygrid->AxisComm,0,data->mygrid->nproc[KDIR]/2, &procRecv,&procSend )); - MPI_SAFE_CALL(MPI_Send_init(bufferSend.data(), bufferSize, realMPI, procSend, + MPI_SAFE_CALL(MPI_Send_init(bufferSend.data(), bufferSend.Size(), realMPI, procSend, 650, data->mygrid->AxisComm, &sendRequest)); - MPI_SAFE_CALL(MPI_Recv_init(bufferRecv.data(), bufferSize, realMPI, procRecv, + MPI_SAFE_CALL(MPI_Recv_init(bufferRecv.data(), bufferRecv.Size(), realMPI, procRecv, 650, data->mygrid->AxisComm, &recvRequest)); #endif diff --git a/src/fluid/boundary/axis.hpp b/src/fluid/boundary/axis.hpp index 229c772f..390a628c 100644 --- a/src/fluid/boundary/axis.hpp +++ b/src/fluid/boundary/axis.hpp @@ -11,6 +11,7 @@ #include #include "idefix.hpp" #include "grid.hpp" +#include "buffer.hpp" // Forward class hydro declaration #include "physics.hpp" @@ -59,8 +60,8 @@ class Axis { MPI_Request sendRequest; MPI_Request recvRequest; - IdefixArray1D bufferSend; - IdefixArray1D bufferRecv; + Buffer bufferSend; + Buffer bufferRecv; int bufferSize; diff --git a/src/mpi/buffer.hpp b/src/mpi/buffer.hpp index ff4dff14..35164355 100644 --- a/src/mpi/buffer.hpp +++ b/src/mpi/buffer.hpp @@ -31,6 +31,10 @@ class Buffer { return(array.data()); } + const IdefixArray1D & getArray(void) const { + return this->array; + } + int Size() { return(array.size()); } @@ -159,6 +163,33 @@ class Buffer { this->pointer += ninjnk; } + void UnpackJDirSymetric(IdefixArray4D& out, + const int var, + const int symMultiplier, + BoundingBox box) { + const int ni = box[IDIR][1]-box[IDIR][0]; + const int ninj = (box[JDIR][1]-box[JDIR][0])*ni; + const int ninjnk = (box[KDIR][1]-box[KDIR][0])*ninj; + const int ibeg = box[IDIR][0]; + const int jbeg = box[JDIR][0]; + const int kbeg = box[KDIR][0]; + const int iend = box[IDIR][1]; + const int jend = box[JDIR][1]; + const int kend = box[KDIR][1]; + const int offset = this->pointer; + + auto arr = this->array; + idefix_for("UnLoadBuffer4D_var_sym",kbeg,kend,jbeg,jend,ibeg,iend, + KOKKOS_LAMBDA (int k, int j, int i) { + const int jinverted = jend-(j-jbeg)-1; + const int arrIndex = i-ibeg + (j-jbeg)*ni + (k-kbeg)*ninj + offset; + out(var,k,jinverted,i) = symMultiplier * arr(arrIndex); + }); + + // Update pointer + this->pointer += ninjnk; + } + void Unpack(IdefixArray4D& out, IdefixArray1D& map, BoundingBox box) { @@ -186,6 +217,36 @@ class Buffer { this->pointer += ninjnk*map.size(); } + void UnpackJDirSymetric(IdefixArray4D& out, + IdefixArray1D& map, + IdefixArray1D& SymMap, + BoundingBox box) { + const int ni = box[IDIR][1]-box[IDIR][0]; + const int ninj = (box[JDIR][1]-box[JDIR][0])*ni; + const int ninjnk = (box[KDIR][1]-box[KDIR][0])*ninj; + const int ibeg = box[IDIR][0]; + const int jbeg = box[JDIR][0]; + const int kbeg = box[KDIR][0]; + const int iend = box[IDIR][1]; + const int jend = box[JDIR][1]; + const int kend = box[KDIR][1]; + const int offset = this->pointer; + + auto arr = this->array; + idefix_for("UnLoadBuffer4D_map_sym",0,map.size(), + kbeg,kend, + jbeg,jend, + ibeg,iend, + KOKKOS_LAMBDA (int n, int k, int j, int i) { + const int jinverted = jend-(j-jbeg)-1; + const int arrIndex = i-ibeg + (j-jbeg)*ni + (k-kbeg)*ninj + n*ninjnk + offset; + out(map(n),k,jinverted,i) = SymMap(map(n)) * arr(arrIndex); + }); + + // Update pointer + this->pointer += ninjnk*map.size(); + } + private: size_t pointer;