Matrices For 1.1

A= Array[a,{4,4}];
MatrixForm[A]

M1= {
{2,0,0,0},
{0,1,0,0},
{0,0,1,0},
{0,0,0,1}
};
MatrixForm[A.M1]

M2= {
{1,0,0,0},
{0,1,0,0},
{0,0,1/2,0},
{0,0,0,1}
};
MatrixForm[M2.A]

M3= {
{1,0,1,0},
{0,1,0,0},
{0,0,1,0},
{0,0,0,1}
};
MatrixForm[M3.A]

M4= {
{0,0,0,1},
{0,1,0,0},
{0,0,1,0},
{1,0,0,0}
};
MatrixForm[A.M4]

M5= {
{1,-1,0,0},
{0,1,0,0},
{0,-1,1,0},
{0,-1,0,1}
};
MatrixForm[M5.A]

M6= {
{1,0,0,0},
{0,1,0,0},
{0,0,1,1},
{0,0,0,0}
};
MatrixForm[A.M6]

M7= {
{0,0,0},
{1,0,0},
{0,1,0},
{0,0,1}
};
MatrixForm[A.M7]

The order is (rows at the front and cols at the back stacked outwards)
		M5.M3.M2.B.M1.M4.M6.M7

MatrixForm[M5.M3.M2]
MatrixForm[M1.M4.M6.M7]

Unit Balls

FirstSome defintitions.

UnitBall[t_,1]:= {
{1,t},
-{1,t},
{t,1},
-{t,1}
}
UnitBall[t_,∞]:= {
(1/2){1+t,1-t},
-(1/2){1+t,1-t},
(1/2){-1+t,1+t},
-(1/2){-1+t,1+t}
}
UnitBall[t_,2]:= {{t,(Sqrt[1-(t^2)])},{t,-(Sqrt[1-(t^2)])}}
UnitBall[t_,p_]:= {{t,(((1-((Abs[t])^p)))^(1/p))},{t,-(((1-((Abs[t])^p)))^(1/p))}}
ParametricPlot[
{UnitBall[t,1],UnitBall[t,∞],UnitBall[t,2],UnitBall[t,4]},
{t,-1,1}]

Induced Norms

AUnitBall[A_,t_,1]:= {
A.{1,t},
-A.{1,t},
A.{t,1},
-A.{t,1}
}
AUnitBall[A_,t_,∞]:= {
(1/2)A.{1+t,1-t},
-(1/2)A.{1+t,1-t},
(1/2)A.{-1+t,1+t},
-(1/2)A.{-1+t,1+t}
}
AUnitBall[A_,t_,2]:= {A.{t,(Sqrt[1-(t^2)])},A.{t,-(Sqrt[1-(t^2)])}}
AUnitBall[A_,t_,p_]:= {A.{t,(((1-((Abs[t])^p)))^(1/p))},A.{t,-(((1-((Abs[t])^p)))^(1/p))}}

A= { {1, 3},{2, 1}};
Show[
ContourPlot[
Norm[{x,y},1],
{x, -6,6},
{y,-6,6}
],
ParametricPlot[
{AUnitBall[A,t,1],AUnitBall[A,t,∞],AUnitBall[A,t,2]},
{t,-1,1},
AspectRatio->Automatic,
PlotStyle->Thick,
PlotLabel->A]
]

Singular Values

A= {
{1.2, 3.4},
{4.5, 2.3}
};
SingularValueList[A]
SingularValueDecomposition[A]

?SingularValueDecomposition

Householder Qs and a QR factorization

Here is a simple example of the procedure for a 4x4 matrix.

QFun[v_,i_]:= Module[
{m,w,Q},
m=Length[v];
w=v;
w[[i]] +=  Norm[v]; 
w =w/Norm[v];
Q = IdentityMatrix[m]-(2/(w.Conjugate[w]))Outer[Times,w, Conjugate[w]]
]

QRFactor[A_]:=Module[
{m,R,Q,QTemp,v,i,j},
m = Length[A];
Q=IdentityMatrix[m];
R=A;
Do[
v=R[[All,i]];
Do[ v[[j]]=0.0,
{j,1,i-1}];
QTemp=QFun[v,i];
R=QTemp.R;
Q=Q.QTemp,
{i,1,m-1}];
{Q,R}
]

A= RandomReal[ {-1.0,1.0},{4,4}];
{Q,R}=QRFactor[A];
MatrixForm[Q.Transpose[Q]]
MatrixForm[Q.R-A]

Back Substitution

The plan is to solve an m×m upper triangular system
	U x = b
for the m vector x.

Example:
	(({{(Subscript[u, 1,1]), (Subscript[u, 1,2]), (Subscript[u, 1,3]), (Subscript[u, 1,4])}, {0, (Subscript[u, 2,2]), (Subscript[u, 2,3]), (Subscript[u, 2,4])}, {0, 0, (Subscript[u, 3,3]), (Subscript[u, 3,4])}, {0, 0, 0, (Subscript[u, 4,4])}}))(({{(Subscript[x, 1])}, {(Subscript[x, 2])}, {(Subscript[x, 3])}, {(Subscript[x, 4])}}))=(({{(Subscript[b, 1])}, {(Subscript[b, 2])}, {(Subscript[b, 3])}, {(Subscript[b, 4])}}))
The equations are
	({{(Subscript[b, 4])==(Subscript[x, 4]) (Subscript[u, 4,4])}, {(Subscript[b, 3])==(Subscript[x, 3]) (Subscript[u, 3,3])+(Subscript[x, 4]) (Subscript[u, 3,4])}, {(Subscript[b, 2])==(Subscript[x, 2]) (Subscript[u, 2,2])+(Subscript[x, 3]) (Subscript[u, 2,3])+(Subscript[x, 4]) (Subscript[u, 2,4])}, {(Subscript[b, 1])==(Subscript[x, 1]) (Subscript[u, 1,1])+(Subscript[x, 2]) (Subscript[u, 1,2])+(Subscript[x, 3]) (Subscript[u, 1,3])+(Subscript[x, 4]) (Subscript[u, 1,4])}})
Solveing them sequentially gives
	({{({{(Subscript[x, 4])->((Subscript[b, 4])/(Subscript[u, 4,4]))}})}, {({{(Subscript[x, 3])->(((Subscript[b, 3])-(Subscript[x, 4]) (Subscript[u, 3,4]))/(Subscript[u, 3,3]))}})}, {({{(Subscript[x, 2])->(((Subscript[b, 2])-(Subscript[x, 3]) (Subscript[u, 2,3])-(Subscript[x, 4]) (Subscript[u, 2,4]))/(Subscript[u, 2,2]))}})}, {({{(Subscript[x, 1])->(((Subscript[b, 1])-(Subscript[x, 2]) (Subscript[u, 1,2])-(Subscript[x, 3]) (Subscript[u, 1,3])-(Subscript[x, 4]) (Subscript[u, 1,4]))/(Subscript[u, 1,1]))}})}})
Here is the start (last two entries) of a sample computation that might help you write the backsolve function.

Clear[U]
m=4;
U=RandomReal[{-1,1},{m,m}];
Do[ U[[i,j]]=0,{j,1,m},{i,j+1,m}]
b=RandomReal[{-1,1},{m}];
{MatrixForm[U],MatrixForm[b]}

x=Table[0,{m}];
x[[4]]=b[[4]]/U[[4,4]];
x[[3]]=((b[[3]]-U[[3,4]] x[[4]])/(U[[3,3]]));
MatrixForm[{LinearSolve[U,b],x}]

Your backsolve function should take in a matrix U and a vector b, check that they are the correct size, and output the solution vector.  You should check your code with randomly generated matrices and vectors against a built-in solver for various problem sizes.  You do not need to check that the matrix U is actually upper triangular!

SVD

{m,n}={5,4};
A=RandomReal[{-1,1},{m,n}];
{U,W,V}=SingularValueDecomposition[A];
{MatrixForm[U],MatrixForm[W],MatrixForm[V]}

The matrices U and V are Unitary

{MatrixForm[U.ConjugateTranspose[U]],MatrixForm[V.ConjugateTranspose[V]]}

This decomposes the matrix as a product of three matrices.

{MatrixForm[A], MatrixForm[U.W.ConjugateTranspose[V]]}

Of course, a good question is "SVD, what is it good for?" A defensible answer is "Almost everything!"

The matlab command is svd(A) producing the same format output. Both matlab and Mathematica produce the unhatted singular value decomposition discussed in the text!