- Email: [email protected]
Two-satellite formation flying control by cell-structured solar sail
Journal Pre-proof Two-satellite formation flying control by cell-structured solar sail Yaroslav Mashtakov, Mikhail Ovchinnikov, Tatyana Petrova, Stepan Tkachev PII:
S0094-5765(20)30086-2
DOI:
https://doi.org/10.1016/j.actaastro.2020.02.024
Reference:
AA 7890
To appear in:
Acta Astronautica
Received Date: 23 September 2019 Revised Date:
7 February 2020
Accepted Date: 12 February 2020
Please cite this article as: Y. Mashtakov, M. Ovchinnikov, T. Petrova, S. Tkachev, Two-satellite formation flying control by cell-structured solar sail, Acta Astronautica, https://doi.org/10.1016/ j.actaastro.2020.02.024. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier Ltd on behalf of IAA.
1 2
Two-satellite formation flying control by cell-structured solar sail 1
3
Yaroslav Mashtakov, Mikhail Ovchinnikov, Tatyana Petrova, Stepan Tkachev ,
4 5 6 7 8 9 10 11 12 13 14 15
Keldysh Institute of Applied Mathematics of RAS, Miusskaya Square, 4, 125047 Moscow, Russia Abstract The method of simultaneous relative motion and attitude control by a solar radiation pressure is proposed. The control goal is to stabilize the given closed relative orbit. The principle idea is to use special materials for solar sail that able to change their optical properties. The solar sail consists of a number of rectangular cells. Each of them can change its optical properties. The necessary control force is produced by varying the average reflectivity of solar sail, and the control torque is achieved by the appropriate pattern of black and white cells. Key words: formation flying, solar radiation pressure, variable radiation reflectivity, attitude control, relative motion control
16
Introduction
17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39
Equation Section (Next) Utilization of a group of satellites, for example, formation flying brings new opportunities in space missions. In addition, group of satellites is more reliable because even if one satellite fails, others can continue the mission operation. The main problem of formation flying utilization is deployment and maintenance of the particular group configuration. The simplest solution for this problem is to use thrusters installed on board of all or several satellites. On the other hand, thrusters require propellant, which can greatly affect the satellite lifetime or the payload mass. To overcome this problem environmental forces for formation flying motion control are used [1]. This approach can be applied relatively easily by installing a special high area-to-mass ratio equipment like a flat sail. There are two sorts of forces that can be used: aerodynamic drag [2–6] and solar radiation pressure (SRP) [7– 12]. The principal idea here is to use a difference in environmental forces acting on each satellite in formation. This difference usually appears when sail changes its attitude with respect to the sun direction, although the effective size variation is also considered in literature [13]. In general, solar sails can be used either for attitude control [14–18] or for orbit control [19,20]. In the paper the case when both attitude and relative motion are controlled via solar sail with variable optical properties is considered. Simultaneous orbit and attitude control using solar sail with no additional devices seems impossible. However, there are number of approaches that permit both orbit and attitude control using either sail shape changing [21] or reflectivity properties variation [13], [22–24]. The latter approach is of interest is the paper. There are several techniques to use variable optical properties. In [13] it is considered that sail consists of two part with different optical properties each and the size of these parts can change continuously. Usually, the more realistic model is used when there is only a set of cells with variable optical properties are place one the sail. This approach was used during well1
Corresponding author E-mail address: [email protected]
40 41 42 43 44 45 46 47 48 49 50 51
known mission IKAROS [25]. In [22] the problem of appropriate pattern determination for the required torque and force implementation is aroused in case of sail with cell structure and the table of possible patterns is created. However this approach can’t be scaled for the large amount of cells since it requires too much memory to save all patterns. Another approach is described in [24]. The pattern is built step-by-step: for achieving attitude control and then – orbit control. Thanks to the symmetry of the sail, it is divided into four quadrants. In each quadrant the square domain which will provide necessary in-plane torque is selected. These domains are the same ones and are located symmetrical with respect to the sail symmetry center. The signs of the inplane torque define the ‘on-off’ state of each domain and the values – the coordinates of the center of the domain and its square. Finally, for orbit control symmetrical cells are added to the pattern. In this paper, similar approach is used but it differs in details, control synthesis scheme and application.
52
1. Problem statement and reference frames
53 54 55 56
Deployment and maintenance of two satellites in the required relative orbit are considered. Each satellite is assumed to be equipped with a solar sail. The initial orbit of leadersatellite is circular. The follower-satellite is moving along the orbit which is close to the first one. Satellites move under the solar radiation pressure and are both disturbed by J 2 harmonic of
57 58 59
the Earth gravitational field. The following right-handed reference frames are used: – O1 XYZ is the inertial frame (IF) with the origin at the Earth center of mass, O1Z is
60 61 62 63 64 65 66 67
perpendicular to the equator plane, O1 X is directed to the vernal equinox; Оxyz is the orbital frame (OF), its origin is in the leader satellite’s center of mass, Оz is directed along its radius vector, Оy is perpendicular to the orbit plane; Оξηζ is the body-fixed frame (BF), its axes are the principal axes of inertia (it is also assumed that Оζ is perpendicular to the sail plain); Оxs y s z s is the solar frame (SF), Оzs is directed to the Sun, Оys is perpendicular to the ecliptic plane. Transition between IF and OF is performed by the following matrix
68 69 70
71
– – –
S = ( e1 e 2
e 3 ) , e3 =
r r×v , e2 = , e1 = e 2 × e 3 , r |r × v |
where r is the radius vector and v is the velocity of the leader satellite. Transition between IF and SF is determined by the matrix Ssun = ( l1
l2
cos Λ 0 l 3 ) , l 3 = sin Λ cos ε , l 2 = − sin ε , l1 = l 2 × l 3 . sin Λ sin ε cos ε
72
Here Λ is the ecliptic longitude, ε is the obliquity of the ecliptic.
73
Problem statement and reference frames
74
There are three types of motion equations that are used in this paper.
75
Orbital dynamics
76
78
Orbital dynamics is described by the following vector equation r & r&= − µE 3 + g , r where µ E is the Earth gravity constant and g is the result vector of the external disturbances. As
79
it was mentioned before, only the effects of J 2 and SRP force are taken into account. The first
80
has a form
77
81
3sin 2 i sin 2 u − 1 3J µ R f J 2 = 2 E4 − sin 2 i sin 2u . 2r − sin 2i sin u
82
Here J 2 = 1.082 × 10−3 , R⊕ is the mean Earth radius, i is the orbit inclination and u is the
83
argument of latitude. SRP force acting on the elementary area dS can be written as follows [26] Φ dFs = − 0 ( rs , n ) × c 2 × (1 − α ) rs + 2αµ ( rs , n ) n + α (1 − µ ) rs + n dS , 3
84
2 ⊕
85
where Φ 0 = 1357 W m is the solar flux constant, rs is the unit vector from the Sun to the
86
satellite, n is the solar sail normal (SSN), α is the reflection coefficient, µ is the specularity
87
coefficient. Further the case of µ = 1 is considered, i.e. there is no diffuse reflection, so
2
dFs = −
88 89 90 91
Φ0 (rs , n ) ( (1 − α ) rs + 2α ( rs , n ) n ) dS . c
Since the parameter α varies from point to point, the total SRP force becomes Φ Fs = − 0 ( rs , n ) S − ∫ α dS rs + 2 ( rs , n ) n ∫ α dS . c
((
Denote f =
∫ α dS ( 0 ≤ f S
≤ 1) and A = −
)
)
Φ0S , then c
Fs = A ( rs , n ) ( (1 − f ) rs + 2 f ( rs , n ) n ) .
92 93
These equations are written for both satellites and are used in numerical simulation.
94
Attitude dynamics
95
Angular dynamics is described in the BF by the Euler equations
96
&+ ω × Jω = M control + M grav , Jω
97
where J is the satellite inertia tensor, ω is the satellite angular velocity with respect to the IF,
98
M control is the control torque and M grav = 3
99 100
µE
(1)
r × Jr is the gravity gradient torque. r5 Attitude kinematics (BF w.r.t. IF) is defined by the quaternion ( λ0 , λ ) , λ02 + λ 2 = 1 and described by the corresponding equations
λ&0 = −0.5 ( λ , ω ) ,
101 102
λ&= 0.5 ( λ0ω + λ × ω) . These equations are used for the numerical simulation and control torque synthesis.
103 104
Relative motion dynamics
105 106
The control synthesis for orbital relative motion is based on the Hill-Clohessy-Wiltshire equations written in the curvilinear parameters ( a0θ r , a0ϕ r , ρ ) (see Fig. 1). These parameters are
107 108 109
more preferable than the regular ones since it is assumed that the leader satellite is moving along circular orbit while the relative orbit is small with respect to the size of the orbit. Motion equations in the OF can be written as follows [27]
110
d2 d a ϕ + 2ω ρ = 0, 2 ( 0 r) dt dt 2 d a θ + ω 2 ( a0θ r ) = 0, 2 ( 0 r) dt d2 d ρ − 2ω ( a0ϕ r ) − 3ω 2 ρ = 0 2 dt dt
111
where ω is the orbital angular velocity of the leader satellite, ρ = r2 − r1 , a0 = r1 . Index “1”
112
corresponds to the leader satellite and “2” to the follower.
113 114 115
116
Figure 1. Curvilinear parameters If control is taken into account then Eq.(2) can be rewritten in the following form
d2 d a ϕ + 2ω ρ = u x , 2 ( 0 r) dt dt 2 d a θ + ω 2 ( a0θ r ) = u y , 2 ( 0 r) dt d2 d ρ − 2ω ( a0ϕ r ) − 3ω 2 ρ = uz . 2 dt dt
(2)
117 118
Here u x , uy , u z are the components of the control vector u =
Fs ,2 − Fs ,1 . m
Solution of (2) is a0ϕ r = −3C1ωt + 2C2 cos ωt − 2C3 sin ωt + C4 ; a0θ r = C5 cos ωt + C6 sin ωt;
119
ρ = 2C1 + C2 sin ωt + C3 cos ωt. 120
One can introduce new variables B j , ψ j based on this solution a0ϕ r = 2 B2 cosψ 1 + B3 ,
ρ = B2 sinψ 1 + 2 B1 , d ( a0ϕ r ) = −2 B2ω sinψ 1 − 3B1ω, dt d ρ = B2ω cosψ 1 , dt a0θ r = B4 cosψ 2 ,
121
d ( a0θ r ) = − B4ω sinψ 2 . dt 122
Below the explicit definition of parameters B j and ψ j is given B1 = C1 , B2 = C22 + C32 , B3 = −3C1ωt + C4 ,
123
B4 = C52 + C62 , C3
ψ 1 = ωt + β1 , where sin β1 =
C22 + C32
ψ 2 = ωt + β 2 , where sin β 2 = − 124
C6 C52 + C62
C2 C22 + C32
, cos β 2 =
,
C5 C52 + C62
.
The equations corresponding to these variables are .
B1 =
1
ω
ux ,
.
B3 = −3B1ω − .
B2 = 125
, cos β1 =
1
ω
2
ω
uz ,
( uz cosψ 1 − 2ux sinψ 1 ) ,
1 ψ1 = ω − ( uz sinψ 1 + 2ux cosψ 1 ) , B2ω .
.
B4 = − .
1
ω
ψ2 = ω −
u y sinψ 2 , 1 u cosψ 2 . ω B4 y
(3)
126 127
Here B1 corresponds to the drift velocity of the follower satellite along axis Ox of the OF, B2 is the size of the in-plane ellipse, B3 is the shift of this ellipse along Ox . Finally, B4 is
128 129
the out-of-plane motion amplitude. The control purpose is to achieve the required Bi , i = 1, 2, 3, 4 parameter, i.e. to have the required relative orbit.
130
Control synthesis
131
The general control synthesis scheme is presented in Fig. 2.
132 133 134 135
Figure 2. Control synthesis scheme First of all, the ideal control that provides required relative motion is found. Then corresponding integral reflection coefficient f i and angles of SSN θi , ϕ i are determined.
136 137
Normal directions define the reference angular motion of each satellite. After that the control torque Mi , control (in Оξη plane) is calculated. Finally, Mi , control and f i determine the solar sail
138
optical pattern. Further in this section each step is discussed.
139
Relative motion control
140 141
Purpose of the control is to deploy the formation and maintain the required relative orbit. This orbit is defined by Bi ( i = 1, 2,3, 4 ). In the paper the relative orbit defined by
142
B1 = 0, B2 = B0 , B3 = 0, B4 = 0
143
is considered. This means that relative orbit is the ellipse with major and minor semiaxes 2B0
144
and B0 respectively lying in the orbital plane. Its center coincides with the origin of the OF. The
145
relative orbit stabilization is divided into two stages: firstly B1 = 0 and B3 = 0 are provided, then
146
B2 = B0 is achieved. The out-of-plane motion control process is separated, so B4 = 0 can be
147 148
guaranteed independently. On the first stage the following Lyapunov control function (LCF) is used
149
V=
1 2 1 2 B1 + B3 . 2 2
150
Its time derivative in accordance with (3) is . . 1 2 V&= B1 B1 + B3 B3 = B1ux + B3 −3B1ω − uz . ω ω
151 152
So the control that ensures global asymptotic stability of B1 = 0 and B3 = 0 is the following u x = −k1 B1 , k1 > 0,
153
uz =
1 ( −3B1ω 2 + k2ω B3 ) , k2 > 0 2
(4)
154
by the Barbashin-Krassovskiy theorem [28]. It should be noted when B1 is large the value of
155
required u z might exceed the accessible control effort umax (in the paper this value is 10−6 N / kg
156
). Hence, one should wait until B1 becomes small which is possible since the first equation of (3)
157 158
is independent. On the second stage the following LCF is used
V=
159 160 161 162
and its time derivative
1 1 V&= ( B1 − 2 ( B2 − B0 ) sinψ 1 ) ux + ( ( B2 − B0 ) cosψ 1 − 2 B3 ) uz − 3B1B3ω.
ω
ω
Then the control is u x = − k x ( B1 − 2 ( B2 − B0 ) sinψ 1 ) , k x > 0
163 164
1 2 1 1 2 B1 + ( B2 − B0 ) + B32 2 2 2
u z = − k z ( ( B2 − B0 ) cosψ 1 − 2 B3 ) , k z > 0
(5)
which leads to
165
2 2 1 1 V&= − ( B1 − 2 ( B2 − B0 ) sinψ 1 ) − ( ( B2 − B0 ) cosψ 1 − 2 B3 ) − 3B1B3ω.
166 167 168 169 170
The last term does not depend on control and its sign is undefined. However, its value after the first stage would be smaller than the sum of the first and second terms. This is why the first stage is required. In case the magnitude of the control force received from the Barbashin-Krassovskiy theorem is larger the maximum allowed control force is set to the umax . Practically, in both cases
171
(4) and (5) the following scheme was used
ω
ω
172
−u sign( a ), a > umax u = max −a, | a |≤ umax
173 174
where a can be either first or second expression in (5) or first expression in (4). Additionally, since B1 is the drift velocity, it could be used to control the convergence speed of B3 to zero. So
175
with the help of nonzero drift the shift might be reduced faster. If B3 is large, then instead of the
176
first expression of (4) the following control is used
177
ux = −k1 ( B1 − γ sign ( B3 ) )
178
where γ is a some positive parameter. The out-of-plane motion control is
u y = −k y B4 cosψ 2 , k y > 0 .
179 180
The control u x , uy , u z is an ideal one. The control is called ideal in terms of being
181 182
received analytically. It should be implemented through the solar sails rotation and integral reflectivity coefficients f i .
183
Relative motion control implementation
184 185
Let θ be the angle between the SSN and the Sun direction, ϕ is the angle between normal projection to the plane Оx s y s of the SF and axis Оxs . Then in the SF the SSN is sin θ cos ϕ n = sin θ sin ϕ . cos θ
186
187
Relative motion control force u in the SF is
uxs = 2 Af 2 cos2 θ2 sin θ 2 cos ϕ 2 − 2 Af1 cos2 θ1 sin θ1 cos ϕ1 , 188
u ys = 2 Af 2 cos2 θ2 sin θ2 sin ϕ 2 − 2 Af1 cos2 θ1 sin θ1 sin ϕ1 ,
(6)
uzs = A (1 − f 2 ) cos θ 2 − A (1 − f1 ) cos θ1 + Af 2 cos3 θ 2 − 2 Af1 cos3 θ1. 189
These expression are the non-linear equations with respect to f i , ϕ i θi . Moreover, there are
190 191
only three equations for six unknown variables. This freedom is used for the maximization of the possible values of the u components. First of all, as the SRP force is decreasing when θi tends
192
to 90 degrees, it is reasonable to consider the case of small θi , so (6) transforms to
uxs = 2 Af 2θ 2 cos ϕ 2 − 2 Af1θ1 cos ϕ1 , 193
u ys = 2 Af 2θ 2 sin ϕ 2 − 2 Af1θ1 sin ϕ1 ,
(7)
uzs = Af 2 − Af1. 194
From the last equation one can see, that f i define u zs completely. The maximum torque
195
is achieved when f = 0.5 , while for f = 0 and f = 1 the torque is zero. It is reasonable to
196
demand f i to be as close to 0.5 as possible, so min J = ( f1 − 0.5 ) + ( f 2 − 0.5 ) 2
2
f1 , f 2
subject to
197 f 2 − f1 = 0 < f min
198
u zs
A ≤ fi ≤ f max < 1
The solution of this problem in the inner domain is
u zs , 2A u f 2 = 0.5 + zs . 2A f1 = 0.5 −
199
(8)
200
It exists when
2 f min − 1 ≤
201
uzs ≤ 2 f max − 1 A
(9)
202
In the further discussion it is supposed that f min = 0.25 and f max = 0.75 . So (9) might be
203
rewritten as
−0.5 ≤
204 205
uzs ≤ 0.5 . A
When f i is known one can find ϕ i from the maximization (with θi fixed) of
206
max L = ( f 2θ 2 cos ϕ 2 − f1θ1 cos ϕ1 ) + ( f 2θ 2 sin ϕ 2 − f1θ1 sin ϕ1 ) .
207
It means that the result values of ϕ i give the largest domain of the possible control force. This
208
problem has two groups of solutions
2
2
ϕ1 , ϕ 2
209
ϕ1 = ϕ 2 , θ1θ 2 < 0 ,
210
ϕ1 = ϕ 2 + π , θ1θ2 > 0 .
211 212
The actual position of the SSN is the same for both solutions. Hence, further the case ϕ1 = ϕ 2 = ϕ is considered. The first and the second equations in (7) one can rewrite as follows
uxs , 2A u ( f 2θ2 − f1θ1 ) sin ϕ = ys . 2A
( f 2θ2 − f1θ1 ) cosϕ = 213
214
Then tgϕ =
215
f 2θ 2 − f1θ1 =
216 217
218 219
u ys u xs
,
u xs2 + u 2ys 2A
.
If uxs cos ϕ > 0 f 2θ 2 − f1θ1 =
u xs2 + u 2ys 2A
.
To find θi one can solve the following optimization problem min L = θ12 + θ 22 θ1 ,θ 2
subject to 220 f 2θ 2 − f1θ1 =
u xs2 + u ys2 2A
−θ max ≤ θ i ≤ θ max
(10)
221
The solution of this problem gives the minimum possible values of θi that satisfy constraints.
222
This keeps θi in the vicinity of zero.
223
The solution in the inner domain is
θ1 = − 224
θ2 =
uxs2 + u ys2
f1 , f + f 22 2 1
2A uxs2 + u 2ys
f2 . f + f 22 2 1
2A
225
Thus, once u is known, attitude of the SSN can be found.
226
Attitude control
227 228
The expression for the SRP torque for the square sail and given sun direction in the BF
rs = ( sin θ cos β
cos θ ) is T
Mξ P cos θ Φ0 −Q cos θ M η = − c cos θ − Q sin θ sin β + P sin θ cos β M ζ
229
230
sin θ sin β
where P=
231
a 2
a 2
∫ ∫ ηα dξ dη ,
Q=
a a − − 2 2
232 233 234 235
a 2
a 2
∫ ∫ ξα dξ dη. −
a a − 2 2
There are only two components that can be defined independently by the optical properties variation and SSN attitude. So, there can be only two independent control torques to be produced. The following LCF is considered
V=
( (
))
1 T 2 2 J ξ ωrel,1 + Jηωrel,2 + ka 1 − ( 0 0 1) , Bn . ( ) 2
(11)
236
Here J ξ , Jη are the in-plane moments of inertia, ωrel ,1 , ωrel ,1 are the corresponding relative
237
angular velocity components ( ωrel = ω − ωref ), ωref = n × n&, n is the required SSN attitude in the
238
IF, ka > 0 and B is the transition matrix between the IF and the BF. The goal of the control is to
239
guarantee asymptotic stability of the motion when the axis Oζ of the BF and n coincide.
240
Derivative of (11) is
241
T d V&= J ξ ωrel,1ω&rel,1 + Jηωrel,2ω&rel,2 − ka ( 0 0 1) , ( Bn ) . dt
242
As there is no need to control the third component of the angular velocity let us take relative
243
angular velocity vector as ωrel = (ωrel ,1 ωrel ,2
T
d ( Bn ) = −ωrel × Bn dt
244 245
0 ) . In this case
and (12) one can rewrite as follows
(12)
. T V&= ωTrel J ωrel + ka Bn × ( 0 0 1) .
246 247
To guarantee V&< 0 it is sufficient to satisfy the following equation
248
&rel + ka Bn × ( 0 0 1) = − kω ω rel . Jω
249 250 251
T
It guarantees by the following control
&ref − ka Bn × ( 0 0 1) . Mcontrol = −kωωrel − Mext + ω × Jω − Jω × Bωref + JBω T
(13)
Here only two first components of M control are taken. The last component is defined by
252
the other ones.
253
Solar sail pattern
254 255 256
The final step of the control synthesis is the determination of the solar sail pattern for each satellite. This pattern must give the required integral reflection coefficients, as well as control torque components M ξ and M η . In general, the reachability domain is the inner area of
257
rhombus LMNK which is presented in Fig. 3.
258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273
Figure 3. The torque reachability domain The components of control torque, which based on Lyapunov theory, are independent. However, the domain between rhombus and square the components of control torque could not be chosen independently that means that the less component of torque could not be increased without decreasing the larger one and vice versa. To overcome this, the control torque components are limited by the square. Each solar sail is considered to be divided into n × n cells. Each cell can either reflect radiation totally or absorb it. The problem is to find the fully reflective cells that give the required torque components and integral reflectivity coefficient. There are various approaches. Here one of them is presented, which requires neither complex calculations nor large amount of memory. Let initially each cell have the reflectivity coefficient α = 0 . The required torque will be produced by the rectangle with α = 1 . Its center is located at the edge of the square, which is described in the following way: its center coincides with sail center, and the edge length equals to a half of sail's edge length (Fig. 4). The SRP torque which is produced by this rectangle is
ξ Fs ,ξ η Fs ,ζ M ξ Ms = r × Fs = η × Fs ,η = −ξ Fs ,ζ = Mη , 0 F ξ F −ηF M s ,ζ η ξ ζ
274 275
where r = (ξ η 0) is the center position vector. Therefore, T
M ξ =− η, η Mξ 2
a ξ + η = , 4 2
sign (η ) = −sign ( M ξ ) .
276
277 278 279 280
2
Figure 4. Solar sail pattern
The area of the rectangle is determined by the following expressions
281
S rec =
282
Fs ,ζ =
cFs ,ζ Φ 0 cos2 θ
.
M ξ2 + Mη2
ξ 2 +η2
283 284
The rectangle with a given center and area is not unique, but it should fit cell pattern on the sail surface. First of all it should be noted that
285
a S rec 2
286 287
since the control parameters in (13) are chosen to guarantee control torque to be in the square in Fig. 3. There are several cases: 2
288 289
2
a 1. ≤ S rec which means that the area corresponds to, at least, one cell and the sides of the n rectangle are taken as the side of the square rounded to the nearest cell.
2
2
290
a a 2. ≥ S rec then the area is taken as S rec = , but the center of rectangle is moving to the n n
291
ξ point k
η k
T
0 . Here
k=
292
M ξ2 + Mη2
ξ 2 + η 2 Fζ
.
293 294
When k ≥ n the distance from the center of sail to the center of the rectangle is less than the half-size of the cell. In this case Sпр = 0 and no torque is produced.
295
Once the cell pattern for torque implementation is determined, f should be realized. The
296 297
number of cells with α = 1 is determined from expression N req = f ⋅ n 2 .
298 299 300
However, there are already N cells with α = 1 which realize the control torque. In general, there are several options: 1. N = N req which means that there is no need in cell pattern modifications and required f is
301 302
achieved. 2. N < N req then the symmetric pairs with zero torque are added. If N − N req is an odd number
303 304
then one cell is left (this is an error in force realization). 3. N > N req . Torque producing cells give pattern with f ≤ 0.25 since no more than quarter of
305
sail has α = 1 . On the other hand in (9) minimum force producing f min = 0.25 . So this situation
306 307 308 309
is avoided. Thus, the control synthesis is finished. The presented scheme allows determining the solar sail pattern which gives the required control force and torque under given constraints.
310
NUMERICAL EXAMPLE
311 312 313 314
Numerical simulation is carried out using the parameters presented in Table 1. Integration method is Runge-Kutta with constant step of one second. The control force and torque are calculated in each integration step. After that the control is substituted to the orbital and angular models for each satellites. And, finally, the relative motion parameters ( Bi ) are calculated. These
315 316
parameters are presented in Figures 5-7. Table 1. Numerical simulation parameters Parameter Value Leader orbit radii Rorb = 9060 km Inclination 110 deg Initial relative orbit rrel = ( 200 100 50 ) m, Mass of each satellite Solar sail size
v rel = ( 0.05 0.5 1) m/s; m = 10 kg 5×5 m
Inertia tensor
J = diag ( 2.1 2.1 3.8 ) kg ⋅ m 2
Initial angular velocity
ω1 = ( 0.005 0.003 0.001) s−1 , ω2 = ( 0.001 0.003 0.002 ) s −1;
Attitude of the SSN
θ1 = 10 deg, θ 2 = 10 deg;
Control law parameters
k1 = 1 s −2 , k2 = 2 s−2 , k x = 0.5 s −2 , k y = 0.1 s−2 , k z = 0.05 s −2 , kω = 0.02 N ⋅ m ⋅ s, ka = 10−4 N ⋅ m
Maximum allowed control force
umax = 10−6 N
Maximum allowed control torque
M = 3 ⋅ 10−5 N ⋅ m
317 318 319
Results of the numerical simulation are presented in Figs 5-12. Figs 5-8 show the evolution of the parameters Bi which determine the relative orbit (drift, shift, orbit size and
320 321 322 323
amplitude of the out-of-plane motion). In Fig. 8 the relative orbit is shown. The blue lines here show relative in-plane motion over the whole time interval. The red ones illustrate resulting steady state relative in-plane motion. The numerical simulation is carried out for 4x4 (left column) and 50x50 (right column) cells sail.
Figure 5a. Drift and shift evolution for 4x4 cells sail
Figure 5b. Drift and shift evolution for 50x50 cells sail
324
Figure 6a. Relative orbit size evolution for 4x4 cells sail
Figure 6b. Relative orbit size evolution for 50x50 cells sail
Figure 7a. Out-of-plane motion amplitude for 4x4 cells sail
Figure 7b. Out-of-plane motion amplitude for 50x50 cells sail
Figure 8a. Relative in-plane motion evolution for 4x4 cells sail
Figure 8b. Relative in-plane motion evolution for 50x50 cells sail
Figure 9a. Angle between SSN and sun direction evolution for 4x4 cells sail
Figure 9b. Angle between SSN and sun direction evolution for 50x50 cells sail
325 326 327
The control successfully solve the task and formation reach the required relative orbit (Figs 5-9). It also can be seen in Fig. 5 that B1 ≈ −20 m during first : 15 revolutions. This
328
allows reaching required B3 value much faster. The Figures 9 show that θ is bounded in the
329 330
desired region. Figures 10a and 10b contain the switching rate of the solar cells for the 4x4 cell sails for both satellites.
Figure 10a. Solar sail cells switching rate (the black one is the highest, the white one is the lowest) for the chief satellite
Figure 10a. Solar sail cells switching rate (the black one is the highest, the white one is the lowest) for the deputy satellite
331 332 333 334 335 336 337
From Figure 10 one can see that there are four cells which almost all the time have constant reflectivity. So the ordinary solar cell material can replace these cells. The maximum usage is approximately 40% (black cell in Firuges). Which means it switches every 0.4 control iteration. The simulation results for 2x2 cells sail are shown in Fig. 11. This configuration fails to solve the stabilization problem. The similar situation is in case of 3x3 cells sail. So, it seems that 4x4 is the minimal configuration when the proposed control scheme solves the problem.
338 339
Figure 11. Drift and shift for 2x2 cells sail
340 341 342 343 344
The main problem here is that the 2x2 configuration does not allow producing the different torque components with the same sign and torque value is of no importance (only sign have an effect). The 3x3 meets difficulties on the stage when symmetric pairs with zero torque must be added. The 50x50 cells configuration provide in several times better accuracy but this is obviously more complex from the practical point of view.
345
Conclusion
346 347 348
In the paper the scheme of the two satellites formation flying control using the solar sail is proposed. It allows deploying and maintaining the required relative motion using the solar sail with variable optical properties only. In continuous case this scheme guaranties asymptotic
349 350 351 352 353 354 355 356 357 358
stability for both angular and relative motion. When the sails consist of discrete cells with variable optical properties the minimum number of cells appears to be sixteen i.e. four along each side. Note that the proposed control scheme is suitable for the time intervals when the Clohessy-Wiltshire approximation is valid. Under long-term SRP effects the eccentricity and other orbital parameters of the chief satellite evolve. For example, for the case considered in the paper the maximum eccentricity value is about 10−3 which still can be considered as perturbation (see Appendix). The long-term modification of the proposed scheme can be the subject of the future studies.
359
Acknowledgments
360 361 362
The work is supported by the RFBR Grant No. 18-31-20014. Authors are grateful to Prof. Anna Guerman for productive discussions.
363
References
364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393
1. 2.
3.
4.
5. 6. 7. 8. 9. 10. 11.
12. 13.
Kumar K.D. et al. Maintenance of satellite formations using environmental forces // Acta Astronaut. Elsevier, 2014. Vol. 102. P. 341–354. Leonard C.L., Hollister W.M., Bergmann E. V. Orbital formation keeping with differential drag // J. Guid. Control. Dyn. American Institute of Aeronautics and Astronautics, 1989. Vol. 12, № 1. P. 108–113. Kumar B.S., Ng A. A bang-bang control approach to maneuver spacecraft in a formation with differential drag // Proc. AIAA Guid. Navig. Control Conf. Exhib. 2008. № August. P. 1–11. Zeng G., Ru M., Yao R. Relative orbit estimation and formation keeping control of satellite formations in low Earth orbits // Acta Astronaut. Elsevier, 2012. Vol. 76. P. 164– 175. Varma S., Kumar K.D. Multiple Satellite Formation Flying Using Differential Aerodynamic Drag // J. Spacecr. Rockets. 2012. Vol. 49, № 2. P. 325–336. Omar S.R., Wersinger J.M. Satellite Formation Control using Differential Drag // 53rd AIAA Aerosp. Sci. Meet. 2015. № January. P. 1–11. Gong S., Yunfeng G., Li J. Solar sail formation flying on an inclined Earth orbit // Acta Astronaut. Elsevier, 2011. Vol. 68, № 1–2. P. 226–239. Shahid K., Kumar K.D. Multiple spacecraft formation reconfiguration using solar radiation pressure // Acta Astronaut. Elsevier, 2014. Vol. 103. P. 269–281. Shahid K., Kumar K.D. Formation Control at the Sun-Earth L(2) Libration Point Using Solar Radiation Pressure // J. Spacecr. Rockets. 2010. Vol. 47, № 4. P. 614–626. Gong S., Baoyin H., Li J. Solar Sail Formation Flying Around Displaced Solar Orbits // J. Guid. Control. Dyn. 2007. Vol. 30, № 4. P. 1148–1152. Williams T., Wang Z.-S. Solar Radiation Pressure and Formation-Keeping in Highly Elliptical Orbits // AIAA/AAS Astrodynamics Specialist Conference and Exhibit. 2002. № August. P. 1–14. Parsay K., Schaub H. Preliminary Results on Optimal Establishment of Solar Sail Formations // J. Astronaut. Sci. 2019. Vol. 66, № 1. P. 32–45. Mu J., Gong S., Li J. Reflectivity-controlled solar sail formation flying for magnetosphere mission // Aerosp. Sci. Technol. Elsevier Masson SAS, 2013. Vol. 30, № 1. P. 339–348.
394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430
14. 15. 16.
17. 18. 19. 20. 21. 22.
23.
24.
25. 26. 27. 28.
Wie B. Solar Sail Attitude Control and Dynamics, Part 1 // J. Guid. Control. Dyn. 2004. Vol. 27, № 4. P. 526–535. Wie B., Murphy D. Solar-sail attitude control design for a sail flight validation mission // J. Spacecr. Rockets. 2007. Vol. 44, № 4. P. 809–821. Wie B. et al. Robust attitude control systems design for solar sails, part 2: MicroPPTbased secondary ACS // Collection of Technical Papers - AIAA Guidance, Navigation, and Control Conference. 2004. Vol. 2, № August. P. 1325–1340. Bolle A., Circi C. Solar sail attitude control through in-plane moving masses // Proc. Inst. Mech. Eng. Part G J. Aerosp. Eng. 2008. Vol. 222, № 1. P. 81–94. Scholz C. et al. Performance analysis of an attitude control system for solar sails using sliding masses // Adv. Sp. Res. COSPAR, 2011. Vol. 48, № 11. P. 1822–1835. Bookless J., McInnes C. Dynamics and control of displaced periodic orbits using solar-sail propulsion // J. Guid. Control. Dyn. 2006. Vol. 29, № 3. P. 527–537. McInnes C.R. Passive control of displaced solar sail orbits // J. Guid. Control. Dyn. 1998. Vol. 21, № 6. P. 975–982. Wu R. et al. Heliogyro Solar Sail with Self-Regulated Centrifugal Deployment Enabled by an Origami-Inspired Morphing Reflector // Acta Astronaut. IAA, 2018. Borggrafe A. et al. Attitude Control of Large Gossamer Spacecraft using Surface Reflectivity Modulation // Proceedings of the 65th International Astronautical Congress. 2014. P. 8. Funase R. et al. On-orbit verification of fuel-free attitude control system for spinning solar sail utilizing solar radiation pressure // Adv. Sp. Res. COSPAR, 2011. Vol. 48, № 11. P. 1740–1746. Biggs J.D., Negri A. Orbit-Attitude Control in a Circular Restricted Three-Body Problem Using Distributed Reflectivity Devices // J. Guid. Control. Dyn. 2019. Vol. 42, № 12. P. 2712–2721. Tsuda Y. et al. Flight status of IKAROS deep space solar sail demonstrator // Acta Astronaut. Pergamon, 2011. Vol. 69, № 9–10. P. 833–840. Fliegel H.F., Gallini T.E., Swift E.R. Global Positioning System Radiation Force Model for geodetic applications // J. Geophys. Res. 1992. Vol. 97, № B1. P. 559. Alfriend K. et al. Spacecraft Formation Flying: Dynamics, Control and Navigation. 1 Edition. Butterworth-Heinemann, 2010. 402 p. David R. Merkin. Introduction to the Theory of Stability / ed. Afagh F.F., Smirnov A.L. Springer-Verlag New York, 1997. 320 p.
Appendix The perturbed equations (3) are
.
B1 =
1
ω
( ux + g x ) ,
.
B3 = −3B1ω − .
B2 =
1
ω
((u
.
ψ1 = ω − .
B4 = − .
1
ω
ψ2 = ω − 431 432 433
z
2
ω
( uz + g x ) ,
+ g x ) cosψ 1 − 2 ( u x + g x ) sinψ 1 ) ,
1 ( ( uz + g x ) sinψ 1 + 2 ( ux + g x ) cosψ 1 ) , B2ω
(u
y
+ g x ) sinψ 2 ,
1 ( u + g x ) cosψ 2 . ω B4 y
In case of small eccentricity (only the first order of e is held) ω can be replaced by
ω0 = µ p3 and the perturbation terms become
434
g x = eω02 ( 2sin ϑ ( B2 sinψ 1 + 2 B1 ) + 4cosϑ ( B2 cosψ 1 + B3 ) ) ,
435
g y = eω02 B4 cosψ 2 ,
436
g z = eω02 ( 2sin ϑ ( 2 B2 cosψ 1 + B3 ) + 4cosϑ ( B2 sinψ 1 + B1 ) ) .
437
The perturbation is proportional to eω02 . Since umax = 10−6 N these terms are several orders
438
smaller than the control. Indeed, during transient motion the maximum value of g is about
439
5 × 10−10 N , while during the steady state motion this value is about 5 × 10−8 N .
-
The control of the two-satellite formation flight is considered The cell-structured solar sail is used for the relative motion control The angular motion control is performed by the sail optical properties variation The control scheme is able to stabilize the required relative orbit
Declaration of interests ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:
S0094-5765(20)30086-2
DOI:
https://doi.org/10.1016/j.actaastro.2020.02.024
Reference:
AA 7890
To appear in:
Acta Astronautica
Received Date: 23 September 2019 Revised Date:
7 February 2020
Accepted Date: 12 February 2020
Please cite this article as: Y. Mashtakov, M. Ovchinnikov, T. Petrova, S. Tkachev, Two-satellite formation flying control by cell-structured solar sail, Acta Astronautica, https://doi.org/10.1016/ j.actaastro.2020.02.024. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier Ltd on behalf of IAA.
1 2
Two-satellite formation flying control by cell-structured solar sail 1
3
Yaroslav Mashtakov, Mikhail Ovchinnikov, Tatyana Petrova, Stepan Tkachev ,
4 5 6 7 8 9 10 11 12 13 14 15
Keldysh Institute of Applied Mathematics of RAS, Miusskaya Square, 4, 125047 Moscow, Russia Abstract The method of simultaneous relative motion and attitude control by a solar radiation pressure is proposed. The control goal is to stabilize the given closed relative orbit. The principle idea is to use special materials for solar sail that able to change their optical properties. The solar sail consists of a number of rectangular cells. Each of them can change its optical properties. The necessary control force is produced by varying the average reflectivity of solar sail, and the control torque is achieved by the appropriate pattern of black and white cells. Key words: formation flying, solar radiation pressure, variable radiation reflectivity, attitude control, relative motion control
16
Introduction
17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39
Equation Section (Next) Utilization of a group of satellites, for example, formation flying brings new opportunities in space missions. In addition, group of satellites is more reliable because even if one satellite fails, others can continue the mission operation. The main problem of formation flying utilization is deployment and maintenance of the particular group configuration. The simplest solution for this problem is to use thrusters installed on board of all or several satellites. On the other hand, thrusters require propellant, which can greatly affect the satellite lifetime or the payload mass. To overcome this problem environmental forces for formation flying motion control are used [1]. This approach can be applied relatively easily by installing a special high area-to-mass ratio equipment like a flat sail. There are two sorts of forces that can be used: aerodynamic drag [2–6] and solar radiation pressure (SRP) [7– 12]. The principal idea here is to use a difference in environmental forces acting on each satellite in formation. This difference usually appears when sail changes its attitude with respect to the sun direction, although the effective size variation is also considered in literature [13]. In general, solar sails can be used either for attitude control [14–18] or for orbit control [19,20]. In the paper the case when both attitude and relative motion are controlled via solar sail with variable optical properties is considered. Simultaneous orbit and attitude control using solar sail with no additional devices seems impossible. However, there are number of approaches that permit both orbit and attitude control using either sail shape changing [21] or reflectivity properties variation [13], [22–24]. The latter approach is of interest is the paper. There are several techniques to use variable optical properties. In [13] it is considered that sail consists of two part with different optical properties each and the size of these parts can change continuously. Usually, the more realistic model is used when there is only a set of cells with variable optical properties are place one the sail. This approach was used during well1
Corresponding author E-mail address: [email protected]
40 41 42 43 44 45 46 47 48 49 50 51
known mission IKAROS [25]. In [22] the problem of appropriate pattern determination for the required torque and force implementation is aroused in case of sail with cell structure and the table of possible patterns is created. However this approach can’t be scaled for the large amount of cells since it requires too much memory to save all patterns. Another approach is described in [24]. The pattern is built step-by-step: for achieving attitude control and then – orbit control. Thanks to the symmetry of the sail, it is divided into four quadrants. In each quadrant the square domain which will provide necessary in-plane torque is selected. These domains are the same ones and are located symmetrical with respect to the sail symmetry center. The signs of the inplane torque define the ‘on-off’ state of each domain and the values – the coordinates of the center of the domain and its square. Finally, for orbit control symmetrical cells are added to the pattern. In this paper, similar approach is used but it differs in details, control synthesis scheme and application.
52
1. Problem statement and reference frames
53 54 55 56
Deployment and maintenance of two satellites in the required relative orbit are considered. Each satellite is assumed to be equipped with a solar sail. The initial orbit of leadersatellite is circular. The follower-satellite is moving along the orbit which is close to the first one. Satellites move under the solar radiation pressure and are both disturbed by J 2 harmonic of
57 58 59
the Earth gravitational field. The following right-handed reference frames are used: – O1 XYZ is the inertial frame (IF) with the origin at the Earth center of mass, O1Z is
60 61 62 63 64 65 66 67
perpendicular to the equator plane, O1 X is directed to the vernal equinox; Оxyz is the orbital frame (OF), its origin is in the leader satellite’s center of mass, Оz is directed along its radius vector, Оy is perpendicular to the orbit plane; Оξηζ is the body-fixed frame (BF), its axes are the principal axes of inertia (it is also assumed that Оζ is perpendicular to the sail plain); Оxs y s z s is the solar frame (SF), Оzs is directed to the Sun, Оys is perpendicular to the ecliptic plane. Transition between IF and OF is performed by the following matrix
68 69 70
71
– – –
S = ( e1 e 2
e 3 ) , e3 =
r r×v , e2 = , e1 = e 2 × e 3 , r |r × v |
where r is the radius vector and v is the velocity of the leader satellite. Transition between IF and SF is determined by the matrix Ssun = ( l1
l2
cos Λ 0 l 3 ) , l 3 = sin Λ cos ε , l 2 = − sin ε , l1 = l 2 × l 3 . sin Λ sin ε cos ε
72
Here Λ is the ecliptic longitude, ε is the obliquity of the ecliptic.
73
Problem statement and reference frames
74
There are three types of motion equations that are used in this paper.
75
Orbital dynamics
76
78
Orbital dynamics is described by the following vector equation r & r&= − µE 3 + g , r where µ E is the Earth gravity constant and g is the result vector of the external disturbances. As
79
it was mentioned before, only the effects of J 2 and SRP force are taken into account. The first
80
has a form
77
81
3sin 2 i sin 2 u − 1 3J µ R f J 2 = 2 E4 − sin 2 i sin 2u . 2r − sin 2i sin u
82
Here J 2 = 1.082 × 10−3 , R⊕ is the mean Earth radius, i is the orbit inclination and u is the
83
argument of latitude. SRP force acting on the elementary area dS can be written as follows [26] Φ dFs = − 0 ( rs , n ) × c 2 × (1 − α ) rs + 2αµ ( rs , n ) n + α (1 − µ ) rs + n dS , 3
84
2 ⊕
85
where Φ 0 = 1357 W m is the solar flux constant, rs is the unit vector from the Sun to the
86
satellite, n is the solar sail normal (SSN), α is the reflection coefficient, µ is the specularity
87
coefficient. Further the case of µ = 1 is considered, i.e. there is no diffuse reflection, so
2
dFs = −
88 89 90 91
Φ0 (rs , n ) ( (1 − α ) rs + 2α ( rs , n ) n ) dS . c
Since the parameter α varies from point to point, the total SRP force becomes Φ Fs = − 0 ( rs , n ) S − ∫ α dS rs + 2 ( rs , n ) n ∫ α dS . c
((
Denote f =
∫ α dS ( 0 ≤ f S
≤ 1) and A = −
)
)
Φ0S , then c
Fs = A ( rs , n ) ( (1 − f ) rs + 2 f ( rs , n ) n ) .
92 93
These equations are written for both satellites and are used in numerical simulation.
94
Attitude dynamics
95
Angular dynamics is described in the BF by the Euler equations
96
&+ ω × Jω = M control + M grav , Jω
97
where J is the satellite inertia tensor, ω is the satellite angular velocity with respect to the IF,
98
M control is the control torque and M grav = 3
99 100
µE
(1)
r × Jr is the gravity gradient torque. r5 Attitude kinematics (BF w.r.t. IF) is defined by the quaternion ( λ0 , λ ) , λ02 + λ 2 = 1 and described by the corresponding equations
λ&0 = −0.5 ( λ , ω ) ,
101 102
λ&= 0.5 ( λ0ω + λ × ω) . These equations are used for the numerical simulation and control torque synthesis.
103 104
Relative motion dynamics
105 106
The control synthesis for orbital relative motion is based on the Hill-Clohessy-Wiltshire equations written in the curvilinear parameters ( a0θ r , a0ϕ r , ρ ) (see Fig. 1). These parameters are
107 108 109
more preferable than the regular ones since it is assumed that the leader satellite is moving along circular orbit while the relative orbit is small with respect to the size of the orbit. Motion equations in the OF can be written as follows [27]
110
d2 d a ϕ + 2ω ρ = 0, 2 ( 0 r) dt dt 2 d a θ + ω 2 ( a0θ r ) = 0, 2 ( 0 r) dt d2 d ρ − 2ω ( a0ϕ r ) − 3ω 2 ρ = 0 2 dt dt
111
where ω is the orbital angular velocity of the leader satellite, ρ = r2 − r1 , a0 = r1 . Index “1”
112
corresponds to the leader satellite and “2” to the follower.
113 114 115
116
Figure 1. Curvilinear parameters If control is taken into account then Eq.(2) can be rewritten in the following form
d2 d a ϕ + 2ω ρ = u x , 2 ( 0 r) dt dt 2 d a θ + ω 2 ( a0θ r ) = u y , 2 ( 0 r) dt d2 d ρ − 2ω ( a0ϕ r ) − 3ω 2 ρ = uz . 2 dt dt
(2)
117 118
Here u x , uy , u z are the components of the control vector u =
Fs ,2 − Fs ,1 . m
Solution of (2) is a0ϕ r = −3C1ωt + 2C2 cos ωt − 2C3 sin ωt + C4 ; a0θ r = C5 cos ωt + C6 sin ωt;
119
ρ = 2C1 + C2 sin ωt + C3 cos ωt. 120
One can introduce new variables B j , ψ j based on this solution a0ϕ r = 2 B2 cosψ 1 + B3 ,
ρ = B2 sinψ 1 + 2 B1 , d ( a0ϕ r ) = −2 B2ω sinψ 1 − 3B1ω, dt d ρ = B2ω cosψ 1 , dt a0θ r = B4 cosψ 2 ,
121
d ( a0θ r ) = − B4ω sinψ 2 . dt 122
Below the explicit definition of parameters B j and ψ j is given B1 = C1 , B2 = C22 + C32 , B3 = −3C1ωt + C4 ,
123
B4 = C52 + C62 , C3
ψ 1 = ωt + β1 , where sin β1 =
C22 + C32
ψ 2 = ωt + β 2 , where sin β 2 = − 124
C6 C52 + C62
C2 C22 + C32
, cos β 2 =
,
C5 C52 + C62
.
The equations corresponding to these variables are .
B1 =
1
ω
ux ,
.
B3 = −3B1ω − .
B2 = 125
, cos β1 =
1
ω
2
ω
uz ,
( uz cosψ 1 − 2ux sinψ 1 ) ,
1 ψ1 = ω − ( uz sinψ 1 + 2ux cosψ 1 ) , B2ω .
.
B4 = − .
1
ω
ψ2 = ω −
u y sinψ 2 , 1 u cosψ 2 . ω B4 y
(3)
126 127
Here B1 corresponds to the drift velocity of the follower satellite along axis Ox of the OF, B2 is the size of the in-plane ellipse, B3 is the shift of this ellipse along Ox . Finally, B4 is
128 129
the out-of-plane motion amplitude. The control purpose is to achieve the required Bi , i = 1, 2, 3, 4 parameter, i.e. to have the required relative orbit.
130
Control synthesis
131
The general control synthesis scheme is presented in Fig. 2.
132 133 134 135
Figure 2. Control synthesis scheme First of all, the ideal control that provides required relative motion is found. Then corresponding integral reflection coefficient f i and angles of SSN θi , ϕ i are determined.
136 137
Normal directions define the reference angular motion of each satellite. After that the control torque Mi , control (in Оξη plane) is calculated. Finally, Mi , control and f i determine the solar sail
138
optical pattern. Further in this section each step is discussed.
139
Relative motion control
140 141
Purpose of the control is to deploy the formation and maintain the required relative orbit. This orbit is defined by Bi ( i = 1, 2,3, 4 ). In the paper the relative orbit defined by
142
B1 = 0, B2 = B0 , B3 = 0, B4 = 0
143
is considered. This means that relative orbit is the ellipse with major and minor semiaxes 2B0
144
and B0 respectively lying in the orbital plane. Its center coincides with the origin of the OF. The
145
relative orbit stabilization is divided into two stages: firstly B1 = 0 and B3 = 0 are provided, then
146
B2 = B0 is achieved. The out-of-plane motion control process is separated, so B4 = 0 can be
147 148
guaranteed independently. On the first stage the following Lyapunov control function (LCF) is used
149
V=
1 2 1 2 B1 + B3 . 2 2
150
Its time derivative in accordance with (3) is . . 1 2 V&= B1 B1 + B3 B3 = B1ux + B3 −3B1ω − uz . ω ω
151 152
So the control that ensures global asymptotic stability of B1 = 0 and B3 = 0 is the following u x = −k1 B1 , k1 > 0,
153
uz =
1 ( −3B1ω 2 + k2ω B3 ) , k2 > 0 2
(4)
154
by the Barbashin-Krassovskiy theorem [28]. It should be noted when B1 is large the value of
155
required u z might exceed the accessible control effort umax (in the paper this value is 10−6 N / kg
156
). Hence, one should wait until B1 becomes small which is possible since the first equation of (3)
157 158
is independent. On the second stage the following LCF is used
V=
159 160 161 162
and its time derivative
1 1 V&= ( B1 − 2 ( B2 − B0 ) sinψ 1 ) ux + ( ( B2 − B0 ) cosψ 1 − 2 B3 ) uz − 3B1B3ω.
ω
ω
Then the control is u x = − k x ( B1 − 2 ( B2 − B0 ) sinψ 1 ) , k x > 0
163 164
1 2 1 1 2 B1 + ( B2 − B0 ) + B32 2 2 2
u z = − k z ( ( B2 − B0 ) cosψ 1 − 2 B3 ) , k z > 0
(5)
which leads to
165
2 2 1 1 V&= − ( B1 − 2 ( B2 − B0 ) sinψ 1 ) − ( ( B2 − B0 ) cosψ 1 − 2 B3 ) − 3B1B3ω.
166 167 168 169 170
The last term does not depend on control and its sign is undefined. However, its value after the first stage would be smaller than the sum of the first and second terms. This is why the first stage is required. In case the magnitude of the control force received from the Barbashin-Krassovskiy theorem is larger the maximum allowed control force is set to the umax . Practically, in both cases
171
(4) and (5) the following scheme was used
ω
ω
172
−u sign( a ), a > umax u = max −a, | a |≤ umax
173 174
where a can be either first or second expression in (5) or first expression in (4). Additionally, since B1 is the drift velocity, it could be used to control the convergence speed of B3 to zero. So
175
with the help of nonzero drift the shift might be reduced faster. If B3 is large, then instead of the
176
first expression of (4) the following control is used
177
ux = −k1 ( B1 − γ sign ( B3 ) )
178
where γ is a some positive parameter. The out-of-plane motion control is
u y = −k y B4 cosψ 2 , k y > 0 .
179 180
The control u x , uy , u z is an ideal one. The control is called ideal in terms of being
181 182
received analytically. It should be implemented through the solar sails rotation and integral reflectivity coefficients f i .
183
Relative motion control implementation
184 185
Let θ be the angle between the SSN and the Sun direction, ϕ is the angle between normal projection to the plane Оx s y s of the SF and axis Оxs . Then in the SF the SSN is sin θ cos ϕ n = sin θ sin ϕ . cos θ
186
187
Relative motion control force u in the SF is
uxs = 2 Af 2 cos2 θ2 sin θ 2 cos ϕ 2 − 2 Af1 cos2 θ1 sin θ1 cos ϕ1 , 188
u ys = 2 Af 2 cos2 θ2 sin θ2 sin ϕ 2 − 2 Af1 cos2 θ1 sin θ1 sin ϕ1 ,
(6)
uzs = A (1 − f 2 ) cos θ 2 − A (1 − f1 ) cos θ1 + Af 2 cos3 θ 2 − 2 Af1 cos3 θ1. 189
These expression are the non-linear equations with respect to f i , ϕ i θi . Moreover, there are
190 191
only three equations for six unknown variables. This freedom is used for the maximization of the possible values of the u components. First of all, as the SRP force is decreasing when θi tends
192
to 90 degrees, it is reasonable to consider the case of small θi , so (6) transforms to
uxs = 2 Af 2θ 2 cos ϕ 2 − 2 Af1θ1 cos ϕ1 , 193
u ys = 2 Af 2θ 2 sin ϕ 2 − 2 Af1θ1 sin ϕ1 ,
(7)
uzs = Af 2 − Af1. 194
From the last equation one can see, that f i define u zs completely. The maximum torque
195
is achieved when f = 0.5 , while for f = 0 and f = 1 the torque is zero. It is reasonable to
196
demand f i to be as close to 0.5 as possible, so min J = ( f1 − 0.5 ) + ( f 2 − 0.5 ) 2
2
f1 , f 2
subject to
197 f 2 − f1 = 0 < f min
198
u zs
A ≤ fi ≤ f max < 1
The solution of this problem in the inner domain is
u zs , 2A u f 2 = 0.5 + zs . 2A f1 = 0.5 −
199
(8)
200
It exists when
2 f min − 1 ≤
201
uzs ≤ 2 f max − 1 A
(9)
202
In the further discussion it is supposed that f min = 0.25 and f max = 0.75 . So (9) might be
203
rewritten as
−0.5 ≤
204 205
uzs ≤ 0.5 . A
When f i is known one can find ϕ i from the maximization (with θi fixed) of
206
max L = ( f 2θ 2 cos ϕ 2 − f1θ1 cos ϕ1 ) + ( f 2θ 2 sin ϕ 2 − f1θ1 sin ϕ1 ) .
207
It means that the result values of ϕ i give the largest domain of the possible control force. This
208
problem has two groups of solutions
2
2
ϕ1 , ϕ 2
209
ϕ1 = ϕ 2 , θ1θ 2 < 0 ,
210
ϕ1 = ϕ 2 + π , θ1θ2 > 0 .
211 212
The actual position of the SSN is the same for both solutions. Hence, further the case ϕ1 = ϕ 2 = ϕ is considered. The first and the second equations in (7) one can rewrite as follows
uxs , 2A u ( f 2θ2 − f1θ1 ) sin ϕ = ys . 2A
( f 2θ2 − f1θ1 ) cosϕ = 213
214
Then tgϕ =
215
f 2θ 2 − f1θ1 =
216 217
218 219
u ys u xs
,
u xs2 + u 2ys 2A
.
If uxs cos ϕ > 0 f 2θ 2 − f1θ1 =
u xs2 + u 2ys 2A
.
To find θi one can solve the following optimization problem min L = θ12 + θ 22 θ1 ,θ 2
subject to 220 f 2θ 2 − f1θ1 =
u xs2 + u ys2 2A
−θ max ≤ θ i ≤ θ max
(10)
221
The solution of this problem gives the minimum possible values of θi that satisfy constraints.
222
This keeps θi in the vicinity of zero.
223
The solution in the inner domain is
θ1 = − 224
θ2 =
uxs2 + u ys2
f1 , f + f 22 2 1
2A uxs2 + u 2ys
f2 . f + f 22 2 1
2A
225
Thus, once u is known, attitude of the SSN can be found.
226
Attitude control
227 228
The expression for the SRP torque for the square sail and given sun direction in the BF
rs = ( sin θ cos β
cos θ ) is T
Mξ P cos θ Φ0 −Q cos θ M η = − c cos θ − Q sin θ sin β + P sin θ cos β M ζ
229
230
sin θ sin β
where P=
231
a 2
a 2
∫ ∫ ηα dξ dη ,
Q=
a a − − 2 2
232 233 234 235
a 2
a 2
∫ ∫ ξα dξ dη. −
a a − 2 2
There are only two components that can be defined independently by the optical properties variation and SSN attitude. So, there can be only two independent control torques to be produced. The following LCF is considered
V=
( (
))
1 T 2 2 J ξ ωrel,1 + Jηωrel,2 + ka 1 − ( 0 0 1) , Bn . ( ) 2
(11)
236
Here J ξ , Jη are the in-plane moments of inertia, ωrel ,1 , ωrel ,1 are the corresponding relative
237
angular velocity components ( ωrel = ω − ωref ), ωref = n × n&, n is the required SSN attitude in the
238
IF, ka > 0 and B is the transition matrix between the IF and the BF. The goal of the control is to
239
guarantee asymptotic stability of the motion when the axis Oζ of the BF and n coincide.
240
Derivative of (11) is
241
T d V&= J ξ ωrel,1ω&rel,1 + Jηωrel,2ω&rel,2 − ka ( 0 0 1) , ( Bn ) . dt
242
As there is no need to control the third component of the angular velocity let us take relative
243
angular velocity vector as ωrel = (ωrel ,1 ωrel ,2
T
d ( Bn ) = −ωrel × Bn dt
244 245
0 ) . In this case
and (12) one can rewrite as follows
(12)
. T V&= ωTrel J ωrel + ka Bn × ( 0 0 1) .
246 247
To guarantee V&< 0 it is sufficient to satisfy the following equation
248
&rel + ka Bn × ( 0 0 1) = − kω ω rel . Jω
249 250 251
T
It guarantees by the following control
&ref − ka Bn × ( 0 0 1) . Mcontrol = −kωωrel − Mext + ω × Jω − Jω × Bωref + JBω T
(13)
Here only two first components of M control are taken. The last component is defined by
252
the other ones.
253
Solar sail pattern
254 255 256
The final step of the control synthesis is the determination of the solar sail pattern for each satellite. This pattern must give the required integral reflection coefficients, as well as control torque components M ξ and M η . In general, the reachability domain is the inner area of
257
rhombus LMNK which is presented in Fig. 3.
258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273
Figure 3. The torque reachability domain The components of control torque, which based on Lyapunov theory, are independent. However, the domain between rhombus and square the components of control torque could not be chosen independently that means that the less component of torque could not be increased without decreasing the larger one and vice versa. To overcome this, the control torque components are limited by the square. Each solar sail is considered to be divided into n × n cells. Each cell can either reflect radiation totally or absorb it. The problem is to find the fully reflective cells that give the required torque components and integral reflectivity coefficient. There are various approaches. Here one of them is presented, which requires neither complex calculations nor large amount of memory. Let initially each cell have the reflectivity coefficient α = 0 . The required torque will be produced by the rectangle with α = 1 . Its center is located at the edge of the square, which is described in the following way: its center coincides with sail center, and the edge length equals to a half of sail's edge length (Fig. 4). The SRP torque which is produced by this rectangle is
ξ Fs ,ξ η Fs ,ζ M ξ Ms = r × Fs = η × Fs ,η = −ξ Fs ,ζ = Mη , 0 F ξ F −ηF M s ,ζ η ξ ζ
274 275
where r = (ξ η 0) is the center position vector. Therefore, T
M ξ =− η, η Mξ 2
a ξ + η = , 4 2
sign (η ) = −sign ( M ξ ) .
276
277 278 279 280
2
Figure 4. Solar sail pattern
The area of the rectangle is determined by the following expressions
281
S rec =
282
Fs ,ζ =
cFs ,ζ Φ 0 cos2 θ
.
M ξ2 + Mη2
ξ 2 +η2
283 284
The rectangle with a given center and area is not unique, but it should fit cell pattern on the sail surface. First of all it should be noted that
285
a S rec 2
286 287
since the control parameters in (13) are chosen to guarantee control torque to be in the square in Fig. 3. There are several cases: 2
288 289
2
a 1. ≤ S rec which means that the area corresponds to, at least, one cell and the sides of the n rectangle are taken as the side of the square rounded to the nearest cell.
2
2
290
a a 2. ≥ S rec then the area is taken as S rec = , but the center of rectangle is moving to the n n
291
ξ point k
η k
T
0 . Here
k=
292
M ξ2 + Mη2
ξ 2 + η 2 Fζ
.
293 294
When k ≥ n the distance from the center of sail to the center of the rectangle is less than the half-size of the cell. In this case Sпр = 0 and no torque is produced.
295
Once the cell pattern for torque implementation is determined, f should be realized. The
296 297
number of cells with α = 1 is determined from expression N req = f ⋅ n 2 .
298 299 300
However, there are already N cells with α = 1 which realize the control torque. In general, there are several options: 1. N = N req which means that there is no need in cell pattern modifications and required f is
301 302
achieved. 2. N < N req then the symmetric pairs with zero torque are added. If N − N req is an odd number
303 304
then one cell is left (this is an error in force realization). 3. N > N req . Torque producing cells give pattern with f ≤ 0.25 since no more than quarter of
305
sail has α = 1 . On the other hand in (9) minimum force producing f min = 0.25 . So this situation
306 307 308 309
is avoided. Thus, the control synthesis is finished. The presented scheme allows determining the solar sail pattern which gives the required control force and torque under given constraints.
310
NUMERICAL EXAMPLE
311 312 313 314
Numerical simulation is carried out using the parameters presented in Table 1. Integration method is Runge-Kutta with constant step of one second. The control force and torque are calculated in each integration step. After that the control is substituted to the orbital and angular models for each satellites. And, finally, the relative motion parameters ( Bi ) are calculated. These
315 316
parameters are presented in Figures 5-7. Table 1. Numerical simulation parameters Parameter Value Leader orbit radii Rorb = 9060 km Inclination 110 deg Initial relative orbit rrel = ( 200 100 50 ) m, Mass of each satellite Solar sail size
v rel = ( 0.05 0.5 1) m/s; m = 10 kg 5×5 m
Inertia tensor
J = diag ( 2.1 2.1 3.8 ) kg ⋅ m 2
Initial angular velocity
ω1 = ( 0.005 0.003 0.001) s−1 , ω2 = ( 0.001 0.003 0.002 ) s −1;
Attitude of the SSN
θ1 = 10 deg, θ 2 = 10 deg;
Control law parameters
k1 = 1 s −2 , k2 = 2 s−2 , k x = 0.5 s −2 , k y = 0.1 s−2 , k z = 0.05 s −2 , kω = 0.02 N ⋅ m ⋅ s, ka = 10−4 N ⋅ m
Maximum allowed control force
umax = 10−6 N
Maximum allowed control torque
M = 3 ⋅ 10−5 N ⋅ m
317 318 319
Results of the numerical simulation are presented in Figs 5-12. Figs 5-8 show the evolution of the parameters Bi which determine the relative orbit (drift, shift, orbit size and
320 321 322 323
amplitude of the out-of-plane motion). In Fig. 8 the relative orbit is shown. The blue lines here show relative in-plane motion over the whole time interval. The red ones illustrate resulting steady state relative in-plane motion. The numerical simulation is carried out for 4x4 (left column) and 50x50 (right column) cells sail.
Figure 5a. Drift and shift evolution for 4x4 cells sail
Figure 5b. Drift and shift evolution for 50x50 cells sail
324
Figure 6a. Relative orbit size evolution for 4x4 cells sail
Figure 6b. Relative orbit size evolution for 50x50 cells sail
Figure 7a. Out-of-plane motion amplitude for 4x4 cells sail
Figure 7b. Out-of-plane motion amplitude for 50x50 cells sail
Figure 8a. Relative in-plane motion evolution for 4x4 cells sail
Figure 8b. Relative in-plane motion evolution for 50x50 cells sail
Figure 9a. Angle between SSN and sun direction evolution for 4x4 cells sail
Figure 9b. Angle between SSN and sun direction evolution for 50x50 cells sail
325 326 327
The control successfully solve the task and formation reach the required relative orbit (Figs 5-9). It also can be seen in Fig. 5 that B1 ≈ −20 m during first : 15 revolutions. This
328
allows reaching required B3 value much faster. The Figures 9 show that θ is bounded in the
329 330
desired region. Figures 10a and 10b contain the switching rate of the solar cells for the 4x4 cell sails for both satellites.
Figure 10a. Solar sail cells switching rate (the black one is the highest, the white one is the lowest) for the chief satellite
Figure 10a. Solar sail cells switching rate (the black one is the highest, the white one is the lowest) for the deputy satellite
331 332 333 334 335 336 337
From Figure 10 one can see that there are four cells which almost all the time have constant reflectivity. So the ordinary solar cell material can replace these cells. The maximum usage is approximately 40% (black cell in Firuges). Which means it switches every 0.4 control iteration. The simulation results for 2x2 cells sail are shown in Fig. 11. This configuration fails to solve the stabilization problem. The similar situation is in case of 3x3 cells sail. So, it seems that 4x4 is the minimal configuration when the proposed control scheme solves the problem.
338 339
Figure 11. Drift and shift for 2x2 cells sail
340 341 342 343 344
The main problem here is that the 2x2 configuration does not allow producing the different torque components with the same sign and torque value is of no importance (only sign have an effect). The 3x3 meets difficulties on the stage when symmetric pairs with zero torque must be added. The 50x50 cells configuration provide in several times better accuracy but this is obviously more complex from the practical point of view.
345
Conclusion
346 347 348
In the paper the scheme of the two satellites formation flying control using the solar sail is proposed. It allows deploying and maintaining the required relative motion using the solar sail with variable optical properties only. In continuous case this scheme guaranties asymptotic
349 350 351 352 353 354 355 356 357 358
stability for both angular and relative motion. When the sails consist of discrete cells with variable optical properties the minimum number of cells appears to be sixteen i.e. four along each side. Note that the proposed control scheme is suitable for the time intervals when the Clohessy-Wiltshire approximation is valid. Under long-term SRP effects the eccentricity and other orbital parameters of the chief satellite evolve. For example, for the case considered in the paper the maximum eccentricity value is about 10−3 which still can be considered as perturbation (see Appendix). The long-term modification of the proposed scheme can be the subject of the future studies.
359
Acknowledgments
360 361 362
The work is supported by the RFBR Grant No. 18-31-20014. Authors are grateful to Prof. Anna Guerman for productive discussions.
363
References
364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393
1. 2.
3.
4.
5. 6. 7. 8. 9. 10. 11.
12. 13.
Kumar K.D. et al. Maintenance of satellite formations using environmental forces // Acta Astronaut. Elsevier, 2014. Vol. 102. P. 341–354. Leonard C.L., Hollister W.M., Bergmann E. V. Orbital formation keeping with differential drag // J. Guid. Control. Dyn. American Institute of Aeronautics and Astronautics, 1989. Vol. 12, № 1. P. 108–113. Kumar B.S., Ng A. A bang-bang control approach to maneuver spacecraft in a formation with differential drag // Proc. AIAA Guid. Navig. Control Conf. Exhib. 2008. № August. P. 1–11. Zeng G., Ru M., Yao R. Relative orbit estimation and formation keeping control of satellite formations in low Earth orbits // Acta Astronaut. Elsevier, 2012. Vol. 76. P. 164– 175. Varma S., Kumar K.D. Multiple Satellite Formation Flying Using Differential Aerodynamic Drag // J. Spacecr. Rockets. 2012. Vol. 49, № 2. P. 325–336. Omar S.R., Wersinger J.M. Satellite Formation Control using Differential Drag // 53rd AIAA Aerosp. Sci. Meet. 2015. № January. P. 1–11. Gong S., Yunfeng G., Li J. Solar sail formation flying on an inclined Earth orbit // Acta Astronaut. Elsevier, 2011. Vol. 68, № 1–2. P. 226–239. Shahid K., Kumar K.D. Multiple spacecraft formation reconfiguration using solar radiation pressure // Acta Astronaut. Elsevier, 2014. Vol. 103. P. 269–281. Shahid K., Kumar K.D. Formation Control at the Sun-Earth L(2) Libration Point Using Solar Radiation Pressure // J. Spacecr. Rockets. 2010. Vol. 47, № 4. P. 614–626. Gong S., Baoyin H., Li J. Solar Sail Formation Flying Around Displaced Solar Orbits // J. Guid. Control. Dyn. 2007. Vol. 30, № 4. P. 1148–1152. Williams T., Wang Z.-S. Solar Radiation Pressure and Formation-Keeping in Highly Elliptical Orbits // AIAA/AAS Astrodynamics Specialist Conference and Exhibit. 2002. № August. P. 1–14. Parsay K., Schaub H. Preliminary Results on Optimal Establishment of Solar Sail Formations // J. Astronaut. Sci. 2019. Vol. 66, № 1. P. 32–45. Mu J., Gong S., Li J. Reflectivity-controlled solar sail formation flying for magnetosphere mission // Aerosp. Sci. Technol. Elsevier Masson SAS, 2013. Vol. 30, № 1. P. 339–348.
394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430
14. 15. 16.
17. 18. 19. 20. 21. 22.
23.
24.
25. 26. 27. 28.
Wie B. Solar Sail Attitude Control and Dynamics, Part 1 // J. Guid. Control. Dyn. 2004. Vol. 27, № 4. P. 526–535. Wie B., Murphy D. Solar-sail attitude control design for a sail flight validation mission // J. Spacecr. Rockets. 2007. Vol. 44, № 4. P. 809–821. Wie B. et al. Robust attitude control systems design for solar sails, part 2: MicroPPTbased secondary ACS // Collection of Technical Papers - AIAA Guidance, Navigation, and Control Conference. 2004. Vol. 2, № August. P. 1325–1340. Bolle A., Circi C. Solar sail attitude control through in-plane moving masses // Proc. Inst. Mech. Eng. Part G J. Aerosp. Eng. 2008. Vol. 222, № 1. P. 81–94. Scholz C. et al. Performance analysis of an attitude control system for solar sails using sliding masses // Adv. Sp. Res. COSPAR, 2011. Vol. 48, № 11. P. 1822–1835. Bookless J., McInnes C. Dynamics and control of displaced periodic orbits using solar-sail propulsion // J. Guid. Control. Dyn. 2006. Vol. 29, № 3. P. 527–537. McInnes C.R. Passive control of displaced solar sail orbits // J. Guid. Control. Dyn. 1998. Vol. 21, № 6. P. 975–982. Wu R. et al. Heliogyro Solar Sail with Self-Regulated Centrifugal Deployment Enabled by an Origami-Inspired Morphing Reflector // Acta Astronaut. IAA, 2018. Borggrafe A. et al. Attitude Control of Large Gossamer Spacecraft using Surface Reflectivity Modulation // Proceedings of the 65th International Astronautical Congress. 2014. P. 8. Funase R. et al. On-orbit verification of fuel-free attitude control system for spinning solar sail utilizing solar radiation pressure // Adv. Sp. Res. COSPAR, 2011. Vol. 48, № 11. P. 1740–1746. Biggs J.D., Negri A. Orbit-Attitude Control in a Circular Restricted Three-Body Problem Using Distributed Reflectivity Devices // J. Guid. Control. Dyn. 2019. Vol. 42, № 12. P. 2712–2721. Tsuda Y. et al. Flight status of IKAROS deep space solar sail demonstrator // Acta Astronaut. Pergamon, 2011. Vol. 69, № 9–10. P. 833–840. Fliegel H.F., Gallini T.E., Swift E.R. Global Positioning System Radiation Force Model for geodetic applications // J. Geophys. Res. 1992. Vol. 97, № B1. P. 559. Alfriend K. et al. Spacecraft Formation Flying: Dynamics, Control and Navigation. 1 Edition. Butterworth-Heinemann, 2010. 402 p. David R. Merkin. Introduction to the Theory of Stability / ed. Afagh F.F., Smirnov A.L. Springer-Verlag New York, 1997. 320 p.
Appendix The perturbed equations (3) are
.
B1 =
1
ω
( ux + g x ) ,
.
B3 = −3B1ω − .
B2 =
1
ω
((u
.
ψ1 = ω − .
B4 = − .
1
ω
ψ2 = ω − 431 432 433
z
2
ω
( uz + g x ) ,
+ g x ) cosψ 1 − 2 ( u x + g x ) sinψ 1 ) ,
1 ( ( uz + g x ) sinψ 1 + 2 ( ux + g x ) cosψ 1 ) , B2ω
(u
y
+ g x ) sinψ 2 ,
1 ( u + g x ) cosψ 2 . ω B4 y
In case of small eccentricity (only the first order of e is held) ω can be replaced by
ω0 = µ p3 and the perturbation terms become
434
g x = eω02 ( 2sin ϑ ( B2 sinψ 1 + 2 B1 ) + 4cosϑ ( B2 cosψ 1 + B3 ) ) ,
435
g y = eω02 B4 cosψ 2 ,
436
g z = eω02 ( 2sin ϑ ( 2 B2 cosψ 1 + B3 ) + 4cosϑ ( B2 sinψ 1 + B1 ) ) .
437
The perturbation is proportional to eω02 . Since umax = 10−6 N these terms are several orders
438
smaller than the control. Indeed, during transient motion the maximum value of g is about
439
5 × 10−10 N , while during the steady state motion this value is about 5 × 10−8 N .
-
The control of the two-satellite formation flight is considered The cell-structured solar sail is used for the relative motion control The angular motion control is performed by the sail optical properties variation The control scheme is able to stabilize the required relative orbit
Declaration of interests ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: