# The Hamiltonian Structure-Preserving Control and Some Applications to Nonlinear Astrodynamics

Authors: Ming Xu; Yan Wei; Shengli Liu
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

A systematic research on the structure-preserving controller is investigated in this paper, including its applications to the second-order, first-order, time-periodic, or degenerated astrodynamics, respectively. The general form of the controller is deduced for the typical Hamiltonian system in full feedback and position-only feedback modes, which is successful in changing the hyperbolic equilibrium to an elliptic one. With the poles assigned at any different positions on imaginary axis, the controlled Hamiltonian system is Lyapunov stable. The Floquet multiplier is employed to measure the stability of time-dependent Hamiltonian system, because the equilibrium of periodic system may be unstable even though the equilibrium is always elliptic. One type of periodic orbits is achieved by the resonant conditions of control gains, and another type is making judicious choice in the foundational motions with different frequencies. The control gains are selected from the viewpoint of both the local and global optimizations on fuel cost. This controller is applied to some astrodynamics to achieve some interesting conclusions, including stable lissajous orbits in solar sail’s three-body problem and degenerated two-body problem, quasiperiodic formation flying on aJ2-perturbed mean circular orbit, and controlled frozen orbits for a spacecraft with a high area-to-mass ratio.


---

## Body

## 1. Introduction

It is well known that most of the astrodynamical problems could be classified as hyperbolic Hamiltonian systems, for example, circular restricted 3-body problem (CR3BP). There are the hyperbolic equilibria which have stable, unstable, and centre manifolds for hyperbolic Hamiltonian system, just like the collinear libration pointsL1,L2, and L3 in CR3BP [1], of which the unstable manifolds may cause instability in the Lyapunov sense.Some astronautical missions, requiring increasing the coverage rate to the object (ground station or target spacecraft) or avoiding communication signals being lost in the sun, are preferring bounded motions near hyperbolic equilibria, such as lissajous (or halo) orbit generated in CR3BP [1, 2]. There presented two different approaches to generate bounded relative trajectories, of which one is obtaining the bounded orbits by making judicious choice in initial conditions. Certain zero-steady controllers could be applied to stabilize the nominal trajectories, owing to the hyperbolicity inheriting from the topological type of the equilibrium.However, for the other case without a natural trajectory to track, some station keeping maneuvers are necessary to implement modification of the topology of the equilibrium. This paper proposes a structure-preserving controller to generate bounded trajectories around the equilibrium, which has many potential applications in astrodynamics. Firstly, a simple structure-preserving controller was constructed by Scheeres et al. [3], in which involved only the instantaneous stable and unstable manifolds, and the manifold’s gains are constrained to be equal. By choosing the proper manifolds’ gains, M. Xu and S. Xu [4] have demonstrated that the poles of the typical second-order Hamiltonian system can be allocated to any position along the imaginary axis and then obtained a new quasiperiodic orbit type referred to as a stable lissajous orbit. Therefore, the controller was powerful enough to work out some famous and difficult problems in sail’s restricted three-body problem that there exist the bounded trajectories near the off-axis solar sail equilibrium [4], and in maintaining the eccentricity in the neighborhood of the unstable fixed point. Besides, the controller implemented the stabilization on a saddle-node bifurcation point for the degenerated case, which makes the modified elliptic equilibrium unique for the controlled system to change the instability of the equilibrium [5]. Further work was implemented on the stabilization of the equilibrium for time-periodic system, which has time varying topological types and no fixed-dimensional unstable/stable/center manifolds [6]. The 2-dimensions temporally independent systems possess fixed-dimensions of unstable/stable/center manifolds and a permanent pair of hyperbolic eigenvalues. Tabuada and Pappas [7] focused on local accessibility preserving abstractions, and provided conditions under which local accessibility properties of the abstracted Hamiltonian system are equivalent to the local accessibility properties of the original Hamiltonian control system.A systematic investigation on the structure-preserving controller is presented in this paper, including some applications in second-order, first-order, time-periodic, and degenerated Hamiltonian astrodynamics, respectively. The general form of the controller is designed for the typical second-order Hamiltonian system, which can work as full feedback and position-only modes. The controller is successful in changing the hyperbolic equilibrium (saddle) to an elliptic one (center) with the poles on the imaginary axis (marginal stability), and then the controlled Hamiltonian system can achieve the Lyapunov stability by means of the Morse lemma. In contrast to the time-independent Hamiltonian system, the equilibrium of periodic system may be unstable even though the equilibrium is always elliptic during its period. Therefore, it is necessary to check whether the Floquet multipliers lie on the unit circle in the complex plane, or whether the moduli of all the multipliers are equal to 1. Two different approaches to generate the periodic orbits are discussed, one of which is achieved by the resonant conditions of control gains for the resonant periodic orbits, and another is making judicious choice in the foundational motions with different frequencies. Furthermore, some selection techniques on the gains are considered from the point view of both the local and global optimizations.This controller is applied to some astrodynamics to achieve some interesting conclusions, including stable lissajous orbits in solar sail’s three-body problem and degenerated two-body problem, quasiperiodic formation flying on aJ2-perturbed mean circular orbit, controlled frozen orbits for a spacecraft with high area-to-mass ratio.


## 2. Hamiltonian System and Hyperbolic Equilibrium

### 2.1. Typical Hamiltonian System in Astrodynamics

For the typical Hamiltonian system in astrodynamics, the Hamiltonian function has the following form as [4]:
(1)H=12pTp+ωpTJq-V(q,t),
where V is the general form of pseudopotential function relying only on the coordinate q and time; ω is the characteristic parameter introduced by the Coriolis acceleration, such as the mean angular motion of the synodic frame used in CR3BP, and Local Vertical-Local Horizontal (LVLH) frame used in formation flying.In general, the variables(q,p) can be expressed by the physical variables (r,r˙) from the Legendre transformation:
(2)[qp]=[I0-ωJI][rr˙],
where I is the identical matrix, and J is the symplectic matrix. The astrodynamics in physical variables can be expressed as
(3)x¨-2ωy˙=∂V∂x,y¨+2ωx˙=∂V∂y.
For a hyperbolic Hamiltonian system, there exist hyperbolic equilibria that have stable, unstable, and center manifolds [1], which are similar to the collinear libration points L1, L2, and L3 in CR3BP. Hence, the equilibrium of the typical 2-dimensional Hamiltonian system in astrodynamics can be obtained from the following equations:
(4)r=r0,r˙=0∂V∂r|r0=0.
Denote Vrr as the second derivative matrix of the pseudopotential function V about the position vector r=[x,y], and the elements of the matrix have the general expression as Vmn=∂2V/∂m∂n,(m,n)=(x,y). And then a hyperbolic (saddle) equilibrium results from Vxx·Vyy-Vxy2<0.Generally, unstable manifolds of the hyperbolic equilibrium may cause instability. The spacecraft needs to be located in the unstable region around the hyperbolic equilibrium in order to meet the requirements of some missions. Thus, it is necessary to establish a stabilizing strategy for the hyperbolic equilibrium.


### 2.2. Typical Hyperbolic Equilibrium in Astrodynamics

For the typical second-order Hamiltonian system, differentV and ω can be used to define different astrodynamical problems, illustrated in the following section.ω
=
1 and V=(1/2)(x2+y2)+μ/r1+(1-μ)/r2 define the circular restricted three-body problem (CR3BP) [8], where μ is the mass ratio between the primaries, r1 and r2, are, respectively, the distances between the spacecraft and the two primaries, that is, m1 and m2. There are five libration points, three collinear ones of which are hyperbolic, labeled as L1,L2, and L3 shown in Figure 1.Figure 1
Hyperbolic equilibria in circular restricted three-body problem: The three collinear equilibria are hyperbolic, labeled asL1,L2, and L3. The equilibria in CR3BP are renamed as libration points.ω
=
1 and V=(1/2)(x2+y2)+μ/r1+((1-μ)/r2)+a·(nx·x+ny·y) define the planar solar sail three body problem [4], where a is the solar pressure acceleration and n=[nxny]T is the sail surface normal vector. There are infinite equilibria existing in this system, which can be parameterized by the sail lightness number β:
(5)β=r141-μ·∥-∇V∥3(r1·-∇V)2,
where β is the sail lightness number which is the ratio of the radiation pressure acceleration to gravitational acceleration. This parameterization for Sun-Earth/Moon system generates level surfaces as demonstrated in Figure 2.Figure 2
Section of the level surfaces in Sun-Earth/Moon system.β0=0.5059.ω
=
0 and U=hz2/2ρ2-1/r-κ·z define the solar sail’s two-body problem [1], where κ is the constant determined by solar pressure; the polar coordinates (z,ρ) are parallel to the sun line and the minimal radius from the z axis, respectively; hz=hzmax is the constant angular momentum directed along the sun line; r=ρ2+z2. ω=0 indicates that this Hamiltonian system is degenerated for the absence of the Coriolis acceleration and a unique double equilibrium due to hz=hzmax.Considering the formation flying on aJ2-perturbed mean circular orbit, the linearized relative dynamics can be deduced as [6, 9]
(6)Δr¨+2F-1F˙Δr˙+F-1F¨Δr=-(∇Δr·∇)UF|Δr=0·Δr,
where Δr=[xyz]T is the relative position in the chief’s LVLH rotating frame, U is the gravitational potential function of the deputy including J2 perturbation, F is the coordinate transformation matrix from the LVLH to inertial frames, and ∇ is the gradient vector described in orthogonal or spherical coordinates. Due to the osculating orbital elements, the relative dynamics on a J2-perturbed mean circular orbit is time-periodic Hamiltonian system. This time-periodic system has the fixed equilibrium Δr=[0,0,0]T, that is, the location of chief; however, the equilibrium of the system has time-varying topological types and no fixed dimensional unstable/stable/center manifolds, which are quite different from the two-dimensional time-independent system with a permanent pair of hyperbolic eigenvalues and fixed dimensions of unstable/stable/center manifolds, as shown in Figure 3.Topological type of the equilibrium in the modified uncontrolled system: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period.
(a)
(b)Some averaging techniques can be used to reduce the typical astrodynamics into a first-order Hamiltonian system. For a spacecraft with high area-to-mass ratio orbiting the Earth, its dynamics is strongly perturbed by the term of the gravitational field due to the Earth’s oblateness and by the effect of solar radiation pressure (SRP). The secular rate of the orbital elements due to SRP andJ2 term can be written as [10]
(7)e˙=-Csinϕ1-e2ϕ˙=-Ccosϕ1-e2e+W(1-e2)2-1,
where C and W are the constant parameters measuring the solar radiation pressure and the oblateness. The reduced system exists three equilibria, that is, (e1,0),(e2,0) and (e3,π)(e1<e2<e3), where (e2,0) is hyperbolic and the others are elliptic, as shown in Figure 4.Figure 4
Hyperbolic regime in the eccentricity-ϕ phase space evolution under the effects of J2+SRP: the bold lines represent the separatrices in correspondence of the saddle point and zero eccentricity Hamiltonian function.


## 2.1. Typical Hamiltonian System in Astrodynamics

For the typical Hamiltonian system in astrodynamics, the Hamiltonian function has the following form as [4]:
(1)H=12pTp+ωpTJq-V(q,t),
where V is the general form of pseudopotential function relying only on the coordinate q and time; ω is the characteristic parameter introduced by the Coriolis acceleration, such as the mean angular motion of the synodic frame used in CR3BP, and Local Vertical-Local Horizontal (LVLH) frame used in formation flying.In general, the variables(q,p) can be expressed by the physical variables (r,r˙) from the Legendre transformation:
(2)[qp]=[I0-ωJI][rr˙],
where I is the identical matrix, and J is the symplectic matrix. The astrodynamics in physical variables can be expressed as
(3)x¨-2ωy˙=∂V∂x,y¨+2ωx˙=∂V∂y.
For a hyperbolic Hamiltonian system, there exist hyperbolic equilibria that have stable, unstable, and center manifolds [1], which are similar to the collinear libration points L1, L2, and L3 in CR3BP. Hence, the equilibrium of the typical 2-dimensional Hamiltonian system in astrodynamics can be obtained from the following equations:
(4)r=r0,r˙=0∂V∂r|r0=0.
Denote Vrr as the second derivative matrix of the pseudopotential function V about the position vector r=[x,y], and the elements of the matrix have the general expression as Vmn=∂2V/∂m∂n,(m,n)=(x,y). And then a hyperbolic (saddle) equilibrium results from Vxx·Vyy-Vxy2<0.Generally, unstable manifolds of the hyperbolic equilibrium may cause instability. The spacecraft needs to be located in the unstable region around the hyperbolic equilibrium in order to meet the requirements of some missions. Thus, it is necessary to establish a stabilizing strategy for the hyperbolic equilibrium.


## 2.2. Typical Hyperbolic Equilibrium in Astrodynamics

For the typical second-order Hamiltonian system, differentV and ω can be used to define different astrodynamical problems, illustrated in the following section.ω
=
1 and V=(1/2)(x2+y2)+μ/r1+(1-μ)/r2 define the circular restricted three-body problem (CR3BP) [8], where μ is the mass ratio between the primaries, r1 and r2, are, respectively, the distances between the spacecraft and the two primaries, that is, m1 and m2. There are five libration points, three collinear ones of which are hyperbolic, labeled as L1,L2, and L3 shown in Figure 1.Figure 1
Hyperbolic equilibria in circular restricted three-body problem: The three collinear equilibria are hyperbolic, labeled asL1,L2, and L3. The equilibria in CR3BP are renamed as libration points.ω
=
1 and V=(1/2)(x2+y2)+μ/r1+((1-μ)/r2)+a·(nx·x+ny·y) define the planar solar sail three body problem [4], where a is the solar pressure acceleration and n=[nxny]T is the sail surface normal vector. There are infinite equilibria existing in this system, which can be parameterized by the sail lightness number β:
(5)β=r141-μ·∥-∇V∥3(r1·-∇V)2,
where β is the sail lightness number which is the ratio of the radiation pressure acceleration to gravitational acceleration. This parameterization for Sun-Earth/Moon system generates level surfaces as demonstrated in Figure 2.Figure 2
Section of the level surfaces in Sun-Earth/Moon system.β0=0.5059.ω
=
0 and U=hz2/2ρ2-1/r-κ·z define the solar sail’s two-body problem [1], where κ is the constant determined by solar pressure; the polar coordinates (z,ρ) are parallel to the sun line and the minimal radius from the z axis, respectively; hz=hzmax is the constant angular momentum directed along the sun line; r=ρ2+z2. ω=0 indicates that this Hamiltonian system is degenerated for the absence of the Coriolis acceleration and a unique double equilibrium due to hz=hzmax.Considering the formation flying on aJ2-perturbed mean circular orbit, the linearized relative dynamics can be deduced as [6, 9]
(6)Δr¨+2F-1F˙Δr˙+F-1F¨Δr=-(∇Δr·∇)UF|Δr=0·Δr,
where Δr=[xyz]T is the relative position in the chief’s LVLH rotating frame, U is the gravitational potential function of the deputy including J2 perturbation, F is the coordinate transformation matrix from the LVLH to inertial frames, and ∇ is the gradient vector described in orthogonal or spherical coordinates. Due to the osculating orbital elements, the relative dynamics on a J2-perturbed mean circular orbit is time-periodic Hamiltonian system. This time-periodic system has the fixed equilibrium Δr=[0,0,0]T, that is, the location of chief; however, the equilibrium of the system has time-varying topological types and no fixed dimensional unstable/stable/center manifolds, which are quite different from the two-dimensional time-independent system with a permanent pair of hyperbolic eigenvalues and fixed dimensions of unstable/stable/center manifolds, as shown in Figure 3.Topological type of the equilibrium in the modified uncontrolled system: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period.
(a)
(b)Some averaging techniques can be used to reduce the typical astrodynamics into a first-order Hamiltonian system. For a spacecraft with high area-to-mass ratio orbiting the Earth, its dynamics is strongly perturbed by the term of the gravitational field due to the Earth’s oblateness and by the effect of solar radiation pressure (SRP). The secular rate of the orbital elements due to SRP andJ2 term can be written as [10]
(7)e˙=-Csinϕ1-e2ϕ˙=-Ccosϕ1-e2e+W(1-e2)2-1,
where C and W are the constant parameters measuring the solar radiation pressure and the oblateness. The reduced system exists three equilibria, that is, (e1,0),(e2,0) and (e3,π)(e1<e2<e3), where (e2,0) is hyperbolic and the others are elliptic, as shown in Figure 4.Figure 4
Hyperbolic regime in the eccentricity-ϕ phase space evolution under the effects of J2+SRP: the bold lines represent the separatrices in correspondence of the saddle point and zero eccentricity Hamiltonian function.


## 3. Structure-Preserving Stabilization for Hamiltonian System

### 3.1. Controller Design Preserving Hamiltonian Structure

For the typical second-order Hamiltonian system, the variation equation around the equilibrium is:(8)δr¨-2ωJδr˙-Urrδr=0.
For this unstable system, there exist hyperbolic eigenvalues ±σ associated with the stable and unstable manifolds u±, and elliptic eigenvalues ±γi associated with the centre manifolds u and u-.The Hamiltonian structure-preserving controller then is constructed as(9)T=-σ2[G1u+u+T+G2u-u-T]-γ2G3[uuT+u-u-T],(10)K=2ΔJ,TC=Tδr+Kδr˙,
where Δ is used to change the Coriolis acceleration, G1,G2 and G3 are, respectively, the gains of unstable, stable, and centre manifolds. The symmetry of matrix T and the skew symmetry of matrix K guarantee the linear feedback controller preserves the Hamiltonian structure [3].


### 3.2. Stability Analysis

For the 2-dimensional and time-independent Hamiltonian system, some theorem and proposition about the stability of the controller can be deduced as follows.Theorem 1.
The poles can be assigned at any different positions in imaginary axis, andG1,G2, and G3 required are not unique.Proof.
Define the matrix as(11)A=[σ21+u+2σ21+u-22γ21+uu-σ2u+21+u+2σ2u-21+u-22γ2uu-1+uu-σ2u+1+u+2σ2u-1+u-2γ2(u+u-)1+uu-],
where u+ and u- are composing the stable and unstable manifolds of the hyperbolic eigenvalues ±σ as
(12)u±=11+u±2[1u±]
and u and u- are composing the center manifolds of the elliptic eigenvalues ±γi as
(13)u=11+uu-[1u].
Thus Vxx,Vyy, and Vxy will be transformed by the controller:
(14)[V~xxV~yyV~xy]=[VxxVyyVxy]-A[G1G2G3],
where V~xx,V~yy, and V~xy can be assigned arbitrarily to satisfy the linearly stable conditions, so long as matrix A is not regular. Obviously, A cannot be regular due to the fact that u+ is different from u- for the different eigenvalues +σ and -σ.
Denoteλ12 and λ22 as the solutions to the characteristic equation of the second-order dynamics, and then B and C can be expressed as
(15)B=-(λ12+λ22),C=λ12·λ22.
Consequently, V~xx,V~yy, and V~xy can be solved from B=4ω2-Vxx-Vyy and C=Vxx·Vyy-Vxy2.
Furthermore,λ12 and λ22 are restricted by
(16)0>λ12≥-12B≥λ22,
but if one sets
(17)λ12=λ22=-12B
then (λ+iB/2)2(λ-iB/2)2=0. Since this system has second-order elementary factors, the Jordan form of [0IVrr2ωJ] will have long-term dispersions because the matrix cannot be diagonalized.
Because of the nonuniqueness ofV~xx,V~yy, and V~xy solved from (9), G1,G2, and G3 are nonunique. Suppose that sets of V~xx,V~yy,V~xy and V-xx,V-yy,V-xy generate the same values of B and C in (9), and the two sets have the following relationships as
(18)V-xx+V-yy=V~xx+V~yy,V-xx·V-yy-V-xy2=V~xx·V~yy-V~xy2.
It is obtained from (7) that
(19)V-xx2-(V~xx+V~yy)V-xx+(V-xy2+V~xxV~yy-V~xy2)=0
then the real solution is:
(20)(V~xx-V~yy)2+4(V~xy2-V-xy2)≥0
we can fix V~xy2<V-xy2 to get two different sets generated from B and C. Besides, the fact that there exist different sets generating the same values of B and C indicate G1,G2, and G3 are nonunique.Hence, the controlled Hamiltonian can be expressed as(21)H=12pTp+ω~pTJq+12ω~2qTq-V(q)-12δqTTδq.
Expand H around the equilibrium region to get
(22)H=H0+12[δqTδpT][ω2I-Vrr-ω~Jω~JI][δqδp]-12δqTTδq+O(3),
where H0 represents the Hamiltonian value at the equilibrium and the first-order polynomial disappears because of the equilibrium.Proposition 2.
The Coriolis acceleration can be modified byω~=ω+Δ, which cannot stabilize the system independently for 2Δ·Jδr˙ does not change Uxx,Uyy, and Uxy.The controller developed here can transform the hyperbolic equilibrium (saddle) to an elliptic one (center) according to the theorems and propositions. The elliptic equilibrium has a linear symplectic transformation, which transformsH to the following form:
(23)H~(q~,p~)=H0+12λ1(q~12+p~12)+12λ2(q~22+p~22)+O(3).
According to Morse lemma, there exists an analytic diffeomorphism around the equilibrium, which transforms H~ to the following form:
(24)K(q-,p-)=K(q-0,p-0)+q-Tq-+p-Tp-;
thus the energy surfaces are locally diffeomorphic to a family of spheres which shrink down to the equilibria: K→K(q-0,p-0). The Lyapunov stability follows because the trajectories are tangent to the energy surfaces [11]. Therefore, we have obtained that the nonlinear full dynamics is also Lyapunov stable.Theorem 3.
The gains can be chosen large enough to guarantee the modified elliptic equilibrium is the only equilibrium for the controlled system.Proof.
Since it is easy to demonstrate thatT is negative definite, the geometric structure of the Hill’s region for the controlled system will be changed by adding -(1/2)(q-q0)T·T·(q-q0) to the Hill’s region which depends on V.
Another equilibriumLc of the controlled system (if exist) will stay in Island or Mainland [1] because it must stay in the region where ∂xV<0 and ∂yV<0 for the negative definition of T.
Here it is just illustrated that the equilibriumLc cannot stay in Mainland, and the negative inclusion for Island can be obtained by the similar way. We can verify whether ∂rV=0 is true by increasing x and y, since Lc is the bottom point of Mainland; moreover, we need to test which one descends faster between (1/2)(q-q0)T·T·(q-q0) and V. The comparison is very tedious but simply, thus it will be ignored here.
For the controlled system, there exists only the modified elliptic equilibrium, which means the equilibrium has the global minimum pseudopotential. So the controller’s stabilizing region is the whole phase space.Furthermore, the equilibrium of time-periodic system may be unstable even if the equilibrium is always elliptic during its period, compared with the independent Hamiltonian system. Consequently, it is necessary to verify if the modules of all the multipliers are 1, or if the Floquet multipliers lie on the unit circle in the complex plane. Considering the relative dynamics on aJ2-perturbed mean circular orbit, the Floquet stability of the periodic Hamiltonian system is obtained by the larger controller gain rather than the critical gain (G1=G2=G,G3=Δ=0), as shown in Figure 5.Floquet multipliers of the Periodic Hamiltonian System: (a) the relationship between the moduli of Floquet multipliers and the controller gainG; (b) the Floquet multipliers and the unit circle in the complex plane.
(a)
(b)


### 3.3. Resonant Condition

A resonance happens when the eigenvalues of a controlled system have relationλ1:λ2=m:n, where m<n and m and n are reduced with each other. Thus the trajectories will be periodic for any initial condition of spacecraft. We will give a particular solution for the resonance in the following section.Supposeλ1:λ2=m:n, and then:
(25)m2(B+B2-4C)=n2(B-B2-4C)
hence,
(26)mn(4ω~2-V~xx-V~yy)=(m2+n2)V~xxV~yy-V~xy2.
A particular solution can be constructed as
(27)V~xx=V~yy=ω~2,V~xy=±ω~2·n2-m2m2+n2.
If ω~=0, a particular solution can be constructed as
(28)V~xx=V~yy=-ω^2V~xy=±ω^2·n2-m2m2+n2,
where ω^ is an arbitrary parameter.Compared with the osculating frequencies of the controlled elliptic equilibrium for the global investigation during its period, the characteristic frequencies defined by the Floquet multipliers have more applications for the time-periodic Hamiltonian system. Floquet multiplierse±θki(k=1,2,3) are the eigenvalues of the linear Poincaré mapping [Δr,Δr˙]T(t)→[Δr,Δr˙]T(t+Td) with the unit eigenvectors ν1ν-1,ν2ν-2, and ν3ν-3. Characteristic frequencies θ2 and θ3 are approximately equal to zeros inheriting from two of osculating frequencies around ϖ, while θ1 depends on the controller gain. The linear combinations of eigenvectors, (ν1+ν-1),(ν1-ν-1)i,(ν2+ν-2),(ν2-ν-2)i,(ν3+ν-3), and (ν3-ν-3)i, are utilized to generate six different types of near-periodic relative orbits with their periods of (2π/θ1)Td and Td, respectively. The six different types of bounded relative orbits generated by the linear combinations can be then regarded as the foundational motions around the equilibrium.Particularly, the linear combinations(ν1+ν-1) and (ν1-ν-1)i can be used as the initial values to generate the periodic relative orbits with their periods of N·Td when θ1 satisfies the resonant condition N·θ1=2π. To be stricter, all the relative periodic orbits, which produced by the judicious choice combined of the eigenvectors, are near periodic only for the nonlinear term that ignored by the linearized differential equations describing the relative dynamics.It is worth mentioning that the six eigenvectors of the Floquet multipliers span the entire space of the relative position and velocity[ΔrT,Δr˙T]T. Therefore, any initial relative positions and velocity vectors can be decomposed as [ΔrT,Δr˙T]T=∑i=16αiνi, where αi is constant coefficients, and νi is chosen among (ν1+ν-1),(ν1-ν-1)i,(ν2+ν-2),(ν2-ν-2)i,(ν3+ν-3), and (ν3-ν-3)i. Hence, all the initial relative position and velocity can generate a bounded trajectory around the controlled elliptic equilibrium, which depends on the topological type of the equilibrium. Moreover, all the general trajectories are quasiperiodic for the different frequencies between the foundational motions and all the bounded trajectories are involved in the foundational motions caused by the eigenvectors.


### 3.4. Cost and Optimization

BecauseG1,G2,G3, and Δ are nonunique, it is necessary to investigate the rule of choosing G1,G2,G3, and Δ to stabilize the system. Denote U-λ1,λ2 as the collection of G1,G2,G3, and Δ which will allocate the controlled system with the expected poles λ1 and λ2.The controller’s output can be then determined as(29)TC=T·δr+K·δr˙=[T00K][δrδr˙].
Define z=[δrTδr˙T]T to specify the sensitivity of the controller as κ=∥TC∥2/∥z∥2.The Frobenius norm which is consistent to the Euclidean norm of the vector, can be used to measure the controller’s sensitivity as(30)κ=∥TC∥2∥z∥2≤∥[T00K]∥F=tr(TTT)+8Δ2.Define(31)Gm=[G1G2G3G1G2G3G1G2G3G1G2G3]T,Am=[A11A12A13A31A32A33A31A32A33A21A22A23]
and then deduce the matrix T from Gm and Am as(32)T=[G1·σ21+u+2+G2·σ21+u-2+G3·2γ21+uu-G1·σ2u+1+u+2+G2·σ2u-1+u-2+G3·γ2(u+u-)1+uu-G1·σ2u+1+u+2+G2·σ2u-1+u-2+G3·γ2(u+u-)1+uu-G1·σ2u+21+u+2+G2·σ2u-21+u-2+G3·2γ2uu-1+uu-]=Am·Gmwhich has the following relationship as:(33)∥T∥F2=[tr(TTT)]2≤∥Am∥F2·∥Gm∥F2=∥Am∥F2·4(G12+G22+G32).
Thus, put it into (33) and (30) to obtain:
(34)κ≤∥Am∥F2·4(G12+G22+G32)+8Δ2,
where ∥Am∥F can be calculated directly from the constant matrix Am according to the invariant manifolds of the equilibrium.For full feedback,L=∥Am∥F2·4(G12+G22+G32)+8Δ2 can be treated as the optimization index to select more suitable values for the manifolds gains; however, for the position feedback, the optimization index can be set as L=G12+G22+G32. With the constraint equation (15), the selection for controller’s gains can be considered as a nonlinear constrained optimization, as:
(35)min(G1,G2,G3,Δ)∈U-λ1,λ2L,s.t.B+(λ12+λ22)=0,C-λ12·λ22=0.
The programming functions “fmincon” and “confuneq” in MATLAB can solve this optimization problem.As confirmed by zero real part of the eigenvalues in Figure6, the equilibrium is elliptic during the period for time-periodic Hamiltonian system. Therefore, there exist three osculating imaginary eigenvalues ±ω1(t)i,±ω2(t)i, and ±ω3(t)i for the controlled equilibrium. The quasiperiodic relative trajectories around the equilibrium comprise three different oscillating motions with the osculating frequencies of ω1,ω2, and ω3, which can be verified by the second-order term in the new Hamiltonian function yielded by a linear symplectic transformation:
(36)H~(q~,p~,t)=H~0(0,0,t)+12ω1(t)(q~12+p~12)+12ω2(t)(q~22+p~22)+12ω3(t)(q~32+p~32)+O(3).Topological type of the uncontrolled equilibrium: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period.
(a)
(b)According to the local optimization, less cost is obtained by small controller gain illustrated by (34); however, the smaller gain cannot maintain the relative position Δr in a domain closer to the equilibrium, which may consume more fuels.According to the global optimization, an averaging quadratic cost function is defined to measure the fuel consumption during[0,tf] for the continuous thrust engine [12]:
(37)J=12tf∫0tfTcTTcdt=12tf∫0tf[ΔrTGTCT(t)C(t)GΔr]dt=12tf∑i=16∑j=16αiαj∫0tf[ζi(t)TGTCT(t)C(t)Gζi(t)]dt,
where ζi(t) is the position component of the periodic orbit which is developed by the eigenvector νi. The component in the ith row and jth column of the measuring matrix M(G) is Mij=(1/tf)∫0tf[ζi(t)TGTCT(t)C(t)Gζi(t)]dt, which inherits the near-periodicity from the characteristic frequencies and depends only on the controlled elliptic equilibrium. Therefore, M inherited from the longest period of the characteristic frequency θ1, that is, tf=(2π/θ1)Td, can serve as the foundational measuring matrix that can be obtained offline (or on the ground) from the linearized dynamics for the fast and approximative estimation. The consumption function during [0,(2π/θ1)Td] has the form as:
(38)J(G,α)=12αTM(G)α,
where M is symmetric and positive-definite and α=[α1α2α3α4α5α6]T. Generally, there exists a unitary matrix U from Singular Value Decomposition (SVD) to transform the function as J=(1/2)∑i=16λiβi2, where λi represents the eigenvalues of M and βi is the ith column element of the intermediate variable β=Uα.The identified gain can be derived to minimizeJ if we know the initial relative position and velocity. The specific gain satisfying the minimal trace of M(G), that is, trM=∑i=16λi, is however preferable to minimize J for the maximum likelihood, for an unknown case on the initial conditions.


## 3.1. Controller Design Preserving Hamiltonian Structure

For the typical second-order Hamiltonian system, the variation equation around the equilibrium is:(8)δr¨-2ωJδr˙-Urrδr=0.
For this unstable system, there exist hyperbolic eigenvalues ±σ associated with the stable and unstable manifolds u±, and elliptic eigenvalues ±γi associated with the centre manifolds u and u-.The Hamiltonian structure-preserving controller then is constructed as(9)T=-σ2[G1u+u+T+G2u-u-T]-γ2G3[uuT+u-u-T],(10)K=2ΔJ,TC=Tδr+Kδr˙,
where Δ is used to change the Coriolis acceleration, G1,G2 and G3 are, respectively, the gains of unstable, stable, and centre manifolds. The symmetry of matrix T and the skew symmetry of matrix K guarantee the linear feedback controller preserves the Hamiltonian structure [3].


## 3.2. Stability Analysis

For the 2-dimensional and time-independent Hamiltonian system, some theorem and proposition about the stability of the controller can be deduced as follows.Theorem 1.
The poles can be assigned at any different positions in imaginary axis, andG1,G2, and G3 required are not unique.Proof.
Define the matrix as(11)A=[σ21+u+2σ21+u-22γ21+uu-σ2u+21+u+2σ2u-21+u-22γ2uu-1+uu-σ2u+1+u+2σ2u-1+u-2γ2(u+u-)1+uu-],
where u+ and u- are composing the stable and unstable manifolds of the hyperbolic eigenvalues ±σ as
(12)u±=11+u±2[1u±]
and u and u- are composing the center manifolds of the elliptic eigenvalues ±γi as
(13)u=11+uu-[1u].
Thus Vxx,Vyy, and Vxy will be transformed by the controller:
(14)[V~xxV~yyV~xy]=[VxxVyyVxy]-A[G1G2G3],
where V~xx,V~yy, and V~xy can be assigned arbitrarily to satisfy the linearly stable conditions, so long as matrix A is not regular. Obviously, A cannot be regular due to the fact that u+ is different from u- for the different eigenvalues +σ and -σ.
Denoteλ12 and λ22 as the solutions to the characteristic equation of the second-order dynamics, and then B and C can be expressed as
(15)B=-(λ12+λ22),C=λ12·λ22.
Consequently, V~xx,V~yy, and V~xy can be solved from B=4ω2-Vxx-Vyy and C=Vxx·Vyy-Vxy2.
Furthermore,λ12 and λ22 are restricted by
(16)0>λ12≥-12B≥λ22,
but if one sets
(17)λ12=λ22=-12B
then (λ+iB/2)2(λ-iB/2)2=0. Since this system has second-order elementary factors, the Jordan form of [0IVrr2ωJ] will have long-term dispersions because the matrix cannot be diagonalized.
Because of the nonuniqueness ofV~xx,V~yy, and V~xy solved from (9), G1,G2, and G3 are nonunique. Suppose that sets of V~xx,V~yy,V~xy and V-xx,V-yy,V-xy generate the same values of B and C in (9), and the two sets have the following relationships as
(18)V-xx+V-yy=V~xx+V~yy,V-xx·V-yy-V-xy2=V~xx·V~yy-V~xy2.
It is obtained from (7) that
(19)V-xx2-(V~xx+V~yy)V-xx+(V-xy2+V~xxV~yy-V~xy2)=0
then the real solution is:
(20)(V~xx-V~yy)2+4(V~xy2-V-xy2)≥0
we can fix V~xy2<V-xy2 to get two different sets generated from B and C. Besides, the fact that there exist different sets generating the same values of B and C indicate G1,G2, and G3 are nonunique.Hence, the controlled Hamiltonian can be expressed as(21)H=12pTp+ω~pTJq+12ω~2qTq-V(q)-12δqTTδq.
Expand H around the equilibrium region to get
(22)H=H0+12[δqTδpT][ω2I-Vrr-ω~Jω~JI][δqδp]-12δqTTδq+O(3),
where H0 represents the Hamiltonian value at the equilibrium and the first-order polynomial disappears because of the equilibrium.Proposition 2.
The Coriolis acceleration can be modified byω~=ω+Δ, which cannot stabilize the system independently for 2Δ·Jδr˙ does not change Uxx,Uyy, and Uxy.The controller developed here can transform the hyperbolic equilibrium (saddle) to an elliptic one (center) according to the theorems and propositions. The elliptic equilibrium has a linear symplectic transformation, which transformsH to the following form:
(23)H~(q~,p~)=H0+12λ1(q~12+p~12)+12λ2(q~22+p~22)+O(3).
According to Morse lemma, there exists an analytic diffeomorphism around the equilibrium, which transforms H~ to the following form:
(24)K(q-,p-)=K(q-0,p-0)+q-Tq-+p-Tp-;
thus the energy surfaces are locally diffeomorphic to a family of spheres which shrink down to the equilibria: K→K(q-0,p-0). The Lyapunov stability follows because the trajectories are tangent to the energy surfaces [11]. Therefore, we have obtained that the nonlinear full dynamics is also Lyapunov stable.Theorem 3.
The gains can be chosen large enough to guarantee the modified elliptic equilibrium is the only equilibrium for the controlled system.Proof.
Since it is easy to demonstrate thatT is negative definite, the geometric structure of the Hill’s region for the controlled system will be changed by adding -(1/2)(q-q0)T·T·(q-q0) to the Hill’s region which depends on V.
Another equilibriumLc of the controlled system (if exist) will stay in Island or Mainland [1] because it must stay in the region where ∂xV<0 and ∂yV<0 for the negative definition of T.
Here it is just illustrated that the equilibriumLc cannot stay in Mainland, and the negative inclusion for Island can be obtained by the similar way. We can verify whether ∂rV=0 is true by increasing x and y, since Lc is the bottom point of Mainland; moreover, we need to test which one descends faster between (1/2)(q-q0)T·T·(q-q0) and V. The comparison is very tedious but simply, thus it will be ignored here.
For the controlled system, there exists only the modified elliptic equilibrium, which means the equilibrium has the global minimum pseudopotential. So the controller’s stabilizing region is the whole phase space.Furthermore, the equilibrium of time-periodic system may be unstable even if the equilibrium is always elliptic during its period, compared with the independent Hamiltonian system. Consequently, it is necessary to verify if the modules of all the multipliers are 1, or if the Floquet multipliers lie on the unit circle in the complex plane. Considering the relative dynamics on aJ2-perturbed mean circular orbit, the Floquet stability of the periodic Hamiltonian system is obtained by the larger controller gain rather than the critical gain (G1=G2=G,G3=Δ=0), as shown in Figure 5.Floquet multipliers of the Periodic Hamiltonian System: (a) the relationship between the moduli of Floquet multipliers and the controller gainG; (b) the Floquet multipliers and the unit circle in the complex plane.
(a)
(b)


## 3.3. Resonant Condition

A resonance happens when the eigenvalues of a controlled system have relationλ1:λ2=m:n, where m<n and m and n are reduced with each other. Thus the trajectories will be periodic for any initial condition of spacecraft. We will give a particular solution for the resonance in the following section.Supposeλ1:λ2=m:n, and then:
(25)m2(B+B2-4C)=n2(B-B2-4C)
hence,
(26)mn(4ω~2-V~xx-V~yy)=(m2+n2)V~xxV~yy-V~xy2.
A particular solution can be constructed as
(27)V~xx=V~yy=ω~2,V~xy=±ω~2·n2-m2m2+n2.
If ω~=0, a particular solution can be constructed as
(28)V~xx=V~yy=-ω^2V~xy=±ω^2·n2-m2m2+n2,
where ω^ is an arbitrary parameter.Compared with the osculating frequencies of the controlled elliptic equilibrium for the global investigation during its period, the characteristic frequencies defined by the Floquet multipliers have more applications for the time-periodic Hamiltonian system. Floquet multiplierse±θki(k=1,2,3) are the eigenvalues of the linear Poincaré mapping [Δr,Δr˙]T(t)→[Δr,Δr˙]T(t+Td) with the unit eigenvectors ν1ν-1,ν2ν-2, and ν3ν-3. Characteristic frequencies θ2 and θ3 are approximately equal to zeros inheriting from two of osculating frequencies around ϖ, while θ1 depends on the controller gain. The linear combinations of eigenvectors, (ν1+ν-1),(ν1-ν-1)i,(ν2+ν-2),(ν2-ν-2)i,(ν3+ν-3), and (ν3-ν-3)i, are utilized to generate six different types of near-periodic relative orbits with their periods of (2π/θ1)Td and Td, respectively. The six different types of bounded relative orbits generated by the linear combinations can be then regarded as the foundational motions around the equilibrium.Particularly, the linear combinations(ν1+ν-1) and (ν1-ν-1)i can be used as the initial values to generate the periodic relative orbits with their periods of N·Td when θ1 satisfies the resonant condition N·θ1=2π. To be stricter, all the relative periodic orbits, which produced by the judicious choice combined of the eigenvectors, are near periodic only for the nonlinear term that ignored by the linearized differential equations describing the relative dynamics.It is worth mentioning that the six eigenvectors of the Floquet multipliers span the entire space of the relative position and velocity[ΔrT,Δr˙T]T. Therefore, any initial relative positions and velocity vectors can be decomposed as [ΔrT,Δr˙T]T=∑i=16αiνi, where αi is constant coefficients, and νi is chosen among (ν1+ν-1),(ν1-ν-1)i,(ν2+ν-2),(ν2-ν-2)i,(ν3+ν-3), and (ν3-ν-3)i. Hence, all the initial relative position and velocity can generate a bounded trajectory around the controlled elliptic equilibrium, which depends on the topological type of the equilibrium. Moreover, all the general trajectories are quasiperiodic for the different frequencies between the foundational motions and all the bounded trajectories are involved in the foundational motions caused by the eigenvectors.


## 3.4. Cost and Optimization

BecauseG1,G2,G3, and Δ are nonunique, it is necessary to investigate the rule of choosing G1,G2,G3, and Δ to stabilize the system. Denote U-λ1,λ2 as the collection of G1,G2,G3, and Δ which will allocate the controlled system with the expected poles λ1 and λ2.The controller’s output can be then determined as(29)TC=T·δr+K·δr˙=[T00K][δrδr˙].
Define z=[δrTδr˙T]T to specify the sensitivity of the controller as κ=∥TC∥2/∥z∥2.The Frobenius norm which is consistent to the Euclidean norm of the vector, can be used to measure the controller’s sensitivity as(30)κ=∥TC∥2∥z∥2≤∥[T00K]∥F=tr(TTT)+8Δ2.Define(31)Gm=[G1G2G3G1G2G3G1G2G3G1G2G3]T,Am=[A11A12A13A31A32A33A31A32A33A21A22A23]
and then deduce the matrix T from Gm and Am as(32)T=[G1·σ21+u+2+G2·σ21+u-2+G3·2γ21+uu-G1·σ2u+1+u+2+G2·σ2u-1+u-2+G3·γ2(u+u-)1+uu-G1·σ2u+1+u+2+G2·σ2u-1+u-2+G3·γ2(u+u-)1+uu-G1·σ2u+21+u+2+G2·σ2u-21+u-2+G3·2γ2uu-1+uu-]=Am·Gmwhich has the following relationship as:(33)∥T∥F2=[tr(TTT)]2≤∥Am∥F2·∥Gm∥F2=∥Am∥F2·4(G12+G22+G32).
Thus, put it into (33) and (30) to obtain:
(34)κ≤∥Am∥F2·4(G12+G22+G32)+8Δ2,
where ∥Am∥F can be calculated directly from the constant matrix Am according to the invariant manifolds of the equilibrium.For full feedback,L=∥Am∥F2·4(G12+G22+G32)+8Δ2 can be treated as the optimization index to select more suitable values for the manifolds gains; however, for the position feedback, the optimization index can be set as L=G12+G22+G32. With the constraint equation (15), the selection for controller’s gains can be considered as a nonlinear constrained optimization, as:
(35)min(G1,G2,G3,Δ)∈U-λ1,λ2L,s.t.B+(λ12+λ22)=0,C-λ12·λ22=0.
The programming functions “fmincon” and “confuneq” in MATLAB can solve this optimization problem.As confirmed by zero real part of the eigenvalues in Figure6, the equilibrium is elliptic during the period for time-periodic Hamiltonian system. Therefore, there exist three osculating imaginary eigenvalues ±ω1(t)i,±ω2(t)i, and ±ω3(t)i for the controlled equilibrium. The quasiperiodic relative trajectories around the equilibrium comprise three different oscillating motions with the osculating frequencies of ω1,ω2, and ω3, which can be verified by the second-order term in the new Hamiltonian function yielded by a linear symplectic transformation:
(36)H~(q~,p~,t)=H~0(0,0,t)+12ω1(t)(q~12+p~12)+12ω2(t)(q~22+p~22)+12ω3(t)(q~32+p~32)+O(3).Topological type of the uncontrolled equilibrium: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period.
(a)
(b)According to the local optimization, less cost is obtained by small controller gain illustrated by (34); however, the smaller gain cannot maintain the relative position Δr in a domain closer to the equilibrium, which may consume more fuels.According to the global optimization, an averaging quadratic cost function is defined to measure the fuel consumption during[0,tf] for the continuous thrust engine [12]:
(37)J=12tf∫0tfTcTTcdt=12tf∫0tf[ΔrTGTCT(t)C(t)GΔr]dt=12tf∑i=16∑j=16αiαj∫0tf[ζi(t)TGTCT(t)C(t)Gζi(t)]dt,
where ζi(t) is the position component of the periodic orbit which is developed by the eigenvector νi. The component in the ith row and jth column of the measuring matrix M(G) is Mij=(1/tf)∫0tf[ζi(t)TGTCT(t)C(t)Gζi(t)]dt, which inherits the near-periodicity from the characteristic frequencies and depends only on the controlled elliptic equilibrium. Therefore, M inherited from the longest period of the characteristic frequency θ1, that is, tf=(2π/θ1)Td, can serve as the foundational measuring matrix that can be obtained offline (or on the ground) from the linearized dynamics for the fast and approximative estimation. The consumption function during [0,(2π/θ1)Td] has the form as:
(38)J(G,α)=12αTM(G)α,
where M is symmetric and positive-definite and α=[α1α2α3α4α5α6]T. Generally, there exists a unitary matrix U from Singular Value Decomposition (SVD) to transform the function as J=(1/2)∑i=16λiβi2, where λi represents the eigenvalues of M and βi is the ith column element of the intermediate variable β=Uα.The identified gain can be derived to minimizeJ if we know the initial relative position and velocity. The specific gain satisfying the minimal trace of M(G), that is, trM=∑i=16λi, is however preferable to minimize J for the maximum likelihood, for an unknown case on the initial conditions.


## 4. Application to Nonlinear Astrodynamics

### 4.1. Stable Lissajous Orbits of Solar Sail: Application to Second-Order Hamiltonian System

A solar sail is a new type of spacecraft with no fuel, and its orbital control is realized by the solar radiation pressure by means of orientating its attitude relative to the Sun. Therefore, we can neglect the control consumption and implement the complicated control law, which is incomparable to the low-trust ion propulsion or the impulse propulsion.We noticed the fact that there exists quasiperiodic or periodic orbits near the sail’s equilibria, and several periodic orbits have been proposed in the solar sail restricted three-body problem. A survey [13] aimed to summarize the combined wealth of literature concerned with the dynamics, stability and control of highly non-Keplerian orbits for various low thrust propulsion devices, and to demonstrate some of these potential applications, like Capturing near earth objects [14]. McInnes [15] is the first person who applied the classical libration point orbit theory to generate the halo orbit around the on-axis sail’s equilibrium. Baoyin and McInnes [16] achieved two different halo orbits types around the on-axis sail’s equilibrium. McInnes [15] has tried to generate lissajous orbits near the equilibrium as well, but failed to keep a bounded trajectory for a long time. Waters and McInnes [17] applied the classical theory to generate some halo orbits that is only available for the off-axis sail equilibria in some specified regions. What is more, all the interesting orbits presented by researchers are unstable and have one-dim unstable manifolds, hence the station-keeping strategies becomes quite necessary.The sail associated with the Sun and Earth-Moon system, (where the Earth and Moon are regarded as a whole celestial body, denoted as Earth/Moon), is regarded as the solar sail restricted three-body system. In this section we suppose that the Sun-Earth/Moon system revolves around a circular orbit.The structure-preserving controller is proposed to generate stable lissajous orbits and Lyapunov orbits for nonresonance shown in Figure7, and resonant periodic orbits shown in Figure 8 (1 : 2, 1 : 3, 1 : 4, 2 : 3, resp.).Figure 7
Stable Lissajous orbit and Lyapunov orbits near the sail equilibrium:G1=20,G2=10,G3=10.Figure 8
Resonant periodic orbits near the sail equilibrium.The feedback gains are initially chosen withG1=20,G2=10,G3=10 and Δ=0 in view of the local optimization, and then the gains will be optimized so as to obtain the same poles. The gains refined by optimization are G1=13.0147,G2=13.0147,G3=12.0063, and Δ=0, with the optimization index improved from 600 to 483 (80.5%) for the position feedback. While for full feedback, the gains refined by optimization are G1=11.3253,G2=11.3254,G3=14.4533, and Δ=-0.5956, with the optimization index improved from 2.0058×103 to 1.5587×103 (77.71%).


### 4.2. Bounded Orbits by Degenerated Solar Sail: Application to Degenerated Hamiltonian System

McInnes [18] has researched the displaced non-Keplerian orbits above the Sun. Similar steps have been taken by Bookless and McInnes [19, 20] to the displaced non-Keplerian orbits above the planet. Nevertheless, just as the stability analysis for motions around the equilibrium, Bookless’ analysis for the dynamics and control is linear and local. There presented the nonlinear dynamical analysis for a displaced orbit above a planet, and there also investigated the motions around the equilibria for the nonresonance case with the help of the Birkhoff normal form and dynamical system techniques. M. Xu and S. Xu [1] has obtained two important contributions to derive necessary and sufficiency conditions for the motion stability around the equilibria, and to illustrate numerically that the non-transition critical KAM torus is filled with the Lyapunov (1,1)-homoclinic orbits. hz=hzmax defines the degradation case for the two body problem of solar sails, which is the saddle-node bifurcation point.For the degenerated system, there exists only one double equilibrium, and the frequencies for the liberalized motion are, respectively,ω1=0 and ω2>0. What is more, it is demonstrated from the variation equation around the equilibrium that the equilibrium has no stable or unstable manifolds but centre manifolds. Owing to Theorem 1, it guarantees that the stabilization for the degenerated hzmax can be achieved by the equilibrium centre manifolds.Similar to the nondegradation case, the suitable gains can be chosen to guarantee the modified elliptic equilibrium is the only equilibrium for the controlled system. So the controller’s stabilizing region is all the phase space. The quasiperiodic orbits generated by the controller are shown in Figure9, and their positions in 3D physical space are shown in Figure 10.Figure 9
Quasiperiodic orbits by the controller near the degenerated equilibrium:G1=0,G2=0,G3=0.3,Δ=0.0005.Quasiperiodic orbits by the controller in the 3D physical space:G1=0,G2=0,G3=0.3,Δ=0.0005.
(a)
(b)
(c)
(d)


### 4.3. Quasiperiodic Formation Flying: Application to Time-Periodic Hamiltonian System

Roughly, we can classify the relative motion control approaches into two branches: continuous control that using low-thrust electric propulsion and impulsive control that relying on chemical thrusters. While the former approaches used mostly full-state feedback of the Cartesian relative position and velocity to develop high-accuracy tracking control laws, thus capable to maintain small steady-state errors under a myriad of orbital perturbations [21–25]; however, the maintenance of the later assumes that the tracking reference orbital elements are known and attempts to generate control commands that will match the osculating orbital elements to the expected values [12, 26, 27].The initial relative position ofx0=y0=z0=1km and the relative velocity generated from the classical C-W equations y˙0=-2ϖx0,x˙0=ϖy0/2,z˙0=0 for the formation flying on a J2-perturbed mean circular orbit, are applied to generate the quasiperiodic relative trajectories with the preferred gain 0.31 illustrated in Figure 11. More rapid changes in configuration are obtained by the controller developed in this paper, rather than by the change in J2 invariant relative orbits caused by the secular term of perigee argument, from the projection of the configuration on the x-z or y-z plane. Therefore, the mission that requires rapid changes in formation configuration, like on-orbit inspection and repair could potentially use the addressed controller. Furthermore, the quasiperiodic relative trajectories which comprise six foundational motions with different frequencies have no reference relative configuration to track for the continuous thrust, make the trajectory prediction difficult by other nonallied spacecrafts.Quasiperiodic relative trajectories generated by the optimal gain of 0.31 during 30 days: presented is the relative orbit during the first orbital period in every day; the black thicker line indicates the initial formation configuration.
(a)
(b)
(c)
(d)Figure12 illustrates the topological type of the equilibrium of the controlled system, which is elliptic during the period, as confirmed by zero real part of the eigenvalues, for the mean circular orbit whose mean inclination is 80°, mean semi-major axis is 7178.137 km, and initial argument of latitude is 0°.Topological type of the equilibrium in the controlled system: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period; on the mean circular orbit witha-=7178.137km,i- = 80°, and an initial argument of latitude of 0°.
(a)
(b)G
+ should share its value with G- to minimize the controller outputs from view of local optimization, according to the contour curves of κ illustrated in Figure 13. The specific gain owning the minimal trace of M(G), that is, trM=∑i=16λi, is preferable to minimize J for the maximum likelihood if the initial conditions are unknown. The optimizing result obtained from local optimization is illustrated in Figure 14 that the optimal gain of 0.31 has a smaller trace than the critical gain of 0.27.Figure 13
Stable regions of the controller in the(G+-G-) space: the hyperbolic-like curve C1 is the boundary keeping the topological type of the equilibrium elliptic; the line-like curve C2 is the boundary keeping the Floquet multipliers lying on the unit circle in the complex plane; the closed curves labeled by numbers are contour curves of the averaging values of κ; the circle C3 is the approximation of the contour curve.Figure 14
Trace of the measuring matrix: compared with the trace of the critical gain of 0.27 obtained from local optimization, the gain of 0.31 can achieve the minimal trace.


### 4.4. Controlled Frozen Orbits: Application to First-Order Hamiltonian System

The high area-to-mass spacecraft’s orbital evolution behaves strangely under the influence of the solar radiation pressure and the perturbation due to the oblateness of the Earth. Since 1960, there have been observed long-term oscillations of the eccentricity in the satellites’ orbital behaviours, like Vanguard and the ECHO balloon. The two variables, the osculating orbital eccentricity and the solar angleϕ between the orbit pericentre and the Sun-line [10, 28], can be used to describe this dynamical system with low inclinations through a Hamiltonian function.All spacecrafts around the hyperbolic eccentricity cannot be maintained on such an orbit; on the contrary, they will move on the unstable manifold towards higher eccentricities, or on the unstable manifold towards lower eccentricities based on the quite small change of the orbital elements. Controlled frozen trajectories can be obtained in the phase space around the hyperbolic eccentricity by means of applying the Hamiltonian structure-preserving controller, such as the one illustrated by red dotted line in Figure15.Controlled frozen orbits generated by the controller ine-ϕ space.
(a)
(b)The mission illustrated in Figure15(a) comprises three phases: (i) transferring (blue line) from a circular orbit to the hyperbolic equilibrium, (ii) stabilizing near the hyperbolic equilibrium (red line), and (iii) transferring (dotted blue line) beyond the critical eccentricity. While the mission illustrated in Figure 15(b) comprises three phases as well: (i) transferring from a circular orbit to the hyperbolic equilibrium (blue line), (ii) stabilizing (red line) near the hyperbolic equilibrium, and (iii) transferring (the dotted blue line) back to the circular orbit.


## 4.1. Stable Lissajous Orbits of Solar Sail: Application to Second-Order Hamiltonian System

A solar sail is a new type of spacecraft with no fuel, and its orbital control is realized by the solar radiation pressure by means of orientating its attitude relative to the Sun. Therefore, we can neglect the control consumption and implement the complicated control law, which is incomparable to the low-trust ion propulsion or the impulse propulsion.We noticed the fact that there exists quasiperiodic or periodic orbits near the sail’s equilibria, and several periodic orbits have been proposed in the solar sail restricted three-body problem. A survey [13] aimed to summarize the combined wealth of literature concerned with the dynamics, stability and control of highly non-Keplerian orbits for various low thrust propulsion devices, and to demonstrate some of these potential applications, like Capturing near earth objects [14]. McInnes [15] is the first person who applied the classical libration point orbit theory to generate the halo orbit around the on-axis sail’s equilibrium. Baoyin and McInnes [16] achieved two different halo orbits types around the on-axis sail’s equilibrium. McInnes [15] has tried to generate lissajous orbits near the equilibrium as well, but failed to keep a bounded trajectory for a long time. Waters and McInnes [17] applied the classical theory to generate some halo orbits that is only available for the off-axis sail equilibria in some specified regions. What is more, all the interesting orbits presented by researchers are unstable and have one-dim unstable manifolds, hence the station-keeping strategies becomes quite necessary.The sail associated with the Sun and Earth-Moon system, (where the Earth and Moon are regarded as a whole celestial body, denoted as Earth/Moon), is regarded as the solar sail restricted three-body system. In this section we suppose that the Sun-Earth/Moon system revolves around a circular orbit.The structure-preserving controller is proposed to generate stable lissajous orbits and Lyapunov orbits for nonresonance shown in Figure7, and resonant periodic orbits shown in Figure 8 (1 : 2, 1 : 3, 1 : 4, 2 : 3, resp.).Figure 7
Stable Lissajous orbit and Lyapunov orbits near the sail equilibrium:G1=20,G2=10,G3=10.Figure 8
Resonant periodic orbits near the sail equilibrium.The feedback gains are initially chosen withG1=20,G2=10,G3=10 and Δ=0 in view of the local optimization, and then the gains will be optimized so as to obtain the same poles. The gains refined by optimization are G1=13.0147,G2=13.0147,G3=12.0063, and Δ=0, with the optimization index improved from 600 to 483 (80.5%) for the position feedback. While for full feedback, the gains refined by optimization are G1=11.3253,G2=11.3254,G3=14.4533, and Δ=-0.5956, with the optimization index improved from 2.0058×103 to 1.5587×103 (77.71%).


## 4.2. Bounded Orbits by Degenerated Solar Sail: Application to Degenerated Hamiltonian System

McInnes [18] has researched the displaced non-Keplerian orbits above the Sun. Similar steps have been taken by Bookless and McInnes [19, 20] to the displaced non-Keplerian orbits above the planet. Nevertheless, just as the stability analysis for motions around the equilibrium, Bookless’ analysis for the dynamics and control is linear and local. There presented the nonlinear dynamical analysis for a displaced orbit above a planet, and there also investigated the motions around the equilibria for the nonresonance case with the help of the Birkhoff normal form and dynamical system techniques. M. Xu and S. Xu [1] has obtained two important contributions to derive necessary and sufficiency conditions for the motion stability around the equilibria, and to illustrate numerically that the non-transition critical KAM torus is filled with the Lyapunov (1,1)-homoclinic orbits. hz=hzmax defines the degradation case for the two body problem of solar sails, which is the saddle-node bifurcation point.For the degenerated system, there exists only one double equilibrium, and the frequencies for the liberalized motion are, respectively,ω1=0 and ω2>0. What is more, it is demonstrated from the variation equation around the equilibrium that the equilibrium has no stable or unstable manifolds but centre manifolds. Owing to Theorem 1, it guarantees that the stabilization for the degenerated hzmax can be achieved by the equilibrium centre manifolds.Similar to the nondegradation case, the suitable gains can be chosen to guarantee the modified elliptic equilibrium is the only equilibrium for the controlled system. So the controller’s stabilizing region is all the phase space. The quasiperiodic orbits generated by the controller are shown in Figure9, and their positions in 3D physical space are shown in Figure 10.Figure 9
Quasiperiodic orbits by the controller near the degenerated equilibrium:G1=0,G2=0,G3=0.3,Δ=0.0005.Quasiperiodic orbits by the controller in the 3D physical space:G1=0,G2=0,G3=0.3,Δ=0.0005.
(a)
(b)
(c)
(d)


## 4.3. Quasiperiodic Formation Flying: Application to Time-Periodic Hamiltonian System

Roughly, we can classify the relative motion control approaches into two branches: continuous control that using low-thrust electric propulsion and impulsive control that relying on chemical thrusters. While the former approaches used mostly full-state feedback of the Cartesian relative position and velocity to develop high-accuracy tracking control laws, thus capable to maintain small steady-state errors under a myriad of orbital perturbations [21–25]; however, the maintenance of the later assumes that the tracking reference orbital elements are known and attempts to generate control commands that will match the osculating orbital elements to the expected values [12, 26, 27].The initial relative position ofx0=y0=z0=1km and the relative velocity generated from the classical C-W equations y˙0=-2ϖx0,x˙0=ϖy0/2,z˙0=0 for the formation flying on a J2-perturbed mean circular orbit, are applied to generate the quasiperiodic relative trajectories with the preferred gain 0.31 illustrated in Figure 11. More rapid changes in configuration are obtained by the controller developed in this paper, rather than by the change in J2 invariant relative orbits caused by the secular term of perigee argument, from the projection of the configuration on the x-z or y-z plane. Therefore, the mission that requires rapid changes in formation configuration, like on-orbit inspection and repair could potentially use the addressed controller. Furthermore, the quasiperiodic relative trajectories which comprise six foundational motions with different frequencies have no reference relative configuration to track for the continuous thrust, make the trajectory prediction difficult by other nonallied spacecrafts.Quasiperiodic relative trajectories generated by the optimal gain of 0.31 during 30 days: presented is the relative orbit during the first orbital period in every day; the black thicker line indicates the initial formation configuration.
(a)
(b)
(c)
(d)Figure12 illustrates the topological type of the equilibrium of the controlled system, which is elliptic during the period, as confirmed by zero real part of the eigenvalues, for the mean circular orbit whose mean inclination is 80°, mean semi-major axis is 7178.137 km, and initial argument of latitude is 0°.Topological type of the equilibrium in the controlled system: (a) the time history of real part of the equilibrium during the orbital period; (b) the time history of imaginary part of the equilibrium during the orbital period; on the mean circular orbit witha-=7178.137km,i- = 80°, and an initial argument of latitude of 0°.
(a)
(b)G
+ should share its value with G- to minimize the controller outputs from view of local optimization, according to the contour curves of κ illustrated in Figure 13. The specific gain owning the minimal trace of M(G), that is, trM=∑i=16λi, is preferable to minimize J for the maximum likelihood if the initial conditions are unknown. The optimizing result obtained from local optimization is illustrated in Figure 14 that the optimal gain of 0.31 has a smaller trace than the critical gain of 0.27.Figure 13
Stable regions of the controller in the(G+-G-) space: the hyperbolic-like curve C1 is the boundary keeping the topological type of the equilibrium elliptic; the line-like curve C2 is the boundary keeping the Floquet multipliers lying on the unit circle in the complex plane; the closed curves labeled by numbers are contour curves of the averaging values of κ; the circle C3 is the approximation of the contour curve.Figure 14
Trace of the measuring matrix: compared with the trace of the critical gain of 0.27 obtained from local optimization, the gain of 0.31 can achieve the minimal trace.


## 4.4. Controlled Frozen Orbits: Application to First-Order Hamiltonian System

The high area-to-mass spacecraft’s orbital evolution behaves strangely under the influence of the solar radiation pressure and the perturbation due to the oblateness of the Earth. Since 1960, there have been observed long-term oscillations of the eccentricity in the satellites’ orbital behaviours, like Vanguard and the ECHO balloon. The two variables, the osculating orbital eccentricity and the solar angleϕ between the orbit pericentre and the Sun-line [10, 28], can be used to describe this dynamical system with low inclinations through a Hamiltonian function.All spacecrafts around the hyperbolic eccentricity cannot be maintained on such an orbit; on the contrary, they will move on the unstable manifold towards higher eccentricities, or on the unstable manifold towards lower eccentricities based on the quite small change of the orbital elements. Controlled frozen trajectories can be obtained in the phase space around the hyperbolic eccentricity by means of applying the Hamiltonian structure-preserving controller, such as the one illustrated by red dotted line in Figure15.Controlled frozen orbits generated by the controller ine-ϕ space.
(a)
(b)The mission illustrated in Figure15(a) comprises three phases: (i) transferring (blue line) from a circular orbit to the hyperbolic equilibrium, (ii) stabilizing near the hyperbolic equilibrium (red line), and (iii) transferring (dotted blue line) beyond the critical eccentricity. While the mission illustrated in Figure 15(b) comprises three phases as well: (i) transferring from a circular orbit to the hyperbolic equilibrium (blue line), (ii) stabilizing (red line) near the hyperbolic equilibrium, and (iii) transferring (the dotted blue line) back to the circular orbit.


## 5. Conclusion

Bounded motions near hyperbolic equilibria have been broadly applied to various astronautical missions. The Hamiltonian structure-preserving controller can change the hyperbolic equilibrium (saddle) to an elliptic one (center), and then some quasiperiodic bounded trajectories will emerge naturally from the KAM tori generated by the controller.In this paper, a systematic investigation is presented on the general form for full feedback and position-only feedback modes, Lyapunov and Floquet stability analysis, and control gain selection of the structure-preserving controller, in the context of second-order, first-order, time-periodic, and degenerated Hamiltonian system, respectively. One type of periodic trajectories is achieved by the resonant conditions of control gains for the resonant periodic orbits, and another type is making judicious choice in the foundational motions with different frequencies. Therefore, the controller has potential applications in stable lissajous orbits for solar sail’s three-body problem and degenerated two-body problem, quasiperiodic formation flying on aJ2-perturbed mean circular orbit, or the controlled frozen orbits for a spacecraft with high area-to-mass ratio.Furthermore, some further researches can be implemented, including how to use unfixed dimensional center manifolds to feedback the controller, the further investigation for the global Floquet stable controller gain for quasiperiodic Hamiltonian system, and the relationships between the different frequencies and the critical gain for Floquet stability or the optimal gain.



---

*Source: 107674-2013-03-21.xml*