【IDL】讨论浮点运算精度问题（IDL与Excel和Matlab相比）

PRO aaaaaaaa
a=[19.426952,9.1238899, 21.078484, 24.719807  ,27.361130 , 22.608578 ,     24.159088 , 21.957047 , 19.305515$
      ,19.103474 , 17.451942  ,18.552963  ,21.653984 , 18.305515,  21.406536,  21.159088,  18.810619,  15.911640$
     , 17.563171 ,17.866234 ,15.618786 ,15.022869 ,20.876442, 20.977463 ,23.078484 ,16.527973 ,13.876442$
,14.527973 ,11.977463, 7.3259315, 5.9774628, 4.8764420, 1.9774628 ,6.3815460 ,6.0330772 ,7.2351189$
   ,1.1897125,8.3917542,6.6948166 ,9.5483894, 3.6494102 ,4.5029831,-2.2041626, -238.76700 ,-117.44246$
    ,-18.193954 , -6.1633301 , -11.456184  ,-11.491383  ,-19.875050]
  b=[  437.0,443.000,434.000  ,451.000  , 468.000  ,473.000  ,472.000 , 476.000,  479.000  , 483.000  ,$
486.000,484.000 $
     ,482.000  ,  479.000  ,  477.000  ,   472.000 ,    469.000 ,    467.000 $
     ,464.000  ,  458.000  ,  453.000  ,   445.000,     438.000 ,    436.000 $
     ,434.000 ,    435.000  ,   438.000 ,    435.000 ,    436.000,     439.000 $
     ,436.000 ,    438.000  ,   436.000 ,    428.000 ,    425.000 ,    421.000 $
    , 412.000 ,    408.000  ,   402.000 ,    395.000 ,    393.000 ,    386.000 $
    , 400.000 ,    817.000 ,    573.000 ,    380.000,     320.000 ,    306.000 $
    , 277.000 ,    245.000]     
     w=REFORM(b,50,1)##REFORM(a,1,50)
END

PS：该有续行符的地方都是有的~编辑好了发表出来确实乱的~矩阵赋值没有问题~
不知道怎么编辑，显得有点乱，其实就是一个50*1的矩阵乘以一个1*50的矩阵，结果保存在w中。
IDL的计算结果是：-0.13714981，但是MATLAB中的结果是：-3.492459654808044e-010。对计算机语言而言，浮点数计算肯定丢失精度的。首先：在excel下进行了验证，结果如下：前面从1*1开始一直算到50*50，二者对比如下。

IDL下转换为精度输出的程序段与结果

openw,lun,'c:\temp\ab.txt',/get
  for i=0,49 do begin
    nA = a[0:i]*1d
    nB = b[0:i]*1d
   printf,lun,i+1,'个数:',(REFORM(na,i+1,1))##(REFORM(nb,1,i+1))
  endfor

      1个数:      8489.5782
      2个数:      12531.461
      3个数:      21679.523
      4个数:      32828.156
      5个数:      45633.165
      6个数:      56327.022
      7个数:      67730.112
      8个数:      78181.666
      9个数:      87429.008
     10个数:      96655.985
     11个数:      105137.63
     12个数:      114117.26
     13个数:      124554.48
     14个数:      133322.83
     15个数:      143533.74
     16个数:      153520.83
     17个数:      162343.01
     18个数:      169773.75
     19个数:      177923.06
     20个数:      186105.80
     21个数:      193181.11
     22个数:      199866.28
     23个数:      209010.16
     24个数:      218156.34
     25个数:      228172.40
     26个数:      235362.07
     27个数:      241439.95
     28个数:      247759.62
     29个数:      252981.79
     30个数:      256197.88
     31个数:      258804.05
     32个数:      260939.93
     33个数:      261802.11
     34个数:      264533.41
     35个数:      267097.46
     36个数:      270143.45
     37个数:      270633.61
     38个数:      274057.45
     39个数:      276748.76
     40个数:      280520.38
     41个数:      281954.60
     42个数:      283692.75
     43个数:      282811.08
     44个数:      87738.444
     45个数:      20443.915
     46个数:      13530.212
     47个数:      11557.947
     48个数:      8052.3544
     49个数:      4869.2414
     50个数:    -0.14574051

后来按照下面流程处理了：先转换为64位长整型（*10^7），全部按照长整型计算，最后结果除以10^7，这么计算精度应该是最高了吧。

程序代码如下：

OPENW,lun,'c:\temp\ab-long.txt',/get
  FOR i=0,49 DO BEGIN
    nA = LONG64(a[0:i]*10LL^7)
    nB = LONG64(b[0:i])
   PRINTF,lun,i+1,'个数:',((REFORM(na,i+1,1))##(REFORM(nb,1,i+1)))/(10D^7)
  ENDFOR
  FREE_LUN,lun

结果如下：

      1个数:      8489.5780
      2个数:      12531.461
      3个数:      21679.523
      4个数:      32828.156
      5个数:      45633.164
      6个数:      56327.021
      7个数:      67730.111
      8个数:      78181.665
      9个数:      87429.007
     10个数:      96655.985
     11个数:      105137.63
     12个数:      114117.26
     13个数:      124554.48
     14个数:      133322.82
     15个数:      143533.74
     16个数:      153520.83
     17个数:      162343.01
     18个数:      169773.75
     19个数:      177923.06
     20个数:      186105.79
     21个数:      193181.10
     22个数:      199866.28
     23个数:      209010.16
     24个数:      218156.34
     25个数:      228172.40
     26个数:      235362.07
     27个数:      241439.95
     28个数:      247759.62
     29个数:      252981.79
     30个数:      256197.87
     31个数:      258804.05
     32个数:      260939.93
     33个数:      261802.10
     34个数:      264533.40
     35个数:      267097.46
     36个数:      270143.45
     37个数:      270633.61
     38个数:      274057.44
     39个数:      276748.76
     40个数:      280520.37
     41个数:      281954.59
     42个数:      283692.74
     43个数:      282811.08
     44个数:      87738.439
     45个数:      20443.911
     46个数:      13530.208
     47个数:      11557.942
     48个数:      8052.3497
     49个数:      4869.2368
     50个数:    -0.15034480