Re: [PATCH] Improve the accuracy of tgamma for negative inputs (BZ 26983)

[Joseph Myers <joseph@codesourcery.com>]

On Fri, 26 Feb 2021, Paul Zimmermann wrote:

> ---
>  sysdeps/ieee754/dbl-64/e_gamma_r.c | 12 ++++++++++--
>  1 file changed, 10 insertions(+), 2 deletions(-)

Please include a longer commit message explaining the change made and why 
it's an appropriate fix for the issue seen.

> @@ -188,8 +189,15 @@ __ieee754_gamma_r (double x, int *signgamp)
>  			       ? __sin (M_PI * frac)
>  			       : __cos (M_PI * (0.5 - frac)));
>  	      int exp2_adj;
> -	      double tret = M_PI / (-x * sinpix
> -				    * gamma_positive (-x, &exp2_adj));
> +              double h1, l1, h2, l2;
> +              mul_split (&h1, &l1, sinpix, gamma_positive (-x, &exp2_adj));
> +              /* sinpix*gamma_positive(.) = h1 + l1 */
> +              mul_split (&h2, &l2, h1, x);
> +              /* h1*x = h2 + l2 */
> +              /* (h1 + l1) * x = h1*x + l1*x = h2 + l2 + l1*x */
> +              l2 += l1 * x;
> +              h2 += l2;
> +	      double tret = M_PI / -h2;

Note that glibc uses tabs for the initial multiple-of-8-columns of 
indentation; this is inserting code using spaces instead.

I'd expect similar changes for ldbl-96 and ldbl-128 at least, in the 
absence of some reason it makes sense for the implementations to diverge.

The test input that justifies this change as being needed should also be 
added to auto-libm-test-in by the patch (or XFAILs for that test input 
removed if already present), with any corresponding libm-test-ulps updates 
also being included for at least one architecture on which the patch was 
tested.  Likewise the test input you found giving too-large ulps for 
binary128, unless the binary128 problem is still present after making a 
corresponding fix to the implementation for ldbl-128.

-- 
Joseph S. Myers
joseph@codesourcery.com