sampledTriSurfaceMesh.H
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) 2011-2019 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 Class
25  Foam::sampledSurfaces::triSurfaceMesh
26 
27 Description
28  A sampledSurface from a triSurfaceMesh. It samples on the points/triangles
29  of the triSurface.
30 
31  - it either samples cells or (non-coupled) boundary faces
32 
33  - 6 different modes:
34  - source=cells, interpolate=false:
35  finds per triangle centre the nearest cell centre and uses its value
36  - source=cells, interpolate=true
37  finds per triangle centre the nearest cell centre.
38  Per surface point checks if this nearest cell is the one containing
39  point; otherwise projects the point onto the nearest point on
40  the boundary of the cell (to make sure interpolateCellPoint
41  gets a valid location)
42 
43  - source=insideCells, interpolate=false:
44  finds per triangle centre the cell containing it and uses its value.
45  Trims triangles outside mesh.
46  - source=insideCells, interpolate=true
47  Per surface point interpolate cell containing it.
48 
49  - source=boundaryFaces, interpolate=false:
50  finds per triangle centre the nearest point on the boundary
51  (uncoupled faces only) and uses the value (or 0 if the nearest
52  is on an empty boundary)
53  - source=boundaryFaces, interpolate=true:
54  finds per triangle centre the nearest point on the boundary
55  (uncoupled faces only).
56  Per surface point projects the point onto this boundary face
57  (to make sure interpolateCellPoint gets a valid location)
58 
59  - since it finds a nearest per triangle each triangle is guaranteed
60  to be on one processor only. So after stitching (by sampledSurfaces)
61  the original surface should be complete.
62 
63 SourceFiles
64  sampledTriSurfaceMesh.C
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #ifndef sampledTriSurfaceMesh_H
69 #define sampledTriSurfaceMesh_H
70 
71 #include "sampledSurface.H"
72 #include "triSurfaceMesh.H"
73 #include "MeshedSurface.H"
74 
75 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76 
77 namespace Foam
78 {
79 
80  class treeDataFace;
81  class meshSearch;
82 
83 namespace sampledSurfaces
84 {
85 
86 /*---------------------------------------------------------------------------*\
87  Class triSurfaceMesh Declaration
88 \*---------------------------------------------------------------------------*/
89 
90 class triSurfaceMesh
91 :
92  public sampledSurface,
93  public MeshedSurface<face>
94 {
95 public:
96  //- Types of communications
97  enum samplingSource
98  {
102  };
103 
104 private:
105 
106  // Private Typedefs
107 
109 
110 
111  // Private Data
112 
113  static const NamedEnum<samplingSource, 3> samplingSourceNames_;
114 
115  //- Surface to sample on
116  const Foam::triSurfaceMesh surface_;
117 
118  //- Whether to sample internal cell values or boundary values
119  const samplingSource sampleSource_;
120 
121  //- Track if the surface needs an update
122  mutable bool needsUpdate_;
123 
124  //- Search tree for all non-coupled boundary faces
125  mutable autoPtr<indexedOctree<treeDataFace>> boundaryTreePtr_;
126 
127  //- From local surface triangle to mesh cell/face.
128  labelList sampleElements_;
129 
130  //- Local points to sample per point
131  pointField samplePoints_;
132 
133 
134  // Private Member Functions
135 
136  //- Get tree of all non-coupled boundary faces
137  const indexedOctree<treeDataFace>& nonCoupledboundaryTree() const;
138 
139  //- Sample field on faces
140  template<class Type>
141  tmp<Field<Type>> sampleField
142  (
144  ) const;
145 
146 
147  template<class Type>
149  interpolateField(const interpolation<Type>&) const;
150 
151  bool update(const meshSearch& meshSearcher);
152 
153 public:
154 
155  //- Runtime type information
156  TypeName("triSurfaceMesh");
157 
158 
159  // Constructors
160 
161  //- Construct from components
163  (
164  const word& name,
165  const polyMesh& mesh,
166  const word& surfaceName,
167  const samplingSource sampleSource
168  );
169 
170  //- Construct from dictionary
172  (
173  const word& name,
174  const polyMesh& mesh,
175  const dictionary& dict
176  );
177 
178  //- Construct from triSurface
180  (
181  const word& name,
182  const polyMesh& mesh,
183  const triSurface& surface,
184  const word& sampleSourceName
185  );
186 
187 
188  //- Destructor
189  virtual ~triSurfaceMesh();
190 
191 
192  // Member Functions
193 
194  //- Does the surface need an update?
195  virtual bool needsUpdate() const;
196 
197  //- Mark the surface as needing an update.
198  // May also free up unneeded data.
199  // Return false if surface was already marked as expired.
200  virtual bool expire();
201 
202  //- Update the surface as required.
203  // Do nothing (and return false) if no update was needed
204  virtual bool update();
205 
206  //- Update the surface using a bound box to limit the searching.
207  // For direct use, i.e. not through sample.
208  // Do nothing (and return false) if no update was needed
209  bool update(const treeBoundBox&);
210 
211  //- Points of surface
212  virtual const pointField& points() const
213  {
214  return MeshStorage::points();
215  }
216 
217  //- Faces of surface
218  virtual const faceList& faces() const
219  {
220  return MeshStorage::faces();
221  }
222 
223 
224  //- Sample field on surface
225  virtual tmp<scalarField> sample
226  (
227  const volScalarField&
228  ) const;
229 
230  //- Sample field on surface
231  virtual tmp<vectorField> sample
232  (
233  const volVectorField&
234  ) const;
235 
236  //- Sample field on surface
238  (
240  ) const;
241 
242  //- Sample field on surface
244  (
245  const volSymmTensorField&
246  ) const;
247 
248  //- Sample field on surface
249  virtual tmp<tensorField> sample
250  (
251  const volTensorField&
252  ) const;
253 
254 
255  //- Interpolate field on surface
257  (
258  const interpolation<scalar>&
259  ) const;
260 
261 
262  //- Interpolate field on surface
264  (
265  const interpolation<vector>&
266  ) const;
267 
268  //- Interpolate field on surface
270  (
272  ) const;
273 
274  //- Interpolate field on surface
276  (
278  ) const;
279 
280  //- Interpolate field on surface
282  (
283  const interpolation<tensor>&
284  ) const;
285 
286  //- Write
287  virtual void print(Ostream&) const;
288 };
289 
290 
291 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
292 
293 } // End namespace sampledSurfaces
294 } // End namespace Foam
295 
296 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
297 
298 #ifdef NoRepository
300 #endif
301 
302 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
303 
304 #endif
305 
306 // ************************************************************************* //
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
Definition: meshSearch.H:57
dictionary dict
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
An abstract class for surfaces with sampling.
const word & name() const
Name of surface.
triSurfaceMesh(const word &name, const polyMesh &mesh, const word &surfaceName, const samplingSource sampleSource)
Construct from components.
A surface geometry mesh with zone information, not to be confused with the similarly named surfaceMes...
Definition: MeshedSurface.H:72
bool interpolate() const
Interpolation requested for surface.
virtual bool update()
Update the surface as required.
Generic GeometricField class.
IOoject and searching on triSurface.
virtual void print(Ostream &) const
Write.
A class for handling words, derived from string.
Definition: word.H:59
const Field< PointType > & points() const
Return reference to global points.
virtual const faceList & faces() const
Faces of surface.
const polyMesh & mesh() const
Access to the underlying mesh.
virtual bool expire()
Mark the surface as needing an update.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
const List< face > & faces() const
Return const access to the faces.
TypeName("triSurfaceMesh")
Runtime type information.
virtual const pointField & points() const
Points of surface.
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
Abstract base class for interpolation.
virtual bool needsUpdate() const
Does the surface need an update?
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
Triangulated surface description with patch information.
Definition: triSurface.H:66
A sampledSurface from a triSurfaceMesh. It samples on the points/triangles of the triSurface...
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
Namespace for OpenFOAM.