Tue Feb 10 12:26:19 2004 run on Linux
% Tests of REDUCE Arithmetic.
% Authors: Anthony C. Hearn and Stanley L. Kameny.
% Copyright (c) 1987, 1988, 1989, 1991, Stanley L. Kameny.
% Copyright (c) 1998, Anthony C. Hearn.
% All Rights Reserved.
% This test file is a combination of three tests files from versions
% 3.6 or earlier: math, rounded and complex.
% Simple tests of rounded arithmetic.
% Tests in the exact mode.
x := 1/2;
1
x := ---
2
y := x + 0.7;
6
y := ---
5
% Tests in approximate mode.
on rounded;
y;
1.2
% As expected not converted to approximate form.
z := y+1.2;
z := 2.4
z/3;
0.8
% Let's raise this to a high power.
ws^24;
0.00472236648287
% Now a high exponent value.
% 10.2^821;
% Elementary function evaluation.
cos(pi);
- 1
symbolic ws;
(!*sq ((!:rd!: . -1.0) . 1) t)
z := sin(pi);
z := 1.22464679915e-16
symbolic ws;
(!*sq ((!:rd!: . 1.2246467991474e-16) . 1) t)
% Handling very small quantities.
% With normal defaults, underflows are converted to 0.
exp(-100000.1**2);
0
% However, if you really want that small number, roundbf can be used.
on roundbf;
exp(-100000.1**2);
1.18441281937e-4342953505
off roundbf;
% Now let us evaluate pi.
pi;
3.14159265359
% Let us try a higher precision.
prec0 := precision 50;
prec0 := 12
pi;
3.1415926535897932384626433832795028841971693993751
% Now find the cosine of pi/6.
cos(ws/6);
0.86602540378443864676372317075293618347140262690519
% This should be the sqrt(3)/2.
ws**2;
0.75
% Here are some well known examples which show the power of this system.
precision 10;
50
% This should give the usual default again.
let xx=e**(pi*sqrt(163));
let yy=1-2*cos((6*log(2)+log(10005))/sqrt(163));
% First notice that xx looks like an integer.
xx;
2.625374126e+17
% And that yy looks like zero.
yy;
0
% But of course it's an illusion.
precision 50;
10
xx;
2.6253741264076874399999999999925007259719818568888e+17
yy;
- 1.2815256559456092775159749532170513334408547400481e-16
% Now let's look at an unusual way of finding an old friend.
procedure agm;
begin scalar a,b,u,x,y,p,pn;
a := 1; b := 1/sqrt 2; u:= 1/4; x := 1$ pn := 4;
repeat <<p := pn;
y := a; a := (a+b)/2; b := sqrt(y*b); % Arith-geom mean.
u := u-x*(a-y)**2; x := 2*x; pn := a**2/u;
write "pn=",pn>>
until pn>=p;
return p
end;
agm
agm();
pn=3.1876726427121086272019299705253692326510535718594
pn=3.1416802932976532939180704245600093827957194388154
pn=3.1415926538954464960029147588180434861088792372613
pn=3.1415926535897932384663606027066313217577024113424
pn=3.1415926535897932384626433832795028841971699491647
pn=3.1415926535897932384626433832795028841971693993751
pn=3.1415926535897932384626433832795028841971693993751
3.1415926535897932384626433832795028841971693993751
% The limit is obviously.
pi;
3.1415926535897932384626433832795028841971693993751
off rounded;
clear x;
precision prec0;
50
% Tests of Complex arithmetic.
on complex;
(31+i)/74;
31 + i
--------
74
ws/(b+1);
31 + i
------------
74*(b + 1)
% This now comes out right!
w:=(x+3*i)**2;
2
w := x + 6*i*x - 9
on gcd;
(x**3-7*x**2+x-7)/(x**2+(3+i)*x+3*i);
2
x - (7 + i)*x + 7*i
----------------------
x + 3
off gcd;
sqrt(x**4+14*i*x**3-51*x**2-14*i*x+1);
2
x + 7*i*x - 1
% All rounded tests are done twice: first, they are done at the default
% precision, in which all rounded operations use standard floating point
% logic. Then precision is increased, causing all rounded operations to
% use extended precision bigfloat arithmetic. This is necessary to
% exercise and test the bigfloat-based arithmetic functions.
prec0 := precision 0;
prec0 := 12
% To determine the nominal default precision.
% Tests using default precision.
on rounded;
*** Domain mode complex changed to complex-rounded
(3.25 + 8.5i) + (6.75 - 8.5i);
10.0
(3.25 + 8.5i) - (6.0 - 9.5i);
- 2.75 + 18.0*i
(1.0 + 10.0*i)*(-6.5 + 2.5*i);
- 31.5 - 62.5*i
(1.2 - 3.4*i)*(-5.6 + 7.8*i);
19.8 + 28.4*i
(19.8 + 28.4*i)/(-5.6 + 7.8*i);
1.2 - 3.4*i
e;
2.71828182846
pi;
3.14159265359
17*i**2;
-17
(-7.0 + 24.0*i)**(1/2);
3.0 + 4.0*i
sqrt(-7.0 + 24.0*i);
3.0 + 4.0*i
sqrt(-10.12 - 8.16*i);
1.2 - 3.4*i
sin(0.0 + 0.0*i);
0
sin(1.0 + 0.0*i);
0.841470984808
sin(1.0 + 1.0*i);
1.29845758142 + 0.634963914785*i
cos(0.0 + 0.0*i);
1
cos(1.0 - 0.0*i);
0.540302305868
cos(1.0 + 1.0*i);
0.833730025131 - 0.988897705763*i
tan(0.0 + 0.0*i);
0
tan(1.0 + 0.0*i);
1.55740772465
tan(1.0 + 1.0*i);
0.27175258532 + 1.08392332734*i
asin(1.0 + 1.0*i);
0.666239432493 + 1.06127506191*i
acos(1.0 + 1.0*i);
0.904556894302 - 1.06127506191*i
atan(1.0 + 1.0*i);
1.0172219679 + 0.402359478109*i
log(1.0 + 1.0*i);
0.34657359028 + 0.785398163397*i
asin 2;
1.57079632679 - 1.31695789692*i
sin ws;
2.0 - 1.06057523872e-16*i
acos 2;
1.31695789692*i
cos ws;
2.0
atan(1+i);
1.0172219679 + 0.402359478109*i
tan ws;
1 + i
log(2+i);
0.804718956217 + 0.463647609001*i
exp ws;
2.0 + i
e**(i*pi);
- 1 + 1.22464679915e-16*i
e**i;
0.540302305868 + 0.841470984808*i
z := sqrt i;
z := 0.707106781187 + 0.707106781187*i
z**2;
i
off rounded;
*** Domain mode complex-rounded changed to complex
%-----------------end of normal floating point tests--------------------
precision(prec0+6);
12
% Arbitrary precision increase -> bigfloat functions
%----------------------start of bigfloat tests--------------------------
on rounded;
*** Domain mode complex changed to complex-rounded
(3.25 + 8.5i) + (6.75 - 8.5i);
10.0
(3.25 + 8.5i) - (6.0 - 9.5i);
- 2.75 + 18.0*i
(1.0 + 10.0*i)*(-6.5 + 2.5*i);
- 31.5 - 62.5*i
(1.2 - 3.4*i)*(-5.6 + 7.8*i);
19.8 + 28.4*i
(19.8 + 28.4*i)/(-5.6 + 7.8*i);
1.2 - 3.4*i
e;
2.71828182845904524
pi;
3.14159265358979324
17*i**2;
-17
(-7.0 + 24.0*i)**(1/2);
3.0 + 4.0*i
sqrt(-7.0 + 24.0*i);
3.0 + 4.0*i
sqrt(-10.12 - 8.16*i);
1.2 - 3.4*i
sin(0.0 + 0.0*i);
0
sin(1.0 + 0.0*i);
0.841470984807896507
sin(1.0 + 1.0*i);
1.29845758141597729 + 0.634963914784736108*i
cos(0.0 + 0.0*i);
1
cos(1.0 - 0.0*i);
0.540302305868139717
cos(1.0 + 1.0*i);
0.833730025131149049 - 0.988897705762865096*i
tan(0.0 + 0.0*i);
0
tan(1.0 + 0.0*i);
1.55740772465490223
tan(1.0 + 1.0*i);
0.271752585319511717 + 1.08392332733869454*i
asin(1.0 + 1.0*i);
0.666239432492515255 + 1.06127506190503565*i
acos(1.0 + 1.0*i);
0.904556894302381364 - 1.06127506190503565*i
atan(1.0 + 1.0*i);
1.01722196789785137 + 0.402359478108525094*i
log(1.0 + 1.0*i);
0.346573590279972655 + 0.78539816339744831*i
asin 2;
1.57079632679489662 - 1.31695789692481671*i
sin ws;
2.0
acos 2;
1.31695789692481671*i
cos ws;
2.0
atan(1+i);
1.01722196789785137 + 0.402359478108525094*i
tan ws;
1 + i
log(2+i);
0.804718956217050187 + 0.463647609000806116*i
exp ws;
2.0 + i
e**(i*pi);
- 1
e**i;
0.540302305868139717 + 0.841470984807896507*i
z := sqrt i;
z := 0.707106781186547524 + 0.707106781186547524*i
z**2;
i
off rounded;
*** Domain mode complex-rounded changed to complex
% ---------------------------------------------------------------------
% The following examples are independent of precision.
precision prec0;
18
% Restores default precision.
s:= 1.1+2.3i;
11 + 23*i
s := -----------
10
s/4;
11 + 23*i
-----------
40
% This would have had a common factor of 4.
x:= a+1.1+2.3i;
10*a + 11 + 23*i
x := ------------------
10
y:= b+1.2+1.3i;
10*b + 12 + 13*i
y := ------------------
10
z:= x/y;
100*a*b + (120 - 130*i)*a + (110 + 230*i)*b + 431 + 133*i
z := -----------------------------------------------------------
2
100*b + 240*b + 313
z/4;
100*a*b + (120 - 130*i)*a + (110 + 230*i)*b + 431 + 133*i
-----------------------------------------------------------
2
4*(100*b + 240*b + 313)
% This would have had a common polynomial factor b^2 + ...
z*7/4;
7*(100*a*b + (120 - 130*i)*a + (110 + 230*i)*b + 431 + 133*i)
---------------------------------------------------------------
2
4*(100*b + 240*b + 313)
s/(c^2+c+1);
11 + 23*i
-----------------
2
10*(c + c + 1)
% This would have had a common factor of c^2+c+1,
clear x;
zz:= x^2+(1.1+2.3i)*x+1.2+1.3i;
2
10*x + (11 + 23*i)*x + 12 + 13*i
zz := -----------------------------------
10
ss:=1.23456789x^2+1.3579i*x+5.6789;
2
123456789*x + 135790000*i*x + 567890000
ss := ------------------------------------------
100000000
z:= x+1.1+2.3i;
10*x + 11 + 23*i
z := ------------------
10
on rationalize;
z;
10*x + 11 + 23*i
------------------
10
% Same as previous answer.
off rationalize;
1.23456789x^2+2.3456i*x+7.89;
2
123456789*x + 234560000*i*x + 789000000
------------------------------------------
100000000
on factor;
x**2+1;
(x + i)*(x - i)
x**4-1;
(x + i)*(x - i)*(x + 1)*(x - 1)
x**4+(i+2)*x**3+(2*i+5)*x**2+(2*i+6)*x+6;
2
(x + i*x + 3)*(x + 1 + i)*(x + 1 - i)
(2*i+3)*x**4+(3*i-2)*x**3-2*(i+1)*x**2+i*x-1;
2 2
i*((2 - 3*i)*x - i)*(x + i*x - 1)
% Multivariate examples.
x**2+y**2;
(10*b + 10*i*x + 12 + 13*i)*(10*b - 10*i*x + 12 + 13*i)
---------------------------------------------------------
100
off factor;
factorize(x**2+1);
{{x + i,1},{x - i,1}}
off complex;
% Tests of some elementary functions.
comment Integer functions that work in all domain modes, independent of
switch NUMVAL, so long as their arguments evaluate to real numbers.
Functions of one argument:
FIX, SGN, ROUND, CEILING, FLOOR
(The following functions are available only in symbolic mode, so they
are not tested here: ISQRT, ICBRT, ILOG2, IROOTN);
fix a;
fix(a)
% Will be evaluated only if a evaluates to a real number.
a := 27/4;
27
a := ----
4
fix a;
6
fix 12.345;
12
sign (-15/2);
-1
round 12.5;
13
ceiling 12.5;
13
floor 12.5;
12
% isqrt 12.5;
% icbrt 12.5;
% ilog2 130.7;
% irootn(72,4);
% irootn(72,3/2); % This will not evaluate.
comment Functions which require arguments which evaluate to integers.
Function of one argument: FACTORIAL
Fumction of two arguments: PERM, CHOOSE;
$
factorial 10;
3628800
perm(5,10);
30240
% Permutations of 5 out of 10.
choose(5,10);
252
% Choose 5 out of 10.
comment
These functions are evaluated in dmodes ROUNDED and COMPLEX-ROUNDED
(ON ROUNDED,COMPLEX) so long as their arguments and values evaluate
to real numbers and NUMVAL (normally ON) is ON.
Variable treated as function of no arguments: E, PI.
Functions of one argument:
EXP, LOG, LN, LOG10, NORM, ARG, SQRT,
RAD2DEG, RAD2DMS, DEG2RAD, DEG2DMS, DMS2DEG, DMS2RAD,
SIN, ASIN, COS, ACOS, TAN, ATAN, COT, ACOT, SEC, ASEC, CSC, ACSC,
SINH, ASINH, COSH, ACOSH, TANH, ATANH, COTH, ACOTH, SECH, ASECH,
CSCH, ACSCH.
Functions of two arguments:
EXPT, LOGB, HYPOT, ATAN2.
Function evaluation is carried out to the precision specified in the
latest PRECISION statement.
(The following functions are available only in symbolic mode, so they
are not tested here:
SIND, ASIND, COSD, ACOSD, TAND, ATAND, COTD, ACOTD, SECD, ASECD,
CSCD, ACSCD, ATAN2D, CBRT);
on rounded;
precision 6;
12
a := exp 3;
a := 20.0855
log a;
3.0
ln a;
3.0
log10 1000;
3.0
norm (-12.345);
12.345
% For real x, this is equivalent to ABS x.
arg (-12.345);
3.14159
% For real x, this -> if x<0 then pi else 0.0.
sqrt 3;
1.73205
ws**2;
3.0
deg2rad 30;
0.523599
rad2deg ws;
30.0
a := deg2dms 12.345;
a := {12,20,42.0}
% a will be a list.
dms2deg ws;
12.345
dms2rad a;
0.215461
rad2deg ws;
12.345
asin 0.5;
0.523599
sin ws;
0.5
acos 0.5;
1.0472
cos ws;
0.5
atan 0.5;
0.463648
tan ws;
0.5
acot 0.5;
1.10715
cot ws;
0.5
asec 3;
1.23096
sec ws;
3.0
acsc 3;
0.339837
csc ws;
3.0
asinh 0.5;
0.481212
sinh ws;
0.5
acosh 2;
1.31696
cosh ws;
2.0
atanh 0.5;
0.549306
tanh ws;
0.5
acoth 2;
0.549306
coth ws;
2.0
sech 1;
0.648054
asech ws;
1
csch 1;
0.850918
acsch ws;
1
expt(2,1.234);
2.35218
logb(ws,2);
1.234
hypot(3,4);
5.0
a := -3*pi/4;
a := - 2.35619
% Any -pi<a<=pi should work.
atan2(sin a,cos a);
- 2.35619
ws - a;
0
% Should be 0.
precision 20;
6
% Functions will be computed to 20 places.
sin 1.5;
0.99749498660405443094
asin ws;
1.5
precision 50;
20
% Functions computed to 50 places.
sin 1.5;
0.99749498660405443094172337114148732270665142592212
asin ws;
1.5
precision 6;
50
comment If argument or value are complex, functions are not computed
when dmode is ROUNDED;
$
sin(1+i);
sin(i + 1)
% Complex argument.
asin 2;
asin(2)
% Value would be complex.
on complex;
*** Domain mode rounded changed to complex-rounded
% Now complex arguments and complex results will be handled.
comment Complex functions of one argument:
EXP, LOG, NORM, ARG, SQRT,
SIN, ASIN, COS, ACOS, TAN, ATAN, COT, ACOT, SEC, ASEC, CSC, ACSC,
SINH, ASINH, COSH. ACOSH, TANH, ATANH, COTH, ACOTH, SECH, ASECH,
CSCH, ACSCH.
(The following functions are available only in symbolic mode, so they
are not tested here:
SIND, ASIND, COSD, ACOSD, TAND, ATAND, COTD, ACOTD, SECD, ASECD,
CSCD, ACSCD.)
Complex function of two variables: EXPT, LOGB, ATAN2;
e**(pi*i);
- 1 + 1.22465e-16*i
% Should be -1 (except for computational error.)
log(1+i);
0.346574 + 0.785398*i
exp ws;
1 + i
norm(5*exp(2i));
5.0
arg(5*exp(2i));
2.0
sqrt(1+i);
1.09868 + 0.45509*i
ws**2;
1 + i
asin 2;
1.5708 - 1.31696*i
sin ws;
2.0 - 1.06058e-16*i
acos 2;
1.31696*i
cos ws;
2.0
atan(1+i);
1.01722 + 0.402359*i
tan ws;
1 + i
acot(1+i);
0.553574 - 0.402359*i
cot ws;
1 + i
asec 0.1;
2.99322*i
sec ws;
0.1
acsc 0.1;
1.5708 - 2.99322*i
csc ws;
0.1 + 6.09254e-18*i
sinh(1+i);
0.634964 + 1.29846*i
asinh ws;
1 + i
cosh(1+i);
0.83373 + 0.988898*i
acosh ws;
1 + i
atanh 2;
0.549306 + 1.5708*i
tanh ws;
2.0 + 1.83697e-16*i
acoth 0.3;
0.30952 + 1.5708*i
coth ws;
0.3 - 5.57214e-17*i
asech(1-i);
0.530638 + 1.11852*i
sech ws;
1 - i
acsch(1-i);
0.530638 + 0.452278*i
csch ws;
1 - i
expt(1+i,1-i);
2.80788 + 1.31787*i
logb(ws,1+i);
1 - i
a := 1+i;
a := 1 + i
% Any a such that - pi < repart a <= pi should work.
atan2(sin a,cos a);
1 + i
ws - a;
0
% Should be 0.
clear a;
% Further math package tests.
%*********************************************************************
%**
%** The math package will compute the floating point values of **
%** the usual elementary functions, namely: **
%** sin asin sind asind sinh asinh **
%** cos acos cosd acosd cosh acosh **
%** tan atan tand atand tanh atanh **
%** cot acot cotd acotd coth acoth **
%** sec asec secd asecd sech asech **
%** csc acsc cscd acscd csch acsch **
%** atan2 atan2d **
%** exp ln sqrt **
%** expt log cbrt **
%** logb hypot **
%** log10 floor **
%** ceiling **
%** round **
%** **
%** All functions are computed to the accuracy of the floating- **
%** point precision of the system set up at the time. **
%** **
%*********************************************************************
%** File #1===Trig Function Tests===
%** Trig functions are tested in both degrees and radians modes.
%*********************************************************************
symbolic;
nil
math!!label;
"Math package mod 1.7, 1 May 93"
symbolic procedure terpr(i,j); if remainder(i,j)=0 then terpri()$
on rounded;
nil
% We need !!plumin, etc.
% #1: sind**2+cosd**2 test: ideal answers 1.0 for all args.
for i:=0:45 do <<write " ",i,"->",sind float i**2+cosd float i**2;
terpr(i,4)>>;
0->1.0
1->1.0 2->1.0 3->1.0 4->1.0
5->1.0 6->1.0 7->1.0 8->1.0
9->1.0 10->1.0 11->1.0 12->1.0
13->1.0 14->1.0 15->1.0 16->1.0
17->1.0 18->1.0 19->1.0 20->1.0
21->1.0 22->1.0 23->1.0 24->1.0
25->1.0 26->1.0 27->1.0 28->1.0
29->1.0 30->1.0 31->1.0 32->1.0
33->1.0 34->1.0 35->1.0 36->1.0
37->1.0 38->1.0 39->1.0 40->1.0
41->1.0 42->1.0 43->1.0 44->1.0
45->1.0
nil
% #2: Quadrant test of sind, cosd: proper answers + +,+ -,- -,- +.
begin scalar a;
a:= sind 45.0;
for i:= 0.0:3.0 do
<<write " ",sind(i*90+45)/a," ", cosd (i*90+45)/a;terpri()>>
end$
1.0 1.0
1.0 -1.0
-1.0 -1.0
-1.0 1.0
% #3: Scaling test: all values should be 1 exactly.
begin scalar a; a:= cosd 60.0;
% for i:= -10.0:10.0 do write fix(cosd(60+i*360)/a)," "
for i:= -10.0:10.0 do write round(cosd(60+i*360)/a)," "
end$
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
% #4: Test of radians -> degrees evaluation: ideal values 1.0.
array a(6)$
begin
for i:=1:6 do a(i):=sind(15.0*i);
for i:=1:6 do <<write sin(!!pii2*i/6.0)/a(i)," "; terpr(i,3)>>
end$
1.0 1.0 1.0
1.0 1.0 1.0
% #5: Test of tand*cotd: ideal values 1.0.
begin
for i:=5 step 5 until 85 do
<<write tand float i*cotd float i," "; terpr(i,25)>>;
terpri()
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0
% #6: Test of secd**2-tand**2: ideal values 1.0.
begin
for i:=5 step 5 until 85 do
<<write secd float i**2-tand float i**2," "; terpr(i,25)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0
% #7: Test of cscd**2-cotd**2: ideal values 1.0.
begin
for i:=5 step 5 until 85 do
<<write cscd float i**2-cotd float i**2," "; terpr(i,25)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0
% #8: Test of asind+acosd: ideal values 1.0.
begin write "sind and cosd"; terpri();
for i:=-10:10 do
<<write (asind(0.1*i)+acosd(0.1*i))/90," "; terpr(i,5)>>;
write "sin and cos";terpri();
for i:=-10:10 do
<<write (acos(0.1*i)+asin(0.1*i))/!!pii2," "; terpr(i,5)>>
end$
sind and cosd
1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
sin and cos
1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #9: Test of atand+acotd: ideal values 1.0.
begin scalar x; write "tand, atand and acotd"; terpri();
for i:=-80 step 10 until 80 do
<<x:=tand float i; write (atand x+acotd x)/90," "; terpr(i,50)>>;
terpri();
write "tan, atan and acot";terpri();
for i:=-80 step 10 until 80 do
<<x:= tan (!!pii2*i/90.0); write (atan x+acot x)/!!pii2," ";
terpr(i,50)>>
end$
tand, atand and acotd
1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0
tan, atan and acot
1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0
% #10: Test of atand tand: ideal values i for i:=-9:89.
begin
for i:=-9:89 do
<<write " ",i,"->",if i=0 then 1.0 else atand tand float i;
terpr(i,4)>>
end$
-9->-9.0 -8->-8.0
-7->-7.0 -6->-6.0 -5->-5.0 -4->-4.0
-3->-3.0 -2->-2.0 -1->-1.0 0->1.0
1->1.0 2->2.0 3->3.0 4->4.0
5->5.0 6->6.0 7->7.0 8->8.0
9->9.0 10->10.0 11->11.0 12->12.0
13->13.0 14->14.0 15->15.0 16->16.0
17->17.0 18->18.0 19->19.0 20->20.0
21->21.0 22->22.0 23->23.0 24->24.0
25->25.0 26->26.0 27->27.0 28->28.0
29->29.0 30->30.0 31->31.0 32->32.0
33->33.0 34->34.0 35->35.0 36->36.0
37->37.0 38->38.0 39->39.0 40->40.0
41->41.0 42->42.0 43->43.0 44->44.0
45->45.0 46->46.0 47->47.0 48->48.0
49->49.0 50->50.0 51->51.0 52->52.0
53->53.0 54->54.0 55->55.0 56->56.0
57->57.0 58->58.0 59->59.0 60->60.0
61->61.0 62->62.0 63->63.0 64->64.0
65->65.0 66->66.0 67->67.0 68->68.0
69->69.0 70->70.0 71->71.0 72->72.0
73->73.0 74->74.0 75->75.0 76->76.0
77->77.0 78->78.0 79->79.0 80->80.0
81->81.0 82->82.0 83->83.0 84->84.0
85->85.0 86->86.0 87->87.0 88->88.0
89->89.0
% #11: Test of acot cotd: ideal values 10*i for i:=1:17.
begin
for i:=10 step 10 until 170 do
<<write " ",i,"->",acotd cotd i; terpr(i,40)>>; terpri();terpri() end$
10->10.0 20->20.0 30->30.0 40->40.0
50->50.0 60->60.0 70->70.0 80->80.0
90->90.0 100->100.0 110->110.0 120->120.0
130->130.0 140->140.0 150->150.0 160->160.0
170->170.0
% #12: Test of asind sind: ideal values 10*i for i:=-9:9.
begin
for i:=-90 step 10 until 90 do
<<write " ",i,"->",asind sind float i; terpr(i,40)>>
end$
-90->-90.0 -80->-80.0
-70->-70.0 -60->-60.0 -50->-50.0 -40->-40.0
-30->-30.0 -20->-20.0 -10->-10.0 0->0.0
10->10.0 20->20.0 30->30.0 40->40.0
50->50.0 60->60.0 70->70.0 80->80.0
90->90.0
% #13: Test of acosd cosd: ideal values 10*i for i:=1:18.
begin
for i:=10 step 10 until 180 do
<<write " ",i,"->",acosd cosd float i; terpr(i,40)>>
end$
10->10.0 20->20.0 30->30.0 40->40.0
50->50.0 60->60.0 70->70.0 80->80.0
90->90.0 100->100.0 110->110.0 120->120.0
130->130.0 140->140.0 150->150.0 160->160.0
170->170.0 180->180.0
% #14: Test of acscd cscd: ideal values 10*i for i:=-9:9, except
% error for i=0.
begin
for i:=-90 step 10 until 90 do
<<write " ",i,"->",if i=0 then "error" else acscd cscd float i;
terpr(i,40)>>
end$
-90->-90.0 -80->-80.0
-70->-70.0 -60->-60.0 -50->-50.0 -40->-40.0
-30->-30.0 -20->-20.0 -10->-10.0 0->error
10->10.0 20->20.0 30->30.0 40->40.0
50->50.0 60->60.0 70->70.0 80->80.0
90->90.0
% #15: Test of asecd secd: ideal values 10*i for i :=0:18. except
% error for i=9.
begin
for i:=0 step 10 until 180 do
<<write" ",i,"->",if i=90 then "error" else asecd secd float i;
terpr(i,40)>>
end$
0->0.0
10->10.0 20->20.0 30->30.0 40->40.0
50->50.0 60->60.0 70->70.0 80->80.0
90->error 100->100.0 110->110.0 120->120.0
130->130.0 140->140.0 150->150.0 160->160.0
170->170.0 180->180.0
%*********************************************************************
%** ===Exp,Log,Sqrt,Cbrt, and Expt Function tests===
%*********************************************************************
% #16: Test of properties of exp function: ideal results 1.0.
array b(5)$
begin scalar x; x:=0;
write "multiplicative property";terpri();
for i:=0:5 do b(i):=1+i/6.0; for i:=0:5 do for j:=i:5 do
<<write " ",exp (b(i)+b(j))/exp(b(i))/exp(b(j));
terpr(x:=x+1,5)>>
end$
multiplicative property
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0
% #17: Various properties of exp: ideal results 1.0.
begin write "inverse property"$ terpri()$
for i:=1:5 do write " ",exp(b(i))*exp(-b(i));terpri();
write "squares"; terpri();
for i:=-10:10 do
<<write " ",sqrt(exp(0.2*i))/exp(0.1*i); terpr(i,5)>>;
write "cubes"; terpri();
for i:=-10:10 do
<<write " ",cbrt(exp(0.3*i))/exp(0.1*i); terpr(i,5)>>
end$
inverse property
1.0 1.0 1.0 1.0 1.0
squares
1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
cubes
1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #18: Test of log exp: ideal results 1.0.
begin for i:=-5:5 do
<<write if i=0 then "0/0" else (log exp float i)/i," "; terpr(i,5)>>
end$
1.0
1.0 1.0 1.0 1.0 0/0
1.0 1.0 1.0 1.0 1.0
% #19: Test of log10 expt(10.0,i): ideal results 1.0.
begin scalar i; write "small values i:=-5:5"; terpri();
for j:=-5:5 do
<<write if j neq 0 then log10 float expt(10.0,j)/j
else "zero"," ";
terpr(j,5)>>;
write "large i=2**j where j:=0:6"; terpri();
for j:=0:5 do
<<write (log10 float expt(10.0,2**j))/2**j," "; terpr(j,5)>>;
terpri();
write "noninteger values of i=j/10.0 where j:=1:20";terpri();
for j:=1:20 do
<<i:=j/10.0; write (log10 float expt(10,i))/i," ";
terpr(j,5)>>
end$
small values i:=-5:5
1.0
1.0 1.0 1.0 1.0 zero
1.0 1.0 1.0 1.0 1.0
large i=2**j where j:=0:6
1.0
1.0 1.0 1.0 1.0 1.0
noninteger values of i=j/10.0 where j:=1:20
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #20: Test of properties of expt(x,i)*(expt(x,-i). ideal result 1.0.
begin integer j;
for x:=2:6 do for i:=2:6 do
<<write expt(float x,i)*expt(float x,-i)," "; terpr(j:=j+1,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #21: Test of expt(-x,i)/expt(x,i) for fractional i.
begin integer j,k; write "odd numerator. ideal result -1.0"; terpri();
for i:=1:10 do
<<k:=(2*i-1.0)/(8*i+1); write rexpt(-8,k)/rexpt(8,k)," ";
terpr(j:=j+1,5)>>;
write "even numerator. ideal result 1.0"; terpri();
for i:=1:10 do
<<k:=(2.0*i)/(8*i+1); write rexpt(-8,k)/rexpt(8,k)," ";
terpr(j:=j+1,5)>>
end$
odd numerator. ideal result -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
even numerator. ideal result 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #22: Test of properties of ln or log or logb:
% inverse argument: ideal result -1.0.
begin integer x;
for i:=2:5 do for j:= 2:10 do
<<x:=x+1; write logb(float i,float j)/logb(float i,1.0/j)," ";
terpr(x,5)>>
end$
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0 -1.0 -1.0 -1.0 -1.0
-1.0
% #23: Test of log(a*b) = log a+log b: ideal result 1.0.
begin integer x;
for i:=1:5 do for j:=i:5 do
<<write log (i*j*0.01)/(log (i*0.1)+log(j*0.1))," ";
terpr(x:=x+1,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #24: Test of sqrt x*sqrt x/x for x:=5i*(5i/3)**i where i:=1:20
% (test values strictly arbitrary): ideal results 1.0.
begin scalar x,s;
for i:=1:20 do
<<x:= 5.0*i;s:=sqrt(x:=x*(expt(x/3,i))); write s*s/x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #25: Test of cbrt x**3/x for x:=5i*(5i/3)**i where i:=-9:10
% (test values strictly arbitrary):ideal results 1.0.
begin scalar x,s;
for i:=-9:10 do
<<x:= 5.0*i; if i neq 0 then s:= cbrt(x:=x*(expt(x/3,i)));
write if i=0 then 1 else s*s*s/x," "; terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
%*********************************************************************
%** ===Hyperbolic Function Tests===
%*********************************************************************
% #26: Test of sinh x+ cosh x= exp x: ideal results 1.0.
begin scalar x;
for i:=1:10 do
<<x:=ln float i$ write (sinh x+cosh x)/exp x," "$ terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #27: Test of cosh x-sinh x= exp(-x): ideal results 1.0.
begin scalar x;
for i:=1:10 do
<<x:=ln float i$ write(cosh x-sinh x)*exp x," "$ terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #28: Test of (cosh x)**2-(sinh x)**2: ideal results 1.0.
begin scalar x$
for i:=1:10 do
<<x:=ln float i$write(cosh x)**2-(sinh x)**2," "; terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #29: Test of tanh*cosh/sinh: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1);
write if i=10 then 1 else tanh x*cosh x/sinh x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #30: Test of tanh*coth: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else tanh x*coth x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #31: Test of sech*cosh: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write sech x*cosh x," "; terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #32: Test of csch*sinh: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else csch x*sinh x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #33: Test of asinh sinh: ideal results 1.0.
begin scalar x; for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else asinh sinh x/x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #34: Test of acosh cosh: ideal results 1.0. However, acosh x
% loses accuracy as x -> 1 since d/dx cosh x -> 0.
begin scalar x;
for i:=1:20 do
<<x:=ln(1+i*0.05); write acosh cosh x/x," "; terpr(i,5)>>
end$
1.0 1.0 0.99999999999999 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #35: Test of cosh acosh:ideal results 1.0.
begin scalar x;
for i:=1:50 do
<<x:=1+i/25.0; write (cosh acosh x)/x," "; terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #36: Test of atanh tanh: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else (atanh tanh x)/x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #37: Test of acoth coth: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else (acoth coth x)/x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #38: Test of asech sech: ideal results 1.0. However, asech x
% loses accuracy as x -> 1 since d/dx sech x -> 0.
begin scalar x;
for i:=1:20 do
<<x:=ln(1+i*0.05); write (asech sech x)/x," "; terpr(i,5)>>
end$
1.0 1.0 0.99999999999999 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
% #39: Test of acsch csch: ideal results 1.0.
begin scalar x;
for i:=1:20 do
<<x:=ln(i*0.1); write if i=10 then 1 else (acsch csch x)/x," ";
terpr(i,5)>>
end$
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1
1.0 1.0 1.0 1.0 1.0
1.0 1.0 1.0 1.0 1.0
end;
Time for test: 50 ms, plus GC time: 10 ms