Fortran: fix compile-time simplification of SET_EXPONENT [PR109511]
Checks
Commit Message
[Now with subject added...]
Dear all,
the compile-time simplification of intrinsic SET_EXPONENT was
broken since the early days of gfortran for argument X < 1
(including negative X) and for I < 0. I identified and fixed
several issues in the implementation. The testcase explores
argument space comparing compile-time and runtime results and
is checked against Intel.
Regtested on x86_64-pc-linux-gnu. OK for mainline?
This is not a regression, but can lead to wrong code.
Would it be OK to backport to open branches?
Thanks,
Harald
Comments
On Fri, Apr 14, 2023 at 08:59:24PM +0200, Harald Anlauf via Fortran wrote:
>
> the compile-time simplification of intrinsic SET_EXPONENT was
> broken since the early days of gfortran for argument X < 1
> (including negative X) and for I < 0. I identified and fixed
> several issues in the implementation. The testcase explores
> argument space comparing compile-time and runtime results and
> is checked against Intel.
>
> Regtested on x86_64-pc-linux-gnu. OK for mainline?
>
Yes, it is certainly better than the current situation.
> This is not a regression, but can lead to wrong code.
> Would it be OK to backport to open branches?
Sure. Looks simply and fairly specific.
I was wondering about the difference between set_exponent()
and scale(), and found that set_exponent() talks about IEEE
values while scale() doesn't. I'm wondering if we should
add the IEEE special cases to the testsuite. Of particular
note, I doubt that this is true:
If X is an IEEE NaN, the result is the same NaN.
program foo
real x, y
x = 1
y = x - x
x = (x - x) / y
print '(F4.0,1X,Z8.8)', x, x
y = set_exponent(x,1)
print '(F4.0,1X,Z8.8)', y, y
end program foo
gfcx -o z a.f90 && ./z
NaN FFC00000
NaN 7FC00000
Those are not the same NaN. The second is a qNaN.
The first looks like a qNaN with the sign bit set.
On Fri, Apr 14, 2023 at 12:33:18PM -0700, Steve Kargl via Fortran wrote:
>
> If X is an IEEE NaN, the result is the same NaN.
>
A better testcase as gfortran will quiet a NaN on assignment.
program foo
integer i
equivalence(i,y)
i = int(z'7FC0BEEF',4) ! Add payload to NaN.
print '(F4.0,2(1X,Z8.8))', y, y, set_exponent(y,1)
end program foo
gfortran -o z a.f90 && ./z
N 7FC0BEEF 7FC00000
Hi Steve,
On 4/14/23 21:33, Steve Kargl via Gcc-patches wrote:
> On Fri, Apr 14, 2023 at 08:59:24PM +0200, Harald Anlauf via Fortran wrote:
>>
>> the compile-time simplification of intrinsic SET_EXPONENT was
>> broken since the early days of gfortran for argument X < 1
>> (including negative X) and for I < 0. I identified and fixed
>> several issues in the implementation. The testcase explores
>> argument space comparing compile-time and runtime results and
>> is checked against Intel.
>>
>> Regtested on x86_64-pc-linux-gnu. OK for mainline?
>>
>
> Yes, it is certainly better than the current situation.
>
>> This is not a regression, but can lead to wrong code.
>> Would it be OK to backport to open branches?
>
> Sure. Looks simply and fairly specific.
OK, thanks, will proceed.
> I was wondering about the difference between set_exponent()
> and scale(), and found that set_exponent() talks about IEEE
> values while scale() doesn't. I'm wondering if we should
> add the IEEE special cases to the testsuite. Of particular
> note, I doubt that this is true:
>
> If X is an IEEE NaN, the result is the same NaN.
>
> program foo
> real x, y
> x = 1
> y = x - x
> x = (x - x) / y
> print '(F4.0,1X,Z8.8)', x, x
> y = set_exponent(x,1)
> print '(F4.0,1X,Z8.8)', y, y
> end program foo
>
> gfcx -o z a.f90 && ./z
> NaN FFC00000
> NaN 7FC00000
>
> Those are not the same NaN. The second is a qNaN.
> The first looks like a qNaN with the sign bit set.
Until now there was no testing at all of SET_EXPONENT in the testsuite.
It would be really good to have better coverage of compile-time and
runtime behavior of the intrinsics and checking consistency ... ;-)
I think you have much more experience in that area. (Hint!)
Cheers,
Harald
On Fri, Apr 14, 2023 at 09:49:22PM +0200, Harald Anlauf wrote:
>
> On 4/14/23 21:33, Steve Kargl via Gcc-patches wrote:
> > I was wondering about the difference between set_exponent()
> > and scale(), and found that set_exponent() talks about IEEE
> > values while scale() doesn't. I'm wondering if we should
> > add the IEEE special cases to the testsuite. Of particular
> > note, I doubt that this is true:
> >
> > If X is an IEEE NaN, the result is the same NaN.
> >
> > program foo
> > real x, y
> > x = 1
> > y = x - x
> > x = (x - x) / y
> > print '(F4.0,1X,Z8.8)', x, x
> > y = set_exponent(x,1)
> > print '(F4.0,1X,Z8.8)', y, y
> > end program foo
> >
> > gfcx -o z a.f90 && ./z
> > NaN FFC00000
> > NaN 7FC00000
> >
> > Those are not the same NaN. The second is a qNaN.
> > The first looks like a qNaN with the sign bit set.
>
> Until now there was no testing at all of SET_EXPONENT in the testsuite.
> It would be really good to have better coverage of compile-time and
> runtime behavior of the intrinsics and checking consistency ... ;-)
>
> I think you have much more experience in that area. (Hint!)
>
I might have some experience. :-) Unfortunately, I have no time
(for at least 4-5 years) to jump down that rabbit hole (as I would
try to fix things ;-).
However, writing testcases to exercise the intrinsic subprograms
is something that the LURKERS here in the mailing list might
comtemplate. Any takers?
Note returning the 'same Nan" is not special to set_exponent().
At least, fraction() and rrspacing() have "If X is an IEEE NaN,
the result is that NaN."
Oh, and thanks for your relentless assault on bugs.
On Fri, 14 Apr 2023 20:59:24 +0200
Harald Anlauf via Fortran <fortran@gcc.gnu.org> wrote:
> - mpfr_clears (absv, log2, pow2, frac, NULL);
> + mpfr_clears (absv, log2, exp, pow2, frac, NULL);
Thanks for taking care of that leak.
Consider that hunk dropped from
https://inbox.sourceware.org/fortran/D9328809-A533-4D6B-A474-DE4C582D9C02@gmail.com/
although i had it in the order of allocation to ease reading..
cheers,
From fa4cb42870df60deb8888dbd51e2ddc6d6ab9e6a Mon Sep 17 00:00:00 2001
From: Harald Anlauf <anlauf@gmx.de>
Date: Fri, 14 Apr 2023 20:45:19 +0200
Subject: [PATCH] Fortran: fix compile-time simplification of SET_EXPONENT
[PR109511]
gcc/fortran/ChangeLog:
PR fortran/109511
* simplify.cc (gfc_simplify_set_exponent): Fix implementation of
compile-time simplification of intrinsic SET_EXPONENT for argument
X < 1 and for I < 0.
gcc/testsuite/ChangeLog:
PR fortran/109511
* gfortran.dg/set_exponent_1.f90: New test.
---
gcc/fortran/simplify.cc | 12 +++----
gcc/testsuite/gfortran.dg/set_exponent_1.f90 | 36 ++++++++++++++++++++
2 files changed, 42 insertions(+), 6 deletions(-)
create mode 100644 gcc/testsuite/gfortran.dg/set_exponent_1.f90
@@ -7364,7 +7364,7 @@ gfc_simplify_set_exponent (gfc_expr *x, gfc_expr *i)
{
gfc_expr *result;
mpfr_t exp, absv, log2, pow2, frac;
- unsigned long exp2;
+ long exp2;
if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
return NULL;
@@ -7396,19 +7396,19 @@ gfc_simplify_set_exponent (gfc_expr *x, gfc_expr *i)
mpfr_abs (absv, x->value.real, GFC_RND_MODE);
mpfr_log2 (log2, absv, GFC_RND_MODE);
- mpfr_trunc (log2, log2);
+ mpfr_floor (log2, log2);
mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
/* Old exponent value, and fraction. */
mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
- mpfr_div (frac, absv, pow2, GFC_RND_MODE);
+ mpfr_div (frac, x->value.real, pow2, GFC_RND_MODE);
/* New exponent. */
- exp2 = (unsigned long) mpz_get_d (i->value.integer);
- mpfr_mul_2exp (result->value.real, frac, exp2, GFC_RND_MODE);
+ exp2 = mpz_get_si (i->value.integer);
+ mpfr_mul_2si (result->value.real, frac, exp2, GFC_RND_MODE);
- mpfr_clears (absv, log2, pow2, frac, NULL);
+ mpfr_clears (absv, log2, exp, pow2, frac, NULL);
return range_check (result, "SET_EXPONENT");
}
new file mode 100644
@@ -0,0 +1,36 @@
+! { dg-do run }
+! PR fortran/109511
+! Check compile-time simplification of SET_EXPONENT against runtime
+
+program exponent
+ implicit none
+ integer :: i
+ i = 0
+ print *, i, set_exponent(1., 0), set_exponent(1., i)
+ if (set_exponent(1., 0) /= set_exponent(1., i)) stop 1
+ i = 1
+ print *, i, set_exponent(1., 1), set_exponent(1., i)
+ if (set_exponent(1., 1) /= set_exponent(1., i)) stop 2
+ i = 2
+ print *, i, set_exponent(-1.75, 2), set_exponent(-1.75, i)
+ if (set_exponent(-1.75, 2) /= set_exponent(-1.75, i)) stop 3
+ print *, i, set_exponent(0.1875, 2), set_exponent(0.1875, i)
+ if (set_exponent(0.1875, 2) /= set_exponent(0.1875, i)) stop 4
+ i = 3
+ print *, i, set_exponent(0.75, 3), set_exponent(0.75, i)
+ if (set_exponent(0.75, 3) /= set_exponent(0.75, i)) stop 5
+ i = 4
+ print *, i, set_exponent(-2.5, 4), set_exponent(-2.5, i)
+ if (set_exponent(-2.5, 4) /= set_exponent(-2.5, i)) stop 6
+ i = -1
+ print *, i, set_exponent(1., -1), set_exponent(1., i)
+ if (set_exponent(1., -1) /= set_exponent(1., i)) stop 7
+ i = -2
+ print *, i, set_exponent(1.125, -2), set_exponent(1.125, i)
+ if (set_exponent(1.125, -2) /= set_exponent(1.125, i)) stop 8
+ print *, i, set_exponent(-0.25, -2), set_exponent(-0.25, i)
+ if (set_exponent(-0.25, -2) /= set_exponent(-0.25, i)) stop 9
+ i = -3
+ print *, i, set_exponent(0.75, -3), set_exponent(0.75, i)
+ if (set_exponent(0.75, -3) /= set_exponent(0.75, i)) stop 10
+end program exponent
--
2.35.3