Use rational numbers for arbitrary precision

This Item was suggested by Shawn Corey as part of our Free eBook give-away. He’s our September winner!

Perl is mostly a platform-agnostic and portable language, but there are a few corners where you notice that Perl isn’t in control of everything. In Item 14: Handle big numbers with bignum, you saw that Perl’s normally limited in the maximum integer value it can represent, but the bignum pragma takes that power away from the system and handles it all within Perl, giving you virtually unlimited integral magnitudes. You can do the same thing for rational numbers, too.

The first answer in perlfaq4 is for “Why am I getting long decimals (eg, 19.9499999999999) instead of the numbers I should be getting (eg, 19.95)?”. It’s a common problem in computer programming because programs typically rely on the underlying architecture to store its numbers (there are many languages where this is not true, though). In short, given a certain number of bits, you cannot exactly store every number. Fractional numbers experience a bit of drift, and Perl performs arithmetic immediately.

One way you can get around this problem is to not have a floating point number. That is, keep it as a ratio as long as possible. In that way, you’re really storing two integers.

Before you get to far, don’t get too angry that you’re about to see a transcendental number in examples about rational numbers. It’s just a more interesting number because people recognize it on the street, and it provides an interesting performance case at the end of the Item.

Here’s an example where you have a subroutine that returns π up to 50 decimal places. That might not be exact enough for you, but it fits nicely on a line. Notice that you can actually call the subroutine π if you use the utf8 pragma (Item 72: Use Unicode in your source code):

use utf8;
use 5.010;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

say π;

The output truncates the value:

3.14159265358979

This truncation isn’t the problem though. The number that you display isn’t what Perl stores. It’s a formatted version based on Perl’s best guess on what it should show (that’s the %g format) based on what it stored. You can change that to display two decimal places:

use utf8;
use 5.010;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

printf "%.2f\n", π;
say π;

Now you see that you’re merely displaying it differently without changing the value:

3.14
3.14159265358979

So far, that’s the “right” number insofar as it shows the correct digit in each position even if it is truncated (and not rounded).

Perl’s default number format is %g, which has its own idea about how many decimal places to display, which might be vary for different systems. If you don’t want to use that, you can specify your own format. Try displaying all 50 decimal places of π:

use utf8;
use 5.010;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

printf "%.50f\n", π;

Now the number goes wonky around the 16th decimal point (your precision may vary):

3.14159265358979311599796346854418516159057617187500

What if that’s really important to you though? What if you need the number to be exactly right, not just close to right? What if you need to land your spacecraft on the asteroid, not just miss it by millimeters?

You might think that adding bigrat, but at first blush it might seem like it doesn’t do anything:

use utf8;
use 5.010;
use bigrat;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

printf "%.50f\n", π;

You get the same output:

3.14159265358979311599796346854418516159057617187500

The module behind bigrat, Math::BigRat, represents numbers as objects, completely skipping Perl’s built-in number handling. However, you reinvoke the built-in number handling with a printf number format, you’re back to where you started. When you treat it as a string, you get the right value:

use utf8;
use 5.010;
use bigrat;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

say π

The output is only missing the trailing 0:

3.1415926535897932384626433832795028841971693993751

It’s easier to see what’s going if you do some math with the number. In this case, you just add 0:

use utf8;
use 5.010;
use bigrat;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }

say π + 0

Now you see the number as a rational number, represented as two big integers:

31415926535897932384626433832795028841971693993751/10000000000000000000000000000000000000000000000000

Let’s ignore that π is an irrational number. For now it’s a number that we know has a lot of decimal places.

Since bigrat doesn’t really deal with floating point numbers, it doesn’t have a chance to lose precision, at least until you decide to do the division. You can add fractions:

use 5.010;
use bigrat;

say 1/4 + 1/3 + 1/5;

Your result is still a fraction and is always exactly right no matter what your computer thinks it can store:

47/60

This might be important if you need to model the universe. It’s important to get certain values right:

use utf8;
use 5.010;
use bigrat;

sub π () { 3.1415926535_8979323846_2643383279_5028841971_6939937510 }
sub e () { 4.80320427e-10 }
sub c () { 29979245800 }
sub h () { 6.62606896e-27 }

my $a = (e**2) / ( h/(2*π) * c );

say $a;

That looks like a big number because its string representation takes up a lot of space on the page:

7247896550101266775674757583548458886480233486241129977318779906079/993222750197951840000000000000000000000000000000000000000000000000000

It’s really not so large, and when you are ready to use it, you have as much precision as you need.

A note on performance

When you want to extreme accuracy you’d get through bigrat, you might lose on performance, depending on what you are doing. All of the numbers are objects now, so everytime you want to do something you’re calling methods.

There’s a series you can use to get progressively closer to π. Here’s the Perl implementation:

use utf8;
# use bigrat;

my $π = 0;
foreach my $k ( 0 .. $ARGV[0] ) {
	$π += 4 * ((-1)**$k) / ( 2*$k + 1 );
	}

printf "%f\n", $π;

And here are the rough timings for a single run on a Mac Pro running OS X.6.4 with a processor dedicated to the just that Perl program, using Perl 5.12.2:

k Built-in math, seconds bigrat, seconds
10 0.016 0.090
100 0.017 0.359
1,000 0.022 83.128
10,000 0.064
100,000 0.479
1,000,000 4.696
10,000,000 47.284

This doesn’t mean that every use of bigrat is going to cause this problem, but even languages that have rationals as first-class data types eventually suffer from creating larger and larger denominators.

And, although it’s outside of the scope of this Item, you can switch out the implementation of rational numbers. Here you used the default implementation, which is Math::BigInt::Calc.

Things to remember

  • You can’t exactly store every number in a limited number of bits.
  • Perl relies on the underlying architecture to store numbers.
  • You can use bigrat to store rational number exactly.
Leave a comment

1 Comments.

  1. Do you know if PDL does have any similar module? Would it be faster?

Leave a Reply


[ Ctrl + Enter ]