probesTemplates.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "probes.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "IOmanip.H"
30 #include "interpolation.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 template<class T>
38 class isNotEqOp
39 {
40 public:
41 
42  void operator()(T& x, const T& y) const
43  {
44  const T unsetVal(-vGreat*pTraits<T>::one);
45 
46  if (x != unsetVal)
47  {
48  // Keep x.
49 
50  // Note: should check for y != unsetVal but multiple sample cells
51  // already handled in read().
52  }
53  else
54  {
55  // x is not set. y might be.
56  x = y;
57  }
58  }
59 };
60 
61 }
62 
63 
64 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
65 
66 template<class Type>
67 void Foam::probes::sampleAndWrite
68 (
70 )
71 {
72  Field<Type> values(sample(vField));
73 
74  if (Pstream::master())
75  {
76  const unsigned int w = IOstream::defaultPrecision() + 7;
77  OFstream& os = *probeFilePtrs_[vField.name()];
78 
79  os << setw(w) << vField.time().timeToUserTime(vField.time().value());
80 
81  forAll(values, probei)
82  {
83  OStringStream buf;
84  buf << values[probei];
85  os << ' ' << setw(w) << buf.str().c_str();
86  }
87  os << endl;
88  }
89 }
90 
91 
92 template<class Type>
93 void Foam::probes::sampleAndWrite
94 (
96 )
97 {
98  Field<Type> values(sample(sField));
99 
100  if (Pstream::master())
101  {
102  const unsigned int w = IOstream::defaultPrecision() + 7;
103  OFstream& os = *probeFilePtrs_[sField.name()];
104 
105  os << sField.time().timeToUserTime(sField.time().value());
106 
107  forAll(values, probei)
108  {
109  OStringStream buf;
110  buf << values[probei];
111  os << ' ' << setw(w) << buf.str().c_str();
112  }
113  os << endl;
114  }
115 }
116 
117 
118 template<class Type>
119 void Foam::probes::sampleAndWrite(const fieldGroup<Type>& fields)
120 {
121  forAll(fields, fieldi)
122  {
123  if (loadFromFiles_)
124  {
125  sampleAndWrite
126  (
128  (
129  IOobject
130  (
131  fields[fieldi],
132  mesh_.time().timeName(),
133  mesh_,
136  false
137  ),
138  mesh_
139  )
140  );
141  }
142  else
143  {
144  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
145 
146  if
147  (
148  iter != objectRegistry::end()
149  && iter()->type()
151  )
152  {
153  sampleAndWrite
154  (
155  mesh_.lookupObject
157  (
158  fields[fieldi]
159  )
160  );
161  }
162  }
163  }
164 }
165 
166 
167 template<class Type>
168 void Foam::probes::sampleAndWriteSurfaceFields(const fieldGroup<Type>& fields)
169 {
170  forAll(fields, fieldi)
171  {
172  if (loadFromFiles_)
173  {
174  sampleAndWrite
175  (
177  (
178  IOobject
179  (
180  fields[fieldi],
181  mesh_.time().timeName(),
182  mesh_,
185  false
186  ),
187  mesh_
188  )
189  );
190  }
191  else
192  {
193  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
194 
195  if
196  (
197  iter != objectRegistry::end()
198  && iter()->type()
200  )
201  {
202  sampleAndWrite
203  (
204  mesh_.lookupObject
206  (
207  fields[fieldi]
208  )
209  );
210  }
211  }
212  }
213 }
214 
215 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
216 
217 template<class Type>
220 (
222 ) const
223 {
224  const Type unsetVal(-vGreat*pTraits<Type>::one);
225 
226  tmp<Field<Type>> tValues
227  (
228  new Field<Type>(this->size(), unsetVal)
229  );
230 
231  Field<Type>& values = tValues.ref();
232 
233  if (fixedLocations_)
234  {
235  autoPtr<interpolation<Type>> interpolator
236  (
237  interpolation<Type>::New(interpolationScheme_, vField)
238  );
239 
240  forAll(*this, probei)
241  {
242  if (elementList_[probei] >= 0)
243  {
244  const vector& position = operator[](probei);
245 
246  values[probei] = interpolator().interpolate
247  (
248  position,
249  elementList_[probei],
250  -1
251  );
252  }
253  }
254  }
255  else
256  {
257  forAll(*this, probei)
258  {
259  if (elementList_[probei] >= 0)
260  {
261  values[probei] = vField[elementList_[probei]];
262  }
263  }
264  }
265 
268 
269  return tValues;
270 }
271 
272 
273 template<class Type>
275 Foam::probes::sample(const word& fieldName) const
276 {
277  return sample
278  (
280  (
281  fieldName
282  )
283  );
284 }
285 
286 
287 template<class Type>
290 (
292 ) const
293 {
294  const Type unsetVal(-vGreat*pTraits<Type>::one);
295 
296  tmp<Field<Type>> tValues
297  (
298  new Field<Type>(this->size(), unsetVal)
299  );
300 
301  Field<Type>& values = tValues.ref();
302 
303  forAll(*this, probei)
304  {
305  if (faceList_[probei] >= 0)
306  {
307  values[probei] = sField[faceList_[probei]];
308  }
309  }
310 
313 
314  return tValues;
315 }
316 
317 
318 template<class Type>
321 {
322  return sample
323  (
325  (
326  fieldName
327  )
328  );
329 }
330 
331 // ************************************************************************* //
Foam::surfaceFields.
void operator()(T &x, const T &y) const
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const word & name() const
Return name.
Definition: IOobject.H:303
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
Output to file stream.
Definition: OFstream.H:82
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Traits class for primitives.
Definition: pTraits.H:50
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
Generic GeometricField class.
tmp< Field< Type > > sampleSurfaceFields(const word &fieldName) const
Sample a single scalar field on all sample locations.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:103
tmp< Field< Type > > sample(const GeometricField< Type, fvPatchField, volMesh > &) const
Sample a volume field at all locations.
scalar y
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
friend class const_iterator
Declare friendship with the const_iterator.
Definition: HashTable.H:197
Istream and Ostream manipulators taking arguments.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const Time & time() const
Return time.
Definition: IOobject.C:328
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
string str() const
Return the string.
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
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Output to memory buffer stream.
Definition: OStringStream.H:49
Namespace for OpenFOAM.