src/dynamicMesh/meshCut/meshModifiers/boundaryCutter/boundaryCutter.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 "boundaryCutter.H"
  27 #include "polyMesh.H"
  28 #include "polyTopoChange.H"
  29 #include "polyAddCell.H"
  30 #include "polyAddFace.H"
  31 #include "polyAddPoint.H"
  32 #include "polyModifyFace.H"
  33 #include "polyModifyPoint.H"
  34 #include "mapPolyMesh.H"
  35 #include "meshTools.H"
  36
  37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  38
  39 defineTypeNameAndDebug(Foam::boundaryCutter, 0);
  40
  41
  42 // * * * * * * * * * * * * * Private Static Functions  * * * * * * * * * * * //
  43
  44
  45 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  46
  47 void Foam::boundaryCutter::getFaceInfo
  48 (
  49     const label faceI,
  50     label& patchID,
  51     label& zoneID,
  52     label& zoneFlip
  53 ) const
  54 {
  55     patchID = -1;
  56
  57     if (!mesh_.isInternalFace(faceI))
  58     {
  59         patchID = mesh_.boundaryMesh().whichPatch(faceI);
  60     }
  61
  62     zoneID = mesh_.faceZones().whichZone(faceI);
  63
  64     zoneFlip = false;
  65
  66     if (zoneID >= 0)
  67     {
  68         const faceZone& fZone = mesh_.faceZones()[zoneID];
  69
  70         zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
  71     }
  72 }
  73
  74
  75 // Adds additional vertices (from edge cutting) to face. Used for faces which
  76 // are not split but still might use edge that has been cut.
  77 Foam::face Foam::boundaryCutter::addEdgeCutsToFace
  78 (
  79     const label faceI,
  80     const Map<labelList>& edgeToAddedPoints
  81 ) const
  82 {
  83     const edgeList& edges = mesh_.edges();
  84     const face& f = mesh_.faces()[faceI];
  85     const labelList& fEdges = mesh_.faceEdges()[faceI];
  86
  87     // Storage for face
  88     DynamicList<label> newFace(2 * f.size());
  89
  90     forAll(f, fp)
  91     {
  92         // Duplicate face vertex .
  93         newFace.append(f[fp]);
  94
  95         // Check if edge has been cut.
  96         label v1 = f.nextLabel(fp);
  97
  98         label edgeI = meshTools::findEdge(edges, fEdges, f[fp], v1);
  99
 100         Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
 101
 102         if (fnd != edgeToAddedPoints.end())
 103         {
 104             // edge has been cut. Introduce new vertices. Check order.
 105             const labelList& addedPoints = fnd();
 106
 107             if (edges[edgeI].start() == f[fp])
 108             {
 109                 // Introduce in same order.
 110                 forAll(addedPoints, i)
 111                 {
 112                     newFace.append(addedPoints[i]);
 113                 }
 114             }
 115             else
 116             {
 117                 // Introduce in opposite order.
 118                 forAllReverse(addedPoints, i)
 119                 {
 120                     newFace.append(addedPoints[i]);
 121                 }
 122             }
 123         }
 124     }
 125
 126     face returnFace;
 127     returnFace.transfer(newFace);
 128
 129     if (debug)
 130     {
 131         Pout<< "addEdgeCutsToFace:" << nl
 132             << "    from : " << f << nl
 133             << "    to   : " << returnFace << endl;
 134     }
 135
 136     return returnFace;
 137 }
 138
 139
 140 void Foam::boundaryCutter::addFace
 141 (
 142     const label faceI,
 143     const face& newFace,
 144
 145     bool& modifiedFace,     // have we already 'used' faceI
 146     polyTopoChange& meshMod
 147 ) const
 148 {
 149     // Information about old face
 150     label patchID, zoneID, zoneFlip;
 151     getFaceInfo(faceI, patchID, zoneID, zoneFlip);
 152     label own = mesh_.faceOwner()[faceI];
 153     label masterPoint = mesh_.faces()[faceI][0];
 154
 155     if (!modifiedFace)
 156     {
 157         meshMod.setAction
 158         (
 159             polyModifyFace
 160             (
 161                 newFace,       // face
 162                 faceI,
 163                 own,           // owner
 164                 -1,            // neighbour
 165                 false,         // flux flip
 166                 patchID,       // patch for face
 167                 false,         // remove from zone
 168                 zoneID,        // zone for face
 169                 zoneFlip       // face zone flip
 170             )
 171         );
 172
 173         modifiedFace = true;
 174     }
 175     else
 176     {
 177         meshMod.setAction
 178         (
 179             polyAddFace
 180             (
 181                 newFace,       // face
 182                 own,           // owner
 183                 -1,            // neighbour
 184                 masterPoint,   // master point
 185                 -1,            // master edge
 186                 -1,            // master face for addition
 187                 false,         // flux flip
 188                 patchID,       // patch for face
 189                 zoneID,        // zone for face
 190                 zoneFlip       // face zone flip
 191             )
 192         );
 193     }
 194 }
 195
 196
 197
 198 // Splits a face using the cut edges and modified points
 199 bool Foam::boundaryCutter::splitFace
 200 (
 201     const label faceI,
 202     const Map<point>& pointToPos,
 203     const Map<labelList>& edgeToAddedPoints,
 204     polyTopoChange& meshMod
 205 ) const
 206 {
 207     const edgeList& edges = mesh_.edges();
 208     const face& f = mesh_.faces()[faceI];
 209     const labelList& fEdges = mesh_.faceEdges()[faceI];
 210
 211     // Count number of split edges and total number of splits.
 212     label nSplitEdges = 0;
 213     label nModPoints = 0;
 214     label nTotalSplits = 0;
 215
 216     forAll(f, fp)
 217     {
 218         if (pointToPos.found(f[fp]))
 219         {
 220             nModPoints++;
 221             nTotalSplits++;
 222         }
 223
 224         // Check if edge has been cut.
 225         label nextV = f.nextLabel(fp);
 226
 227         label edgeI = meshTools::findEdge(edges, fEdges, f[fp], nextV);
 228
 229         Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
 230
 231         if (fnd != edgeToAddedPoints.end())
 232         {
 233             nSplitEdges++;
 234             nTotalSplits += fnd().size();
 235         }
 236     }
 237
 238     if (debug)
 239     {
 240         Pout<< "Face:" << faceI
 241             << " nModPoints:" << nModPoints
 242             << " nSplitEdges:" << nSplitEdges
 243             << " nTotalSplits:" << nTotalSplits << endl;
 244     }
 245
 246     if (nSplitEdges == 0 && nModPoints == 0)
 247     {
 248         FatalErrorIn("boundaryCutter::splitFace") << "Problem : face:" << faceI
 249             << " nSplitEdges:" << nSplitEdges
 250             << " nTotalSplits:" << nTotalSplits
 251             << abort(FatalError);
 252         return false;
 253     }
 254     else if (nSplitEdges + nModPoints == 1)
 255     {
 256         // single or multiple cuts on a single edge or single modified point
 257         // Dont cut and let caller handle this.
 258         Warning << "Face " << faceI << " has only one edge cut " << endl;
 259         return false;
 260     }
 261     else
 262     {
 263         // So guaranteed to have two edges cut or points modified. Split face:
 264         // - find starting cut
 265         // - walk to next cut. Make face
 266         // - loop until face done.
 267
 268         // Information about old face
 269         label patchID, zoneID, zoneFlip;
 270         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
 271
 272         // Get face with new points on cut edges for ease of looping
 273         face extendedFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
 274
 275         // Find first added point. This is the starting vertex for splitting.
 276         label startFp = -1;
 277
 278         forAll(extendedFace, fp)
 279         {
 280             if (extendedFace[fp] >= mesh_.nPoints())
 281             {
 282                 startFp = fp;
 283                 break;
 284             }
 285         }
 286
 287         if (startFp == -1)
 288         {
 289             // No added point. Maybe there is a modified point?
 290             forAll(extendedFace, fp)
 291             {
 292                 if (pointToPos.found(extendedFace[fp]))
 293                 {
 294                     startFp = fp;
 295                     break;
 296                 }
 297             }
 298         }
 299
 300         if (startFp == -1)
 301         {
 302             FatalErrorIn("boundaryCutter::splitFace")
 303                 << "Problem" << abort(FatalError);
 304         }
 305
 306         // Have we already modified existing face (first face gets done
 307         // as modification; all following ones as polyAddFace)
 308         bool modifiedFace = false;
 309
 310         // Example face:
 311         //    +--+
 312         //   /   |
 313         //  /    |
 314         // +     +
 315         //  \    |
 316         //   \   |
 317         //    +--+
 318         //
 319         // Needs to get split into:
 320         // - three 'side' faces a,b,c
 321         // - one middle face d
 322         //    +--+
 323         //   /|\A|
 324         //  / | \|
 325         // + C|D +
 326         //  \ | /|
 327         //   \|/B|
 328         //    +--+
 329
 330
 331         // Storage for new face
 332         DynamicList<label> newFace(extendedFace.size());
 333
 334         label fp = startFp;
 335
 336         forAll(extendedFace, i)
 337         {
 338             label pointI = extendedFace[fp];
 339
 340             newFace.append(pointI);
 341
 342             if
 343             (
 344                 newFace.size() > 2
 345              && (
 346                     pointI >= mesh_.nPoints()
 347                  || pointToPos.found(pointI)
 348                 )
 349             )
 350             {
 351                 // Enough vertices to create a face from.
 352                 face tmpFace;
 353                 tmpFace.transfer(newFace);
 354
 355                 // Add face tmpFace
 356                 addFace(faceI, tmpFace, modifiedFace, meshMod);
 357
 358                 // Starting point is also the starting point for the new face
 359                 newFace.append(extendedFace[startFp]);
 360                 newFace.append(extendedFace[fp]);
 361             }
 362
 363             fp = (fp+1) % extendedFace.size();
 364         }
 365
 366         // Check final face.
 367         if (newFace.size() > 2)
 368         {
 369             // Enough vertices to create a face from.
 370             face tmpFace;
 371             tmpFace.transfer(newFace);
 372
 373             // Add face tmpFace
 374             addFace(faceI, tmpFace, modifiedFace, meshMod);
 375         }
 376
 377         // Split something
 378         return true;
 379     }
 380 }
 381
 382
 383 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 384
 385 // Construct from components
 386 Foam::boundaryCutter::boundaryCutter(const polyMesh& mesh)
 387 :
 388     mesh_(mesh),
 389     edgeAddedPoints_(),
 390     faceAddedPoint_()
 391 {}
 392
 393
 394 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 395
 396 Foam::boundaryCutter::~boundaryCutter()
 397 {}
 398
 399
 400 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 401
 402 void Foam::boundaryCutter::setRefinement
 403 (
 404     const Map<point>& pointToPos,
 405     const Map<List<point> >& edgeToCuts,
 406     const Map<labelPair>& faceToSplit,
 407     const Map<point>& faceToFeaturePoint,
 408     polyTopoChange& meshMod
 409 )
 410 {
 411     // Clear and size maps here since mesh size will change.
 412     edgeAddedPoints_.clear();
 413
 414     faceAddedPoint_.clear();
 415     faceAddedPoint_.resize(faceToFeaturePoint.size());
 416
 417
 418     //
 419     // Points that just need to be moved
 420     // Note: could just as well be handled outside of setRefinement.
 421     //
 422
 423     forAllConstIter(Map<point>, pointToPos, iter)
 424     {
 425         meshMod.setAction
 426         (
 427             polyModifyPoint
 428             (
 429                 iter.key(), // point
 430                 iter(),     // position
 431                 false,      // no zone
 432                 -1,         // zone for point
 433                 true        // supports a cell
 434             )
 435         );
 436     }
 437
 438
 439     //
 440     // Add new points along cut edges.
 441     //
 442
 443     // Map from edge label to sorted list of points
 444     Map<labelList> edgeToAddedPoints(edgeToCuts.size());
 445
 446     forAllConstIter(Map<List<point> >, edgeToCuts, iter)
 447     {
 448         label edgeI = iter.key();
 449
 450         const edge& e = mesh_.edges()[edgeI];
 451
 452         // Sorted (from start to end) list of cuts on edge
 453         const List<point>& cuts = iter();
 454
 455         forAll(cuts, cutI)
 456         {
 457             // point on feature to move to
 458             const point& featurePoint = cuts[cutI];
 459
 460             label addedPointI =
 461                 meshMod.setAction
 462                 (
 463                     polyAddPoint
 464                     (
 465                         featurePoint,               // point
 466                         e.start(),                  // master point
 467                         -1,                         // zone for point
 468                         true                        // supports a cell
 469                     )
 470                 );
 471
 472             Map<labelList>::iterator fnd = edgeToAddedPoints.find(edgeI);
 473
 474             if (fnd != edgeToAddedPoints.end())
 475             {
 476                 labelList& addedPoints = fnd();
 477
 478                 label sz = addedPoints.size();
 479                 addedPoints.setSize(sz+1);
 480                 addedPoints[sz] = addedPointI;
 481             }
 482             else
 483             {
 484                 edgeToAddedPoints.insert(edgeI, labelList(1, addedPointI));
 485             }
 486
 487             if (debug)
 488             {
 489                 Pout<< "Added point " << addedPointI << " for edge " << edgeI
 490                     << " with cuts:" << edgeToAddedPoints[edgeI] << endl;
 491             }
 492         }
 493     }
 494
 495
 496     //
 497     // Introduce feature points.
 498     //
 499
 500     forAllConstIter(Map<point>, faceToFeaturePoint, iter)
 501     {
 502         label faceI = iter.key();
 503
 504         const face& f = mesh_.faces()[faceI];
 505
 506         if (faceToSplit.found(faceI))
 507         {
 508             FatalErrorIn("boundaryCutter::setRefinement")
 509                 << "Face " << faceI << " vertices " << f
 510                 << " is both marked for face-centre decomposition and"
 511                 << " diagonal splitting."
 512                 << abort(FatalError);
 513         }
 514
 515         if (mesh_.isInternalFace(faceI))
 516         {
 517             FatalErrorIn("boundaryCutter::setRefinement")
 518                 << "Face " << faceI << " vertices " << f
 519                 << " is not an external face. Cannot split it"
 520                 << abort(FatalError);
 521         }
 522
 523         label addedPointI =
 524             meshMod.setAction
 525             (
 526                 polyAddPoint
 527                 (
 528                     iter(), // point
 529                     f[0],   // master point
 530                     -1,     // zone for point
 531                     true    // supports a cell
 532                 )
 533             );
 534         faceAddedPoint_.insert(faceI, addedPointI);
 535
 536         if (debug)
 537         {
 538             Pout<< "Added point " << addedPointI << " for feature point "
 539                 << iter() << " on face " << faceI << " with centre "
 540                 << mesh_.faceCentres()[faceI] << endl;
 541         }
 542     }
 543
 544
 545     //
 546     // Split or retriangulate faces
 547     //
 548
 549
 550     // Maintain whether face has been updated (for -split edges
 551     // -new owner/neighbour)
 552     boolList faceUptodate(mesh_.nFaces(), false);
 553
 554
 555     // Triangulate faces containing feature points
 556     forAllConstIter(Map<label>, faceAddedPoint_, iter)
 557     {
 558         label faceI = iter.key();
 559
 560         // Get face with new points on cut edges.
 561         face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
 562
 563         label addedPointI = iter();
 564
 565         // Information about old face
 566         label patchID, zoneID, zoneFlip;
 567         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
 568         label own = mesh_.faceOwner()[faceI];
 569         label masterPoint = mesh_.faces()[faceI][0];
 570
 571         // Triangulate face around mid point
 572
 573         face tri(3);
 574
 575         forAll(newFace, fp)
 576         {
 577             label nextV = newFace.nextLabel(fp);
 578
 579             tri[0] = newFace[fp];
 580             tri[1] = nextV;
 581             tri[2] = addedPointI;
 582
 583             if (fp == 0)
 584             {
 585                 // Modify the existing face.
 586                 meshMod.setAction
 587                 (
 588                     polyModifyFace
 589                     (
 590                         tri,                        // face
 591                         faceI,
 592                         own,                        // owner
 593                         -1,                         // neighbour
 594                         false,                      // flux flip
 595                         patchID,                    // patch for face
 596                         false,                      // remove from zone
 597                         zoneID,                     // zone for face
 598                         zoneFlip                    // face zone flip
 599                     )
 600                 );
 601             }
 602             else
 603             {
 604                 // Add additional faces
 605                 meshMod.setAction
 606                 (
 607                     polyAddFace
 608                     (
 609                         tri,                        // face
 610                         own,                        // owner
 611                         -1,                         // neighbour
 612                         masterPoint,                // master point
 613                         -1,                         // master edge
 614                         -1,                         // master face for addition
 615                         false,                      // flux flip
 616                         patchID,                    // patch for face
 617                         zoneID,                     // zone for face
 618                         zoneFlip                    // face zone flip
 619                     )
 620                 );
 621             }
 622         }
 623
 624         faceUptodate[faceI] = true;
 625     }
 626
 627
 628     // Diagonally split faces
 629     forAllConstIter(Map<labelPair>, faceToSplit, iter)
 630     {
 631         label faceI = iter.key();
 632
 633         const face& f = mesh_.faces()[faceI];
 634
 635         if (faceAddedPoint_.found(faceI))
 636         {
 637             FatalErrorIn("boundaryCutter::setRefinement")
 638                 << "Face " << faceI << " vertices " << f
 639                 << " is both marked for face-centre decomposition and"
 640                 << " diagonal splitting."
 641                 << abort(FatalError);
 642         }
 643
 644
 645         // Get face with new points on cut edges.
 646         face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
 647
 648         // Information about old face
 649         label patchID, zoneID, zoneFlip;
 650         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
 651         label own = mesh_.faceOwner()[faceI];
 652         label masterPoint = mesh_.faces()[faceI][0];
 653
 654         // Split face from one side of diagonal to other.
 655         const labelPair& diag = iter();
 656
 657         label fp0 = findIndex(newFace, f[diag[0]]);
 658         label fp1 = findIndex(newFace, f[diag[1]]);
 659
 660         if (fp0 == -1 || fp1 == -1 || fp0 == fp1)
 661         {
 662             FatalErrorIn("boundaryCutter::setRefinement")
 663                 << "Problem : Face " << faceI << " vertices " << f
 664                 << " newFace:" << newFace << " diagonal:" << f[diag[0]]
 665                 << ' ' << f[diag[1]]
 666                 << abort(FatalError);
 667         }
 668
 669         // Replace existing face by newFace from fp0 to fp1 and add new one
 670         // from fp1 to fp0.
 671
 672         DynamicList<label> newVerts(newFace.size());
 673
 674         // Get vertices from fp0 to (and including) fp1
 675         label fp = fp0;
 676
 677         do
 678         {
 679             newVerts.append(newFace[fp]);
 680
 681             fp = (fp == newFace.size()-1 ? 0 : fp+1);
 682
 683         } while (fp != fp1);
 684
 685         newVerts.append(newFace[fp1]);
 686
 687
 688         // Modify the existing face.
 689         meshMod.setAction
 690         (
 691             polyModifyFace
 692             (
 693                 face(newVerts.shrink()),    // face
 694                 faceI,
 695                 own,                        // owner
 696                 -1,                         // neighbour
 697                 false,                      // flux flip
 698                 patchID,                    // patch for face
 699                 false,                      // remove from zone
 700                 zoneID,                     // zone for face
 701                 zoneFlip                    // face zone flip
 702             )
 703         );
 704
 705
 706         newVerts.clear();
 707
 708         // Get vertices from fp1 to (and including) fp0
 709
 710         do
 711         {
 712             newVerts.append(newFace[fp]);
 713
 714             fp = (fp == newFace.size()-1 ? 0 : fp+1);
 715
 716         } while (fp != fp0);
 717
 718         newVerts.append(newFace[fp0]);
 719
 720         // Add additional face
 721         meshMod.setAction
 722         (
 723             polyAddFace
 724             (
 725                 face(newVerts.shrink()),    // face
 726                 own,                        // owner
 727                 -1,                         // neighbour
 728                 masterPoint,                // master point
 729                 -1,                         // master edge
 730                 -1,                         // master face for addition
 731                 false,                      // flux flip
 732                 patchID,                    // patch for face
 733                 zoneID,                     // zone for face
 734                 zoneFlip                    // face zone flip
 735             )
 736         );
 737
 738         faceUptodate[faceI] = true;
 739     }
 740
 741
 742     // Split external faces without feature point but using cut edges.
 743     // Does right handed walk but not really.
 744     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
 745     {
 746         label edgeI = iter.key();
 747
 748         const labelList& eFaces = mesh_.edgeFaces()[edgeI];
 749
 750         forAll(eFaces, i)
 751         {
 752             label faceI = eFaces[i];
 753
 754             if (!faceUptodate[faceI] && !mesh_.isInternalFace(faceI))
 755             {
 756                 // Is external face so split
 757                 if (splitFace(faceI, pointToPos, edgeToAddedPoints, meshMod))
 758                 {
 759                     // Successfull split
 760                     faceUptodate[faceI] = true;
 761                 }
 762             }
 763         }
 764     }
 765
 766
 767     // Add cut edges (but don't split) any other faces using any cut edge.
 768     // These can be external faces where splitFace hasn't cut them or
 769     // internal faces.
 770     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
 771     {
 772         label edgeI = iter.key();
 773
 774         const labelList& eFaces = mesh_.edgeFaces()[edgeI];
 775
 776         forAll(eFaces, i)
 777         {
 778             label faceI = eFaces[i];
 779
 780             if (!faceUptodate[faceI])
 781             {
 782                 // Renumber face to include split edges.
 783                 face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
 784
 785                 label own = mesh_.faceOwner()[faceI];
 786
 787                 label nei = -1;
 788
 789                 if (mesh_.isInternalFace(faceI))
 790                 {
 791                     nei = mesh_.faceNeighbour()[faceI];
 792                 }
 793
 794                 label patchID, zoneID, zoneFlip;
 795                 getFaceInfo(faceI, patchID, zoneID, zoneFlip);
 796
 797                 meshMod.setAction
 798                 (
 799                     polyModifyFace
 800                     (
 801                         newFace,            // modified face
 802                         faceI,              // label of face being modified
 803                         own,                // owner
 804                         nei,                // neighbour
 805                         false,              // face flip
 806                         patchID,            // patch for face
 807                         false,              // remove from zone
 808                         zoneID,             // zone for face
 809                         zoneFlip            // face flip in zone
 810                     )
 811                 );
 812
 813                 faceUptodate[faceI] = true;
 814             }
 815         }
 816     }
 817
 818     // Convert edge to points storage from edge labels (not preserved)
 819     // to point labels
 820     edgeAddedPoints_.resize(edgeToCuts.size());
 821
 822     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
 823     {
 824         edgeAddedPoints_.insert(mesh_.edges()[iter.key()], iter());
 825     }
 826 }
 827
 828
 829 void Foam::boundaryCutter::updateMesh(const mapPolyMesh& morphMap)
 830 {
 831     // Update stored labels for mesh change.
 832
 833     //
 834     // Do faceToAddedPoint
 835     //
 836
 837     {
 838         // Create copy since we're deleting entries.
 839         Map<label> newAddedPoints(faceAddedPoint_.size());
 840
 841         forAllConstIter(Map<label>, faceAddedPoint_, iter)
 842         {
 843             label oldFaceI = iter.key();
 844
 845             label newFaceI = morphMap.reverseFaceMap()[oldFaceI];
 846
 847             label oldPointI = iter();
 848
 849             label newPointI = morphMap.reversePointMap()[oldPointI];
 850
 851             if (newFaceI >= 0 && newPointI >= 0)
 852             {
 853                 newAddedPoints.insert(newFaceI, newPointI);
 854             }
 855         }
 856
 857         // Copy
 858         faceAddedPoint_.transfer(newAddedPoints);
 859     }
 860
 861
 862     //
 863     // Do edgeToAddedPoints
 864     //
 865
 866
 867     {
 868         // Create copy since we're deleting entries
 869         HashTable<labelList, edge, Hash<edge> >
 870             newEdgeAddedPoints(edgeAddedPoints_.size());
 871
 872         for
 873         (
 874             HashTable<labelList, edge, Hash<edge> >::const_iterator iter =
 875                 edgeAddedPoints_.begin();
 876             iter != edgeAddedPoints_.end();
 877             ++iter
 878         )
 879         {
 880             const edge& e = iter.key();
 881
 882             label newStart = morphMap.reversePointMap()[e.start()];
 883
 884             label newEnd = morphMap.reversePointMap()[e.end()];
 885
 886             if (newStart >= 0 && newEnd >= 0)
 887             {
 888                 const labelList& addedPoints = iter();
 889
 890                 labelList newAddedPoints(addedPoints.size());
 891                 label newI = 0;
 892
 893                 forAll(addedPoints, i)
 894                 {
 895                     label newAddedPointI =
 896                         morphMap.reversePointMap()[addedPoints[i]];
 897
 898                     if (newAddedPointI >= 0)
 899                     {
 900                         newAddedPoints[newI++] = newAddedPointI;
 901                     }
 902                 }
 903                 if (newI > 0)
 904                 {
 905                     newAddedPoints.setSize(newI);
 906
 907                     edge newE = edge(newStart, newEnd);
 908
 909                     newEdgeAddedPoints.insert(newE, newAddedPoints);
 910                 }
 911             }
 912         }
 913
 914         // Copy
 915         edgeAddedPoints_.transfer(newEdgeAddedPoints);
 916     }
 917 }
 918
 919
 920 // ************************************************************************* //