Files
wp_xxyyzz 21484fc4bd LazStats: Add tests.
git-svn-id: https://svn.code.sf.net/p/lazarus-ccr/svn@7925 8e941d3f-bd1b-0410-a28a-d453659cc2b4
2020-12-05 19:03:07 +00:00

65 lines
1.4 KiB
ObjectPascal

program Project1;
uses
SysUtils,
spe,
MathUnit;
// Numerical recipes
function GammaLn(x: double): Double;
{ gibt den Log der vollständigen Gamma-Funktion zurück (x>0).
Gamma(x) = integral ( t^(x-1) * exp(-t) dt ) (von 0 bis Unendlich)
(Log, um Floating Point Underflow zu vermeiden). }
const
stp = 2.50662827465;
var
xx,tmp,ser : extended;
begin
if x<=0 then
raise Exception.Create('Argument für GammaLn ist negativ.');
if (x > 1) then begin
xx := x - 1.0;
tmp := xx + 5.5;
tmp := (xx+0.5) * ln(tmp) - tmp;
ser := 1.0 + 76.18009173 /(xx+1.0) - 86.50532033/(xx+2.0)
+ 24.01409822 /(xx+3.0) - 1.231739516/(xx+4.0)
+ 0.120858003E-2/(xx+5.0) - 0.536382E-5/(xx+6.0);
result := tmp + ln(stp*ser);
end else
if (x < 1) then
result := GammaLn(x+1.0) - ln(x)
else
if (x=1) then
result := 0.0;
end;
const
xmin = 1e-3;
xmax = 5.0;
procedure Test;
var
i: Integer;
x, y_lazStats, y_numlib, y_numrecip: Double;
begin
WriteLn('x':20, 'y(lazstats)':20, 'y(numlib)':20, 'y(Num.Recipies)':20);
x := xmin;
while (x <= xmax) do begin
y_numlib := 1.0 - spelga(x);
y_lazstats := MathUnit.gammaln(x);
y_numrecip := gammaln(x);
WriteLn(x:20:6, y_lazstats:20:5, y_numlib:20:5, y_numrecip:20:5);
x := x*1.25
// x := x + dx;
end;
end;
begin
WriteLn('GammaLn function');
WriteLn;
Test;
ReadLn;
end.