MISN-0-104

THE NUMEROV ALGORITHM FOR SATELLITE ORBITS

THE NUMEROV ALGORITHM FOR SATELLITE ORBITS by Peter Signell

x(t)

x4
1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 2. The Equation of Motion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 a. Newtons Second Law Plus the Law of Gravitation . . . . . . .1 b. The Equation in Cartesian Coordinates . . . . . . . . . . . . . . . . . .2

4D

3. Second-Order Numerov Solution . . . . . . . . . . . . . . . . . . . . . . . 3 a. Net-Point Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .3 b. Finite Dierence Approximation . . . . . . . . . . . . . . . . . . . . . . . . . 4 c. The Recurrence Relations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 d. Specifying Initial Position and Velocity . . . . . . . . . . . . . . . . . . 4 4. Choice of Step Size . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 5. Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 Acknowledgments. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .5

x n+1 = a n x n + b n x n-1

Project PHYSNET Physics Bldg. Michigan State University East Lansing, MI

ID Sheet: MISN-0-104 THIS IS A DEVELOPMENTAL-STAGE PUBLICATION OF PROJECT PHYSNET The goal of our project is to assist a network of educators and scientists in transferring physics from one person to another. We support manuscript processing and distribution, along with communication and information systems. We also work with employers to identify basic scientic skills as well as physics topics that are needed in science and technology. A number of our publications are aimed at assisting users in acquiring such skills. Our publications are designed: (i) to be updated quickly in response to eld tests and new scientic developments; (ii) to be used in both classroom and professional settings; (iii) to show the prerequisite dependencies existing among the various chunks of physics knowledge and skill, as a guide both to mental organization and to use of the materials; and (iv) to be adapted quickly to specic user needs ranging from single-skill instruction to complete custom textbooks. New authors, reviewers and eld testers are welcome. PROJECT STAFF Andrew Schnepp Eugene Kales Peter Signell Webmaster Graphics Project Director

Title: The Numerov Algorithm for Satellite Orbits Author: Peter Signell, Michigan State University Version: 2/1/2000 Length: 1 hr; 12 pages Input Skills: 1. State the 3-point nite dierence approximation to a second derivative (MISN-0-4). 2. State Newtons Second Law (MISN-0-14). 3. State Newtons Law of Gravitation (MISN-0-101). Output Skills (Knowledge): K1. Derive the recurrence relation for the Numerov Algorithm, to second order and in two dimensions, in a form suitable for use in obtaining satellite trajectories numerically and showing all steps in the derivation. K2. Derive equations for insertion of initial position and velocity in the Numerov Algorithm and communicate a method of obtaining a particular desired accuracy. Post-Options: 1. Orbits in an Inverse Square Law Force Field: A Computer Project (MISN-0-105). 2. Computer Algorithm For the Damped Driven Oscillator (MISN0-39). Evaluation: Stage 0

ADVISORY COMMITTEE D. Alan Bromley E. Leonard Jossem A. A. Strassenburg Yale University The Ohio State University S. U. N. Y., Stony Brook

Views expressed in a module are those of the module author(s) and are not necessarily those of other project participants. c 2001, Peter Signell for Project PHYSNET, Physics-Astronomy Bldg., Mich. State Univ., E. Lansing, MI 48824; (517) 355-3784. For our liberal use policies see: http://www.physnet.org/home/modules/license.html.

MISN-0-104

MISN-0-104

THE NUMEROV ALGORITHM FOR SATELLITE ORBITS by Peter Signell

y-axis ^ r y ` r q x x-axis
Figure 1. The vectors used to dene a particles position. time t.3 The origin of the coordinate system has been placed at the force center.4 2b. The Equation in Cartesian Coordinates. We can write the above vector equation in Cartesian component form by writing r in terms of Cartesian unit vectors x and y : r(t) = x(t) x + y(t) y . The reason for going over to Cartesian coordinates is that x and y are independent of time whereas r(t) is not.That is, x and y form a time independent reference system. Then taking components and writing the magnitude of r in terms of x and y , |r(t)|2 = x2 (t) + y 2 (t), we get the x equation: x = G mE [x2 (t) + y 2 (t)]
3/2

1. Introduction
Many problems in science and engineering cannot be solved in terms of known functions,even when the underlying equation is known. Such a problem is the position as a function of time along a satellite orbit or, for that matter, along the earths orbit around the sun. For such cases, one must resort to approximate numerical techniques; one such technique is examined in this module.

2. The Equation of Motion

2a. Newtons Second Law Plus the Law of Gravitation. We will numerically determine the trajectory of a satellite in the eld of a single gravitational source, such as the earth. The equation governing such motion can be obtained by combining Newtons Second Law,1 F =m d r(t) mr (t), dt2
2

with Newtons Law of Gravitation,2 F = We get: r (t) = GmmE r(t). r2 (t)

x(t) f (t) x(t).

GmE GmE r(t) = 3 r(t). r2 (t) r (t)

Similarly for the y equation: y (t) = f (t) y(t).

3 See

Here r(t) is the time changing position vector of the mass m and r(t) is the radially outward pointing unit vector at the position of the mass at
4 For

1 The number of primes will indicate the order of the derivative with respect to time in this module. 2 See Newtons Law of Gravitation (MISN-O-101).

Kinematics: Circular Motion (MISN-O-9). most satellite problems the center of the Earth can be taken to be a xed force center. For the motion of the moon, however, that is insucient due to its large mass. The methods presented here still apply, but the Earths mass becomes replaced by the mass of the total system. Other quantities require careful interpretation. For further details, see Two Body Kinematics and Dynamics (MISN-O-45).

MISN-0-104

MISN-0-104

x(t 4) x(t)

x(t)

3b. Finite Dierence Approximation. We will now connect the consecutive values of x and y by using the nite dierence approximation to any functions gs second derivative: gn (gn+1 2gn + gn1 ) /2 . (2)

Equation (2) is the result of truncating the following series after the second term:6 iv gn+1 + gn1 = 2gn + 2 gn + (4 /12)gn + . . . (3)

D 2D 3D t 4

7D

t 12

Figure 2. A function x(t), specied at equally-spaced values of t.

The truncation is justiable to the extent that numerical values of Delta will be small so that 4 will be insignicant compared to 2 . Keeping that in mind, we apply Eq. (2) to the x and y functions of Eq. (1) and get, with a little rearrangement: xn+1 = (2 + 2 fn )xn xn1 . Similarly, the y-coordinate relation is: yn+1 = (2 + 2 fn )yn yn1 . (5) (4)

3. Second-Order Numerov Solution

3a. Net-Point Notation. We will here use the Numerov method of solving equations involving derivatives. This method has the advantage of being easily understandable, although it is less appropriate for orbit problems than are some other less understandable methods.5 In the Numerov method we deal with the solution functions, x(t) and y(t), as a series of numbers at net-point times that are integrally spaced: tn n . This is illustrated in Figure 2 for x(t). We now introduce a more succinct notation for the coordinate positions at the net-point times: xn x(tn ) x(n), yn y(tn ) y(n), and our equations become: xn yn where:
2 fn G mE x2 + yn n
5 The

Such truncation of the series at the second term results in an algorithm that is referred to as being of second order. Derive Eq. (5). 3c. The Recurrence Relations. Remembering that fn = 2 GmE (x2 + yn )3/2 , we see that knowledge of x0 , y0 , x1 , and y1 will n enable us to compute x2 and y2 . Then, using x1 , y1 , x2 , and y2 , we can compute x3 , and y3 . This process can be repeated until x and y are known at any time you wish. The relations in Eqs. (4)-(5) are thus called recurrence relations. The only remaining problem is the specication of x0 , y0 , xl and y1. 3d. Specifying Initial Position and Velocity. In order to specify the trajectory, we usually nd it most practical to specify position and velocity at some particular time. This means that we specify: position: x0 , y0 . In order to convert x0 , y0 , x0 , y0 , into x0, y0, x1 , velocity: x0 , y0 . and y1 , we subtract Eq.(2) at t = 0 to obtain: x1 x1 = 2x0 ,
6 See

= f n xn , = f n yn ,
3/2

(1)

Numerov method tends to allow errors to accumulate. This is all right if the errors are both positive and negative so they tend to cancel out. However, for orbit problems the errors tend to all have the same sign so a method like the Runge-Kutta is better.

(6)

Taylors Series for the Expansion of a Function About a Point (MISN-O-4).

MISN-0-104 while the t = 0 sum equations (4)-(5) are: x1 = (2 + 2 f0 ) x0 x1 . The unknown x1 can be eliminated from Eq. (6) to give: x1 = (1 + 22 f0 ) x0 + x0 . Similarly, y1 = (1 + 22 f0 ) y0 + y0 .

MISN-0-104

ME-1

(7)

MODEL EXAM
1. See Output Skills K1-K2 in this modules ID Sheet. One or both of these skills may be on the actual exam.

(8) (9)

To summarize, we specify the initial position and velocity components, x0 , y0 , x0 , y0 , and then use Eqs. (8) to obtain x1 and y1 . These are combined with x0 and y0 to start the recurrence relations (4)-(5).7

4. Choice of Step Size

Finally, how do you know what size time interval to use? You could attempt to assess the importance of successive terms in Eq. (3), but a more reliable method is to decrease until the predicted trajectory stabilizes. That is, until it does not change signicantly when is made even smaller.

Brief Answers: 1. See this modules text.

5. Summary
The algorithm, then, consists of: 1. a recurrence relation; 2. two initial conditions; and 3. a method of assuring desired accuracy.

Acknowledgments
I would like to thank Dr. N. R. Yoder for the many insights I gained in this area while collaborating with him. Preparation of this module was supported in part by the National Science Foundation, Division of Science Education Development and Research, through Grant #SED 74-20088 to Michigan State University.
7 See Orbits in an Inverse Square Law Force Field: A Computer Project MISN0-105), wherein the presently derived algorithm is used to calculate orbits and trajectories.

10

11

12