The reason he was asking that question is that he wanted to give
various objects a random but physically plausible rotation. The
question I intend to answer is **how do I give an object in blender a
continuous rotation about a random (but stable) axis?**

Blender uses a quaternion representation based on a rotation θ
about a normalized
axis $a$

$w=\mathrm{cos}(\theta /2)$

$x={a}_{x}\mathrm{sin}(\theta /2)$

$y={a}_{y}\mathrm{sin}(\theta /2)$

$z={a}_{z}\mathrm{sin}(\theta /2)$

$q=[w,x,y,z]$

We will express our random rotation as the composition of one fixed orientation quaternion ${q}_{1}$ composed with another quaternion ${q}_{2}$ representing the spin which varies over time.

To pick ${q}_{1}$ we can use the Hypersphere Point Picking procedure. To pick a random axis we can use the Sphere Point Picking procedure, although I like to modify the procedure so that instead of calculating φ I pick a random $z$ in the range [-1,1] and then calculate the $x$ and $y$ based on the random $\theta $ .

Once we have the random orientation ${q}_{1}$ and the varying rotation ${q}_{2}$ we want to compose them ${q}_{3}={q}_{1}{q}_{2}$ . The formula for composing two quaternions is

${w}_{3}={w}_{1}{w}_{2}-{x}_{1}{x}_{2}-{y}_{1}{y}_{2}-{z}_{1}{z}_{2}$

${x}_{3}={w}_{1}{x}_{2}+{x}_{1}{w}_{2}+{y}_{1}{z}_{2}-{z}_{1}{y}_{2}$

${y}_{3}={w}_{1}{y}_{2}-{x}_{1}{z}_{2}+{y}_{1}{w}_{2}+{z}_{1}{x}_{2}$

${z}_{3}={w}_{1}{z}_{2}+{x}_{1}{y}_{2}-{y}_{1}{x}_{2}+{z}_{1}{w}_{2}$

When you substitute the axis+angle formulation for ${q}_{2}$ you end up with

$$\begin{array}{c}\hfill {w}_{3}={w}_{1}\mathrm{cos}(\theta /2)-{x}_{1}{a}_{x}\mathrm{sin}(\theta /2)-{y}_{1}{a}_{y}\mathrm{sin}(\theta /2)-{z}_{1}{a}_{z}\mathrm{sin}(\theta /2)\\ \hfill {x}_{3}={w}_{1}{a}_{x}\mathrm{sin}(\theta /2)+{x}_{1}\mathrm{cos}(\theta /2)+{y}_{1}{a}_{z}\mathrm{sin}(\theta /2)-{z}_{1}{a}_{y}\mathrm{sin}(\theta /2)\\ \hfill {y}_{3}={w}_{1}{a}_{y}\mathrm{sin}(\theta /2)-{x}_{1}{a}_{z}\mathrm{sin}(\theta /2)+{y}_{1}\mathrm{cos}(\theta /2)+{z}_{1}{a}_{x}\mathrm{sin}(\theta /2)\\ \hfill {z}_{3}={w}_{1}{a}_{z}\mathrm{sin}(\theta /2)+{x}_{1}{a}_{y}\mathrm{sin}(\theta /2)-{y}_{1}{a}_{x}\mathrm{sin}(\theta /2)+{z}_{1}\mathrm{cos}(\theta /2)\end{array}$$

If we regroup the equation around the cos and sin terms we get $$\begin{array}{c}\hfill {w}_{3}={w}_{1}\mathrm{cos}(\theta /2)+(-{x}_{1}{a}_{x}-{y}_{1}{a}_{y}-{z}_{1}{a}_{z})\mathrm{sin}(\theta /2)\\ \hfill {x}_{3}={x}_{1}\mathrm{cos}(\theta /2)+({w}_{1}{a}_{x}+{y}_{1}{a}_{z}-{z}_{1}{a}_{y})\mathrm{sin}(\theta /2)\\ \hfill {y}_{3}={y}_{1}\mathrm{cos}(\theta /2)+({w}_{1}{a}_{y}-{x}_{1}{a}_{z}+{z}_{1}{a}_{x})\mathrm{sin}(\theta /2)\\ \hfill {z}_{3}={z}_{1}\mathrm{cos}(\theta /2)+({w}_{1}{a}_{z}+{x}_{1}{a}_{y}-{y}_{1}{a}_{x})\mathrm{sin}(\theta /2)\end{array}$$

With recent versions of blender that support sinusoidal easing on fcurves it is possible to create an fcurve in an action that accurately represents a sine wave or a cosine wave of a fixed magnitude by creating keyframe points where the value of the sine/cosine equals 0, 1, and -1 and setting the ease in/out for each keyframe point appropriately. Even functions of the form $g\mathrm{cos}\left(\theta \right)$ are easy to express by just scaling the y element of the fcurve in the action. Unfortunately the formulas we have derived above combine both sine and cosine terms. In order to express this using blender fcurves, we would have to use NLA strips with blend mode addition.

$\mathrm{sin}(\varphi +\theta \text{'})=\mathrm{sin}\left(\varphi \right)\mathrm{cos}(\theta \text{'})+\mathrm{cos}\left(\varphi \right)\mathrm{sin}(\theta \text{'})$

Let's multiply by r and replace $\theta \text{'}$ with $\theta /2$ :

$r\mathrm{sin}(\varphi +\theta /2)=r\mathrm{sin}\left(\varphi \right)\mathrm{cos}(\theta /2)+r\mathrm{cos}\left(\varphi \right)\mathrm{sin}(\theta /2)$

Now those terms involving r and φ look a lot like polar coordinates, and we can use that technique to transform

$b\mathrm{cos}(\theta /2)+c\mathrm{sin}(\theta /2)$

into

$r\mathrm{sin}(\varphi +\theta /2)$

using the following relationships familiar to anyone versed in polar coordinates:

$$\begin{array}{c}c=r\mathrm{cos}\left(\varphi \right)\\ b=r\mathrm{sin}\left(\varphi \right)\\ r=\sqrt{{b}^{2}+{c}^{2}}\\ \varphi =\mathrm{atan2(c,b)}\end{array}$$

To actually replicate the sine wave in the fcurve we need 5 keyframes corresponding to the following expressions:

$$\begin{array}{ccc}\hfill 0& =& r\mathrm{sin}\left(0\right)\hfill \\ \hfill \mathrm{r}& =& r\mathrm{sin}\left(\genfrac{}{}{0.1ex}{}{\pi}{2}\right)\hfill \\ \hfill 0& =& r\mathrm{sin}\left(\pi \right)\hfill \\ \hfill -\mathrm{r}& =& r\mathrm{sin}\left(\genfrac{}{}{0.1ex}{}{3\pi}{2}\right)\hfill \\ \hfill 0& =& r\mathrm{sin}\left(2\pi \right)\hfill \end{array}$$

Let us use $\alpha $ to represent the values $0$, $\genfrac{}{}{0.1ex}{}{\pi}{2}$, $\pi $, $\genfrac{}{}{0.1ex}{}{3\pi}{2}$, and $2\pi $. We will also designate how long it takes for the object to complete a single rotation by defining $p$ as the period of the rotation as a number of video frames. Next, we designate $f$ to represent a frame number in the following formulas. $$\begin{array}{ccc}\hfill f& =& {f}_{0}+\frac{p}{2\pi}\theta \hfill \\ \hfill \varphi +\theta /2& =& \alpha \hfill \\ \hfill \theta & =& 2(\alpha -\varphi )\hfill \\ \hfill f& =& {f}_{0}+\frac{p}{2\pi}2(\alpha -\varphi )\hfill \\ \hfill f\left(\alpha \right)& =& {f}_{0}+\frac{p}{\pi}(\alpha -\varphi )\hfill \end{array}$$ Now we have the formulas to compute the keyframe and value for the five control points of the fcurve for each quaternion channel. $$\begin{array}{c}\left(f\right(0\left),0\right)\\ \left(f\right(\genfrac{}{}{0.1ex}{}{\pi}{2}\left),r\right)\\ \left(f\right(\pi \left),0\right)\\ \left(f\right(\genfrac{}{}{0.1ex}{}{3\pi}{2}),-r)\\ \left(f\right(2\pi \left),0\right)\end{array}$$

We will also want to use a CYCLES modifier on the fcurve so that it repeats infinitely and the object will keep spinning no matter how long the animation runs.

One oddity that you might notice is that the frame span on the fcurve is actually $2p$ . This is because the quaternions $q$ and $-q$ actually represent the same orientation, and since the sine waves we have crafted make each element vary over both negative and positive values, one cycle of the fcurve actually has one half that can be considered $q$ and the second half can be considered $-q$ .

When we put this into practice with the animate-random-spin.py script we end up with fcurves that look like this:

- How to make a true linear quaternion rotation? (blender stackexchange)
- on-line converter from LaTeX to MathML