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-2026 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"
32 
33 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
34 
35 template<class Type>
36 inline void Foam::functionObjects::fieldValues::surfaceFieldValue::combineField
37 (
38  Field<Type>& field
39 )
40 {
41  List<Field<Type>> allValues(Pstream::nProcs());
42 
43  allValues[Pstream::myProcNo()] = field;
44 
45  Pstream::gatherList(allValues);
46 
47  if (Pstream::master())
48  {
49  field = ListListOps::combine<Field<Type>>
50  (
51  allValues,
52  accessOp<Field<Type>>()
53  );
54  }
55 }
56 
57 
58 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
59 
60 template<class Type>
62 (
63  const word& fieldName
64 ) const
65 {
66  if
67  (
68  selectionType_ != selectionTypes::sampledSurface
69  && obr_.foundObject<SurfaceField<Type>>(fieldName)
70  )
71  {
72  return true;
73  }
74  else if (obr_.foundObject<VolField<Type>>(fieldName))
75  {
76  return true;
77  }
78 
79  return false;
80 }
81 
82 
83 template<class Type>
86 (
87  const word& fieldName
88 ) const
89 {
90  if (selectionType_ == selectionTypes::sampledSurface)
91  {
92  if (obr_.foundObject<VolField<Type>>(fieldName))
93  {
94  const VolField<Type>& fld =
95  obr_.lookupObject<VolField<Type>>(fieldName);
96 
97  if (surfacePtr_().interpolate())
98  {
99  // Interpolate the field to the surface points
100  const interpolations::cellPoint<Type> interp(fld);
101  tmp<Field<Type>> tintFld(surfacePtr_().interpolate(interp));
102  const Field<Type>& intFld = tintFld();
103 
104  // Average the interpolated field onto the surface faces
105  const faceList& faces = surfacePtr_().faces();
106  tmp<Field<Type>> tavg(new Field<Type>(faces.size(), Zero));
107  Field<Type>& avg = tavg.ref();
108  forAll(faces, facei)
109  {
110  const face& f = faces[facei];
111  forAll(f, fp)
112  {
113  avg[facei] += intFld[f[fp]];
114  }
115  avg[facei] /= f.size();
116  }
117 
118  return tavg;
119  }
120  else
121  {
122  return surfacePtr_().sample(fld);
123  }
124  }
125  else if (obr_.foundObject<SurfaceField<Type>>(fieldName))
126  {
128  << "Surface field " << fieldName
129  << " cannot be sampled onto surface " << surfacePtr_().name()
130  << ". Only vol fields can be sampled onto surfaces."
131  << abort(FatalError);
132  }
133  }
134  else
135  {
136  if (obr_.foundObject<VolField<Type>>(fieldName))
137  {
138  const VolField<Type>& fld =
139  obr_.lookupObject<VolField<Type>>(fieldName);
140  return filterField(fld);
141  }
142  else if (obr_.foundObject<SurfaceField<Type>>(fieldName))
143  {
144  const SurfaceField<Type>& fld =
145  obr_.lookupObject<SurfaceField<Type>>(fieldName);
146  return filterField(fld);
147  }
148  }
149 
151  << "Field " << fieldName << " not found in database"
152  << abort(FatalError);
153 
154  return tmp<Field<Type>>(nullptr);
155 }
156 
157 
158 template<class Type, class ResultType>
160 (
161  const Field<Type>& values,
162  const scalarField& signs,
163  const scalarField& weights,
164  const vectorField& Sf,
165  ResultType& result
166 ) const
167 {
168  return false;
169 }
170 
171 
172 template<class Type>
174 (
175  const Field<Type>& values,
176  const scalarField& signs,
177  const scalarField& weights,
178  const vectorField& Sf,
179  Type& result
180 ) const
181 {
182  return processValuesTypeType(values, signs, weights, Sf, result);
183 }
184 
185 
186 template<class Type>
188 (
189  const Field<Type>& values,
190  const scalarField& signs,
191  const scalarField& weights,
192  const vectorField& Sf,
193  scalar& result
194 ) const
195 {
196  switch (operation_)
197  {
198  case operationType::minMag:
199  {
200  result = gMin(mag(values));
201  return true;
202  }
203  case operationType::maxMag:
204  {
205  result = gMax(mag(values));
206  return true;
207  }
208  default:
209  {
210  // No fall through
211  return false;
212  }
213  }
214 }
215 
216 
217 template<class Type>
220 (
221  const Field<Type>& values,
222  const scalarField& signs,
223  const scalarField& weights,
224  const vectorField& Sf,
225  Type& result
226 ) const
227 {
228  switch (operation_)
229  {
230  case operationType::sum:
231  {
232  result = gSum(weights*values);
233  return true;
234  }
236  {
237  result = gSum(weights*cmptMag(values));
238  return true;
239  }
240  case operationType::orientedSum:
241  {
242  result = gSum(signs*weights*values);
243  return true;
244  }
246  {
247  result =
248  gSum(weights*values)
249  /stabilise(gSum(weights), vSmall);
250  return true;
251  }
252  case operationType::areaAverage:
253  {
254  const scalarField magSf(mag(Sf));
255  result =
256  gSum(weights*magSf*values)
257  /stabilise(gSum(weights*magSf), vSmall);
258  return true;
259  }
260  case operationType::areaIntegrate:
261  {
262  const scalarField magSf(mag(Sf));
263  result = gSum(weights*magSf*values);
264  return true;
265  }
266  case operationType::min:
267  {
268  result = gMin(values);
269  return true;
270  }
271  case operationType::max:
272  {
273  result = gMax(values);
274  return true;
275  }
276  case operationType::CoV:
277  {
278  const scalarField magSf(mag(Sf));
279 
280  const Type meanValue = gSum(values*magSf)/gSum(magSf);
281 
282  const label nComp = pTraits<Type>::nComponents;
283 
284  for (direction d=0; d<nComp; ++d)
285  {
286  const scalarField vals(values.component(d));
287  const scalar mean = component(meanValue, d);
288 
289  setComponent(result, d) =
291  (
292  sqrt(gSum(magSf*sqr(vals - mean))/gSum(magSf)),
293  mean
294  );
295  }
296 
297  return true;
298  }
299  case operationType::UI:
300  {
301  const scalarField magSf(mag(Sf));
302 
303  const Type meanValue = gSum(values*magSf)/gSum(magSf);
304 
305  const label nComp = pTraits<Type>::nComponents;
306 
307  for (direction d=0; d<nComp; ++d)
308  {
309  const scalarField vals(values.component(d));
310  const scalar mean = component(meanValue, d);
311 
312  setComponent(result, d) =
313  1 - 0.5*protectedDivide
314  (
315  gSum(magSf*mag(vals - mean))/gSum(magSf),
316  mean
317  );
318  }
319 
320  return true;
321  }
322  case operationType::none:
323  {
324  return true;
325  }
326  default:
327  {
328  return false;
329  }
330  }
331 }
332 
333 
334 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
335 
336 template<class Type>
338 (
339  const word& fieldName,
340  const Field<Type>& values,
341  const scalarField& signs,
342  const scalarField& weights,
343  const vectorField& Sf
344 )
345 {
346  // Do the operation
347  if (operation_ != operationType::none)
348  {
349  bool ok = false;
350 
351  #define writeValuesFieldType(fieldType, none) \
352  { \
353  fieldType result; \
354  \
355  const bool typeOk = \
356  processValues(values, signs, weights, Sf, result); \
357  \
358  if (typeOk) \
359  { \
360  /* Add to result dictionary, over-writing any previous entry */\
361  resultDict_.add(fieldName, result, true); \
362  \
363  /* Write into the file and the log */ \
364  if (Pstream::master()) \
365  { \
366  file() << tab << result; \
367  \
368  Log << " " << operationTypeNames_[operation_] \
369  << "(" << selectionName_.c_str() << ") of " \
370  << fieldName << " = " << result << endl; \
371  } \
372  } \
373  \
374  ok = ok || typeOk; \
375  }
377  #undef writeValuesFieldType
378 
379  if (!ok)
380  {
382  << "Operation " << operationTypeNames_[operation_]
383  << " not available for values of type "
385  << exit(FatalError);
386  }
387  }
388 }
389 
390 
391 template<class Type>
394 (
395  const VolField<Type>& field
396 ) const
397 {
398  tmp<Field<Type>> tvalues(new Field<Type>(faceId_.size()));
399  Field<Type>& values = tvalues.ref();
400 
401  forAll(values, i)
402  {
403  const label facei = faceId_[i];
404  const label patchi = facePatchId_[i];
405 
406  if (patchi >= 0)
407  {
408  values[i] = field.boundaryField()[patchi][facei];
409  }
410  else
411  {
413  << type() << " " << name() << ": "
414  << selectionTypeNames[selectionType_]
415  << "(" << selectionName_.c_str() << "):"
416  << nl
417  << " Unable to process internal faces for volume field "
418  << field.name() << nl << abort(FatalError);
419  }
420  }
421 
422  return tvalues;
423 }
424 
425 
426 template<class Type>
429 (
430  const SurfaceField<Type>& field
431 ) const
432 {
433  tmp<Field<Type>> tvalues(new Field<Type>(faceId_.size()));
434  Field<Type>& values = tvalues.ref();
435 
436  forAll(values, i)
437  {
438  const label facei = faceId_[i];
439  const label patchi = facePatchId_[i];
440 
441  if (patchi >= 0)
442  {
443  values[i] = field.boundaryField()[patchi][facei];
444  }
445  else
446  {
447  values[i] = field[facei];
448  }
449  }
450 
451  return tvalues;
452 }
453 
454 
455 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:449
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:486
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const
Return name.
Definition: IOobject.H:307
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
void writeValues(const word &fieldName, const Field< Type > &values, const scalarField &signs, const scalarField &weights, const vectorField &Sf)
Templated helper function to output field values.
tmp< Field< Type > > getFieldValues(const word &fieldName) const
Return field values by looking up field name.
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.
Piecewise-linear interpolation method. Uses volPointInterpolation to create values on the points....
Traits class for primitives.
Definition: pTraits.H:53
Motion of the mesh specified as a list of pointMeshMovers.
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:197
A class for handling words, derived from string.
Definition: word.H:63
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(lagrangian::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(), lagrangian::cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const unitSet none
Type gMin(const UList< Type > &f, const label comm)
Scalar protectedDivide(const Scalar a, const Scalar b)
Divide two numbers and protect the result from overflowing.
Definition: Scalar.H:381
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const zero Zero
Definition: zero.H:97
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
label & setComponent(label &l, const direction)
Definition: label.H:86
errorManip< error > abort(error &err)
Definition: errorManip.H:131
dimensioned< Type > average(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
FOR_ALL_FIELD_TYPES(makeDimensionedPointFieldFunctions)
dimensioned< scalar > sumMag(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
tmp< DimensionedField< scalar, GeoMesh, Field > > stabilise(const DimensionedField< scalar, GeoMesh, PrimitiveField > &dsf, const dimensioned< scalar > &ds)
Type gSum(const UList< Type > &f, const label comm)
tmp< DimensionedField< typename outerProduct< Type, Type >::type, GeoMesh, Field >> sqr(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
void component(GeometricField< typename GeometricField< Type, GeoMesh, PrimitiveField1 >::cmptType, GeoMesh, PrimitiveField1 > &gcf, const GeometricField< Type, GeoMesh, PrimitiveField2 > &gf, const direction d)
dimensioned< Type > min(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
Type gMax(const UList< Type > &f, const label comm)
word name(const LagrangianState state)
Return a string representation of a Lagrangian state enumeration.
tmp< DimensionedField< scalar, GeoMesh, Field > > mag(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
error FatalError
tmp< DimensionedField< Type, GeoMesh, Field > > cmptMag(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
void sqrt(LagrangianPatchField< scalar > &f, const LagrangianPatchField< scalar > &f1)
static const char nl
Definition: Ostream.H:297
uint8_t direction
Definition: direction.H:45
dimensioned< Type > max(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
labelList f(nPoints)
#define writeValuesFieldType(fieldType, none)
Foam::surfaceFields.