sampledCuttingPlane.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2016-2022 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
27Class
28 Foam::sampledCuttingPlane
29
30Description
31 A sampledSurface defined by a plane using an iso-surface algorithm
32 to \a cut the mesh.
33
34 This is often embedded as part of a sampled surfaces function object.
35
36Usage
37 Example of function object partial specification:
38 \verbatim
39 surfaces
40 {
41 surface1
42 {
43 type cuttingPlane;
44 point ...;
45 normal ...;
46 }
47 }
48 \endverbatim
49
50 Where the sub-entries comprise:
51 \table
52 Property | Description | Required | Default
53 type | cuttingPlane | yes |
54 planeType | Plane description (pointAndNormal etc) | no |
55 offsets | Offsets of the origin in the normal direction | no | (0)
56 isoMethod | Iso-algorithm (cell/topo/point) | no | topo
57 bounds | limit with bounding box | no |
58 zone | limit to cell zone (name or regex) | no |
59 zones | limit to cell zones (names, regexs) | no |
60 exposedPatchName | name for zone subset | optional |
61 regularise | Face simplification (enum or bool) | no | true
62 mergeTol | tolerance for merging points | no | 1e-6
63 coordinateSystem | define plane within given cartesian system | no |
64 transform | define plane within given cartesian system | no |
65 \endtable
66
67Note
68 The keyword \c zones has priority over \c zone.
69
70SeeAlso
71 Foam::plane
72
73SourceFiles
74 sampledCuttingPlane.C
75
76\*---------------------------------------------------------------------------*/
77
78#ifndef Foam_sampledCuttingPlane_H
79#define Foam_sampledCuttingPlane_H
80
81#include "sampledSurface.H"
82#include "plane.H"
83#include "fvMeshSubset.H"
84#include "isoSurfaceBase.H"
85
86// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
87
88namespace Foam
89{
90
91/*---------------------------------------------------------------------------*\
92 Class sampledCuttingPlane Declaration
93\*---------------------------------------------------------------------------*/
94
95class sampledCuttingPlane
96:
97 public sampledSurface
98{
99 // Private Data
100
101 //- Plane
102 const plane plane_;
103
104 //- The offsets to the plane - defaults to (0).
105 List<scalar> offsets_;
106
107 //- Parameters (filtering etc) for iso-surface
108 isoSurfaceParams isoParams_;
109
110 //- Whether to recalculate cell values as average of point values
111 bool average_;
112
113 //- Use simple sub-meshing in algorithm itself
114 bool simpleSubMesh_;
115
116 //- The zone or zones in which cutting is to occur
117 wordRes zoneNames_;
118
119 //- For zones: patch to put exposed faces into
120 mutable word exposedPatchName_;
121
122 //- Track if the surface needs an update
123 mutable bool needsUpdate_;
124
125
126 // Sampling geometry. Directly stored or via an iso-surface (ALGO_POINT)
127
128 //- The extracted surface (direct storage)
129 mutable meshedSurface surface_;
130
131 //- For every face the original cell in mesh (direct storage)
132 mutable labelList meshCells_;
133
134 //- Constructed iso-surface (ALGO_POINT), for interpolators
135 autoPtr<isoSurfaceBase> isoSurfacePtr_;
136
137
138 // Mesh Subsetting
139
140 //- Cached subMesh for (pre-)subset of cell zones
141 mutable autoPtr<fvMeshSubset> subMeshPtr_;
142
143 //- Cached ignore cells for (post-)subset of cell zones
144 mutable autoPtr<bitSet> ignoreCellsPtr_;
145
146
147 // Fields
148
149 //- Distance to cell centres
150 autoPtr<volScalarField> cellDistancePtr_;
151
152 //- Distance to points
153 scalarField pointDistance_;
154
155
156 // Private Member Functions
157
158 //- Define plane from dictionary entry and apply any coordinate
159 //- transformations
160 static plane definePlane(const polyMesh& mesh, const dictionary& dict);
161
162 //- Check and warn if bounding box does not intersect mesh or plane
163 void checkBoundsIntersection
164 (
165 const plane& pln,
166 const boundBox& meshBb
167 ) const;
168
169 //- Fill cellDistance, pointDistance fields for the specified plane
170 void setDistanceFields(const plane& pln);
171
172 //- Collect iso-surfaces into a single surface (No point merging)
173 void combineSurfaces(PtrList<isoSurfaceBase>& isoSurfPtrs);
174
175 //- Create iso surface
176 void createGeometry();
177
178 //- Sample volume field onto surface faces
179 template<class Type>
180 tmp<Field<Type>> sampleOnFaces
181 (
182 const interpolation<Type>& sampler
183 ) const;
184
185 //- Interpolate volume field onto surface points
186 template<class Type>
187 tmp<Field<Type>> sampleOnPoints
188 (
189 const interpolation<Type>& interpolator
190 ) const;
191
192 //- Use isoSurfacePtr_ for point interpolation
193 template<class Type>
194 tmp<Field<Type>> sampleOnIsoSurfacePoints
195 (
196 const interpolation<Type>& interpolator
197 ) const;
198
199
200protected:
201
202 // Protected Member Functions
203
204 //- Is currently backed by an isoSurfacePtr_
205 bool hasIsoSurface() const
206 {
207 return bool(isoSurfacePtr_);
208 }
209
210public:
211
212 //- Runtime type information
213 TypeName("sampledCuttingPlane");
214
215
216 // Constructors
217
218 //- Construct from dictionary
220 (
221 const word& name,
222 const polyMesh& mesh,
223 const dictionary& dict
224 );
225
226
227 //- Destructor
228 virtual ~sampledCuttingPlane() = default;
229
230
231 // Member Functions
232
233 //- Does the surface need an update?
234 virtual bool needsUpdate() const;
235
236 //- Mark the surface as needing an update.
237 // May also free up unneeded data.
238 // Return false if surface was already marked as expired.
239 virtual bool expire();
240
241 //- Update the surface as required.
242 // Do nothing (and return false) if no update was needed
243 virtual bool update();
244
245 //- The current surface geometry
246 const meshedSurface& surface() const
247 {
248 if (isoSurfacePtr_)
249 {
250 return *isoSurfacePtr_;
251 }
252 return surface_;
253 }
254
255 //- For each face, the original cell in mesh
256 const labelList& meshCells() const
257 {
258 if (isoSurfacePtr_)
259 {
260 return isoSurfacePtr_->meshCells();
261 }
262 return meshCells_;
263 }
264
265 //- Points of surface
266 virtual const pointField& points() const
267 {
268 return surface().points();
270
271 //- Faces of surface
272 virtual const faceList& faces() const
273 {
274 return surface().surfFaces();
275 }
276
277 //- Per-face zone/region information
278 virtual const labelList& zoneIds() const
279 {
280 return labelList::null();
281 }
282
283 //- Face area magnitudes
284 virtual const vectorField& Sf() const
285 {
286 return surface().Sf();
287 }
288
289 //- Face area magnitudes
290 virtual const scalarField& magSf() const
291 {
292 return surface().magSf();
293 }
294
295 //- Face centres
296 virtual const vectorField& Cf() const
297 {
298 return surface().Cf();
299 }
300
301
302 // Sample
303
304 //- Sample volume field onto surface faces
305 virtual tmp<scalarField> sample
306 (
307 const interpolation<scalar>& sampler
308 ) const;
309
310 //- Sample volume field onto surface faces
312 (
313 const interpolation<vector>& sampler
314 ) const;
315
316 //- Sample volume field onto surface faces
318 (
319 const interpolation<sphericalTensor>& sampler
320 ) const;
321
322 //- Sample volume field onto surface faces
324 (
325 const interpolation<symmTensor>& sampler
326 ) const;
327
328 //- Sample volume field onto surface faces
331 const interpolation<tensor>& sampler
332 ) const;
333
334
335 // Interpolate
337 //- Interpolate volume field onto surface points
339 (
340 const interpolation<scalar>& interpolator
341 ) const;
343 //- Interpolate volume field onto surface points
345 (
346 const interpolation<vector>& interpolator
347 ) const;
349 //- Interpolate volume field onto surface points
351 (
352 const interpolation<sphericalTensor>& interpolator
353 ) const;
355 //- Interpolate volume field onto surface points
357 (
358 const interpolation<symmTensor>& interpolator
359 ) const;
361 //- Interpolate volume field onto surface points
363 (
364 const interpolation<tensor>& interpolator
365 ) const;
366
367
368 // Output
369
370 //- Print information
371 virtual void print(Ostream& os, int level=0) const;
372};
373
374
375// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
376
377} // End namespace Foam
378
379// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
380
381#ifdef NoRepository
383#endif
384
385// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
386
387#endif
388
389// ************************************************************************* //
Minimal example by using system/controlDict.functions:
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
static const List< label > & null()
Return a null List.
Definition: ListI.H:109
const vectorField & Sf() const
Face area vectors (normals)
const scalarField & magSf() const
Face area magnitudes.
const vectorField & Cf() const
Face centres.
const List< Face > & surfFaces() const
Return const access to the faces.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
const Field< point_type > & points() const noexcept
Return reference to global points.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:64
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base class for volume field interpolation.
Definition: interpolation.H:60
Preferences for controlling iso-surface algorithms.
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
Definition: plane.H:95
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
A sampledSurface defined by a plane using an iso-surface algorithm to cut the mesh.
virtual const pointField & points() const
Points of surface.
virtual ~sampledCuttingPlane()=default
Destructor.
virtual void print(Ostream &os, int level=0) const
Print information.
const meshedSurface & surface() const
The current surface geometry.
virtual const faceList & faces() const
Faces of surface.
virtual const vectorField & Cf() const
Face centres.
virtual const labelList & zoneIds() const
Per-face zone/region information.
virtual const scalarField & magSf() const
Face area magnitudes.
virtual bool expire()
Mark the surface as needing an update.
virtual bool needsUpdate() const
Does the surface need an update?
virtual bool update()
Update the surface as required.
virtual const vectorField & Sf() const
Face area magnitudes.
const labelList & meshCells() const
For each face, the original cell in mesh.
TypeName("sampledCuttingPlane")
Runtime type information.
bool hasIsoSurface() const
Is currently backed by an isoSurfacePtr_.
An abstract class for surfaces with sampling.
const word & name() const noexcept
Name of surface.
const polyMesh & mesh() const noexcept
Access to the underlying mesh.
bool interpolate() const noexcept
Same as isPointData()
A class for managing temporary objects.
Definition: tmp.H:65
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
A class for handling words, derived from Foam::string.
Definition: word.H:68
bool
Definition: EEqn.H:20
OBJstream os(runTime.globalPath()/outputName)
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: List.H:66
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
MeshedSurface< face > meshedSurface
List< treeBoundBox > meshBb(1, treeBoundBox(boundBox(coarseMesh.points(), false)).extend(rndGen, 1e-3))
dictionary dict
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73