1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*----------------------------------------------------------------------------*/
27 #include "meshRefinement.H"
28 #include <dynamicMesh/combineFaces.H>
29 #include <dynamicMesh/polyTopoChange.H>
30 #include <dynamicMesh/removePoints.H>
32 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
34 // Merge faces that are in-line.
35 Foam::label Foam::meshRefinement::mergePatchFaces
38 const scalar concaveCos,
39 const labelList& patchIDs
42 // Patch face merging engine
43 combineFaces faceCombiner(mesh_);
45 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
47 // Pick up all candidate cells on boundary
48 labelHashSet boundaryCells(mesh_.nFaces()-mesh_.nInternalFaces());
52 label patchI = patchIDs[i];
54 const polyPatch& patch = patches[patchI];
60 boundaryCells.insert(mesh_.faceOwner()[patch.start()+i]);
65 // Get all sets of faces that can be merged
66 labelListList mergeSets
68 faceCombiner.getMergeSets
76 label nFaceSets = returnReduce(mergeSets.size(), sumOp<label>());
78 Info<< "mergePatchFaces : Merging " << nFaceSets
79 << " sets of faces." << endl;
83 // Topology changes container
84 polyTopoChange meshMod(mesh_);
86 // Merge all faces of a set into the first face of the set. Remove
88 faceCombiner.setRefinement(mergeSets, meshMod);
90 // Change the mesh (no inflation)
91 autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, true);
94 mesh_.updateMesh(map);
96 // Move mesh (since morphing does not do this)
97 if (map().hasMotionPoints())
99 mesh_.movePoints(map().preMotionPoints());
103 // Delete mesh volumes. No other way to do this?
109 mesh_.setInstance(oldInstance());
112 faceCombiner.updateMesh(map);
114 // Get the kept faces that need to be recalculated.
115 // Merging two boundary faces might shift the cell centre
116 // (unless the faces are absolutely planar)
117 labelHashSet retestFaces(6*mergeSets.size());
119 forAll(mergeSets, setI)
121 label oldMasterI = mergeSets[setI][0];
123 label faceI = map().reverseFaceMap()[oldMasterI];
125 // faceI is always uncoupled boundary face
126 const cell& cFaces = mesh_.cells()[mesh_.faceOwner()[faceI]];
130 retestFaces.insert(cFaces[i]);
133 updateMesh(map, retestFaces.toc());
141 // Remove points not used by any face or points used by only two faces where
142 // the edges are in line
143 Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::mergeEdges
148 // Point removal analysis engine
149 removePoints pointRemover(mesh_);
151 // Count usage of points
152 boolList pointCanBeDeleted;
153 label nRemove = pointRemover.countPointUsage(minCos, pointCanBeDeleted);
155 Info<< "Removing " << nRemove
156 << " straight edge points." << endl;
158 autoPtr<mapPolyMesh> map;
162 // Save my local faces that will change. These changed faces might
163 // cause a shift in the cell centre which needs to be retested.
164 // Have to do this before changing mesh since point will be removed.
165 labelHashSet retestOldFaces(nRemove / Pstream::nProcs());
168 const faceList& faces = mesh_.faces();
172 const face& f = faces[faceI];
176 if (pointCanBeDeleted[f[fp]])
178 retestOldFaces.insert(faceI);
185 // Topology changes container
186 polyTopoChange meshMod(mesh_);
188 pointRemover.setRefinement(pointCanBeDeleted, meshMod);
190 // Change the mesh (no inflation)
191 map = meshMod.changeMesh(mesh_, false, true);
194 mesh_.updateMesh(map);
196 // Move mesh (since morphing does not do this)
197 if (map().hasMotionPoints())
199 mesh_.movePoints(map().preMotionPoints());
203 // Delete mesh volumes. No other way to do this?
209 mesh_.setInstance(oldInstance());
212 pointRemover.updateMesh(map);
214 // Get the kept faces that need to be recalculated.
215 labelHashSet retestFaces(6*retestOldFaces.size());
217 const cellList& cells = mesh_.cells();
219 forAllConstIter(labelHashSet, retestOldFaces, iter)
221 label faceI = map().reverseFaceMap()[iter.key()];
223 const cell& ownFaces = cells[mesh_.faceOwner()[faceI]];
227 retestFaces.insert(ownFaces[i]);
230 if (mesh_.isInternalFace(faceI))
232 const cell& neiFaces = cells[mesh_.faceNeighbour()[faceI]];
236 retestFaces.insert(neiFaces[i]);
240 updateMesh(map, retestFaces.toc());
247 // ************************ vim: set sw=4 sts=4 et: ************************ //