Opened 9 years ago

Closed 9 years ago

#9663 closed enhancement (fixed)

Fast computation of Stirling numbers of 2nd kind

Reported by: fredrik.johansson Owned by: sage-combinat
Priority: major Milestone: sage-4.6.1
Component: combinatorics Keywords:
Cc: was Merged in: sage-4.6.1.alpha0
Authors: Fredrik Johansson, Nathann Cohen Reviewers: Nathann Cohen, Nicolas Borie, Jeroen Demeyer
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Description

Currently, Stirling numbers are computed by calling GAP. The patch provides fast Cython code for Stirling numbers of the second kind, and allows using GAP or Maxima as an optional algorithm.

By having less overhead, the Cython code is about 10000 times faster than GAP or Maxima for tiny inputs, and it remains much faster than GAP for larger inputs as well. Apparently Maxima uses a fast algorithm unlike GAP, but my code is still about twice as fast as Maxima for huge n due to an algorithmic optimization.

sage: %timeit stirling_number2(10,5)
625 loops, best of 3: 2.33 µs per loop
sage: %timeit stirling_number2(10,5,algorithm='gap')
25 loops, best of 3: 20 ms per loop
sage: %timeit stirling_number2(10,5,algorithm='maxima')
5 loops, best of 3: 40 ms per loop

625 loops, best of 3: 16.2 µs per loop
sage: %timeit stirling_number2(100,50,algorithm='gap')
25 loops, best of 3: 20 ms per loop
sage: %timeit stirling_number2(100,50,algorithm='maxima')
5 loops, best of 3: 40 ms per loop

sage: %timeit stirling_number2(2000,1500)
25 loops, best of 3: 35.9 ms per loop
sage: %timeit stirling_number2(2000,1500,algorithm='gap')
5 loops, best of 3: 348 ms per loop
sage: %timeit stirling_number2(2000,1500,algorithm='maxima')
5 loops, best of 3: 210 ms per loop

sage: %timeit stirling_number2(4000,3000)
5 loops, best of 3: 249 ms per loop
sage: %timeit stirling_number2(4000,3000,algorithm='gap')
5 loops, best of 3: 2.96 s per loop
sage: %timeit stirling_number2(4000,3000,algorithm='maxima')
5 loops, best of 3: 948 ms per loop

sage: %time stirling_number2(20000,15000);
CPU times: user 20.30 s, sys: 0.23 s, total: 20.53 s
Wall time: 21.82 s
sage: %time stirling_number2(20000,15000,algorithm='maxima');
CPU times: user 0.00 s, sys: 0.01 s, total: 0.01 s
Wall time: 51.90 s

Mathematica seems to be about as slow as GAP (warning: timed on a different system):

In[1]:= Timing[StirlingS2[4000,3000];]
Out[1]= {27.1809, Null}

Attachments (3)

stirling2.patch (11.5 KB) - added by fredrik.johansson 9 years ago.
fast implementation of stirling_number2 -- updated patch
trac_9663 - additional test.patch (1.7 KB) - added by ncohen 9 years ago.
trac_9663-additional_tests.patch (2.0 KB) - added by ncohen 9 years ago.

Download all attachments as: .zip

Change History (20)

comment:1 Changed 9 years ago by fredrik.johansson

  • Status changed from new to needs_review

comment:2 Changed 9 years ago by was

Please explain the *massive* number of changes to module_list.py of the form:

153	 	                [[blank looking line]]
 	153	                [[another blank looking line]]

comment:3 Changed 9 years ago by fredrik.johansson

Oops, my editor was set to "strip trailing whitespace when saving files".

Changed 9 years ago by fredrik.johansson

fast implementation of stirling_number2 -- updated patch

comment:4 Changed 9 years ago by fredrik.johansson

  • Cc was added

Please see the new version of the patch.

comment:5 Changed 9 years ago by ncohen

  • Authors changed from fredrik.johansson to Fredrik Johansson
  • Milestone set to sage-4.5.3
  • Reviewers set to Nathann Cohen

Nice one !!

I learned many things while reviewing this patch :-)

Would you mind adding this small doctest in the patch I attach ? If not, you can set this ticket to "positive_review" !

Thanksssssssssssss !!!

Nathann

Changed 9 years ago by ncohen

comment:6 Changed 9 years ago by ncohen

Are you around ? There's basically nothing to do on this patch :-)

Nathann

comment:7 Changed 9 years ago by nborie

  • Reviewers changed from Nathann Cohen to Nathann Cohen, Nicolas Borie
  • Status changed from needs_review to positive_review

The two patch applied on 4.5.3 All tests pass, no warning in docbuild... Nice documentation, powerful implantation... Good job!

I give the two patch a positive review.

For the release manager, please apply :

  • stirling2.patch
  • trac_9663 - additional test.patch

comment:8 follow-up: Changed 9 years ago by jdemeyer

  • Authors changed from Fredrik Johansson to Fredrik Johansson, Nathann Cohen
  • Milestone changed from sage-4.6 to sage-4.6.1

I think you should add a patch with a test for "unknown algorithm".

comment:9 Changed 9 years ago by jdemeyer

  • Status changed from positive_review to needs_work

comment:10 Changed 9 years ago by jdemeyer

Also, please do not put spaces in patch filenames.

comment:11 in reply to: ↑ 8 ; follow-up: Changed 9 years ago by ncohen

Replying to jdemeyer:

I think you should add a patch with a test for "unknown algorithm".

What do you mean ?

Nathann

comment:12 in reply to: ↑ 11 Changed 9 years ago by jdemeyer

Replying to ncohen:

Replying to jdemeyer:

I think you should add a patch with a test for "unknown algorithm".

What do you mean ?

A test which does something like

sage: n = stirling_number2(20,11,algorithm='foobar')

to check the "unknown algorithm" code.

comment:13 follow-up: Changed 9 years ago by ncohen

  • Status changed from needs_work to needs_review

Here is a new version of my patch with the requested doctest.

Nathann

comment:14 in reply to: ↑ 13 Changed 9 years ago by jdemeyer

  • Reviewers changed from Nathann Cohen, Nicolas Borie to Nathann Cohen, Nicolas Borie, Jeroen Demeyer
  • Status changed from needs_review to needs_work

Replying to ncohen:

Here is a new version of my patch with the requested doctest.

Nathann

On line 670, TESTS:: should be TESTS: (the :: should precede a block of code, which is not the case here).

comment:15 Changed 9 years ago by ncohen

  • Status changed from needs_work to needs_review

Done.

Nathann

Changed 9 years ago by ncohen

comment:16 Changed 9 years ago by jdemeyer

  • Status changed from needs_review to positive_review

comment:17 Changed 9 years ago by jdemeyer

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