Pathological floating point problems

Most programmers are familiar with the inexactness of floating point calculations in a binary processor. The classic example being:

0.1 + 0.2 =  0.30000000000000004

In many situations the amount of error in such calculations is very small and can be overlooked or eliminated with rounding.

There are pathological problems however, where seemingly simple, straight-forward calculations are extremely sensitive to even tiny amounts of imprecision.

Show how your language deals with such classes of problems.

A sequence that seems to converge to a wrong limit. Consider the sequence:

   v₁ = 2
   v₂ = -4
   v_n = 111 - 1130/v_n-1 + 3000/(v_n-1 * v_n-2)

As n grows larger, the series should converge to 6 but small amounts of error will cause it to approach 100.

Display the values of the sequence where n = 3, 4, 5, 6, 7, 8, 20, 30, 50 & 100 to at least 16 decimal places.

   n = 3   18.5
   n = 4   9.378378
   n = 5   7.801153
   n = 6   7.154414
   n = 7   6.806785
   n = 8   6.5926328
   n = 20  6.0435521101892689
   n = 30  6.006786093031205758530554
   n = 50  6.0001758466271871889456140207471954695237
   n = 100 6.000000019319477929104086803403585715024350675436952458072592750856521767230266

The Chaotic Bank Society is offering a new investment account to their customers. You first deposit $e-1 where e is 2.7182818... the base of natural logarithms.

After each year, your account balance will be multiplied by the number of years that have passed, and $1 in service charges will be removed. So after 1 year, your balance will be multiplied by 1 and $1 will be removed for service charges, after 2 years your balance will be doubled and $1 removed, after 3 years your balance will be tripled and $1 removed, after 10 years, multiplied by 10 and $1 removed and so on. What will your balance be after 25 years?

   Starting balance: $e-1
   Balance = (Balance * year) - 1 for 25 years
   Balance after 25 years: $0.0399387296732302

Siegfried Rump's example. Consider the following function, designed by Siegfried Rump in 1988.

   f(a,b) = 333.75b⁶ + a²( 11a²b² - b⁶ - 121b⁴ - 2) + 5.5b⁸ + a/2b
   compute f(a,b) where a = 77617.0 and b = 33096.0.
   f(77617.0, 33096.0) = -0.827396059946821

Demonstrate how to solve at least one of the first two problems, or both, and the third if you're feeling particularly jaunty.

See Floating-Point Arithmetic Section 1.3.2 Difficult problems.

Perl 6

Works with: Rakudo version 2016-01

The simple solution to doing calculations where floating point numbers exhibit pathological behavior is: don't do floating point calculations. :-) Perl 6 is just as susceptible to floating point error as any other C based language, however, it offers built-in rational Types; where numbers are represented as a ratio of two integers. For normal precision it uses Rats - accurate to 1/2^64, and for arbitrary precision, FatRats, which can grow as large as available memory. Rats don't require any special special setup to use. Any decimal number is automatically stored as a Rat. FatRats require explicit coercion and are "sticky". Any FatRat operand in a calculation will cause all further results to be stored as FatRats. <lang perl6>say '1st: Convergent series'; my @series = 2.FatRat, -4, { 111 - 1130 / $^v + 3000 / ( $^v * $^u ) } ... *; for flat 3..8, 20, 30, 50, 100 -> $n {say "n = {$n.fmt("%3d")} @series[$n-1]"};

say "\n2nd: Chaotic banking society"; sub postfix:<!> (Int $n) { [*] 2..$n } # factorial operator my $years = 25; my $balance = sum map { 1 / FatRat.new($_!) }, 1 .. $years + 15; # Generate e-1 to sufficient precision with a Taylor series put "Starting balance, \$(e-1): \$$balance"; for 1..$years -> $i { $balance = $i * $balance - 1 } printf("After year %d, you will have \$%1.16g in your account.\n", $years, $balance);

print "\n3rd: Rump's example: f(77617.0, 33096.0) = "; sub f (\a, \b) { 333.75*b⁶ + a²*( 11*a²*b² - b⁶ - 121*b⁴ - 2) + 5.5*b⁸ + a/(2*b) } say f(77617.0, 33096.0).fmt("%0.16g");</lang>

Output:

1st: Convergent series
n =   3 18.5
n =   4 9.378378
n =   5 7.801153
n =   6 7.154414
n =   7 6.806785
n =   8 6.5926328
n =  20 6.0435521101892689
n =  30 6.006786093031205758530554
n =  50 6.0001758466271871889456140207471954695237
n = 100 6.000000019319477929104086803403585715024350675436952458072592750856521767230266

2nd: Chaotic banking society
Starting balance, $(e-1): $1.7182818284590452353602874713526624977572470936999
After year 25, you will have $0.0399387296732302 in your account.

3rd: Rump's example: f(77617.0, 33096.0) = -0.827396059946821