# Dynamics and Control of Tethered Satellite System in Elliptical Orbits under Resonances

Authors: Zhaojun Pang; Zhonghua Du; Chun Cheng; Qingtao Wang
Journal: International Journal of Aerospace Engineering
 (2020)

---

## Abstract

This paper studies resonance motions of a tethered satellite system (TSS) in elliptical orbits. A perturbation analysis is carried out to obtain all possible resonance types and corresponding parameter relations, including internal resonances and parametrically excited resonances. Besides, a resonance parametric domain is given to provide a reference for the parameter design of the system. The bifurcation behaviors of the system under resonances are studied numerically. The results show that resonant cases more easily enter chaotic motion than nonresonant cases. The extended time-delay autosynchronization (ETDAS) method is applied to stabilize the chaotic motion to a periodic one. Stability analysis shows that the stable domains become smaller in resonance cases than in the nonresonance case. Finally, it is shown that the large amplitudes of periodic solutions under resonances are the main reason why the system is difficult to control.


---

## Body

## 1. Introduction

Tethered satellite system is a promising new type of spacecraft [1–3]. It has great potential in applications such as space debris capture [4], space elevator [5], and orbital transfer [6]. For different applications, there will be different combinations of system parameters, which will increase the possibility of resonance. Especially for the system in elliptical orbit, the resonance of the system becomes more complex due to the periodic excitation. Moreover, the dynamic characteristics of the resonance system are very different from those of the nonresonance system, so it is necessary to study all kinds of resonance motion of the tethered satellite system.The researches on the dynamic behavior of the tethered satellite system under nonresonance mainly focus on the steady-state solutions and stabilities. Takeichi et al. [7] studied the periodic solution of the librational motion of a tethered system in elliptic orbit via the Lindstedt perturbation method. Their analysis showed that the periodic solution is the minimum energy solution, and from the mechanical point of view, the periodic solution in an elliptic orbit has the same significance as the equilibrium state in a circular orbit. Sidorenko and Celletti [8] obtained the family of planar periodic motions for the “spring-mass” model of tethered satellite systems and studied the bifurcations and the stability of these periodic motions concerning in-plane and out-of-plane perturbations. Nakanishi et al. [9] investigated the in-plane periodic solutions of a dumbbell satellite system in elliptic orbits by bifurcations with respect to the orbital eccentricity. The periodic solutions are projected on the van der Pol plane to observe the directional change points of the trajectories which can become a useful tool to develop new control schemes. Peláez and Lara [10] employed the Poincaré method of continuation of periodic orbits to obtain other periodic solutions in electrodynamic tethers on inclined orbits that cannot be obtained by using asymptotic techniques. Burov et al. [11] studied the families of periodic motions analytically in a cabin-elevator system, and the condition for the existence of periodic solutions was determined.Usually, the periodic motions of the system are the desired target motions, but most of these periodic motions are unstable; many control methods are proposed to stabilize the unstable periodic motions. Fujii et al. [12] used a time-delay feedback control (DFC) to control chaotic librational motion of a tethered satellite in an elliptic orbit; the result shows that the DFC scheme is particularly effective in controlling the chaotic motion of the system. The same control method is employed by Peláez and Lorenzini [13] to control an electrodynamic tethered system in inclined orbit. But this control law does not stabilize the unstable periodic orbit for reasonable values of the control parameters. Thus, for the librational control of the electrodynamic tethered system, Iñarrea and Peláez [14] adopted the extended time-delay autosynchronization (ETDAS), while Williams [15] used the predictive control with time-delayed feedback scheme; all achieved good control performances. Additional, Kojima et al. [16] applied a model-following, decoupling-control method combining with the delayed feedback control method in a new approach to control the librational motion of the tethered satellite system in elliptic orbits. The numerical simulations showed that the control method has good performance such as a short settling time for convergence and small control forces.However, there are few studies on the nonlinear resonances of a tethered satellite system, and the limited researches mainly focus on the system in circular orbit [17, 18], and the more complex elliptical orbit system which may have combination resonance needs to be studied. This paper studies the nonlinear resonance of TSS in elliptical orbits. The mathematical model of the in-plane tethered satellite system in elliptical orbits in which the main satellite is treated as a rigid body is derived in Section 2. Then, perturbation analysis of the system is executed to get the parametric relations of every kind of resonance in Section 3. In Section 4, the bifurcation analyses of the TSS under resonance are studied numerically. Section 5 implements the ETDAS technique to stabilize the unstable motion under resonance to a periodic motion. Periodic solutions under different resonances are compared in Section 6. The conclusions are given in Section 7.


## 2. Mathematical Model

Consider an in-plane tethered satellite system as shown in Figure1, where the main satellite is treated as a rigid body of mass M and the subsatellite is envisioned to be a point of mass m. The subsatellite is connected with an inelastic massless tether l at a joint point of distance ρ to the mass center C of the main satellite. The mass of the main satellite is assumed to be much greater than the mass of the subsatellite, such that the center of mass of the system can be assumed to coincide with the main satellite that moves in an unperturbed Kepler elliptical orbit of semimajor axis a and eccentricity ratio e.Figure 1
Schematic view of the in-plane tethered satellite system.Two right-oriented reference frames are used to describe the system motion. The first is the earth-centered inertial frame denoted byO‐XYZ, the origin of which is located at the center of the Earth, where i, j, and k are unit vectors in the directions of OX, OY, and OZ, respectively. The second is that the body frame established with the directions of the axes coincided with the principal axes of the main satellite. Ix, Iy, and Iz are the principal moments of inertia of the main satellite expressed in the body frame. In this paper, only a specific relation Iy=Ix is considered.According to Figure1, the position vectors of the main satellite and the subsatellite in the inertial frame are, respectively,
(1)R1=Rcosvi+Rsinvj,(2)R2=R1−ρcosν+α+lcosν+θi−ρsinν+α+lsinν+θj,where ν is the true anomaly and R is a distance between the main satellite and the center of the earth. Because the main satellite moves in an unperturbed Kepler elliptical orbit, it is easy to get R=a1−e2/1+ecosν.The potential energy of the system is in the following form [19]:
(3)V=−μemR2−μeMR1−12μeR13Ix+Iy+Iz+32μeR13Ixcos2α+Iysin2α+Iz,where μe is the gravitational constant.The kinetic energy of the system is(4)T=12Izdαdt+dνdt2+12mdR2dt2+12MdR1dt2.The equation of motion of the system can be derived by Lagrange’s equations,(5)ddt∂L∂q̇i−∂L∂qi=0,i=1,2,where the Lagrangian function L=T−V and the generalized coordinates are q1,q2=θ,α.Considering that the tether length is much smaller than the orbit radius, the joint point of distanceρ is smaller than the tether length, and only a small eccentricity ratio e is studied here, the equations of motion of the system are obtained as follows:
(6)θ¨=−3sinθcosθ−β1α̇2+2α̇sinθ−α+3cosαsinθ+e21+θ̇sinν+3sinθcosθcosν,α¨=β2θ̇2+2θ̇+3cos2θsinθ−α+2e1+α̇sinν,where β1=ρ/l, β2=Mρl/Iz, and the overdot denotes the derivative with respect to the true anomaly ν.


## 3. Perturbation Analysis

The resonance motions of the system are analyzed by a multiscale method. It is assumed that the eccentricity is small, and a small scaling factorε is introduced to rescale the eccentricity as follows:
(7)e→εe.The first-order approximate analytic solutions of Equation (6) can be expressed as
(8)θ=εθ1T0,T1+ε2θ2T0,T1,α=εα1T0,T1+ε2α2T0,T1,where Ti=εiτ, i=0,1, represents different time scales.By substituting Equation (8) into Equation (6) and equating the same power of ε, one can get the following:
(i)
Order (ε1)(9)D02θ1+3β1+1θ1=2esinν,D02α1−3β2θ1+3β2α1=2esinν(ii)
Order (ε2)(10)D02θ2+2β1θ1−α1D0α1+3β1+1θ2=−2D1D0θ1+2eD0θ1sinν+3eθ1cosν,D02α2−2β2θ1−α1D0θ1+3β2α2−θ2=−2D1D0α1+2eD0α1sinνwhere Dr=∂/∂Trr=0,1 represents the partial differential operators.The first-order approximation solutions of Equation (9) take the complex form
(11)θ1=A1T1expjω1T0+jB1expjT0+cc,α1=Γ1A1T1expjω1T0+A2T1expjω2T0+jB1expjT0+cc,where j is an imaginary unit and the coefficients A1T1 and A2T1 are functions of time T1. And cc denotes the complex conjugate of the preceding term. The coefficients B1, B2, and Γ1 are time-independent, their expressions are as follows:
(12)B1=−eω12−1,B2=−eω12+ω22−1ω12−1ω22−1,Γ1=β2β2−β1−1.Andjω1 and jω2 are the two mutually different pure imaginary roots of the following characteristic equation of λ:
(13)det31+β1+λ20−3β23β2+λ2=λ4+3β1+β2+1λ2+9β1β2+9β2=0.So it is easy to getω12=3+3β1 and ω22=3β2.Substituting Equation (11) into Equation (10) yields
(14)D02θ2+3β1+1θ2=2jD1A1ω1expjω1T0+2jβ1Γ1A12ω11−Γ1exp2jω1T0−2jβ1A22ω2exp2jω2T0+R1exp2jT0+R2expjω2T0+jT0+R3expjω2T0−jT0+R4expjω1T0+jT0+R5expjω1T0−jT0+R6expjω1T0+jω2T0+R7expjω1T0−jω2T0+cc(15)D02α2−3β2θ2+3β2α2=2jω2D1A2expjω2T0+2jω1Γ1D1A1expjω1T0+2jβ2A12ω1Γ1−1exp2jω1T0+P1exp2jT0+P2expjω2T0+jT0+P3expjω2T0−jT0+P4expjω1T0+jT0+P5expjω1T0−jT0+P6expjω1T0+jω2T0+P7expjω1T0−jω2T0+ccwhere the coefficients Ri and Pii=1,2,⋯7 are shown in the appendix.It is easy to see that the exponential terms in the right-hand side of Equations (14) and (15) are ±jω1T0, ±jω2T0, ±2jω1T0, ±2jω2T0, ±2jT0, ±jω1+1T0, ±jω1−1T0, ±jω2+1T0, ±jω2−1T0, ±jω1+ω2T0, and ±jω1−ω2T0. So, it can be seen that internal resonances and parametrically excited resonances may occur in the system. All possible resonance types of the system and their corresponding parameter relationships are given in Table 1. In the table, cases 1–3 are internal resonances and cases 4–8 are parametrically excited resonances.Table 1
Resonance types and their parametric relationships.
CaseTypes of resonanceFrequency relationshipParametric relationshipsResonance terms1Primary resonance (PR)ω1=ω2β2=1+β1jω1T0, jω2T02Internal resonance 1 (IR1)ω2=2ω1β2=41+β12jω1T0, jω2−ω1T03Internal resonance 2 (IR2)ω1=2ω2β2=1+β1/42jω2T0, jω1−ω2T04Combination resonance 1 (CR1)ω1−ω2=1β2=3+3β1−12/3jω1−1T0, jω2+1T05Combination resonance 2 (CR2)ω2−ω1=1β2=3+3β1+12/3jω2−1T0, jω1+1T06Combination resonance 3 (CR3)ω1+ω2=1β2=1−3+3β12/3j1−ω1T0, j1−ω2T07Superharmonic resonance 1 (SR1)ω2=2β2=4/32jT08Superharmonic resonance 2 (SR2)ω1=2β1=1/32jT0As can be seen from Table1, when a certain parametric relationship is satisfied, the system may have the internal resonance or parametrically excited resonance. β1 is always positive, so one can get ω1>3. Therefore, case 6 does not exist, that is, combination resonance 3 (CR3) will never occur in the system.From Table1, in cases 1–5, the relationship between parameters β1 and β2 is positively correlated. Figure 2 shows the relation curves between parameters β1 and β2 under every resonance type. Considering that most applications of TSS have ρ<l, the range of β1 in the figure is 0,1. It can be seen from the figure that, except for SR2, that is β1=1/3; the value range of parameter β2 is always in [0,8]. That means, in most applications of space tethered systems (0<β1<1), including elliptic orbits, the resonance region can be avoided when β2>8. It can provide a reference for the parameter design of TSS.Figure 2
The relation curves betweenβ1 and β2 under every resonance type.


## 4. Bifurcation Analysis

### 4.1. Bifurcation Diagram Changing with Orbital Eccentricitye

In Section3, the parametric relations of every kind of resonance are obtained by the analytic method. In this section, the bifurcation behaviors of the system under nonresonance (NR) are compared with those under resonances by numerical simulations. Consider that β1 is usually a small quantity; β1=0.001 is chosen in the simulations. So ω1≈3 can be obtained. Table 2 shows the values of parameter β2 corresponding to every resonance type. It can be seen that when β1=0.001, six kinds of resonances can occur in the system by choosing the value of parameter β2, where SR2 and CR3 will not occur. And notice that these six kinds of resonance will not occur at the same time.Case 1.
NR
Firstly, the bifurcation behavior of the system under NR with orbital eccentricitye is investigated. The parameters β1=0.001 and β2=10 are selected. No resonance occurs in the system. Figure 3 shows the bifurcation diagrams of the system changing with orbital eccentricity e. As can be seen from the figure, the bifurcation diagrams of θ and α are similar. In e∈0,0.284, they experience bounded motions of periods, period-doubling, and quasiperiods. When e≥0.284, the motions of θ and α are all chaotic. As we all know, for a dumbbell model, the pitch motion θ enters a tumbling state around the eccentricity of 0.313. So when the mother satellite is treated as a rigid body, this system enters a chaotic motion easier than a dumbbell model. But the critical orbital eccentricity does not decrease much. The coupling of the motions of θ and α makes it easier for the system to enter chaotic motions. But the coupling is not very strong in the case of nonresonance.Case 2.
CR1
Next, the bifurcation diagrams ofθ and α are observed under combination resonance. A set of parameters is selected as β1=0.001 and β2=0.18; one can get ω1≈3 and ω2≈3−1, so ω1−ω2≈1. The bifurcation diagrams under CR1 are shown in Figure 4. As can be seen from the figure, the bifurcation diagrams of θ and α show significant differences. The bifurcation diagram of pitch angle θ is similar to that in the case of nonresonance, and the critical value of e is 0.279. However, the bifurcation diagram of α is quite different from that of the nonresonant case. The motion of α enters into chaotic motion at a very small value of e, and the critical value of e is just 0.08. And it is worth noting that the amplitudes of the periodic, doubling-periodic, and quasiperiodic motion of α also increase significantly compared with those of the nonresonance. This type of parametrically excited resonance has a great impact on the motion of α, while it has a small impact on the motion of θ.Case 3.
PR
In this case, the bifurcation diagrams under primary resonance are obtained. Here,β1=0.001 and β2=1 are selected; one can get ω1≈3 and ω2=3, so ω1≈ω2. Figure 5 shows the bifurcation diagrams under the primary resonance. As can be seen from the figure, the bifurcation diagrams of θ and α are similar to those in the case of combination resonance. The critical values of e are 0.282 and 0.133 in the bifurcation diagrams of θ and α, respectively.Table 2
The values of parameterβ2 corresponding to every resonance type at β1=0.001.
The value ofβ2The value ofω2Frequency relationshipTypes of resonanceβ2≈0.18ω2≈3−1ω1−ω2≈1CR1β2≈0.25ω2≈3/2ω1≈2ω2IR2β2≈1ω2≈3ω1≈ω2PRβ2≈1.33ω2≈2ω2≈2SR1β2≈2.49ω2≈1+3ω2−ω1=1CR2β2≈4ω2≈23ω2=2ω1IR1Figure 3
Bifurcation diagrams changing with orbital eccentricitye under nonresonance.Figure 4
Bifurcation diagrams changing with orbital eccentricitye under combination resonance.Figure 5
Bifurcation diagrams changing with orbital eccentricitye under primary resonance.After simulating the bifurcation diagrams of the other four types of resonances atβ1=0.001, the results similar to the above two kinds of resonances can be obtained. That is, the resonances have a great impact on the motion of α, while they have a small impact on the motion of θ. The main effects of resonance on the motion of α are that the amplitudes of its periodic, doubling-periodic, and quasiperiodic motions increase compared with those of nonresonance. The critical value of orbital eccentricity that causes it to move into a chaotic motion is much smaller than that of nonresonance. The main reason for this situation can be obtained by observing the system’s first-order approximate differential Equation (9). In the equation, the motion of θ decouples with α. So internal resonances and parametrically excited resonances related toαhave very little effect on it. The first-order approximation differential equation of α in Equation (9) is coupled with θ, so the motion of α is greatly affected by these kinds of resonances.Case 4.
SR2
By observing Equation (9), it can be seen that the motion of θ is affected by the parametric excitation term, so it can be suspected that the SR2 may have great effects on the motion of θ. A set of parameters is selected as β1=0.33 and β2=10 in bifurcation simulations; one can get ω1≈2 and ω2=30, so the system is under SR2. The bifurcation diagrams under superharmonic resonance are shown in Figure 6. As can be seen from the figure, the bifurcation diagrams of θ and α are similar. The bifurcation diagram of θ is significantly different from those of nonresonance and other types of resonances. The main manifestation is that the amplitudes of the periodic, doubling-periodic, and quasiperiodic motions of θ are significantly larger than those in other cases when the eccentricity is small. Moreover, the critical value of e that makes the motion of θ be chaotic is also reduced to 0.204. As the motion of α is coupled with that of θ, their bifurcation diagrams are similar. The critical value of e is also 0.204 in the bifurcation diagram of α.Figure 6
Bifurcation diagrams changing with orbital eccentricitye under superharmonic resonance.


### 4.2. Bifurcation Diagram Changing with Parameterβ2

The bifurcation behaviors of the system changing with parameterβ2 are investigated next. Given that β1=0.001 and e=0.15, Figure 7 shows the bifurcation diagrams changing with parameter β2. There are six vertical solid lines in the figure; they are the values of parameter β2 under various kinds of resonance when β1=0.001 as shown in Table 2. It can be seen in the figure that the change of parameter β2 has little effects on the motion of θ, that is, the six types of resonance listed in Table 2 have little effect on the motion ofθ. The motion of α varies significantly with the change of parameter β2, and the regions with a large amplitude of α correspond to the positions of resonance regions given in Table 2. So the resonance regions obtained by analytical analysis are in good consistency with the results obtained by numerical simulations. This is also consistent with the above results of bifurcation behaviors with the change of parameter e.Figure 7
Bifurcation diagrams changing with parameterβ2 at β1=0.001.In Figure7, it can also be seen that there are two areas with an obvious sudden increase of amplitude in the bifurcation diagram of α, which are not corresponding to the above resonance regions. The two areas are located in the center of the region β2=3 and β2≈4.64, respectively. They correspond to superharmonic resonance which is ω2=3 and combination resonance which is ω2−ω1=2. The two resonances are ignored in the theoretical analysis because the dynamic equations are truncated by order 2 in the previous perturbation analysis.The above analysis shows that SR2 has a greater effect on the motion ofθ. Here β1=1/3 is selected, so ω1=2, and SR2 occurs in the system. Table 3 shows the values of parameter β2 corresponding to every resonance type at β1=1/3. It can be seen that when β2 takes the values in Table 3, multiple resonances occur simultaneously in the system. Figure 8 shows the bifurcation diagrams changing with parameter β2 at β1=1/3. It can be seen that the bifurcation diagrams of θ and α are similar, and both of them have a sharp increase in amplitude in the resonance regions. It also confirms the above results of bifurcation behaviors with the change of parameter e.Table 3
The values of parameterβ2 corresponding to every resonance type at β1=1/3.
The value ofβ2The value ofω2Frequency relationshipTypes of resonanceβ2=1/3ω2=1ω1−ω2=1, ω1=2ω2CR1, IR2β2=4/3ω2=2ω1=ω2, ω2=2PR, SR1β2=3ω2=3ω2−ω1=1CR2β2=16/3ω2=4ω2=2ω1IR1Figure 8
Bifurcation diagrams changing with parameterβ2 at β1=1/3.From the above, since the first-order approximate equation ofθ is not coupled with α, only SR2 has a great effect on its motion. However, the amplitude of the motion of α is significantly increased under various kinds of resonances and can even be chaotic. In short, resonances have a huge effect on the motion of the system.


## 4.1. Bifurcation Diagram Changing with Orbital Eccentricitye

In Section3, the parametric relations of every kind of resonance are obtained by the analytic method. In this section, the bifurcation behaviors of the system under nonresonance (NR) are compared with those under resonances by numerical simulations. Consider that β1 is usually a small quantity; β1=0.001 is chosen in the simulations. So ω1≈3 can be obtained. Table 2 shows the values of parameter β2 corresponding to every resonance type. It can be seen that when β1=0.001, six kinds of resonances can occur in the system by choosing the value of parameter β2, where SR2 and CR3 will not occur. And notice that these six kinds of resonance will not occur at the same time.Case 1.
NR
Firstly, the bifurcation behavior of the system under NR with orbital eccentricitye is investigated. The parameters β1=0.001 and β2=10 are selected. No resonance occurs in the system. Figure 3 shows the bifurcation diagrams of the system changing with orbital eccentricity e. As can be seen from the figure, the bifurcation diagrams of θ and α are similar. In e∈0,0.284, they experience bounded motions of periods, period-doubling, and quasiperiods. When e≥0.284, the motions of θ and α are all chaotic. As we all know, for a dumbbell model, the pitch motion θ enters a tumbling state around the eccentricity of 0.313. So when the mother satellite is treated as a rigid body, this system enters a chaotic motion easier than a dumbbell model. But the critical orbital eccentricity does not decrease much. The coupling of the motions of θ and α makes it easier for the system to enter chaotic motions. But the coupling is not very strong in the case of nonresonance.Case 2.
CR1
Next, the bifurcation diagrams ofθ and α are observed under combination resonance. A set of parameters is selected as β1=0.001 and β2=0.18; one can get ω1≈3 and ω2≈3−1, so ω1−ω2≈1. The bifurcation diagrams under CR1 are shown in Figure 4. As can be seen from the figure, the bifurcation diagrams of θ and α show significant differences. The bifurcation diagram of pitch angle θ is similar to that in the case of nonresonance, and the critical value of e is 0.279. However, the bifurcation diagram of α is quite different from that of the nonresonant case. The motion of α enters into chaotic motion at a very small value of e, and the critical value of e is just 0.08. And it is worth noting that the amplitudes of the periodic, doubling-periodic, and quasiperiodic motion of α also increase significantly compared with those of the nonresonance. This type of parametrically excited resonance has a great impact on the motion of α, while it has a small impact on the motion of θ.Case 3.
PR
In this case, the bifurcation diagrams under primary resonance are obtained. Here,β1=0.001 and β2=1 are selected; one can get ω1≈3 and ω2=3, so ω1≈ω2. Figure 5 shows the bifurcation diagrams under the primary resonance. As can be seen from the figure, the bifurcation diagrams of θ and α are similar to those in the case of combination resonance. The critical values of e are 0.282 and 0.133 in the bifurcation diagrams of θ and α, respectively.Table 2
The values of parameterβ2 corresponding to every resonance type at β1=0.001.
The value ofβ2The value ofω2Frequency relationshipTypes of resonanceβ2≈0.18ω2≈3−1ω1−ω2≈1CR1β2≈0.25ω2≈3/2ω1≈2ω2IR2β2≈1ω2≈3ω1≈ω2PRβ2≈1.33ω2≈2ω2≈2SR1β2≈2.49ω2≈1+3ω2−ω1=1CR2β2≈4ω2≈23ω2=2ω1IR1Figure 3
Bifurcation diagrams changing with orbital eccentricitye under nonresonance.Figure 4
Bifurcation diagrams changing with orbital eccentricitye under combination resonance.Figure 5
Bifurcation diagrams changing with orbital eccentricitye under primary resonance.After simulating the bifurcation diagrams of the other four types of resonances atβ1=0.001, the results similar to the above two kinds of resonances can be obtained. That is, the resonances have a great impact on the motion of α, while they have a small impact on the motion of θ. The main effects of resonance on the motion of α are that the amplitudes of its periodic, doubling-periodic, and quasiperiodic motions increase compared with those of nonresonance. The critical value of orbital eccentricity that causes it to move into a chaotic motion is much smaller than that of nonresonance. The main reason for this situation can be obtained by observing the system’s first-order approximate differential Equation (9). In the equation, the motion of θ decouples with α. So internal resonances and parametrically excited resonances related toαhave very little effect on it. The first-order approximation differential equation of α in Equation (9) is coupled with θ, so the motion of α is greatly affected by these kinds of resonances.Case 4.
SR2
By observing Equation (9), it can be seen that the motion of θ is affected by the parametric excitation term, so it can be suspected that the SR2 may have great effects on the motion of θ. A set of parameters is selected as β1=0.33 and β2=10 in bifurcation simulations; one can get ω1≈2 and ω2=30, so the system is under SR2. The bifurcation diagrams under superharmonic resonance are shown in Figure 6. As can be seen from the figure, the bifurcation diagrams of θ and α are similar. The bifurcation diagram of θ is significantly different from those of nonresonance and other types of resonances. The main manifestation is that the amplitudes of the periodic, doubling-periodic, and quasiperiodic motions of θ are significantly larger than those in other cases when the eccentricity is small. Moreover, the critical value of e that makes the motion of θ be chaotic is also reduced to 0.204. As the motion of α is coupled with that of θ, their bifurcation diagrams are similar. The critical value of e is also 0.204 in the bifurcation diagram of α.Figure 6
Bifurcation diagrams changing with orbital eccentricitye under superharmonic resonance.


## 4.2. Bifurcation Diagram Changing with Parameterβ2

The bifurcation behaviors of the system changing with parameterβ2 are investigated next. Given that β1=0.001 and e=0.15, Figure 7 shows the bifurcation diagrams changing with parameter β2. There are six vertical solid lines in the figure; they are the values of parameter β2 under various kinds of resonance when β1=0.001 as shown in Table 2. It can be seen in the figure that the change of parameter β2 has little effects on the motion of θ, that is, the six types of resonance listed in Table 2 have little effect on the motion ofθ. The motion of α varies significantly with the change of parameter β2, and the regions with a large amplitude of α correspond to the positions of resonance regions given in Table 2. So the resonance regions obtained by analytical analysis are in good consistency with the results obtained by numerical simulations. This is also consistent with the above results of bifurcation behaviors with the change of parameter e.Figure 7
Bifurcation diagrams changing with parameterβ2 at β1=0.001.In Figure7, it can also be seen that there are two areas with an obvious sudden increase of amplitude in the bifurcation diagram of α, which are not corresponding to the above resonance regions. The two areas are located in the center of the region β2=3 and β2≈4.64, respectively. They correspond to superharmonic resonance which is ω2=3 and combination resonance which is ω2−ω1=2. The two resonances are ignored in the theoretical analysis because the dynamic equations are truncated by order 2 in the previous perturbation analysis.The above analysis shows that SR2 has a greater effect on the motion ofθ. Here β1=1/3 is selected, so ω1=2, and SR2 occurs in the system. Table 3 shows the values of parameter β2 corresponding to every resonance type at β1=1/3. It can be seen that when β2 takes the values in Table 3, multiple resonances occur simultaneously in the system. Figure 8 shows the bifurcation diagrams changing with parameter β2 at β1=1/3. It can be seen that the bifurcation diagrams of θ and α are similar, and both of them have a sharp increase in amplitude in the resonance regions. It also confirms the above results of bifurcation behaviors with the change of parameter e.Table 3
The values of parameterβ2 corresponding to every resonance type at β1=1/3.
The value ofβ2The value ofω2Frequency relationshipTypes of resonanceβ2=1/3ω2=1ω1−ω2=1, ω1=2ω2CR1, IR2β2=4/3ω2=2ω1=ω2, ω2=2PR, SR1β2=3ω2=3ω2−ω1=1CR2β2=16/3ω2=4ω2=2ω1IR1Figure 8
Bifurcation diagrams changing with parameterβ2 at β1=1/3.From the above, since the first-order approximate equation ofθ is not coupled with α, only SR2 has a great effect on its motion. However, the amplitude of the motion of α is significantly increased under various kinds of resonances and can even be chaotic. In short, resonances have a huge effect on the motion of the system.


## 5. Periodic Control

According to the analysis in the previous section, resonances will make the system more prone to chaotic motion, which will bring harm to the operation of the system. In general, we can make the system avoid the resonance region by selecting appropriate system parameters, but sometimes, when the cost of changing system parameters is higher, control schemes are needed to make the system operate normally. Therefore, this section discusses a control strategy to make the chaotic motion of the system to the periodic orbit in the resonance region.The extended time-delay autosynchronization (ETDAS) has been applied successfully to stabilize the chaotic motion of tethered satellite systems. This control technique has some great advantages that make it widely used. Firstly, it only requires the knowledge of the period of the desired periodic orbit instead of a reference signal of the desired regular motion. Secondly, when the system operates in the neighborhood of the desired periodic solution, the control inputs will take small values. And this method does not require fast switching or sampling that makes it easy to apply in practice. Based on these excellent advantages, the ETDAS will be employed to stabilize the unstable periodic motions in this paper. Assuming that the system is acted upon by additional forces regardless of the actual actuators, the controlled governing equation becomes(16)θ¨=−3sinθcosθ−β1α̇2+2α̇sinθ−α+3cosαsinθ+e21+θ̇sinν+3sinθcosθcosν+F1ν,α¨=β2θ̇2+2θ̇+3cos2θsinθ−α+2e1+α̇sinν+F2ν,where the feedback control inputs F1ν and F2ν are given by
(17)F1ν=k1θ̇ν−1−R1∑j=1∞R1j−1θ̇ν−jτ,F2ν=k2α̇ν−1−R2∑j=1∞R2j−1α̇ν−jτ,where k1 and k2 are the feedback gains and R1 and R2 are the memory parameters. The delayed time τ is periodic of the desired periodic orbit.The basic block chart of the ETDAS control scheme is shown in Figure9.Figure 9
Basic block chart of the ETDAS control method.Pang and Jin [20] have verified the validity of ETDAS for the TSS studied in this paper when it is in a nonresonant region by numerical and experimental methods. Therefore, this section mainly investigates whether the ETDAS method is still valid in resonance regions. The stability analysis method proposed by Bleich and Socolar [21] is adopted, the stability domains of the controlled tethered satellite system have been calculated. For simplicity, the study of the stability is limited to the case of R1=R2=R and k1=k2=k. As previous studies have shown, the system becomes harder to control as the orbital eccentricity increases. Therefore, there is a threshold value e∗ of orbital eccentricity. When the threshold value is exceeded, the method cannot control the system to a desired periodic orbit. The stability domains will be shown as functions of the control parameters R, k, and e∗. Without loss of generality, three cases analyzed in the previous section are selected for comparison, namely, NR, CR1, and SR2.Figure10 gives the stability domains provided by the ETDAS method in the three-dimensional space (R, k, and e∗) for the three cases. In the figure, the upper bound surfaces of the stability domains are given for various cases, and the stability domains are the regions between each surface and the (R, k) plane of e∗=0. The top surface is NR, the middle one is CR1, and the bottom one is SR2. It can be seen from the figure that the stability domains of resonance cases are smaller than that of nonresonance one. Among them, the largest threshold value of e is e∗=0.47 for NR and e∗=0.41 for CR1, and the smallest one is e∗=0.09 for SR2. Besides, it can be seen from the figure that the control stability domains are not particularly sensitive to the change of parameter R. Therefore, in order to show the changes of the stable domains more clearly, Figure 11 shows the stable domains in two-dimensional space (k, e∗) when R=0.5. As can be seen from the figure, the stable domain becomes smaller in resonance cases than in nonresonance cases, among which the SR2 is reduced very much. It can also be seen that when the control parameter k is small, it is of great help to increase the threshold value of eccentricity. When k is large, its increase can hardly increase e∗. So the control purposes can be achieved without setting the value of k too large.Figure 10
The stability domains of the ETDAS method in space (k, R, and e∗) for the three cases.Figure 11
The stability domains of the ETDAS method in space (k, e∗) at R=0.5.


## 6. Periodic Solutions under Resonances

The multiscale method described in Section3 of this paper can be used to obtain the periodic solution of the system in the case of small orbital eccentricity. But for every case of resonances, the method has to be recalculated to find every analytic solution. To simplify the analysis, a straightforward expansion method is used to obtain a 2π basic periodic solution in this section. This method has been successfully used by Peláez et al. [22] to obtain the 2π periodic solution of an electrodynamic tethered system. The basic periodic solutions of the system studied in this paper have been given by Pang et al. [23]. The correctness of the analytical solutions in the nonresonant case has been verified. To prove that the solutions are still reliable in resonance, the numerical and analytical solutions of CR1 and SR2 are compared.Figure12 shows the time histories of θ and α in two orbits corresponding to different orbital eccentricities under CR1. Three cases are given in the figure corresponding to the following values of the orbital eccentricity: e=0.01,0.03,0.05. In the figure, the curve with a larger amplitude corresponds to higher orbital eccentricity, and the straight line corresponds to the numerical solutions and the circle to the analytical solutions. It can be seen from Figure 12 that, within two orbital periods, the analytical solutions of periodic motion are in good agreement with the numerical solutions, and the smaller the eccentricity value is, the better the coincidence is. This shows that the approximate periodic solutions obtained by the perturbation method are correct and reliable in this kind of resonance.Figure 12
Periodic solutions given by asymptotic expansions and numerical results for different values ofe under CR1.Then, the time histories ofθ and α in two orbits corresponding to different orbital eccentricities under SR2 are shown in Figure 13. The meanings of the curves in the figure are consistent with those in Figure 12. It can be seen from Figure 13 that, within two orbital periods, the analytical solutions of periodic motion are in good agreement with the numerical solutions, indicating that the approximate periodic solutions obtained by the perturbation method are correct and reliable under such resonance conditions. It is worth noting that the amplitudes of the periodic solutions corresponding to different orbital eccentricity are larger than those under CR1.Figure 13
Periodic solutions given by asymptotic expansions and numerical results for different values ofe under SR2.To observe the differences in periodic solutions under different resonance conditions, Figure14 shows the phase diagram of basic periodic solutions in phase space in the same orbital eccentricity e=0.05 under NR, CR1, and SR2. It can be seen from the figure that SR2 has the largest amplitude of periodic solution, CR1 has the second-largest amplitude, and NR has the smallest amplitude. The amplitude of the periodic solution in SR2 is about 0.4,0.46, CR1 is about 0.05,0.28, and NR is about 0.05,0.06. The amplitude is much larger in resonance than in nonresonance. To further illustrate, Figure 15 shows the phase diagram of the basic periodic solution in phase space under three cases of NR, CR1, and SR2 at the same orbital eccentricity e=0.1. As can be seen from the figure, the amplitude of SR2 is about 0.71,1.45, CR1 is about 0.1,0.56, and NR is about 0.11,0.12. Obviously, SR2 has a much larger amplitude than NR, especially the motion of α.Figure 14
The phase diagram of basic periodic solutions ine=0.05 under NR, CR1, and SR2.Figure 15
The phase diagram of basic periodic solutions ine=0.1 under NR, CR1, and SR2.The results from the above phase diagram simulations show the major cause for the decrease of the control ability of ETDAS under resonance, because the purpose of ETDAS is to control the motion of the system to a periodic orbit. As can be seen from the above phase diagram simulations, the amplitudes of periodic motion of the system in resonance are larger than those in nonresonance. Even in some resonance cases, the amplitudes increase very much, which directly makes the system more difficult to control in resonance than in nonresonance.


## 7. Conclusion

In this study, resonance motions of tethered satellite systems in elliptical orbits are studied. The perturbation analysis of the system is carried out by using the multiscale method. All possible resonance types and corresponding parameter relations of the system are given, including internal resonances and combination resonances. The resonance parameter domain which is in0<β2≤8 is given to provide a reference for parameter design of a space tethered system. The bifurcation behaviors of the system with orbital eccentricity e and parameter β2 are studied by the numerical method. The results show that the threshold values of e in resonance cases are much smaller than those in nonresonance cases, that is to say, it is easier to enter into chaotic motion under resonances. And the resonance regions obtained by analytical analysis are verified by bifurcation diagram changing with parameter β2.The ETDAS method is applied to stabilize the chaotic motion of the system. Stability analysis shows that stable domains become smaller in resonance cases than in nonresonance cases. The2π periodic solutions are obtained analytically. By comparing the phase diagrams, it can be observed that the amplitudes of the periodic solutions in resonances are much larger than those in nonresonances. This is the main reason for the decrease of the control ability of ETDAS under resonances. In future work, some new control strategies should be designed so that the system can still be controlled to periodic motion from a chaotic one in high orbital eccentricity under resonances.



---

*Source: 8844139-2020-09-21.xml*