The Miller-Rabin Primality Test in Clojure

I am new to clojure and I wanted to test out my skills by taking an existing implementation written in Python and writing it in Clojure.

Concerns

My main concerns are where I use cond and put an upto parameter in the function definitions of isWitness and probablyPrime. I am not sure if that's the best way to loop (I tried using the actual loop construct but found out, you couldn't break out of it).

Furthermore, on line 41 where I call the isWitness function from probablyPrime seems very messy. I am calling decompose twice when calling it once would be best (I don't know how to simulate local "variables" or use let properly).

Also I have separated a part of my isWitness into another function called checkWitness. I am not sure if that's the best way to do it.

(ns jeremy-kun.miller-rabin
(:require [clojure.math.numeric-tower :as math]))

;; This is the clojure version of the original python code done on jeremykun website
;; https://jeremykun.com/2013/06/16/miller-rabin-primality-test/
(defn pow [x y & [z]]
(if z
(mod (math/expt x y) z)
(math/expt x y)))

(defn rand-between [n]
(rand-nth (range 2 n)))

(defn decompose [exponent-of-two n]
(if (not= 0 (mod n 2))
{:exponent exponent-of-two :n n}
(decompose (+ 1 exponent-of-two) (/ n 2))))

(defn checkWitness [possibleWitness p exponent remainder upto]
(cond (> upto exponent) true
(= (- p 1) (pow possibleWitness 2 p)) false
:else (checkWitness (pow possibleWitness 2 p) p exponent remainder (inc upto)))
)

(defn isWitness [possibleWitness p exponent remainder upto]
(cond
(or (= (pow possibleWitness remainder p) 1) (= (pow possibleWitness remainder p) (- p 1))) false
(zero? exponent) true
:else (checkWitness (pow possibleWitness remainder p) p exponent remainder upto)))

(defn probablyPrime [p accuracy upto]
(if (or (= p 2) (= p 3))
true)

(if (< p 2)
false)

(cond
(> upto accuracy) true
(isWitness (rand-between (- p 2)) p (:exponent (decompose 0 (- p 1))) (:n (decompose 0 (- p 1))) 0) false
:else (probablyPrime p accuracy (inc upto))
))

;;(def test-witness [[10 5] [11 9] [12 5] [13 10] [14 11] [15 3] [16 13] [17 3] [18 7] [19 11]])

;;
;(prn
;  (map (fn [x] (isWitness (last x) (first x)
;                         (:exponent (decompose 0 (- (first x) 1))) (:n (decompose 0 (- (first x) 1)))
;                         0))
;       test-witness))

(prn (probablyPrime 25 100 0))                              ;; false
(prn (probablyPrime 100 100 0))                             ;; false
(prn (probablyPrime 151 100 0))                             ;; true
(prn (probablyPrime 97 100 0))                              ;; true

Test Cases

I have already included four test cases for primes. If you want to test more then just replace the first parameter and leave the rest as is, especially upto.

The block of code that is commented out tests the isWitness function. You can test it out by uncommenting the code and the value test-witness. The value test-witness is a list of lists. Each element in the list of list is the value to check and a random number between 2 and the value you want to test for primality minus 2.

For instance, [10 5], 10 is the value to check for primality and 5 is a random number between 2 and 10.

If you uncomment the block and run the code, it should return

(true false true false true true true false true false)

which is what the python code runs if you use the same values.

There is an error in the translation of the probablyPrime function. In Clojure, a function body is an implicit do form. The two if forms and the cond form:

(defn probablyPrime [p accuracy upto]
(if...
)

(if ...
)

(cond
...
))

... will be evaluated in succession. The values of the (if ...)s will be discarded; the result of the (cond ...) will be returned. The Python returns the associated value if either if condition is satisfied.

For example, the Clojure gives

(probablyPrime -2 0 1)
;true

... whereas the Python returns false.

I suggest you bring all the conditions under the cond:

(defn probablyPrime [p accuracy upto]
(cond
(or (= p 2) (= p 3)) true
(< p 2) false
(> upto accuracy) true
(isWitness (rand-between (- p 2)) p (:exponent (decompose 0 (- p 1))) (:n (decompose 0 (- p 1))) 0) false
:else (probablyPrime p accuracy (inc upto))))

There are other problems that make your Clojure unworkable for testing large numbers:

• You use proper recursion to express Python for loops. This may cause stack overflow and will certainly slow things down. Use Clojure's recur instead.
• Your rand-between function is $\mathcal{\Theta}(n)$. This is a performance disaster when $n$ is huge.

There are also some constructions that would improve your code in small ways:

• Use let forms to avoid recalculating expressions.
• Use inc and dec instead of (+ ... 1) and (- ... 1).
• Prefer explicit arities to rest parameters in function pow.

With these improvements, the rest of the code now looks as follows:

(ns miller-rabin.core
(:gen-class)
(:require [clojure.math.numeric-tower :as math]))

(defn pow
([x y] (math/expt x y))
([x y z] (mod (pow x y) z)))

(defn rand-between [n]
(+ (rand-int (- n 2)) 2))

(defn decompose [exponent-of-two n]
(if (odd? n)
{:exponent exponent-of-two :n n}
(recur (inc exponent-of-two) (quot n 2))))

(defn checkWitness [possibleWitness p exponent remainder upto]
(let [possibleWitness (pow possibleWitness 2 p)]
(or
(> upto exponent)
(and
(not= (dec p) possibleWitness)
(recur possibleWitness p exponent remainder (inc upto))))))

(defn isWitness [possibleWitness p exponent remainder upto]
(let [possibleWitness (pow possibleWitness remainder p)]
(cond
(or (= possibleWitness 1) (= possibleWitness (dec p))) false
(zero? exponent) true
:else (checkWitness possibleWitness p exponent remainder upto))))

I've noticed that, in translation, several functions have acquired an extra argument that is used internally to control repetition. I'll just restore decompose:

(defn decompose [n]
(loop [exponent-of-two 0, n n]
(if (odd? n)
{:exponent exponent-of-two :n n}
(recur (inc exponent-of-two) (quot n 2)))))

The others - isWitness and probablyPrime - are up to you. If you follow the Python, you don't need to factor out checkWitness.

To do:

• I have not seriously tested the above.
• Factoring out checkWitness from isWitness can be bettered using the sequence library.
• @ferada Thanks. I would prefer big theta to big O. Is this possible. Sep 4 '16 at 13:16
• Sure, it's just LaTeX markup, take a look at the source code now. Sep 4 '16 at 13:25
• @ferada Thanks. I've used MathJax in the Mathematics site - didn't know it was available here. Sep 4 '16 at 13:30
• Just a quick question for: [exponent-of-two 0, n n] does the exponent-of-two 0, set the value of exponent-of-two to 0? Sep 5 '16 at 13:31
• @JeelShah Yes. In the approved lingo, it binds the local identifier exponent-of-two to the value 0. And it binds thee local identifier n to the global value n, which is got from the argument. Sep 5 '16 at 22:07