N88-BASICでクォータニオン(quaternion) (1回目)
2021/6/25(金)
N88-BASICでクォータニオン(quaternion) (1回目)
(by ULproject for N88-BASIC,NL-BASIC)
行列(matrix)の代わりに四元数(しげんすう)
(quaternion)を使用して見ました。
以後、四元数と行列を大文字、それ以外を
小文字で表します
四元数は複素数(2元)の四元版です。
ii = jj = kk = -1
ij = k, jk = i, ki = j
ji = -k, kj = -i, ik = -j
となる、i,j,kを定義します。
複素数 c = a + bi
共役複素数 c~ = a - bi
(cにオーバーバーですがここではc~
とします。)
四元数 Q = a + bi + cj + dk
共役四元数 Q~ = a - bi - cj - dk
複素数は複素平面上の点cを表せます。
四元数では1つの実数wと三次元空間
の点q=(x,y,z)を表すことにします。
Q = w + ix + jy + kz
= (w, x,y,z)
= (w, q)
wはスカラー、qはベクトル
i,j,kがx,y,z軸の単位ベクトル
になっているのに似ています。
q = (x,y,z) = ix + iy + kz
四元数の積
Q = w + ix + jy + kz
Q' = w' + ix' + jy' + kz'
QQ' = ww'+ wix'+ wjy'+ wkz'
+ ixw'+iixx'+ijxy'+ikxz'
+ jyw'+jiyx'+jjyy'+jkyz'
+ kzw'+kizx'+kjzy'+kkzz'
= ww'+ wix'+ wjy'+ wkz'
+ ixw'- xx'+ kxy'- jxz'
+ jyw'- kyx'- yy'+ iyz'
+ kzw'+ jzx'- izy'- zz'
QQ'= (ww'- xx'- yy'- zz',
wx'+ xw'+ yz'- zy',
wy'+ yw'+ zx'- xz',
wz'+ zw'+ xy'- yx')
となります。また、
QQ' = (w , q )(w', q')
= (ww'- q･q',wq'+ qw'+ q×q')
となります。
ここで、ベクトルa,b,c,dの公式を
まとめておきます。
(a×b)･c = (c×a)･b
= (b×c)･a
(a×b)･(c×d) = (a･c)(b･d)
- (a･d)(b･c)
(a×b)×c = (a･c)b - (b･c)a
a×(b×c) = (a･c)b - (a･b)c
次に点pと四元数の積ですが、
P = (0,p)、P' = (0, p')
Q = (w,q)、Q~ = (w,-q )
として、
P'= QPQ~を求めます。
P'= QPQ~
= Q(p･q,pw-p×q)
= (w(p･q)-q･(pw-p×q),w(pw-p×q)+q(p･q)+q×(pw-p×q))
= (-q･(p×q),w2p-w(p×q)+q(p･q)+(q×p)w-q×(p×q))
= (-p･(q×q),w2p+2w(q×p)+q(p･q)-(q･q)p+(q･p)q)
= (0,(w2-q･q)p+2(q･p)q+2w(q×p))
よって、
p'= (w2-q･q)p+2(q･p)q+2w(q×p)
q = nuと置いて、
p'= (w2-u2)p+2u2(n･p)n+2wu(n×p)
ここで、nベクトル回りの回転の式
p'= pcosθ+(1-cosθ)(n･p)n+(n×p)sinθ
と比べると、
半角公式より
cos2(θ/2)-sin2(θ/2)=(1+cosθ)/2-(1-cosθ)/2
=cosθ
2sin(θ/2)sin(θ/2) = 1-cosθ
2cos(θ/2)sin(θ/2) = sinθ
に
w = cos(θ/2)
u = sin(θ/2)
と置いて左辺に代入すると
w2-u2 = cosθ
2u2 = 1-cosθ
2wu = sinθ
となり
p'= (w2-u2)p+2u2(n･p)n+2wu(n×p)
p'= pcosθ+(1-cosθ)(n･p)n+(n×p)sinθ
が同じ式だと分かります
w = cos(θ/2)
q = nsin(θ/2) (q = nuより)
として、
P'= QPQ~
(0,p') = (w,q)(0,p)(w,-q)
でpをnベクトル回りにθ回転させる
事ができます
平行移動(並進)を追加するために
回転用と並進用2つの四元数が必要なため、
四元数の複合体(四元数2つ分)を定義します。
回転四元数をQr、並進四元数をQtとして、
Q=(Qr, Qt)を定義し、ここでは複合四元数
と呼ぶ事にします。
四元数の回転と並進は
p'=Qr p Qr~ + Qt
これを、行列と比較すると、
行列Mの回転行列をMr並進行列をMt
として、M=(Mr, Mt)と表すと、
|a d g j| |a d g| |j|
M =|b e h k|, Mr=|b e h|, Mt=|k|
|c f i l| |c f i| |l|
|0 0 0 1|
p'= Mp は、
|a d g j||x| |a d g||x| |j|
p'=|b e h k||y|=|b e h||y|+|k|
|c f i l||z| |c f i||z| |l|
|0 0 0 1||1| つまり、
p'= Mr p + Mtとなり、
p'= Qr p Qr~ + Qtと似たような
式になります。
M = MB
|? ? ? l||･ ･ ･ a|
M =|? ? ? m||･ ･ ･ b|
|? ? ? n||･ ･ ･ c|
|0 0 0 1||0 0 0 1|
|? ? ? ?a+?b+?c+l|
=|? ? ? ?a+?b+?c+m|
|? ? ? ?a+?b+?c+n|
|0 0 0 1|
M = (MrBr, Mt + MrBt )、よって、
Q = (QrAr, Qr + QrAtQr~)
M = BM
|? ? ? a||･ ･ ･ l|
M =|? ? ? b||･ ･ ･ m|
|? ? ? c||･ ･ ･ n|
|0 0 0 1||0 0 0 1|
|? ? ? ?l+?m+?n+a|
=|? ? ? ?l+?m+?n+b|
|? ? ? ?l+?m+?n+c|
|0 0 0 1|
M = (BrMr, Bt BrMt )、よって、
Q = (ArQr, At + ArQtAr~)
となります。
M,Bは行列、Q,Aは複合四元数です。
まとめ、
p'= Mp = Mr p + Mt
p'= QpQ~ = Qr p Qr~ + Qt
M = MB = (MrBr, Mt + MrBt )
Q = QA = (QrAr, Qt + QrAtQr~)
M = BM = (BrMr, Bt + BrMt )
Q = AQ = (ArQr, At + ArQtAr~)
BASICでは、DIM Q(7)で
Q(0)～Q(3)でQr
Q(4)～Q(7)でQt
を表現します。
単位複合四元数は
Q(0)～Q(7)=
(1, 0, 0, 0, 0, 0, 0, 0)
としています。
サブルーチン*MULTQQを使わなければ
、0での掛け算などの無駄を省いて、
もっと短いプログラムになります
今回は、上記式の通りプログラムを
書いている事が分かり易いように、
四元数の積は*MULTQQで計算させて
います。
*MULTQQは四元数の積
(q0, q1, q2, q3)
= (w1, x1, y1, z1)(w2, x2, y2, z2)
の計算をします。
QQ'= (ww'- xx'- yy'- zz',
wx'+ xw'+ yz'- zy',
wy'+ yw'+ zx'- xz',
wz'+ zw'+ xy'- yx')
の通りに計算させています
quad001.basを参照して下さい。
NL-BASICとblg～.zip(quad001.bas)は
以下のリンクからダウンロードできます
Readme.txtを読んで遊んで下さい
動作は行列と同じなので行列の動画を貼っておきます
