The determination of orbits by using Picard iteration is reported. This is a direct extension of the classical method of Picard that has been used in finding approximate solutions of nonlinear differential equations for a variety of problems. The application of the Picard method of successive approximations to the initial value and the two point boundary value problems is given

Feagin, T.

Mikkilineni, R. P.

English

NASA Technical Reports Server

THE DETERMINATION OF ORBITS USING 
PICARD ITERATION 
Rajendra Prasad Mikkilineni and Terry Feagin 
University of Tennessee Space Institute 
Tullahoma, Tennessee 
The topic of this presentation is the determination of orbits using Picard iteration, which is 
a direct extension of the classical method of Picard that has been used in finding approxi- 
mate solutions of nonlinear differential equations for a variety of problems. The application 
of the Picard method of successive approximations to the initial value and the two-point 
boundary value problems is given below. 
The initial value problem, 
i = F (z, t) 
n = 1,2, .  . . 
is solved ' ~ y  means of the iterative scheme, 
i n = F ( z n - 1 4  
zn (t,) = c ,  
whereby the nth approximation to the state vector z is computed from the (n - 1) approxi- 
mation (beginning with some initial approximation, z,(t) j. The computation involves per- 
forming a simple integration or quadrature at  edch stage of the iterative process. The con- 
stant of integration is chosen such that each iterate satisfies the initial condition, z (t ,  ) = c. 
Similarly, the two-point boundary value problem, 
i = F (z, t) 
n = 1 , 2 , .  . . 
I 
https://ntrs.nasa.gov/search.jsp?R=19760003080 2020-03-22T19:28:22+00:00Z
where g denotes certain constraints or boundary conditions that must be satisfied by the 
solution at t, and t,. 
This paper presents an investigation of the suitability of this type of iterative scheme for tile 
problem of estimating orbits. In the estimation problem, we again have to solve the differ. 
entia1 equation for the state. However, we are now given a set of imprecise observatiuns, 
y(t, ) iz, which are, in general, nonlinear functions of the state. The problem is to find 
a z(t) which comes close (in some sense) to satisfying the observations. The estimation 
problem can thus be written as 
t = F (z, t) 
subject to the condition that z(t) minimizes some function Q, where 
For instance, Q could be the sum of the ;+ares of the residuals (the residual being the 
difference between the observed apd me computed values of the observations). Again, !he 
iteration is set up. 
now choosing the constant oi integration for z,(t) such that 
is minimized. As before, each stage of the iteration requires that a quadrature is performed. 
In order to perform the integration, we have chosen to approximate the right-hand side, 
F(Z,-~, t), by 2 series of Chebyshev polynomials, Tj(7), and to integrate term-5y-term. 
Thus, the approximation 
F iz,., , t) 2 ' bjTj(7) 
j= 0 
is made where 7 = 1 -. 2 t/(t, - t, is normalized time. The coefficients of the series are 
determined readily using the orthqgonality properties from the relation 
where i = cos (in/N). The integration can be performed by manipulating these coefficients, 
4 
giving 
i=o 
2 
where 
I 1 
-bj+l)/2j forj=1,2,. ..N 
and a, = the arbitrary constant of integration that is determined from the requirement that 
Q be minirnized. 
We will now consider the application of this technique to an orbit estimation proh!.:rr., 
namely, the determination of the orbit of an earth satellite from observations o! the : ange 
and range rate by several tracking stations. 
There are no external forces except the gravity field of the earth. The equations o! motion 
for the satellite are written as a set of six, first-order equations for the six state variables: 
x, y, z, u, v, and w. For instance, %(ti) represents the range observed at ti from ihe fifth 
tracking station, and range rate is the time derivative of the range denoted by R,(ti): these 
are given. For the computer simulation presented here, we have taken the nominal values of 
the range and the range rate and added random noise. The standard deviations assumed are 
30 meters in range and about 50 cm per second in range rate. 
For this particular problem, the Q function should be considered. 
weighted squares of the residuals: 
where 
g, ($1 = { [. ($1 - x, (ti) I 
+ 
+ 
[Y (ti) - Y, ($1 I 
1 I2 
[z ($1 - zs ($1 I 1 
and 
x,(ti), y,(ti), z,(t,) = the fifth station coordinate1 a! ti, 
R,( ti) = the range from the particular tracking station, 
g, = the range computed using the present iterate, 
ks( t i )  = observed range rate, and 
W, = weight for the particuld observation. 
This is the sum of the 
3 
1 
18.0 
0.0 
10 
168.0 
I 
! 2.0 10.0 
28.0 14.0 
20 30 
79.0 5 2.0 
I I 
We can see that gs is the function of the state at t i ,  and x, represents the coordinates of the 
station itself. The minimization of this function Q with respect to the integratior. constant 
is carried out using Newton’s method. 
Table 1 provides the data for the three tracking stations we have assumed. The longitude 
and latitude, number of observations, and the tctal interval considered-dhout 1 500 sec- 
onds-are listed. The elliptic orbit conside.-ed has ar- :centricity of 0.0557, with a semi- 
inajor axis of about 7178 km. This corresponds to  a perigee at 400 km altitude, an apagee 
at 1200 km altitude, and an inclination of 20” to the equatorial plane. 
Tab:; 1 
Statior, ant1 Orbit Data 
Parameters 
Latitude (deg) 
Longitude (deg) 
Number of Observations 
Interval Between Observations (s) 
Semimajor Axis (km: 
Eccentricity 
Inclination (deg) 
Longitude of Ascending Node 
Argument of Perigee 
- 
Station 
1 1 2 1 3  
a =  0.0 
o =  0.0 
Some of the results obtained are depicted in figure 1,  whicb shows the error as a function 
of time for different iterations. The initial guess is off by about 70 km from the true solu- 
tion, and with !wo iterations the error is brcught dowli to something like 15 km. 
In order to  sc.e the convergence properly, the log of absolute error is plotted as a function 
of the iterhtion number in figure 2. The error is plotted for SevPral p i n t s  along the trajec- 
tory. The first point is T = 0, and the last point is T = 15 13.1. I t  is seen that the conver- 
gence is linear as was expected because of its relation t o  the cldqsii-al Picard method. 
One of the disadvantages of this methad is that the convergence is obtaied only for arcs 
of length less than a t  -it one-third of a revolution. Another problem is that of aeveloping 
a pocedure for ‘Ig t-:io successive arcs. 
4 
0-  
6 0 -  
a- 
FIRST ITERATION 
0 500 lo00 1500 
INITIAL GUESS 
figure 1. Error as a function of time for different iterations 
(elliptic orbit). 
On the other hand, there are many advantages: For example, the method is quite simple and 
does not require the linear perturbation equations. The solution is in the form of a poly- 
nomial and is convenient t o  store. There is no interpolation required, should the solution 
be needed at some intermediate point. The error in representing this by a polynomial can, 
of course, be estimated by observing the last few terms of the polynomial. Also, the method 
is not sensitive to a poor initial guess. As we have seen here, the initial guess was off by 
about 70 km and even then the prcjcess will converge without difficulty. 
The question is, is this method really more efficient than a more traditional approach to 
solving a problem? Actually, we have not compared it with the existing methods. However, 
we have some estimates on the time taken for this particular problem-it takcs about 30 
seconds on the IBM 360/6S. 
5 
a 
m 
U 
LL 
0 
(1 -1 
-2 
ITERATION NUMBER 
Figure 2. Convergence of the method for an elliptic orbit. 
6 
A FINlTE ELEMENT GRAVITY FIELLl? 
J.  L. Junkins 
University of Vitginia 
Charlottesville, Virginia 
The paper details an approach for constructing a globally valid, piecewise continuous family 
of locally valid potential functions. The thesis put forth here is that the higher frequency 
terms of the geopotential can be more efficiently computed from such locally convergent 
approximations (typically three variable power series of order less than 5 )  than from any 
Y dobally convergent gravity representation. This approach appears to  be a step in the direc- 
tion of voiding the recent trend that the better model we have of the geopotential, the more 
expensive it is t o  integrate orbits with it! Numerical experiments conducted thus far confirm 
the validity of the approach and that acceleration errors of 0 ( km/sz) are achievable. 
The trade-offs involved in selecting the finite element shape and size, and the order of the 
loch'  approximations versus resulting accuracy, computational speed, and storage require- 
ment, are currently under study. 
7 
MEAN RATES OF THE ORBITAL ELEMENTS OF A SATELLITE 
PERTURBED BY A LENS SHAPED MASCON 
M. t .  Ananda 
Jet Propulsion Labomtories 
Pasadena, Culifontia 
A set of mean orbital rates are computed for a satellite perturbed by a lens shaped mascon. 
A disturbing potential in terms of the orbital elements of the satellite and the mascon 
parameters is developed. The partial derivatives of the potential with respect to  the orbital 
elements are formed. These partials are averaged over the pcriod of the satellite orbit to 
eliminate the short periodic terms. The averaged partials are substituted into the variation 
of parameters equations t o  give the mean orbital rates. In the limiting case, when the radius 
of the lens shaped mascon reaches zero the mean orbital rate due to  a point mass is ab- 
tained. The orbita! rates developed by the method described here are compared against the 
rates obtained by numerical differencing. The method developed here is used to  reduce the 
Apollo-15 and -16 subsatellite data for lunar farside gravity determination. 
8 


The determination of orbits using Picard iteration

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server