septernionI.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-2019 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
29 {}
30 
32 :
33  t_(t),
34  r_(r)
35 {}
36 
38 :
39  t_(t),
40  r_(quaternion::I)
41 {}
42 
44 :
45  t_(Zero),
46  r_(r)
47 {}
48 
50 :
51  t_(st.r()),
52  r_(st.E())
53 {}
54 
55 
56 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
57 
58 inline const Foam::vector& Foam::septernion::t() const
59 {
60  return t_;
61 }
62 
63 
65 {
66  return r_;
67 }
68 
69 
71 {
72  return t_;
73 }
74 
75 
77 {
78  return r_;
79 }
80 
81 
83 {
84  return r().transform(v - t());
85 }
86 
87 
89 {
90  return t() + r().invTransform(v);
91 }
92 
93 
94 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
95 
97 {
98  t_ = tr.t() + tr.r().invTransform(t_);
99  r_ *= tr.r();
100 }
101 
102 
103 inline void Foam::septernion::operator=(const vector& t)
104 {
105  t_ = t;
106  r_ = quaternion::I;
107 }
108 
110 {
111  t_ += t;
112 }
113 
115 {
116  t_ -= t;
117 }
118 
119 
121 {
122  t_ = Zero;
123  r_ = r;
124 }
125 
127 {
128  t_ = r.invTransform(t_);
129  r_ *= r;
130 }
131 
133 {
134  t_ = r.transform(t_);
135  r_ /= r;
136 }
137 
138 
139 inline void Foam::septernion::operator*=(const scalar s)
140 {
141  t_ *= s;
142  r_ *= s;
143 }
144 
145 inline void Foam::septernion::operator/=(const scalar s)
146 {
147  t_ /= s;
148  r_ /= s;
149 }
150 
151 
152 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
153 
155 {
156  return septernion(-tr.r().transform(tr.t()), conjugate(tr.r()));
157 }
158 
159 
160 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
161 
162 inline bool Foam::operator==(const septernion& tr1, const septernion& tr2)
163 {
164  return (tr1.t() == tr2.t() && tr1.r() == tr2.r());
165 }
166 
167 
168 inline bool Foam::operator!=(const septernion& tr1, const septernion& tr2)
169 {
170  return !operator==(tr1, tr2);
171 }
172 
173 
174 inline Foam::septernion Foam::operator+
175 (
176  const septernion& tr,
177  const vector& t
178 )
179 {
180  return septernion(tr.t() + t, tr.r());
181 }
182 
183 
184 inline Foam::septernion Foam::operator+
185 (
186  const vector& t,
187  const septernion& tr
188 )
189 {
190  return septernion(t + tr.t(), tr.r());
191 }
192 
193 
194 inline Foam::septernion Foam::operator-
195 (
196  const septernion& tr,
197  const vector& t
198 )
199 {
200  return septernion(tr.t() - t, tr.r());
201 }
202 
203 
204 inline Foam::septernion Foam::operator*
205 (
206  const quaternion& r,
207  const septernion& tr
208 )
209 {
210  return septernion(tr.t(), r*tr.r());
211 }
212 
213 
214 inline Foam::septernion Foam::operator*
215 (
216  const septernion& tr,
217  const quaternion& r
218 )
219 {
220  return septernion(r.invTransform(tr.t()), tr.r()*r);
221 }
222 
223 
224 inline Foam::septernion Foam::operator/
225 (
226  const septernion& tr,
227  const quaternion& r
228 )
229 {
230  return septernion(r.transform(tr.t()), tr.r()/r);
231 }
232 
233 
234 inline Foam::septernion Foam::operator*
235 (
236  const septernion& tr1,
237  const septernion& tr2
238 )
239 {
240  return septernion
241  (
242  tr2.r().invTransform(tr1.t()) + tr2.t(),
243  tr1.r().transform(tr2.r())
244  );
245 }
246 
247 
248 inline Foam::septernion Foam::operator/
249 (
250  const septernion& tr1,
251  const septernion& tr2
252 )
253 {
254  return tr1*inv(tr2);
255 }
256 
257 
258 inline Foam::septernion Foam::operator*(const scalar s, const septernion& tr)
259 {
260  return septernion(s*tr.t(), s*tr.r());
261 }
262 
263 
264 inline Foam::septernion Foam::operator*(const septernion& tr, const scalar s)
265 {
266  return septernion(s*tr.t(), s*tr.r());
267 }
268 
269 
270 inline Foam::septernion Foam::operator/(const septernion& tr, const scalar s)
271 {
272  return septernion(tr.t()/s, tr.r()/s);
273 }
274 
275 
276 // ************************************************************************* //
tmp< fvMatrix< Type > > operator*(const volScalarField::Internal &, const fvMatrix< Type > &)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
Septernion class used to perform translations and rotations in 3D space.
Definition: septernion.H:65
tmp< fvMatrix< Type > > operator/(const fvMatrix< Type > &, const volScalarField::Internal &)
septernion()
Construct null.
Definition: septernionI.H:28
void operator*=(const septernion &)
Definition: septernionI.H:96
void operator=(const vector &)
Definition: septernionI.H:103
void operator-=(const vector &)
Definition: septernionI.H:114
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.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
quaternion conjugate(const quaternion &q)
Return the conjugate of the given quaternion.
Definition: quaternionI.H:590
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Quaternion class used to perform rotations in 3D space.
Definition: quaternion.H:60
static const zero Zero
Definition: zero.H:97
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
static const septernion I
Definition: septernion.H:83
void operator/=(const quaternion &)
Definition: septernionI.H:132
static const quaternion I
Definition: quaternion.H:118
void operator+=(const vector &)
Definition: septernionI.H:109
vector invTransformPoint(const vector &v) const
Inverse Transform the given coordinate point.
Definition: septernionI.H:88
const vector & t() const
Definition: septernionI.H:58
Compact representation of the Plücker spatial transformation tensor in terms of the rotation tensor E...
const quaternion & r() const
Definition: septernionI.H:64
vector invTransform(const vector &v) const
Rotate the given vector anti-clockwise.
Definition: quaternionI.H:302
vector transformPoint(const vector &v) const
Transform the given coordinate point.
Definition: septernionI.H:82
bool operator!=(const particle &, const particle &)
Definition: particle.C:1282
vector transform(const vector &v) const
Rotate the given vector.
Definition: quaternionI.H:296