urv_qrit.m

UTV工具包提供46个Matlab函数

function [R,V,U] = urv_qrit(p,num_ref,R,V,U)
%  urv_qrit --> Refinement of R in the URV decomposition via QR-iterations.
%
%  <Synopsis>
%    [R,V,U] = urv_qrit(p,num_ref,R,V,U)
%    [R,V]   = urv_qrit(p,num_ref,R,V)
%    [R]     = urv_qrit(p,num_ref,R)
%
%  <Description>
%    Given the URV decomposition U*R*V' with numerical rank p, the
%    function refines the rank-revealing decomposition via num_ref
%    steps of block QR iterations.
%
%  <Algorithm>
%    Refinement is identical to block QR iteration, in which the
%    off-diagonal block of the upper triangular matrix R is "flipped"
%    to the (2,1)-position and then back again.
%
%  <See Also>
%    urv_ref --> Refine one column of R in the URV decomposition.

%  <References>
%  [1] R. Mathias and G.W. Stewart, "A Block QR Algorithm and the Singular
%      Value Decomposition", Lin. Alg. Appl., 182 (1993), pp. 91--100.
%
%  <Revision>
%    Ricardo D. Fierro, California State University San Marcos
%    Per Christian Hansen, IMM, Technical University of Denmark
%    Peter S.K. Hansen, IMM, Technical University of Denmark
%
%    Last revised: June 22, 1999
%-----------------------------------------------------------------------

% Check the number of input and output arguments.
vflag = 1;
uflag = 1;
if (nargout == 1 & nargin == 3)
  vflag = 0;
  uflag = 0;
elseif (nargout == 2 & nargin == 4)
  uflag = 0;
elseif ~(nargout == 3 & nargin == 5)
  error('Wrong number of input or output arguments.')
end

% Initialize.
[n,n] = size(R);

if (p ~= abs(round(p))) | (p > n)
  error('Requires p to be an integer between 0 and n.')
end
if (num_ref ~= abs(round(num_ref)))
  error('Requires positive integer value for num_ref.')
end

% Refinement loop.
for (i = 1:num_ref)
  % Flip the off-diagonal block of R using Givens rotations on the right.
  for (r = p+1:n)
    for (i = p:-1:1)
      % Apply rotation to R on the right.
      [c,s,R(i,i)] = gen_giv(R(i,i),R(i,r));
      R(i,r) = 0;                               % Eliminate R(i,r).
      [R(1:i-1,i),R(1:i-1,r)] = app_giv(R(1:i-1,i),R(1:i-1,r),c,s);
      [R(p+1:r,i),R(p+1:r,r)] = app_giv(R(p+1:r,i),R(p+1:r,r),c,s);

      % Apply rotation to V on the right.
      if (vflag)
        [V(1:n,i),V(1:n,r)] = app_giv(V(1:n,i),V(1:n,r),c,s);
      end
    end
  end

  % Restore upper triangular form with Givens rotations on the left.
  for (r = n:-1:p+1)
    for (i = 1:p)
      % Apply rotation to R on the left.
      [c,s,R(i,i)] = gen_giv(R(i,i),R(r,i));
      R(r,i) = 0;                               % Eliminate R(r,i).
      [R(i,i+1:p),R(r,i+1:p)] = app_giv(R(i,i+1:p),R(r,i+1:p),c,s);
      [R(i,r:n),R(r,r:n)]     = app_giv(R(i,r:n),R(r,r:n),c,s);

      % Apply rotation to U on the right.
      if (uflag)
        [U(:,i),U(:,r)] = app_giv(U(:,i),U(:,r),c,s);
      end
    end
  end
end

%-----------------------------------------------------------------------
% End of function urv_qrit
%-----------------------------------------------------------------------