NUMAL Section 3.1.1.3.1.2
BEGIN SECTION : 3.1.1.3.1.2 (December, 1975)
AUTHOR : D.T.WINTER
INSTITUTE : MATHEMATICAL CENTRE
RECEIVED : 731217
BRIEF DESCRIPTION :
THIS SECTION CONTAINS 2 PROCEDURES FOR THE SOLUTION OF AN
UNDERDETERMINED SYSTEM OF LINEAR EQUATIONS:
SOLUND EXPECTS AS INPUT THE MATRIX OF THE SYSTEM OF EQUATIONS,
CALCULATES THE SINGULAR VALUES DECOMPOSITION BY MEANS OF THE
PROCEDURE QRISNGVALDEC, AND SOLVES THE SYSTEM BY MEANS OF THE
PROCEDURE SOLSVDUND.
SOLSVDUND ASSUMES THAT THE MATRIX IS
ALREADY DECOMPOSED AND SOLVES THE SYSTEM OF EQUATIONS, MULTIPLYING
THE RIGHT-HAND SIDE BY THE PSEUDO-INVERSE OF THE GIVEN MATRIX.
KEYWORDS :
BEST LEAST-SQUARES SOLUTION
SINGULAR VALUES
PSEUDO-INVERSE
SUBSECTION : SOLSVDUND
CALLING SEQUENCE :
THE HEADING OF THE PROCEDURE IS :
"PROCEDURE" SOLSVDUND(U, VAL, V, M, N, X, EM);
"VALUE" M, N; "INTEGER" M, N; "ARRAY" U, VAL, V, X, EM;
"CODE" 34282;
THE MEANING OF THE FORMAL PARAMETERS IS :
U: <ARRAY IDENTIFIER>;
"ARRAY" U[1:M,1:N];
ENTRY:THE MATRIX U IN THE SINGULAR VALUES DECOMPOSITION V*S*U'.
VAL: <ARRAY IDENTIFIER>;
"ARRAY" VAL[1:N];
ENTRY:THE SINGULAR VALUES;
V: <ARRAY IDENTIFIER>;
"ARRAY" V[1:N,1:N];
ENTRY:THE MATRIX V IN THE SINGULAR VALUES DECOMPOSITION.
N: <ARITHMETIC EXPRESSION>;
THE LENGTH OF THE RIGHT-HAND SIDE VECTOR;
M: <ARITHMETIC EXPRESSION>;
THE NUMBER OF UNKNOWNS, N SHOULD SATISFY N <= M;
X: <ARRAY IDENTIFIER>;
"ARRAY" X[1:M];
ENTRY: THE RIGHT-HAND SIDE VECTOR IN X[1:N];
EXIT: THE SOLUTION VECTOR IN X[1:M];
EM: <ARRAY IDENTIFIER>;
"ARRAY" EM[6:6];
ENTRY: EM[6]: THE MINIMAL NON-NEGLECTABLE SINGULAR VALUE.
PROCEDURES USED :
MATVEC = CP34011
TAMVEC = CP34012
REQUIRED CENTRAL MEMORY : AN AUXILIARY ARRAY OF N REALS IS DECLARED
RUNNING TIME : ROUGHLY PROPORTIONAL TO (M + N) * N
METHOD AND PERFORMANCE :
THE SOLUTION IS FOUND IN THREE STEPS :
1. V' * X = X1 IS CALCULATED,
2. VAL+ * X1 = X2 IS CALCULATED, HERE VAL+ DENOTES THE DIAGONAL
MATRIX OBTAINED FROM VAL BY SETTING VAL+[I,I] = 1/VAL[I] IF
VAL[I] GREATER THAN OR EQUAL TO EM[6], AND 0 OTHERWISE,
3. THE SOLUTION U * X2 IS CALCULATED.
LANGUAGE : ALGOL 60
SUBSECTION : SOLUND
CALLING SEQUENCE :
THE HEADING OF THE PROCEDURE IS :
"INTEGER" "PROCEDURE" SOLUND(A, M, N, X, EM);
"VALUE" M, N; "INTEGER" M, N; "ARRAY" A, X, EM;
"CODE" 34283;
SOLUND:= THE NUMBER OF SINGULAR VALUES NOT FOUND, I.E. ZERO IF ALL
SINGULAR VALUES ARE CALCULATED.
THE MEANING OF THE FORMAL PAREMETERS IS :
A: <ARRAY IDENTIFIER>;
"ARRAY" A[1:M,1:N];
ENTRY: THE TRANSPOSE OF THE MATRIX;
M: <ARITHMETIC EXPRESSION>;
THE NUMBER OF ROWS OF A;
N: <ARITHMETIC EXPRESSION>;
THE NUMBER OF COLUMNS OF A, N <= M;
X: <ARRAY IDENTIFIER>;
"ARRAY" X[1:M];
ENTRY: THE RIGHT-HAND SIDE VECTOR IN X[1:N];
EXIT: THE SOLUTION VECTOR.
EM: <ARRAY IDENTIFIER>;
"ARRAY" EM[0:7];
ENTRY: EM[0]: THE MACHINE PRECISION;
EM[2]: THE RELATIVE PRECISION FOR THE SINGULAR VALUES;
EM[4]: THE MAXIMAL NUMBER OF ITERATIONS TO BE PERFORMED IN
THE SINGULAR VALUES DECOMPOSITION;
EM[6]: THE MINIMAL NON-NEGLECTABLE SINGULAR VALUE;
EXIT: EM[1]: THE INFINITY NORM OF THE MATRIX;
EM[3]: THE MAXIMAL NEGLECTED SUPERDIAGONAL ELEMENT;
EM[5]: THE NUMBER OF ITERATIONS PERFORMED IN THE SINGULAR
VALUES DECOMPOSITION;
EM[7]: THE NUMERICAL RANK OF THE MATRIX, I.E. THE NUMBER OF
SINGULAR VALUES GREATER THAN OR EQUAL TO EM[6].
PROCEDURES USED :
QRISNGVALDEC = CP34273
SOLSVDUND = CP34282
REQUIRED CENTRAL MEMORY :
AUXILIARY ARRAYS ARE DECLARED TO A TOTAL OF (N + 1) * N REALS
METHOD AND PERFORMANCE :
THE SOLUTION IS FOUND IN TWO STEPS :
1. THE SINGULAR VALUES DECOMPOSITION IS CALCULATED BY MEANS OF THE
PROCEDURE QRISNGVALDEC;
2. THE SOLUTION IS CALCULATED BY MEANS OF THE PROCEDURE SOLSVDUND.
RUNNING TIME : ROUGHLY PROPORTIONAL TO (M + N) * N * N
LANGUAGE : ALGOL-60
REFERENCES :
WILKINSON, J.H. AND C.REINSCH
HANDBOOK OF AUTOMATIC COMPUTATION, VOL. 2
LINEAR ALGEBRA
HEIDELBERG (1971)
EXAMPLE OF USE :
FIRST WE GIVE A PROGRAM, AND THAN THE RESULTS OF THIS PROGRAM :
"BEGIN" "ARRAY" A[1:8,1:5], B[1:8], EM[0:7];
"INTEGER" I;
A[1,1]:=22; A[1,2]:= A[2,3]:=10; A[1,3]:= A[7,1]:= A[8,5]:=2;
A[1,4]:= A[3,5]:=3; A[1,5]:= A[2,2]:=7; A[2,1]:=14; A[2,5]:=8;
A[2,4]:= A[8,3]:=0; A[3,1]:= A[3,3]:= A[6,5]:=-1; A[3,2]:=13;
A[3,4]:=-11; A[4,1]:=-3; A[4,2]:= A[4,4]:= A[5,4]:= A[8,4]:=-2;
A[4,3]:=13; A[4,5]:= A[5,5]:= A[8,1]:=4; A[5,1]:= A[6,1]:=9;
A[5,2]:=8; A[5,3]:= A[6,2]:= A[7,5]:=1; A[6,3]:=-7;
A[6,4]:= A[7,4]:= A[8,2]:=5; A[7,2]:=-6; A[7,3]:=6;
B[1]:=-1; B[2]:=2; B[3]:=1; B[4]:=4; B[5]:=0;
EM[0]:="-14; EM[2]:="-12; EM[4]:=80; EM[6]:="-10;
I:= SOLUND(A, 8, 5, B, EM);
OUTPUT(61, "("4B, "("NUMBER SINGULAR VALUES NOT FOUND : ")",
3ZD,/, 4B, "("NORM : ")",N, /, 4B, "("MAX NEGL SUBD ELEM : ")",
N,/, 4B, "("NUMBER ITERATIONS : ")", 3ZD, /, 4B, "("RANK : ")",
3ZD, /")", I, EM[1], EM[3], EM[5], EM[7]);
OUTPUT(61, "("/, 4B, "("SOLUTION VECTOR")", /,/, 8(4B, N,/)")",
B[1], B[2], B[3], B[4], B[5], B[6], B[7], B[8])
"END"
NUMBER SINGULAR VALUES NOT FOUND : 0
NORM : +4.4000000000000"+001
MAX NEGL SUBD ELEM : +4.3977072741076"-014
NUMBER ITERATIONS : 6
RANK : 3
SOLUTION VECTOR
+1.6410256410255"-002
+1.4807692307694"-002
-4.8397435897438"-002
+1.0000000000002"-002
-6.7948717948740"-003
+1.1602564102565"-002
+2.9999999999996"-002
-8.3974358974328"-003
SOURCE TEXT(S):
"CODE" 34282;
"CODE" 34283;