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