A Robust Control Algorithm of a Descent Vehicle Angular Motion in the Earth’s Atmosphere

Abstract

A new approach to synthesize a robust controller for the angular motion of the Earth lander by decomposition method of output modal control is proposed. A universal analytical solution for the problem of stabilizing the angular position of the lander is obtained. A comparative analysis of the presented algorithm with the currently used onboard algorithm for descent control of the manned spacecraft Soyuz is carried out. The advantages of the new algorithm relative to the existing algorithm are presented, both in terms of stabilization accuracy and the consumption of the working fluid of the control motors.

1. Introduction

The descent and landing of spacecraft are one of the most important and crucial stages of their flight [1,2]. The control methods at these stages of flight are significantly different due to the design of the spacecraft [3,4,5]. The control methods of the shuttle-type spacecraft reentry into the atmosphere are based on well-developed algorithms for aerospaceplane. Control methods for capsule-type spacecraft have been begun to be developed in Soviet Union era, and are continue to be actively improved (transport manned spacecraft “Soyuz”), to which this article is devoted.

The works consider both uncontrolled movements (disturbed rotations in rarefied atmosphere [6], the search for stability conditions [6], predict resonance [7]), and trajectory following [8] and orientation of suborbital ships [9,10], reusable launch vehicles [11], Earth landers [12,13], Martian scientific laboratories [14,15]. The applied management methods differ significantly from each other. These are pulse modulation [11], and sliding mode control [12], and control with forward and feedback [15,16].

The problem of increasing the landing accuracy requires research to improve the on-board lander algorithms. One of the ways to improve the accuracy of landing on the Earth is the development of new control algorithms for the stabilization of the angular position of the spacecraft when moving in the atmosphere, which should have advantages over existing algorithms that have actually found application at present [15], both in accuracy and in fuel consumption for control [17,18,19,20]. This article is devoted to such research.

It should be noted that the standard control of the capsule-type spacecraft orientation in the atmosphere [12,13] is aimed at damping the angular velocities and tracking the programmed roll angle. In this case, the balancing position of the spacecraft at the other two angles (attack and glide) is maintained only due to the static stability of the spacecraft—the accuracy of such stabilization is low. It is required to improve accuracy without increasing fuel consumption.

The disturbing forces and moments during the motion of the spacecraft in the atmosphere are stochastic. They largely depend on the shape and configuration of the spacecraft and are difficult to define. Therefore, if the control law depends on constantly changing aerodynamic parameters, the fuel consumption for maintaining a given orientation of the spacecraft will increase significantly. The aim of this work is to synthesize a universal (suitable for any type of spacecraft) controller that would ensure the fulfillment of higher requirements for stability and stabilization accuracy without increasing control costs, and at the same time would be robust, i.e., would not depend on the aerodynamic parameters of the spacecraft.

We propose an approach to the synthesis of a controller based on the J.W. Van der Woude’s modal output control method, modified for dynamic systems with multiple inputs and outputs [21]. The novelty of the approach lies in the fact that due to the parameterization of the matrices with the desired spectra and the proper choice of the assigned poles, it is possible to achieve independence of the control channels for pitch, roll and yaw, both from each other and from the aerodynamic parameters of the spacecraft.

2. Statement of the Research Problem

The orientation of the capsule-type lander and the associated velocity coordinate system (CS) $E_v=o{x}_v{y}_v{z}_v$. (rotated relative to the associated geometric CS $E=oxyz$ around the axis z by the calculated “balancing” angle of attack ${\alpha }_{bal}^{*}$) relative to the reference velocity CS $Q_v=o{X}_v{Y}_v{Z}_v$ during landing from the near-earth orbit is considered [9,14,15].

For the study, a section of the trajectory is selected, on which the values of the atmosphere density $\rho$ and the linear velocity v of the lander contribute to the most effective control of the lander motion. Such a section approximately corresponds to heights $h\le 80$ and Mach numbers [22] $M\ge 6$, and the average time of movement along it is $T=280$. In this section, the aerodynamic coefficients of the SA are practically independent of the Mach number and altitude, and the balancing position, characterized by the balancing angles of attack ${\alpha }_{bal}$ and side slip ${\beta }_{bal}$ [23], can be considered constant.

The task for the study consists in high-precision stabilization of the lander in the programmed balancing position (with tracking the programmed roll angle ${\gamma }_{pr}$), i.e., in maintaining the state vector of angular motion.

$$\begin{array}{cc}{\mathit{x}}_{att}=\left[\begin{array}{c}{\theta }^{E÷Q_v}\\ {\omega }_{E_v}^{E_v}\end{array}\right]& (\begin{array}{cc}{\mathit{\theta }}^{E÷Q_v}={\left[\begin{array}{ccc}\gamma & \beta & \alpha \end{array}\right]}^T,& {\mathit{\omega }}_{E_v}^{E_v}={\left[\begin{array}{ccc}{\omega }_{x_v}^{E_v}& {\omega }_{y_v}^{E_v}& {\omega }_{z_v}^{E_v}\end{array}\right]}^T),\end{array}\end{array}$$

where $\gamma$ – speed roll angle, $\beta$ – angles of side slip, $\alpha$ – angles of attack (GOST 20058-80), near its programmed value

$$\begin{array}{cc}{\mathit{x}}_{att,pr}=\left[\begin{array}{c}{\mathit{\theta }}_{pr}^{E÷Q_v}\\ {\mathit{\omega }}_{E_{v,pr}}^{E_v}\end{array}\right]& \left({\mathit{\theta }}_{pr}^{E÷Q_v}\begin{array}{cc}={\left[\begin{array}{ccc}{\gamma }_{pr}& {\beta }_{bal}& {\alpha }_{bal}\end{array}\right]}^T,& {\mathit{\omega }}_{E_{v,pr}}^{E_v}={\left[\begin{array}{ccc}{\omega }_{x_v}^{Q_v}& {\omega }_{y_v}^{Q_v}& {\omega }_{z_v}^{Q_v}\end{array}\right]}^T\end{array}\right)\end{array}$$

(1)

with the help of reaction engines of the system of executive organs of landing of constant thrust with variable pulse durations according to information from the gyroscopic angle measurement system (GAMS) (inertial roll angle) and angular velocity sensor (vector of angular velocity ${\mathit{\omega }}_{E_v}^{E_v}$) about the observation vector

$${\mathit{y}}_{att}=\left[\begin{array}{c}\tilde{\gamma }\\ {\mathit{\omega }}_{E_v}^{E_v}\end{array}\right].$$

(2)

Hereinafter, the subscript with the name of the CS or its axes denotes the mappings of vectors (matrix columns) and matrices (tensors of inertia, kinematic equations) to the corresponding bases or projections of vectors on the indicated axes. The superscript indicates the basis, the absolute movement of which characterizes the given vector. If relative movement is considered, the symbol “÷” is added in the superscript followed by the designation of the basis relative to which the movement occurs.

It is required, using analytical methods of modal control, to increase, as far as possible, the accuracy of stabilization of the lander without increasing fuel consumption in comparison with the standard algorithm.

3. Mathematical Model

In the modeling, we will consider the lander as an absolutely rigid body, symmetric about the longitudinal axis, the center of mass (CM) does not displace along the transverse axis [24]. We introduce the generalized state vector

$${\mathit{x}}^T=\left[\begin{array}{ccc}t& {\mathit{x}}_{trj}^T& {\mathit{x}}_{att}^T\end{array}\right],$$

where t – time from the moment of powering on the GAMS; ${\mathit{x}}_{trj}$ – the state vector during the movement of the CM (three coordinates of position and velocity).

The model of motion of the CM of the lander [25] has the general form

$${\dot{\mathit{x}}}_{trj}={\mathit{f}}_{trj}\left({\mathit{x}}_{trj},{\mathit{F}}_{Q_v}\left({\mathit{x}}_{trj},{\mathit{x}}_{att}\right)\right),$$

(3)

where ${\mathbf{F}}_{Q_v}$ – resultant of aerodynamic, gravitational and inertial forces. The model of the angular motion of the lander [9,10] (between the CS $E_v$ and $Q_v$) consists of the kinematic and dynamic equations [26,27,28], as well as the measurement model:

$$\begin{array}{c}{\dot{\mathit{x}}}_{att}={\mathit{f}}_{att}\left({\mathit{x}}_{trj},{\mathit{x}}_{att},{\mathit{M}}_{E_v}^{ctr}\right)=\left[\begin{array}{c}{\mathit{G}}_{E_v}^{E_v÷Q_v}\left({\mathit{\theta }}^{E÷Q_v}\right)\left({\mathit{\omega }}_{E_v}^{E_v}-{\mathit{\omega }}_{E_v}^{Q_v}\left({\mathit{x}}_{trj},{\mathit{\theta }}^{E÷Q_v}\right)\right)\\ {\mathit{J}}_{E_v}^{-1}\left({\mathit{M}}_{E_v}^{aer}\left({\mathit{x}}_{trj},{\mathit{\theta }}^{E÷Q_v}\right)+{\mathit{M}}_{E_v}^{gyr}\left({\mathit{\omega }}_{E_v}^{E_v}\right)+{\mathit{M}}_{E_v}^{dst}+{\mathit{M}}_{E_v}^{ctr}\right)\end{array}\right],\\ {\mathit{y}}_{att}={\mathit{g}}_{att}\left(t,{\mathit{x}}_{trj},{\mathit{x}}_{att}\right)=\left[\begin{array}{c}\tilde{\gamma }\left(t,{\mathit{x}}_{trj},{\mathit{\theta }}^{E÷Q_v}\right)\\ {\mathit{\omega }}_{E_v}^{E_v}\end{array}\right],\end{array}$$

(4)

here ${\mathit{G}}_{E_v}^{E_v÷Q_v}$ – the matrix of the kinematic equations of motion of the CS $E_v$ relative to the CS $Q_v$ in the Krylov angles; ${\mathit{J}}_{E_v}$ – tensor of inertia of lander relative to CM in basis $E_v$; ${\mathit{M}}_{E_v}^{aer}$, ${\mathit{M}}_{E_v}^{gyr}$, ${\mathit{M}}_{E_v}^{dst}$ and ${\mathit{M}}_{E_v}^{ctr}$ – respectively, aerodynamic, gyroscopic, disturbing unbalanced and control moments relative to the CM of lander.

4. Standard Algorithm and Conditions for Its Comparison with the New Algorithm

The standard motion control of a capsule-type lander from a constant thrust provided by reaction engine is formed in two stages. First, the vector of control signals is found

$$\mathit{u}={\mathit{u}}_0-\mathit{F}\Delta {\mathit{y}}_{att},$$

(5)

where ${\mathit{u}}_0$ – feedforward control to control, $\mathit{F}$ – output regulator matrix, $\Delta {\mathit{y}}_{att}={\mathit{y}}_{att}-{\mathit{y}}_{att,pr}$, and ${\mathit{y}}_{att,pr}$ – the programmed value of the output vector (2), consisting of the programmed values of the inertial roll angle ${\tilde{\gamma }}_{pr}$ and the vector of the angular velocity of the associated CS ${\mathit{\omega }}_{E_v,pr}^{E_v}$. Then the signals (5) are converted into the control torque of constant thrust of reaction engine (the duration of switching on the reaction engine) according to the piecewise linear law with a dead zone and saturation.

The standard stabilization algorithm is empirical. In it, control (5) turns out to be autonomous in the roll $\left(\begin{array}{cc}\gamma ,& {\omega }_{x_v}^{E_v}\end{array}\right)$, yaw $\begin{array}{cc}(\beta ,& {\omega }_{y_v}^{E_v}\end{array})$ and pitch $\left(\begin{array}{cc}\alpha ,& {\omega }_{z_v}^{E_v}\end{array}\right)$ channels. In the atmospheric section, it is aimed at damping the angular velocities ${\mathit{\omega }}_{E_v}^{E_v}$ and tracking the programmed roll angle $(\gamma \to {\gamma }_{pr})$. The balancing position $\left(\begin{array}{cc}\alpha \to {\alpha }_{bal},& \beta \to {\beta }_{bal}\end{array}\right)$ is maintained due to the static stability of the lander [23].

In the simplest version of law (5)

$$\begin{array}{cc}{\mathit{u}}_0={\mathit{\omega }}_{E_v,pr}^{E_v}={\mathbf{0}}_{3\times 1},& {\tilde{\gamma }}_{pr}={\gamma }_{pr}.\end{array}$$

Hereinafter, ${\mathbf{0}}_{n\times m}$ is the zero matrix of dimension $n\times m$. To improve the orientation accuracy and create the same conditions when comparing the standard and new algorithms (according to the influence of the regulator matrices), it is advisable to use the conversion of the speed roll angle into the inertial roll angle

$${\tilde{\gamma }}_{pr}=\tilde{\gamma }\left(t,{\mathit{x}}_{trj},{\mathit{\theta }}_{pr}^{E÷Q_v}\right)$$

and taking into account the angular velocity of the reference CS $Q_v$

$${\mathit{\omega }}_{E_v,pr}^{E_v}={\mathit{\omega }}_{E_v}^{Q_v}\left({\mathit{x}}_{trj},{\mathit{\theta }}_{pr}^{E÷Q_v}\right).$$

Since in model (4) in the balancing position

$${\mathit{M}}_{E_v}^{aer}\left({\mathit{x}}_{trj},{\mathit{\theta }}_{pr}^{E÷Q_v}\right)\equiv {\mathbf{0}}_{3\times 1},$$

the feedforward control without taking into account the current kinematics of the lander motion is

$${\underset{-}{\mathit{u}}}_0=-{\mathit{M}}_{E_v}^{gyr}\left({\mathit{\omega }}_{E_v,pr}^{E_v}\right).$$

(6)

In real time, the angular velocity ${\mathit{\omega }}_{E_v}^{Q_v}$, which depends, in particular, on the aerodynamic force (a strongly varying stochastic vector), is difficult to calculate. But as the programmed angular velocity with some approximation, you can use its part

$${\mathit{\omega }}_{E_v,pr}^{E_v}={\mathit{\omega }}_{E_v}^{Q_{grn}}\left({\mathit{x}}_{trj},{\mathit{\theta }}_{pr}^{E÷Q_v}\right),$$

(7)

i.e., take into account only the angular velocity of the Earth’s rotation around its axis.

5. Linearization of Angular Motion Model of the Lander

We linearize the system of Equations (4) at each computational step of the on-board computer by expanding the right-hand sides into a Taylor series in terms of the coordinates of the vector ${\mathit{x}}_{att}$ near their programmed values

$${\mathit{x}}_{att,pr}=\left[\begin{array}{c}{\mathit{\theta }}_{pr}^{E÷Q_v}\\ {\mathit{\omega }}_{E_v}^{Q_{grn}}\left({\mathit{x}}_{trj},{\mathit{\theta }}_{pr}^{E÷Q_v}\right)\end{array}\right],$$

which are constant per cycle and written taking into account equalities (1) and (7) at the current values of time t and coordinates of the CM ${\mathit{x}}_{trj}$ from model (3). This linearization is possible because in the considered range of heights (above 40 km), the parameters of the CM motion of the lander change more slowly than the parameters of the angular motion.

After linearization, a system of approximate equations is formed

$$\left\{\begin{array}{cc}\hfill & {\dot{\mathit{x}}}_{att}=\underset{\xi }{\underbrace{{\mathit{f}}_{att}\left({\mathit{x}}_{trj},{\mathit{x}}_{att},{\mathbf{0}}_{3\times 1}\right)}}+\underset{\mathit{A}}{\underbrace{{\mathit{f}}_{att}^{\prime }{\mathit{x}}_{att}\left({\mathit{x}}_{trj},{\mathit{x}}_{att},{\mathbf{0}}_{3\times 1}\right)}}\underset{\Delta {\mathit{x}}_{att}}{\underbrace{\left({\mathit{x}}_{att}-{\mathit{x}}_{att,pr}\right)}}+\underset{\mathit{B}}{\underbrace{\left[\begin{array}{c}{\mathbf{0}}_{3\times 3}\\ {\mathit{J}}_{E_v}^{-1}\end{array}\right]}}\mathit{u},\hfill \\ \hfill & {\mathit{y}}_{att}=\underset{{\mathit{y}}_{att,pr}={\mathit{y}}_{att}-\Delta {\mathit{y}}_{att}}{\underbrace{{\mathit{g}}_{att}\left(t,{\mathit{x}}_{trj},{\mathit{x}}_{att,pr}\right)}}+\underset{\mathit{C}}{\underbrace{{\mathit{g}}_{att}^{\prime }{\mathit{x}}_{att}\left(t,{\mathit{x}}_{trj},{\mathit{x}}_{att,pr}\right)}}\underset{\Delta {\mathit{x}}_{att}}{\underbrace{\left({\mathit{x}}_{att}-{\mathit{x}}_{att,pr}\right)}}.\hfill \end{array}\right.$$

(8)

In deviations from the programmed values, this system takes the classical form with a perturbation

$$\left\{\begin{array}{cc}\hfill & \Delta {\dot{\mathit{x}}}_{att}=\mathit{\xi }+\mathit{A}\Delta {\mathit{x}}_{att}+\mathit{Bu},\hfill \\ \hfill & \Delta {\mathit{y}}_{att}=\mathit{C}\Delta {\mathit{x}}_{att},\hfill \end{array}\right.$$

(9)

where $\mathit{\xi }$ – vector of disturbances, $\mathit{A}$, $\mathit{B}$ and $\mathit{C}$ – the matrices of state, control and observation, respectively.

According to the theory presented in [29], the control for such a model that ensures the fulfillment of the condition

$$eig\left(\mathit{A}-\mathit{BFC}\right)={\Lambda }^{*}$$

(10)

for the desired spectrum

$${\Lambda }^{*}=\left\{\begin{array}{cccccc}{\varphi }_1,& {\varphi }_2,& {\varphi }_3,& {\varphi }_4,& {\varphi }_5,& {\varphi }_6\end{array}\right\}$$

has the form (5). The static addition caused by the presence of a perturbation (zero Taylor term) with the matrix of controlled parameters

$$\mathit{N}=\left[\begin{array}{cc}{\mathit{I}}_3& {\mathbf{0}}_{3\times 3}\end{array}\right]$$

(the angular position ${\mathit{\theta }}^{E÷Q_v}$ is regulated) is equal to

$${\mathit{u}}_0=-\left[\begin{array}{cc}\mathit{FC}& {\mathit{I}}_3\end{array}\right]{\left[\begin{array}{cc}\mathit{A}& \mathit{B}\\ \mathit{N}& {\mathbf{0}}_{3\times 3}\end{array}\right]}^{-1}\left[\begin{array}{c}{\mathit{I}}_6\\ {\mathbf{0}}_{3\times 6}\end{array}\right]\mathit{\xi },$$

(11)

here and below, ${\mathit{I}}_n$ is the identity matrix of order n.

Thus, the control problem is reduced to finding the output controller matrix $\mathit{F}$ for a triple of matrices $\mathit{A}\left(t\right)$, $\mathit{B}$, and $\mathit{C}\left(t\right)$ (matrices $\mathit{A}$ and $\mathit{C}$, due to linearization (8), change from cycle to cycle), written in block form

$$\begin{array}{ccc}\mathit{A}\left(t\right)=\left[\begin{array}{cc}{\mathit{A}}_{[1,1]}\left(t\right)& {\mathit{A}}_{[1,2]}\left(t\right)\\ {\mathit{A}}_{[2,1]}\left(t\right)& {\mathit{A}}_{[2,2]}\left(t\right)\end{array}\right],& \mathit{B}=\left[\begin{array}{c}{\mathbf{0}}_{3\times 3}\\ {\mathit{J}}^{-1}\end{array}\right],& \mathit{C}\left(t\right)=\left[\begin{array}{cc}{\mathit{c}}^T\left(t\right)& {\mathbf{0}}_{1\times 3}\\ {\mathbf{0}}_{3\times 3}& {\mathit{I}}_3\end{array}\right],\end{array}$$

(12)

where ${\mathit{A}}_{[1,1]}\left(t\right),{\mathit{A}}_{[1,2]}\left(t\right),{\mathit{A}}_{[2,1]}\left(t\right),{\mathit{A}}_{[2,2]}\left(t\right)\in {\mathit{R}}^{3\times 3}$ – blocks of the state matrix (in the general case, not zero), $\mathit{c}\left(t\right)={\left[\begin{array}{ccc}{c}_{\gamma }\left(t\right)\hfill & c_{\beta }\left(t\right)& c_{\alpha }\left(t\right)\end{array}\right]}^T$ – vector of measured combinations of kinematic parameters, $\mathit{J}={\mathit{J}}_{E_v}$.

Based on the simulation results in MATLAB for the full linear model (9), the values of the variable coefficients of the matrices $\mathit{A}\left(t\right)$ and $\mathit{C}\left(t\right)$ from the record (12) were estimated. It turned out that many coefficients change insignificantly during the motion of the lander, and are close to zero or one in magnitude. The most significant changes are the coefficients $a_{4,2}\left(t\right)$, $a_{5,2}\left(t\right)$, and $a_{6,3}\left(t\right)$ of the state matrix $\mathit{A}\left(t\right)=\left[a_{i,j}\left(t\right)\right]$.

Let us form a simplified linear model (Figure 1) by changing the record of the state and observation matrices in comparison with the record (12) and introducing an underscore for the simplified matrices:

$$\begin{array}{ccc}\underset{-}{\mathit{A}}\left(t\right)=\left[\begin{array}{cc}{\mathbf{0}}_{3\times 3}& {\mathit{I}}_3\\ {\underset{-}{\mathit{A}}}_{[2,1]}\left(t\right)& {\mathbf{0}}_{3\times 3}\end{array}\right],& \mathit{B}=\left[\begin{array}{c}{\mathbf{0}}_{3\times 3}\\ {\mathit{J}}^{-1}\end{array}\right],& \underset{-}{\mathit{C}}=\left[\begin{array}{cc}{\mathit{e}}_1^T& {\mathbf{0}}_{1\times 3}\\ {\mathbf{0}}_{3\times 3}& {\mathit{I}}_3\end{array}\right],\end{array}$$

(13)

where

$$\begin{array}{c}{\underset{-}{\mathit{A}}}_{[2,1]}\left(t\right)=\left[\begin{array}{ccc}0& -a_{4,2}\left(t\right)& 0\\ 0& -a_{5,2}\left(t\right)& 0\\ 0& 0& -a_{6,3}\left(t\right)\end{array}\right]=-\tilde{q}\left(t\right)\left[\begin{array}{ccc}0& {\overline{a}}_{4,2}& 0\\ 0& {\overline{a}}_{5,2}& 0\\ 0& 0& {\overline{a}}_{6,3}\end{array}\right];\end{array}$$

${\mathit{e}}_1={\left[\begin{array}{ccc}1& 0& 0\end{array}\right]}^T$; $\tilde{q}\left(t\right)>0$ – scaled value of the velocity head; ${\overline{a}}_{4,2}>0$, ${\overline{a}}_{5,2}>0$, ${\overline{a}}_{6,3}>0$ – constant of linearization of the scaled aerodynamic moment.

The static addition ${\mathit{u}}_0$ in the control law (5) for the simplified model instead of formula (11), is calculated in a simplified way: ${\mathit{u}}_0={\underset{-}{\mathit{u}}}_0$, using formulas (6) and (7).

Figure 1. Onboard model of the spacecraft angular motion.

Figure 1. Onboard model of the spacecraft angular motion.

Since the coefficients ${\overline{a}}_{4,2}$, ${\overline{a}}_{5,2}$, ${\overline{a}}_{6,3}$, and $\tilde{q}\left(t\right)$ are stochastic and difficult to determine, the problem arises for object (9) with matrices $\mathit{A}=\underset{-}{\mathit{A}}$ and $\mathit{C}=\underset{-}{\mathit{C}}$ to synthesize a robust output controller (5) with a stationary matrix $\mathit{F}$ independent of the matrix components $\mathit{A}$, which nevertheless ensures that condition (10) is satisfied.

6. Robust Output Regulator

The simplified linear model described by the triple matrices (13) can be split into two components:

• autonomous model in the Roll-Yaw channel

$$\begin{array}{ccc}{\mathit{A}}_{RY}=\left[\begin{array}{cc}{\mathbf{0}}_{2\times 2}& {\mathit{I}}_2\\ {\mathit{A}}_{RY[2,1]}& {\mathbf{0}}_{2\times 2}\end{array}\right],& {\mathit{B}}_{RY}=\left[\begin{array}{c}{\mathbf{0}}_{2\times 2}\\ {\mathit{J}}_{RY}^{-1}\end{array}\right],& {\mathit{C}}_{RY}=\left[\begin{array}{cc}{\mathit{e}}_R^T& {\mathbf{0}}_{1\times 2}\\ {\mathbf{0}}_{2\times 2}& {\mathit{I}}_2\end{array}\right]\end{array}$$

(14)

where

$$\begin{array}{ccc}\begin{array}{cc}{\mathit{A}}_{RY[2,1]}=\left[\begin{array}{cc}0& -a_{4,2}\\ 0& -a_{5,2}\end{array}\right],& \begin{array}{c}{a}_{4,2}>0,\hfill \\ a_{5,2}>0,\hfill \end{array}\end{array}& {\mathit{J}}_{RY}=\left[\begin{array}{cc}{J}_{x_v}& -J_{x_v{y}_v}\\ -J_{x_v{y}_v}& J_{y_v}\end{array}\right],& {\mathit{e}}_R=\left[\begin{array}{c}1\\ 0\end{array}\right]\end{array},$$

with the desired spectrum eigenvectors

$$\begin{array}{c}{\Lambda }_{RY}=\left\{\begin{array}{cccc}{\varphi }_1,& {\varphi }_2,& {\varphi }_3,& {\varphi }_4\end{array}\right\};\end{array}$$

(15)

• autonomous model in the Pitch channel

  $$\begin{array}{ccc}{\mathit{A}}_P=\left[\begin{array}{cc}0& 1\\ -a_{6,3}& 0\end{array}\right],& {\mathit{B}}_P=\left[\begin{array}{c}0\\ J_{z_v}^{-1}\end{array}\right],& {\mathit{C}}_P=\left[\begin{array}{cc}0& 1\end{array}\right],\end{array}$$

  (16)

  where $a_{63}>0$, with the desired spectrum eigenvectors

  $${\Lambda }_P=\left\{\begin{array}{cc}{\varphi }_5,& {\varphi }_6\end{array}\right\}.$$

  (17)

In records (14) and (16), the axial moments of inertia $J_{x_v}$, $J_{y_v}$, $J_{z_v}$ and the centrifugal moment of inertia $J_{x_v{y}_v}$ in the SC $E_v$ were used.

Let us consider an autonomous problem of modal output control in the Roll-Yaw channel, described by a completely controllable and completely observable (by state) triple of matrices (14) and spectrum (15). We obtain a parameterized set of its solutions based on a modification of the direct van der Wood approach [21,30].

At the zero level of decomposition

$$\begin{array}{cc}{\mathit{A}}_{RY0}={\mathit{A}}_{RY},& {\mathit{B}}_{RY0}={\mathit{B}}_{RY}\end{array}$$

the left annihilator and its pseudoinverse matrix [31] are respectively equal to

$$\begin{array}{cc}{\mathit{B}}_{RY0}^{\perp L}=\left[\begin{array}{cc}{\mathit{I}}_2& {\mathbf{0}}_{2\times 2}\end{array}\right],& {\mathit{B}}_{RY0}^{\perp L+}={\mathit{B}}_{RY0}^{\perp LT}{\left({\mathit{B}}_{RY0}^{\perp L}{\mathit{B}}_{RY0}^{\perp LT}\right)}^{-1}={\mathit{B}}_{RY0}^{\perp LT}.\end{array}$$

The first level of decomposition

$$\begin{array}{cc}{\mathit{A}}_{RY1}={\mathit{B}}_{RY0}^{\perp L}{\mathit{A}}_{RY0}{\mathit{B}}_{RY0}^{\perp L+}={\mathbf{0}}_{2\times 2},& {\mathit{B}}_{RY1}={\mathit{B}}_{RY0}^{\perp L}{\mathit{A}}_{RY0}{\mathit{B}}_{RY0}={\mathit{J}}_{RY}^{-1}\end{array}$$

is finite due to the invertibility of matrix ${\mathit{B}}_{RY1}$.

The controller matrix at the first decomposition level is

$${\mathit{K}}_{RY1}={\mathit{B}}_{RY1}^{-1}{\mathit{A}}_{RY1}-{\mathbf{\Phi }}_{RY1}{\mathit{B}}_{RY1}^{-1}=-{\mathbf{\Phi }}_{RY1}{\mathit{J}}_{RY}=-{\mathit{J}}_{RY}{\tilde{\mathbf{\Phi }}}_{RY1},$$

where ${\tilde{\mathbf{\Phi }}}_{RY1}$ and ${\mathbf{\Phi }}_{RY1}={\mathit{J}}_{RY}{\tilde{\mathbf{\Phi }}}_{RY1}{\mathit{J}}_{RY}^{-1}$ are mutually similar matrices with the spectrum

$$eig{\tilde{\mathbf{\Phi }}}_{RY1}=eig{\mathbf{\Phi }}_{RY1}=\left\{\begin{array}{cc}{\varphi }_1,& {\varphi }_2\end{array}\right\}.$$

(18)

The pseudoinverse matrix and the auxiliary matrix at the zero decomposition level are, respectively, equal

$$\begin{array}{c}{\mathit{B}}_{RY0}^{+}={\left({\mathit{B}}_{RY0}^T{\mathit{B}}_{RY0}\right)}^{-1}{\mathit{B}}_{RY0}^T=\left[\begin{array}{cc}{\mathbf{0}}_{2\times 2}& {\mathit{J}}_{RY}\end{array}\right],\\ {\mathit{B}}_{RY0}^{-}={\mathit{B}}_{RY0}^{+}+{\mathit{K}}_{RY1}{\mathit{B}}_{RY0}^{\perp L}={\mathit{J}}_{RY}\underset{{\tilde{\mathit{B}}}_{RY0}^{-}}{\underbrace{\left[\begin{array}{cc}-{\tilde{\mathbf{\Phi }}}_1& {\mathit{I}}_2\end{array}\right]}},\end{array}$$

and the controller matrix at this level is

$$\begin{array}{c}{\mathit{K}}_{RY}={\mathit{K}}_{RY0}={\mathit{B}}_{RY0}^{-}{\mathit{A}}_{RY0}-{\mathbf{\Phi }}_{RY0}{\mathit{B}}_{RY0}^{-}=\\ ={\mathit{J}}_{RY}\left[\begin{array}{cc}{\mathit{A}}_{RY[2,1]}+{\tilde{\mathbf{\Phi }}}_{RY0}{\tilde{\mathbf{\Phi }}}_{RY1}& -{\tilde{\mathbf{\Phi }}}_{RY0}-{\tilde{\mathbf{\Phi }}}_{RY1}\end{array}\right],\end{array}$$

(19)

where ${\tilde{\mathbf{\Phi }}}_{RY0}$ and ${\mathbf{\Phi }}_{RY0}=\mathit{J}{\tilde{\mathbf{\Phi }}}_{RY0}{\mathit{J}}^{-1}$ are mutually similar matrices with the spectrum

$$eig{\tilde{\mathbf{\Phi }}}_{RY0}=eig{\mathbf{\Phi }}_{RY0}=\left\{\begin{array}{cc}{\varphi }_3,& {\varphi }_4\end{array}\right\}.$$

(20)

Matrix (19) characterizes the modal state controller for the pair of matrices $\left({\mathit{A}}_{RY},{\mathit{B}}_{RY}\right)$ from the recorder (14) and the spectrum (15).

To calculate the modal output controller, we write the right annihilator of the matrix ${\mathit{C}}_{RY0}={\mathit{C}}_{RY}$

$${\mathit{C}}_{RY0}^{\perp R}=-{\left[\begin{array}{cccc}0& 1& 0& 0\end{array}\right]}^T$$

and in the equation:

$${\tilde{\mathbf{\Phi }}}_{RY0}\underset{{\tilde{\mathit{G}}}_{RY0}}{\underbrace{{\tilde{\mathit{B}}}_{RY0}^{-}{\mathit{C}}_{RY0}^{\perp R}}}=\underset{{\tilde{\mathit{H}}}_{RY0}}{\underbrace{{\tilde{\mathit{B}}}_{RY0}^{-}{\mathit{A}}_{RY0}{\mathit{C}}_{RY0}^{\perp R}}},$$

(21)

calculate the matrix coefficients

$$\begin{array}{cc}{\tilde{\mathit{G}}}_{RY0}={\tilde{\mathbf{\Phi }}}_{RY1}\left[\begin{array}{c}0\\ 1\end{array}\right],& {\tilde{\mathit{H}}}_{RY0}=\left[\begin{array}{c}{a}_{4,2}\\ a_{5,2}\end{array}\right].\end{array}$$

(22)

From equations (22) it can be seen that the right side of Equation (21) is unchanged, and the matrix coefficient for the calculated matrix ${\tilde{\mathbf{\Phi }}}_{RY0}$ on the left side can be changed depending on the value of the matrix ${\tilde{\mathbf{\Phi }}}_{RY1}$. Let us write the value of this coefficient in a general parameterized form $({\mu }_1,{\mu }_2\in \mathbb{C})$:

$${\tilde{\mathit{G}}}_{RY0}^{*}=\left[\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right].$$

(23)

Let us find the matrix ${\mathbf{\Phi }}_1$ with the desired spectrum (18) at the first level of decomposition, at which the equality $G_1=G_2$ is fulfilled. To do this, let us consider the expression for the matrix G1 from the record (22) as the equation

$${\tilde{\mathbf{\Phi }}}_{RY1}\underset{{\tilde{\mathit{G}}}_{RY1}}{\underbrace{\left[\begin{array}{c}0\\ 1\end{array}\right]}}=\underset{{\tilde{\mathit{H}}}_{RY1}}{\underbrace{\left[\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right]}},$$

(24)

solvable with respect to the matrix ${\tilde{\mathbf{\Phi }}}_{RY1}$. To provide the spectrum (5.5) to the matrix ${\tilde{\mathbf{\Phi }}}_{RY1}$, we form a pair of state and observation matrices ${\tilde{\mathbf{\Phi }}}_{RY1}$.

$$\begin{array}{cc}{\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}={\tilde{\mathit{H}}}_{RY1}\underset{{\tilde{\mathit{G}}}_{RY1}^{+}}{\underbrace{{\left({\tilde{\mathit{G}}}_{RY1}^T{\tilde{\mathit{G}}}_{RY1}\right)}^{-1}{\tilde{\mathit{G}}}_{RY1}^T}}=\left[\begin{array}{cc}0& {\mu }_1\\ 0& {\mu }_2\end{array}\right],& {\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}={\tilde{\mathit{G}}}_{RY1}^{\perp L}=\left[\begin{array}{cc}1& 0\end{array}\right].\end{array}$$

(25)

For this pair, the observability matrix and its determinant are

$$\begin{array}{cc}{\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}=\left[\begin{array}{c}{\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}\\ {\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}{\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}\end{array}\right]=\left[\begin{array}{cc}1& 0\\ 0& {\mu }_1\end{array}\right],& det{\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}={\mu }_1,\end{array}$$

i.e., full observability takes place if

$${\mu }_1\ne 0.$$

(26)

Taking into account the Ackerman formula [30], we calculate the state observer matrix for a pair of matrices (25) and spectrum (18):

$${\mathit{L}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}=\left({\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}-{\varphi }_1{\mathit{I}}_2\right)\left({\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}-{\varphi }_2{\mathit{I}}_2\right){\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}^{-1}\left[\begin{array}{c}0\\ 1\end{array}\right]=\left[\begin{array}{c}{\tilde{\mu }}_2\\ {\tilde{\mu }}_1\end{array}\right],$$

where

$$\begin{array}{cc}{\tilde{\mu }}_1=\frac{\left({\mu }_2-{\varphi }_1\right)\left({\mu }_2-{\varphi }_2\right)}{{\mu }_1},& {\tilde{\mu }}_2={\mu }_2-{\varphi }_1-{\varphi }_2.\end{array}$$

The matrix with the desired spectrum at the first level of decomposition generally has a form that is neither diagonal nor triangular:

$${\tilde{\mathbf{\Phi }}}_{RY1}={\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}-{\mathit{L}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}{\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY1}}=\left[\begin{array}{cc}-{\tilde{\mu }}_2& {\mu }_1\\ -{\tilde{\mu }}_1& {\mu }_2\end{array}\right].$$

(27)

Next, from Equation (21), we find the matrix ${\tilde{\mathbf{\Phi }}}_{RY0}$ with the desired spectrum (20) at the zero level of decomposition. Let the matrix ${\tilde{\mathbf{\Phi }}}_{RY1}$ of the form (27) be assigned at the first level of decomposition. Then the matrix coefficients (22) are

$$\begin{array}{cc}{\tilde{\mathit{G}}}_{RY0}={\tilde{\mathit{G}}}_{RY0}^{*}=\left[\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right],& {\tilde{\mathit{H}}}_{RY0}=\left[\begin{array}{c}{a}_{42}\\ a_{52}\end{array}\right].\end{array}$$

Hence, when satisfying inequality (26), Equation (21) is solvable with respect to the matrix ${\tilde{\mathbf{\Phi }}}_{RY0}$. To provide the spectrum (20) to the matrix ${\tilde{\mathbf{\Phi }}}_{RY0}$, we form a pair of state and observation matrices

$$\begin{array}{c}{\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}={\tilde{\mathit{H}}}_{RY0}\underset{{\tilde{\mathit{G}}}_{RY0}^{+}}{\underbrace{{\left({\tilde{\mathit{G}}}_{RY0}^T{\tilde{\mathit{G}}}_{RY0}\right)}^{-1}{\tilde{\mathit{G}}}_{RY0}^T}}=\frac{1}{\mu }\left[\begin{array}{cc}{a}_{4,2}{\mu }_1& a_{4,2}{\mu }_2\\ a_{5,2}{\mu }_1& a_{5,2}{\mu }_2\end{array}\right],\\ {\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}={\tilde{\mathit{G}}}_{RY0}^{\perp L}=\left[\begin{array}{cc}-{\mu }_2& {\mu }_1\end{array}\right],\end{array}$$

(28)

where $\mu ={\mu }_1^2+{\mu }_2^2$. For this pair, the observability matrix and its determinant are

$$\begin{array}{cc}{\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}=\left[\begin{array}{c}{\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}\\ {\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}{\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}\end{array}\right]=-\left[\begin{array}{cc}{\mu }_2& -{\mu }_1\\ \frac{{\mu }_1}{\mu }\delta & \frac{{\mu }_2}{\mu }\delta \end{array}\right],& det\left({\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}\right)=\delta ,\end{array}$$

where $\delta =a_{42}{\mu }_2-a_{52}{\mu }_1$, i.e., full observability takes place if

$${\mu }_2\ne \frac{a_{5,2}}{a_{4,2}}{\mu }_1.$$

(29)

By Ackerman formula [30], we calculate the state observer matrix for a pair of matrices (28) and spectrum (20):

$${\mathit{L}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}=\left({\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}-{\varphi }_3{\mathit{I}}_2\right)\left({\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}-{\varphi }_4{\mathit{I}}_2\right){\mathit{N}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}^{-1}\left[\begin{array}{c}0\\ 1\end{array}\right]=-\frac{1}{\mu }\left[\begin{array}{c}{\kappa }_{11}{\mu }_1+{\kappa }_{12}{\mu }_2\\ {\kappa }_{21}{\mu }_1+{\kappa }_{22}{\mu }_2\end{array}\right],$$

where

$$\begin{array}{cc}& {\kappa }_{11}={\delta }^{-1}\left(a_{4,2}-{\varphi }_3{\mu }_1\right)\left(a_{4,2}-{\varphi }_4{\mu }_1\right),\hfill \\ & {\kappa }_{12}={\delta }^{-1}\left(a_{4,2}\left(a_{5,2}-{\varphi }_3{\mu }_2\right)-{\varphi }_4{\mu }_2\left(a_{4,2}-{\varphi }_3{\mu }_1\right)\right),\hfill \\ & {\kappa }_{21}={\delta }^{-1}\left(a_{5,2}\left(a_{4,2}-{\varphi }_3{\mu }_1\right)-{\varphi }_4{\mu }_1\left(a_{5,2}-{\varphi }_3{\mu }_2\right)\right),\hfill \\ & {\kappa }_{22}={\delta }^{-1}\left(a_{5,2}-{\varphi }_3{\mu }_2\right)\left(a_{5,2}-{\varphi }_4{\mu }_2\right).\hfill \end{array}$$

The matrix with the desired spectrum at the zero decomposition level will take the form

$${\tilde{\mathbf{\Phi }}}_{RY0}={\mathit{A}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}-{\mathit{L}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}{\mathit{C}}_{{\tilde{\mathbf{\Phi }}}_{RY0}}=\left[\begin{array}{cc}-{\kappa }_{12}& {\kappa }_{11}\\ -{\kappa }_{22}& {\kappa }_{21}\end{array}\right].$$

(30)

Next, we substitute matrices ${\tilde{\mathbf{\Phi }}}_{RY0}$ (30) and ${\tilde{\mathbf{\Phi }}}_{RY1}$ (27) into the calculation formula of the state regulator matrix (19) and calculate the output regulator matrix

$${\mathit{F}}_{RY}={\mathit{K}}_{RY}\underset{{\mathit{C}}_{RY}^{+}}{\underbrace{{\mathit{C}}_{RY}^T{\left({\mathit{C}}_{RY}{\mathit{C}}_{RY}^T\right)}^{-1}}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}\underset{{\tilde{f}}_{RY1,1}}{\underbrace{{\kappa }_{12}{\tilde{\mu }}_2-{\kappa }_{11}{\tilde{\mu }}_1}}& \underset{{\tilde{f}}_{RY1,2}}{\underbrace{{\tilde{\mu }}_2+{\kappa }_{12}}}& \underset{{\tilde{f}}_{RY1,3}}{\underbrace{-{\mu }_1-{\kappa }_{11}}}\\ \underset{{\tilde{f}}_{RY2,1}}{\underbrace{{\kappa }_{22}{\tilde{\mu }}_2-{\kappa }_{21}{\tilde{\mu }}_1}}& \underset{{\tilde{f}}_{RY2,2}}{\underbrace{{\tilde{\mu }}_1+{\kappa }_{22}}}& \underset{{\tilde{f}}_{RY3,3}}{\underbrace{-{\mu }_2-{\kappa }_{21}}}\end{array}\right].$$

(31)

This matrix describes the set of solutions to the modal output control problem (14), (15), characterized by the parameters ${\mu }_1,{\mu }_2$ and the poles $f_1,f_2,f_3,f_4$. Symbolic calculations in MATLAB confirm the validity of the equation

$$eig\left({\mathit{A}}_{RY}-{\mathit{B}}_{RY}{\mathit{F}}_{RY}{\mathit{C}}_{RY}\right)={\Lambda }_{RY}$$

written on the basis of expressions (14), (15) and (31).

In order to reduce the mutual influence of control channels, let us set the problem of zeroing the cross coefficients ${\tilde{f}}_{RY1,3}$, ${\tilde{f}}_{RY2,1}$ and ${\tilde{f}}_{RY2,2}$ between the “Roll” and “Yaw” channels from the record (31), having at our disposal the parameters ${\mu }_1$, ${\mu }_2$ obeying conditions (26) and (29), as well as any ratios of the poles ${\varphi }_1$, ${\varphi }_2$, ${\varphi }_3$, ${\varphi }_4$, that do not violate the location of these poles in the left complex half-plane. This problem is described by the following system of equations and inequalities:

$$\begin{array}{cc}\left.\left\{\begin{array}{ccc}{\tilde{f}}_{RY1,3}=0\hfill & \iff & {\kappa }_{11}+{\mu }_1=0,\\ {\left.{\tilde{f}}_{RY2,1}=0\right|}_{{\tilde{f}}_{RY2,2}=0}\hfill & \iff & \left[\begin{array}{cc}\hfill & {\tilde{\mu }}_1=0,\hfill \\ \hfill & {\kappa }_{21}+{\tilde{\mu }}_2=0,\hfill \end{array}\right.\\ {\tilde{f}}_{RY2,2}=0\hfill & \iff & {\kappa }_{22}+{\tilde{\mu }}_1=0,\end{array}\right.\right|& \begin{array}{cc}\begin{array}{cc}\hfill & {\mu }_1\ne 0,\hfill \\ \hfill & {\mu }_2\ne \frac{a_{5,2}}{a_{4,2}}{\mu }_1,\hfill \end{array}& \begin{array}{c}Re{\varphi }_1<0,\hfill \\ Re{\varphi }_2<0,\hfill \\ Re{\varphi }_3<0,\hfill \\ Re{\varphi }_4<0.\hfill \end{array}\end{array}\end{array}$$

(32)

Since $a_{52}>0$, the solution to problem (32) is one of the systems

$$\begin{array}{cc}\left\{\begin{array}{cc}\hfill & \left[\begin{array}{cc}\hfill & {\mu }_2={\varphi }_1,\hfill \\ \hfill & {\mu }_2={\varphi }_2,\hfill \end{array}\right.\hfill \\ \hfill & \left[\begin{array}{cc}\hfill & \left\{\begin{array}{cc}\hfill & {\mu }_2\notin \left\{\begin{array}{cc}{\varphi }_4,& {\varphi }_4-{\varphi }_3\end{array}\right\},\hfill \\ \hfill & {\mu }_1=\frac{a_{4,2}}{{\varphi }_4-{\mu }_2},\hfill \\ \hfill & {\mu }_2{\varphi }_3=a_{5,2},\hfill \end{array}\right.\hfill \\ \hfill & \left\{\begin{array}{cc}\hfill & {\mu }_2\notin \left\{\begin{array}{cc}{\varphi }_3,& {\varphi }_3-{\varphi }_4\end{array}\right\},\hfill \\ \hfill & {\mu }_1=\frac{a_{4,2}}{{\varphi }_3-{\mu }_2},\hfill \\ \hfill & {\mu }_2{\varphi }_4=a_{5,2},\hfill \end{array}\right.\hfill \end{array}\right.\hfill \end{array}\right.& \left\{\begin{array}{cc}\hfill & \left({\varphi }_1+{\varphi }_2\right)\notin \left\{\begin{array}{ccc}{\varphi }_3,& {\varphi }_4,& {\varphi }_3+{\varphi }_4\end{array}\right\},\hfill \\ \hfill & {\mu }_1=-\frac{a_{4,2}}{{\varphi }_1+{\varphi }_2-{\varphi }_3-{\varphi }_4},\hfill \\ \hfill & {\mu }_2=\frac{{\varphi }_3{\varphi }_4-a_{5,2}}{{\varphi }_1+{\varphi }_2-{\varphi }_3-{\varphi }_4}+{\varphi }_1+{\varphi }_2,\hfill \\ \hfill & {\varphi }_1{\varphi }_2=a_{5,2}.\hfill \end{array}\right.\end{array}$$

(33)

Here the first system corresponds to the first equation ${\tilde{\mu }}_1=0$ in the set from the record (32), and the second system corresponds to the second equation ${\kappa }_{21}+{\tilde{\mu }}_2=0$ in the same set. Thus, system (32) has five qualitatively different solutions.

The first solution

$$\begin{array}{cccc}{\varphi }_1\notin \left\{\begin{array}{cc}{\varphi }_4,& {\varphi }_4-{\varphi }_3\end{array}\right\},& {\varphi }_1{\varphi }_3=a_{5,2},& {\mu }_1=\frac{a_{4,2}}{{\varphi }_4-{\varphi }_1},& {\mu }_2={\varphi }_1\end{array}$$

corresponds to matrices (30) and (27) with spectra (20) and (18) in the form

$$\begin{array}{cc}{\left.{\tilde{\mathbf{\Phi }}}_{RY0}\right|}_{{\varphi }_1\ne {\varphi }_4}=\left[\begin{array}{cc}{\varphi }_4& -\frac{a_{4,2}}{{\varphi }_4-{\varphi }_1}\\ 0& {\varphi }_3\end{array}\right],& {\left.{\tilde{\mathbf{\Phi }}}_{RY1}\right|}_{{\varphi }_1\ne {\varphi }_4}=\left[\begin{array}{cc}{\varphi }_2& \frac{a_{4,2}}{{\varphi }_4-{\varphi }_1}\\ 0& {\varphi }_1\end{array}\right]\end{array}$$

(34)

and the output regulator matrix (31) equal to

$${\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_1{\varphi }_3=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{\varphi }_2{\varphi }_4& -{\varphi }_2-{\varphi }_4& 0\\ 0& 0& -{\varphi }_1-{\varphi }_3\end{array}\right].$$

The second solution

$$\begin{array}{cccc}{\varphi }_1\notin \left\{\begin{array}{cc}{\varphi }_3,& {\varphi }_3-{\varphi }_4\end{array}\right\},& {\varphi }_1{\varphi }_4=a_{5,2},& {\mu }_1=\frac{a_{4,2}}{{\varphi }_3-{\varphi }_1},& {\mu }_2={\varphi }_1\end{array}$$

corresponds to matrices (30) and (27) with spectra (20) and (18) in the form

$$\begin{array}{cc}{\left.{\tilde{\mathbf{\Phi }}}_{RY0}\right|}_{{\varphi }_1\ne {\varphi }_3}=\left[\begin{array}{cc}{\varphi }_3& -\frac{a_{4,2}}{{\varphi }_3-{\varphi }_1}\\ 0& {\varphi }_4\end{array}\right],& {\left.{\tilde{\mathbf{\Phi }}}_{RY1}\right|}_{{\varphi }_1\ne {\varphi }_3}=\left[\begin{array}{cc}{\varphi }_2& \frac{a_{4,2}}{{\varphi }_3-{\varphi }_1}\\ 0& {\varphi }_1\end{array}\right]\end{array}$$

(35)

and the output regulator matrix (31) equal to

$${\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_1{\varphi }_4=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{\varphi }_2{\varphi }_3& -{\varphi }_2-{\varphi }_3& 0\\ 0& 0& -{\varphi }_1-{\varphi }_4\end{array}\right].$$

The third solution

$$\begin{array}{cccc}{\varphi }_2\notin \left\{\begin{array}{cc}{\varphi }_4,& {\varphi }_4-{\varphi }_3\end{array}\right\},& {\varphi }_2{\varphi }_3=a_{5,2},& {\mu }_1=\frac{a_{4,2}}{{\varphi }_4-{\varphi }_2},& {\mu }_2={\varphi }_2\end{array}$$

corresponds to matrices (30) and (27) with spectra (20) and (18) in the form

$$\begin{array}{cc}{\left.{\tilde{\mathbf{\Phi }}}_{RY0}\right|}_{{\varphi }_2\ne {\varphi }_4}=\left[\begin{array}{cc}{\varphi }_4& -\frac{a_{42}}{{\varphi }_4-{\varphi }_2}\\ 0& {\varphi }_3\end{array}\right],& {\left.{\tilde{\mathbf{\Phi }}}_{RY1}\right|}_{{\varphi }_2\ne {\varphi }_4}=\left[\begin{array}{cc}{\varphi }_1& \frac{a_{42}}{{\varphi }_4-{\varphi }_2}\\ 0& {\varphi }_2\end{array}\right]\end{array}$$

(36)

and the output regulator matrix (31) equal to

$${\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_2{\varphi }_3=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{\varphi }_1{\varphi }_4& -{\varphi }_1-{\varphi }_4& 0\\ 0& 0& -{\varphi }_2-{\varphi }_3\end{array}\right].$$

The forth solution

$$\begin{array}{cccc}{\varphi }_2\notin \left\{\begin{array}{cc}{\varphi }_3,& {\varphi }_3-{\varphi }_4\end{array}\right\},& {\varphi }_2{\varphi }_4=a_{5,2},& {\mu }_1=\frac{a_{4,2}}{{\varphi }_3-{\varphi }_2},& {\mu }_2={\varphi }_2\end{array}$$

corresponds to matrices (30) and (27) with spectra (20) and (18) in the form

$$\begin{array}{cc}{\left.{\tilde{\mathbf{\Phi }}}_{RY0}\right|}_{{\varphi }_2\ne {\varphi }_3}=\left[\begin{array}{cc}{\varphi }_3& -\frac{a_{42}}{{\varphi }_3-{\varphi }_2}\\ 0& {\varphi }_4\end{array}\right],& {\left.{\tilde{\mathbf{\Phi }}}_{RY1}\right|}_{{\varphi }_2\ne {\varphi }_3}=\left[\begin{array}{cc}{\varphi }_1& \frac{a_{42}}{{\varphi }_3-{\varphi }_2}\\ 0& {\varphi }_2\end{array}\right]\end{array}$$

(37)

and the output regulator matrix (31) equal to

$${\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_2{\varphi }_4=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{\varphi }_1{\varphi }_3& -{\varphi }_1-{\varphi }_3& 0\\ 0& 0& -{\varphi }_2-{\varphi }_4\end{array}\right].$$

The fifth solution

$$\begin{array}{cccc}{s}_{12}\notin \left\{\begin{array}{ccc}-{\varphi }_3,& -{\varphi }_4,& s_{34}\end{array}\right\},& m_{12}=a_{5,2},& {\mu }_1=\frac{a_{4,2}}{s_{12}-s_{34}},& {\mu }_2=\frac{a_{5,2}-m_{34}}{s_{12}-s_{34}}-s_{12},\end{array}$$

where $s_{12}=-{\varphi }_1-{\varphi }_2$, $s_{34}=-{\varphi }_3-{\varphi }_4$, $m_{12}={\varphi }_1{\varphi }_2$, $m_{34}={\varphi }_3{\varphi }_4$, corresponds to matrices (30) and (27) with spectra (20) and (18) in the form

$$\begin{array}{c}{\left.{\tilde{\mathbf{\Phi }}}_{RY0}\right|}_{s_{12}\ne s_{34}}=\left[\begin{array}{cc}\frac{m_{12}-m_{34}}{s_{12}-s_{34}}-s_{34}& -\frac{a_{4,2}}{s_{12}-s_{34}}\\ \frac{1}{a_{4,2}}\left(\frac{{\left(m_{12}-m_{34}\right)}^2}{s_{12}-s_{34}}+s_{12}{m}_{34}-s_{34}{m}_{12}\right)& -\frac{m_{12}-m_{34}}{s_{12}-s_{34}}\end{array}\right],\\ {\left.{\tilde{\mathbf{\Phi }}}_{RY1}\right|}_{s_{12}\ne s_{34}}=\left[\begin{array}{cc}-\frac{m_{12}-m_{34}}{s_{12}-s_{34}}& \frac{a_{4,2}}{s_{12}-s_{34}}\\ -\frac{1}{a_{4,2}}\left(\frac{{\left(m_{12}-m_{34}\right)}^2}{s_{12}-s_{34}}+s_{12}{m}_{34}-s_{34}{m}_{12}\right)& \frac{m_{12}-m_{34}}{s_{12}-s_{34}}-s_{12}\end{array}\right],\end{array}$$

(38)

and the output regulator matrix (31) equal to

$${\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_1{\varphi }_2=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{\varphi }_3{\varphi }_4& -{\varphi }_3-{\varphi }_4& 0\\ 0& 0& -{\varphi }_1-{\varphi }_2\end{array}\right].$$

Having compared the results of the five presented solutions, and also taking into account the fact that the spectra (18) and (20) can be swapped between the matrices ${\tilde{\mathbf{\Phi }}}_{RY0}$ and ${\tilde{\mathbf{\Phi }}}_{RY1}$, we draw the following conclusion. If the product of any two poles in a given spectrum (15) (we denote these poles by the symbols ${\varphi }_{y1}$ and ${\varphi }_{y2}$, and the other two poles by the symbols ${\varphi }_{x1}$ and ${\varphi }_{x2}$) is equal to a positive number $a_{52}$, then the output regulator matrix

$${\left.{\mathit{F}}_{RY}\right|}_{m_y=a_{5,2}}={\mathit{J}}_{RY}\left[\begin{array}{ccc}{m}_x& s_x& 0\\ 0& 0& s_y\end{array}\right],$$

(39)

where $s_x=-{\varphi }_{x1}-{\varphi }_{x2}$, $s_y=-{\varphi }_{y1}-{\varphi }_{y2}$, $m_x={\varphi }_{x1}{\varphi }_{x2}$ – positive constants and $m_y={\varphi }_{y1}{\varphi }_{y2}=a_{5,2}>0$, provides a spectrum

$$eig\left({\mathit{A}}_{RY}-{\mathit{B}}_{RY}{\left.{\mathit{F}}_{RY}\right|}_{{\varphi }_{y1}{\varphi }_{y2}=a_{5,2}}{\mathit{C}}_{RY}\right)=\left\{\begin{array}{cccc}{\varphi }_{x1},& {\varphi }_{x2},& {\varphi }_{y1},& {\varphi }_{y2}\end{array}\right\},$$

which, according to the Hurwitz criterion [29], corresponds to a stable system.

Matrix (39) does not contain cross coefficients between the “Roll” and “Yaw” channels. Moreover, it is robust because does not depend on the variable parameters $a_{4,2}$ and $a_{5,2}$ of the state matrix ${\mathit{A}}_{RY}$.

Next, we will consider the autonomous problem of modal output control in the “Pitch” channel, described by a completely controllable and completely observable (by state) triple of matrices (16) and spectrum (17). This is a problem with one control input, which means that the state regulator matrix is uniquely found using the Ackerman formula [30]:

$${\mathit{K}}_P=\left[\begin{array}{cc}0& 1\end{array}\right]{\left[\begin{array}{cc}{\mathit{B}}_P& {\mathit{A}}_P{\mathit{B}}_P\end{array}\right]}^{-1}\left({\mathit{A}}_P-{\varphi }_5{\mathit{I}}_2\right)\left({\mathit{A}}_P-{\varphi }_6{\mathit{I}}_2\right)=J_z\left[\begin{array}{cc}{m}_z-a_{6,3}& s_z\end{array}\right],$$

where $s_z=-{\varphi }_{z1}-{\varphi }_{z2}$, $m_z={\varphi }_{z1}{\varphi }_{z2}$, ${\varphi }_{z1}={\varphi }_5$, ${\varphi }_{z2}={\varphi }_6$. Output control is possible if the equation as below is fulfilled.

$$\begin{array}{ccc}{\mathit{K}}_P\underset{{\mathit{C}}_P^{\perp R}}{\underbrace{\left[\begin{array}{c}1\\ 0\end{array}\right]}}=0& \iff & m_z=a_{6,3}.\end{array}$$

The output regulator matrix will be

$${\left.{\mathit{F}}_P\right|}_{m_z=a_{6,3}}={\mathit{K}}_P\underset{{\mathit{C}}_P^{+}}{\underbrace{{\mathit{C}}_P^T{\left({\mathit{C}}_P{\mathit{C}}_P^T\right)}^{-1}}}=J_z{s}_z.$$

(40)

If $s_z$ is a positive constant, then, since $m_z={\varphi }_{z1}{\varphi }_{z2}=a_{6,3}>0$, matrix (40) will provide the spectrum

$$eig\left({\mathit{A}}_P-{\mathit{B}}_P{\left.{\mathit{F}}_P\right|}_{{\varphi }_{z1}{\varphi }_{z2}=a_{6,3}}{\mathit{C}}_P\right)=\left\{\begin{array}{cc}{\varphi }_{z1},& {\varphi }_{z2}\end{array}\right\},$$

which, according to the Hurwitz criterion [29], corresponds to a stable system. Moreover, matrix (40) is robust, since does not depend on the variable parameter $a_{6,3}$ of the state matrix ${\mathit{A}}_P$.

Combining the results (39) and (40), we write down the robust matrix of the output regulator, which does not contain cross coefficients between the “Roll”, “Yaw” and “Pitch” channels:

$${\left.\mathit{F}\right|}_{\begin{array}{c}{m}_y=a_{5,2}\\ m_z=a_{6,3}\end{array}}=\mathit{J}\left[\begin{array}{cccc}{m}_x& s_x& 0& 0\\ 0& 0& s_y& 0\\ 0& 0& 0& s_z\end{array}\right].$$

(41)

Thus, a robust output regulator (Figure 2) has been synthesized, in which there are no cross-connections between control channels, and the remaining coefficients do not depend on the variable parameters of the state matrix. The resulting regulator with matrix (41) is universal in the sense that it is determined only by the inertial characteristics of the lander and the desired poles, which allows its use on vehicles with any aerodynamic characteristics, including any aerodynamic quality [23,25].

It is shown that a robust regulator can be synthesized on the basis of a modified van der Wood approach using parameterization and methodically different calculations of matrices with the desired spectra at the zero and first decomposition levels (at the zero level, the output regulator is calculated, and due to the first level, its desired properties are provided). Two separate parametrized modal control problems have been solved: in the “Roll-Yaw” and “Pitch” channels. In the first task, the goal was set to zero cross-connections between the “Roll” and “Yaw” channels. To achieve it, complex sets and systems of equations and inequalities were compiled, and 5 cumbersome solutions were obtained. The complexity is caused by the presence of decomposition, and the robustness of the output regulator for a specific linear model is obtained as a concomitant factor in zeroing cross-connections. In the second task (“Pitch”), a robust solution is found by means of a special designation of the poles.

7. Selecting the Desired Poles

Let the control signals (5) be recalculated into the duration of switching on of the relay reaction engines of the system of executive organs according to the standard on-board logic.

Let us consider the choice of the desired poles used in the formation of the robust regulator output matrix (41). From the given constraints written in its subscript,

$$\begin{array}{cc}& m_y={\varphi }_{y1}{\varphi }_{y2}=a_{5,2}\left(t\right)=\tilde{q}\left(t\right){\overline{a}}_{5,2}>0,\hfill \\ & m_z={\varphi }_{z1}{\varphi }_{z2}=a_{6,3}\left(t\right)=\tilde{q}\left(t\right){\overline{a}}_{6,3}>0\hfill \end{array}$$

it can be seen that in the “Roll” channel the constant poles f1 and f2 can be assigned arbitrarily if the stability conditions are met

$$\begin{array}{cc}{s}_x=-\left({\varphi }_{x1}+{\varphi }_{x2}\right)>0,& m_x={\varphi }_{x1}{\varphi }_{x2}>0,\end{array}$$

and in the channels “Yaw” and “Pitch” the variable poles have the form

$$\begin{array}{cc}{\varphi }_{y1}\left(t\right)=-\frac{s_y}{2}-\sqrt{\frac{s_y^2}{4}-\tilde{q}\left(t\right){\overline{a}}_{5,2}},\hfill & {\varphi }_{y2}\left(t\right)=-\frac{s_y}{2}+\sqrt{\frac{s_y^2}{4}-\tilde{q}\left(t\right){\overline{a}}_{5,2}};\\ {\varphi }_{z1}\left(t\right)=-\frac{s_z}{2}-\sqrt{\frac{s_z^2}{4}-\tilde{q}\left(t\right){\overline{a}}_{6,3}},\hfill & {\varphi }_{z2}\left(t\right)=-\frac{s_z}{2}+\sqrt{\frac{s_z^2}{4}-\tilde{q}\left(t\right){\overline{a}}_{6,3}}.\end{array}$$

(42)

It is known from model (13) that the coefficients ${\overline{a}}_{5,2}$, ${\overline{a}}_{6,3}$ and $\tilde{q}\left(t\right)$ are positive. Therefore, regardless of the specific values of the constants $s_y=-\left({\varphi }_{y1}+{\varphi }_{y2}\right)>0$ and $s_z=-\left({\varphi }_{z1}+{\varphi }_{z2}\right)>0$,the poles ${\varphi }_{y1}$, ${\varphi }_{y2}$ and ${\varphi }_{z1}$, ${\varphi }_{z2}$ will be located in the left complex half-plane, ensuring, together with the poles, the stability of the closed-loop system (13) and (41).

With an appropriate choice of the values of $s_y$ and $s_z$, at the beginning of the flight section under consideration at a low velocity head q, there is a time interval where

$$\tilde{q}\left(t\right)\le \frac{1}{4}max\left(\frac{s_y^2}{{\overline{a}}_{5,2}},\frac{s_z^2}{{\overline{a}}_{6,3}}\right),$$

that is, one or both pairs of poles ${\varphi }_{y1}$, ${\varphi }_{y2}$ and ${\varphi }_{z1}$, ${\varphi }_{z2}$ consist of real numbers. This allows the process of bringing the orientation in the corresponding control channel to be aperiodic (with less overshoot and fuel consumption than during the oscillatory process). Further, with a descend of the lander and an increase in the velocity head q, a moment of time begins, starting from which

$$\tilde{q}\left(t\right)>\frac{1}{4}max\left(\frac{s_y^2}{{\overline{a}}_{5,2}},\frac{s_z^2}{{\overline{a}}_{6,3}}\right),$$

and the parameters $s_y$ and $s_z$, defining constant real parts in pairs of complex conjugate poles (42) with variable imaginary parts, provide a fixed stability margin in the Yaw and Pitch channels.

Specific positive values of the constants $s_x$, $s_y$, $s_z$ and $m_x$ in formula (41) are selected based on the results of mathematical modeling for a specific object and initial conditions (IC). The selection criterion is the accuracy of maintaining the orientation of the vehicle at a given restriction on the total fuel consumption $Q_{max}$. Accuracy ($\delta \gamma ,\delta \beta ,\delta \alpha$ in angles and $\delta w_{x_v},\delta w_{y_v},\delta w_{z_v}$ in angular velocities) is understood as the maximum modulus deviations of the state parameters from their programmed values over the time interval $[t_0+T_{trans};t_0+T]$, where $T_{trans}$ is the duration of the PP of alignment, T is the total simulation time.

As a working variant of the IC of angular motion, one of the most fuel-intensive options is used in terms of initial deviations and signs of angles and angular velocities. Further, for this option, a certain sample of stochastic descent processes (with scatter of the parameters of the Earth’s atmosphere, aerodynamics of the lander and measurement errors) is modeled without roll-overs with control according to the existing algorithm. As a result, the most probable total fuel consumption is determined, and its value is assigned to the variable $Q_{max}$. After that, one similar descent process is modeled using a new algorithm for various combinations of the constants $s_x$, $s_y$, $s_z$ and $m_x$. Under the restriction $Q_{max}$, from the simulated variants of the new algorithm, the “optimal” variant with the values $s_x^{*}$, $s_y^{*}$, $s_z^{*}$ and $m_x^{*}$ corresponding constants is selected according to the criterion

$$\Delta \left(Q_{max}\right)={\chi }_{\gamma }\delta \gamma +{\chi }_{\beta }\delta \beta +{\chi }_{\alpha }\delta \alpha +{\chi }_{w_{x_v}}\delta w_{x_v}+{\chi }_{w_{y_v}}\delta w_{y_v}+{\chi }_{w_{z_v}}\delta w_{z_v}\to \underset{m_x,s_x,s_y,s_z}{min}$$

(43)

where ${\chi }_{\gamma }=4$, ${\chi }_{\beta }={\chi }_{\alpha }=2$, ${\chi }_{w_{x_v}}={\chi }_{w_{y_v}}={\chi }_{w_{z_v}}=1$ – weight coefficients of the corresponding state parameters.

8. Numerical Example

Let us investigate the stabilization of the lander by the example of an object with typical characteristics. Consider the model (3), (4) with the parameters indicated in

Table 1. Simulation parameters.

Parameter	Value
total simulation time	$T=280\mathrm{s}$
numerical integration method	Runge-Kutta 4th Order Method
size and board steps	$\begin{array}{cc}{h}_m=0.005,& h_b=0.020\end{array}$

and the IC presented in

Table 2. Initial Simulation Conditions.

Parameter	Value
the initial moment of time (from the time of powering on the GAMS)	$t_0=280\mathrm{s}$
height above the surface of the earth	$h_0=80\mathrm{km}$
longitude (east)	$\begin{array}{cc}{h}_m=0.005,& h_b=0.020\end{array}$
longitude (east)	${\lambda }_0=51{}^{\circ }$
latitude (north)	${\phi }_0=38{}^{\circ }$
the angle of inclination of the trajectory	${\theta }_0=-1.5{}^{\circ }$
linear velocity of the CM of the lander (relative to the Earth)	$v_0=7.5\mathrm{km}/\mathrm{s}$
heading angle	${\eta }_0=0$
speed roll angle	${\gamma }_0={\gamma }_{pr}\pm 5{}^{\circ }$
angles of side slip	${\beta }_0={\beta }_{bal}\pm 5{}^{\circ }$
angles of attack	${\alpha }_0={\alpha }_{bal}\pm 5{}^{\circ }$
angular velocity of channel $x_v$	${\omega }_{x_{v}0}^{E_v}=\pm 2{}^{\circ }/\mathrm{s}$
angular velocity of channel $y_v$	${\omega }_{y_{v}0}^{E_v}=\pm 2{}^{\circ }/\mathrm{s}$
angular velocity of channel $z_v$	${\omega }_{z_{v}0}^{E_v}=\pm 2{}^{\circ }/\mathrm{s}$

. The nominal values of the characteristics of the lander, reaction engines, measuring instruments and the Earth’s atmosphere (GOST 4401-81), as well as their typical spreads [31] are used (

Table 3. Parameters of the spacecraft and engines (test version).

Channel	Moment of Inertia, kg·m${}^2$	Engine Thrust, N
$x_v$	1500	150
$y_v$	1700	75
$z_v$	1800	75
$x_v$$y_v$	50	–

).

The balancing position of the lander and the programmed value of the roll angle (for controlling the trajectory of the lander) are determined by the components of the vector

$${\mathit{\theta }}_{pr}^{E÷Q_v}={\left[\begin{array}{ccc}{\gamma }_{pr}& {\beta }_{bal}& {\alpha }_{bal}\end{array}\right]}^T={\left[\begin{array}{ccc}\pm 60{}^{\circ }\hfill & 0& -21.5{}^{\circ }\end{array}\right]}^T.$$

The position of the measuring base of the GAMS (“frozen” orbital CS at the moment of powering on the GAMS) relative to the reference base (“frozen” Greenwich CS) is set, respectively, by the angles of course, longitude and latitude

$$\begin{array}{ccc}{\eta }_{gyr}=15{}^{\circ },& {\lambda }_{gyr}=31{}^{\circ },& {\phi }_{gyr}=35{}^{\circ }.\end{array}$$

The initial moment of time, counted from the moment of powering on the GAMS, as well as the IC for the motion of the CM in the study of various processes of angular motion are taken to be the same. Combinations of the IC of angular motion (64 variants) are used in statistics when testing control algorithms.

We will call the set of ICs characteristic,

$$\begin{array}{cccc}{\gamma }_0={\gamma }_{pr}-5{}^{\circ },& {\beta }_0={\beta }_{bal}-5{}^{\circ },& {\alpha }_0={\alpha }_{bal}-5{}^{\circ },& {\omega }_{x_{v}0}^{E_v}={\omega }_{y_{v}0}^{E_v}={\omega }_{z_{v}0}^{E_v}=-2{}^{\circ }/\mathrm{s}\end{array}$$

(44)

and further we will use it for a visual graphical comparison of control processes with the same ICs on the existing and new algorithms.

MATLAB simulates 64 samples (Table 2) out of 100 lander stabilization processes without roll overturns according to the standard algorithm. The worst statistical characteristics of the processes for samples are presented in

Table 4. Statistics of lander stabilization processes without roll overturns according to the standard algorithm.

Channel		Stabilization Accuracy				Consumption (Units ${\overline{\mathit{Q}}}_{\mathit{alg}0}^{\mathit{BT}0}$)
		By Angles (°)		By Velocities (°/s)
		${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	In Channel	Total
$x_v$	ME	2.489	1.347	1.250	0.487	0.289	1.000 0.163
	MSE	0.147	0.100	0.104	0.008	0.046
$y_v$	ME	0.931	1.000	2.922	1.760	0.231
	MSE	0.018	0.006	0.030	0.007	0.009
$z_v$	ME	1.805	1.360	6.661	2.518	0.495
	MSE	0.197	0.053	0.756	0.172	0.158

, and the graphs of the corresponding process under typical IC (44) are shown in Figure 3.

The statistics contains mathematical expectations (ME) and mean square deviations (MSE) for stabilization accuracy: absolute ${\delta }_{max}^{trn+}$ (maximum deviations from programmed values after PP in 30 s) and average ${\delta }_{avr}^{all}$ (average deviations from programmed values for the entire control time). The consumptions are given in the units of the highest ME total flow rate ${\overline{Q}}_{alg0}^{BT0}$. The upper index is the number of bank turns, the lower one is the number of the algorithm (0—standard, 1—new).

Further, according to criterion (43) with fuel limitation

$$Q_{max}={\overline{Q}}_{alg0}^{BT0}-2\sigma \left(Q_{alg0}^{BT0}\right),$$

where $Q_{alg0}^{BT0}$ – fuel consumption for the existing algorithm (its ME and MSE are used) under the characteristic IC (44), the “optimal” values of the parameters of the controller matrix (41) were determined

$$\begin{array}{cccc}{s}_x^{*}=0.3,& s_y^{*}=1.4,& s_z^{*}=0.9,& m_x^{*}=0.2.\end{array}$$

(45)

MATLAB simulates 64 samples (Table 2) of 100 processes of the lander stabilization without roll overturns using a new robust algorithm with controller matrix (41) and parameters (45). The worst statistical characteristics of the processes for samples are presented in

Table 5. Statistics of lander stabilization processes without roll overturns according to the new algorithm.

Channel		Stabilization Accuracy				Consumption (Units ${\overline{\mathit{Q}}}_{\mathit{alg}0}^{\mathit{BT}0}$)
		By Angles (°)		By Velocities (°/s)
		${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	In Channel	Total
$x_v$	ME	1.438	0.727	0.453	0.147	0.235	0.944 0.148
	MSE	0.162	0.065	0.062	0.011	0.026
$y_v$	ME	0.078	0.120	0.264	0.108	0.151
	MSE	0.014	0.007	0.044	0.011	0.010
$z_v$	ME	1.412	0.771	5.492	1.427	0.578
	MSE	0.249	0.036	1.030	0.109	0.143

, and the graphs of the corresponding process under typical IC (44) are shown in Figure 4.

As an additional example of using the proposed robust algorithm, a similar simulation of the spacercraft motion with increased moments of inertia and proportionally increased thrust of control motors was carried out (

Table 6. Characteristics of the spacecraft and engines (test version).

Channel	Moment of Inertia, kg·m${}^2$	Engine Thrust, N
$x_v$	1600…3000	160…3000
$y_v$	1800…3400	90…160
$z_v$	1900…3600	90…160
$x_v$$y_v$	60…100	–

and

Table 7. Averaged statistics of spacecraft stabilization processes without roll overturns according to the new algorithm for the range of tensors of inertia and thrust of engines.

Channel		Stabilization Accuracy				Consumption (Units ${\overline{\mathit{Q}}}_{\mathit{alg}0}^{\mathit{BT}0}$)
		By Angles (°)		By Velocities (°/s)
		${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	${\mathit{\delta }}_{max}^{\mathit{trn}+}$	${\mathit{\delta }}_{\mathit{avr}}^{\mathit{all}}$	In Channel	Total
$x_v$	ME	1.843	0.882	0.558	0.175	0.259	0.985 0.151
	MSE	0.212	0.114	0.131	0.072	0.115
$y_v$	ME	0.117	0.212	0.378	0.218	0.201
	MSE	0.024	0.091	0.107	0.061	0.054
$z_v$	ME	1.534	0.934	5.552	1.529	0.525
	MSE	0.258	0.092	1.196	0.162	0.148

).

Thus, according to the results of statistical modeling in MATLAB for a new robust algorithm (Table 5 and Table 7, Figure 4) with the output feedback matrix (41), such desired poles with characteristics (45) were found that, in comparison with the standard damping algorithm (Table 4, Figure 3), the accuracy of stabilization of the lander is doubled at approximately the same consumption.

9. Conclusions

As a result of applying a new approach to the synthesis of output control, a robust regulator is analytically synthesized to stabilize the angular position of the lander when it moves in the Earth’s atmosphere. A comparative analysis of this algorithm with the corresponding algorithm currently used on board the lander of the Soyuz transport manned spacecraft is carried out. As can be seen from Table 4 and Table 5, and also shown in Figure 3 and Figure 4, the new robust control algorithm for the angular motion of the spacecraft allows doubling the spacecraft stabilization accuracy in comparison with the existing standard algorithm without increasing the fuel consumption for control.