Question/request/bug(?) re. floating-point in dmd

[qznc]

On Thursday, 24 October 2013 at 06:12:03 UTC, qznc wrote:
>> import std.stdio;
>> struct pair { double hi, lo; }
>> pair normalize(pair q)
>> {
>>  double h = q.hi + q.lo;
>>  double l = q.lo + (q.hi - h);
>>  return pair(h,l);
>> }
>> void main()
>> {
>>  immutable static pair spn = normalize(pair(1.0,0.1));
>>  writeln(spn);
>>  writeln(pair(1.0,0.1));
>>  writeln(normalize(pair(1.0,0.1)));
>> }
> I cannot see any significant difference. The fadd-fsub-fadd 
> sequence seems to be the same in both cases.

More observations. It requires -m32 to reproduce. On 64bit is 
gives the desired output even in optimized form.

I wrote a small C program for comparison:

------
#include <stdio.h>

typedef struct pair {
   double hi, lo;
} pair;

pair normalize(pair q) {
   double h = q.hi + q.lo;
   double l = q.lo + (q.hi - h);
   //double l = q.lo + (q.hi - (q.hi + q.lo));
   pair res = {h,l};
   return res;
}

int main() {
   pair x = {1.0, 0.1};
   pair y = normalize(x);
   printf("%g %g\n", y.hi, y.lo);
   return 0;
}
------
gcc -m32 normtest_replicate.c; ./a.out
Output: 1.1 -8.32667e-17

Note the commented line, if you use that instead of the line 
above, then
Output: 1.1 0

The C semantic says that h must be rounded to double. In the 
second case l is computed with full hardware precision instead.

For a shorter example, two versions:

    double x = a - b;     double y = c + (a - b);
    double y = c + x;

In C those have different semantics in terms of the precision, 
but in D they are the equivalent.