QuaternionNumber：比 GLM 代码可读性更强、更有学习价值的四元数类
Andy16673328152 · · 科技·工程
前言
你也许比较讨厌 GLM 的四元数，因为它不仅不是标准库实现，使用时还要包含三个头文件，而且精度不高。为了弥补这些缺点，我编写出了这个类。这个类不仅有题目所说的优点，还更为轻量。
代码
下载链接：link
云剪贴板地址：https://www.luogu.com.cn/paste/6xjdxuzb
:::info[也可以在这里查看]
#include <cmath>
#include <iostream>
#include <string>
#include <sstream>
#include <vector>
#include <algorithm>
#include <type_traits>
#include <tuple>
#include <limits>
#include <array>
using namespace std;
#ifndef QUATERNION_HPP
#define QUATERNION_HPP
string TrimTrailingZeros(string a){
while(a[a.size()-1]=='0'){
a.erase(a.size()-1,1);
}
return a;
}
template <typename T>
struct Quaternion{
static_assert(is_same<T,float>::value||is_same<T,double>::value||is_same<T,long double>::value,"T must be 'float','double' or 'long double',not 'int'");
T r,i,j,k;
Quaternion conj()const{
return {r,-i,-j,-k};
}
constexpr T real()const{
return r;
}
void real(T val){
r=val;
}
constexpr T imagi()const{
return i;
}
void imagi(T val){
i=val;
}
constexpr T imagj()const{
return j;
}
void imagj(T val){
j=val;
}
constexpr T imagk()const{
return k;
}
void imagk(T val){
k=val;
}
constexpr Quaternion():r(0.0),i(0.0),j(0.0),k(0.0){}
constexpr Quaternion(T a,T b,T c,T d):r(a),i(b),j(c),k(d){}
constexpr Quaternion(T a):r(a){}
template <typename U>
Quaternion(const Quaternion<U> &a):r(a.r),i(a.i),j(a.j),k(a.k){}
Quaternion operator =(Quaternion a){
r=a.r;
i=a.i;
j=a.j;
k=a.k;
return *this;
}
template <typename U>
Quaternion operator =(const Quaternion<U> &a){
r=a.r;
i=a.i;
j=a.j;
k=a.k;
return *this;
}
Quaternion& operator +=(T b){
r+=b;
return *this;
}
Quaternion& operator -=(T b){
r-=b;
return *this;
}
Quaternion& operator *=(T b){
r*=b;
i*=b;
j*=b;
k*=b;
return *this;
}
Quaternion& operator /=(T b){
r/=b;
i/=b;
j/=b;
k/=b;
return *this;
}
bool is_real()const{
return (abs(i)<1e-10)&&(abs(j)<1e-10)&&(abs(k)<1e-10);
}
bool is_imaginary()const{
return (abs(r)<1e-10);
}
template<size_t N>
decltype(auto) get()const{
static_assert(N<4,"Index out of range");
if(N==0){
return r;
}
if(N==1){
return i;
}
if(N==2){
return j;
}
if(N==3){
return k;
}
}
static Quaternion identity(){
return Quaternion(1.0);
}
static Quaternion unit(T r,T i,T j,T k){
Quaternion q(r,i,j,k);
return q.normalized();
}
Quaternion& normalize(){
T al=sqrt(r*r+i*i+j*j+k*k);
if(al>1e-10){
r/=al;
i/=al;
j/=al;
k/=al;
}
return *this;
}
Quaternion normalized()const{
Quaternion re(*this);
return re.normalize();
}
bool isunit(T eps=1e-10)const{
T alsq=r*r+i*i+j*j+k*k;
return abs(alsq-1)<eps;
}
static Quaternion FromAxisAngle(T angle,T x,T y,T z){
return Quaternion(cos(angle/2),x*sin(angle/2),y*sin(angle/2),z*sin(angle/2));
}
array<T,9> ToRotationMatrix()const{
T xx=r*r,yy=i*i,zz=j*j,ww=k*k,xy=r*i,xz=r*j,xw=r*k,yz=i*j,yw=i*k,zw=j*k;
return {
1-2*(yy+zz),2*(xy-zw),2*(xz+yw),
2*(xy+zw),1-2*(xx+zz),2*(yz-xw),
2*(xz-yw),2*(yz+xw),1-2*(xx+yy)
};
}
static Quaternion FromRotationMatrix(const array<T,9>& mat){
T t=mat[0]+mat[4]+mat[8];
Quaternion q;
if(t>0){
T s=0.5*(sqrt(t+1));
q.r=0.25/s;
q.i=(mat[5]-mat[7])*s;
q.j=(mat[6]-mat[2])*s;
q.k=(mat[1]-mat[3])*s;
}else{
if(mat[0]>mat[4]&&mat[0]>mat[8]){
T s=2*sqrt(1+mat[0]-mat[4]-mat[8]);
q.r=(mat[5]-mat[7])/s;
q.i=0.25*s;
q.j=(mat[3]+mat[1])/s;
q.k=(mat[6]+mat[2])/s;
}else{
if(mat[4]>mat[8]){
T s=2*sqrt(1+mat[4]-mat[0]-mat[8]);
q.r=(mat[6]-mat[2])/s;
q.i=(mat[3]+mat[1])/s;
q.j=0.25*s;
q.k=(mat[7]+mat[5])/s;
}else{
T s=2*sqrt(1+mat[8]-mat[4]-mat[0]);
q.r=(mat[1]-mat[3])/s;
q.i=(mat[6]+mat[2])/s;
q.j=(mat[7]+mat[5])/s;
q.k=0.25*s;
}
}
}
return q;
}
enum class EulerOrder{
XYZ,XZY,YXZ,YZX,ZXY,ZYX
};
array<T,3> ToEulerAngles(EulerOrder order=EulerOrder::XYZ)const{
Quaternion q=this->normalized();
T rr=q.r,ii=q.i,jj=q.j,kk=q.k;
T roll,pitch,yaw;
switch(order){
case EulerOrder::XYZ:
pitch=asin(2*(rr*jj-kk*ii));
if(abs(pitch-T(3.1415926535/2))<1e-10){
roll=0;
yaw=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(pitch+T(3.1415926535/2))<1e-10){
roll=0;
yaw=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
roll=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
yaw=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
case EulerOrder::XZY:
yaw=asin(2*(rr*jj-kk*ii));
if(abs(yaw-T(3.1415926535/2))<1e-10){
roll=0;
pitch=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(yaw+T(3.1415926535/2))<1e-10){
roll=0;
pitch=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
roll=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
pitch=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
case EulerOrder::YXZ:
roll=asin(2*(rr*jj-kk*ii));
if(abs(roll-T(3.1415926535/2))<1e-10){
pitch=0;
yaw=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(roll+T(3.1415926535/2))<1e-10){
pitch=0;
yaw=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
pitch=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
yaw=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
case EulerOrder::YZX:
yaw=asin(2*(rr*jj-kk*ii));
if(abs(yaw-T(3.1415926535/2))<1e-10){
pitch=0;
roll=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(yaw+T(3.1415926535/2))<1e-10){
pitch=0;
roll=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
pitch=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
roll=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
case EulerOrder::ZXY:
roll=asin(2*(rr*jj-kk*ii));
if(abs(roll-T(3.1415926535/2))<1e-10){
yaw=0;
pitch=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(roll+T(3.1415926535/2))<1e-10){
yaw=0;
pitch=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
yaw=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
pitch=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
case EulerOrder::ZYX:
pitch=asin(2*(rr*jj-kk*ii));
if(abs(pitch-T(3.1415926535/2))<1e-10){
yaw=0;
roll=atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
if(abs(pitch+T(3.1415926535/2))<1e-10){
yaw=0;
roll=-atan2(2*(rr*kk-ii*jj),1-2*(jj*jj+kk*kk));
}else{
yaw=atan2(2*(rr*ii+jj*kk),1-2*(ii*ii+jj*jj));
roll=atan2(2*(rr*kk+ii*jj),1-2*(jj*jj+kk*kk));
}
}
return {roll,pitch,yaw};
default:
return {0,0,0};
}
}
static Quaternion FromEulerAngles(T roll,T pitch,T yaw,EulerOrder order=EulerOrder::XYZ){
T c1=cos(roll/2);
T s1=sin(roll/2);
T c2=cos(pitch/2);
T s2=sin(pitch/2);
T c3=cos(yaw/2);
T s3=sin(yaw/2);
switch(order){
case EulerOrder::XYZ:
Quaternion q(
c1*c2*c3-s1*s2*s3,
s1*c2*c3+c1*s2*s3,
c1*s2*c3-s1*c2*s3,
c1*c2*s3+s1*s2*c3
);
return q.normalized();
case EulerOrder::XZY:
Quaternion q1(
c1*c2*c3+s1*s2*s3,
s1*c2*c3-c1*s2*s3,
c1*c2*s3-s1*s2*c3,
c1*s2*c3+s1*c2*s3
);
return q1.normalized();
case EulerOrder::YXZ:
Quaternion q2(
c1*c2*c3+s1*s2*s3,
s1*c2*c3+c1*s2*s3,
c1*s2*c3-s1*c2*s3,
c1*c2*s3-s1*s2*c3
);
return q2.normalized();
case EulerOrder::YZX:
Quaternion q3(
c1*c2*c3-s1*s2*s3,
c1*s2*c3+s1*c2*s3,
s1*c2*c3+c1*s2*s3,
c1*c2*s3-s1*s2*c3
);
return q3.normalized();
case EulerOrder::ZXY:
Quaternion q4(
c1*c2*c3-s1*s2*s3,
c1*s2*c3-s1*c2*s3,
c1*c2*s3+s1*s2*c3,
s1*c2*c3+c1*s2*s3
);
return q4.normalized();
case EulerOrder::ZYX:
Quaternion q5(
c1*c2*c3+s1*s2*s3,
c1*c2*s3-s1*s2*c3,
c1*s2*c3+s1*c2*s3,
s1*c2*c3-c1*s2*s3
);
return q5.normalized();
default:
return Quaternion(1,0,0,0);
}
}
array<T,3> RotateVector(T x,T y,T z)const{
T rx=2*(j*z-k*y),ry=2*(k*x-i*z),rz=2*(i*y-j*x);
return {x+r*rx+(j*rz-k*ry),y+r*ry+(k*rx-i*rz),z+r*rz+(i*ry-j*rx)};
}
array<T,3> RotateVector(const array<T,3> &v)const{
return RotateVector(v[0],v[1],v[2]);
}
array<T,3> axis()const{
T a=sqrt(1-r*r);
if(a<1e-10){
return {1,0,0};
}
return {i/a,j/a,k/a};
}
T angle()const{
return 2*acos(r);
}
T angleDegrees()const{
return angle()*180/3.14159265358979323846;
}
Quaternion inverse()const{
return (this->conj())/(r*r+i*i+j*j+k*k);
}
};
namespace std{
template<typename T>
struct tuple_size<Quaternion<T> >:integral_constant<size_t,4>{};
template<typename T>
struct tuple_element<0,Quaternion<T> >{using type=T;};
template<typename T>
struct tuple_element<1,Quaternion<T> >{using type=T;};
template<typename T>
struct tuple_element<2,Quaternion<T> >{using type=T;};
template<typename T>
struct tuple_element<3,Quaternion<T> >{using type=T;};
}
template<typename T>
Quaternion<T> operator +(Quaternion<T> a,T b){
a.r+=b;
return a;
}
template<typename T>
Quaternion<T> operator +(T a,Quaternion<T> b){
b.r+=a;
return b;
}
template<typename T>
Quaternion<T> operator -(Quaternion<T> a,T b){
a.r-=b;
return a;
}
template<typename T>
Quaternion<T> operator -(Quaternion<T> a){
return {-a.r,-a.i,-a.j,-a.k};
}
template<typename T>
Quaternion<T> operator -(T a,Quaternion<T> b){
b.r-=a;
b=-b;
return b;
}
template<typename T>
Quaternion<T> operator *(Quaternion<T> a,T b){
a.r*=b;
a.i*=b;
a.j*=b;
a.k*=b;
return a;
}
template<typename T>
Quaternion<T> operator *(T a,Quaternion<T> b){
b.r*=a;
b.i*=a;
b.j*=a;
b.k*=a;
return b;
}
template<typename T>
Quaternion<T> operator /(Quaternion<T> a,T b){
a.r/=b;
a.i/=b;
a.j/=b;
a.k/=b;
return a;
}
template<typename T>
T abs(Quaternion<T> a){
return sqrt(a.r*a.r+a.i*a.i+a.j*a.j+a.k*a.k);
}
template<typename T>
Quaternion<T> operator /(T a,Quaternion<T> b){
return a*(b.conj()/(abs(b)*abs(b)));
}
template<typename T>
Quaternion<T> operator +(Quaternion<T> a,Quaternion<T> b){
a.r+=b.r;
a.i+=b.i;
a.j+=b.j;
a.k+=b.k;
return a;
}
template<typename T>
Quaternion<T> operator -(Quaternion<T> a,Quaternion<T> b){
a.r-=b.r;
a.i-=b.i;
a.j-=b.j;
a.k-=b.k;
return a;
}
template<typename T>
Quaternion<T> operator *(Quaternion<T> a,Quaternion<T> b){
return Quaternion<T>(a.r*b.r-a.i*b.i-a.j*b.j-a.k*b.k,a.r*b.i+a.i*b.r+a.j*b.k-a.k*b.j,a.r*b.j+a.j*b.r+a.k*b.i-a.i*b.k,a.r*b.k+a.k*b.r+a.i*b.j-a.j*b.i);
}
template<typename T>
Quaternion<T> operator /(Quaternion<T> a,Quaternion<T> b){
return a*(1.0/b);
}
template<typename T>
bool operator ==(Quaternion<T> a,Quaternion<T> b){
if(isnan(a.r)||isnan(a.i)||isnan(a.j)||isnan(a.k)||isnan(b.r)||isnan(b.i)||isnan(b.j)||isnan(b.k)){
return false;
}
if(isinf(a.r)||isinf(a.i)||isinf(a.j)||isinf(a.k)||isinf(b.r)||isinf(b.i)||isinf(b.j)||isinf(b.k)){
return (a.r==b.r)&&(a.i==b.i)&&(a.j==b.j)&&(a.k==b.k);
}
return (abs(a.r-b.r)<1e-10)&&(abs(a.i-b.i)<1e-10)&&(abs(a.j-b.j)<1e-10)&&(abs(a.k-b.k)<1e-10);
}
template<typename T>
bool operator !=(Quaternion<T> a,Quaternion<T> b){
return !(a==b);
}
template<typename T>
Quaternion<T>& operator +=(Quaternion<T> &a,const Quaternion<T> &b){
a=a+b;
return a;
}
template<typename T>
Quaternion<T>& operator -=(Quaternion<T> &a,const Quaternion<T> &b){
a=a-b;
return a;
}
template<typename T>
Quaternion<T>& operator *=(Quaternion<T> &a,const Quaternion<T> &b){
a=a*b;
return a;
}
template<typename T>
Quaternion<T>& operator /=(Quaternion<T> &a,const Quaternion<T> &b){
a=a/b;
return a;
}
template<typename T>
Quaternion<T> log(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(log(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return log(abs(a))+u*acos(a.r/abs(a));
}
template<typename T>
Quaternion<T> log10(Quaternion<T> a){
return log(a)/log(10);
}
template<typename T>
Quaternion<T> log2(Quaternion<T> a){
return log(a)/log(2);
}
template<typename T>
Quaternion<T> logp(Quaternion<T> a,Quaternion<T> b){
return log(a)/log(b);
}
template<typename T>
Quaternion<T> exp(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(exp(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return exp(a.r)*(cos(abs(a-a.r))+u*sin(a-a.r));
}
template<typename T>
Quaternion<T> pow(Quaternion<T> a,Quaternion<T> b){
return exp(log(a)*b);
}
template<typename T>
Quaternion<T> nthrt(Quaternion<T> a,Quaternion<T> b){
return pow(a,1.0/b);
}
template<typename T>
Quaternion<T> sqrt(Quaternion<T> a){
return nthrt(a,Quaternion<T>(2.0));
}
template<typename T>
Quaternion<T> cbrt(Quaternion<T> a){
return nthrt(a,Quaternion<T>(3.0));
}
template<typename T>
T dot(Quaternion<T> a,Quaternion<T> b){
return a.r*b.r+a.i*b.i+a.j*b.j+a.k*b.k;
}
template<typename T>
Quaternion<T> pure(Quaternion<T> a){
return Quaternion<T>(0.0,a.i,a.j,a.k);
}
template <typename T>
Quaternion<T> unit(T b,T c,T d){
return Quaternion<T>(0.0,b/sqrt(b*b+c*c+d*d),c/sqrt(b*b+c*c+d*d),d/sqrt(b*b+c*c+d*d));
}
template<typename T>
Quaternion<T> slerp(const Quaternion<T> &a,const Quaternion<T> &b,T t){
T d=dot(a,b),sign=1;
if(d<0){
d=-d;
sign=-1;
}
const T dt=T(0.9995);
if(d>dt){
Quaternion<T> q=a+(b*sign-a)*t;
return q.normalized();
}
T theta=acos(d)*t;
T sin_theta=sin(theta);
T sin_theta0=sin(acos(d));
T s2=sin_theta/sin_theta0;
T s1=cos(theta)-d*s2;
return (a*s1)+(b*sign*s2);
}
template<typename T>
Quaternion<T> sin(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(sin(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return sin(a.r)*cosh(abs(a-a.r))+u*cos(a.r)*sinh(abs(a-a.r));
}
template<typename T>
Quaternion<T> cos(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(cos(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return cos(a.r)*cosh(abs(a-a.r))-u*sin(a.r)*sinh(abs(a-a.r));
}
template<typename T>
Quaternion<T> tan(Quaternion<T> a){
return sin(a)/cos(a);
}
template<typename T>
Quaternion<T> asin(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(asin(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return -u*log(u*a+sqrt(1.0));
}
template<typename T>
Quaternion<T> acos(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(acos(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return -u*log(a+u*sqrt(1.0-a*a));
}
template<typename T>
Quaternion<T> atan(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(atan(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return u/2*log((u+a)/(u-a));
}
template<typename T>
Quaternion<T> sinh(Quaternion<T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(sinh(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return sinh(a.r)*cos(abs(a-a.r))+u*cosh(a.r)*sin(abs(a-a.r));
}
template<typename T>
Quaternion<T> cosh(Quaternion <T> a){
if(abs(a-a.r)<1e-10){
return Quaternion<T>(cosh(a.r),0,0,0);
}
Quaternion<T> u(0.0,a.i/abs(a-a.r),a.j/abs(a-a.r),a.k/abs(a-a.r));
return cos(a.r)*cos(abs(a-a.r))+u*sinh(a.r)*sin(abs(a-a.r));
}
template<typename T>
Quaternion<T> tanh(Quaternion<T> a){
return sinh(a)/cosh(a);
}
template<typename T>
Quaternion<T> asinh(Quaternion<T> a){
return log(a+sqrt(a*a+1));
}
template<typename T>
Quaternion<T> acosh(Quaternion<T> a){
return log(a+sqrt(a-1.0)*sqrt(a+1.0));
}
template<typename T>
Quaternion<T> atanh(Quaternion<T> a){
return 0.5*log((1.0+a)/(1.0-a));
}
template<typename T>
string to_string(Quaternion<T> q){
stringstream ss;
ss<<q;
return ss.str();
}
template<typename T>
istream &operator >>(istream &in,Quaternion<T> &q){
string s;
in>>s;
s.erase(remove_if(s.begin(),s.end(),::isspace),s.end());
if(s.empty()){
q=Quaternion<T>(0.0,0.0,0.0,0.0);
return in;
}
if(s.find("<")!=-1){
string s1=s.substr(0,s.find("<"));
string s2=s.substr(s.find("<")+1,s.size()-s.find("<")-1);
s1.erase(0,1);
s1.erase(s1.size()-1,1);
string s3=s1.substr(0,s1.find(","));
s1.erase(0,s3.size()+1);
string s4=s1.substr(0,s1.find(","));
s1.erase(0,s4.size()+1);
Quaternion<T> u(0.0,stold(s3),stold(s4),stold(s1));
T theta=stold(s2);
q=cos(theta/2)+u*sin(theta/2);
return in;
}
if(s[0]!='+'&&s[0]!='-'){
s='+'+s;
}
vector<string>v;
int start=0;
for(int i=1;i<s.size();i++){
if(s[i]=='+'||s[i]=='-'){
v.push_back(s.substr(start,i-start));
start=i;
}
}
v.push_back(s.substr(start));
for(const auto& t:v){
if(t.empty()){
continue;
}
char sign=t[0];
string coe=t.substr(1);
if(coe.find("i")!=-1){
string nn=coe.substr(0,coe.size()-1);
T val;
try{
val=nn.empty()?1.0L:stold(nn);
}catch(const exception& e){
in.setstate(ios::failbit);
return in;
}
q.i=(sign=='+')?val:-val;
}else{
if(coe.find("j")!=-1){
string nn=coe.substr(0,coe.size()-1);
T val;
try{
val=nn.empty()?1.0L:stold(nn);
}catch(const exception& e){
in.setstate(ios::failbit);
return in;
}
q.j=(sign=='+')?val:-val;
}else{
if(coe.find("k")!=-1){
string nn=coe.substr(0,coe.size()-1);
T val;
try{
val=nn.empty()?1.0L:stold(nn);
}catch(const exception& e){
in.setstate(ios::failbit);
return in;
}
q.k=(sign=='+')?val:-val;
}else{
T val=stold(coe);
q.r=(sign=='+')?val:-val;
}
}
}
}
return in;
}
template<typename T>
ostream &operator <<(ostream &o,const Quaternion<T> &q){
string s1=to_string(q.r),s2=to_string(q.i),s3=to_string(q.j),s4=to_string(q.k);
s1=TrimTrailingZeros(s1);
s2=TrimTrailingZeros(s2);
s3=TrimTrailingZeros(s3);
s4=TrimTrailingZeros(s4);
bool flag=0;
if(s1.find(".")==s1.size()-1){
s1.erase(s1.size()-1,1);
}
if(s2.find(".")==s2.size()-1){
s2.erase(s2.size()-1,1);
}
if(s3.find(".")==s3.size()-1){
s3.erase(s3.size()-1,1);
}
if(s4.find(".")==s4.size()-1){
s4.erase(s4.size()-1,1);
}
if(abs(q.r)>1e-10){
o<<s1;
flag=1;
}
if(q.i){
if(abs(q.i)>1e-10){
if((!flag)||q.i<0){
flag=1;
}else{
o<<"+";
}
}
if(abs(q.i-1)<1e-10){
o<<"i";
}else{
if(abs(q.i+1)<1e-10){
o<<"-i";
}else{
o<<s2<<"i";
}
}
}
if(q.j){
if(abs(q.j)>1e-10){
if((!flag)||q.j<0){
flag=1;
}else{
o<<"+";
}
}
if(abs(q.j-1)<1e-10){
o<<"j";
}else{
if(abs(q.j+1)<1e-10){
o<<"-j";
}else{
o<<s3<<"j";
}
}
}
if(q.k){
if(abs(q.k)>1e-10){
if((!flag)||q.k<0){
flag=1;
}else{
o<<"+";
}
}
if(abs(q.k-1)<1e-10){
o<<"k";
}else{
if(abs(q.k+1)<1e-10){
o<<"-k";
}else{
o<<s4<<"k";
}
}
}
return o;
}
#if __cplusplus>201103L
namespace QLiterals{
Quaternion<float> operator"" _if(unsigned long long x){
return Quaternion<float>(0,(float)x,0,0);
}
Quaternion<float> operator"" _if(long double x){
return Quaternion<float>(0,(float)x,0,0);
}
Quaternion<float> operator"" _jf(unsigned long long x){
return Quaternion<float>(0,0,(float)x,0);
}
Quaternion<float> operator"" _jf(long double x){
return Quaternion<float>(0,0,(float)x,0);
}
Quaternion<float> operator"" _kf(unsigned long long x){
return Quaternion<float>(0,0,0,(float)x);
}
Quaternion<float> operator"" _kf(long double x){
return Quaternion<float>(0,0,0,(float)x);
}
Quaternion<float> operator"" _rf(unsigned long long x){
return Quaternion<float>((float)x,0,0,0);
}
Quaternion<float> operator"" _rf(long double x){
return Quaternion<float>((float)x,0,0,0);
}
Quaternion<float> operator"" _qf(const char* str,size_t len){
string s(str,len);
stringstream ss(s);
Quaternion<float> c;
ss>>c;
return c;
}
Quaternion<float> operator"" _pf(const char* str,size_t len){
string s(str,len);
if(s.find("<")==s.npos){
return Quaternion<float>(stof(s));
}
string s1=s.substr(0,s.find("<"));
string s2=s.substr(s.find("<")+1,s.size()-s.find("<")-1);
s1.erase(0,1);
s1.erase(s1.size()-1,1);
string s3=s1.substr(0,s1.find(","));
s1.erase(0,s3.size()+1);
string s4=s1.substr(0,s1.find(","));
s1.erase(0,s4.size()+1);
Quaternion<float> u(0.0,stof(s3),stof(s4),stof(s1));
float theta=stof(s2);
Quaternion<float> q=cos(theta/2)+u*sin(theta/2);
return q;
}
Quaternion<double> operator"" _id(unsigned long long x){
return Quaternion<double>(0,(double)x,0,0);
}
Quaternion<double> operator"" _id(long double x){
return Quaternion<double>(0,(double)x,0,0);
}
Quaternion<double> operator"" _jd(unsigned long long x){
return Quaternion<double>(0,0,(double)x,0);
}
Quaternion<double> operator"" _jd(long double x){
return Quaternion<double>(0,0,(double)x,0);
}
Quaternion<double> operator"" _kd(unsigned long long x){
return Quaternion<double>(0,0,0,(double)x);
}
Quaternion<double> operator"" _kd(long double x){
return Quaternion<double>(0,0,0,(double)x);
}
Quaternion<double> operator"" _rd(unsigned long long x){
return Quaternion<double>((double)x,0,0,0);
}
Quaternion<double> operator"" _rd(long double x){
return Quaternion<double>((double)x,0,0,0);
}
Quaternion<double> operator"" _qd(const char* str,size_t len){
string s(str,len);
stringstream ss(s);
Quaternion<double> c;
ss>>c;
return c;
}
Quaternion<double> operator"" _pd(const char* str,size_t len){
string s(str,len);
if(s.find("<")==s.npos){
return Quaternion<double>(stod(s));
}
string s1=s.substr(0,s.find("<"));
string s2=s.substr(s.find("<")+1,s.size()-s.find("<")-1);
s1.erase(0,1);
s1.erase(s1.size()-1,1);
string s3=s1.substr(0,s1.find(","));
s1.erase(0,s3.size()+1);
string s4=s1.substr(0,s1.find(","));
s1.erase(0,s4.size()+1);
Quaternion<double> u(0.0,stod(s3),stod(s4),stod(s1));
double theta=stod(s2);
Quaternion<double> q=cos(theta/2)+u*sin(theta/2);
return q;
}
Quaternion<long double> operator"" _ild(unsigned long long x){
return Quaternion<long double>(0,(long double)x,0,0);
}
Quaternion<long double> operator"" _ild(long double x){
return Quaternion<long double>(0,(long double)x,0,0);
}
Quaternion<long double> operator"" _jld(unsigned long long x){
return Quaternion<long double>(0,0,(long double)x,0);
}
Quaternion<long double> operator"" _jld(long double x){
return Quaternion<long double>(0,0,(long double)x,0);
}
Quaternion<long double> operator"" _kld(unsigned long long x){
return Quaternion<long double>(0,0,0,(long double)x);
}
Quaternion<long double> operator"" _kld(long double x){
return Quaternion<long double>(0,0,0,(long double)x);
}
Quaternion<long double> operator"" _rld(unsigned long long x){
return Quaternion<long double>((long double)x,0,0,0);
}
Quaternion<long double> operator"" _rld(long double x){
return Quaternion<long double>((long double)x,0,0,0);
}
Quaternion<long double> operator"" _qld(const char* str,size_t len){
string s(str,len);
stringstream ss(s);
Quaternion<long double> c;
ss>>c;
return c;
}
Quaternion<long double> operator"" _pld(const char* str,size_t len){
string s(str,len);
if(s.find("<")==s.npos){
return Quaternion<long double>(stold(s));
}
string s1=s.substr(0,s.find("<"));
string s2=s.substr(s.find("<")+1,s.size()-s.find("<")-1);
s1.erase(0,1);
s1.erase(s1.size()-1,1);
string s3=s1.substr(0,s1.find(","));
s1.erase(0,s3.size()+1);
string s4=s1.substr(0,s1.find(","));
s1.erase(0,s4.size()+1);
Quaternion<long double> u(0.0,stold(s3),stold(s4),stold(s1));
long double theta=stold(s2);
Quaternion<long double> q=cos(theta/2)+u*sin(theta/2);
return q;
}
}
namespace std{
using namespace QLiterals;
}
#endif//C++14
#endif/*QUATERNION_HPP*/:::