splitMesh.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2016 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 splitMesh
29
30Group
31 grpMeshManipulationUtilities
32
33Description
34 Splits mesh by making internal faces external. Uses attachDetach.
35
36 Generates a meshModifier of the form:
37
38 Splitter
39 {
40 type attachDetach;
41 faceZoneName membraneFaces;
42 masterPatchName masterPatch;
43 slavePatchName slavePatch;
44 triggerTimes runTime.value();
45 }
46
47 so will detach at the current time and split all faces in membraneFaces
48 into masterPatch and slavePatch (which have to be present but of 0 size)
49
50\*---------------------------------------------------------------------------*/
51
52#include "argList.H"
53#include "polyMesh.H"
54#include "Time.H"
55#include "polyTopoChange.H"
56#include "mapPolyMesh.H"
57#include "faceSet.H"
58#include "attachDetach.H"
60#include "regionSide.H"
61#include "primitivePatch.H"
62#include "processorMeshes.H"
63
64using namespace Foam;
65
66// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67
68// Find edge between points v0 and v1.
69label findEdge(const primitiveMesh& mesh, const label v0, const label v1)
70{
71 const labelList& pEdges = mesh.pointEdges()[v0];
72
73 forAll(pEdges, pEdgeI)
74 {
75 label edgeI = pEdges[pEdgeI];
76
77 const edge& e = mesh.edges()[edgeI];
78
79 if (e.otherVertex(v0) == v1)
80 {
81 return edgeI;
82 }
83 }
84
86 << "Cannot find edge between mesh points " << v0 << " and " << v1
87 << abort(FatalError);
88
89 return -1;
90}
91
92
93// Checks whether patch present
94void checkPatch(const polyBoundaryMesh& bMesh, const word& name)
95{
96 const label patchi = bMesh.findPatchID(name);
97
98 if (patchi == -1)
99 {
101 << "Cannot find patch " << name << nl
102 << "It should be present but of zero size" << endl
103 << "Valid patches are " << bMesh.names()
104 << exit(FatalError);
105 }
106
107 if (bMesh[patchi].size())
108 {
110 << "Patch " << name << " is present but non-zero size"
111 << exit(FatalError);
112 }
113}
114
115
116
117int main(int argc, char *argv[])
118{
119 argList::addNote
120 (
121 "Splits mesh by making internal faces external at defined faceSet"
122 );
123
124 argList::noParallel();
125 argList::noFunctionObjects(); // Never use function objects
126
127 #include "addOverwriteOption.H"
128
129 argList::addArgument("faceSet", "The faces used for splitting");
130 argList::addArgument("master", "The master patch name");
131 argList::addArgument("slave", "The slave patch name");
132
133 #include "setRootCase.H"
134 #include "createTime.H"
135 #include "createPolyMesh.H"
136
137 const word oldInstance = mesh.pointsInstance();
138
139 const word setName = args[1];
140 const word masterPatch = args[2];
141 const word slavePatch = args[3];
142 const bool overwrite = args.found("overwrite");
143
144 // List of faces to split
145 faceSet facesSet(mesh, setName);
146
147 Info<< "Read " << facesSet.size() << " faces to split" << endl << endl;
148
149
150 // Convert into labelList and check
151
152 labelList faces(facesSet.toc());
153
154 forAll(faces, i)
155 {
156 if (!mesh.isInternalFace(faces[i]))
157 {
159 << "Face " << faces[i] << " in faceSet " << setName
160 << " is not an internal face."
161 << exit(FatalError);
162 }
163 }
164
165
166 // Check for empty master and slave patches
167 checkPatch(mesh.boundaryMesh(), masterPatch);
168 checkPatch(mesh.boundaryMesh(), slavePatch);
169
170
171 //
172 // Find 'side' of all faces on splitregion. Uses regionSide which needs
173 // set of edges on side of this region. Use PrimitivePatch to find these.
174 //
175
176 // Addressing on faces only in mesh vertices.
177 primitiveFacePatch fPatch
178 (
180 (
182 (
183 mesh.faces(),
184 faces
185 )
186 ),
187 mesh.points()
188 );
189
190 const labelList& meshPoints = fPatch.meshPoints();
191
192 // Mark all fence edges : edges on boundary of fPatch but not on boundary
193 // of polyMesh
194 labelHashSet fenceEdges(fPatch.size());
195
196 const labelListList& allEdgeFaces = fPatch.edgeFaces();
197
198 forAll(allEdgeFaces, patchEdgeI)
199 {
200 if (allEdgeFaces[patchEdgeI].size() == 1)
201 {
202 const edge& e = fPatch.edges()[patchEdgeI];
203
204 label edgeI =
206 (
207 mesh,
208 meshPoints[e.start()],
209 meshPoints[e.end()]
210 );
211
212 fenceEdges.insert(edgeI);
213 }
214 }
215
216 // Find sides reachable from 0th face of faceSet
217 label startFacei = faces[0];
218
219 regionSide regionInfo
220 (
221 mesh,
222 facesSet,
223 fenceEdges,
224 mesh.faceOwner()[startFacei],
225 startFacei
226 );
227
228 // Determine flip state for all faces in faceSet
229 boolList zoneFlip(faces.size());
230
231 forAll(faces, i)
232 {
233 zoneFlip[i] = !regionInfo.sideOwner().found(faces[i]);
234 }
235
236
237 // Create and add face zones and mesh modifiers
238 List<pointZone*> pz(0);
239 List<faceZone*> fz(1);
240 List<cellZone*> cz(0);
241
242 fz[0] =
243 new faceZone
244 (
245 "membraneFaces",
246 std::move(faces),
247 std::move(zoneFlip),
248 0,
249 mesh.faceZones()
250 );
251
252 Info<< "Adding point and face zones" << endl;
253 mesh.addZones(pz, fz, cz);
254
255 attachPolyTopoChanger splitter(mesh);
256 splitter.setSize(1);
257
258 // Add the sliding interface mesh modifier to start working at current
259 // time
260 splitter.set
261 (
262 0,
263 new attachDetach
264 (
265 "Splitter",
266 0,
267 splitter,
268 "membraneFaces",
269 masterPatch,
270 slavePatch,
271 scalarField(1, runTime.value())
272 )
273 );
274
275 Info<< nl << "Constructed topologyModifier:" << endl;
276 splitter[0].writeDict(Info);
277
278 if (!overwrite)
279 {
280 ++runTime;
281 }
282
283 splitter.attach();
284
285 if (overwrite)
286 {
287 mesh.setInstance(oldInstance);
288 }
289 else
290 {
291 mesh.setInstance(runTime.timeName());
292 }
293
294 Info<< "Writing mesh to " << runTime.timeName() << endl;
295 if (!mesh.write())
296 {
298 << "Failed writing polyMesh."
299 << exit(FatalError);
300 }
301 topoSet::removeFiles(mesh);
302 processorMeshes::removeFiles(mesh);
303
304
305 Info<< nl << "End" << nl << endl;
306
307 return 0;
308}
309
310
311// ************************************************************************* //
A list of faces which address into the list of points.
A List with indirect addressing. Like IndirectList but does not store addressing.
Definition: IndirectList.H:79
bool found(const word &optName) const
Return true if the named option is found.
Definition: argListI.H:178
Attach/detach boundary mesh modifier. This modifier takes a set of internal faces and converts them i...
Definition: attachDetach.H:65
This class is derived from polyMesh and serves as a tool for statically connecting pieces of a mesh b...
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
A list of face labels.
Definition: faceSet.H:54
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:67
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:79
Determines the 'side' for every face and connected to a singly-connected (through edges) region of fa...
Definition: regionSide.H:64
A class for handling words, derived from Foam::string.
Definition: word.H:68
dynamicFvMesh & mesh
engineTime & runTime
Required Variables.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
label findEdge(const edgeList &edges, const labelList &candidates, const label v0, const label v1)
Return edge among candidates that uses the two vertices.
Definition: meshTools.C:359
Namespace for OpenFOAM.
void checkPatch(const bool allGeometry, const word &name, const PatchType &pp, pointSet &points)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManip< error > abort(error &err)
Definition: errorManip.H:144
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Foam::argList args(argc, argv)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333