src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26 #include "mapDistribute.H"
  27 #include "commSchedule.H"
  28 #include "HashSet.H"
  29 #include "globalIndex.H"
  30 #include "globalIndexAndTransform.H"
  31 #include "transformField.H"
  32 #include "ListOps.H"
  33
  34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  35
  36 defineTypeNameAndDebug(Foam::mapDistribute, 0);
  37
  38 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  39
  40 template<>
  41 void Foam::mapDistribute::transform::operator()
  42 (
  43     const vectorTensorTransform&,
  44     const bool,
  45     List<label>&
  46 ) const
  47 {}
  48 template<>
  49 void Foam::mapDistribute::transform::operator()
  50 (
  51     const coupledPolyPatch&,
  52     UList<label>&
  53 ) const
  54 {}
  55 template<>
  56 void Foam::mapDistribute::transform::operator()
  57 (
  58     const coupledPolyPatch&,
  59     Map<label>&
  60 ) const
  61 {}
  62 template<>
  63 void Foam::mapDistribute::transform::operator()
  64 (
  65     const coupledPolyPatch&,
  66     EdgeMap<label>&
  67 ) const
  68 {}
  69
  70
  71 template<>
  72 void Foam::mapDistribute::transform::operator()
  73 (
  74     const vectorTensorTransform&,
  75     const bool,
  76     List<scalar>&
  77 ) const
  78 {}
  79 template<>
  80 void Foam::mapDistribute::transform::operator()
  81 (
  82     const coupledPolyPatch&,
  83     UList<scalar>&
  84 ) const
  85 {}
  86 template<>
  87 void Foam::mapDistribute::transform::operator()
  88 (
  89     const coupledPolyPatch&,
  90     Map<scalar>&
  91 ) const
  92 {}
  93 template<>
  94 void Foam::mapDistribute::transform::operator()
  95 (
  96     const coupledPolyPatch&,
  97     EdgeMap<scalar>&
  98 ) const
  99 {}
 100
 101
 102 template<>
 103 void Foam::mapDistribute::transform::operator()
 104 (
 105     const vectorTensorTransform&,
 106     const bool,
 107     List<bool>&
 108 ) const
 109 {}
 110 template<>
 111 void Foam::mapDistribute::transform::operator()
 112 (
 113     const coupledPolyPatch&,
 114     UList<bool>&
 115 ) const
 116 {}
 117 template<>
 118 void Foam::mapDistribute::transform::operator()
 119 (
 120     const coupledPolyPatch&,
 121     Map<bool>&
 122 ) const
 123 {}
 124 template<>
 125 void Foam::mapDistribute::transform::operator()
 126 (
 127     const coupledPolyPatch&,
 128     EdgeMap<bool>&
 129 ) const
 130 {}
 131
 132
 133 Foam::List<Foam::labelPair> Foam::mapDistribute::schedule
 134 (
 135     const labelListList& subMap,
 136     const labelListList& constructMap
 137 )
 138 {
 139     // Communications: send and receive processor
 140     List<labelPair> allComms;
 141
 142     {
 143         HashSet<labelPair, labelPair::Hash<> > commsSet(Pstream::nProcs());
 144
 145         // Find what communication is required
 146         forAll(subMap, procI)
 147         {
 148             if (procI != Pstream::myProcNo())
 149             {
 150                 if (subMap[procI].size())
 151                 {
 152                     // I need to send to procI
 153                     commsSet.insert(labelPair(Pstream::myProcNo(), procI));
 154                 }
 155                 if (constructMap[procI].size())
 156                 {
 157                     // I need to receive from procI
 158                     commsSet.insert(labelPair(procI, Pstream::myProcNo()));
 159                 }
 160             }
 161         }
 162         allComms = commsSet.toc();
 163     }
 164
 165
 166     // Reduce
 167     if (Pstream::master())
 168     {
 169         // Receive and merge
 170         for
 171         (
 172             int slave=Pstream::firstSlave();
 173             slave<=Pstream::lastSlave();
 174             slave++
 175         )
 176         {
 177             IPstream fromSlave(Pstream::scheduled, slave);
 178             List<labelPair> nbrData(fromSlave);
 179
 180             forAll(nbrData, i)
 181             {
 182                 if (findIndex(allComms, nbrData[i]) == -1)
 183                 {
 184                     label sz = allComms.size();
 185                     allComms.setSize(sz+1);
 186                     allComms[sz] = nbrData[i];
 187                 }
 188             }
 189         }
 190         // Send back
 191         for
 192         (
 193             int slave=Pstream::firstSlave();
 194             slave<=Pstream::lastSlave();
 195             slave++
 196         )
 197         {
 198             OPstream toSlave(Pstream::scheduled, slave);
 199             toSlave << allComms;
 200         }
 201     }
 202     else
 203     {
 204         {
 205             OPstream toMaster(Pstream::scheduled, Pstream::masterNo());
 206             toMaster << allComms;
 207         }
 208         {
 209             IPstream fromMaster(Pstream::scheduled, Pstream::masterNo());
 210             fromMaster >> allComms;
 211         }
 212     }
 213
 214
 215     // Determine my schedule.
 216     labelList mySchedule
 217     (
 218         commSchedule
 219         (
 220             Pstream::nProcs(),
 221             allComms
 222         ).procSchedule()[Pstream::myProcNo()]
 223     );
 224
 225     // Processors involved in my schedule
 226     return List<labelPair>(UIndirectList<labelPair>(allComms, mySchedule));
 227
 228
 229     //if (debug)
 230     //{
 231     //    Pout<< "I need to:" << endl;
 232     //    const List<labelPair>& comms = schedule();
 233     //    forAll(comms, i)
 234     //    {
 235     //        const labelPair& twoProcs = comms[i];
 236     //        label sendProc = twoProcs[0];
 237     //        label recvProc = twoProcs[1];
 238     //
 239     //        if (recvProc == Pstream::myProcNo())
 240     //        {
 241     //            Pout<< "    receive from " << sendProc << endl;
 242     //        }
 243     //        else
 244     //        {
 245     //            Pout<< "    send to " << recvProc << endl;
 246     //        }
 247     //    }
 248     //}
 249 }
 250
 251
 252 const Foam::List<Foam::labelPair>& Foam::mapDistribute::schedule() const
 253 {
 254     if (schedulePtr_.empty())
 255     {
 256         schedulePtr_.reset
 257         (
 258             new List<labelPair>
 259             (
 260                 schedule(subMap_, constructMap_)
 261             )
 262         );
 263     }
 264     return schedulePtr_();
 265 }
 266
 267
 268 void Foam::mapDistribute::checkReceivedSize
 269 (
 270     const label procI,
 271     const label expectedSize,
 272     const label receivedSize
 273 )
 274 {
 275     if (receivedSize != expectedSize)
 276     {
 277         FatalErrorIn
 278         (
 279             "template<class T>\n"
 280             "void mapDistribute::distribute\n"
 281             "(\n"
 282             "    const Pstream::commsTypes commsType,\n"
 283             "    const List<labelPair>& schedule,\n"
 284             "    const label constructSize,\n"
 285             "    const labelListList& subMap,\n"
 286             "    const labelListList& constructMap,\n"
 287             "    List<T>& field\n"
 288             ")\n"
 289         )   << "Expected from processor " << procI
 290             << " " << expectedSize << " but received "
 291             << receivedSize << " elements."
 292             << abort(FatalError);
 293     }
 294 }
 295
 296
 297 void Foam::mapDistribute::printLayout(Ostream& os) const
 298 {
 299     // Determine offsets of remote data.
 300     labelList minIndex(Pstream::nProcs(), labelMax);
 301     labelList maxIndex(Pstream::nProcs(), labelMin);
 302     forAll(constructMap_, procI)
 303     {
 304         const labelList& construct = constructMap_[procI];
 305         minIndex[procI] = min(minIndex[procI], min(construct));
 306         maxIndex[procI] = max(maxIndex[procI], max(construct));
 307     }
 308
 309     label localSize;
 310     if (maxIndex[Pstream::myProcNo()] == labelMin)
 311     {
 312         localSize = 0;
 313     }
 314     else
 315     {
 316         localSize = maxIndex[Pstream::myProcNo()]+1;
 317     }
 318
 319     os  << "Layout: (constructSize:" << constructSize_ << ")" << endl
 320         << "local (processor " << Pstream::myProcNo() << "):" << endl
 321         << "    start : 0" << endl
 322         << "    size  : " << localSize << endl;
 323
 324     label offset = localSize;
 325     forAll(minIndex, procI)
 326     {
 327         if (procI != Pstream::myProcNo())
 328         {
 329             if (constructMap_[procI].size() > 0)
 330             {
 331                 if (minIndex[procI] != offset)
 332                 {
 333                     FatalErrorIn("mapDistribute::printLayout(..)")
 334                         << "offset:" << offset
 335                         << " procI:" << procI
 336                         << " minIndex:" << minIndex[procI]
 337                         << abort(FatalError);
 338                 }
 339
 340                 label size = maxIndex[procI]-minIndex[procI]+1;
 341                 os  << "processor " << procI << ':' << endl
 342                     << "    start : " << offset << endl
 343                     << "    size  : " << size << endl;
 344
 345                 offset += size;
 346             }
 347         }
 348     }
 349     forAll(transformElements_, trafoI)
 350     {
 351         if (transformElements_[trafoI].size() > 0)
 352         {
 353             os  << "transform " << trafoI << ':' << endl
 354                 << "    start : " << transformStart_[trafoI] << endl
 355                 << "    size  : " << transformElements_[trafoI].size() << endl;
 356         }
 357     }
 358 }
 359
 360
 361 // Construct per processor compact addressing of the global elements
 362 // needed. The ones from the local processor are not included since
 363 // these are always all needed.
 364 void Foam::mapDistribute::calcCompactAddressing
 365 (
 366     const globalIndex& globalNumbering,
 367     const labelList& elements,
 368     List<Map<label> >& compactMap
 369 ) const
 370 {
 371     compactMap.setSize(Pstream::nProcs());
 372
 373     // Count all (non-local) elements needed. Just for presizing map.
 374     labelList nNonLocal(Pstream::nProcs(), 0);
 375
 376     forAll(elements, i)
 377     {
 378         label globalIndex = elements[i];
 379
 380         if (globalIndex != -1 && !globalNumbering.isLocal(globalIndex))
 381         {
 382             label procI = globalNumbering.whichProcID(globalIndex);
 383             nNonLocal[procI]++;
 384         }
 385     }
 386
 387     forAll(compactMap, procI)
 388     {
 389         compactMap[procI].clear();
 390         if (procI != Pstream::myProcNo())
 391         {
 392             compactMap[procI].resize(2*nNonLocal[procI]);
 393         }
 394     }
 395
 396
 397     // Collect all (non-local) elements needed.
 398     forAll(elements, i)
 399     {
 400         label globalIndex = elements[i];
 401
 402         if (globalIndex != -1 && !globalNumbering.isLocal(globalIndex))
 403         {
 404             label procI = globalNumbering.whichProcID(globalIndex);
 405             label index = globalNumbering.toLocal(procI, globalIndex);
 406             label nCompact = compactMap[procI].size();
 407             compactMap[procI].insert(index, nCompact);
 408         }
 409     }
 410 }
 411
 412
 413 void Foam::mapDistribute::calcCompactAddressing
 414 (
 415     const globalIndex& globalNumbering,
 416     const labelListList& cellCells,
 417     List<Map<label> >& compactMap
 418 ) const
 419 {
 420     compactMap.setSize(Pstream::nProcs());
 421
 422     // Count all (non-local) elements needed. Just for presizing map.
 423     labelList nNonLocal(Pstream::nProcs(), 0);
 424
 425     forAll(cellCells, cellI)
 426     {
 427         const labelList& cCells = cellCells[cellI];
 428
 429         forAll(cCells, i)
 430         {
 431             label globalIndex = cCells[i];
 432
 433             if (globalIndex != -1 && !globalNumbering.isLocal(globalIndex))
 434             {
 435                 label procI = globalNumbering.whichProcID(globalIndex);
 436                 nNonLocal[procI]++;
 437             }
 438         }
 439     }
 440
 441     forAll(compactMap, procI)
 442     {
 443         compactMap[procI].clear();
 444         if (procI != Pstream::myProcNo())
 445         {
 446             compactMap[procI].resize(2*nNonLocal[procI]);
 447         }
 448     }
 449
 450
 451     // Collect all (non-local) elements needed.
 452     forAll(cellCells, cellI)
 453     {
 454         const labelList& cCells = cellCells[cellI];
 455
 456         forAll(cCells, i)
 457         {
 458             label globalIndex = cCells[i];
 459
 460             if (globalIndex != -1 && !globalNumbering.isLocal(globalIndex))
 461             {
 462                 label procI = globalNumbering.whichProcID(globalIndex);
 463                 label index = globalNumbering.toLocal(procI, globalIndex);
 464                 label nCompact = compactMap[procI].size();
 465                 compactMap[procI].insert(index, nCompact);
 466             }
 467         }
 468     }
 469 }
 470
 471
 472 void Foam::mapDistribute::exchangeAddressing
 473 (
 474     const globalIndex& globalNumbering,
 475     labelList& elements,
 476     List<Map<label> >& compactMap,
 477     labelList& compactStart
 478 )
 479 {
 480     // The overall compact addressing is
 481     // - myProcNo data first (uncompacted)
 482     // - all other processors consecutively
 483
 484     compactStart.setSize(Pstream::nProcs());
 485     compactStart[Pstream::myProcNo()] = 0;
 486     constructSize_ = globalNumbering.localSize();
 487     forAll(compactStart, procI)
 488     {
 489         if (procI != Pstream::myProcNo())
 490         {
 491             compactStart[procI] = constructSize_;
 492             constructSize_ += compactMap[procI].size();
 493         }
 494     }
 495
 496
 497
 498     // Find out what to receive/send in compact addressing.
 499
 500     // What I want to receive is what others have to send
 501     labelListList wantedRemoteElements(Pstream::nProcs());
 502     // Compact addressing for received data
 503     constructMap_.setSize(Pstream::nProcs());
 504     forAll(compactMap, procI)
 505     {
 506         if (procI == Pstream::myProcNo())
 507         {
 508             // All my own elements are used
 509             label nLocal = globalNumbering.localSize();
 510             wantedRemoteElements[procI] = identity(nLocal);
 511             constructMap_[procI] = identity(nLocal);
 512         }
 513         else
 514         {
 515             // Remote elements wanted from processor procI
 516             labelList& remoteElem = wantedRemoteElements[procI];
 517             labelList& localElem = constructMap_[procI];
 518             remoteElem.setSize(compactMap[procI].size());
 519             localElem.setSize(compactMap[procI].size());
 520             label i = 0;
 521             forAllIter(Map<label>, compactMap[procI], iter)
 522             {
 523                 const label compactI = compactStart[procI] + iter();
 524                 remoteElem[i] = iter.key();
 525                 localElem[i]  = compactI;
 526                 iter() = compactI;
 527                 i++;
 528             }
 529         }
 530     }
 531
 532     subMap_.setSize(Pstream::nProcs());
 533     labelListList sendSizes;
 534     Pstream::exchange<labelList, label>
 535     (
 536         wantedRemoteElements,
 537         subMap_,
 538         sendSizes
 539     );
 540
 541     // Renumber elements
 542     forAll(elements, i)
 543     {
 544         elements[i] = renumber(globalNumbering, compactMap, elements[i]);
 545     }
 546 }
 547
 548
 549 void Foam::mapDistribute::exchangeAddressing
 550 (
 551     const globalIndex& globalNumbering,
 552     labelListList& cellCells,
 553     List<Map<label> >& compactMap,
 554     labelList& compactStart
 555 )
 556 {
 557     // The overall compact addressing is
 558     // - myProcNo data first (uncompacted)
 559     // - all other processors consecutively
 560
 561     compactStart.setSize(Pstream::nProcs());
 562     compactStart[Pstream::myProcNo()] = 0;
 563     constructSize_ = globalNumbering.localSize();
 564     forAll(compactStart, procI)
 565     {
 566         if (procI != Pstream::myProcNo())
 567         {
 568             compactStart[procI] = constructSize_;
 569             constructSize_ += compactMap[procI].size();
 570         }
 571     }
 572
 573
 574
 575     // Find out what to receive/send in compact addressing.
 576
 577     // What I want to receive is what others have to send
 578     labelListList wantedRemoteElements(Pstream::nProcs());
 579     // Compact addressing for received data
 580     constructMap_.setSize(Pstream::nProcs());
 581     forAll(compactMap, procI)
 582     {
 583         if (procI == Pstream::myProcNo())
 584         {
 585             // All my own elements are used
 586             label nLocal = globalNumbering.localSize();
 587             wantedRemoteElements[procI] = identity(nLocal);
 588             constructMap_[procI] = identity(nLocal);
 589         }
 590         else
 591         {
 592             // Remote elements wanted from processor procI
 593             labelList& remoteElem = wantedRemoteElements[procI];
 594             labelList& localElem = constructMap_[procI];
 595             remoteElem.setSize(compactMap[procI].size());
 596             localElem.setSize(compactMap[procI].size());
 597             label i = 0;
 598             forAllIter(Map<label>, compactMap[procI], iter)
 599             {
 600                 const label compactI = compactStart[procI] + iter();
 601                 remoteElem[i] = iter.key();
 602                 localElem[i]  = compactI;
 603                 iter() = compactI;
 604                 i++;
 605             }
 606         }
 607     }
 608
 609     subMap_.setSize(Pstream::nProcs());
 610     labelListList sendSizes;
 611     Pstream::exchange<labelList, label>
 612     (
 613         wantedRemoteElements,
 614         subMap_,
 615         sendSizes
 616     );
 617
 618     // Renumber elements
 619     forAll(cellCells, cellI)
 620     {
 621         labelList& cCells = cellCells[cellI];
 622
 623         forAll(cCells, i)
 624         {
 625             cCells[i] = renumber(globalNumbering, compactMap, cCells[i]);
 626         }
 627     }
 628 }
 629
 630
 631 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 632
 633 //- Construct null
 634 Foam::mapDistribute::mapDistribute()
 635 :
 636     constructSize_(0),
 637     schedulePtr_()
 638 {}
 639
 640
 641 //- Construct from components
 642 Foam::mapDistribute::mapDistribute
 643 (
 644     const label constructSize,
 645     const Xfer<labelListList>& subMap,
 646     const Xfer<labelListList>& constructMap
 647 )
 648 :
 649     constructSize_(constructSize),
 650     subMap_(subMap),
 651     constructMap_(constructMap),
 652     schedulePtr_()
 653 {}
 654
 655
 656 //- Construct from components
 657 Foam::mapDistribute::mapDistribute
 658 (
 659     const label constructSize,
 660     const Xfer<labelListList>& subMap,
 661     const Xfer<labelListList>& constructMap,
 662     const Xfer<labelListList>& transformElements,
 663     const Xfer<labelList>& transformStart
 664 )
 665 :
 666     constructSize_(constructSize),
 667     subMap_(subMap),
 668     constructMap_(constructMap),
 669     transformElements_(transformElements),
 670     transformStart_(transformStart),
 671     schedulePtr_()
 672 {}
 673
 674
 675 Foam::mapDistribute::mapDistribute
 676 (
 677     const labelList& sendProcs,
 678     const labelList& recvProcs
 679 )
 680 :
 681     constructSize_(0),
 682     schedulePtr_()
 683 {
 684     if (sendProcs.size() != recvProcs.size())
 685     {
 686         FatalErrorIn
 687         (
 688             "mapDistribute::mapDistribute(const labelList&, const labelList&)"
 689         )   << "The send and receive data is not the same length. sendProcs:"
 690             << sendProcs.size() << " recvProcs:" << recvProcs.size()
 691             << abort(FatalError);
 692     }
 693
 694     // Per processor the number of samples we have to send/receive.
 695     labelList nSend(Pstream::nProcs(), 0);
 696     labelList nRecv(Pstream::nProcs(), 0);
 697
 698     forAll(sendProcs, sampleI)
 699     {
 700         label sendProc = sendProcs[sampleI];
 701         label recvProc = recvProcs[sampleI];
 702
 703         // Note that also need to include local communication (both
 704         // RecvProc and sendProc on local processor)
 705
 706         if (Pstream::myProcNo() == sendProc)
 707         {
 708             // I am the sender. Count destination processor.
 709             nSend[recvProc]++;
 710         }
 711         if (Pstream::myProcNo() == recvProc)
 712         {
 713             // I am the receiver.
 714             nRecv[sendProc]++;
 715         }
 716     }
 717
 718     subMap_.setSize(Pstream::nProcs());
 719     constructMap_.setSize(Pstream::nProcs());
 720     forAll(nSend, procI)
 721     {
 722         subMap_[procI].setSize(nSend[procI]);
 723         constructMap_[procI].setSize(nRecv[procI]);
 724     }
 725     nSend = 0;
 726     nRecv = 0;
 727
 728     forAll(sendProcs, sampleI)
 729     {
 730         label sendProc = sendProcs[sampleI];
 731         label recvProc = recvProcs[sampleI];
 732
 733         if (Pstream::myProcNo() == sendProc)
 734         {
 735             // I am the sender. Store index I need to send.
 736             subMap_[recvProc][nSend[recvProc]++] = sampleI;
 737         }
 738         if (Pstream::myProcNo() == recvProc)
 739         {
 740             // I am the receiver.
 741             constructMap_[sendProc][nRecv[sendProc]++] = sampleI;
 742             // Largest entry inside constructMap
 743             constructSize_ = sampleI+1;
 744         }
 745     }
 746 }
 747
 748
 749 Foam::mapDistribute::mapDistribute
 750 (
 751     const globalIndex& globalNumbering,
 752     labelList& elements,
 753     List<Map<label> >& compactMap
 754 )
 755 :
 756     constructSize_(0),
 757     schedulePtr_()
 758 {
 759     // Construct per processor compact addressing of the global elements
 760     // needed. The ones from the local processor are not included since
 761     // these are always all needed.
 762     calcCompactAddressing
 763     (
 764         globalNumbering,
 765         elements,
 766         compactMap
 767     );
 768
 769     //// Sort remote elements needed (not really necessary)
 770     //forAll(compactMap, procI)
 771     //{
 772     //    if (procI != Pstream::myProcNo())
 773     //    {
 774     //        Map<label>& globalMap = compactMap[procI];
 775     //
 776     //        SortableList<label> sorted(globalMap.toc().xfer());
 777     //
 778     //        forAll(sorted, i)
 779     //        {
 780     //            Map<label>::iterator iter = globalMap.find(sorted[i]);
 781     //            iter() = i;
 782     //        }
 783     //    }
 784     //}
 785
 786
 787     // Exchange what I need with processor that supplies it. Renumber elements
 788     // into compact numbering
 789     labelList compactStart;
 790     exchangeAddressing
 791     (
 792         globalNumbering,
 793         elements,
 794         compactMap,
 795         compactStart
 796     );
 797
 798     if (debug)
 799     {
 800         printLayout(Pout);
 801     }
 802 }
 803
 804
 805 Foam::mapDistribute::mapDistribute
 806 (
 807     const globalIndex& globalNumbering,
 808     labelListList& cellCells,
 809     List<Map<label> >& compactMap
 810 )
 811 :
 812     constructSize_(0),
 813     schedulePtr_()
 814 {
 815     // Construct per processor compact addressing of the global elements
 816     // needed. The ones from the local processor are not included since
 817     // these are always all needed.
 818     calcCompactAddressing
 819     (
 820         globalNumbering,
 821         cellCells,
 822         compactMap
 823     );
 824
 825     //// Sort remote elements needed (not really necessary)
 826     //forAll(compactMap, procI)
 827     //{
 828     //    if (procI != Pstream::myProcNo())
 829     //    {
 830     //        Map<label>& globalMap = compactMap[procI];
 831     //
 832     //        SortableList<label> sorted(globalMap.toc().xfer());
 833     //
 834     //        forAll(sorted, i)
 835     //        {
 836     //            Map<label>::iterator iter = globalMap.find(sorted[i]);
 837     //            iter() = i;
 838     //        }
 839     //    }
 840     //}
 841
 842
 843     // Exchange what I need with processor that supplies it. Renumber elements
 844     // into compact numbering
 845     labelList compactStart;
 846     exchangeAddressing
 847     (
 848         globalNumbering,
 849         cellCells,
 850         compactMap,
 851         compactStart
 852     );
 853
 854     if (debug)
 855     {
 856         printLayout(Pout);
 857     }
 858 }
 859
 860
 861 Foam::mapDistribute::mapDistribute
 862 (
 863     const globalIndex& globalNumbering,
 864     labelList& elements,
 865     const globalIndexAndTransform& globalTransforms,
 866     const labelPairList& transformedElements,
 867     labelList& transformedIndices,
 868     List<Map<label> >& compactMap
 869 )
 870 :
 871     constructSize_(0),
 872     schedulePtr_()
 873 {
 874     // Construct per processor compact addressing of the global elements
 875     // needed. The ones from the local processor are not included since
 876     // these are always all needed.
 877     calcCompactAddressing
 878     (
 879         globalNumbering,
 880         elements,
 881         compactMap
 882     );
 883
 884     // Add all (non-local) transformed elements needed.
 885     forAll(transformedElements, i)
 886     {
 887         labelPair elem = transformedElements[i];
 888         label procI = globalIndexAndTransform::processor(elem);
 889         if (procI != Pstream::myProcNo())
 890         {
 891             label index = globalIndexAndTransform::index(elem);
 892             label nCompact = compactMap[procI].size();
 893             compactMap[procI].insert(index, nCompact);
 894         }
 895     }
 896
 897
 898     // Exchange what I need with processor that supplies it. Renumber elements
 899     // into compact numbering
 900     labelList compactStart;
 901     exchangeAddressing
 902     (
 903         globalNumbering,
 904         elements,
 905         compactMap,
 906         compactStart
 907     );
 908
 909
 910     // Renumber the transformed elements
 911     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 912     // Count per transformIndex
 913     label nTrafo = globalTransforms.transformPermutations().size();
 914     labelList nPerTransform(nTrafo, 0);
 915     forAll(transformedElements, i)
 916     {
 917         labelPair elem = transformedElements[i];
 918         label trafoI = globalIndexAndTransform::transformIndex(elem);
 919         nPerTransform[trafoI]++;
 920     }
 921     // Offset per transformIndex
 922     transformStart_.setSize(nTrafo);
 923     transformElements_.setSize(nTrafo);
 924     forAll(transformStart_, trafoI)
 925     {
 926         transformStart_[trafoI] = constructSize_;
 927         constructSize_ += nPerTransform[trafoI];
 928         transformElements_[trafoI].setSize(nPerTransform[trafoI]);
 929     }
 930
 931     // Sort transformed elements into their new slot.
 932     nPerTransform = 0;
 933
 934     transformedIndices.setSize(transformedElements.size());
 935     forAll(transformedElements, i)
 936     {
 937         labelPair elem = transformedElements[i];
 938         label procI = globalIndexAndTransform::processor(elem);
 939         label index = globalIndexAndTransform::index(elem);
 940         label trafoI = globalIndexAndTransform::transformIndex(elem);
 941
 942         // Get compact index for untransformed element
 943         label rawElemI =
 944         (
 945             procI == Pstream::myProcNo()
 946           ? index
 947           : compactMap[procI][index]
 948         );
 949
 950         label& n = nPerTransform[trafoI];
 951         // index of element to transform
 952         transformElements_[trafoI][n] = rawElemI;
 953         // destination of transformed element
 954         transformedIndices[i] = transformStart_[trafoI]+n;
 955         n++;
 956     }
 957
 958     if (debug)
 959     {
 960         printLayout(Pout);
 961     }
 962 }
 963
 964
 965 Foam::mapDistribute::mapDistribute
 966 (
 967     const globalIndex& globalNumbering,
 968     labelListList& cellCells,
 969     const globalIndexAndTransform& globalTransforms,
 970     const List<labelPairList>& transformedElements,
 971     labelListList& transformedIndices,
 972     List<Map<label> >& compactMap
 973 )
 974 :
 975     constructSize_(0),
 976     schedulePtr_()
 977 {
 978     // Construct per processor compact addressing of the global elements
 979     // needed. The ones from the local processor are not included since
 980     // these are always all needed.
 981     calcCompactAddressing
 982     (
 983         globalNumbering,
 984         cellCells,
 985         compactMap
 986     );
 987
 988     // Add all (non-local) transformed elements needed.
 989     forAll(transformedElements, cellI)
 990     {
 991         const labelPairList& elems = transformedElements[cellI];
 992
 993         forAll(elems, i)
 994         {
 995             label procI = globalIndexAndTransform::processor(elems[i]);
 996             if (procI != Pstream::myProcNo())
 997             {
 998                 label index = globalIndexAndTransform::index(elems[i]);
 999                 label nCompact = compactMap[procI].size();
1000                 compactMap[procI].insert(index, nCompact);
1001             }
1002         }
1003     }
1004
1005
1006     // Exchange what I need with processor that supplies it. Renumber elements
1007     // into compact numbering
1008     labelList compactStart;
1009     exchangeAddressing
1010     (
1011         globalNumbering,
1012         cellCells,
1013         compactMap,
1014         compactStart
1015     );
1016
1017
1018     // Renumber the transformed elements
1019     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1020     // Count per transformIndex
1021     label nTrafo = globalTransforms.transformPermutations().size();
1022     labelList nPerTransform(nTrafo, 0);
1023     forAll(transformedElements, cellI)
1024     {
1025         const labelPairList& elems = transformedElements[cellI];
1026
1027         forAll(elems, i)
1028         {
1029             label trafoI = globalIndexAndTransform::transformIndex(elems[i]);
1030             nPerTransform[trafoI]++;
1031         }
1032     }
1033     // Offset per transformIndex
1034     transformStart_.setSize(nTrafo);
1035     transformElements_.setSize(nTrafo);
1036     forAll(transformStart_, trafoI)
1037     {
1038         transformStart_[trafoI] = constructSize_;
1039         constructSize_ += nPerTransform[trafoI];
1040         transformElements_[trafoI].setSize(nPerTransform[trafoI]);
1041     }
1042
1043     // Sort transformed elements into their new slot.
1044     nPerTransform = 0;
1045
1046     transformedIndices.setSize(transformedElements.size());
1047     forAll(transformedElements, cellI)
1048     {
1049         const labelPairList& elems = transformedElements[cellI];
1050         transformedIndices[cellI].setSize(elems.size());
1051
1052         forAll(elems, i)
1053         {
1054             label procI = globalIndexAndTransform::processor(elems[i]);
1055             label index = globalIndexAndTransform::index(elems[i]);
1056             label trafoI = globalIndexAndTransform::transformIndex(elems[i]);
1057
1058             // Get compact index for untransformed element
1059             label rawElemI =
1060             (
1061                 procI == Pstream::myProcNo()
1062               ? index
1063               : compactMap[procI][index]
1064             );
1065
1066             label& n = nPerTransform[trafoI];
1067             // index of element to transform
1068             transformElements_[trafoI][n] = rawElemI;
1069             // destination of transformed element
1070             transformedIndices[cellI][i] = transformStart_[trafoI]+n;
1071             n++;
1072         }
1073     }
1074
1075     if (debug)
1076     {
1077         printLayout(Pout);
1078     }
1079 }
1080
1081
1082 Foam::mapDistribute::mapDistribute(const mapDistribute& map)
1083 :
1084     constructSize_(map.constructSize_),
1085     subMap_(map.subMap_),
1086     constructMap_(map.constructMap_),
1087     transformElements_(map.transformElements_),
1088     transformStart_(map.transformStart_),
1089     schedulePtr_()
1090 {}
1091
1092
1093 Foam::mapDistribute::mapDistribute(const Xfer<mapDistribute>& map)
1094 :
1095     constructSize_(map().constructSize_),
1096     subMap_(map().subMap_.xfer()),
1097     constructMap_(map().constructMap_.xfer()),
1098     transformElements_(map().transformElements_.xfer()),
1099     transformStart_(map().transformStart_.xfer()),
1100     schedulePtr_()
1101 {}
1102
1103
1104 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
1105
1106 Foam::label Foam::mapDistribute::whichTransform(const label index) const
1107 {
1108     return findLower(transformStart_, index+1);
1109 }
1110
1111
1112 void Foam::mapDistribute::transfer(mapDistribute& rhs)
1113 {
1114     constructSize_ = rhs.constructSize_;
1115     subMap_.transfer(rhs.subMap_);
1116     constructMap_.transfer(rhs.constructMap_);
1117     transformElements_.transfer(rhs.transformElements_);
1118     transformStart_.transfer(rhs.transformStart_);
1119     schedulePtr_.clear();
1120 }
1121
1122
1123 Foam::Xfer<Foam::mapDistribute> Foam::mapDistribute::xfer()
1124 {
1125     return xferMove(*this);
1126 }
1127
1128
1129 Foam::label Foam::mapDistribute::renumber
1130 (
1131     const globalIndex& globalNumbering,
1132     const List<Map<label> >& compactMap,
1133     const label globalI
1134 )
1135 {
1136     if (globalI == -1)
1137     {
1138         return globalI;
1139     }
1140     if (globalNumbering.isLocal(globalI))
1141     {
1142         return globalNumbering.toLocal(globalI);
1143     }
1144     else
1145     {
1146         label procI = globalNumbering.whichProcID(globalI);
1147         label index = globalNumbering.toLocal(procI, globalI);
1148         return compactMap[procI][index];
1149     }
1150 }
1151
1152
1153 void Foam::mapDistribute::compact(const boolList& elemIsUsed)
1154 {
1155     // 1. send back to sender. Have sender delete the corresponding element
1156     //    from the submap and do the same to the constructMap locally
1157     //    (and in same order).
1158
1159     // Send elemIsUsed field to neighbour. Use nonblocking code from
1160     // mapDistribute but in reverse order.
1161     if (Pstream::parRun())
1162     {
1163         List<boolList> sendFields(Pstream::nProcs());
1164
1165         for (label domain = 0; domain < Pstream::nProcs(); domain++)
1166         {
1167             const labelList& map = constructMap_[domain];
1168
1169             if (domain != Pstream::myProcNo() && map.size())
1170             {
1171                 boolList& subField = sendFields[domain];
1172                 subField.setSize(map.size());
1173                 forAll(map, i)
1174                 {
1175                     subField[i] = elemIsUsed[map[i]];
1176                 }
1177
1178                 OPstream::write
1179                 (
1180                     Pstream::nonBlocking,
1181                     domain,
1182                     reinterpret_cast<const char*>(subField.begin()),
1183                     subField.size()*sizeof(bool)
1184                 );
1185             }
1186         }
1187
1188         // Set up receives from neighbours
1189
1190         List<boolList> recvFields(Pstream::nProcs());
1191
1192         for (label domain = 0; domain < Pstream::nProcs(); domain++)
1193         {
1194             const labelList& map = subMap_[domain];
1195
1196             if (domain != Pstream::myProcNo() && map.size())
1197             {
1198                 recvFields[domain].setSize(map.size());
1199                 IPstream::read
1200                 (
1201                     Pstream::nonBlocking,
1202                     domain,
1203                     reinterpret_cast<char*>(recvFields[domain].begin()),
1204                     recvFields[domain].size()*sizeof(bool)
1205                 );
1206             }
1207         }
1208
1209
1210         // Set up 'send' to myself - write directly into recvFields
1211
1212         {
1213             const labelList& map = constructMap_[Pstream::myProcNo()];
1214
1215             recvFields[Pstream::myProcNo()].setSize(map.size());
1216             forAll(map, i)
1217             {
1218                 recvFields[Pstream::myProcNo()][i] = elemIsUsed[map[i]];
1219             }
1220         }
1221
1222
1223         // Wait for all to finish
1224
1225         Pstream::waitRequests();
1226
1227
1228         // Compact out all submap entries that are referring to unused elements
1229         for (label domain = 0; domain < Pstream::nProcs(); domain++)
1230         {
1231             const labelList& map = subMap_[domain];
1232
1233             labelList newMap(map.size());
1234             label newI = 0;
1235
1236             forAll(map, i)
1237             {
1238                 if (recvFields[domain][i])
1239                 {
1240                     // So element is used on destination side
1241                     newMap[newI++] = map[i];
1242                 }
1243             }
1244             if (newI < map.size())
1245             {
1246                 newMap.setSize(newI);
1247                 subMap_[domain].transfer(newMap);
1248             }
1249         }
1250     }
1251
1252
1253     // 2. remove from construct map - since end-result (element in elemIsUsed)
1254     //    not used.
1255
1256     label maxConstructIndex = -1;
1257
1258     for (label domain = 0; domain < Pstream::nProcs(); domain++)
1259     {
1260         const labelList& map = constructMap_[domain];
1261
1262         labelList newMap(map.size());
1263         label newI = 0;
1264
1265         forAll(map, i)
1266         {
1267             label destinationI = map[i];
1268
1269             // Is element is used on destination side
1270             if (elemIsUsed[destinationI])
1271             {
1272                 maxConstructIndex = max(maxConstructIndex, destinationI);
1273
1274                 newMap[newI++] = destinationI;
1275             }
1276         }
1277         if (newI < map.size())
1278         {
1279             newMap.setSize(newI);
1280             constructMap_[domain].transfer(newMap);
1281         }
1282     }
1283
1284     constructSize_ = maxConstructIndex+1;
1285
1286     // Clear the schedule (note:not necessary if nothing changed)
1287     schedulePtr_.clear();
1288 }
1289
1290
1291 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
1292
1293 void Foam::mapDistribute::operator=(const mapDistribute& rhs)
1294 {
1295     // Check for assignment to self
1296     if (this == &rhs)
1297     {
1298         FatalErrorIn
1299         (
1300             "Foam::mapDistribute::operator=(const Foam::mapDistribute&)"
1301         )   << "Attempted assignment to self"
1302             << abort(FatalError);
1303     }
1304     constructSize_ = rhs.constructSize_;
1305     subMap_ = rhs.subMap_;
1306     constructMap_ = rhs.constructMap_;
1307     transformElements_ = rhs.transformElements_;
1308     transformStart_ = rhs.transformStart_;
1309     schedulePtr_.clear();
1310 }
1311
1312
1313 // ************************************************************************* //