The Divisor Function Sigma

January 10, 2014

The divisor function first factors n, by any convenient method, then sorts the factors and applies the uniq-c function of the Standard Prelude, which collects the factors and their multiplicities; for instance, (uniq-c = (sort < (factors 60))) returns the list ((2 . 2) (3 . 1) (5 . 1)). Then the appropriate definition of the divisor function is applied, using map instead of explicit loops:

(define (divisor x n)
  (define (prod xs) (apply * xs))
  (if (= n 1) 1
    (let ((fs (uniq-c = (sort < (factors n)))))
      (if (zero? x)
          (prod (map add1 (map cdr fs)))
          (prod (map (lambda (p a)
                       (/ (- (expt p (* (+ a 1) x)) 1)
                          (- (expt p x) 1)))
                     (map car fs) (map cdr fs)))))))

The summatory function calculates the number of times each divisor appears and applies the appropriate powering operation:

(define (summatory-divisor x n)
  (do ((d 1 (+ d 1))
       (s 0 (+ s (* (expt d x) (quotient n d)))))
      ((< n d) s)))

It may help to see visually what is happening. Consider what happens when n = 10; the divisors of the numbers from 1 to 10 are:

1: 1
2: 1 2
3: 1 3
4: 1 2 4
5: 1 5
6: 1 2 3 6
7: 1 7
8: 1 2 4 8
9: 1 3 9
10: 1 2 5 10

If you rearrange the divisors, you can see how many times each divisor appears:

1:  (10) 1 2 3 4 5 6 7 8 9 10
2:  (5)  2 4 6 8 10
3:  (3)  3 6 9
4:  (2)  4 8
5:  (2)  5 10
6:  (1)  6
7:  (1)  7
8:  (1)  8
9:  (1)  9
10: (1) 10

Thus, each divisor d appears (quotient n d) times, which the function then raises to the appropriate power and sums over all divisors. Here are some examples:

> (divisor 0 60)
12
> (divisor 1 60)
168
> (divisor 2 60)
5460
> (summatory-divisor 0 60)
261
> (summatory-divisor 1 60)
3014
> (summatory-divisor 2 60)
89286

You can run the program at http://programmingpraxis.codepad.org/gpTcbMQr.

Today’s exercise was inspired by Project Euler 401, which you are now equipped to solve. But if you want to solve it within the one-minute rule, you will need to improve the summatory-divisor function, because even though it operates in linear time, doing a quadrillion of anything will take a while. Out of respect for Project Euler, I won’t give the code, but the idea is to combine adjacent divisors with identical counts and compute the sum of their squares using the formula for square pyramidal numbers.

Advertisement

Pages: 1 2

One Response to “The Divisor Function Sigma”

  1. Paul said

    In Python.

    from collections import Counter
    from ma.primee import td_factors
    
    def divisor_sigma(n, x):
        res = 1
        if x  == 0 :
            for pi, ai in Counter(td_factors(n)).iteritems():
                if ai:
                    res *= (ai + 1)
        else:
            for pi, ai in Counter(td_factors(n)).iteritems():
                if ai:
                    res *= (pi ** ((ai + 1) * x) - 1) // (pi ** x - 1)
        return res
        
    def sum_divisor_sigma(n, x):
        if x == 0:
            return sum(n // d for d in xrange(1, n + 1))
        else:
            return sum(n // d * d ** x for d in xrange(1, n + 1))
    

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: