Changeset 27381


Ignore:
Timestamp:
11/10/22 17:04:28 (2 years ago)
Author:
Mathieu Morlighem
Message:

CHG: propagating deadlocking fix to other parts of the code

File:
1 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/classes/Nodes.cpp

    r27380 r27381  
    181181         * object that is not a clone, tell them to show their dofs, so that later on, they can get picked
    182182         * 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
    185189        for(int rank=0;rank<num_procs;rank++){
    186190                if(this->common_send[rank]){
    187191                        int  numids = this->common_send[rank];
     192                        send_truedofs[rank] = xNew<int>(numids*maxdofspernode);
    188193                        for(int i=0;i<numids;i++){
    189194                                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]);
    193198                }
    194199        }
     
    196201                if(this->common_recv[rank]){
    197202                        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);
    199204                        for(int i=0;i<numids;i++){
    200205                                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);
    206217
    207218        /*Update indexingupdateflag*/
Note: See TracChangeset for help on using the changeset viewer.