Twiddling with 64-bit values as 2 ints;

["Stefan Kanthak" <stefan.kanthak@nexgo.de>]

Hi,

32 years ago, C89 introduced 64-bit integers: [un]signed long long
IEEE 754 defined the 64-bit double-precision floating-point format,
now called binary64. in 1985.

Especially SunSoft's [fd]libm, which (to my knowledge) started around
this time, and also IBM's APMathLib/libultim, which followed a little
later, and also quite some ACM TOMS routines, but use (pairs of) 32-bit
integers for bit-twiddling on the representation of double/binary64:
additions/subtractions/shifts on the 52-bit mantissa/fraction, and
operations on the full 64-bit double, involve both ints, and need to
take care of the carry/borrow -- explicitly, and quite ugly!
It's also generally unknown whether a compiler will recognize this
sort of carry/borrow/overflow handling and generate proper machine
code using "add with carry"/"subtract with borrow" instructions.

JFTR: while sticking with 32-bit integers MAY give better performance
      on 32-bit processors, especially when an operations only involves
      either low or high part, the explicit carry/borrow handling can
      have negative performance impact.

See for example <http://www.netlib.no/netlib/toms/722>, written by
William J. Cody (known from Cody/Waite range reduction):

|    W. J. Cody, J. T. Coonen, March 30, 1992
...
|       /* Otherwise, use integer arithmetic to increment or      */
|       /* decrement least significant half of z, being careful   */
|       /* with carries and borrows involving most significant    */
|       /* half.                                                  */
|          else if (((argx < Zero) && (argx < argy)) ||
|                   ((argx > Zero) && (argx > argy))) {
|                   --lowpart(z);
|                   if (lowpart(z) == -1)
|                      --highpart(z);
|                   }
|                else {
|                   ++lowpart(z);
|                   if (lowpart(z) == 0)
|                      ++highpart(z);
|                   }
|

Compare this with the REALLY UGLY
<https://sourceware.org/git/?p=glibc.git;a=blob_plain;f=math/s_nextafter.c;hb=HEAD>

|  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
...
|        if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) ||   /* x is nan */
|           ((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0))     /* y is nan */
|           return x+y;
...
|        if(hx>=0) {                               /* x > 0 */
|            if(hx>hy||((hx==hy)&&(lx>ly))) {      /* x > y, x -= ulp */
|                if(lx==0) hx -= 1;
|                lx -= 1;
|            } else {                              /* x < y, x += ulp */
|                lx += 1;
|                if(lx==0) hx += 1;
|            }
|        } else {                                  /* x < 0 */
|            if(hy>=0||hx>hy||((hx==hy)&&(lx>ly))){/* x < y, x -= ulp */
|                if(lx==0) hx -= 1;
|                lx -= 1;
|            } else {                              /* x > y, x += ulp */
|                lx += 1;
|                if(lx==0) hx += 1;
|            }
|        }

(Heretic.-) questions:
- why does glibc still employ such ugly code?
- Why doesn't glibc take advantage of 64-bit integers in such code?

JFTR: on 64-bit processors, when the compiler does not recognize
      that hx:lx and hy:ly are in fact a single 64-bit integer it
      can hold in a SINGLE register, but smears it over 2 registers,
      such cruft kills performance.

For 32-bit processors, the JFTR from above still holds: using 64-bit
integers with a C89 compiler should give better machine code.

Stefan