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