Home > Archive > Scheme > February 2007 > go drscheme going on freebsd!! drscheme, making i calculate pi to 30 places?
You are viewing an archived Text-only version of the thread.
To view this thread in it's original format and/or if you want to reply to
this thread please [click here]
| Author |
go drscheme going on freebsd!! drscheme, making i calculate pi to 30 places?
|
|
| gavino 2007-02-07, 7:09 pm |
| How do I set dr scheme to calculate pi to 30 places?
Anyone?
| |
| Prabhakar Ragde 2007-02-07, 10:06 pm |
| "gavino" <gavcomedy@gmail.com> writes:
> How do I set dr scheme to calculate pi to 30 places?
> Anyone?
(+ 3
(foldr (lambda (x y) (/ 1 (+ x y))) 0
'(7 15 1 292 1 1 1 2 1 3 1 14 2 1 1 2 2 2 2 1 84 2 1 1 15 3 13)))
--
Prabhakar Ragde, Professor plragde at uwaterloo dot ca
Cheriton School of Computer Science http://www.cs.uwaterloo.ca/~plragde
Faculty of Mathematics DC 1314, (519)888-4567,x4660
University of Waterloo Waterloo, Ontario CANADA N2L 3G1
| |
| gavino 2007-02-07, 10:06 pm |
| hm
well
hm
okayyy
heh
uh
I dont this this is a generizable equaiton...
| |
| Abdulaziz Ghuloum 2007-02-09, 7:10 pm |
| On Feb 9, 3:34 pm, Nathan Thern <nth...@yahoo.com> wrote:
> 6162494174639481259093616275458736114027
149893421880624492184769763/
> 1961582819337766313385796154546284449822
433760466765902190005452800
>
> which is pi, accurate to 32 decimal places.
That pi is not even close to 3.
Aziz,,,
| |
| Abdulaziz Ghuloum 2007-02-09, 7:10 pm |
| On Feb 9, 4:47 pm, "Abdulaziz Ghuloum" <aghul...@gmail.com> wrote:
> On Feb 9, 3:34 pm, Nathan Thern <nth...@yahoo.com> wrote:
>
>
>
> That pi is not even close to 3.
Silly me. I didn't see the slash at the end :-)
| |
| Jan Skibinski 2007-02-10, 7:09 pm |
| On Feb 7, 7:17 pm, "gavino" <gavcom...@gmail.com> wrote:
> How do I set dr scheme to calculate pi to 30 places?
> Anyone?
If you sign up to Google account you will be able to browse the
"Ramanujan's Lost Notebook" By George Eyre Andrews, Bruce C. Berndt.
I am not providing any link, because the one I have seems to be a
Google account specific. But you can easily trace the book in Google.
Anyway, here is an excerpt from my haskell module fractions.hs, which
calculates the pi with any desired accuracy eps, based on Ramanujah's
formula:
==================
pi' eps =
--
-- pi with accuracy eps
--
-- Based on Ramanujan formula, as described in Ref. 3
-- Accuracy: extremely good, 10^-19 for one term of continued
-- fraction
--
(sqrt' eps d) / (approxCF eps (fromTaylorToCF s x))
where
x = 1//(640320^3)::Fraction
s = [((-1)^k*(fac (6*k))//((fac k)^3*(fac (3*k))))*((a*k+b)//c) |
k<-[0..]]
a = 545140134
b = 13591409
c = 426880
d = 10005
=======================
I posted it in 2000 to Haskell mailing list. I have no time now to
translate it to Scheme though, but you are welcome to this and other
pages of mine:
http://www.haskell.org/haskellwiki/...ental_functions
Regards,
Jan
| |
| Michel Salim 2007-02-11, 4:12 am |
| On Feb 10, 10:42 am, "Jan Skibinski" <jan.skibin...@gmail.com> wrote:
> On Feb 7, 7:17 pm, "gavino" <gavcom...@gmail.com> wrote:
>
>
> If you sign up to Google account you will be able to browse the
> "Ramanujan's Lost Notebook" By George Eyre Andrews, Bruce C. Berndt.
> I am not providing any link, because the one I have seems to be a
> Google account specific. But you can easily trace the book in Google.
>
Which page is that reference on?
> Anyway, here is an excerpt from my haskell module fractions.hs, which
> calculates the pi with any desired accuracy eps, based on Ramanujah's
> formula:
>
I'm using the Bailey-Borwein "One Billion Digits of Pi" paper:
http://www.cecm.sfu.ca/organics/pap...wein/index.html
and, using algorithm1 (rather than sum1), surprisingly the
approximation of pi does not improve beyond 2 iterations.
> (nth-a 1)
0.31830988693116063
> (/ 1 (nth-a 1))
3.1415926462135473
> (/ 1 (nth-a 2))
3.1415926535897927
> (/ 1 (nth-a 10))
3.1415926535897927
Is there a way to compute square roots with arbitrary precision in
Scheme? Otherwise the sum method, while slower, is the easiest to
implement (algorithm2 involves taking the fifth root)
(define a0 (- 6 (* 4 (sqrt 2))))
(define y0 (- (sqrt 2) 1))
(define (next-y cur-y)
(let ([t (sqrt (sqrt (- 1 (pow4 cur-y))))])
(/ (- 1 t) (+ 1 t))))
(define (next-a a_n y_n+1 n)
(- (* (pow4 (+ 1 y_n+1)) a_n)
(* (fast-pow 2 (+ (* 2 n) 3))
y_n+1
(+ 1 y_n+1 (square y_n+1)))))
(define (nth-a n)
(let loop ([an a0]
[yn y0]
[i 0])
(if (= i n)
an
(let* ([yn++ (next-y yn)]
[an++ (next-a an yn++ i)])
(loop an++ yn++ (add1 i))))))
Regards,
--
Michel
| |
| Michel Salim 2007-02-11, 4:12 am |
| On Feb 11, 1:12 am, "Michel Salim" <michel.sa...@gmail.com> wrote:
> On Feb 10, 10:42 am, "Jan Skibinski" <jan.skibin...@gmail.com> wrote:> On Feb 7, 7:17 pm, "gavino" <gavcom...@gmail.com> wrote:
>
>
>
> Which page is that reference on?
>
>
> I'm using the Bailey-Borwein "One Billion Digits of Pi" paper:http://www.cecm.sfu.ca/organics/pap...wein/index.html
>
> and, using algorithm1 (rather than sum1), surprisingly the
> approximation of pi does not improve beyond 2 iterations.
>
> 0.31830988693116063
> 3.1415926462135473
> 3.1415926535897927
>
> 3.1415926535897927
>
> Is there a way to compute square roots with arbitrary precision in
> Scheme? Otherwise the sum method, while slower, is the easiest to
> implement (algorithm2 involves taking the fifth root)
>
Partly answering myself: by iteration one could get arbitrarily close
to the real square root of a number. But how close is close enough?
Presumably the error bound needs to get tighter as the number of
iterations increases.
--
Michel
| |
| Jan Skibinski 2007-02-11, 4:12 am |
| Hi Michel,
> Which page is that reference on?
Sorry for the sloppy reference. Go to Google, switch to category
"search books" and enter "Ramanujan's Lost Notebook". You can browse
some of its pages on line, including interesting introduction. I cited
this book because the story of Ramanujan sounds like from Hollywood
movies. The Wikipedia article says that someone is making a movie
about him this month. I vaguely remember some reference to him as a
genius who use to imagine numbers not as digits - like most of us do,
but as solutions to some equations instead.
> I'm using the Bailey-Borwein "One Billion Digits of Pi" paper:http://www.cecm.sfu.ca/organics/pap...wein/index.html
>
> and, using algorithm1 (rather than sum1), surprisingly the
> approximation of pi does not improve beyond 2 iterations.
I was using a sum-based formula, but not the one cited in Bailey-
Borwein.
Since two references in my Haskell module are no longer valid, I
looked around and found this link, which relates to my approach.
B503 Algorithms Determining Pi Write up 3, Jeremy Engle, April 21,
2004
http://www.cs.indiana.edu/~jtengle/...b503writeup.pdf
We both seem to basing our comptations on the same Ramanujan sum, but
Engel is summing directly, while I am first converting the series to
continued fraction and only then proceed with computations. Since this
happens in Haskell the continued fraction is computed lazily
(theoretically to infinity :-). Last time I tested the results (6
years ago) I was amazed to see that claims of extreme good accuracy
for this formula were not exaggerated. Just the first term leads to
the accuracy of 10^-19.
Why is it so? If I may use the analogy to zero-pole analysis in DSP:
changing positions of zeros does not have such a dramatic effects on
transfer functions as changing positions of poles. Direct summing is
like manipulating zeros, continued-fraction approach relates to poles.
I can only take a partial credit for this, since the original idea is
Ramanujan's, and continued-fraction approach was either borrowed from
Gosper or from Potts (both cited in the fraction.hs module). I might
have changed something in the approach and I definitely coded it
myself in Haskell, but that's far from being original. Unfortunately,
Pott's both articles were based and are no longer available and I
cannot verify who is to be credited.
Jan
P.S. On the light side, this is what I found at:
http://www.boo.net/~jasonp/pi-ref.txt
"1897
Indiana (U.S.) state legislature almost passes a law declaring that
pi is equal to 16/5 (the value varies depending on the reference.)
They do it because the author of the 'discovery' offers to let them
use it for free, while other states would have had to pay royalties
on it. This may also be the source of a rumor about the government
doing it for religious reasons."
Below is a more detailed comparison of the two approaches.
Engle (2004):
============
k6 sqrt(k3)
pi = -------------
S
where:
(-1)^n (6n)! (k2 + nk1)
S = SUM(n=0 .. n=inf): --------------------------------
(n!)^3(3n)! (8k4*k5)^n
k1 = 545140134
k2 = 13591409
k3 = 10005
k4 = 100100025
k5 = 327843840
k6 = 426880
Skibinski (1998-2000):
=====================
Using Engle's constants, and recovering formulas from the Haskell
code,
shown in the previous post, here is some quasi-Scheme code mixed with
mathematical formula for the sum, to show similarity with the sum used
by Engle. Format of continued fractions is shown at the bottom of this
post.
(sqrt eps k3)
pi =
----------------------------------------------------------------------------
(cont-fraction->fraction eps (taylor->cont-fraction S x))
(-1)^n (6n)! (k2 + nk1)
S = SUM(n=0 .. n=inf): --------------------------------
(n!)^3(3n)! k6
1
x = ----------
640320^3
eps = desired accuracy, given as fraction, such as 1/10^40
procedure: (sqrt eps x)
Compute square root of x with accuracy eps
using continued fractions CF.
The CF pattern is: [(m x-m^2) (2m x-m^2) (2m x-m^3)....]
where m is the biggest integer, such that x-m^2 >= 0
procedure: (cont-fraction->fraction eps cont-fraction)
Approximate infinite continued fraction x by a fraction,
evaluating from left to right, and stopping when
accuracy eps is achieved, or when a partial numerator
is zero -- as it indicates the end of CF.
procedure: (taylor->cont-fraction taylor x)
Convert infinite number of terms of Taylor expansion of
a function f(x) to an infinite continued fraction,
where s = [s0 s1 s2 s3....] is a list of Taylor
series coefficients, such that f(x)=s0 + s1*x + s2*x^2....
Require: No Taylor coefficient is zero
The following representation of continued fractions has been used:
Continued fraction: Its Scheme representation:
================== ==========================
b0 + a0
------- ==> [(b0 a0) (b1 a1) (b2
a2).....]
b1 + a1
-------
b2 + ...
where "a's" and "b's" are Fractions.
Note:
=====
Many continued fractions could be represented simply as lists of
integers
[b1 b2 b3 b4..], with coefficients a's assumed 1.
However, there are some useful continued fractions that are
given with fractional coefficients, where "a" or "b" are either
fractions or integers.
A fractional form can always be converted to an integer form, but
a conversion process is not always simple and such an effort is not
always worth of the achieved savings in the storage space or the
computational efficiency.
| |
| Michel Salim 2007-02-11, 7:12 pm |
| Hi Jan,
On 11 f=E9v, 05:04, "Jan Skibinski" <jan.skibin...@gmail.com> wrote:
> Hi Michel,
>
>
> Sorry for the sloppy reference. Go to Google, switch to category
> "search books" and enter "Ramanujan's Lost Notebook". You can browse
> some of its pages on line, including interesting introduction. I cited
> this book because the story of Ramanujan sounds like from Hollywood
> movies. The Wikipedia article says that someone is making a movie
> about him this month. I vaguely remember some reference to him as a
> genius who use to imagine numbers not as digits - like most of us do,
> but as solutions to some equations instead.
>
Ah, OK. I found the book alright, but was trying to find the equation
there. Incidentally, the first Ramanujan anecdote I heard was that of
Hardy visiting him in the hospital and remarking on the pedestrian
nature of the cab license plate that took him there, upon which
Ramanujan leaped up from bed and showed otherwise:
http://en.wikipedia.org/wiki/1729_(number)
> Skibinski (1998-2000):
> =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
3D=3D=3D=3D=3D=3D=3D=3D
> Using Engle's constants, and recovering formulas from the Haskell
> code,
> shown in the previous post, here is some quasi-Scheme code mixed with
> mathematical formula for the sum, to show similarity with the sum used
> by Engle. Format of continued fractions is shown at the bottom of this
> post.
>
> (sqrt eps k3)
> pi =3D
> -------------------------------------------------------------------------=
---
> (cont-fraction->fraction eps (taylor->cont-fraction S x))
>
> (-1)^n (6n)! (k2 + nk1)
> S =3D SUM(n=3D0 .. n=3Dinf): --------------------------------
> (n!)^3(3n)! k6
> 1
> x =3D ----------
> 640320^3
>
> eps =3D desired accuracy, given as fraction, such as 1/10^40
>
Aha! Thanks. Hmm, an interesting extension would be to make every
function that deals with approximation returns another function that
can be used to resume the computation, given a new accuracy target;
that way one does not have to commit to a specific accuracy target
upfront.
--
Michel
| |
| Jan Skibinski 2007-02-11, 7:12 pm |
| Hi Michel,
> Aha! Thanks. Hmm, an interesting extension would be to make every
> function that deals with approximation returns another function that
> can be used to resume the computation, given a new accuracy target;
> that way one does not have to commit to a specific accuracy target
> upfront.
Good thought!
One might also take a shot at combining laziness with continuations.
I like the former mainly for a code clarity since I have no illusions
about its efficiency. But I know that both languages, Haskell and
Scheme, can have both features if they want to.
Disclaimer: When writing fraction.hs I was not interested in pushing
the envelope to "one billion digits" for sqrt-2, pi, e and such, but
to have a workable module of some basic functions (sin cos tan arctan
log sqrt etc.) providing results with a desired accuracy - as is
sometimes needed in certain calculations. I originally intended to
extend it to a domain of complex fractions as well, but never found
enough motivation to pursue this topic.
Jan
| |
| Jan Skibinski 2007-02-13, 4:16 am |
| On Feb 7, 8:28=C2=A0pm, "gavino" <gavcom...@gmail.com> wrote:
> hm
> well
> hm
> okayyy
> heh
> uh
> I dont this this is a generizable equaiton...
Take a look at the beginning of section "Continued fraction expansion
of pi" in wikipedia article, http://en.wikipedia.org/wiki/Continued_fractio=
n,
and you should see the explanation for the list '(7 15 1 292 1 1 1 2 1
3 1 14 2 1 1 2 2 2 2 1 84 2 1 1 15 3 13) used by Prabhakar Ragde in
his representation of pi -- in the post you were responding to.
I must admit, I also did not recognize it at first as a continued
fraction.
But what I like in this wikipedia article is even simpler
representation of pi, which uses *generalized* continued fraction - a
format I have described somewhere else in this thread:
((3 1^2) (6 3^2) (6 5^2) (6 7^2) (6 9^2) (6 11^2) ....)
Cool! As opposed to my original implementation, this formula does not
depend on (sqrt eps 1005), and could be used directly in fraction.hs,
like this:
---------------------------------------------
pi eps =3D approxCF eps s
where
s =3D [3,1]:[[6, n^2] | n->[3,5..]]
----------------------------------------------
In case you wonder about syntax of s: this is a so-called list
comprehension - a powerful construct, which can be easily defined in
Scheme via macros.
Swindle has it, and I defined it too somewhere on this list. Using my
macro, the s above could be defined in Scheme, more or less, as:
(define s (cons '(3 1)
(list-of (list 6 (* n n)) =C2=A6 (=E2=88=88 n (enumerate 3 5 infinity)))=
))
There is one remaining problem: enumeration to infinity.There are two
choices: either to replace lists with delayed streams, or to switch to
a lazy version of Scheme. If I remember correctly, DrScheme has a
"plugguble module" with lazy Scheme.
Jan
|
|
|
|
|