This version is about 2x faster :
(import 'java.lang.Math)
(import 'java.math.MathContext)
(import 'java.math.RoundingMode)
(import 'java.math.BigInteger)
(import 'java.math.BigDecimal)
(defn sb-pi [places]
"Calculates PI digits using the Salamin-Brent algorithm
and Java's BigDecimal class."
(let [digits (+ 10 places) ;; add some guard digits
round-mode BigDecimal/ROUND_DOWN
zero (new BigDecimal 0)
one (new BigDecimal 1)
two (new BigDecimal 2)
four (new BigDecimal 4)]
(defn big-sqrt[#^BigDecimal num]
"Calculates square root using Newton's method."
(defn big-sqrt-int
[#^BigDecimal num
#^BigDecimal x0
#^BigDecimal x1]
(let [x0new x1
x1new (. num divide x0new digits round-mode)
x1tot (. (. x1new add x0new)
divide two digits round-mode)]
(if (. x0 equals x1)
x1tot
(recur num x1 x1tot))))
(big-sqrt-int
num zero (. BigDecimal valueOf (Math/sqrt (. num
doubleValue)))))
(defn sb-pi-int
[#^BigDecimal a #^BigDecimal b #^BigDecimal t #^BigDecimal x
#^BigDecimal y]
(let
[#^BigDecimal y1 a
#^BigDecimal a1 (. (. a add b) divide two digits round-
mode)
#^BigDecimal b1 (big-sqrt (. b multiply y1))
#^BigDecimal t1 (. t subtract
(. x multiply
(. (. y1 subtract a1)
multiply (. y1 subtract a1))))
#^BigDecimal x1 (. x multiply two)]
(if (. a equals b)
(. (. (. (. a1 add b1)
multiply (. a1 add b1))
divide (. t1 multiply four) digits round-mode)
setScale places round-mode)
(recur a1 b1 t1 x1 y1))))
(sb-pi-int one
(. one divide (big-sqrt two) digits round-mode)
(. one divide four)
one
nil)))
;(big-sqrt (new BigDecimal 2))))
(println (sb-pi (Integer/parseInt (second *command-line-args*))))
$ time clj pi.clj 1 --> 0.977s
$ time clj pi.clj 10 --> 0.975s
$ time clj pi.clj 100 --> 0.979s
$ time clj pi.clj 1000 --> 1.089s
$ time clj pi.clj 10000 --> 11.759s
The same algorithm in Java is faster than the Clojure version by a
rough factor of 3x :
import java.math.BigDecimal;
import static java.math.BigDecimal.*;
class Pi {
private static final BigDecimal TWO = new BigDecimal(2);
private static final BigDecimal FOUR = new BigDecimal(4);
private static int ROUND_MODE = ROUND_DOWN;
public static void main(String[] args) {
System.out.println(pi(Integer.parseInt(args[0])));
//System.out.println(sqrt(TWO, Integer.parseInt(args[0])));
}
// Salamin-Brent Algorithm
public static BigDecimal pi(final int digits) {
final int SCALE = 10 + digits;
BigDecimal a = ONE;
BigDecimal b = ONE.divide(sqrt(TWO, SCALE), SCALE, ROUND_MODE);
BigDecimal t = new BigDecimal(0.25);
BigDecimal x = ONE;
BigDecimal y;
while (!a.equals(b)) {
y = a;
a = a.add(b).divide(TWO, SCALE, ROUND_MODE);
b = sqrt(b.multiply(y), SCALE);
t = t.subtract(x.multiply(y.subtract(a).multiply(y.subtract
(a))));
x = x.multiply(TWO);
}
return a.add(b)
.multiply(a.add(b))
.divide(t.multiply(FOUR), SCALE, ROUND_MODE)
.setScale(digits, ROUND_MODE);
}
// square root method (Newton's)
public static BigDecimal sqrt(BigDecimal A, final int SCALE) {
BigDecimal x0 = new BigDecimal("0");
BigDecimal x1 = new BigDecimal(Math.sqrt(A.doubleValue()));
while (!x0.equals(x1)) {
x0 = x1;
x1 = A.divide(x0, SCALE, ROUND_MODE);
x1 = x1.add(x0);
x1 = x1.divide(TWO, SCALE, ROUND_MODE);
}
return x1;
}
}
$ time java Pi 1 ----> 0.367s
$ time java Pi 10 ----> 0.366s
$ time java Pi 100 ----> 0.370s
$ time java Pi 1000 ----> 0.518s
$ time java Pi 10000 ----> 3.425s
--~--~---------~--~----~------------~-------~--~----~
You received this message because you are subscribed to the Google
Groups "Clojure" group.
To post to this group, send email to clojure@googlegroups.com
Note that posts from new members are moderated - please be patient with your
first post.
To unsubscribe from this group, send email to
clojure+unsubscr...@googlegroups.com
For more options, visit this group at
http://groups.google.com/group/clojure?hl=en
-~----------~----~----~----~------~----~------~--~---
