AveragingMethod.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) 2013-2019 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 "AveragingMethod.H"
27 #include "runTimeSelectionTables.H"
28 #include "pointMesh.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class Type>
34 {}
35 
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 template<class Type>
41 (
42  const IOobject& io,
43  const dictionary& dict,
44  const fvMesh& mesh,
45  const labelList& size
46 )
47 :
48  regIOobject(io),
50  dict_(dict),
51  mesh_(mesh)
52 {
53  forAll(size, i)
54  {
56  (
57  new Field<Type>(size[i], Zero)
58  );
59  }
60 }
61 
62 
63 template<class Type>
65 (
66  const AveragingMethod<Type>& am
67 )
68 :
69  regIOobject(am),
71  dict_(am.dict_),
72  mesh_(am.mesh_)
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
77 
78 template<class Type>
81 (
82  const IOobject& io,
83  const dictionary& dict,
84  const fvMesh& mesh
85 )
86 {
87  word averageType(dict.lookup(typeName));
88 
89  // Info<< "Selecting averaging method "
90  // << averageType << endl;
91 
92  typename dictionaryConstructorTable::iterator cstrIter =
93  dictionaryConstructorTablePtr_->find(averageType);
94 
95  if (cstrIter == dictionaryConstructorTablePtr_->end())
96  {
98  << "Unknown averaging method " << averageType
99  << ", constructor not in hash table" << nl << nl
100  << " Valid averaging methods are:" << nl
101  << dictionaryConstructorTablePtr_->sortedToc()
102  << abort(FatalError);
103  }
104 
105  return autoPtr<AveragingMethod<Type>>(cstrIter()(io, dict, mesh));
106 }
107 
108 
109 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
110 
111 template<class Type>
113 {}
114 
115 
116 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
117 
118 template<class Type>
120 {
121  updateGrad();
122 }
123 
124 
125 template<class Type>
127 (
128  const AveragingMethod<scalar>& weight
129 )
130 {
131  updateGrad();
132 
133  *this /= max(weight, small);
134 }
135 
136 
137 template<class Type>
139 {
140  return os.good();
141 }
142 
143 
144 template<class Type>
146 {
147  const pointMesh pointMesh_(mesh_);
148 
149  // point volumes
150  Field<scalar> pointVolume(mesh_.nPoints(), 0);
151 
152  // output fields
154  (
155  IOobject
156  (
157  this->name() + ":cellValue",
158  this->time().timeName(),
159  mesh_
160  ),
161  mesh_,
162  dimensioned<Type>("zero", dimless, Zero)
163  );
165  (
166  IOobject
167  (
168  this->name() + ":cellGrad",
169  this->time().timeName(),
170  mesh_
171  ),
172  mesh_,
174  );
176  (
177  IOobject
178  (
179  this->name() + ":pointValue",
180  this->time().timeName(),
181  mesh_
182  ),
183  pointMesh_,
184  dimensioned<Type>("zero", dimless, Zero)
185  );
187  (
188  IOobject
189  (
190  this->name() + ":pointGrad",
191  this->time().timeName(),
192  mesh_
193  ),
194  pointMesh_,
196  );
197 
198  // Barycentric coordinates of the tet vertices
200  tetCrds
201  ({
202  barycentric(1, 0, 0, 0),
203  barycentric(0, 1, 0, 0),
204  barycentric(0, 0, 1, 0),
205  barycentric(0, 0, 0, 1)
206  });
207 
208  // tet-volume weighted sums
209  forAll(mesh_.C(), celli)
210  {
211  const List<tetIndices> cellTets =
212  polyMeshTetDecomposition::cellTetIndices(mesh_, celli);
213 
214  forAll(cellTets, tetI)
215  {
216  const tetIndices& tetIs = cellTets[tetI];
217  const triFace triIs = tetIs.faceTriIs(mesh_);
218  const scalar v = tetIs.tet(mesh_).mag();
219 
220  cellValue[celli] += v*interpolate(tetCrds[0], tetIs);
221  cellGrad[celli] += v*interpolateGrad(tetCrds[0], tetIs);
222 
223  forAll(triIs, vertexI)
224  {
225  const label pointi = triIs[vertexI];
226 
227  pointVolume[pointi] += v;
228  pointValue[pointi] += v*interpolate(tetCrds[vertexI], tetIs);
229  pointGrad[pointi] += v*interpolateGrad(tetCrds[vertexI], tetIs);
230  }
231  }
232  }
233 
234  // average
235  cellValue.primitiveFieldRef() /= mesh_.V();
236  cellGrad.primitiveFieldRef() /= mesh_.V();
237  pointValue.primitiveFieldRef() /= pointVolume;
238  pointGrad.primitiveFieldRef() /= pointVolume;
239 
240  // write
241  if (!cellValue.write(write)) return false;
242  if (!cellGrad.write(write)) return false;
243  if (!pointValue.write(write)) return false;
244  if (!pointGrad.write(write)) return false;
245 
246  return true;
247 }
248 
249 
250 // ************************************************************************* //
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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 autoPtr< AveragingMethod< Type > > New(const IOobject &io, const dictionary &dict, const fvMesh &mesh)
Selector.
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:54
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Generic GeometricField class.
Generic dimensioned Type class.
Barycentric< scalar > barycentric
A scalar version of the templated Barycentric.
Definition: barycentric.H:45
triFace faceTriIs(const polyMesh &mesh) const
Return the indices corresponding to the tri on the face for.
Definition: tetIndicesI.H:64
tetPointRef tet(const polyMesh &mesh) const
Return the geometry corresponding to this tet.
Definition: tetIndicesI.H:98
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
Generic field type.
Definition: FieldField.H:51
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
dynamicFvMesh & mesh
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
scalar mag() const
Return volume.
Definition: tetrahedronI.H:170
Base class for lagrangian averaging methods.
Definition: MPPICParcel.H:55
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const fvMesh & mesh_
The mesh on which the averaging is to be done.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:260
const dictionary & dict_
Protected data.
virtual bool write(const bool write=true) const
Write using setting from DB.
virtual void updateGrad()
Protected member functions.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
virtual bool writeData(Ostream &) const
Dummy write.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual void average()
Calculate the average.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Macros to ease declaration of run-time selection tables.
virtual bool write(const bool write=true) const
Write using setting from DB.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
AveragingMethod(const IOobject &io, const dictionary &dict, const fvMesh &mesh, const labelList &size)
Constructors.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
virtual ~AveragingMethod()
Destructor.