数值分析实验报告
题目
一、 问题提出
地球卫星轨道是一个椭圆，椭圆周长的计算公式是

，这里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)(khThTT
kmmkmmmk
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
0

即加速3次求得：
k
0
1
2
计算得：I =
所以卫星轨道的周长S = 4aI = 48708 km

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