Opened 12 years ago

Closed 10 years ago

#10170 closed enhancement (fixed)

Speed up the computation of Bell numbers

Reported by: Robert Gerbicz Owned by: Sage Combinat CC user
Priority: major Milestone: sage-5.10
Component: combinatorics Keywords: bell number
Cc: Sage Combinat CC user Merged in: sage-5.10.beta1
Authors: Robert Gerbicz, Travis Scrimshaw Reviewers: Ben Salisbury
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Status badges

Description (last modified by Jeroen Demeyer)

For large n values use the formula: bell(n) = exp(-1)*sum(k=0, infinity, k**n/k!)

Currently my code beats all known implementations, also sage's code which wraps GAP's Bell is slow and using lots of memory, and unable to compute bell_number(n) if n is large.

Some timings on Amd 2.1 GHz:

n time (Wall time)

300 0.01 sec.

1000 0.03 sec.

3000 0.18 sec.

10000 2.88 sec.

30000 46.79 sec.

100000 819.37 sec.

Compare these with another programs: http://fredrik-j.blogspot.com/2009/03/computing-generalized-bell-numbers.html

Also exposes mpmath's bell() as an alternative algorithm.


Apply: trac_10170-bell_number_improvements-ts.patch

Attachments (3)

trac10170-bell_number_improvement.patch (1.5 KB) - added by Robert Gerbicz 12 years ago.
bell_proof.txt (1.5 KB) - added by Robert Gerbicz 12 years ago.
trac_10170-bell_number_improvements-ts.patch (7.7 KB) - added by Travis Scrimshaw 10 years ago.
Folded in mpmath patch

Download all attachments as: .zip

Change History (19)

Changed 12 years ago by Robert Gerbicz

Changed 12 years ago by Robert Gerbicz

Attachment: bell_proof.txt added

comment:1 Changed 12 years ago by Jason Bandlow

Status: newneeds_info

Hello and thanks for the patch! I'm curious, have you compared your solution with wrapping mpmath directly? For example,

sage: import mpmath
sage: mpmath.bell?

Also, this ticket is basically a dup of #7812, so one of these (maybe that one, since this has a patch) should be closed.

Finally, don't forget to set the status as 'needs_review' once you've added a patch!

Thanks, Jason

comment:2 in reply to:  1 ; Changed 12 years ago by Robert Gerbicz

Status: needs_infoneeds_review

Replying to jbandlow:

Hello and thanks for the patch! I'm curious, have you compared your solution with wrapping mpmath directly? For example,

sage: import mpmath
sage: mpmath.bell?

Yes, but times are not saved I will rerun them: here are the timings for mpmath:

n time (Wall time)

300 0.02 sec.

1000 0.06 sec.

3000 0.41 sec.

10000 10.26 sec.

30000 168.27 sec.

100000 (still computing)

So here mpmath is about 2-4 times slower, and seems that using more memory than my code.

comment:3 in reply to:  2 Changed 12 years ago by Robert Gerbicz

Replying to gerbicz:

100000 (still computing)

it took 2383 sec. for mpmath.

comment:4 Changed 12 years ago by martin

Status: needs_reviewneeds_work

Patch will not apply to sage 4.5.3

comment:5 Changed 10 years ago by Travis Scrimshaw

Authors: Robert GerbiczRobert Gerbicz, Travis Scrimshaw
Cc: Sage Combinat CC user added
Description: modified (diff)
Milestone: sage-5.9
Status: needs_workneeds_review

I've uploaded a new version of the patch which also adds mpmath (so #7812 can be closed as a duplicate after this patch is done). I've also come across an error in calling with the 'mpmath' agrument and started #14247.

For patchbot:

Apply: trac_10170-bell_number_improvements-ts.patch

comment:6 Changed 10 years ago by Travis Scrimshaw

I've folded in the changes from #14247 and marked that as a duplicate since there's no need for a separate ticket.

For patchbot:

Apply: trac_10170-bell_number_improvements-ts.patch

comment:7 Changed 10 years ago by Nathann Cohen

Hmmm... Looks nice, but shouldn't this be in that module ?

http://www.sagemath.org/doc/reference/combinat/sage/combinat/expnums.html

Nathann

comment:8 Changed 10 years ago by Travis Scrimshaw

Description: modified (diff)

From my 10 minutes of googling, no, it should not be. The exponential numbers seem to be a generalization of the Bell numbers and I don't know that this implementation can be generalized to that (and hence go into that module). My guess is that it's probably yes, but I'll leave that to someone with more expertise to answer.

comment:9 Changed 10 years ago by Nathann Cohen

At the line where you define E_1 and E_2 they appear to be perfectly equal. Is that intended ? O_o

Nathann

comment:10 Changed 10 years ago by Robert Gerbicz

That is not true, in fact: E_2=exp(-1)*E_1 as you can read in my original bell_proof.txt

About exponential numbers: that is a different problem. Here we compute only a single Bell number not the first n Bell numbers and there the algorithm is optimal (O(n*n) time), the same question for Bell numbers is an open problem.

Last edited 10 years ago by Robert Gerbicz (previous) (diff)

comment:11 in reply to:  10 Changed 10 years ago by Nathann Cohen

That is not true, in fact:

I have no idea, but this is what is written in the patch.

Nathann

comment:12 Changed 10 years ago by Robert Gerbicz

Yes, there is a missing parentheses in line 255 in the patch.

Changed 10 years ago by Travis Scrimshaw

Folded in mpmath patch

comment:13 Changed 10 years ago by Travis Scrimshaw

Fixed.

comment:14 Changed 10 years ago by Ben Salisbury

Milestone: sage-5.9sage-5.10
Reviewers: Ben Salisbury
Status: needs_reviewpositive_review

comment:15 Changed 10 years ago by Jeroen Demeyer

Description: modified (diff)

comment:16 Changed 10 years ago by Jeroen Demeyer

Merged in: sage-5.10.beta1
Resolution: fixed
Status: positive_reviewclosed
Note: See TracTickets for help on using tickets.