home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
C!T ROM 2
/
ctrom_ii_b.zip
/
ctrom_ii_b
/
PROGRAM
/
PASCAL
/
TPL60N19
/
TESTPRGS
/
UNIT1.PAS
< prev
next >
Wrap
Pascal/Delphi Source File
|
1992-02-28
|
23KB
|
701 lines
{$a+,n-,x-,s-,i-,r-,b-,v-}
unit Unit1;
interface
uses mainvars;
procedure start;
procedure mile2060;
implementation
procedure start;
begin (* PARA *)
{First two assignments use integer right-hand sides.}
Zero := 0;
One := 1;
Two := One + One;
Three := Two + One;
Four := Three + One;
Five := Four + One;
Eight := Four + Four;
Nine := Three * Three;
TwentySeven := Nine * Three;
ThirtyTwo := Four * Eight;
TwoForty := Four * Five * Three * Four;
MinusOne := -One;
Half := One / Two;
OneAndHalf := One + Half;
NoErrors [Failure] := 0;
NoErrors [SeriousDefect] := 0;
NoErrors [Defect] := 0;
NoErrors [Flaw] := 0;
PageNo := 0;
{=============================================}
Milestone := 0;
{=============================================}
writeln ('Type any character to start the program.');
{ assign(input,'con:');} { for TURBO Pascal version 2 }
{ reset (input); } { for old Cray Pascal }
while not eoln (input) do
read (input, ch);
Instructions;
Pause;
Heading;
Pause;
Characteristics;
Pause;
History;
{=============================================}
Milestone := 7;
{=============================================}
Pause;
writeln ('Program is now RUNNING tests on small integers:');
TestCondition (Failure, (Zero + Zero = Zero) and (One - One = Zero)
and (One > Zero)
and (One + One = Two), ' 0+0<>0 or 1-1<>0 or 1<=0 or 1+1<>2 '
);
Z := - Zero;
if Z <> 0.0 then
begin
NoErrors [Failure] := NoErrors [Failure] + 1;
writeln ('Comparison alleges that -0.0 is Non-zero!');
U2 := 0.001;
Radix := 1;
TestPartialUnderflow;
end;
TestCondition (Failure, (Three = Two + One) and (Four = Three + One)
and (Four + Two * (- Two) = Zero)
and (Four - Three - One = Zero),
' 3<>2+1, 4<>3+1, 4+2*(-2)<>0 or 4-3-1<>0');
TestCondition (Failure, (MinusOne = - One)
and (MinusOne + One = Zero ) and (One + MinusOne = Zero)
and (MinusOne + abs (One) = Zero)
and (MinusOne + MinusOne * MinusOne = Zero),
'-1+1<>0, -1+abs(1)<>0 or -1+(-1)*(-1)<>0');
TestCondition (Failure, Half + MinusOne + Half = Zero,
' 1/2 + (-1) + 1/2 <> 0 ');
{=============================================}
Milestone := 10;
{=============================================}
TestCondition (Failure, (Nine = Three * Three)
and (TwentySeven = Nine * Three) and (Eight = Four + Four)
and (ThirtyTwo = Eight * Four)
and (ThirtyTwo - TwentySeven - Four - One = Zero),
'9<>3*3, 27<>9*3, 32<>8*4 or 32-27-4-1<>0');
TestCondition (Failure, (Five = Four + One)
and (TwoForty = Four * Five * Three * Four)
and (TwoForty / Three - Four * Four * Five = Zero)
and ( TwoForty / Four - Five * Three * Four = Zero)
and ( TwoForty / Five - Four * Three * Four = Zero),
'5<>4+1,240/3<>80,240/4<>60, or 240/5<>48');
if NoErrors [Failure] = 0 then
begin
writeln (' -1, 0, 1/2, 1, 2, 3, 4, 5, 9, 27, 32 & 240 are O.K.');
writeln
end;
writeln ('Searching for Radix and Precision.');
W := One;
repeat
W := W + W;
Y := W + One;
Z := Y - W;
Y := Z - One;
until (MinusOne + abs (Y) >= Zero);
{.. now W is just big enough that |((W+1)-W)-1| >= 1 ...}
Precision := 0;
Y := One;
repeat
Radix := W + Y;
Y := Y + Y;
Radix := Radix - W;
until (Radix <> Zero);
if Radix < Two then
Radix := One;
writeln ('Radix = ', Radix);
if Radix <> 1 then
begin
W := One;
repeat
Precision := Precision + One;
W := W * Radix;
Y := W + One;
until (Y - W) <> One;
{... now W = Radix^Precision is barely too big to satisfy (W+1)-W = 1
...}
end;
U1 := One / W;
U2 := Radix * U1;
writeln ('Closest relative separation found is U1 = ', U1);
writeln;
writeln ('Recalculating radix and precision');
E0 := Radix;
E1 := U1;
E9 := U2;
{save old values}
X := Four / Three;
Third := X - One;
F6 := Half - Third;
X := F6 + F6;
X := abs (X - Third);
if X < U2 then
X := U2;
{... now X = (unknown no.) ulps of 1+...}
repeat
U2 := X;
Y := Half * U2 + ThirtyTwo * U2 * U2;
Y := One + Y;
X := Y - One;
until (U2 <= X) or (X <= Zero);
{... now U2 = 1 ulp of 1 + ... }
X := Two / Three;
F6 := X - Half;
Third := F6 + F6;
X := Third - Half;
X := abs (X + F6);
if X < U1 then
X := U1;
{... now X = (unknown no.) ulps of 1 -... }
repeat
U1 := X;
Y := Half * U1 + ThirtyTwo * U1 * U1;
Y := Half - Y;
X := Half + Y;
Y := Half - X;
X := Half + Y;
until (U1 <= X) or (X <= Zero);
{... now U1 = 1 ulp of 1 - ... }
if U1 = E1 then
writeln (' confirms closest relative separation U1 .')
else
writeln (' gets better closest relative separation U1 = ', U1);
W := One / U1;
F9 := (Half - U1) + Half;
Radix := Int (0.01 + U2 / U1);
if Radix = E0 then
writeln ('Radix confirmed.')
else
writeln ('MYSTERY: recalculated Radix = ', Radix);
TestCondition (Defect, Radix <= Eight + Eight,
'Radix is too big: roundoff problems ');
TestCondition (Flaw, (Radix = Two) or (Radix = 10)
or (Radix = One), 'Radix is not as good as 2 or 10. ');
end (*start*);
procedure mile2060;
begin
{=============================================}
Milestone := 20;
{=============================================}
TestCondition (Failure, F9 - Half < Half,
' (1-U1)-1/2 < 1/2 is FALSE, prog. fails?');
X := F9;
I := 1;
Y := X - Half;
Z := Y - Half;
TestCondition (Failure, (X <> One)
or (Z = Zero), 'Comparison is fuzzy,X=1 but X-1/2-1/2<>1');
X := One + U2;
I := 0;
{=============================================}
Milestone := 25;
{=============================================}
BMinusU2 := Radix - One;
BMinusU2 := (BMinusU2 - U2) + One;
if Radix <> One then
begin {... BMinusU2 = nextafter(Radix, 0) }
X := - TwoForty * ln (U1) / ln (Radix);
Y := Int (Half + X);
if abs (X - Y) * Four < One then
X := Y;
Precision := X / TwoForty;
Y := Int (Half + Precision);
if abs (Precision - Y) * TwoForty < Half then
Precision := Y;
{ Purify integers }
end;
if (Precision <> Int (Precision)) or (Radix = One) then
begin
writeln ('Precision cannot be characterized by an integer',
' number of sig. digits,');
writeln ('but, by itself, this is a minor flaw.');
end;
if Radix = One then
writeln ('logarithmic encoding has precision characterized',
'solely by U1.')
else
writeln ('The number of significant digits of the Radix is ',
Precision);
TestCondition (SeriousDefect, U2 * Nine * Nine * TwoForty < One,
' Precision worse than 5 decimal figures ');
{=============================================}
Milestone := 30;
{=============================================}
{ Test for extra-precise subepressions }
X := abs (((Four / Three - One) - One / Four) * Three - One / Four);
repeat
Z2 := X;
X := (One + (Half * Z2 + ThirtyTwo * Z2 * Z2)) - One;
until (Z2 <= X) or (X <= Zero);
Y := abs ((Three / Four - Two / Three) * Three - One / Four);
Z := Y;
X := Y;
repeat
Z1 := Z;
Z := (One / Two - ((One / Two - (Half * Z1 + ThirtyTwo * Z1 * Z1))
+ One / Two)) + One / Two;
until (Z1 <= Z) or (Z <= Zero);
repeat
repeat
Y1 := Y;
Y := (Half - ((Half - (Half * Y1 + ThirtyTwo * Y1 * Y1)) + Half
)) + Half;
until (Y1 <= Y) or (Y <= Zero);
X1 := X;
X := ((Half * X1 + ThirtyTwo * X1 * X1) - F9) + F9;
until (X1 <= X) or (X <= Zero);
if (X1 <> Y1) or (X1 <> Z1) then
begin
NoErrors [SeriousDefect] := NoErrors [SeriousDefect] + 1;
writeln ('SERIOUS DEFECT: Disagreements among the values X1, Y1, Z1');
writeln ('resp. ', X1, Y1, Z1);
writeln ('are symptoms of inconsistencies introduced');
writeln ('by extra-precise evaluation of allegedly');
writeln ('"optimized" arithmetic subexpressions.');
writeln ('Possibly some part of this test is inconsistent.');
if (X1 = U1) or (Y1 = U1) or (Z1 = U1) then
writeln ('That feature is not tested further by this program.');
end
else if (Z1 <> U1) or (Z2 <> U2) then
begin
if (Z1 >= U1) or (Z2 >= U2) then
begin
NoErrors [Failure] := NoErrors [Failure] + 1;
writeln ('FAILURE: Precision ', Precision);
writeln ('U1 = ', U1, ' Z1 - U1 = ', Z1 - U1);
writeln ('U2 = ', U2, ' Z2 - U2 = ', Z2 - U2);
end
else begin
if (Z1 <= Zero) or (Z2 <= Zero) then begin
writeln ('Because of unusual Radix = ', Radix);
writeln (' or exact rational arithmetic a result');
writeln (' Z1 = ', Z1, ' or Z2 = ', Z2);
writeln (' of an extra precision test is inconsistent.');
if Z1 = Z2 then
end;
if (Z1 <> Z2) or (Z1 > Zero) then begin
X := Z1 / U1;
Y := Z2 / U2;
if Y > X then X := Y;
Q := - ln (X);
writeln ('Some subexpressions appear to be calculated');
writeln ('extra precisely with about ');
writeln (Q / ln (Radix), 'extra B-digits i.e. ');
writeln ('roughly ', Q / ln (10),
' extra significant decimals.');
end;
writeln ('That feature is not tested further by this program.')
end
end;
Pause;
{=============================================}
Milestone := 35;
{=============================================}
if Radix >= Two then
begin
X := W / (Radix * Radix);
Y := X + One;
Z := Y - X;
T := Z + U2;
X := T - Z;
TestCondition (Failure, X = U2,
'Subtraction is not normlzd X=Y,X+Z<>Y+Z!');
if X = U2 then
writeln ('Subtraction appears to be normalized,',
' as it should be.');
end;
writeln;
writeln ('Checking for guard digit on *, /, and -.');
Y := F9 * One;
Z := One * F9;
X := F9 - Half;
Y := (Y - Half) - X;
Z := (Z - Half) - X;
X := One + U2;
T := X * Radix;
R := Radix * X;
X := T - Radix;
X := X - Radix * U2;
T := R - Radix;
T := T - Radix * U2;
X := X * (Radix - One);
T := T * (Radix - One);
if (X = Zero) and (Y = Zero) and (Z = Zero) and (T = Zero) then
GMult := Yes
else
begin
GMult := No;
TestCondition (SeriousDefect, false,
' * lacks guard digit, 1*X <> X ');
end;
Z := Radix * U2;
X := One + Z;
Y := abs ((X + Z) - X * X) - U2;
X := One - U2;
Z := abs ((X - U2) - X * X) - U1;
TestCondition (Failure, (Y <= Zero)
and (Z <= Zero), ' * gets too many final digits wrong. ');
Y := One - U2;
X := One + U2;
Z := One / Y;
Y := Z - X;
X := One / Three;
Z := Three / Nine;
X := X - Z;
T := Nine / TwentySeven;
Z := Z - T;
TestCondition (Defect, (X = Zero) and (Y = Zero)
and (Z = Zero), 'Division error > ulp, 1/3 <> 3/9 <> 9/27');
Y := F9 / One;
X := F9 - Half;
Y := (Y - Half) - X;
X := One + U2;
T := X / One;
X := T - X;
if (X = Zero) and (Y = Zero) and (Z = Zero) then
GDiv := Yes
else
begin
GDiv := No;
TestCondition (SeriousDefect, false,
' Division lacks guard digit so X/1 <> X');
end;
X := One / (One + U2);
Y := X - Half - Half;
TestCondition (SeriousDefect, Y < Zero,
' Computed value of 1/1.000..1 >= 1. ');
X := One - U2;
Y := One + Radix * U2;
Z := X * Radix;
T := Y * Radix;
R := Z / Radix;
StickyBit := T / Radix;
X := R - X;
Y := StickyBit - Y;
TestCondition (Failure, (X = Zero) and (Y = Zero),
' * &or / gets too many last digits wrong');
Y := One - U1;
X := One - F9;
Y := One - Y;
T := Radix - U2;
Z := Radix - BMinusU2;
T := Radix - T;
if (X = U1) and (Y = U1) and (Z = U2) and (T = U2) then
GAddSub := Yes
else
begin
GAddSub := No;
TestCondition (SeriousDefect, false,
'- lacks guard dig.,cancellation obscured');
end;
if (F9 <> One) and (F9 - One >= Zero) then begin
TestCondition (SeriousDefect, false,
'comparison alleges (1-U1) < 1 although');
writeln(' subtration yields (1-U1) - 1 = 0 , thereby vitiating');
writeln(' such precautions against division by zero as');
writeln(' ... if (X=1.0) then ..... else .../(X-1.0)...');
end;
if (GMult = Yes) and (GDiv = Yes) and (GAddSub = Yes) then
writeln (' *, /, and - have guard digits, as they should.');
{=============================================}
Milestone := 40;
{=============================================}
Pause;
writeln ('Checking rounding on multiply, divide and add/subtract.');
RMult := Other;
RDiv := Other;
RAddSub := Other;
RadixD2 := Radix / Two;
A1 := Two;
Done := false;
repeat
AInverse := Radix;
repeat
X := AInverse;
AInverse := AInverse / A1;
until Int (AInverse) <> AInverse;
Done := (X = One) or (A1 > Three);
if not Done then
A1 := Nine + One;
until Done;
if X = One then
A1 := Radix;
AInverse := One / A1;
X := A1;
Y := AInverse;
Done := false;
repeat
Z := X * Y - Half;
TestCondition (Failure, Z = Half,
' X * (1/X) differs from 1. ');
Done := X = Radix;
X := Radix;
Y := One / X;
until Done;
Y2 := One + U2;
Y1 := One - U2;
X := OneAndHalf - U2;
Y := OneAndHalf + U2;
Z := (X - U2) * Y2;
T := Y * Y1;
Z := Z - X;
T := T - X;
X := X * Y2;
Y := (Y + U2) * Y1;
X := X - OneAndHalf;
Y := Y - OneAndHalf;
if (X = Zero) and (Y = Zero) and (Z = Zero) and (T <= Zero) then
begin
X := (OneAndHalf + U2) * Y2;
Y := OneAndHalf - U2 - U2;
Z := OneAndHalf + U2 + U2;
T := (OneAndHalf - U2) * Y1;
X := X - (Z + U2);
StickyBit := Y * Y1;
S := Z * Y2;
T := T - Y;
Y := (U2 - Y) + StickyBit;
Z := S - (Z + U2 + U2);
StickyBit := (Y2 + U2) * Y1;
Y1 := Y2 * Y1;
StickyBit := StickyBit - Y2;
Y1 := Y1 - Half;
if (X = Zero) and (Y = Zero) and (Z = Zero) and (T = Zero)
and ( StickyBit = Zero) and (Y1 = Half) then
begin
RMult := Rounded;
writeln ('Multiplication appears to round correctly.');
end
else if (X + U2 = Zero) and (Y < Zero) and (Z + U2 = Zero)
and (T < Zero) and (StickyBit + U2 = Zero)
and (Y1 < Half) then
begin
RMult := Chopped;
writeln ('Multiplication appears to chop.');
end
else
writeln ('* is neither chopped nor correctly rounded.');
if (RMult = Rounded) and (GMult = No) then
notify('multiplication');
end
else
writeln ('* is neither chopped nor correctly rounded.');
{=============================================}
Milestone := 45;
{=============================================}
Y2 := One + U2;
Y1 := One - U2;
Z := OneAndHalf + U2 + U2;
X := Z / Y2;
T := OneAndHalf - U2 - U2;
Y := (T - U2) / Y1;
Z := (Z + U2) / Y2;
X := X - OneAndHalf;
Y := Y - T;
T := T / Y1;
Z := Z - (OneAndHalf + U2);
T := (U2 - OneAndHalf) + T;
if not ((X > Zero) or (Y > Zero) or (Z > Zero) or (T > Zero)) then
begin
X := OneAndHalf / Y2;
Y := OneAndHalf - U2;
Z := OneAndHalf + U2;
X := X - Y;
T := OneAndHalf / Y1;
Y := Y / Y1;
T := T - (Z + U2);
Y := Y - Z;
Z := Z / Y2;
Y1 := (Y2 + U2) / Y2;
Z := Z - OneAndHalf;
Y2 := Y1 - Y2;
Y1 := (F9 - U1) / F9;
if (X = Zero) and (Y = Zero) and (Z = Zero) and (T = Zero)
and (Y2 = Zero) and (Y2 = Zero)
and (Y1 - Half = F9 - Half ) then
begin
RDiv := Rounded;
writeln ('Division appears to round correctly.');
if GDiv = No then notify(' division ');
end
else if (X < Zero) and (Y < Zero) and (Z < Zero) and (T < Zero)
and (Y2 < Zero) and (Y1 - Half < F9 - Half) then
begin
RDiv := Chopped;
writeln ('Division appears to chop.');
end;
end;
if RDiv = Other then
writeln ('/ is neither chopped nor correctly rounded.');
BInverse := One / Radix;
TestCondition (Failure, (BInverse * Radix - Half = Half),
' Radix * ( 1 / Radix ) differs from 1. ');
{=============================================}
Milestone := 50;
{=============================================}
TestCondition (Failure, ((F9 + U1) - Half = Half)
and ((BMinusU2 + U2 ) - One = Radix - One),
'Incomplete carry-propagation in Addition');
X := One - U1 * U1;
Y := One + U2 * (One - U2);
Z := F9 - Half;
X := (X - Half) - Z;
Y := Y - One;
if (X = Zero) and (Y = Zero) then
begin
RAddSub := Chopped;
writeln ('Add/Subtract appears to be chopped.');
end;
if GAddSub = Yes then
begin
X := (Half + U2) * U2;
Y := (Half - U2) * U2;
X := One + X;
Y := One + Y;
X := (One + U2) - X;
Y := One - Y;
if (X = Zero) and (Y = Zero) then
begin
X := (Half + U2) * U1;
Y := (Half - U2) * U1;
X := One - X;
Y := One - Y;
X := F9 - X;
Y := One - Y;
if (X = Zero) and (Y = Zero) then
begin
RAddSub := Rounded;
writeln ('Addition/Subtraction appears to round correctly.');
if GAddSub = No then notify(' add/subtract ');
end
else
writeln ('Addition/Subtraction neither rounds nor chops.');
end
else
writeln ('Addition/Subtraction neither rounds nor chops.');
end
else
writeln ('Addition/Subtraction neither rounds nor chops.');
S := One;
X := One + Half * (One + Half);
Y := (One + U2) * Half;
Z := X - Y;
T := Y - X;
StickyBit := Z + T;
if StickyBit <> 0 then
begin
S := 0;
NoErrors [Flaw] := NoErrors [Flaw] + 1;
write('FLAW: (X - Y) + (Y - X) is non zero!');
end;
StickyBit := Zero;
if (GMult = Yes) and (GDiv = Yes) and (GAddSub = Yes)
and (RMult = Rounded) and (RDiv = Rounded)
and (RAddSub = Rounded) and (Int (RadixD2) = RadixD2) then
begin
writeln (' Checking for sticky bit.');
X := (Half + U1) * U2;
Y := Half * U2;
Z := One + Y;
T := One + X;
if (Z - One <= Zero) and (T - One >= U2) then
begin
Z := T + Y;
Y := Z - X;
if (Z - T >= U2) and (Y - T = Zero) then
begin
X := (Half + U1) * U1;
Y := Half * U1;
Z := One - Y;
T := One - X;
if (Z - One = Zero) and (T - F9 = Zero) then
begin
inline($fa/$fb);
Z := (Half - U1) * U1;
T := F9 - Z;
Q := F9 - Y;
if (T - F9 = Zero) and (F9 - U1 - Q = Zero) then
begin
Z := (One + U2) * OneAndHalf;
T := (OneAndHalf + U2) - Z + U2;
X := One + Half / Radix;
Y := One + Radix * U2;
Z := X * Y;
if (T = Zero) and (X + Radix * U2 - Z = Zero) then
begin
if Radix <> Two then
begin
X := Two + U2;
Y := X / Two;
if (Y - One = Zero) then
StickyBit := S;
end
else StickyBit := S;
end;
end;
end;
end;
end;
end;
if StickyBit = One then
writeln ('Sticky bit apparently used correctly.')
else writeln ('Sticky bit used incorrectly or not at all.');
if (GMult = No) or (GDiv = No) or (GAddSub = No) or (RMult = Other)
or (RDiv = Other) or (RAddSub = Other) then begin
TestCondition (Flaw, false,
'lack(s) of guard digits or failure(s) to');
writeln('correctly round or chop (noted above) count as one');
writeln('flaw in the final tally below.')
end;
{=============================================}
Milestone := 60;
{=============================================}
writeln;
writeln ('Does Multiplication commute? Testing on ', NoTrials,
' random pairs.');
R9 := sqrt (3.0);
RandomNumber1 := Third;
I := 1;
repeat
X := Random;
Y := Random;
Z9 := Y * X;
Z := X * Y;
Z9 := Z - Z9;
I := I + 1;
until (I > NoTrials) or (Z9 <> Zero);
if I = NoTrials then
begin
RandomNumber1 := One + Half / Three;
RandomNumber2 := (U2 + U1) + One;
Z := RandomNumber1 * RandomNumber2;
Y := RandomNumber2 * RandomNumber1;
Z9 := (One + Half / Three) * ((U2 + U1) + One) - (One + Half /
Three) * ((U2 + U1) + One);
end;
if not ((I = NoTrials) or (Z9 = Zero)) then
begin
NoErrors [Defect] := NoErrors [Defect] + 1;
writeln ('DEFECT: X * Y = Y * X trail fails.');
end
else
writeln ('No failures found in ', NoTrials, ' integer pairs.');
end;
end.