sampledCuttingPlane.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-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 Class
25  Foam::sampledCuttingPlane
26 
27 Description
28  A sampledSurface defined by a plane using the iso-surface algorithm
29  to 'cut' the mesh.
30 
31 SourceFiles
32  sampledCuttingPlane.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef sampledCuttingPlane_H
37 #define sampledCuttingPlane_H
38 
39 #include "sampledSurface.H"
40 #include "isoSurface.H"
41 //#include "isoSurfaceCell.H"
42 #include "plane.H"
43 #include "ZoneIDs.H"
44 #include "fvMeshSubset.H"
45 
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48 namespace Foam
49 {
50 
51 /*---------------------------------------------------------------------------*\
52  Class sampledCuttingPlane Declaration
53 \*---------------------------------------------------------------------------*/
54 
56 :
57  public sampledSurface
58 {
59  // Private data
60 
61  //- Plane
62  const plane plane_;
63 
64  //- Merge tolerance
65  const scalar mergeTol_;
66 
67  //- Whether to coarsen
68  const Switch regularise_;
69 
70  //- Whether to recalculate cell values as average of point values
71  const Switch average_;
72 
73  //- Zone name/index (if restricted to zones)
74  mutable cellZoneID zoneID_;
75 
76  //- For zones: patch to put exposed faces into
77  mutable word exposedPatchName_;
78 
79  //- Track if the surface needs an update
80  mutable bool needsUpdate_;
81 
82 
83  //- Optional subsetted mesh
84  autoPtr<fvMeshSubset> subMeshPtr_;
85 
86  //- Distance to cell centres
87  autoPtr<volScalarField> cellDistancePtr_;
88 
89  //- Distance to points
90  scalarField pointDistance_;
91 
92  //- Constructed iso surface
93  // autoPtr<isoSurfaceCell> isoSurfPtr_;
94  autoPtr<isoSurface> isoSurfPtr_;
95 
96  //- Triangles converted to faceList
97  mutable autoPtr<faceList> facesPtr_;
98 
99 
100  // Private Member Functions
101 
102  //- Create iso surface
103  void createGeometry();
104 
105  //- Sample field on faces
106  template<class Type>
107  tmp<Field<Type>> sampleField
108  (
110  ) const;
111 
112 
113  template<class Type>
115  interpolateField(const interpolation<Type>&) const;
116 
117 
118 public:
119 
120  //- Runtime type information
121  TypeName("sampledCuttingPlane");
122 
123 
124  // Constructors
125 
126  //- Construct from dictionary
128  (
129  const word& name,
130  const polyMesh& mesh,
131  const dictionary& dict
132  );
133 
134 
135  //- Destructor
136  virtual ~sampledCuttingPlane();
137 
138 
139  // Member Functions
140 
141  //- Does the surface need an update?
142  virtual bool needsUpdate() const;
143 
144  //- Mark the surface as needing an update.
145  // May also free up unneeded data.
146  // Return false if surface was already marked as expired.
147  virtual bool expire();
148 
149  //- Update the surface as required.
150  // Do nothing (and return false) if no update was needed
151  virtual bool update();
152 
153  //- Points of surface
154  virtual const pointField& points() const
155  {
156  return surface().points();
157  }
158 
159  //- Faces of surface
160  virtual const faceList& faces() const
161  {
162  if (facesPtr_.empty())
163  {
164  const triSurface& s = surface();
165 
166  facesPtr_.reset(new faceList(s.size()));
167 
168  forAll(s, i)
169  {
170  facesPtr_()[i] = s[i].triFaceFace();
171  }
172  }
173  return facesPtr_;
174  }
175 
176 
177  // const isoSurfaceCell& surface() const
178  const isoSurface& surface() const
179  {
180  return isoSurfPtr_();
181  }
182 
183  //- Sample field on surface
184  virtual tmp<scalarField> sample
185  (
186  const volScalarField&
187  ) const;
188 
189  //- Sample field on surface
190  virtual tmp<vectorField> sample
191  (
192  const volVectorField&
193  ) const;
194 
195  //- Sample field on surface
197  (
199  ) const;
200 
201  //- Sample field on surface
203  (
204  const volSymmTensorField&
205  ) const;
206 
207  //- Sample field on surface
208  virtual tmp<tensorField> sample
209  (
210  const volTensorField&
211  ) const;
212 
213 
214  //- Interpolate field on surface
216  (
217  const interpolation<scalar>&
218  ) const;
219 
220  //- Interpolate field on surface
222  (
223  const interpolation<vector>&
224  ) const;
225 
226  //- Interpolate field on surface
228  (
230  ) const;
231 
232  //- Interpolate field on surface
234  (
236  ) const;
237 
238  //- Interpolate field on surface
240  (
241  const interpolation<tensor>&
242  ) const;
243 
244  //- Write
245  virtual void print(Ostream&) const;
246 };
247 
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 } // End namespace Foam
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 #ifdef NoRepository
257 #endif
258 
259 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
260 
261 #endif
262 
263 // ************************************************************************* //
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
virtual ~sampledCuttingPlane()
Destructor.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
void reset(T *=nullptr)
If object pointer already set, delete object and set to given.
Definition: autoPtrI.H:114
virtual bool expire()
Mark the surface as needing an update.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
An abstract class for surfaces with sampling.
const word & name() const
Name of surface.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
bool interpolate() const
Interpolation requested for surface.
List< face > faceList
Definition: faceListFwd.H:43
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
Generic GeometricField class.
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
Definition: plane.H:60
bool empty() const
Return true if the autoPtr is empty (ie, no pointer set)
Definition: autoPtrI.H:76
virtual const pointField & points() const
Points of surface.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
virtual bool needsUpdate() const
Does the surface need an update?
A class for handling words, derived from string.
Definition: word.H:59
const isoSurface & surface() const
TypeName("sampledCuttingPlane")
Runtime type information.
const Field< PointType > & points() const
Return reference to global points.
const polyMesh & mesh() const
Access to the underlying mesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual void print(Ostream &) const
Write.
virtual const faceList & faces() const
Faces of surface.
sampledCuttingPlane(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
Abstract base class for interpolation.
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 surface formed by the iso value. After "Regularised Marching Tetrahedra: improved iso-surface extra...
Definition: isoSurface.H:83
A class for managing temporary objects.
Definition: PtrList.H:53
Triangulated surface description with patch information.
Definition: triSurface.H:66
virtual bool update()
Update the surface as required.
Namespace for OpenFOAM.
A sampledSurface defined by a plane using the iso-surface algorithm to &#39;cut&#39; the mesh.