edgeStats.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-2012 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 "edgeStats.H"
27 #include "Time.H"
28 #include "polyMesh.H"
29 #include "Ostream.H"
30 #include "twoDPointCorrector.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 const Foam::scalar Foam::edgeStats::edgeTol_ = 1e-3;
35 
36 
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
40 Foam::direction Foam::edgeStats::getNormalDir
41 (
42  const twoDPointCorrector* correct2DPtr
43 ) const
44 {
45  direction dir = 3;
46 
47  if (correct2DPtr)
48  {
49  const vector& normal = correct2DPtr->planeNormal();
50 
51  if (mag(normal & vector(1, 0, 0)) > 1-edgeTol_)
52  {
53  dir = 0;
54  }
55  else if (mag(normal & vector(0, 1, 0)) > 1-edgeTol_)
56  {
57  dir = 1;
58  }
59  else if (mag(normal & vector(0, 0, 1)) > 1-edgeTol_)
60  {
61  dir = 2;
62  }
63  }
64  return dir;
65 }
66 
67 
68 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
69 
70 // Construct from mesh
71 Foam::edgeStats::edgeStats(const polyMesh& mesh)
72 :
73  mesh_(mesh),
74  normalDir_(3)
75 {
76  IOobject motionObj
77  (
78  "motionProperties",
79  mesh.time().constant(),
80  mesh,
83  );
84 
85  if (motionObj.headerOk())
86  {
87  Info<< "Reading " << mesh.time().constant() / "motionProperties"
88  << endl << endl;
89 
90  IOdictionary motionProperties(motionObj);
91 
92  Switch twoDMotion(motionProperties.lookup("twoDMotion"));
93 
94  if (twoDMotion)
95  {
96  Info<< "Correcting for 2D motion" << endl << endl;
97 
98  autoPtr<twoDPointCorrector> correct2DPtr
99  (
100  new twoDPointCorrector(mesh)
101  );
102 
103  normalDir_ = getNormalDir(&correct2DPtr());
104  }
105  }
106 }
107 
108 
109 // Construct from components
110 Foam::edgeStats::edgeStats
111 (
112  const polyMesh& mesh,
113  const twoDPointCorrector* correct2DPtr
114 )
115 :
116  mesh_(mesh),
117  normalDir_(getNormalDir(correct2DPtr))
118 {}
119 
120 
121 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
122 
123 Foam::scalar Foam::edgeStats::minLen(Ostream& os) const
124 {
125  label nX = 0;
126  label nY = 0;
127  label nZ = 0;
128 
129  scalar minX = GREAT;
130  scalar maxX = -GREAT;
131  vector x(1, 0, 0);
132 
133  scalar minY = GREAT;
134  scalar maxY = -GREAT;
135  vector y(0, 1, 0);
136 
137  scalar minZ = GREAT;
138  scalar maxZ = -GREAT;
139  vector z(0, 0, 1);
140 
141  scalar minOther = GREAT;
142  scalar maxOther = -GREAT;
143 
144  const edgeList& edges = mesh_.edges();
145 
146  forAll(edges, edgeI)
147  {
148  const edge& e = edges[edgeI];
149 
150  vector eVec(e.vec(mesh_.points()));
151 
152  scalar eMag = mag(eVec);
153 
154  eVec /= eMag;
155 
156  if (mag(eVec & x) > 1-edgeTol_)
157  {
158  minX = min(minX, eMag);
159  maxX = max(maxX, eMag);
160  nX++;
161  }
162  else if (mag(eVec & y) > 1-edgeTol_)
163  {
164  minY = min(minY, eMag);
165  maxY = max(maxY, eMag);
166  nY++;
167  }
168  else if (mag(eVec & z) > 1-edgeTol_)
169  {
170  minZ = min(minZ, eMag);
171  maxZ = max(maxZ, eMag);
172  nZ++;
173  }
174  else
175  {
176  minOther = min(minOther, eMag);
177  maxOther = max(maxOther, eMag);
178  }
179  }
180 
181  os << "Mesh bounding box:" << boundBox(mesh_.points()) << nl << nl
182  << "Mesh edge statistics:" << nl
183  << " x aligned : number:" << nX << "\tminLen:" << minX
184  << "\tmaxLen:" << maxX << nl
185  << " y aligned : number:" << nY << "\tminLen:" << minY
186  << "\tmaxLen:" << maxY << nl
187  << " z aligned : number:" << nZ << "\tminLen:" << minZ
188  << "\tmaxLen:" << maxZ << nl
189  << " other : number:" << mesh_.nEdges() - nX - nY - nZ
190  << "\tminLen:" << minOther
191  << "\tmaxLen:" << maxOther << nl << endl;
192 
193  if (normalDir_ == 0)
194  {
195  return min(minY, min(minZ, minOther));
196  }
197  else if (normalDir_ == 1)
198  {
199  return min(minX, min(minZ, minOther));
200  }
201  else if (normalDir_ == 2)
202  {
203  return min(minX, min(minY, minOther));
204  }
205  else
206  {
207  return min(minX, min(minY, min(minZ, minOther)));
208  }
209 }
210 
211 
212 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
213 
214 
215 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
216 
217 
218 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
219 
220 
221 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
uint8_t direction
Definition: direction.H:46
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
static const scalar edgeTol_
Definition: edgeStats.H:82
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:979
scalar y
dynamicFvMesh & mesh
List< edge > edgeList
Definition: edgeList.H:38
static const char nl
Definition: Ostream.H:262
label nEdges() const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
scalar minLen(Ostream &os) const
Calculate minimum edge length and print.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)