compiler option -pc80 for better precision

[breil1]

Hello

I had strange precision problem with ifort and I recently discover -pc80
I just want to have some advice about other option people use to have real double precision result. I found -prec-div for division and -prec-sqrt for square root. Do some people recommand thoose options or other ?

At least I don't understand why -pc80 is not by default.

[TimP]

-pc80 increases only the precision of intermediate results, only when using the 32-bit compiler to generate code without an SSE option. There is a potential reduction in performance, as well as confusion when results don't match results obtained with SSE2 options or on Windows platforms. If your accuracy problem is due to not setting -assume protect_parens, -pc80 is probably not the best solution.
-prec-div and -prec-sqrt require the compiler to generate IEEE standard instructions for those operations, rather than performing optimizations such as invert and multiply, or successive approximation. Those optimizations are not so likely to gain performance, when using the latest CPU families (Intel Penryn or AMD Opteron). The successive approximation is not possible with code where -pc80 would take effect.

[Ron_Green]

And adding to what Tim said, if you have a fairly recent compiler, you will want to explore the -fp-model compiler option. Please read the documentation on this option as it's too lengthy to discuss here. Roughly, it sets up the floating point environment and allows varying degrees of IEEE compliance, from IEEE strict to a fully optimized env with loose compliance but maximum performance.

It's often a balancing act between performance and accuracy. The key is to find the sweet spot for your particular code, and this does vary by application.

[TimP]

Ron is correct to point out that you should consider whether -fp-model will solve your issues. Fortran doesn't support as many options for -fp-model as C, and there is no fp-model option for ifort which deals with extra precision promotion of intermediates. The Intel default is -fp-model fast (unlike other brands of compilers, which default to -fp-model precise or similar equivalent). -fp-model precise for ifort includes the options -no-ftz -prec-sqrt -prec-sqrt -assume protect_parens, as well as removing vector optimizations on sum reductions and many math intrinsics. Only the protect_parens option would apply to the default x87 code generation, where -pc80 has an effect. As far as I can see, -fp-model is more attuned to operation with at least SSE2 code generation (-xW).

[breil1]

Thank you to tim and ron for your replies but I tried the option -assume protect_parens and -fp-model but it didn't work on my problem.
I'm dealing with geometric problem and I need to compute some volume
using green formula. I give you a test of my problem. It's just a test of what I need to do. At the end I need all the error less than at least 1.e-14. For compiler I tried intel fortran ifort 9.0 and 10.1 I also tried G95. g95 give good result with at least -O1. ifort only work with -pc80.

program volume

implicit none

integer,parameter :: prec = selected_real_kind(p=15,r=99)

integer :: i
real(prec) :: vol,vol_plan,vol_poly,vol1,vol2,x1,x2,y1,y2,I00
real(prec) :: xp(4),yp(4),xt(3),yt(3),vt(6),v1(3),v2(3)

print*,'precision',prec

vol=0._prec

vol_plan=0._prec


xp(1)=0.881981658900060750_prec
yp(1)=0.560270236940128785_prec
xp(2)=0.898937686133987057_prec
yp(2)=0.561694546915867732_prec
xp(3)=0.899037081398744609_prec
yp(3)=0.581752936897710904_prec
xp(4)=0.880176007237810309_prec
yp(4)=0.578007759690138179_prec

do i=1,4
x1=xp(i)
y1=yp(i)
write(*, '(A,I,2(E28.20,1X))') 'COORDONNEES ',i,x1,y1
enddo

vol_poly=0._prec

do i=1,4
x1=xp(i)
y1=yp(i)
x2=xp(mod(i,4)+1)
y2=yp(mod(i,4)+1)
I00 = 0.5_prec*(x1+x2)*(y2-y1)
vol_poly = vol_poly + I00
write(*, '(A,I,2(E28.20,1X))') 'CONTRIBUTION VOLUME PLAN',i,I00,vol_poly
enddo

write(*,'(A,1(E28.20,1X))') 'vol_poly',vol_poly

vol=0._prec
vol_plan=0._prec

xt(1)=xp(1)
yt(1)=yp(1)
xt(2)=xp(2)
yt(2)=yp(2)
xt(3)=xp(3)
yt(3)=yp(3)

vol1=0._prec
do i=1,3
x1=xt(i)
y1=yt(i)
x2=xt(mod(i,3)+1)
y2=yt(mod(i,3)+1)
I00 = 0.5_prec*(x1+x2)*(y2-y1)
v1(i)= I00
vt(i)= I00
vol1 = vol1 + I00
write(*, '(A,I,2(E28.20,1X))') 'CONTRIBUTION VOLUME PLAN',i,I00,vol1
enddo

write(*,'(A,1(E28.20,1X))') 'vol1 ',vol1,sum(v1(:))

xt(1)=xp(1)
yt(1)=yp(1)
xt(2)=xp(3)
yt(2)=yp(3)
xt(3)=xp(4)
yt(3)=yp(4)

vol2=0._prec

vol2=0._prec
do i=1,3
x1=xt(i)
y1=yt(i)
x2=xt(mod(i,3)+1)
y2=yt(mod(i,3)+1)
I00 = 0.5_prec*(x1+x2)*(y2-y1)
vt(i+3)= I00
v2(i)=I00
&nbs p; vol2 = vol2 + I00
write(*, '(A,I,2(E28.20,1X))') 'CONTRIBUTION VOLUME PLAN',i,I00,vol2
enddo

write(*,'(A,1(E28.20,1X))') 'vol2 ',vol2,sum(v2(:))


vol=sum(vt(:))
write(*,'(A,1(E28.20,1X))') 'vol',vol
write(*,'(A,1(E28.20,1X))') 'error 1',abs(vol_poly-vol)/vol


vol=vol1+vol2
write(*,'(A,1(E28.20,1X))') 'vol',vol
write(*,'(A,1(E28.20,1X))') 'error 2',abs(vol_poly-vol)/vol

vol=sum(v1(:)+v2(:))
write(*,'(A,1(E28.20,1X))') 'vol',vol
write(*,'(A,1(E28.20,1X))') 'error 3',abs(vol_poly-vol)/vol



end program volume

[TimP]

The (y2-y1) calculations produce loss of precision sufficient to prevent reaching your goal, unless they are done in x87 extended precision. In your example, you could set those differences as your initial constants, and achieve more precision without extended precision than you get now with extended precision.
Microsoft decided to eliminate extended precision as a Windows "standard" back around Visual Studio 5, as they considered it a source of confusion and didn't make the run-time library consistent. You may lament the fading away of extra precision, but the deed was set in stone years ago. Since the design of Pentium 4, you give away a great deal of performance to enable extra precision.
g95 presumably has x87 code with 64-bit rounding mode as a default only in the 32-bit compiler. It may support a real(10) data type, as gfortran does. If you intend to write code which works only with implicit promotion to real(10), you may have better results by using such a compiler with explicit data typing.

[Steven_L_Intel1]

You can also use the REAL(16) type for critical calculations.

[TimP]

In line with Steve's suggestion, there are portable ways to get the next bigger selected_real_kind after double precision, so as to take care of such situations.