bugGNU Octave - Bugs: bug #35611, Bug in bincoeff with large n and k

 
 

bug #35611: Bug in bincoeff with large n and k

Submitter:  None
Submitted:  Thu 23 Feb 2012 03:28:28 PM UTC
   
 
Category:  None Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Inaccurate Result
Status:  Wont Fix Assigned to:  None
Originator Name:  Marcus Meyer Originator Email:  -email is unavailable-
Open/Closed:  * Closed Release:  * dev
Operating System:  * Microsoft Windows Fixed Release:  None
Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Fri 24 Feb 2012 04:31:33 PM UTC, comment #4: 

Marco,

This is just what I suspected.  The underlying routine for gammaln probably differs by a few eps depending on the backend library.  For your system, the result is slightly under the true result and rounding brings it up to the correct value.  For the bug reporter's system the result is slightly higher than the true result and rounding causes an error.

To recap, there isn't any consistency that Octave can exploit and so this bug can't be "fixed".  On the other hand, there is a clear workaround which is to use the appropriate function, nchoosek for integers, and bincoeff for real values.

Rik <rik5>
Group administrator
Fri 24 Feb 2012 08:16:51 AM UTC, comment #3: 

Dear all,

just to inform that with 3.6.0 in linux, I get, for bincoeff,  1.73103094564398e+13 before rounding, 1.73103094564400e+13 after rounding and for nchoosek 1.73103094564400e+13. The same is true for 3.2.4.
This time I'm lucky that the results are right, but should I think that my system is different from all the others?

Cheers,

Marco

Marco Caliari <caliari>
Group Member
Thu 23 Feb 2012 05:58:19 PM UTC, comment #2: 

Check the documentation for nchoosek and bincoeff.  They compute the same mathematical value, but the underlying algorithms are very different.  In particular, nchoosek is optimized for integer inputs and gives the correct answer in both Octave and Matlab.  bincoeff is a much broader function and can handle real inputs as well as negative inputs.

Quoting from the documentation for bincoeff:


In most cases, the `nchoosek' function is faster for small scalar
integer arguments.  It also warns about loss of precision for big
arguments.


Quoting from the documentation for nchoosek:


`nchoosek' works only for non-negative, integer arguments.  Use
`bincoeff' for non-integer and negative scalar arguments, or for
computing many binomial coefficients at once with vector inputs
for N or K.


nchoosek (100,10) = 1.7e13 which is representable exactly as an integer in doubele precision arithmetic (you can check bitmax() which is 9e15 for doubles).  Thus, nchoosek is really the correct function to use for this example.

The line to take a look at is 98 of bincoeff.m which I reproduce below.


## Clean up rounding errors.
b(n_int) = round (b(n_int));


Before the rounding, the value of b is 17310309456440.8.  In this case, if we used truncation, fix(), rather than round we would get the correct answer.  However, I'm note sure this is a general result.  I believe round was chosen so that 50% of the time we would guess low and 50% of the time we would guess high and on average the algorithm would be more accurate.  It would be valuable if you could prove any relation that might exist between the input and the fractional part of the bincoeff algorithm.  For example, is it always the case for integers less than bitmax that the result should use fix() rather than round()?


Rik <rik5>
Group administrator
Thu 23 Feb 2012 03:39:08 PM UTC, comment #1: 

The code for bincoeff isn't exceptionally complicated and it's fairly self-contained:

   http://hg.savannah.gnu.org/hgweb/octave/file/d17237256856/scripts/miscellaneous/bincoeff.m#l55

Do you want to attempt a fix yourself? There's an obvious rounding error.

Jordi GutiƩrrez Hermoso <jordigh>
Group Member
Thu 23 Feb 2012 03:28:28 PM UTC, original submission:  

octave-3.2.4.exe:10> bincoeff(100,10)
ans =  17310309456441

The correct result is:

octave-3.2.4.exe:13> 10*11*98*97*19*31*23*13*94
ans =  17310309456440

Matlab gives the correct result with nchoosek(100,10).
Maybe bincoeff can be improved for large input?

Best regards and thanks
Marcus

Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by caliari (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by jordigh (Posted a comment)
  • -email is unavailable- added by None (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 4 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2012-02-24 rik5 StatusConfirmed Wont Fix
        Open/ClosedOpen Closed
    2012-02-23 jordigh StatusNone Confirmed
        Release3.2.4 dev

    Back to the top

    Powered by Savane 3.13-758e.
    Corresponding source code