function XL0_calc(xt,zn,t,n){ //xt星体,zn坐标号,t儒略世纪数,n计算项数
t/=10; //转为儒略千年数
var i,j,v=0,tn=1,c;
var F=XL0[xt],n1,n2,N;
var n0, pn=zn*6+1, N0 = F[pn+1]-F[pn]; //N0序列总数
for(i=0;i<6;i++,tn*=t){
n1=F[pn+i], n2=F[pn+1+i], n0=n2-n1;
if(!n0) continue;
if(n<0) N=n2; //确定项数
else { N = int2(3*n*n0/N0+0.5)+n1; if(i) N+=3; if(N > n2) N=n2; }
for(j=n1,c=0;j<N;j+=3) c += F[j]*Math.cos(F[j+1] +t*F[j+2]);
v += c*tn;
}
v/=F[0];
if(xt==0){ //地球
var t2=t*t, t3=t2*t; //千年数的各次方
if(zn==0) v += (-0.0728 -2.7702*t -1.1019*t2 -0.0996*t3) / rad;
if(zn==1) v += (+0.0000 +0.0004*t +0.0004*t2 -0.0026*t3) / rad;
if(zn==2) v += (-0.0020 +0.0044*t +0.0213*t2 -0.0250*t3) / 1000000;
}else{ //其它行星
var dv = XL0_xzb[ (xt-1)*3+zn ];
if(zn==0) v += -3*t/rad;
if(zn==2) v += dv/1000000;
else v += dv/rad;
}
return v;
}
我计算 XL0_calc(1,0,1.114532,10)
当计算到 v/=F[0]; 的上一行时 v=2912419233242.7744 经验证，结果 正确
再往下算，由于 F[0]=1000000000,所以 v/=F[0] = 2912.4192332428 也正确
由于 xt=1 zn=0
所以后面的代码中，只需计算 v += -3*t/rad;
而t=t/10=1.114532/10=0.1114532 ，-3*t/rad= -0.0000016210
因此，v += -3*t/rad=2912.4192332428 -0.0000016210 = 2912.4192316218
从程序的运行规律规律来看，这样计算的答案 应该是正确的
但是，通过运行整个程序，得到的最终结果 却有差异，实际结果为：2912.4192312033106
请高人指点迷津
加好友 
发短信
Post By：2011/6/19 21:32:00
