previous chapter  


Chapter 4

Possible Operational Orbits

    In this chapter I will discuss the different orbits that were considered for ALP-SAT. To calculate them I introduce the theoretical background of orbit calculation. Mathematical formulas and procedures are evaluated by using different software tools which are described in appendix A.

4.1     Celestial Mechanics

    The theoretical framework of the theory of celestial mechanics is Newton's theory of gravity. Newton's theory was derived from Kepler's laws for the planetary motion. These are:
  1. The orbit of each planet is an ellipse with the Sun occupying one focus.
  2. The line joining the Sun to a planet sweeps out equal areas in equal intervals of time.
  3. A planet's orbital period is proportional to the mean distance between Sun and planet, raised to the power 3/2.
Johannes Kepler discovered them during the first twenty years of the seventeenth century. Only hundred years later did Newton develop his theory of gravity.
This theory states that the force of gravitational attraction between two point-like bodies of mass m1, m2 with a distance r from each other is equal to
         FG=G m1 m2 / r^2(4.1)
G is the universal constant of gravitation, having the value of 6.67×10-11 Nm^2/kg^2. The movement of a spacecraft in the vicinity of a central body follows this law. Further the gravitational field of a homogenous sphere is the same as the field of a point-like mass. Thus computation of the orbit can be restricted to the motion of two point-like objects in their mutual gravitational fields.
The solution for this so-called two-body problem, is the Keplerian orbit.
    Influences like atmospheric drag and radiation pressure play a minor role in celestial mechanics. The references for sections 4.1.1 to 4.1.4 are Nussbaumer 1989 chapter 2 and Fortescue and Stark 1995 chapter 4. In the following x = ¦x¦ for any vector x.

4.1.1     Solution of the two body problem

    The differential equation for the two-body problem is
where vec r1 and vec r2 are the cartesian coordinates of the bodies, m1 and m2 their masses and r their distance. Subtracting these equations brings the differential equation for the relative motion
         ddot(vec r)=-Gm vec r/r^3(4.4)
in which vec r = vec r2-vec r1 is the relative location and m = m1+m2 the total mass.
For the restricted two body problem where a small mass orbits a heavy mass equation (4.4) can be used with vec r being the position of the smaller body from the center of the heavier one and m the mass of the heavier one. This is the case for a spacecraft orbiting a celestial body. In this case the gravitational parameter u = Gm is introduced as given in table 4.1.
Table 4.1: Gravitational parameters of selected bodies.
    (from Fortescue and Stark 1995)
Bodyu ( km/s )

Jupiter    1.267×108

    Taking the cross product of (4.4) with vec r yields to
vec r times ddot(vec r)=0
And as d/dt(vec r x dot(vec r))=dot(vec r) x dot(vec r) + vec r x ddot(vec r) = 0 follows
         vec r x dot(vec r) = constant(4.5)
Therefore the angular momentum vec L = vec r×mdot vec r is constant. As vec r and dot vec r are perpendicular to the angular momentum the orbit is in a plane normal to the angular momentum. Be A the area over which the line joining the two bodies has swept up to a specific time. The change of this area is then equal to
         dA/dt=1/2 abs(vec r x dot vec r) = 1/2m abs(vec L) = const(4.6)
This implies Kepler's second law. The solution to equation (4.4) may be obtained by taking the cross-product of (4.4) with vec L
vec L x ddot(vec r) = Gm/r^3(vec L x vec r) = -Gm^2 d/dt ((vec r)/r)
and integrating once with respect to time.
This yields
         vec L x dot(vec r) = -Gm^2((vec r)/r + vec e)(4.7)
where vec e is the vector constant of integration called the eccentricity vector. The final equation is obtained by taking the dot product of (4.7) with vec r and using

vec r . (vec L x dot(vec r)) = vec L . (dot(vec r) x vec r) = - 1/m vec L^2
         r = L^2/Gm^3 . 1/(1+e cos theta) = p/(1+e cos theta)(4.8)
where & is the angle between vec e and vec r and is called true anomaly. Equation (4.8) is the formula of a conic section with the semi-latus rectum
and is therefore in agreement with Kepler's first law.
e is the eccentricity and determines the type of the conic. It's a circle when e=0, an ellipse for 0<e<1, a parabola when e=1 and a hyperbola when e>1. The eccentricity vector has the direction from the focus to the perifocus. The semi-latus rectum is the distance from the focal point to the curve measured perpendicular to the eccentricity vector (see graph below).

    The general form of an ellipse is shown in figure 4.1. The shape of the ellipse is characterized by the semi-major axis a and the eccentricity e. The symmetry axis parallel to the major axis is called line of apsides. Any point on the ellipse can be characterized by the true anomaly & or by the eccentric anomaly E. The radius of a specific point is given by equation (4.8) or by
         r=a(1-e cosE)(4.10)
For the semi-latus rectum p follows:

Graph of Ellipse
a=semi-mayor axis &=true anomaly
b=semi-minor axi E=eccentric anomaly
p=semi-latus rectum        rp=pericenter radius
ra=apocenter radius

Figure 4.1: Ellipses

    For the restricted problem of two bodies the heavier body occupies one focal point. The point closest to the main body is called periapsis or perifocus, the opposite point apoapsis or apofocus. The focal points normally carry the name of the central body; thus the perifocus is called perigee if the orbit is around Earth, perihelion if round the Sun etc. The distances from the focal point to the peri- and apofocus are given by
         rp=p/(1+e)=a(1-e)and   r_a=p/(1-e)=a(1+e)(4.12)
    Equating the energy and the angular momentum at peri- and apofocus leads to a total energy of
         Etot=-(G m1 m2)/(2a)(4.13)
The orbital period T is calculated by the formula A = pia2sqrt(1-e^2), for the area of the ellipse, and (4.6)
from which follows
         T=2\pi sqrt(a^3/Gm)(4.14)
which implies Kepler's third law.
    To find the position versus time relationship the differentiation of (4.10) and (4.8) gives
ae sinE dE/dt = L^2/Gm^3 e sin theta)/((e cos theta+1)^2) dtheta/dt

Using equation (4.6) and r  sin & = asqrt(1-e^2)  sinE yields

dE/dt=1/r sqrt(Gm/a)
Separating the variables and integration gives Kepler's equation
which is:
M=tp sqrt(Gm/a^3)
is the mean anomaly and tp is the elapsed time since any previous perifocus. Kepler's equation makes it possible to compute the movement of the two bodies in time.

    For the restricted problem of two bodies the energy per unit mass of the smaller body e is
To evaluate the velocity of the smaller body V one can use the following equation
         varepsilon=1/2 V^2-mu/r(4.17)
which is also valid for a parabolic or hyperbolic curve. The velocity for a circular orbit (e=0, r=a) follows as
The apofocus and perifocus velocity are given by
         Vp=sqrt(mu/a (1+e)/(1-e))and  Va=sqrt(mu/a (1-e)/(1+e))(4.19)

Figure 4.2: ARIANE 5 GTO Orbit around Earth: Orbit trajectory with velocity vectors and absolute velocity to the x-axis.

    The evolution of the velocity along an orbit is shown in Figure 4.2. The figure shows the feature of equation (4.17). The velocity rises to a maximum at the apogee while being low for most of the trajectory.

    In the case of a parabolic trajectory the perigee height rp is
As the energy in this case is zero the parabolic trajectory represents a boundary condition between a captive, orbiting satellite and an escaping one. The energy for a parabolic trajectory is the minimum energy to escape the central body. Consequently the minimal velocity required to escape the central body follows from equation (4.17)
Vsc is called escape velocity. The escape velocity from the Earth surface is 11.12  km/s, from the solar surface 617  km/s and from the Sun at the Earth orbit 42.1  km/s.
Figure 4.3: Hyperbola

    From equation (4.11) we see that hyperbolic curves have a negative semi-major axis. The hyperbolic orbit is open and applies to swing-by manoeuvres and escape trajectories. The shape of a hyperbola is given in figure 4.3. The time evolution is calculated in a similar way than for the ellipse. A hyperbolic eccentric anomaly
F is introduced which satisfies
         r=a (1-e cosh F)(4.21)
An analogous to Kepler's equation is found
         M=e sinh F - F(4.22)
where in this case
M=tp sqrt(mu/((-a)^3))
The energy for a hyperbolic orbit is positive. Therefore the motion at very large distances is dominated by kinetic energy. As r becomes very large, equation (4.17) shows that the speed becomes
At very large distances the orbit approaches the asymptotes which are characterized by their direction angle
and the distance to the focal point
d=-a sqrt(e^2-1)
which are obtained by equation (4.8) and figure 4.3. The angle through which the trajectory is deflected d is given by
         delta=2 sin^(-1)(1/e)(4.24)
The perifocus radius follows rp = a (1-e) as for the other curves.

4.1.2     Specifying the orbit

At a given time a satellite's movement can be characterized by the location and the velocity in a reference frame. For ALP-SAT two different reference frames are used, shown in figure 4.4. Both of them are fixed relative to the stars. For the frames the vernal equinox is used, the point where the Sun crosses the Earth's equatorial plane from south to north (21 March). Due to lunisolar perturbations of the Earth's rotation axis this point moves along the equator at 50,4" per year. The Earth equator frame has the X-axis pointed towards the vernal equinox and the Z-axis parallel to the Earth's spin axis towards north. The Earth ecliptic frame has the same X-axis direction but the Z-axis towards the ecliptical north pole. The Y-axis makes up a 'right-handed' orthogonal set with the X- and Z-axes. For both frames the origin is in the center of the Earth.
Coordinate frames
Figure 4.4: Coordinate reference frames

Celestial Coordinates
Figure 4.5: Position of an object on the celestial sphere

The celestial sphere is the sphere with infinite radius centered on Earth. As shown in figure 4.5 the position of any object on the celestial sphere can be given by its right ascension a and its declination d.
The orbit of the satellite is generally defined by the six orbital elements. The orbital plane intersects with the equatorial plane in the line of nodes. The ascending and descending nodes (point where the spacecraft crosses the equatorial plane) are on this line. The orbital plane is given by its inclination i to the equatorial plane and by the angle between the vernal equinox and the line of nodes called Right ascension of the ascending node O. The location of the curve in the orbital plane is given by the argument of perigee w which is angle between the line of nodes and the line of apsides. The shape of the curve is characterized by the eccentricity e and by the semi-major axis a. The sixth element gives the location of the spacecraft on the orbit. It can be the true anomaly & or the time since the last perigee tp. Figure 4.6 shows these angles and the location of the orbital plane. The orbital parameters O, i, and w can be interpreted as Euler angles for the positioning of the orbit in space.
Orbital Parameters Figure 4.6: Orbital parameters

4.1.3     Orbit calculation

The first approach to perform orbital computations is by supposing that every orbit is a keplerian one. If a spacecraft passes the influence sphere of different central bodies the method of patched conics can be used. There the orbit is a keplerian one around the body with the strongest gravitational attraction at the current satellite position. The final flight path is obtained by patching together the trajectories resulting from each conical section. The zone of strongest gravitational attraction is called sphere of influence. In the case of a smaller body circling a bigger one the sphere of influence of the smaller body can be computed with
         rs=rb((m_small/(m_big))^2/5 (1+3 cos^2 beta)^-1/10(4.25)
where rs is the distance to the smaller body and b the angle between this radius vector and the direction to the bigger body. The zone of influence of Earth has a radius of 805000 km to 925000 km in the gravitational field of the Sun. The Moon has a zone of influence with a radius from 58000 km to 66000 km in the terrestrial field. For patched conics, orbital manoeuvres are considered to happen in one point, changing the spacecraft velocity and the orbit instantly.
    There are two fundamentally different techniques to calculate an orbit more precisely. Special perturbation techniques deal with the direct numerical integration of the equation of motion including all relevant perturbing accelerations. These techniques are straightforward and applicable to any number of bodies and perturbing accelerations. However they require considerable computational effort and don't lead to generally usable formulas. General perturbation techniques on the other hand provide an analytical description of the variation with time of the orbit. Although being more complicated, these techniques provide a better physical insight and allow the development of analytical expressions that require a minimum of computational effort.
One of them is the method of the variation of orbital elements. The idea is to compute the keplerian orbital elements from the present location and velocity for a keplerian orbit. These elements are called osculating elements. Then calculate the perturbing forces upon the satellite and consider the variation of the osculating elements to these forces. This method assumes that the perturbation forces are much smaller than the force from the central body.
Spacecrafts are under influence of many forces disturbing the 1/r potential. These include gravitation effects, atmospheric drag and radiation pressure.
The terrestrial field differs from the 1/r potential. Earth's gravitational field is best described by the following equation
         U(r,theta,phi)=mu/r {-1+sum(n=2 to infty) [(RE/r)^n Jn Pn0 (cos theta)         
         +sum(m=1,n) (RE/r)^n (Cnm cos(m phi)+Snm sin(m phi) ) Pnm (cos theta) ] }(4.26)
    where U(r,&,phi) is the gravitational potential at distance r from the center of the Earth and & and phi are the latitude and longitude. Pnm are Legendre polynomials.
Jn,  Cnm,  and  Snm are the harmonic coefficients that describe the mass distribution of the Earth. The terms in equation (4.26) decrease with (RE/r)n. Table 4.2 gives the terrestrial harmonic coefficients of lower order.
Table 4.2: Harmonic coefficients for Earth (from Fortescue and Stark 1995)
J21082.6×10-6  C210  S210
J3-2.53×10-6  C221.57×10-6  S22-0.904×106
J4-1.62×10-6  C312.19×10-6  S310.27×10-6

The J2 term is larger than the others. It results from the Earth's oblateness and dominates the gravitational perturbation of the Earth. The Earth's oblateness has two main effects on an orbit around Earth. The regression of the line of nodes is caused by the torque induced on the orbit. Per orbit the right ascension of ascending node O decreases by
         DeltaOmega=-(3 pi J2 RE^2)/p^2 cos i(4.27)
where p is the semi-latus rectum.
The other effect of the Earth's oblateness is the precession of the line of apsides. Caused by the stronger attraction above the equator the increase in argument of perigee w is given by
         Delta omega=3 pi (J2 RE^2)/p^2 (2-5/2 sin(i)^2)(4.28)
This variation is zero at an inclination of ~ 63.4°, a fact used by Russian communication satellites.
    Atmospheric drag decreases the satellite energy and thus the orbital height. This decrease has to be watched carefully as it may eventually cause atmospheric reentry and consequent loss of the supervised satellite. The force parallel to the velocity vector is best given by
         vec(F)D=1/2  rho S CD Vr^2 -(vec(Vr))/Vr (4.29)
where Vr is the velocity vector of the spacecraft relative to the atmosphere, rho the atmospheric density, S a reference area for the satellite and CD the satellite's coefficient of drag referred to the reference area. The drag is the strongest at the perigee when the altitude is the lowest and the satellite has the biggest velocity.
    A satellite is also influenced by the gravitational pull from other bodies.
Orbits around Earth are generally equally influenced by the Moon and the Sun. These perturbations are termed lunisolar perturbations. The force of Jupiter, the next important body, is five orders of magnitude lower. The disturbing acceleration is given by
         ad= mu(d) sqrt(vec(R) . vec(R))(4.30)
vec(R)= vec(r)sd/rsd^3 - vec(r)d/rd^3
Here ud is the gravitational parameter of the disturbing body, vec rsd is the vector from the satellite to the disturbing body and vec rd is the vector from the satellites central body to the disturbing body.
    A spacecraft moving through the solar system is exposed to a force from the impact of solar radiation on it's surface. Electromagnetic radiation carries a momentum that exerts a pressure
on any surface where F is the energy flux and c the speed of light. The irradiance pressure from other bodies than the Sun is insignificant.
The solar radiation pressure (SRP) for the environment of the Earth is PE = 4.7×106 N/m^2. The effect of the SRP on a spacecraft can be expressed by the force per unit spacecraft mass that is
         f(srp)=s A/m P(E) (a(o)/r(o))^2(4.31)
in which A/m is the spacecraft's area to mass ratio, ro and ao the spacecraft's and the Earth's mean distance from the Sun. s is a constant value between 0 and 2 depending on the reflective properties of the spacecraft's surface.
graph comparing disturbances
Figure 4.6: Comparison of the disturbing accelerations for Earth orbits
    (from Fortescue and Stark 1995)

    Figure 4.6 gives the a comparison of the disturbing accelerations in the vicinity of Earth. The logarithm of the forces normalized with 1 g is shown. The dominant force is the primary inverse square law gravity field, followed by the J2 perturbation of the Earth's oblateness. The higher order terms are only important for low orbit operations, for example rendezvous strategies. For the surface forces a A/m ratio of 0.005 m^2/kg was taken. The atmospheric density depends on the solar activity (see section 3.1.3). Therefore the curve for the drag is very uncertain and can vary up to one order of magnitude.

4.1.4     Restricted Three Body Problem

The study of a small mass in the gravitational field of two primary bodies is called restricted three-body problem or RTBP. The two primary bodies with masses m1 > m2 move along circular orbits around their mutual center of mass. A third body with negligible mass moves in space under the influence of the two primaries, without affecting their motion. A synodical system of coordinates is introduced as a rotating system with the origin at the center of mass. The x-axis is in the direction from m2 to m1, the z-axis in the direction of the momentum vector of the two bodies and the y-axis to complete the right handed system.
    The equation of motion for the third body is
         ddot(x)-2 dot(y)= dOmega/dx         
         ddot(y)+2 dot(x)dOmega/dy(4.32)
         ddot(z) = dOmega/dz         
                   Omega=(x^2+y^2\right)/2 + (1-mu)/r1 + mu/r2         
         (r1)^2 =  (x-mu)^2+y^2+z^2         
         (r2)^2 =  (x-mu+1)^2+y^2+z^2         
         mu = m1/(m1+m2)         
The unit of length is chosen as the distance between the two primaries where the unit of time is selected such that the angular velocity is one. The unit of mass satisfies G (M1+M2) = 1.
By setting velocity and acceleration to zero five points are obtained within the frame of reference at which stationary body will be at equilibrium. These Lagrangian or Liberation points, denoted L1 to L5, are shown in figure 4.7. The co-linear points L1, L2 and L3 lie on the line joining the primaries and are positions of unstable equilibrium. The trialgular liberation points L4 and L5 are on the x-y plane and form equilateral triangles with the primaries. The equilibrium at these points is stable.
Location of the Liberation Points
Figure 4.7: RTBP and Liberation Points (u=0.2)

    In the restricted tree body problem for the Earth Sun system the L1 point is 1.5 million km towards and L2 is 1.5 million km from Earth away from the Sun. The L3 point is behind the Sun at 150 million km which is the Sun Earth distance.

4.1.5     Orbital manoeuvres and propulsion

    Orbital transfer manoeuvres aim at reaching an orbit from a given one with a minimum of additional energy. The characteristical value for a propulsion manoeuvre is
         D vec V = vec Vend - vec Vbeg(4.33)
where vec Vbeg is the velocity at the beginning and vec Vend at the end of the manoeuvre. Immediate transfer from one to another orbit can only happen where the two orbits intersect and at least two manoeuvres are required for a transfer between two nonintersecting orbits. The increase of the spacecraft energy is best at the perifocus where the spacecraft has the biggest velocity. A change in the orbital plane is best performed at the apofocus, which is the point of minimum velocity. Figure 4.8 shows the variation of the eccentricity de for a small thrust dv in the direction of the current velocity versus the eccentric anomaly E. The shape of the curves does not vary with the length of the semi-major axis. The curves were computed with the formulas presented in appendix A using Matlab. It can be seen that the eccentricity increases most for a perigee manoeuvre and decreases most for a apogee manoeuvre (only if the velocity increase has the same direction than the current velocity).
infinitesimal change in eccentricity
Figure 4.8: change in eccentricity de for a minimum thrust dv along an orbit.

    To get a notion for the energy required for a plane change let's consider the following: Injection into geostationary Earth orbit from a geostationary transfer orbit with inclination 7° (launch from Kourou) needs a DV of 1.47km/s where the same manoeuvre from a geostationary transfer orbit with inclination 28° (launch from Kennedy Space Center) needs a DV of 1.80km/s. This calculation was made with (4.18), (4.19) and (4.27) and the vector composition shown below.
    If a change in orbital shape and a plane change have to be performed it is better to do them in one manoeuvre.
The transfer between two coplanar circular orbits that requires a minimum of energy is the Hohmann transfer as shown in figure 4.9. It uses two manoeuvres the first to bring the spacecraft on a elliptical orbit with perifocus in the smaller orbit and apofocus on the larger one. After half a orbit period the second manoeuvre is performed to place the spacecraft on the final orbit as shown in figure . The orbit is the most effective as the eccentricity decreases most for an apofocus manoeuvre and the semi-major axis increases most for a perifocus manoeuvre.
Hohmann transfer
Figure 4.9: Hohmann transfer

The idea behind a swing-by manoeuvre is to gather speed relative to a primary body by passing through the gravitational field of a second one that orbits the primary. It is a powerful technique to change the spacecraft velocity relative to the primary body. The trajectory in the field of the secondary body is a hyperbole. By passing behind the secondary body (relative to its motion) increases the spacecraft velocity whereas passing in front of the secondary body decreases the velocity. A swing-by manoeuvre can also be used to change the orbital plane. Swing-by's for velocity increase were performed by the Voyager 1 and 2 probes. The Ulysses satellite was inserted into its polar orbit around the Sun by a swing-by at Jupiter. (see also the Ulysses homepage at ESA)
    Another way to change an orbit is by progressive change by perturbations. But this method was never used as it requires a lot of numerical computation.

    During its life a spacecraft needs different types of propulsion. The biggest thrust (DV~10km/s) is required for the launch vehicle to lift the spacecraft from the Earth's surface into a low Earth orbit just above the atmosphere. Manoeuvres follow to bring the spacecraft into the desired operational orbit. These manoeuvres which require a DV between 2 and 8km/s are performed by the launch vehicle, by the on board engine or by both. Finally the orbit and the spacecraft attitude needs to be maintained which needs a DV from 0.03 to 0.3km/s.
    A rocket engine is characterized by its thrust
F=dm/dt Ve
which is the force exerted upon the spacecraft by the engine. Ve is the effective exhaust velocity which combines the thrust by the outflow velocity and by the difference between the exhaust pressure and the outside pressure. Neglecting external forces the differential equation for the spacecraft motion follows
M dV/dt=dm/dt Ve
with dot m being the exhaust mass flow rate. The solution is the Tsiolkovsky equation
         DV=Ve ln(Mbeg/Mend)(4.34)
where Mbeg is the mass at the beginning and Mend at the end. The specific impulse
Isp is defined by the total impulse per unit propellant weight consumed and is given by
         Isp=F/(dot m go)=Ve/go(4.35)
for a constant thrust and exhaust mass flow rate. Here go=9.81m/s2 is the gravitational acceleration at the Earth's surface.
    Different types of thrust engines were conceived but up to now only chemical rockets are widely used. Solid propellant engines have a specific impulse from 200 to 260 s. They have a simple design, but once ignited the combustion proceeds at a given thrust until all the propellant is consumed. The liquid propellant engine on the other hand is more complex in design but flexible in burn duration and thrust level. Monopropellant engines have a specific impulse about 200 s whereas bi-propellant engines have a Isp of 300 to 400 s. Bi-propellant engines are nearly twice as complex as monopropellant engines.
    For orbital manoeuvres the minimum energy consumption is achieved by giving the spacecraft the necessary DV within an instant at a certain point.
As engine burns take a certain time a bigger burn is required that in this theoretical case. This loss of thrust energy is called gravity loss and is given by the formula
         gravity loss=(DVreal-DVideal)/DVideal 100    [%](4.36)


next section