forums.silverfrost.com

[ljfarrugia]

The trivial code below produces what I consider inconsistent results. Specifically there is a difference between results when the mod intrinsic function is presented with a plain argument (c) or with an argument which requires evaluation (a+b). The code works exactly as expected when compiled with gfortran on Linux, but gives confusing results with Silverfrost Fortran. I am using version 8.82.0

program main
double precision :: a,b,c,d

logical :: f1, f2
a = 1.0d0/3.0d0
b = 2.0d0/3.0d0
write(*,'(a,e15.10e1)') 'a = ',a
write(*,'(a,e15.10e1)') 'b = ',b
c = a + b
write(*,'(a,e15.10e1)') 'a + b = ',a+b
write(*,'(a,e15.10e1)') 'c = ',c
f1 = (c > 1.0d0)
f2 = (c < 1.0d0)
write(*,*) 'c is smaller than 1.0 = ',f2
write(*,*) 'c is larger than 1.0 = ',f1
f1 = (a+b > 1.0d0)
f2 = (a+b < 1.0d0)
write(*,*) 'a+b is smaller than 1.0 = ',f2
write(*,*) 'a+b is larger than 1.0 = ',f1
c = abs(mod(c,1.0d0))
d = abs(mod(a+b,1.0d0))
write(*,'(a,e15.10e1)') 'abs(mod(c,1.0)) = ',c
write(*,'(a,e15.10e1)') 'abs(mod(a+b,1.0)) = ',d
end

[JohnCampbell]

I changed your program to:

[Kenneth_Smith]

Running John's program with FTN95 v9.03.0.

The 32 bit complier returns:

[JohnCampbell]

OK,

A new program that demonstrates 32-bit has different round-off.

[PaulLaidler]

For maximum precision you can use real*10 with Win32 rather than x64.

[ljfarrugia]

Many thanks for all the comments. I kind of guessed that rounding errors were the cause of the discrepancy, but I didn't realise that 32bit executables were so subject to this issue.

I came across this issue when trying to understand why I was getting seriously differing results from a rather large program compiled under Windows & Linux. I normally have used 32 bit executables and in the past I havent noticed any great issues between the platforms

Any suggestions as the the best way forward? Just use 64 bit executables?
Louis

[PaulLaidler]

64 bit executables are not generally the answer to issues relating to floating point precision.

For example 32 bit FTN95 provides for 80 bit floating point values as well as 32 and 64. Most Fortran compilers only provide 32 bit and 64 bit floating point values.

It is the programmer's responsibility to be aware/beware of calculations that are prone to round-off error.

[Kenneth_Smith]

This is a contrived example, but it illustrates Paul's point about 64 bit executables not being the "fix all" answer to floating point issues.

[DanRRight]

I am curious why 32bit FTN95 still retain 80bit accuracy while 64bit one has only 64bit one? So both Intel and AMD processors have hardware for 80bit while 64bit FTN95 ignores it (as well as other Fortran compilers?)

Also how about SSE2 and AVX instructions with 64bit compilers. Do they still support 80bit accuracy? I've heard this opinion on the net. Also have heard that it's Windows which does not support double precision 80bit accuracy, Linux supports it (is this true?)

If anyone has some free time to investigate all that mess please try to make short confirming tests

[PaulLaidler]

64 bit applications can in theory access the the x87 FPU coprocessor which has 80 bit floating point precision. This coprocessor is separate from the the CPU that uses SSE etc. instructions and which basically does not extend to 80 bits.

The problem is that, at least under Windows, a massive time penalty is incurred when mixing x87 instructions with SSE etc. instructions. So compilers normally avoid this mix.

[JohnCampbell]