next up previous
Next: Eigenvalue Problems. Up: Maple Implementation Previous: Help File

Test File for the Share Library.

To successfully submit a program for the share library, you must provide a help file (as above) for the users' benefit, and a test file to enable the Maple group to detect when your code begins to decay (a program is said to decay when the underlying operating system or language evolves so that the program no longer runs in the current version. With a language committed to upward compatibility this is not a problem, but Maple has always preferred to advance rather than be committed to the mistakes of the past).

The test file must perform some simple tests of the program's performance and on successful completion of each subtask print out `okay' and nothing else. If the subtask fails then some diagnostic information should be printed out instead of the word `okay'.

The tests now run in less than 1 minute on a 33Mhz IBM- compatible 486 machine.

This file sets interface(echo=0).

It also defines the utility operators ones and zeros. In Maple V Release 2 there are better ways to define these utility operators. (namely ones := (m,n) -> linalg[matrix](m,n,1);, etc.)

interface(echo=0):
read `Echelon.m`;
ones := <array(1..m,1..n,[seq([seq(1,i=1..n)],j=1..m)])|m,n|i,j>:
zeros := <array(1..m,1..n,[seq([seq(0,i=1..n)],j=1..m)])|m,n|i,j>:
Strang := array(1..3,1..4,
  [ [ 1, 3, 3, 2],
    [ 2, 6, 9, 5],
    [-1,-3, 3, 0] ]):
R := RowEchelon(Strang,'dt'):
if (op(op(R)) = op(array(1..3,1..4,
  [ [ 1, 3, 0, 1],
    [ 0, 0, 1, 1/3],
    [ 0, 0, 0, 0] ]))  )        and 
   ( dt = 3) then print(okay) else print(R,dt) fi;

Goldberg1 := array(1..3,1..5,
  [ [ 0, 1, 1, 1, 1],
    [ 0,-1,-1, 1, 1],
    [ 0, 1, 1,-1, 2] ]):
R := RowEchelon(Goldberg1,'dt','rank','P','L','U'):
if (op(op(R)) = op(array(1..3,1..5,
  [ [ 0, 1, 1, 0, 0],
    [ 0, 0, 0, 1, 0],
    [ 0, 0, 0, 0, 1] ])) )   and
    (op(evalm( Goldberg1 - P &* L &* U &* R)) =
     op(zeros(3,5))) then print(okay) else print(R,dt) fi;

Goldberg2 := array(1..3,1..5,
  [ [ 1, 1, 1, 1, a],
    [ 0,-1,-1, 1, b],
    [ 0, 1, 1,-1, c]]):
R := RowEchelon(Goldberg2, 'dt','rank','P','L','U'):
if (op(op(R)) = op(array(
    1 .. 3, 1 .. 5,[(3, 1)=0,(3, 2)=0,(1, 1)=1,(3, 3)=0,(1, 2)=0,(3, 4
    )=0,(2, 1)=0,(1, 3)=0,(3, 5)=1,(2, 2)=1,(1, 4)=2,(2, 3)=1,(1, 5)=0,
    (2, 4)=-1,(2, 5)=0]))) and
   (op(evalm( Goldberg2 - P &* L &* U &* R)) = op(zeros(3,5))) then
   print(okay) else print(R,dt) fi;

Goldberg3 := subs(b=-c,op(Goldberg2)):
R := RowEchelon(Goldberg3,'dt','rank','P','L','U'):
if (op(evalm(Goldberg3 - P &* L &* U &* R)) = op(zeros(3,5))) then
   print(okay) else print(R,dt) fi;

Noble1 := array(1..3,1..4,[[1,-2,3,1],[2,k,6,6],[-1,3,k-3,0]]):
R := RowEchelon(Noble1,'dt','rank','P','L','U'):
if (dt=k*(k+4)) and 
   (op(evalm( Noble1 - P&*L&*U&*R)) = op(zeros(3,4))) 
then print(okay) else print(R,dt) fi;
Test two degenerate cases
Degenerate1 := zeros(5,8):
R := RowEchelon(Degenerate1,'dt','rank','P','L','U'):
if ( op(evalm(P&*L&*U&*R-Degenerate1)) = op(zeros(5,8)) )
then print(okay) else print(R,dt,rank) fi;

Degenerate2 := ones(1,1):
R := RowEchelon(Degenerate2,'dt','rank','P','L','U'):
if (dt = 1) and (rank=1) and 
   ( op(evalm(P&*L&*U&*R-Degenerate2)) = op(zeros(1,1)) ) 
then print(okay) else print(R,dt,rank) fi;
Test complex numbers
Normalizer := evalc:
Complex1 := array(1..2,1..2,[[ I, 1],[-I,1]]):
R := RowEchelon(Complex1,'dt','rank','P','L','U'):
if (rank = 2) and (dt = 2*I) and
   (op(evalm(P&*L&*U&*R-Complex1)) = op(zeros(2,2)))
then print(okay) else print(R,dt,rank) fi;

Complex2 := array(1..2,1..2,[[ 1, I],[-I,1]]):
R := RowEchelon(Complex2,'dt','rank','P','L','U'):
if (rank = 1) and
   (op(evalm(P&*L&*U&*R-Complex2)) = op(zeros(2,2)))
then print(okay) else print(R,dt,rank) fi;
Now test the algebraic numbers (in a simple way)
 alias(I=RootOf(z^2+1,z)):
Normalizer := x -> evala(Normal(x)):
Algebraic1 := array(1..2,1..2,[[ I, 1],[-I,1]]):
R := RowEchelon(Algebraic1,'dt','rank','P','L','U'):
if (rank = 2) and
   (op(map(Normalizer,evalm(P&*L&*U&*R-Algebraic1)))
        = op(zeros(2,2)))
then print(okay) else print(R,dt,rank) fi;

Algebraic2 := array(1..2,1..2,[[ 1, I],[-I,1]]):
R := RowEchelon(Algebraic2,'dt','rank','P','L','U'):
if (rank = 1) and
   (op(map(Normalizer,evalm(P&*L&*U&*R-Algebraic2)))
        = op(zeros(2,2)))
then print(okay) else print(R,dt,rank) fi;
One that is not convenient for Hermite normal form:
Normalizer := normal:
Transcendental1 := array(1..2,1..2,[[x,tan(x)],[tan(x),x]]):
R := RowEchelon(Transcendental1,'dt'):
if (dt = x^2-tan(x)^2) then print(okay) else print(R,dt) fi;

quit



Robert Corless
Wed Jan 31 11:33:59 EST 1996