From help-octave-request at bevo dot che dot wisc dot edu  Thu Jan 29 18:36:47 1998
Subject: Re: potential bug in betai; any hints at a fix or workaround?
From: king at cogsci dot ucsd dot edu (Jonathan King)
To: jwe at bevo dot che dot wisc dot edu
Cc: king at cogsci dot ucsd dot edu, help-octave@bevo.che.wisc.edu,         bug-octave@bevo.che.wisc.edu
Date: Thu, 29 Jan 1998 16:36:06 -0800

jwe at bevo dot che dot wisc dot edu wrote:
>On 29-Jan-1998, Jonathan King <king at cogsci dot ucsd dot edu> wrote:
>

[etc., code, etc. snipped]

>| Betai seems to be the problem here; both Matlab 5.1 and a
>| hand-written C version essentially yanked from Numerical Recipes
>| provided what seem to be the correct answers.

[more snipped]

>| Is there a fix for this?
>| An obvious work-around?
>
>For Octave 2.0.10, the M-file implemntation for betai is replaced by a
>function from the Slatec library.  Here's what I get with my current
>sources:
>
>  octave:4> rtop2 ([.56141,.70710], 144)
>  ans =
>
>     1.6902e-13   1.9771e-23
>
>Is that about right?

A-yup.  Or, rather, that's about what Matlab and the hand-coded C
gives you.  Octave 2.0.9 gave:

    octave.bin:3> rtop2 ([.56141, .70710], 144)
    ans =
    
        1.3371e-10   -1.8682e-01

Note that the second answer is off by rather a lot... Those of you
you want to see a clear picture of how wrong this was might try:

    plot(rtop2 ([linspace(-.99,.99,200)], 144))

But this gets much worse for larger second arguments in rtop2, so:

    plot(rtop2 ([linspace(-.99,.99,200)], 200))

is truly scary looking.

By this point, Octave's answers are at least one point off by not
just the sign, but over 40 orders of magnitude.  But the answers
outside of this range appear to be completely correct, as far as I
can tell.  It's really spooky...

jking