Re: Confused with frexp and atangent (mathfp branch)

[Roman Belenov <rbelenov at yandex dot ru>]

Thanks, the patch fixes the problem.

Jeff Johnston <jjohnstn@redhat.com> writes:

> Roman Belenov wrote:
>> I'm working with a toolchain that uses newlib with mathfp branch used for
>> libm. I've encountered some problems - it seems that newlib functions give
>> incorrect results.
>> First, frexp(0., &exp) writes -1022 to exp (instead of 0 required by
>> standard)
>> as a result of calculating
>>     *exp = ((hd & 0x7ff00000) >> 20) - 1022;
>> (since 0. is represented as all-zero-bits-value). Probably explicit test for
>> zero (together with NAN and INF) is required.
>> 
>
> Yes, a fix is required.  The check for special values is done on the result
> instead of the input.
>
>> Second, there is confusing usage of branch variable in atangent() function,
>> resulting in incorrect results from atan2. As a result of the following code
>>   if (arctan2)
>>     {
>>       if (u < 0.0 || branch == 2)
>>         res = __PI - res;
>>       if (v < 0.0 || branch == 1)
>>         res = -res;
>>     }
>> result is always set to PI for very small values of y/x (including y==0. due
>> to aforementioned behaviour of frexp) and to -PI/2 for very high numbers
>> ("very" means that preceding exponent test sets branch to nonzero value),
>> since signs of u and v are effectively ignored. Seems like checks for branch
>> here are erroneous. I don't have the textbook mention in s_atangent.c
>> available, so can't check whether this snippet comes from its text and whether
>> there is any rationale.
>> Can anybody elaborate on these issues ?
>> 
>
> Another bug.  Yes, it appears the check for branch is extraneous and causing
> the problem.  Please try the accompanying patch.
>
> -- Jeff J.
>
>
>
> ? mathfp.patch
> Index: s_atangent.c
> ===================================================================
> RCS file: /cvs/src/src/newlib/libm/mathfp/s_atangent.c,v
> retrieving revision 1.2
> diff -u -p -r1.2 s_atangent.c
> --- s_atangent.c	20 Oct 2003 18:46:38 -0000	1.2
> +++ s_atangent.c	9 Jan 2006 22:29:18 -0000
> @@ -197,9 +197,9 @@ _DEFUN (atangent, (double, double, doubl
>  
>    if (arctan2)
>      {
> -      if (u < 0.0 || branch == 2)
> +      if (u < 0.0)
>          res = __PI - res;
> -      if (v < 0.0 || branch == 1)
> +      if (v < 0.0)
>          res = -res;
>      }
>    else if (x < 0.0)
> Index: s_frexp.c
> ===================================================================
> RCS file: /cvs/src/src/newlib/libm/mathfp/s_frexp.c,v
> retrieving revision 1.2
> diff -u -p -r1.2 s_frexp.c
> --- s_frexp.c	20 Oct 2003 18:46:38 -0000	1.2
> +++ s_frexp.c	9 Jan 2006 22:29:18 -0000
> @@ -82,6 +82,17 @@ double frexp (double d, int *exp)
>    double f;
>    __uint32_t hd, ld, hf, lf;
>  
> +  /* Check for special values. */
> +  switch (numtest (d))
> +    {
> +      case NAN:
> +      case INF:
> +        errno = EDOM;
> +      case 0:
> +        *exp = 0;
> +        return (d);
> +    }
> +
>    EXTRACT_WORDS (hd, ld, d);
>  
>    /* Get the exponent. */
> @@ -94,16 +105,6 @@ double frexp (double d, int *exp)
>  
>    INSERT_WORDS (f, hf, lf);
>  
> -  /* Check for special values. */
> -  switch (numtest (f))
> -    {
> -      case NAN:
> -      case INF:
> -        errno = EDOM;
> -        *exp = 0;
> -        return (f);
> -    }
> -
>    return (f);
>  }
>  
> Index: sf_frexp.c
> ===================================================================
> RCS file: /cvs/src/src/newlib/libm/mathfp/sf_frexp.c,v
> retrieving revision 1.1.1.1
> diff -u -p -r1.1.1.1 sf_frexp.c
> --- sf_frexp.c	17 Feb 2000 19:39:52 -0000	1.1.1.1
> +++ sf_frexp.c	9 Jan 2006 22:29:18 -0000
> @@ -24,6 +24,17 @@ float frexpf (float d, int *exp)
>    float f;
>    __int32_t wf, wd;
>  
> +  /* Check for special values. */
> +  switch (numtestf (d))
> +    {
> +      case NAN:
> +      case INF:
> +        errno = EDOM;
> +      case 0:
> +        *exp = 0;
> +        return (d);
> +    }
> +
>    GET_FLOAT_WORD (wd, d);
>  
>    /* Get the exponent. */
> @@ -35,16 +46,6 @@ float frexpf (float d, int *exp)
>  
>    SET_FLOAT_WORD (f, wf);
>  
> -  /* Check for special values. */
> -  switch (numtestf (f))
> -    {
> -      case NAN:
> -      case INF:
> -        errno = EDOM;
> -        *exp = 0;
> -        return (f);
> -    }
> -
>    return (f);
>  }
>  

-- 
 							With regards, Roman.