NAG Toolbox: nag_lapack_zgeqp3 (f08bt)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – complex array

The first dimension of the array a must be at least $\max \left(1,\mathbf{m}\right)$.

The second dimension of the array a must be at least $\max \left(1,\mathbf{n}\right)$.

The $m$ by $n$ matrix $A$.

2:     $\mathrm{jpvt}\left(:\right)$ – int64int32nag_int array

The dimension of the array jpvt must be at least $\max \left(1,\mathbf{n}\right)$

If $\mathbf{jpvt}\left(j\right)\ne 0$, then the $j$ th column of $A$ is moved to the beginning of $AP$ before the decomposition is computed and is fixed in place during the computation. Otherwise, the $j$ th column of $A$ is a free column (i.e., one which may be interchanged during the computation with any other free column).

Optional Input Parameters

1:     $\mathrm{m}$ – int64int32nag_int scalar

Default: the first dimension of the array a.

$m$, the number of rows of the matrix $A$.

Constraint: $\mathbf{m}\ge 0$.

2:     $\mathrm{n}$ – int64int32nag_int scalar

Default: the second dimension of the array a.

$n$, the number of columns of the matrix $A$.

Constraint: $\mathbf{n}\ge 0$.

Output Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – complex array

The first dimension of the array a will be $\max \left(1,\mathbf{m}\right)$.

The second dimension of the array a will be $\max \left(1,\mathbf{n}\right)$.

If $m\ge n$, the elements below the diagonal store details of the unitary matrix $Q$ and the upper triangle stores the corresponding elements of the $n$ by $n$ upper triangular matrix $R$.

If $m<n$, the strictly lower triangular part stores details of the unitary matrix $Q$ and the remaining elements store the corresponding elements of the $m$ by $n$ upper trapezoidal matrix $R$.

The diagonal elements of $R$ are real.

2:     $\mathrm{jpvt}\left(:\right)$ – int64int32nag_int array

The dimension of the array jpvt will be $\max \left(1,\mathbf{n}\right)$

Details of the permutation matrix $P$. More precisely, if $\mathbf{jpvt}\left(j\right)=k$, then the $k$th column of $A$ is moved to become the $j$ th column of $AP$; in other words, the columns of $AP$ are the columns of $A$ in the order $\mathbf{jpvt}\left(1\right),\mathbf{jpvt}\left(2\right),\dots ,\mathbf{jpvt}\left(n\right)$.

3:     $\mathrm{tau}\left(:\right)$ – complex array

The dimension of the array tau will be $\max \left(1,\min \left(\mathbf{m},\mathbf{n}\right)\right)$

Further details of the unitary matrix $Q$.

4:     $\mathrm{info}$ – int64int32nag_int scalar

$\mathbf{info}=0$ unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

   $\mathbf{info}=-i$

If $\mathbf{info}=-i$, parameter $i$ had an illegal value on entry. The parameters are numbered as follows:

1: m, 2: n, 3: a, 4: lda, 5: jpvt, 6: tau, 7: work, 8: lwork, 9: rwork, 10: info.

It is possible that info refers to a parameter that is omitted from the MATLAB interface. This usually indicates that an error in one of the other input parameters has caused an incorrect value to be inferred.

Accuracy

Further Comments

[a, info] = f08at(a(:,1:m), tau);

[a, info] = f08at(a, tau);

[c, info] = f08au('Left','Conjugate Transpose', a(:,min(m,n)), tau, c);

Example

Open in the MATLAB editor: f08bt_example

function f08bt_example


fprintf('f08bt example results\n\n');

a = [ 0.47 - 0.34i, -0.40 + 0.54i,  0.60 + 0.01i,  0.80 - 1.02i;
     -0.32 - 0.23i, -0.05 + 0.20i, -0.26 - 0.44i, -0.43 + 0.17i;
      0.35 - 0.60i, -0.52 - 0.34i,  0.87 - 0.11i, -0.34 - 0.09i;
      0.89 + 0.71i, -0.45 - 0.45i, -0.02 - 0.57i,  1.14 - 0.78i;
     -0.19 + 0.06i,  0.11 - 0.85i,  1.44 + 0.80i,  0.07 + 1.14i];

b = [-1.08 - 2.59i,  2.22 + 2.35i;
     -2.61 - 1.49i,  1.62 - 1.48i;
      3.13 - 3.61i,  1.65 + 3.43i;
      7.33 - 8.01i, -0.98 + 3.08i;
      9.12 + 7.63i, -2.84 + 2.78i];
[m,n] = size(a);
nrhs = size(b,2);
jpvt = zeros(n,1,'int64');

% Compute the QR factorization of a
[qr, jpvt, tau, info] = f08bt( ...
			      a, jpvt);

% Compute C = (C1) = (Q^H)*b, storing the result in c
%             (C2)
[c, info] = f08au( ...
		   'Left', 'Conjugate Transpose', qr, tau, b);

% Choose tol to reflect the relative accuracy of the input data
tol = 0.01;

% Determine and print the rank, k, of r relative to tol
k = find(abs(diag(qr)) <= tol*abs(qr(1,1)));
if numel(k) == 0
  k = numel(diag(qr));
else
  k = k(1)-1;
end

fprintf('\nTolerance used to estimate the rank of a\n     %11.2e\n', tol);
fprintf('Estimated rank of a\n        %d\n', k);

% Compute least-squares solution by backsubstitution in r(1:k, 1:k)*c = c1
c1 = zeros(m, nrhs);
c1(1:k, :) = inv(triu(qr(1:k,1:k)))*c(1:k,:);

% Compute estimates of the square roots of the residual sums of
% squares (2-norm of each of the columns of C2)
for j = 1:nrhs
  rnorm(j) = norm(c(k+1:m,j));
end

% Permute the least-squares solutions stored in c1 to give x = p*y
x = zeros(n, nrhs);
for i=1:n
   x(jpvt(i), :) = c1(i, :);
end
fprintf('\nLeast-squares solution(s)\n');
disp(x);
fprintf('Square root(s) of the residual sum(s) of squares\n');
disp(rnorm);

f08bt example results


Tolerance used to estimate the rank of a
        1.00e-02
Estimated rank of a
        3

Least-squares solution(s)
   0.0000 + 0.0000i   0.0000 + 0.0000i
   2.7020 + 8.0911i  -2.2682 - 2.9884i
   2.8888 + 2.5012i   0.9779 + 1.3565i
   2.7100 + 0.4791i  -1.3734 + 0.2212i

Square root(s) of the residual sum(s) of squares
    0.2513    0.0810