Quaternions – If they can’t see you, they can’t eat you

I did not know quaternions can manage the art of camouflage – until today. Camouflage and metamorphosis, they come together. People are wearing masks, to disguise, to cheat, to the extent that they cheat themselves. Quaternions took me by surprise. As a kid I was told about metamorphosis.

I never witnessed it myself. Until, as I said, today, until I looked at these two pictures:

Quaternion path a’la van Zon and Schofield

Quaternion path a’la van Celledoni and Zanna (Kosenko)

In theory, if I am not making errors here, both images are supposed to represent the same story, except that it is being told by two different witnesses. It looks like a butterfly turning into caterpillar. How can it be? It is, after all, precise mathematics, not a psychology.
The first image represents the path in the group \mathrm{SU}(2) that we have discussed before in Introducing geodesics. Then, in Standing on the shoulders of giants – Reboot, examining the image obtained using my old code based on the paper of Van Zon and Schofield, I realized that:

We have a nice curly trajectory, but there are also strange straight line spikes. On animation they are, perhaps, harmless. But if something like this happens with the software controlling the flight of airplanes, space rockets or satellites – people may die. There is a BUG in the algorithm.

Bugs are no good. As I have written it elsewhere:

I believe that the Universe has Purpose, that it is much like a computer program of great complexity, and that “we” – the IGUS-es – have a role in its evolution. For a while our role can be described simply as “debugging units.” In short, my present answer to the question “why are we here?” reads: DEBUGGING THE UNIVERSE.

With the somewhat less ambitious task of removing the zigzag bug from the spinning quaternion history we decided to follow the giants and we have ended in the new algorithm based on the paper and the code by Celledoni and Zanna, the code described in The final answer for the Universe in which m<1. I will repeat this code here, with a small change that makes it more universal, as discussed in Attitude matrix for m<1.

The case m<1.
The solution \mathbf{L}(t) of the Euler’s equations has two trajectories. For one with L_3(t)>0 we take

(1)   \begin{eqnarray*} A_1&=&\sqrt{\frac{I_1 (d I_3-1)}{I_3-I_1}},\\ A_2&=&\sqrt{\frac{I_2 (d I_3-1)}{I_3-I_2}},\\ A_3&=&\sqrt{\frac{I_3 (1-d I_1)}{I_3-I_1}},\\ \end{eqnarray*}

while for the other one, with L_3(t)<0 we take

(2)   \begin{eqnarray*} A_1&=&-\sqrt{\frac{I_1 (d I_3-1)}{I_3-I_1}},\\ A_2&=&\sqrt{\frac{I_2 (d I_3-1)}{I_3-I_2}},\\ A_3&=&-\sqrt{\frac{I_3 (1-d I_1)}{I_3-I_1}}.\\ \end{eqnarray*}

Let

(3)   \begin{eqnarray*} B&=&\sqrt{\frac{(1-d I_1) (I_3-I_2)}{I_1 I_2 I_3}},\\ m&=&\frac{(d I_3-1) (I_2-I_1)}{(1-d I_1) (I_3-I_2)}. \end{eqnarray*}

The solution \mathbf{L}(t) of the Euler’s equations is given by

(4)   \begin{eqnarray*} L_1(t)&=&A_1\, \cn (Bt,m),\\ L_2(t)&=&A_2\, \sn (Bt,m),\\ L_3(t)&=&A_3\, \dn (Bt,m), \end{eqnarray*}

With constants \alpha and \nu defined as

(5)   \begin{equation*} \alpha=\frac{I_3-I_1}{\sqrt{\frac{I_1(1-d I_1)(I_3-I_2)I_3}{I_2}}}, \end{equation*}

(6)   \begin{equation*} \nu=\frac{I_1-dI_1I_3}{I_3-dI_1I_3} \end{equation*}

we set the phase variable \psi(t) as

(7)   \begin{equation*} \psi(t)=\frac{t}I_1+\arctan\left((A_2/A_3)\mathrm{sd}(Bt,m)  \right)-\alpha \Pi(\nu,\am(Bt,m),m),\ \end{equation*}

where the Jacobi function \mathrm{sd} is defined as \mathrm{sd}(u,m)=\sn(u,m)/\dn(u,m).
Then the quaternionic attitude solution is given by q(t)=W(t)+\mathbf{i}X(t)+\mathbf{j}Y(t)+\mathbf{k}Z(t), with

(8)   \begin{eqnarray*} W(t)&=&\cos \frac{\psi(t)}{2}\sqrt{\frac{1+L_1(t)}{2}},\\ X(t)&=&\sin \frac{\psi(t)}{2}\sqrt{\frac{1+L_1(t)}{2}},\\ Y(t)&=&\frac{L_3(t)\cos \frac{\psi(t)}{2}+L_2(t)\sin \frac{\psi(t)}{2}}{\sqrt{2(1+L_1(t))}},\\ Z(t)&=&\frac{-L_2(t)\cos \frac{\psi(t)}{2}+L_3(t)\sin \frac{\psi(t)}{2}}{\sqrt{2(1+L_1(t))}} \end{eqnarray*}

Why is this better than the previous formulas that were producing BUGS? The secret is in the denominator!

In the formulas for Y(t) and Z(t) we have 1+L_1(t) in the denominator. That is 1+A_1\, \cn (Bt,m). The constant A_1
is given by

    \[A_1=\sqrt{\frac{I_1 (d I_3-1)}{I_3-I_1}}.\]

The smallest possible value of d is 1/I_3. then A_1=0.. The largest possible value of d, while m\leq 1, is 1/I_2, when m=1. Then, with I_1<I_2<I_3,

    \[|A_1|=\sqrt{I_1/I_2}\sqrt{\frac{I_3-I_2}{I_3-I_1}}<1.\]

Therefore A_1<1, therefore 1+A_1\, \cn (Bt,m)>0 (Remember that \cn is a cosinus of the amplitude!). So, we never have a problem with the denominator! That is the advantage of the new algorithm! No infinities, all is finite, what a relieve!

So, we take our flagship example with I_1=1, I_2=2, I_3=3, we choose d=0.34, just a little bit over the minimal value of 0.33.., and we draw the stereographic projection of q(t) from S^3 to \mathbf{R}^3.
What we get is the second image.
No resemblance whatsoever to the first one.
And yet, they are supposed to represent the same reality! How can it be possible? They are like day and night …

So, I created the METAMORPHOSIS – a continuous group translation.

metamorphosis

There is still something strange there…. And it needs to understood. Life is exciting.