Relative Dynamics and Control of Satellite Formation Flying Representing the Synthetic Aperture Telescope on Geostationary Orbit

1. Introduction

The need to monitor such emergencies as fires and floods in real time for the large territory of Kazakhstan requires the use of remote sensing satellites in high orbits. In this paper, to cover the largest possible area and provide high temporal and spatial resolution, it is proposed to use a satellite formation flying at GEO and representing a synthetic aperture telescope in the form of a Fizeau interferometer that provides a high spatial resolution corresponding to the resolution of an individual telescope with a high quality optical system [1]. The use of multi-aperture telescope makes it possible to obtain a high quality image of the Earth's surface with the help of 3 and more satellites representing its optical elements [2,3,4,5,6,7,8].

The considered in this paper satellite formation flying forms a telescope operating in infrared range and consisting of three mirror satellites located in one plain and imaging satellite located above this plane on the distance equal to the focal length of telescope. In other words, satellite formation flying have the form of tetrahedron or regular pyramid. High-resolution monitoring requires to keep the geometrical configuration of formation with high accuracy [9].

Several authors study the problem of developing control for the satellite formation flying representing the multi-aperture telescope. In [10] the result of control for satellite formation as three satellite interferometer is obtained. The formation is kept rigid with position error tolerance 0.2 m. Developed adaptive controller provides tracking the desired trajectories for every satellite in formation taking into account actuator saturation constraints.

In [11] the role of relative dynamics models for formation flying control system analysis was established, and the need for simple and high-fidelity dynamics models was identified. Authors considered the problem of keeping the satellite formation flying for a DARWIN-type mission in LEO and near L2 point. At first the relative dynamics models for formation flying were established taking into account appropriate perturbations. Further LQR controller was represented for formation keeping with the accuracy of 0.2 - 0.4 m.

In [12] the problem of control of satellite formation for the mission of synthetic aperture telescope at GEO is considered. The precise relative motion dynamic model of satellite formation flying is represented using Relative Orbital Elements (ROEs) with account of J22 effect. Control system is developed based on LQR technique and provided the relative position accuracy 0.15-0.2m.

Thus, to provide the operation of a synthetic aperture telescope it is necessary to have a control system that would ensure the conservation of the satellite formation configuration taking into account perturbations typical for GEO, such as the inhomogeneity of the Earth’s gravity field, perturbations arising from the attraction of the Moon and the Sun, Solar Radiation Pressure and etc. It means that much attention must be paid to the development of relative dynamic model and controller for satellite formation flying.

Many relative dynamic models are used for satellite formation flying. The most of the satellite relative motion dynamics models are represented by direct ordinary differential equation models [13,14,15,16,17] formed relative the Local-Vertical-Local-Horizontal (LVLH) frame attached to the reference (chief) satellite. Hill, Clohessy and Wiltshire [18,19] developed a linear relative motion dynamic model by neglecting the perturbation forces and assuming that the Earth gravitational field is uniform. Tschauner and Hempel [20] developed time varying linearized dynamic model of the relative motion of deputy satellite with respect to an elliptical reference orbit of the chief satellite. To improve the accuracy of relative motion dynamic model the third-body effects [21], J2 perturbations [22,23,24,25] and atmospheric drag [26,27] were taken into account. Evidently the models based on the ordinary differential equation are convenient in controller design. In this work the Sedwick-Schweighart equations with J2 effect were chosen to describe the tetrahedron satellite formation dynamics in a geostationary orbit. The mathematical model of satellite relative dynamics and problem formulation is given in the Section 2. The simulation results and analysis of uncontrolled perturbed motion of satellite formation flying is given in the Section 3.

Many works are devoted to the control of the motion of satellites in a formation based on various approaches and control techniques. In particular the authors of [28] presented a solution to the problem of maintaining a formation of nanosatellites based on the leader-follower approach using time optimal controller, LQR and H∞ controller. They made a comparison of convergence time to the required position and control effort for every controller. In [29] formation flying dynamics under consideration of J2 perturbation and nonlinear dynamics using distributed game strategy is considered. In this work control strategy to obtain a desired formation configuration, minimal energy consumption and minimal impact of disturbance on formation system is determined.

A fuzzy control system was applied to the formation. In [30] a fuzzy controller is used together with a PD-controller to motion control of the satellites formation moving in the elliptical orbits. In [31], a low-thrust fuzzy control was presented for keeping the satellite formation flying with relative position accuracy 0.02 m.

A robust controller for satellite formation flying subject to nonlinearity, parametric uncertainties, and external disturbances is presented in [32]. The proposed robust formation controller yields a position controller to form desired formation trajectories with the accuracy 0.05 m.

In this work we focus on the problem of maintaining the configuration of satellite formation flying representing the optical elements of synthetic aperture telescope for monitoring the fire in the forest. The quality characteristics of the telescope, in particular, the resolution and point spread function (PSF) of telescope, as well as quality of the image obtained with its help, is greatly influenced by the mutual arrangement of the main optical elements of the telescope. It means that to make high-quality and high-resolution observations, it is necessary to develop the control system for keeping the configuration of formation with ultrahigh accuracy. Suzumoto [1] and Rousset [3] states that to achieve a high-resolution telescope, it is required that the relative positioning accuracy of its elements reaches 1/10 of the observation wavelength. In the case of the infrared range, the accuracy of the relative positioning of satellites in the formation must be 0.1 µm. Currently such accuracy is available when using deformable mirrors [1]. In this paper we focus on the telescope using conventional mirrors and µm positioning accuracy. Thus, the main task that is solved in the process of designing a formation control system is to keep the formation configuration with an µm-class accuracy. In particular, relative position error of mirror satellites in tetrahedral formation should be kept within a range of 0.1-1µm. To provide this conditions the root method (RLM), LQR, $H_2$ optimal control, $H_{\infty }$ optimal control and mixed $H_2/H_{\infty }$ control methods are applied for synthesizing the control system for satellite formation. The mathematical model of these control methods is given in the Section 4. The simulation results for the controlled motion are given in the Section 5.

2. Problem Formulation and Mathematical Model of the Motion of Satellite in a Tetrahedral Formation

The satellite formation flying represents the synthetic aperture telescope for forest fire monitoring and consists of four satellites forming a tetrahedron or regular pyramid. The imaging satellite S₀ is located above the three mirror satellites S₁,S₂,S₃ located in the plane of the tetrahedron base and moves in an orbit called the reference orbit (Figure 1). The mirror satellites S₁,S₂,S₃ form a regular triangle, and the height of the tetrahedron passes through the center of this triangle.

For monitoring the fire area in this work it was assumed that the ground sample distance (GSD) should be equal to 30 m. To achieve this value of GSD the distance between the plane with mirror satellites and imaging satellite that corresponds to the focal length of the telescope is assumed to be equal to 21.6 m, distance between mirror satellites affecting the aperture of telescope is assumed to be equal to 26.5 m.

The utilization of distributed satellite systems for aperture synthesis presents new opportunities for Earth monitoring from space. Nevertheless, this approach also introduces novel hurdles in maintaining the satellite formation flying geometry under action of disturbances. In this work the gravitational perturbations due to non-sphericity and inhomogeneity of Earth represented by $J_2$ zonal harmonic is considered in the process of derivation of relative dynamic model for the satellite formation.

In the Figure 1 are shown following coordinate systems for describing the motion of satellite formation:

OXYZ - inertial coordinate system (ICS), the center of which is located at the Earth's center of mass, the OZ axis is directed along the Earth's rotation axis, the OX is directed to the point of the Spring equinox of the J2000 epoch;

Axyz – local vertical local horizontal (LVLH) coordinate system, the center of which is in the reference satellite represented by imaging satellite, the Ax axis is directed along the radius vector of the reference satellite from the center of the Earth, the Az axis is normal to the orbit plane in the direction of the orbital momentum, the Ay axis completes the system to the right handed system.

As it is known [33] to take into account the inhomogeneity of the Earth's gravitational field represented by $J_2$ zonal harmonic we need to write the gravitational potential in the following form:

$$U_2=-\frac{1}{2}\mu J_2\frac{r_3^2}{R^3}(3{sin}^2\phi -1),$$

(1)

where $r_{З}$ is the radius of the Earth, $R$ is the radius vector of the reference satellite in the ICS, $\phi$ is the latitude of the reference satellite in the ICS.

Further using (1) the equations of satellite motion in formation relative to the reference satellite in LVLH frame was derived in [33] in the form:

$$\begin{gathered} {\ddot{x}}_i-2nc{\dot{y}}_i-\left(5c^2-2\right)n^2{x}_i=0, \\ {\ddot{y}}_i+2nc{\dot{x}}_i=0, \\ {\ddot{z}}_i+q_i^2{z}_i=2f_i{q}_{i}cos\left(q_{i}t+\phi \right), \end{gathered}$$

(2)

$$\begin{gathered} n=\sqrt{\frac{\mu }{r_E^3}},c=\sqrt{1+s},s=\frac{3J_2{r}_E^2}{8R^2}\left(1+3\mathrm{c}\mathrm{o}\mathrm{s}2I\right), \\ {q}_i=nc-\left(\cos {\gamma }_i\sin {\gamma }_i\cot {\mathsf{\Omega }}_i-{\sin }^2{\gamma }_i\cos I_{1i}\right)\left({\dot{\mathsf{\Omega }}}_{1i}-{\dot{\mathsf{\Omega }}}_{2i}\right)-{\dot{\mathsf{\Omega }}}_{1i}\cos I_{1i}, \end{gathered}$$

(3)

$${\gamma }_i={\cot }^{-1}\left(\frac{\cot I_{2i}\sin I_{1i}-\cos I_{1i}\cos {\mathsf{\Omega }}_i}{\sin {\mathsf{\Omega }}_i}\right),{\mathsf{\Omega }}_i=\frac{z_i\left(t_0\right)}{R\sin I}$$

(4)

$$\begin{gathered} I_{1i}=\frac{v_{iz}\left(t_0\right)}{kR}+I_{2i},I_{2i}\approx I,k=nc+\frac{3n{J}_2{r}_E^2}{2R^2}{\mathit{cos}}^{2}I, \\ {\dot{\Omega }}_{1i}=-\frac{3n{J}_2{r}_E^2}{2R^2}\mathit{cos}{I}_{1i},{\dot{\Omega }}_{2i}=-\frac{3n{J}_2{r}_E^2}{2R^2}\mathit{cos}{I}_{2i}, \end{gathered}$$

(5)

$$\begin{gathered} f_i=-\frac{\mathit{sin}{I}_{1i}\mathit{sin}{I}_{2i}\mathit{sin}{\Omega }_i}{\mathit{sin}{\Phi }_i}\left({\dot{\Omega }}_{1i}-{\dot{\Omega }}_{2i}\right)R, \\ {\Phi }_i={\mathit{cos}}^{-1}\left(\mathit{cos}{I}_{1i}\mathit{cos}{I}_{2i}+\mathit{sin}{I}_{1i}\mathit{sin}{I}_{2i}\mathit{cos}{\Omega }_i\right),i=1\dots 3. \end{gathered}$$

(6)

where $x,y,z$ is the position of satellite in LVLH (position of imaging satellite), I is the inclination of orbit of reference satellite, $I_2$ is the inclination of the deputy satellites in formation represented by mirror satellites, $z_i\left(t_0\right)$ is the initial value for coordinate z in relative position of satellites in formation, $v_{iz}\left(t_0\right)$ the initial value of the velocity on the z axis, R is the radius-vector of reference satellite in ICS, i is the deputy satellite number.

To determine the deviation of the current geometrical configuration from the tetrahedral one we introduced the volume of the tetrahedron [34]:

$$V_T=\frac{1}{6}det\left|\overrightarrow{L_1},\overrightarrow{L_2},\overrightarrow{L_3}\right|,$$

(7)

where $\overrightarrow{L_1}$= $[x_1,y_1,z_1]$, $\overrightarrow{L_2}$= $[x_2,y_2,z_2]$, $\overrightarrow{L_3}$= $[x_3,y_3,z_3]$ are relative distances between the imaging satellite and mirror satellites in LVLH.

3. Analysis of Uncontrolled Motion of Tetrahedral Satellite Formation

At first the study of dynamics of passive uncontrolled motion of formation flying was performed. To assess the "run-up" (or "drift") of satellites in a formation during flight without control, graphs of changes in the distances between the reference imaging satellite - S₀ and mirror satellites S₁, S₂ and S₃ (Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9), and a graph of changes in the configuration volume (Figure 10) were obtained in the process of numerical simulation. All of these graphs represent the accuracy of relative positioning of satellites in the formation in case of perturbed uncontrolled motion.

The initial position and initial velocity of satellites in formation are given in Table 1. The required geometrical configuration of satellite formation corresponds to the initial position of satellites given in the table.

It is assumed that the reference unperturbed orbit of the imaging satellite has following parameters: the radius of the circular orbit is 42164 km, the inclination $I=0^0$,

$$n=0.000072939,c=1.000017129,q=0.000072943.$$

(8)

For simplicity of analysis of numerical simulation the following parameters was introduced:

- the magnitude of distance between the imaging satellite and the mirror satellite $L_i$, $i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of distance between the mirror satellites $l_i$, $i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of required distance between the imaging satellite and the mirror satellite $L_i\left(t_0\right)=L_{i0}$, $i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of required distance between the mirror satellites $l_i\left(t_0\right)=l_{i0}$, $i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of relative position error ${\Delta L}_i$ of the mirror satellite relative to imaging satellite, ${\Delta L}_i=L_i-L_{i0}$, $i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of relative position error ${\Delta l}_i$ of the mirror satellites, ${\Delta l}_i=l_i-l_{i0},i=\underset{\_}{\mathrm{1,3}}$;

- the magnitude of the formation volume deviation $\Delta V$ from the required value, $\Delta V=\left|\frac{V-V_0}{V_0}\right|;$

- the magnitude of the tetrahedral formation height deviation $\Delta H$ from the required value, $\Delta H=H-H_0$.

Changing of the relative position error of the mirror satellite relative to imaging satellite for time intervals, i.e. up to 100 seconds, is ambiguous for each mirror satellite (Figure 2). For example, the most intense change in distance is observed for S₂ and S₃, i.e. for the mirror satellites S₂ and S₃, and at 100 seconds it reaches the values of 26.421 m and 26.579 m, respectively. The distance between satellites S₀ and S₁ changes much more slowly; for example, at t = 100 s it equals 26.501 m. Thus, the relative position error for S₂ and S₃ is 7.9 cm (or ΔL₂ = -7.9 cm, ΔL₃ = +7.9 cm, where the minus sign means that the satellite is moving towards reference point of the local coordinate system, see Figure 2), while for S₁ it is 1 mm (or ΔL₁ = +1 mm, see Figure 2).

Figure 2. The relative position error ${\Delta L}_i,i=1,2,3$ for the time interval 100 s.

Figure 2. The relative position error ${\Delta L}_i,i=1,2,3$ for the time interval 100 s.

In the case of a passive flight of the formation without control for an hour (see Figure 3), the displacement of S₂ continues up to a value of 2.1 m. ΔL₁ increases with lesser intensity, for example, at t = 3600 with ΔL₁ = +1.3 m. Whereas within an hour ΔL₃ reaches 3.2 m, i.e. in this case, in this case ΔL₃ is already greater than ΔL₁ by almost 2.5 times. Based on this, it can be assumed that at longer time intervals, the mirror satellite S₁ will move away more intensively, while the satellite S₂ at some point in time will leave a stable orbit while keeping its constant value, or will reapproach the imaging satellite S₀.

Figure 3. The relative position error ${\Delta L}_i,i=1,2,3$ for the time interval 3600 s.

Figure 3. The relative position error ${\Delta L}_i,i=1,2,3$ for the time interval 3600 s.

In the case of passive flight, the change in the distance between mirror satellites S₁, S₂ and S₃ at short time intervals, i.e. at t ≤ 3600 s is much less than the case of changing the distances between the imaging satellite and mirror satellites (Figure 4). For example, at t = 100 s, the maximum values for satellites S₁, S₂ and S₃ are Δl₁ = Δl₂ = –1.4 mm, Δl₃ = + 0.7 mm, which is almost 56 and 113 times less than the cases of ΔL₂ and ΔL₃, while the change in Δl₁ and ΔL₁ comparable.

Figure 4. The relative position error ${\Delta l}_i,i=1,2,3$ for the time interval 100 s.

Figure 4. The relative position error ${\Delta l}_i,i=1,2,3$ for the time interval 100 s.

At t = 3600 s, the magnitude of the distances between the imaging and mirror satellites, and between the mirror satellites become almost the same order (Figure 5). For example, Δl₂ = –1.9 m, which is only 1.1 times less than ΔL₂, respectively, Δl₃ = + 0.9 m, which is 3.56 times less than ΔL₃. For the distances between the satellites S₀ and S₁, and between S₁ and S₂, the opposite picture is observed, in this case Δl₁ = –1.9 m, which is already 1.46 times greater than ΔL₁.

Figure 5. The relative position error ${\Delta l}_i,i=1,2,3$ for the time interval 3600 s.

Figure 5. The relative position error ${\Delta l}_i,i=1,2,3$ for the time interval 3600 s.

As mentioned above, to assess the “run-up” of satellites in a formation during flight without control, or in other words, to assess the keeping of the formation volume, graphs of volume deviation $\mathsf{\Delta }V$ from the required value were plotted for different time intervals ranging from 100 seconds to 1 hour. For example, with ΔV = 1, we can conclude that the volume has not changed from the initial value.

For short time intervals, a slow but monotonous increase in volume is observed (Figure 6, cases (a) and (b)). For example, at 100 seconds, the volume increases by only 0.005% compared to the required value (Figure 6, case (a)), while with an hour of passive flight, this value already reaches 6.5% (Figure 6, case (b)). For most tasks, these volume fluctuations are insignificant, and should in no way affect to the successful completion of the mission, which cannot be admissible for the formation representing the synthetic aperture telescope. This mission requires relative position accuracy of about 0.1-1μm, which leads to active control of the formation in the early stages of the motion.

Figure 6. The magnitude of volume deviation $\Delta V$ from the required value for different time intervals, cases: (a) t = 100 s, (b) t = 3600 s.

Figure 6. The magnitude of volume deviation $\Delta V$ from the required value for different time intervals, cases: (a) t = 100 s, (b) t = 3600 s.

Table 2. Table of maximum values for relative position errors of satellites.

t	max(ΔL1)	max(ΔL2)	max(ΔL3)	max(Δl1)	max(Δl2)	max(Δl3)
0-100 s	+1 mm	-7.9 сm	+7.9 сm	-1.4 mm	-1.4 mm	+ 0.7 mm
0-3600 s	+1.3 m	-2.1 m	+3.2 m	-1.9 m	-1.9 m	+ 0.9 m

Table 2. Table of maximum values for relative position errors of satellites.

t	max(ΔL1)	max(ΔL2)	max(ΔL3)	max(Δl1)	max(Δl2)	max(Δl3)
0-100 s	+1 mm	-7.9 сm	+7.9 сm	-1.4 mm	-1.4 mm	+ 0.7 mm
0-3600 s	+1.3 m	-2.1 m	+3.2 m	-1.9 m	-1.9 m	+ 0.9 m

4. Control System for Tetrahedron Satellite Formation

As it was stated earlier, in this paper we consider the problem of developing a control system for a tetrahedral satellite formation flying in the geostationary orbit for high temporal and spatial resolution monitoring of fire in the forest area of Kazakhstan.

The satellite formation flying forms a synthetic aperture telescope operating in infrared range. The main task that is solved by this formation control system is to keep the formation geometrical configuration with a µm-class accuracy. In particular, relative position error of mirror satellites in tetrahedral formation should be kept within a range of 0.1-1 µm.

The equations for the controlled motion of a satellite in a formation has the form (2) with the control forces assumed as the linear function

$$\begin{gathered} u_x=-k_x\left(x_T-x\right)-k_{Vx}{v}_x, \\ {u}_y=-k_y\left(y_T-y\right)-k_{Vy}{v}_y, \\ {u}_z=-k_z\left(z_T-z\right)-k_{Vz}{v}_z, \end{gathered}$$

(9)

where $x_T,y_T,z_T$ are required position of satellite in formation relative the reference coordinate system, $v_x,v_y,v_z$ are relative velocity of satellites in formation, $k_x$, $k_{Vx}$, $k_y,k_{Vy},k_z,k_{Vz}$ are constant coefficients.

To determine the unknown feedback coefficients $k_x$, $k_{Vx}$, $k_y,k_{Vy},k_z,k_{Vz}$ in expressions for the control acceleration (9) we applied several approaches of linear control theory:

1) root method (RLM) providing the required response speed, damping and stability of the system;

2) linear quadratic regulator (LQR) minimizing the state error and control effort of the system;

3) $H_2$ optimal control minimizing the system's sensitivity to disturbances;

4) $H_{\infty }$ optimal control providing the robustness of the system to disturbances;

5) mixed $H_2/H_{\mathrm{\infty }}$ control allowing the flexible trade-off between performance and robustness.

Let’s consider every approach in detail.

Root method (RLM)

At the first stage we considered the classical approach for the root method providing the special location of the closed-loop poles of the system in the complex plain.

As can be seen from (9) the third equation of the system does not contain the variables included in the first and second equations, then it can be considered individually. We will look for a solution to the equation in the form:

$$z=C{e}^{\lambda t}.$$

(10)

From here we find the corresponding $\dot{z}$,$\ddot{z}$and introduce the corresponding values into the third equation of system (2). Making transformations we obtain an algebraic equation:

$$C\left({\lambda }^2+k_{Vy}\lambda +q^2+k_y\right)=0.$$

(11)

Obviously, equation (12) is a polynomial of the second degree of the form:

$${a_0\lambda }^2+a_1\lambda +a_2=0.$$

(12)

Consider the first and second equations of system (2). We will look for their solutions in the form:

$$x=A{e}^{\lambda t},y=B{e}^{\lambda t}.$$

(13)

From here we find the corresponding time derivatives and put them into the equations of system (2). Making transformations, we obtain the corresponding system of algebraic equations:

$$\begin{gathered} A\left({\lambda }^2+k_{Vx}\lambda +k_x-\left(5n^2{c}^2-2n^2\right)\right)-2nc\lambda B=0, \\ B\left({\lambda }^2+k_{Vy}\lambda +k_y\right)+2nc\lambda A=0. \end{gathered}$$

(14)

Since the system of equations (14) must have a non-zero solution for A and B, the determinant of this system is equal to zero:

$$\left|\begin{array}{cc}{\lambda }^2+k_{Vx}\lambda +k_x-\left(5n^2{c}^2-2n^2\right)& -2nc\lambda \\ 2nc\lambda & {\lambda }^2+k_{Vy}\lambda +k_y\end{array}\right|=0$$

(15)

or

$$\begin{gathered} {\lambda }^4+{\lambda }^3\left(k_{Vx}+k_{Vy}\right)+{\lambda }^2\left(k_x+k_y+k_{Vx}{k}_{Vy}-\left(5n^2{c}^2-2n^2\right)-4n^2{c}^2\right) \\ +\lambda (k_x{k}_{Vy}+k_{Vx}{k}_y-k_{Vy}(5n^2{c}^2-2n^2\left)\right)+k_x{k}_y \\ -(5n^2{c}^2-2n^2)k_y=0 \end{gathered}$$

(16)

$${a_0\lambda }^4+a_1{\lambda }^3+a_2{\lambda }^2+a_3\lambda +a_4=0.$$

(17)

We define the required roots of the second-order characteristic equation as the roots of the Butterworth polynomial since it provides maximally flat and monotonic transition process [35] in the form:

$${\lambda }^2+1.414{\Omega }_r\lambda +{{\Omega }_r}^2=0,$$

(18)

where ${\Omega }_r=\frac{t_{н}}{t_p}$, $t_{н}$ is the normalized time of the transient process, $t_p$ is the real time of the transient process.

And the roots of the characteristic equation of the fourth order as the roots of the Butterworth polynomial in the form:

$${\lambda }^4+2.613{{\Omega }_r\lambda }^3+3.4141{{\Omega }_r}^2{\lambda }^2+2.613{{\Omega }_r}^3\lambda +{{\Omega }_r}^4=0$$

(19)

The roots of Butterworth polynomial locate on the circle of radius ${\Omega }_r$ in the left half of complex plain.

Further equating in formulas (12) and (18) terms at the similar power of $\lambda$ we obtain expressions for determining the unknown feedback coefficients:

$$k_{Vy}=1.414{\Omega }_r,k_y={{\Omega }_r}^2-q^2.$$

(20)

Equating in formulas (17) and (19) terms at the similar power of $\lambda$ we obtain a system of algebraic equations for determining the unknown feedback coefficients:

$$\begin{gathered} k_{Vx}+k_{Vy}=2.613{\Omega }_r, \\ {k}_x+k_y+k_{Vx}{k}_{Vy}-\left(5n^2{c}^2-2n^2\right)-4n^2{c}^2=3.4141{{\Omega }_r}^2 \\ {k}_x{k}_{Vy}+k_{Vx}{k}_y-k_{Vy}(5n^2{c}^2-2n^2)=2.613{{\Omega }_r}^3, \\ {k}_x{k}_y-\left(5n^2{c}^2-2n^2\right)k_y={{\Omega }_r}^4. \end{gathered}$$

(21)

At the following stage we considered another approach for the root method proving the location of the closed-loop poles of the system in LMI regions.

LMI regions are convex subsets D of the complex plane characterized by [36]:

$$D=\{z\in C:L+Mz+M^T\overline{z}<0\},$$

(22)

where $M={\left\{{\mu }_{ij}\right\}}_{1\le i,j\le m}$ and $L=L^T={\left\{{\lambda }_{ij}\right\}}_{1\le i,j\le m}$ are fixed real matrices. The matrix-valued function

$$f_D\left(z\right):=L+Mz+M^T\overline{z}$$

(23)

is called characteristic function of the region $D$.

The approach of control synthesis as the function $\overrightarrow{u}=K\overrightarrow{X}$ that provides the robust pole placement in LMI regions usually is considered for the linear systems of the type:

$$\begin{array}{c}\dot{\overrightarrow{X}}=A\overrightarrow{X}+B_1\overrightarrow{w}+B_2\overrightarrow{u},\\ \overrightarrow{Z}=C_1\overrightarrow{X}+D_{11}\overrightarrow{w}+D_{12}\overrightarrow{u},\\ \overrightarrow{Y}=C_2\overrightarrow{X}+D_{12}\overrightarrow{w}+D_{22}\overrightarrow{u},\end{array}$$

(24)

where $\overrightarrow{X}$ is the state vector of the system, $\overrightarrow{w}$ is the exogenous input, $\overrightarrow{u}$ is control, $\overrightarrow{Z}$ is the regulated output, $\overrightarrow{Y}$ is the measured output.

The closed-loop poles of (24) lie in the LMI region (22) if and only if there exists a symmetric positive definite matrix $P$ satisfying [36]:

$$\begin{array}{c}{\left[{\lambda }_{ij}P+{\mu }_{ij}\left(A+B_{2}K\right)P+{\mu }_{ij}P+{\mu }_{ji}P{(A+B_{2}K)}^T\right]}_{1\le i,j\le m}<0\\ P>0\end{array}$$

(25)

with the notation ${\left[S_{ij}\right]}_{1\le i,j\le m}=\left(\begin{array}{ccc}{S}_{11}& \cdots & S_{i1}\\ ⋮& \ddots & ⋮\\ S_{m1}& \cdots & S_{mm}\end{array}\right)$.

LMI region can be represented individually or as intersection of several LMI regions. In this article this region is given as the intersection of a disk of radius $r$ with origin at the point $(-h,0)$ and a conical sector with inner angle $\theta$:

$$f_D\left(z\right):=\left(\begin{array}{cc}-r& \overline{z}+h\\ z+h& -r\end{array}\right),$$

(26)

$$f_D\left(z\right):=\left(\begin{array}{cc}sin\frac{\theta }{2}(z+\overline{z})& -cos\frac{\theta }{2}(z-\overline{z})\\ cos\frac{\theta }{2}(z-\overline{z})& sin\frac{\theta }{2}(z-\overline{z})\end{array}\right).$$

(27)

The matrices of equation (24) used for synthesis of controller has the following form:

$$A=\left[\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ \left(5c^2-2\right)n^2\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ {-q}^2\end{array}\end{array}\end{array}\begin{array}{c}1\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}-2nc\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 1\\ \begin{array}{c}0\\ 2nc\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}1\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\right],B_1=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\right],B_2=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ \begin{array}{c}0\\ 1\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}1\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 1\end{array}\end{array}\end{array}\right]$$

(28)

$$C_1=diag[\mathrm{1,1},\mathrm{1,1},\mathrm{1,1}],C_2=diag[\mathrm{1,1},\mathrm{1,1},\mathrm{1,1}],$$

(29)

$$D_{11}=0_{6\times 3},D_{12}=[0_{3\times 3};diag[\mathrm{1,1},1\left]\right],D_{21}=0_{6\times 3},D_{22}=0_{6\times 3}$$

(30)

Linear quadratic regulator (LQR)

The second method to synthesize the controller for keeping the satellite formation with high accuracy is a linear quadratic regulator (LQR), which is obtained by minimizing the quality criterion of the form [37]:

$$J=\frac{1}{2}\int \left({\overrightarrow{∆r}}^{T}Q\overrightarrow{∆r}+{\overrightarrow{u}}^{T}R\overrightarrow{u}\right)dt$$

(31)

where $\overrightarrow{∆r}=[x_T-x,y_T-y,z_T-z,v_x,v_{y,}{v}_z]$, Q, R are positive matrices with constant components, $x_T$, $y_T$, $z_T$ - is the required position of the satellite in a formation.

The control acceleration $\overrightarrow{u}=[u_x,u_y,u_z$] for the linear system in the form (24) obtained as result of minimizing the performance criterion (31) has the form [37]:

$$\overrightarrow{u}=-R^{-1}{B}^{T}P\overrightarrow{∆r}=K\overrightarrow{∆r},$$

(32)

where $B=B_2$ introduced in (28), the matrix R has the form:

$$R=0.0001·diag[1,1,1].$$

(33)

And the matrix P can be determined from the equation:

$$A^{T}P+PA-PB{R}^{-1}{B}^{T}P+Q=0,$$

(34)

where $A$, $B=B_2$ are the matrices introduced in (28), the matrix Q has the form:

$$Q=0.0001·diag[\mathrm{1,1},1,1,1,1].$$

(35)

$H_2$ optimal control

$H_2$ optimal controller is considered for the linear systems with a state space representation similar to (24):

$$\begin{array}{c}\dot{\overrightarrow{X}}=A\overrightarrow{X}+B_1\overrightarrow{w}+B_2\overrightarrow{u},\\ \begin{array}{c}\overrightarrow{Z}=C_1\overrightarrow{X}+D_{12}\overrightarrow{u},\\ \overrightarrow{Y}=C_2\overrightarrow{X}+D_{12}\overrightarrow{w}+D_{22}\overrightarrow{u},\end{array}\end{array}$$

(36)

where $D_{12}^T{D}_{12}>0$.

The state feedback synthesis control problem reduces to developing a controller $\overrightarrow{u}=K\overrightarrow{X}$ that minimizes the $H_2$ norm of the transfer function $H_{\mathrm{z}\mathrm{w}}\left(s\right)$ [38]:

$${‖H_{\mathrm{z}\mathrm{w}}\left(s\right)‖}_2^2=\frac{1}{\pi }{\int }_0^{\infty }trace\left({H_{\mathrm{z}\mathrm{w}}\left(j\omega \right)}^T{H}_{\mathrm{z}\mathrm{w}}\left(j\omega \right)\right)d\omega$$

(37)

To calculate this norm we need to solve the following problem [38]:

$$\begin{array}{c}\min trace\left[\left(C_1+D_{12}K\right)P{(C_1+D_{12}K)}^T\right]\\ s.t.\left(A+B_{2}K\right)P+P{(A+B_{2}K)}^T+B_1{B}_1^T\le 0,\\ P=P^T>0\end{array}$$

(38)

Introducing $X_2=X_2^T=P$, $L=KP$, $Q=Q^T$ and using appropriate transformations (38) is equivalent to LMI problem:

$$\begin{array}{c}{\min }_{Q,X_2,L}trace\left(Q\right)\\ s.t.A{X}_2+X_2{A}^T+B_{2}L+L^T{B}_2^T+B_1{B}_1^T\le 0,\\ \left[\begin{array}{cc}{X}_2& X_2{C}_1^T+L^T{D}_{12}^T\\ C_1{X}_2+D_{12}L& Q\end{array}\right]\ge 0.\end{array}$$

(39)

Using the matrices $A$, $B_1$, $B_2$, $C_1$, $D_{12}$ introduced in (28), (29), (30) and finding $(X_2,L)$ satisfying (39) the state feedback gain matrix is calculated as $K=L{X}_2^{-1}$.

$H_{\mathrm{\infty }}$ optimal control

$H_{\infty }$ optimal controller is considered for the linear systems with a state space representation similar to (24). The state feedback synthesis control problem reduces to developing a controller $\overrightarrow{u}=K\overrightarrow{X}$ that minimizes the $H_{\mathrm{\infty }}$ norm of the transfer function $H_{\mathrm{z}\mathrm{w}}\left(s\right)$ [39]:

$${‖H_{\mathrm{z}\mathrm{w}}\left(s\right)‖}_{\infty }={\sup }_{\omega \in R}{\sigma }_{max}\left(F\right(j\omega \left)\right)$$

(40)

To calculate this norm we need to solve the following problem [39]:

$$\begin{array}{c}{\min }_{P,K,\gamma }\gamma \\ s.t.\left[\begin{array}{ccc}{(A+B_{2}K)}^{T}P+P(A+B_{2}K)& {PB}_1& {(C_1+D_{12}K)}^T\\ B_1^{T}P& -\gamma I& D_{11}^T\\ C_1+D_{12}K& D_{11}& -\gamma I\end{array}\right]<0,\\ P>0.\end{array}$$

(41)

Introducing $X_{\infty }=X_{\infty }^T=P$, $L=KP$ and using appropriate transformations (41) is equivalent to:

$$\begin{array}{c}{\min }_{X_{\infty },L,\gamma }\gamma \\ s.t.\left[\begin{array}{ccc}{(A{X}_{\infty }+B_{2}L)}^T+(A{X}_{\infty }+B_{2}L)& B_1& {(C_1+D_{12}L)}^T\\ B_1^T& -\gamma I& D_{11}^T\\ C_1+D_{12}L& D_{11}& -\gamma I\end{array}\right]<0,\\ X_{\infty }>0.\end{array}$$

(42)

Using the matrices $A$, $B_1$, $B_2$, $C_1$, $D_{12}$ introduced in (28), (29), (30) and finding $(X_{\infty },L)$ satisfying (42) the state feedback gain matrix is calculated as $K=L{X}_{\infty }^{-1}$.

Mixed$H_2/H_{\infty }$control

Mixed $H_2/H_{\mathrm{\infty }}$ control is considered for the linear systems with a state space representation:

$$\begin{array}{c}\dot{\overrightarrow{X}}=A\overrightarrow{X}+B_1\overrightarrow{w}+B_2\overrightarrow{u},\\ \begin{array}{c}\overrightarrow{Z_{\infty }}=C_{\infty }\overrightarrow{X}+D_{\infty 1}\overrightarrow{w}+D_{\infty 2}\overrightarrow{u},\\ \overrightarrow{Z_2}=C_2\overrightarrow{X}+D_{21}\overrightarrow{w}+D_{22}\overrightarrow{u},\end{array}\\ \overrightarrow{Y}=C_y\overrightarrow{X}+D_{y1}\overrightarrow{w},\end{array}$$

(43)

where $\overrightarrow{X}$ is the state, $\overrightarrow{u}$ is the control, $\overrightarrow{Y}$ is the output, where $\overrightarrow{w}\to \overrightarrow{Z_{\infty }}$ is the $H_{\mathrm{\infty }}$ performance channel, $\overrightarrow{w}\to \overrightarrow{Z_2}$ is the $H_2$ performance channel.

Mixed $H_2/H_{\mathrm{\infty }}$ output feedback controller for the plant (43) can be obtained in the form [40]:

$$\begin{array}{c}\dot{\overrightarrow{x_k}}=A_k\overrightarrow{x_k}+B_k\overrightarrow{Y},\\ \overrightarrow{u}=C_k\overrightarrow{x_k}+D_k\overrightarrow{Y}.\end{array}$$

(44)

This controller stabilizes the system exponentially, the $H_{\mathrm{\infty }}$ performance channel has a performance level ${‖H_{\overrightarrow{w}\to \overrightarrow{Z_{\infty }}}\left(K\right)‖}_{\infty }\le \gamma$, the $H_2$ performance ${‖H_{\overrightarrow{w}\to \overrightarrow{Z_2}}\left(K\right)‖}_2$ is minimized among all $K$ stabilizing the system exponentially and having prescribed $H_{\mathrm{\infty }}$ performance level $\gamma$.

Thus, the mixed $H_2/H_{\mathrm{\infty }}$-synthesis problem is the optimization program [40]:

$$\begin{array}{c}\min {‖H_{\overrightarrow{w}\to \overrightarrow{Z_2}}\left(K\right)‖}_2\\ s.t.{‖H_{\overrightarrow{w}\to \overrightarrow{Z_{\infty }}}\left(K\right)‖}_{\infty }\le \gamma ,\\ Kstabilizessysteminternally,\end{array}$$

(45)

where $H_{\overrightarrow{w}\to \overrightarrow{Z_2}}(K,s)$ denotes the transfer function of the $H_2$ closed-loop performance channel, $H_{\overrightarrow{w}\to \overrightarrow{Z_{\infty }}}\left(K\right)$ denotes the transfer function of the $H_{\mathrm{\infty }}$ closed-loop performance channel, $\gamma$ is a threshold.

The solution of presented problem is presented in [40].

The matrices $A$, $B_1$, $B_2$, $C_{\infty }$, $D_{\infty 1}$, $D_{\infty 2}$, $C_2$, $D_{21}$, $D_{22}$, $C_y$, $D_{y1}$ used for synthesis of $H_2/H_{\mathrm{\infty }}$ controller has the following form correspondingly:

$$A=\left[\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ \left(5c^2-2\right)n^2\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ {-q}^2\end{array}\end{array}\end{array}\begin{array}{c}1\\ 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}-2nc\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 1\\ \begin{array}{c}0\\ 2nc\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\begin{array}{c}0\\ 0\\ \begin{array}{c}1\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\right],B_1=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 0\end{array}\end{array}\end{array}\right],B_2=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ \begin{array}{c}0\\ 1\\ \begin{array}{c}0\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}1\\ 0\end{array}\end{array}& \begin{array}{c}0\\ 0\\ \begin{array}{c}0\\ 1\end{array}\end{array}\end{array}\right]$$

(46)

$$C_{\infty }=diag[\mathrm{1,1},\mathrm{1,1},\mathrm{1,1}],C_2=diag[\mathrm{1,1},\mathrm{1,1},\mathrm{1,1}],C_y=diag[\mathrm{1,1},\mathrm{1,1},\mathrm{1,1}],$$

(47)

$$D_{\infty 1}=0_{6\times 3},D_{\infty 2}=[0_{3\times 3};diag[\mathrm{1,1},1\left]\right],D_{21}=0_{6\times 3},D_{22}=0_{6\times 3},D_{y1}=0_{6\times 3}$$

(48)

5. Results and Discussion of Numerical Simulation of Tetrahedron Satellite Formation under Control

Numerical modeling of the satellite formation dynamics was carried under action of RLM, LQR, $H_2$, $H_{\infty }$ controllers, mixed $H_2/H_{\infty }$ controllers. Initial conditions for numerical integration and initial geometrical configuration are given in Table 1. The problem of keeping the actual positions of satellites in this initial geometrical configuration is considered in this section.

For numerical simulation it was assumed that the reference unperturbed orbit of the imaging satellite has the following parameters: the radius of the circular orbit is 42000 km, the inclination $I=0^0$.

As in the section 3, for simplicity the following parameters are used in this section for the description of numerical results:

- the relative position error ${\mathsf{\Delta }L}_i$ of the mirror satellite relative to imaging satellite, ${\mathsf{\Delta }L}_i=L_i-L_{i0}$, $i=\overline{\mathrm{1,3}}$;

- the relative position error ${\mathsf{\Delta }l}_i$ of the mirror satellites, ${\mathsf{\Delta }l}_i=l_i-l_{i0},i=\overline{\mathrm{1,3}}$;

- magnitude of the formation volume deviation $\mathsf{\Delta }V$ from the required value, $\mathsf{\Delta }V=\left|\frac{V-V_0}{V_0}\right|$;

- magnitude of the tetrahedral formation height deviation $\mathsf{\Delta }H$ from the required value, $\mathsf{\Delta }H=H-H_0$.

The plots of relative position error ${\Delta L}_i$ of the mirror satellite $S_i$ relative to imaging satellite $S_0$ are given in the Figure 7, Figure 8 and Figure 9. The required distance between the imaging satellite and the mirror satellite $L_{i0}$ is taken as $L_{i0}=26.5m$.

The plots of the relative position error ${\Delta l}_i$ of the mirror satellites are given in the Figure 10, Figure 11 and Figure 12. The required distance between the mirror satellites $l_{i0}$ is taken as $l_{i0}=26.5m$.

The plot of magnitude of the formation volume deviation $\Delta V$ from the initial required value $V_0$ is given in the Figure 13. The magnitude of control force for imaging satellites $S_0,S_1,S_2$ are given in the Figure 15, Figure 16 and Figure 17. Below in Table 3 are given the results of analysis of controllers:

- the quality analysis of these controllers on the base of the parameters of accuracy, transient time, overshoot, damping character;

- quality analysis of the satellite formation configuration based on the volume character of the formation;

- analysis of control effort for keeping the configuration with required accuracy.

In accordance with the results of simulation modeling mixed $H_2/H_{\infty }$ controller provides best results in terms of transient time and accuracy, as the maximum relative position error for mirror satellites is 0.17 $\mu m$ and maximum relative position error for mirror satellites relative to imaging satellite is 1129 $\mu m$ under the maximum control effort of 0.231 N. This control provides the minimal changes in height and volume of satellite formation (Figure 13 and Figure 14) that guarantee the keeping of formation configuration in tetrahedral form with high accuracy. $H_{\infty }$ controller shows the similar results in accuracy, transient time, control effort, character of volume and height of satellite formation. It provides the maximum relative position error for mirror satellites of 0.20 $\mu m$ and maximum relative position error for mirror satellites relative to imaging satellite of 1366 $\mu m$. Mixed $H_2/H_{\infty }$ controller gives better accuracy and performance in the process of keeping the configuration as it combines the advantages of $H_2$ and $H_{\infty }$ control techniques. In contrast the pure $H_2$ controller results in oscillatory character of relative distances variation between the satellites and doubles the transient time (Figure 7, Figure 8, Figure 9, Figure 10, Figure 11 and Figure 12), whereas the maximum relative position error for mirror satellites is 0.21 $\mu m$ and maximum relative position error for mirror satellites relative to imaging satellite is 1477 $\mu m$ under the maximum control effort of 0.271 N. This accuracy is close to $H_{\infty }$ controller, but damping character of transient process is worse.

LQR controller shows good results in stabilization of transient process, stability of volume and height of the satellite formation, but in contrast to $H_{\infty }$ controller and mixed $H_2/H_{\infty }$ controller it provides worse accuracy: the maximum relative position error for mirror satellites is 0.27 $\mu m$ and maximum relative position error for mirror satellites relative to imaging satellite is 1932 $\mu m$.

For the root method (RLM) two approaches of poles placement in complex plane is considered in this paper and corresponding RLM controller and RLM2 controller was obtained. RLM controller provides the placement of poles of control system on the circle of radius 1.2 in left part of complex plane. RLM2 controller provides the placement of poles of control system in the region bounded by a disk of radius 2 with origin at the point (-2, 0) and a conical sector with apex at the point (-0.1,0) and inner angle $\frac{3\pi }{4}$. This technique allows to achieve good transient process quality such as transient time and damping characteristics. RLM2 controller provides the accuracy similar in values to $H_2$ and $H_{\infty }$ controllers: the maximum relative position error for mirror satellites is 0.22 $\mu m$ and maximum relative position error for mirror satellites relative to imaging satellite is 1544 $\mu m$. RLM is worse in accuracy. Both RLM and RLM2 controllers provide almost the similar transient time 9-10 seconds, but in contrast to other considered here controllers it is more sensitive to disturbances, that results in non-permanent character of satellite formation volume and height (Figure 13 and Figure 14) and it can lead to the violation of tetrahedral configuration which is not valid for a mission presented in this article.

The analysis made it possible to conclude that $H_{\infty }$ controller and mixed $H_2/H_{\infty }$ controller is the most appropriate for providing the satellite formation mission that forms a telescope for remote sensing, since it provides the necessary accuracy of satellite relative positioning, the quality of transient processes and keep the configuration with the high accuracy. Regarding the results of simulation and thrust calculation (Figure 15, Figure 16 and Figure 17) in case of using cold gas thrusters with specific impulse 70 seconds the necessary fuel consumption for mirror satellite is 9.86 kg per year. It indicates the necessity to improve the controller for decreasing the amount of necessary fuel for satellites in formation.

Figure 7. The relative position error ${\Delta L}_1$.

Figure 7. The relative position error ${\Delta L}_1$.

Figure 8. The relative position error ${\Delta L}_2$.

Figure 8. The relative position error ${\Delta L}_2$.

Figure 9. The relative position error ${\Delta L}_3$.

Figure 9. The relative position error ${\Delta L}_3$.

Figure 10. The relative position error ${\mathsf{\Delta }l}_1$.

Figure 10. The relative position error ${\mathsf{\Delta }l}_1$.

Figure 11. The relative position error ${\mathsf{\Delta }l}_2$.

Figure 11. The relative position error ${\mathsf{\Delta }l}_2$.

Figure 12. The relative position error ${\mathsf{\Delta }l}_3$.

Figure 12. The relative position error ${\mathsf{\Delta }l}_3$.

Figure 13. Magnitude of formation volume deviation $\mathsf{\Delta }V$.

Figure 13. Magnitude of formation volume deviation $\mathsf{\Delta }V$.

Figure 14. Magnitude of the tetrahedral formation height deviation $\mathsf{\Delta }H$.

Figure 14. Magnitude of the tetrahedral formation height deviation $\mathsf{\Delta }H$.

Figure 15. Magnitude of control force for satellite S1.

Figure 15. Magnitude of control force for satellite S1.

Figure 16. Magnitude of control force for satellite S2.

Figure 16. Magnitude of control force for satellite S2.

Figure 17. Magnitude of control force for satellite S3.

Figure 17. Magnitude of control force for satellite S3.

6. Conclusion

The problem of development the control system for a satellite formation flying in the form of synthetic aperture telescope on the geostationary orbit for monitoring of fire in the forest area of Kazakhstan is considered in this paper. To achieve a high spatial resolution it was assumed that the formation has definite tetrahedral configuration in the form of synthetic aperture telescope operating in infrared range and consists of three mirror satellites and one imaging satellite. Synthesizing the image with good quality requires to keep the configuration of formation with µm-class accuracy and provide minimum volume deviation of the configuration from its required value. It has been shown that in uncontrolled flight volume conservation is not a sufficient condition for successful missions, since for a single configuration volume value there are many configurations. To obtain the best quality of images, it is necessary to keep the relative position of satellites with 0.1-1µm-class accuracy. Currently adaptive optics and MEMS deformable mirror can provide such accuracy. In this paper we focused on the synthetic aperture telescope using conventional mirrors and µm-class positioning accuracy.

To keep the configuration in the required form with µm accuracy the analysis of efficiency of various control methods, such as the root method (RLM), LQR, $H_2$ optimal control, $H_{\infty }$ optimal control and mixed $H_2/H_{\infty }$ control method, was carried out in the process of numerical simulation. The simulation results made it possible to highlight the following features of various approaches to the development of control system for the satellite formation flying representing a synthetic aperture telescope:

- RLM controllers is sensitive to disturbances, that results in non-permanent character of satellite formation volume and height and it can lead to the violation of satellite configuration which is not valid for a mission presented in this article;

- $H_{\infty }$ controller and mixed $H_2/H_{\infty }$ controller provide best results in terms of transient time and accuracy;

- LQR controller shows good results in stabilization of transient process, stability of volume and height of the satellite formation, but in contrast to controller and mixed $H_2/H_{\infty }$ controller it provides worse accuracy;

- $H_2$ controller results in oscillatory character of relative distances variation between the satellites and doubles the transient time.

Thus, in accordance with the results of analysis, it possible to conclude that $H_{\infty }$ controller and mixed $H_2/H_{\infty }$ controller is the most appropriate for providing the satellite formation mission that forms a synthetic aperture telescope for remote sensing, since it provides the necessary accuracy of satellite relative positioning, the quality of transient processes and keep the configuration geometry with the high accuracy.