/[escript]/branches/trilinos_from_5897/dudley/src/ElementFile_distributeByRankOfDOF.cpp
ViewVC logotype

Diff of /branches/trilinos_from_5897/dudley/src/ElementFile_distributeByRankOfDOF.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 6078 by caltinay, Wed Mar 2 04:13:26 2016 UTC revision 6079 by caltinay, Mon Mar 21 12:22:38 2016 UTC
# Line 14  Line 14 
14  *  *
15  *****************************************************************************/  *****************************************************************************/
16    
 /****************************************************************************  
   
  Dudley: ElementFile: this will redistribute the Elements including overlap by  
   
 *****************************************************************************/  
   
17  #include "ElementFile.h"  #include "ElementFile.h"
18    
19  namespace dudley {  namespace dudley {
20    
21  void Dudley_ElementFile_distributeByRankOfDOF(Dudley_ElementFile* self, int* mpiRankOfDOF, index_t* Id)  void ElementFile::distributeByRankOfDOF(const int* mpiRankOfDOF,
22                                            const index_t* nodesId)
23  {  {
24      if (self == NULL)      const int size = MPIInfo->size;
25          return;      if (size > 1) {
     dim_t e, i;  
     int myRank = self->MPIInfo->rank;  
     dim_t NN = self->numNodes;  
     dim_t size = self->MPIInfo->size;  
     if (size > 1)  
     {  
26  #ifdef ESYS_MPI  #ifdef ESYS_MPI
27          int p, *Owner_buffer = NULL, loc_proc_mask_max;          const int myRank = MPIInfo->rank;
28          dim_t j, *send_count = NULL, *recv_count = NULL, *newOwner = NULL;          int numRequests = 0;
29          dim_t *loc_proc_mask = NULL, *loc_send_count = NULL;          std::vector<MPI_Request> mpi_requests(8 * size);
30          dim_t newNumElements, numElementsInBuffer;          std::vector<MPI_Status> mpi_stati(8 * size);
31          index_t *send_offset = NULL, *recv_offset = NULL, *Id_buffer = NULL, *Tag_buffer = NULL, *Nodes_buffer = NULL, k;  
         bool *proc_mask = NULL;  
         size_t size_size = size * sizeof(dim_t);  
         dim_t numRequests = 0;  
         MPI_Request *mpi_requests = NULL;  
         MPI_Status *mpi_stati = NULL;  
         mpi_requests = new  MPI_Request[8 * size];  
         mpi_stati = new  MPI_Status[8 * size];  
32          // count the number elements that have to be sent to each processor          // count the number elements that have to be sent to each processor
33          // (send_count) and define a new element owner as the processor with          // (send_count) and define a new element owner as the processor with
34          // the largest number of DOFs and the smallest id          // the largest number of DOFs and the smallest id
35          send_count = new dim_t[size];          std::vector<dim_t> send_count(size);
36          recv_count = new dim_t[size];          std::vector<dim_t> recv_count(size);
37          newOwner = new int[self->numElements];          int* newOwner = new int[numElements];
38          memset(send_count, 0, size_size);  #pragma omp parallel
 #pragma omp parallel private(p,loc_proc_mask,loc_send_count)  
39          {          {
40              loc_proc_mask = new dim_t[size];              std::vector<dim_t> loc_proc_mask(size);
41              loc_send_count = new dim_t[size];              std::vector<dim_t> loc_send_count(size);
42              memset(loc_send_count, 0, size_size);  #pragma omp for
43  #pragma omp for private(e,j,loc_proc_mask_max) schedule(static)              for (index_t e = 0; e < numElements; e++) {
44              for (e = 0; e < self->numElements; e++)                  if (Owner[e] == myRank) {
             {  
                 if (self->Owner[e] == myRank)  
                 {  
45                      newOwner[e] = myRank;                      newOwner[e] = myRank;
46                      memset(loc_proc_mask, 0, size_size);                      loc_proc_mask.assign(size, 0);
47                      for (j = 0; j < NN; j++)                      for (int j = 0; j < numNodes; j++) {
48                      {                          const int p = mpiRankOfDOF[Nodes[INDEX2(j, e, numNodes)]];
                         p = mpiRankOfDOF[self->Nodes[INDEX2(j, e, NN)]];  
49                          loc_proc_mask[p]++;                          loc_proc_mask[p]++;
50                      }                      }
51                      loc_proc_mask_max = 0;                      dim_t loc_proc_mask_max = 0;
52                      for (p = 0; p < size; ++p)                      for (int p = 0; p < size; ++p) {
                     {  
53                          if (loc_proc_mask[p] > 0)                          if (loc_proc_mask[p] > 0)
54                              loc_send_count[p]++;                              loc_send_count[p]++;
55                          if (loc_proc_mask[p] > loc_proc_mask_max)                          if (loc_proc_mask[p] > loc_proc_mask_max) {
                         {  
56                              newOwner[e] = p;                              newOwner[e] = p;
57                              loc_proc_mask_max = loc_proc_mask[p];                              loc_proc_mask_max = loc_proc_mask[p];
58                          }                          }
59                      }                      }
60                  }                  } else {
                 else  
                 {  
61                      newOwner[e] = -1;                      newOwner[e] = -1;
62                  }                  }
63              }              }
64  #pragma omp critical  #pragma omp critical
65              {              {
66                  for (p = 0; p < size; ++p)                  for (int p = 0; p < size; ++p)
67                      send_count[p] += loc_send_count[p];                      send_count[p] += loc_send_count[p];
68              }              }
69              delete[] loc_proc_mask;          } // end parallel section
70              delete[] loc_send_count;          MPI_Alltoall(&send_count[0], 1, MPI_DIM_T, &recv_count[0], 1,
71          }                       MPI_DIM_T, MPIInfo->comm);
72          MPI_Alltoall(send_count, 1, MPI_INT, recv_count, 1, MPI_INT, self->MPIInfo->comm);          // get the new number of elements for this processor
73          /* get the new number of elements for this processor */          dim_t newNumElements = 0;
74          newNumElements = 0;          dim_t numElementsInBuffer = 0;
75          for (p = 0; p < size; ++p)          for (int p = 0; p < size; ++p) {
76              newNumElements += recv_count[p];              newNumElements += recv_count[p];
   
         /* get the new number of elements for this processor */  
         numElementsInBuffer = 0;  
         for (p = 0; p < size; ++p)  
77              numElementsInBuffer += send_count[p];              numElementsInBuffer += send_count[p];
78          /* allocate buffers */          }
79          Id_buffer = new  index_t[numElementsInBuffer];  
80          Tag_buffer = new  index_t[numElementsInBuffer];          std::vector<index_t> Id_buffer(numElementsInBuffer);
81          Owner_buffer = new  int[numElementsInBuffer];          std::vector<int> Tag_buffer(numElementsInBuffer);
82          Nodes_buffer = new  index_t[numElementsInBuffer * NN];          std::vector<int> Owner_buffer(numElementsInBuffer);
83          send_offset = new  index_t[size];          std::vector<index_t> Nodes_buffer(numElementsInBuffer * numNodes);
84          recv_offset = new  index_t[size];          std::vector<index_t> send_offset(size);
85          proc_mask = new  bool[size];          std::vector<index_t> recv_offset(size);
86            std::vector<unsigned char> proc_mask(size);
87          /* calculate the offsets for the processor buffers */  
88          recv_offset[0] = 0;          // calculate the offsets for the processor buffers
89          for (p = 0; p < size - 1; ++p)          for (int p = 0; p < size - 1; ++p) {
90              recv_offset[p + 1] = recv_offset[p] + recv_count[p];              recv_offset[p + 1] = recv_offset[p] + recv_count[p];
         send_offset[0] = 0;  
         for (p = 0; p < size - 1; ++p)  
91              send_offset[p + 1] = send_offset[p] + send_count[p];              send_offset[p + 1] = send_offset[p] + send_count[p];
92            }
93    
94          memset(send_count, 0, size_size);          send_count.assign(size, 0);
95          /* copy element into buffers. proc_mask makes sure that an          // copy element into buffers. proc_mask makes sure that an element is
96           * element is copied once only for each processor */          // copied once only for each processor
97          for (e = 0; e < self->numElements; e++)          for (index_t e = 0; e < numElements; e++) {
98          {              if (Owner[e] == myRank) {
99              if (self->Owner[e] == myRank)                  proc_mask.assign(size, 1);
100              {                  for (int j = 0; j < numNodes; j++) {
101                  memset(proc_mask, true, size*sizeof(bool));                      const int p = mpiRankOfDOF[Nodes[INDEX2(j, e, numNodes)]];
102                  for (j = 0; j < NN; j++)                      if (proc_mask[p]) {
103                  {                          const index_t k = send_offset[p] + send_count[p];
104                      p = mpiRankOfDOF[self->Nodes[INDEX2(j, e, NN)]];                          Id_buffer[k] = Id[e];
105                      if (proc_mask[p])                          Tag_buffer[k] = Tag[e];
                     {  
                         k = send_offset[p] + send_count[p];  
                         Id_buffer[k] = self->Id[e];  
                         Tag_buffer[k] = self->Tag[e];  
106                          Owner_buffer[k] = newOwner[e];                          Owner_buffer[k] = newOwner[e];
107                          for (i = 0; i < NN; i++)                          for (int i = 0; i < numNodes; i++)
108                              Nodes_buffer[INDEX2(i, k, NN)] = Id[self->Nodes[INDEX2(i, e, NN)]];                              Nodes_buffer[INDEX2(i, k, numNodes)] =
109                                             nodesId[Nodes[INDEX2(i, e, numNodes)]];
110                          send_count[p]++;                          send_count[p]++;
111                          proc_mask[p] = false;                          proc_mask[p] = 0;
112                      }                      }
113                  }                  }
114              }              }
115          }          }
116          /* allocate new tables */          // allocate new tables
117          Dudley_ElementFile_allocTable(self, newNumElements);          allocTable(newNumElements);
118    
119          /* start to receive new elements */          // start to receive new elements
120          for (p = 0; p < size; ++p)          for (int p = 0; p < size; ++p) {
121          {              if (recv_count[p] > 0) {
122              if (recv_count[p] > 0)                  MPI_Irecv(&Id[recv_offset[p]], recv_count[p], MPI_DIM_T, p,
123              {                            MPIInfo->counter() + myRank, MPIInfo->comm,
124                  MPI_Irecv(&(self->Id[recv_offset[p]]), recv_count[p],                            &mpi_requests[numRequests]);
125                            MPI_INT, p, self->MPIInfo->counter() + myRank,                  numRequests++;
126                            self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Irecv(&Tag[recv_offset[p]], recv_count[p], MPI_INT, p,
127                  numRequests++;                            MPIInfo->counter() + size + myRank, MPIInfo->comm,
128                  MPI_Irecv(&(self->Tag[recv_offset[p]]), recv_count[p],                            &mpi_requests[numRequests]);
129                            MPI_INT, p, self->MPIInfo->counter() + size + myRank,                  numRequests++;
130                            self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Irecv(&Owner[recv_offset[p]], recv_count[p], MPI_INT, p,
131                  numRequests++;                            MPIInfo->counter() + 2 * size + myRank,
132                  MPI_Irecv(&(self->Owner[recv_offset[p]]), recv_count[p],                            MPIInfo->comm, &mpi_requests[numRequests]);
133                            MPI_INT, p, self->MPIInfo->counter() + 2 * size + myRank,                  numRequests++;
134                            self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Irecv(&Nodes[recv_offset[p] * numNodes],
135                  numRequests++;                            recv_count[p] * numNodes, MPI_DIM_T, p,
136                  MPI_Irecv(&(self->Nodes[recv_offset[p] * NN]), recv_count[p] * NN,                            MPIInfo->counter() + 3 * size + myRank,
137                            MPI_INT, p, self->MPIInfo->counter() + 3 * size + myRank,                            MPIInfo->comm, &mpi_requests[numRequests]);
                           self->MPIInfo->comm, &mpi_requests[numRequests]);  
138                  numRequests++;                  numRequests++;
139              }              }
140          }          }
141          /* now the buffers can be send away */          // now the buffers can be sent away
142          for (p = 0; p < size; ++p)          for (int p = 0; p < size; ++p) {
143          {              if (send_count[p] > 0) {
144              if (send_count[p] > 0)                  MPI_Issend(&Id_buffer[send_offset[p]], send_count[p],
145              {                             MPI_DIM_T, p, MPIInfo->counter() + p,
146                  MPI_Issend(&(Id_buffer[send_offset[p]]), send_count[p],                             MPIInfo->comm, &mpi_requests[numRequests]);
147                             MPI_INT, p, self->MPIInfo->counter() + p,                  numRequests++;
148                             self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Issend(&Tag_buffer[send_offset[p]], send_count[p],
149                  numRequests++;                             MPI_INT, p, MPIInfo->counter() + size + p,
150                  MPI_Issend(&(Tag_buffer[send_offset[p]]), send_count[p],                             MPIInfo->comm, &mpi_requests[numRequests]);
151                             MPI_INT, p, self->MPIInfo->counter() + size + p,                  numRequests++;
152                             self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Issend(&Owner_buffer[send_offset[p]], send_count[p],
153                  numRequests++;                             MPI_INT, p, MPIInfo->counter() + 2 * size + p,
154                  MPI_Issend(&(Owner_buffer[send_offset[p]]), send_count[p],                             MPIInfo->comm, &mpi_requests[numRequests]);
155                             MPI_INT, p, self->MPIInfo->counter() + 2 * size + p,                  numRequests++;
156                             self->MPIInfo->comm, &mpi_requests[numRequests]);                  MPI_Issend(&Nodes_buffer[send_offset[p] * numNodes],
157                  numRequests++;                             send_count[p] * numNodes, MPI_DIM_T, p,
158                  MPI_Issend(&(Nodes_buffer[send_offset[p] * NN]), send_count[p] * NN,                             MPIInfo->counter() + 3 * size + p,
159                             MPI_INT, p, self->MPIInfo->counter() + 3 * size + p,                             MPIInfo->comm, &mpi_requests[numRequests]);
                            self->MPIInfo->comm, &mpi_requests[numRequests]);  
160                  numRequests++;                  numRequests++;
161              }              }
162          }          }
163          /* wait for the requests to be finalized */          MPIInfo->incCounter(4 * size);
164          self->MPIInfo->incCounter(4 * size);          // wait for the requests to be finalized
165          MPI_Waitall(numRequests, mpi_requests, mpi_stati);          MPI_Waitall(numRequests, &mpi_requests[0], &mpi_stati[0]);
         /* clear buffer */  
         delete[] Id_buffer;  
         delete[] Tag_buffer;  
         delete[] Owner_buffer;  
         delete[] Nodes_buffer;  
         delete[] send_offset;  
         delete[] recv_offset;  
         delete[] proc_mask;  
         delete[] mpi_requests;  
         delete[] mpi_stati;  
         delete[] send_count;  
         delete[] recv_count;  
166          delete[] newOwner;          delete[] newOwner;
167  #endif  #endif
168      } else { // single rank      } else { // single rank
169  #pragma omp for private(e,i) schedule(static)  #pragma omp parallel for
170          for (e = 0; e < self->numElements; e++)          for (index_t e = 0; e < numElements; e++) {
171          {              Owner[e] = 0;
172              self->Owner[e] = myRank;              for (int i = 0; i < numNodes; i++)
173              for (i = 0; i < NN; i++)                  Nodes[INDEX2(i, e, numNodes)] =
174                  self->Nodes[INDEX2(i, e, NN)] = Id[self->Nodes[INDEX2(i, e, NN)]];                                       nodesId[Nodes[INDEX2(i, e, numNodes)]];
175          }          }
176      }      }
177  }  }

Legend:
Removed from v.6078  
changed lines
  Added in v.6079

  ViewVC Help
Powered by ViewVC 1.1.26