Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
sampling
sampledSurface
isoSurface
isoSurface.H
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) 1991-2010 OpenCFD Ltd.
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
Class
25
Foam::isoSurface
26
27
Description
28
A surface formed by the iso value.
29
After "Regularised Marching Tetrahedra: improved iso-surface extraction",
30
G.M. Treece, R.W. Prager and A.H. Gee.
31
32
Note:
33
- does tets without using cell centres/cell values. Not tested.
34
- regularisation can create duplicate triangles/non-manifold surfaces.
35
Current handling of those is bit ad-hoc for now and not perfect.
36
- regularisation does not do boundary points so as to preserve the
37
boundary perfectly.
38
- uses geometric merge with fraction of bounding box as distance.
39
- triangles can be between two cell centres so constant sampling
40
does not make sense.
41
- on empty patches behaves like zero gradient.
42
- does not do 2D correctly, creates non-flat iso surface.
43
- on processor boundaries might two overlapping (identical) triangles
44
(one from either side)
45
46
The handling on coupled patches is a bit complex. All fields
47
(values and coordinates) get rewritten so
48
- empty patches get zerogradient (value) and facecentre (coordinates)
49
- separated processor patch faces get interpolate (value) and facecentre
50
(coordinates). (this is already the default for cyclics)
51
- non-separated processor patch faces keep opposite value and cell centre
52
53
Now the triangle generation on non-separated processor patch faces
54
can use the neighbouring value. Any separated processor face or cyclic
55
face gets handled just like any boundary face.
56
57
58
SourceFiles
59
isoSurface.C
60
61
\*---------------------------------------------------------------------------*/
62
63
#ifndef isoSurface_H
64
#define isoSurface_H
65
66
#include <
triSurface/triSurface.H
>
67
#include <
OpenFOAM/labelPair.H
>
68
#include <
meshTools/pointIndexHit.H
>
69
#include <
OpenFOAM/PackedBoolList.H
>
70
#include <
finiteVolume/volFields.H
>
71
#include <
finiteVolume/slicedVolFields.H
>
72
73
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
74
75
namespace
Foam
76
{
77
78
class
fvMesh;
79
80
/*---------------------------------------------------------------------------*\
81
Class isoSurface Declaration
82
\*---------------------------------------------------------------------------*/
83
84
class
isoSurface
85
:
86
public
triSurface
87
{
88
// Private data
89
90
enum
segmentCutType
91
{
92
NEARFIRST,
// intersection close to e.first()
93
NEARSECOND,
// ,, e.second()
94
NOTNEAR
// no intersection
95
};
96
97
enum
cellCutType
98
{
99
NOTCUT,
// not cut
100
SPHERE,
// all edges to cell centre cut
101
CUT
// normal cut
102
};
103
104
105
//- Reference to mesh
106
const
fvMesh
& mesh_;
107
108
const
scalarField
& pVals_;
109
110
//- Input volScalarField with separated coupled patches rewritten
111
autoPtr<slicedVolScalarField>
cValsPtr_;
112
113
//- Isosurface value
114
const
scalar iso_;
115
116
//- Regularise?
117
const
Switch
regularise_;
118
119
//- When to merge points
120
const
scalar mergeDistance_;
121
122
123
//- Whether face might be cut
124
List<cellCutType>
faceCutType_;
125
126
//- Whether cell might be cut
127
List<cellCutType>
cellCutType_;
128
129
//- Estimated number of cut cells
130
label nCutCells_;
131
132
//- For every triangle the original cell in mesh
133
labelList
meshCells_;
134
135
//- For every unmerged triangle point the point in the triSurface
136
labelList
triPointMergeMap_;
137
138
139
// Private Member Functions
140
141
// Point synchronisation
142
143
//- Does tensor differ (to within mergeTolerance) from identity
144
bool
noTransform(
const
tensor
& tt)
const
;
145
146
//- Is patch a collocated (i.e. non-separated) coupled patch?
147
static
bool
collocatedPatch(
const
polyPatch
&);
148
149
//- Per face whether is collocated
150
PackedBoolList
collocatedFaces(
const
coupledPolyPatch
&)
const
;
151
152
//- Take value at local point patchPointI and assign it to its
153
// correct place in patchValues (for transferral) and sharedValues
154
// (for reduction)
155
void
insertPointData
156
(
157
const
processorPolyPatch
& pp,
158
const
Map<label>
& meshToShared,
159
const
pointField
& pointValues,
160
const
label patchPointI,
161
pointField
& patchValues,
162
pointField
& sharedValues
163
)
const
;
164
165
//- Synchonise points on all non-separated coupled patches
166
void
syncUnseparatedPoints
167
(
168
pointField
& collapsedPoint,
169
const
point
& nullValue
170
)
const
;
171
172
173
//- Get location of iso value as fraction inbetween s0,s1
174
scalar isoFraction
175
(
176
const
scalar s0,
177
const
scalar s1
178
)
const
;
179
180
//- Check if any edge of a face is cut
181
bool
isEdgeOfFaceCut
182
(
183
const
scalarField
& pVals,
184
const
face
&
f
,
185
const
bool
ownLower,
186
const
bool
neiLower
187
)
const
;
188
189
void
getNeighbour
190
(
191
const
labelList
&
boundaryRegion
,
192
const
volVectorField
& meshC,
193
const
volScalarField
& cVals,
194
const
label cellI,
195
const
label faceI,
196
scalar& nbrValue,
197
point
& nbrPoint
198
)
const
;
199
200
//- Set faceCutType,cellCutType.
201
void
calcCutTypes
202
(
203
const
labelList
& boundaryRegion,
204
const
volVectorField
& meshC,
205
const
volScalarField
& cVals,
206
const
scalarField
& pVals
207
);
208
209
static
labelPair
findCommonPoints
210
(
211
const
labelledTri
&,
212
const
labelledTri
&
213
);
214
215
static
point
calcCentre(
const
triSurface
&);
216
217
static
pointIndexHit
collapseSurface
218
(
219
pointField
&
localPoints
,
220
DynamicList<labelledTri, 64>
& localTris
221
);
222
223
//- Determine per cc whether all near cuts can be snapped to single
224
// point.
225
void
calcSnappedCc
226
(
227
const
labelList
& boundaryRegion,
228
const
volVectorField
& meshC,
229
const
volScalarField
& cVals,
230
const
scalarField
& pVals,
231
DynamicList<point>
& snappedPoints,
232
labelList
& snappedCc
233
)
const
;
234
235
//- Determine per point whether all near cuts can be snapped to single
236
// point.
237
void
calcSnappedPoint
238
(
239
const
PackedBoolList
& isBoundaryPoint,
240
const
labelList
& boundaryRegion,
241
const
volVectorField
& meshC,
242
const
volScalarField
& cVals,
243
const
scalarField
& pVals,
244
DynamicList<point>
& snappedPoints,
245
labelList
& snappedPoint
246
)
const
;
247
248
249
//- Return input field with coupled (and empty) patch values rewritten
250
template
<
class
Type>
251
tmp
<
SlicedGeometricField
252
<Type,
fvPatchField
,
slicedFvPatchField
,
volMesh
> >
253
adaptPatchFields
254
(
255
const
GeometricField<Type, fvPatchField, volMesh>
& fld
256
)
const
;
257
258
//- Generate single point by interpolation or snapping
259
template
<
class
Type>
260
Type generatePoint
261
(
262
const
scalar s0,
263
const
Type& p0,
264
const
bool
hasSnap0,
265
const
Type& snapP0,
266
267
const
scalar s1,
268
const
Type& p1,
269
const
bool
hasSnap1,
270
const
Type& snapP1
271
)
const
;
272
273
template
<
class
Type>
274
void
generateTriPoints
275
(
276
const
scalar s0,
277
const
Type& p0,
278
const
bool
hasSnap0,
279
const
Type& snapP0,
280
281
const
scalar s1,
282
const
Type& p1,
283
const
bool
hasSnap1,
284
const
Type& snapP1,
285
286
const
scalar s2,
287
const
Type& p2,
288
const
bool
hasSnap2,
289
const
Type& snapP2,
290
291
const
scalar s3,
292
const
Type& p3,
293
const
bool
hasSnap3,
294
const
Type& snapP3,
295
296
DynamicList<Type>
&
points
297
)
const
;
298
299
template
<
class
Type>
300
label generateFaceTriPoints
301
(
302
const
volScalarField
& cVals,
303
const
scalarField
& pVals,
304
305
const
GeometricField<Type, fvPatchField, volMesh>
& cCoords,
306
const
Field<Type>
& pCoords,
307
308
const
DynamicList<Type>
& snappedPoints,
309
const
labelList
& snappedCc,
310
const
labelList
& snappedPoint,
311
const
label faceI,
312
313
const
scalar neiVal,
314
const
Type& neiPt,
315
const
bool
hasNeiSnap,
316
const
Type& neiSnapPt,
317
318
DynamicList<Type>
& triPoints,
319
DynamicList<label>
& triMeshCells
320
)
const
;
321
322
template
<
class
Type>
323
void
generateTriPoints
324
(
325
const
volScalarField
& cVals,
326
const
scalarField
& pVals,
327
328
const
GeometricField<Type, fvPatchField, volMesh>
& cCoords,
329
const
Field<Type>
& pCoords,
330
331
const
DynamicList<Type>
& snappedPoints,
332
const
labelList
& snappedCc,
333
const
labelList
& snappedPoint,
334
335
DynamicList<Type>
& triPoints,
336
DynamicList<label>
& triMeshCells
337
)
const
;
338
339
triSurface
stitchTriPoints
340
(
341
const
bool
checkDuplicates,
342
const
List<point>
& triPoints,
343
labelList
& triPointReverseMap,
// unmerged to merged point
344
labelList
& triMap
// merged to unmerged triangle
345
)
const
;
346
347
//- Check single triangle for (topological) validity
348
static
bool
validTri(
const
triSurface
&,
const
label);
349
350
//- Determine edge-face addressing
351
void
calcAddressing
352
(
353
const
triSurface
& surf,
354
List
<
FixedList<label, 3>
>&
faceEdges
,
355
labelList
& edgeFace0,
356
labelList
& edgeFace1,
357
Map<labelList>
& edgeFacesRest
358
)
const
;
359
360
//- Determine orientation
361
static
void
walkOrientation
362
(
363
const
triSurface
& surf,
364
const
List
<
FixedList<label, 3>
>&
faceEdges
,
365
const
labelList
& edgeFace0,
366
const
labelList
& edgeFace1,
367
const
label seedTriI,
368
labelList
& flipState
369
);
370
371
//- Orient surface
372
static
void
orientSurface
373
(
374
triSurface
&,
375
const
List
<
FixedList<label, 3>
>&
faceEdges
,
376
const
labelList
& edgeFace0,
377
const
labelList
& edgeFace1,
378
const
Map<labelList>
& edgeFacesRest
379
);
380
381
//- Is triangle (given by 3 edges) not fully connected?
382
static
bool
danglingTriangle
383
(
384
const
FixedList<label, 3>
& fEdges,
385
const
labelList
& edgeFace1
386
);
387
388
//- Mark all non-fully connected triangles
389
static
label markDanglingTriangles
390
(
391
const
List
<
FixedList<label, 3>
>&
faceEdges
,
392
const
labelList
& edgeFace0,
393
const
labelList
& edgeFace1,
394
const
Map<labelList>
& edgeFacesRest,
395
boolList
& keepTriangles
396
);
397
398
static
triSurface
subsetMesh
399
(
400
const
triSurface
& s,
401
const
labelList
& newToOldFaces,
402
labelList
& oldToNewPoints,
403
labelList
& newToOldPoints
404
);
405
406
public
:
407
408
//- Runtime type information
409
TypeName
(
"isoSurface"
);
410
411
412
// Constructors
413
414
//- Construct from cell values and point values. Uses boundaryField
415
// for boundary values. Holds reference to cellIsoVals and
416
// pointIsoVals.
417
isoSurface
418
(
419
const
volScalarField
& cellIsoVals,
420
const
scalarField
& pointIsoVals,
421
const
scalar iso,
422
const
bool
regularise,
423
const
scalar mergeTol = 1
E
-6
// fraction of bounding box
424
);
425
426
427
// Member Functions
428
429
//- For every face original cell in mesh
430
const
labelList
&
meshCells
()
const
431
{
432
return
meshCells_;
433
}
434
435
//- For every unmerged triangle point the point in the triSurface
436
const
labelList
&
triPointMergeMap
()
const
437
{
438
return
triPointMergeMap_;
439
}
440
441
//- Interpolates cCoords,pCoords. Uses the references to the original
442
// fields used to create the iso surface.
443
template
<
class
Type>
444
tmp<Field<Type>
>
interpolate
445
(
446
const
GeometricField<Type, fvPatchField, volMesh>
& cCoords,
447
const
Field<Type>
& pCoords
448
)
const
;
449
450
};
451
452
453
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
454
455
}
// End namespace Foam
456
457
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
458
459
#ifdef NoRepository
460
# include "
isoSurfaceTemplates.C
"
461
#endif
462
463
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
464
465
#endif
466
467
// ************************ vim: set sw=4 sts=4 et: ************************ //