faceSourceTemplates.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 | Copyright (C) 2011-2015 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 "faceSource.H"
27 #include "surfaceFields.H"
28 #include "volFields.H"
29 #include "sampledSurface.H"
30 #include "interpolationCellPoint.H"
31 
32 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
33 
34 template<class Type>
35 bool Foam::fieldValues::faceSource::validField(const word& fieldName) const
36 {
39 
40  if (source_ != stSampledSurface && obr_.foundObject<sf>(fieldName))
41  {
42  return true;
43  }
44  else if (obr_.foundObject<vf>(fieldName))
45  {
46  return true;
47  }
48 
49  return false;
50 }
51 
52 
53 template<class Type>
55 (
56  const word& fieldName,
57  const bool mustGet,
58  const bool applyOrientation
59 ) const
60 {
63 
64  if (source_ != stSampledSurface && obr_.foundObject<sf>(fieldName))
65  {
66  return filterField(obr_.lookupObject<sf>(fieldName), applyOrientation);
67  }
68  else if (obr_.foundObject<vf>(fieldName))
69  {
70  const vf& fld = obr_.lookupObject<vf>(fieldName);
71 
72  if (surfacePtr_.valid())
73  {
74  if (surfacePtr_().interpolate())
75  {
76  const interpolationCellPoint<Type> interp(fld);
77  tmp<Field<Type> > tintFld(surfacePtr_().interpolate(interp));
78  const Field<Type>& intFld = tintFld();
79 
80  // Average
81  const faceList& faces = surfacePtr_().faces();
82  tmp<Field<Type> > tavg
83  (
85  );
86  Field<Type>& avg = tavg();
87 
88  forAll(faces, faceI)
89  {
90  const face& f = faces[faceI];
91  forAll(f, fp)
92  {
93  avg[faceI] += intFld[f[fp]];
94  }
95  avg[faceI] /= f.size();
96  }
97 
98  return tavg;
99  }
100  else
101  {
102  return surfacePtr_().sample(fld);
103  }
104  }
105  else
106  {
107  return filterField(fld, applyOrientation);
108  }
109  }
110 
111  if (mustGet)
112  {
114  (
115  "Foam::tmp<Foam::Field<Type> > "
116  "Foam::fieldValues::faceSource::getFieldValues"
117  "("
118  "const word&, "
119  "const bool, "
120  "const bool"
121  ") const"
122  ) << "Field " << fieldName << " not found in database"
123  << abort(FatalError);
124  }
125 
126  return tmp<Field<Type> >(new Field<Type>(0));
127 }
128 
129 
130 template<class Type>
132 (
133  const Field<Type>& values,
134  const vectorField& Sf,
135  const scalarField& weightField
136 ) const
137 {
138  Type result = pTraits<Type>::zero;
139  switch (operation_)
140  {
141  case opSum:
142  {
143  result = sum(values);
144  break;
145  }
146  case opSumMag:
147  {
148  result = sum(cmptMag(values));
149  break;
150  }
151  case opSumDirection:
152  {
154  (
155  "template<class Type>"
156  "Type Foam::fieldValues::faceSource::processSameTypeValues"
157  "("
158  "const Field<Type>&, "
159  "const vectorField&, "
160  "const scalarField&"
161  ") const"
162  )
163  << "Operation " << operationTypeNames_[operation_]
164  << " not available for values of type "
166  << exit(FatalError);
167 
168  result = pTraits<Type>::zero;
169  break;
170  }
172  {
174  (
175  "template<class Type>"
176  "Type Foam::fieldValues::faceSource::processSameTypeValues"
177  "("
178  "const Field<Type>&, "
179  "const vectorField&, "
180  "const scalarField&"
181  ") const"
182  )
183  << "Operation " << operationTypeNames_[operation_]
184  << " not available for values of type "
186  << exit(FatalError);
187 
188  result = pTraits<Type>::zero;
189  break;
190  }
191  case opAverage:
192  {
193  result = sum(values)/values.size();
194  break;
195  }
196  case opWeightedAverage:
197  {
198  if (weightField.size())
199  {
200  result = sum(weightField*values)/sum(weightField);
201  }
202  else
203  {
204  result = sum(values)/values.size();
205  }
206  break;
207  }
208  case opAreaAverage:
209  {
210  const scalarField magSf(mag(Sf));
211 
212  result = sum(magSf*values)/sum(magSf);
213  break;
214  }
216  {
217  const scalarField magSf(mag(Sf));
218 
219  if (weightField.size())
220  {
221  result = sum(weightField*magSf*values)/sum(magSf*weightField);
222  }
223  else
224  {
225  result = sum(magSf*values)/sum(magSf);
226  }
227  break;
228  }
229  case opAreaIntegrate:
230  {
231  const scalarField magSf(mag(Sf));
232 
233  result = sum(magSf*values);
234  break;
235  }
236  case opMin:
237  {
238  result = min(values);
239  break;
240  }
241  case opMax:
242  {
243  result = max(values);
244  break;
245  }
246  case opCoV:
247  {
248  const scalarField magSf(mag(Sf));
249 
250  Type meanValue = sum(values*magSf)/sum(magSf);
251 
252  const label nComp = pTraits<Type>::nComponents;
253 
254  for (direction d=0; d<nComp; ++d)
255  {
256  scalarField vals(values.component(d));
257  scalar mean = component(meanValue, d);
258  scalar& res = setComponent(result, d);
259 
260  res = sqrt(sum(magSf*sqr(vals - mean))/sum(magSf))/mean;
261  }
262 
263  break;
264  }
265  default:
266  {
267  // Do nothing
268  }
269  }
270 
271  return result;
272 }
273 
274 
275 template<class Type>
277 (
278  const Field<Type>& values,
279  const vectorField& Sf,
280  const scalarField& weightField
281 ) const
282 {
283  return processSameTypeValues(values, Sf, weightField);
284 }
285 
286 
287 
288 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
289 
290 template<class Type>
292 (
293  const word& fieldName,
294  const scalarField& weightField,
295  const bool orient
296 )
297 {
298  const bool ok = validField<Type>(fieldName);
299 
300  if (ok)
301  {
302  Field<Type> values(getFieldValues<Type>(fieldName, true, orient));
303 
304  vectorField Sf;
305  if (surfacePtr_.valid())
306  {
307  // Get oriented Sf
308  Sf = surfacePtr_().Sf();
309  }
310  else
311  {
312  // Get oriented Sf
313  Sf = filterField(mesh().Sf(), true);
314  }
315 
316  // Combine onto master
317  combineFields(values);
318  combineFields(Sf);
319 
320  // Write raw values on surface if specified
321  if (surfaceWriterPtr_.valid())
322  {
323  faceList faces;
325 
326  if (surfacePtr_.valid())
327  {
328  combineSurfaceGeometry(faces, points);
329  }
330  else
331  {
332  combineMeshGeometry(faces, points);
333  }
334 
335  if (Pstream::master())
336  {
337  fileName outputDir =
338  baseFileDir()/name_/"surface"/obr_.time().timeName();
339 
340  surfaceWriterPtr_->write
341  (
342  outputDir,
344  points,
345  faces,
346  fieldName,
347  values,
348  false
349  );
350  }
351  }
352 
353 
354  // Apply scale factor
355  values *= scaleFactor_;
356 
357  if (Pstream::master())
358  {
359  Type result = processValues(values, Sf, weightField);
360 
361  // Add to result dictionary, over-writing any previous entry
362  resultDict_.add(fieldName, result, true);
363 
364  file()<< tab << result;
365 
366  if (log_) Info<< " " << operationTypeNames_[operation_]
367  << "(" << sourceName_ << ") of " << fieldName
368  << " = " << result << endl;
369  }
370  }
371 
372  return ok;
373 }
374 
375 
376 template<class Type>
378 (
380  const bool applyOrientation
381 ) const
382 {
383  tmp<Field<Type> > tvalues(new Field<Type>(faceId_.size()));
384  Field<Type>& values = tvalues();
385 
386  forAll(values, i)
387  {
388  label faceI = faceId_[i];
389  label patchI = facePatchId_[i];
390  if (patchI >= 0)
391  {
392  values[i] = field.boundaryField()[patchI][faceI];
393  }
394  else
395  {
397  (
398  "fieldValues::faceSource::filterField"
399  "("
400  "const GeometricField<Type, fvPatchField, volMesh>&, "
401  "const bool"
402  ") const"
403  ) << type() << " " << name_ << ": "
404  << sourceTypeNames_[source_] << "(" << sourceName_ << "):"
405  << nl
406  << " Unable to process internal faces for volume field "
407  << field.name() << nl << abort(FatalError);
408  }
409  }
410 
411  if (applyOrientation)
412  {
413  forAll(values, i)
414  {
415  values[i] *= faceSign_[i];
416  }
417  }
418 
419  return tvalues;
420 }
421 
422 
423 template<class Type>
425 (
427  const bool applyOrientation
428 ) const
429 {
430  tmp<Field<Type> > tvalues(new Field<Type>(faceId_.size()));
431  Field<Type>& values = tvalues();
432 
433  forAll(values, i)
434  {
435  label faceI = faceId_[i];
436  label patchI = facePatchId_[i];
437  if (patchI >= 0)
438  {
439  values[i] = field.boundaryField()[patchI][faceI];
440  }
441  else
442  {
443  values[i] = field[faceI];
444  }
445  }
446 
447  if (applyOrientation)
448  {
449  forAll(values, i)
450  {
451  values[i] *= faceSign_[i];
452  }
453  }
454 
455  return tvalues;
456 }
457 
458 
459 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
OFstream & file()
Return access to the file (if only 1)
labelList faceSign_
List of +1/-1 representing face flip map.
Definition: faceSource.H:424
const pointField & points
static const NamedEnum< sourceType, 3 > sourceTypeNames_
Source type names.
Definition: faceSource.H:321
unsigned char direction
Definition: direction.H:43
word name_
Name of this fieldValue object.
Definition: fieldValue.H:72
word sourceName_
Name of source object.
Definition: fieldValue.H:87
Foam::surfaceFields.
labelList faceId_
Local list of face IDs.
Definition: faceSource.H:417
bool writeValues(const word &fieldName, const scalarField &weightField, const bool orient)
Templated helper function to output field values.
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
static const NamedEnum< operationType, 15 > operationTypeNames_
Operation type names.
Definition: faceSource.H:345
labelList f(nPoints)
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:578
labelList facePatchId_
Local list of patch ID per face.
Definition: faceSource.H:420
bool validField(const word &fieldName) const
Return true if the field name is valid.
tmp< Field< Type > > filterField(const GeometricField< Type, fvsPatchField, surfaceMesh > &field, const bool applyOrientation) const
Filter a surface field according to faceIds.
bool foundObject(const word &name) const
Is the named Type found?
sourceType source_
Source type.
Definition: faceSource.H:387
void combineFields(Field< Type > &field)
Combine fields from all processor domains into single field.
label & setComponent(label &l, const direction)
Definition: label.H:79
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
A class for handling words, derived from string.
Definition: word.H:59
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
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
autoPtr< surfaceWriter > surfaceWriterPtr_
Surface writer.
Definition: faceSource.H:384
Switch log_
Switch to send output to Info as well as to file.
Definition: fieldValue.H:84
messageStream Info
tmp< Field< Type > > getFieldValues(const word &fieldName, const bool mustGet=false, const bool applyOrientation=false) const
Return field values by looking up field name.
volScalarField sf(fieldObject, mesh)
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const fvMesh & mesh() const
Helper function to return the reference to the mesh.
Definition: fieldValueI.H:79
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const objectRegistry & obr_
Database this class is registered to.
Definition: fieldValue.H:75
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:739
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Generic GeometricField class.
static const char tab
Definition: Ostream.H:259
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
A class for handling file names.
Definition: fileName.H:69
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
Type processValues(const Field< Type > &values, const vectorField &Sf, const scalarField &weightField) const
Apply the &#39;operation&#39; to the values. Wrapper around.
dictionary resultDict_
Results dictionary for external access of results.
Definition: fieldValue.H:96
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
const Time & time() const
Return time.
Type processSameTypeValues(const Field< Type > &values, const vectorField &Sf, const scalarField &weightField) const
Apply the &#39;operation&#39; to the values. Operation has to.
virtual fileName baseFileDir() const
Return the base directory for output.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
A class for managing temporary objects.
Definition: PtrList.H:118
scalar scaleFactor_
Scale factor - optional.
Definition: faceSource.H:402
operationType operation_
Operation to apply to values.
Definition: faceSource.H:390
autoPtr< sampledSurface > surfacePtr_
Underlying sampledSurface.
Definition: faceSource.H:430