首页 文章

Clojure中快速素数生成

提问于
浏览 次
47

我一直在努力解决Clojure中的问题，以便变得更好，而且我已经遇到了几次素数 . 我的问题是它只是花了太长时间 . 我希望有人可以帮我找到一种以Clojure-y方式做到这一点的有效方法 .

当我拳头做到这一点时，我粗暴地强迫它 . 这很容易做到 . 但是计算10001个素数在Xeon 2.33GHz上用了2分钟，对规则来说太长了，一般来说太长了 . 这是算法：

(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中筛选出Eratosthenes . 我不认为我得到了所有的错误，但我放弃了，因为它太慢了（我认为甚至比上面更糟糕） .

(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较小，它还会导致堆栈溢出 . 尽管事实上我正在使用复发 . 那可能是我的错 .

接下来，我尝试使用Java ArrayList上的Java方法筛选 . 这需要相当多的时间和记忆 .

我最近的尝试是使用Clojure哈希映射的筛子，在筛子中插入所有数字然后分解不是素数的数字 . 最后，它采用密钥列表，它是它找到的素数 . 找到10000个素数需要大约10-12秒 . 我不确定它是否已经完全调试过了 . 它也是递归的（使用recur和loop），因为我试图成为Lispy .

因此，对于这些问题，问题10（总计2000000以下的所有素数）正在扼杀我 . 我最快的代码提出了正确的答案，但它需要105秒才能完成，并需要相当多的内存（我给它512 MB只是所以我不必大惊小怪） . 我的其他算法花了这么长时间我总是先停止它们 .

我可以使用筛子来快速计算Java或C中的许多质数，而不需要使用太多的内存 . 我知道我必须在我的Clojure / Lisp风格中遗漏导致问题的东西 .

有什么我做错了吗？ Clojure对大序列有点慢吗？阅读一些项目的欧拉讨论，人们已经在不到100毫秒的时间内计算了其他Lisps中的前10000个素数 . 我意识到JVM可能会减慢速度并且Clojure相对年轻，但我不希望它有100倍的差异 .

有人可以通过快速的方式启发我来计算Clojure中的素数吗？

clojure lisp primes

14 回答

  • 1

    这是另一种庆祝 Clojure's Java interop 的方法 . 这需要在2.4 Ghz Core 2 Duo上运行374ms（运行单线程） . 我让Java的 BigInteger#isProbablePrime 中的高效 Miller-Rabin 实现处理primality检查 .

    (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)))))

    对于比这大得多的数字，5的确定性可能不是很好 . 这个确定性等于 96.875% 确定它是黄金时期（ 1 - .5^certainty ）

    回复于
  • 0

    我意识到这是一个非常古老的问题，但我最近最终寻找相同的东西，这里的链接不是我正在寻找的东西（尽可能限制功能类型，懒洋洋地生成〜每个我想要的素数） .

    我偶然发现了一个很好的F# implementation，所以所有的积分都是他的 . 我只把它移植到Clojure：

    (defn gen-primes "Generates an infinite, lazy sequence of prime numbers"
  []
  (let [reinsert (fn [table x prime]
                   (update-in table [(+ prime x)] conj prime))]
    (defn 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)))

    用法很简单

    (->>
  (gen-primes)
  (filter #(< % 2000000)
  ; do what you need with the stuff
)
    回复于
  • 11

    聚会很晚，但我会举一个例子，使用Java BitSet：

    (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))

    在2014 Macbook Pro（2.3GHz Core i7）上运行，我得到：

    user=> (time (do (sieve 1e6) nil))
"Elapsed time: 64.936 msecs"
    回复于
  • 3

    请参见此处的最后一个示例：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)
    回复于
  • 2

    这是一个很好而简单的实现：

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

    ...但它是为一些1.0之前版本的Clojure编写的 . 请参阅Clojure Contrib中的lazy_seqs，了解与当前语言版本一起使用的版本 .

    回复于
  • 0
    (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)))

    在10G堆栈大小的情况下，我在2.1Gz的macbook上在~33秒内获得第1001个素数 .

    回复于
  • 0

    所以我刚刚开始使用Clojure，是的，这在项目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

    我对速度非常满意：它已经耗尽了2009 MBP上的REPL . 它主要是快速的，因为我完全避开惯用的Clojure，而是像猴子一样环绕 . 它也快了4倍，因为我使用瞬态矢量在筛子上工作而不是完全停留不可改变的 .

    Edit: 经过Will Ness的一些建议/错误修复后，它现在运行得更快 .

    回复于
  • 29

    这是Scheme中的简单筛子：

    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...
    回复于
  • 0

    根据Will的评论，这是我对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))))))))

    初次提交比较：

    这是我尝试将this prime number generator从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)
    回复于
  • 0

    来自：http://steloflute.tistory.com/entry/Clojure-%ED%94%84%EB%A1%9C%EA%B7%B8%EB%9E%A8-%EC%B5%9C%EC%A0%81%ED%99%94

    使用Java数组

    (defmacro loopwhile [init-symbol init whilep step & body]
  `(loop [~init-symbol ~init]
     (when ~whilep ~@body (recur (+ ~init-symbol ~step)))))

(defn primesUnderb [limit]
  (let [p (boolean-array limit true)]
    (loopwhile i 2 (< i (Math/sqrt limit)) 1
               (when (aget p i)
                 (loopwhile j (* i 2) (< j limit) i (aset p j false))))
    (filter #(aget p %) (range 2 limit))))

    用法和速度：

    user=> (time (def p (primesUnderb 1e6)))
"Elapsed time: 104.065891 msecs"
    回复于
  • 21

    我刚开始使用Clojure，所以我不知道它是否好，但这是我的解决方案：

    (defn divides? [x i]
  (zero? (mod x i)))

(defn factors [x]
    (flatten (map #(list % (/ x %)) (filter #(divides? x %) (range 1 (inc (Math/floor (Math/sqrt x))))))))

(defn prime? [x]
  (empty? (filter #(and divides? (not= x %) (not= 1 %)) (factors x))))

(def primes 
  (filter prime? (range 2 java.lang.Integer/MAX_VALUE)))

(defn sum-of-primes-below [n]
  (reduce + (take-while #(< % n) primes)))
    回复于
  • 3

    来到这个线程并寻找更快的替代已经在这里的那些，我很惊讶没有人链接到以下article by Christophe Grand：

    (defn primes3 [max]
  (let [enqueue (fn [sieve n factor]
                  (let [m (+ n (+ factor factor))]
                    (if (sieve m)
                      (recur sieve m factor)
                      (assoc sieve m factor))))
        next-sieve (fn [sieve candidate]
                     (if-let [factor (sieve candidate)]
                       (-> sieve
                         (dissoc candidate)
                         (enqueue candidate factor))
                       (enqueue sieve candidate candidate)))]
    (cons 2 (vals (reduce next-sieve {} (range 3 max 2))))))

    除了懒惰版本：

    (defn lazy-primes3 []
  (letfn [(enqueue [sieve n step]
            (let [m (+ n step)]
              (if (sieve m)
                (recur sieve m step)
                (assoc sieve m step))))
          (next-sieve [sieve candidate]
            (if-let [step (sieve candidate)]
              (-> sieve
                (dissoc candidate)
                (enqueue candidate step))
              (enqueue sieve candidate (+ candidate candidate))))
          (next-primes [sieve candidate]
            (if (sieve candidate)
              (recur (next-sieve sieve candidate) (+ candidate 2))
              (cons candidate 
                (lazy-seq (next-primes (next-sieve sieve candidate) 
                            (+ candidate 2))))))]
    (cons 2 (lazy-seq (next-primes {} 3)))))
    回复于
  • 10

    惯用语，而不是太糟糕

    (def primes
  (cons 1 (lazy-seq
            (filter (fn [i]
                      (not-any? (fn [p] (zero? (rem i p)))
                                (take-while #(<= % (Math/sqrt i))
                                            (rest primes))))
                    (drop 2 (range))))))
=> #'user/primes
(first (time (drop 10000 primes)))
"Elapsed time: 0.023135 msecs"
=> 104729
    回复于
  • 4

    已经有很多答案，但我有一个替代解决方案，可以生成无限的素数序列 . 我也对一些解决方案的标记感兴趣 .

    首先是一些Java互操作 . 以供参考：

    (defn prime-fn-1 [accuracy]
  (cons 2
    (for [i (range)
          :let [prime-candidate (-> i (* 2) (+ 3))]
          :when (.isProbablePrime (BigInteger/valueOf prime-candidate) accuracy)]
      prime-candidate)))

    本杰明@ https://stackoverflow.com/a/7625207/3731823是 primes-fn-2

    nha @ https://stackoverflow.com/a/36432061/3731823是 primes-fn-3

    我的实现是 primes-fn-4 ：

    (defn primes-fn-4 []
  (let [primes-with-duplicates
         (->> (for [i (range)] (-> i (* 2) (+ 5))) ; 5, 7, 9, 11, ...
              (reductions
                (fn [known-primes candidate]
                  (if (->> known-primes
                           (take-while #(<= (* % %) candidate))
                           (not-any?   #(-> candidate (mod %) zero?)))
                   (conj known-primes candidate)
                   known-primes))
                [3])     ; Our initial list of known odd primes
              (cons [2]) ; Put in the non-odd one
              (map (comp first rseq)))] ; O(1) lookup of the last element of the vec "known-primes"

    ; Ugh, ugly de-duplication :(
    (->> (map #(when (not= % %2) %) primes-with-duplicates (rest primes-with-duplicates))
         (remove nil?))))

    报告的数字（计算前N个素数的时间以毫秒计）是5次运行中最快的，实验之间没有JVM重启，因此您的里程可能会有所不同：

    1e6      3e6

(primes-fn-1  5)     808     2664
(primes-fn-1 10)     952     3198
(primes-fn-1 20)    1440     4742
(primes-fn-1 30)    1881     6030
(primes-fn-2)       1868     5922
(primes-fn-3)        489     1755  <-- WOW!
(primes-fn-4)       2024     8185
    回复于

相关问题