-
Notifications
You must be signed in to change notification settings - Fork 0
/
catalan.lisp
31 lines (29 loc) · 1.11 KB
/
catalan.lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
;; Catalan
(defun fpcatalan ()
(let ((catalan (comp-catalan (+ fpprec 12))))
(bcons (list (fpround (car catalan)) (cadr catalan)))))
(defun comp-catalan (prec)
(let ((fpprec prec)
nr tt bb qq)
(setq nr (ceiling fpprec (log (/ 729 4) 2d0))) ;; nr of summands
(multiple-value-setq (tt bb qq) (split-catalan 1 (1+ nr)))
(setq qq (intofp qq))
(fpquotient (intofp tt)
(fptimes* (intofp bb)
(list (car qq) (+ 6 (cadr qq)))))))
(defun split-catalan (i j)
(let (aa bb pp qq tt)
(if (= (- j i) 1)
(setq aa (+ 15 (* i (- (* 580 i) 184)))
bb (* (expt i 3) (- (* 2 i) 1))
pp (ash bb 5)
qq (expt (* (- (* 6 i) 1) (- (* 18 i) 15)) 2)
tt (* aa pp))
(let ((m (ash (+ i j) -1)))
(multiple-value-bind (tl bl ql pl) (split-catalan i m)
(multiple-value-bind (tr br qr pr) (split-catalan m j)
(setq pp (* pl pr)
qq (* ql qr)
bb (* bl br)
tt (+ (* br qr tl) (* bl pl tr)))))))
(values tt bb qq pp)))