Solution to a Linear Equality Constrained Least Squares Problem

using a generalized RQ factorization of matrices A and B, where A is M×N, B is P×N, and P N M+P. It is assumed that rank(B) = p and that the null spaces of A and B intersect only trivially; these conditions ensure that the problem has a unique solution.

Calling Sequence

CALL DGGLSE	(M, N, P, DA, LDA, DB, LDB, DC, DD, DX, DWORK, LDWORK, INFO)
CALL SGGLSE	(M, N, P, SA, LDA, SB, LDB, SC, SD, SX, SWORK, LDWORK, INFO)
CALL ZGGLSE	(M, N, P, ZA, LDA, ZB, LDB, ZC, ZD, ZX, ZWORK, LDWORK, INFO)
CALL CGGLSE	(M, N, P, CA, LDA, CB, LDB, CC, CD, CX, CWORK, LDWORK, INFO)
void dgglse	(int m, int n, int p, double *da, int lda, double *db, int ldb, double *dc, double *dd, double *dx, int *info)
void sgglse	(int m, int n, int p, float *sa, int lda, float *sb, int ldb, float *sc, float *sd, float *sx, int *info)
void zgglse	(int m, int n, int p, doublecomplex *za, int lda, doublecomplex *zb, int ldb, doublecomplex *zc, doublecomplex *zd, doublecomplex *zx, int *info)
void cgglse	(int m, int n, int p, complex *ca, int lda, complex *cb, int ldb, complex *cc, complex *cd, complex *cx, int *info)

Arguments

M	Number of rows of the matrix A. M 0.
N	Number of columns of the matrices A and B. N 0.
P	Number of rows of the matrix B. M+P N P 0.
xA	On entry, the M×N matrix A. On exit, A is overwritten.
LDA	Leading dimension of the array A as specified in a dimension or type statement. LDA max(1, M).
xB	On entry, the P×N matrix B. On exit, B is overwritten.
LDB	Leading dimension of the array B as specified in a dimension or type statement. LDB max(1, P).
xC	On entry, the right-hand side vector C of length M for the least squares part of the LSE problem. On exit, the residual sum of squares for the solution is given by the sum of squares of elements N-P+1 to M of vector C.
xD	On entry, the right-hand side vector D of length P for the constrained equation. On exit, D is overwritten.
xX	On exit, the N-element array X contains the solution of the LSE problem.
xWORK	Scratch array with a dimension of LDWORK.
LDWORK	Leading dimension of the array WORK as specified in a dimension or type statement. LDWORK N + P + max(M, N, P).
INFO	On exit:
	INFO = 0	Subroutine completed normally.
	INFO < 0	The ith argument, where i = |INFO|, had an illegal value.

Sample Program

PROGRAM TEST

IMPLICIT NONE

C

INTEGER LDA, LDB, LDWORK

INTEGER N

PARAMETER (N = 3)

PARAMETER (LDA = N)

PARAMETER (LDB = N)

PARAMETER (LDWORK = 3 * N)

C

DOUBLE PRECISION A(LDA,N), B(LDB,N), C(N), D(N)

DOUBLE PRECISION WORK(LDWORK), X(N)

INTEGER ICOL, INFO, IROW

C

EXTERNAL DGGLSE

C

C Initialize the array A to store the matrix A shown below.

C Initialize the array B to store the matrix B shown below.

C Initialize the array C to store the right hand side vector

C c for the least squares problem shown below. Initialize

C the array D to store the right hand side vector d for the

C constrained equation as shown below.

C

C 1 4 6 11 -8 2 100 1

C A = 2 9 14 B = -4 5 -2 c = 10 d = 10

C 3 14 23 1 -2 1 1 100

C

DATA A / 1.0D0, 2.0D0, 3.0D0,

$ 4.0D0, 9.0D0, 1.4D1,

$ 6.0D0, 1.4D1, 2.3D1 /

DATA B / 1.1D1, -4.0D0, 1.0D0,

$ -8.0D0, 5.0D0, 2.0D0,

$ 2.0D0, -2.0D0, 1.0D0 /

DATA C / 1.0D2, 1.0D1, 1.0D0 /

DATA D / 1.0D0, 1.0D1, 1.0D2 /

C

C Print the initial values of A and B.

C

PRINT 1000

PRINT 1010, ((A(IROW,ICOL), ICOL = 1, N), IROW = 1, N)

PRINT 1020

PRINT 1010, ((B(IROW,ICOL), ICOL = 1, N), IROW = 1, N)

C

C Compute and print the solution.

C

CALL DGGLSE (N, N, N, A, LDA, B, LDB, C, D, X, WORK, LDWORK,

$ INFO)

IF (INFO .LT. 0) THEN

PRINT 1030, ABS(INFO)

STOP 1

END IF

PRINT 1040

PRINT 1050, X

C

1000 FORMAT (1X, 'A:')

1010 FORMAT (3(3X, F4.1))

1020 FORMAT (/1X, 'B:')

1030 FORMAT (1X, 'Illegal value for argument #', I1,

$ ' in DGGLSE.')

1040 FORMAT (/1X, ' X')

1050 FORMAT (1X, F6.1)

C

END

Sample Output

A:
1.0 4.0 6.0
2.0 9.0 14.0
3.0 14.0 23.0
B:
11.0 -8.0 2.0
-4.0 5.0 -2.0
1.0 2.0 1.0
X
12.8
26.2
34.9