Conjugated 3D Virtual Reality Worlds in Spacecraft Attitude Control

Abstract

The chapter aims to contribute to the application of virtual reality (VR) in spacecraft attitude control for systems, possessing “Inertial Morphing (IM)” capabilities.The concept of IM, proposed/applied by authors in 2017, is attractive in design of the fully autonomous future space missions, as swift control, requiring minimized energy and computations, can be achieved with exiguous/paltry morphings. To assist rapid planning of the optimized maneuver scenarios, we propose collocated merging of the various VR simulation worlds for the same spacecraft or merging of the VR worlds for the same system (but in different configurations). This enabled concurrent utilization of different methods of modeling, including Poinsot’s and Bine’s construction and Euler’s equations. Therefore, superimposed VR worlds are called “conjugated.” We present the classical methods of modeling of the torque-free systems and then show their implementation in the developed VR-integrated interactive package. Effectiveness of the VR-conjugated environment is illustrated with its use for planning of spacecraft de-tumbling, spacecraft 180-degrees inversion, and 90-degrees inversion. The developed VR environment enables utilization of both the body-axes spacecraft coordinate system and inertial coordinate system with instant transition from one into another, switching on/off various virtual reality worlds for multiple supports in the process of mission design.

1. Introduction

With the onset of computers and development of advanced computer graphics, virtual reality proved to be a very useful and powerful assistant in various areas of human life, including education, professional activities, and entertainment [1]. It may involve virtual and hybrid reality systems, where plain virtual reality (VR) is an environment that is fully digitized and contains no physical elements and hybrid reality (HR) is mostly a virtualized environment; however, some objects have physical existence. By placing trackers on real-world objects, they can be represented and manipulated inside the hybrid environment while providing tactile cues that are not present in virtual reality simulations [2].

Often, VR software, packages, and hardware aim to focus on achieving high graphical and physical realism or aim to add different levels of immersion, in some cases creating a feeling of actual presence in the virtual world [3]. VR and animations are becoming popular in education, enabling explanation of complex principles, simulated scenarios, and interpretation of results. In medicine, VR can provide effective skill transfer into the operating room or remotely [4]. In space flight contexts, “real-life training”, which is often expensive, demanding for large facilities or even impossible, VR is used to efficiently extend mission control specialists and astronauts training possibilities [5, 6].

Despite wide embrace of VR application areas, literature on the specialized VR tools in relation to the spacecraft attitude dynamics is quite limited. In the most of these specialized application cases, VR is used as a graphical visualization tool for two main purposes: (1) to assist the mission professionals to visualize orbits and/or spacecraft systems externally (illustrated in [7]) or (2) to assist in space crew mission training, helping to understand the local spacecraft or station environment [8]. These two main areas, essentially, are enhancing understanding of the inertial or non-inertial simulated environments.

The current work aims to contribute to the application of VR in spacecraft attitude dynamics, in particular, to the systems, enabling “Inertial Morphing (IM)”: IM with exiguous/paltry number of IM control actions is a new concept, proposed by the authors in 2017 and applied to the spacecraft attitude control [9, 10]. The exiguous control enables swift control, requiring minimized energy and computations, and allows design of the fully autonomous future space missions [11, 12].

This work involves collocated merging of the various virtual reality simulation worlds for the same spacecraft system (where VR worlds are used to illustrate different methods of modeling) or involves merging of the virtual reality worlds for the same system but in different configurations (in this case, the corresponding VR worlds are to be simultaneously used in the development of the integrated spacecraft attitude maneuver scenario). Therefore, these superimposed VR worlds are called “conjugated”. We present the theory behind the developed VR-integrated interactive environment and then the applications of the programmed VR environment, based on the presented theory. Use of the VR-conjugated environment is illustrated with its use for planning of a few illustration attitude maneuvers, including spacecraft de-tumbling, spacecraft 180-degrees inversion, and 90-degrees inversion [13]. The developed VR environment enables utilization of both, the body-axes spacecraft coordinate system, and inertial coordinate system with instant transition from one into another, switching on/off various virtual reality worlds for multiple supports in the design process [10].

2. Main conventions and abbreviations in this chapter

It would be appropriate to declare some of the main conventions to be used in this work. Firstly, SI units will be employed; therefore, all area moments of inertia will be in kg×m², angular velocity and its components will be in rad/s, kinetic energy will be in Joules (kg×m²×s⁻²), angular momentum will be in kg×m²×s⁻¹ (which is the same as N·m·s), and these units for brevity will be omitted everywhere in the text after the numerical values.

Also, we will often use abbreviation KEE, AME, and AMS for the “kinetic energy ellipsoid”, “angular momentum ellipsoid”, and “angular momentum sphere”, VR for the “virtual reality”. In terms of notations, for brevity, we will use [I] and [ω] to express components of the principal area moments of inertia and angular velocity of the system.

3. “Poinsot’s Construction”: graphical method in attitude kinematics

3.1 Poinsot’s Construction

Various spacecraft systems can be modeled as torque-free systems. For simulation of their motions, including regular spin, flipping, or tumbling, Euler’s equations can be used. This method, involving analytics and numerical methods, will be employed in the later sections. However, basic explanation of the motions can be explained even without solving the Euler’s equations, but just employing the geometrical interpretation of the torque-free motion. This elegant method of representing kinematically the motion of a body is known in classical mechanics, as “Poinsot's construction” (after Louis Poinsot, who published the method in 1834 [14]).

Paying tribute to his amazing discoveries, we wish to present the basic information on this wonderful scientist. Louis Poinsot (January 3, 1777–December 5, 1859) was a French mathematician and physicist. Being passionate about abstract mathematics, he invented the geometrical mechanics. For his numerous contributions, he was elected Fellow of the Royal Society of London in 1858. He died in 1859 at the age of 82. His is one of the 72 names inscribed on the Eiffel Tower: where names of the French scientists, engineers, and mathematicians are engraved in recognition of their contributions [15]. A related fragment from the Eiffel Tower is shown in Figure 1. Symbolically that Poinsot’s name is next to Focault’s name: Poinsot was a contemporary of Léon Foucault, who invented the gyroscope and whose pendulum experiments provided incontrovertible evidence that the Earth rotates.

Figure 1.

Name of Louis Poinsot on Eiffel Tower: (a) the location of the names on the Tower [16]; (b) a fragment of plaque on the South-East side of the Tower with Poinsot’s name [17].

Poinsot’s construction will be shortly reviewed in this section. Torque-free motion is subject to four constants: the kinetic energy of the body and the three components of the angular momentum, expressed in the inertial coordinate system. The angular velocity vector $\mathbf{\omega }$ of the rigid rotor is not constant, but satisfies Euler's equations. Without explicitly solving these equations, Louis Poinsot was able to visualize the motion of the endpoint of the angular velocity vector. He used the conservation of kinetic energy and angular momentum as two constraints on the motion of the angular velocity vector $\mathbf{\omega }$.

3.2 Poinsot’s Kinetic Energy Ellipsoid

Let us first explore the kinetic energy K constraint, writing its analytical expression:

$$K=\frac{1}{2}{I}_{\mathit{xx}}{\omega }_x^2+\frac{1}{2}{I}_{\mathit{yy}}{\omega }_y^2+\frac{1}{2}{I}_{\mathit{zz}}{\omega }_z^2=\mathit{\text{const}}=K_0$$E1

This can be written as

$$I_{\mathit{xx}}{\omega }_x^2+I_{\mathit{yy}}{\omega }_y^2+I_{\mathit{zz}}{\omega }_z^2=2K_0$$E2

which represents an ellipsoid, often called “Poinsot’s ellipsoid (PE)”. It should not be mixed up with the “standard inertial ellipsoid”, which differs from the Poinsot's ellipsoid by the scale factor $\sqrt{2K}$ [18]. However, we further re-write (2) in more convenient, non-dimensional form:

$${\left({\overline{\omega }}_x/a\right)}^2+{\left({\overline{\omega }}_y/b\right)}^2+{\left({\overline{\omega }}_z/c\right)}^2=1$$E3

where

$${\overline{\omega }}_x={\omega }_x/\sqrt{2K_0},{\overline{\omega }}_y={\omega }_y/\sqrt{2K_0},{\overline{\omega }}_z={\omega }_z/\sqrt{2K_0}\text{and}$$E4

$$a=1/\sqrt{I_{\mathit{xx}}},b=1/\sqrt{I_{\mathit{yy}}},c=1/\sqrt{I_{\mathit{zz}}}$$E5

Eq. (3) is a standard form representation of a 3D ellipsoid, called “Poinsot’s kinetic energy ellipsoid (KEE)”, which is shown with its notations in Figure 2a. Each quantity a, b, c is equal to the half of the length of the corresponding principal axis of the ellipsoid. If two of the axes have the same length, then the ellipsoid is an ellipsoid of revolution, also called a spheroid. In this case, the ellipsoid is invariant under a rotation around the third axis, and there are thus infinitely many ways of choosing the two perpendicular axes of the same length. If the third axis is shorter, the ellipsoid is an oblate spheroid; if it is longer, it is a prolate spheroid. If the three axes have the same length, the ellipsoid is a sphere.

Figure 2.

(a) Poinsot’s kinetic energy ellipsoid; (b) Poinsot’s angular momentum ellipsoid.

3.3 Poinsot’s Angular Momentum Ellipsoid

Similar to K, the angular momentum H is another constraint of the system, as for the torque-free motion case has constant squared length, which can be expressed in terms of angular velocity components as follows:

$${\left|\mathbf{H}\right|}^2={\left(I_{\mathit{xx}}{\omega }_x\right)}^2+{\left(I_{\mathit{yy}}{\omega }_y\right)}^2+{\left(I_{\mathit{zz}}{\omega }_z\right)}^2=\mathit{\text{const}}=H_0$$E6

This relationship can be also re-written in terms of the same non-dimensional quantities ${\overline{\omega }}_x,{\overline{\omega }}_y,\text{and}{\overline{\omega }}_z$as follows:

$${\left({\overline{\omega }}_x/A\right)}^2+{\left({\overline{\omega }}_y/B\right)}^2+{\left({\overline{\omega }}_z/C\right)}^2=1$$E7

where

$$A=H_0/\left(I_{\mathit{xx}}\sqrt{2K}\right),B=H_0/\left(I_{\mathit{yy}}\sqrt{2K}\right),C=H_0/\left(I_{\mathit{zz}}\sqrt{2K}\right)$$E8

Eq. (7) represents another 3D ellipsoid, called “angular momentum ellipsoid” (AME), or “momental ellipsoid”, or “ellipsoid of inertia”. With its notations, it is shown in Figure 2b. Note, that the polar reciprocal of AME with regard to its center is another ellipsoid, which is sometimes called the ellipsoid of gyration [19].

3.4 Poinsot’s angular velocity polhodes

After the main notations were introduced, let us consider a particular illustration case of the system with the following parameters: $\left[I_{\mathit{xx}},I_{\mathit{yy}},I_{\mathit{zz}}\right]=\left[2,3,6\right]$ and $\left[{\omega }_x,{\omega }_y,{\omega }_z\right]=\left[0.1,8,1\right]$. Then, system’s constants are $K=99.01$ and $H_0=24.74$ and the key characteristics of the KEE and AME are:

$\left[a,b,c\right]=\left[0.7071,0.5774,0.4082\right]$;$\left[A,B,C\right]=\left[0.8790,0.5860,0.2930\right]$ .

We use the VR tools and show the corresponding KEE and AME surfaces separately in Figure 3a and 3b. In the figure, we also display $\overline{\mathbf{\omega }}$ (scaled vector of the angular velocity) and $\overline{\mathbf{H}}$ (non-dimensional vector of the angular momentum). The scaling factor of 2 is applied to $\overline{\mathbf{\omega }}$ to ensure that the vector is visible, otherwise it would be completely hidden inside of the KEE surface. So, instead of $\overline{\mathbf{\omega }}$ to enable visualization we show $2\times \overline{\mathbf{\omega }}$.

Figure 3.

Illustration of polhodes: (a) KEE; (b) AME; (c) merging KEE and AME via co-centering reveals polhodes; (d) for better visualization, AME is shown semi-transparent.

In Figure 3c we show KEE and AME ellipsoids as merged and co-centered. Figure 3c then instantly reveals that after completing collocation, two ellipsoids intersect along two lines. For better observation of the intersecting surfaces, we also show in Figure 3d the AME as semi-transparent object. This last figure graphically illustrates data in Eq. (3) and Eq. (6), in particular, the difference in the values of a and A, that is, values of the semi-major axes of the KEE and AME.

It can be observed from Figure 3c that the tip of the vector $\overline{\mathbf{\omega }}$ resides on one of the polhodes. In fact, the angular velocity vector $\mathbf{\omega }$ is tracing out the path, called the polhode, which is generally circular or taco-shaped. In the fashion of the day, Poinsot coined the terms polhode and its counterpart, herpolhode, to describe this wobble in the motion of rotating rigid bodies. Poinsot derived these terms from the ancient Greek πόλος (pólos—pivot or end of an axis) + ὁδός (hodós—path or way)—thus, polhode is the “path of the pole” [20].

If the rotating rigid body is symmetric (has two equal moments of inertia), the vector $\mathbf{\omega }$ slides along the sides of a cone (and its endpoint draws a circle). This is known as the torque-free precession of the rotation axis of the rotor.

An observer, attached to the body frame coordinate system, sees the angular velocity vector $\overline{\mathbf{\omega }}$ as sliding along the surface of a cone, called the body cone, whose intersection with the inertia ellipsoid is the polhode. The observer also sees that the vector of the angular velocity is generally circling around one of the principal axes, as illustrated in Figure 4 with three representative cases. In a particular special case of a symmetrical body, the inertia ellipsoid is an ellipsoid of revolution, so that the polhode becomes a circle perpendicular to the spin axis for the coning vector $\overline{\mathbf{\omega }}$.

Figure 4.

Motion of $\overline{\mathbf{\omega }}$ (vector of the scaled angular velocity) in the body axes: (a) around x-axis; (b) around y-axis; (c) around z-axis.

In this work, a special attention will be given to another class of other special cases, when $H^2\approx 2K_0{I}_{\mathit{int}}$, where I_int denotes the intermediate value of the moment of inertia. In these cases, the Poinsot’s polhodes shapes start getting “taco” shapes (shown in Figure 5) while approaching the separatrices and the system motion would represent unstable flipping along the intermediate body axis, if the system is initially provided with significant spin about the intermediate axis.

Figure 5.

“Taco”-shapes of the Poinsot’s polhodes for the cases in which $H^2\approx 2K_0{I}_{\mathit{int}}$.

Polhode can be described analytically and without loss of generality, we perform derivation of the associated equation, assuming that a < b < c. We consider the following parametrization for the xy plane and substitute it into Eq. (7):

$$\left\{\begin{array}{c}{\overline{\omega }}_x=\mathit{A\rho }cos\theta \\ {\overline{\omega }}_y=\mathit{B\rho }sin\theta \end{array}\Rightarrow {\overline{\omega }}_z^2=C^2\left(1-{\rho }^2\right)\right.$$E9

Substitution into Eq. (3) gives:

$${\left(\frac{A\rho cos\theta }{a}\right)}^2+{\left(\frac{B\rho cos\theta }{b}\right)}^2+\frac{C^2\left(1-{\rho }^2\right)}{c^2}=1$$E10

which results in the parametric solution for $\rho$:

$$\rho =\pm \sqrt{\left(1-\frac{C^2}{c^2}\right)/\left(\frac{A^2}{a^2}{cos}^2\theta +\frac{B^2}{b^2}{sin}^2\theta -\frac{C^2}{c^2}\right)}$$E11

With this equation, for each value of $\theta$, the corresponding value of $\rho$ can be calculated and then the values of x, y, and z can be finally determined, using Eqs. (9). Alternatively, all Eqs. (9) can be expressed in terms of the single angular parameter $\theta$:

$$\begin{array}{l}{\overline{\omega }}_x=\pm \frac{Acos\theta }{a}\sqrt{\left(1-\frac{C^2}{c^2}\right)/\left(\frac{A^2}{a^2}{cos}^2\theta +\frac{B^2}{b^2}{sin}^2\theta -\frac{C^2}{c^2}\right)};\\ {\overline{\omega }}_y=\pm \frac{Bsin\theta }{b}\sqrt{\left(1-\frac{C^2}{c^2}\right)/\left(\frac{A^2}{a^2}{cos}^2\theta +\frac{B^2}{b^2}{sin}^2\theta -\frac{C^2}{c^2}\right)};\\ {\overline{\omega }}_z=\pm \frac{C}{c}\sqrt{\left(1-\frac{C^2}{c^2}\right)/\left(\frac{A^2}{a^2}{cos}^2\theta +\frac{B^2}{b^2}{sin}^2\theta -\frac{C^2}{c^2}\right)}.\end{array}$$E12

3.5 Poinsot’s Herpolhode and Invariable Plane

The angular momentum vector $\mathbf{H}$ can be expressed in terms of the moment of inertia tensor $\mathbf{I}$ and the angular velocity scaled vector $\overline{\mathbf{\omega }}$:

$$\mathbf{H}=\left\{\begin{array}{c}{H}_x\\ H_y\\ H_z\end{array}\right\}=\left[\begin{array}{ccc}{I}_{\mathit{xx}}& 0& 0\\ 0& I_{\mathit{yy}}& 0\\ 0& 0& I_{\mathit{zz}}\end{array}\right]\left\{\begin{array}{c}{\omega }_x\\ {\omega }_y\\ {\omega }_z\end{array}\right\}=\mathbf{I}\cdot \mathbf{\omega }=\sqrt{2K_0}\left(\mathbf{I}\cdot \overline{\mathbf{\omega }}\right)$$E13

Similarly, kinetic energy can be also expressed as follows:

$$K_0=\frac{1}{2}\left\{\begin{array}{c}{\omega }_x\\ {\omega }_y\\ {\omega }_z\end{array}\right\}\left[\begin{array}{ccc}{I}_{\mathit{xx}}& 0& 0\\ 0& I_{\mathit{yy}}& 0\\ 0& 0& I_{\mathit{zz}}\end{array}\right]\left\{\begin{array}{c}{\omega }_x\\ {\omega }_y\\ {\omega }_z\end{array}\right\}=\frac{1}{2}\mathbf{\omega }\cdot \mathbf{I}\cdot \mathbf{\omega }=\frac{\sqrt{2K_0}}{2}\overline{\mathbf{\omega }}\cdot \mathbf{H}$$E14

Eq. (14) then leads to the following:

$$\overline{\mathbf{\omega }}\cdot \mathbf{H}=\sqrt{2K_0}$$E15

The dot product of two vectors can be calculated, using (a) projection of $\overline{\mathbf{\omega }}$ on $\mathbf{H}$, giving result ${\overline{\mathbf{\omega }}}_{\mathbf{H}}\times \mid \mathbf{H}\mid =H_0{\overline{\mathbf{\omega }}}_{\mathbf{H}}$ or (b) projection of $\mathbf{H}$ on $\overline{\mathbf{\omega }}$, giving result $\left|\overline{\mathbf{\omega }}\right|\times {\mathbf{H}}_{\overline{\mathbf{\omega }}}$. Both methods lead to the constant value of $\sqrt{2K_0}$; however, method (a) is more preferred, as involves a product of two constants, value of the angular momentum (which is constant due to the “Law of Conservation of angular momentum” for the torque-free system), and, hence, the constant projection $d={\overline{\mathbf{\omega }}}_{\mathbf{H}}$ of the angular velocity scaled vector $\overline{\mathbf{\omega }}$ on $\mathbf{H}$, shown as OD in Figure 6. This choice enabled Poinsot to establish a remarkable plane p, the invariable plane, being perpendicular to the angular momentum $\mathbf{H}$ and at the distance OD from the origin O:

Figure 6.

Invariable plane, normal to the inertia ellipsoid at point P (being a tip of the vector $\overline{\mathbf{\omega }}$) and at the invariable (constant) distance “d” from O: $d=\mid \mathit{OD}\mid =\sqrt{2K}/H_0$.

$$d=\mid \mathit{OD}\mid ={\overline{\mathbf{\omega }}}_{\mathbf{H}}=\sqrt{2K_0}/H_0$$E16

As the rigid body moves, its inertia ellipsoid (which is locked in the body axes system) revolves around the fixed point O in a peculiar way: It rolls on the invariable plane without slipping, with every new point of contact “painting” two continuous curves: polhode on the inertia ellipsoid and herpolhode on the invariable plane, as if it was a double-sided carbon paper between the invariable plane and the ellipsoid.

The herpolhode is always concave to the origin, described by Goldstein [18] as “snakelike”. The herpolhode is generally an open curve, which means that the rotation does not perfectly repeat, but the polhode is a closed curve.

In summary of the presentation of the Poinsot’s construction, we see that in an absolute reference frame, the instantaneous angular velocity vector is the point of intersection between a fixed invariable plane and a kinetic-energy ellipsoid that is tangent to it and rolls around on it without slipping.

Vector of the angular momentum for the torque-free systems is fixed in space; however, angular momentum vector is not. Poinsot’s construction is useful in understanding on how the vector $\overline{\mathbf{\omega }}$ moves. In order to understand the motion of the vector of the angular momentum $\overline{\mathbf{H}}$, another graphical method is to be used. It is based on the Binet’s ellipsoid, named after the French mathematician, physicist and astronomer Jacques Philippe Marie Binet (1786–1856).

4. Analysis of the attitude kinematics of the Inertially Morphed Systems

4.1 Concept of inertial morphing

In the previous publications by the authors [11, 21, 22], we proposed for wide range variations of the attitude motions of the spacecraft to deliberately change its principle moments of inertia during the attitude maneuvering flight. We called this concept as “Inertial Morphing (IM)” and envisaged that it can be achieved in many different ways via variety of methods/principles, including reposition of control masses with smart mechanisms or via components docking/undocking; redistribution of the magnetic control liquids with magnetic fields; rejection of controlled liquids via ejection/ablation, etc. Efforts of the IM systems can be also enhanced via specifically designed motions of the spacecraft appendages to assist in attitude control. We believe that the novel concept of IM enables design and construction of the inertially morphed spacecraft, possessing acrobatic capabilities, and may allow design of new class of gyroscopic systems with a “sense” of time.

The concept of IM was prompted by the observation of rigid bodies in space performing the intriguing flipping motion. The phenomenon was first observed in space by famous US scientist-astronaut Owen Kay Garriott on-board Skylab in 1973, but later was publicized as “Dzhanibekov’s effect” named after V. Dzhanibekov, who observed the flipping of a spinning wing nut during his space flight in 1985 [23]. Interestingly, the resulting flipping motions of rigid bodies with various shapes can be explained by the L.Euler’s equations, derived in 1785 for rigid bodies with fixed area moments of inertia.

We wish to utilize and expand the functionality of the Poinsot’s graphical method. As per the IM concept, we will be specifically emphasizing on the variation of the principal moments of inertia (i.e., implying that $\left[\mathbf{I}\right]=\left[\mathbf{I}\left(t\right)\right]$), not initially envisaged and explored by L. Poinsot and L. Euler in their works, we will need to treat moments of inertia as variables, expand classical Euler’s equation of motion, and involve another graphical construction, based on Binet’s construction, which is mostly suitable for the analysis of morphed systems, as is centered around the angular momentum, involving $\left[\mathbf{I}\right]$.

4.2 Definition of the AMS and Binet’s KEE

For the future analysis and interpretation of results, let us introduce the following non-dimensional coordinates ${\overline{H}}_x$ , ${\overline{H}}_y$ and ${\overline{H}}_z$, associated with x, y, and z body axes:

$$\begin{array}{l}{\overline{H}}_x\left(t\right)=H_x/H_0;{\overline{H}}_y\left(t\right)=H_y/H_0;{\overline{H}}_z\left(t\right)=H_z/H_0,\\ \text{where}{H}_x=I_{\mathit{xx}}{\omega }_x;H_x=I_{\mathit{xx}}{\omega }_x;H_x=I_{\mathit{xx}}{\omega }_x\\ \text{and}{H}_0=\sqrt{H_x^2+H_y^2+H_z^2}=\sqrt{{\left(I_{\mathit{xx}}{\omega }_x\right)}^2+{\left(I_{\mathit{yy}}{\omega }_y\right)}^2+{\left(I_{\mathit{zz}}{\omega }_z\right)}^2}\end{array}$$E17

We assume that the spacecraft system is in free flight. Then the fundamental “Law of Conservation of angular momentum” can be applied to the system and this law can expressed in a compact form:

$${{\overline{H}}_x}^2+{{\overline{H}}_y}^2+{{\overline{H}}_z}^2=1$$E18

This equation can be interpreted graphically as a unit sphere (i.e., simply a sphere of radius one), shown in Figure 7a.

Figure 7.

Surfaces for attitude dynamics: (a) unit radius angular momentum sphere (AMS) in the non-dimensional coordinates; (b) example of kinetic energy ellipsoid (KEEs) for the system with a_x = 1.2; a_y = 1.4; a_z = 0.8; (c) KEE for a_x = 0.8; a_y = 1; a_z = 1.2 case.

Kinetic energy of the system can be written in terms of the area moments of inertia and components of the angular velocity:

$$K\left(t\right)=\frac{1}{2}{\left(I_{\mathit{xx}}{\omega }_x\right)}^2+\frac{1}{2}{\left(I_{\mathit{yy}}{\omega }_y\right)}^2+\frac{1}{2}{\left(I_{\mathit{zz}}{\omega }_z\right)}^2$$E19

This expression can be further re-written as follows:

$$K\left(t\right)={\left[\frac{H_x\left(t\right)}{\sqrt{2I_{\mathit{xx}}\left(t\right)}}\right]}^2+{\left[\frac{H_y\left(t\right)}{\sqrt{2I_{\mathit{yy}}\left(t\right)}}\right]}^2+{\left[\frac{H_z\left(t\right)}{\sqrt{2I_{\mathit{zz}}\left(t\right)}}\right]}^2.$$

In view of the introduced non-dimensional coordinates ${\overline{H}}_x$, ${\overline{H}}_y$, and ${\overline{H}}_z$, expression for $K\left(t\right)$, similar to $H_0$, can be also expressed in the concise form:

$${\left(\frac{{\overline{H}}_x}{a_x}\right)}^2+{\left(\frac{{\overline{H}}_y}{a_y}\right)}^2+{\left(\frac{{\overline{H}}_z}{a_z}\right)}^2=1$$E20

where

$$a_x=\frac{\sqrt{2K\left(t\right)I_{\mathit{xx}}\left(t\right)}}{H_0};a_y=\frac{\sqrt{2K\left(t\right)I_{\mathit{yy}}\left(t\right)}}{H_0};a_z=\frac{\sqrt{2K\left(t\right)I_{\mathit{zz}}\left(t\right)}}{H_0}$$E21

Eq. (20) corresponds to the so-called Binet’s KEE [18], also fixed to the body axes non-inertial coordinate system, as Poinsot’s AME and KEE ellipsoids.

However, Binet’s KEE, given by Eq. (20) is constructed in the non-dimensional coordinates ${\overline{H}}_x$, ${\overline{H}}_y$, and ${\overline{H}}_z$, and has the semi-major axes, equal to a_x, a_y, and a_z. In the applications of the IM concept, any combinations of a_x, a_y, and a_z may occur. In Figure 7b and 7c, we present two contrast cases of the kinetic energy ellipsoids (KEEs). In the first example, a_z < a_x < a_y, however in the second example, the ascending order of the values of the semi-major axes is different: a_x < a_y < a_z. The axis, associated with the smallest (out of a_x, a_y and a_z) value, will be called “minimum axis” of inertia. The axis, associated with the largest (out of a_x, a_y and a_z) value, will be called “maximum axis” of inertia and the third axis will be called “intermediate axis” of inertia. Special attention should be paid to the identification of the intermediate axis, as if the system is provided with the spin about this axis, the resultant motion will be unstable and the system would start flipping motion. In relation to two cases in Figure 7b and 7c, for initiation of the unstable flips, the system (b) should be provided with the predominant spin about the x axis, whereas the second system would start flips if the predominant spin is provided about the y axis.

4.3 Co-centering AMS and KEE in a numerical example

Let us consider an illustrative numerical example of the system with the following moments of inertia and initial components of the angular velocity:

$$I_{\mathit{xx}}=2;I_{\mathit{yy}}=4;I_{\mathit{zz}}=6;{\omega }_{x0}=0.1;{\omega }_{y0}=1;{\omega }_{z0}=1$$E22

These parameters enable to use Eq. (17), and Eq. (19) to determine $H_0=7.21\mathrm{kg}\times {\mathrm{m}}^2/\mathrm{s}$ and $K=5.01\mathrm{J}$, then use Eq. (21) to calculate a_x = 0.62; a_y = 0.88; a_z = 1.07 and therefore identify “y” axis as an axis, corresponding to the intermediate area moment of inertia. The corresponding KEE can be constructed and plotted together with the co-centered, collocated AMS (shown as semi-transparent surface), as shown in Figure 8. From the figure, we can clearly see that the largest portion of the KEE is inside the AMS; however, top and bottom portions of the KEE are extruding outside the AMS. Figure 8 enables us to observe two symmetrical intersection lines between KEE and AMS, painted with white color. Note that only the top intersection line, polhode, is seen in the static Figure 8, whereas the bottom polhode is obscured with the KEE. Their size and shape reflect the nature of the attitude motion of the system; therefore, polhodes are very important visual integrated identifiers of the rotation of the system. By just viewing the shape, the designer, for example, could instantly determine if the system is in complex tumbling, unstable flips or regular spin, could also determine the “distribution” of rotation between the body axes.

Figure 8.

Co-centered, collocated angular momentum sphere and kinetic energy ellipsoid presented in the virtual reality interface for the illustration case.

Figure 8 also presents with a red bottom dot an initial state of the analyzed system. Indeed, Eq. (17) can be used to calculate non-dimensional initial coordinates of the system:

$${\overline{H}}_{x0}=0.0277;{\overline{H}}_{y0}=0.5545;{\overline{H}}_{z0}=0.8317.$$E23

It would be important to note that there may be infinite number of combinations of other initial conditions, resulting in their initial state points to reside on the same polhode. As an example, we show with the left red dot in Figure 8 another set of system’s initial conditions, placing the state point on the same polhode. Similar to (23), we present the other set of initial conditions with the following numbers:

$${\overline{H}}_{x0}=0.2217;{\overline{H}}_{y0}=0.3374;{\overline{H}}_{z0}=0.9149.$$E24

For illustration purpose, in Figure 9 we represent two different systems in a single combined plot with their co-centered KEEs and identical AMSs. The feature of this case is that they both have identical separatrices, shown with red- and yellow-coinciding lines. The condition for the separatrices is: H² = 2KI_yy; however, selection of the initial conditions on two sides of the separatrix does not necessary mean that the switch from H² > 2KI_yy to H² < 2KI_yy (or vice versa) would occur. In this respect, for both, “A” and “B” systems, the same relationship $H^2>2K{I}_{\mathit{yy}}$ is satisfied.

Figure 9.

Two co-centered conjugated systems with I_A = [3,4,9], ${\omega }_A=\left[2,1,2\right]$ and I_B = [8, 6.53, 5], ${\omega }_B=\left[-15,-15,-12\right]$, which have the same separatrices: (a) the identical AMS shapes are shown as opaque (non-transparent surfaces); (b) the AMSs are shown as semi-transparent, revealing red KEE for the first system and yellow KEE for the second system.

4.4 Torque-free Euler’s equations of motion of the rigid body

The famous Euler’s equations, describing attitude dynamics of the rigid body, can be written on the matrix form:

$$\left[\begin{array}{ccc}{I}_{\mathit{xx}}& 0& 0\\ 0& I_{\mathit{yy}}& 0\\ 0& 0& I_{\mathit{zz}}\end{array}\right]\left\{\begin{array}{c}{\dot{\omega }}_x\\ {\dot{\omega }}_y\\ {\dot{\omega }}_z\end{array}\right\}=\left\{\begin{array}{c}\left(I_{\mathit{yy}}-I_{\mathit{zz}}\right){\omega }_y{\omega }_z\\ \left(I_{\mathit{zz}}-I_{\mathit{xx}}\right){\omega }_z{\omega }_x\\ \left(I_{\mathit{xx}}-I_{\mathit{yy}}\right){\omega }_x{\omega }_y\end{array}\right\}$$E25

These equations can be solved numerically for any combination of I_xx, I_yy, I_zz, and initial conditions, including, for example, Eq. (23) and/or Eq. (24). Then, the time histories for ${\omega }_x$, ${\omega }_y$, and ${\omega }_z$ can be determined, results presented in the non-dimensional coordinates ${\overline{H}}_x$, ${\overline{H}}_y$, and ${\overline{H}}_z$ and eventually plotted as trajectories of the tip of the angular momentum vector $\overline{\mathbf{H}}$. It is remarkable that the solution trajectory would exactly correspond to the polhode, that is, intersection line between the KEE and AMS surfaces. Therefore, one of the aims of the dynamics analyst may be determination of the shape of the polhode of the system, which would enable fast classification of the attitude motion by visual inspection.

By the way, Eq. (25) can be used to determine direction of motion of the vector $\overline{\mathbf{H}}$. For example, considering a red point in Figure 8, symbolizing the tip of the angular momentum vector, located in the vicinity of x = 0 plane, we set ω_x = 0 and reduce Euler’s equations to one equation only:

$$I_{\mathit{xx}}{\dot{\omega }}_x=\left(I_{\mathit{yy}}-I_{\mathit{zz}}\right){\omega }_y{\omega }_z$$E26

Then, sliding slightly the point from the x = 0 plane along polhode, to the position, where ω_y > 0 and ω_z > 0 (i.e., to the xyz-positive octant), we determine that ${\dot{\omega }}_x<0$ , as I_yy<I_zz for the selected example. Therefore, red dot point must move in anti-clock direction around z axis (determined with the right-hand rule), to ensure reduction in the ω_x value with time.

4.5 Analytical expressions for Polhodes

Without the loss of generality, let us assume that a < b < c. For this case, the xz cross section of the co-centric AMS and KEE is presented in Figure 10.

Figure 10.

Cross section of the KEE and AMS with the xz plane for two particular cases of the systems with: (a)-(b) a_x = 0.8; a_y = 1; a_z = 1.2; (c)-(d) a_x = 1.4; a_y = 1.02; a_z = 0.8.

For the analysis of the intersection of the AMS and KEE along polhodes, we introduce the following parametrization in the xy plane, being an orthogonal plane to the plane of the xz cross section:

$${\overline{H}}_x=rcos\theta ;{\overline{H}}_y=rsin\theta$$E27

In Eq. (27), the angle $\theta$ and radius $r$ are the polar coordinates in the xy plane. In view of Eq. (27) and Eq. (18), describing the AMS, we can deduct that

$${\overline{H}}_z=\sqrt{1-r^2}$$E28

Substitution of Eqs. (27) and (28) into Eq. (20) for the KEE enables us to re-write it in terms of the polar coordinates as follows:

$${\left(rcos\theta /a_x\right)}^2+{\left(rsin\theta /a_y\right)}^2+{\left(\sqrt{1-r^2}/a_z\right)}^2=1$$E29

It can be solved for the radius $r$ becoming a function of only angle $\theta$, assuming that the system’s dynamics characteristics are all known:

$$r=\pm \sqrt{\frac{1-{\left(1/a_z\right)}^2}{{\left(cos\theta /a_x\right)}^2+{\left(sin\theta /a_y\right)}^2-{\left(1/a_z\right)}^2}}$$E30

Running a cycle for all polar angle $\theta$ in the range between 0 and 360 degrees, corresponding values of $r$ can be determined, using Eq. (30). And then, for each matching pair of angle $\theta$ and radius $r$, the triplets of corresponding coordinates ${\overline{H}}_x$, ${\overline{H}}_y$, and ${\overline{H}}_z$ of the polhodes can be determined, using Eqs. (27-28).

Eqs. (27, 28, 30) can be used to plot polhodes in the upper and lower quarters of the AMS for the illustrative system in Figure 11a. To plot polhodes in the left and right quarters, in Eq. (30)x- and z-related parameters should be swapped.

Figure 11.

(a) Plotted polhodes and separatrices on the AMS; (b) notations for separatrices.

Polhodes are seen on orthogonal projections as ellipses and hyperbolas. As for two separatrices, subdividing polhodes into four groups, they are seen as ellipses or as X-shaped two lines, illustrated in Figure 11b with the view from the intermediate axis y.

Dihedral angle for the separatrix can be determined, using the following Eq. (31):

$${\gamma }_{\mathit{yz}}=arctan\frac{{\overline{H}}_x}{{\overline{H}}_z},\text{where}\frac{{\overline{H}}_x}{{\overline{H}}_z}=\sqrt{\frac{I_{\mathit{xx}}\left(I_{\mathit{zz}}-I_{\mathit{yy}}\right)}{I_{\mathit{zz}}\left(I_{\mathit{yy}}-I_{\mathit{xx}}\right)}}$$E31

In the literature, it is a common approach to present Poinsot’s and Binet’s constructions separately. However, as we wish to extract maximum advantages from both methods, we produce a co-centered combination of the conjugated virtual reality worlds (see Figure 12a and 12b), corresponding to these two graphical methods and constructed for the system with $\left[I\right]=\left[2,4,6\right]$ and $\left[{\omega }_0\right]=\left[0.1,1,1\right]$. In Figure 12c, we also show the xz-cross section of these worlds. It is interesting to visually observe orthogonal directions of the semi-major axes in IE and KEE and also relative scaling of the vectors of angular velocity and angular momentum in two constructions. With the scaling adopted, the Poinsot’s construction is at the lower level, compared to the Poinsot’s construction, and, figuratively speaking, is “under the bonnet” in the combined methods setup.

Figure 12.

Combined Poinsot’s and Binet’s constructions: (a) layout for the conjugated VR worlds with semi-transparent KEE, AMS, AME; (b) obscure KEE, AMS; (c) xz-cross section.

5. Analysis of the unstable flipping motion

5.1 Calculation of the period of the tumbling and flipping motions

There are analytical expressions for the period of the periodic motion of the torque-free rotor; however, their notations depend upon the order of the principal moments of inertia. Let us consider first them for the $I_{\mathit{xx}}<I_{\mathit{yy}}<I_{\mathit{zz}}$ case, which, in turn, has two sub-cases:

• If $H^2>2K_0{I}_{\mathit{yy}}$ , which is equivalent to $a_y<1$ , then

$$T=4\mathcal{K}\sqrt{\frac{I_{\mathit{xx}}{I}_{\mathit{yy}}{I}_{\mathit{zz}}}{\left(I_{\mathit{zz}}-I_{\mathit{yy}}\right)\left(H^2-2K_0{I}_{\mathit{xx}}\right)}}$$E32

• If $H^2<2K_0{I}_{\mathit{yy}}$ , which is equivalent to $a_y>1$ , then

$$T=4\mathcal{K}\sqrt{\frac{I_{\mathit{xx}}{I}_{\mathit{yy}}{I}_{\mathit{zz}}}{\left(I_{\mathit{xx}}-I_{\mathit{yy}}\right)\left(H^2-2K_0{I}_{\mathit{zz}}\right)}}$$E33

where $\mathcal{K}=\mathcal{K}\left(k\right)$ is complete elliptic integral of the first kind:

$$\mathcal{K}\left(k\right)=\underset{0}{\overset{1}{\int }}\frac{\mathit{ds}}{\sqrt{\left(1-s^2\right)\left(1-k^2{s}^2\right)}}=\underset{0}{\overset{\pi /2}{\int }}\frac{\mathit{du}}{\sqrt{1-k^2{sin}^{2}u}}$$E34

being a function of the parameter $k\left(-1\le k\le 1\right)$, often called “elliptical modulus” and in the current task formulation is to be calculated as follows:

$$k=\sqrt{\frac{\left(I_{\mathit{yy}}-I_{\mathit{xx}}\right)\left(2K_0{I}_{\mathit{zz}}-H^2\right)}{\left(I_{\mathit{xx}}-I_{\mathit{yy}}\right)\left(H^2-2K_0{I}_{\mathit{xx}}\right)}}$$E35

The elliptic integral Eq. (34) can be calculated, using various on-line calculators, sometimes providing amazing accuracy of up to 50 digits, if needed [24]. In many other cases, where this high precision is not a requirement, elliptic integral $\mathcal{K}$ can be calculated, using MATLAB® and Wolfram MATHEMATICA®; however, a special care should be given to the specifications of the relevant specialized functions to be employed. For example, $\mathcal{K}$ can be calculated, using MATLAB command ellipke or ellipticK command [25]; however, instead of $k$ to be used in KEISAN calculator, $k^2$ (i.e., the squared value of $k$) should be provided as argument in these commands in MATLAB: this is due to the specific definition of these commands in MATLAB, presented for completeness below:

$${\left.\mathcal{K}\left(m\right)\right|}_{\mathit{\text{MATLAB}}}=\underset{0}{\overset{\pi /2}{\int }}\frac{\mathit{d\theta }}{\sqrt{1-m{sin}^2\theta }},$$E36

Note that the elliptical modulus k or, in other terminology, the modular angle α are related to the parameter m as follows:

$$m=k^2=sin2\mathrm{\alpha }$$E37

In view of the different formats for the arguments in elliptic integral calculation, a series of simple verification tests may be advisable. Reconciliation of the specification differences is illustrated in Figure 13 and the user can reproduce conforming results, using the following:

• KEISAN, by entering test input of k=0.5 and requesting 22-digit accuracy: 1.685750354812596042871

• MATLAB: by calling any of the commands ellipticK(0.25) or ellipke(0.25), where 0.25=0.5² (with switch to the “long format”): 1.685750354812596

Figure 13.

Calculation of the complete elliptic integral of the first kind for the illustration example: (a) KEISAN definition with “k” argument; (b) MATLAB definition with “m” argument.

In our case, we have high precision calculation of the period programmed in Java and embedded in the VR, enabling to get this information interactively straight in the virtual reality environment and display it in real time with any inertial morphings, applied to the conjugated concurrent configurations of the system.

The VR interface system, programmed by authors, is illustrated in Figure 14.

Figure 14.

Two conjugated configurations of the spacecraft, corresponding to the red and yellow kinetic energy ellipsoids KEEs.

Figure 14 displays two initial configurations of the spacecraft, using AMS and KEEs. The first, corresponding to the red KEE, has the following characteristics: I_xx = 3.24, ω_x0 = 0.49, I_yy = 4.63, ω_y0 = 1, I_zz = 5.34, ω_z0 = 1. Instead of showing traditional vector of the angular momentum $\overline{\mathbf{H}}$, for visualization purposes, a magnified vector $2\overline{\mathbf{H}}$ is shown with the light blue color. Magnification should be applied, otherwise vector $\overline{\mathbf{H}}$, without magnification, would not be seen, as would be completely residing inside the AMS.

Parameters of the first configuration of the system are in compliance with $H^2>2K_0{I}_{\mathit{yy}}$ condition.. The period of its tumbling motion, therefore, can be calculated, using Eq. (32), and is equal to T=19.01 s. Another configuration of the spacecraft, corresponding to the yellow KEE (shown as semi-transparent surface), has the following characteristics: I_xx =3.36, ω_x0 = 10, I_yy =6, ω_y0 = 0.1, I_zz = 2, ω_z0 = 0.1. This configuration is in compliance with $H^2<2K_0{I}_{\mathit{yy}}$ condition; therefore, for the calculation of the period of the motion, Eq. (34) should be used. Its application results in T=3.84 s. It may be of interest to notice that both KEEs have no intersections, but have four touching points, one of which is shown with the red dot, corresponding to the common set of the angular velocity components: [ω] = [0, 1.22, 0.85].

In our case, we have high precision calculation of the period programmed in Java and embedded in the VR, enabling to get this information interactively straight in the virtual reality environment and display it in real time with any inertial morphings, applied to the conjugated concurrent configurations of the system.

5.2 Differential equations of the torque-free spacecraft with Inertial Morphing

In order to be able to simulate the spacecraft with IM capabilities, the Euler’s equations should be expanded, as follows, to allow change of the inertial properties of the system:

$$\left[\begin{array}{cccccc}{I}_{\mathit{xx}}& 0& 0& 0& 0& 0\\ 0& I_{\mathit{yy}}& 0& 0& 0& 0\\ 0& 0& I_{\mathit{zz}}& 0& 0& 0\\ 0& 0& 0& sin\theta sin\varphi & cos\varphi & 0\\ 0& 0& 0& sin\theta cos\varphi & -sin\varphi & 0\\ 0& 0& 0& cos\theta & 0& 1\end{array}\right]\left\{\begin{array}{c}{\dot{\omega }}_x\\ {\dot{\omega }}_y\\ {\dot{\omega }}_z\\ \dot{\psi }\\ \dot{\theta }\\ \dot{\varphi }\end{array}\right\}=\left\{\begin{array}{c}\left(I_{\mathit{yy}}-I_{\mathit{zz}}\right){\omega }_y{\omega }_z-{\dot{I}}_{\mathit{xx}}{\omega }_x\\ \left(I_{\mathit{zz}}-I_{\mathit{xx}}\right){\omega }_z{\omega }_x-{\dot{I}}_{\mathit{yy}}{\omega }_y\\ \left(I_{\mathit{xx}}-I_{\mathit{yy}}\right){\omega }_x{\omega }_y-{\dot{I}}_{\mathit{zz}}{\omega }_z\\ {\omega }_x\\ {\omega }_y\\ {\omega }_z\end{array}\right\}$$

where ψ, θ, and ϕ are 313 Euler angles [26].

5.3 Application of IM for stopping flipping motions or their initiation

In our earlier publications [9, 11], we suggested a method of switching between stable and unstable motions of the spinning systems. Unstable motion is possible only when the predominant rotation of the system is about its intermediate axis. Therefore, for the flipping systems, in order to stop flips, we suggested to use a method of deliberate change of the system’s moments of inertia in a way to transform the axis of the predominant spin into the minimum or maximum inertia axis. Therefore, for stopping the flips, if they are not desired types of motion, there could be two strategies, illustrated in Figure 15.

Figure 15.

Conceptual strategies (with multiple solutions) of: (a) stopping unstable flipping motion of the system about its intermediate axis of inertia, y-axis; (b) initiation of unstable flips about y-axis for two cases, when I_yy is initially the minimum or the maximum axis.

In a similar way, the flipping motion can be initiated on a system, being initially in a regular stable spin. This could be only possible if the axis of the predominant rotation is initiated about the minimum or maximum principal axis. For initiation of the spin, when needed, IM can be applied to transform the axis of predominant rotation into the intermediate axis.

Pre-planning of the IM maneuver can be done with the VR interface, which has numerous features in one. It enables to include Poinsot’s and Binet’s constructions, display several designs concurrently, interactively review them in real time, and select the solution, satisfying key requirements. The starting configuration is set with the red sliders, evoking calculation of the periods of the motions, visualization of the IEs, AMEs, KEEs, AMS, patterns of the motion of the angular momentum and angular velocity vectors, polhodes, herpolhodes, separatrices, alternative separatrices, etc. Two additional conjugate configurations are also available for activation, allowing to foresee attitude motion of the system after the application of various inertial morphings.

One of the interface windows is displayed in Figure 16, where just for the illustration purposes, the KEEs for the conjugated variants are presented as semi-transparent.

Figure 16.

Conjugated VR worlds, combining Poinsot’s and Bine’s constructions and provisions for three concurrent morphed variants of the system: (a) interface; (b-e) review stages.

The interface illustrates possibility to design a scenario of stopping the flipping motion. We consider the system which has the following initial conditions: [I]=[2, 2.5, 3] and [ω]=[.1, 10, .1], corresponding to the flipping motion with the period of T=12 s.

There could be multiple solutions to the stabilization of the spacecraft. As illustration of the IM, we firstly use strategy-1 in Figure 15a and apply at the instant t=12.4 s the only one single prompt morphing, enabling change of the _iI_yy of the system to its new value of _fI_yy=3.2. The results of the numerical simulations, partially presented in Figure 17a, confirm that with the applied IM, the system was transferred into the almost regular spin about the y-axis. This simulation is in agreement with the Law of Conservation of angular momentum, requiring and achieved increase of the dominant spin rate as per the (_iI_yy/ _fI_yy) ration: 10 × (_iI_yy/ _fI_yy) = 10 × (2.5/3.2) ≈ 8. Furthermore, as illustration of strategy-2 we apply only one single prompt morphing, enabling change of the _iI_yy of the system to its new value of _fI_yy = 1.8. The results of the numerical simulations, partially presented in Figure 17b, confirm that with the applied IM, the system was transferred into the almost regular spin about the y-axis. This simulation is in perfect agreement with the Law of Conservation of angular momentum, requiring and achieved increase of the dominant spin rate as per the (_iI_yy/ _fI_yy) ratio: 10 × (_iI_yy/_fI_yy) = 10 × (2.5/1.8) ≈ 14. With the first solution, the KEE of the system “shrinks” getting smaller than the co-centered AMS with only two small touching spots on the y-axis (see Figure 17c). With the second solution, the KEE “expands” to completely embrace the AMS, with two small touching spots on the y-axis (see Figure 17d).

Figure 17.

Stopping flipping motions of the system with only one inertial morphing: (a) Strategy-1 with y becoming minimal axis of inertia; (b) Strategy-2 with y becoming maximum axis; (c) and (d) transformations of the KEEs, corresponding to Strategies-1 and 2.

Various IM maneuvers can be combined together. An example of the possible practical application for the spacecraft is presented in Figure 18, showing 180-degrees inversion of the spacecraft. This enables utilization of the same thruster for boost and for braking. After conceptual design of the sequence of the IMs, mainly performed in the body axes, the VR interface can be further used to observe the attitude maneuver of the spacecraft in the inertial coordinates.

Figure 18.

Representation of the 180° inversion of the spacecraft, allowing use of the same thruster for boost and braking: (a–f) stages of the flight; (g) performing flip in VR.

5.4 Conjugated Virtual Reality assist de-tumbling of the spacecraft

For the system with three distinct principal area moments of inertia, there could be six different cases: I_xx < I_yy < I_zz; I_xx < I_zz < I_yy; I_yy < I_xx < I_zz; I_yy < I_zz < I_xx; I_zz < I_xx< I_yy; I_zz < I_yy < I_xx and the format of the equations for polhodes and separatrices would depend upon the order of the moments of inertia. Instead of re-writing relevant equations for various cases, we will introduce a universal statement, employing the subscripts below, which distinguish minimal, intermediate, and maximum moments of inertia. In these notations, the dihedral angle $\gamma$ of the separatrix, measured from the axis with maximum moment of inertia to the axis with minimal moment of inertia (note: order${}^{"}max\to {min}^{"}$ is important), can be expressed as follows [21]:

$$tan\gamma ={tan}^{-1}\left(\frac{{\overline{H}}_{\min }}{{\overline{H}}_{\max }}\right),\text{where}\frac{{\overline{H}}_{\min }}{{\overline{H}}_{\max }}=\sqrt{\frac{I_{\min }\left(I_{\max }-I_{\text{intermed}}\right)}{I_{\max }\left(I_{\text{ntermed}}-I_{\min }\right)}}$$E38

Therefore, if at the particular instant of time, we have a tip of the vector $\overline{\mathbf{H}}$, marked as dot on the AMS with current coordinates ${\overline{H}}_{\mathit{xC}}$, ${\overline{H}}_{\mathit{yC}}$, and ${\overline{H}}_{\mathit{zC}}$, then we can use Eq. (38) to calculate $tan\gamma$ and then can apply inertia morphing in such a way that the point would become on one of the separatrices of the morphed system. With instantly applied morphing, the dot would be “intercepted” by the separatrix and would start moving along the newly “constructed railway”, instead of moving along its previous path. We called this method “insertion into separatrice” [23]. The easiest way to implement intercepts would be to assume them occurring when the dot is at one of the points of intersection of the dot’s polhode with one of the xy, yz, or xz planes. To get the required paramters for the “insertion”, we can assign any new values of the minimal and maximum moments of inertia (I_min and I_max) and then use Eq. (38) to calculate the value of the intermediate moment of inertia, required for the transfer to occur:

$$I_{\text{intermed}}=I_{\min }{I}_{\max }\left(1+{tan}^2\gamma \right)/\left(I_{\min }+I_{\max }{tan}^2\gamma \right)$$E39

The developed VR interface is very efficient in creating and “interviewing” various solutions and then selection of the most suitable maneuver scenario. Let us consider a particular numerical example. In Figure 19a we show the system with [I] = [3.44, 6, 3], [ω₀] = [2.46, 1.44, 0.96], described with red KEE. Let us “intercept” the point, moving along its polhode, at the xz-plane and initiate its flipping motion along the “y” axis. For the point on its original polhode we have

Figure 19.

Two competing variants of the de-tumbling of the spacecraft presented side-to-side: (a–b) maneuver planning interface; (c–d) simulation raw results; (e–f) animation in VR.

$$tan\gamma =H_{\mathit{zC}}/H_{\mathit{yC}}=0.6313/0.7756=0.8139.$$

Therefore, further assuming that the hardware could ensure new morphed values of the minimum and maximum moments of inertia I_yy = 6; I_zz = 4, we can calculate the required new value of the intermediate moment of inertia I_xx = 5.0032. So, the set of the moments of inertia after the first IM is: [I]_1B = [5, 6, 4]. However, the VR interface would immediately produce the same parameters, with quick interactive selection of the I_yy and I_zz. As an alternative to this solution, if the new values of the minimum and maximum moments of inertia would be selected as I_yy = 6; I_zz = 3, then we would get I_xx = 4.2904: [I]_1B = [4.29, 6, 3]. Also, the third illustration solution would be selection of the following values for the first IM: [I] _1C = [4.58, 4.33, 5]. Interestingly, cases A and C have the same separatrices, but have significant difference in execution time, which can be important for the planning of the maneuvre. In planning of the case A, red and yellow KEEs conjugates were employed; in planning of the case C, red and green conjugates were used. In the Figure 19, some of the surfaces were interactively made semi-transparent. In Figure 19e and 19f we show the spacecraft VR world, merged together with the AMS+KEE VR world. For convenience, their motion can be synchronised or disconnected, enabling for the user to observe AMS-KEE in both body and inertial coordinate systems. In addition, change of the viewpoints enables observations of the spacecrat from its body axes or from inertial reference set.

To complete stabilization of the spacecraft, application of the second IM is needed: it should be applied at the moment when the vector $\overline{\mathbf{H}}$ is passing one of the ‘y’ axis AMS poles. As illustration, we select new values after the second IM applied: [I]_A,B,C = [6.5, 6, 4].

The VR planning maneuver interface, resultant time histories of the angular velocities for the two-morphing meneuvers (A and C), and the screen snapshots from the VR animations are shown in Figure 19.

5.5 Conjugated Virtual Reality assist 90° inversion of the spacecraft

In this subsection we will illustrate the application of the IM for the 90-degrees spacecraft inversion with the sequence shown in Figure 20. For the planning, the VR system is used. As our system with conjugated VR worlds allows superposition of the various morphed configurations for the spacecraft, we can select two separatrices: the first—coming from the initial “destination” of the tip of the vector $\overline{\mathbf{H}}$ and the second—coming to the axis, to which the transfer of the rotation is to be passed. Their intersection constitutes the point, where instantaneous IM should be applied. In the interface in Figure 21, the first separatrix corresponds to the red KEE, and the second—to the yellow KEE.

Figure 20.

Illustration of the 90-degrees inversion of the spacecraft, aiming to transfer its rotation from one axis to another axis at 90 degrees: (a–e) stages of the maneuver.

Figure 21.

Illustration of the utilization of the conjugated VR worlds for planning spacecraft maneuver: (a) selection of the “separatrix to separatrix” trajectory; (b) making AMS and KEEs highly transparent for the visualization of the transfer; (c) time history of components of the angular velocity; (d) time history of 313 Euler angles.

Remarkably, only three IMs would be required for implementation of such complex maneuver. One possible example of the parameters for the 90-degrees inversion is presented in Table 1.

IM index	t, s	Ixx	Iyy	Izz
0	0	3	2.8	5
1	10	3	3.1	5
2	18.25	4.5	3.1	5
3	25.15	5.5	3.1	5
f	35.14	5.5	3.1	5

Table 1.

IM parameters for the 90-degrees spacecraft inversion.

Figure 21 presents results of the simulation of the selected case. It can be seen that the spacecraft, which initially was in stable spin about its longitudinal axis, after three IMs was transferred to the rotation about its lateral axis. This is confirmed not only with the plots of the components of the angular velocity in Figure 21c, but also with the plots of the Euler angles Figure 21d.

Using similar principle, it is possible to design cascaded maneuvers. For example, in [12] we designed the sequence, where the spacecraft performed all-axes inversions “parade”. With the simplicity of the required control, involving only exiguous/paltry number of IM control actions, it is believed that the systems with IM could find wide application in autonomous small spacecraft.

6. Conclusions

This book chapter is further exploring application of the “inertial morphing” control principles and is dedicated to the design of the integrated VR platform, enabling design, test, optimization, and visualization of various missions for the autonomous spacecraft.

The key feature of this work is introduction of several VR worlds for effective evolutionary visualization of the attitude dynamics of the spacecraft. The demonstrated set of the conjugated VR worlds include the following: “Inertial VR World”, “Body-Axes Non-Inertial VR World”, “Inertial Geometric Interpretation VR Worlds”, “Body-Axes Non-Inertial Geometric Interpretation VR World”. This integrated virtual reality environment of several conjugated VR worlds proved to be an efficient tool for the design of the autonomous spacecraft missions with specific requirements of re-orientation, acrobatic attitude maneuvering, for simulation of the competing scenarios, their optimization, and animated visualization. Most significantly, for the inertially morphed systems, it enables design of the astonishingly simple, limited number discrete control inputs. For example, it enables design of the inversion of the spacecraft with two only instantaneous control adjustments. And the VR environment also enables design of the spacecraft de-tumbling, as well as 90-degrees inversion of the spacecraft with only three instantaneous control adjustments.