1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright held by original author
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
29 Quaternion class used to perform rotations in 3D space.
35 \*---------------------------------------------------------------------------*/
44 #include "contiguous.H"
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 // Forward declaration of friend functions and operators
54 Istream& operator>>(Istream& is, quaternion&);
55 Ostream& operator<<(Ostream& os, const quaternion& C);
58 /*---------------------------------------------------------------------------*\
59 Class quaternion Declaration
60 \*---------------------------------------------------------------------------*/
66 //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
69 //- Vector part of the quaternion ( = axis of rotation)
73 //- Multiply vector v by quaternion as if v is a pure quaternion
74 inline quaternion mulq0v(const vector& v) const;
79 // Static data members
81 static const char* const typeName;
83 static const quaternion zero;
84 static const quaternion I;
92 //- Construct given scalar and vector parts
93 inline quaternion(const scalar w, const vector& v);
95 //- Construct a rotation quaternion given the direction d
97 inline quaternion(const vector& d, const scalar theta);
99 //- Construct given scalar part, the vector part = vector::zero
100 inline explicit quaternion(const scalar w);
102 //- Construct a pure quaternion given the vector part, scalar part = 0
103 inline explicit quaternion(const vector& v);
105 //- Construct a quaternion given the three Euler angles
113 //- Construct from Istream
114 quaternion(Istream&);
121 //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
122 inline scalar w() const;
124 //- Vector part of the quaternion ( = axis of rotation)
125 inline const vector& v() const;
127 //- The rotation tensor corresponding the quaternion
128 inline tensor R() const;
133 //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
136 //- Vector part of the quaternion ( = axis of rotation)
139 inline void normalize();
144 //- Rotate the given vector
145 inline vector transform(const vector& v) const;
147 //- Rotate the given vector anti-clockwise
148 inline vector invTransform(const vector& v) const;
150 //- Rotate the given quaternion (and normalize)
151 inline quaternion transform(const quaternion& q) const;
153 //- Rotate the given quaternion anti-clockwise (and normalize)
154 inline quaternion invTransform(const quaternion& q) const;
159 inline void operator=(const quaternion&);
160 inline void operator+=(const quaternion&);
161 inline void operator-=(const quaternion&);
162 inline void operator*=(const quaternion&);
163 inline void operator/=(const quaternion&);
165 inline void operator=(const scalar);
167 inline void operator=(const vector&);
169 inline void operator*=(const scalar);
170 inline void operator/=(const scalar);
173 // IOstream operators
175 friend Istream& operator>>(Istream& is, quaternion&);
176 friend Ostream& operator<<(Ostream& os, const quaternion& C);
180 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
182 inline scalar magSqr(const quaternion& q);
183 inline scalar mag(const quaternion& q);
185 //- Return the conjugate of the given quaternion
186 inline quaternion conjugate(const quaternion& q);
188 //- Return the normailzed (unit) quaternion of the given quaternion
189 inline quaternion normalize(const quaternion& q);
191 //- Return the inverse of the given quaternion
192 inline quaternion inv(const quaternion& q);
194 //- Return a string representation of a quaternion
195 word name(const quaternion&);
197 //- Data associated with quaternion type are contiguous
199 inline bool contiguous<quaternion>() {return true;}
202 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
204 inline bool operator==(const quaternion& q1, const quaternion& q2);
205 inline bool operator!=(const quaternion& q1, const quaternion& q2);
206 inline quaternion operator+(const quaternion& q1, const quaternion& q2);
207 inline quaternion operator-(const quaternion& q);
208 inline quaternion operator-(const quaternion& q1, const quaternion& q2);
209 inline scalar operator&(const quaternion& q1, const quaternion& q2);
210 inline quaternion operator*(const quaternion& q1, const quaternion& q2);
211 inline quaternion operator/(const quaternion& q1, const quaternion& q2);
212 inline quaternion operator*(const scalar s, const quaternion& q);
213 inline quaternion operator*(const quaternion& q, const scalar s);
214 inline quaternion operator/(const quaternion& q, const scalar s);
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
219 } // End namespace Foam
221 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
223 #include "quaternionI.H"
225 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
229 // ************************************************************************* //