vectorTools.H
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-2016 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
26Class
27 Foam::vectorTools
28
29Description
30 Functions for analysing the relationships between vectors
31
32\*---------------------------------------------------------------------------*/
33
34#ifndef vectorTools_H
35#define vectorTools_H
36
37#include "vector.H"
38#include "unitConversion.H"
39
40// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41
42namespace Foam
43{
44
45/*---------------------------------------------------------------------------*\
46 Namespace vectorTools Declaration
47\*---------------------------------------------------------------------------*/
48
49//- Collection of functions for testing relationships between two vectors.
50namespace vectorTools
51{
52 //- Test if a and b are parallel: a^b = 0
53 // Uses the cross product, so the tolerance is proportional to
54 // the sine of the angle between a and b in radians
55 template<class T>
56 bool areParallel
57 (
58 const Vector<T>& a,
59 const Vector<T>& b,
60 const T& tolerance = SMALL
61 )
62 {
63 return (mag(a ^ b) < tolerance) ? true : false;
64// return ( mag( mag(a & b)/(mag(a)*mag(b)) - 1.0 ) < tolerance )
65// ? true
66// : false;
67 }
68
69 //- Test if a and b are orthogonal: a.b = 0
70 // Uses the dot product, so the tolerance is proportional to
71 // the cosine of the angle between a and b in radians
72 template<class T>
73 bool areOrthogonal
74 (
75 const Vector<T>& a,
76 const Vector<T>& b,
77 const T& tolerance = SMALL
78 )
79 {
80 return (mag(a & b) < tolerance) ? true : false;
81 }
82
83 //- Test if angle between a and b is acute: a.b > 0
84 template<class T>
85 bool areAcute
86 (
87 const Vector<T>& a,
88 const Vector<T>& b
89 )
90 {
91 return ((a & b) > 0) ? true : false;
92 }
93
94 //- Test if angle between a and b is obtuse: a.b < 0
95 template<class T>
96 bool areObtuse
97 (
98 const Vector<T>& a,
99 const Vector<T>& b
100 )
101 {
102 return ((a & b) < 0) ? true : false;
103 }
104
105 //- Calculate angle between a and b in radians
106 template<class T>
107 T cosPhi
108 (
109 const Vector<T>& a,
110 const Vector<T>& b,
111 const T& tolerance = SMALL
112 )
113 {
114 scalar cosPhi = (a & b)/(mag(a)*mag(b) + tolerance);
115
116 // Enforce bounding between -1 and 1
117 return min(max(cosPhi, -1), 1);
118 }
119
120 //- Calculate angle between a and b in radians
121 template<class T>
123 (
124 const Vector<T>& a,
125 const Vector<T>& b,
126 const T& tolerance = SMALL
127 )
128 {
129 scalar cosPhi = (a & b)/(mag(a)*mag(b) + tolerance);
130
131 // Enforce bounding between -1 and 1
132 return acos( min(max(cosPhi, -1), 1) );
133 }
134
135 //- Calculate angle between a and b in degrees
136 template<class T>
138 (
139 const Vector<T>& a,
140 const Vector<T>& b,
141 const T& tolerance = SMALL
142 )
143 {
144 return radToDeg(radAngleBetween(a, b, tolerance));
145 }
146
147} // End namespace vectorTools
148
149
150// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
151
152} // End namespace Foam
153
154// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
155
156#endif
157
158// ************************************************************************* //
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
Definition: Vector.H:65
Functions for analysing the relationships between vectors.
const volScalarField & T
bool areOrthogonal(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Test if a and b are orthogonal: a.b = 0.
Definition: vectorTools.H:73
T degAngleBetween(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Calculate angle between a and b in degrees.
Definition: vectorTools.H:137
T cosPhi(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Calculate angle between a and b in radians.
Definition: vectorTools.H:107
bool areObtuse(const Vector< T > &a, const Vector< T > &b)
Test if angle between a and b is obtuse: a.b < 0.
Definition: vectorTools.H:96
bool areParallel(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Test if a and b are parallel: a^b = 0.
Definition: vectorTools.H:56
T radAngleBetween(const Vector< T > &a, const Vector< T > &b, const T &tolerance=SMALL)
Calculate angle between a and b in radians.
Definition: vectorTools.H:122
bool areAcute(const Vector< T > &a, const Vector< T > &b)
Test if angle between a and b is acute: a.b > 0.
Definition: vectorTools.H:85
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
constexpr scalar radToDeg() noexcept
Multiplication factor for radians to degrees conversion.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
dimensionedScalar acos(const dimensionedScalar &ds)
volScalarField & b
Definition: createFields.H:27
Unit conversion functions.