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.

trans is 'N'. If m is greater than or equal to n, gels() solves the least-squares problem

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert AX-B\Vert _F. \end{array} \end{displaymath}$

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

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert X\Vert _F \\ \mbox{subject to} & AX = B. \end{array}\end{displaymath}$
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

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert X\Vert _F \\ \mbox{subject to} & A^TX=B. \end{array}\end{displaymath}$

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

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert A^TX-B\Vert _F. \end{array}\end{displaymath}$
trans is 'C' and A and B are complex. If m is greater than or equal to n, gels() solves the least-norm problem

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert X\Vert _F \\ \mbox{subject to} & A^HX=B. \end{array}\end{displaymath}$

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

$\begin{displaymath} \begin{array}{ll} \mbox{minimize} & \Vert A^HX-B\Vert _F. \end{array}\end{displaymath}$

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

is full rank.

geqrf( A, tau)

QR factorization of a real or complex matrix A:

$\begin{displaymath} A = Q R. \end{displaymath}$

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:

$\begin{displaymath} C := \mathop{\mathrm{op}}(Q)C \quad (\mathrm{side} = \mathr... ...} \\ Q^T & \mathrm{trans} = \mathrm{'T'}, \end{array}\right. \end{displaymath}$

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:

$\begin{displaymath} C := \mathop{\mathrm{op}}(Q)C \quad (\mathrm{side} = \mathr... ...} \\ Q^H & \mathrm{trans} = \mathrm{'C'}, \end{array}\right. \end{displaymath}$

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
>>> from cvxopt.base import 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