Post Reply 
Viète's Formula for PI
06-17-2020, 05:06 PM (This post was last modified: 06-17-2020 08:55 PM by pinkman.)
Post: #1
Viète's Formula for PI
Hi,
Once again, I had another discussion with my cousin about PI (it's either PI or prime numbers...).
He told me about the beauty of the Viète's formula :
   

History here: link

This video shows the formula entered in CAS, and yes! I find it beautiful: https://youtu.be/BEYDBl94UpM
Regards,
Thibault
Find all posts by this user
Quote this message in a reply
06-17-2020, 09:37 PM
Post: #2
RE: Viète's Formula for PI
Cool! Really beautiful to see in the CAS the pi approximation in action! :-)

Ramón
Valladolid, Spain
TI-50, Casio fx-180P, HP48GX, HP50g, HP Prime G2
Visit this user's website Find all posts by this user
Quote this message in a reply
06-18-2020, 12:51 PM (This post was last modified: 06-18-2020 02:49 PM by pinkman.)
Post: #3
RE: Viète's Formula for PI
And the convergence speed is acceptable : for p(18) we have 10 digits.
p(18) means Σ(x,x,1,18) = 171 iterations.
Compared to Wallis formula (https://www.hpmuseum.org/forum/thread-14...ght=wallis), it's absolutely fast!
Find all posts by this user
Quote this message in a reply
06-19-2020, 09:00 PM
Post: #4
RE: Viète's Formula for PI
very nice, thank you!


tiny tip: in CAS, for algebraic or fractional input, the a.b/c key also runs approx(); one less key-stroke than shift-enter.

Cambridge, UK
41CL/DM41X 12/15C/16C DM15/16 17B/II/II+ 28S 42S/DM42 32SII 48GX 50g 35s WP34S PrimeG2 WP43S/pilot
Casio, Rockwell 18R
Find all posts by this user
Quote this message in a reply
06-19-2020, 11:12 PM
Post: #5
RE: Viète's Formula for PI
(06-18-2020 12:51 PM)pinkman Wrote:  Compared to Wallis formula (https://www.hpmuseum.org/forum/thread-14...ght=wallis), it's absolutely fast!

Yes, but it pales in comparison to the Wallis-Wasicki formula :-)

+---+---------------------+---------------------+
| N |         2*W         |        2*W*C        |
+---+---------------------+---------------------+
| 2 | 2.84444444444444444 | 3.14385964912280701 |
| 4 | 2.97215419501133786 | 3.14158816337302932 |
| 6 | 3.02317019200136082 | 3.14159266276745771 |
| 8 | 3.05058999605551092 | 3.14159265357083669 |
|10 | 3.06770380664349896 | 3.14159265358983256 |
|12 | 3.07940134316788626 | 3.14159265358979314 |
|14 | 3.08790206983111306 | 3.14159265358979321 |
+---+---------------------+---------------------+


Only 7 iterations (or 14, depending on how you implement the algorithm) for 18 correct decimal digits
(21697209162666264236130304/6906436179074198667175275 = 3.141592653589793238[633]...).
The Pascal code is available here. It should translate easily into the Prime programming language, but probably no joy with only 12 significant digits...
Find all posts by this user
Quote this message in a reply
06-23-2020, 06:39 PM
Post: #6
RE: Viète's Formula for PI
(06-18-2020 12:51 PM)pinkman Wrote:  And the convergence speed is acceptable : for p(18) we have 10 digits.
p(18) means Σ(x,x,1,18) = 171 iterations.

For faster convergence, you might want to try this algorithm:
Code:

Program  Viete;                                                                                 
                                                                                                                    
Uses Crt;                                                                                                            
                                                                                                                    
var         k:  Byte;                                                                                                
            n: LongInt;                                                                                             
            c: Char;                                                                                                
   d, p, t, v: Extended;                                                                                             
                                                                                                                     
begin                                                                                                               
  ClrScr;                                                                                                           
  WriteLn('+---+----------------------+----------------------+');                                                    
  WriteLn('| k |         2*v          |         ~ pi         |');                                                    
  WriteLn('+---+----------------------+----------------------+');                                                    
  v := 1;                                                                                                            
  n := 1;                                                                                                           
  d := 0;                                                                                                           
  for k := 1 to 15 do                                                                                                
    begin                                                                                                           
      n := n + n;                                                                                                   
      d := Sqrt(d + 2);                                                                                             
      v := v*d;                                                                                                      
      t := 6*v*v - 2/3;                                                                                             
      p := 2*n/v*(t + 1)/t;                                                                                         
      WriteLn('|',k:2,' |',2*n/v:21:18,' |',p:21:18,' |')
    end;                                                                                                             
  WriteLn('+---+----------------------+----------------------+');                                                   
  c := ReadKey                                                                                                       
end.

+---+----------------------+----------------------+
| k |         2*v          |         ~ pi         |
+---+----------------------+----------------------+
| 1 | 2.828427124746190098 | 3.077994223988500989 |
| 2 | 3.061467458920718174 | 3.137427049842311165 |
| 3 | 3.121445152258052286 | 3.141329731135014592 |
| 4 | 3.136548490545939264 | 3.141576182507746574 |
| 5 | 3.140331156954752913 | 3.141591623553754287 |
| 6 | 3.141277250932772868 | 3.141592589203296386 |
| 7 | 3.141513801144301077 | 3.141592649565492850 |
| 8 | 3.141572940367091385 | 3.141592653338272210 |
| 9 | 3.141587725277159701 | 3.141592653574073140 |
|10 | 3.141591421511199975 | 3.141592653588810732 |
|11 | 3.141592345570117743 | 3.141592653589731833 |
|12 | 3.141592576584872667 | 3.141592653589789401 |
|13 | 3.141592634338562990 | 3.141592653589793000 |
|14 | 3.141592648776985670 | 3.141592653589793224 |
|15 | 3.141592652386591346 | 3.141592653589793238 |
+---+----------------------+----------------------+
Find all posts by this user
Quote this message in a reply
06-23-2020, 09:58 PM
Post: #7
RE: Viète's Formula for PI
Quote:Yes, but it pales in comparison to the Wallis-Wasicki formula :-)

Here is a quick PPL port of your Wallis-Wasicki implementation:

Code:


EXPORT WallisWasicki()
BEGIN
LOCAL i, k, m, n;
LOCAL c, d, t, w;
  PRINT("+---+-------------+------------+");
  PRINT("| N |     2*W     |    2*W*C   |");
  PRINT("+---+-------------+------------+");
  FOR k := 1 TO 7 DO
      n := 2*k;
      w := 1;
      c := 0;
      d := 4*n;
      m := n + 1;
      FOR i := 1 TO k DO
          t := i*(16*i - 8);
          w := w*t*t/(t*(t - 2) - 3);
          c := (m - 4)/(d + m/(2 - c));
          m := m - 2
        END;
      c := 1 - c;
      PRINT("|" + n + " |" + 2*w + " |" + 2*w*c + " |");
    END;
  PRINT("+---+-------------+------------+");
  RETURN 2*w*c;
END;

The terminal output is not well formatted, but still easy to read.

I’m also porting your 314 pi digits algorithm, I had to stop (work...) but it will be ready in a few hours (if I find the time).
Find all posts by this user
Quote this message in a reply
06-23-2020, 10:04 PM
Post: #8
RE: Viète's Formula for PI
(06-23-2020 06:39 PM)Gerson W. Barbosa Wrote:  
(06-18-2020 12:51 PM)pinkman Wrote:  And the convergence speed is acceptable : for p(18) we have 10 digits.
p(18) means Σ(x,x,1,18) = 171 iterations.

For faster convergence, you might want to try this algorithm:
[...]

Yes! But... I really love to see the CAS in action, even if it reaches its limit in terms of recursivity.
Find all posts by this user
Quote this message in a reply
06-23-2020, 10:52 PM (This post was last modified: 07-15-2020 05:07 PM by Gerson W. Barbosa.)
Post: #9
RE: Viète's Formula for PI
(06-23-2020 06:39 PM)Gerson W. Barbosa Wrote:  ...
For faster convergence, you might want to try this algorithm:
...

The following appears to be better (more tests required) . Only one line has been changed. Anyway, here is it again:

Code:
Program  Viete;                                                                                 
                                                                                                                    
Uses Crt;                                                                                                            
                                                                                                                    
var         k: Byte;                                                                                                
            n: LongInt;                                                                                             
            c: Char;                                                                                                
   d, p, t, v: Extended;                                                                                             
                                                                                                                     
begin                                                                                                               
  ClrScr;                                                                                                           
  WriteLn('+---+----------------------+----------------------+');                                                    
  WriteLn('| k |         2*v          |         ~ pi         |');                                                    
  WriteLn('+---+----------------------+----------------------+');                                                    
  v := 1;                                                                                                            
  n := 1;                                                                                                           
  d := 0;                                                                                                           
  for k := 1 to 15 do                                                                                                
    begin                                                                                                           
      n := n + n;                                                                                                   
      d := Sqrt(d + 2);                                                                                             
      v := v*d;                                                                                                      
      t := 6*v*v - 27/10 - 1/(n*n);                                                                                            
      p := 2*n/v*(t + 1)/t;                                                                                         
      WriteLn('|',k:2,' |',2*n/v:21:18,' |',p:21:18,' |')
    end;                                                                                                             
  WriteLn('+---+----------------------+----------------------+');                                                   
  c := ReadKey                                                                                                       
end.

+---+----------------------+----------------------+
| k |         2*v          |         ~ pi         |
+---+----------------------+----------------------+
| 1 | 2.828427124746190098 | 3.140960508696045357 |
| 2 | 3.061467458920718174 | 3.141593674140638978 |
| 3 | 3.121445152258052286 | 3.141592707164788547 |
| 4 | 3.136548490545939264 | 3.141592654569488290 |
| 5 | 3.140331156954752913 | 3.141592653605653769 |
| 6 | 3.141277250932772868 | 3.141592653590043215 |
| 7 | 3.141513801144301077 | 3.141592653589797153 |
| 8 | 3.141572940367091385 | 3.141592653589793300 |
| 9 | 3.141587725277159701 | 3.141592653589793240 |
|10 | 3.141591421511199975 | 3.141592653589793239 |
+---+----------------------+----------------------+



P.S.:

This yields 1.8 digits per iteration, three times as much when compared to the plain Viète's formula.

P.P.S.:

This new formula is easier to program and will yield slightly more than two digits per term:

\(\pi \approx 2\left ( \frac{4}{3} \times \frac{16}{15}\times \frac{36}{35}\times\frac{64}{63} \times \cdots \times \frac{ 4n ^{2}}{ 4n ^{2}-1}\right )\cdot \left ( 1+\frac{2}{8n+3+\frac{3}{8n+4+\frac{15}{8n+4+ \frac{35}{8n+4 + \frac{63}{\dots\frac{\ddots }{8n+4+\frac{4n^{2}-1}{8n+4}}} }}} } \right )\)

TurboBCD program:

Code:

Program Wallis_Wasicki;

var d, n, i: Byte;
          t: Integer;
       c, w: Real;

begin
  ClrScr;
  WriteLn('+---+---------------------+---------------------+');
  WriteLn('| N |         2*W         |     2*W*C/(C-2)     |');
  WriteLn('+---+---------------------+---------------------+');
  for n := 1 to 8 do
    begin
      c := 0;
      d := 8*n + 4;
      w := 1;
      for i := n downto 1 do
        begin
          t := 4*i*i;
          w := w*t/(t - 1);
          c := (t - 1)/(c + d)
        end;
      c := c + d + 1;
      WriteLn('|',n:2,' |',2*w:20:17,' |',2*c*w/(c - 2):20:17,' |')
    end;
    WriteLn('+---+---------------------+---------------------+')
end.

+---+---------------------+---------------------+
| N |         2*W         |     2*W*C/(C-2)     |
+---+---------------------+---------------------+
| 1 | 2.66666666666666666 | 3.14074074074074073 |
| 2 | 2.84444444444444446 | 3.14159848961611078 |
| 3 | 2.92571428571428570 | 3.14159260997123044 |
| 4 | 2.97215419501133786 | 3.14159265392705764 |
| 5 | 3.00217595455690690 | 3.14159265358714120 |
| 6 | 3.02317019200136082 | 3.14159265358981426 |
| 7 | 3.03867362888341912 | 3.14159265358979309 |
| 8 | 3.05058999605551094 | 3.14159265358979325 |
+---+---------------------+---------------------+


In this table N is the number of terms, W is the Wallis Product evaluated to N terms and C/(C - 2) is the correction factor.
Find all posts by this user
Quote this message in a reply
06-23-2020, 11:00 PM
Post: #10
RE: Viète's Formula for PI
(06-23-2020 09:58 PM)pinkman Wrote:  
Quote:Yes, but it pales in comparison to the Wallis-Wasicki formula :-)

Here is a quick PPL port of your Wallis-Wasicki implementation:

...

I’m also porting your 314 pi digits algorithm, I had to stop (work...) but it will be ready in a few hours (if I find the time).

Thank you very much for the PPL port!

Perhaps it's time I should get myself a Prime. But I think I will wait until a good arbitrary precision package is available, either built-in or third-party.

Gerson.
Find all posts by this user
Quote this message in a reply
06-24-2020, 01:15 PM
Post: #11
RE: Viète's Formula for PI
The port of the 200-314-997 decimals will need the CAS, which has capabilities for manipulating big integers but not arbitrary precision decimal numbers.
I’ll check if I can avoid using decimals by calculating on a base of estimated_pi * 10^200 (or less, depending on CAS limit).
Need a few hours more!
Find all posts by this user
Quote this message in a reply
06-29-2020, 05:52 AM
Post: #12
RE: Viète's Formula for PI
I wish the extremely fast HP Prime had the LongFloat Library integrated.
- -
VP
Find all posts by this user
Quote this message in a reply
06-29-2020, 10:54 PM
Post: #13
RE: Viète's Formula for PI
Me too!

Thibault - not collector but in love with the few HP models I own - Also musician : http://walruspark.co
Find all posts by this user
Quote this message in a reply
06-30-2020, 03:05 PM (This post was last modified: 07-09-2020 04:07 AM by compsystems.)
Post: #14
RE: Viète's Formula for PI
There must be some open source LongFloat library, if so Xcas could include it.
[Image: 01-38.jpg]
https://keisan.casio.com/
Find all posts by this user
Quote this message in a reply
07-16-2020, 04:42 PM
Post: #15
RE: Viète's Formula for PI
(06-23-2020 09:58 PM)pinkman Wrote:  Here is a quick PPL port of your Wallis-Wasicki implementation:

[...]

Gerson sent me a PM to thank me for having posted this code to hpcalc.org, but in fact I did not post anything, I guess Eric did.
Funny, and good idea, but the credits come to Gerson and his continuous fraction quick convergence for Wallis product.

Thibault - not collector but in love with the few HP models I own - Also musician : http://walruspark.co
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 1 Guest(s)