FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

boundaryCutter.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 Description
00025 
00026 \*---------------------------------------------------------------------------*/
00027 
00028 #include "boundaryCutter.H"
00029 #include <OpenFOAM/polyMesh.H>
00030 #include <dynamicMesh/polyTopoChange.H>
00031 #include <dynamicMesh/polyAddCell.H>
00032 #include <dynamicMesh/polyAddFace.H>
00033 #include <dynamicMesh/polyAddPoint.H>
00034 #include <dynamicMesh/polyModifyFace.H>
00035 #include <dynamicMesh/polyModifyPoint.H>
00036 #include <OpenFOAM/mapPolyMesh.H>
00037 #include <meshTools/meshTools.H>
00038 
00039 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00040 
00041 namespace Foam
00042 {
00043 
00044 defineTypeNameAndDebug(boundaryCutter, 0);
00045 
00046 }
00047 
00048 
00049 // * * * * * * * * * * * * * Private Static Functions  * * * * * * * * * * * //
00050 
00051 
00052 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00053 
00054 void Foam::boundaryCutter::getFaceInfo
00055 (
00056     const label faceI,
00057     label& patchID,
00058     label& zoneID,
00059     label& zoneFlip
00060 ) const
00061 {
00062     patchID = -1;
00063 
00064     if (!mesh_.isInternalFace(faceI))
00065     {
00066         patchID = mesh_.boundaryMesh().whichPatch(faceI);
00067     }
00068 
00069     zoneID = mesh_.faceZones().whichZone(faceI);
00070 
00071     zoneFlip = false;
00072 
00073     if (zoneID >= 0)
00074     {
00075         const faceZone& fZone = mesh_.faceZones()[zoneID];
00076 
00077         zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
00078     }
00079 }
00080 
00081 
00082 // Adds additional vertices (from edge cutting) to face. Used for faces which
00083 // are not split but still might use edge that has been cut.
00084 Foam::face Foam::boundaryCutter::addEdgeCutsToFace
00085 (
00086     const label faceI,
00087     const Map<labelList>& edgeToAddedPoints
00088 ) const
00089 {
00090     const edgeList& edges = mesh_.edges();
00091     const face& f = mesh_.faces()[faceI];
00092     const labelList& fEdges = mesh_.faceEdges()[faceI];
00093 
00094     // Storage for face
00095     DynamicList<label> newFace(2 * f.size());
00096 
00097     forAll(f, fp)
00098     {
00099         // Duplicate face vertex .
00100         newFace.append(f[fp]);
00101 
00102         // Check if edge has been cut.
00103         label v1 = f.nextLabel(fp);
00104 
00105         label edgeI = meshTools::findEdge(edges, fEdges, f[fp], v1);
00106 
00107         Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
00108 
00109         if (fnd != edgeToAddedPoints.end())
00110         {
00111             // edge has been cut. Introduce new vertices. Check order.
00112             const labelList& addedPoints = fnd();
00113 
00114             if (edges[edgeI].start() == f[fp])
00115             {
00116                 // Introduce in same order.
00117                 forAll(addedPoints, i)
00118                 {
00119                     newFace.append(addedPoints[i]);
00120                 }
00121             }
00122             else
00123             {
00124                 // Introduce in opposite order.
00125                 forAllReverse(addedPoints, i)
00126                 {
00127                     newFace.append(addedPoints[i]);
00128                 }
00129             }
00130         }
00131     }
00132 
00133     face returnFace;
00134     returnFace.transfer(newFace);
00135 
00136     if (debug)
00137     {
00138         Pout<< "addEdgeCutsToFace:" << nl
00139             << "    from : " << f << nl
00140             << "    to   : " << returnFace << endl;
00141     }
00142 
00143     return returnFace;
00144 }
00145 
00146 
00147 void Foam::boundaryCutter::addFace
00148 (
00149     const label faceI,
00150     const face& newFace,
00151 
00152     bool& modifiedFace,     // have we already 'used' faceI
00153     polyTopoChange& meshMod
00154 ) const
00155 {
00156     // Information about old face
00157     label patchID, zoneID, zoneFlip;
00158     getFaceInfo(faceI, patchID, zoneID, zoneFlip);
00159     label own = mesh_.faceOwner()[faceI];
00160     label masterPoint = mesh_.faces()[faceI][0];
00161     
00162     if (!modifiedFace)
00163     {
00164         meshMod.setAction
00165         (
00166             polyModifyFace
00167             (
00168                 newFace,       // face
00169                 faceI,
00170                 own,           // owner
00171                 -1,            // neighbour
00172                 false,         // flux flip
00173                 patchID,       // patch for face
00174                 false,         // remove from zone
00175                 zoneID,        // zone for face
00176                 zoneFlip       // face zone flip
00177             )
00178         );
00179 
00180         modifiedFace = true;
00181     }
00182     else
00183     {
00184         meshMod.setAction
00185         (
00186             polyAddFace
00187             (
00188                 newFace,       // face
00189                 own,           // owner
00190                 -1,            // neighbour
00191                 masterPoint,   // master point
00192                 -1,            // master edge
00193                 -1,            // master face for addition
00194                 false,         // flux flip
00195                 patchID,       // patch for face
00196                 zoneID,        // zone for face
00197                 zoneFlip       // face zone flip
00198             )
00199         );
00200     }
00201 }
00202 
00203 
00204 
00205 // Splits a face using the cut edges and modified points 
00206 bool Foam::boundaryCutter::splitFace
00207 (
00208     const label faceI,
00209     const Map<point>& pointToPos,
00210     const Map<labelList>& edgeToAddedPoints,
00211     polyTopoChange& meshMod
00212 ) const
00213 {
00214     const edgeList& edges = mesh_.edges();
00215     const face& f = mesh_.faces()[faceI];
00216     const labelList& fEdges = mesh_.faceEdges()[faceI];
00217 
00218     // Count number of split edges and total number of splits.
00219     label nSplitEdges = 0;
00220     label nModPoints = 0;
00221     label nTotalSplits = 0;
00222 
00223     forAll(f, fp)
00224     {
00225         if (pointToPos.found(f[fp]))
00226         {
00227             nModPoints++;
00228             nTotalSplits++;
00229         }
00230 
00231         // Check if edge has been cut.
00232         label nextV = f.nextLabel(fp);
00233 
00234         label edgeI = meshTools::findEdge(edges, fEdges, f[fp], nextV);
00235 
00236         Map<labelList>::const_iterator fnd = edgeToAddedPoints.find(edgeI);
00237 
00238         if (fnd != edgeToAddedPoints.end())
00239         {
00240             nSplitEdges++;
00241             nTotalSplits += fnd().size();
00242         }
00243     }
00244 
00245     if (debug)
00246     {
00247         Pout<< "Face:" << faceI
00248             << " nModPoints:" << nModPoints
00249             << " nSplitEdges:" << nSplitEdges
00250             << " nTotalSplits:" << nTotalSplits << endl;
00251     }
00252 
00253     if (nSplitEdges == 0 && nModPoints == 0)
00254     {
00255         FatalErrorIn("boundaryCutter::splitFace") << "Problem : face:" << faceI
00256             << " nSplitEdges:" << nSplitEdges
00257             << " nTotalSplits:" << nTotalSplits
00258             << abort(FatalError);
00259         return false;
00260     }
00261     else if (nSplitEdges + nModPoints == 1)
00262     {
00263         // single or multiple cuts on a single edge or single modified point
00264         // Dont cut and let caller handle this.
00265         Warning << "Face " << faceI << " has only one edge cut " << endl;
00266         return false;
00267     }
00268     else
00269     {
00270         // So guaranteed to have two edges cut or points modified. Split face:
00271         // - find starting cut
00272         // - walk to next cut. Make face
00273         // - loop until face done.
00274 
00275         // Information about old face
00276         label patchID, zoneID, zoneFlip;
00277         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
00278 
00279         // Get face with new points on cut edges for ease of looping
00280         face extendedFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
00281 
00282         // Find first added point. This is the starting vertex for splitting.
00283         label startFp = -1;
00284 
00285         forAll(extendedFace, fp)
00286         {
00287             if (extendedFace[fp] >= mesh_.nPoints())
00288             {
00289                 startFp = fp;
00290                 break;
00291             }
00292         }
00293 
00294         if (startFp == -1)
00295         {
00296             // No added point. Maybe there is a modified point?
00297             forAll(extendedFace, fp)
00298             {
00299                 if (pointToPos.found(extendedFace[fp]))
00300                 {
00301                     startFp = fp;
00302                     break;
00303                 }
00304             }
00305         }
00306 
00307         if (startFp == -1)
00308         {
00309             FatalErrorIn("boundaryCutter::splitFace")
00310                 << "Problem" << abort(FatalError);
00311         }
00312 
00313         // Have we already modified existing face (first face gets done
00314         // as modification; all following ones as polyAddFace)
00315         bool modifiedFace = false;
00316 
00317         // Example face:
00318         //    +--+
00319         //   /   |
00320         //  /    |
00321         // +     +
00322         //  \    |
00323         //   \   |
00324         //    +--+
00325         //
00326         // Needs to get split into:
00327         // - three 'side' faces a,b,c
00328         // - one middle face d
00329         //    +--+
00330         //   /|\A|
00331         //  / | \|
00332         // + C|D +
00333         //  \ | /|
00334         //   \|/B|
00335         //    +--+
00336 
00337 
00338         // Storage for new face
00339         DynamicList<label> newFace(extendedFace.size());
00340 
00341         label fp = startFp;
00342 
00343         forAll(extendedFace, i)
00344         {
00345             label pointI = extendedFace[fp];
00346 
00347             newFace.append(pointI);
00348 
00349             if
00350             (
00351                 newFace.size() > 2
00352              && (
00353                     pointI >= mesh_.nPoints()
00354                  || pointToPos.found(pointI)
00355                 )
00356             )
00357             {
00358                 // Enough vertices to create a face from.
00359                 face tmpFace;
00360                 tmpFace.transfer(newFace);
00361 
00362                 // Add face tmpFace
00363                 addFace(faceI, tmpFace, modifiedFace, meshMod);
00364 
00365                 // Starting point is also the starting point for the new face
00366                 newFace.append(extendedFace[startFp]);
00367                 newFace.append(extendedFace[fp]);
00368             }
00369 
00370             fp = (fp+1) % extendedFace.size();
00371         }
00372 
00373         // Check final face.
00374         if (newFace.size() > 2)
00375         {
00376             // Enough vertices to create a face from.
00377             face tmpFace;
00378             tmpFace.transfer(newFace);
00379 
00380             // Add face tmpFace
00381             addFace(faceI, tmpFace, modifiedFace, meshMod);
00382         }
00383 
00384         // Split something
00385         return true;
00386     }
00387 }
00388 
00389 
00390 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00391 
00392 // Construct from components
00393 Foam::boundaryCutter::boundaryCutter(const polyMesh& mesh)
00394 :
00395     mesh_(mesh),
00396     edgeAddedPoints_(),
00397     faceAddedPoint_()
00398 {}
00399 
00400 
00401 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00402 
00403 Foam::boundaryCutter::~boundaryCutter()
00404 {}
00405 
00406 
00407 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00408 
00409 void Foam::boundaryCutter::setRefinement
00410 (
00411     const Map<point>& pointToPos,
00412     const Map<List<point> >& edgeToCuts,
00413     const Map<labelPair>& faceToSplit,
00414     const Map<point>& faceToFeaturePoint,
00415     polyTopoChange& meshMod
00416 )
00417 {
00418     // Clear and size maps here since mesh size will change.
00419     edgeAddedPoints_.clear();
00420 
00421     faceAddedPoint_.clear();
00422     faceAddedPoint_.resize(faceToFeaturePoint.size());
00423 
00424 
00425     //
00426     // Points that just need to be moved
00427     // Note: could just as well be handled outside of setRefinement.
00428     //
00429 
00430     forAllConstIter(Map<point>, pointToPos, iter)
00431     {
00432         meshMod.setAction
00433         (
00434             polyModifyPoint
00435             (
00436                 iter.key(), // point
00437                 iter(),     // position
00438                 false,      // no zone
00439                 -1,         // zone for point
00440                 true        // supports a cell
00441             )
00442         );
00443     }
00444 
00445 
00446     //
00447     // Add new points along cut edges.
00448     //
00449 
00450     // Map from edge label to sorted list of points
00451     Map<labelList> edgeToAddedPoints(edgeToCuts.size());
00452 
00453     forAllConstIter(Map<List<point> >, edgeToCuts, iter)
00454     {
00455         label edgeI = iter.key();
00456 
00457         const edge& e = mesh_.edges()[edgeI];
00458 
00459         // Sorted (from start to end) list of cuts on edge
00460         const List<point>& cuts = iter();
00461 
00462         forAll(cuts, cutI)
00463         {
00464             // point on feature to move to
00465             const point& featurePoint = cuts[cutI];
00466 
00467             label addedPointI =
00468                 meshMod.setAction
00469                 (
00470                     polyAddPoint
00471                     (
00472                         featurePoint,               // point
00473                         e.start(),                  // master point
00474                         -1,                         // zone for point
00475                         true                        // supports a cell
00476                     )
00477                 );
00478 
00479             Map<labelList>::iterator fnd = edgeToAddedPoints.find(edgeI);
00480 
00481             if (fnd != edgeToAddedPoints.end())
00482             {
00483                 labelList& addedPoints = fnd();
00484 
00485                 label sz = addedPoints.size();
00486                 addedPoints.setSize(sz+1);
00487                 addedPoints[sz] = addedPointI;
00488             }
00489             else
00490             {
00491                 edgeToAddedPoints.insert(edgeI, labelList(1, addedPointI));
00492             }
00493 
00494             if (debug)
00495             {
00496                 Pout<< "Added point " << addedPointI << " for edge " << edgeI
00497                     << " with cuts:" << edgeToAddedPoints[edgeI] << endl;
00498             }
00499         }
00500     }
00501 
00502 
00503     //
00504     // Introduce feature points.
00505     //
00506 
00507     forAllConstIter(Map<point>, faceToFeaturePoint, iter)
00508     {
00509         label faceI = iter.key();
00510 
00511         const face& f = mesh_.faces()[faceI];
00512 
00513         if (faceToSplit.found(faceI))
00514         {
00515             FatalErrorIn("boundaryCutter::setRefinement")
00516                 << "Face " << faceI << " vertices " << f
00517                 << " is both marked for face-centre decomposition and"
00518                 << " diagonal splitting."
00519                 << abort(FatalError);
00520         }
00521 
00522         if (mesh_.isInternalFace(faceI))
00523         {
00524             FatalErrorIn("boundaryCutter::setRefinement")
00525                 << "Face " << faceI << " vertices " << f
00526                 << " is not an external face. Cannot split it"
00527                 << abort(FatalError);
00528         }
00529 
00530         label addedPointI =
00531             meshMod.setAction
00532             (
00533                 polyAddPoint
00534                 (
00535                     iter(), // point
00536                     f[0],   // master point
00537                     -1,     // zone for point
00538                     true    // supports a cell
00539                 )
00540             );
00541         faceAddedPoint_.insert(faceI, addedPointI);
00542 
00543         if (debug)
00544         {
00545             Pout<< "Added point " << addedPointI << " for feature point "
00546                 << iter() << " on face " << faceI << " with centre "
00547                 << mesh_.faceCentres()[faceI] << endl;
00548         }
00549     }
00550 
00551 
00552     //
00553     // Split or retriangulate faces
00554     //
00555 
00556 
00557     // Maintain whether face has been updated (for -split edges
00558     // -new owner/neighbour)
00559     boolList faceUptodate(mesh_.nFaces(), false);
00560 
00561 
00562     // Triangulate faces containing feature points
00563     forAllConstIter(Map<label>, faceAddedPoint_, iter)
00564     {
00565         label faceI = iter.key();
00566 
00567         // Get face with new points on cut edges.
00568         face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
00569 
00570         label addedPointI = iter();
00571 
00572         // Information about old face
00573         label patchID, zoneID, zoneFlip;
00574         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
00575         label own = mesh_.faceOwner()[faceI];
00576         label masterPoint = mesh_.faces()[faceI][0];
00577 
00578         // Triangulate face around mid point
00579 
00580         face tri(3);
00581 
00582         forAll(newFace, fp)
00583         {
00584             label nextV = newFace.nextLabel(fp);
00585 
00586             tri[0] = newFace[fp];
00587             tri[1] = nextV;
00588             tri[2] = addedPointI;
00589 
00590             if (fp == 0)
00591             {
00592                 // Modify the existing face.
00593                 meshMod.setAction
00594                 (
00595                     polyModifyFace
00596                     (
00597                         tri,                        // face
00598                         faceI,
00599                         own,                        // owner
00600                         -1,                         // neighbour
00601                         false,                      // flux flip
00602                         patchID,                    // patch for face
00603                         false,                      // remove from zone
00604                         zoneID,                     // zone for face
00605                         zoneFlip                    // face zone flip
00606                     )
00607                 );
00608             }
00609             else
00610             {
00611                 // Add additional faces
00612                 meshMod.setAction
00613                 (
00614                     polyAddFace
00615                     (
00616                         tri,                        // face
00617                         own,                        // owner
00618                         -1,                         // neighbour
00619                         masterPoint,                // master point
00620                         -1,                         // master edge
00621                         -1,                         // master face for addition
00622                         false,                      // flux flip
00623                         patchID,                    // patch for face
00624                         zoneID,                     // zone for face
00625                         zoneFlip                    // face zone flip
00626                     )
00627                 );
00628             }
00629         }
00630 
00631         faceUptodate[faceI] = true;
00632     }
00633 
00634 
00635     // Diagonally split faces
00636     forAllConstIter(Map<labelPair>, faceToSplit, iter)
00637     {
00638         label faceI = iter.key();
00639 
00640         const face& f = mesh_.faces()[faceI];
00641 
00642         if (faceAddedPoint_.found(faceI))
00643         {
00644             FatalErrorIn("boundaryCutter::setRefinement")
00645                 << "Face " << faceI << " vertices " << f
00646                 << " is both marked for face-centre decomposition and"
00647                 << " diagonal splitting."
00648                 << abort(FatalError);
00649         }
00650 
00651 
00652         // Get face with new points on cut edges.
00653         face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
00654 
00655         // Information about old face
00656         label patchID, zoneID, zoneFlip;
00657         getFaceInfo(faceI, patchID, zoneID, zoneFlip);
00658         label own = mesh_.faceOwner()[faceI];
00659         label masterPoint = mesh_.faces()[faceI][0];
00660 
00661         // Split face from one side of diagonal to other.
00662         const labelPair& diag = iter();
00663 
00664         label fp0 = findIndex(newFace, f[diag[0]]);
00665         label fp1 = findIndex(newFace, f[diag[1]]);
00666 
00667         if (fp0 == -1 || fp1 == -1 || fp0 == fp1)
00668         {
00669             FatalErrorIn("boundaryCutter::setRefinement")
00670                 << "Problem : Face " << faceI << " vertices " << f
00671                 << " newFace:" << newFace << " diagonal:" << f[diag[0]]
00672                 << ' ' << f[diag[1]]
00673                 << abort(FatalError);
00674         }
00675 
00676         // Replace existing face by newFace from fp0 to fp1 and add new one
00677         // from fp1 to fp0.
00678 
00679         DynamicList<label> newVerts(newFace.size());
00680 
00681         // Get vertices from fp0 to (and including) fp1
00682         label fp = fp0;
00683 
00684         do
00685         {
00686             newVerts.append(newFace[fp]);
00687 
00688             fp = (fp == newFace.size()-1 ? 0 : fp+1);
00689 
00690         } while (fp != fp1);
00691 
00692         newVerts.append(newFace[fp1]);
00693 
00694 
00695         // Modify the existing face.
00696         meshMod.setAction
00697         (
00698             polyModifyFace
00699             (
00700                 face(newVerts.shrink()),    // face
00701                 faceI,
00702                 own,                        // owner
00703                 -1,                         // neighbour
00704                 false,                      // flux flip
00705                 patchID,                    // patch for face
00706                 false,                      // remove from zone
00707                 zoneID,                     // zone for face
00708                 zoneFlip                    // face zone flip
00709             )
00710         );
00711 
00712 
00713         newVerts.clear();
00714 
00715         // Get vertices from fp1 to (and including) fp0
00716 
00717         do
00718         {
00719             newVerts.append(newFace[fp]);
00720 
00721             fp = (fp == newFace.size()-1 ? 0 : fp+1);
00722 
00723         } while (fp != fp0);
00724 
00725         newVerts.append(newFace[fp0]);
00726 
00727         // Add additional face
00728         meshMod.setAction
00729         (
00730             polyAddFace
00731             (
00732                 face(newVerts.shrink()),    // face
00733                 own,                        // owner
00734                 -1,                         // neighbour
00735                 masterPoint,                // master point
00736                 -1,                         // master edge
00737                 -1,                         // master face for addition
00738                 false,                      // flux flip
00739                 patchID,                    // patch for face
00740                 zoneID,                     // zone for face
00741                 zoneFlip                    // face zone flip
00742             )
00743         );
00744 
00745         faceUptodate[faceI] = true;
00746     }
00747 
00748 
00749     // Split external faces without feature point but using cut edges.
00750     // Does right handed walk but not really.
00751     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
00752     {
00753         label edgeI = iter.key();
00754 
00755         const labelList& eFaces = mesh_.edgeFaces()[edgeI];
00756 
00757         forAll(eFaces, i)    
00758         {
00759             label faceI = eFaces[i];
00760 
00761             if (!faceUptodate[faceI] && !mesh_.isInternalFace(faceI))
00762             {
00763                 // Is external face so split
00764                 if (splitFace(faceI, pointToPos, edgeToAddedPoints, meshMod))
00765                 {
00766                     // Successfull split
00767                     faceUptodate[faceI] = true;
00768                 }
00769             }
00770         }
00771     }
00772 
00773 
00774     // Add cut edges (but don't split) any other faces using any cut edge.
00775     // These can be external faces where splitFace hasn't cut them or
00776     // internal faces.
00777     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
00778     {
00779         label edgeI = iter.key();
00780 
00781         const labelList& eFaces = mesh_.edgeFaces()[edgeI];
00782 
00783         forAll(eFaces, i)    
00784         {
00785             label faceI = eFaces[i];
00786 
00787             if (!faceUptodate[faceI])
00788             {
00789                 // Renumber face to include split edges.
00790                 face newFace(addEdgeCutsToFace(faceI, edgeToAddedPoints));
00791 
00792                 label own = mesh_.faceOwner()[faceI];
00793 
00794                 label nei = -1;
00795 
00796                 if (mesh_.isInternalFace(faceI))
00797                 {
00798                     nei = mesh_.faceNeighbour()[faceI];
00799                 }
00800 
00801                 label patchID, zoneID, zoneFlip;
00802                 getFaceInfo(faceI, patchID, zoneID, zoneFlip);
00803 
00804                 meshMod.setAction
00805                 (
00806                     polyModifyFace
00807                     (
00808                         newFace,            // modified face
00809                         faceI,              // label of face being modified
00810                         own,                // owner
00811                         nei,                // neighbour
00812                         false,              // face flip
00813                         patchID,            // patch for face
00814                         false,              // remove from zone
00815                         zoneID,             // zone for face
00816                         zoneFlip            // face flip in zone
00817                     )
00818                 );
00819 
00820                 faceUptodate[faceI] = true;
00821             }
00822         }
00823     }
00824 
00825     // Convert edge to points storage from edge labels (not preserved)
00826     // to point labels
00827     edgeAddedPoints_.resize(edgeToCuts.size());
00828 
00829     forAllConstIter(Map<labelList>, edgeToAddedPoints, iter)
00830     {
00831         edgeAddedPoints_.insert(mesh_.edges()[iter.key()], iter());
00832     }
00833 }
00834 
00835 
00836 void Foam::boundaryCutter::updateMesh(const mapPolyMesh& morphMap)
00837 {
00838     // Update stored labels for mesh change.
00839 
00840     //
00841     // Do faceToAddedPoint
00842     //
00843 
00844     {
00845         // Create copy since we're deleting entries.
00846         Map<label> newAddedPoints(faceAddedPoint_.size());
00847 
00848         forAllConstIter(Map<label>, faceAddedPoint_, iter)
00849         {
00850             label oldFaceI = iter.key();
00851 
00852             label newFaceI = morphMap.reverseFaceMap()[oldFaceI];
00853 
00854             label oldPointI = iter();
00855 
00856             label newPointI = morphMap.reversePointMap()[oldPointI];
00857 
00858             if (newFaceI >= 0 && newPointI >= 0)
00859             {
00860                 newAddedPoints.insert(newFaceI, newPointI);
00861             }
00862         }
00863 
00864         // Copy
00865         faceAddedPoint_.transfer(newAddedPoints);
00866     }
00867 
00868 
00869     //
00870     // Do edgeToAddedPoints
00871     //
00872 
00873 
00874     {
00875         // Create copy since we're deleting entries
00876         HashTable<labelList, edge, Hash<edge> >
00877             newEdgeAddedPoints(edgeAddedPoints_.size());
00878 
00879         for
00880         (
00881             HashTable<labelList, edge, Hash<edge> >::const_iterator iter =
00882                 edgeAddedPoints_.begin();
00883             iter != edgeAddedPoints_.end();
00884             ++iter
00885         )
00886         {
00887             const edge& e = iter.key();
00888 
00889             label newStart = morphMap.reversePointMap()[e.start()];
00890 
00891             label newEnd = morphMap.reversePointMap()[e.end()];
00892 
00893             if (newStart >= 0 && newEnd >= 0)
00894             {
00895                 const labelList& addedPoints = iter();
00896 
00897                 labelList newAddedPoints(addedPoints.size());
00898                 label newI = 0;
00899 
00900                 forAll(addedPoints, i)
00901                 {
00902                     label newAddedPointI =
00903                         morphMap.reversePointMap()[addedPoints[i]];
00904 
00905                     if (newAddedPointI >= 0)
00906                     {
00907                         newAddedPoints[newI++] = newAddedPointI;
00908                     }
00909                 }
00910                 if (newI > 0)
00911                 {
00912                     newAddedPoints.setSize(newI);
00913 
00914                     edge newE = edge(newStart, newEnd);
00915 
00916                     newEdgeAddedPoints.insert(newE, newAddedPoints);
00917                 }
00918             }
00919         }
00920 
00921         // Copy
00922         edgeAddedPoints_.transfer(newEdgeAddedPoints);
00923     }
00924 }
00925 
00926 
00927 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines