Electrostatic Forces for Swarm Navigation and Reconfiguration

Full text


Electrostatic Forces for Swarm Navigation and Reconfiguration

Muhammad Kamran Saleem (791127-P133)

This thesis is presented as part of Degree of Master of Science in Electrical Engineering

Blekinge Institute of Technology February 2007

Blekinge Institute of Technology ZARM Center of Applied Space School of Engineering Technology and Microgravity Department of Applied Signal Processing Department of Space Technology


This thesis is based on the research work done for a project sponsored by European Space Agency Advanced Concepts Team (ACT). I would like to thank ESA for providing the financial support for this project. I would like to express my sincere gratitude to Dr. Stephan Theil for giving me the opportunity to work on my master thesis at ZARM. It is an honor to work in his team under his kind supervision. I am grateful to Lorenzo Pettazzi and Hans Kruger for the guidance and really appreciate both for helping me by all possible means and being patient with me during my thesis. I am very thankful to Prof. Hans -J¨urgen Zepernick for being my internal supervisor and for his support and encouragement for my master thesis. I would also like to thank all my friends in Sweden and Germany for their academic/non-academic help.

My outmost thanks to my family for supporting me during my studies, in particular to my parents Saleem and Zahida, without whose well wishes and efforts, I could not have been able to pursue my masters.

Finally, I would like to express all my love to my fiancee Sidra for all her everlasting love and support.


1 Orbital Mechanics 3

1.1 Relative Two-Body Equation of Motion . . . 3

1.2 Constant of Motion . . . 4

1.3 Geometry of Conic Section . . . 6

1.4 Orbital Elements . . . 8

1.5 Euler’s Angles . . . 10

1.6 Clohessy Wiltshire Equations . . . 12

1.6.1 Particle Kinematics with Moving Frames . . . 13

1.6.2 General Relative Orbit Description . . . 15

1.6.3 Derivation of Hills Equation . . . 16

1.7 Lagrangian Points and Three Body Problem . . . 20

2 Space Environment 23 2.1 Space Plasma Environment . . . 25

2.2 Debye Length . . . 27

2.3 Charged Particle Exchange in Space Environment . . . 28

2.4 Differential Charging . . . 29

2.5 Current Voltage Characteristics . . . 29

2.6 Spacecraft Charging . . . 30

2.6.1 Eclipse Charging . . . 30

2.6.2 Sunlight Charging . . . 31

2.6.3 Controlling Spacecraft Potential . . . 31

2.7 Environmental and Artificial Currents . . . 33

2.8 Spacecraft Charging Model . . . 34

3 Hybrid Propulsion System 37 3.1 Natural Spacecraft Charging at GEO . . . 37

3.2 Review of Charge Mitigation Methods . . . 39

3.2.1 Active Charge Control Requirements . . . 41

3.3 Micro-Propulsion Devices . . . 42

3.3.1 Thruster Technologies Overview . . . 42

3.3.2 Field Emission Electric Propulsion (FEEP) . . . 43

3.3.3 Colloid Thruster . . . 43

3.3.4 Pulsed Plasma Thruster (PPT) . . . 44

3.3.5 Radio Frequency Ion Thruster (RIT) . . . 45

3.4 Electron Gun . . . 47

3.5 Spacecraft Charge Control Strategy . . . 48

3.5.1 Charging Current . . . 48

3.5.2 Combination of Charging and Stabilizing Current . . . 50

3.6 Evaluation of Duration of the Charging Process . . . 52


4.6 Compatibility between ES and CS Concepts . . . 66

4.7 Cost Function Definition . . . 67

4.8 Differential Evolution Optimization Algorithm . . . 69

4.9 ES and CS Compatibility Test Outcomes . . . 70

4.9.1 Results for GEO Environment . . . 70

4.9.2 Results for L1 Environment . . . 73

4.10 Stand-by Formation . . . 76

5 Summary 79


With regard to the benefits of using multiple vehicles, the concept of formation control has been studied extensively in the literature with application to the coordination of multiple robots, unmanned air vehicles (UAVs), autonomous underwater vehicles (AUVs), satellites, aircraft and spacecraft. There are many advantages of using formations of multiple vehicles, such as increased feasibility, robustness, accuracy, flexibility, cost efficiency, energy efficiency, and probability of success.

The propulsion system on a spacecraft is of key importance in formation flying missions and life time of mission depends on the available fuel on the spacecraft. Furthermore, in the missions where spacecrafts are flying in close proximity the exhaust plums from the thrusters of a satellite can damage the near by spacecraft which might result in total failure to the mission. In this study a recently developed concept of Coulomb satellite equipped with hybrid propulsion system will be presented. The coulomb satellites utilizes µN thrusters to eject charges in certain range in order to drive spacecraft to some specific potential level through which controlled attractive and repulsive forces are generated within the spacecrafts in the formations. Through utilization of Electrostatic forces between the swarm members the fuel consumption is reduced sufficiently.

Furthermore, formation flying missions requires a lot of inter spacecraft communication result- ing in complex systems with excessive payloads. Utilizing Equilibrium shaping technique this complexity and inter spacecraft communications can be minimized to a very low level. Further- more, the introduction of Equilibrium shaping technique in formation flying missions produces autonomous decision making capabilities within the agents of the swarm through which final target position are preassigned to each agent in the formation. With the combination of Equilib- rium shaping technique and Coulomb satellite concepts certain formations can be produced for which minimum communication between swarm agents and maximum fuel saving (up to 80%) are possible.


1.1 Gravity and disturbance forces between two bodies [28]. . . 3

1.2 Geometry of an elliptic conic section [21]. . . 6

1.3 Conic section geometry [21]. . . 7

1.4 Inclination angle [21]. . . 9

1.5 Right ascension of ascending node [21]. . . 10

1.6 Argument of perigee [21]. . . 10

1.7 Graphical representation for Euler angle [11]. . . 12

1.8 Two coordinate frames with moving origins [28]. . . 14

1.9 Illustration of leader follower type of two spacecraft formation [28]. . . 15

1.10 General type of spacecraft formation with out-of-orbit plane relative motion [2]. . 16

1.11 Lagrangian points [12]. . . 20

1.12 Circular restricted three body problem [37]. . . 21

2.1 Value of point charge q1 needed to apply specified electrostatic forces in vacuum, q2= 2µC, q1max= 2µC [23]. . . 24

2.2 Plot of altitude Vs electron density for the LEO [10]. . . 26

2.3 plot of altitude Vs ion composition for the LEO [10]. . . 26

2.4 Current exchange between a surface and its space environment [8]. . . 28

2.5 Current flowing to and from the spacecraft. . . 29

2.6 J-V graph for spacecraft [19]. . . 30

2.7 Illustration of the exploiting the net torques due to non-uniform spacecraft charg- ing [29]. . . 32

2.8 Simple illustration of the spacecraft charging in a plasma environment [29]. . . . 32

2.9 Mitigation of negative surface potential by means of ion emission from SCATHA [17]. 33 2.10 Emission of electrons and ions from a negatively charged surface [17]. . . 33

2.11 Artificial and Environmental currents. . . 34

3.1 Natural space plasma properties [18]. . . 37

3.2 Spacecraft floating potential at eclipse [23]. . . 38

3.3 V-I characteristics for sunlight charging [23]. . . 39

3.4 Discharge of SC surface by sharp spike or a hot filament [17]. . . 40

3.5 GEO worst case I-V characteristics [23]. . . 42

3.6 Basic principle design of FEEP [27]. . . 44

3.7 Principle operation of Colloid Thruster [6]. . . 45

3.8 Principle design of PPT [25]. . . 46

3.9 Principle design of RIT [30]. . . 47

3.10 Electron guns [14]. . . 47

3.11 Charge control using Ich [23]. . . 49

3.12 Charge control with combination of the charging and the stabilizing current [23]. 52 3.13 Charge control with combination of the charging and the stabilizing current in sunlight; plasma Case 2 [23]. . . 53


4.8 Illustration of crossover and mutation process. . . 70

4.9 Polygonal equilibrium formations [23]. . . 71

4.10 Solids equilibrium formations [23]. . . 73

4.11 Formations without Debye length [23]. . . 75

4.12 Formations with Debye length [23]. . . 76

4.13 Switch maneuver between the design and stand-by formation [23]. . . 77

4.14 Stand-by and design formation for 8 spacecrafts [23]. . . 78


1.1 Orbital Elements . . . 11

2.1 Average GEO environment [19]. . . 27

2.2 Worst case GEO environment [19]. . . 27

2.3 Range of Debye length in various plasma environment [16]. . . 28

3.1 Plasma parameters for different GEO environment during eclipse [33]. . . 38

3.2 Some important parameter values for FEEP [30]. . . 43

3.3 Some important parameter values for colloid thruster [30]. . . 44

3.4 Some important parameter values for PPT [30]. . . 45

3.5 Some important parameter values for RIT thruster [30]. . . 46

3.6 Plasma parameters for low energy GEO environments. . . 48

3.7 control current parameters for using RIT-4 as ion emitting device for a given Vdes; ∆tcharge=1 s, I(∆tcharge)= 8 µA, Inat=[-65µA, +65µA] . . . 51

3.8 Worst case residual forces due to EA utilizing RIT-4 engine; Control current parameters from Table 3.7. . . 55

4.1 Equations and unknowns for different formations. . . 64

4.2 List of the cost function values for some formations that can be acquired by means of the ES technique. . . 72

4.3 Results Without Debye length. . . 74

4.4 Results With Debye length. . . 74

4.5 Characteristics for stand-by and design formation. . . 78


The relative dynamics of two spacecrafts flying in close proximity has been studied since the beginning of the manned space program. Spacecraft flying in close proximity is not a new concept. Since 1960’s numerous manned and unmanned spacecraft have performed operations in close proximity, mostly for the purpose of docking for re-supply or crew transfer, and earth observation. The recent research focus of spacecraft formation flying has now extended to maintain a formation of various spacecrafts with autonomous control. For example, U.S Air Force is studying concepts of having identical satellites form a sparse aperture radar dish in space. Use of multiple satellites flying in a specific geometry avoids the significant technical and financial challenge of attempting to build a single radar dish of equivalent size. These satellite formations can have variable formation diameters ranging from several dozens of meters to several kilometers [28]. Once the formation flying technology gets fully matured, then might result in swarms of space vehicles flying as virtual platforms or distributed space systems and sensor webs which gather significantly more and better science data than is possible today [4].

The first big achievement toward formation flying in space was made on December 15, 1965, when Gemini VI accomplished the first space rendezvous and performed station-keeping with Gemini VII over three orbits at distances of 0.3 to 90 meters. Then Gemini VIII followed with the first docking of two spacecraft in orbit on March 16, 1966. The first automated rendezvous was performed by Soviet spacecraft Cosmos 186 and Cosmos 188 on October 27, 1967, when the two unmanned spacecraft performed a preprogrammed closure and docking maneuver [3, 20].

Another of the big achievement in the field of autonomous formation flight was made in May 17, 2001, when Earth Observing-1 (EO-1) performed its first autonomous maneuver to maintain a one minute in-track separation with Landsat-7 [5]. The future space missions involving multiple spacecrafts flying in deferent formations such as Darwin, Terrestrial Planet Finder (TPF), LISA etc will make use of the useful data that has been collected from the successful formation flight of the EO-1 and Landsat7.

The knowledge of orbital mechanics is essential for a full understanding of space operations.

Chapter 1 gives a brief theoretical and mathematical introduction to the field of orbital me- chanics. The different commonly used terms, laws and equations are explained briefly to give a solid background and to easily understand the work that is done in the later stages of this project. One of the main objective in this study is to integrate the Coulomb satellite concept into a swarm navigation technique in order to enhance the performance of swarm of satellite in terms of fuel consumption and formation control. The interaction of spacecraft with ambient plasma environment in space is of prime concern since it plays an important role in spacecraft charging.

Chapter 2 contains a brief introduction to space environment and natural phenomenas like spacecraft charging in eclipse and in sunlight are explained in detail. Furthermore, the different currents to and from (natural and artificial) the spacecraft are described and the spacecraft charging model is given at the end of this chapter.

The Coulomb satellites are equipped with hybrid propulsions system through which controlled emission of charges is carried out in order to build desired potential level on the spacecraft. In Chapter 3 the active charge control requirements are presented for such a propulsion system


in the optimization problem are discussed with the Differential Evolution (DE) algorithm and compatibility results of both concepts are shown at GEO and L1(Lagrangian Point1 between Earth and Sun) environment. Furthermore, a very interesting and useful application of combi- nation of equilibrium shaping and Coulomb satellite concept i.e. stand by formation is discussed at the end of this chapter.

Finally, Chapter 5 contains a brief summary of the work being done in this study.


Orbital mechanics is the study of the motions of artificial satellites and space vehicles moving under the influence of forces such as gravity, atmospheric drag, thrust, etc. The root of orbital mechanics can be traced back to the 17th century when mathematician Isaac Newton put forward his laws of motion and formulated his law of universal gravitation.

1.1 Relative Two-Body Equation of Motion

Consider two particles of mass m1 and m2 moving in space as shown in Figure 1.1. The only forces acting on them is the mutual gravitational attraction and some disturbance forces fd1and fd2 which might be present due to various reasons (drag, gravitational force etc). The magni- tude of the gravitational attraction between the masses is given by Newton’s law of universal gravitation as

F12= −F21= Gm1m2

|r12|2 r12

|r12| (1.1)

Figure 1.1: Gravity and disturbance forces between two bodies [28].

The position of m2 relative to mass m1 is given by the vector r12= R2− R1, where the position vectors R1and R2are measured relative to an inertial reference frame N. According to Newton’s second law of motion we know that

F = ma (1.2)

Where, m and a denotes the mass and acceleration of the body respectively. From now on the notation ˙r will be used for first derivative with respect to time (i.e. velocity, drdt = ˙r = v) and ¨r


µ = G(m1+ m2) (1.5) where G is the universal gravity constant. Note that for many systems m1 ≫ m2 and µ can therefore be approximated as

µ ≈ Gm1

An example of this situation would be a satellite in earth’s orbit. Since, the mass m2 of satellite would be negligible compared to earth. Now taking the difference between Eq(1.3) and Eq(1.4), we get the motion of m2 relative to m1 as

¨r = −µ

r3r + ad (1.6)

where ad is the disturbance acceleration defined as ad= 1

m2fd2− 1

m1fd1 (1.7)

The vector differential equation shown in Eq(1.6) is the most important result in celestial me- chanics. It forms the basis for various developments. The two disturbance accelerations are often near a cancellation. Consider earth-moon system with the sun’s gravitational attraction modeled as external influence. Let m3 be the sun’s mass then we can express the disturbance acceleration ad as

ad= 1 m2


|r23|3 r23− 1 m1


|r13|3 r13≈ 0 (1.8)

Since r23 ≈ r13 (i.e. distance from sun to moon and sun to earth respectively). Thus, even if sun’s gravitational force is very large, its effect on the relative two-body motion is often negligible. Therefore, the relative disturbance acceleration vector ad is typically considered to be small or actually set to zero to obtain a good approximate solution. Hence, we attain the relative two body equation of motion as

¨r + µ

r3r = 0 (1.9)

1.2 Constant of Motion

The gravitational field is “conservative”, which means that an object moving under the influence of gravity alone does not loose or gain mechanical energy but only exchanges one form of energy,

“kinetic” for another form called “potential energy”.


We also know that it takes a tangential component of force to change the angular momentum of a system in rotational motion about some center of rotation. Since the gravitational force is always directed radially toward the center of the large mass, we would expect that the angular momentum of the satellite about the center of our reference frame (large mass) does not change.

The proof of above mentioned statements will be shown in the sequel. There are two important constants of motion i.e.

ˆ specific mechanical energy

ˆ specific angular momentum

The energy constant of motion can be derived as follows:

Dot multiplying Eq(1.9) by ˙r gives

˙r · ¨r + ˙r · µ

r3r = 0 (1.10)

We can write the above equation as

v · ˙v + ˙r · µ

r3r = 0 (1.11)

Noticing that dtd(v22) = v · ˙v and dtd(−µr) = rµ2˙r we can write Eq(1.11) as d


 v2 2 −µ


= 0 (1.12)

If any expression has time derivative equals to zero then that expression must be a constant which we call ξ and in our case is the specific mechanical energy

ξ = v2 2 −µ


(1.13) The first term of ξ is the kinetic energy per unit mass of the satellite and second term is the potential energy per unit mass. Thus, we can conclude that the specific mechanical energy ξ of a satellite is the sum of kinetic energy per unit mass and its potential energy per unit mass remains constant along its orbit. This means that ξ of a satellite neither increases nor decreases as a result of its motion in its orbit. The angular momentum constant of the motion is obtained as follows: Cross multiplying Eq(1.9) by r will give us

r × ¨r + r × µ

r3r = 0 (1.14)

We know that a × a = 0 so the second term disappears and we get r × ¨r = 0. Noticing that


dt(r × ˙r) = ˙r × ˙r + r × ¨r the equation above becomes d

dt(r × ˙r) = 0 (1.15)

As ˙r = v, we can write above equation as d

dt(r × v) = 0

The expression r × v is called specific angular momentum and expressed as h. Since the time derivative of h is also zero, we can say that the specific angular momentum h of a satellite remains constant along its orbit and is shown as


1.3 Geometry of Conic Section

The basic knowledge of geometry of conic sections is of fundamental importance for under- standing orbital mechanics. This section provides a review of some of the important aspects of geometry of circular, elliptic, parabolic and hyperbolic orbits.

Kepler discovered that the orbit of one body about another body is an ellipse and Newton proved that all of the conic sections are feasible orbits i.e. circle, ellipse, parabola, hyperbola. A sample elliptical orbit is illustrated in Figure 1.2.

Figure 1.2: Geometry of an elliptic conic section [21].

The shape of ellipse is defined through its semi-major axis a and semi-minor axis b where a ≥ b.

For the circle 2a is simply the diameter, for the parabola 2a is infinite and for the hyperbola 2a is taken as negative. Every ellipse has two focal points F1 and F2. A special case of ellipse is a circle where the two focal points occupy the same point (i.e., a = b). The distance between the foci is represented as 2c. For the circle 2c is zero, for parabola 2c is infinite and for the hyperbola 2c is taken as a negative. The eccentricity e is an important parameter and describes the shape of orbit. It indicates whether the conic intersection is elliptic, parabolic or hyperbolic.

For ellipses the eccentricity ranges between 0 ≤ e ≤ 1. Parabolas always have e = 1, hyperbolas


always have eccentricities greater than 1 and circles have e = 0. The extreme end points of the major axis of an orbit are referred to as apses. The point nearest to the prime focus is called

“periapsis” and the point farthest from the prime focus is called “apoapsis”. The geometry of circle, ellipse, parabola and hyperbola are shown in Figure 1.3.

Figure 1.3: Conic section geometry [21].

The conic section can be described mathematically as following. Let the vector r point from the focus F to the current orbit position with r being its magnitude. The distance p is the perpendicular distance (to the major axis) between the focus and the orbit and is called semilatus rectum as shown in Figure 1.2. The angle f measures the heading of position vector r relative to the semimajor axis and is called anomaly. The vector components (x,y) of r are expressed as

x = rcosf (1.17)

y = rsinf (1.18)

it follows directly from the definition of a conic section given above that for any conic except a parabola

e = c

a (1.19)


p = a(1 − e2) (1.20)

The general polar equation of conic section is expressed as

r = p

1 + ecosf (1.21)

The distance from the prime focus to either periapsis or apoapsis can be expressed by simply inserting f = 0 or f = 180 (radians) in Eq(1.21). Thus for any conic, we have


rmax= rapoapsis= p

1 + ecos180 (1.24)


ra= p

1 − e = a(1 + e) (1.25)

1.4 Orbital Elements

The inertial reference frame is the first requirement for describing an orbit. There are differ- ent types of coordinate systems. The most important are heliocentric coordinate frame and geocentric coordinate frame as described below [1]:

ˆ Heliocentric coordinate system: The heliocentric coordinate system is mostly used for the orbits around the sun such as planets, asteroids, comets. As the name implies, the center of this coordinate system is the center of the sun. The fundamental plane is “ecliptic” which is the plane of the earth’s revolution around the sun. The line of intersection of the ecliptic plane and the earth’s equatorial plane defines the direction of the X-axis. On first day of spring a line joining the center of the earth and the center of the sun points in the direction of positive X-axis. This is called vernal equinox direction and points toward the direction of the constellation Aries.

ˆ Geocentric coordinate system: The geocentric equatorial system has its origin at the earth’s center. The fundamental plane is the equator and the positive X-axis points in the direction of vernal equinox. The Z-axis points in the direction of north pole.

To know the state of motion of the body, it is not necessary to know the position and the velocity of the object at every instant. The knowledge of the numerical values of the 6 constant parameters (orbital elements) that characterize the orbit is sufficient. The orbital elements are the parameters needed to specify the orbit uniquely. The commonly used set of orbit elements are

{a, e, i, Ω, ω, ν} (1.26)

The first two invariants a and e determine the orbit size and shape. The next three i, Ω and ω are the Euler angles that define the orbit plane orientation. Finally, the true anomaly ν specifies where the object is within the orbit trajectory at a specified time.

Semi-major axis (a): The semi-major axis describes an orbit’s size and is half of the distance between apogee and perigee on the ellipse. This is a significant measurement since it also equals the average radius, and thus is a measure of the mechanical energy of the orbiting object.


Eccentricity (e): Eccentricity measures the shape of an orbit and determines the positional relationship to the central body which occupies one of the foci. It is a ratio of the foci separation to the size (semi-major axis). which is calculated as

e = c

a (1.27)

where c is the distance between any one of foci to the center and a is semi major axis.

Inclination (i): The angle used to orient the orbital plane is called inclination angle. It is a measurement of the orbital plane’s tilt. This is an angular measurement from the equatorial plane to the orbital plane (0o ≤ i ≤ 180o), measured counter clockwise at the ascending node as shown in Figure 1.4. Inclination is utilized to define several general classes of orbits. Orbits with inclinations equal to 0o or 180o are equatorial orbits. If an orbit has an inclination of 90o, it is a polar orbit. Inclination orients the orbital plane with respect to the equatorial plane.

Figure 1.4: Inclination angle [21].

Right ascension of the ascending node (Ω): It is the measurement of the orbital planes’s rotation around the Earth. It is an angular measurement with in the equatorial plane from the first point of Aries eastward to the ascending node. The line made by the intersection of orbital plan and equatorial plan is called line of nodes. The orbital plan cuts the equator at two points called nodes. The node where satellite comes up is called ascending node and where satellite goes down is called descending node The first point of Aries is simply a fixed point in space.

Argument of perigee (ω): Inclination and right ascension fix the orbital plane in inertial space. The argument of perigee orients the orbit within the orbital plane. It is an angular measurement within the orbital plane from the ascending node to perigee in direction of satellite motion (0o ≤ ω ≤ 360o) as shown in Figure 1.6.

True anomaly (ν): Now that we have the size, shape and orientation of orbit established, the only thing left is to specify where exactly the satellite is on this orbit at some particular time.

Anomaly is simply another astronomer word for angle. True anomaly is an angular measurement that describes where the satellite is in its orbit at a specified time. It is simply an angle that marches in time from 0 to 360 degrees during one revolution. It is defined to be 0 degree at perigee, and 180 degrees at apogee. The summary of orbital elements is shown in Table 1.1.


Figure 1.5: Right ascension of ascending node [21].

Figure 1.6: Argument of perigee [21].

1.5 Euler’s Angles

The formulation of spacecraft attitude dynamics and control problems involves considerations of kinematics. In kinematics the orientation of a body which is in rotational motion is of primary interest. One scheme for orienting a rigid body to a desired attitude is called a body-axis rotation, which involves three successive rotations about the axes of the rotated body fixed reference. The Euler angles relate two coordinate systems having a common origin. The transformation from one coordinate system to the other is achieved by a series of two-dimensional rotations, which are performed about the coordinate system axes generated by the previous rotation step explained as following.

Consider Figure 1.7(a), where (X,Y,Z) is initial (fixed) coordinate system related to a rotating reference frame (x,y,z) by the Euler angle. At the start both coordinate systems (X,Y,Z) and (x(1),y(1),z(1)) are coincident.

In Figure 1.7(b) the first rotation is shown which involves Euler angle α. The rotating reference frame is rotated about the Z-axis of fixed reference frame through an angel α counterclockwise relative to X,Y,Z to give the new system (x(2),y(2),z(2)). It can be seen clearly from the Figure


Table 1.1: Orbital Elements

Element Name Description Definition Remarks

a Semi-major


orbit size half of the long axis of the ellipse

orbital period and energy de- pend on orbit size

e eccentricity orbit shape ratio of half the foci separation (c) to the semimajor axis

closed orbit: 0 ≤ e ≤ 1, open orbits:1 ≤ e

i inclination orbital plane’s tilt

angle between the orbital plane and equatorial plane, measured counterclockwise at the ascending node

equatorial: i = 0o or 180o, prograde: 0o ≤ i ≤ 90o, polar: i = 90o, ret- rograde: 90o ≤ i ≤ 180o

Ω right ascen-

sion of the ascending node

orbital plane’s rotation about the Earth

angle, mea-

sured eastward, from the vernal equinox to the ascending node

0o ≤ Ω ≤ 360o, undefined when i = 0o or 180o (equato- rial orbit)

ω argument of


orbit’s ori- entation in the orbital plane

angle, measured in the direction of satellite mo- tion from the as- cending node to perigee

0o ≤ ω ≤ 360o, undefined when i = 0o or 180o, or e = 0 (circular orbit)

ν true


satellite’s location in it’s orbit

angle, measured in the direction of satellite mo- tion, from perigee to the satellite’s location

0o ≤ ν ≤ 360o, undefined when e = 0 (circular orbit)

1.7(b) that this rotation mixes the coordinates along X and Y, while the coordinates along Z-axis remains unaffected.

The second rotation involves the Euler angle β. The (x(2),y(2),z(2)) axis system is rotated about the y(2) axis through an angle β counterclockwise to produce new coordinate system (x(3),y(3),z(3)) which results in mixing the coordinates along x(2) and z(2) while y(2) is unaf- fected as shown in Figure 1.7(c). This operation also generates a line of nodes parallel to the direction of y(2).

The final rotation involves the Euler angle γ. The (x(3),y(3),z(3)) axis system is rotated about the z(3) axis to generate the final coordinate system x,y,z. This mixes the coordinates along x(3) and y(3), while the coordinates along z(3) remains unaffected shown in Figure 1.7(d).

The usual ranges for α, β, γ are


(a) Initial position. (b) 1st rotation about Z-axis.

(c) 2nd rotation about Y-axis (d) 3rd rotation about Z-axis

Figure 1.7: Graphical representation for Euler angle [11].

ˆ 0 ≤ α ≤ 360

ˆ 0 ≤ β ≤ 180

ˆ 0 ≤ γ ≤ 360

As a final conclusion we can say that the angle β is simply the angle between the Z-axis of both coordinate systems. The angle α is the angle between the X-axis of the reference coordinate system and the projection of Z onto the X,Y plane and the γ is the angle between the Y-axis and the line of nodes.

1.6 Clohessy Wiltshire Equations

With the desire to place spacecrafts in formations comes the need to predict accurately and un- derstand the relative motion between the satellites. To describe this relative motion, researchers initially turned to Hill’s equations, also known as Cohessy-Wiltshire equations. Hill’s equations are a set of linearized differential equations that describe the relative motion of two spacecrafts in similar near circular orbits [31].

Models for spacecraft relative motion can be developed using different choices of relative motion coordinates and reference frames. The two main choices are geocentric inertial frame of reference and a chief centered local vertical local horizontal (LVLH) frame. The choices for the coordinates are inertial cartesian coordinates and LVLH cartesian coordinates. Relative orbits of interest


for formation flying are described in a rotating frame of reference attached to the chief satellite.

Hence, it is desirable to write down the relative motion dynamics in this frame of reference. The Hill-Clohessy-Wiltshire (HCW) equations are the simplest set of equations in LVLH coordinates.

These equations are commonly used for spacecraft formation studies and express the linearized motion of one satellite relative to a circularly orbiting reference point or chief location.

1.6.1 Particle Kinematics with Moving Frames

We know that it is simpler to define a particle position in terms of cylindrical or spherical coordinate system. The transport theorem allows to take the derivative of one vector with respect to one coordinate system, even though the vector itself has its components taken in another, possibly rotating coordinate system.

Transport theorem: Let N and O be two frames with a relative angular velocity vector ωO/N, and r be a generic vector then derivative of r in the N frame can be related to the derivative of r in the O frame as

Ndr dt =


dt + ωO/N× r (1.28)

The Eq(1.28) is a very fundamental and important equation and it is used every time kinematic equations are derived. The notation Ndtdr means the derivative of r taken with respect to the N reference frame. The derivation of transport theorem is not given here and can be found in many books related to orbital mechanics [28].

Now considering a more general problem in which the coordinate frames are free to translate, while frame orientation (defined through the three respective unit direction vectors) might be rotating. Let p be a generic particle in a three dimensional space. Assume two different frames A = {O, ˆa1, ˆa2, ˆa3} and B = n

O, ˆb1, ˆb2, ˆb3o

as shown in Figure 1.8. The position of O relative to O is given by vector R. Let r and ρ be the position vectors of particle p in the A and B reference frames, respectively. The angular velocity vector of frame B relative to frame A is given by ωB/A. The position vector r can be related as

r = R + ρ (1.29)

For simplicity we can use notation (vp)B to express the velocity vector of particle P with derivative taken relative to frame B



dt (ρ) (1.30)

The velocity vector (vp)A is given by (vp)A

Ad dt (r) =


dt (R + ρ) =

Ad dt (R) +


dt (ρ) (1.31)

Similarly, we can write


)A =


dt (R) (1.32)

where (vO)A is the velocity vector of origin O in frame A. Using transport theorem shown in Eq(1.28) and Eq(1.32) and writing Eq(1.31) as


Figure 1.8: Two coordinate frames with moving origins [28].

(vp)A= (vO

)A+ (vp)B+ ωB/A× ρ (1.33)

to find the acceleration (ap)A of particle p in the frame A. The derivative of Eq(1.33) is taken in the frame A and we get

(ap)A =


dt ((vp)A) =

Ad dt


)A+ (vp)B+ ωB/A× ρ

(1.34) Applying differentiation and using the transport theorem we get

(ap)A =

Ad dt(vO





dt ωB/A×ρ+ωB/A×(


dt (ρ)+ωB/A×ρ) (1.35) The acceleration of origin O in the frame A is defined to be



Ad dt(vO

)A (1.36)

and the acceleration of particle p in the frame B is



dt(vp)B (1.37)

Using the definitions in Eq(1.30) and Eq(1.36) and Eq(1.37) the particle p acceleration vector (ap)A can be written in more simplified form as

(ap)A = (aO

)A+ (ap)B+ ˙ωB/A× ρ + 2ωB/A× (vp)B+ ωB/A× (ωB/A× ρ) (1.38) Note that Eq(1.38) holds between any two reference frames. It is not necessary that A or B be inertially fixed. The vector components used in the various terms on the right hand side of Eq(1.38) can be taken along any choice of unit vectors.


1.6.2 General Relative Orbit Description

The relative distance between the spacecrafts is the key element in formation flying missions. The simplest type of spacecraft formation flying geometry is the leader follower type of formation flying. In leader-follower type of formation flying two spacecrafts are essentially in identical orbits, but are separated only by having different anomalies as shown in Figure 1.9. In case of circular orbit the spacecraft separation will remain fixed because both are always moving at the same orbital speed but if the orbit is elliptic, then the spacecraft separation will contract and expand, depending on whether the formation is approaching the orbit apoapsis or periapsis [28].

Figure 1.9: Illustration of leader follower type of two spacecraft formation [28].

The satellite about which all other satellite motions are referenced is called chief satellite. The remaining satellites are referred to as deputy satellites flying in formation with the chief satellite.

Note that it is not necessary that the chief position actually be occupied by a physical satellite.

Sometimes this chief position is used as an orbiting reference point about which the deputy satellites orbit.

The inertial chief position is expressed through the vector rc, while the deputy satellite position is given by rd. Now Hill’s coordinate frame is introduced to express how the relative orbit geometry is seen by the chief. The origin of Hill’s coordinate frame is at chief satellite position and its orientation is given by the vector triad {ˆor, ˆoθ, ˆoh} shown in Figure 1.9. The vector ˆor

is in the orbit radius direction, while ˆoh is parallel to the orbit momentum vector in the orbit normal direction. The vector ˆoθ completes the right hand coordinates system. Mathematically, these vectors are expressed as:



 ˆ or= rrc



oθ= ˆoh× ˆor ˆ

oh= hh


Another general type of spacecraft formation flying is shown in Figure 1.10 The relative orbit position vector ρ is expressed in O frame components as

ρ = (x, y, z)T (1.40)


Figure 1.10: General type of spacecraft formation with out-of-orbit plane relative motion [2].

Here, the various spacecrafts are on slightly different orbits that will satisfy some specific con- straints. These constraints ensure that the relative orbit is bounded and that the spacecraft will not drift apart.

1.6.3 Derivation of Hills Equation

To derive the relative equation of motion using cartesian coordinates in the rotating Hill frame, we write the deputy satellite position vector as:

rd= rc+ ρ = (rc+ x)ˆor+ yˆoθ+ ˆoh (1.41) where rc is the current orbit radius of the chief satellite as shown in Figure 1.10. The angular velocity vector of the rotating Hills frame (O) relative to the inertial frame (N) is given by

ωO/N = ˙f ˆoh (1.42)

where, ωO/N defines the angular velocity of the frame O relative to frame N and ˙f is the chief frame true anomaly. Now making use of Eq(1.38) and taking two derivatives with respect to the inertial frame, we get the deputy satellite acceleration vector as

¨rd = 


rc+ ¨x − 2 ˙y ˙f − ¨f y − ˙f2(rc+ x) ˆ or+


y + 2 ˙f ( ˙rc+ ˙x) + ¨f (rc+ x) − ˙f2y ˆ

oθ+ ¨zˆoh (1.43) The above kinematic expression can be simplified by making use of following identities. The chief orbit angular momentum magnitude is given by h = rc2f as h is constant for Keplerian˙ motion. Taking first derivative of h yields

˙h = 0 = 2rc˙rcf + r˙ c2f¨ (1.44) This orbit element constraint can be used to solve the true anomaly acceleration. From Eq(1.44) we get


f = −2¨ ˙rc

rcf˙ (1.45)

Further, we write the chief satellite position as rc = rcˆor. Taking two time derivatives with respect to the inertial frame and using the orbit equations of motion ¨r = −rµ3r, and knowing that ˆo = rrc as shown in Eq(1.39), we get the chief acceleration vector expressed as

¨rc =


rc− rc2 ˆ

or= −µ

rc3rc = −µ

rc2ˆor (1.46)

Equating vector components in Eq(1.46), the chief orbit radius acceleration is expressed as

¨rc = rc2− µ

rc2 = rc2

 1 − rc


(1.47) Putting Eq(1.47) and Eq(1.45) into Eq(1.43) gives

¨rd =

rc2− µ

rc2 + ¨x − 2 ˙y ˙f −

−2˙rc rf˙

y − ˙f2(rc+ x)

 ˆ or+


y + 2 ˙f ˙rc+ 2 ˙f ˙x +



(rc+ x) − ˙f2y


oθ+ ¨zˆoh (1.48) Simplifying the above equation, the deputy acceleration vector expression is reduced to


 x − 2 ˙¨ f

˙y − y˙rc


− x ˙f2− µ rc2

 ˆ or+

¨ y + 2 ˙f

˙x − x˙rc


− y ˙f2


oθ+ ¨zˆoh (1.49) The deputy satellite orbital equations of motion are given by

¨rd= −µ

rd3rd= −µ rd3 =

 rc+ x

y z

 (1.50)

with rd = p(rc+ x)2+ y2+ z2. Equating Eq(1.49) and Eq(1.50), we get the exact nonlinear relative equations of motion as







x − 2 ˙¨ f

˙y − yr˙rcc

− x ˙f2rµ2c = −rµ3


(rc+ x)

¨ y + 2 ˙f

˙x − xr˙rcc

− y ˙f2= −rµ3


y z = −¨ rµ3




The relative equations of motion in Eq(1.51) are valid for arbitrarily large relative orbits, and the chief orbit may be eccentric. If the relative orbit coordinates (x, y, z) are small compared to the chief orbit radius rc, then Eq(1.51) can be further simplified. As we know that

rd=p(rc+ x)2+ y2+ z2 =pr2c + x2+ y2+ z2+ 2rcx and multiplying and dividing this expression by pr2c, we obtain



r3d = µ

rcq 1 + 2xr


3 = µ r3c

 1 + 2x




Using binomial expansion to expand the term 1 + 2xr




 1 + 2x



= 1 +

−3 2

 2x rc


 1 − 3x


then Eq(1.53) becomes

µ rd3 = µ


 1 − 3x



We know from our knowledge of conic section that p = hµ2 and h = ˙f r2 and we can write f˙2= h2

r4 Using the relationship h2 = pµ, we get

2= pµ r4 The above equation can also be written as


rc3 = rc2

p (1.55)

We also know that

rc = p 1 + ecosf Multiplying both sides by ˙f2 and dividing by p, we get


p = f˙2

1 + ecosf (1.56)

Comparing Eq(1.55) and Eq(1.56) the term rµ3

c can be written in more useful form as µ

r3c = rc2

p = f˙2

1 + ecosf (1.57)


Note that the orbit elements shown in Eq(1.57) are chief orbit elements. Now neglecting higher order terms, and simplifying the right hand side of Eq(1.50), we obtain

−µ r3d =

 rc+ x

y z

≈ µ r3c

 1 − 3x


 rc+ x

y z

= µ r3c

rc− 2x y z

 (1.58)

Substituting Eq(1.58) into Eq(1.51) gives







x − 2 ˙¨ f

˙y − yr˙rcc

− x ˙f2rµ2c +rpcrc2+ rpc22x = 0

¨ y + 2 ˙f

˙x − xr˙rcc

− y ˙f2+rpc2y = 0


z +rpc2z = 0

Simplifying above equation yields relative orbit equation of motion assuming that x, y, z are small compared to chief radius rc, that is







x − 2 ˙¨ f

˙y − yr˙rcc

− x ˙f2

1 + 2rpc

rµ2c +rp2c2= 0

¨ y + 2 ˙f

˙x − xr˙rcc

− y ˙f2

1 − rpc

= 0


z +rpc2z = 0 From Eq(1.57), we deduce that rµ2

c = rp2cf . Therefore, we can further simplify as˙







x − 2 ˙¨ f

˙y − yr˙rcc

− x ˙f2

1 + 2rpc

= 0

¨ y + 2 ˙f

˙x − xr˙rcc

− y ˙f2

1 − rpc

= 0


z +rpc2z = 0


Using Eq(1.45) and Eq(1.57) and knowing that θ = ω + f , the general relative equation of motion are rewritten in the common form as







x − x¨  ˙θ2+ 2rµ3 c

− y ¨θ − 2 ˙y ˙θ = 0

y − x¨¨ θ + 2 ˙x ˙θ − y ˙θ2rµc3

= 0

¨ z +rµ3

cz = 0


If the chief satellite orbit is assumed to be circular, then e = 0, p = rc, and the chief orbit radius rc is constant. As we know that for a circular orbit the mean orbital rate n is equal to the true anomaly rate f , the relative equations of motion reduce to the simple form known as the Clohessy-Wiltshire (CW) equations shown as


and the direction that completes the right handed coordinate frame respectively. The motion along x, y, and z is also referred to as radial, along-track, and out-of-plane motion, respectively, and n is the mean motion of the reference spacecraft.

1.7 Lagrangian Points and Three Body Problem

The five equilibrium points, in a system of two large masses on which all the forces acting on the small mass placed between two large masses would cancel out are known as Lagrangian points as shown in Figure 1.11. According to Kepler’s laws we know that the closer a planet is to the sun, the faster it will move. So, if any spacecraft going around the sun in an orbit smaller than earth’s orbit will also soon overtake and move away, and will not keep a fixed position relative to earth.

Figure 1.11: Lagrangian points [12].

However, there is a loophole. Considering the system where the two large masses are sun and the earth and a spacecraft is placed between. Then, the earth’s gravitational force will pull the spacecraft in the opposite direction and cancels some of the pull of the sun. This will cause a weaker pull toward the sun and the spacecraft then needs less speed to maintain its orbit. If the distance is just right i.e. about a hundredth of the distance to the sun. Then spacecraft will keep its position between the sun and the earth and will need just one year to go around the sun. This point is called L1 as shown in Figure 1.11. This location is very good for monitoring


the sun. The solar wind reaches it about one hour before reaching earth. The spacecrafts at this point must have thrusters with large quantity of fuel because the position is unstable. If a spacecraft slips off L1, it will slowly drift away, and some correcting action are needed.

The L1, L2 and L3 point are unstable equilibrium positions, just like a ball on top of a hill. A little push and it starts moving away. A spacecraft at one of these points has to use thrusters frequently to remain in the same place. As seen from the sun the L4 and L5 points lie at 60 degrees ahead or behind earth in its orbit. Unlike the other Lagrange points, L4 and L5 are resistant to gravitational perturbations. Because of this stability, objects tend to accumulate in these points, such as dust and some asteroid-type objects. At L4 or L5, a spacecraft is truly stable, like a ball in a bowl i.e. when gently pushed away, it orbits the Lagrange point without drifting farther and farther, and will not require frequent use of thrusters.

Lets consider the circularly restricted three body problem where the two massive bodies are the sun and the earth. Moreover let’s consider a reference frame hˆiI, ˆjI, ˆkIi

in the center of mass of the earth-sun system. The ˆiI vector is always directed along the line connecting the center of earth and sun. As shown in Figure 1.12. M represents the earth and M represents the sun.

The third mass msat is a satellite which is very small as compared to the other two bodies. It is assumed that the two primary bodies rotate about their barycenter (composite center of mass) in circular orbit with a constant angular velocity, i.e.

n =pG (M1+ M2/D3) (1.62)

where M1and M2 are the masses of the sun and earth respectively and D is the constant distance between them. The position vector of the spacecraft relative to the barycenter is expressed in terms of basis vector n~iI,~jI, ~kIo

of a rotating reference frame with an angular velocity of n~k and with its origin at the barycenter as follows

R = X~i + Y ~j + Z~k~ (1.63)

With the given assumption the equation of motion of the third mass is shown in Eq(1.64).

Figure 1.12: Circular restricted three body problem [37].


ˆ µ = GM = 398601 km /sec

ˆ µ = GM = 1.3471 ×1011 km3/sec2

ˆ D is the distance between the earth and the sun (approx: 93,000,000 miles)

ˆ R = [X, Y, Z] represents the position of spacecraft in [iI, jI, kI] reference frame.

ˆ n = p(µ+ µ) /D3 = 4.06188 × 10−7 rad/s

The terms 2n ˙Y and 2n ˙X are the Coriolis accelerations (acceleration corresponding to Coriolis forces) and n2X and n2Y are centrifugal acceleration. The Coriolis force is a fictitious force exerted on a body when it moves in a rotating reference frame.

The Lagrangian points can be found out by setting all derivatives in Eq(1.64) to zero and solving for X,Y,Z. The concerned point for the analysis in this study is L1 point as shown in Figure 1.11. for simplicity lets define the position of the L1 point in the [iI, jI, kI] reference frame with [Xo, Yo, Zo] and expressing the position of m3 body with respect to the reference frame [iL1, jL1, kL1] centered in L1 with r = [x, y, z], then we have

X = Xo+ x, Y = Yo+ y, Z = Zo+ z (1.65) Substituting above definition in Eq(1.64) and expanding the right hand side of the equation with binomial theorem it is possible to derive the linearized equation of motion of a body around the L1 point. For detailed derivation see [37]. In particular, we have





x − 2n ˙y − (2σ + 1) x = 0¨ y + 2n ˙x + (σ − 1) y = 0¨


z + σz = 0


where σ is defined as

σ = µ

D3L1−S + µ DL1−E3


and DL1−S and DL1−E are the distances from the L1 point to the sun and the earth respectively.


Recently the concept of controlling the relative motion of spacecraft using electrostatic charging has been proposed. For tight spacecraft formations with separation distances ranging in the order of tens of meters, the Coulomb forces between the spacecraft promise to provide a very fuel and power efficient way of propulsion. The rocket principle is used in most of the spacecraft propulsion systems which have been built up to now and i.e. mass is ejected from a vehicle to affect momentum transfer and propulsive force. Varieties of this principle utilize chemical reactions to accelerate the mass as well as electromagnetic forces, however the thruster lifetime is fundamentally constrained by the amount of propellant available on board.

The spacecraft formation flying using Coulomb forces is a relatively new technology for spacecraft control, and may have application for a wide variety of mission objectives including collision avoidance, attitude control, and orbit perturbation correction. Coulomb-controlled formations appear ideally suited for close formation-flying in high earth orbits to perform wide field of view imaging missions using separated spacecraft interferometry. Formation flying on the order of tens of meters is very difficult using conventional ion propulsion methods, because the exhaust plumes will quickly interfere with the delicate on-board sensors. The Coulomb forces would allow the relative motion of satellites to be controlled with reduced contaminations. Furthermore, the fuel efficiency in Coulomb formation flying missions makes very long duration missions possible [29].

The Coulomb control principle is most easily explained by examining the interaction between two neighboring bodies capable of transferring electric charge. More detailed analysis of the physical processes will be discussed later. Consider two vehicles separated by a distance d in space. Initially, both spacecraft are electrically neutral, i.e., the amount of negative charge is equal to the amount of positive charge producing a net vehicle charge of zero and no interaction between the spacecraft. Now, allow one spacecraft to change its charge state through the emission of electrons. This is done via utilizing a cathode device or electron gun. If the electron beam is used to transfer an amount of negative charge, qsc, from spacecraft 1 (SC1) to spacecraft 2 (SC2), the net negative charge of SC2 will equal the net positive charge remaining on SC1, producing an attractive force between the spacecraft given by [16]

Fo= 1 4πǫo


d2 (2.1)

where ǫo is permittivity of free space (8.85418782 × 10−12m−3kg−1s4A2). The charge required to produce a 10 µN attractive force at a spacecraft separation of d = 10 m is qsc = 3.3 × 10−7C.

According to definition of ampere I = Q/t, with a 1 mA electron beam current this charge can be transferred in only 330 µsec.

The potential of the charged spherical spacecraft can be evaluated from Gauss’ law as Vsc = 1

4πǫo qsc

rsc (2.2)

where Vsc is the spacecraft potential in volts and rsc is the spacecraft radius. For example, if we consider a system where two spacecraft are driven to different potential level by transferring charge from one spacecraft to the other spacecraft then for a charge of qsc = 3.3 × 10−7 C and


the charge on second particle q2 is held fixed at 2µC.

Figure 2.1: Value of point charge q1 needed to apply specified electrostatic forces in vacuum, q2= 2µC, q1max= 2µC [23].

Spacecraft charging has been associated with negative impacts on satellite payloads. Arcs and other breakdown phenomena arising from differential charging are a danger to the sensitive electronics. Differential charging results when parts of one spacecraft are charged to different potentials relative to each other. Breakdowns can cause spontaneous interruption of payload functions. If spacecraft potential with respect to space is adjusted uniformly over the space- craft, then Vsc can be driven to large values (such as many kilo-volts), which might affect the experiments related to plasma measurements.

The Coulomb spacecraft formation flying concept envisions ejecting ions or electrons from the spacecraft. That thrusting phenomena is inherently different from that of a traditional ion thruster. Coulomb forces of tens to hundreds of micro-Newtons can be generated (at spacecraft separations of tens of meters) with a few milli-Watts of spacecraft power [22]. The Coulomb propulsion concept development requires no inherently new devices or technology. The tech- nological revolutionary nature of the Coulomb Satellite Formation (CSF) concept relies on an innovative integration of existing technologies and simple physical principles to provide a fuel





Relaterade ämnen :