3 ms·
The following implementation of the Chudnovsky algorithm takes about the same time in Gambit-Scheme on M1 Mac for 20M digits (and allocates 16GB): (define (p
by bakul 4y ago
The following implementation of the Chudnovsky algorithm takes about the same time in Gambit-Scheme on M1 Mac for 20M digits (and allocates 16GB):
(define (pi digits)
(let* ((A 13591409)
(B 545140134)
(C 640320)
(C^3/24 (quotient (expt 640320 3) 24))
(D 12))
(define (-1^n*g n g) (if (odd? n) (- g) g))
(define (split m n)
(if (= 1 (- n m))
(let* ((6n (* 6 n)) (g (* (- 6n 5) (- (+ n n) 1) (- 6n 1))))
(list g (* C^3/24 (expt n 3)) (* (-1^n*g n g) (+ (* n B) A))))
(let* ((mid (quotient (+ m n) 2))
(gpq1 (split m mid))
(gpq2 (split mid n))
(g1 (car gpq1)) (p1 (cadr gpq1)) (q1 (caddr gpq1))
(g2 (car gpq2)) (p2 (cadr gpq2)) (q2 (caddr gpq2)))
(list (* g1 g2) (* p1 p2) (+ (* q1 p2) (* q2 g1))))))
(let* ((num-terms (inexact->exact (floor (+ 2 (/ digits 14.181647462)))))
(sqrt-C (integer-sqrt (* C (expt 100 digits))))
(gpq (split 0 num-terms))
(g (car gpq)) (p (cadr gpq)) (q (caddr gpq)))
(quotient (* p C sqrt-C) (* D (+ q (* p A)))))))
I am surprised it is much slower in Julia (as per what is noted in the gist).
- juliusgeo 4y agoInteresting! What does Gauss Legendre look like in Gambit-Scheme? I also wonder if perhaps there is a way to get the Julia implementation up to par…
- bakul 4y agoI will give it a try, time permitting. The scheme code above is a variation of something I converted from CommonLisp (orignal code by Robert Smith) in 2013. Its speed depends on bignum multiplication speed. I imagine most serious bignum implementations use the Schönhage–Strassen algorithm for very large numbers.
- juliusgeo 4y agoThanks for looking into this! That is interesting, I wonder if it varies a great deal based on the implementation of the BigFloat rather than the language that it is implemented in.