math - need divisors algorithm


P

Philp Smith

Hi

Does anyone have suggested code for a compact, efficient, elegant, most of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)
 
Ad

Advertisements

T

Tim Peters

[Philp Smith]
Does anyone have suggested code for a compact, efficient, elegant, most of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)

If the canonical factorization of N is the product of p_i**e_i, the
number of divisors is the product of e_i+1. This can be
astronomically large, so it's important to minimize the number of
multiplications performed, and to write a generator to produce the
divisors (it may, e.g., be impossible to materialize a list of all of
them).

This one's easy and efficient:

def gen_divisors(prime_list):
elts = sorted(set(prime_list))
numelts = len(elts)

# Split on the number of copies of elts.
def gen_inner(i):
if i >= numelts:
yield 1
return
thiselt = elts
thismax = prime_list.count(thiselt)

# powers <- [thiselt**i for i in range(thismax+1)],
# but using only thismax multiplications.
powers = [1]
for j in xrange(thismax):
powers.append(powers[-1] * thiselt)

for d in gen_inner(i+1):
for prime_power in powers:
yield prime_power * d

for d in gen_inner(0):
yield d

For example, do:

for d in gen_divisors([2, 3, 2, 3, 2, 3, 5, 2, 3]):
print d

You're on your own for the "proper" business; try not to use more than
one line to do it <wink>.

I find recursion clearest here. An iterative version is easy enough
to do by incrementing from 0 thru product(e_i) (and just skipping the
endpoints if you really want "proper") in a mixed-radix system defined
by the e_i. Exercise for the reader.

Harder: generate divisors in increasing order.
 
P

Philp Smith

excellent - thanks
in the meantime I managed to write something myself which
a) works
b) looks reasonably elegant and pythonic
but
c) is recursive
d) requires a mainroutine and co-routine

I have a feeling if I analysed it it would prove to be relatively
inefficient as I'm sure it creates a lot of temporary objects (partial
lists)

So I will gratefully try your version - and for interest's sake post mine at
some point (its on other computer) for comment.

Phil

Tim Peters said:
[Philp Smith]
Does anyone have suggested code for a compact, efficient, elegant, most
of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)

If the canonical factorization of N is the product of p_i**e_i, the
number of divisors is the product of e_i+1. This can be
astronomically large, so it's important to minimize the number of
multiplications performed, and to write a generator to produce the
divisors (it may, e.g., be impossible to materialize a list of all of
them).

This one's easy and efficient:

def gen_divisors(prime_list):
elts = sorted(set(prime_list))
numelts = len(elts)

# Split on the number of copies of elts.
def gen_inner(i):
if i >= numelts:
yield 1
return
thiselt = elts
thismax = prime_list.count(thiselt)

# powers <- [thiselt**i for i in range(thismax+1)],
# but using only thismax multiplications.
powers = [1]
for j in xrange(thismax):
powers.append(powers[-1] * thiselt)

for d in gen_inner(i+1):
for prime_power in powers:
yield prime_power * d

for d in gen_inner(0):
yield d

For example, do:

for d in gen_divisors([2, 3, 2, 3, 2, 3, 5, 2, 3]):
print d

You're on your own for the "proper" business; try not to use more than
one line to do it <wink>.

I find recursion clearest here. An iterative version is easy enough
to do by incrementing from 0 thru product(e_i) (and just skipping the
endpoints if you really want "proper") in a mixed-radix system defined
by the e_i. Exercise for the reader.

Harder: generate divisors in increasing order.
 
E

Ed Suominen

Philp said:
Hi

Does anyone have suggested code for a compact, efficient, elegant, most of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)

Is this compact enough? :)

def properDivisors(N):
return [x for x in range(1,N-2) if not divmod(N,x)[1]]
 
J

Joal Heagney

Ed said:
Philp Smith wrote:

Hi

Does anyone have suggested code for a compact, efficient, elegant, most of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)


Is this compact enough? :)

def properDivisors(N):
return [x for x in range(1,N-2) if not divmod(N,x)[1]]

---
Ed Suominen
Registered Patent Agent
Open-Source Software Author (yes, both...)
Web Site: http://www.eepatents.com

I tried this out, (Python 2.3) and ran into the situation where neither
xrange or range won't accept anything other than an int. (To quote the
"bards", "53487861112345 is RIGHT out.") I was wondering how others
would deal with this situation?

Joal
 
D

David Eppstein

Does anyone have suggested code for a compact, efficient, elegant, most of
I have code for this in
<http://www.ics.uci.edu/~eppstein/numth/egypt/egypt.py>
(look for the function called, surprisingly, "divisors").

It's not very compact -- 155 lines from the start of the Pollard Rho
factorization (thanks to Kirby Urner) to the end of the divisors
function, but it's reasonably efficient even for largish numbers.
 
Ad

Advertisements

P

Peter Schorn

Philp said:
Hi

Does anyone have suggested code for a compact, efficient, elegant, most of
all pythonic routine to produce a list of all the proper divisors of an
integer (given a list of prime factors/powers)

What about

# Returns a list of all divisors of n = a1^e1*a2^e2* ... *an^en where
# input parameter l = [(a1, e1), (a2, e2), ..., (an, en)]
def divisors(l):
if l: return [i*j for i in [l[0][0]**k for k in range(l[0][1] + 1)]
for j in divisors(l[1:])]
else: return [1]


# divisors([(2,3),(3,2)]) == [1, 3, 9, 2, 6, 18, 4, 12, 36, 8, 24, 72]


Regards, Peter
 

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments. After that, you can post your question and our members will help you out.

Ask a Question

Top