Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Matrices/standard_complex_least_squares.ads, Revision 1.1.1.1
1.1 maekawa 1: with Standard_Complex_Vectors;
2: with Standard_Complex_Matrices;
3:
4: package Standard_Complex_Least_Squares is
5:
6: -- DESCRIPTION :
7: -- This package provides an implementation of Least Squares algorithm
8: -- for matrices of standard complex numbers.
9:
10: procedure QRLS ( x : in out Standard_Complex_Matrices.Matrix;
11: ldx,n,k : in natural;
12: qraux,y : in Standard_Complex_Vectors.Vector;
13: qy,qty,b,rsd,xb : out Standard_Complex_Vectors.Vector;
14: job : in natural; info : out natural );
15:
16: -- DESCRIPTION :
17: -- Applies the output of Standard_Complex_QR_Decomposition.QRD to
18: -- compute coordinate transformations, projections and least quares
19: -- solutions. For k <= min(n,p), let xk be the matrix
20: --
21: -- xk = (x(jpvt(1)),x(jpvt(2)), ... ,x(jpvt(k)))
22: --
23: -- formed from columnns jpvt(1), ... ,jpvt(k) of the original
24: -- n x p matrix x that was input to dqrdc (if no pivoting was done,
25: -- xk consists of the first k columns of x in their original order).
26: -- QRD produces a factored orthogonal matrix q and an upper triangular
27: -- matrix r such that
28: -- xk = q * (r)
29: -- (0)
30: -- this information is contained in coded form in x and qraux.
31:
32: -- ON ENTRY :
33: -- x contains the output of QRD, of size ldx times p;
34: -- ldx leading dimension of x;
35: -- n number of rows in the matrix xk, must be same as in QRD
36: -- k number of columns of the matrix xk, k <= min(n,p)
37: -- qraux contains p entries, auxiliary output from QRD;
38: -- y n-vector to be manipulated by QRLS;
39: -- job specifies what is to be computed, job has the decimal
40: -- expansion abcde, with the following meaning :
41: -- if a /= 0, compute qy,
42: -- if b,c,d, or e /= 0, compute qty,
43: -- if c /= 0, compute b,
44: -- if d /= 0, compute rsd,
45: -- if e /= 0, compute xb.
46: -- Note that a request to compute b, rsd, or xb
47: -- automatically triggers the computation of qty, for
48: -- which an array must be provided in the calling sequence.
49:
50: -- ON RETURN :
51: -- x may be altered when used as work space;
52: -- qy qy'length = n,
53: -- qy conntains q*y, if its computation has been requested;
54: -- qty qty'length = n,
55: -- qty contains trans(q)*y, if its computation has been
56: -- requested; here trans(q) is the transpose of the matrix q;
57: -- b b'length = k,
58: -- b contains the solution of the least squares problem
59: --
60: -- minimize norm2(y - xk*b),
61: --
62: -- if its computation has been requested. (Note that
63: -- if pivoting was requested in dqrdc, the j-th
64: -- component of b will be associated with column jpvt(j)
65: -- of the original matrix x that was input into dqrdc.)
66: -- rsd rsd'length = n,
67: -- rsd contains the least squares residual y - xk*b,
68: -- if its computation has been requested; rsd is also the
69: -- orthogonal projection of y onto the orthogonal complement
70: -- of the column space of xk;
71: -- xb x'length=n, xb contains the least squares approximation xk*b,
72: -- if its computation has been requested; xb is also the
73: -- orthogonal projection of y onto the column space of x;
74: -- info is zero unless the computation of b has been requested
75: -- and r is exactly singular.
76: -- In this case, info is the index of the first zero
77: -- diagonal element of r and b is left unaltered.
78:
79: -- The parameters qy, qty, b, rsd, and xb are not referenced
80: -- if their computation is not requested and in this case can be
81: -- replaced by dummy variables in the calling program.
82: -- To save storage, the user may in some cases use the same array
83: -- for different parameters in the calling sequence.
84: -- A frequently occuring example is when one wishes to compute
85: -- any of b, rsd, or xb and does not need y or qty. In this
86: -- case one may identify y, qty, and one of b, rsd, or xb, while
87: -- providing separate arrays for anything else that is to be computed.
88: -- Thus the calling sequence
89: --
90: -- QRLS(x,ldx,n,k,qraux,y,dum,y,b,y,dum,110,info)
91: --
92: -- will result in the computation of b and rsd, with rsd overwriting y.
93: -- More generally, each item in the following list contains groups of
94: -- permissible identifications for a single callinng sequence.
95: --
96: -- 1. (y,qty,b) (rsd) (xb) (qy)
97: -- 2. (y,qty,rsd) (b) (xb) (qy)
98: -- 3. (y,qty,xb) (b) (rsd) (qy)
99: -- 4. (y,qy) (qty,b) (rsd) (xb)
100: -- 5. (y,qy) (qty,rsd) (b) (xb)
101: -- 6. (y,qy) (qty,xb) (b) (rsd)
102: --
103: -- In any group the value returned in the array allocated to the group
104: -- corresponds to the last member of the group.
105:
106: -- ACKNOWLEDGMENT :
107: -- This Ada version is a translation of the LINPACK version, datet 08/14/78,
108: -- written by G.W. Stewart, University of Maryland, Argonne National Lab.
109:
110: end Standard_Complex_Least_Squares;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>