It seems preposterous at first glance. I just want to be sure. Even in 3D the behaviour of rotating objects can be surprising (see the Dzhanibekov effect); in 4D it could be more surprising.

A 2D or 3D spacecraft (with no reaction wheels or gimbaling etc.) needs at least two thrusters to control its spin. See my answer on space.SE.

For any two $4\times4$ matrices $X$ and $Y$, define the commutator $[X,Y]=XY-YX$, the anticommutator $\{X,Y\}=XY+YX$, and the Frobenius inner product $\langle X,Y\rangle=\text{tr}(X^\top Y)$, where $\text{tr}$ is the trace and $\top$ is the transpose.

Let $M$ be a symmetric positive-definite $4\times4$ matrix, and $T=\mathbf f\,\mathbf r^\top-\mathbf r\,\mathbf f^\top$ an antisymmetric $4\times4$ matrix. ($M$ describes the distribution of mass in the spacecraft, $\mathbf r$ is a vector locating the thruster relative to the centre of mass, $\mathbf f$ is the force produced by the thruster, and $T$ is the torque produced by the thruster. More details here. Everything is described in the rotating reference frame.)

The angular velocity $\Omega(t)$, an antisymmetric $4\times4$ matrix, changes with time $t$ according to Euler's equation


where $f(t)\geq0$ is a function (continuous, piecewise-constant, or just integrable) describing when and how strongly the thruster is used.

Question: Can $M$ and $T$ be chosen such that, for any two antisymmetric $4\times4$ matrices $\Omega_0$ and $\Omega_1$, there exist $t_1>0$ and $f$ such that the solution $\Omega$ to Euler's equation with initial value $\Omega(0)=\Omega_0$ has final value $\Omega(t_1)=\Omega_1$?

To simplify things, we may take $f$ to be a combination of Dirac deltas instead of an ordinary function. Then the angular velocity satisfies $\{M,\Omega'(t)\}+[\Omega(t),\{M,\Omega(t)\}]=0$ except at a finite set of times when $\{M,\Omega(t)\}$ changes by a positive multiple of $T$.

Here is Euler's equation in terms of the components $\omega_{ij}$ of the angular velocity. Assume that $M$ is diagonal, with components $m_i>0$.


(And permute the indices to get 6 equations like this.)

If $m_1=m_2$, then $\omega_{12}'$ has constant sign; $\omega_{12}$ is either always non-increasing, or always non-decreasing. So, if the difference between initial and final values of $\omega_{12}$ has the wrong sign compared to $\tau_{12}$, then there is no solution. Thus, we must take $m_1\neq m_2$, and similarly $m_i\neq m_j$ for $i\neq j$.

The Frobenius inner product of two antisymmetric matrices has every term duplicated: $\sum_{i,j}x_{ij}y_{ij}=2\sum_{i<j}x_{ij}y_{ij}$ (since $x_{ji}=-x_{ij}$ and $y_{ji}=-y_{ij}$). So it's natural to take half of this.

The angular momentum is $L(t)=\{M,\Omega(t)\}$. Its squared magnitude is $\tfrac12\langle L(t),L(t)\rangle$; the derivative of this is $\langle L(t),f(t)T\rangle$. So the angular momentum has constant magnitude whenever $f(t)=0$. (The angular momentum itself would be constant in an inertial reference frame, but here we're using a rotating reference frame.)

The rotational energy is $\tfrac14\langle\Omega(t),L(t)\rangle$; the derivative of this is $\tfrac12\langle\Omega(t),f(t)T\rangle$. So the energy is constant whenever $f(t)=0$.

All of that applies in any dimension. But in 4D there's another constant, the angular momentum bivector's exterior square, or equivalently its inner product with its Hodge dual.

Here are those constants in component form. Squared magnitude of angular momentum:


Doubled energy:


Halved exterior square of angular momentum:


(Its derivative is the exterior product of angular momentum and torque: $(m_1+m_2)\omega_{12}\,f\,\tau_{34}+\cdots$ .)

These three constants define two ellipsoids and another quadratic surface (with signature $+^3-^3$) in the 6D space of antisymmetric matrices. When $f=0$, $\Omega$ must remain on the intersection of these three quadratic surfaces. When $f>0$, the above expressions are just quadratic forms, not necessarily constant.

It seems more convenient here to work with $L$ instead of $\Omega$. Of course $l_{ij}=(m_i+m_j)\omega_{ij}$. Let's write a quadratic form as $q(L)=q(L,L)$, using the same symbol for the corresponding bilinear form. All three of the above forms have derivative $2q(L,L')=2q(L,f\,T)$; in other words $q(L,[\Omega,L])=0$.

If we can find some combination of those forms such that $q(L,T)\geq0$ for all relevant $L$, then $q(L)$ would be always non-decreasing in time; there would be no solution when the initial and final values have $q(L_1)<q(L_0)$. But I'm not sure what "relevant" means here; I haven't thought this through.

Notice that Euler's equation has a time symmetry: Given a solution $\Omega(t)$, we can construct another solution $\tilde\Omega(t)=-\Omega(t_1-t)$, so that $\{M,\tilde\Omega'(t)\}+[\tilde\Omega(t),\{M,\tilde\Omega(t)\}]=f(t_1-t)\,T$, and the initial and final values are swapped (and negated): $\tilde\Omega(0)=-\Omega_1$ and $\tilde\Omega(t_1)=-\Omega_0$. Therefore, any angular velocity can be reached from any other, if and only if any angular velocity can be reached from $0$. (We can reverse a path from $0$ to $-\Omega_0$ to get a path from $\Omega_0$ to $0$, and concatenate that with a path from $0$ to $\Omega_1$, to get a path from $\Omega_0$ to $\Omega_1$.)

Also, Euler's equation has a scale symmetry: Given a solution $\Omega(t)$, for any $k>0$ we can construct another solution $\tilde\Omega(t)=k\,\Omega(kt)$, so that $\{M,\tilde\Omega'(t)\}+[\tilde\Omega(t),\{M,\tilde\Omega(t)\}]=k^2f(kt)\,T$, and the final value is $\tilde\Omega(t_1/k)=k\,\Omega(t_1)=k\,\Omega_1$. Therefore, any angular velocity can be reached from $0$, if and only if any angular velocity in a small neighbourhood of $0$ can be reached from $0$.

(It looks like these two symmetries are related by $k=-1$, but we should keep a distinction between positive time and negative time.)

This has the form of a quadratic differential equation, $\mathbf x'(t)=\mathbf x(t)\odot\mathbf x(t)$ where $\odot$ is some bilinear function. (Specifically, for antisymmetric matrices $X$ and $Y$, define $X\odot Y=-\{M,\}^{-1}([X,\{M,Y\}])$, so Euler's equation is $\Omega'=\Omega\odot\Omega$, as long as $f=0$. Alternatively, define $X\odot Y=-[\{M,\}^{-1}(X),Y]$, so Euler's equation is $L'=L\odot L$.)

Fix a norm on the space, and find some constant $\lVert\odot\rVert>0$ such that $\lVert\mathbf x\odot\mathbf y\rVert\leq\lVert\odot\rVert\lVert\mathbf x\rVert\lVert\mathbf y\rVert$ for all $\mathbf x,\mathbf y$.

Given initial value $\mathbf x(0)=\mathbf a$, the equation $\mathbf x'=\mathbf x\odot\mathbf x$ has the power series solution

$$\mathbf x(t)=\sum_{n=0}^\infty t^n\frac{\sum^{n!}\mathbf a^{n+1}}{\sum^{n!}1}$$

where the coefficient of $t^n$ is the average of all possible ways of evaluating $\mathbf a^{n+1}$ using the product $\odot$. For example, the coefficient of $t^3$ is $\tfrac16$ of

$$\sum^6\mathbf a^4=((\mathbf a\mathbf a)\mathbf a)\mathbf a+(\mathbf a(\mathbf a\mathbf a))\mathbf a+2(\mathbf a\mathbf a)(\mathbf a\mathbf a)+\mathbf a((\mathbf a\mathbf a)\mathbf a)+\mathbf a(\mathbf a(\mathbf a\mathbf a))$$ $$\newcommand{\aaaa}[3]{\mathbf a\underset{#1}\odot\mathbf a\underset{#2}\odot\mathbf a\underset{#3}\odot\mathbf a} =\aaaa{1}{2}{3}+\aaaa{2}{1}{3}\begin{matrix}{}+\aaaa{1}{3}{2} \\ {}+\aaaa{2}{3}{1}\end{matrix}+\aaaa{3}{1}{2}+\aaaa{3}{2}{1}.$$

(In the $1$-dimensional algebra, the product is associative, and the power series simplifies to $\sum_nt^na^{n+1}=a/(1-at)$.)

This converges absolutely, for small $t$:

$$\sum_{n=0}^\infty\left\lVert t^n\frac{\sum^{n!}\mathbf a^{n+1}}{\sum^{n!}1}\right\rVert\leq\sum_{n=0}^\infty|t|^n\lVert\odot\rVert^n\lVert\mathbf a\rVert^{n+1}=\frac{\lVert\mathbf a\rVert}{1-|t|\lVert\odot\rVert\lVert\mathbf a\rVert},$$ $$|t|<\frac{1}{\lVert\odot\rVert\lVert\mathbf a\rVert}.$$

And at certain times (when the thruster is used) an impulse may be applied to $\mathbf x$, with a fixed direction $\mathbf b$ but an arbitrary magnitude $c>0$, thus: $\mathbf x(t+)=\mathbf x(t-)+c\mathbf b$. These impulses, and the time intervals between them, are many variables that we can control. If the number of variables is at least $6$ (or the dimension of the space $\mathbf x$ is in), then there is hope of surrounding $0$ in an open set, and thus reaching everywhere in the space (according to the previous section).

$$\mathbf x(0+)=0+c_0\mathbf b$$ $$\mathbf x(t_1-)=\sum_{n=0}^\infty\frac{t_1^n}{n!}\sum^{n!}\mathbf x(0+)^{n+1}$$ $$\mathbf x(t_1+)=\mathbf x(t_1-)+c_1\mathbf b$$ $$\mathbf x(t_1+t_2-)=\sum_{n=0}^\infty\frac{t_2^n}{n!}\sum^{n!}\mathbf x(t_1+)^{n+1}$$ $$\mathbf x(t_1+t_2+)=\mathbf x(t_1+t_2-)+c_2\mathbf b$$ $$\mathbf x(t_1+t_2+t_3-)=\sum_{n=0}^\infty\frac{t_3^n}{n!}\sum^{n!}\mathbf x(t_1+t_2+)^{n+1}$$ $$c_0,t_1,c_1,t_2,c_2,t_3\geq0$$ $$\mathbf x(t_1+t_2+t_3-)\approx0\quad?$$

  • 1
    $\begingroup$ What's the relationship of $M$ with the inertia tensor? My understanding was that the inertia tensor is a linear transformation on bivectors, which in this case are 6-dimensional, so it would be a symmetric 6x6 matrix. Is that not correct? $\endgroup$ Jun 26, 2022 at 21:38
  • $\begingroup$ The inertia tensor is the function $\Omega\mapsto L=\{M,\Omega\}=M\Omega+\Omega M$. $\endgroup$
    – mr_e_man
    Jun 27, 2022 at 12:47
  • $\begingroup$ And of course I'm considering bivectors as antisymmetric matrices. $\endgroup$
    – mr_e_man
    Jun 27, 2022 at 13:01
  • $\begingroup$ Given that $M$ is diagonal (which is no loss of generality, since any symmetric matrix has an orthogonal eigenbasis), the inertia tensor is also diagonal, with components $m_1+m_2$, $m_1+m_3$, $m_2+m_3$, $m_1+m_4$, $m_2+m_4$, $m_3+m_4$. $\endgroup$
    – mr_e_man
    Jul 4, 2022 at 0:06
  • $\begingroup$ Oh, I didn't explain what $M$ is (though the linked physics.SE post did explain it). It's the matrix of second moments of mass: $M=\int\mathbf r\,\mathbf r^T\,\rho\,dV$, where $\mathbf r$ is position relative to the centre of mass, $\rho$ is density, and $V$ is volume. But that's not necessary for understanding this Question. $\endgroup$
    – mr_e_man
    Jul 4, 2022 at 17:53


You must log in to answer this question.