@@ -34,10 +34,197 @@ class MPMesh{
3434 std::vector<std::vector<int >> ownerHaloLocalIDs;
3535
3636 void startCommunication ();
37- void communicateFields (const std::vector<std::vector<double >>& fieldData, const int numEntities, const int numEntries, int mode,
38- std::vector<std::vector<int >>& recvIDVec, std::vector<std::vector<double >>& recvDataVec);
37+
3938 void communicate_and_take_halo_contributions (const Kokkos::View<double **>& meshField, int nEntities, int numEntries, int mode, int op);
39+ // Now Kokkos views are made 1D
40+ template <typename ViewType>
41+ void communicate_and_take_halo_contributions1 (
42+ const ViewType& meshField,
43+ int nEntities,
44+ int numEntries,
45+ int mode ,
46+ int op){
47+
48+ int self;
49+ MPI_Comm comm = p_MPs->getMPIComm ();
50+ MPI_Comm_rank (comm, &self);
51+
52+ Kokkos::Timer timer;
53+ auto reconVals_host = Kokkos::create_mirror_view_and_copy (Kokkos::HostSpace (), meshField);
54+ Kokkos::fence ();
55+ pumipic::RecordTime (" Communication-GPU to CPU-E-" + std::to_string (numEntries) + " -" + std::to_string (self), timer.seconds ());
56+
57+ timer.reset ();
58+ std::vector<double > fieldData1 (nEntities*numEntries);
59+ std::memcpy (fieldData1.data (), reconVals_host.data (), nEntities*numEntries*sizeof (double ));
60+ pumipic::RecordTime (" New allocation + copy" + std::to_string (numEntries) + " -" + std::to_string (self), timer.seconds ());
61+
62+ timer.reset ();
63+ std::vector<std::vector<int >> recvIDVec;
64+ std::vector<std::vector<double >> recvDataVec;
65+ communicateFields1 (fieldData1, nEntities, numEntries, mode, recvIDVec, recvDataVec);
66+ pumipic::RecordTime (" Communication-InterProcess-E-" + std::to_string (numEntries) + " -" + std::to_string (self), timer.seconds ());
67+
68+ timer.reset ();
69+ int numProcsTot = recvIDVec.size ();
70+ // Flatten IDs
71+ int totalSize = 0 ;
72+ std::vector<int > offsets (numProcsTot, 0 );
73+ for (int i=0 ; i<numProcsTot; i++) {
74+ offsets[i] = totalSize;
75+ totalSize += recvIDVec[i].size ();
76+ }
77+ std::vector<int > flatIDVec (totalSize, 0 );
78+ for (int i=0 ; i<numProcsTot; i++) {
79+ std::copy (recvIDVec[i].begin (), recvIDVec[i].end (), flatIDVec.begin () + offsets[i]);
80+ }
4081
82+ Kokkos::View<int *> recvIDGPU (" recvIDGPU" , totalSize);
83+ auto hostView = Kokkos::View<int *, Kokkos::HostSpace>(" recvIDCPU" , totalSize);
84+ std::copy (flatIDVec.begin (), flatIDVec.end (), hostView.data ());
85+ Kokkos::deep_copy (recvIDGPU, hostView);
86+
87+ // Flatten Data
88+ int totalSize_data=0 ;
89+ std::vector<int > offsets_data (numProcsTot, 0 );
90+ for (int i=0 ; i<numProcsTot; i++) {
91+ offsets_data[i] = totalSize_data;
92+ totalSize_data += recvDataVec[i].size ();
93+ }
94+ std::vector<double > flatDataVec (totalSize_data, 0 );
95+ for (int i=0 ; i<numProcsTot; i++) {
96+ std::copy (recvDataVec[i].begin (), recvDataVec[i].end (), flatDataVec.begin () + offsets_data[i]);
97+ }
98+ Kokkos::View<double *> recvDataGPU (" recvDataGPU" , totalSize_data);
99+ auto hostView_data= Kokkos::View<double *, Kokkos::HostSpace>(" recvDataCPU" , totalSize_data);
100+ std::copy (flatDataVec.begin (), flatDataVec.end (), hostView_data.data ());
101+ Kokkos::deep_copy (recvDataGPU, hostView_data);
102+ Kokkos::fence ();
103+ // Assertions
104+ assert (totalSize_data == totalSize*numEntries);
105+ for (int i=0 ; i<numProcsTot; i++){
106+ assert (recvDataVec[i].size () == recvIDVec[i].size () * numEntries);
107+ }
108+ pumipic::RecordTime (" Communication-CPU to GPU-E-" + std::to_string (numEntries) + " -" + std::to_string (self), timer.seconds ());
109+
110+ // Take contributions from other procs
111+ timer.reset ();
112+ Kokkos::parallel_for (" halo contribution" , recvIDGPU.size (), KOKKOS_LAMBDA (const int i){
113+ int vertex = recvIDGPU (i);
114+ for (int k=0 ; k<numEntries; k++){
115+ if (op==0 ) Kokkos::atomic_add (&meshField (vertex,k), recvDataGPU (i*numEntries+k));
116+ if (op==1 ) meshField (vertex, k) = recvDataGPU (i * numEntries + k);
117+ }
118+ });
119+ Kokkos::fence ();
120+ pumipic::RecordTime (" Communication-GPU reduction-E-" + std::to_string (numEntries) + " -" + std::to_string (self), timer.seconds ());
121+ }
122+
123+
124+ void communicateFields (const std::vector<std::vector<double >>& fieldData, const int numEntities, const int numEntries, int mode,
125+ std::vector<std::vector<int >>& recvIDVec, std::vector<std::vector<double >>& recvDataVec);
126+ // 1D Kokkos view will be copied direcctly to 1D std::vector
127+ void communicateFields1 (
128+ const std::vector<double >& fieldData,
129+ const int numEntities, const int numEntries, int mode,
130+ std::vector<std::vector<int >>& recvIDVec,
131+ std::vector<std::vector<double >>& recvDataVec){
132+
133+ int self, numProcsTot;
134+ MPI_Comm comm = p_MPs->getMPIComm ();
135+ MPI_Comm_rank (comm, &self);
136+ MPI_Comm_size (comm, &numProcsTot);
137+
138+ assert (numEntities == numOwnersTot + numHalosTot);
139+
140+ std::vector<std::vector<double >> sendDataVec (numProcsTot);
141+
142+ recvIDVec.resize (numProcsTot);
143+ recvDataVec.resize (numProcsTot);
144+
145+ for (int i = 0 ; i < numProcsTot; i++){
146+ if (i==self) continue ;
147+
148+ int numToSend = 0 , numToRecv = 0 ;
149+ if (mode == 0 ) {
150+ // gather (halos send to owners)
151+ numToSend = numOwnersOnOtherProcs[i];
152+ numToRecv = numHalosOnOtherProcs[i];
153+ }
154+ else {
155+ // scatter (owners send to halos)
156+ numToSend = numHalosOnOtherProcs[i];
157+ numToRecv = numOwnersOnOtherProcs[i];
158+ }
159+
160+ if (numToSend > 0 ){
161+ sendDataVec[i].reserve (numToSend*numEntries);
162+ }
163+ if (numToRecv > 0 ){
164+ recvDataVec[i].resize (numToRecv*numEntries);
165+ recvIDVec[i].resize (numToRecv);
166+ }
167+ }
168+
169+ if (mode == 0 ){
170+ // Halos sends to owners
171+ for (int iEnt = 0 ; iEnt < numHalosTot; iEnt++){
172+ auto ownerProc = haloOwnerProcs[iEnt];
173+ for (int iDouble = 0 ; iDouble < numEntries; iDouble++)
174+ sendDataVec[ownerProc].push_back (fieldData[(numOwnersTot+iEnt)*numEntries + iDouble]);
175+ }
176+ }
177+ else if (mode == 1 ){
178+ // Owner sends to halos
179+ for (size_t iProc=0 ; iProc<ownerOwnerLocalIDs.size (); iProc++) {
180+ for (auto & ownerID : ownerOwnerLocalIDs[iProc]) {
181+ for (int iDouble = 0 ; iDouble < numEntries; iDouble++)
182+ sendDataVec[iProc].push_back (fieldData[ownerID*numEntries + iDouble]);
183+ }
184+ }
185+ }
186+
187+ std::vector<MPI_Request> requests;
188+ requests.reserve (4 *numProcsTot);
189+ for (int proc = 0 ; proc < numProcsTot; proc++){
190+ if (proc == self) continue ;
191+ if (mode == 0 && numHalosOnOtherProcs[proc]){
192+ assert (recvIDVec[proc].size () == (size_t )numHalosOnOtherProcs[proc]);
193+ assert (recvDataVec[proc].size () == recvIDVec[proc].size () * (size_t )numEntries);
194+ MPI_Request req3, req4;
195+ MPI_Irecv (recvIDVec[proc].data (), recvIDVec[proc].size (), MPI_INT , proc, 1 , comm, &req3);
196+ MPI_Irecv (recvDataVec[proc].data (), recvDataVec[proc].size (), MPI_DOUBLE , proc, 2 , comm, &req4);
197+ requests.push_back (req3);
198+ requests.push_back (req4);
199+ }
200+ if (mode == 0 && numOwnersOnOtherProcs[proc]) {
201+ assert (haloOwnerLocalIDs[proc].size () == (size_t )numOwnersOnOtherProcs[proc]);
202+ assert (sendDataVec[proc].size () == haloOwnerLocalIDs[proc].size () * (size_t )numEntries);
203+ MPI_Request req1, req2;
204+ MPI_Isend (haloOwnerLocalIDs[proc].data (), haloOwnerLocalIDs[proc].size (), MPI_INT , proc, 1 , comm, &req1);
205+ MPI_Isend (sendDataVec[proc].data (), sendDataVec[proc].size (), MPI_DOUBLE , proc, 2 , comm, &req2);
206+ requests.push_back (req1);
207+ requests.push_back (req2);
208+ }
209+
210+ if (mode == 1 && numOwnersOnOtherProcs[proc]){
211+ MPI_Request req3, req4;
212+ MPI_Irecv (recvIDVec[proc].data (), recvIDVec[proc].size (), MPI_INT , proc, 1 , comm, &req3);
213+ MPI_Irecv (recvDataVec[proc].data (), recvDataVec[proc].size (), MPI_DOUBLE , proc, 2 , comm, &req4);
214+ requests.push_back (req3);
215+ requests.push_back (req4);
216+ }
217+ if (mode == 1 && numHalosOnOtherProcs[proc]) {
218+ MPI_Request req1, req2;
219+ MPI_Isend (ownerHaloLocalIDs[proc].data (), ownerHaloLocalIDs[proc].size (), MPI_INT , proc, 1 , comm, &req1);
220+ MPI_Isend (sendDataVec[proc].data (), sendDataVec[proc].size (), MPI_DOUBLE , proc, 2 , comm, &req2);
221+ requests.push_back (req1);
222+ requests.push_back (req2);
223+ }
224+ }
225+ MPI_Waitall (requests.size (), requests.data (), MPI_STATUSES_IGNORE );
226+ }
227+
41228 MPMesh (Mesh* inMesh, MaterialPoints* inMPs):
42229 p_mesh (inMesh), p_MPs(inMPs) {
43230 };
0 commit comments