55 const scalar cosTheta,
56 const bool isNormalised
59 const scalar cosHalfTheta2 = 0.5*(cosTheta + 1);
60 w_ =
sqrt(cosHalfTheta2);
64 v_ =
sqrt(1 - cosHalfTheta2)*d;
95 const eulerOrder order,
111 operator=(quaternion(
vector(0, 0, 1), angles.x()));
112 operator*=(quaternion(
vector(0, 1, 0), angles.y()));
113 operator*=(quaternion(
vector(0, 0, 1), angles.z()));
119 operator=(quaternion(
vector(0, 0, 1), angles.x()));
120 operator*=(quaternion(
vector(1, 0, 0), angles.y()));
121 operator*=(quaternion(
vector(0, 1, 0), angles.z()));
127 operator=(quaternion(
vector(0, 0, 1), angles.x()));
128 operator*=(quaternion(
vector(1, 0, 0), angles.y()));
129 operator*=(quaternion(
vector(0, 0, 1), angles.z()));
135 operator=(quaternion(
vector(0, 1, 0), angles.x()));
136 operator*=(quaternion(
vector(1, 0, 0), angles.y()));
137 operator*=(quaternion(
vector(0, 0, 1), angles.z()));
143 operator=(quaternion(
vector(0, 1, 0), angles.x()));
144 operator*=(quaternion(
vector(1, 0, 0), angles.y()));
145 operator*=(quaternion(
vector(0, 1, 0), angles.z()));
151 operator=(quaternion(
vector(0, 1, 0), angles.x()));
152 operator*=(quaternion(
vector(0, 0, 1), angles.y()));
153 operator*=(quaternion(
vector(1, 0, 0), angles.z()));
159 operator=(quaternion(
vector(0, 1, 0), angles.x()));
160 operator*=(quaternion(
vector(0, 0, 1), angles.y()));
161 operator*=(quaternion(
vector(0, 1, 0), angles.z()));
167 operator=(quaternion(
vector(1, 0, 0), angles.x()));
168 operator*=(quaternion(
vector(0, 1, 0), angles.y()));
169 operator*=(quaternion(
vector(0, 0, 1), angles.z()));
175 operator=(quaternion(
vector(1, 0, 0), angles.x()));
176 operator*=(quaternion(
vector(0, 1, 0), angles.y()));
177 operator*=(quaternion(
vector(1, 0, 0), angles.z()));
183 operator=(quaternion(
vector(1, 0, 0), angles.x()));
184 operator*=(quaternion(
vector(0, 0, 1), angles.y()));
185 operator*=(quaternion(
vector(0, 1, 0), angles.z()));
191 operator=(quaternion(
vector(1, 0, 0), angles.x()));
192 operator*=(quaternion(
vector(0, 0, 1), angles.y()));
193 operator*=(quaternion(
vector(1, 0, 0), angles.z()));
199 <<
"Unknown euler rotation order " 327 return quaternion(-(v() & u), w()*u + (v() ^ u));
339 return (
conjugate(*this).mulq0v(u)*(*this)).v();
360 const scalar
w2 =
sqr(w());
361 const scalar x2 =
sqr(v().
x());
362 const scalar y2 =
sqr(v().
y());
363 const scalar z2 =
sqr(v().z());
365 const scalar txy = 2*v().x()*v().y();
366 const scalar twz = 2*w()*v().z();
367 const scalar txz = 2*v().x()*v().z();
368 const scalar twy = 2*w()*v().y();
369 const scalar tyz = 2*v().y()*v().z();
370 const scalar twx = 2*w()*v().x();
374 w2 + x2 - y2 - z2, txy - twz, txz + twy,
375 txy + twz,
w2 - x2 + y2 - z2, tyz - twx,
376 txz - twy, tyz + twx,
w2 - x2 - y2 + z2
409 const eulerOrder order
412 const scalar
w2 =
sqr(w());
413 const scalar x2 =
sqr(v().
x());
414 const scalar y2 =
sqr(v().
y());
415 const scalar z2 =
sqr(v().z());
423 2*(v().
x()*v().
y() + w()*v().z()),
425 2*(w()*v().
y() - v().
x()*v().z()),
426 2*(v().
y()*v().z() + w()*v().
x()),
436 2*(v().
y()*v().z() - w()*v().
x()),
437 2*(v().
x()*v().z() + w()*v().
y()),
439 2*(v().
y()*v().z() + w()*v().
x()),
440 2*(w()*v().
y() - v().
x()*v().z())
449 2*(w()*v().z() - v().
x()*v().
y()),
451 2*(v().
y()*v().z() + w()*v().
x()),
452 2*(w()*v().
y() - v().
x()*v().z()),
462 2*(v().
x()*v().z() + w()*v().
y()),
463 2*(w()*v().
x() - v().
y()*v().z()),
465 2*(v().
x()*v().z() - w()*v().
y()),
466 2*(v().
y()*v().z() + w()*v().
x())
475 2*(v().
x()*v().z() + w()*v().
y()),
477 2*(w()*v().
x() - v().
y()*v().z()),
478 2*(v().
x()*v().
y() + w()*v().z()),
488 2*(v().
x()*v().
y() - w()*v().z()),
489 2*(v().
y()*v().z() + w()*v().
x()),
491 2*(v().
x()*v().
y() + w()*v().z()),
492 2*(w()*v().
x() - v().
y()*v().z())
501 2*(w()*v().
y() - v().
x()*v().z()),
503 2*(v().
x()*v().
y() + w()*v().z()),
504 2*(w()*v().
x() - v().
y()*v().z()),
514 2*(v().
y()*v().z() + w()*v().
x()),
515 2*(w()*v().z() - v().
x()*v().
y()),
517 2*(v().
y()*v().z() - w()*v().
x()),
518 2*(v().
x()*v().
y() + w()*v().z())
527 2*(w()*v().
x() - v().
y()*v().z()),
529 2*(v().
x()*v().z() + w()*v().
y()),
530 2*(w()*v().z() - v().
x()*v().
y()),
540 2*(v().
x()*v().
y() + w()*v().z()),
541 2*(w()*v().
y() - v().
x()*v().z()),
543 2*(v().
x()*v().
y() - w()*v().z()),
544 2*(v().
x()*v().z() + w()*v().
y())
553 2*(v().
y()*v().z() + w()*v().
x()),
555 2*(w()*v().z() - v().
x()*v().
y()),
556 2*(v().
x()*v().z() + w()*v().
y()),
566 2*(v().
x()*v().z() - w()*v().
y()),
567 2*(v().
x()*v().
y() + w()*v().z()),
569 2*(v().
x()*v().z() + w()*v().
y()),
570 2*(w()*v().z() - v().
x()*v().
y())
577 <<
"Unknown euler rotation order " 603 w() = w()*q.
w() - (v() & q.
v());
604 v() =
w0*q.
v() + q.
w()*v() + (v() ^ q.
v());
609 return operator*=(
inv(q));
711 const quaternion& q1,
715 return quaternion(q1.w() + q2.w(), q1.v() + q2.v());
727 const quaternion& q1,
731 return quaternion(q1.w() - q2.w(), q1.v() - q2.v());
737 return q1.
w()*q2.
w() + (q1.
v() & q2.
v());
743 const quaternion& q1,
749 q1.w()*q2.w() - (q1.v() & q2.v()),
750 q1.w()*q2.v() + q2.w()*q1.v() + (q1.v() ^ q2.v())
757 const quaternion& q1,
779 return quaternion(q.w()/
s, q.v()/
s);
dimensionedScalar acos(const dimensionedScalar &ds)
tensor R() const
The rotation tensor corresponding to the quaternion.
void operator/=(const quaternion &q)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
vector eulerAngles(const eulerOrder order) const
Return the Euler rotation angles corresponding to the specified rotation order.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
bool equal(const T &a, const T &b)
Compare two values for equality.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const vector & v() const noexcept
Vector part of the quaternion ( = axis of rotation)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
static quaternion unit(const vector &v)
Return the unit quaternion (versor) from the given vector (w = sqrt(1 - |sqr(v)|)) ...
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
dimensionedScalar asin(const dimensionedScalar &ds)
void operator*=(const quaternion &q)
scalar w() const noexcept
Scalar part of the quaternion ( = cos(theta/2) for rotation)
quaternion conjugate(const quaternion &q)
Return the conjugate of the given quaternion.
dimensionedScalar cos(const dimensionedScalar &ds)
tmp< faMatrix< Type > > operator*(const areaScalarField::Internal &, const faMatrix< Type > &)
quaternion & normalise()
Inplace normalise the quaternion by its magnitude.
Quaternion class used to perform rotations in 3D space.
errorManip< error > abort(error &err)
void operator-=(const quaternion &q)
tmp< faMatrix< Type > > operator-(const faMatrix< Type > &)
Unary negation.
dimensionedScalar sin(const dimensionedScalar &ds)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
quaternion()=default
Default construct.
quaternion & operator=(const quaternion &)=default
Copy assignment.
tmp< GeometricField< Type, faPatchField, areaMesh > > operator &(const faMatrix< Type > &, const DimensionedField< Type, areaMesh > &)
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
bool operator!=(const eddy &a, const eddy &b)
vector invTransform(const vector &v) const
Rotate the given vector anti-clockwise.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
Tensor of scalars, i.e. Tensor<scalar>.
vector transform(const vector &v) const
Rotate the given vector.
void operator+=(const quaternion &q)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)