Быстрая генерация простых чисел в Clojure

Question

Я работаю над решением Project Euler проблем в Clojure, чтобы стать лучше, и я уже несколько раз сталкивался с генерацией простых чисел. Моя проблема в том, что это занимает слишком много времени. Я надеялся, что кто-нибудь может помочь мне найти эффективный способ сделать это Clojure-y.

Когда я сделал это кулаком, я грубо заставил это. Это было легко сделать. Но вычисление 10001 простых чисел заняло 2 минуты на Xeon 2,33 ГГц, слишком долго для правил и вообще слишком долго. Вот алгоритм:

(defn next-prime-slow
    "Find the next prime number, checking against our already existing list"
    ([sofar guess]
        (if (not-any? #(zero? (mod guess %)) sofar)
            guess                         ; Then we have a prime
            (recur sofar (+ guess 2)))))  ; Try again                               

(defn find-primes-slow
    "Finds prime numbers, slowly"
    ([]
        (find-primes-slow 10001 [2 3]))   ; How many we need, initial prime seeds
    ([needed sofar]
        (if (<= needed (count sofar))
            sofar                         ; Found enough, we're done
            (recur needed (concat sofar [(next-prime-slow sofar (last sofar))])))))

Заменив next-prime-slow на более новую процедуру, которая учитывала некоторые дополнительные правила (например, свойство 6n +/- 1), я смог ускорить процесс примерно до 70 секунд.

Затем я попытался сделать сито из Эратосфена в чистом Clojure. Я не думаю, что вытащил все ошибки, но я сдался, потому что это было просто слишком медленно (я думаю, даже хуже, чем выше).

(defn clean-sieve
    "Clean the sieve of what we know isn't prime based"
    [seeds-left sieve]
    (if (zero? (count seeds-left))
        sieve              ; Nothing left to filter the list against
        (recur
            (rest seeds-left)    ; The numbers we haven't checked against
            (filter #(> (mod % (first seeds-left)) 0) sieve)))) ; Filter out multiples

(defn self-clean-sieve  ; This seems to be REALLY slow
    "Remove the stuff in the sieve that isn't prime based on it's self"
    ([sieve]
        (self-clean-sieve (rest sieve) (take 1 sieve)))
    ([sieve clean]
        (if (zero? (count sieve))
            clean
            (let [cleaned (filter #(> (mod % (last clean)) 0) sieve)]
                (recur (rest cleaned) (into clean [(first cleaned)]))))))

(defn find-primes
    "Finds prime numbers, hopefully faster"
    ([]
        (find-primes 10001 [2]))
    ([needed seeds]
        (if (>= (count seeds) needed)
            seeds        ; We have enough
            (recur       ; Recalculate
                needed
                (into
                    seeds    ; Stuff we've already found
                    (let [start (last seeds)
                            end-range (+ start 150000)]   ; NOTE HERE
                        (reverse                                                
                            (self-clean-sieve
                            (clean-sieve seeds (range (inc start) end-range))))))))))

Это плохо. Это также вызывает переполнение стека, если число 150000 меньше. Это несмотря на то, что я использую recur. Это может быть моя вина.

Далее я попробовал сито, используя методы Java на Java ArrayList. Это заняло немало времени и памяти.

Моя последняя попытка - использовать сито с использованием хеш-карты Clojure, вставляя все числа в сито, а затем разделяя числа, которые не являются простыми. В конце он берет список ключей, которые являются простыми числами, которые он нашел. Требуется приблизительно 10-12 секунд, чтобы найти 10000 простых чисел. Я не уверен, что он полностью отлажен. Это тоже рекурсивно (с использованием recur и loop), так как я пытаюсь быть Lispy.

Так что с такими проблемами проблема 10 (суммируя все простые числа до 2000000) убивает меня. Мой самый быстрый код нашел правильный ответ, но это заняло 105 секунд, и мне потребовалось немало памяти (я выделил ему 512 МБ, чтобы мне не пришлось с этим возиться). Другие мои алгоритмы занимают так много времени, что я всегда заканчивал тем, что останавливал их первым.

Я мог бы использовать сито для вычисления такого количества простых чисел в Java или C довольно быстро и без использования большого количества памяти. Я знаю, что в моем стиле Clojure / Lisp я что-то упускаю, что вызывает проблему.

Есть ли что-то, что я делаю действительно неправильно? Clojure просто медленно работает с большими последовательностями? Читая некоторые из обсуждений проекта Эйлера, люди вычислили первые 10000 простых чисел в других Лиспах менее чем за 100 миллисекунд. Я понимаю, что JVM может все замедлить, и Clojure относительно молод, но я не ожидал бы разницы в 100 раз.

Может ли кто-нибудь объяснить мне быстрый способ вычисления простых чисел в Clojure?

Answer 1

Вот еще один подход, который празднует Clojure's Java interop. Это занимает 374 мс на 2,4 ГГц Core 2 Duo (работает однопоточным). Я позволил эффективной Miller-Rabin реализации в Java BigInteger#isProbablePrime иметь дело с проверкой первичности.

(def certainty 5)

(defn prime? [n]
      (.isProbablePrime (BigInteger/valueOf n) certainty))

(concat [2] (take 10001 
   (filter prime? 
      (take-nth 2 
         (range 1 Integer/MAX_VALUE)))))

Вероятность Miller-Rabin 5, вероятно, не очень хороша для чисел, намного превышающих это. Эта уверенность равна 96.875%, что она проста (1 - .5^certainty)

Answer 2

Я понимаю, что это очень старый вопрос, но недавно я закончил тем, что искал то же самое, и ссылки здесь не были тем, что я искал (максимально ограничен функциональными типами, лениво генерируя ~ каждое ~ простое Я хочу).

Я наткнулся на красивую F # реализацию , так что все кредиты принадлежат ему. Я просто перенес это на Clojure:

(defn gen-primes "Generates an infinite, lazy sequence of prime numbers"
  []
  (letfn [(reinsert [table x prime]
            (update-in table [(+ prime x)] conj prime))
          (primes-step [table d]
                       (if-let [factors (get table d)]
                         (recur (reduce #(reinsert %1 d %2) (dissoc table d) factors)
                                (inc d))
                         (lazy-seq (cons d (primes-step (assoc table (* d d) (list d))
                                                        (inc d))))))]
    (primes-step {} 2)))

Использование просто

(take 5 (gen-primes))    

Answer 3

Очень поздно для вечеринки, но я приведу пример, используя Java BitSets:

(defn sieve [n]
  "Returns a BitSet with bits set for each prime up to n"
  (let [bs (new java.util.BitSet n)]
    (.flip bs 2 n)
    (doseq [i (range 4 n 2)] (.clear bs i))
    (doseq [p (range 3 (Math/sqrt n))]
      (if (.get bs p)
        (doseq [q (range (* p p) n (* 2 p))] (.clear bs q))))
    bs))

Запустив это на Macbook Pro 2014 года (2,3 ГГц Core i7), я получаю:

user=> (time (do (sieve 1e6) nil))
"Elapsed time: 64.936 msecs"

Answer 4

Смотрите последний пример здесь: http://clojuredocs.org/clojure_core/clojure.core/lazy-seq

;; An example combining lazy sequences with higher order functions
;; Generate prime numbers using Eratosthenes Sieve
;; See http://en.wikipedia.org/wiki/Sieve_of_Eratosthenes
;; Note that the starting set of sieved numbers should be
;; the set of integers starting with 2 i.e., (iterate inc 2) 
(defn sieve [s]
  (cons (first s)
        (lazy-seq (sieve (filter #(not= 0 (mod % (first s)))
                                 (rest s))))))

user=> (take 20 (sieve (iterate inc 2)))
(2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71)

Answer 5

Вот хорошая и простая реализация:

http://clj -me.blogspot.com / 2008/06 / primes.html

... но написано для какой-то версии Clojure до 1.0. Смотрите lazy_seqs в Clojure Contrib для того, чтобы работать с текущей версией языка.

Answer 6

Итак, я только что начал с Clojure, и да, это часто встречается в Project Euler, не так ли? Я написал довольно быстрый алгоритм пробного деления, но он не слишком масштабируется до того, как каждый цикл делений становится слишком медленным.

Итак, я снова начал, на этот раз используя метод сита:

(defn clense
  "Walks through the sieve and nils out multiples of step"
  [primes step i]
  (if (<= i (count primes))
    (recur 
      (assoc! primes i nil)
      step
      (+ i step))
    primes))

(defn sieve-step
  "Only works if i is >= 3"
  [primes i]
  (if (< i (count primes))
    (recur
      (if (nil? (primes i)) primes (clense primes (* 2 i) (* i i)))
      (+ 2 i))
    primes))

(defn prime-sieve
  "Returns a lazy list of all primes smaller than x"
  [x]
  (drop 2 
    (filter (complement nil?)
    (persistent! (sieve-step 
      (clense (transient (vec (range x))) 2 4) 3)))))

Использование и скорость:

user=> (time (do (prime-sieve 1E6) nil))
"Elapsed time: 930.881 msecs

Я очень доволен скоростью: у нее заканчивается REPL на MBP 2009 года. Это в основном быстро, потому что я полностью избегаю идиоматического Clojure и вместо этого зацикливаюсь, как обезьяна. Это также в 4 раза быстрее, потому что я использую переходный вектор для работы с ситом вместо того, чтобы оставаться полностью неизменным.

Редактировать: После нескольких предложений / исправлений ошибок от Will Ness теперь он работает намного быстрее.

Answer 7

(defn sieve
  [[p & rst]]
  ;; make sure the stack size is sufficiently large!
  (lazy-seq (cons p (sieve (remove #(= 0 (mod % p)) rst)))))

(def primes (sieve (iterate inc 2)))

при размере стека 10M я получаю 1001-е простое число за ~ 33 секунды на MacBook 2.1 Гц.

Answer 8

Вот простое сито в схеме:

http://telegraphics.com.au/svn/puzzles/trunk/programming-in-scheme/primes-up-to.scm

Вот пробег для простых чисел до 10 000:

#;1> (include "primes-up-to.scm")
; including primes-up-to.scm ...
#;2> ,t (primes-up-to 10000)
0.238s CPU time, 0.062s GC time (major), 180013 mutations, 130/4758 GCs (major/minor)
(2 3 5 7 11 13...

Answer 9

На основании комментария Уилла, вот мое мнение о postponed-primes:

(defn postponed-primes-recursive
  ([]
     (concat (list 2 3 5 7)
             (lazy-seq (postponed-primes-recursive
                        {}
                        3
                        9
                        (rest (rest (postponed-primes-recursive)))
                        9))))
  ([D p q ps c]
     (letfn [(add-composites
               [D x s]
               (loop [a x]
                 (if (contains? D a)
                   (recur (+ a s))
                   (persistent! (assoc! (transient D) a s)))))]
       (loop [D D
              p p
              q q
              ps ps
              c c]
         (if (not (contains? D c))
           (if (< c q)
             (cons c (lazy-seq (postponed-primes-recursive D p q ps (+ 2 c))))
             (recur (add-composites D
                                    (+ c (* 2 p))
                                    (* 2 p))
                    (first ps)
                    (* (first ps) (first ps))
                    (rest ps)
                    (+ c 2)))
           (let [s (get D c)]
             (recur (add-composites
                     (persistent! (dissoc! (transient D) c))
                     (+ c s)
                     s)
                    p
                    q
                    ps
                    (+ c 2))))))))

Первоначальное представление для сравнения:

Вот моя попытка перенести этот генератор простых чисел с Python на Clojure. Ниже возвращается бесконечная ленивая последовательность.

(defn primes
  []
  (letfn [(prime-help
            [foo bar]
            (loop [D foo
                   q bar]
              (if (nil? (get D q))
                (cons q (lazy-seq
                         (prime-help
                          (persistent! (assoc! (transient D) (* q q) (list q)))
                          (inc q))))
                (let [factors-of-q (get D q)
                      key-val (interleave
                               (map #(+ % q) factors-of-q)
                               (map #(cons % (get D (+ % q) (list)))
                                    factors-of-q))]
                  (recur (persistent!
                          (dissoc!
                           (apply assoc! (transient D) key-val)
                           q))
                         (inc q))))))]
    (prime-help {} 2)))

Использование:

user=> (first (primes))
2
user=> (second (primes))
3
user=> (nth (primes) 100)
547
user=> (take 5 (primes))
(2 3 5 7 11)
user=> (time (nth (primes) 10000))
"Elapsed time: 409.052221 msecs"
104743

редактирование:

Сравнение производительности, где postponed-primes использует очередь простых чисел, замеченных до сих пор, вместо рекурсивного вызова postponed-primes:

user=> (def counts (list 200000 400000 600000 800000))
#'user/counts
user=> (map #(time (nth (postponed-primes) %)) counts)
("Elapsed time: 1822.882 msecs"
 "Elapsed time: 3985.299 msecs"
 "Elapsed time: 6916.98 msecs"
 "Elapsed time: 8710.791 msecs"
2750161 5800139 8960467 12195263)
user=> (map #(time (nth (postponed-primes-recursive) %)) counts)
("Elapsed time: 1776.843 msecs"
 "Elapsed time: 3874.125 msecs"
 "Elapsed time: 6092.79 msecs"
 "Elapsed time: 8453.017 msecs"
2750161 5800139 8960467 12195263)

Answer 10

Если вам не нужно ленивое решение и вы просто хотите, чтобы последовательность простых чисел была ниже определенного предела, прямая реализация Sieve of Eratosthenes довольно быстрая. Вот моя версия с использованием переходных процессов:

(defn classic-sieve
  "Returns sequence of primes less than N"
  [n]
  (loop [nums (transient (vec (range n))) i 2]
    (cond
     (> (* i i) n) (remove nil? (nnext (persistent! nums)))
     (nums i) (recur (loop [nums nums j (* i i)]
                       (if (< j n)
                         (recur (assoc! nums j nil) (+ j i))
                         nums))
                     (inc i))
     :else (recur nums (inc i)))))