Skip site navigation (1)Skip section navigation (2)
Date:      Thu, 5 Dec 2013 11:52:25 -0800
From:      Steve Kargl <sgk@troutmask.apl.washington.edu>
To:        Filipe Maia <filipe.c.maia@gmail.com>
Cc:        freebsd-numerics@freebsd.org
Subject:   Re: dead code in lgamma_r[f]?
Message-ID:  <20131205195225.GA65732@troutmask.apl.washington.edu>
In-Reply-To: <CAN5hRiV6N1arc4RDv=9JbRiXp-J9o3WzAbeZSOpAxks2ZeG%2B_w@mail.gmail.com>
References:  <20131205173700.GA64575@troutmask.apl.washington.edu> <20131205182324.GA65135@troutmask.apl.washington.edu> <CAN5hRiV6N1arc4RDv=9JbRiXp-J9o3WzAbeZSOpAxks2ZeG%2B_w@mail.gmail.com>

next in thread | previous in thread | raw e-mail | index | archive | help
On Thu, Dec 05, 2013 at 08:06:22PM +0100, Filipe Maia wrote:
> What if ix == 0x433xxxxx?

See below.

> Wouldn't the if clause be false?

See below.

> > Index: src/e_lgamma_r.c
> > ===================================================================
> > --- src/e_lgamma_r.c    (revision 1427)
> > +++ src/e_lgamma_r.c    (working copy)
> > @@ -173,19 +173,15 @@
> >       */
> >         z = floor(y);
> >         if(z!=y) {                              /* inexact anyway */
> > -           y  *= 0.5;
> > -           y   = 2.0*(y - floor(y));           /* y = |x| mod 2.0 */
> > -           n   = (int) (y*4.0);
> > +           y /= 2;
> > +           y  = 2*(y - floor(y));              /* y = |x| mod 2.0 */
> > +           n  = (int)(y*4);

If you have a y value that corresponds to ix=0x433xxxxx where 
one of the x is nonzero, I believe that you end up going through
this portion of code.

> >         } else {

You end up here if and only if z==y, which means y is integral.
As noted in the original email, the 'if(ix>=0x43400000)' below
is dead code because this code in __ieee754_lgamma_r():

   if(ix>=0x43300000) 	/* |x|>=2**52, must be -integer */
       return one/zero;
   t = sin_pi(x);

prevents the call to sin_pi(x).

> > -            if(ix>=0x43400000) {
> > -                y = zero; n = 0;                 /* y must be even */
> > -            } else {
> > -                if(ix<0x43300000) z = y+two52; /* exact */

If we again look at the code from __ieee754_lgamma_r(), we see
that sin_pi() is called if ix < 0x43300000, so by the time we
arrive at the 'if(ix<0x43300000)' statement we already know that
the condition is true.

> > -               GET_LOW_WORD(n,z);
> > -               n &= 1;
> > -                y  = n;
> > -                n<<= 2;
> > -            }
> > +           z  = y+two52;                       /* exact */
> > +           GET_LOW_WORD(n,z);
> > +           n &= 1;
> > +           y  = n;
> > +           n<<= 2;
> >          }
> >         switch (n) {
> >             case 0:   y =  __kernel_sin(pi*y,zero,0); break;

-- 
Steve



Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20131205195225.GA65732>