求卫星轨道的周长

数值分析实验报告

题目

一、问题提出

地球卫星轨道是一个椭圆，椭圆周长的计算公式是

，这里a是椭圆的半长轴，c是地球中心（椭圆中心）的距离，记h为近地点距离，H为远地点距离，R= 6371（km）为地球半径，则a=（2R+H+h）/2,c=(H-h)/2.我国第一颗人找地球卫星近地点距离h=439（km），远地点距离H=2384（km），试求卫星轨道的周长.

二、模型建立

龙贝格求积算法公式为：

,2,1 , )(141)2(144 )

(1)1(1)( =---=-+-k h T h T T

k m m

k m m m k m

椭圆周长的计算公式：

R= 6371（km ），则a=（2R+H+h ）/2,c=(H-h)/2. R= 6371（km ）， h=439（km ），H=2384（km ）

三、 求解方法

Matlab M 文件：

function R = romberg(f,a,b,n) format long

R = zeros([n + 1, n + 1]);

R(0+1, 0+1) = (b - a) / 2 * (feval(f, a) + feval(f, b)); for i = 1 : n, h = (b - a) / 2^i; s = 0; for k = 1 : 2^(i-1),

s = s + feval(f, a + (2*k - 1)*h); end

R(i+1, 0+1) = R(i-1+1, 0+1)/2 + h*s;

end

for j = 1 : n, fac = 1 / (4^j - 1);

for m = j : n,

R(m+1, j+1) = R(m+1, j-1+1) + fac*(R(m+1, j-1+1) - R(m-1+1, j-1+1));

end

end

function I=f(x)

R=6371;h=439;H=2384;

a=(2*R+H+h)/2;c=(H-h)/2;

I=sqrt(1-(c/a)^2*(sin(x)^2));

四、输出结果

积分I输出结果：

ans =

0 0

即加速3次求得：

k

1

2

计算得：I =

所以卫星轨道的周长S = 4aI = 48708 km

五、结果分析

由计算结果可知，利用龙贝格算法计算积分，利用外推法，提高了计算精度，加快了收敛速度，求得的结果比较精确。