Friday, May 1, 2009

python and divisor summing

There are lots of areas of number theory concerned with exact divisors of integers - perfect numbers, amicable numbers, sociable numbers and much more.

A good place to begin understanding the behaviour of the divisor sum functions
[ divsum() in maxima, sigma() in pari/gp ] is to consider the integers 1 to 100.

A simple python program to show you what happens for
those integers is given below.

#!/usr/bin/python
## Aliquot sequence demonstrator
## This program was written as a very simple demonstrator of what
## happens for numbers in the range 1 through 100 as, for each,
## you sum the 'proper' divisors and repeat.
## In order to make viewing on a terminal output practical the
## routines work in batches of 10 (or 20) at a time as
## required.
## A crude approach to knowing when to stop printing relies in
## part on a list of known terminators.
## A more general program would instead maybe look at how many
## repetitions of the same value had occurred instead. Try that
## once you have decided how many repetitions you think are
## appropriate.
## Further reading on aliquot sequences might be had by searching:
##  The Lehmer Five are the first unterminated aliquot sequences
##
## (c) 2009 Gary Wright
##
##This program is free software: you can redistribute it and/or modify
##it under the terms of the GNU General Public License as published by
##the Free Software Foundation, either version 3 of the License, or
##(at your option) any later version.
##
##This program is distributed in the hope that it will be useful,
##but WITHOUT ANY WARRANTY; without even the implied warranty of
##MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
##GNU General Public License for more details.
##
##You should have received a copy of the GNU General Public License
##along with this program.  If not, see http://www.gnu.org/licenses/
##
## http://www.wrightsolutions.co.uk/guide/code
## V1.0 20090501GW
##
#import divsum_out as ds
#ds.divsum_out20(range(20,40))# Readable when do divsum_out20(range(20)) ...
#...but not much higher
#ds.divsum_out10(0);divsum_out10(10)# More reabable and presentation suitable

known_terminators = ( 0, 6, 28 )     # tuple as immutable data

def divsum_out(n1):
print "%s%3i%s" % ('#',n1,'#')   # For first entries/line use
                                 # bang rather than pipe (|)
                                 # so easy to distinguish the first entries
                                 # which are just the numbers themselves.
#printf( ,'|',n1,'|')
while(n1 not in known_terminators):
    n1 = divsum(n1)
    print "%s%3i%s" % ('|',n1,'|')
return

def divsum(n):
# Sum all those values in the range 1..n-1 that are aliquot divisors of n
return sum(i for i in range(1, n-1) if n % i == 0)

def divsum_out20(nlist):
""" Run through divisor sums for 20 integers given as an input list
"""
delim = '#'
still_running_counter = 20
while(still_running_counter > 0):
    line = ''
    for i in nlist:
        line += (len(line) == 0)*delim
        line += str(i).rjust(3)+delim
    print line
    still_running_counter = sum(i for i in nlist if i not in known_terminators)
    delim = '|'
    nlist = map(divsum, nlist)
return

def divsum_out10(sval=0):
""" Run through divisor sums for 10 integers beginning with the supplied integer sval
sval is an abbreviation for Start value and defaults to 0.
You might try 10 or 20.
"""
delim = '#'
nlist = range(sval,sval+10)
still_running_counter = 10
while(still_running_counter > 0):
    line = ''
    for i in nlist:
        line += (len(line) == 0)*delim
        line += str(i).rjust(3)+delim
    print line
    still_running_counter = sum(i for i in nlist if i not in known_terminators)
    delim = '|'
    nlist = map(divsum, nlist)
return


I ran this as follows:
python -i divsum_out.py
>>> divsum_out10(0)

If you are new to python then you should refrain from typing the '>>>' which I included just to indicate (as python does) that you are then within an interactive python shell.

If you prefer to do your maths in maple then you might websearch for aliquotparts and aliquotsum functions.

Note: The web url and mail address given in the code are yet to be operational and will be usable when the server setup is complete in late summer 2009. I include them for future reference only.

1 comment:

  1. Being pretty new to using Python for mathematics back in 2009, the original script was a good first attempt with just one minor* mistake.

    Can you spot the mistake?

    I know I am opening myself up to a lot of finger poking by asking that, but if you have intermediate python knowledge then please comment, and I will name the reader who first spots the mistake I refer to.

    *I say minor mistake because it only affects a single result, I think, so the vast majority of what is output is correct.

    ReplyDelete