Programming Praxis – Divisors And Totatives

In today’s Programming Praxis exercise,our goal is to make a small library to compute the following of a number: the divisors, the sum of the divisors, the number of divisors, the totatives and the totient. Let’s get started, shall we?

Some imports:

import Data.List
import Data.Numbers.Primes
import Data.Ratio

Calculating the divisors of a number could be done the compact but naive way (divisors n = filter ((== 0) . mod n) [1..n]), but let’s go with the slightly longer but more efficient version of multiplying every combination of prime factors.

divisors :: Integral a => a -> [a]
divisors = nub . sort . map product . subsequences . primeFactors

Getting the sum of divisors is trivial.

divisorSum :: Integral a => a -> a
divisorSum = sum . divisors

For the amount of divisors we’ll ignore the trivial method (length . divisors) and go with the more efficient algorithm.

divisorCount :: Integral a => a -> Int
divisorCount = product . map (succ . length) . group . primeFactors

The totatives can be be trivially calculated using trial division.

totatives :: Integral a => a -> [a]
totatives n = filter ((== 1) . gcd n) [1..n]

For the totient there is also a better way than the trivial one (length . totatives):

totient :: (Integral a, Integral b) => a -> b
totient n = floor . product $ n % 1 : [1 - 1 % f | f <- nub $ primeFactors n]

As usual, some tests to see if everything is working properly:

main :: IO ()
main = do print $ divisors 60 == [1,2,3,4,5,6,10,12,15,20,30,60]
          print $ divisorSum 60 == 168
          print $ divisorCount 60 == 12
          print $ totatives 30 == [1,7,11,13,17,19,23,29]
          print $ totient 30 == 8
          print $ totient 60 == 16

Looks like it is, and at five lines of code I’d call this a pretty small library :) Fortunately, the non-trivial algorithms don’t add much in the way of length. Sure, the lines are a little longer, but they’re all still one-liners. I’d say it’s worth it for the improved speed.

About these ads

Tags: , , , , , , , ,

11 Responses to “Programming Praxis – Divisors And Totatives”

  1. programmingpraxis Says:

    I am gradually learning Haskell by reading your solutions to my exercises. Your totient function confuses me. nub $ primeFactors n gives the list of distinct factors of n. I assume % is an operator that performs division on fractions (it’s not in my list of Haskell operators), so the result of the list comprehension, for n=60, is the list 1/2, 2/3, 4/5. I understand why you cons n onto the list to take the product, but don’t understand why you take n % 1; is that some kind of type coercion so n is the same type as the rest of the list? Can’t Haskell work out the types so the coercion can be done automatically? Then I don’t understand the floor function. Assuming all the arithmetic is done as fractions, the calculation is guaranteed to produce an exact integer (all the denominators are factors of n, so they all cancel), so why do you need floor? And one last question: why use the Integral type instead of the Fractional type (since you are working with fractions)? If you used the Fractional type, could Haskell work out the type coercion of n % 1 automatically?

    Phil

  2. Remco Niemeijer Says:

    > I assume % is an operator that performs division on fractions (it’s not in my list of Haskell operators)

    Correct. % is found in Data.Ratio. It’s a way of doing mathematically exact division so you don’t have to worry about floating point inaccuracies.

    > I understand why you cons n onto the list to take the product, but don’t understand why you take n % 1; is that some kind of type coercion so n is the same type as the rest of the list?

    Yes. n has type Integral a => a, so I need to coerce it to Ratio a before I can cons it onto a [Ratio a].

    > Can’t Haskell work out the types so the coercion can be done automatically?

    Haskell can work out the types just fine. However, the only implicit coercion in Haskell is done on literals. You can do something like 2 : [1 % 2] and the 2 will automatically be converted to a Ratio. But n already has a type, so you will have to explicitly coerce it. The reasoning behind this is that implicit type coercion can lead to subtle bugs.

    > Then I don’t understand the floor function. Assuming all the arithmetic is done as fractions, the calculation is guaranteed to produce an exact integer (all the denominators are factors of n, so they all cancel), so why do you need floor?

    Technically, I don’t. Like you say, the result is mathematically exact. However, it makes more sense to return an Integer than to return a Ratio Integer, since the result will always be a whole number. floor is one way of doing so.

    > And one last question: why use the Integral type instead of the Fractional type (since you are working with fractions)? If you used the Fractional type, could Haskell work out the type coercion of n % 1 automatically?

    I’m using Integral for the arguments because the functions only work on whole numbers; there’s no such thing as the divisors of 2.4125. If you were to use Fractional then you wouldn’t need Ratios or coercion, but the function would not accept integers, which makes little sense. Besides, primeFactors returns a list of Integrals, so the conversion of n to a Ratio would be replaced by having to coerce the factors to Fractionals.

  3. programmingpraxis Says:

    That’s clear. Thank you.

    Scheme’s numeric tower (integer, rational, real, complex) is simpler than Haskell’s, which I find confusing. Do you gain any real advantage from the added complexity?

  4. Remco Niemeijer Says:

    One benefit of Haskell’s more finely grained typeclasses is that when defining a new datatype you can make it an instance of only the ones for which your datatype actually has an implementation rather than having to do a whole bunch of foo x = undefined. There are probably more reasons, since there’s a numeric-prelude package that splits everything up into dozens of type classes, but I haven’t looked at the issue in enough detail to give you any. Personally though, I do find myself wishing for a simpler stack on occasion, which would cut down on some of the necessary coercion. Then again, I’d probably find myself bemoaning the lack of granularity in other situations if that were the case.

  5. programmingpraxis Says:

    I wrote a set of one-liners equivalent to Remco’s version. You can run the program at http://programmingpraxis.codepad.org/poR5avz2. Follow the link to see all the prelude functions that are required; here is the meat of the library:

    (define (divisors n)
      (unique = (sort < (map product (subsets (factors n))))))
    
    (define (sumdiv n) (sum (divisors n)))
    
    (define (numdiv n)
      (apply * (map add1 (map cdr (uniq-c = (factors n))))))
    
    (define (totatives n)
      (filter (lambda (x) (= (gcd n x) 1)) (range 1 n)))
    
    (define (totient n)
      (product (cons n (map (lambda (x) (- 1 (/ x))) (unique = (factors n))))))
    
  6. programmingpraxis Says:

    Remco: please fix my indentation. Sorry.

  7. Remco Niemeijer Says:

    Fixed the indentation. A small tip: if you use square brackets instead of angle brackets for the code tag it will format your code as intended, even adding line numbers.

  8. Grégory LEOCADIE Says:

    I don’t understand why you use prime factors of a number to determine it divisors ?(some math properties I’m missing here)

  9. Remco Niemeijer Says:

    It’s an optimization. Instead of trying every number below n to see if it’s a divisor, which is rather slow, we use the following property: every divisor of n is either prime or can be obtained by multiplying two or more primes. These primes are necessarily also divisors of n. Therefore, by multiplying every possible combination of prime divisors of n, we obtain every divisor of n, possibly more than once. Once we have these products, we eliminate the duplicates and sort the list, producing the list of all divisors.

  10. Grégory LEOCADIE Says:

    Thank you for your explanation.

  11. Programming Praxis – Polite Numbers « Bonsai Code Says:

    [...] We need a function to calculate the divisors of a number, which we recycle from a previous exercise: [...]

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 )

Google+ photo

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

Connecting to %s


Follow

Get every new post delivered to your Inbox.

Join 35 other followers

%d bloggers like this: