[Rd] pcauchy precision (PR#6756)
maechler at stat.math.ethz.ch
maechler at stat.math.ethz.ch
Tue Apr 13 09:07:59 CEST 2004
>>>>> "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 <><
More information about the R-devel
mailing list