sampledSurfaceTemplates.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 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 "sampledSurface.H"
27 
28 template<class Type>
29 bool Foam::sampledSurface::checkFieldSize(const Field<Type>& field) const
30 {
31  if (faces().empty() || field.empty())
32  {
33  return false;
34  }
35 
36  if (field.size() != faces().size())
37  {
39  (
40  "sampledSurface::checkFieldSize(const Field<Type>&) const"
41  )
42  << "size mismatch: "
43  << "field (" << field.size()
44  << ") != surface (" << faces().size() << ")"
45  << exit(FatalError);
46  }
47 
48  return true;
49 }
50 
51 
52 template<class Type>
54 {
55  Type value = pTraits<Type>::zero;
56 
57  if (checkFieldSize(field))
58  {
59  value = sum(field*magSf());
60  }
61 
62  reduce(value, sumOp<Type>());
63  return value;
64 }
65 
66 
67 template<class Type>
69 {
70  Type value = integrate(field());
71  field.clear();
72  return value;
73 }
74 
75 
76 template<class Type>
78 {
79  Type value = pTraits<Type>::zero;
80 
81  if (checkFieldSize(field))
82  {
83  value = sum(field*magSf());
84  }
85 
86  reduce(value, sumOp<Type>());
87 
88  // avoid divide-by-zero
89  if (area())
90  {
91  return value/area();
92  }
93  else
94  {
95  return pTraits<Type>::zero;
96  }
97 }
98 
99 
100 template<class Type>
102 {
103  Type value = average(field());
104  field.clear();
105  return value;
106 }
107 
108 
109 template<class ReturnType, class Type>
110 void Foam::sampledSurface::project
111 (
112  Field<ReturnType>& res,
113  const Field<Type>& field
114 ) const
115 {
116  if (checkFieldSize(field))
117  {
118  const vectorField& norm = Sf();
119 
120  forAll(norm, faceI)
121  {
122  res[faceI] = field[faceI] & (norm[faceI]/mag(norm[faceI]));
123  }
124  }
125  else
126  {
127  res.clear();
128  }
129 }
130 
131 
132 template<class ReturnType, class Type>
133 void Foam::sampledSurface::project
134 (
135  Field<ReturnType>& res,
136  const tmp<Field<Type> >& field
137 ) const
138 {
139  project(res, field());
140  field.clear();
141 }
142 
143 
144 template<class ReturnType, class Type>
146 Foam::sampledSurface::project
147 (
148  const tmp<Field<Type> >& field
149 ) const
150 {
151  tmp<Field<ReturnType> > tRes(new Field<ReturnType>(faces().size()));
152  project(tRes(), field);
153  return tRes;
154 }
155 
156 
157 template<class Type>
160 (
162 ) const
163 {
164  const fvMesh& mesh = dynamic_cast<const fvMesh&>(pfld.mesh()());
165 
167  (
169  (
170  IOobject
171  (
172  "cellAvg",
173  mesh.time().timeName(),
174  pfld.db(),
177  false
178  ),
179  mesh,
181  )
182  );
183  GeometricField<Type, fvPatchField, volMesh>& cellAvg = tcellAvg();
184 
185  labelField nPointCells(mesh.nCells(), 0);
186  {
187  for (label pointI = 0; pointI < mesh.nPoints(); pointI++)
188  {
189  const labelList& pCells = mesh.pointCells(pointI);
190 
191  forAll(pCells, i)
192  {
193  label cellI = pCells[i];
194 
195  cellAvg[cellI] += pfld[pointI];
196  nPointCells[cellI]++;
197  }
198  }
199  }
200  forAll(cellAvg, cellI)
201  {
202  cellAvg[cellI] /= nPointCells[cellI];
203  }
204  // Give value to calculatedFvPatchFields
205  cellAvg.correctBoundaryConditions();
206 
207  return tcellAvg;
208 }
209 
210 
211 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensioned< scalar > mag(const dimensioned< Type > &)
bool empty() const
Return true if the UList is empty (ie, size() is zero).
Definition: UListI.H:313
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
Type integrate(const Field< Type > &) const
Integration of a field across the surface.
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
label nCells() const
const labelListList & pointCells() const
const polyMesh & mesh() const
Access to the underlying mesh.
virtual const scalarField & magSf() const
Return face area magnitudes.
const Mesh & mesh() const
Return mesh.
Type average(const Field< Type > &) const
Area-averaged value of a field across the surface.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
void clear()
Clear the list, i.e. set size to zero.
virtual const vectorField & Sf() const
Return face area vectors.
Generic dimensioned Type class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define forAll(list, i)
Definition: UList.H:421
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Pre-declare SubField and related Field type.
Definition: Field.H:57
scalar area() const
The total surface area.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
void correctBoundaryConditions()
Correct boundary field.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual const faceList & faces() const =0
Faces of surface.
tmp< GeometricField< Type, fvPatchField, volMesh > > pointAverage(const GeometricField< Type, pointPatchField, pointMesh > &pfld) const
Interpolate from points to cell centre.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
A class for managing temporary objects.
Definition: PtrList.H:118
label nPoints() const