CV2DIO.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2013-2018 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 "CV2D.H"
27 #include "OFstream.H"
28 
29 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
30 
31 void Foam::CV2D::writePoints(const fileName& fName, bool internalOnly) const
32 {
33  Info<< "Writing points to " << fName << nl << endl;
34  OFstream str(fName);
35 
36  for
37  (
38  Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
39  vit != finite_vertices_end();
40  ++vit
41  )
42  {
43  if (!internalOnly || vit->internalOrBoundaryPoint())
44  {
45  meshTools::writeOBJ(str, toPoint3D(vit->point()));
46  }
47  }
48 }
49 
50 
51 void Foam::CV2D::writeTriangles(const fileName& fName, bool internalOnly) const
52 {
53  Info<< "Writing triangles to " << fName << nl << endl;
54  OFstream str(fName);
55 
56  labelList vertexMap(number_of_vertices(), -2);
57  label verti = 0;
58 
59  for
60  (
61  Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
62  vit != finite_vertices_end();
63  ++vit
64  )
65  {
66  if (!internalOnly || !vit->farPoint())
67  {
68  vertexMap[vit->index()] = verti++;
69  meshTools::writeOBJ(str, toPoint3D(vit->point()));
70  }
71  }
72 
73  for
74  (
75  Triangulation::Finite_faces_iterator fit = finite_faces_begin();
76  fit != finite_faces_end();
77  ++fit
78  )
79  {
80  if
81  (
82  !internalOnly
83  || (
84  fit->vertex(0)->internalOrBoundaryPoint()
85  || fit->vertex(1)->internalOrBoundaryPoint()
86  || fit->vertex(2)->internalOrBoundaryPoint()
87  )
88  )
89  {
90  str << "f";
91  for (label i = 0; i < 3; ++i)
92  {
93  str << " " << vertexMap[fit->vertex(i)->index()] + 1;
94  }
95  str << nl;
96  }
97  }
98 }
99 
100 
101 void Foam::CV2D::writeFaces(const fileName& fName, bool internalOnly) const
102 {
103  Info<< "Writing dual faces to " << fName << nl << endl;
104  OFstream str(fName);
105 
106  label dualVerti = 0;
107 
108  for
109  (
110  Triangulation::Finite_faces_iterator fit = finite_faces_begin();
111  fit != finite_faces_end();
112  ++fit
113  )
114  {
115  if
116  (
117  !internalOnly
118  || (
119  fit->vertex(0)->internalOrBoundaryPoint()
120  || fit->vertex(1)->internalOrBoundaryPoint()
121  || fit->vertex(2)->internalOrBoundaryPoint()
122  )
123  )
124  {
125  fit->faceIndex() = dualVerti++;
126  meshTools::writeOBJ(str, toPoint3D(circumcenter(fit)));
127  }
128  else
129  {
130  fit->faceIndex() = -1;
131  }
132  }
133 
134  for
135  (
136  Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
137  vit != finite_vertices_end();
138  ++vit
139  )
140  {
141  if (!internalOnly || vit->internalOrBoundaryPoint())
142  {
143  Face_circulator fcStart = incident_faces(vit);
144  Face_circulator fc = fcStart;
145 
146  str<< 'f';
147 
148  do
149  {
150  if (!is_infinite(fc))
151  {
152  if (fc->faceIndex() < 0)
153  {
155  << "Dual face uses vertex defined by a triangle"
156  " defined by an external point"
157  << exit(FatalError);
158  }
159 
160  str<< ' ' << fc->faceIndex() + 1;
161  }
162  } while (++fc != fcStart);
163 
164  str<< nl;
165  }
166  }
167 }
168 
169 
171 (
172  wordList& patchNames,
173  labelList& patchSizes,
174  EdgeMap<label>& mapEdgesRegion,
175  EdgeMap<label>& indirectPatchEdge
176 ) const
177 {
178  label nPatches = qSurf_.patchNames().size() + 1;
179  label defaultPatchIndex = qSurf_.patchNames().size();
180 
181  patchNames.setSize(nPatches);
182  patchSizes.setSize(nPatches, 0);
183  mapEdgesRegion.clear();
184 
185  const wordList& existingPatches = qSurf_.patchNames();
186 
187  forAll(existingPatches, sP)
188  {
189  patchNames[sP] = existingPatches[sP];
190  }
191 
192  patchNames[defaultPatchIndex] = "CV2D_default_patch";
193 
194  for
195  (
196  Triangulation::Finite_edges_iterator eit = finite_edges_begin();
197  eit != finite_edges_end();
198  ++eit
199  )
200  {
201  Face_handle fOwner = eit->first;
202  Face_handle fNeighbor = fOwner->neighbor(eit->second);
203 
204  Vertex_handle vA = fOwner->vertex(cw(eit->second));
205  Vertex_handle vB = fOwner->vertex(ccw(eit->second));
206 
207  if
208  (
209  (vA->internalOrBoundaryPoint() && !vB->internalOrBoundaryPoint())
210  || (vB->internalOrBoundaryPoint() && !vA->internalOrBoundaryPoint())
211  )
212  {
213  point ptA = toPoint3D(vA->point());
214  point ptB = toPoint3D(vB->point());
215 
216  label patchIndex = qSurf_.findPatch(ptA, ptB);
217 
218  if (patchIndex == -1)
219  {
220  patchIndex = defaultPatchIndex;
221 
223  << "Dual face found that is not on a surface "
224  << "patch. Adding to CV2D_default_patch."
225  << endl;
226  }
227 
228  edge e(fOwner->faceIndex(), fNeighbor->faceIndex());
229  patchSizes[patchIndex]++;
230  mapEdgesRegion.insert(e, patchIndex);
231 
232  if (!pointPair(*vA, *vB))
233  {
234  indirectPatchEdge.insert(e, 1);
235  }
236  }
237  }
238 }
239 
240 
242 (
243  point2DField& dualPoints,
244  faceList& dualFaces,
245  wordList& patchNames,
246  labelList& patchSizes,
247  EdgeMap<label>& mapEdgesRegion,
248  EdgeMap<label>& indirectPatchEdge
249 ) const
250 {
251  // Dual points stored in triangle order.
252  dualPoints.setSize(number_of_faces());
253  label dualVerti = 0;
254 
255  for
256  (
257  Triangulation::Finite_faces_iterator fit = finite_faces_begin();
258  fit != finite_faces_end();
259  ++fit
260  )
261  {
262  if
263  (
264  fit->vertex(0)->internalOrBoundaryPoint()
265  || fit->vertex(1)->internalOrBoundaryPoint()
266  || fit->vertex(2)->internalOrBoundaryPoint()
267  )
268  {
269  fit->faceIndex() = dualVerti;
270 
271  dualPoints[dualVerti++] = toPoint2D(circumcenter(fit));
272  }
273  else
274  {
275  fit->faceIndex() = -1;
276  }
277  }
278 
279  dualPoints.setSize(dualVerti);
280 
281  extractPatches(patchNames, patchSizes, mapEdgesRegion, indirectPatchEdge);
282 
283  forAll(patchNames, patchi)
284  {
285  Info<< "Patch " << patchNames[patchi]
286  << " has size " << patchSizes[patchi] << endl;
287  }
288 
289  // Create dual faces
290  // ~~~~~~~~~~~~~~~~~
291 
292  dualFaces.setSize(number_of_vertices());
293  label dualFacei = 0;
294  labelList faceVerts(maxNvert);
295 
296  for
297  (
298  Triangulation::Finite_vertices_iterator vit = finite_vertices_begin();
299  vit != finite_vertices_end();
300  ++vit
301  )
302  {
303  if (vit->internalOrBoundaryPoint())
304  {
305  Face_circulator fcStart = incident_faces(vit);
306  Face_circulator fc = fcStart;
307  label verti = 0;
308 
309  do
310  {
311  if (!is_infinite(fc))
312  {
313  if (fc->faceIndex() < 0)
314  {
316  << "Dual face uses vertex defined by a triangle"
317  " defined by an external point"
318  << exit(FatalError);
319  }
320 
321  // Look up the index of the triangle
322  faceVerts[verti++] = fc->faceIndex();
323  }
324  } while (++fc != fcStart);
325 
326  if (faceVerts.size() > 2)
327  {
328  dualFaces[dualFacei++] =
329  face(labelList::subList(faceVerts, verti));
330  }
331  else
332  {
333  Info<< "From triangle point:" << vit->index()
334  << " coord:" << toPoint2D(vit->point())
335  << " generated illegal dualFace:" << faceVerts
336  << endl;
337  }
338  }
339  }
340 
341  dualFaces.setSize(dualFacei);
342 }
343 
344 
345 void Foam::CV2D::writePatch(const fileName& fName) const
346 {
347  point2DField dual2DPoints;
348  faceList dualFaces;
350  labelList patchSizes;
351  EdgeMap<label> mapEdgesRegion;
352  EdgeMap<label> indirectPatchEdge;
353 
354  calcDual
355  (
356  dual2DPoints,
357  dualFaces,
358  patchNames,
359  patchSizes,
360  mapEdgesRegion,
361  indirectPatchEdge
362  );
363 
364  pointField dualPoints(dual2DPoints.size());
365  forAll(dualPoints, ip)
366  {
367  dualPoints[ip] = toPoint3D(dual2DPoints[ip]);
368  }
369 
370  // Dump as primitive patch to be read by extrudeMesh.
371  OFstream str(fName);
372 
373  Info<< "Writing patch to be used with extrudeMesh to file " << fName
374  << endl;
375 
376  str << dualPoints << nl << dualFaces << nl;
377 }
378 
379 
380 // ************************************************************************* //
#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
void writePatch(const fileName &fName) const
Write patch.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
List< face > faceList
Definition: faceListFwd.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
point toPoint3D(const point2D &) const
Definition: CV2DI.H:141
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
SubList< label > subList
Declare type of subList.
Definition: List.H:192
bool pointPair(const indexedVertex< Gt, Vb > &v0, const indexedVertex< Gt, Vb > &v1)
const point2D & toPoint2D(const point &) const
Definition: CV2DI.H:124
void writeFaces(const fileName &fName, bool internalOnly) const
Write dual faces as .obj file.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
wordList patchNames(nPatches)
List< label > labelList
A List of labels.
Definition: labelList.H:56
void calcDual(point2DField &dualPoints, faceList &dualFaces, wordList &patchNames, labelList &patchSizes, EdgeMap< label > &mapEdgesRegion, EdgeMap< label > &indirectPatchEdge) const
Calculates dual points (circumcentres of tets) and faces.
static const char nl
Definition: Ostream.H:265
const List< word > & patchNames() const
Return the patch names.
vector2DField point2DField
point2DField is a vector2DField.
void extractPatches(wordList &patchNames, labelList &patchSizes, EdgeMap< label > &mapEdgesRegion, EdgeMap< label > &indirectPatchEdge) const
Extract patch names and sizes.
List< word > wordList
A List of words.
Definition: fileName.H:54
void writeTriangles(const fileName &fName, bool internalOnly) const
Write triangles as .obj file.
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
void writePoints(const fileName &fName, bool internalOnly) const
Write internal points to .obj file.
label findPatch(const point &ptA, const point &ptB) const
Find which patch is intersected by the line from one point to.