On Apr 12, 2024, at 4:33 PM, Rich Felker <dalias@libc.org> wrote:

On Fri, Apr 12, 2024 at 03:57:28PM -0400, Rich Felker wrote:

On Thu, Apr 11, 2024 at 06:17:25PM -0700, Peter Ammon wrote:

if (p>=0 && p<(LDBL_MANT_DIG-1+3)/4) {
   int re = LDBL_MANT_DIG-1-(p*4);
   long double round = 1ULL<<re;

This expression overflows. re can be as large as 108 but 1ULL has
64-bit type.

There's probably some reasonable way to write it that works for
reasonable sizes of long double (I think it's safe to assume IEEE quad
is the max that will ever be supported), but I presume I wrote it with
the original inefficient loop to be fully general to arbitrary
precision.

It could just be written to use scalbn, I think. At the time I
probably was trying to avoid dependency on a libm that could have been
separate from libc (bad historical thing). IIRC there was a time when
frexp was not used either.

Does the attached look ok?

Rich
<pct_a.diff>