Brad Gilbert

I was looking through the code of the Rakudo implementation of Perl6 where I noticed that it defines pi as my constant pi = 3.14159_26535_89793_238e0; ( with an alias my constant π := pi; )
I immediately remembered that for a subset of fractional numbers Perl6 has a type that stores them without the loss of precision that generally accompanies floating point math. That type is of course the Rat (and FatRat) type. So of course I type pi.Rat into the REPL, it then prints 3.141593 which is obviously nowhere near as precise as the result of just typing pi into the REPL 3.14159265358979.
I wanted to see the numerator and denominator values that Perl chose to use so I typed pi.Rat.perl and got <355/113>, which is nowhere near as precise as the Rat type is capable of handling. That was just the entrance to the rabbit hole.

So I went searching online and found the Wikipedia article on approximations of pi. Which lead to the article of continued fractions, that has the sequence [3;7,15,1,292,1,1,1,2,1,3,1,…]. I won’t go into too much detail about continued fractions in this post, for more information I would recommend the Wikipedia article.

That bit about ... had me curious so I looked for a longer list, which lead to a pdf of the first 250 elements.

I think that is enough backstory, time to look at the code that I wrote to play around with this new bit of information.

I started by creating an array of the elements of the continued fraction.

#!/usr/bin/env perl6
use v6;

my @continued-fraction = <
  3 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 1 4 2 6 6 99 1
  2 2 6 3 5 1 1 6 8 1 7 1 2 3 7 1 2 1 1 12 1 1 1 3 1 1 8 1 1 2 1 6 1 1 5 2 2 3
  1 2 4 4 16 1 161 45 1 22 1 2 2 1 4 1 2 24 1 2 1 3 1 2 1 1 10 2 5 4 1 2 2 8 1
  5 2 2 26 1 4 1 1 8 2 42 2 1 7 3 3 1 1 7 2 4 9 7 2 3 1 57 1 18 1 9 19 1 2 18 1
  3 7 30 1 1 1 3 3 3 1 2 8 1 1 2 1 15 1 2 13 1 2 1 4 1 12 1 1 3 3 28 1 10 3 2
  20 1 1 1 1 4 1 1 1 5 3 2 1 6 1 4 1 120 2 1 1 3 1 23 1 15 1 3 7 1 16 1 2 1 21
  2 1 1 2 9 1 6 4 127 14 5 1 3 13 7 9 1 1 1 1 1 5 4 1 1 3 1 1 29 3 1 1 2 2 1 3 1
>\ >>.Int; # converts each element into an integer

To make the rest of the code simpler I created a helper function to get the reciprocal of a number.

multi sub reciprocal ( Rational ::R $rat ){ # Rat or FatRat
  R.new( $rat.denominator, $rat.numerator );
}
multi sub reciprocal ( Int $integer ){
  FatRat.new( 1, $integer );
}

Now that I think about it I could have created a new operator by naming it prefix:<1/>, but that would be a new rabbit hole into precedence and associativity of operators.

Anyway the basic formula is start from the right adding the reciprocal of the last number to the second to last, then adding the reciprocal of the previous result to the third to last, etc until you have added the first number.

@continued-fraction.reverse.reduce( ->$a,$b{ reciprocal($a) + $b } );

The result of which is:

11470607307161132716588859789437271840000798468219508525085448960062402821241970285202775001996784319885772780138414238670402064758
/
3651207706401417758577438077383177726074790491037819709434096621864877255825019944180015629313589782838738505101104441533461516175

A rather difficult to comprehend number. It is also way more precise than is necessary. Its accurate to more than 100 decimal places.

3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482534211706798214808651328230664709384460955

To play around with the precision it would be nice if that algorithm was packaged up in a sub that took a count of elements as the only parameter. ( zero based )

multi sub pi-fraction ( 0 ){ Rat.new( 3 ) }
multi sub pi-fraction ( Int $index ){
  @continued-fraction[0..$index].reverse.reduce(->$a,$b{ reciprocal($a) + $b });
}

Let’s print out the first four fractions

for ^4 { pi-fraction($_).Rat.perl.say }

3.0
<22/7>
<333/106>
<355/113>

At this point it would be a good idea to create some MAIN subroutines. I also want to control how the numbers are printed, so I add some helper subs.

sub as-fraction ( Rational $rat ){
  $rat.numerator ~ ' / ' ~ $rat.denominator
}

sub display ( Rational $number ){
  my $as-fraction = as-fraction $number;
  my $display-width = $as-fraction.chars; # good enough for the chars we are using
  my $padding = ' ' x 8 - ($display-width % 4);
  $padding ~= ' ' x 4 if $display-width < 12;
  $padding ~= ' ' x 4 if $display-width < 8;
  return $as-fraction ~ $padding ~ $number;
}

#| list the first <count> fractions
multi sub MAIN ( 'list', $count = 0 ){
  my $total-count = @continued-fraction.elems;
  my @index = do given $count {
    when 0 { ^$total-count }
    when $_ < 0 { ($total-count + $count) ..^ $total-count }
    default { ^$count }
  }
  for @index -> $index {
    my $pi = pi-fraction $index;
    say display $pi;
  }
}

So now we can run it like this ./pi.p6 list 4

3 / 1               3
22 / 7              3.142857
333 / 106           3.141509
355 / 113           3.141593

I wanted to know what the largest one that could be stored as a Rat so I added this:

#| Display the most precise Rat
multi sub MAIN ( 'max', 'Rat' ) {
  my $max-rat;
  my $max-index;

  for @continued-fraction.keys -> $index {
    my $current = pi-fraction $index;

    last unless $current.Rat; # stop if it can't be represented as a Rat

    $max-rat = $current;
    $max-index = $index;
  }
  say '#', $max-index;
  say display $max-rat;
}

$ ./pi.p6 max Rat
#32
2646693125139304345 / 842468587426513207        3.1415926535897932385

If you are new to Perl6 you might not know why I put comments before the MAIN subs that start with #|. The reason should become apparent if you attempt to run it with the wrong arguments, or in this case, no arguments.

$ ./pi.p6
Usage:
  ./pi.p6 list [<count>] -- list the first <count> fractions
  ./pi.p6 max Rat -- Display the most precise Rat

Without those comments it will just print:

$ ./pi.p6
Usage:
  ./pi.p6 list [<count>]
  ./pi.p6 max Rat

I will probably continue to play around with this some more, as programming in Perl6 is a lot of fun.