[688] | 1 | #include <cassert> |
---|
| 2 | #include "node.hpp" |
---|
| 3 | #include "timerRemap.hpp" |
---|
| 4 | #include "circle.hpp" |
---|
| 5 | #include "meshutil.hpp" |
---|
| 6 | #include "polyg.hpp" |
---|
| 7 | #include "gridRemap.hpp" |
---|
| 8 | #include "intersect.hpp" |
---|
| 9 | #include "errhandle.hpp" |
---|
| 10 | #include "mpi_routing.hpp" |
---|
| 11 | #include "misc.hpp" |
---|
| 12 | |
---|
| 13 | #include "parallel_tree.hpp" |
---|
| 14 | |
---|
| 15 | namespace sphereRemap { |
---|
| 16 | |
---|
| 17 | static const int assignLevel = 2; |
---|
| 18 | |
---|
| 19 | // only the circle is packed, rest of node will be initialized on unpacking |
---|
[2269] | 20 | static void packNode(NodePtr node, char *buffer, int& index) |
---|
[688] | 21 | { |
---|
| 22 | if (buffer == NULL) // compute size only |
---|
| 23 | index += 4 * sizeof(double); |
---|
| 24 | else |
---|
| 25 | { |
---|
[2269] | 26 | *(Coord *)(&buffer[index]) = node->centre; |
---|
[688] | 27 | index += sizeof(Coord); |
---|
| 28 | |
---|
[2269] | 29 | *(double *)(&buffer[index]) = node->radius; |
---|
[688] | 30 | index += sizeof(double); |
---|
| 31 | } |
---|
| 32 | } |
---|
| 33 | |
---|
[2269] | 34 | static void unpackNode(NodePtr node, char *buffer, int& index) |
---|
[688] | 35 | { |
---|
| 36 | Coord centre; |
---|
| 37 | double r; |
---|
| 38 | |
---|
| 39 | centre = *(Coord *)(&buffer[index]); |
---|
| 40 | index += sizeof(Coord); |
---|
| 41 | |
---|
| 42 | r = *(double *)(&buffer[index]); |
---|
| 43 | index += sizeof(double); |
---|
| 44 | |
---|
[2269] | 45 | node->centre = centre; |
---|
| 46 | node->radius = r; |
---|
[688] | 47 | } |
---|
| 48 | |
---|
| 49 | |
---|
| 50 | static void packElement(Elt *ptElement, char *buffer, int& index) |
---|
| 51 | { |
---|
| 52 | if (buffer == NULL) index += packedPolygonSize(*ptElement); |
---|
| 53 | else packPolygon(*ptElement, buffer, index); |
---|
| 54 | } |
---|
| 55 | |
---|
| 56 | static void unpackElement(Elt *ptElement, char *buffer, int& index) |
---|
| 57 | { |
---|
| 58 | unpackPolygon(*ptElement, buffer, index); |
---|
| 59 | } |
---|
| 60 | |
---|
| 61 | static void packVector(const vector<int>& vec, char *buf, int& pos) |
---|
| 62 | { |
---|
| 63 | if (buf == NULL) |
---|
| 64 | pos += sizeof(int) + vec.size()*sizeof(int); |
---|
| 65 | else |
---|
| 66 | { |
---|
| 67 | *((int *) &(buf[pos])) = vec.size(); |
---|
| 68 | pos += sizeof(int); |
---|
| 69 | for (int i = 0; i < vec.size(); i++) |
---|
| 70 | { |
---|
| 71 | *((int *) &(buf[pos])) = vec[i]; |
---|
| 72 | pos += sizeof(int); |
---|
| 73 | } |
---|
| 74 | } |
---|
| 75 | } |
---|
| 76 | |
---|
| 77 | static void unpackVector(vector<int>& vec, const char *buf, int& pos) |
---|
| 78 | { |
---|
| 79 | vec.resize(*((int *) &(buf[pos]))); |
---|
| 80 | pos += sizeof(int); |
---|
| 81 | for (int i = 0; i < vec.size(); i++) |
---|
| 82 | { |
---|
| 83 | vec[i] = *((int *) &(buf[pos])); |
---|
| 84 | pos += sizeof(int); |
---|
| 85 | } |
---|
| 86 | } |
---|
| 87 | |
---|
| 88 | static void assignRoute(CSampleTree& tree, const CCascadeLevel& cl) // newroot <- root |
---|
| 89 | { |
---|
| 90 | vector<int> routeRank(cl.group_size); |
---|
| 91 | for (int i = 0; i < cl.group_size; i++) |
---|
| 92 | routeRank[i] = i; //(cl.group_size + cl.polour() - i) % cl.group_size; |
---|
| 93 | std::vector<int>::iterator rank = routeRank.begin(); |
---|
| 94 | tree.root->assignRoute(rank, assignLevel); |
---|
| 95 | } |
---|
| 96 | |
---|
[2269] | 97 | static void transferRoutedNodes(CMPIRouting& MPIRoute, /*const*/ vector<NodePtr>& nodeSend, const vector<int>& route, vector<NodePtr>& nodeRecv) |
---|
[688] | 98 | { |
---|
| 99 | /* `route` knows what goes where */ |
---|
| 100 | MPIRoute.init(route); |
---|
| 101 | int nRecv = MPIRoute.getTotalSourceElement(); |
---|
| 102 | nodeRecv.resize(nRecv); |
---|
[2534] | 103 | for (int i=0;i<nRecv;i++) |
---|
| 104 | { |
---|
| 105 | nodeRecv[i] = nodeRecv[i]->create(); |
---|
| 106 | } |
---|
[688] | 107 | MPIRoute.transferToTarget(&nodeSend[0], &nodeRecv[0], packNode, unpackNode); |
---|
| 108 | } |
---|
| 109 | |
---|
[2269] | 110 | static void transferRoutedIntersections(CMPIRouting& MPIRoute, const vector<NodePtr>& nodeSend, const vector<vector<int> >& route, vector<NodePtr>& nodeRecv) |
---|
[688] | 111 | { |
---|
| 112 | // `route` knows what goes where |
---|
| 113 | MPIRoute.init(route); |
---|
| 114 | int nRecv = MPIRoute.getTotalSourceElement(); |
---|
| 115 | nodeRecv.resize(nRecv); |
---|
[2534] | 116 | for (int i=0;i<nRecv;i++) |
---|
| 117 | { |
---|
| 118 | nodeRecv[i] = nodeRecv[i]->create(); |
---|
| 119 | } |
---|
[2269] | 120 | MPIRoute.transferToTarget((NodePtr*) &nodeSend[0], &nodeRecv[0], packNode, unpackNode); |
---|
[688] | 121 | //cout << MPIRoute.mpiRank << " ROUTE " << nRecv << ": " << nodeSend.size() << " " << nodeRecv.size() << " "; |
---|
| 122 | } |
---|
| 123 | |
---|
[923] | 124 | //CParallelTree::CParallelTree(MPI_Comm comm) : communicator(comm), cascade(MIN_NODE_SZ*MIN_NODE_SZ, comm) |
---|
[1639] | 125 | CParallelTree::CParallelTree(MPI_Comm comm) : communicator(comm), cascade(MAX_NODE_SZ*MAX_NODE_SZ*2, comm) |
---|
[688] | 126 | { |
---|
| 127 | treeCascade.reserve(cascade.num_levels); |
---|
| 128 | for (int lev = 0; lev < cascade.num_levels; lev++) |
---|
| 129 | treeCascade.push_back(CSampleTree(cascade.level[lev].group_size, assignLevel)); |
---|
| 130 | } |
---|
| 131 | |
---|
[2269] | 132 | void CParallelTree::buildSampleTreeCascade(vector<NodePtr>& sampleNodes /*route field will be modified*/, int level) |
---|
[688] | 133 | { |
---|
| 134 | buildSampleTree(treeCascade[level], sampleNodes, cascade.level[level]); |
---|
| 135 | assignRoute(treeCascade[level] /*out*/, cascade.level[level] /*in*/); |
---|
| 136 | |
---|
| 137 | if (level+1 < cascade.num_levels) |
---|
| 138 | { |
---|
| 139 | vector<int> route(sampleNodes.size()); |
---|
| 140 | treeCascade[level].routeNodes(route, sampleNodes, assignLevel); |
---|
| 141 | |
---|
[2269] | 142 | vector<NodePtr> routedNodes; |
---|
[688] | 143 | CMPIRouting mpiRoute(cascade.level[level].pg_comm); |
---|
| 144 | transferRoutedNodes(mpiRoute, sampleNodes, route, routedNodes); |
---|
| 145 | buildSampleTreeCascade(routedNodes, level+1); |
---|
| 146 | } |
---|
| 147 | } |
---|
| 148 | |
---|
[2269] | 149 | void buildSampleTree(CSampleTree& tree, const vector<NodePtr>& node, const CCascadeLevel& comm) |
---|
[688] | 150 | /* |
---|
| 151 | In the beginning all the sample elements are distributed |
---|
| 152 | -> communicate to make available at each rank |
---|
| 153 | so that each rank can build the same sample tree |
---|
| 154 | */ |
---|
| 155 | { |
---|
| 156 | int n = node.size(); // number of samples intially on this rank (local) |
---|
| 157 | |
---|
| 158 | int blocSize = comm.group_size * ipow(MAX_NODE_SZ, assignLevel); |
---|
| 159 | |
---|
| 160 | int nrecv; // global number of samples THIS WILL BE THE NUMBER OF LEAFS IN THE SAMPLE TREE |
---|
[1639] | 161 | MPI_Allreduce(&n, &nrecv, 1, MPI_INT, MPI_SUM, comm.comm); // => size of sample tree does not depend on keepNodes! |
---|
[688] | 162 | double ratio = blocSize / (1.0 * nrecv); |
---|
| 163 | int nsend = ratio * n + 1; // nsend = n_local_samples / n_global_samples * blocksize + 1 = blocksize/comm.size |
---|
| 164 | if (nsend > n) nsend = n; |
---|
| 165 | |
---|
| 166 | int *counts = new int[comm.size]; |
---|
[1639] | 167 | MPI_Allgather(&nsend, 1, MPI_INT, counts, 1, MPI_INT, comm.comm); |
---|
[688] | 168 | |
---|
| 169 | nrecv = 0; |
---|
| 170 | int *displs = new int[comm.size]; |
---|
| 171 | for (int i = 0; i < comm.size; i++) |
---|
| 172 | { |
---|
| 173 | displs[i] = 4 * nrecv; |
---|
| 174 | nrecv += counts[i]; |
---|
| 175 | counts[i] = 4 * counts[i]; |
---|
| 176 | } |
---|
| 177 | |
---|
| 178 | /* pack circles around sample elements */ |
---|
| 179 | double *sendBuffer = new double[nsend*4]; |
---|
| 180 | int index = 0; |
---|
| 181 | vector<int> randomArray(n); |
---|
| 182 | randomizeArray(randomArray); |
---|
| 183 | for (int i = 0; i < nsend; i++) |
---|
| 184 | { |
---|
[2269] | 185 | const Node& no = *node[randomArray[i]]; |
---|
[688] | 186 | *((Coord *) (sendBuffer + index)) = no.centre; |
---|
| 187 | index += sizeof(Coord)/sizeof(*sendBuffer); |
---|
| 188 | sendBuffer[index++] = no.radius; |
---|
| 189 | } |
---|
| 190 | |
---|
| 191 | /* each process needs the sample elements from all processes */ |
---|
| 192 | double *recvBuffer = new double[nrecv*4]; |
---|
[1639] | 193 | MPI_Allgatherv(sendBuffer, 4 * nsend, MPI_DOUBLE, recvBuffer, counts, displs, MPI_DOUBLE, comm.comm); |
---|
[688] | 194 | delete[] sendBuffer; |
---|
| 195 | delete[] counts; |
---|
| 196 | delete[] displs; |
---|
| 197 | |
---|
| 198 | /* unpack */ |
---|
[923] | 199 | /* |
---|
[688] | 200 | randomArray.resize(nrecv); |
---|
| 201 | randomizeArray(randomArray); |
---|
| 202 | tree.leafs.resize(nrecv); |
---|
| 203 | index = 0; |
---|
[923] | 204 | |
---|
| 205 | |
---|
| 206 | for (int i = 0; i < nrecv; i++) |
---|
[688] | 207 | { |
---|
| 208 | Coord x = *(Coord *)(&recvBuffer[index]); |
---|
| 209 | index += sizeof(Coord)/sizeof(*recvBuffer); |
---|
| 210 | double radius = recvBuffer[index++]; |
---|
| 211 | tree.leafs[randomArray[i]].centre = x; |
---|
| 212 | tree.leafs[randomArray[i]].radius = radius; |
---|
| 213 | |
---|
| 214 | } |
---|
[923] | 215 | */ |
---|
[688] | 216 | |
---|
[923] | 217 | randomArray.resize(blocSize); |
---|
| 218 | randomizeArray(randomArray); |
---|
| 219 | tree.leafs.resize(blocSize); |
---|
[2269] | 220 | for(auto& it : tree.leafs) it=make_shared<Node>() ; |
---|
| 221 | |
---|
[923] | 222 | index = 0; |
---|
| 223 | |
---|
| 224 | size_t s=(sizeof(Coord)/sizeof(*recvBuffer)+1)*nrecv ; |
---|
| 225 | |
---|
| 226 | for (int i = 0; i < blocSize; i++) |
---|
| 227 | { |
---|
| 228 | Coord x = *(Coord *)(&recvBuffer[index%s]); |
---|
| 229 | index += sizeof(Coord)/sizeof(*recvBuffer); |
---|
| 230 | double radius = recvBuffer[index%s]; |
---|
| 231 | index++ ; |
---|
[2269] | 232 | tree.leafs[randomArray[i]]->centre = x; |
---|
| 233 | tree.leafs[randomArray[i]]->radius = radius; |
---|
[923] | 234 | |
---|
| 235 | } |
---|
| 236 | |
---|
| 237 | |
---|
[688] | 238 | delete [] recvBuffer; |
---|
| 239 | |
---|
| 240 | CTimer::get("buildSampleTree(local)").resume(); |
---|
| 241 | tree.build(tree.leafs); |
---|
[923] | 242 | // cout << "SampleTree build : assign Level " << assignLevel << " nb Nodes : " << tree.levelSize[assignLevel] << endl; |
---|
[688] | 243 | CTimer::get("buildSampleTree(local)").suspend(); |
---|
| 244 | CTimer::get("buildSampleTree(local)").print(); |
---|
| 245 | /* End gracefully if sample tree could not be constructed with desired number of nodes on assignLevel */ |
---|
| 246 | int allok, ok = (tree.levelSize[assignLevel] == comm.group_size); |
---|
| 247 | if (!ok) |
---|
| 248 | { |
---|
[923] | 249 | cerr << comm.rank << ": PROBLEM: (node assign)" << tree.levelSize[assignLevel] << " != " << comm.group_size << " (keepNodes)" |
---|
| 250 | << " node size : "<<node.size()<<" bloc size : "<<blocSize<<" total number of leaf : "<<tree.leafs.size()<<endl ; |
---|
[688] | 251 | /* |
---|
[1639] | 252 | MPI_Allreduce(&ok, &allok, 1, MPI_INT, MPI_PROD, communicator); |
---|
[688] | 253 | if (!allok) { |
---|
| 254 | MPI_Finalize(); |
---|
| 255 | exit(1); |
---|
| 256 | } |
---|
| 257 | */ |
---|
[1639] | 258 | MPI_Abort(MPI_COMM_WORLD,-1) ; |
---|
[688] | 259 | } |
---|
[923] | 260 | /* |
---|
[688] | 261 | cout<<"*******************************************"<<endl ; |
---|
| 262 | cout<<"******* Sample Tree output *********"<<endl ; |
---|
| 263 | cout<<"*******************************************"<<endl ; |
---|
| 264 | tree.output(cout,1) ; |
---|
| 265 | |
---|
| 266 | cout<<"*******************************************"<<endl ; |
---|
[923] | 267 | */ |
---|
[688] | 268 | |
---|
| 269 | assert(tree.root->incluCheck() == 0); |
---|
| 270 | } |
---|
| 271 | |
---|
| 272 | |
---|
[2269] | 273 | void CParallelTree::buildLocalTree(const vector<NodePtr>& node, const vector<int>& route) |
---|
[688] | 274 | { |
---|
| 275 | CMPIRouting MPIRoute(communicator); |
---|
[1639] | 276 | MPI_Barrier(communicator); |
---|
[688] | 277 | CTimer::get("buildLocalTree(initRoute)").resume(); |
---|
| 278 | MPIRoute.init(route); |
---|
| 279 | CTimer::get("buildLocalTree(initRoute)").suspend(); |
---|
| 280 | CTimer::get("buildLocalTree(initRoute)").print(); |
---|
| 281 | |
---|
| 282 | nbLocalElements = MPIRoute.getTotalSourceElement(); |
---|
| 283 | localElements = new Elt[nbLocalElements]; |
---|
| 284 | |
---|
| 285 | vector<Elt*> ptElement(node.size()); |
---|
| 286 | for (int i = 0; i < node.size(); i++) |
---|
[2269] | 287 | ptElement[i] = (Elt *) (node[i]->data); |
---|
[688] | 288 | |
---|
| 289 | vector<Elt*> ptLocalElement(nbLocalElements); |
---|
| 290 | for (int i = 0; i < nbLocalElements; i++) |
---|
| 291 | ptLocalElement[i] = &localElements[i]; |
---|
| 292 | |
---|
| 293 | CTimer::get("buildLocalTree(transfer)").resume(); |
---|
| 294 | MPIRoute.transferToTarget(&ptElement[0], &ptLocalElement[0], packElement, unpackElement); |
---|
| 295 | CTimer::get("buildLocalTree(transfer)").suspend(); |
---|
| 296 | CTimer::get("buildLocalTree(transfer)").print(); |
---|
| 297 | |
---|
| 298 | CTimer::get("buildLocalTree(local)").resume(); |
---|
| 299 | |
---|
| 300 | int mpiRank; |
---|
[1639] | 301 | MPI_Comm_rank(communicator, &mpiRank); |
---|
[688] | 302 | localTree.leafs.reserve(nbLocalElements); |
---|
| 303 | for (int i = 0; i < nbLocalElements; i++) |
---|
| 304 | { |
---|
| 305 | Elt& elt = localElements[i]; |
---|
| 306 | elt.id.ind = i; |
---|
| 307 | elt.id.rank = mpiRank; |
---|
[2269] | 308 | localTree.leafs.push_back(make_shared<Node>(elt.x, cptRadius(elt), &localElements[i])); |
---|
[688] | 309 | } |
---|
| 310 | localTree.build(localTree.leafs); |
---|
| 311 | |
---|
| 312 | cptAllEltsGeom(localElements, nbLocalElements, srcGrid.pole); |
---|
| 313 | CTimer::get("buildLocalTree(local)").suspend(); |
---|
| 314 | CTimer::get("buildLocalTree(local)").print(); |
---|
| 315 | } |
---|
| 316 | |
---|
[2269] | 317 | void CParallelTree::build(vector<NodePtr>& node, vector<NodePtr>& node2) |
---|
[688] | 318 | { |
---|
| 319 | |
---|
| 320 | int assignLevel = 2; |
---|
| 321 | int nbSampleNodes = 2*ipow(MAX_NODE_SZ + 1, assignLevel); |
---|
| 322 | |
---|
[852] | 323 | |
---|
| 324 | long int nb1, nb2, nb, nbTot ; |
---|
[898] | 325 | nb1=node.size() ; nb2=node2.size() ; |
---|
[852] | 326 | nb=nb1+nb2 ; |
---|
[1639] | 327 | MPI_Allreduce(&nb, &nbTot, 1, MPI_LONG, MPI_SUM, communicator) ; |
---|
[852] | 328 | int commSize ; |
---|
[1639] | 329 | MPI_Comm_size(communicator,&commSize) ; |
---|
[852] | 330 | |
---|
[688] | 331 | // make multiple of two |
---|
| 332 | nbSampleNodes /= 2; |
---|
| 333 | nbSampleNodes *= 2; |
---|
[898] | 334 | // assert( nbTot > nbSampleNodes*commSize) ; |
---|
[852] | 335 | |
---|
[853] | 336 | int nbSampleNodes1 = nbSampleNodes * (nb1*commSize)/(1.*nbTot) ; |
---|
| 337 | int nbSampleNodes2 = nbSampleNodes * (nb2*commSize)/(1.*nbTot) ; |
---|
[852] | 338 | |
---|
[688] | 339 | |
---|
[812] | 340 | // assert(node.size() > nbSampleNodes); |
---|
| 341 | // assert(node2.size() > nbSampleNodes); |
---|
[852] | 342 | // assert(node.size() + node2.size() > nbSampleNodes); |
---|
[2269] | 343 | vector<NodePtr> sampleNodes; sampleNodes.reserve(nbSampleNodes1+nbSampleNodes2); |
---|
[688] | 344 | |
---|
| 345 | vector<int> randomArray1(node.size()); |
---|
| 346 | randomizeArray(randomArray1); |
---|
| 347 | vector<int> randomArray2(node2.size()); |
---|
| 348 | randomizeArray(randomArray2); |
---|
| 349 | |
---|
[852] | 350 | /* |
---|
[812] | 351 | int s1,s2 ; |
---|
| 352 | |
---|
| 353 | if (node.size()< nbSampleNodes/2) |
---|
| 354 | { |
---|
| 355 | s1 = node.size() ; |
---|
| 356 | s2 = nbSampleNodes-s1 ; |
---|
| 357 | } |
---|
| 358 | else if (node2.size()< nbSampleNodes/2) |
---|
| 359 | { |
---|
| 360 | s2 = node.size() ; |
---|
| 361 | s1 = nbSampleNodes-s2 ; |
---|
| 362 | } |
---|
| 363 | else |
---|
| 364 | { |
---|
| 365 | s1=nbSampleNodes/2 ; |
---|
| 366 | s2=nbSampleNodes/2 ; |
---|
| 367 | } |
---|
[852] | 368 | */ |
---|
[2269] | 369 | for (int i = 0; i <nbSampleNodes1; i++) sampleNodes.push_back(make_shared<Node>(node[randomArray1[i%nb1]]->centre, node[randomArray1[i%nb1]]->radius, nullptr)); |
---|
| 370 | for (int i = 0; i <nbSampleNodes2; i++) sampleNodes.push_back(make_shared<Node>(node2[randomArray2[i%nb2]]->centre, node2[randomArray2[i%nb2]]->radius, nullptr)); |
---|
[812] | 371 | |
---|
| 372 | /* |
---|
| 373 | for (int i = 0; i < nbSampleNodes/2; i++) |
---|
[688] | 374 | { |
---|
[812] | 375 | sampleNodes.push_back(Node(node[randomArray1[i]].centre, node[randomArray1[i]].radius, NULL)); |
---|
| 376 | sampleNodes.push_back(Node(node2[randomArray2[i]].centre, node2[randomArray2[i]].radius, NULL)); |
---|
[688] | 377 | } |
---|
[812] | 378 | */ |
---|
[688] | 379 | CTimer::get("buildParallelSampleTree").resume(); |
---|
| 380 | //sampleTree.buildParallelSampleTree(sampleNodes, cascade); |
---|
| 381 | buildSampleTreeCascade(sampleNodes); |
---|
| 382 | CTimer::get("buildParallelSampleTree").suspend(); |
---|
| 383 | CTimer::get("buildParallelSampleTree").print(); |
---|
| 384 | |
---|
| 385 | //route source mesh |
---|
| 386 | CTimer::get("parallelRouteNode").resume(); |
---|
| 387 | vector<int> route(node.size()); |
---|
| 388 | routeNodes(route /*out*/, node); |
---|
| 389 | CTimer::get("parallelRouteNode").suspend(); |
---|
| 390 | CTimer::get("parallelRouteNode").print(); |
---|
| 391 | |
---|
| 392 | CTimer::get("buildLocalTree").resume(); |
---|
| 393 | buildLocalTree(node, route); |
---|
| 394 | CTimer::get("buildLocalTree").suspend(); |
---|
| 395 | CTimer::get("buildLocalTree").print(); |
---|
| 396 | |
---|
| 397 | CTimer::get("buildRouteTree").resume(); |
---|
| 398 | /* update circles of tree cascade so it can be used for routing */ |
---|
| 399 | updateCirclesForRouting(localTree.root->centre, localTree.root->radius); |
---|
| 400 | CTimer::get("buildRouteTree").suspend(); |
---|
| 401 | CTimer::get("buildRouteTree").print(); |
---|
| 402 | } |
---|
| 403 | |
---|
[2269] | 404 | void CParallelTree::routeNodes(vector<int>& route, vector<NodePtr>& nodes /*route field used*/, int level) |
---|
[688] | 405 | { |
---|
| 406 | treeCascade[level].routeNodes(route /*assign*/, nodes, assignLevel); |
---|
| 407 | |
---|
| 408 | if (level+1 < cascade.num_levels) |
---|
| 409 | { |
---|
[2269] | 410 | vector<NodePtr> routedNodes; |
---|
[688] | 411 | CMPIRouting MPIRoute(cascade.level[level].pg_comm); |
---|
| 412 | transferRoutedNodes(MPIRoute, nodes, route /*use*/, routedNodes); |
---|
| 413 | vector<int> globalRank(routedNodes.size()); |
---|
| 414 | routeNodes(globalRank, routedNodes, level + 1); |
---|
| 415 | MPIRoute.transferFromSource(&route[0] /*override*/, &globalRank[0]); |
---|
| 416 | } |
---|
| 417 | else |
---|
| 418 | { |
---|
| 419 | CMPIRouting MPIRoute(cascade.level[level].comm); // or use pg_comm, on last cascade level identical |
---|
| 420 | MPIRoute.init(route); |
---|
| 421 | int nbRecvNode = MPIRoute.getTotalSourceElement(); |
---|
| 422 | vector<int> globalRank(nbRecvNode); |
---|
| 423 | for (int i = 0; i < globalRank.size(); i++) |
---|
| 424 | globalRank[i] = cascade.level[0].rank; |
---|
| 425 | MPIRoute.transferFromSource(&route[0] /*override*/, &globalRank[0]); |
---|
| 426 | } |
---|
| 427 | } |
---|
| 428 | |
---|
| 429 | /* assume `to` to be empty vector at entry */ |
---|
| 430 | void linearize(const vector<vector<int> >& from, vector<int>& to) |
---|
| 431 | { |
---|
| 432 | int cnt = 0; |
---|
| 433 | for (int i = 0; i < from.size(); ++i) |
---|
| 434 | cnt += from[i].size(); |
---|
| 435 | to.resize(cnt); |
---|
| 436 | vector<int>::iterator dest = to.begin(); |
---|
| 437 | for (int i = 0; i < from.size(); ++i) |
---|
| 438 | dest = copy(from[i].begin(), from[i].end(), dest); |
---|
| 439 | } |
---|
| 440 | |
---|
| 441 | /* at entry `to` must already have it's final shape and only values are overriden */ |
---|
| 442 | void delinearize(const vector<int>& from, vector<vector<int> >& to) |
---|
| 443 | { |
---|
| 444 | vector<int>::const_iterator end, src = from.begin(); |
---|
| 445 | for (int i = 0; i < to.size(); ++i, src=end) |
---|
| 446 | copy(src, end = src + to[i].size(), to[i].begin()); |
---|
| 447 | } |
---|
| 448 | |
---|
[2269] | 449 | void CParallelTree::routeIntersections(vector<vector<int> >& routes, vector<NodePtr>& nodes, int level) |
---|
[688] | 450 | { |
---|
| 451 | treeCascade[level].routeIntersections(routes /*assign*/, nodes); |
---|
| 452 | |
---|
| 453 | if (level+1 < cascade.num_levels) |
---|
| 454 | { |
---|
[2269] | 455 | vector<NodePtr> routedNodes; |
---|
[688] | 456 | CMPIRouting MPIRoute(cascade.level[level].pg_comm); |
---|
| 457 | |
---|
| 458 | vector<int> flattenedRoutes1; |
---|
| 459 | linearize(routes, flattenedRoutes1); |
---|
[2269] | 460 | vector<NodePtr> double_nodes(flattenedRoutes1.size()); |
---|
[688] | 461 | int j = 0; |
---|
| 462 | for (int i = 0; i < routes.size(); ++i) |
---|
| 463 | for (int k = 0; k < routes[i].size(); ++k, ++j) |
---|
[2269] | 464 | { |
---|
| 465 | double_nodes[j] = make_shared<Node>() ; |
---|
| 466 | *double_nodes[j] = *nodes[i]; |
---|
| 467 | } |
---|
[688] | 468 | transferRoutedNodes(MPIRoute, double_nodes, flattenedRoutes1 /*use*/, routedNodes); |
---|
| 469 | vector<vector<int> > globalRanks(routedNodes.size()); |
---|
| 470 | routeIntersections(globalRanks /*out*/, routedNodes /*in*/, level + 1); |
---|
| 471 | vector<vector<int> > flattenedRoutes(flattenedRoutes1.size()); |
---|
| 472 | // transferFromSource expects sizes (nbSendNode=flattenedRoutes, nbRecvNode=routedNodes.size()) |
---|
| 473 | MPIRoute.transferFromSource(&flattenedRoutes[0], &globalRanks[0], packVector, unpackVector); |
---|
| 474 | for (int i = 0, j = 0; i < routes.size(); ++i) |
---|
| 475 | { |
---|
| 476 | int old_size = routes[i].size(); |
---|
| 477 | routes[i].resize(0); |
---|
| 478 | for (int k = 0; k < old_size; ++k, ++j) |
---|
| 479 | for (int l = 0; l < flattenedRoutes[j].size(); ++l) |
---|
| 480 | routes[i].push_back(flattenedRoutes[j][l]); |
---|
| 481 | } |
---|
| 482 | assert(j == flattenedRoutes1.size()); |
---|
| 483 | |
---|
| 484 | } |
---|
| 485 | else |
---|
| 486 | { |
---|
| 487 | CMPIRouting MPIRoute(cascade.level[level].comm); |
---|
| 488 | MPIRoute.init(routes); |
---|
| 489 | int nbRecvNode = MPIRoute.getTotalSourceElement(); |
---|
| 490 | vector<int> globalRanks(nbRecvNode, cascade.level[0].rank); |
---|
| 491 | vector<int> flattened_routes; |
---|
| 492 | linearize(routes, flattened_routes); |
---|
| 493 | MPIRoute.transferFromSource(&flattened_routes[0], &globalRanks[0]); |
---|
| 494 | delinearize(flattened_routes, routes); |
---|
| 495 | } |
---|
| 496 | if (level!=level) |
---|
| 497 | { |
---|
| 498 | for (int i = 0; i < routes.size(); ++i) |
---|
| 499 | for (int k = 0; k < routes[i].size(); ++k) |
---|
| 500 | if (routes[i][k] == cascade.level[0].rank) routes[i].erase(routes[i].begin() + (k--)); |
---|
| 501 | } |
---|
| 502 | } |
---|
| 503 | |
---|
| 504 | void CParallelTree::updateCirclesForRouting(Coord rootCentre, double rootRadius, int level) |
---|
| 505 | { |
---|
| 506 | if (level + 1 < cascade.num_levels) // children in the cascade have to update first |
---|
| 507 | { |
---|
| 508 | updateCirclesForRouting(rootCentre, rootRadius, level + 1); |
---|
| 509 | rootCentre = treeCascade[level+1].root->centre; |
---|
| 510 | rootRadius = treeCascade[level+1].root->radius; |
---|
| 511 | } |
---|
| 512 | |
---|
| 513 | // gather circles on this level of the cascade |
---|
| 514 | int pg_size; |
---|
[1639] | 515 | MPI_Comm_size(cascade.level[level].pg_comm, &pg_size); |
---|
[688] | 516 | vector<Coord> allRootCentres(pg_size); |
---|
| 517 | vector<double> allRootRadia(pg_size); |
---|
[1639] | 518 | MPI_Allgather(&rootCentre, 3, MPI_DOUBLE, &allRootCentres[0], 3, MPI_DOUBLE, cascade.level[level].pg_comm); |
---|
| 519 | MPI_Allgather(&rootRadius, 1, MPI_DOUBLE, &allRootRadia[0], 1, MPI_DOUBLE, cascade.level[level].pg_comm); |
---|
[688] | 520 | |
---|
| 521 | // now allRootsRadia and allRootCentres must be inserted into second levels of us and propagated to root |
---|
| 522 | treeCascade[level].root->assignCircleAndPropagateUp(&allRootCentres[0], &allRootRadia[0], assignLevel); |
---|
| 523 | } |
---|
| 524 | |
---|
| 525 | CParallelTree::~CParallelTree() |
---|
| 526 | { |
---|
| 527 | delete [] localElements; |
---|
| 528 | } |
---|
| 529 | |
---|
| 530 | } |
---|