Changeset 27381
- Timestamp:
- 11/10/22 17:04:28 (2 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/c/classes/Nodes.cpp
r27380 r27381 181 181 * object that is not a clone, tell them to show their dofs, so that later on, they can get picked 182 182 * up by their clones: */ 183 int maxdofspernode = this->MaxNumDofs(setenum); 184 int* truedofs = xNewZeroInit<int>(this->Size()*maxdofspernode); //only one alloc 183 int maxdofspernode = this->MaxNumDofs(GsetEnum); 184 int **send_truedofs = xNewZeroInit<int*>(num_procs); 185 int *recv_truedofs = xNewZeroInit<int>(this->Size()*maxdofspernode); 186 ISSM_MPI_Request *send_requests = xNew<ISSM_MPI_Request>(num_procs); 187 for (int rank = 0;rank<num_procs;rank++) send_requests[rank] = ISSM_MPI_REQUEST_NULL; 188 185 189 for(int rank=0;rank<num_procs;rank++){ 186 190 if(this->common_send[rank]){ 187 191 int numids = this->common_send[rank]; 192 send_truedofs[rank] = xNew<int>(numids*maxdofspernode); 188 193 for(int i=0;i<numids;i++){ 189 194 Node* node=xDynamicCast<Node*>(this->GetObjectByOffset(this->common_send_ids[rank][i])); 190 node->ShowMasterDofs(& truedofs[i*maxdofspernode+0],setenum);191 } 192 ISSM_MPI_ Send(truedofs,numids*maxdofspernode,ISSM_MPI_INT,rank,0,IssmComm::GetComm());195 node->ShowMasterDofs(&send_truedofs[rank][i*maxdofspernode+0],setenum); 196 } 197 ISSM_MPI_Isend(send_truedofs[rank],numids*maxdofspernode,ISSM_MPI_INT,rank,0,IssmComm::GetComm(),&send_requests[rank]); 193 198 } 194 199 } … … 196 201 if(this->common_recv[rank]){ 197 202 int numids = this->common_recv[rank]; 198 ISSM_MPI_Recv( truedofs,numids*maxdofspernode,ISSM_MPI_INT,rank,0,IssmComm::GetComm(),&status);203 ISSM_MPI_Recv(recv_truedofs,numids*maxdofspernode,ISSM_MPI_INT,rank,0,IssmComm::GetComm(),&status); 199 204 for(int i=0;i<numids;i++){ 200 205 Node* node=xDynamicCast<Node*>(this->GetObjectByOffset(this->common_recv_ids[rank][i])); 201 node->UpdateCloneDofs(&truedofs[i*maxdofspernode+0],setenum); 202 } 203 } 204 } 205 xDelete<int>(truedofs); 206 node->UpdateCloneDofs(&recv_truedofs[i*maxdofspernode+0],setenum); 207 } 208 } 209 } 210 xDelete<int>(recv_truedofs); 211 for(int rank=0;rank<num_procs;rank++){ 212 if(rank!=my_rank) ISSM_MPI_Wait(&send_requests[rank],&status); 213 xDelete<int>(send_truedofs[rank]); 214 } 215 xDelete<int*>(send_truedofs); 216 xDelete<ISSM_MPI_Request>(send_requests); 206 217 207 218 /*Update indexingupdateflag*/
Note:
See TracChangeset
for help on using the changeset viewer.