Opened 11 years ago

Closed 8 years ago

#10170 closed enhancement (fixed)

Speed up the computation of Bell numbers

Reported by: gerbicz Owned by: sage-combinat
Priority: major Milestone: sage-5.10
Component: combinatorics Keywords: bell number
Cc: sage-combinat 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 jdemeyer)

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 gerbicz 11 years ago.
bell_proof.txt (1.5 KB) - added by gerbicz 11 years ago.
trac_10170-bell_number_improvements-ts.patch (7.7 KB) - added by tscrim 8 years ago.
Folded in mpmath patch

Download all attachments as: .zip

Change History (19)

Changed 11 years ago by gerbicz

Changed 11 years ago by gerbicz

comment:1 follow-up: Changed 11 years ago by jbandlow

  • Status changed from new to needs_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 ; follow-up: Changed 11 years ago by gerbicz

  • Status changed from needs_info to needs_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 11 years ago by gerbicz

Replying to gerbicz:

100000 (still computing)

it took 2383 sec. for mpmath.

comment:4 Changed 11 years ago by allmar

  • Status changed from needs_review to needs_work

Patch will not apply to sage 4.5.3

comment:5 Changed 8 years ago by tscrim

  • Authors changed from Robert Gerbicz to Robert Gerbicz, Travis Scrimshaw
  • Cc sage-combinat added
  • Description modified (diff)
  • Milestone set to sage-5.9
  • Status changed from needs_work to needs_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 8 years ago by tscrim

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 8 years ago by ncohen

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 8 years ago by tscrim

  • 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 8 years ago by ncohen

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 follow-up: Changed 8 years ago by 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 8 years ago by gerbicz (previous) (diff)

comment:11 in reply to: ↑ 10 Changed 8 years ago by ncohen

That is not true, in fact:

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

Nathann

comment:12 Changed 8 years ago by gerbicz

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

Changed 8 years ago by tscrim

Folded in mpmath patch

comment:13 Changed 8 years ago by tscrim

Fixed.

comment:14 Changed 8 years ago by bsalisbury1

  • Milestone changed from sage-5.9 to sage-5.10
  • Reviewers set to Ben Salisbury
  • Status changed from needs_review to positive_review

comment:15 Changed 8 years ago by jdemeyer

  • Description modified (diff)

comment:16 Changed 8 years ago by jdemeyer

  • Merged in set to sage-5.10.beta1
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.