damn numerical error in Fortran

if we run the following simple code:

F90:

program foo
real*8 :: a = 0
integer :: n
do n=0,100,1
print*,a
a = a + 0.01
end do
end program foo

#include <stdio.h>
void main(void)
{
double a = 0;
int n;
for(n=0; n<=100; n++)
{
printf(” %.10f n”,a);
a += 0.01;
}
}

we may find different outputs from different compilers:

Intel Fortran	g95	Gcc
0.909999979659915 0.919999979436398 0.929999979212880 0.939999978989363 0.949999978765845 0.959999978542328 0.969999978318810 0.979999978095293 0.989999977871776 0.999999977648258	0.909999979659915 0.9199999794363976 0.9299999792128801 0.9399999789893627 0.9499999787658453 0.9599999785423279 0.9699999783188105 0.979999978095293 0.9899999778717756 0.9999999776482582	0.9100000000 0.9200000000 0.9300000000 0.9400000000 0.9500000000 0.9600000000 0.9700000000 0.9800000000 0.9900000000 1.0000000000

if we change that +0.01 to +0.01d0, it will improve the accuracy but still unsatisfactory

Ifort real*16	Ifort real*8	G95 real*8
0.910000000000000018943180357666733 0.920000000000000019151347174783950 0.930000000000000019359513991901167 0.940000000000000019567680809018384 0.950000000000000019775847626135601 0.960000000000000019984014443252818 0.970000000000000020192181260370035 0.980000000000000020400348077487251 0.990000000000000020608514894604468 1.00000000000000002081668171172169	0.910000000000001 0.920000000000001 0.930000000000001 0.940000000000001 0.950000000000001 0.960000000000001 0.970000000000001 0.980000000000001 0.990000000000001 1.00000000000000	0.9100000000000006 0.9200000000000006 0.9300000000000006 0.9400000000000006 0.9500000000000006 0.9600000000000006 0.9700000000000006 0.9800000000000006 0.9900000000000007 1.0000000000000007

so.. god bless fortran…

2 thoughts on “damn numerical error in Fortran”

Published by 犊犊