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 | Copyright (C) 2011-2016 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 Application
25  splitMesh
26 
27 Description
28  Splits mesh by making internal faces external. Uses attachDetach.
29 
30  Generates a meshModifier of the form:
31 
32  Splitter
33  {
34  type attachDetach;
35  faceZoneName membraneFaces;
36  masterPatchName masterPatch;
37  slavePatchName slavePatch;
38  triggerTimes runTime.value();
39  }
40 
41  so will detach at the current time and split all faces in membraneFaces
42  into masterPatch and slavePatch (which have to be present but of 0 size)
43 
44 \*---------------------------------------------------------------------------*/
45 
46 #include "argList.H"
47 #include "polyMesh.H"
48 #include "Time.H"
49 #include "polyTopoChange.H"
50 #include "mapPolyMesh.H"
51 #include "faceSet.H"
52 #include "attachDetach.H"
53 #include "attachPolyTopoChanger.H"
54 #include "regionSide.H"
55 #include "primitiveFacePatch.H"
56 
57 using namespace Foam;
58 
59 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60 
61 // Find edge between points v0 and v1.
62 label findEdge(const primitiveMesh& mesh, const label v0, const label v1)
63 {
64  const labelList& pEdges = mesh.pointEdges()[v0];
65 
66  forAll(pEdges, pEdgeI)
67  {
68  label edgeI = pEdges[pEdgeI];
69 
70  const edge& e = mesh.edges()[edgeI];
71 
72  if (e.otherVertex(v0) == v1)
73  {
74  return edgeI;
75  }
76  }
77 
79  << "Cannot find edge between mesh points " << v0 << " and " << v1
80  << abort(FatalError);
81 
82  return -1;
83 }
84 
85 
86 // Checks whether patch present
87 void checkPatch(const polyBoundaryMesh& bMesh, const word& name)
88 {
89  const label patchi = bMesh.findPatchID(name);
90 
91  if (patchi == -1)
92  {
94  << "Cannot find patch " << name << nl
95  << "It should be present but of zero size" << endl
96  << "Valid patches are " << bMesh.names()
97  << exit(FatalError);
98  }
99 
100  if (bMesh[patchi].size())
101  {
103  << "Patch " << name << " is present but non-zero size"
104  << exit(FatalError);
105  }
106 }
107 
108 
109 
110 int main(int argc, char *argv[])
111 {
113  #include "addOverwriteOption.H"
114 
115  argList::validArgs.append("faceSet");
116  argList::validArgs.append("masterPatch");
117  argList::validArgs.append("slavePatch");
118 
119  #include "setRootCase.H"
120  #include "createTime.H"
121  runTime.functionObjects().off();
122  #include "createPolyMesh.H"
123  const word oldInstance = mesh.pointsInstance();
124 
125  const word setName = args[1];
126  const word masterPatch = args[2];
127  const word slavePatch = args[3];
128  const bool overwrite = args.optionFound("overwrite");
129 
130  // List of faces to split
131  faceSet facesSet(mesh, setName);
132 
133  Info<< "Read " << facesSet.size() << " faces to split" << endl << endl;
134 
135 
136  // Convert into labelList and check
137 
138  labelList faces(facesSet.toc());
139 
140  forAll(faces, i)
141  {
142  if (!mesh.isInternalFace(faces[i]))
143  {
145  << "Face " << faces[i] << " in faceSet " << setName
146  << " is not an internal face."
147  << exit(FatalError);
148  }
149  }
150 
151 
152  // Check for empty master and slave patches
153  checkPatch(mesh.boundaryMesh(), masterPatch);
154  checkPatch(mesh.boundaryMesh(), slavePatch);
155 
156 
157  //
158  // Find 'side' of all faces on splitregion. Uses regionSide which needs
159  // set of edges on side of this region. Use PrimitivePatch to find these.
160  //
161 
162  // Addressing on faces only in mesh vertices.
163  primitiveFacePatch fPatch
164  (
165  faceList
166  (
168  (
169  mesh.faces(),
170  faces
171  )
172  ),
173  mesh.points()
174  );
175 
176  const labelList& meshPoints = fPatch.meshPoints();
177 
178  // Mark all fence edges : edges on boundary of fPatch but not on boundary
179  // of polyMesh
180  labelHashSet fenceEdges(fPatch.size());
181 
182  const labelListList& allEdgeFaces = fPatch.edgeFaces();
183 
184  forAll(allEdgeFaces, patchEdgeI)
185  {
186  if (allEdgeFaces[patchEdgeI].size() == 1)
187  {
188  const edge& e = fPatch.edges()[patchEdgeI];
189 
190  label edgeI =
191  findEdge
192  (
193  mesh,
194  meshPoints[e.start()],
195  meshPoints[e.end()]
196  );
197 
198  fenceEdges.insert(edgeI);
199  }
200  }
201 
202  // Find sides reachable from 0th face of faceSet
203  label startFacei = faces[0];
204 
205  regionSide regionInfo
206  (
207  mesh,
208  facesSet,
209  fenceEdges,
210  mesh.faceOwner()[startFacei],
211  startFacei
212  );
213 
214  // Determine flip state for all faces in faceSet
215  boolList zoneFlip(faces.size());
216 
217  forAll(faces, i)
218  {
219  zoneFlip[i] = !regionInfo.sideOwner().found(faces[i]);
220  }
221 
222 
223  // Create and add face zones and mesh modifiers
224  List<pointZone*> pz(0);
225  List<faceZone*> fz(1);
226  List<cellZone*> cz(0);
227 
228  fz[0] =
229  new faceZone
230  (
231  "membraneFaces",
232  faces,
233  zoneFlip,
234  0,
235  mesh.faceZones()
236  );
237 
238  Info<< "Adding point and face zones" << endl;
239  mesh.addZones(pz, fz, cz);
240 
241  attachPolyTopoChanger splitter(mesh);
242  splitter.setSize(1);
243 
244  // Add the sliding interface mesh modifier to start working at current
245  // time
246  splitter.set
247  (
248  0,
249  new attachDetach
250  (
251  "Splitter",
252  0,
253  splitter,
254  "membraneFaces",
255  masterPatch,
256  slavePatch,
257  scalarField(1, runTime.value())
258  )
259  );
260 
261  Info<< nl << "Constructed topologyModifier:" << endl;
262  splitter[0].writeDict(Info);
263 
264  if (!overwrite)
265  {
266  runTime++;
267  }
268 
269  splitter.attach();
270 
271  if (overwrite)
272  {
273  mesh.setInstance(oldInstance);
274  }
275  else
276  {
277  mesh.setInstance(runTime.timeName());
278  }
279 
280  Info<< "Writing mesh to " << runTime.timeName() << endl;
281  if (!mesh.write())
282  {
284  << "Failed writing polyMesh."
285  << exit(FatalError);
286  }
287 
288  Info<< "End\n" << endl;
289  return 0;
290 }
291 
292 
293 // ************************************************************************* //
label end() const
Return end vertex label.
Definition: edgeI.H:92
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A list of face labels.
Definition: faceSet.H:48
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
wordList names() const
Return a list of patch names.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual const pointField & points() const =0
Return mesh points.
static void noParallel()
Remove the parallel options.
Definition: argList.C:146
static SLList< string > validArgs
A list of valid (mandatory) arguments.
Definition: argList.H:154
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:337
label otherVertex(const label a) const
Given one vertex, return the other.
Definition: edgeI.H:103
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
A list of faces which address into the list of points.
Determines the &#39;side&#39; for every face and connected to a singly-connected (through edges) region of fa...
Definition: regionSide.H:61
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const labelListList & pointEdges() const
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::polyBoundaryMesh.
This class is derived from polyMesh and serves as a tool for statically connecting pieces of a mesh b...
label start() const
Return start vertex label.
Definition: edgeI.H:81
static const char nl
Definition: Ostream.H:262
Attach/detach boundary mesh modifier. This modifier takes a set of internal faces and converts them i...
Definition: attachDetach.H:59
A List with indirect addressing.
Definition: fvMatrix.H:106
virtual const faceList & faces() const =0
Return faces.
messageStream Info
virtual const labelList & faceOwner() const =0
Face face-owner addresing.
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
Foam::argList args(argc, argv)
label findPatchID(const word &patchName) const
Find patch index given a name.
Namespace for OpenFOAM.