src/mesh/cfMesh/utilities/octrees/meshOctree/meshOctreeInsideOutside/meshOctreeInsideOutside.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | cfMesh: A library for mesh generation
   4    \\    /   O peration     |
   5     \\  /    A nd           | Author: Franjo Juretic (franjo.juretic@c-fields.com)
   6      \\/     M anipulation  | Copyright (C) Creative Fields, Ltd.
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of cfMesh.
  10
  11     cfMesh is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     cfMesh 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 cfMesh.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25
  26 \*---------------------------------------------------------------------------*/
  27
  28 #include "meshOctreeInsideOutside.H"
  29 #include "triSurf.H"
  30 #include "boundBox.H"
  31 #include "labelLongList.H"
  32
  33 # ifdef USE_OMP
  34 #include <omp.h>
  35 # endif
  36
  37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  38
  39 namespace Foam
  40 {
  41
  42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  43
  44 meshOctreeInsideOutside::meshOctreeInsideOutside
  45 (
  46     meshOctree& octree
  47 )
  48 :
  49     octreeModifier_(octree),
  50     cubeGroup_(octree.numberOfLeaves(), -1),
  51     cubesInGroup_(),
  52     groupType_(),
  53     boundaryDATACubes_(),
  54     hasOutsideNeighbour_(octree.numberOfLeaves(), false),
  55     communicationCubes_(),
  56     neighbouringGroups_()
  57 {
  58     initialiseBoxes();
  59
  60     frontalMarking();
  61
  62     markOutsideCubes();
  63
  64     reviseDataBoxes();
  65
  66     markInsideCubes();
  67
  68     label nInternal(0), nUnknown(0), nData(0), nOutside(0);
  69
  70     const label nLeaves = octree.numberOfLeaves();
  71     for(label leafI=0;leafI<nLeaves;++leafI)
  72     {
  73         const meshOctreeCubeBasic& oc = octree.returnLeaf(leafI);
  74
  75         if( oc.cubeType() & meshOctreeCube::INSIDE )
  76         {
  77             ++nInternal;
  78         }
  79         else if( oc.cubeType() & meshOctreeCube::UNKNOWN )
  80         {
  81             ++nUnknown;
  82         }
  83         else if( oc.cubeType() & meshOctreeCube::DATA )
  84         {
  85             ++nData;
  86         }
  87         else if( oc.cubeType() & meshOctreeCube::OUTSIDE )
  88         {
  89             ++nOutside;
  90         }
  91     }
  92
  93     if( octree.neiProcs().size() )
  94     {
  95         reduce(nInternal, sumOp<label>());
  96         reduce(nUnknown, sumOp<label>());
  97         reduce(nData, sumOp<label>());
  98         reduce(nOutside, sumOp<label>());
  99     }
 100
 101     Info << "Number of internal boxes is " << nInternal << endl;
 102     Info << "Number of outside boxes is " << nOutside << endl;
 103     Info << "Number of data boxes is " << nData << endl;
 104     Info << "Number of unknown boxes is " << nUnknown << endl;
 105 }
 106
 107 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 108
 109 meshOctreeInsideOutside::~meshOctreeInsideOutside()
 110 {
 111 }
 112
 113 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 114
 115 void meshOctreeInsideOutside::initialiseBoxes()
 116 {
 117     const LongList<meshOctreeCube*>& leaves = octreeModifier_.leavesAccess();
 118
 119     # ifdef USE_OMP
 120     # pragma omp parallel for if( leaves.size() > 1000 )
 121     # endif
 122     forAll(leaves, leafI)
 123     {
 124         if( leaves[leafI]->hasContainedElements() )
 125         {
 126             leaves[leafI]->setCubeType(meshOctreeCubeBasic::DATA);
 127         }
 128         else
 129         {
 130             leaves[leafI]->setCubeType(meshOctreeCubeBasic::UNKNOWN);
 131         }
 132     }
 133 }
 134
 135 void meshOctreeInsideOutside::frontalMarking()
 136 {
 137     communicationCubes_.clear();
 138     neighbouringGroups_.clear();
 139
 140     labelLongList frontCubes;
 141     DynList<label> neighbours;
 142
 143     label nGroup(0), nThreads(1);
 144
 145     const LongList<meshOctreeCube*>& leaves = octreeModifier_.leavesAccess();
 146     const meshOctree& octree = octreeModifier_.octree();
 147
 148     boolList commCubes(leaves.size(), false);
 149
 150     # ifdef USE_OMP
 151     if( leaves.size() > 1000 )
 152         nThreads = 3 * omp_get_num_procs();
 153
 154     # pragma omp parallel num_threads(nThreads) \
 155     private(frontCubes, neighbours)
 156     # endif
 157     {
 158         LongList<std::pair<label, label> > threadCommPairs;
 159
 160         # ifdef USE_OMP
 161         const label threadI = omp_get_thread_num();
 162         # else
 163         const label threadI(0);
 164         # endif
 165
 166         const label chunkSize = leaves.size() / nThreads + 1;
 167
 168         const label minLeaf = threadI * chunkSize;
 169
 170         const label maxLeaf = Foam::min(leaves.size(), minLeaf + chunkSize);
 171
 172         for(label leafI=minLeaf;leafI<maxLeaf;++leafI)
 173         {
 174             if( leaves[leafI]->hasContainedElements() )
 175                 continue;
 176             if( cubeGroup_[leafI] != -1 )
 177                 continue;
 178
 179             label groupI;
 180             # ifdef USE_OMP
 181             # pragma omp critical
 182             # endif
 183             groupI = nGroup++;
 184
 185             direction cType(meshOctreeCubeBasic::UNKNOWN);
 186             frontCubes.clear();
 187             frontCubes.append(leafI);
 188             cubeGroup_[leafI] = groupI;
 189
 190             labelLongList neiDATACubes;
 191
 192             while( frontCubes.size() )
 193             {
 194                 const label fLabel = frontCubes.removeLastElement();
 195                 octree.findNeighboursForLeaf(fLabel, neighbours);
 196
 197                 forAll(neighbours, neiI)
 198                 {
 199                     const label nei = neighbours[neiI];
 200                     if( (nei >= minLeaf) && (nei < maxLeaf) )
 201                     {
 202                         if( cubeGroup_[nei] != -1 )
 203                             continue;
 204
 205                         if( leaves[nei]->hasContainedElements() )
 206                         {
 207                             neiDATACubes.append(nei);
 208                         }
 209                         else
 210                         {
 211                             frontCubes.append(nei);
 212                             cubeGroup_[nei] = groupI;
 213                         }
 214                     }
 215                     else if( nei == -1 )
 216                     {
 217                         cType = meshOctreeCubeBasic::OUTSIDE;
 218                     }
 219                     else if( nei == meshOctreeCubeBasic::OTHERPROC )
 220                     {
 221                         commCubes[fLabel] = true;
 222                     }
 223                     else
 224                     {
 225                         if( leaves[nei]->hasContainedElements() )
 226                         {
 227                             neiDATACubes.append(nei);
 228                         }
 229                         else
 230                         {
 231                             threadCommPairs.append
 232                             (
 233                                 std::make_pair(fLabel, nei)
 234                             );
 235                         }
 236                     }
 237                 }
 238             }
 239
 240             # ifdef USE_OMP
 241             # pragma omp critical
 242             # endif
 243             {
 244                 if( groupI >= boundaryDATACubes_.size() )
 245                     boundaryDATACubes_.setSize(groupI+1);
 246
 247                 boundaryDATACubes_.setRow(groupI, neiDATACubes);
 248                 groupType_[groupI] = cType;
 249             }
 250         }
 251
 252         # ifdef USE_OMP
 253         # pragma omp barrier
 254         # endif
 255
 256         //- find group to neighbouring groups addressing
 257         List<DynList<label> > localNeiGroups(nGroup);
 258         forAll(threadCommPairs, cfI)
 259         {
 260             const std::pair<label, label>& lp = threadCommPairs[cfI];
 261             const label groupI = cubeGroup_[lp.first];
 262             const label neiGroup = cubeGroup_[lp.second];
 263
 264             if( (neiGroup >= nGroup) || (groupI >= nGroup) )
 265                 FatalError << "neiGroup " << neiGroup
 266                     << " groupI " << groupI << " are >= than "
 267                     << "nGroups " << nGroup << abort(FatalError);
 268
 269             if( neiGroup != -1 )
 270             {
 271                 localNeiGroups[groupI].appendIfNotIn(neiGroup);
 272                 localNeiGroups[neiGroup].appendIfNotIn(groupI);
 273             }
 274         }
 275
 276         # ifdef USE_OMP
 277         # pragma omp critical
 278         # endif
 279         {
 280             neighbouringGroups_.setSize(nGroup);
 281
 282             forAll(localNeiGroups, groupI)
 283             {
 284                 const DynList<label>& lGroups = localNeiGroups[groupI];
 285
 286                 neighbouringGroups_.appendIfNotIn(groupI, groupI);
 287
 288                 forAll(lGroups, i)
 289                 {
 290                     neighbouringGroups_.append(groupI, lGroups[i]);
 291                 }
 292             }
 293         }
 294     }
 295
 296     //- create cubesInGroup_ addressing
 297     labelList nCubesInGroup(nGroup, 0);
 298     forAll(cubeGroup_, leafI)
 299     {
 300         if( cubeGroup_[leafI] < 0 )
 301             continue;
 302
 303         ++nCubesInGroup[cubeGroup_[leafI]];
 304     }
 305
 306     cubesInGroup_.setSizeAndRowSize(nCubesInGroup);
 307
 308     forAllReverse(cubeGroup_, leafI)
 309     {
 310         const label groupI = cubeGroup_[leafI];
 311
 312         if( groupI < 0 )
 313             continue;
 314
 315         cubesInGroup_(groupI, --nCubesInGroup[groupI]) = leafI;
 316     }
 317
 318     //- mark cubes at inter-processor boundaries
 319     forAll(commCubes, leafI)
 320     {
 321         if( commCubes[leafI] )
 322             communicationCubes_.append(leafI);
 323     }
 324
 325     # ifdef DEBUGSearch
 326     label nMarked(0);
 327     forAll(cubeGroup_, leafI)
 328     {
 329         if( cubeGroup_[leafI] != -1 )
 330             ++nMarked;
 331     }
 332     reduce(nMarked, sumOp<label>());
 333     const label totalLeaves = returnReduce(leaves_.size(), sumOp<label>());
 334     Info << "Total number of leaves " << totalLeaves << endl;
 335     Info << "Number of marked leaves " << nMarked << endl;
 336     # endif
 337 }
 338
 339 void meshOctreeInsideOutside::markOutsideCubes()
 340 {
 341     const LongList<meshOctreeCube*>& leaves = octreeModifier_.leavesAccess();
 342     const meshOctree& octree = octreeModifier_.octree();
 343
 344     DynList<label> neighbours;
 345     label nChanged;
 346     bool keepUpdating;
 347
 348     do
 349     {
 350         keepUpdating = false;
 351
 352         do
 353         {
 354             nChanged = 0;
 355
 356             //- make sure that groups created by different threads
 357             //- have the same information
 358             forAll(neighbouringGroups_, groupI)
 359             {
 360                 if( groupType_[groupI] & meshOctreeCubeBasic::OUTSIDE )
 361                 {
 362                     forAllRow(neighbouringGroups_, groupI, i)
 363                     {
 364                         const label neiGroup = neighbouringGroups_(groupI, i);
 365                         if( groupType_[neiGroup] & meshOctreeCube::UNKNOWN )
 366                         {
 367                             ++nChanged;
 368                             groupType_[neiGroup] = meshOctreeCube::OUTSIDE;
 369                         }
 370                     }
 371                 }
 372             }
 373
 374             if( nChanged != 0 )
 375                 keepUpdating = true;
 376
 377         } while( nChanged != 0 );
 378
 379         do
 380         {
 381             nChanged = 0;
 382             LongList<meshOctreeCubeCoordinates> dataToSend;
 383
 384             //- go through the list of communicationCubes and send the ones
 385             //- which are marked as outside
 386             forAll(communicationCubes_, i)
 387             {
 388                 const label groupI = cubeGroup_[communicationCubes_[i]];
 389
 390                 if( groupI < 0 )
 391                     continue;
 392
 393                 if( groupType_[groupI] & meshOctreeCube::OUTSIDE )
 394                     dataToSend.append(*leaves[communicationCubes_[i]]);
 395             }
 396
 397             LongList<meshOctreeCubeCoordinates> receivedCoords;
 398             octree.exchangeRequestsWithNeighbourProcessors
 399             (
 400                 dataToSend,
 401                 receivedCoords
 402             );
 403
 404             //- go through the list of received coordinates and check if any
 405             //- local boxes are their neighbours. If a local neighbour is
 406             //- a DATA box set the hasOutsideNeighbour_ flag to true. If the
 407             //- local neighbour is of UNKNOWN type set it to OUTSIDE.
 408             # ifdef USE_OMP
 409             # pragma omp parallel for if( receivedCoords.size() > 100 ) \
 410             private(neighbours) schedule(dynamic, 20)
 411             # endif
 412             forAll(receivedCoords, i)
 413             {
 414                 octree.findNeighboursForLeaf(receivedCoords[i], neighbours);
 415
 416                 forAll(neighbours, neiI)
 417                 {
 418                     const label nei = neighbours[neiI];
 419
 420                     if( nei < 0 )
 421                         continue;
 422
 423                     if( leaves[nei]->hasContainedElements() )
 424                     {
 425                         hasOutsideNeighbour_[nei] = true;
 426                         continue;
 427                     }
 428                     if( groupType_[cubeGroup_[nei]] & meshOctreeCube::UNKNOWN )
 429                     {
 430                         groupType_[cubeGroup_[nei]] = meshOctreeCube::OUTSIDE;
 431                         ++nChanged;
 432                     }
 433                 }
 434             }
 435
 436             if( nChanged != 0 )
 437                 keepUpdating = true;
 438
 439             reduce(nChanged, sumOp<label>());
 440         } while( nChanged != 0 );
 441
 442         reduce(keepUpdating, maxOp<bool>());
 443
 444     } while( keepUpdating );
 445
 446     //- set OUTSIDE type to the cubes in OUTSIDE groups
 447     for
 448     (
 449         std::map<label, direction>::const_iterator it=groupType_.begin();
 450         it!=groupType_.end();
 451         ++it
 452     )
 453     {
 454         if( it->first < 0 )
 455             continue;
 456
 457         if( it->second & meshOctreeCubeBasic::OUTSIDE )
 458         {
 459             const label groupI = it->first;
 460
 461             //- set the cube type to OUTSIDE
 462             forAllRow(cubesInGroup_, groupI, i)
 463                 leaves[cubesInGroup_(groupI, i)]->setCubeType
 464                 (
 465                     meshOctreeCube::OUTSIDE
 466                 );
 467
 468             //- set true to the collected DATA boxes
 469             forAllRow(boundaryDATACubes_, groupI, neiI)
 470                 hasOutsideNeighbour_[boundaryDATACubes_(groupI, neiI)] = true;
 471         }
 472     }
 473 }
 474
 475 void meshOctreeInsideOutside::reviseDataBoxes()
 476 {
 477     //- remove DATA flag from boxes which do not have an OUTSIDE neighbour
 478     //- and are not surrounded with DATA boxes containing different surface
 479     //- triangles in different patches
 480     const LongList<meshOctreeCube*>& leaves = octreeModifier_.leavesAccess();
 481     const meshOctree& octree = octreeModifier_.octree();
 482     const triSurf& surface = octree.surface();
 483     DynList<label> neighbours;
 484
 485     boolList checkedPatches(leaves.size(), false);
 486
 487     label nMarked;
 488
 489     do
 490     {
 491         nMarked = 0;
 492
 493         LongList<meshOctreeCubeCoordinates> checkCoordinates;
 494         labelHashSet transferCoordinates;
 495
 496         # ifdef USE_OMP
 497         # pragma omp parallel for if( leaves.size() > 1000 ) \
 498         private(neighbours) schedule(dynamic, 20) reduction(+ : nMarked)
 499         # endif
 500         forAll(leaves, leafI)
 501             if( Pstream::parRun() && hasOutsideNeighbour_[leafI] )
 502             {
 503                 octree.findAllLeafNeighbours(leafI, neighbours);
 504                 forAll(neighbours, neiI)
 505                     if( neighbours[neiI] == meshOctreeCubeBasic::OTHERPROC )
 506                     {
 507                         # ifdef USE_OMP
 508                         # pragma omp critical
 509                         # endif
 510                         {
 511                             if( !transferCoordinates.found(leafI) )
 512                             {
 513                                 checkCoordinates.append
 514                                 (
 515                                     leaves[leafI]->coordinates()
 516                                 );
 517                                 transferCoordinates.insert(leafI);
 518                             }
 519                         }
 520
 521                         break;
 522                     }
 523             }
 524             else if(
 525                 (leaves[leafI]->cubeType() & meshOctreeCube::DATA) &&
 526                 !hasOutsideNeighbour_[leafI]
 527             )
 528             {
 529                 meshOctreeCube* oc = leaves[leafI];
 530
 531                 # ifdef DEBUGSearch
 532                 Info << "Box " << leafI << " may not be a DATA box" << endl;
 533                 # endif
 534
 535                 DynList<label> patches;
 536                 const VRWGraph& ct =
 537                     oc->slotPtr()->containedTriangles_;
 538                 const constRow el = ct[oc->containedElements()];
 539                 forAll(el, elI)
 540                     patches.appendIfNotIn(surface[el[elI]].region());
 541
 542                 if( patches.size() > 1 )
 543                     continue;
 544
 545                 checkedPatches[leafI] = true;
 546
 547                 //- check if there exist neighbours
 548                 //- which have some DATA neighbours
 549                 octree.findAllLeafNeighbours(leafI, neighbours);
 550                 forAll(neighbours, neiI)
 551                 {
 552                     const label nei = neighbours[neiI];
 553
 554                     if( nei < 0 )
 555                         continue;
 556
 557                     if( hasOutsideNeighbour_[nei] )
 558                     {
 559                         oc->setCubeType(meshOctreeCube::INSIDE);
 560
 561                         ++nMarked;
 562                         break;
 563                     }
 564                 }
 565             }
 566
 567         if( octree.neiProcs().size() )
 568         {
 569             LongList<meshOctreeCubeCoordinates> receivedCoords;
 570             octree.exchangeRequestsWithNeighbourProcessors
 571             (
 572                 checkCoordinates,
 573                 receivedCoords
 574             );
 575
 576             //- check if any of the local neighbours is a data box with
 577             //- no OUTSIDE neighbours
 578             # ifdef USE_OMP
 579             # pragma omp parallel for if( receivedCoords.size() > 100 ) \
 580             private(neighbours) schedule(dynamic, 20) reduction(+ : nMarked)
 581             # endif
 582             forAll(receivedCoords, i)
 583             {
 584                 octree.findAllLeafNeighbours(receivedCoords[i], neighbours);
 585
 586                 forAll(neighbours, neiI)
 587                 {
 588                     const label nei = neighbours[neiI];
 589
 590                     if( nei < 0 )
 591                         continue;
 592
 593                     if
 594                     (
 595                         (leaves[nei]->cubeType() & meshOctreeCube::DATA) &&
 596                         !hasOutsideNeighbour_[nei] && checkedPatches[nei]
 597                     )
 598                     {
 599                         leaves[nei]->setCubeType(meshOctreeCube::INSIDE);
 600
 601                         ++nMarked;
 602                     }
 603                 }
 604             }
 605
 606             reduce(nMarked, sumOp<label>());
 607         }
 608     } while( nMarked != 0 );
 609
 610     # ifdef DEBUGSearch
 611     label nOutside(0), nData(0), hasOutNei(0);
 612     forAll(leaves, leafI)
 613     {
 614         const direction cType = leaves[leafI]->cubeType();
 615         if( cType & meshOctreeCubeBasic::OUTSIDE )
 616             ++nOutside;
 617         else if( cType & meshOctreeCubeBasic::DATA )
 618             ++nData;
 619
 620         if( hasOutsideNeighbour_[leafI] )
 621             ++hasOutNei;
 622     }
 623
 624     reduce(hasOutNei, sumOp<label>());
 625     reduce(nData, sumOp<label>());
 626     reduce(nOutside, sumOp<label>());
 627     Info << "Number of outside boxes " << nOutside << endl;
 628     Info << "Number of data boxes " << nData << " real data "
 629         << hasOutNei << endl;
 630     returnReduce(1, sumOp<label>());
 631     //::exit(1);
 632     # endif
 633 }
 634
 635 void meshOctreeInsideOutside::markInsideCubes()
 636 {
 637     const LongList<meshOctreeCube*>& leaves = octreeModifier_.leavesAccess();
 638     const meshOctree& octree = octreeModifier_.octree();
 639     label nChanged;
 640     bool keepUpdating;
 641     DynList<label> neighbours;
 642
 643     //- make INSIDE groups for which it is possible
 644     for
 645     (
 646         std::map<label, direction>::iterator it=groupType_.begin();
 647         it!=groupType_.end();
 648         ++it
 649     )
 650     {
 651         const label groupI = it->first;
 652
 653         if( groupI < 0 )
 654             continue;
 655
 656         if( it->second & meshOctreeCubeBasic::UNKNOWN )
 657         {
 658             forAllRow(boundaryDATACubes_, groupI, neiI)
 659             {
 660                 const label cLabel = boundaryDATACubes_(groupI, neiI);
 661                 if(
 662                     hasOutsideNeighbour_[cLabel] ||
 663                     (
 664                         leaves[cLabel]->cubeType() & meshOctreeCube::INSIDE
 665                     )
 666                 )
 667                 {
 668                     it->second = meshOctreeCube::INSIDE;
 669                     break;
 670                 }
 671             }
 672         }
 673     }
 674
 675     do
 676     {
 677         keepUpdating = false;
 678
 679         //- mark INSIDE groups created by different threads
 680         do
 681         {
 682             nChanged = 0;
 683
 684             forAll(neighbouringGroups_, groupI)
 685             {
 686                 if( groupType_[groupI] & meshOctreeCube::INSIDE )
 687                 {
 688                     forAllRow(neighbouringGroups_, groupI, i)
 689                     {
 690                         const label neiGroup = neighbouringGroups_(groupI, i);
 691
 692                         if( groupType_[neiGroup] & meshOctreeCube::UNKNOWN )
 693                         {
 694                             ++nChanged;
 695                             groupType_[neiGroup] = meshOctreeCube::INSIDE;
 696                         }
 697                     }
 698                 }
 699             }
 700
 701             if( nChanged != 0 )
 702                 keepUpdating = true;
 703
 704         } while( nChanged != 0 );
 705
 706         if( octree.neiProcs().size() == 0 )
 707             continue;
 708
 709         //- the code for exchanging data between different processes
 710         LongList<meshOctreeCubeCoordinates> dataToSend, receivedCoords;
 711
 712         //- send coordinates of boxes with hasOutsideNeighbour_ flag and
 713         //- the boxes which have been marked as INSIDE to the neighbouring procs
 714         forAll(hasOutsideNeighbour_, leafI)
 715         {
 716             if(
 717                 hasOutsideNeighbour_[leafI] ||
 718                 (leaves[leafI]->cubeType() & meshOctreeCubeBasic::INSIDE)
 719             )
 720             {
 721                 octree.findNeighboursForLeaf(leafI, neighbours);
 722
 723                 forAll(neighbours, neiI)
 724                     if( neighbours[neiI] == meshOctreeCube::OTHERPROC )
 725                     {
 726                         dataToSend.append(leaves[leafI]->coordinates());
 727                         break;
 728                     }
 729             }
 730         }
 731
 732         octree.exchangeRequestsWithNeighbourProcessors
 733         (
 734             dataToSend,
 735             receivedCoords
 736         );
 737
 738         # ifdef USE_OMP
 739         # pragma omp parallel for if( receivedCoords.size() > 100 ) \
 740         private(neighbours) schedule(dynamic, 20)
 741         # endif
 742         forAll(receivedCoords, i)
 743         {
 744             octree.findNeighboursForLeaf(receivedCoords[i], neighbours);
 745
 746             forAll(neighbours, neiI)
 747             {
 748                 const label nei = neighbours[neiI];
 749
 750                 if( nei < 0 )
 751                     continue;
 752
 753                 const label groupI = cubeGroup_[nei];
 754
 755                 if( groupI < 0 )
 756                     continue;
 757
 758                 if( groupType_[groupI] & meshOctreeCube::UNKNOWN )
 759                     groupType_[groupI] = meshOctreeCubeBasic::INSIDE;
 760             }
 761         }
 762
 763         do
 764         {
 765             nChanged = 0;
 766             dataToSend.clear();
 767
 768             //- go through the list of communicationCubes and send the ones
 769             //- which are marked as outside
 770             forAll(communicationCubes_, i)
 771             {
 772                 if(
 773                     groupType_[cubeGroup_[communicationCubes_[i]]] &
 774                     meshOctreeCubeBasic::INSIDE
 775                 )
 776                     dataToSend.append
 777                     (
 778                         leaves[communicationCubes_[i]]->coordinates()
 779                     );
 780             }
 781
 782             receivedCoords.clear();
 783             octree.exchangeRequestsWithNeighbourProcessors
 784             (
 785                 dataToSend,
 786                 receivedCoords
 787             );
 788
 789             //- go through the list of received coordinates and check if any
 790             //- local boxes are their neighbours. If a local neighbour is
 791             //- a DATA box set the hasOutsideNeighbour_ flag to true. If the
 792             //- local neighbour is of UNKNOWN type set it to OUTSIDE.
 793             # ifdef USE_OMP
 794             # pragma omp parallel for if( receivedCoords.size() > 100 ) \
 795             private(neighbours) schedule(dynamic, 20) reduction(+ : nChanged)
 796             # endif
 797             forAll(receivedCoords, i)
 798             {
 799                 octree.findNeighboursForLeaf(receivedCoords[i], neighbours);
 800
 801                 forAll(neighbours, neiI)
 802                 {
 803                     const label nei = neighbours[neiI];
 804
 805                     if( nei < 0 )
 806                         continue;
 807
 808                     const label groupI = cubeGroup_[nei];
 809
 810                     if( groupI < 0 )
 811                         continue;
 812
 813                     if( groupType_[groupI] & meshOctreeCube::UNKNOWN )
 814                     {
 815                         groupType_[groupI] = meshOctreeCube::INSIDE;
 816                         ++nChanged;
 817                     }
 818                 }
 819             }
 820
 821             reduce(nChanged, sumOp<label>());
 822
 823             if( nChanged != 0 )
 824                 keepUpdating = true;
 825
 826         } while( nChanged != 0 );
 827
 828         reduce(keepUpdating, maxOp<bool>());
 829
 830     } while( keepUpdating );
 831
 832     //- set INSIDE type to the cubes in INSIDE groups
 833     for
 834     (
 835         std::map<label, direction>::const_iterator it=groupType_.begin();
 836         it!=groupType_.end();
 837         ++it
 838     )
 839     {
 840         if( it->first < 0 )
 841             continue;
 842
 843         if( it->second & meshOctreeCubeBasic::INSIDE )
 844         {
 845             const label groupI = it->first;
 846
 847             //- set the cube type to OUTSIDE
 848             forAllRow(cubesInGroup_, groupI, i)
 849                 leaves[cubesInGroup_(groupI, i)]->setCubeType
 850                 (
 851                     meshOctreeCube::INSIDE
 852                 );
 853         }
 854     }
 855 }
 856
 857 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 858
 859 } // End namespace Foam
 860
 861 // ************************************************************************* //