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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "edgeStats.H"
29#include "Time.H"
30#include "polyMesh.H"
31#include "Ostream.H"
32#include "twoDPointCorrector.H"
33#include "IOdictionary.H"
34
35// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36
37const Foam::scalar Foam::edgeStats::edgeTol_ = 1e-3;
38
39
40
41// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42
43Foam::direction Foam::edgeStats::getNormalDir
44(
45 const twoDPointCorrector* correct2DPtr
46) const
47{
48 direction dir = 3;
49
50 if (correct2DPtr)
51 {
52 const vector& normal = correct2DPtr->planeNormal();
53
54 if (mag(normal & vector(1, 0, 0)) > 1-edgeTol_)
55 {
56 dir = 0;
57 }
58 else if (mag(normal & vector(0, 1, 0)) > 1-edgeTol_)
59 {
60 dir = 1;
61 }
62 else if (mag(normal & vector(0, 0, 1)) > 1-edgeTol_)
63 {
64 dir = 2;
65 }
66 }
67 return dir;
68}
69
70
71// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
72
73Foam::edgeStats::edgeStats(const polyMesh& mesh)
74:
75 mesh_(mesh),
76 normalDir_(3)
77{
78 IOobject motionObj
79 (
80 "motionProperties",
81 mesh.time().constant(),
82 mesh,
83 IOobject::MUST_READ_IF_MODIFIED,
84 IOobject::NO_WRITE
85 );
86
87 if (motionObj.typeHeaderOk<IOdictionary>(true))
88 {
89 Info<< "Reading " << mesh.time().constant() / "motionProperties"
90 << endl << endl;
91
92 IOdictionary motionProperties(motionObj);
93
94 if (motionProperties.get<bool>("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
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
123Foam::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 for (const edge& e : edges)
147 {
148 vector eVec(e.vec(mesh_.points()));
149
150 scalar eMag = mag(eVec);
151
152 eVec /= eMag;
153
154 if (mag(eVec & x) > 1-edgeTol_)
155 {
156 minX = min(minX, eMag);
157 maxX = max(maxX, eMag);
158 nX++;
159 }
160 else if (mag(eVec & y) > 1-edgeTol_)
161 {
162 minY = min(minY, eMag);
163 maxY = max(maxY, eMag);
164 nY++;
165 }
166 else if (mag(eVec & z) > 1-edgeTol_)
167 {
168 minZ = min(minZ, eMag);
169 maxZ = max(maxZ, eMag);
170 nZ++;
171 }
172 else
173 {
174 minOther = min(minOther, eMag);
175 maxOther = max(maxOther, eMag);
176 }
177 }
178
179 os << "Mesh bounding box:" << boundBox(mesh_.points()) << nl << nl
180 << "Mesh edge statistics:" << nl
181 << " x aligned : number:" << nX << "\tminLen:" << minX
182 << "\tmaxLen:" << maxX << nl
183 << " y aligned : number:" << nY << "\tminLen:" << minY
184 << "\tmaxLen:" << maxY << nl
185 << " z aligned : number:" << nZ << "\tminLen:" << minZ
186 << "\tmaxLen:" << maxZ << nl
187 << " other : number:" << mesh_.nEdges() - nX - nY - nZ
188 << "\tminLen:" << minOther
189 << "\tmaxLen:" << maxOther << nl << endl;
190
191 if (normalDir_ == 0)
192 {
193 return min(minY, min(minZ, minOther));
194 }
195 else if (normalDir_ == 1)
196 {
197 return min(minX, min(minZ, minOther));
198 }
199 else if (normalDir_ == 2)
200 {
201 return min(minX, min(minY, minOther));
202 }
203 else
204 {
205 return min(minX, min(minY, min(minZ, minOther)));
206 }
207}
208
209
210// ************************************************************************* //
scalar y
Helper class to calculate minimum edge length on mesh.
Definition: edgeStats.H:57
static const scalar edgeTol_
Definition: edgeStats.H:84
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
label nEdges() const
Number of mesh edges.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
uint8_t direction
Definition: direction.H:56
List< edge > edgeList
A List of edges.
Definition: edgeList.H:63
Vector< scalar > vector
Definition: vector.H:61
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
volScalarField & e
Definition: createFields.H:11