Документ взят из кэша поисковой машины. Адрес оригинального документа : http://www.adass.org/adass/proceedings/adass96/reprints/taidakovat.pdf Дата изменения: Thu Jan 15 01:05:57 1998 Дата индексирования: Tue Oct 2 16:00:52 2012 Кодировка: Поисковые слова: п п п п р п р п р п р п р п р п р п р п р п

Astronomical Data Analysis Software and Systems VI ASP Conference Series, Vol. 125, 1997 Gareth Hunt and H. E. Payne, eds.

A New Stable Method for Long-Time Integration in an N-Body Problem
Tanya Taidakova Crimean Astrophysical Observatory, Simeiz, 334242, Ukraine Abstract. The most serious error in numerical simulations is the accumulation of discretization error due to the finite stepsize. Traditional integrators such as Runge-Kutta methods cause linear secular errors to the energy, the semi-ma jor axis, and the eccentricity of orbiting ob jects. Potter (1973) describ ed an implicit second-order integrator for particles in a plasma with a magnetic field. We have used this integrator for an investigation of the dynamics of particles around a planet (or star) in a co-rotating coordinate system. A big advantage of this numerical integrator is its stability: the error in the semi-ma jor axis and the eccentricity dep ends only on the step size and does not grow with an increasing numb er of time steps. The argument of p ericenter changes linearly with time and more slowly than in the case of the Runge-Kutta integrator. In addition, this implicit integrator takes much less computing time than the second-order Runge-Kutta method. We tested this method for several astronomical systems and for motion of an asteroid in a 1:1 Jupiter resonance during 200 million time steps (ab out 5 million years or 800 thousand p eriods of asteroid resonance motion).

1.

Introduction

In the last few years there has b een great interest in the numerical study of long term evolution of b odies of the Solar system. As the integration time increases, the numerical results b ecome more contaminated by various errors. The most serious error is the accumulation of the discretization (truncation) error due to a finite stepsize (or the replacement of continuous differential equations by finite difference equations). The conventional integrators such as Runge-Kutta, multistep and Taylor methods, generate linear secular errors in orbital energy and angular momentum. This means that the semi-ma jor axis and the eccentricity change linearly with time and the linear secular error in the semi-ma jor axis produces a quadratic secular error in the planetary longitude. A new symplectic integrator produces no secular truncation errors in the actions of a Hamiltonian system. 2. Implicit Integrator

In this pap er we briefly discuss an implicit numerical integrator. The discretization errors in the energy, the semi-ma jor axis, and the eccentricity by the implicit second-order integrator show only p eriodic changes. The truncation error in the 174

© Copyright 1997 Astronomical Society of the Pacific. All rights reserved.


Stable Long-Time Integration in an N-Body Problem

175

argument of p ericenter grows linearly in time. The equations of motion of a particle in the gravitational field of the Sun and the planet with mass mpl in the corotating coordinate system take the form: x= Ё y= Ё z= Ё where: - - - x (x2 + y 2 + z y (x + y + z 2 ) z (x2 + y 2 + z 2 )
2 2 3/2 2 )3/2

2y + x + Fx -2 + y + Fy x Fz ,

(1)

Fx Fy Fz

= = =

- - -

mpl (x - xpl ) [(x - xpl )2 +(y - ypl )2 +(z - zpl )2 )] mpl (y - ypl ) [(x - xpl )2 +(y - ypl )2 +(z - zpl )2 )] mpl z [(x - xpl )2 +(y - ypl )2 +(z - zpl )2 )]
3/2

3/2

3/2

3/2

.

Here the total mass of the Sun and the planet is taken as the unit of mass, and the distance b etween the planet and the Sun is taken as the unit of length. The unit of time is chosen in such a way that the angular velocity of orbital motion of the planet is equal to unity, and, hence, its orbital p eriod is 2 . Let v = vx + ivy , x = x + iy , and F = Fx + iFy . We obtain rather than (1): dz dt dvz dt dx dt dv dt = = = vz Fz v


(2)

= -2 iv + x + F

.

We may solve the equation dU (t) + R(U (t),t) = 0 with initial conditions dt U (t0 ) = U0 by the implicit second-order integrator describ ed in Potter (1973): v
[n+1]

=

v

[ n]

-

1 (R 2

[ n]

+R

[n+1]

)t

.

(3)

In our equations (2), R is a function of x, y , z . We will calculate this function R in space-time p oints n + 1/2. From (2) by the use of Eq. (3) we derive the equations for new integrator (Taidakova 1990, 1995): vx (1 - 2 t)+ (2vy + x
[ n] [ n] [ n+ 1 ] 2

[ vxn+1]

=

+ Fx 2 )t +(y 1+2 t

[ n+ 1 ]

[ n+ 1 ] 2

+ Fy

[ n+ 1 ] 2

)2 t


176

Taidakova

Figure 1. Numerical errors in the energy with second order RungeKutta, fourth order Runge-Kutta and second order Potter integrator. vy (1 - 2 t) - (2vx - y
[ vzn] + Fz [ n] [ n+ 1 ] 2 [ n] [ n] [ n+ 1 ] 2

[ vyn+1] [ vzn +1]

= = = = =

- Fy 2 )t +(x 1+2 t

[ n+ 1 ]

[ n+ 1 ] 2

+ Fx

[ n+ 1 ] 2

)2 t

t
[ + vxn] )t/2 [ + vyn] )t/2 [ + vzn] )t/2 [ n]

x y z

[n+1] [n+1] [n+1]

[ x +(vxn+1] [ y [n] +(vyn+1]

z

[ n]

[ +(vzn

+1] [ n]

,
[ n+ 1 ] 2

(4) =y
[ n]

where: x[ z 3.
[ n+ 1 ] 2

n+ 1 ] 2 [ n]

=x

+ vx t/2 ; y
[ n+ ] y, z
1 2

+ vy t/2 ;
[ n+ 1 ] 2

[ n]

=z

+ vz t/2 ; Fx,

[ n]

=F x

[ n+ 1 ] 2

,y

,z

[ n+ 1 ] 2

,t

[ n+ 1 ] 2

.

Numerical Examples

We tested this method for several astronomical systems (Taidakova 1995). In order to see the prop erties of the implicit integrator, we first choose the 2-b ody problem. Figure 1 shows the numerical errors in the energy log(E/E ) in the numerical integrations of the motion of the particle in circular orbit around the Sun with different integrators. Figure 2 shows the errors in the parameter i = (x0 - xi )2 +(vx0 - vxi )2 +(vy0 - vyi )2 , where i is the numb er of revolutions, x is the coordinate (y = 0), and vx ,vy are the velocities in the numerical integration of the orbital motion of asteroid in the 1:1 Jupiter resonance during 200 million steps. The computer time with the implicit second-order integrator is ab out 1.52­1.94 times faster than that with second-order Runge-Kutta integrator. The calculation were carried out in a PC with DX4/100.


Stable Long-Time Integration in an N-Body Problem

177

Figure 2. Numerical errors in the parameter i with second order Runge-Kutta method and second order Potter integrator. Acknowledgments. The author thanks the Local Organizing Committee and U.S. Civilian Research & Development Foundation (Grant # 96023) for supp ort. References Potter, D. 1973, Computational Physics (New York: Wiley) Taidakova, T. 1990, Nauch.Inform.Astrosoveta Akademii Nauk SSSR (Riga: Zinatne), 68, 72 Taidakova, T. 1995, Ph.D. Thesis, Moscow State University