volFieldValueTemplates.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-2017 OpenFOAM Foundation
9  Copyright (C) 2015-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 "volFieldValue.H"
30 #include "volFields.H"
31 
32 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
33 
34 template<class Type>
36 (
37  const word& fieldName
38 ) const
39 {
41  typedef typename VolFieldType::Internal IntVolFieldType;
42 
43  return
44  (
45  obr_.foundObject<VolFieldType>(fieldName)
46  || obr_.foundObject<IntVolFieldType>(fieldName)
47  );
48 }
49 
50 
51 template<class Type>
54 (
55  const word& fieldName,
56  const bool mandatory
57 ) const
58 {
60  typedef typename VolFieldType::Internal IntVolFieldType;
61 
62  if (obr_.foundObject<VolFieldType>(fieldName))
63  {
64  return filterField(obr_.lookupObject<VolFieldType>(fieldName));
65  }
66  else if (obr_.foundObject<IntVolFieldType>(fieldName))
67  {
68  return filterField(obr_.lookupObject<IntVolFieldType>(fieldName));
69  }
70 
71  if (mandatory)
72  {
74  << "Field " << fieldName << " not found in database" << nl
75  << abort(FatalError);
76  }
77 
78  return tmp<Field<Type>>::New();
79 }
80 
81 
82 template<class Type>
84 (
85  const Field<Type>& values,
86  const scalarField& V,
87  const scalarField& weightField
88 ) const
89 {
90  Type result = Zero;
91  switch (operation_)
92  {
93  case opNone:
94  {
95  break;
96  }
97  case opMin:
98  {
99  result = gMin(values);
100  break;
101  }
102  case opMax:
103  {
104  result = gMax(values);
105  break;
106  }
107  case opSumMag:
108  {
109  result = gSum(cmptMag(values));
110  break;
111  }
112  case opSum:
113  case opWeightedSum:
114  {
115  if (is_weightedOp() && canWeight(weightField))
116  {
117  result = gSum(weightField*values);
118  }
119  else
120  {
121  // Unweighted form
122  result = gSum(values);
123  }
124  break;
125  }
126  case opAverage:
127  case opWeightedAverage:
128  {
129  if (is_weightedOp() && canWeight(weightField))
130  {
131  result =
132  gSum(weightField*values)/(gSum(weightField) + ROOTVSMALL);
133  }
134  else
135  {
136  // Unweighted form
137  const label n = returnReduce(values.size(), sumOp<label>());
138  result = gSum(values)/(scalar(n) + ROOTVSMALL);
139  }
140  break;
141  }
142  case opVolAverage:
143  case opWeightedVolAverage:
144  {
145  if (is_weightedOp() && canWeight(weightField))
146  {
147  result = gSum(weightField*V*values)
148  /(gSum(weightField*V) + ROOTVSMALL);
149  }
150  else
151  {
152  // Unweighted form
153  result = gSum(V*values)/(gSum(V) + ROOTVSMALL);
154  }
155  break;
156  }
157  case opVolIntegrate:
158  case opWeightedVolIntegrate:
159  {
160  if (is_weightedOp() && canWeight(weightField))
161  {
162  result = gSum(weightField*V*values);
163  }
164  else
165  {
166  // Unweighted form
167  result = gSum(V*values);
168  }
169  break;
170  }
171  case opCoV:
172  {
173  const scalar sumV = gSum(V);
174 
175  Type meanValue = gSum(V*values)/sumV;
176 
177  for (direction d=0; d < pTraits<Type>::nComponents; ++d)
178  {
179  tmp<scalarField> vals(values.component(d));
180  const scalar mean = component(meanValue, d);
181  scalar& res = setComponent(result, d);
182 
183  res = sqrt(gSum(V*sqr(vals - mean))/sumV)/(mean + ROOTVSMALL);
184  }
185 
186  break;
187  }
188  }
189 
190  return result;
191 }
192 
193 
194 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
195 
196 template<class Type>
198 (
199  const word& fieldName,
200  const scalarField& V,
201  const scalarField& weightField
202 )
203 {
204  const bool ok = validField<Type>(fieldName);
205 
206  if (ok)
207  {
208  Field<Type> values(getFieldValues<Type>(fieldName));
209 
210  if (writeFields_)
211  {
212  word outName = fieldName + '_' + regionTypeNames_[regionType_];
213  if (this->volRegion::regionName_ != polyMesh::defaultRegion)
214  {
215  outName = outName + '-' + this->volRegion::regionName_;
216  }
217 
219  (
220  IOobject
221  (
222  outName,
223  obr_.time().timeName(),
224  obr_,
225  IOobject::NO_READ,
226  IOobject::NO_WRITE
227  ),
228  weightField.empty()
229  ? scaleFactor_*values
230  : scaleFactor_*weightField*values
231  ).write();
232  }
233 
234  if (operation_ != opNone)
235  {
236  // Apply scale factor
237  values *= scaleFactor_;
238 
239  Type result = processValues(values, V, weightField);
240 
241  switch (postOperation_)
242  {
243  case postOpSqrt:
244  {
245  // sqrt: component-wise - does not change the type
246  for (direction d=0; d < pTraits<Type>::nComponents; ++d)
247  {
248  setComponent(result, d)
249  = sqrt(mag(component(result, d)));
250  }
251  break;
252  }
253  default:
254  {
255  break;
256  }
257  }
258 
259  // Write state/results information
260  word prefix, suffix;
261  {
262  if (postOperation_ != postOpNone)
263  {
264  // Adjust result name to include post-operation
265  prefix += postOperationTypeNames_[postOperation_];
266  prefix += '(';
267  suffix += ')';
268  }
269 
270  prefix += operationTypeNames_[operation_];
271  prefix += '(';
272  suffix += ')';
273  }
274 
275  word regionPrefix;
276  if (this->volRegion::regionName_ != polyMesh::defaultRegion)
277  {
278  regionPrefix = this->volRegion::regionName_ + ',';
279  }
280 
281  word resultName = prefix + regionPrefix + fieldName + suffix;
282 
283  Log << " " << prefix << this->volRegion::regionName_ << suffix
284  << " of " << fieldName << " = ";
285 
286 
287  // Operation or post-operation returns scalar?
288 
289  scalar sresult{0};
290 
291  bool alwaysScalar(operation_ & typeScalar);
292 
293  if (alwaysScalar)
294  {
295  sresult = component(result, 0);
296 
297  if (postOperation_ == postOpMag)
298  {
299  sresult = mag(sresult);
300  }
301  }
302  else if (postOperation_ == postOpMag)
303  {
304  sresult = mag(result);
305  alwaysScalar = true;
306  }
307 
308 
309  if (alwaysScalar)
310  {
311  file()<< tab << sresult;
312 
313  Log << sresult << endl;
314 
315  this->setResult(resultName, sresult);
316  }
317  else
318  {
319  file()<< tab << result;
320 
321  Log << result << endl;
322 
323  this->setResult(resultName, result);
324  }
325  }
326  }
327 
328  return ok;
329 }
330 
331 
332 template<class Type>
335 (
336  const Field<Type>& field
337 ) const
338 {
339  if (this->volRegion::useAllCells())
340  {
341  return field;
342  }
343 
344  return tmp<Field<Type>>::New(field, cellIDs());
345 }
346 
347 
348 // ************************************************************************* //
volFields.H
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::setComponent
label & setComponent(label &l, const direction)
Definition: label.H:123
Log
#define Log
Definition: PDRblock.C:35
Foam::component
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
Definition: FieldFieldFunctions.C:44
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::IOField
A primitive field of type <T> with automated input and output.
Definition: foamVtkLagrangianWriter.H:61
Foam::returnReduce
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Definition: PstreamReduceOps.H:94
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::HashTableOps::values
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
Definition: HashOps.H:149
volFieldValue.H
Foam::functionObjects::fieldValues::volFieldValue::processValues
Type processValues(const Field< Type > &values, const scalarField &V, const scalarField &weightField) const
Apply the 'operation' to the values.
Definition: volFieldValueTemplates.C:84
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::sumOp
Definition: ops.H:213
Foam::cmptMag
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:400
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::functionObjects::fieldValues::volFieldValue::validField
bool validField(const word &fieldName) const
Return true if the field name is valid.
Definition: volFieldValueTemplates.C:36
Foam::Field
Generic templated field type.
Definition: Field.H:63
field
rDeltaTY field()
Foam::functionObjects::fieldValues::volFieldValue::writeValues
bool writeValues(const word &fieldName, const scalarField &V, const scalarField &weightField)
Templated helper function to output field values.
Definition: volFieldValueTemplates.C:198
Foam::functionObjects::fieldValues::volFieldValue::filterField
tmp< Field< Type > > filterField(const Field< Type > &field) const
Filter a field according to cellIds.
Foam::FatalError
error FatalError
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::tab
constexpr char tab
Definition: Ostream.H:403
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::functionObjects::fieldValues::volFieldValue::getFieldValues
tmp< Field< Type > > getFieldValues(const word &fieldName, const bool mandatory=false) const
Insert field values into values list.
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::direction
uint8_t direction
Definition: direction.H:52
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:36
Foam::gMin
Type gMin(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:593
Foam::GeometricField< Type, fvPatchField, volMesh >
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:592