questionnable nint() behavior

Timothy C Prince tprince@myrealbox.com
Tue Mar 21 00:52:00 GMT 2006



-----Original Message-----
From: Steve Kargl <sgk@troutmask.apl.washington.edu>
To: Dominique Dhumieres <dominiq@lps.ens.fr>
Date: Mon, 20 Mar 2006 16:35:03 -0800
Subject: Re: questionnable nint() behavior

On Mon, Mar 20, 2006 at 10:24:12PM +0100, Dominique Dhumieres wrote:
> I have found in my file the following program
> 
>       program main
>       real x, y
>       x = 8388609.0
>       y = 0.4999999701976776123046875
>       print '(A6,F9.1,A3,I8)', 'nint (', x, ') =', nint (x)
>       print '(A6,F10.8,A3,I2,A18,L1)', 'nint (', y, ') =', nint (y),    &
>      &      ', where y < 0.5 = ', y < 0.5
>       end
> 
> xlf, ifc, and g95 return
> 
> nint (8388609.0) = 8388609
> nint (0.49999997) = 0, where y < 0.5 = T
> 
> IMHO the correct result, g77, pgf, and gfc return
> 
> nint (8388609.0) = 8388610
> nint (0.49999997) = 1, where y < 0.5 = T
> 
> as if nint(a) was implemented as int(a+0.5) without
> taking care of rounding mode (rounding to nearest
> rounds to even on tie). I don't remember where I borrowed
> the code, but I thought the problem was fixed in gfortran.
>

In trans-intrinsics.c, you'll find

/* This is needed because the gcc backend only implements FIX_TRUNC_EXPR
   NINT(x) = INT(x + ((x > 0) ? 0.5 : -0.5)).  */

static tree
build_round_expr (stmtblock_t * pblock, tree arg, tree type)

A much better algorithm may be found in FreeBSD's round[fl].c. :-)

--
______________________________________
We could get away with this, by promoting float to double, and double to long double, as happens already in gfortran i386 for some -march= values. g77 (libf2c) had a more careful and slow version, using floor().
Tim Prince



More information about the Fortran mailing list