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

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

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

branches/trilinos_from_5897/dudley/src/Mesh_createNodeFileMappings.cpp revision 6078 by caltinay, Tue Mar 8 02:07:16 2016 UTC branches/trilinos_from_5897/dudley/src/NodeFile_createMappings.cpp revision 6079 by caltinay, Mon Mar 21 12:22:38 2016 UTC
# Line 25  Line 25 
25    
26  namespace dudley {  namespace dudley {
27    
28  #define UNUSED -1  void NodeFile::createDOFMappingAndCoupling()
   
 void Dudley_Mesh_createDOFMappingAndCoupling(Dudley_Mesh* in, bool use_reduced_elements)  
29  {  {
30      index_t min_DOF, max_DOF, *shared = NULL, *offsetInShared = NULL, *locDOFMask =      const index_t myFirstDOF = dofDistribution->getFirstComponent();
31          NULL, i, k, myFirstDOF, myLastDOF, *nodeMask = NULL, firstDOF, lastDOF, *globalDOFIndex, *wanted_DOFs = NULL;      const index_t myLastDOF = dofDistribution->getLastComponent();
32      dim_t mpiSize, len_loc_dof, numNeighbors, n, lastn, numNodes, *rcv_len = NULL, *snd_len = NULL, count;      const int mpiSize = MPIInfo->size;
33      int myRank, p, p_min, p_max, *neighbor = NULL;      const int myRank = MPIInfo->rank;
34      paso::SharedComponents_ptr rcv_shcomp, snd_shcomp;  
35      Dudley_NodeMapping *this_mapping = NULL;      index_t min_DOF, max_DOF;
36      paso::Connector_ptr this_connector;      std::pair<index_t,index_t> DOF_range(util::getFlaggedMinMaxInt(
37      paso::Distribution_ptr dof_distribution;                                              numNodes, globalDegreesOfFreedom, -1));
     escript::JMPI mpi_info(in->MPIInfo);  
 #ifdef ESYS_MPI  
     MPI_Request *mpi_requests = NULL;  
     MPI_Status *mpi_stati = NULL;  
 #else  
     int *mpi_requests = NULL, *mpi_stati = NULL;  
 #endif  
   
     numNodes = in->Nodes->numNodes;  
     if (use_reduced_elements)  
     {  
         dof_distribution = in->Nodes->reducedDegreesOfFreedomDistribution;  
         globalDOFIndex = in->Nodes->globalReducedDOFIndex;  
     }  
     else  
     {  
         dof_distribution = in->Nodes->degreesOfFreedomDistribution;  
         globalDOFIndex = in->Nodes->globalDegreesOfFreedom;  
     }  
     myFirstDOF = dof_distribution->getFirstComponent();  
     myLastDOF = dof_distribution->getLastComponent();  
   
     mpiSize = mpi_info->size;  
     myRank = mpi_info->rank;  
   
     min_DOF = Dudley_Util_getFlaggedMinInt(1, numNodes, globalDOFIndex, -1);  
     max_DOF = Dudley_Util_getFlaggedMaxInt(1, numNodes, globalDOFIndex, -1);  
38    
39      if (max_DOF < min_DOF)      if (DOF_range.second < DOF_range.first) {
     {  
40          min_DOF = myFirstDOF;          min_DOF = myFirstDOF;
41          max_DOF = myLastDOF - 1;          max_DOF = myLastDOF - 1;
42        } else {
43            min_DOF = DOF_range.first;
44            max_DOF = DOF_range.second;
45      }      }
46    
47      p_min = mpiSize;      int p_min = mpiSize;
48      p_max = -1;      int p_max = -1;
49      if (max_DOF >= min_DOF)      if (max_DOF >= min_DOF) {
50      {          for (int p = 0; p < mpiSize; ++p) {
51          for (p = 0; p < mpiSize; ++p)              if (dofDistribution->first_component[p] <= min_DOF)
         {  
             if (dof_distribution->first_component[p] <= min_DOF)  
52                  p_min = p;                  p_min = p;
53              if (dof_distribution->first_component[p] <= max_DOF)              if (dofDistribution->first_component[p] <= max_DOF)
54                  p_max = p;                  p_max = p;
55          }          }
56      }      }
57    
     len_loc_dof = max_DOF - min_DOF + 1;  
58      std::stringstream ss;      std::stringstream ss;
59      if (myFirstDOF<myLastDOF && !(min_DOF <= myFirstDOF && myLastDOF-1 <= max_DOF)) {      if (myFirstDOF<myLastDOF && !(min_DOF <= myFirstDOF && myLastDOF-1 <= max_DOF)) {
60          ss << "createDOFMappingAndCoupling: Local elements do not span local "          ss << "createDOFMappingAndCoupling: Local elements do not span local "
61                "degrees of freedom. min_DOF=" << min_DOF << ", myFirstDOF="                "degrees of freedom. min_DOF=" << min_DOF << ", myFirstDOF="
62             << myFirstDOF << ", myLastDOF-1=" << myLastDOF-1             << myFirstDOF << ", myLastDOF-1=" << myLastDOF-1
63             << ", max_DOF=" << max_DOF << " on rank=" << in->MPIInfo->rank;             << ", max_DOF=" << max_DOF << " on rank=" << MPIInfo->rank;
64      }      }
65      const std::string msg(ss.str());      const std::string msg(ss.str());
66      int error = msg.length();      int error = msg.length();
67      int gerror = error;      int gerror = error;
68      escript::checkResult(error, gerror, in->MPIInfo);      escript::checkResult(error, gerror, MPIInfo);
69      if (gerror > 0) {      if (gerror > 0) {
70          char* gmsg;          char* gmsg;
71          escript::shipString(msg.c_str(), &gmsg, in->MPIInfo->comm);          escript::shipString(msg.c_str(), &gmsg, MPIInfo->comm);
72          throw DudleyException(gmsg);          throw DudleyException(gmsg);
73      }      }
74    
75      rcv_len = new  dim_t[mpiSize];      const index_t UNUSED = -1;
76      snd_len = new  dim_t[mpiSize];      const dim_t len_loc_dof = max_DOF - min_DOF + 1;
77  #ifdef ESYS_MPI      index_t* shared = new index_t[numNodes * (p_max - p_min + 1)];
78      mpi_requests = new  MPI_Request[mpiSize * 2];      std::vector<index_t> offsetInShared(mpiSize + 1);
79      mpi_stati = new  MPI_Status[mpiSize * 2];      index_t* locDOFMask = new index_t[len_loc_dof];
80  #else      index_t* nodeMask = new index_t[numNodes];
81      mpi_requests = new  int[mpiSize * 2];  
82      mpi_stati = new  int[mpiSize * 2];      ESYS_ASSERT(myLastDOF-min_DOF <= len_loc_dof, "out of bounds!");
 #endif  
     wanted_DOFs = new  index_t[numNodes];  
     nodeMask = new  index_t[numNodes];  
     neighbor = new  int[mpiSize];  
     shared = new  index_t[numNodes * (p_max - p_min + 1)];  
     offsetInShared = new  index_t[mpiSize + 1];  
     locDOFMask = new  index_t[len_loc_dof];  
83    
     memset(rcv_len, 0, sizeof(dim_t) * mpiSize);  
84  #pragma omp parallel  #pragma omp parallel
85      {      {
86  #pragma omp for private(i) schedule(static)  #pragma omp for
87          for (i = 0; i < len_loc_dof; ++i)          for (index_t i = 0; i < len_loc_dof; ++i)
88              locDOFMask[i] = UNUSED;              locDOFMask[i] = UNUSED;
89  #pragma omp for private(i) schedule(static)  #pragma omp for
90          for (i = 0; i < numNodes; ++i)          for (index_t i = 0; i < numNodes; ++i)
91              nodeMask[i] = UNUSED;              nodeMask[i] = UNUSED;
92  #pragma omp for private(i,k) schedule(static)  #pragma omp for
93          for (i = 0; i < numNodes; ++i)          for (index_t i = 0; i < numNodes; ++i) {
94          {              const index_t k = globalDegreesOfFreedom[i];
95              k = globalDOFIndex[i];              if (k > -1) {
             if (k > -1)  
             {  
                 locDOFMask[k - min_DOF] = UNUSED - 1;  
96  #ifdef BOUNDS_CHECK  #ifdef BOUNDS_CHECK
97                  if ((k - min_DOF) >= len_loc_dof)                  if ((k - min_DOF) >= len_loc_dof) {
                 {  
98                      printf("BOUNDS_CHECK %s %d i=%d k=%d min_DOF=%d\n", __FILE__, __LINE__, i, k, min_DOF);                      printf("BOUNDS_CHECK %s %d i=%d k=%d min_DOF=%d\n", __FILE__, __LINE__, i, k, min_DOF);
99                      exit(1);                      exit(1);
100                  }                  }
101  #endif  #endif
102                    locDOFMask[k - min_DOF] = UNUSED - 1;
103              }              }
104          }          }
105    #pragma omp for
106  #pragma omp for private(i) schedule(static)          for (index_t i = myFirstDOF - min_DOF; i < myLastDOF - min_DOF; ++i) {
         for (i = myFirstDOF - min_DOF; i < myLastDOF - min_DOF; ++i)  
         {  
107              locDOFMask[i] = i - myFirstDOF + min_DOF;              locDOFMask[i] = i - myFirstDOF + min_DOF;
 #ifdef BOUNDS_CHECK  
             if (i < 0 || i >= len_loc_dof)  
             {  
                 printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i);  
                 exit(1);  
             }  
 #endif  
108          }          }
109      }      }
110    
111      numNeighbors = 0;      index_t* wanted_DOFs = new index_t[numNodes];
112      n = 0;      std::vector<index_t> rcv_len(mpiSize);
113      lastn = n;      std::vector<index_t> snd_len(mpiSize);
114      for (p = p_min; p <= p_max; ++p)      std::vector<int> neighbor(mpiSize);
115      {      int numNeighbors = 0;
116          firstDOF = std::max(min_DOF, dof_distribution->first_component[p]);      dim_t n = 0;
117          lastDOF = std::min(max_DOF + 1, dof_distribution->first_component[p + 1]);      dim_t lastn = n;
118          if (p != myRank)  
119          {      for (int p = p_min; p <= p_max; ++p) {
120              for (i = firstDOF - min_DOF; i < lastDOF - min_DOF; ++i)          if (p != myRank) {
121              {              const index_t firstDOF = std::max(min_DOF, dofDistribution->first_component[p]);
122                const index_t lastDOF = std::min(max_DOF + 1, dofDistribution->first_component[p + 1]);
123  #ifdef BOUNDS_CHECK  #ifdef BOUNDS_CHECK
124                  if (i < 0 || i >= len_loc_dof)              if (lastDOF-min_DOF > len_loc_dof) {
125                  {                  printf("BOUNDS_CHECK %s %d p=%d\n", __FILE__, __LINE__, p);
126                      printf("BOUNDS_CHECK %s %d p=%d i=%d\n", __FILE__, __LINE__, p, i);                  exit(1);
127                      exit(1);              }
                 }  
128  #endif  #endif
129                  if (locDOFMask[i] == UNUSED - 1)              for (index_t i = firstDOF - min_DOF; i < lastDOF - min_DOF; ++i) {
130                  {                  if (locDOFMask[i] == UNUSED - 1) {
131                      locDOFMask[i] = myLastDOF - myFirstDOF + n;                      locDOFMask[i] = myLastDOF - myFirstDOF + n;
132                      wanted_DOFs[n] = i + min_DOF;                      wanted_DOFs[n] = i + min_DOF;
133                      ++n;                      ++n;
134                  }                  }
135              }              }
136              if (n > lastn)              if (n > lastn) {
             {  
137                  rcv_len[p] = n - lastn;                  rcv_len[p] = n - lastn;
                 neighbor[numNeighbors] = p;  
138  #ifdef BOUNDS_CHECK  #ifdef BOUNDS_CHECK
139                  if (numNeighbors < 0 || numNeighbors >= mpiSize + 1)                  if (numNeighbors >= mpiSize + 1) {
140                  {                      printf("BOUNDS_CHECK %s %d p=%d numNeighbors=%d n=%d\n", __FILE__, __LINE__, p, numNeighbors, n);
                     printf("BOUNDS_CHECK %s %d p=%d numNeighbors=%d n=%d\n", __FILE__, __LINE__, p, numNeighbors,  
                            n);  
141                      exit(1);                      exit(1);
142                  }                  }
143  #endif  #endif
144                    neighbor[numNeighbors] = p;
145                  offsetInShared[numNeighbors] = lastn;                  offsetInShared[numNeighbors] = lastn;
146                  numNeighbors++;                  numNeighbors++;
147                  lastn = n;                  lastn = n;
148              }              }
149          }          } // if p!=myRank
150      }      } // for p
151    
152  #ifdef BOUNDS_CHECK  #ifdef BOUNDS_CHECK
153      if (numNeighbors < 0 || numNeighbors >= mpiSize + 1)      if (numNeighbors >= mpiSize + 1) {
     {  
154          printf("BOUNDS_CHECK %s %d numNeighbors=%d\n", __FILE__, __LINE__, numNeighbors);          printf("BOUNDS_CHECK %s %d numNeighbors=%d\n", __FILE__, __LINE__, numNeighbors);
155          exit(1);          exit(1);
156      }      }
157  #endif  #endif
158      offsetInShared[numNeighbors] = lastn;      offsetInShared[numNeighbors] = lastn;
159    
160      /* assign new DOF labels to nodes */      // assign new DOF labels to nodes
161  #pragma omp parallel for private(i,k) schedule(static)  #pragma omp parallel for
162      for (i = 0; i < numNodes; ++i)      for (index_t i = 0; i < numNodes; ++i) {
163      {          const index_t k = globalDegreesOfFreedom[i];
         k = globalDOFIndex[i];  
164          if (k > -1)          if (k > -1)
165              nodeMask[i] = locDOFMask[k - min_DOF];              nodeMask[i] = locDOFMask[k - min_DOF];
166      }      }
167    
168      /* now we can set the mapping from nodes to local DOFs */      degreesOfFreedomMapping.assign(nodeMask, numNodes, UNUSED);
169      this_mapping = Dudley_NodeMapping_alloc(numNodes, nodeMask, UNUSED);  
170      /* define how to get DOF values for controlled bu other processors */      // define how to get DOF values for controlled but other processors
171  #ifdef BOUNDS_CHECK  #ifdef BOUNDS_CHECK
172      for (i = 0; i < offsetInShared[numNeighbors]; ++i)      if (numNodes && offsetInShared[numNeighbors] >= numNodes * (p_max - p_min + 1)) {
173      {          printf("BOUNDS_CHECK %s %d\n", __FILE__, __LINE__);
174          if (i < 0 || i >= numNodes * (p_max - p_min + 1))          exit(1);
         {  
             printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i);  
             exit(1);  
         }  
175      }      }
176  #endif  #endif
177  #pragma omp parallel for private(i) schedule(static)  #pragma omp parallel for
178      for (i = 0; i < offsetInShared[numNeighbors]; ++i)      for (index_t i = 0; i < lastn; ++i)
179          shared[i] = myLastDOF - myFirstDOF + i;          shared[i] = myLastDOF - myFirstDOF + i;
180    
181      rcv_shcomp.reset(new paso::SharedComponents(myLastDOF - myFirstDOF,      paso::SharedComponents_ptr rcv_shcomp(new paso::SharedComponents(
182              numNeighbors, neighbor, shared, offsetInShared, 1, 0,                  myLastDOF - myFirstDOF, numNeighbors, &neighbor[0], shared,
183              mpi_info));                  &offsetInShared[0], 1, 0, MPIInfo));
184    
185      /*      /////////////////////////////////
186       *    now we build the sender      //   now we build the sender   //
187       */      /////////////////////////////////
188  #ifdef ESYS_MPI  #ifdef ESYS_MPI
189      MPI_Alltoall(rcv_len, 1, MPI_INT, snd_len, 1, MPI_INT, mpi_info->comm);      std::vector<MPI_Request> mpi_requests(mpiSize * 2);
190        std::vector<MPI_Status> mpi_stati(mpiSize * 2);
191        MPI_Alltoall(&rcv_len[0], 1, MPI_DIM_T, &snd_len[0], 1, MPI_DIM_T,
192                     MPIInfo->comm);
193        int count = 0;
194  #else  #else
195      for (p = 0; p < mpiSize; ++p)      snd_len[0] = rcv_len[0];
         snd_len[p] = rcv_len[p];  
196  #endif  #endif
197      count = 0;  
198      for (p = 0; p < rcv_shcomp->numNeighbors; p++)      for (int p = 0; p < rcv_shcomp->numNeighbors; p++) {
     {  
199  #ifdef ESYS_MPI  #ifdef ESYS_MPI
200          MPI_Isend(&(wanted_DOFs[rcv_shcomp->offsetInShared[p]]),          MPI_Isend(&wanted_DOFs[rcv_shcomp->offsetInShared[p]],
201                    rcv_shcomp->offsetInShared[p + 1] - rcv_shcomp->offsetInShared[p], MPI_INT,                  rcv_shcomp->offsetInShared[p+1]-rcv_shcomp->offsetInShared[p],
202                    rcv_shcomp->neighbor[p], mpi_info->counter() + myRank, mpi_info->comm,                  MPI_DIM_T, rcv_shcomp->neighbor[p],
203                    &mpi_requests[count]);                  MPIInfo->counter() + myRank, MPIInfo->comm,
204  #endif                  &mpi_requests[count]);
205          count++;          count++;
206    #endif
207      }      }
208      n = 0;      n = 0;
209      numNeighbors = 0;      numNeighbors = 0;
210      for (p = 0; p < mpiSize; p++)      for (int p = 0; p < mpiSize; p++) {
211      {          if (snd_len[p] > 0) {
         if (snd_len[p] > 0)  
         {  
212  #ifdef ESYS_MPI  #ifdef ESYS_MPI
213              MPI_Irecv(&(shared[n]), snd_len[p],              MPI_Irecv(&shared[n], snd_len[p], MPI_DIM_T, p,
214                        MPI_INT, p, mpi_info->counter() + p, mpi_info->comm, &mpi_requests[count]);                      MPIInfo->counter() + p, MPIInfo->comm,
215  #endif                      &mpi_requests[count]);
216              count++;              count++;
217    #endif
218              neighbor[numNeighbors] = p;              neighbor[numNeighbors] = p;
219              offsetInShared[numNeighbors] = n;              offsetInShared[numNeighbors] = n;
220              numNeighbors++;              numNeighbors++;
# Line 279  void Dudley_Mesh_createDOFMappingAndCoup Line 223  void Dudley_Mesh_createDOFMappingAndCoup
223      }      }
224      offsetInShared[numNeighbors] = n;      offsetInShared[numNeighbors] = n;
225  #ifdef ESYS_MPI  #ifdef ESYS_MPI
226      mpi_info->incCounter(mpi_info->size);      MPIInfo->incCounter(MPIInfo->size);
227      MPI_Waitall(count, mpi_requests, mpi_stati);      MPI_Waitall(count, &mpi_requests[0], &mpi_stati[0]);
228  #endif  #endif
229      /* map global ids to local id's */      // map global IDs to local IDs
230  #pragma omp parallel for private(i) schedule(static)  #pragma omp parallel for
231      for (i = 0; i < offsetInShared[numNeighbors]; ++i)      for (index_t i = 0; i < n; ++i) {
     {  
232          shared[i] = locDOFMask[shared[i] - min_DOF];          shared[i] = locDOFMask[shared[i] - min_DOF];
233      }      }
234    
235      snd_shcomp.reset(new paso::SharedComponents(myLastDOF - myFirstDOF,      paso::SharedComponents_ptr snd_shcomp(new paso::SharedComponents(
236              numNeighbors, neighbor, shared, offsetInShared, 1, 0,                  myLastDOF - myFirstDOF, numNeighbors, &neighbor[0], shared,
237              dof_distribution->mpi_info));                  &offsetInShared[0], 1, 0, MPIInfo));
238    
239      this_connector.reset(new paso::Connector(snd_shcomp, rcv_shcomp));      degreesOfFreedomConnector.reset(new paso::Connector(snd_shcomp, rcv_shcomp));
240    
     /* assign new DOF labels to nodes */  
     delete[] rcv_len;  
     delete[] snd_len;  
     delete[] mpi_requests;  
     delete[] mpi_stati;  
241      delete[] wanted_DOFs;      delete[] wanted_DOFs;
242      delete[] nodeMask;      delete[] nodeMask;
     delete[] neighbor;  
243      delete[] shared;      delete[] shared;
     delete[] offsetInShared;  
244      delete[] locDOFMask;      delete[] locDOFMask;
     if (use_reduced_elements)  
     {  
         in->Nodes->reducedDegreesOfFreedomMapping = this_mapping;  
         in->Nodes->reducedDegreesOfFreedomConnector = this_connector;  
     }  
     else  
     {  
         in->Nodes->degreesOfFreedomMapping = this_mapping;  
         in->Nodes->degreesOfFreedomConnector = this_connector;  
     }  
245  }  }
246    
247  void Dudley_Mesh_createMappings(Dudley_Mesh* mesh, index_t* dof_distribution, index_t* node_distribution)  void NodeFile::createNodeMappings(const std::vector<index_t>& dofDist,
248                                      const std::vector<index_t>& nodeDist)
249  {  {
250      int i;      // ==== distribution of Nodes ====
251      index_t *maskReducedNodes = NULL, *indexReducedNodes = NULL;      nodesDistribution.reset(new paso::Distribution(
252      dim_t numReducedNodes;                                                  MPIInfo, &nodeDist[0], 1, 0));
253    
254      maskReducedNodes = new index_t[mesh->Nodes->numNodes];      // ==== distribution of DOFs ====
255      indexReducedNodes = new index_t[mesh->Nodes->numNodes];      dofDistribution.reset(new paso::Distribution(MPIInfo, &dofDist[0], 1, 0));
256    
257  #pragma omp parallel for private(i) schedule(static)      index_t* nodeMask = new index_t[numNodes];
258      for (i = 0; i < mesh->Nodes->numNodes; ++i)      const index_t UNUSED = -1;
259          maskReducedNodes[i] = -1;  
260      Dudley_Mesh_markNodes(maskReducedNodes, 0, mesh, true);      // ==== nodes mapping (dummy) ====
261    #pragma omp parallel for
262      numReducedNodes = Dudley_Util_packMask(mesh->Nodes->numNodes, maskReducedNodes, indexReducedNodes);      for (index_t i = 0; i < numNodes; ++i)
     Dudley_Mesh_createNodeFileMappings(mesh, numReducedNodes,  
             indexReducedNodes, dof_distribution, node_distribution);  
     delete[] maskReducedNodes;  
     delete[] indexReducedNodes;  
 }  
   
 void Dudley_Mesh_createNodeFileMappings(Dudley_Mesh* in, dim_t numReducedNodes,  
         index_t* indexReducedNodes, index_t* dof_first_component, index_t* nodes_first_component)  
 {  
   
     index_t myFirstDOF, myLastDOF, myFirstNode, myLastNode, *reduced_dof_first_component = NULL, *nodeMask = NULL,  
         *reduced_nodes_first_component = NULL, k, *maskMyReducedDOF = NULL, *indexMyReducedDOF =  
         NULL, *maskMyReducedNodes = NULL, *indexMyReducedNodes = NULL;  
     dim_t myNumDOF, myNumNodes, myNumReducedNodes, myNumReducedDOF, globalNumReducedNodes, globalNumReducedDOF, i,  
         mpiSize;  
     int myRank;  
   
     mpiSize = in->Nodes->MPIInfo->size;  
     myRank = in->Nodes->MPIInfo->rank;  
   
     /* mark the nodes used by the reduced mesh */  
   
     reduced_dof_first_component = new  index_t[mpiSize + 1];  
     reduced_nodes_first_component = new  index_t[mpiSize + 1];  
   
     myFirstDOF = dof_first_component[myRank];  
     myLastDOF = dof_first_component[myRank + 1];  
     myNumDOF = myLastDOF - myFirstDOF;  
   
     myFirstNode = nodes_first_component[myRank];  
     myLastNode = nodes_first_component[myRank + 1];  
     myNumNodes = myLastNode - myFirstNode;  
   
     maskMyReducedDOF = new index_t[myNumDOF];  
     indexMyReducedDOF = new index_t[myNumDOF];  
     maskMyReducedNodes = new index_t[myNumNodes];  
     indexMyReducedNodes = new index_t[myNumNodes];  
   
 #pragma omp parallel private(i)  
     {  
 #pragma omp for schedule(static)  
         for (i = 0; i < myNumNodes; ++i)  
             maskMyReducedNodes[i] = -1;  
 #pragma omp for schedule(static)  
         for (i = 0; i < myNumDOF; ++i)  
             maskMyReducedDOF[i] = -1;  
 #pragma omp for private(k) schedule(static)  
         for (i = 0; i < numReducedNodes; ++i)  
         {  
             k = in->Nodes->globalNodesIndex[indexReducedNodes[i]];  
             if ((k >= myFirstNode) && (myLastNode > k))  
                 maskMyReducedNodes[k - myFirstNode] = i;  
             k = in->Nodes->globalDegreesOfFreedom[indexReducedNodes[i]];  
             if ((k >= myFirstDOF) && (myLastDOF > k))  
             {  
                 maskMyReducedDOF[k - myFirstDOF] = i;  
             }  
         }  
     }  
     myNumReducedNodes = Dudley_Util_packMask(myNumNodes, maskMyReducedNodes, indexMyReducedNodes);  
     myNumReducedDOF = Dudley_Util_packMask(myNumDOF, maskMyReducedDOF, indexMyReducedDOF);  
   
 #ifdef ESYS_MPI  
     MPI_Allgather(&myNumReducedNodes, 1, MPI_INT, reduced_nodes_first_component, 1, MPI_INT,  
                   in->Nodes->MPIInfo->comm);  
     MPI_Allgather(&myNumReducedDOF, 1, MPI_INT, reduced_dof_first_component, 1, MPI_INT,  
                   in->Nodes->MPIInfo->comm);  
 #else  
     reduced_nodes_first_component[0] = myNumReducedNodes;  
     reduced_dof_first_component[0] = myNumReducedDOF;  
 #endif  
     globalNumReducedNodes = 0;  
     globalNumReducedDOF = 0;  
     for (i = 0; i < mpiSize; ++i)  
     {  
         k = reduced_nodes_first_component[i];  
         reduced_nodes_first_component[i] = globalNumReducedNodes;  
         globalNumReducedNodes += k;  
   
         k = reduced_dof_first_component[i];  
         reduced_dof_first_component[i] = globalNumReducedDOF;  
         globalNumReducedDOF += k;  
     }  
     reduced_nodes_first_component[mpiSize] = globalNumReducedNodes;  
     reduced_dof_first_component[mpiSize] = globalNumReducedDOF;  
     /* ==== distribution of Nodes =============================== */  
     in->Nodes->nodesDistribution.reset(new paso::Distribution(in->Nodes->MPIInfo, nodes_first_component, 1, 0));  
   
     /* ==== distribution of DOFs =============================== */  
     in->Nodes->degreesOfFreedomDistribution.reset(  
             new paso::Distribution(in->Nodes->MPIInfo, dof_first_component, 1, 0));  
   
     /* ==== distribution of reduced Nodes =============================== */  
     in->Nodes->reducedNodesDistribution.reset(new paso::Distribution(  
             in->Nodes->MPIInfo, reduced_nodes_first_component, 1, 0));  
   
     /* ==== distribution of reduced DOF =============================== */  
     in->Nodes->reducedDegreesOfFreedomDistribution.reset(  
             new paso::Distribution(in->Nodes->MPIInfo, reduced_dof_first_component, 1, 0));  
   
     delete[] maskMyReducedDOF;  
     delete[] indexMyReducedDOF;  
     delete[] maskMyReducedNodes;  
     delete[] indexMyReducedNodes;  
     delete[] reduced_dof_first_component;  
     delete[] reduced_nodes_first_component;  
   
     nodeMask = new index_t[in->Nodes->numNodes];  
     /* ==== nodes mapping which is a dummy structure ======== */  
 #pragma omp parallel for private(i) schedule(static)  
     for (i = 0; i < in->Nodes->numNodes; ++i)  
263          nodeMask[i] = i;          nodeMask[i] = i;
264      in->Nodes->nodesMapping = Dudley_NodeMapping_alloc(in->Nodes->numNodes, nodeMask, UNUSED);      nodesMapping.assign(nodeMask, numNodes, UNUSED);
265    
266      /* ==== mapping between nodes and reduced nodes ========== */      // ==== mapping between nodes and DOFs + DOF connector ====
267  #pragma omp parallel for private(i) schedule(static)      createDOFMappingAndCoupling();
     for (i = 0; i < in->Nodes->numNodes; ++i)  
         nodeMask[i] = UNUSED;  
 #pragma omp parallel for private(i) schedule(static)  
     for (i = 0; i < numReducedNodes; ++i)  
         nodeMask[indexReducedNodes[i]] = i;  
     in->Nodes->reducedNodesMapping = Dudley_NodeMapping_alloc(in->Nodes->numNodes, nodeMask, UNUSED);  
     delete[] nodeMask;  
     /* ==== mapping between nodes and DOFs + DOF connector ========== */  
     Dudley_Mesh_createDOFMappingAndCoupling(in, false);  
     /* == mapping between nodes and reduced DOFs + reduced DOF connector == */  
     Dudley_Mesh_createDOFMappingAndCoupling(in, true);  
268    
269      /* get the Ids for DOFs and reduced nodes */      // get the IDs for DOFs
270  #pragma omp parallel private(i)  #pragma omp parallel for
271      {      for (index_t i = 0; i < degreesOfFreedomMapping.numTargets; ++i)
272  #pragma omp for          degreesOfFreedomId[i] = Id[degreesOfFreedomMapping.map[i]];
         for (i = 0; i < in->Nodes->reducedNodesMapping->numTargets; ++i)  
             in->Nodes->reducedNodesId[i] = in->Nodes->Id[in->Nodes->reducedNodesMapping->map[i]];  
 #pragma omp for  
         for (i = 0; i < in->Nodes->degreesOfFreedomMapping->numTargets; ++i)  
             in->Nodes->degreesOfFreedomId[i] = in->Nodes->Id[in->Nodes->degreesOfFreedomMapping->map[i]];  
 #pragma omp for  
         for (i = 0; i < in->Nodes->reducedDegreesOfFreedomMapping->numTargets; ++i)  
             in->Nodes->reducedDegreesOfFreedomId[i] =  
                 in->Nodes->Id[in->Nodes->reducedDegreesOfFreedomMapping->map[i]];  
     }  
273  }  }
274    
275  } // namespace dudley  } // namespace dudley

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

  ViewVC Help
Powered by ViewVC 1.1.26