循环函数花费的时间太长

Question

我试图做一个函数来实现n个立方体之和

1^3 + 2^3 + 3^3 +.+ n^3 =和

我的函数应该接收sum，如果不存在n或-1，则返回n或-1。

下面是一些例子：

(find-n 9)   ; should return 2 because 1^3 + 2^3 = 9
(find-n 100) ; should return 4 because 1^3 + 2^3 + 3^3 + 4^3 = 100
(find-n 10)  ; should return -1

在做了一些工作之后，我完成了以下两个功能：

; aux function
(defn exp-3 [base] (apply *' (take 3 (repeat base))))

; main function
(defn find-n [m]
  (loop [sum 0
         actual-base 0]
       (if (= sum m) 
           actual-base
           (if (> sum m)
               -1
               (recur (+' sum (exp-3 (inc actual-base))) (inc actual-base))))))

这些函数正常工作，但是使用BigNumbers评估操作所花费的时间太长，例如：

(def sum 1025247423603083074023000250000N)
(time (find-n sum))
; => "Elapsed time: 42655.138544 msecs"
; => 45001000

我问这个问题是想提出一些建议，如何使这个函数更快。

Answer 1

这都是关于代数的，与Clojure或编程几乎没有关系。因为这个站点不支持数学排版，所以让我们用Clojure来表达它。

定义

(defn sigma [coll] (reduce + coll))

和

(defn sigma-1-to-n [f n]
  (sigma (map f (rest (range (inc n))))))

(或

(defn sigma-1-to-n [f n]
  (->> n inc range rest (map f) sigma))

)

那么问题是，给定n，找到i这样的(= (sigma-1-to-n #(* % % %) i) n)。

快速完成此操作的关键是多维数据集的Faulhaber公式。它告诉我们，对于任何自然数i，以下内容都是相等的：

(#(*' % %) (sigma-1-to-n identity i))

(sigma-1-to-n #(* % % %) i)

(#(*' % %) (/ (*' i (inc i)) 2))

所以，作为立方体的和，数字

• 一定是一个完美的正方形
• 它的平方根是第一个这么多数字的和。

为了确定一个整数是否是一个完美的平方，我们取其近似浮点平方根，并查看最近整数的平方是否恢复了我们的整数：

(defn perfect-square-root [n]
  (let [candidate (-> n double Math/sqrt Math/round)]
    (when (= (*' candidate candidate) n)
      candidate)))

如果参数不是一个完美的平方，则返回nil。

现在我们有了平方根，我们必须确定它是否是自然数范围的总和:在普通代数中，它是(j (j + 1)) / 2吗，对于某些自然数j。

我们可以用类似的方法直接回答这个问题。

j (j + 1) = (j + 1/2)^2 + 1/4

因此，如果有参数，下面的函数返回与参数相加的连续数数：

(defn perfect-sum-of [n]
  (let [j (-> n (*' 2)
                (- 1/4)
                double
                Math/sqrt
                (- 0.5)
                Math/round)]
    (when (= (/ (*' j (inc j)) 2) n)
      j)))

我们可以结合这些来做你想做的事：

(defn find-n [big-i]
  {:pre [(integer? big-i) ((complement neg?) big-i)]}
  (let [sqrt (perfect-square-root big-i)]
    (and sqrt (perfect-sum-of sqrt))))

(def sum 1025247423603083074023000250000N)

(time (find-n sum))
"Elapsed time: 0.043095 msecs"
=> 45001000

(请注意，时间比以前快了大约20倍，可能是因为HotSpot必须在find-n上工作，附加的测试已经对其进行了彻底的执行)

这显然比原来的快多了。

警告

我担心由于浮点数的有限精度，上述程序可能会产生假阴性(永远不会产生假阳性)。然而，测试表明，对于该问题所使用的数字类型，该过程是不可打破的。

一个Java double有52位精度，大约是小数点15.6位。需要注意的是，如果数字比此大得多，则该过程可能会忽略确切的整数解，因为舍入只能与它开始时的浮点数一样精确。

然而，这个过程正确地解决了一个31位整数的例子。和许多测试(一千万！)类似的数字不会导致一次失败。

为了测试解决方案，我们生成了一个(懒惰)序列的[limit cube-sum]对：

(defn generator [limit cube-sum]
  (iterate
    (fn [[l cs]]
      (let [l (inc l)
            cs (+' cs (*' l l l))]
        [limit cs]))
    [limit cube-sum]))

例如,

(take 10 (generator 0 0))
=> ([0 0] [1 1] [2 9] [3 36] [4 100] [5 225] [6 441] [7 784] [8 1296] [9 2025])

现在我们

• 从给定的示例开始，
• 试试接下来的一千万箱
• 移除那些有效的。

所以

(remove (fn [[l cs]] (= (find-n cs) l)) (take 10000000 (generator 45001000 1025247423603083074023000250000N)))
=> () 

他们都在工作。没有失败。为了确保我们的测试是有效的：

(remove (fn [[l cs]] (= (find-n cs) l)) (take 10 (generator 45001001 1025247423603083074023000250000N)))
=>
([45001001 1025247423603083074023000250000N]
 [45001002 1025247514734170359564546262008N]
 [45001003 1025247605865263720376770289035N]
 [45001004 1025247696996363156459942337099N]
 [45001005 1025247788127468667814332412224N]
 [45001006 1025247879258580254440210520440N]
 [45001007 1025247970389697916337846667783N]
 [45001008 1025248061520821653507510860295N]
 [45001009 1025248152651951465949473104024N]
 [45001010 1025248243783087353664003405024N])

一切都应该失败，而他们确实失败了。

Answer 2

只要避免使用apply (在CLJ中不是很快)，就可以加速4倍：

(defn exp-3 [base]
  (*' base base base))

另有10%：

(defn find-n [m]
  (loop [sum 0
         actual-base 0]
    (if (>= sum m)
      (if (= sum m) actual-base -1)
      (let [nb (inc actual-base)]
        (recur (+' sum (*' nb nb nb)) nb)))))

Answer 3

以下基于算法的方法依赖于一个简单的公式，即第一个N个自然数的立方体之和是：(N*(N+1)/2)^2

(defn sum-of-cube
  "(n*(n+1)/2)^2"
  [n]
  (let [n' (/ (*' n (inc n)) 2)]
    (*' n' n')))

(defn find-nth-cube
  [n]
  ((fn [start end prev]
     (let [avg (bigint (/ (+' start end) 2))
           cube (sum-of-cube avg)]
       (cond (== cube n) avg
             (== cube prev) -1
             (> cube n) (recur start avg cube)
             (< cube n) (recur avg end cube))))
    1 n -1))

(time (find-nth-cube 1025247423603083074023000250000N))
"Elapsed time: 0.355177 msecs"
=> 45001000N

我们想要找到数N，使得1.N个立方体之和是某个数字X。如果存在这样的数字，我们可以在某个范围内对它执行二进制搜索，方法是应用上述公式查看公式的结果是否等于X。这种方法有效，因为顶部的函数在增加，因此太大的值f(n)意味着我们必须寻找较低的数字n，而任何值f(n)太小意味着我们必须寻找更大的数字n。

我们选择一个0到X的范围(大于必要，但简单而安全)。如果应用于给定候选数的公式产生X，我们就会知道这个数字是存在的。如果没有，我们继续进行二进制搜索，直到我们找到该数字，或者直到我们尝试了两次相同的数字，这表明这个数字不存在。

在logN的上限下，计算1E100 (1 googol)只需1毫秒，因此对于算法方法来说非常有效。