Hi David,

Thanks for your answer, I should have subscribed to the polyml mailing list 
earlier.
Using the SSE2 instructions and/or setting the intermediate precision to 64 bit 
would give a nice and clear semantics to floating point expressions in Poly/ML, 
so I would be very much looking forward to such a change.

Best regards,
Fabian


> Am 07.10.2015 um 15:56 schrieb David Matthews 
> <[email protected]>:
> 
> 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 <immler at in.tum.de>:
>>> 
>>> 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

Reply via email to