Opened 10 years ago

Closed 10 years ago

#14247 closed defect (duplicate)

bell numbers in mpmath is incorrect

Reported by: Travis Scrimshaw Owned by: Sage Combinat CC user
Priority: major Milestone: sage-duplicate/invalid/wontfix
Component: combinatorics Keywords:
Cc: Fredrik Johansson Merged in:
Authors: Reviewers: Travis Scrimshaw
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: #10170 Stopgaps:

Status badges

Description

With #10170, for at least n >= 30, mpmath's bell(n) is incorrect:

sage: bell_number(30, 'mpmath') 
846749014511809388871680
sage: bell_number(30, 'gap')   
846749014511809332450147

The GAP implementation is the correct one (it agrees with OEIS and the implementation in #10170 for n > 200). I don't know if this is a me (i.e. this is fixed in an upgrade of mpmath not yet merged as of 5.8.beta1), sage, or an mpmath issue.

Attachments (1)

trac_14247-bell_number_mpmath-ts.patch (3.0 KB) - added by Travis Scrimshaw 10 years ago.

Download all attachments as: .zip

Change History (9)

comment:1 Changed 10 years ago by Karl-Dieter Crisman

Cc: Fredrik Johansson added

comment:2 Changed 10 years ago by Fredrik Johansson

It's not enough to call mpmath.bell(n) and convert to an integer since it gives a floating-point approximation; you need to set the precision to slightly more than log_2(B_n) bits.

Changed 10 years ago by Travis Scrimshaw

comment:3 Changed 10 years ago by Travis Scrimshaw

Authors: Travis Scrimshaw
Status: newneeds_review

Ah okay. Thank you.

I've uploaded a patch which gives bell_number() an optional argument to temporarily set mpmath's precision higher (thus bell_number() [with the algorithm set to 'mpmath'] effectively just wraps mpmath's bell).

comment:4 Changed 10 years ago by Fredrik Johansson

It doesn't make sense to return an incorrect value by default.

A better solution would be to set the precision to say 10 bits, compute p = mpmath.mag(mpmath.bell(n)) + 10, set mpmath.prec = p, then compute ret = ZZ(int(bell(n))).

comment:5 Changed 10 years ago by Travis Scrimshaw

Authors: Travis Scrimshaw
Milestone: sage-5.9sage-duplicate/invalid/wontfix

I've done what you suggested and folded this patch into #10170 since #10170 is not yet reviewed and this is basically part of #10170. I've also set this ticket as a duplicate.

comment:6 Changed 10 years ago by Nathann Cohen

Status: needs_reviewpositive_review

(so that Jeroen sees it and closes it)

Nathann

comment:7 Changed 10 years ago by Nathann Cohen

Could you also add this to the description of #10170 and check that the bug is fixed with a doctest there ?

Nathann

comment:8 Changed 10 years ago by Jeroen Demeyer

Resolution: duplicate
Reviewers: Travis Scrimshaw
Status: positive_reviewclosed
Note: See TracTickets for help on using tickets.