The motion of tethered tug–debris system with fuel residuals

The motion of tethered tug–debris system with fuel residuals

Accepted Manuscript The motion of tethered tug-debris system with fuel residuals Vladimir S. Aslanov, Vadim V. Yudintsev PII: DOI: Reference: S0273-1...

1MB Sizes 2 Downloads 49 Views

Accepted Manuscript The motion of tethered tug-debris system with fuel residuals Vladimir S. Aslanov, Vadim V. Yudintsev PII: DOI: Reference:

S0273-1177(15)00469-X http://dx.doi.org/10.1016/j.asr.2015.06.032 JASR 12322

To appear in:

Advances in Space Research

Please cite this article as: Aslanov, V.S., Yudintsev, V.V., The motion of tethered tug-debris system with fuel residuals, Advances in Space Research (2015), doi: http://dx.doi.org/10.1016/j.asr.2015.06.032

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. 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.

The motion of tethered tug-debris system with fuel residuals Vladimir S. Aslanov∗, Vadim V. Yudintsev Samara State Aerospace University (SSAU) 34, Moskovskoye shosse, Samara, 443086, Russia

Abstract Active debris removal using a space tug with a tether is one of the promising techniques to decrease the population of large non-functional satellites and orbital stages in near Earth orbits. Properties of debris should be taken into account in the development of the space tugs. In this paper we consider the motion of a debris objects with fuel residuals that can affect the safety of the debris transportation process. The equations of the attitude motion of the tug-debris system in a central gravitational field are derived. Stationary solutions of the equations are found. The system of linearized equations are introduced that can be used for short term analysis. The numerical simulation results are provided that show good accuracy of the linearized equations. Proposed equations can be used to analyze the attitude motion of the tug-debris system and to determine the conventional parameters for safe tethered transportation of space debris.

Note 2.1

Note 2.2

Keywords: Space debris, Tether, Space tug, Fuel residuals 1. Introduction and problem formulation During the last years several active debris removal methods were developed [1–3]. Active debris removal using tethered space debris is one of the promising techniques to decrease the population of large non-functional satel∗

Corresponding author Email addresses: [email protected] (Vladimir S. Aslanov), [email protected] (Vadim V. Yudintsev) URL: http://aslanov.ssau.ru (Vladimir S. Aslanov)

Preprint submitted to Elsevier

June 27, 2015

Note 2.3

lites and orbital stages on the near Earth orbit [4, 5]. Reference [6] noted that there are two types of large space debris: spacecraft and orbital stages. On the one hand, orbital stages may be more easily deorbited because they don’t have large appendages such as antennas and solar panels. Tethered transportation of spacecraft with flexible appendages is considered in [7, 8]. In the other hand orbital stages may contain fuel residuals that affect to the deorbit process. The effect of liquid fuel slosh on spacecraft has been explored in the literature [9–11]. This literature considers the control of vehicle with fuel slosh dynamics. The aim of the present paper is to develop a simple mathematical model of the tug-debris system with fuel in terms of a multibody system model that can be used to analyze active debris removal missions. The paper is divided into four parts. In the section 2 nonlinear equations of the system are obtained that describe the motion of tug-debris system in a central gravitational field. The attitude motion near the stationary point is considered in the section 3. Simplified linear equations are derived. The numerical simulation results provided in the section 4. These results show that simplified equations give a good approximation to the exact solution.

Note 2.4

2. Mathematical model This section formulates the dynamics of the tug-debris system. The debris is represented as a rigid body with fuel sloshing mass. We suppose that the space tug has an attitude control system that maintains required orientation of the tug, so the tug is represented as a mass point. The tether is considered as a massless rod. We use Lagrange formalism to derive the equations of the relative motion. A scheme of the system is shown in Fig. 1. A solution of the general problem of oscillations of residual fuel in a container is extremely difficult. Here we use the simplest model where the sloshing liquid is modelled as an equivalent pendulum model. Reference [12] demonstrates that equivalent pendulum model can approximate motion of the fuel residuals [13]. This model can be used when the oscillations of liquid are small [14, 15]. 2.1. Kinematics of the system The plain motion of the tug-debris system is considered in orbital frame Cxo yo attached to the center of mass of the entire system C. The motion

2

Note 2.6

Figure 1: Orbital frame

of the system occurs throughout the action of the thrust F and a central gravity force. The thrust F is assumed to act along the axis Cxo . The position of the debris relative to Cxo yo frame is described by the angle θ + α and the vector R2 . The angle α defines the orientation of the tether. The angle θ is an angle between the tether and the longitudinal axis of the debris (Fig. 2). The tether length is l1 . The tether attachment point is determined by the vector ρ1 . We suppose that the tether is attached at the longitudinal axis of debris ρ1 = {x1 , 0}T . The attachment point of the equivalent pendulum is defined by the column vector ρ3 = {x3 , 0}T . The pendulum length is l3 . The angle φ of the pendulum with respect to the debris longitudinal axis representing the fuel slosh (Fig. 2). The column vector R2 denotes the position of the center of mass (C2 ) of the debris relative to the frame Cxo yo . The position column vector of the tug is expressed as R1 = R2 + A(θ + α) · ρ1 + A(α) · ex l1

(1)

where ex = [1, 0]T is an unit column vector, A(α), A(θ + α) are rotation matrices [ ] [ ] cos α − sin α cos(θ + α) − sin(θ + α) A(α) = , A(θ + α) = sin α cos α sin(θ + α) cos(θ + α) Position of the pendulum (fuel lump) relative to the frame Cxo yo can be 3

Note 2.1,2.8

Note 1.1 Note 1.2, 1.5, 2.9

Tug Tether Debris

Fuel tank Figure 2: The space debris and the space tug

computed as R3 = R2 + A(θ + α) · (ρ3 − A(φ) · ex l3 ) where

(2)

[

] cos φ − sin φ A(φ) = . sin φ cos φ

The motion of tug, debris and the fuel is considered relative to the center of mass of the system therefore we can write 3 ∑

Ri mi = 0

(3)

i=1

where mi is a mass of the tug (i = 1), space debris (i = 2) and fuel (i = 3). Using (3) we can eliminate R2 from the expressions (1) and (2).

4

2.2. Kinetic energy and Lagrange equations The kinetic energy of the relative motion of the tug debris and the fuel slug is given by 3 ∑ 2T = mi Vi2 + Jz (θ˙ + α) ˙ 2 (4) i=1

where Jz is the moment of inertia of the debris. We suppose we know all the moments inertia of the debris Jx , Jy , Jz and Jy = Jz ; V i = dRi /dt is the velocity of the body i: [ ] 1 l1 α(m ˙ 1 − M )sα − l3 m3 ω3 sα+θ+φ + ω2 (x1 (m1 − M ) + m3 x3 )sα+θ V1= ˙ − m1 )cα + l3 m3 ω3 cα+θ+φ + ω2 (x1 (M − m1 ) − m3 x3 )cα+θ M l1 α(M [ ] 1 l1 m1 αs ˙ α − l3 m3 ω3 sα+θ+φ + ω2 (m1 x1 + m3 x3 )sα+θ V2= ˙ α + l3 m3 ω3 cα+θ+φ − ω2 (m1 x1 + m3 x3 )cα+θ M −l1 m1 αc [ ] 1 m1 (l1 αs ˙ α + x1 ω2 sα+θ ) − (M − m3 )(x3 ω2 sα+θ − l3 ω3 sα+θ+φ ) V3= ˙ α − l3 ω3 (M − m3 )cα+θ+φ − ω2 cα+θ (x3 (m3 − M ) + m1 x1 ) M −l1 m1 αc where

Note 1.4

˙ ω2 = α˙ + θ,

ω3 = ω2 + φ˙

The letters c and s with a variable or an expression in the subscript denote the cosine and sine of variable or expression respectively. Gravitational force Gi that act on the body i is Gi = −

µmi ri r3i

(5)

where µ is the Earth’s gravitational parameter, ri is the position vector of the body i relative to the Earth’s center r i = r + Ri , i = 1, 2, 3

(6)

r = {0, r}T is the position column vector of the system center of mass. The column vector of the thrust force is F = {F, 0}T . We assume F = const. The orbital frame Cx0 y0 is a non inertial frame so we have to add inertial forces Φi = −mi (ao + ω o × (ω o × Ri ) + εo × Ri + 2ω o × V i ) (7)

5

where ao is acceleration of the system’s center of mass [ ] −r¨ ν − 2r˙ ν˙ ao = r¨ − rν˙ 2

(8)

ω o , εo are the angular velocity and the angular acceleration of the frame Cx0 y0 ω o = {0, 0, ν} ˙ T , εo = {0, 0, ν¨}T , (9) To get derivatives r, ˙ ν, ˙ r¨, u¨ let us write the equations for the center of mass of the system in osculating elements ν, p, e, σ [16]: √ ν˙ = µ(1 + e cos ν), (10) pF p˙ = −2 , (11) mv e + cos ν e˙ = −2 F, (12) mv cos ν σ˙ = −2 F (13) mve where ν = u − σ is a true anomaly angle, p is the focal parameter of the orbit, e is the eccentricity, σ is the angle of periapsis, v is the orbital velocity √ µ(1 + e2 + 2e cos ν) v= , (14) p r is a distance of the center of mass of the system relative to the center of the Earth p r= . (15) 1 + e cos ν The derivatives r, ˙ ν, ˙ r¨, u¨ can be obtained after differentiating the equations (13) and (15). The space debris is considered as a rigid body, therefore a gravitational torque should be taken into account [17] M2z =

3µ (Jy − Jx ) sin 2(θ + α) 2r23

(16)

The space tug and the fuel lump are considered as point masses, so M1z = M3z = 0. 6

Note 1.8

Now we can write generalized forces [18] 3 ∑ ∂Ri ∂R1 Qk = (Gi + Φi ) · +F · + Mkz , ∂qk ∂qk i=1

k = 1, 2, 3.

(17)

where qk is a generalized coordinate q1 = θ, q2 = α, q3 = φ

(18)

Using (4) and (17) we can build the differential equations of the considered system d ∂T ∂T − = Qi , i = 1, 2, 3 (19) dt ∂ q˙i ∂qi 2.3. Simplified equations The obtained equations (19) can be integrated but they are very cumbersome and inconvenient for motion analysis. To determine the stationary solutions and then to study of the motion near the stationary point, let us write simplified equations taking the following assumptions • the orbit of the mass center is not changed r = const, but we take into account the non inertial motion of the orbital frame, • the tether length is small relative to r. The inertial forces can be simplified as follows. The centrifugal force is Φwi = −m[ao + ω o × (ω o × Ri )]

(20)

√ where ωo = µ/r3 = const is the orbital angular velocity, ao is a acceleration of the mass center ao = −ωo2 r (21) The second term of (20) can be rewritten as ω o × (ω o × Ri ) = −ωo2 Ri , so the centrifugal force get the form Φwi = mωo2 (r + Ri ) 7

(22)

For the circular orbit ωo2 =

µ , r3

(23)

so the expression (5) get the form [ ]3 µmi r 2 Gi = − 3 ri = −ωo mi ri ri ri

(24)

The expression in the braces for Ri << r can be written [ ]3 [ ]−3/2 r 2yi yi ≈ 1+ ≈1−3 ri r r This simplification allows to rewrite Gi [ yi ] Gi ≈ −ωo2 mi 1 − 3 ri r

(25)

The generalized forces get the form ∂R1 ∑ F yi ∂Ri Qk = F · + mi (− + 3ωo2 r − 2ω o × v i ) · + Q∗k , ∂qk M r ∂qk i=1 3

k = 1, 2, 3 (26) where ω o = {0, 0, ωo }T , r = {0, r, 0}; Q∗2 = Q∗3 = 0, Q∗1 =

3µ (Jz − Jx ) cos(θ + α) sin(θ + α) − Jz ν¨ r3

(27)

New generalized forces allow to write the simplified equations of the motion, which will be used to study an evolution of the system around stationary points. The equations are given in the Appendix. 3. Motion of the system near the stationary point 3.1. Stationary solution Considered tug-tether-debris system can be represented as a two mass system connected with a massless rod. In central gravitational field this two mass system has two stationary points α01 = 0 (unstable) and α02 = π/2 (stable) [17]. Tug’s thrust F shifts stable stationary point to α02 < π/2. 8

Note 1.9

This stationary point depends on the tug’s thrust, length of the tether and masses of the tug and the debris. To determine the stationary solutions of the equations, the derivatives are set to zero θ˙ = φ˙ = α˙ = 0, θ¨ = φ¨ = α ¨=0 (28)

Note 2.11

In this case, the system of equations (19) converted to a non-linear system of equations for the unknown angles θ0 , α0 , φ0 . To simplify the search for solutions of this system let us obtain an approximation to get the stationary solution of the system. We equate Qα (26) to zero and set θ = 0 and φ = 0. We get (a cos α − F b) sin α = 0

(29)

where a, b are coefficients that depend on the parameters of the system [ a = 3ω02 (x1 + l1 )2 M1 m1 + (x3 − l3 )2 M3 m3 + 2m1 m3 (l3 − x3 )(l1 + x1 )]

Note 2.12

b = M1 (x1 + l1 ) + m3 (l3 − x3 ) where M = m1 + m2 + m3 , Mi = M − mi , i = 1, 2, 3. Fig. 3a demonstrates the generalized force Qα as a function of α. The figure shows two stationary points α01 = 0 and α02 ≈ 0.35. Fig. 3 is plotted for parameters presented in the Table 1. Eq. (29) shows that the system a) Qα (α): φ = 0, l1 = 900 m

b) α0 as a function of l1 1.4

2

1.2

0

0



0.8 0

Q ,N

1.0 -2

0.6

-4

0.4 -6

0.2

-8

0.0 0.0

0.1

0.2

0.3

0.4

500

0.5

1000

1500

2000

l1 , m

, rad

Figure 3: Generalized force Qα as a function of θ and stationary point α0 as a function of the tether length

has two stationary points. The first stationary point is determined by the 9

Note 1.6,1.7

condition sin α = 0. The second stationary solution is determined by the condition a cos α − F b = 0 (30) that exists only if F b/a ≤ 1. Fig. 3b shows stationary solution α0 as a function of tether length l1 for F = 0.3 N and the parameters of the system that are presented in the Table 1. Fig. 3b demonstrates that there is only one stationary solution for the tether length l1 < l1∗ ≈ 400 m α0 = 0,

θ0 = 0.

(31)

There are two stationary solutions for l1 > l1∗ α01 = 0, θ01 = 0 and α02 > 0, θ02 > 0.

(32)

To reduce the influence of the tether force to the motion of the space tug tether length should be taken that results only one stationary solution α = 0. In this case, the disturbing torque of the tether force N , acting to the tug is minimal. Also the motion with a small angle α ≈ 0 reduces risk of tether rupture due exposure to the tug’s jet Fig. 4.

Figure 4: Exhaust blast area of the space tug

3.1.1. Linearized equations After determining the stationary solutions we can construct linearized equations of motion in the neighborhood of the stationary point θ0 = α0 = φ0 = 0. We rewrite the kinetic energy of the system as a quadratic form with constant coefficients T (q, q) ˙ ≈ T (q0 , q) ˙ (33) 10

Note 2.13

where q0 = (θ0 , α0 , φ0 ) is the stationary point. In the expressions for the ˆα generalized forces we set q = q0 + qˆ, were qˆ = (θ, ˆ , φ) ˆ are new variables that describe the deviation from the stationary point q0 . We expand the generalized forces in series of qˆ leaving only terms of first order of qˆ. Using the new expressions for the kinetic energy (33) and the generalized forces we get equations in well know form 3 ∑

(aij q¨ˆj + bij qˆj ) = 0,

i = 1, 2, 3

(34)

  b11 b12 b13 B = b21 b22 b23  b31 b32 b33

(35)

j=1

  a11 a12 a13 A = a21 a22 a23  , a31 a32 a33

The coefficients of the matrices A and B are a11 = 2l3 m3 (m1 x1 − M3 x3 ) + m1 x1 (M1 x1 − 2m3 x3 )+ m3 M3 x23 + Jz M + l32 m3 M3 (36) a12 = a21 = l3 m3 [m1 (2x1 + l1 ) − 2M3 x3 + l3 M3 ] − m1 m3 x3 (2x1 + l1 ) + m1 M1 x1 (x1 + l1 ) + m3 M3 x23 (37) a13 = a31 = l3 m3 [m1 x1 + M3 (l3 − x3 )]

(38)

b11 = (x1 (m1 (6m3 ω02 (x3 − l3 ) − F ) + F M ) − 3m1 M1 ω02 x21 + 3(Jz − Jx )M µ m3 (−F x3 − 3M3 ω02 (l3 − x3 )2 )) + F l3 m3 − (39) r3 3µJ21 = x1 (m1 (−3ω02 (2m3 (l3 − x3 ) + l1 M1 ) − F ) + F M )− r3 3m1 M1 ω02 x21 + m3 (3ω02 (x3 − l3 )(M3 (l3 − x3 ) + l1 m1 ) − F x3 )+ 3(Jz − Jx )µM F l3 m3 − (40) r3

b12 = b21 −

b13 = b31 = l3 m3 (F − 3ω02 (m1 x1 + M3 (l3 − x3 ))) 11

(41)

a22 = −2l1 m1 (−M1 x1 + m3 x3 − l3 m3 )+ 2l3 m3 (m1 x1 − M3 x3 ) + m1 x1 (M1 x1 − 2m3 x3 )+ m3 M3 x23 + l32 m3 M3 + l12 m1 M1 (42) a23 = a32 = l3 m3 [m1 (x1 + l1 ) + M3 (l3 − x3 )]

(43)

b22 = l3 m3 (6ω02 (M3 x3 − m1 (l1 + x1 )) + F )+ F (M1 (x1 + l1 ) − m3 x3 )− 2 3ω0 (m1 M1 (x1 + l1 )2 + m3 M3 x23 − 2m1 m3 x3 (x1 + l1 ))− 3l32 m3 M3 ω02 (44) { } b23 = b32 = l3 m3 F − 3ω02 [m1 (x1 + l1 ) + M3 (l3 − x3 )] (45) a33 = l32 m3 M3

(46)

b33 = l3 m3 (F − 3l3 M3 ω02 )

(47)

Using the linear equations (34) we can get the natural frequencies of the system. The solutions of the equations (34) have the form qj = Cj sin λt,

j = 1, 2, 3

(48)

Substituting (48) into (34) we get det(Aλ2 + B) = 0

(49)

that allows us to find three frequencies λ1 , λ2 , λ3 . Fig. 5a shows the frequencies of the system as functions of tether length. Fig. 5b and as functions of the thrust F . 4. Numerical example Let us compare the solutions of the nonlinear system (19) to the solutions of the linearized system (34) with the following initial conditions θ0 = 0.1, α0 = 0.3, φ0 = 0, θ˙0 = α˙ 0 = φ˙ 0 = 0.

(50)

The parameters of the system are shown in the Table 1. The solutions are obtained for two cases. Case 1 for l1 = 30 m, F = 2 N and case 2 for l1 = 300 m, F = 2 N. The simulation results are shown in Fig. 6. The comparison of results show good accuracy of the approximate linearized model. 12

Table 1: Parameters of the space tug and the debris

Parameter

Value

Parameter

Value

m1 ,

200

m2 , kg

3000

m3 , kg

500

l3 , m

1

x3 , m

−1

x1 , m

5

Jy = Jz , kg · m2

10000 Jx , kg · m2

3000

b) λj (F ), j = 1, 2, 3. l1 = 100 m

a) λj (l1 ), j = 1, 2, 3

0.04 0.012 0.010

0.03

, 1/s

0.006

0.02

i

i

, 1/s

0.008

0.004 0.01 0.002 0.000

0.00 50

100

150

200

250

300

350

400

l1 , m

0.5

1.0

1.5

2.0

2.5

3.0

F, N

Figure 5: The natural frequencies of the system as functions of the tether length l1 for F = 0.3 N and as functions of F for l1 = 100 m, j=1,2,3

5. Conclusion The motion equations of debris with the fuel residuals during tethered transportation are derived. The stationary solutions of the equations are found. The linearized differential equations with constant coefficients governing the motion of the system near the stationary point are derived. It has been demonstrated that the solutions obtained by means of the linearized system are in good agreement with the solutions of the original nonlinear system of equations. The proposed simplified equations can be readily used in practice to investigate the motion of the tug-debris system with fuel residuals for different parameters of the system.

13

a) Case 1

b) Case 2 L

L 0.10

0.10 0.05



, rad

, rad

0.05 0.00



-0.05

0.00 -0.05

-0.10 -0.10 0

200

400

600

800

1000

0

200

400

t, s

600

1000

800

1000

800

1000

L

0.3

0.3

0.2

0.2

0.1

0.1 , rad

, rad

L



800

t, s



0.0 -0.1

0.0 -0.1

-0.2

-0.2

-0.3

-0.3 0

200

400

600

800

1000

0

200

400

t, s

600 t, s

L

L 0.4

0.5



, rad

, rad

0.2 0.0



0.0 -0.2

-0.5 -0.4 0

200

400

600

800

1000

t, s

-0.6

0

200

400

600 t, s

Figure 6: The solutions of the nonlinear θ, α, φ and linearized θL , αL , φL equations for to cases l1 = 30 m (case 1) and l1 = 300 m (case 2)

14

6. Acknowledgments TThis work is funded by the Russian Foundation for Basic Research, Project No. 15-01-01456 A. References [1] X. Dafu, K. Xianren, Tether modeling study on electro-dynamic tether deorbiting system, Acta Aeronautica et Astronautica Sinica 5 (2008) 018. [2] R. L. Forward, R. P. Hoyt, C. W. Uphoff, Terminator Tether : A Spacecraft Deorbit Device, Journal of spacecraft and rockets 37 (2) (2000) 187–196. [3] S. Kitamura, Y. Hayakawa, S. Kawamoto, A reorbiter for large GEO debris objects using ion beam irradiation, Acta Astronautica 94 (2) (2014) 725–735. doi:10.1016/j.actaastro.2013.07.037. [4] L. Jasper, H. Schaub, C. Seubert, V. Trushlyakov, E. Yutkin, Tethered tug for large low earth orbit debris removal, in: AAS/AIAA Astrodynamics Specialists Conference, Charleston, 2012. [5] L. Jasper, H. Schaub, Input shaped large thrust maneuver with a tethered debris object, Acta Astronautica 96 (2014) 128–137. doi:10.1016/j.actaastro.2013.11.005. [6] C. Bonnal, J.-M. Ruault, M.-C. Desjean, Active debris removal: Recent progress and current trends, Acta Astronautica 85 (2013) 51–60. doi:10.1016/j.actaastro.2012.11.009. [7] V. S. Aslanov, V. V. Yudintsev, Behaviour of tethered debris with flexible appendages, Acta Astronautica 104 (1) (2014) 91–98. doi:10.1016/j.actaastro.2014.07.028. [8] V. S. Aslanov, V. V. Yudintsev, Dynamics, Analytical Solutions and Choice of Parameters for Towed Space Debris with Flexible Appendages, Advances in Space Research 55 (2) (2014) 660–667. doi:10.1016/j.asr.2014.10.034.

15

[9] M. Reyhanoglu, J. Rubio Hervas, Nonlinear dynamics and control of space vehicles with multiple fuel slosh modes, Control Engineering Practice 20 (2012) 912–918. doi:10.1016/j.conengprac.2012.05.011. [10] B.-Z. Yue, Study on the Chaotic Dynamics in Attitude Maneuver of Liquid-Filled Flexible Spacecraft, AIAA Journal 49 (10) (2011) 2090– 2099. doi:10.2514/1.J050144. [11] J. Rubio Hervas, M. Reyhanoglu, Thrust-vector control of a three-axis stabilized upper-stage rocket with fuel slosh dynamics, Acta Astronautica 98 (2014) 120–127. doi:10.1016/j.actaastro.2014.01.022. [12] L. D. Peterson, E. F. Crawley, R. J. Hansman, Nonlinear fluid slosh coupled to the dynamics of a spacecraft, AIAA journal 27 (9) (1989) 1230–1240. [13] M. Reyhanoglu, Modeling and control of space vehicles with fuel slosh dynamics, Advances in Spacecraft Technologies 3 (2010) 549–562. [14] S. Abramson, S. Silverman, The dynamic behaviour of liquids in moving containers, Tech. rep., NASA, No. SP-106, Washington (1966). [15] R. Ibrahim, Liquid sloshing dynamics: theory and applications, Cambridge University Press, 2005. [16] D. Okhotsimskii, Study of the Motion in a Central Field under the Action of Constant Tangential Acceleration, Kosmicheskie issledovaniya 2 (6) (1964) 817–842. [17] V. V. Beletskii, Motion of an Artificial Satellite about its Center of Mass, Israel Program for Scientific Translations, 1966. [18] J. R. Taylor, Classical mechanics, University Science Books, 2005.

16

Appendix A. The nonlinear equations of motion ( ) (J2 θ¨ + J2 u¨)M 2 + l3 m3 {x3 m3 (m2 − 2M ) + m23 + M 2 + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + m21 x1 }β˙ 2 sβ−γ + l1 m1 {l3 m3 M cα−γ ( ) + [x1 m1 (m2 + m3 − 2M ) + m21 + M 2 − m3 M x3 ]cα−β }¨ α ( ) + l3 m3 {l3 m3 (m1 + m2 − 2M ) + m23 + M 2 ( ) − [x3 m3 (m2 − 2M ) + m23 + M 2 + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + m21 x1 ]cβ−γ }¨ γ 2 2 + {−l3 m3 [x3 (m3 (m2 − 2M ) + m3 + M ) + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + m21 x1 ]cβ−γ + m1 [2m3 x3 x1 (m2 + m3 − 2M ) + m23 x23 + M 2 x21 ] ( ) + m3 x23 m3 (m2 − 2M ) + m23 + M 2 + m2 x1 (x1 (m2 + m3 − 2M ) + 2m3 x3 ) + m3 x2 }β¨ 1

1 1 2

+ 2l3 m3 {x3 [m3 (m2 − 2M ) + +M ] + m1 [x1 (m2 + m3 − 2M ) + m3 x3 ] + m21 x1 } ω0 φs ˙ φ 2 2 + 2l1 m1 {x1 [m1 (m2 + m3 − 2M ) + m1 + M ]sθ + M l3 m3 sθ+φ − m3 M x3 sθ }ω0 α˙ = l3 m3 {x3 [m3 (m2 − 2M ) + m23 + M 2 ] + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + m21 x1 }γ˙ 2 sβ−γ + α˙ 2 l1 m1 {l3 m3 M sα−γ − m3 M x3 sα−β + x1 [m1 (m2 + m3 − 2M ) + m21 + M 2 ]sα−β } + 3J21 µM 2 p−3 cβ sβ − F l3 m3 M sγ − F [x1 (m1 (m2 + m3 − 2M ) + m21 + M 2 ) − m3 M x3 ]sβ ) 3ω0 2 ( {l3 m3 m3 (m1 + m2 − 2M ) + m23 + M 2 s2γ 2 − 2l3 m3 [(x3 (m21 x1 + m3 (m2 − 2M ) + m23 + M 2 ) + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ))s2α+2θ+φ − l1 m1 M cγ sα ] ( ( ) ) cβ {l1 m1 sα x1 m1 (m2 + m3 − 2M ) + m21 + M 2 − m3 M x3 + sβ {m1 (2m3 x3 x1 (m2 + m3 − 2M ) + m23 x23 + M 2 x21 ) + m3 x23 (m3 (m2 − 2M ) + m23 + M 2 ) + m21 x1 (x1 (m2 + m3 − 2M ) + 2m3 x3 ) + m31 x21 }} (A.1) m23

17

− 3(m1 (M 2 + m21 + m1 (−2M + m2 + m3 ))s2α l12 + 2m1 (M l3 m3 s2α+θ+φ + s2α+θ ((M 2 + m21 + m1 (−2M + m2 + m3 ))x1 − M m3 x3 ))l1 + l32 m3 (M 2 + m23 + (−2M + m1 + m2 )m3 )s2γ − 2l3 m3 s2α+2θ+φ (x1 m21 + ((−2M + m2 + m3 )x1 + m3 x3 )m1 +(M 2 +m23 +(m2 −2M )m3 )x3 )+s2β (x21 m31 +x1 ((−2M +m2 +m3 )x1 +2m3 x3 )m21 + (M 2 x21 + 2m3 (−2M + m2 + m3 )x3 x1 + m23 x23 )m1 + m3 (M 2 + m23 + (m2 − 2M )m3 )x23 ))ω02 ˙ 1 m1 (M l3 m3 sθ+φ + sθ ((M 2 + m2 + m1 (−2M + m2 + m3 ))x1 − M m3 x3 ))ω0 − 4θl 1

+4φl ˙ 3 m3 (sφ (x1 m21 +((−2M +m2 +m3 )x1 +m3 x3 )m1 +(M 2 +m23 +(m2 −2M )m3 )x3 ) − M l1 m1 sθ+φ )ω0 + 2F (l1 (M 2 + m21 + m1 (−2M + m2 + m3 ))sα + M l3 m3 sγ + sβ ((M 2 + m21 + m1 (m2 + m3 − 2M ))x1 − M m3 x3 )) +2(−l1 m1 (M l3 m3 sα−γ +sα−β ((M 2 +m21 +m1 (−2M +m2 +m3 ))x1 −M m3 x3 ))α˙ 2 +α ¨ l1 m1 (M cα−γ l3 m3 + l1 (M 2 + m21 + m1 (−2M + m2 + m3 )) + cα−β ((M 2 + m21 + m1 (−2M + m2 + m3 ))x1 − M m3 x3 )) − γ˙ 2 l3 m3 (sβ−γ (x1 m21 + ((−2M + m2 + m3 )x1 + m3 x3 )m1 + (M 2 + m23 + (m2 − 2M )m3 )x3 ) − M l1 m1 sα−γ ) + β˙ 2 (l1 m1 sα−β ((M 2 + m2 + m1 (−2M + m2 + m3 ))x1 − M m3 x3 ) 1 2 +l3 m3 sβ−γ (x1 m1 +((−2M +m2 +m3 )x1 +m3 x3 )m1 +(M 2 +m23 +(m2 −2M )m3 )x3 )) + γ¨ l3 m3 (M cα−γ l1 m1 + l3 (M 2 + m23 + (−2M + m1 + m2 )m3 ) − cβ−γ (x1 m21 + ((−2M + m2 + m3 )x1 + m3 x3 )m1

¨ β m1 M1 x1 (cγ l3 m3 + (M 2 + m23 + (m2 − 2M )m3 )x3 )) + β(c − cβ x3 m3 + M1 (cα l1 + cβ x1 )) − cβ m1 m3 x3 (cγ l3 m3 − cβ x3 m3 + M1 (cα l1 + cβ x1 )) + cβ m2 (m1 x1 + m3 x3 )(cα l1 m1 − cγ l3 m3 + cβ (m1 x1 + m3 x3 )) + m2 sβ (m1 x1 + m3 x3 )(l1 m1 sα − l3 m3 sγ + sβ (m1 x1 + m3 x3 )) − m1 m3 sβ x1 (M3 (sβ x3 − l3 sγ ) − m1 (l1 sα + sβ x1 )) + m3 M3 sβ x3 (M3 (sβ x3 − l3 sγ ) − m1 (l1 sα + sβ x1 )) + cβ m1 m3 x1 (m1 (cα l1 + cβ x1 ) + M3 (cγ l3 − cβ x3 )) − cβ m3 M3 x3 (m1 (cα l1 + cβ x1 ) + M3 (cγ l3 − cβ x3 )) + m1 sβ (M1 x1 − m3 x3 )(l1 M1 sα + l3 m3 sγ + sβ (M1 x1 − m3 x3 )))) = 0 (A.2) 18

− 3ω02 (2l1 m1 M cγ sα − 2m21 x1 cγ sβ − 2m2 m3 x3 cγ sβ + l3 (m3 (m1 + m2 − 2M ) + m23 + M 2 )s2γ + 2m1 M x1 s2α+2θ+φ + 2m3 M x3 s2α+2θ+φ − 2m1 M x1 sφ − 2m3 M x3 sφ − m1 m2 x1 s2α+2θ+φ − m1 m3 x1 s2α+2θ+φ − m23 x3 s2α+2θ+φ − m1 m3 x3 s2α+2θ+φ + m1 m2 x1 sφ + m1 m3 x1 sφ + m23 x3 sφ + m1 m3 x3 sφ − M 2 x3 s2α+2θ+φ + M 2 x3 sφ ) ¨ β−γ (x3 (m3 (m2 − 2M ) + m2 + M 2 ) + 2(l1 m1 M α ¨ cα−γ − βc

3 2 +m1 (x1 (m2 +m3 −2M )+m3 x3 )+m1 x1 )+l3 γ¨ (m3 (m1 +m2 −2M )+m23 +M 2 ) + α˙ 2 l1 m1 (−M )sα−γ + β˙ 2 (x3 (m3 (m2 − 2M ) + m23 + M 2 ) + m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + m21 x1 )sβ−γ ) + 2F M sγ − 4αω ˙ 0 (sφ (m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + M32 x3 )

˙ 0 sφ (m1 (x1 (m2 + m3 − 2M ) + m3 x3 ) + M 2 x3 ) − l1 m1 M sθ+φ ) − 4θω 3 2 ˙ − 4βω0 (m x1 + m2 m3 x3 )sφ = 0 (A.3) 1

where γ = α + θ + φ, β = α + θ.

19