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