Re: potential bug in betai; any hints at a fix or workaround?

[Jonathan King]

address@hidden wrote:
>On 29-Jan-1998, Jonathan King <address@hidden> 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