第一次真正碰到了数值精度问题。调试了我整整一个下午加一个晚上加半个上午……

integer, parameter :: NTOT3=1024*1024*1024
real :: mean
real, allocatable :: U(:)
integer I

allocate (U(NTOT3))

!... 读取U的部分

mean = 0.0
Do I = 1, NTOT3
  mean = mean + U(I)
Enddo
Write (*,*) mean/Real(NTOT3)

这样的做法是错误的。在叠加进行到一部分的时候，mean的值已经远远大于U(I)，以至于U(I)比mean的最小分辨还小，因此后面的叠加操作全部都没法执行了。

解决办法：对mean定义成real*8，最好是real*16类型