IDL Analyst Reference Guide: Linear Systems

IMSL_QRSOL

The IMSL_QRSOL function solves a real linear least-squares problem Ax = b.

Syntax

Result = IMSL_QRSOL(b[, a] [, AUXQR=variable] [, BASIS=variable]
[, /DOUBLE] [, QR=variable] [, PIVOT=variable] [, RESIDUAL=variable] [, TOLERANCE=value])

Return Value

The solution, x, of the linear least-squares problem Ax = b.

Arguments

b

Matrix containing the right-hand side.

a

(Optional) Two-dimensional matrix containing the coefficient matrix. Element A (i, j) contains the j-th coefficient of the i-th equation.

Keywords

AUXQR

Named variable in which the matrix containing the scalars t_k of the Householder transformations that define the decomposition, as computed in the IMSL_QRFAC procedure, is stored. The AUXQR, PIVOT, and QR keywords must be used together.

BASIS

Named variable containing an integer specifying the number of columns used in the solution. The value BASIS = k, if |r_k,k| < TOLERANCE*|r_0,0| and |r_i,i| ³ TOLERANCE*|r_0,0| for i = 0, 1, ..., k – 1. For more information on the use of this option, see Discussion.

DOUBLE

If present and nonzero, double precision is used.

QR

Named variable which stores the matrix containing Householder transformations that define the decomposition, as computed in the IMSL_QRFAC procedure. The AUXQR, PIVOT, and QR keywords must be used together.

PIVOT

Named variable in which the array containing the desired variable order and usage information is stored. The AUXQR, PIVOT, and QR keywords must be used together.

RESIDUAL

Named variable in which the matrix containing the residual vector b – Ax is stored.

TOLERANCE

Nonnegative tolerance used to determine the subset of columns of A to be included in the solution. Default: TOLERANCE = SQRT(e), where e is machine precision

Discussion

IMSL_QRSOL solves a system of linear least-squares problems Ax = b with column pivoting. It computes a QR factorization of the matrix AP, where P is the permutation matrix defined by the pivoting, and computes the smallest integer k satisfying |r_k,k| < TOLERANCE*|r_0,0| to the output keyword BASIS.

Householder transformations:

Q_k = I – t_ku_ku_k^T, k = 0, ..., min(m – 1, n) – 1

compute the factorization. The decomposition is computed in the form Q_{min (m – 1, n) – 1} ... Q₀ AP = R, so AP = QR where Q = Q₀ ... Q_{min (m – 1, n) – 1}. Since each Householder vector u_k has zeros in the first k + 1 entries, it is stored as part of column k of QR. The upper-trapezoidal matrix R is stored in the upper-trapezoidal part of the first min(m, n) rows of QR. The solution x to the least-squares problem is computed by solving the upper-triangular system of linear equations R (0:k, 0:k) y (0:k) = (Q^Tb) (0:k) with k = Basis – 1. The solution is completed by setting y (k:n – 1) to zero and rearranging the variables, x = Py.

If the QR and AUXQR keywords are specified, then the function computes the least-squares solution to Ax = b given the QR factorization previously defined. There are Basis columns used in the solution. Hence, in the case that all columns are free, x is computed as described in the default case.

Example

This example illustrates the least-squares solution of four linear equations in three unknowns by using column pivoting. This is equivalent to least-squares quadratic polynomial fitting to four data values. The polynomial is written as p(t) = x₀ + tx₁ + t²x₂ and the data pairs (t_i, b_i ), t_i = 2(i + 1), i = 0, 1, 2, 3. The solution to Ax = b is returned by the IMSL_QRSOL function.

RM, a, 4, 3  
; Define the coefficient matrix.  
   row 0:  1 2 4  
   row 1:  1 4 16  
   row 2:  1 6 36  
   row 3:  1 8 64  
RM, b, 4, 1  
; Define the right-hand side.  
   row 0:  4.999  
   row 1:  9.001  
   row 2:  12.999  
   row 3:  17.001  
x = IMSL_QRSOL(b, a)  
; Call IMSL_QRSOL.  
PM, x, Title = 'Solution', Format = '(f8.5)'  
; Output the results.  
Solution  
   0.99900  
   2.00020  
   0.00000  
PM, a # x - b, Title = 'Residual', Format = '(f10.7)'  
Residual  
    0.0004015  
   -0.0011997  
    0.0012007  
   -0.0004005