SVD.mesa
Copyright Ó 1989, 1992 by Xerox Corporation. All rights reserved.
Ken Shoemake, June 30, 1989 2:33:15 am PDT
RowN: TYPE ~ LinearSystem.RowN;
ColumnN: TYPE ~ LinearSystem.ColumnN;
MatrixN: TYPE ~ LinearSystem.MatrixN;
badShape: ERROR [msg: Rope.ROPE];
noConvergence: SIGNAL [msg: Rope.ROPE];
SVDecomp:
PROC [a: MatrixN, m, n:
INT, w: ColumnN, v: MatrixN];
Given matrix A[0..m)[0..n), this routine computes its singular value decomposition,
A = U W VT . The matrix U replaces A on output. The diagonal matrix of singular values W is output as a vector W[0..n). The matrix V (not the transpose VT) is output as V[0..n)[0..n). The number of rows, m, must be greater or equal to the number of columns, n; if it is smaller, then a should be filled up to square with zero rows.
SVBackSub:
PROC [u: MatrixN, w: ColumnN, v: MatrixN, m, n:
INT, b: ColumnN, x: ColumnN];
Solves A x = B for a vector X, where A is specified by the arrays U, W, V, as returned by SVDecomp. M and N are the logical dimensions of A, and will be equal for square matrices. B[0..m) is the input right-hand side. X[0..n) is the output solution vector. No input quantities are destroyed, so the routine may be called sequentially with different B's.