4.5 Least-Squares and Least-Norm Problems

gels(A, B[, trans=’N’])

Solves least-squares and least-norm problems with a full rank m by n matrix A.

  1. trans is ’N’. If m is greater than or equal to n, gels() solves the least-squares problem
    minimize  ∥AX - B ∥ .
                  F

    If m is less than or equal to n, gels() solves the least-norm problem

    minimize   ∥X∥F
subject to AX  = B.

  2. trans is ’T’ or ’C’ and A and B are real. If m is greater than or equal to n, gels() solves the least-norm problem
    minimize  ∥X ∥F
subject to  AT X = B.

    If m is less than or equal to n, gels() solves the least-squares problem

    minimize ∥AT X - B ∥F .

  3. trans is ’C’ and A and B are complex. If m is greater than or equal to n, gels() solves the least-norm problem
    minimize  ∥XH∥F
subject to A  X = B.

    If m is less than or equal to n, gels() solves the least-squares problem

    minimize  ∥AHX  - B ∥F.

A and B must have the same typecode (’d’ or ’z’). trans = ’T’ is not allowed if A is complex. On exit, the solution X is stored as the leading submatrix of B. The array A is overwritten with details of the QR or the LQ factorization of A. Note that gels() does not check whether A is full rank.

geqrf(A, tau)

QR factorization of a real or complex matrix A:

A = QR.

If A is m by n, then Q is m by m and orthogonal/unitary, and R is m by n and upper triangular (if m is greater than or equal to n), or upper trapezoidal (if m is less than or equal to n). tau is a matrix of the same type as A and of length at least min{m,n}. On exit, R is stored in the upper triangular part of A. The matrix Q is stored as a product of min{m,n} elementary reflectors in the first min{m,n} columns of A and in tau.

ormqr(A, tau, C[, side=’L’[, trans=’N’]])

Product with a real orthogonal matrix:

                                                              {
                   ′ ′                        ′ ′               Q    trans = ′N ′
C := op(Q)C   (side = L),    C := C op(Q)  (side = R ),   op(Q ) =  QT   trans = ′T′,

where Q is square and orthogonal. Q is stored in A and tau as a product of min{A.size[0], A.size[1]} elementary reflectors, as computed by geqrf().

unmqr(A, tau, C[, side=’L’[, trans=’N’]])

Product with a real orthogonal or complex unitary matrix:

                                                              ({ Q    trans = ′N′
C := op(Q)C   (side = ′L′),   C := C op(Q)  (side = ′R′),  op(Q ) =  QT   trans = ′T′
                                                              ( QH   trans = ′C′,

Q is square and orthogonal or unitary. Q is stored in A and tau as a product of min{A.size[0], A.size[1]} elementary reflectors, as computed by geqrf(). The arrays A, tau and C must have the same type. trans = ’T’ is only allowed if the typecode is ’d’.

In the following example, we solve a least-squares problem by a direct call to gels(), and by separate calls to geqrf(), ormqr(), and trtrs().

>>> from cvxopt import random, blas, lapack, matrix  
>>> m, n = 10, 5  
>>> A, b = random.normal(m,n), random.normal(m,1)  
>>> x1 = +b  
>>> lapack.gels(+A, x1)                  # x1[:n] minimizes ||A*x1[:n] - b||_2  
>>> tau = matrix(0.0, (n,1))  
>>> lapack.geqrf(A, tau)                 # A = [Q1, Q2] * [R1; 0]  
>>> x2 = +b  
>>> lapack.ormqr(A, tau, x2, trans=’T’)  # x2 := [Q1, Q2]’ * b  
>>> lapack.trtrs(A[:n,:], x2, uplo=’U’)  # x2[:n] := R1^{-1}*x2[:n]  
>>> blas.nrm2(x1[:n] - x2[:n])  
3.0050798580569307e-16