From: Frank Reininghaus <frank78ac@googlemail.com>
To: gsl-discuss@sourceware.org
Subject: Examining bug #21838
Date: Thu, 14 Feb 2008 23:01:00 -0000 [thread overview]
Message-ID: <47B4C827.1040004@googlemail.com> (raw)
[-- Attachment #1: Type: text/plain, Size: 2107 bytes --]
Hi everyone,
I just tried to examine bug #21838 in gsl_cdf_fdist_P () (source file
cdf/fdist.c), reported in
http://lists.gnu.org/archive/html/bug-gsl/2007-09/msg00023.html
with the attached testcase program (testcases, expected values and
tolerance taken from cdf/test.c).
I examined the first testcase: It seems that the origin of the NaN is in
beta_cont_frac () in the source file cdf/beta_inc.c. The loop aborts
after the maximum number of iterations because fabs (delta_frac - 1.0)
does not get small enough, and NaN is returned.
The funny thing is that I could influence this behaviour in a way I
didn't expect:
1. Linking to the current git version of GSL which was compiled with the
default setting -O2 yielded:
(P - eP) [0] = nan, tolerance = 2.3283064365386963e-10
(P - eP) [1] = nan, tolerance = 2.3283064365386963e-10
(P - eP) [2] = nan, tolerance = 2.3283064365386963e-10
(P - eP) [3] = 4.1520314247533996e-05, tolerance = 2.3283064365386963e-10
2. Linking to GSL which was compiled without optimisation (-O0):
(P - eP) [0] = nan, tolerance = 2.3283064365386963e-10
(P - eP) [1] = nan, tolerance = 2.3283064365386963e-10
(P - eP) [2] = 1.8060553053089734e-11, tolerance = 2.3283064365386963e-10
(P - eP) [3] = nan, tolerance = 2.3283064365386963e-10
3. Running it with Valgrind (independent of the optimisation settings):
(P - eP) [0] = 6.9538819147396680e-12, tolerance = 2.3283064365386963e-10
(P - eP) [1] = -7.3852035598065413e-12, tolerance = 2.3283064365386963e-10
(P - eP) [2] = -1.1277645484142340e-12, tolerance = 2.3283064365386963e-10
(P - eP) [3] = 4.1520320904209207e-05, tolerance = 2.3283064365386963e-10
Maybe these (at least for me) surprising observations could be used to
isolate the cause of the bug. Does anyone have a good idea?
Regards,
Frank
P.S.: I compiled and ran the program on two systems I have access to and
got exactly the same results:
1.
* Athlon 64 X2 in 32 bit mode
* Kubuntu 7.10
* gcc (GCC) 4.1.3 20070929 (prerelease) (Ubuntu 4.1.2-16ubuntu2))
2.
* Pentium III
* Opensuse 10.2
* gcc (GCC) 4.1.2 20061115 (prerelease) (SUSE Linux)
[-- Attachment #2: testcase-bug-21838.c --]
[-- Type: text/x-csrc, Size: 843 bytes --]
#include <gsl/gsl_cdf.h>
#include <stdio.h>
#include <gsl/gsl_math.h>
/* tolerance taken from cdf/test.c */
#define TEST_TOL6 (1048576.0*GSL_DBL_EPSILON)
int main () {
double P [4]; /* P values calculated using gsl_cdf_fdist_P () */
double eP [4]; /* expected values taken from cdf/test.c */
P [0] = gsl_cdf_fdist_P (3.479082213465832574, 1, 4040712);
eP [0] = 0.93785072763723411967;
P [1] = gsl_cdf_fdist_P (3.002774644786533109, 1, 4040712);
eP [1] = 0.91687787379476055771;
P [2] = gsl_cdf_fdist_P (3.000854441173130827, 1, 4040712);
eP [2] = 0.91677930719813578619;
P [3] = gsl_cdf_fdist_P (3.000064021622133037, 1, 4040712);
eP [3] = 0.91678021757407962678;
int i;
for (i = 0; i < 4; i++)
printf ("(P - eP) [%i] = %.16e, tolerance = %.16e\n",
i, eP [i] - P [i], TEST_TOL6);
return 0;
}
next reply other threads:[~2008-02-14 23:01 UTC|newest]
Thread overview: 5+ messages / expand[flat|nested] mbox.gz Atom feed top
2008-02-14 23:01 Frank Reininghaus [this message]
2008-02-15 9:07 ` Jason Stover
2008-02-20 18:04 ` Brian Gough
2008-02-16 18:22 ` Brian Gough
2008-02-18 9:56 ` Brian Gough
Reply instructions:
You may reply publicly to this message via plain-text email
using any one of the following methods:
* Save the following mbox file, import it into your mail client,
and reply-to-all from there: mbox
Avoid top-posting and favor interleaved quoting:
https://en.wikipedia.org/wiki/Posting_style#Interleaved_style
* Reply using the --to, --cc, and --in-reply-to
switches of git-send-email(1):
git send-email \
--in-reply-to=47B4C827.1040004@googlemail.com \
--to=frank78ac@googlemail.com \
--cc=gsl-discuss@sourceware.org \
/path/to/YOUR_REPLY
https://kernel.org/pub/software/scm/git/docs/git-send-email.html
* If your mail client supports setting the In-Reply-To header
via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line
before the message body.
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).