OpenMP and floating-point operations

[ftinetti]

Hi Heinz,

I checked the code and I didnt find a bug or an not initialized variable or anything else. No reductions are present. The problems seems to be caused really by the kind floatingpoint-opreations are done.

Do you have a (reduced yet complete) example to play around with (compile, run)? It's good you examined the code, but I think it's unlikely the compiler changes the kind of floating point operations... Maybe something along the lines of what Mark suggests is happening, but in any case, the problem seems to be ill-conditioned (or something similar to), in which case parallel computing exposed it, hopefully in the early stages...

HTH,

Fernando.

[HeinzM]

Hi Mark,

I didnt know that 80bit precision is used internally. May be that there is the problem. I do cfd-computations, and using Kind=8 variables has well worked for many many years. This is the first time such a problem arises. I also asked the Intel-premium-support. They think about it, since over a week

Thanks a lot!
Best regards,
Heinz

[HeinzM]

Hi Fernando

Do you have a (reduced yet complete) example to play around with (compile, run)?

It is a finite-element-code for computing water flow over floodplains. The problem occours while assembling the element-matrices. A few hundred show the problem, but there are hundred-thousand. So reducing the code will likely not show the problem. I dont believe you want to walk through the hole code, right?

It's good you examined the code, but I think it's unlikely the compiler changes the kind of floating point operations... Maybe something along the lines of what Mark suggests is happening, but in any case, the problem seems to be ill-conditioned (or something similar to), in which case parallel computing exposed it, hopefully in the early stages...

If you are right, I have a hidden problem and no idea how to find it, not even where to look for .....

Thank you and best regards,
Heinz

[ftinetti]

Hi again,

A few hundred show the problem, but there are hundred-thousand. So reducing the code will likely not show the problem. I dont believe you want to walk through the hole code, right?

Hmmm... I don't know, feel free to send the code, I would like to see it, anyway. I do not guarantee any result, of course...

ftinetti wrote:
It's good you examined the code, but I think it's unlikely the compiler changes the kind of floating point operations... Maybe something along the lines of what Mark suggests is happening, but in any case, the problem seems to be ill-conditioned (or something similar to), in which case parallel computing exposed it, hopefully in the early stages...

If you are right

It's a big "if" ...

I have a hidden problem and no idea how to find it, not even where to look for .....

You are right, and it is a bug-like behavior.

Fernando.

[Oldboy]

The errors Heinz has found are not small indeed. Here are examples of smaller errors:

http://people.sc.fsu.edu/~jburkardt/f_s ... output.txt
http://people.sc.fsu.edu/~jburkardt/f_s ... output.txt
http://people.sc.fsu.edu/~jburkardt/f_s ... output.txt
http://en.wikipedia.org/wiki/IEEE-754
Observe that gfortran and intelfortran don't show the same errors using omp.
The serial version is only gfortran.
In iEEE-754 5 rounding rules are specified. I have not seen any compiler directives how to use one of them.