Hi Fabian,
This was a bit bizarre particularly as when I tried it on Windows there
didn't seem to be a problem. It is a problem with the precision and
when rounding occurs. Poly/ML uses the old x87 floating point unit on
both X86/32 and X86/64. As far as I'm aware all x86/64 machines support
the SSE2 instructions so other compilers use SSE2 instructions on 64-bit
but, by default use the x87 instructions on 32-bit.
The x87 instructions use 80-bit registers internally and only round to
64-bits when the registers are stored to double-precision values. The
SSE2 instructions are 64-bit so all the results are the same whether
values are stored to memory or not. The x87 has bits in the control
word that can allow it to compute intermediate results to 80-bit, 64-bit
or 32-bit precision. The question is whether this should be set to
80-bit or 64-bits?
Currently Poly/ML sets it to 80-bit precision at the start. As far as
I've been able to ascertain that is what other x86/32 compilers do. The
idea is, I suppose, to compute intermediate results to the best
precision available. It does, though, produce the strange results that
you've encountered. I think the reason for the difference in results is
whether the function can be computed at compile-time or not.
I'm coming round to the view that it should instead set the FPU to
64-bits so that intermediate results are the same whether they remain in
the register or are written out to memory and reloaded. It looks as
though the libraries in Visual C set that anyway which is why I only saw
the problem when I tested on Linux. This would also mean that results
would not change when, as I plan, I change to using the SSE2
instructions on x86/64.
Best regards,
David
On 07/10/2015 14:07, Fabian Immler wrote:
Hi David,
Perhaps I should rephrase my enquiry to: "Do I get any guarantees from the
implementation of type real in PolyML"?
Naively, I would have expected, that every subexpression is evaluated according
to the IEEE-754 standard (in double precision).
I skimmed the implementation in libpolyml/reals.cpp, it looks like all intermediate
results are stored in and retrieved from memory (real_arg, real_result). Is that
sufficient to comply to IEEE-754? I.e., can there not occur "double rounding"
(as in [1, section 3.3.1])? Or isn't there the need to use something like
FLT_EVAL_METHOD=1 [2,3]?
[1] http://www.springer.com/us/book/9780817647049
[2] https://en.wikipedia.org/wiki/C99#IEEE.C2.A0754_floating_point_support
[3] http://www.cplusplus.com/reference/cfloat/
Best regards,
Fabian
Am 07.10.2015 um 12:26 schrieb Fabian Immler <[email protected]>:
Hi David,
I noticed some strange behavior when working with type real in PolyML:
In the following example, where I have "testres = twosum a", depending on
whether I call
"testres" or "twosum a" I get different results.
If I take the definitions out of the structure Test, the results are the same.
They are also the same if e.g. print-statements are inserted in between the
assignments to b, s, sb and r. The results are also the same when evaluating
inside Isabelle with the debugger turned on.
Do you have an idea what is going on there?
Is it possible that some intermediate results are kept in higher precision in
the FPU?
Does the compiler perform some "optimizations" on expressions of type real?
Best regards,
Fabian
structure Test = struct
fun one n = (if n = 0 then 100.0 else one (n - 1));
val a : real = one 1;
fun twosum a =
let
val b = 1.0/1243313.0;
val s = a + b
val sb = s - b
val r = s - sb
in
r
end;
val testres = twosum a;
end;
open Test;
val a = twosum a
val b = testres
_______________________________________________
polyml mailing list
[email protected]
http://lists.inf.ed.ac.uk/mailman/listinfo/polyml