volFieldValueTemplates.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-2018 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 "volFieldValue.H"
27 #include "volFields.H"
28 
29 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
30 
31 template<class Type>
33 (
34  const word& fieldName
35 ) const
36 {
38 
39  if (obr_.foundObject<vf>(fieldName))
40  {
41  return true;
42  }
43 
44  return false;
45 }
46 
47 
48 template<class Type>
51 (
52  const word& fieldName,
53  const bool mustGet
54 ) const
55 {
57 
58  if (obr_.foundObject<vf>(fieldName))
59  {
60  return filterField(obr_.lookupObject<vf>(fieldName));
61  }
62 
63  if (mustGet)
64  {
66  << "Field " << fieldName << " not found in database"
67  << abort(FatalError);
68  }
69 
70  return tmp<Field<Type>>(new Field<Type>(0.0));
71 }
72 
73 
74 template<class Type>
76 (
77  const Field<Type>& values,
78  const scalarField& V,
79  const scalarField& weightField
80 ) const
81 {
82  Type result = Zero;
83  switch (operation_)
84  {
85  case opSum:
86  {
87  result = gSum(values);
88  break;
89  }
90  case opWeightedSum:
91  {
92  result = gSum(weightField*values);
93  break;
94  }
95  case opSumMag:
96  {
97  result = gSum(cmptMag(values));
98  break;
99  }
100  case opAverage:
101  {
102  result = gSum(values)/nCells();
103  break;
104  }
105  case opWeightedAverage:
106  {
107  result = gSum(weightField*values)/gSum(weightField);
108  break;
109  }
110  case opVolAverage:
111  {
112  result = gSum(V*values)/this->V();
113  break;
114  }
115  case opWeightedVolAverage:
116  {
117  result = gSum(weightField*V*values)/gSum(weightField*V);
118  break;
119  }
120  case opVolIntegrate:
121  {
122  result = gSum(V*values);
123  break;
124  }
125  case opWeightedVolIntegrate:
126  {
127  result = gSum(weightField*V*values);
128  break;
129  }
130  case opMin:
131  {
132  result = gMin(values);
133  break;
134  }
135  case opMax:
136  {
137  result = gMax(values);
138  break;
139  }
140  case opCoV:
141  {
142  Type meanValue = gSum(values*V)/this->V();
143 
144  const label nComp = pTraits<Type>::nComponents;
145 
146  for (direction d=0; d<nComp; ++d)
147  {
148  scalarField vals(values.component(d));
149  scalar mean = component(meanValue, d);
150  scalar& res = setComponent(result, d);
151 
152  res = sqrt(gSum(V*sqr(vals - mean))/this->V())/mean;
153  }
154 
155  break;
156  }
157  case opNone:
158  {}
159  }
160 
161  return result;
162 }
163 
164 
165 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
166 
167 template<class Type>
169 (
170  const word& fieldName
171 )
172 {
173  const bool ok = validField<Type>(fieldName);
174 
175  if (ok)
176  {
177  Field<Type> values(setFieldValues<Type>(fieldName));
178  scalarField V(filterField(fieldValue::mesh_.V()));
179  scalarField weightField(values.size(), 1.0);
180 
181  if (weightFieldName_ != "none")
182  {
183  weightField = setFieldValues<scalar>(weightFieldName_, true);
184  }
185 
186  Type result = processValues(values, V, weightField);
187 
188  if (Pstream::master())
189  {
190  // Add to result dictionary, over-writing any previous entry
191  resultDict_.add(fieldName, result, true);
192 
193  if (writeFields_)
194  {
196  (
197  IOobject
198  (
199  fieldName + '_' + regionTypeNames_[regionType_]
200  + '-' + volRegion::regionName_,
201  obr_.time().timeName(),
202  obr_,
203  IOobject::NO_READ,
204  IOobject::NO_WRITE
205  ),
206  weightField*values
207  ).write();
208  }
209 
210 
211  file()<< tab << result;
212 
213  Log << " " << operationTypeNames_[operation_]
214  << "(" << volRegion::regionName_ << ") of " << fieldName
215  << " = " << result << endl;
216  }
217  }
218 
219  return ok;
220 }
221 
222 
223 template<class Type>
226 (
227  const Field<Type>& field
228 ) const
229 {
230  if (isNull(cellIDs()))
231  {
232  return field;
233  }
234  else
235  {
236  return tmp<Field<Type>>(new Field<Type>(field, cellIDs()));
237  }
238 }
239 
240 
241 // ************************************************************************* //
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 const char tab
Definition: Ostream.H:264
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
tmp< Field< Type > > setFieldValues(const word &fieldName, const bool mustGet=false) const
Insert field values into values list.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
uint8_t direction
Definition: direction.H:45
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Traits class for primitives.
Definition: pTraits.H:50
Generic GeometricField class.
bool isNull(const T &t)
Return true if t is a reference to the nullObject of type T.
Definition: nullObjectI.H:40
Type gSum(const FieldField< Field, Type > &f)
Pre-declare SubField and related Field type.
Definition: Field.H:57
bool validField(const word &fieldName) const
Return true if the field name is valid.
A class for handling words, derived from string.
Definition: word.H:59
tmp< Field< Type > > filterField(const Field< Type > &field) const
Filter a field according to cellIds.
Type processValues(const Field< Type > &values, const scalarField &V, const scalarField &weightField) const
Apply the &#39;operation&#39; to the values.
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:657
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)
#define Log
Report write to Foam::Info if the local log switch is true.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual bool write(const bool valid=true) const
Write using setting from DB.
bool writeValues(const word &fieldName)
Templated helper function to output field values.
label & setComponent(label &l, const direction)
Definition: label.H:79
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50