[Rd] pcauchy precision (PR#6756)

[maechler at stat.math.ethz.ch]

>>>>> "Morten" == Morten Welinder <terra at gnome.org>
>>>>>     on Sun, 11 Apr 2004 19:04:01 +0200 (CEST) writes:

    Morten> Full_Name: Morten Welinder Version: snapshot   OS:..
    Morten> Submission from: (NULL) (65.213.85.129)


    Morten> Two things are wrong.

    Morten> 1. There is nan test outside IEEE_754.

No, that's not wrong.  nmath.h *does* define ISNAN() also "outside IEEE_754"

    Morten> 2. The meat part of the function should really be something
    Morten> like...

    Morten>     if (!lower_tail) x = -x;

good idea for code simplification.  No bug in the current code though.

    Morten>     if (fabs (x) > 1) { double temp = atan (1 / x) /
    Morten> M_PI; return (x > 0) ? R_D_Clog (temp) : R_D_val
    Morten> (-temp); } else return R_D_val (0.5 + atan (x) /
    Morten> M_PI);

    Morten> ...instead of the current heavily truncated series
    Morten> expansion.  The above is much simpler 

"simpler": yes.

    Morten> and more precise.

That's quite a strong claim  which you haven't supported by any evidence!

I don't think there's a precision loss in the current code,
particularly not where the  "heavily truncated series" is used.
I've been very careful to use it only where it's fully accurate.
I would expect that your code may be more accurate (2-3 digits)
when 1 <= |x| <= x_big := 5478.3.  But the current code might well be
more accurate for |x| >= x_big - particularly for cases where
the C "libm" atan() implementation isn't quite perfect.

    Morten>   (Thanks to Ian Smith for pointing the
    Morten> in-hindsight- obvious 1/x trick out.)

I agree that's a cute idea / "trick".

Still, I haven't seen any bug in the current code!

Martin Maechler <maechler at stat.math.ethz.ch>	http://stat.ethz.ch/~maechler/
Seminar fuer Statistik, ETH-Zentrum  LEO C16	Leonhardstr. 27
ETH (Federal Inst. Technology)	8092 Zurich	SWITZERLAND
phone: x-41-1-632-3408		fax: ...-1228			<><