# Low-Thrust Transfer Design of Low-Observable Geostationary Earth Orbit Satellite

Authors: Bing  Hua; Zhujun Shao
Journal: International Journal of Aerospace Engineering
 (2015)

---

## Abstract

With radar and surface-to-air missiles posing an increasing threat to on-orbit spacecraft, low-observable satellites play an important role in low-thrust transfers. This paper presents the design for a low-thrust geostationary earth orbit (GEO) transfer control strategy which takes into consideration the low-observable constraint and discusses Earth shadow and perturbation. A control parameter optimization addresses the orbit transfer problem, and five thrust modes are used. Simulation results show that the method outlined in this paper is simple and feasible and results in reduced transfer time with a small amount of calculation. The method therefore offers a useful reference for low-thrust GEO transfer design.


---

## Body

## 1. Introduction

Creation of a low-observable satellite is accomplished by low-observable technology, which makes it difficult or impossible to avoid detection of satellites in orbit by hostile enemy forces [1–6]. Adding a low-observability module to the spacecraft’s overall design is key to improving operational effectiveness and satellite survivability.A geostationary satellite is usually launched from ground to Low Earth Orbit (LEO) or Middle Earth Orbit (MEO) rather than to GEO. Since the main threat to LEO and MEO satellites comes from radar, a great deal of emphasis in low-observability design is placed on radar analysis. The radar cross section (RCS) is the designers’ only controllable factor in radar detection. Contemporary work on low observability has its roots in efforts at reducing the RCS of a spacecraft, falling into two categories: low-observable shape design and flight attitude planning [1–8]. The former is a simple way of reducing the range at which radar can detect the spacecraft; however, contouring the surface of a spacecraft to reduce the RCS equally in all directions is not possible. As a result, the latter is typically supplemented in designs in order to produce better results.The application of electric propulsion in geostationary orbit platforms is inevitable for the development of the aerospace. Used for spacecraft applications in Earth orbit, such as station-keeping, orbit-raising, and orbit transfer, Boeing-702sp is an all-electric propulsion satellite with a specific impulse of more than 3800 seconds [9–12]. A great deal of research has been directed toward solving the low-thrust transfer problem. Most of that research has been based on optimal control theory. The numerical solution incorporates both direct and indirect methods, due to their high accuracy [13–16]. As is typical with these methods, the solution is often difficult to derive, requiring a complex initial guess and tedious iteration for convergence.Electric propulsion, low thrust, and highly specific impulses have led to greatly improved fuel efficiency at the expense of relatively long transfer times. For these reasons, satellites can be easily caught by detection systems during low-thrust GEO transfer. Existing research on low-thrust transfer emphasizes optimum techniques, giving little consideration to the satellite attitude constraint of control. This paper models an optimization scheme of the low-thrust transfer based on low-observable technology, which can reduce the visibility of satellites. Furthermore, the optimization includes space environment (Earth shadow and perturbation). Offering the advantages of short time, a small amount of calculation, and no complex initial guess, the method designed in this paper can provide a valuable reference for low-thrust GEO transfer design.The paper is organized in five distinct sections, as represented in Figure1. The introduction above has described low-observable satellites and low-thrust transfer optimization techniques. In Section 2, Earth shadow and perturbation in orbit transfer are analyzed. In Section 3, a mathematical model of low-observable constraint is demonstrated, including the radar detection area and the relationship of thrust control angle and attitude angle. Section 4 presents a method based on control parameter analysis to solve minimum-time transfer. Finally, we summarize the paper.Figure 1
Paper structure.


## 2. Dynamic Model

### 2.1. Modified Equinoctial Orbit Model

The dynamical equations of motion for a thrusting spacecraft can be established by the classical orbital elements: the semimajor axisa, the eccentricity e, the inclination i, the right ascension Ω, the argument of perigee ω, and the mean anomaly M: (1)a˙=2a2esin⁡vhfr+2a2phrfte˙=1hpsin⁡vfr+1hp+rcos⁡v+refti˙=rcos⁡ω+vhfnΩ˙=rsin⁡ω+vhsin⁡ifnω˙=-pcos⁡vhefr+p+rsin⁡vheft+rsin⁡vω+vcos⁡ihsin⁡ifnM˙=n+1aehpcos⁡v-2erfr-p+rsin⁡vft,where p=a(1-e2),h=μp, n=μ/a2, r=p/1+ecos⁡v with Earth’s gravitational coefficient μ and the true anomaly v. Performing analyses of transfer orbits using classical orbital elements is a straightforward task, but singularities are exhibited for zero eccentricity and inclinations of 0° and 90°.To eliminate these deficiencies, a modified set of equinoctial orbit elements (p,ex,ey,hy,hx, and L) is frequently used. P is the semilatus rectum, (ex,ey) is the eccentricity vector, (hx,hy) is the inclination vector, and L is the true longitude. The relationship between the modified equinoctial elements and the classical orbital elements is given by(2)p=a1-e2ex=ecos⁡ω+Ωey=esin⁡ω+Ωhx=tan⁡i2cos⁡Ωhy=tani2sin⁡ΩL=ω+Ω+v.The equations of motion, written in terms of the modified equinoctial elements, are(3)p˙=pμ2pωfte˙x=pμfrsin⁡L+1+ωcos⁡L+exftω-hxsin⁡L-hycos⁡Leyfnωe˙y=pμ-frcos⁡L+1+ωsin⁡L+eyftω+hxsin⁡L-hycos⁡Lexfnωh˙x=pμs2fn2ωcos⁡Lh˙y=pμs2fn2ωsin⁡LL˙=μpωp2+1ωpμhxsin⁡L-hycos⁡Lfnwith ω=1+excos⁡L+eysin⁡L, s2=1+hx2+hy2. The acceleration components fr, ft, and fn are denoted in the RTN frame (with fr being the normal direction, ft being the tangential direction, and fn being the direction orthogonal to the orbit plane), and f is the thrust acceleration. The constraint on the control is f≤fmax=Tmax/m.To take into account the true acceleration of the thrust, we must consider the mass flow. Its evolution is given by(4)dmdt=-TIspg0,where T is the thrust modulus, Isp the specific impulse of thruster, and g0 the gravitational acceleration at sea-level and g0=9.80665N/s.As shown in Figure2, thrust acceleration components can also be defined by acceleration f and two control angles (pitch-steering angle α and yaw-steering angle β), which describe the direction of the thrust vector in relation to the velocity vector and the orbital plane, respectively. The acceleration components are expressed by(5)fr=fcos⁡βsin⁡αft=fcos⁡βcos⁡αfn=fsin⁡βwith -π≤α≤π and -π/2≤β≤π/2.Figure 2
Coordinate frame.


### 2.2. Eclipse Effects

To model the electric propulsion system accurately, it is necessary to model the satellite’s trajectory as it passes through the shadow of Earth. In particular, the satellite is at discharge, which leads to high power consumption when in the shadow. In order to ensure the safety of spacecraft, therefore, the thrust is 0 when the vehicle is in the shadow of Earth. A simplified Earth-shadow model, namely, the cylindrical projection model shown in Figure3, is used to estimate the location of the shadow [17–20].Figure 3
Cylindrical projection model.Referring to the geometry illustrated in Figure3, rs is defined as the vector from Earth to the Sun, with norm rs, and r as the vector from Earth to the satellite, with norm r. Re is the radius of Earth, and Rs is the radius of the Sun.The cone angleθ is defined by(6)cos⁡θ=∠rsr=r·rsrrs.The projected spacecraft positionr is given by(7)r=rcos⁡vcos⁡Ωcos⁡ω-cos⁡vsin⁡Ωsin⁡ωcos⁡i-sin⁡vcos⁡Ωcos⁡ω-sin⁡vsin⁡Ωsin⁡ωcos⁡icos⁡vsin⁡Ωcos⁡ω+cos⁡vcos⁡Ωsin⁡ωcos⁡i-sin⁡vsin⁡Ωcos⁡ω+sin⁡vcos⁡Ωsin⁡ωcos⁡icos⁡vsin⁡ωsin⁡i+sin⁡vcos⁡ωsin⁡i.The shadow entry and exit locations are judged by locating the cone terminators at the projected spacecraft location. The shadow can only be found when the angle satisfies(8)cos⁡θ≤0sin⁡θ<Rer.


### 2.3. Perturbation

Earth’s oblateness, atmospheric drag, light pressure, secondary body, and other factors in space can also perturb a satellite’s motion. Among these variables, Earth’s oblateness is of vital importance for predicting the trajectory of the satellite accurately. The oblate Earth perturbation is caused by the reality of Earth’s shape not being perfectly spherical. The impact of Earth’s oblateness due toJ2 is always taken into account in the engineering calculations [20, 21].In terms of orbital elements, the dynamical system is(9)dadt=0dedt=0didt=0dΩdt=-3nJ2Re221-e2a2cos⁡idωdt=3nJ2Re241-e2a25cos2⁡i-1dMdt=n-3nJ2Re241-e23/2a21-3cos2⁡i.The effect of Earth’s oblateness (due toJ2) on the orbital transfer can be included by appending the perturbation to the respective right-hand sides of (5):(10)fr=32J2ur4RE23sin2⁡isin2⁡u-1ft=-32J2ur4RE2sin2⁡isin⁡2ufn=-J2ur4RE2sin⁡2isin⁡u.


## 2.1. Modified Equinoctial Orbit Model

The dynamical equations of motion for a thrusting spacecraft can be established by the classical orbital elements: the semimajor axisa, the eccentricity e, the inclination i, the right ascension Ω, the argument of perigee ω, and the mean anomaly M: (1)a˙=2a2esin⁡vhfr+2a2phrfte˙=1hpsin⁡vfr+1hp+rcos⁡v+refti˙=rcos⁡ω+vhfnΩ˙=rsin⁡ω+vhsin⁡ifnω˙=-pcos⁡vhefr+p+rsin⁡vheft+rsin⁡vω+vcos⁡ihsin⁡ifnM˙=n+1aehpcos⁡v-2erfr-p+rsin⁡vft,where p=a(1-e2),h=μp, n=μ/a2, r=p/1+ecos⁡v with Earth’s gravitational coefficient μ and the true anomaly v. Performing analyses of transfer orbits using classical orbital elements is a straightforward task, but singularities are exhibited for zero eccentricity and inclinations of 0° and 90°.To eliminate these deficiencies, a modified set of equinoctial orbit elements (p,ex,ey,hy,hx, and L) is frequently used. P is the semilatus rectum, (ex,ey) is the eccentricity vector, (hx,hy) is the inclination vector, and L is the true longitude. The relationship between the modified equinoctial elements and the classical orbital elements is given by(2)p=a1-e2ex=ecos⁡ω+Ωey=esin⁡ω+Ωhx=tan⁡i2cos⁡Ωhy=tani2sin⁡ΩL=ω+Ω+v.The equations of motion, written in terms of the modified equinoctial elements, are(3)p˙=pμ2pωfte˙x=pμfrsin⁡L+1+ωcos⁡L+exftω-hxsin⁡L-hycos⁡Leyfnωe˙y=pμ-frcos⁡L+1+ωsin⁡L+eyftω+hxsin⁡L-hycos⁡Lexfnωh˙x=pμs2fn2ωcos⁡Lh˙y=pμs2fn2ωsin⁡LL˙=μpωp2+1ωpμhxsin⁡L-hycos⁡Lfnwith ω=1+excos⁡L+eysin⁡L, s2=1+hx2+hy2. The acceleration components fr, ft, and fn are denoted in the RTN frame (with fr being the normal direction, ft being the tangential direction, and fn being the direction orthogonal to the orbit plane), and f is the thrust acceleration. The constraint on the control is f≤fmax=Tmax/m.To take into account the true acceleration of the thrust, we must consider the mass flow. Its evolution is given by(4)dmdt=-TIspg0,where T is the thrust modulus, Isp the specific impulse of thruster, and g0 the gravitational acceleration at sea-level and g0=9.80665N/s.As shown in Figure2, thrust acceleration components can also be defined by acceleration f and two control angles (pitch-steering angle α and yaw-steering angle β), which describe the direction of the thrust vector in relation to the velocity vector and the orbital plane, respectively. The acceleration components are expressed by(5)fr=fcos⁡βsin⁡αft=fcos⁡βcos⁡αfn=fsin⁡βwith -π≤α≤π and -π/2≤β≤π/2.Figure 2
Coordinate frame.


## 2.2. Eclipse Effects

To model the electric propulsion system accurately, it is necessary to model the satellite’s trajectory as it passes through the shadow of Earth. In particular, the satellite is at discharge, which leads to high power consumption when in the shadow. In order to ensure the safety of spacecraft, therefore, the thrust is 0 when the vehicle is in the shadow of Earth. A simplified Earth-shadow model, namely, the cylindrical projection model shown in Figure3, is used to estimate the location of the shadow [17–20].Figure 3
Cylindrical projection model.Referring to the geometry illustrated in Figure3, rs is defined as the vector from Earth to the Sun, with norm rs, and r as the vector from Earth to the satellite, with norm r. Re is the radius of Earth, and Rs is the radius of the Sun.The cone angleθ is defined by(6)cos⁡θ=∠rsr=r·rsrrs.The projected spacecraft positionr is given by(7)r=rcos⁡vcos⁡Ωcos⁡ω-cos⁡vsin⁡Ωsin⁡ωcos⁡i-sin⁡vcos⁡Ωcos⁡ω-sin⁡vsin⁡Ωsin⁡ωcos⁡icos⁡vsin⁡Ωcos⁡ω+cos⁡vcos⁡Ωsin⁡ωcos⁡i-sin⁡vsin⁡Ωcos⁡ω+sin⁡vcos⁡Ωsin⁡ωcos⁡icos⁡vsin⁡ωsin⁡i+sin⁡vcos⁡ωsin⁡i.The shadow entry and exit locations are judged by locating the cone terminators at the projected spacecraft location. The shadow can only be found when the angle satisfies(8)cos⁡θ≤0sin⁡θ<Rer.


## 2.3. Perturbation

Earth’s oblateness, atmospheric drag, light pressure, secondary body, and other factors in space can also perturb a satellite’s motion. Among these variables, Earth’s oblateness is of vital importance for predicting the trajectory of the satellite accurately. The oblate Earth perturbation is caused by the reality of Earth’s shape not being perfectly spherical. The impact of Earth’s oblateness due toJ2 is always taken into account in the engineering calculations [20, 21].In terms of orbital elements, the dynamical system is(9)dadt=0dedt=0didt=0dΩdt=-3nJ2Re221-e2a2cos⁡idωdt=3nJ2Re241-e2a25cos2⁡i-1dMdt=n-3nJ2Re241-e23/2a21-3cos2⁡i.The effect of Earth’s oblateness (due toJ2) on the orbital transfer can be included by appending the perturbation to the respective right-hand sides of (5):(10)fr=32J2ur4RE23sin2⁡isin2⁡u-1ft=-32J2ur4RE2sin2⁡isin⁡2ufn=-J2ur4RE2sin⁡2isin⁡u.


## 3. Satellite Low Observability

The research activity presented here is focused on the optimum algorithm and does not take into account space environment analysis, especially the low-observable constraint. The low-observable satellite is designed to minimize its frontal RCS, requiring low-observable shape design and flight attitude adjustment [6–8]. The satellite keeps its front toward Earth when it is flying over ground-based radar detection areas.


### 3.1. Radar Detection Area

In order to avoid reflecting radar signals directly, the scanning range of the ground-based radar should be modeled first. Figure4 shows low-thrust GEO satellite transfer with no radar detection. Figure 5 is an illustration of that GEO satellite being detected in low-thrust transfer by a ground-based radar, and the shadow is the radar detection area.Figure 4
Low-thrust transfer with no radar detection.Figure 5
Low-thrust transfer with radar detection.Due to the limited energy and probability of intercept, radar repeats its search for the target in a narrow area, which is modeled as a specific coverage of yaw angle, pitch angle, and operating distance, rather than the entire airspace [22–24].The latitude and longitude of the ground station are defined, respectively, asλp and δp. S is the track of the subsatellite point, whose right ascension and declination are denoted as λ and δ. We is the rotation speed of Earth.Consider(11)λ=artan⁡cos⁡itan⁡u+Ω-Wetδ=arsin⁡sin⁡isin⁡u.The radar pitch angleψh and radar yaw angle ϕh are obtained by the spherical triangle. Respectively,(12)tan⁡ψh=rcos⁡γ-Rercos⁡γsin⁡ϕh=sin⁡λ-λpsin⁡90∘-δsin⁡γ,where r=a(1-e2)/1+ecosv·γ is the geocentric angle between the subsatellite point and observation points. It satisfies(13)cos⁡γ=sin⁡δpsin⁡δ+cos⁡δpcos⁡δcos⁡λ-λp.However their distance is(14)ROP=Re+h2+Re2-2ReRe+hcos⁡γ.Suppose the radar yaw angle[σ1,σ2], pitch angle [τ1,τ2], and the operating distance rop:(15)ψh∈σ1,σ2ϕh∈τ1,τ2ROP≤rop.Only a satellite meeting these conditions has access to the radar detection area.


### 3.2. Thrust Control Angels and Satellite Attitude

Four thrusters are installed on the floor of the satellite, as shown in Figure6. Two work during orbit transfer, and the other two act as a backup system. Thrust is perpendicular to the bottom of the satellite in low-thrust transfer. The satellite has a frontal RCS shape design. Assuming the maximum thrust is F, then thrust components in the satellite coordinates are (16)FxFyFzT=0F0T.Figure 6
Thruster installation.The RTN frame is where the attitude of three axis-stabilized satellites is defined, in which thrust components are described as(17)FxoFyoFzoT=Fcos⁡βcos⁡αFcos⁡βsin⁡αFsin⁡βT.Satellite attitude is related to the order of three rotations [20]. The yaw angle, roll angle, and pitch angle of the satellite are recorded as ϕ, ψ, and θ and are derived in the order of 3-1-2. Define attitude matrix A312=R2(θ)R1(ψ)R3(ϕ) and its transpose matrix A312-1=R3-1ϕR1-1ψR2-1(θ). The relationship between satellite attitude and thrust control angles is(18)cos⁡βcos⁡αsin⁡βcos⁡βsin⁡α=A312-1010=R3-1ϕR1-1ψR2-1θ010=cos⁡ϕ-sin⁡ϕ0sin⁡ϕcos⁡ϕ00011000cos⁡ψ-sin⁡ψ0sin⁡ψcos⁡ψcos⁡θ0sin⁡θ010-sin⁡θ0cos⁡θ010=cos⁡ϕcos⁡θ-sin⁡ϕsin⁡ψsin⁡θ-sin⁡ϕcos⁡ψcos⁡ϕsin⁡θ+sin⁡ϕsin⁡ψcos⁡θsin⁡ϕcos⁡θ+cos⁡ϕsin⁡ψsin⁡θcos⁡ϕcos⁡ψsin⁡ϕsin⁡θ-cos⁡ϕsin⁡ψcos⁡θ-sin⁡ϕsin⁡θsin⁡ψcos⁡ψcos⁡θ010=-sin⁡ϕcos⁡ψcos⁡ϕcos⁡ψsin⁡ψ.It is ultimately expressed as(19)cos⁡βcos⁡α=-sin⁡ϕcos⁡ψsin⁡β=cos⁡ϕcos⁡ψcos⁡βsin⁡α=sin⁡ψ.The attitude adjustment range in the radar irradiation area is assumed to be(20)ϕmin<ϕ<ϕmaxψmin<ψ<ψmaxθmin<θ<θmax.The corresponding acceleration thrust components can be obtained by(21)fr=fcos⁡βsin⁡α=-fsin⁡ϕcos⁡ψft=fcos⁡βcos⁡α=fcos⁡ϕcos⁡ψfn=fsin⁡β=fsin⁡ψ.To realize low RCS toward Earth, our paper sets the optimal satellite attitude at a range of±5°, which can be adjusted according to the real simulation and test. Consider(22)ϕ∈-5∘,5∘ψ∈-5∘,5∘θ∈175∘,185∘.


## 3.1. Radar Detection Area

In order to avoid reflecting radar signals directly, the scanning range of the ground-based radar should be modeled first. Figure4 shows low-thrust GEO satellite transfer with no radar detection. Figure 5 is an illustration of that GEO satellite being detected in low-thrust transfer by a ground-based radar, and the shadow is the radar detection area.Figure 4
Low-thrust transfer with no radar detection.Figure 5
Low-thrust transfer with radar detection.Due to the limited energy and probability of intercept, radar repeats its search for the target in a narrow area, which is modeled as a specific coverage of yaw angle, pitch angle, and operating distance, rather than the entire airspace [22–24].The latitude and longitude of the ground station are defined, respectively, asλp and δp. S is the track of the subsatellite point, whose right ascension and declination are denoted as λ and δ. We is the rotation speed of Earth.Consider(11)λ=artan⁡cos⁡itan⁡u+Ω-Wetδ=arsin⁡sin⁡isin⁡u.The radar pitch angleψh and radar yaw angle ϕh are obtained by the spherical triangle. Respectively,(12)tan⁡ψh=rcos⁡γ-Rercos⁡γsin⁡ϕh=sin⁡λ-λpsin⁡90∘-δsin⁡γ,where r=a(1-e2)/1+ecosv·γ is the geocentric angle between the subsatellite point and observation points. It satisfies(13)cos⁡γ=sin⁡δpsin⁡δ+cos⁡δpcos⁡δcos⁡λ-λp.However their distance is(14)ROP=Re+h2+Re2-2ReRe+hcos⁡γ.Suppose the radar yaw angle[σ1,σ2], pitch angle [τ1,τ2], and the operating distance rop:(15)ψh∈σ1,σ2ϕh∈τ1,τ2ROP≤rop.Only a satellite meeting these conditions has access to the radar detection area.


## 3.2. Thrust Control Angels and Satellite Attitude

Four thrusters are installed on the floor of the satellite, as shown in Figure6. Two work during orbit transfer, and the other two act as a backup system. Thrust is perpendicular to the bottom of the satellite in low-thrust transfer. The satellite has a frontal RCS shape design. Assuming the maximum thrust is F, then thrust components in the satellite coordinates are (16)FxFyFzT=0F0T.Figure 6
Thruster installation.The RTN frame is where the attitude of three axis-stabilized satellites is defined, in which thrust components are described as(17)FxoFyoFzoT=Fcos⁡βcos⁡αFcos⁡βsin⁡αFsin⁡βT.Satellite attitude is related to the order of three rotations [20]. The yaw angle, roll angle, and pitch angle of the satellite are recorded as ϕ, ψ, and θ and are derived in the order of 3-1-2. Define attitude matrix A312=R2(θ)R1(ψ)R3(ϕ) and its transpose matrix A312-1=R3-1ϕR1-1ψR2-1(θ). The relationship between satellite attitude and thrust control angles is(18)cos⁡βcos⁡αsin⁡βcos⁡βsin⁡α=A312-1010=R3-1ϕR1-1ψR2-1θ010=cos⁡ϕ-sin⁡ϕ0sin⁡ϕcos⁡ϕ00011000cos⁡ψ-sin⁡ψ0sin⁡ψcos⁡ψcos⁡θ0sin⁡θ010-sin⁡θ0cos⁡θ010=cos⁡ϕcos⁡θ-sin⁡ϕsin⁡ψsin⁡θ-sin⁡ϕcos⁡ψcos⁡ϕsin⁡θ+sin⁡ϕsin⁡ψcos⁡θsin⁡ϕcos⁡θ+cos⁡ϕsin⁡ψsin⁡θcos⁡ϕcos⁡ψsin⁡ϕsin⁡θ-cos⁡ϕsin⁡ψcos⁡θ-sin⁡ϕsin⁡θsin⁡ψcos⁡ψcos⁡θ010=-sin⁡ϕcos⁡ψcos⁡ϕcos⁡ψsin⁡ψ.It is ultimately expressed as(19)cos⁡βcos⁡α=-sin⁡ϕcos⁡ψsin⁡β=cos⁡ϕcos⁡ψcos⁡βsin⁡α=sin⁡ψ.The attitude adjustment range in the radar irradiation area is assumed to be(20)ϕmin<ϕ<ϕmaxψmin<ψ<ψmaxθmin<θ<θmax.The corresponding acceleration thrust components can be obtained by(21)fr=fcos⁡βsin⁡α=-fsin⁡ϕcos⁡ψft=fcos⁡βcos⁡α=fcos⁡ϕcos⁡ψfn=fsin⁡β=fsin⁡ψ.To realize low RCS toward Earth, our paper sets the optimal satellite attitude at a range of±5°, which can be adjusted according to the real simulation and test. Consider(22)ϕ∈-5∘,5∘ψ∈-5∘,5∘θ∈175∘,185∘.


## 4. Low-Thrust GEO Transfer Design

The goal of this paper is to design a minimum-time transfer for geostationary spacecraft equipped with electric propulsion systems. The transfer problem is thus to find an essentially bound control to reduce eccentricity and inclination and raise the semimajor axis.After simplifying (1), they fall into(23)a˙=2a2hesin⁡vfr+1+ecos⁡vfte˙=phsin⁡vfr+cos⁡v+e+cos⁡v1+ecos⁡vfti˙=rhcos⁡ω+vfnΩ˙=rsin⁡ω+vhsin⁡ifnω˙=-pcos⁡vhefr+p+rsin⁡vheft+rsin⁡vω+vcos⁡ihsin⁡ifnM˙=n+1aehpcos⁡v-2erfr-p+rsin⁡vft.The acceleration componentfn only contributes to a decrease of inclination. The change of the semimajor axis and eccentricity are both related to fr and ft. fr and ft have a greater effect on the semimajor axis than eccentricity.First, a parameterk is introduced to isolate fn from f:(24)frt=fr2+ft2=f·kfn=f1-k2with 0≤k≤1.Then, a control law that reduces eccentricity quickly is derived by observing the time-rate equation for the semimajor axis and eccentricity. Defineep1=sin⁡v and ep2=cos⁡v+e+cos⁡v/1+ecos⁡v. The law is given by(25)fr=f·kft=0ep1>ep2fr=0ft=f·kep1≤ep2.Five thrust modes and both their corresponding modulus and direction are summed up in Table1.Table 1
Thrust modes in low-thrust transfer.
Direction
Modulus
Condition
Mode 1
Mode 2
Mode 3
Mode 4
Mode 5
k
≠
1
k
=
1
k
=
1
Earth shadow
Low observability
i
p
=
cos
⁡
(
ω
+
v
)
ip>0
f
n
=
-
f
1
-
k
2
f
n
=
0
f
n
=
0
f
n
=
0
ip≤0
f
n
=
f
1
-
k
2
e
p
1
>
e
p
2
f
n
=
0
ep1≥0
f
r
=
-
f
⋅
k
,
f
t
=
0
f
r
=
-
f, ft=0
f
r
=
f ft=0
f
r
=
0
f
r
=
0, ft=0
ep1<0
f
r
=
f
⋅
k
,
f
t
=
0
f
r
=
f, ft=0
f
t
=
0
f
r
=
0, ft=0
e
p
1
≤
e
p
2
ep2≥0
f
r
=
0, ft=-f⋅k
f
r
=
0, ft=-f
f
r
=
0, ft=0
ep2<0
f
r
=
0, ft=f⋅k
f
r
=
0, ft=f
f
r
=
0, ft=fAs stated, our transfer problem is parameterized, and the control is given according toip, ep1, and ep2, which are determined by the argument of perigee and the true anomaly. Because the parameter k is between 0 and 1, it is convenient to optimize it to minimize the transfer time by a simple traversal, which is a selection of k from 0 to 1 by varying it in a prefixed step. As Table 1 details fully, our method for achieving terminal orbit is performed in three steps.Step 1 (reduction of inclination).
The thrust mode relies on mode 1 to reduce inclination by rationalk, along with the simultaneous targeting of the semimajor axis and eccentricity. The thrust form is shown in Figure 7.Figure 7
Reduction of inclination.Step 2 (raising the semimajor axis).
Mode 2 represents a law used to control the semimajor axis and the eccentricity when the inclination reaches the target value. The thrust operates as shown in Figure8.Figure 8
Raising semimajor axis.Step 3 (reducing eccentricity).
While the semimajor axis reaches its target,ft remains at 0 to ensure the stability of the semimajor axis. fr is used to reduce eccentricity in mode 3, as illustrated in Figure 9.Figure 9
Reducing eccentricity.Step 4.
Mode 4 is the thrust model for Earth shadow, according to formula (9) to realize real-time judgment of shadow, in which thrust is 0.Step 5.
Mode 5 is used for the radar detection area to realize satellite low observability. The satellite adjusts to its best attitude, corresponding to thrust componentsfr=0, fn=0, and, if ep2<0, ft=f; otherwise, ft=0, because ft=f contributes to both the decrease of eccentricity and the increase of the semimajor axis when ep2<0 (Figure 10).Figure 10
Overall scheme.In short, the method described in this paper determines a control depending on a parameter and the state of current orbit. Our control law is simple, with only one uncertain control parameter, whose optimization is so convenient that our minimum-time transfer problem is greatly simplified.


## 5. Simulation

### 5.1. Initial Parameter

The ideal attitude of the satellite was assumed to have the range of its yaw angle, roll angle, and pitch angle as[-5∘,5∘], [-5∘,5∘], and [175∘,185∘]. The yaw angle of each radar ranges from −60° to 60°, pitch angle ranges from 0 to 90°, and operating distance was 15,000 km. They were located in the following places (longitude/latitude): (65°E, 35°N), (30°E, 75°N), (75°E, 75°N), (120°E, 75°N), (165°E, 75°N), (15°E, 75°N), (60°W, 75°N), (105°W, 75°N), and (150°W, 75°N); (30°E, 75°S), (75°E, 75°S), (120°E, 75°S), (165°E, 75°S), (15°E, 75°S), (60°W, 75°S), (105°W, 75°S), and (150°W, 75°S); (30°E, 45°N), (75°E, 45°N), (120°E, 45°N), (165°E, 45°N), (15°E, 45°N), (60°W, 45°N), (105°W, 45°N), and (150°W, 45°N); (30°E, 45°S), (75°E, 45°N), (120°E, 45°S), (165°E, 45°S), (15°E, 45°S), (60°W, 45°S), (105°W, 45°S), and (150°W, 45°S); (30°E, 15°N), (75°E, 15°N), (120°E, 15°N), (165°E, 15°N), (15°E, 15°N), (60°W, 15°N), (105°W, 15°N), and (150°W, 15°N); (30°E, 15°S), (75°E, 15°S), (120°E, 15°S), (165°E, 15°S), (15°E, 15°S), (60°W, 15°S), (105°W, 15°S), and (150°W, 15°S).The satellite’s initial massm0=2600kg was assumed to be equipped with four thrusters with a specific impulse of Isp=1600s and a maximum thrust of F=80mN. The GTO departure date was fixed as 1 June 2008 for shadow calculations. Table 2 displays the initial GTO and GEO boundary conditions.Table 2
Orbit parameters.
Orbit parameters
Initial orbit
Final orbit
GTO
GEO
Semimajor axisa (km)
24731
42164
Inclinationi (deg)
28.5
0
Eccentricitye
0.7194
0
Right ascensionΩ (deg)
1.7104
Free
Argument of perigeeω (deg)
179.6
Free
Mean anomalyM (deg)
2.9044
Free


### 5.2. Thrust and Control Angles

Figures11 and 12 show the time-evolution of thrust and thrust control angles in the orbit transfer, in which Earth shadow and oblateness (J2) effects were considered with no low-observable constraint. In contrast, Figures 13 and 14 show results taking into consideration all of the aforementioned factors.Figure 11
Time-evolution of thrust with no radar detection.Figure 12
Time-evolution of thrust control angles with no radar detection.Figure 13
Time-evolution of thrust with radar detection.Figure 14
Time-evolution of thrust control angles with radar detection.As can be seen in Figures11–14, the satellite flew over the radar detection area three times in these six orbit periods, respectively, 6.263×106~6.294×106s,6.363×106~6.377×106s, and 6.44×106~6.465×106s. The thrust components were fr=0, ft=f, and fn=0. The thruster operated in mode 5 to maintain low observability of the satellite (Table 3).Table 3
Low-thrust transfer with low observability.
Start/stop time (s)
Total time (s)
Thrust componets N
Control angels
1
6.263
×
1
0
6~6.294×106 s
3.1
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0
2
6.363
×
1
0
6~6.377×106 s
1.4
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0
3
6.44
×
1
0
6~6.465×106 s
2.5
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0


### 5.3. Numerical Results

To demonstrate the validity of our proposed method, four scenarios are presented in Table4. Case 1 represents a transfer based on the optimal control method l, using Pontryagin’s maximum principle to find the optimum control law [14]. Cases 2 and 3 represent transfers based on our method. Case 3 involves Earth shadow and oblateness (J2) effects. Cases 4 and 5 involve Earth shadow, oblateness (J2) effects, and low observability. Case 4 considers only one radar detection area (65°E, 35°N), while Case 5 includes other 48 areas.Table 4
Fuel consumption and transfer days.
Mission
GTO
→
GEO
Case1
Case2
Case3
Case4
Case5
Fuel consumptionΔm (kg)
404
431
443
507
723
Transfer daysΔt (d)
230
245
251
288
410The results of these four cases are summarized in detail in Table4. The time-evolution of state for Case 3 is presented in Figures 15–17. The time-evolution of state for Case 4 is shown in Figures 18–20. The time-evolution of state for Case 5 is presented in Figures 21–23. GTO-GEO transfer times were 230 days and 245 days for the optimal method and our method, respectively. Case 4 led to a 288-day transfer, requiring 507 kg of fuel. Case 5 was a 410-day transfer consuming 723 kg of fuel. Compared with the results of the optimal control, our method clearly offered performance near the results of optimal method with no ground estimation, was particularly flexible when orbit elements changed, and presented no two-point boundary problem and a low computational burden.Figure 15
Time-evolution of the semimajor axis of Case 3.Figure 16
Time-evolution of the eccentricity of Case 3.Figure 17
Time-evolution of the inclination of Case 3.Figure 18
Time-evolution of the semimajor axis of Case 4.Figure 19
Time-evolution of the eccentricity of Case 4.Figure 20
Time-evolution of the inclination of Case 4.Figure 21
Time-evolution of the semimajor axis of Case 5.Figure 22
Time-evolution of the eccentricity of Case 5.Figure 23
Time-evolution of the inclination of Case 5.


## 5.1. Initial Parameter

The ideal attitude of the satellite was assumed to have the range of its yaw angle, roll angle, and pitch angle as[-5∘,5∘], [-5∘,5∘], and [175∘,185∘]. The yaw angle of each radar ranges from −60° to 60°, pitch angle ranges from 0 to 90°, and operating distance was 15,000 km. They were located in the following places (longitude/latitude): (65°E, 35°N), (30°E, 75°N), (75°E, 75°N), (120°E, 75°N), (165°E, 75°N), (15°E, 75°N), (60°W, 75°N), (105°W, 75°N), and (150°W, 75°N); (30°E, 75°S), (75°E, 75°S), (120°E, 75°S), (165°E, 75°S), (15°E, 75°S), (60°W, 75°S), (105°W, 75°S), and (150°W, 75°S); (30°E, 45°N), (75°E, 45°N), (120°E, 45°N), (165°E, 45°N), (15°E, 45°N), (60°W, 45°N), (105°W, 45°N), and (150°W, 45°N); (30°E, 45°S), (75°E, 45°N), (120°E, 45°S), (165°E, 45°S), (15°E, 45°S), (60°W, 45°S), (105°W, 45°S), and (150°W, 45°S); (30°E, 15°N), (75°E, 15°N), (120°E, 15°N), (165°E, 15°N), (15°E, 15°N), (60°W, 15°N), (105°W, 15°N), and (150°W, 15°N); (30°E, 15°S), (75°E, 15°S), (120°E, 15°S), (165°E, 15°S), (15°E, 15°S), (60°W, 15°S), (105°W, 15°S), and (150°W, 15°S).The satellite’s initial massm0=2600kg was assumed to be equipped with four thrusters with a specific impulse of Isp=1600s and a maximum thrust of F=80mN. The GTO departure date was fixed as 1 June 2008 for shadow calculations. Table 2 displays the initial GTO and GEO boundary conditions.Table 2
Orbit parameters.
Orbit parameters
Initial orbit
Final orbit
GTO
GEO
Semimajor axisa (km)
24731
42164
Inclinationi (deg)
28.5
0
Eccentricitye
0.7194
0
Right ascensionΩ (deg)
1.7104
Free
Argument of perigeeω (deg)
179.6
Free
Mean anomalyM (deg)
2.9044
Free


## 5.2. Thrust and Control Angles

Figures11 and 12 show the time-evolution of thrust and thrust control angles in the orbit transfer, in which Earth shadow and oblateness (J2) effects were considered with no low-observable constraint. In contrast, Figures 13 and 14 show results taking into consideration all of the aforementioned factors.Figure 11
Time-evolution of thrust with no radar detection.Figure 12
Time-evolution of thrust control angles with no radar detection.Figure 13
Time-evolution of thrust with radar detection.Figure 14
Time-evolution of thrust control angles with radar detection.As can be seen in Figures11–14, the satellite flew over the radar detection area three times in these six orbit periods, respectively, 6.263×106~6.294×106s,6.363×106~6.377×106s, and 6.44×106~6.465×106s. The thrust components were fr=0, ft=f, and fn=0. The thruster operated in mode 5 to maintain low observability of the satellite (Table 3).Table 3
Low-thrust transfer with low observability.
Start/stop time (s)
Total time (s)
Thrust componets N
Control angels
1
6.263
×
1
0
6~6.294×106 s
3.1
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0
2
6.363
×
1
0
6~6.377×106 s
1.4
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0
3
6.44
×
1
0
6~6.465×106 s
2.5
×
1
0
4
f
r
=
0, ft=f, and fn=0
α
=
0, β=0


## 5.3. Numerical Results

To demonstrate the validity of our proposed method, four scenarios are presented in Table4. Case 1 represents a transfer based on the optimal control method l, using Pontryagin’s maximum principle to find the optimum control law [14]. Cases 2 and 3 represent transfers based on our method. Case 3 involves Earth shadow and oblateness (J2) effects. Cases 4 and 5 involve Earth shadow, oblateness (J2) effects, and low observability. Case 4 considers only one radar detection area (65°E, 35°N), while Case 5 includes other 48 areas.Table 4
Fuel consumption and transfer days.
Mission
GTO
→
GEO
Case1
Case2
Case3
Case4
Case5
Fuel consumptionΔm (kg)
404
431
443
507
723
Transfer daysΔt (d)
230
245
251
288
410The results of these four cases are summarized in detail in Table4. The time-evolution of state for Case 3 is presented in Figures 15–17. The time-evolution of state for Case 4 is shown in Figures 18–20. The time-evolution of state for Case 5 is presented in Figures 21–23. GTO-GEO transfer times were 230 days and 245 days for the optimal method and our method, respectively. Case 4 led to a 288-day transfer, requiring 507 kg of fuel. Case 5 was a 410-day transfer consuming 723 kg of fuel. Compared with the results of the optimal control, our method clearly offered performance near the results of optimal method with no ground estimation, was particularly flexible when orbit elements changed, and presented no two-point boundary problem and a low computational burden.Figure 15
Time-evolution of the semimajor axis of Case 3.Figure 16
Time-evolution of the eccentricity of Case 3.Figure 17
Time-evolution of the inclination of Case 3.Figure 18
Time-evolution of the semimajor axis of Case 4.Figure 19
Time-evolution of the eccentricity of Case 4.Figure 20
Time-evolution of the inclination of Case 4.Figure 21
Time-evolution of the semimajor axis of Case 5.Figure 22
Time-evolution of the eccentricity of Case 5.Figure 23
Time-evolution of the inclination of Case 5.


## 6. Conclusion

A new method based on control parameter analysis has been used to design a low-thrust orbit transfer. The optimization of the parameter introduced here is utilized to solve the optimal control problem. Satellite low observability was considered in GEO transfers, and the effects of Earth shadow and oblateness were also included in our method.This approach possesses three important features: (1) The approach is characterized by simplicity and feasibility, on-board real-time solution, less calculation, and no initial guess for convergence, (2) the influence of Earth shadow and perturbation is taken into account, proving its strong fault tolerance, and (3) a simple but effective approach is used to achieve satellite low observability and reduce the threat of radar. In general, this trajectory method would serve as a useful preliminary design for GEO mission designers.



---

*Source: 439815-2015-12-22.xml*