夏が来れば思い出す、Euler角、四元数(1): ちょっと計算してみた。

はじめに

タイトルに書いておいていきなりこんなことを書くのもなんですが、Euler角は別の機会に取り上げる(ジンバルロックとかあまり考えたくないので… (´・ω・｀))かもしれないことにして、この記事では四元数について思い出したことを書いていきます。

四元数とは?

とりあえず、説明。

数学的に厳密な説明ではないかもしれませんが、説明しよう 🙂

四元数とは簡単に言うと複素数を拡張した数体系であり、複素数では1個であった虚数単位が四元数では3個になっています。

一般に四元数$x$は$p,q,r,s \in \mathbb{R}$を用いて$x = p+qi+rj+sk$と表すことができます。

四元数の3個の虚数単位(1以外の基底元)を$i,j,k$と置くと、以下の関係が成り立ちます。

\begin{eqnarray}
i^2 &=& j^2 = k^2 = ijk = -1 \label{eq:quarterniondef}
\end{eqnarray}

って、ハミルトンが1843年に最初に考えた式そのままですが(汗)、とにかく成り立ちます。

それで、何に使えるの?

物体の姿勢をEuler角を使って表現すると、姿勢に対する角度が決まらない特異点が必ず存在しますが、四元数で姿勢を表すと、そのような問題を回避することができます。

そこで、人工衛星や宇宙探査機、さらには3次元CGにおける物体の運動などといった、どういう姿勢を取るかがわからなかったり、姿勢の変化が微小であるとは限らない物体の姿勢の表現や制御などに使えます。

取扱上の注意

四元数を扱う上で最も注意しなければならないのは、

「乗法に関する交換則が成立しない。」

ことでしょう。行列の演算みたいですね。

そこで、(\ref{eq:quarterniondef})式の関係式を使って乗法に関する交換則が成り立つかどうか調べてみます。

(\ref{eq:quarterniondef})式の等式のうち、$ijk = -1$の両辺に右から$-k$をかけると、(\ref{eq:ijeqk})式が求まります。

\begin{eqnarray}
ij &=& k \label{eq:ijeqk}
\end{eqnarray}

また、$ijk = -1$の両辺に左から$-i$をかけると

\begin{eqnarray}
jk &=& i \label{eq:jkeqi}
\end{eqnarray}

になったりします。さらに(\ref{eq:jkeqi})式の両辺の左から$j$をかけると、(\ref{eq:keqji})式が得られます。

\begin{eqnarray}
-k &=& ji \label{eq:keqji}
\end{eqnarray}

(\ref{eq:ijeqk})式及び(\ref{eq:keqji})式を比較すると、$i,j$の順序を入れ替えると計算結果の符号が入れ替わることがわかります。よって、$i,j$については乗法に関する交換則が成り立たないことがわかります。

また、(\ref{eq:ijeqk})式の両辺に左から$j$をかけると、

\begin{eqnarray}
-i &=& kj \label{eq:ieqkj}
\end{eqnarray}

となるので、(\ref{eq:jkeqi})式及び(\ref{eq:ieqkj})式より$j,k$についても乗法に関する交換則が成り立たないことがわかります。

$k,i$の場合についても(\ref{eq:keqji})式の両辺に$-i$を左からかけると

\begin{eqnarray}
ki &=& j \label{eq:kieqj}
\end{eqnarray}

が、(\ref{eq:jkeqi})式の両辺に右から$k$をかけることで、

\begin{eqnarray}
-j &=& ik \label{eq:jeqik}
\end{eqnarray}

がそれぞれ得られるので、(\ref{eq:kieqj})式及び(\ref{eq:jeqik})式より乗法に関する交換則が成り立たないことがわかります。

したがって、1以外の基底元$i,j,k$については乗法に関する交換則が成り立たないことがわかります。

(\ref{eq:quarterniondef})式から(\ref{eq:jeqik})式までの結果及び1と$i,j,k$の積は以下の表のようにまとめることができます。

	$1$	$i$	$j$	$k$
$1$	$1$	$i$	$j$	$k$
$i$	$i$	$-1$	$k$	$-j$
$j$	$j$	$-k$	$-1$	$i$
$k$	$k$	$j$	$-i$	$-1$

乗算の演算。

ここまでの説明で四元数の真髄を理解していただけたと思いますので、四元数どうしの乗算を計算してみます。
四元数$q_1,q_2$の乗算は、

\begin{eqnarray}
q_1&=&a_1+b_1i+c_1j+d_1k \nonumber \\
q_2&=&a_2+b_2i+c_2j+d_2k \nonumber
\end{eqnarray}

とおくと、

\begin{eqnarray}
q_1q_2 &=& a_1a_2-b_1b_2-c_1c_2-d_1d_2 \nonumber \\
&{}& +a_1b_2i+a_1c_2j+a_1d_2k \nonumber \\
&{}& +b_1a_2i+b_1c_2k-b_1d_2j \nonumber \\
&{}& +c_1a_2j-c_1b_2k+c_1d_2i \nonumber \\
&{}& +d_1a_2k+d_1b_2j-d_1c_2i \nonumber \\
&=& a_1a_2-b_1b_2-c_1c_2-d_1d_2 \nonumber \\
&{}& +(a_1b_2+b_1a_2+c_1d_2-d_1c_2)i \nonumber \\
&{}& +(a_1c_2-b_1d_2+c_1a_2+d_1b_2)j \nonumber \\
&{}& +(a_1d_2+b_1c_2-c_1b_2+d_1a_2)k \label{eq:qmul}
\end{eqnarray}

と計算できます。

…といっても、(\ref{eq:qmul})式のままだと$i,j,k$について整理したにもかかわらず、文字数が多すぎて繁雑な感じがします。

そこで、応用面も見越してもう少し見通しの良さそうな書き方ができないか考えてみます。$\boldsymbol{b}_1 = (b_1,c_1,d_1), \boldsymbol{i}_1 = \displaystyle\left(
\begin{array}{c}
i_1 \\
j_1 \\
k_1
\end{array}
\right)$とおくと、$q_1 = a_1+\boldsymbol{b}_1\cdot\boldsymbol{i}_1$と書けます。$q_2$についても同様に$q_2 = a_2+\boldsymbol{b}_2\cdot\boldsymbol{i}_2$と書けます。

(\ref{eq:qmul})式の右辺を$a_1, a_2$及び$i,j,k$が入らない(=実数になる)項に着目して再度整理し直すと、

\begin{eqnarray}
q_1q_2 &=& a_1a_2-b_1b_2-c_1c_2-d_1d_2 \nonumber \\
&{}& +a_1(b_2i+c_2j+d_2k) \nonumber \\
&{}& +a_2(b_1i+c_1j+d_1k) \nonumber \\
&{}& +(c_1d_2-d_1c_2)i \nonumber \\
&{}& +(-b_1d_2+d_1b_2)j \nonumber \\
&{}& +(b_1c_2-c_1b_2)k \nonumber \\
&=& a_1a_2-\boldsymbol{b}_1\cdot\boldsymbol{b}_2 +
a_1(\boldsymbol{b}_2\cdot\boldsymbol{i}) + a_2(\boldsymbol{b}_1\cdot\boldsymbol{i}) \nonumber \\
&{}& +(c_1d_2-d_1c_2)i \nonumber \\
&{}& +(-b_1d_2+d_1b_2)j \nonumber \\
&{}& +(b_1c_2-c_1b_2)k \nonumber \\
&=& a_1a_2-\boldsymbol{b}_1\cdot\boldsymbol{b}_2 +
a_1(\boldsymbol{b}_2\cdot\boldsymbol{i}) + a_2(\boldsymbol{b}_1\cdot\boldsymbol{i}) + (\boldsymbol{b}_1 \times \boldsymbol{b}_2) \cdot \boldsymbol{i} \label{eq:qmulimproved}
\end{eqnarray}

となります。内積や外積を用いて何とか表現できるので、少しわかりやすい式になったかと思います。（｀・ω・´）

(\ref{eq:qmulimproved})式の右辺第3式の最後の項が虚数成分の外積の式になっているため、一般の四元数$q_1,q_2$については積の交換則が成り立たないことがわかります。

四元数を用いた三次元座標系の回転変換の表現

前項で何とか四元数の乗算ができたので、これを使って四元数を用いた三次元座標系の回転変換を行うことを考えます。

ここで突然ですが、単位ベクトル$\boldsymbol{u} = \displaystyle\left(
\begin{array}{c}
u_x \\
u_y \\
u_z
\end{array}
\right), ||\boldsymbol{u}|| = 1$($||\cdots||$は$\cdots$のノルムを表します。)の周りを右ねじが進む方向に$\theta$だけ回転させる回転移動を考え、さらに、$q = \cos\displaystyle\frac{\theta}{2} + \boldsymbol {u}\boldsymbol{i}\sin\displaystyle\frac{\theta}{2}$とおきます。また$q$の共役四元数を$\overline{q} = \cos\displaystyle\frac{\theta}{2} – \boldsymbol{u}\boldsymbol{i}\sin\displaystyle\frac{\theta}{2}$とします。

ここで、ある四元数$x = w + \boldsymbol{r}\boldsymbol{i}, w \in \mathbb{R}, \boldsymbol{r} \in \mathbb{R}^3, \boldsymbol{i} = \displaystyle\left(
\begin{array}{c}
i \\
j \\
k
\end{array}
\right)$について、以下の式により計算できる$x^{\prime}$を求めることを考えます。

\begin{eqnarray}
x^{\prime} &=& qx\overline{q} \label{eq:qxq}
\end{eqnarray}

少々計算が面倒ではありますが、(\ref{eq:qmulimproved})式の関係を用いて(\ref{eq:qxq})式の右辺を頑張って計算すると、

\begin{eqnarray}
\def\coshalftheta{\cos\frac{\theta}{2}}
\def\sinhalftheta{\sin\frac{\theta}{2}}
\def\bmu{\boldsymbol{u}}
\def\bmi{\boldsymbol{i}}
\def\bmr{\boldsymbol{r}}
x^{\prime} &=& \left[\coshalftheta + \bmu\bmi\sinhalftheta\right]\left[w +\bmr\bmi\right]\left[\coshalftheta – \bmu\bmi\sinhalftheta\right] \nonumber \\
&=& \left[w\coshalftheta – \bmu\cdot\bmr\sinhalftheta+\left[\bmr\coshalftheta+w\bmu\sinhalftheta+(\bmu\times\bmr)\sinhalftheta\right]\bmi\right] \cdot \nonumber \\
&{}& \left[\coshalftheta – \bmu\bmi\sinhalftheta\right] \nonumber \\
&=& w\cos^2\frac{\theta}{2}-\bmu\cdot\bmr\sinhalftheta\coshalftheta \nonumber \\
&{}& +\left[\bmr\coshalftheta+w\bmu\sinhalftheta+(\bmu\times\bmr)\sinhalftheta\right]\cdot\bmu\sinhalftheta \nonumber \\
&{}& -\left[w\coshalftheta-\bmu\cdot\bmr\sinhalftheta\right]\cdot\bmu\bmi\sinhalftheta \nonumber \\
&{}& +\left[\bmr\coshalftheta+w\bmu\sinhalftheta+(\bmu\times\bmr)\sinhalftheta\right]\bmi\coshalftheta \nonumber \\
&{}& -\left[\left[\bmr\coshalftheta+w\bmu\sinhalftheta+(\bmu\times\bmr)\sinhalftheta\right]\times\bmu\right]\bmi\sinhalftheta \label{eq:xprime}
\end{eqnarray}

(\ref{eq:xprime})式を見る限り、この計算の先行きにかなり不安を覚える感じもしますが、まず実部の係数をまとめることを考えます。$\boldsymbol{u}\times\boldsymbol{r}$と$\boldsymbol{u}$は直交することから、$(\boldsymbol{u}\times\boldsymbol{r})\cdot\boldsymbol{u} = 0$となりますので、

\begin{eqnarray}
\def\coshalftheta{\cos\frac{\theta}{2}}
\def\sinhalftheta{\sin\frac{\theta}{2}}
\def\bmu{\boldsymbol{u}}
\def\bmi{\boldsymbol{i}}
\def\bmr{\boldsymbol{r}}
&{}&w\cos^2\frac{\theta}{2}-\bmu\cdot\bmr\sinhalftheta\coshalftheta \nonumber \\
&{}& +\left[\bmr\coshalftheta+w\bmu\sinhalftheta+(\bmu\times\bmr)\sinhalftheta\right]\cdot\bmu\sinhalftheta \nonumber \\
&=& w\cos^2\frac{\theta}{2}-\bmu\cdot\bmr\sinhalftheta\coshalftheta+\bmr\cdot\bmu\coshalftheta\sinhalftheta+w\sin^2\frac{\theta}{2} \nonumber \\
&=& w \label{eq:w}
\end{eqnarray}

となるので、実部の値($w$)は(\ref{eq:qxq})式の計算では変化しないことがわかります。

(\ref{eq:w})式で実部は計算できたので、次は虚部をまとめることを考えます。虚部の項のうち$w$が現れる項は3項ありますが、これらの和は$0$となります。したがって、残りの項を実部まで含めて表すと、(\ref{eq:im})式が成り立ちます。

\begin{eqnarray}
\def\coshalftheta{\cos\frac{\theta}{2}}
\def\sinhalftheta{\sin\frac{\theta}{2}}
\def\bmu{\boldsymbol{u}}
\def\bmi{\boldsymbol{i}}
\def\bmr{\boldsymbol{r}}
x^{\prime} &=& w+\left[(\bmu\cdot\bmr)\cdot\bmu\sin^2\frac{\theta}{2}+\bmr\cos^2\frac{\theta}{2}+(\bmu\times\bmr)\sinhalftheta\coshalftheta\right. \nonumber \\
&{}& \left.-(\bmr\times\bmu)\coshalftheta\sinhalftheta-\left[(\bmu\times\bmr)\times\bmu\right]\sin^2\frac{\theta}{2}\right]\bmi \nonumber \\
&=& w+\left[(\bmu\times\bmr)\sin\theta+\bmr\cos^2\frac{\theta}{2}+(\bmu\cdot\bmr)\cdot\bmu\sin^2\frac{\theta}{2}\right. \nonumber \\
&{}& \left.-\left[(\bmu\times\bmr)\times\bmu\right]\sin^2\frac{\theta}{2}\right]\bmi\label{eq:im}
\end{eqnarray}

ここで、$\boldsymbol{r}$を$\boldsymbol{u}$に平行な成分$\boldsymbol{r}_u$と$\boldsymbol{u}$と直交する成分$\boldsymbol{r}_l$に分解します。すると、

と書けるので、(\ref{eq:bmru})及び(\ref{eq:bmrl})式が成り立ちます。

\begin{eqnarray}
\def\coshalftheta{\cos\frac{\theta}{2}}
\def\sinhalftheta{\sin\frac{\theta}{2}}
\def\bmu{\boldsymbol{u}}
\def\bmi{\boldsymbol{i}}
\def\bmr{\boldsymbol{r}}
(\bmu\cdot\bmr)\bmu &=& \left[\bmu\cdot(\bmr_u+\bmr_l)\right]\cdot\bmu \nonumber \\
&=& || \bmr_u || \bmu \nonumber \\
&=& \bmr_u \label{eq:bmru} \\
(\bmu\times\bmr)\times\bmu &=& (\bmu\times\bmr_l)\times\bmu \nonumber \\
&=& \bmr_l \label{eq:bmrl}
\end{eqnarray}

(\ref{eq:decomposition})式(\ref{eq:bmru})及び(\ref{eq:bmrl})の各式を(\ref{eq:im})式に代入すると、

\begin{eqnarray}
\def\coshalftheta{\cos\frac{\theta}{2}}
\def\sinhalftheta{\sin\frac{\theta}{2}}
\def\bmu{\boldsymbol{u}}
\def\bmi{\boldsymbol{i}}
\def\bmr{\boldsymbol{r}}
x^{\prime} &=& w+\left[(\bmu\times\bmr_l)\sin\theta+(\bmr_u+\bmr_l)\cos^2\frac{\theta}{2}+\bmr_u\sin^2\frac{\theta}{2}\right. \nonumber \\
&{}& \left.-\bmr_l\sin^2\frac{\theta}{2}\right]\bmi \nonumber \\
&=& w+\left[(\bmu\times\bmr_l)\sin\theta+\bmr_u+\bmr_l\cos\theta\right]\bmi \label{eq:result}
\end{eqnarray}

と変形できます。

ここで$\boldsymbol{r}_l$は$\boldsymbol{u}$と直交するので、$||\boldsymbol {u}\times\boldsymbol{r}_l|| = || \boldsymbol{r}_l ||$となり、かつ$\boldsymbol {u}\times\boldsymbol{r}_l$と$\boldsymbol{r}_l$は直交するので、(\ref{eq:result})式右辺第2式の虚数部の係数(ベクトル)は以下の1.及び2.を満たすことがわかります。

$\boldsymbol{r}_u$を中心とし、$\boldsymbol{u}\times\boldsymbol{r}_l$と$\boldsymbol{r}_l$が張る平面上にある円上の点に原点から向かうベクトルであること。
1.のベクトルの終点は$\boldsymbol{r}$を$\boldsymbol{u}$のまわりに$\theta$だけ回転させた点であること。

よって、(10)式が三次元座標系における回転変換を表現していることを示すことができました。

まとめ

実はこの記事を書いている間に四元数については(プログラミング等も含めて)すでに詳しい記事がいくつかある(「参考文献」参照。)ことがわかったので、できるだけ別の視点というか、計算の詳細な過程を書いてみました。

この記事では$w \in \mathbb{R}$とおいてから計算を行いました。実際に(10)式の計算を行うプログラムを作成する場合には、$w=0$とおくと$w$が現れる項の計算はしなくてもよくなりますので、少しだけ計算量を減らすことができるかもしれません。

ご覧いただいた方のお役に立つかどうかはわかりませんが、参考にしていただけると幸いです。

この記事は以上ですが、Scalaで実装してみましたので、ついでに読んでいっていただけると幸いです(リンクはこちら)。