Hi Segher,

Quick version: Jump to the new patch, which I like much more.
Longer version:
On 19.05.21 17:15, Segher Boessenkool wrote:
>>     real(16)               :: y   ! 128bit REAL
>>     integer(16), parameter :: k2 = nint (2 / epsilon (y), kind(k2))
>>     integer(16), parameter :: m2 = 10384593717069655257060992658440192_16
>>     !2**113
>>     if (k2 /= m2) stop 3
>>
>> On x86_64-linux-gnu, k2 == m2 — but on powerpc64le-linux-gnu,
>> k2 == 2**106 instead of 2**113.
>>
>> My solution is to permit also 2**106 besides 2**113.
> I do not understand Fortran well enough, could you explain what the code
> is supposed to do?

First, 2_16 means the integer '2' of the integer kind '16', i.e. int128_t type.

The original bug report (PR96983) was that 'nint'
with the 16byte floating point/quad-precision was giving an ICE
and the complaint was that there was no testcase for the value.

And I think this testcase tries to ensure that the result of 'nint'
both at compile time and at runtime matches what should be the result.

(Quotes from the Fortran standard, augumenty by what the source code does.)

'nint' does "__builtin_round" is available: "The result is the integer
nearest A, or if there are two integers equally near A, the result
is whichever such integer has the greater magnitude" – and in this
testcase, the argument is a quad-precision/16byte/128bit floating-point
number and the result is a 128bit integer.

('a**b' is the Fortran syntax for: 'a' raised to the power of 'b'.)

This testcase does:

nint(2/epsilon(y)). Here, 'epsilon' is the
"Model number that is small compared to 1."
Namely: b**(p-1) = '(radix)**(1-digits)'
alias 'real_format *fmt = REAL_MODE_FORMAT (mode)'
with radix = fmt->b  and digits = fmt->p;

[b**(p-1) is from the Fortran standard but 'b' and 'p' also match the
ME/target names, while radix/digits matches the FE names and also the
Fortran intrinsic inquiry function names.]

This is for radix = 2 equivalent to:

2/2**(1-digits) = 2*2**(digits-1) = 2**(digits)

On x86-64, digits == fmt->p == 113.

Our powerpc64le gives digits == 106.

  * * *

Having written all this, I wonder why we don't just
rely on the assumption that '2**digit(x)' works – and use this
to generate the valid.

Namely, as the attached updated patch does.

As I like that patch and believe it is obvious, I intent to
commit it as such – unless there are further comments.

It passes on both x86-64-gnu-linux and powerpc64le-none-linux-gnu.
I think the radix == 2 is a good bet, but if we ever run into issues,
it can also be changed to use radix(...) as well ...

Tobias

-----------------
Mentor Graphics (Deutschland) GmbH, Arnulfstrasse 201, 80634 München Registergericht München HRB 106955, Geschäftsführer: Thomas Heurung, Frank Thürauf