quaternion.H
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-2021 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 Class
25  Foam::quaternion
26 
27 Description
28  Quaternion class used to perform rotations in 3D space.
29 
30 SourceFiles
31  quaternionI.H
32  quaternion.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef quaternion_H
37 #define quaternion_H
38 
39 #include "scalar.H"
40 #include "vector.H"
41 #include "tensor.H"
42 #include "word.H"
43 #include "contiguous.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
50 // Forward declaration of friend functions and operators
51 
52 class quaternion;
53 Istream& operator>>(Istream& is, quaternion&);
54 Ostream& operator<<(Ostream& os, const quaternion& C);
55 
56 
57 /*---------------------------------------------------------------------------*\
58  Class quaternion Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 class quaternion
62 {
63  // private data
64 
65  //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
66  scalar w_;
67 
68  //- Vector part of the quaternion ( = axis of rotation)
69  vector v_;
70 
71 
72  //- Multiply vector v by quaternion as if v is a pure quaternion
73  inline quaternion mulq0v(const vector& v) const;
74 
75  //- Conversion of two-axis rotation components into Euler-angles
76  inline static vector twoAxes
77  (
78  const scalar t11,
79  const scalar t12,
80  const scalar c2,
81  const scalar t31,
82  const scalar t32
83  );
84 
85  //- Conversion of three-axis rotation components into Euler-angles
86  inline static vector threeAxes
87  (
88  const scalar t11,
89  const scalar t12,
90  const scalar s2,
91  const scalar t31,
92  const scalar t32
93  );
94 
95 
96 public:
97 
98  //- Component type
99  typedef scalar cmptType;
100 
101  //- Euler-angle rotation sequence
102  enum rotationSequence
103  {
105  };
106 
107 
108  // Member constants
109 
110  //- Rank of quaternion is 1
111  static const direction rank = 1;
112 
113 
114  // Static Data Members
115 
116  static const char* const typeName;
117 
118  static const quaternion zero;
119  static const quaternion I;
120 
121 
122  // Constructors
123 
124  //- Construct null
125  inline quaternion();
126 
127  //- Construct given scalar and vector parts
128  inline quaternion(const scalar w, const vector& v);
129 
130  //- Construct a rotation quaternion given the direction d
131  // and angle theta
132  inline quaternion(const vector& d, const scalar theta);
133 
134  //- Construct a rotation quaternion given the direction d
135  // and cosine angle cosTheta and a if d is normalised
136  inline quaternion
137  (
138  const vector& d,
139  const scalar cosTheta,
140  const bool normalised
141  );
142 
143  //- Construct a real from the given scalar part, the vector part = zero
144  inline explicit quaternion(const scalar w);
145 
146  //- Construct a pure imaginary quaternion given the vector part,
147  // the scalar part = 0
148  inline explicit quaternion(const vector& v);
149 
150  //- Return the unit quaternion (versor) from the given vector
151  // (w = sqrt(1 - |sqr(v)|))
152  static inline quaternion unit(const vector& v);
153 
154  //- Construct a quaternion given the three Euler angles
155  inline quaternion
156  (
157  const rotationSequence rs,
158  const vector& angles
159  );
160 
161  //- Construct a quaternion from a rotation tensor
162  inline explicit quaternion(const tensor& rotationTensor);
163 
164  //- Construct from Istream
166 
167 
168  // Member Functions
169 
170  // Access
171 
172  //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
173  inline scalar w() const;
174 
175  //- Vector part of the quaternion ( = axis of rotation)
176  inline const vector& v() const;
177 
178  //- The rotation tensor corresponding the quaternion
179  inline tensor R() const;
180 
181  //- Return a vector of euler angles corresponding to the
182  // specified rotation sequence
183  inline vector eulerAngles(const rotationSequence rs) const;
184 
185  inline quaternion normalised() const;
186 
187 
188  // Edit
189 
190  //- Scalar part of the quaternion ( = cos(theta/2) for rotation)
191  inline scalar& w();
192 
193  //- Vector part of the quaternion ( = axis of rotation)
194  inline vector& v();
195 
196  inline void normalise();
197 
198 
199  // Transform
200 
201  //- Rotate the given vector
202  inline vector transform(const vector& v) const;
203 
204  //- Rotate the given vector anti-clockwise
205  inline vector invTransform(const vector& v) const;
206 
207  //- Rotate the given quaternion (and normalise)
208  inline quaternion transform(const quaternion& q) const;
209 
210  //- Rotate the given quaternion anti-clockwise (and normalise)
211  inline quaternion invTransform(const quaternion& q) const;
212 
213 
214  // Member Operators
215 
216  inline void operator+=(const quaternion&);
217  inline void operator-=(const quaternion&);
218  inline void operator*=(const quaternion&);
219  inline void operator/=(const quaternion&);
220 
221  inline void operator=(const scalar);
222 
223  inline void operator=(const vector&);
224 
225  inline void operator*=(const scalar);
226  inline void operator/=(const scalar);
227 
228 
229  // IOstream Operators
230 
232  friend Ostream& operator<<(Ostream& os, const quaternion& C);
233 };
234 
235 
236 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
237 
238 inline scalar magSqr(const quaternion& q);
239 inline scalar mag(const quaternion& q);
240 
241 //- Return the conjugate of the given quaternion
242 inline quaternion conjugate(const quaternion& q);
243 
244 //- Return the normalised (unit) quaternion of the given quaternion
245 inline quaternion normalise(const quaternion& q);
246 
247 //- Return the inverse of the given quaternion
248 inline quaternion inv(const quaternion& q);
249 
250 //- Return a string representation of a quaternion
251 word name(const quaternion&);
252 
253 //- Spherical linear interpolation of quaternions
255 (
256  const quaternion& qa,
257  const quaternion& qb,
258  const scalar t
259 );
260 
261 //- Simple weighted average with sign change
263 (
264  const UList<quaternion>& qs,
265  const UList<scalar> w
266 );
267 
268 //- Exponent of a quaternion
269 quaternion exp(const quaternion& q);
270 
271 //- Power of a quaternion
272 quaternion pow(const quaternion& q, const label power);
273 
274 //- Power of a quaternion
275 quaternion pow(const quaternion& q, const scalar power);
276 
277 //- Data associated with quaternion type are contiguous
278 template<>
279 inline bool contiguous<quaternion>() {return true;}
280 
281 
282 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
283 
284 inline bool operator==(const quaternion& q1, const quaternion& q2);
285 inline bool operator!=(const quaternion& q1, const quaternion& q2);
286 inline quaternion operator+(const quaternion& q1, const quaternion& q2);
287 inline quaternion operator-(const quaternion& q);
288 inline quaternion operator-(const quaternion& q1, const quaternion& q2);
289 inline scalar operator&(const quaternion& q1, const quaternion& q2);
290 inline quaternion operator*(const quaternion& q1, const quaternion& q2);
291 inline quaternion operator/(const quaternion& q1, const quaternion& q2);
292 inline quaternion operator*(const scalar s, const quaternion& q);
293 inline quaternion operator*(const quaternion& q, const scalar s);
294 inline quaternion operator/(const quaternion& q, const scalar s);
295 
296 
297 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
298 
299 } // End namespace Foam
300 
301 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
302 
303 #include "quaternionI.H"
304 
305 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
306 
307 #endif
308 
309 // ************************************************************************* //
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
Graphite solid properties.
Definition: C.H:51
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:74
Quaternion class used to perform rotations in 3D space.
Definition: quaternion.H:61
const vector & v() const
Vector part of the quaternion ( = axis of rotation)
Definition: quaternionI.H:260
quaternion normalised() const
Definition: quaternionI.H:278
quaternion()
Construct null.
Definition: quaternionI.H:28
vector eulerAngles(const rotationSequence rs) const
Return a vector of euler angles corresponding to the.
Definition: quaternionI.H:373
friend Ostream & operator<<(Ostream &os, const quaternion &C)
void operator-=(const quaternion &)
Definition: quaternionI.H:532
static const direction rank
Rank of quaternion is 1.
Definition: quaternion.H:110
rotationSequence
Euler-angle rotation sequence.
Definition: quaternion.H:102
tensor R() const
The rotation tensor corresponding the quaternion.
Definition: quaternionI.H:323
void operator+=(const quaternion &)
Definition: quaternionI.H:526
scalar w() const
Scalar part of the quaternion ( = cos(theta/2) for rotation)
Definition: quaternionI.H:254
static quaternion unit(const vector &v)
Return the unit quaternion (versor) from the given vector.
Definition: quaternionI.H:81
void operator*=(const quaternion &)
Definition: quaternionI.H:538
void operator=(const scalar)
Definition: quaternionI.H:551
vector transform(const vector &v) const
Rotate the given vector.
Definition: quaternionI.H:296
scalar cmptType
Component type.
Definition: quaternion.H:98
static const quaternion I
Definition: quaternion.H:118
friend Istream & operator>>(Istream &is, quaternion &)
static const quaternion zero
Definition: quaternion.H:117
vector invTransform(const vector &v) const
Rotate the given vector anti-clockwise.
Definition: quaternionI.H:302
static const char *const typeName
Definition: quaternion.H:115
void operator/=(const quaternion &)
Definition: quaternionI.H:545
A class for handling words, derived from string.
Definition: word.H:62
Template function to specify if the data of a type are contiguous.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionedScalar c2
Second radiation constant: default SI units: [m K].
Namespace for OpenFOAM.
bool operator!=(const particle &, const particle &)
Definition: particle.C:1210
dimensionedScalar exp(const dimensionedScalar &ds)
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from unit vector n1 to n2.
Definition: transform.H:47
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
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
bool contiguous< quaternion >()
Data associated with quaternion type are contiguous.
Definition: quaternion.H:278
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Istream & operator>>(Istream &, pistonPointEdgeData &)
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< VolField< Type > > operator&(const fvMatrix< Type > &, const DimensionedField< Type, volMesh > &)
tmp< fvMatrix< Type > > operator*(const volScalarField::Internal &, const fvMatrix< Type > &)
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
tmp< fvMatrix< Type > > operator+(const fvMatrix< Type > &, const fvMatrix< Type > &)
quaternion conjugate(const quaternion &q)
Return the conjugate of the given quaternion.
Definition: quaternionI.H:590
quaternion slerp(const quaternion &qa, const quaternion &qb, const scalar t)
Spherical linear interpolation of quaternions.
Definition: quaternion.C:55
quaternion normalise(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:603
tmp< fvMatrix< Type > > operator-(const fvMatrix< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
uint8_t direction
Definition: direction.H:45
tmp< fvMatrix< Type > > operator/(const fvMatrix< Type > &, const volScalarField::Internal &)