probesTemplates.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2017-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "probes.H"
30#include "volFields.H"
31#include "surfaceFields.H"
32#include "IOmanip.H"
33#include "interpolation.H"
34
35// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37namespace Foam
38{
39
40template<class T>
42{
43 void operator()(T& x, const T& y) const
44 {
45 const T unsetVal(-VGREAT*pTraits<T>::one);
46
47 if (x != unsetVal)
48 {
49 // Keep x.
50
51 // Note: should check for y != unsetVal but multiple sample cells
52 // already handled in read().
53 }
54 else
55 {
56 // x is not set. y might be.
57 x = y;
58 }
59 }
60};
61
62} // End namespace Foam
63
64
65// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
66
67template<class GeoField>
69Foam::probes::getOrLoadField(const word& fieldName) const
70{
71 tmp<GeoField> tfield;
72
74 {
75 tfield.reset
76 (
77 new GeoField
78 (
80 (
81 fieldName,
83 mesh_,
86 false
87 ),
88 mesh_
89 )
90 );
91 }
92 else
93 {
94 // Slightly paranoid here
95 tfield.cref(mesh_.cfindObject<GeoField>(fieldName));
96 }
97
98 return tfield;
99}
100
101
102template<class Type>
104(
105 const word& fieldName,
106 const Field<Type>& values
107)
108{
109 const MinMax<Type> limits(values);
110 const Type avgVal = average(values);
111
112 this->setResult("average(" + fieldName + ")", avgVal);
113 this->setResult("min(" + fieldName + ")", limits.min());
114 this->setResult("max(" + fieldName + ")", limits.max());
115 this->setResult("size(" + fieldName + ")", values.size());
116
117 if (verbose_)
118 {
119 Info<< name() << " : " << fieldName << nl
120 << " avg: " << avgVal << nl
121 << " min: " << limits.min() << nl
122 << " max: " << limits.max() << nl << nl;
123 }
124}
125
126
127// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
128
129template<class Type>
130void Foam::probes::writeValues
131(
132 const word& fieldName,
133 const Field<Type>& values,
134 const scalar timeValue
135)
136{
137 if (Pstream::master())
138 {
139 const unsigned int width(IOstream::defaultPrecision() + 7);
140 OFstream& os = *probeFilePtrs_[fieldName];
141
142 os << setw(width) << timeValue;
143
144 forAll(values, probei)
145 {
146 if (includeOutOfBounds_ || processor_[probei] != -1)
147 {
148 os << ' ' << setw(width) << values[probei];
149 }
150 }
151 os << endl;
152 }
153}
154
155
156template<class GeoField>
157void Foam::probes::performAction
158(
159 const fieldGroup<GeoField>& fieldNames,
160 unsigned request
161)
162{
163 for (const word& fieldName : fieldNames)
164 {
165 tmp<GeoField> tfield = getOrLoadField<GeoField>(fieldName);
166
167 if (tfield)
168 {
169 const auto& field = tfield();
170 const scalar timeValue = field.time().timeOutputValue();
171
172 Field<typename GeoField::value_type> values(sample(field));
173
174 this->storeResults(fieldName, values);
175 if (request & ACTION_WRITE)
176 {
177 this->writeValues(fieldName, values, timeValue);
178 }
179 }
180 }
181}
182
183
184// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
185
186template<class Type>
189{
190 const Type unsetVal(-VGREAT*pTraits<Type>::one);
191
192 auto tvalues = tmp<Field<Type>>::New(Field<Type>(this->size(), unsetVal));
193 auto& values = tvalues.ref();
194
195 if (fixedLocations_)
196 {
198 (
199 interpolation<Type>::New(samplePointScheme_, vField)
200 );
201
202 forAll(*this, probei)
203 {
204 if (elementList_[probei] >= 0)
205 {
206 const vector& position = operator[](probei);
207
208 values[probei] = interpPtr().interpolate
209 (
210 position,
211 elementList_[probei],
212 -1
213 );
214 }
215 }
216 }
217 else
218 {
219 forAll(*this, probei)
220 {
221 if (elementList_[probei] >= 0)
222 {
223 values[probei] = vField[elementList_[probei]];
224 }
225 }
226 }
227
229
230 return tvalues;
231}
232
233
234template<class Type>
237{
238 const Type unsetVal(-VGREAT*pTraits<Type>::one);
239
240 auto tvalues = tmp<Field<Type>>::New(Field<Type>(this->size(), unsetVal));
241 auto& values = tvalues.ref();
242
243 forAll(*this, probei)
244 {
245 if (faceList_[probei] >= 0)
246 {
247 values[probei] = sField[faceList_[probei]];
248 }
249 }
250
252
253 return tvalues;
254}
255
256
257template<class Type>
259Foam::probes::sample(const word& fieldName) const
260{
261 return sample(mesh_.lookupObject<VolumeField<Type>>(fieldName));
262}
263
264
265template<class Type>
268{
269 return sample(mesh_.lookupObject<SurfaceField<Type>>(fieldName));
270}
271
272
273// ************************************************************************* //
Istream and Ostream manipulators taking arguments.
scalar y
Generic templated field type.
Definition: Field.H:82
Minimal example by using system/controlDict.functions:
Generic GeometricField class.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static unsigned int defaultPrecision() noexcept
Return the default precision.
Definition: IOstream.H:342
A min/max value pair with additional methods. In addition to conveniently storing values,...
Definition: MinMax.H:128
const T & max() const noexcept
The max value (second)
Definition: MinMaxI.H:121
const T & min() const noexcept
The min value (first)
Definition: MinMaxI.H:107
Output to file stream, using an OSstream.
Definition: OFstream.H:57
static void listCombineAllGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
After completion all processors have the same data.
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
const fvMesh & mesh_
Reference to the fvMesh.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Abstract base class for volume field interpolation.
Definition: interpolation.H:60
const Type * cfindObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
Definition: one.H:62
void storeResults(const word &fieldName, const Field< Type > &values)
Store results: min/max/average/size.
bool loadFromFiles_
Load fields from files (not from objectRegistry)
Definition: probes.H:191
tmp< GeoField > getOrLoadField(const word &fieldName) const
Get from registry or load from disk.
tmp< Field< Type > > sampleSurfaceField(const word &fieldName) const
Sample a surface field at all locations.
splitCell * master() const
Definition: splitCell.H:113
A class for managing temporary objects.
Definition: tmp.H:65
const T & cref() const
Definition: tmpI.H:213
void reset(tmp< T > &&other) noexcept
Clear existing and transfer ownership.
Definition: tmpI.H:314
A class for handling words, derived from Foam::string.
Definition: word.H:68
const volScalarField & T
rDeltaTY field()
OBJstream os(runTime.globalPath()/outputName)
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
Definition: HashOps.H:149
Namespace for OpenFOAM.
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
void operator()(T &x, const T &y) const
Foam::surfaceFields.