surfaceFieldValueTemplates.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-2023 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 "surfaceFieldValue.H"
27 #include "surfaceFields.H"
28 #include "volFields.H"
29 #include "sampledSurface.H"
30 #include "surfaceWriter.H"
31 #include "interpolationCellPoint.H"
32 
33 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
34 
35 template<class Type>
37 (
38  const word& fieldName
39 ) const
40 {
41  if
42  (
44  && obr_.foundObject<SurfaceField<Type>>(fieldName)
45  )
46  {
47  return true;
48  }
49  else if (obr_.foundObject<VolField<Type>>(fieldName))
50  {
51  return true;
52  }
53 
54  return false;
55 }
56 
57 
58 template<class Type>
61 (
62  const word& fieldName
63 ) const
64 {
65  if (selectionType_ == selectionTypes::sampledSurface)
66  {
67  if (obr_.foundObject<VolField<Type>>(fieldName))
68  {
69  const VolField<Type>& fld =
70  obr_.lookupObject<VolField<Type>>(fieldName);
71 
72  if (surfacePtr_().interpolate())
73  {
74  // Interpolate the field to the surface points
75  const interpolationCellPoint<Type> interp(fld);
76  tmp<Field<Type>> tintFld(surfacePtr_().interpolate(interp));
77  const Field<Type>& intFld = tintFld();
78 
79  // Average the interpolated field onto the surface faces
80  const faceList& faces = surfacePtr_().faces();
81  tmp<Field<Type>> tavg(new Field<Type>(faces.size(), Zero));
82  Field<Type>& avg = tavg.ref();
83  forAll(faces, facei)
84  {
85  const face& f = faces[facei];
86  forAll(f, fp)
87  {
88  avg[facei] += intFld[f[fp]];
89  }
90  avg[facei] /= f.size();
91  }
92 
93  return tavg;
94  }
95  else
96  {
97  return surfacePtr_().sample(fld);
98  }
99  }
100  else if (obr_.foundObject<SurfaceField<Type>>(fieldName))
101  {
103  << "Surface field " << fieldName
104  << " cannot be sampled onto surface " << surfacePtr_().name()
105  << ". Only vol fields can be sampled onto surfaces."
106  << abort(FatalError);
107  }
108  }
109  else
110  {
111  if (obr_.foundObject<VolField<Type>>(fieldName))
112  {
113  const VolField<Type>& fld =
114  obr_.lookupObject<VolField<Type>>(fieldName);
115  return filterField(fld);
116  }
117  else if (obr_.foundObject<SurfaceField<Type>>(fieldName))
118  {
119  const SurfaceField<Type>& fld =
120  obr_.lookupObject<SurfaceField<Type>>(fieldName);
121  return filterField(fld);
122  }
123  }
124 
126  << "Field " << fieldName << " not found in database"
127  << abort(FatalError);
128 
129  return tmp<Field<Type>>(nullptr);
130 }
131 
132 
133 template<class Type, class ResultType>
135 (
136  const Field<Type>& values,
137  const scalarField& signs,
138  const scalarField& weights,
139  const vectorField& Sf,
140  ResultType& result
141 ) const
142 {
143  return false;
144 }
145 
146 
147 template<class Type>
149 (
150  const Field<Type>& values,
151  const scalarField& signs,
152  const scalarField& weights,
153  const vectorField& Sf,
154  Type& result
155 ) const
156 {
157  return processValuesTypeType(values, signs, weights, Sf, result);
158 }
159 
160 
161 template<class Type>
163 (
164  const Field<Type>& values,
165  const scalarField& signs,
166  const scalarField& weights,
167  const vectorField& Sf,
168  scalar& result
169 ) const
170 {
171  switch (operation_)
172  {
173  case operationType::minMag:
174  {
175  result = gMin(mag(values));
176  return true;
177  }
178  case operationType::maxMag:
179  {
180  result = gMax(mag(values));
181  return true;
182  }
183  default:
184  {
185  // No fall through
186  return false;
187  }
188  }
189 }
190 
191 
192 template<class Type>
195 (
196  const Field<Type>& values,
197  const scalarField& signs,
198  const scalarField& weights,
199  const vectorField& Sf,
200  Type& result
201 ) const
202 {
203  switch (operation_)
204  {
205  case operationType::sum:
206  {
207  result = gSum(weights*values);
208  return true;
209  }
211  {
212  result = gSum(weights*cmptMag(values));
213  return true;
214  }
215  case operationType::orientedSum:
216  {
217  result = gSum(signs*weights*values);
218  return true;
219  }
221  {
222  result =
223  gSum(weights*values)
224  /stabilise(gSum(weights), vSmall);
225  return true;
226  }
227  case operationType::areaAverage:
228  {
229  const scalarField magSf(mag(Sf));
230  result =
231  gSum(weights*magSf*values)
232  /stabilise(gSum(weights*magSf), vSmall);
233  return true;
234  }
235  case operationType::areaIntegrate:
236  {
237  const scalarField magSf(mag(Sf));
238  result = gSum(weights*magSf*values);
239  return true;
240  }
241  case operationType::min:
242  {
243  result = gMin(values);
244  return true;
245  }
246  case operationType::max:
247  {
248  result = gMax(values);
249  return true;
250  }
251  case operationType::CoV:
252  {
253  const scalarField magSf(mag(Sf));
254 
255  Type meanValue = gSum(values*magSf)/gSum(magSf);
256 
257  const label nComp = pTraits<Type>::nComponents;
258 
259  for (direction d=0; d<nComp; ++d)
260  {
261  scalarField vals(values.component(d));
262  scalar mean = component(meanValue, d);
263  scalar& res = setComponent(result, d);
264 
265  res = sqrt(gSum(magSf*sqr(vals - mean))/gSum(magSf))/mean;
266  }
267 
268  return true;
269  }
270  case operationType::none:
271  {
272  return true;
273  }
274  default:
275  {
276  return false;
277  }
278  }
279 }
280 
281 
282 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
283 
284 template<class Type>
286 (
287  const word& fieldName,
288  const scalarField& signs,
289  const scalarField& weights,
290  const vectorField& Sf
291 )
292 {
293  const bool ok = validField<Type>(fieldName);
294 
295  if (ok)
296  {
297  // Get the values
298  Field<Type> values(getFieldValues<Type>(fieldName));
299 
300  // Write raw values on surface if specified
301  if (writeFields_)
302  {
303  faceList faces;
305 
306  if (selectionType_ == selectionTypes::sampledSurface)
307  {
308  combineSurfaceGeometry(faces, points);
309  }
310  else
311  {
312  combineMeshGeometry(faces, points);
313  }
314 
315  Field<Type> writeValues(weights*values);
316  combineFields(writeValues);
317 
318  if (Pstream::master())
319  {
320  surfaceWriterPtr_->write
321  (
322  outputDir(),
323  fieldName
324  + '_' + selectionTypeNames[selectionType_]
325  + '_' + selectionName_,
326  points,
327  faces,
328  false,
329  fieldName,
330  writeValues
331  );
332  }
333  }
334 
335  // Do the operation
336  if (operation_ != operationType::none)
337  {
338  // Apply scale factor
339  values *= scaleFactor_;
340 
341  bool ok = false;
342 
343  #define writeValuesFieldType(fieldType, none) \
344  ok = \
345  ok \
346  || writeValues<Type, fieldType> \
347  ( \
348  fieldName, \
349  values, \
350  signs, \
351  weights, \
352  Sf \
353  );
355  #undef writeValuesFieldType
356 
357  if (!ok)
358  {
360  << "Operation " << operationTypeNames_[operation_]
361  << " not available for values of type "
363  << exit(FatalError);
364  }
365  }
366  }
367 
368  return ok;
369 }
370 
371 
372 template<class Type, class ResultType>
374 (
375  const word& fieldName,
376  const Field<Type>& values,
377  const scalarField& signs,
378  const scalarField& weights,
379  const vectorField& Sf
380 )
381 {
382  ResultType result;
383 
384  if (processValues(values, signs, weights, Sf, result))
385  {
386  // Add to result dictionary, over-writing any previous entry
387  resultDict_.add(fieldName, result, true);
388 
389  if (Pstream::master())
390  {
391  file() << tab << result;
392 
393  Log << " " << operationTypeNames_[operation_]
394  << "(" << selectionName_ << ") of " << fieldName
395  << " = " << result << endl;
396  }
397 
398  return true;
399  }
400 
401  return false;
402 }
403 
404 
405 template<class Type>
408 (
409  const VolField<Type>& field
410 ) const
411 {
412  tmp<Field<Type>> tvalues(new Field<Type>(faceId_.size()));
413  Field<Type>& values = tvalues.ref();
414 
415  forAll(values, i)
416  {
417  const label facei = faceId_[i];
418  const label patchi = facePatchId_[i];
419 
420  if (patchi >= 0)
421  {
422  values[i] = field.boundaryField()[patchi][facei];
423  }
424  else
425  {
427  << type() << " " << name() << ": "
428  << selectionTypeNames[selectionType_]
429  << "(" << selectionName_ << "):"
430  << nl
431  << " Unable to process internal faces for volume field "
432  << field.name() << nl << abort(FatalError);
433  }
434  }
435 
436  return tvalues;
437 }
438 
439 
440 template<class Type>
443 (
444  const SurfaceField<Type>& field
445 ) const
446 {
447  tmp<Field<Type>> tvalues(new Field<Type>(faceId_.size()));
448  Field<Type>& values = tvalues.ref();
449 
450  forAll(values, i)
451  {
452  const label facei = faceId_[i];
453  const label patchi = facePatchId_[i];
454 
455  if (patchi >= 0)
456  {
457  values[i] = field.boundaryField()[patchi][facei];
458  }
459  else
460  {
461  values[i] = field[facei];
462  }
463  }
464 
465  return tvalues;
466 }
467 
468 
469 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Pre-declare SubField and related Field type.
Definition: Field.H:83
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:479
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const
Return name.
Definition: IOobject.H:310
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
tmp< Field< Type > > getFieldValues(const word &fieldName) const
Return field values by looking up field name.
bool writeValues(const word &fieldName, const scalarField &signs, const scalarField &weights, const vectorField &Sf)
Templated helper function to output field values.
bool processValuesTypeType(const Field< Type > &values, const scalarField &signs, const scalarField &weights, const vectorField &Sf, Type &result) const
Apply a Type -> Type operation to the values.
bool validField(const word &fieldName) const
Return true if the field name is valid.
bool processValues(const Field< Type > &values, const scalarField &signs, const scalarField &weights, const vectorField &Sf, ResultType &result) const
Apply the operation to the values, and return true if successful.
tmp< Field< Type > > filterField(const SurfaceField< Type > &field) const
Filter a surface field according to faceIds.
const objectRegistry & obr_
Reference to the region objectRegistry.
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
bool foundObject(const word &name) const
Is the named Type in registry.
Traits class for primitives.
Definition: pTraits.H:53
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
const pointField & points
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(Cloud< passiveParticle >, particles, iter) { positions[particlei++]=iter().position(mesh);} for(i=0;i< pTraits< point >::nComponents;i++) { forAll(positions, particlei) { gmvFile<< component(positions[particlei], i)<< ' ';} gmvFile<< nl;}}forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define Log
Report write to Foam::Info if the local log switch is true.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > average(const SurfaceField< Type > &ssf)
Area-weighted average a surfaceField creating a volField.
Definition: fvcAverage.C:46
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const zero Zero
Definition: zero.H:97
Type gSum(const FieldField< Field, Type > &f)
tmp< DimensionedField< scalar, GeoMesh > > stabilise(const DimensionedField< scalar, GeoMesh > &dsf, const dimensioned< scalar > &ds)
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
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
label & setComponent(label &l, const direction)
Definition: label.H:86
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
static const char tab
Definition: Ostream.H:265
errorManip< error > abort(error &err)
Definition: errorManip.H:131
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionSet cmptMag(const dimensionSet &)
Definition: dimensionSet.C:296
error FatalError
FOR_ALL_FIELD_TYPES(makeFieldSourceTypedef)
Type gMin(const FieldField< Field, Type > &f)
static const char nl
Definition: Ostream.H:266
Type gMax(const FieldField< Field, Type > &f)
uint8_t direction
Definition: direction.H:45
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensioned< scalar > sumMag(const DimensionedField< Type, GeoMesh > &df)
labelList f(nPoints)
#define writeValuesFieldType(fieldType, none)
Foam::surfaceFields.