150 lines
2.6 KiB
Plaintext
150 lines
2.6 KiB
Plaintext
program Pathological_floating_point_problems;
|
|
|
|
{$APPTYPE CONSOLE}
|
|
|
|
uses
|
|
System.SysUtils,
|
|
Velthuis.BigIntegers,
|
|
Velthuis.BigRationals,
|
|
Velthuis.BigDecimals;
|
|
|
|
type
|
|
TBalance = record
|
|
e, d: BigInteger;
|
|
end;
|
|
|
|
procedure Sequence(Precision: Integer);
|
|
var
|
|
v, v1, c111, c1130, c3000, t2, t3: BigRational;
|
|
n: integer;
|
|
begin
|
|
BigDecimal.DefaultPrecision := Precision;
|
|
v1 := 2;
|
|
v := -4;
|
|
n := 2;
|
|
c111 := BigRational(111);
|
|
c1130 := BigRational(1130);
|
|
c3000 := BigRational(3000);
|
|
|
|
var r :=
|
|
function(): BigRational
|
|
begin
|
|
t3 := v * v1;
|
|
t3 := BigRational.Create(BigDecimal.Divide(c3000, t3));
|
|
t2 := BigRational.Create(BigDecimal.Divide(c1130, v));
|
|
result := c111 - t2;
|
|
result := result + t3;
|
|
end;
|
|
|
|
writeln(' n sequence value');
|
|
|
|
for var x in [3, 4, 5, 6, 7, 8, 20, 30, 50, 100] do
|
|
begin
|
|
while n < x do
|
|
begin
|
|
var tmp := BigRational.Create(v);
|
|
v := BigRational.Create(r());
|
|
v1 := BigRational.Create(tmp);
|
|
inc(n);
|
|
end;
|
|
|
|
var f := double(v);
|
|
writeln(format('%3d %19.16f', [n, f]));
|
|
end;
|
|
end;
|
|
|
|
procedure Bank();
|
|
var
|
|
balance: TBalance;
|
|
m, one, ef, df: BigInteger;
|
|
e, b: BigDecimal;
|
|
begin
|
|
balance.e := 1;
|
|
balance.d := -1;
|
|
one := BigInteger.One;
|
|
|
|
for var y := 1 to 25 do
|
|
begin
|
|
m := y;
|
|
balance.e := m * balance.e;
|
|
balance.d := m * balance.d;
|
|
balance.d := balance.d - one;
|
|
end;
|
|
|
|
e := BigDecimal.Create('2.71828182845904523536028747135');
|
|
|
|
ef := balance.e;
|
|
df := balance.d;
|
|
|
|
b := e * ef;
|
|
b := b + df;
|
|
|
|
writeln(format('Bank balance after 25 years: $%.2f', [b.AsDouble]));
|
|
end;
|
|
|
|
function f(a, b: Double): Double;
|
|
var
|
|
a1, a2, b1, b2, b4, b6, b8, two, t1, t21, t2, t3, t4, t23: BigDecimal;
|
|
begin
|
|
var fp :=
|
|
function(x: double): BigDecimal
|
|
begin
|
|
Result := BigDecimal.Create(x);
|
|
end;
|
|
|
|
a1 := fp(a);
|
|
b1 := fp(b);
|
|
a2 := a1 * a1;
|
|
b2 := b1 * b1;
|
|
b4 := b2 * b2;
|
|
b6 := b2 * b4;
|
|
b8 := b4 * b4;
|
|
|
|
two := fp(2);
|
|
t1 := fp(333.75);
|
|
t1 := t1 * b6;
|
|
t21 := fp(11);
|
|
t21 := t21 * a2 * b2;
|
|
t23 := fp(121);
|
|
t23 := t23 * b4;
|
|
|
|
t2 := t21 - b6;
|
|
t2 := (t2 - t23) - two;
|
|
t2 := a2 * t2;
|
|
|
|
t3 := fp(5.5);
|
|
t3 := t3 * b8;
|
|
|
|
t4 := two * b1;
|
|
t4 := a1 / t4;
|
|
|
|
var s := t1 + t2;
|
|
s := s + t3;
|
|
s := s + t4;
|
|
result := s.AsDouble;
|
|
end;
|
|
|
|
procedure Rump();
|
|
var
|
|
a, b: Double;
|
|
begin
|
|
a := 77617;
|
|
b := 33096;
|
|
writeln(format('Rump f(%g, %g): %g', [a, b, f(a, b)]));
|
|
end;
|
|
|
|
{ TBigRationalHelper }
|
|
|
|
|
|
begin
|
|
for var Precision in [100, 200] do
|
|
begin
|
|
writeln(#10'Precision: ', Precision);
|
|
sequence(Precision);
|
|
bank();
|
|
rump();
|
|
end;
|
|
|
|
readln;
|
|
end.
|