Hi,
This is a follow on to yesterdays thread about computing
multiplicative partitions of integers...
Today, a co-worker and myself found a (non-obvious) fast algorithm for
doing this. Suprisingly, the tree based methods we began to discuss
yesterday turned out to be quite slow as you end up finding various
partitions multiple times. We haven't worked out how the algorithm
scales, but it seems very fast (as the tests below show) even though
it is pure python as of yet.
The algorithm is recursive - at each step you need to calculate
divisors of integers. The thing which makes the algorithm efficient
is that as it progresses we are able to strongly limit the number of
divisors needed. This is possible as it gets rid of the duplicates
seen in the tree based algorithms.
Question: where would we begin looking to see if this algorithm is
known already?
Thanks everyone for your thoughts on this one.
Here is the code (BSD licensed, authors, Brian Granger and Dan Karipedes):
def divisors_minmax(n, dmin, dmax):
"""Find the divisors of n in the interval (dmin,dmax]."""
i = dmin+1
while i<=dmax:
if n % i == 0:
yield i
i += 1
def list_or_tuple(seq):
return isinstance(seq, (list, tuple))
def flatten(seq, to_expand=list_or_tuple):
"""Flatten a nested sequence"""
for item in seq:
if to_expand(item):
for subitem in flatten(item, to_expand):
yield subitem
else:
yield item
def mult_partitions(n, s):
"""Compute the multiplicative partitions of n of size s"""
return [tuple(flatten(p)) for p in mult_partitions_recurs(n,s)]
def mult_partitions_recurs(n, s, pd=1):
if s == 1:
return [n]
divs = divisors_minmax(n, pd, int(sqrt(n)))
fs = []
for d in divs:
fs.extend([(d,f) for f in mult_partitions_recurs(n/d, s-1, pd)])
pd = d
return fs
A demo:
In [14]: mult_partitions(23452345,3)
Out[14]:
[(5, 7, 670067),
(5, 67, 70007),
(5, 73, 64253),
(5, 137, 34237),
(5, 469, 10001),
(5, 511, 9179),
(5, 959, 4891),
(7, 67, 50005),
(7, 73, 45895),
(7, 137, 24455),
(7, 335, 10001),
(7, 365, 9179),
(7, 685, 4891),
(35, 67, 10001),
(35, 73, 9179),
(35, 137, 4891),
(67, 73, 4795),
(67, 137, 2555),
(67, 365, 959),
(67, 511, 685),
(73, 137, 2345),
(73, 335, 959),
(73, 469, 685),
(137, 335, 511),
(137, 365, 469)]
And some timing:
In [10]: %timeit -n1 -r1 mult_partitions(23452345,2)
1 loops, best of 1: 2.5 ms per loop
In [11]: %timeit -n1 -r1 mult_partitions(23452345,3)
1 loops, best of 1: 5.55 ms per loop
In [12]: %timeit -n1 -r1 mult_partitions(23452345,4)
1 loops, best of 1: 6.79 ms per loop
In [13]: %timeit -n1 -r1 mult_partitions(23452345,5)
1 loops, best of 1: 6.54 ms per loop
--~--~---------~--~----~------------~-------~--~----~
To post to this group, send email to sage-devel@googlegroups.com
To unsubscribe from this group, send email to [EMAIL PROTECTED]
For more options, visit this group at http://groups.google.com/group/sage-devel
URLs: http://www.sagemath.org
-~----------~----~----~----~------~----~------~--~---
