cellModel.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-2016 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 "cellModel.H"
27 #include "pyramid.H"
28 
29 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
30 
32 (
33  const labelList& pointLabels,
34  const pointField& points
35 ) const
36 {
37  // Estimate centre of cell
38  vector cEst = Zero;
39 
40  // Sum the points idicated by the label list
41  forAll(pointLabels, i)
42  {
43  cEst += points[pointLabels[i]];
44  }
45 
46  // Average by dividing by the number summed over.
47  cEst /= scalar(pointLabels.size());
48 
49 
50  // Calculate the centre by breaking the cell into pyramids and
51  // volume-weighted averaging their centres
52  scalar sumV = 0.0;
53  vector sumVc = Zero;
54 
55  const faceList cellFaces = faces(pointLabels);
56 
57  forAll(cellFaces, i)
58  {
59  const face& curFace = cellFaces[i];
60 
61  scalar pyrVol =
63  (
64  curFace,
65  cEst
66  ).mag(points);
67 
68  if (pyrVol > SMALL)
69  {
71  << "zero or negative pyramid volume: " << -pyrVol
72  << " for face " << i
73  << endl;
74  }
75 
76  sumVc -=
77  pyrVol
79  .centre(points);
80 
81  sumV -= pyrVol;
82  }
83 
84  return sumVc/(sumV + VSMALL);
85 }
86 
87 
88 Foam::scalar Foam::cellModel::mag
89 (
90  const labelList& pointLabels,
91  const pointField& points
92 ) const
93 {
94  // Estimate centre of cell
95  vector cEst = Zero;
96 
97  // Sum the points idicated by the label list
98  forAll(pointLabels, i)
99  {
100  cEst += points[pointLabels[i]];
101  }
102 
103  // Average by dividing by the number summed over.
104  cEst /= scalar(pointLabels.size());
105 
106 
107  // Calculate the magnitude by summing the -mags of the pyramids
108  // The sign change is because the faces point outwards
109  // and a pyramid is constructed from an inward pointing face
110  // and the base centre-apex vector
111  scalar v = 0;
112 
113  const faceList cellFaces = faces(pointLabels);
114 
115  forAll(cellFaces, i)
116  {
117  const face& curFace =cellFaces[i];
118 
119  scalar pyrVol =
121  (
122  curFace,
123  cEst
124  ).mag(points);
125 
126  if (pyrVol > SMALL)
127  {
129  << "zero or negative pyramid volume: " << -pyrVol
130  << " for face " << i
131  << endl;
132  }
133 
134  v -= pyrVol;
135  }
136 
137  return v;
138 }
139 
140 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
scalar mag(const labelList &pointLabels, const pointField &points) const
Cell volume.
Definition: cellModel.C:89
vector centre(const labelList &pointLabels, const pointField &points) const
Vector centroid.
Definition: cellModel.C:32
static const zero Zero
Definition: zero.H:91
#define WarningInFunction
Report a warning using Foam::Warning.
dimensioned< scalar > mag(const dimensioned< Type > &)
A geometric pyramid primitive with a base of &#39;n&#39; sides: i.e. a parametric pyramid. A pyramid is constructed from a base polygon and an apex point.
Definition: pyramid.H:47