[PATCH, middle-end, i386]: reciprocal rsqrt pass + full recip x86 backend support

[Richard Guenther]

On 6/15/07, Uros Bizjak <ubizjak@gmail.com> wrote:
> On 6/15/07, Tobias Burnus <burnus@net-b.de> wrote:
>
> > > ! ---       Use sliding step-size to sample nearfield more frequently
> > >            DELN = 2.*FLOAT(NN-NTR)/FLOAT(NTR*(NTR-1))
> > I filled this now as http://gcc.gnu.org/bugzilla/show_bug.cgi?id=32352
> >
> > With -mrecip and -bounds-check, I get an error a bit earlier:
> >
> > At line 35436 of file aermod.f90
> > Fortran runtime error: Array reference out of bounds for array 'xt', lower
> > bound of dimension 1 exceeded
> >
> > This is SUBROUTINE NUMRISE:
> >          DO IN = NBEG , NEND
> > ! ---       Define coordinates of plume relative to bldg.-defined origin
> >             XBB = XT(IN-1) - XBADJ
>
> I have fond the problem in following lines:
>
> ! --- Process new position
>       S = S + DS
>       IF ( FLOAT(NNP/NP).EQ.FLOAT(NNP)/XNP ) THEN          <<<< here
>          NN = NNP/NP
>
> The problem is one ULP of the rcpss insn (+NR step), illustrated by
> the following testcase:
>
> --cut here--
> float recip(float a)
> {
>   float xx = 1.0f/a;
>   printf("%a\n", xx);
>
>   return xx;
> }
>
> int NNP = 1;
> int NP = 1;
> float XNP = 1.0f;
>
> int main()
> {
>   float a, b;
>   int c;
>
>   a = (float)(NNP/NP);
>   b = (float)(NNP)*recip(XNP);
>
>   printf("%a %a %i\n", a, b, a == b);
>
>   return 0;
> }
> --cut here--
>
> gcc -O2 -msse2 -mfpmath=sse -ffast-math equal.c
> ./a.out
> 0x1p+0
> 0x1p+0 0x1p+0 1
>
> gcc -O2 -msse2 -mfpmath=sse -ffast-math -mrecip equal.c
> ./a.out
> 0x1.fffffep-1
> 0x1p+0 0x1.fffffep-1 0
>
> Actually, the recip pass is OK, it is just a precision problem. This
> is why -mrecip should stay separate option, not included in
> -ffast-math...

I agree.  Maybe even split -mrcp and -mrsqrt.  In this particular case
can we detect the division result is used in a comparison?

Richard.