Problem with rounding using lrint

Applications Programming
1

I run into a situation when integer result of 81 * 0.5 is 40 (expected 41) when using lrint() function.

When I use a simple +0.5 approach results are right.

The code is:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
 int main(int argc, char *argv[])
{
        const int       X       = 81;
        double          prc     = (argc == 2) ? atof(argv[1]) : 0.5;
        double          dval    = X * prc;
        long            ival    = lrint(dval);
        printf("%i * %.16f = %.16f  IVAL = %li  ", X, prc, dval, ival);
        ival = (X * prc) + 0.5;
        printf("USING +0.5 IVAL = %li\n", ival);
        return(0);
}

results are:

~$ c 0.5
81 * 0.5000000000000000 = 40.5000000000000000  IVAL = 40  USING +0.5 IVAL = 41

Using gcc ver 5.4.0 on Linux, running on 2 years old Xeon processor, if that matters. I checked floating point rounding defaults in compiler and it is a regular "round to the nearest integer" rule.

My question is - should I abandon rint-family functions and use +0.5 approach?

2

You need to determine what rounding rules you want to use for tie-breaking cases (see Wikipedia's discussion on rounding)
and you need to realize that even in double precision floating point, the result of a floating point calculation is not always exact (even when computing a value using decimal calculations would be exact). For example, the awk program (which uses double precision floating point for its calculations):

printf '.1 .5\n.3 .5\n' | awk '{printf("%.40f\n", $1 * $2)}'

produces the output:

0.0500000000000000027755575615628913510591
0.1499999999999999944488848768742172978818

not the output you might expect:

0.0500000000000000000000000000000000000000
0.1500000000000000000000000000000000000000

Which shows two examples where decimal arithmetic produces exact results, but binary arithmetic produced one result that was a little bit high and one result that was a little bit low.

3

Thanks.
I'd love to hear what others use in C/C++, that is rint() family functions? home grown solutions?

4

Normally you call fegetround() to determine which rounding mode of the four possible ones is active.

fesetround() sets the mode. Only the fpu supported modes are available on a particular box.

For complete details try:
The GNU C Library: Rounding