arcEdge.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-2015 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 "arcEdge.H"
27 #include "unitConversion.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(arcEdge, 0);
35  addToRunTimeSelectionTable(curvedEdge, arcEdge, Istream);
36 }
37 
38 
39 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
40 
41 Foam::cylindricalCS Foam::arcEdge::calcAngle()
42 {
43  vector a = p2_ - p1_;
44  vector b = p3_ - p1_;
45 
46  // find centre of arcEdge
47  scalar asqr = a & a;
48  scalar bsqr = b & b;
49  scalar adotb = a & b;
50 
51  scalar denom = asqr*bsqr - adotb*adotb;
52 
53  if (mag(denom) < VSMALL)
54  {
55  FatalErrorIn("cylindricalCS arcEdge::calcAngle()")
56  << "Invalid arc definition - are the points co-linear? Denom ="
57  << denom
58  << abort(FatalError);
59  }
60 
61  scalar fact = 0.5*(bsqr - adotb)/denom;
62 
63  point centre = 0.5*a + fact*((a ^ b) ^ a);
64 
65  centre += p1_;
66 
67  // find position vectors w.r.t. the arcEdge centre
68  vector r1(p1_ - centre);
69  vector r2(p2_ - centre);
70  vector r3(p3_ - centre);
71 
72  // find angles
73  angle_ = radToDeg(acos((r3 & r1)/(mag(r3) * mag(r1))));
74 
75  // check if the vectors define an exterior or an interior arcEdge
76  if (((r1 ^ r2) & (r1 ^ r3)) < 0.0)
77  {
78  angle_ = 360.0 - angle_;
79  }
80 
81  vector tempAxis;
82 
83  if (angle_ <= 180.0)
84  {
85  tempAxis = r1 ^ r3;
86 
87  if (mag(tempAxis)/(mag(r1)*mag(r3)) < 0.001)
88  {
89  tempAxis = r1 ^ r2;
90  }
91  }
92  else
93  {
94  tempAxis = r3 ^ r1;
95  }
96 
97  radius_ = mag(r3);
98 
99  // set up and return the local coordinate system
100  return cylindricalCS("arcEdgeCS", centre, tempAxis, r1);
101 }
102 
103 
104 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
105 
106 Foam::arcEdge::arcEdge
107 (
108  const pointField& points,
109  const label start,
110  const label end,
111  const point& pMid
112 )
113 :
114  curvedEdge(points, start, end),
115  p1_(points_[start_]),
116  p2_(pMid),
117  p3_(points_[end_]),
118  cs_(calcAngle())
119 {}
120 
121 
122 Foam::arcEdge::arcEdge(const pointField& points, Istream& is)
123 :
124  curvedEdge(points, is),
125  p1_(points_[start_]),
126  p2_(is),
127  p3_(points_[end_]),
128  cs_(calcAngle())
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
133 
135 {
136  if (lambda < -SMALL || lambda > 1 + SMALL)
137  {
138  FatalErrorIn("arcEdge::position(const scalar lambda) const")
139  << "Parameter out of range, lambda = " << lambda
140  << abort(FatalError);
141  }
142 
143  if (lambda < SMALL)
144  {
145  return p1_;
146  }
147  else if (lambda > 1 - SMALL)
148  {
149  return p3_;
150  }
151  else
152  {
153  return cs_.globalPosition(vector(radius_, lambda*angle_, 0.0));
154  }
155 }
156 
157 
158 Foam::scalar Foam::arcEdge::length() const
159 {
160  return degToRad(angle_*radius_);
161 }
162 
163 
164 // ************************************************************************* //
vector point
Point is a vector.
Definition: point.H:41
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
Define a curved edge that is parameterized for 0<lambda<1 between the start and end point...
Definition: curvedEdge.H:53
dimensioned< scalar > mag(const dimensioned< Type > &)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
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
Namespace for OpenFOAM.
scalar length() const
Return the length of the curve.
Definition: arcEdge.C:158
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
scalar radToDeg(const scalar rad)
Conversion from radians to degrees.
Macros for easy insertion into run-time selection tables.
Unit conversion functions.
dimensionedScalar acos(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
point position(const scalar) const
Return the point position corresponding to the curve parameter.
Definition: arcEdge.C:134
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
Cylindrical coordinate system.
Definition: cylindricalCS.H:48
point globalPosition(const point &local) const
Convert from position in local coordinate system to global.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)