Instantly share code, notes, and snippets.

# ericnormand/00 Prime factorization.md

Last active December 27, 2019 00:43
Star You must be signed in to star a gist

Prime factorization

Every natural number can be expressed as a product of primes. For instance, `10 = 2 x 5` and `24 = 2 x 2 x 2 x 3`.

Your task is to write a function that takes a number `n` as an argument and returns a list of its prime factors.

This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 (defn factor ([n] (factor 2 [] n)) ([p ps n] (cond (= 1 n) ps (zero? (rem n p)) (recur p (conj ps p) (quot n p)) :else (recur (inc p) ps n))))
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 ;; My first cut - this seemed to give correct answers, but ;; when I gave it the 50millionth prime to factor, it took about 40 seconds. (defn factors [n] (let [multiple-of? (fn [n div] (zero? (mod n div))) divisors (lazy-seq (cons 2 (iterate (partial + 2) 3)))] (loop [res [] n n divs divisors] (let [d (first divs)] (if (> d n) res (if (multiple-of? n d) (recur (conj res d) (quot n d) divs) (recur res n (rest divs)))))))) ;; First attempt to speed it up - I figured replacing the lazy sequence of divisors ;; with a "next-div" fn would be faster. It was - now did 50millionth prime in 12 secs. (defn factors [n] (let [multiple-of? (fn [n div] (zero? (mod n div))) next-div (fn [d] (if (= 2 d) 3 (+ d 2)))] (loop [res [] n n d 2] (if (> d n) res (if (multiple-of? n d) (recur (conj res d) (quot n d) d) (recur res n (next-div d))))))) ;; Next attempt to speed up further. I figured the 'next-div' logic was simple once ;; we got to odd divisors, so I pulled-out the "odd-ball case" of dividing by 2 and this ;; was an additional win. Factors 50millionth prime in 10 seconds. (defn remove-twos [n] (loop [n n twos []] (if (odd? n) [n twos] (recur (quot n 2) (conj twos 2))))) (defn factors [n] (let [multiple-of? (fn [n div] (zero? (mod n div))) [odd-n twos] (remove-twos n)] (loop [res twos n odd-n d 3] (if (> d n) res (if (multiple-of? n d) (recur (conj res d) (quot n d) d) (recur res n (+ d 2)))))))
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 (defn primes-up-to "Given a positive integer N, returns a vector of the prime numbers up to n." [N] ;; Implementation using an Eratosthenes Sieve (if (< N 2) [] (let [N (int N) sieve (boolean-array (inc N) false)] (loop [n (int 2)] (cond (aget sieve n) ;; not a prime, next (recur (unchecked-inc-int n)) (< N (* n n)) ;; we're done, package the result (into [] (remove #(aget sieve %)) (range 2 (inc N))) :else ;; n is a prime, marking its multiples (do (loop [m (unchecked-add-int n n)] (when (<= m N) (aset sieve m true) (recur (unchecked-add-int m n)))) (recur (unchecked-inc-int n)))))))) (defn divides? [d n] (-> n (rem d) (= 0))) (defn reduce-by-prime "Trims out p in the prime factorization of m, i.e returns the highest divisor of n not divisible by p." [p n] (if (divides? p n) (recur p (/ n p)) n)) (comment (reduce-by-prime 5 100) => 4 (reduce-by-prime 7 100) => 100) (defn prime-factors [n] (if (< n 2) (sorted-set) (into (sorted-set) (let [s (int (Math/ceil (Math/sqrt n))) potential-prime-divisors (primes-up-to s)] (loop [m n ppds potential-prime-divisors prime-divisors []] (if (or (= m 1) (empty? ppds)) (if (empty? prime-divisors) #{n} prime-divisors) (let [pd (first ppds) m1 (reduce-by-prime pd m)] (recur m1 (rest ppds) (cond-> prime-divisors (< m1 m) (conj pd)))))))))) (comment (prime-factors 31) (prime-factors 120) (prime-factors 4) (prime-factors 2) (let [ ;; let's try this on a big random number n (* (rand-int Integer/MAX_VALUE) (rand-int 1000000))] (time ;; takes about 1s on my machine, most of it is spent computing the sieve. [n (prime-factors n)])))
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 (ns prime-factors) (require 'clojure.math.numeric-tower) ; start at 2, then run through every odd number (defn of-1 [n] (loop [m n div 2 acc []] (cond (zero? (mod m div)) (recur (/ m div) div (conj acc div)) (< m 2) acc :else (recur m (if (odd? div) (+ 2 div) (inc div)) acc)))) ; start at 2, then run through the odd numbers that can divide the ; given number evenly (defn of-2 [n] (loop [running-n n divisor 2 factors []] (cond (zero? (mod running-n divisor)) (recur (/ running-n divisor) divisor (conj factors divisor)) (< running-n 2) factors :else (recur running-n (if (odd? divisor) (loop [running-divisor divisor] (cond (zero? (rem running-n running-divisor)) running-divisor :else (+ 2 running-divisor))) (inc divisor)) factors)))) ; run through 2, 3, and every 6n-1, 6n+1 that divides n evenly (defn of [n] (loop [running-n n divisors (conj (filter (fn [x] (or (= 1 (rem x 6)) (= 5 (rem x 6)))) (range 5 (inc n))) 3 2) factors []] (cond (zero? (rem running-n (first divisors))) (recur (/ running-n (first divisors)) divisors (conj factors (first divisors))) (< running-n 2) factors :else (recur running-n (drop-while (fn [x] (pos? (rem running-n x))) divisors) factors))))