Found this while wandering around the internet:
; Copyright (c) 2008 the authors listed at the following URL, and/or
; the authors of referenced articles or incorporated external code:
; http://en.literateprograms.org/Pi_with_Machin's_formula_(Lisp)?action=history&offset=20060307031705
;
; Permission is hereby granted, free of charge, to any person obtaining
; a copy of this software and associated documentation files (the
; "Software"), to deal in the Software without restriction, including
; without limitation the rights to use, copy, modify, merge, publish,
; distribute, sublicense, and/or sell copies of the Software, and to
; permit persons to whom the Software is furnished to do so, subject to
; the following conditions:
;
; The above copyright notice and this permission notice shall be
; included in all copies or substantial portions of the Software.
;
; THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
; EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
; MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
; IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
; CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
; TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
; SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
;
; Retrieved from: http://en.literateprograms.org/Pi_with_Machin's_formula_(Lisp)?oldid=2756
(define (arccot x unity)
(let ((xpower (floor (/ unity x))))
(arccot-plus-helper (* x x) 1 xpower)))
(define (arccot-plus-helper xsq n xpower)
(let ((term (floor (/ xpower n))))
(if (= term 0)
0
(+ (arccot-minus-helper xsq (+ n 2) (floor (/ xpower xsq)))
term))))
(define (arccot-minus-helper xsq n xpower)
(let ((term (floor (/ xpower n))))
(if (= term 0)
0
(- (arccot-plus-helper xsq (+ n 2) (floor (/ xpower xsq)))
term))))
(define (pidigits digits)
(letn (
(unity (pow 10 (+ digits 10)))
(pi (* 4 (- (* 4 (arccot 5 unity)) (arccot 239 unity)))))
(floor (/ pi (pow 10 10)))))
I converted it to newLISP by the obvious means (defun to (define. Didn't have a (expt function so I substituted (pow. On my 32bit machine (pidigits works up to 8 digits(yeah I know--- impressive!) then I get the following:
> (pidigits 6)
3141592
> (pidigits 7)
31415926
> (pidigits 8)
314159265
> (pidigits 9)
-791741031
I'm curious if this is some sort of limit reached with 32bit math? Could someone with a 64bit system give this a run?
--hsm
Quote
I'm curious if this is some sort of limit reached with 32bit math? Could someone with a 64bit system give this a run?
in newLISP all integer math using the +,-,*,/ is 64bit. But you are mixing it with floating point math when using the power and floor functions. Then, when dividing using / you truncate the values to integers again. Convert everthing to pure float using add,sub,mul and div and it will work.
Read here:
http://www.newlisp.org/downloads/newlisp_manual.html#int_float
and look here to see how to format the out-coming float number back to an integer:
http://www.newlisp.org/downloads/newlisp_manual.html#format
Lutz Good catch! Great improvement. Now the limit is understandable:
C:>newlisp H:pi.lsp
newLISP v.9.3.6 on Win32 IPv4, execute 'newlisp -h' for more info.
> (pidigits 18)
"3.141592653589793280"
> (pidigits 19)
"9.223372036854775807"
>
Which I ascribe to the cap on 64 bit unsigned integers being very near 20 digits long at 18,446,744,073,709,551,615, and I presume the next step probably craps out accordingly.
--hsm