probes.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-2022 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::probes
26 
27 Description
28  Set of locations to sample.
29 
30  Call write() to sample and write files.
31 
32 SourceFiles
33  probes.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef probes_H
38 #define probes_H
39 
40 #include "functionObject.H"
41 #include "HashPtrTable.H"
42 #include "OFstream.H"
43 #include "polyMesh.H"
44 #include "pointField.H"
45 #include "volFieldsFwd.H"
46 #include "surfaceFieldsFwd.H"
47 #include "surfaceMesh.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 // Forward declaration of classes
55 class Time;
56 class objectRegistry;
57 class dictionary;
58 class fvMesh;
59 class polyTopoChangeMap;
60 
61 /*---------------------------------------------------------------------------*\
62  Class probes Declaration
63 \*---------------------------------------------------------------------------*/
64 
65 class probes
66 :
67  public functionObject,
68  public pointField
69 {
70 protected:
71 
72  // Protected classes
73 
74  //- Class used for grouping field types
75  template<class Type>
76  class fieldGroup
77  :
78  public DynamicList<word>
79  {
80  public:
81  //- Construct null
82  fieldGroup()
83  :
84  DynamicList<word>(0)
85  {}
86  };
87 
88 
89  // Protected member data
90 
91  //- Const reference to fvMesh
92  const fvMesh& mesh_;
93 
94 
95  // Read from dictionary
96 
97  //- Names of fields to probe
99 
100  //- Fixed locations, default = yes
101  // Note: set to false for moving mesh calculations where locations
102  // should move with the mesh
103  bool fixedLocations_;
104 
105  //- Interpolation scheme name
106  // Note: only possible when fixedLocations_ is true
108 
109 
110  // Calculated
111 
112  //- Categorised scalar/vector/tensor vol fields
118 
119  //- Categorised scalar/vector/tensor surf fields
125 
126  // Cells to be probed (obtained from the locations)
128 
129  // Faces to be probed
131 
132  //- Current open files
134 
135 
136  // Protected Member Functions
137 
138  //- Clear old field groups
139  void clearFieldGroups();
140 
141  //- Append fieldName to the appropriate group
142  label appendFieldGroup(const word& fieldName, const word& fieldType);
143 
144  //- Classify field types, returns the number of fields
146 
147  //- Find cells and faces containing probes
148  virtual void findElements(const fvMesh&);
149 
150  //- Classify field type and Open/close file streams,
151  // returns number of fields to sample
152  label prepare();
153 
154 
155 private:
156 
157  //- Sample and write a particular volume field
158  template<class Type>
159  void sampleAndWrite
160  (
161  const VolField<Type>&
162  );
163 
164 
165  //- Sample and write a particular surface field
166  template<class Type>
167  void sampleAndWrite
168  (
169  const SurfaceField<Type>&
170  );
171 
172  //- Sample and write all the fields of the given type
173  template<class Type>
174  void sampleAndWrite(const fieldGroup<Type>&);
175 
176  //- Sample and write all the surface fields of the given type
177  template<class Type>
178  void sampleAndWriteSurfaceFields(const fieldGroup<Type>&);
179 
180 
181 public:
182 
183  //- Runtime type information
184  TypeName("probes");
185 
186 
187  // Constructors
188 
189  //- Construct from Time and dictionary
190  probes
191  (
192  const word& name,
193  const Time& time,
194  const dictionary& dict
195  );
196 
197  //- Disallow default bitwise copy construction
198  probes(const probes&) = delete;
199 
200 
201  //- Destructor
202  virtual ~probes();
203 
204 
205  // Member Functions
206 
207  //- Return the list of fields required
208  virtual wordList fields() const;
209 
210  //- Return locations to probe
211  virtual const pointField& probeLocations() const
212  {
213  return *this;
214  }
215 
216  //- Return location for probe i
217  virtual const point& probe(const label i) const
218  {
219  return operator[](i);
220  }
221 
222  //- Cells to be probed (obtained from the locations)
223  const labelList& elements() const
224  {
225  return elementList_;
226  }
227 
228  //- Read the probes
229  virtual bool read(const dictionary&);
230 
231  //- Execute, currently does nothing
232  virtual bool execute();
233 
234  //- Sample and write
235  virtual bool write();
236 
237  //- Update topology using the given map
238  virtual void movePoints(const polyMesh&);
239 
240  //- Update topology using the given map
241  virtual void topoChange(const polyTopoChangeMap&);
242 
243  //- Update from another mesh using the given map
244  virtual void mapMesh(const polyMeshMap&);
245 
246  //- Update topology using the given map due to readUpdate
247  virtual void readUpdate(const polyMesh::readUpdateState state)
248  {}
249 
250  //- Sample a volume field at all locations
251  template<class Type>
253  (
254  const VolField<Type>&
255  ) const;
256 
257  //- Sample a single vol field on all sample locations
258  template<class Type>
259  tmp<Field<Type>> sample(const word& fieldName) const;
260 
261  //- Sample a single scalar field on all sample locations
262  template<class Type>
263  tmp<Field<Type>> sampleSurfaceFields(const word& fieldName) const;
264 
265  //- Sample a surface field at all locations
266  template<class Type>
268  (
269  const SurfaceField<Type>&
270  ) const;
271 
272 
273  // Member Operators
274 
275  //- Disallow default bitwise assignment
276  void operator=(const probes&) = delete;
277 };
278 
279 
280 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
281 
282 } // End namespace Foam
283 
284 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
285 
286 #ifdef NoRepository
287  #include "probesTemplates.C"
288 #endif
289 
290 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
291 
292 #endif
293 
294 // ************************************************************************* //
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:78
Generic GeometricField class.
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:67
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
T & operator[](const label)
Return element of UList.
Definition: UListI.H:167
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base-class for Time/database functionObjects.
const word & name() const
Return the name of this functionObject.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:90
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Class used for grouping field types.
Definition: probes.H:78
fieldGroup()
Construct null.
Definition: probes.H:81
Set of locations to sample.
Definition: probes.H:68
virtual const point & probe(const label i) const
Return location for probe i.
Definition: probes.H:216
fieldGroup< symmTensor > surfaceSymmTensorFields_
Definition: probes.H:122
HashPtrTable< OFstream > probeFilePtrs_
Current open files.
Definition: probes.H:132
virtual ~probes()
Destructor.
Definition: probes.C:288
void operator=(const probes &)=delete
Disallow default bitwise assignment.
fieldGroup< tensor > tensorFields_
Definition: probes.H:116
const fvMesh & mesh_
Const reference to fvMesh.
Definition: probes.H:91
void clearFieldGroups()
Clear old field groups.
tmp< Field< Type > > sample(const VolField< Type > &) const
Sample a volume field at all locations.
label classifyFields()
Classify field types, returns the number of fields.
TypeName("probes")
Runtime type information.
virtual wordList fields() const
Return the list of fields required.
Definition: probes.C:320
virtual void readUpdate(const polyMesh::readUpdateState state)
Update topology using the given map due to readUpdate.
Definition: probes.H:246
fieldGroup< tensor > surfaceTensorFields_
Definition: probes.H:123
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: probes.C:364
label prepare()
Classify field type and Open/close file streams,.
Definition: probes.C:167
const labelList & elements() const
Cells to be probed (obtained from the locations)
Definition: probes.H:222
fieldGroup< symmTensor > symmTensorFields_
Definition: probes.H:115
probes(const word &name, const Time &time, const dictionary &dict)
Construct from Time and dictionary.
Definition: probes.C:260
tmp< Field< Type > > sampleSurfaceFields(const word &fieldName) const
Sample a single scalar field on all sample locations.
fieldGroup< sphericalTensor > surfaceSphericalTensorFields_
Definition: probes.H:121
wordList fields_
Names of fields to probe.
Definition: probes.H:97
labelList elementList_
Definition: probes.H:126
fieldGroup< scalar > surfaceScalarFields_
Categorised scalar/vector/tensor surf fields.
Definition: probes.H:119
virtual void findElements(const fvMesh &)
Find cells and faces containing probes.
Definition: probes.C:50
fieldGroup< vector > surfaceVectorFields_
Definition: probes.H:120
bool fixedLocations_
Fixed locations, default = yes.
Definition: probes.H:102
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: probes.C:447
virtual const pointField & probeLocations() const
Return locations to probe.
Definition: probes.H:210
virtual void movePoints(const polyMesh &)
Update topology using the given map.
Definition: probes.C:353
fieldGroup< scalar > scalarFields_
Categorised scalar/vector/tensor vol fields.
Definition: probes.H:112
labelList faceList_
Definition: probes.H:129
word interpolationScheme_
Interpolation scheme name.
Definition: probes.H:106
label appendFieldGroup(const word &fieldName, const word &fieldType)
Append fieldName to the appropriate group.
fieldGroup< sphericalTensor > sphericalTensorFields_
Definition: probes.H:114
virtual bool execute()
Execute, currently does nothing.
Definition: probes.C:326
virtual bool write()
Sample and write.
Definition: probes.C:332
virtual bool read(const dictionary &)
Read the probes.
Definition: probes.C:294
fieldGroup< vector > vectorFields_
Definition: probes.H:113
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
Namespace for OpenFOAM.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dictionary dict