[Proj] Re: Re: Proj4 Bug (rtodms)

Glynn Clements glynn at gclements.plus.com
Thu Nov 9 13:44:54 PST 2006


Gerald I. Evenden wrote:

> > > The numbers have effectively right and exact representations in ten-bytes
> > > format, but the compiler instead of making r/60 makes r * 1/60 so it
> > > changes a little thing. I precise that this is done without any
> > > optimization (force no optimization) and with option "don't correct FDIV
> > > flaw" (for old processors) and "no quick floating points".
> >
> > In that case, the Borland compiler is entirely useless for any code
> > which makes non-trivial use of floating-point arithmetic.
> >
> > Even if this particular issue can be worked around, you are likely to
> > encounter similar issues time and time again.
> >
> > Any compiler written by someone who doesn't understand that x/y and
> > x*(1/y) aren't equivalent simply isn't worth the disk space it
> > occupies.
> 
> That argument does not hold too well.

No? Why?

> I altered the origin problem as 
> follows;
> 
> #include <math.h>
> #include <stdio.h>
> static int res1,res2;
> int main()
> {
>         double r=240.0;
> //      double tmp = 1./60.;
>         double tmp = 0.016666666666666666666666666666666666666;
> //      res1 = r / 60.;
>         res1 = r * tmp;
>         res2 = floor(r / 60.);
>         printf("%s\n",(res1==res2) ? "OK!" : "NOT OK!");
>         printf("res1: %d   res2: %d\n", res1, res2);
> //      getch();
> }
> 
> Both cases (1/60., 0.01666--) run as:
> 
> ./a.out
> OK!
> res1: 4   res2: 4
> 
> using a temporary value for "1/60." a la Borland and in both cases the 
> procedure ran "OK."

The issue isn't that x/y and x*(1/y) definitely won't be equal, it's
that they aren't *necessarily* equal in all cases (there are cases
where they *are* guaranteed to be equal, i.e. if 1/y is exactly
representable).

> I am still suspicious of the setting of the FPU's rounding flag.
> 
> This does raise in my mind: what should the default state of the rounding flag 
> be?

Appendix F of the C99 standard requires round-to-nearest during
translation (F.7.2) and at program startup (F.7.3). However, the
standard doesn't require an implementation to conform to appendix F
unless it indicates conformance by defining the macro
__STDC_IEC_559__.

In any case, 240.0/60.0 should be exactly equal to 4.0 regardless of
the rounding mode; anything else is a bug. The rounding mode is only
relevant if a result cannot be represented exactly.

Contrary to widespread belief, the FPU doesn't get the least
significant bits of a result from /dev/random ;)

-- 
Glynn Clements <glynn at gclements.plus.com>



More information about the Proj mailing list