#
#
# MAPLE SESSION 1
#
>   ?LinearAlgebra
#
# 9.1   Vectors, Arrays and Matrices} 
#
#
#
# MAPLE SESSION 2
#
>     Matrix(3);
>     Matrix(3,4);
>     Matrix(2,3,[[a,b,c],[d,e,f]]);
>     Matrix(2,3,[[a,b],[d,e,f]]);  
>     Matrix(2,3,[[a,b],[c,d,e,f]]);
>     Matrix(2,3,[a,b,c,d,e,f]);
#
#
# MAPLE SESSION 3
#
>    W:=Vector(4);
>    V:=Vector([x,y,z]);
#
#
# MAPLE SESSION 4
#
>    f := (i,j) -> x^(i*j);
>    A := Matrix(2,2,f);
#
#
# MAPLE SESSION 5
#
>    A := Matrix(4,4,f); 
>    factor(LinearAlgebra[Determinant](A)); 
#
#
# MAPLE SESSION 6
#
>    M:=Matrix(5,(i,j)->5*i+j-5);
#
#
# MAPLE SESSION 7
#
>    map(ithprime,M);
#
#
# MAPLE SESSION 8
#
>    Matrix(5,(i,j)->ithprime(5*i+j-5));
#
#
# MAPLE SESSION 9
#
>    Matrix(10,(i,j)->ithprime(10*i+j-10));
#
# 9.1.1   Matrix and Vector entry assignment
#
#
#
# MAPLE SESSION 10
#
>    A:=Matrix(2,3,[[1,2,3],[5,10,16]]);
>    A[2,3];
#
#
# MAPLE SESSION 11
#
>    A[2,3]:=15;
>    A;
#
#
# MAPLE SESSION 12
#
>    A := Matrix(4,(i,j)->(i+j));
>    A[2..3,2..4];
>    B := Matrix(2,3,[[0,1,2],[3,4,5]]);
>    A[2..3,2..4]:=B;
>    A;
#
#
# MAPLE SESSION 13
#
>    B:=Matrix(3,(i,j)->b[i,j]);  
>    B[[3,2,2,1],1..3];
#
#
# MAPLE SESSION 14
#
>    A := Matrix(3,4,[[1,2,3,4],[2,4,6,8],[3,6,9,12]]); 
>    A[[3,2],[4,3,2]]; 
>    V := Vector([a,b,c,d]); 
>    W := V[[3,2]]; 
#
# 9.1.2   The Matrix palette
#
#
#
# MAPLE SESSION 15
#
>   \bold |   
#
#
# MAPLE SESSION 16
#
>   Matrix([[\butbox%?, %?], [%?, %?]]);
#
#
# MAPLE SESSION 17
#
>   Matrix([[23, %?], [%?, %?]]);
#
#
# MAPLE SESSION 18
#
>   Matrix([[23, \butbox%?], [%?, %?]]);
#
#
# MAPLE SESSION 19
#
>   Matrix([[23, int(1/x,x=1..2), [\butbox%?, %?]]);
#
#
# MAPLE SESSION 20
#
>   Matrix([[23, int(1/x,x=1..2), [25, \butbox%?]]);
#
# 9.1.3   Matrix operations
#
#
#
# MAPLE SESSION 21
#
>    A := Matrix(2,[[1,2],[3,4]]);
>    B := Matrix(2,[[-2,3],[-5,1]]);
>    A + B;
>    A - B;
>    5*A;
#
#
# MAPLE SESSION 22
#
>    A := Matrix(2,[[1,2],[3,4]]): 
>    B := Matrix(2,[[-2,3],[-5,1]]): 
>    A . B;
>    AI := 1/A;
>    A . AI;
>    A^3;
#
#
# MAPLE SESSION 23
#
>   with(LinearAlgebra);
#
#
# MAPLE SESSION 24
#
>    with(LinearAlgebra): 
>    A := Matrix(2,[[1,2],[3,4]]): 
>    B := Matrix(2,[[-2,3],[-5,1]]): 
>    Multiply(A , B);
>    Multiply(Multiply(A,A),A);
>    AI := MatrixInverse(A);
>    Transpose(A);
>    Trace(A);
#
#
# MAPLE SESSION 25
#
>    with(LinearAlgebra):  
>    A:=Matrix(2,3,[[1,2,3],[4,5,6]]); 
>    B:=Matrix(3,2,[[2,4],[-7,3],[5,1]]); 
>    C:=Matrix(2,2,[[1,-2],[-3,4]]); 
>    A . B;  
>    Multiply(A,B); 
>    A.B-2*C; 
#
# 9.1.4   Matrix and vector construction shortcuts
#
#
#
# MAPLE SESSION 26
#
>    V := <1,2,3>;
#
#
# MAPLE SESSION 27
#
>    R := <1|2|3>;
#
#
# MAPLE SESSION 28
#
>    U := <a,b,c>;
>    V := <i,j,k>;
>    W := <x,y,z>;
>    M := <U | V | W>;
#
#
# MAPLE SESSION 29
#
>    U := <a|b|c>; 
>    V := <i|j|k>; 
>    W := <x|y|z>; 
>    M := <U , V , W>; 
#
#
# MAPLE SESSION 30
#
>    A:=Matrix(3,(i,j)->a^i*b^j):  
>    B:=Matrix(3,(i,j)->b^i*c^j): 
>    C:=Matrix(3,(i,j)->c^i*a^j): 
>    A,B,C;
#
#
# MAPLE SESSION 31
#
>    <A|B|C>;
#
#
# MAPLE SESSION 32
#
>    <A,B>;
#
#
# MAPLE SESSION 33
#
>    <A,B,C>; 
#
# 9.1.5   Viewing large Matrices and Vectors} 
#
#
#
# MAPLE SESSION 34
#
>    M:=Matrix(50,20,(i,j)->ithprime(20*i+j-20)); 
#
#
# MAPLE SESSION 35
#
>   M := \_rtable[17353556];
#
#
# MAPLE SESSION 36
#
>    M[26..29,10..11];
#
# 9.2   Matrix context menu} 
#
#
#
# MAPLE SESSION 37
#
>    M:=Matrix(4,(i,j)->2^(i*j));
#
#
# MAPLE SESSION 38
#
>   R0 := [LinearAlgebra:-Dimensions(\_rtable[5673280])];
#
#
# MAPLE SESSION 39
#
>   R1 := [LinearAlgebra:-Determinant(\_rtable[5673280])];
#
#
# MAPLE SESSION 40
#
>   ifactor(R1);
> $
> (2)^20\,(3)^2\,(7)
> $
#
# 9.2.1   The ``Export As'' sub-menu
# 9.2.2   The ``Norm'' sub-menu
# 9.2.3   The ``Solvers'' sub-menu
# 9.2.4   The ``Conversions'' sub-menu
#
#
#
# MAPLE SESSION 41
#
>    M:=<<Pi,exp(1)>|<log(2),int(1/(sqrt(1+x^10)),x=0..1)>>;
#
#
# MAPLE SESSION 42
#
>    R0 := evalf(M,10);
#
# 9.2.4   The ``In-place Options'' sub-menu
# 9.3   Elementary row and column operations
#
#
#
# MAPLE SESSION 43
#
>    with(LinearAlgebra): 
>    A:=Matrix([[1, 1, 3, -3], 
>    [5, 5, 13, -7], 
>    [3, 1, 7, -11]]);
>    RowOperation(A,[2,1],-5);
>    A;
#
#
# MAPLE SESSION 44
#
>    RowOperation(A,[2,1],-5,inplace=true);
>    A;
#
#
# MAPLE SESSION 45
#
>    restart: 
>    with(LinearAlgebra): 
>    A:=Matrix([[1, 1, 3, -3], [5, 5, 13, -7],
>   [3, 1, 7, -11]]);
>    A1 := RowOperation(A,[2,1],-5);
>    A2 := RowOperation(A1,[3,1],-3);
>    A3 := RowOperation(A2,[2,3]);
#
#
# MAPLE SESSION 46
#
>    A4:=RowOperation(A3,2,-1/2); 
>    A5:=RowOperation(A4,3,-1/2); 
>    A6:=RowOperation(A5,[2,3],-1); 
>    A7:=RowOperation(A6,[1,3],-3); 
>    A8:=RowOperation(A7,[1,2],-1); 
#
# 9.4   Gaussian elimination
#
#
#
# MAPLE SESSION 47
#
>    with(LinearAlgebra): 
>    A:=Matrix([[1, 1, 3, -3], [5, 5, 13, -7], 
>   [3, 1, 7, -11]]): 
>    LUDecomposition(A,output='U');
#
#
# MAPLE SESSION 48
#
>    LUDecomposition(A,output='R');
#
# 9.5   Inverses, determinants, minors and the adjoint
#
#
#
# MAPLE SESSION 49
#
>    with(LinearAlgebra): 
>    A:= Matrix([[1,1,3],[5,5,13],[3,1,7]]);
>    Determinant(A);
>    B := MatrixInverse(A);
#
#
# MAPLE SESSION 50
#
>    B.A;
#
#
# MAPLE SESSION 51
#
>    with(LinearAlgebra): 
>    A := Matrix([[1,1,3],[5,5,13],[3,1,7]]): 
>    C := Adjoint(A);
#
#
# MAPLE SESSION 52
#
>    C.A;
#
#
# MAPLE SESSION 53
#
>    with(LinearAlgebra): 
>    A := Matrix([[1,1,3],[5,5,13],[3,1,7]]): 
>    Minor(A,2,3);
#
# 9.6   Special matrices and vectors
# 9.6.1   Band matrix
#
#
#
# MAPLE SESSION 54
#
>    with(LinearAlgebra): 
>    BandMatrix([a,b,c,d],1,5,7);
#
#
# MAPLE SESSION 55
#
>   BandMatrix([a,b,c,d],2,5,7);
#
# 9.6.2   Constant matrices and vectors
#
#
#
# MAPLE SESSION 56
#
>    with(LinearAlgebra): 
>    ConstantMatrix(-3,2,4); 
#
#
# MAPLE SESSION 57
#
>    with(LinearAlgebra): 
>    ConstantVector(11,2);   
>    ConstantVector[row](11,3);   
#
# 9.6.3   Diagonal matrices
#
#
#
# MAPLE SESSION 58
#
>    with(LinearAlgebra): 
>    DiagonalMatrix([a,b,c],3);
#
#
# MAPLE SESSION 59
#
>    DiagonalMatrix([a,b,c],4); 
>    DiagonalMatrix([<<a,c>|<b,d>>,e,
>  <<f,i,l>|<g,j,m>|<h,k,n>>]); 
#
# 9.6.4   Givens rotation matrices
#
#
#
# MAPLE SESSION 60
#
>    with(LinearAlgebra): 
>    V := <3,4,5>; 
>    GivensRotationMatrix(V, 1, 2); 
#
# 9.6.5   Hankel matrices
#
#
#
# MAPLE SESSION 61
#
>    with(LinearAlgebra): 
>    HankelMatrix([a,b,c,d,e]);
#
#
# MAPLE SESSION 62
#
>    HankelMatrix([a,b,c,d,e,f,g,h,i,j,k]); 
>    HankelMatrix([a,b,c,d,e,f,g,h,i,j,k],4); 
#
# 9.6.6   Hilbert matrices
#
#
#
# MAPLE SESSION 63
#
>    with(LinearAlgebra): 
>    HilbertMatrix(3); 
#
#
# MAPLE SESSION 64
#
>    HilbertMatrix(3,x);  
>    HilbertMatrix(4,5,x);  
#
# 9.6.7   Householder matrices
#
#
#
# MAPLE SESSION 65
#
>    with(LinearAlgebra): 
>    V:=Vector([1,2,2]); 
>    M := HouseholderMatrix(V);
#
#
# MAPLE SESSION 66
#
>    M.V; 
#
# 9.6.8   Identity matrix
#
#
#
# MAPLE SESSION 67
#
>    with(LinearAlgebra): 
>    IdentityMatrix(3);
#
#
# MAPLE SESSION 68
#
>    IdentityMatrix(4); 
>    IdentityMatrix(4,6); 
#
# 9.6.9   Jordan block matrices
#
#
#
# MAPLE SESSION 69
#
>    with(LinearAlgebra): 
>    JordanBlockMatrix([[-1,2],[5,3]]);
#
#
# MAPLE SESSION 70
#
>    JordanBlockMatrix([[2,2],[3,1],[4,3],[1,2]]); 
>    JordanBlockMatrix([[2,2],[3,1],[4,3],[1,2]],10); 
#
# 9.6.10   Random matrices and vectors
#
#
#
# MAPLE SESSION 71
#
>    with(LinearAlgebra): 
>    RandomMatrix(2,3);
#
#
# MAPLE SESSION 72
#
>    RandomMatrix(2,3,generator=0..9); 
>    RandomMatrix(6,generator=-10.0..20.0); 
#
#
# MAPLE SESSION 73
#
>    RandomMatrix(20,3, density=0.1, generator=0..1.0); 
#
#
# MAPLE SESSION 74
#
>    RandomMatrix(4,generator=1..9,
>  outputoptions=[shape=triangular[upper]]); 
#
#
# MAPLE SESSION 75
#
>   rand1 := 2*rand(0..1)-1: 
#
#
# MAPLE SESSION 76
#
>   RandUniUpMat := proc(n::posint,a::integer,b::integer) 
>    local M1,i: 
>    M1:=LinearAlgebra[RandomMatrix](n,generator=a..b,
>  xx outputoptions=[shape=triangular[upper]]); 
>    for i from 1 to n do 
>     M1[i,i]:=rand1(): 
>    end do: 
>   return M1; 
>   end proc; 
>    RandUniUpMat(3,-5,5);
#
#
# MAPLE SESSION 77
#
>   RandUniMat := proc(n::posint,a::integer,b::integer) 
>    local M1, M2: 
>    M1 := RandUniUpMat(n,a,b): 
>    M2 := RandUniUpMat(n,a,b): 
>    M1.LinearAlgebra[Transpose](M2); 
>   end proc; 
#
#
# MAPLE SESSION 78
#
>    with(LinearAlgebra): 
>    M:=RandUniMat(3,10,50);
>    Determinant(M);
>    MatrixInverse(M);
#
#
# MAPLE SESSION 79
#
>    with(LinearAlgebra): 
>    RandomVector(4); 
>    RandomVector(6,generator=0..1.0);       
>    RandomVector[row](6,generator=0..1.0);       
>    RandomVector[row](6,generator=(2*rand(0..1)-1)); 
#
# 9.6.11   Toeplitz matrices
#
#
#
# MAPLE SESSION 80
#
>    with(LinearAlgebra): 
>    ToeplitzMatrix([a,b,c,d,e]);
#
#
# MAPLE SESSION 81
#
>    ToeplitzMatrix([a,b,c,d,e,f,g]); 
>    ToeplitzMatrix([a,b,c,d,e,f,g],3); 
>    ToeplitzMatrix([a,b,c,d,e,f,g],symmetric); 
#
# 9.6.11   Vandermonde matrices
#
#
#
# MAPLE SESSION 82
#
>    with(LinearAlgebra): 
>    VandermondeMatrix([a,b,c]);
>    Determinant(%);
>    factor(%);
#
#
# MAPLE SESSION 83
#
>    V := VandermondeMatrix([a,b,c,d,e]); 
>    factor(Determinant(V)); 
>    VandermondeMatrix([a,b,c,d,e],3,4); 
#
# 9.6.12   Zero matrices and vectors
#
#
#
# MAPLE SESSION 84
#
>    with(LinearAlgebra): 
>    ZeroMatrix(2,3);
#
#
# MAPLE SESSION 85
#
>    with(LinearAlgebra): 
>    ZeroMatrix(4); 
#
#
# MAPLE SESSION 86
#
>    ZeroVector(4); 
>    ZeroVector[row](4); 
#
# 9.7   Systems of linear equations
#
#
#
# MAPLE SESSION 87
#
>    with(LinearAlgebra): 
>    EqList:= [10*x-27*y+z+r+2*s-11*t = 1,  
>    20*x-62*y+29*z+20*r+11*s-16*t = 1,  
>    -x-8*y+36*z+24*r +9*s+9*t = 1,  
>    -8*x+27*y-19*z-13*r-6*s+5*t = -5];
>    Vars:=[x,y,z,r,s,t];
#
#
# MAPLE SESSION 88
#
>    (A,b) := GenerateMatrix(EqList,Vars);
#
#
# MAPLE SESSION 89
#
>    AM := GenerateMatrix(EqList,Vars,augmented=true);
#
#
# MAPLE SESSION 90
#
>    LinearSolve(AM);
#
#
# MAPLE SESSION 91
#
>    SOL := LinearSolve(AM, free='w');
#
#
# MAPLE SESSION 92
#
>    A . SOL = b;
#
#
# MAPLE SESSION 93
#
>    SOL := LinearSolve(A, b, free='w'); 
#
#
# MAPLE SESSION 94
#
>    with(LinearAlgebra): 
>    ATM := <<1,0,0,0>|<1,0,0,0>|<0,1,0,0>|<-3,0,1,0>
>  |<4,2,-5,0>>;
>    V := BackwardSubstitute(ATM);
#
#
# MAPLE SESSION 95
#
>    A := ATM[1..-1,1..-2];
>    b := ATM[1..-1,-1];
#
#
# MAPLE SESSION 96
#
>    A . V = b;
#
#
# MAPLE SESSION 97
#
>    V := LinearSolve(ATM, method='subs');
#
# 9.8   Row space, column space, nullspace
#
#
#
# MAPLE SESSION 98
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,5,[[1,4,-10,3,-3], 
>                   [10,41,-102, 30,-31], 
>                   [-9,-19,56,-27,10]]);
>    Rank(A);
>    RowSpace(A);
>    ColumnSpace(A);
>    NSA := NullSpace(A);
#
#
# MAPLE SESSION 99
#
>    W:=<NSA[1] | NSA[2] | NSA[3]>;
>    A . W;
#
# 9.9   Eigenvectors and diagonalization
#
#
#
# MAPLE SESSION 100
#
>    with(LinearAlgebra): 
>    A:=<<177,-546,-364>|<77,-236,-154>|<-28,84,51>>;
>    Eigenvalues(A);
#
#
# MAPLE SESSION 101
#
>    Eigenvectors(A);
#
#
# MAPLE SESSION 102
#
>    (EG,P):=Eigenvectors(A); 
>    MatrixInverse(P).A.P; 
#
#
# MAPLE SESSION 103
#
>    JordanForm(A); 
>    JordanForm(A,output='Q'); 
#
#
# MAPLE SESSION 104
#
>    with(LinearAlgebra): 
>    A := Matrix(2,2,[[1.0,2.0],[3.0,4.0]]); 
>    Eigenvalues(A); 
>    Eigenvectors(A); 
>    B := Matrix(2,2,[[1+10*I,-8*I],[12*I, 1-10*I]]); 
>    Eigenvalues(B); 
>    Eigenvectors(B); 
>    P := JordanForm(B,output='Q'); 
#
# 9.10   Jordan form
#
#
#
# MAPLE SESSION 105
#
>    with(LinearAlgebra): 
>    C := Matrix(4,4,[[10,10,-14,15],[0,3,0,0], 
>    	[8,1,-13,8],[1,-8,-2,-4]]); 
>    Q := JordanForm(C,output='Q'); 
>    MatrixInverse(Q).C.Q; 
#
# 9.11   Inner products, and Vector and matrix norms
# 9.11.1   The dot product and bilinear forms} 
#
#
#
# MAPLE SESSION 106
#
>    with(LinearAlgebra): 
>    V:=Vector([seq(v[i],i=1..4)]);
>    W:=Vector([seq(w[i],i=1..4)]);
>    DotProduct(V,W);
>    DotProduct(V,W,conjugate=false);
>    DotProduct(<1,2,3>,<3,2,1>);
#
#
# MAPLE SESSION 107
#
>    with(LinearAlgebra): 
>    L:=Matrix(3,3,[[1,0,0],[1,1,0],[-3,-2,1]]): 
>    DG:=DiagonalMatrix([1,2,3]): 
>    A:=L.DG.Transpose(L);
#
#
# MAPLE SESSION 108
#
>    IsDefinite(A);
#
#
# MAPLE SESSION 109
#
>    BilinearForm(<1,1,1>,<1,3,1>,A);
#
# 9.11.2   Vector Norms} 
#
#
#
# MAPLE SESSION 110
#
>    with(LinearAlgebra): 
>    V:=Vector([seq(v[i],i=1..4)]);
>    VectorNorm(V,3);
>    VectorNorm(V,2);
#
#
# MAPLE SESSION 111
#
>    VectorNorm(V);
>    VectorNorm(<-4,3,-5>,infinity);
#
# 9.11.3   Matrix Norms} 
#
#
#
# MAPLE SESSION 112
#
>    with(LinearAlgebra): 
>    A:=Matrix(2,2,[seq([a[i,1],a[i,2]],i=1..2)]);
>    MatrixNorm(A,Frobenius);
>    B:=<<1,2>|<3,4>>;
>    MatrixNorm(B,Frobenius);
#
#
# MAPLE SESSION 113
#
>    with(LinearAlgebra): 
>    A:=<<177,-546,-364>|<77,-236,-154>|<-28,84,51>>; 
>    MatrixNorm(A,1); 
>    MatrixNorm(A,2); 
>    MatrixNorm(A,infinity); 
#
# 9.12   Least squares problems
#
#
#
# MAPLE SESSION 114
#
>    X:=[0.70, 0.76, 0.37, 0.82, 0.29, 0.56, 0.42, 0.47];
>    Y:=[0.035, 0.025,-0.18, 0.045,-0.16,-0.058,-0.11,
>  -0.085];
>    A:=Matrix([seq([1,X[k]],k=1..8)]); 
>    b:=Vector([seq(Y[k],k=1..8)]); 
#
#
# MAPLE SESSION 115
#
>    with(LinearAlgebra): 
>    c := LeastSquares(A,b);
#
#
# MAPLE SESSION 116
#
>    pts := [seq([X[k],Y[k]],k=1..8)];
>    bestline := c[1] + c[2]*x;
>    with(plots): 
>    PL1 := plot(pts,style=point,symbol=circle): 
>    PL2 := plot(bestline,x=0..1): 
>    display(PL1,PL2); 
#
# 9.13   QR-factorization and the Gram-Schmidt process
#
#
#
# MAPLE SESSION 117
#
>    with(LinearAlgebra): 
>    A:=<<1,2>|<12,9>>;
>    (Q,R):=QRDecomposition(A);
>    Q.R;
#
#
# MAPLE SESSION 118
#
>    Transpose(Q).Q;
#
#
# MAPLE SESSION 119
#
>    V1:=<33,-12,-12,-12>: 
>    V2:=<3,6,-20,-20>: 
>    V3:=<21,29,3,68>: 
>    A := <V1 | V2 | V3>;
#
#
# MAPLE SESSION 120
#
>    with(LinearAlgebra): 
>    (Q,R) := QRDecomposition(A);
#
#
# MAPLE SESSION 121
#
>    Transpose(Q).Q;
#
#
# MAPLE SESSION 122
#
>    M := <A | Q>;
>    Rank(A);
>    Rank(M);
#
#
# MAPLE SESSION 123
#
>    with(LinearAlgebra): 
>    V1 := <33,-12,-12,-12>; 
>    V2 := <3,6,-20,-20>; 
>    V3 := <21,29,3,68>; 
>    OBAS := GramSchmidt([V1,V2,V3]); 
#
#
# MAPLE SESSION 124
#
>    M := convert(OBAS, Matrix); 
>    Transpose(M).M; 
#
#
# MAPLE SESSION 125
#
>    with(LinearAlgebra): 
>    V1 := <33,-12,-12,-12>; 
>    V2 := <3,6,-20,-20>; 
>    V3 := <21,29,3,68>; 
>    ONBAS := GramSchmidt([V1,V2,V3],normalized); 
#
# 9.14   LU-factorization
#
#
#
# MAPLE SESSION 126
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,3,[[2,1,-3],[-4,4,8],[-6,9,3]]);
>    LUDecomposition(A);
#
#
# MAPLE SESSION 127
#
>    with(LinearAlgebra): 
>    A:=Matrix([[1, 1, 3, -3], [5, 5, 13, -7], 
>  [3, 1, 7, -11]]); 
>    LUDecomposition(A); 
>    (p,l,u):=LUDecomposition(A); 
>    p.l.u; 
#
#
# MAPLE SESSION 128
#
>    LUDecomposition(A,method=RREF); 
>    (p,l,u,r):=LUDecomposition(A,method=RREF); 
>    p.l.u.r; 
#
#
# MAPLE SESSION 129
#
>    LUDecomposition(A,output=['P','L','U1','R']); 
>    LUDecomposition(A,output=['R']); 
>    LUDecomposition(A,output=['P','L']); 
#
#
# MAPLE SESSION 130
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,3,[[1,-3,-5],[-3,13,27],[-5,27,62]]);
>    IsDefinite(A);
#
#
# MAPLE SESSION 131
#
>    L:=LUDecomposition(A,method=Cholesky);
#
#
# MAPLE SESSION 132
#
>    L.Transpose(L); 
#
# 9.15   Other {\it LinearAlgebra} functions
#
#
#
# MAPLE SESSION 133
#
>    with(LinearAlgebra): 
>    A := Matrix(3,3,[[0,-1,2],[3,-4,6],[2,-2,3]]); 
>    FrobeniusForm(A); 
>    (F,Q) := FrobeniusForm(A,output=['F','Q']); 
>    MatrixInverse(Q).A.Q; 
#
#
# MAPLE SESSION 134
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,4,[[1,2,3,4],[2,-3,4,5],[3,-7,8,9]]); 
>    GenerateEquations(A,[x,y,z,w]); 
#
#
# MAPLE SESSION 135
#
>    GenerateEquations(A,[x,y,z]);   
#
#
# MAPLE SESSION 136
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,3,[[5-x,5-2*x,2-x],
>  [-x^2+x-4,-2*x^2+2*x-1,-x^2+x], 
>  [-x^3-5,-2*x^3-10,-x^3-5]]); 
>    (H,U) := HermiteForm(A,x,output=['H','U']); 
>    HH:= map(expand,U.A); 
>    Equal(H,HH); 
#
#
# MAPLE SESSION 137
#
>    with(LinearAlgebra): 
>    A:=Matrix(2,[[2,1-I],[1+I,1]]); 
>    U:=1/sqrt(3)*Matrix(2,[[1-I,-1],[1,1+I]]); 
>    HermitianTranspose(U).U; 
>    HermitianTranspose(U).A.U; 
#
#
# MAPLE SESSION 138
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,[[2,1+I,I],[1-I,1,3],[-I,3,1]]); 
>    (H,Q) := HessenbergForm(A, output=['H','Q']); 
>    map(fnormal[6],H); 
>    HH := map(simplify[zero],%); 
>    Q.H.HermitianTranspose(Q); 
>    map(fnormal[6],%); 
>    map(simplify[zero],%); 
#
#
# MAPLE SESSION 139
#
>    with(LinearAlgebra): 
>    V1:=<2|3|5|-1>; 
>    V2:=<3|9|6|-1>; 
>    V3:=<6|32|10|-1>; 
>    V4:=<8|21|17|-3>; 
>    V5:=<19|52|26|-4>; 
>    IntersectionBasis([ [V1,V2,V3], [V4,V5] ]); 
#
#
# MAPLE SESSION 140
#
>    with(LinearAlgebra): 
>    M := Matrix([[1,2],[3,4]]); 
>    Map(x->1/x,M); 
>    M; 
#
#
# MAPLE SESSION 141
#
>    with(LinearAlgebra): 
>    A:=<<1|2|3>,<4|5|6>,<7|8|9>>; 
>    B:=<<a|b|c>,<d|e|f>,<g|h|i>>; 
>    MatrixAdd(A,B,1,-3,inplace); 
>    A; 
#
#
# MAPLE SESSION 142
#
>    with(LinearAlgebra): 
>    A := Matrix(3,3,[[0,-1,2],[3,-4,6],[2,-2,3]]);
>    mpoly := MinimalPolynomial(A,x);
>    charpoly := CharacteristicPolynomial(A,x);
>    normal(charpoly/mpoly);
#
#
# MAPLE SESSION 143
#
>    P := unapply( mpoly, x );
>    P(A);
#
#
# MAPLE SESSION 144
#
>    with(LinearAlgebra): 
>    V := <1 | 2 | 3 | 4>; 
>    W := Normalize(V); 
>    U := Normalize(V,2); 
>    DotProduct(U,U); 
#
#
# MAPLE SESSION 145
#
>    with(LinearAlgebra): 
>    A := Matrix(4,[[1,2,3,4],[0,-2,5,7],[0,3,5,6],
>  [0,11,-9,3]]); 
>    Pivot(A,2,2); 
#
#
# MAPLE SESSION 146
#
>    Pivot(A,2,2,[3,4]); 
#
#
# MAPLE SESSION 147
#
>    with(LinearAlgebra): 
>    ScalarMatrix(lambda,3); 
>    ScalarMatrix(lambda,3,4); 
#
#
# MAPLE SESSION 148
#
>    with(LinearAlgebra): 
>    ScalarVector(x,3,4); 
#
#
# MAPLE SESSION 149
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,[[2,1+I,I],[1-I,1,3],[1,0,1]]); 
>    (T,Z) := SchurForm(A, output=['T','Z']); 
>    map(fnormal[6],T); 
>    TT := map(simplify[zero],%); 
>    Z.T.HermitianTranspose(Z); 
>    map(fnormal[6],%); 
>    map(simplify[zero],%); 
#
#
# MAPLE SESSION 150
#
>    with(LinearAlgebra): 
>    A := RandomMatrix(4,3,outputoptions=[datatype=float]);
>    S := SingularValues(A, output='list');
>    map(evalf[7],%);
#
#
# MAPLE SESSION 151
#
>    Sig := DiagonalMatrix( S[1..3], 4, 3 ); 
>    U, Vt := SingularValues(A, output=['U', 'Vt']); 
>    U.Sig.Vt; 
#
#
# MAPLE SESSION 152
#
>    with(LinearAlgebra): 
>    A:=Matrix(3,3,[[5-x,5-2*x,2-x],
>  [-x^2+x-4,-2*x^2+2*x-1,-x^2+x], 
>  [-x^3-5,-2*x^3-10,-x^3-5]]); 
>    S := SmithForm(A); 
#
#
# MAPLE SESSION 153
#
>    (U,V) := SmithForm(A,x,output=['U','V']); 
>    U.A.V; 
>    map(simplify,%); 
#
#
# MAPLE SESSION 154
#
>    with(LinearAlgebra): 
>    A:=Matrix(6,[seq([seq(a[i,j],j=1..6)],i=1..6)]); 
>    SubMatrix(A,[1,3,5],[1..3,6]); 
#
#
# MAPLE SESSION 155
#
>    with(LinearAlgebra): 
>    V:=Vector([2,4,6,8,10]); 
>    SubVector(V,[1,4,5]);
#
#
# MAPLE SESSION 156
#
>    with(LinearAlgebra): 
>    V1:=<2|3|5|-1>; 
>    V2:=<3|9|6|-1>; 
>    V3:=<6|32|10|-1>; 
>    V4:=<8|21|17|-3>; 
>    V5:=<19|52|26|-4>; 
>    B1:=SumBasis([ [V1,V2,V3], [V4,V5] ]); 
>    B2:=Basis([V1,V2,V3,V4,V5]);
#
#
# MAPLE SESSION 157
#
>    with(LinearAlgebra): 
>    p:=x -> x^2-2: 
>    q:=x -> x^2-3: 
>    SylvesterMatrix(p(t),q(x-t),t); 
>    Determinant(%); 
>    solve(%=0,x);
#
#
# MAPLE SESSION 158
#
>    with(LinearAlgebra): 
>    A := RandomMatrix(3,3,outputoptions=[datatype=float,
>  shape=symmetric]); 
>    TridiagonalForm(A); 
>    TridiagonalForm(A, output=NAG); 
>    (T,Q) := TridiagonalForm(A, output=['T','Q']); 
>    Q.T.Transpose(Q); 
>    map(fnormal[6],%); 
>    map(simplify[zero],%);
#
#
# MAPLE SESSION 159
#
>    with(LinearAlgebra): 
>    UnitVector(3,4);
#
#
# MAPLE SESSION 160
#
>    with(LinearAlgebra): 
>    U := <a | b | c>; 
>    V := <i | j | k>; 
>    VectorAdd(U,V,3,4);
#
#
# MAPLE SESSION 161
#
>    with(LinearAlgebra): 
>    U := <1 | 2 | 3>; 
>    V := <4 | 1 | -2>; 
>    W := <1 | 1 | 1>; 
>    VectorAngle(U,2*U); 
>    VectorAngle(U,V); 
>    VectorAngle(U,W);
#
#
# MAPLE SESSION 162
#
>    with(LinearAlgebra): 
>    A:=<<a | b>, <c | d>>; 
>    B:=<<x | y>, <z | w>>; 
>    Zip(f,A,B); 
>    Zip(`+`,A,B); 
>    Zip(`*`,A,B);
#
# 9.16   The {\it linalg} package
#
#
#
# MAPLE SESSION 163
#
>    with(linalg);
#
# 9.16.1   Matrices and vectors}     
#
#
#
# MAPLE SESSION 164
#
>    with(linalg): 
>    v:=vector([1,2,3]);
>    A := matrix(2,3,[a,b,c,d,e,f]);
>    A := matrix(2,3,[[a,b,c],[d,e,f]]);
>    v;
>    A;
>    print(v);
>    print(A);
#
#
# MAPLE SESSION 165
#
>    op(A); 
>    eval(A); 
>    evalm(A);
#
# 9.16.2   Conversion between {\it linalg} and {\it LinearAlgebra}
#
#
#
# MAPLE SESSION 166
#
>    with(linalg): 
>    with(LinearAlgebra): 
>    A := matrix(3,3,(i,j)->(i+j)); 
>    Determinant(A);
>    det(A);
#
#
# MAPLE SESSION 167
#
>    B := convert(A, Matrix);
>    Determinant(B);
#
#
# MAPLE SESSION 168
#
>    C := convert(B, matrix);
>    det(C);
#
#
# MAPLE SESSION 169
#
>    with(LinearAlgebra): 
>    A := Matrix(3,3,(i,j)->x^(i+j)); 
>    B := A - y*IdentityMatrix(3);
>    simplify(B);
>    C := A - 5*IdentityMatrix(3);
#
#
# MAPLE SESSION 170
#
>    with(linalg):             
>    A := matrix(3,3,(i,j)->x^(i+j)); 
>    I3 := diag(1,1,1);
>    B := A - y*I3;           
>    evalm(B);
#
# 9.16.3   Matrix operations in {\it linalg}
#
#
#
# MAPLE SESSION 171
#
>    with(linalg): 
>    A := matrix(2,2,[1,2,3,4]): 
>    B := matrix(2,2,[-2,3,-5,1]): 
>    A+B;
>    evalm(%);
#
#
# MAPLE SESSION 172
#
>    with(linalg): 
>    A:=matrix(2,3,[1,2,3,4,5,6]); 
>    B:=matrix(3,2,[2,4,-7,3,5,1]); 
>    C:=matrix(2,2,[1,-2,-3,4]); 
>    A\&*B; 
>    evalm(%); 
>    multiply(A,B); 
>    evalm(A\&*B-2*C);