https://gcc.gnu.org/bugzilla/show_bug.cgi?id=107370
--- Comment #20 from Jakub Jelinek <jakub at gcc dot gnu.org> ---
All I'm arguing is that everything I see from gcc seems to be correct.
30 digits precision is not enough for IEEE quad, 36 is enough.
If I try your sqrt (2^255) case on x86_64, I see:
#define _GNU_SOURCE
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
int
main ()
{
volatile _Float128 x = __builtin_sqrtf128(powf128(2, 225));
char buf[256];
strfromf128 (buf, 128, "%a", x);
strfromf128 (buf + 128, 128, "%.36f", x);
printf ("%s\t%s\n", buf, buf + 128);
x = x * x;
strfromf128 (buf, 128, "%a", x);
strfromf128 (buf + 128, 128, "%.36f", x);
printf ("%s\t%s\n", buf, buf + 128);
}
0x1.6a09e667f3bcc908b2fb1366ea95p+112
7343016637207168931428032607349397.000000000000000000000000000000000000
0x1.ffffffffffffffffffffffffffffp+224
53919893334301279589334030174039256154977430310253516431710891278336.000000000000000000000000000000000000
so again, x * x is 1ulp from 0x1.0p+225.
