MatrixError.html

> restart;
Digits := 6;
nicelog := z -> if z < 10.0^(-Digits) then -Digits else log[10](z) fi;

>    with(Matlab):
with(linalg):
nTerms:=150;
i:=array(1..3,1..3,[[1,0,0],[0,1,0],[0,0,1]]);
A:=array(1..3,1..3,[[0.10,0.10,0.10],[0.20,0.75,0.05],[0.10,0.30,0.60]]);
B:=array(1..3,1..3,[[0.90,-0.10,0.00],[-0.20,0.25,-0.05],[-0.10,-0.30,0.40]]);
A_sum:=matrix(3,3,[0,0,0,0,0,0,0,0,0]):
B_sum:=matrix(3,3,[0,0,0,0,0,0,0,0,0]):
A_relerr:=vector(nTerms):
B_relerr:=vector(nTerms):
A_inv:=Matlab[inv](A):
B_inv:=Matlab[inv](B):
for k from 1 to nTerms do
   A_sum:=evalm(A_sum+(i-A)^(k-1)):
   B_sum:=evalm(B_sum+(i-B)^(k-1)):
   A_relerr[k]:=norm(A_sum-A_inv)/norm(A_inv):
   B_relerr[k]:=norm(B_sum-B_inv)/norm(B_inv):
od:
plot( [seq([j,nicelog(A_relerr[j])], j=1..nTerms)], axes=FRAME,
       labels = ["Order of Partial Sum", "Relative Error of Partial Sum"] );
plot( [seq([j,nicelog(B_relerr[j])], j=1..nTerms)], axes=FRAME,
       labels = ["Order of Partial Sum", "Relative Error of Partial Sum"] );

Warning, the names det and transpose have been redefined

Warning, the protected names norm and trace have been redefined and unprotected

>