lurvdemo.m

UTV工具包提供46个Matlab函数

function lurvdemo
%  lurvdemo --> Demonstrates low-rank URV algorithms.

%  <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
%-----------------------------------------------------------------------

fprintf(1,'.                                                           \n');
fprintf(1,'.  The purpose of LURVDEMO is to demonstrate                \n');
fprintf(1,'.  the low-rank URV algorithms in LURV and                  \n');
fprintf(1,'.  LURV_A.                                                  \n');
fprintf(1,'.                                                           \n');
input('.                                 [Press RETURN to continue]');
fprintf(1,'.                                                           \n');

format short e;

% Test matrix generation.
m = 50;
n = 20;
randn('seed',100);
A = randn(m,n);
numrank = 13;

[u,s,v] = svd(A);

s1 = 2*logspace(1,-3,numrank);
s2 = 5*logspace(-4,-6,n-numrank);
s(1:n,:) = diag([s1 s2]);

A = u*s*v';

% Define input parameters.
tol_rank = 0.001;
max_iter = 5;
num_ref  = 0;
est_type = 1;

fprintf(1,'.                                                           \n');
fprintf(1,'.  Input parameters to LURV                                 \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    no. rows m of A                     = %3.0f \n',m);
fprintf(1,'.    no. cols n of A                     = %3.0f \n',n);
fprintf(1,'.    rank tolerance (tol_rank)           = %6.4e \n',tol_rank);
fprintf(1,'.    no. of steps of the principle               \n');
fprintf(1,'.    singular vector estimator (max_iter)= %3.0f \n',max_iter);
fprintf(1,'.    no. refinement steps (num_ref)      = %3.0f \n',num_ref);
fprintf(1,'.    use Lanczos procedure (est_type)    = %3.0f \n',est_type);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  Compute the low-rank-revealing URV decomposition using   \n');
fprintf(1,'.  the warm-started algorithm:                              \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.    [p,R,V,U,vec] = lurv(A,tol_rank,max_iter,num_ref, ...  \n');
fprintf(1,'.                                                est_type)  \n');
fprintf(1,'.                                                           \n');
input('.                                 [Press RETURN to continue]');
fprintf(1,'.                                                           \n');
fprintf(1,'.  wait...                                                  \n');
fprintf(1,'.                                                           \n');

% Compute URV.
flops(0);
t0 = cputime;

[p,R,V,U,vec] = lurv(A,tol_rank,max_iter,num_ref,est_type);

time = cputime - t0;
noflops = flops;

rnorm = norm(R(1:p,p+1:n));
theta = norm(V(:,1:p)'*v(:,p+1:n));
phi   = norm(U(:,1:p)'*u(:,p+1:m));

fprintf(1,'.                                                           \n');
fprintf(1,'.  Output results from LURV                                 \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    numerical rank p of A               = %3.0f \n',p);
fprintf(1,'.    upper bound of ||R(1:p,p+1:n)||     = %6.4e \n',vec(1));
fprintf(1,'.    estimate of pth singular value      = %6.4e \n',vec(2));
fprintf(1,'.    estimate of (p+1)th singular value  = %6.4e \n',vec(3));
fprintf(1,'.    estimate of num. nullspace angle    = %6.4e \n',vec(4));
fprintf(1,'.    estimate of num. range angle        = %6.4e \n',vec(5));
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  SVD based results                                        \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    numerical rank p of A               = %3.0f \n',numrank);
fprintf(1,'.    ||R(1:p,p+1:n)||                    = %6.4e \n',rnorm);
fprintf(1,'.    pth singular value                  = %6.4e \n',s1(numrank));
fprintf(1,'.    (p+1)th singular value              = %6.4e \n',s2(1));
fprintf(1,'.    numerical nullspace angle           = %6.4e \n',theta);
fprintf(1,'.    numerical range angle               = %6.4e \n',phi);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  Computational complexity of LURV                         \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    floating point operations           = %6.4e \n',noflops);
fprintf(1,'.    CPU time                            = %6.4e \n',time);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
input('.                                 [Press RETURN to continue]');
fprintf(1,'.                                                           \n');

% Define input parameters.
tol_rank = 0.001;
max_iter = 5;
num_ref  = 0;
est_type = 1;

fprintf(1,'.                                                           \n');
fprintf(1,'.  Input parameters to LURV_A                               \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    no. rows m of A                     = %3.0f \n',m);
fprintf(1,'.    no. cols n of A                     = %3.0f \n',n);
fprintf(1,'.    rank tolerance (tol_rank)           = %6.4e \n',tol_rank);
fprintf(1,'.    no. of steps of the principle               \n');
fprintf(1,'.    singular vector estimator (max_iter)= %3.0f \n',max_iter);
fprintf(1,'.    no. refinement steps (num_ref)      = %3.0f \n',num_ref);
fprintf(1,'.    use Lanczos procedure (est_type)    = %3.0f \n',est_type);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  Compute the low-rank-revealing URV decomposition using   \n');
fprintf(1,'.  the cold-started algorithm:                              \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.    [p,R,V,U,vec] = lurv_a(A,tol_rank,max_iter,...         \n');
fprintf(1,'.                                        num_ref,est_type)  \n');
fprintf(1,'.                                                           \n');
input('.                                 [Press RETURN to continue]');
fprintf(1,'.                                                           \n');
fprintf(1,'.  wait...                                                  \n');
fprintf(1,'.                                                           \n');

% Compute ULV.
flops(0);
t0 = cputime;

[p,L,V,U,vec] = lurv_a(A,tol_rank,max_iter,num_ref,est_type);

time = cputime - t0;
noflops = flops;

rnorm = norm(R(1:p,p+1:n));
theta = norm(V(:,1:p)'*v(:,p+1:n));
phi   = norm(U(:,1:p)'*u(:,p+1:m));

fprintf(1,'.                                                           \n');
fprintf(1,'.  Output results from LURV_A                               \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    numerical rank p of A               = %3.0f \n',p);
fprintf(1,'.    upper bound of ||R(1:p,p+1:n)||     = %6.4e \n',vec(1));
fprintf(1,'.    estimate of pth singular value      = %6.4e \n',vec(2));
fprintf(1,'.    estimate of (p+1)th singular value  = %6.4e \n',vec(3));
fprintf(1,'.    estimate of num. nullspace angle    = %6.4e \n',vec(4));
fprintf(1,'.    estimate of num. range angle        = %6.4e \n',vec(5));
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  SVD based results                                        \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    numerical rank p of A               = %3.0f \n',numrank);
fprintf(1,'.    ||R(1:p,p+1:n)||                    = %6.4e \n',rnorm);
fprintf(1,'.    pth singular value                  = %6.4e \n',s1(numrank));
fprintf(1,'.    (p+1)th singular value              = %6.4e \n',s2(1));
fprintf(1,'.    numerical nullspace angle           = %6.4e \n',theta);
fprintf(1,'.    numerical range angle               = %6.4e \n',phi);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');
fprintf(1,'.                                                           \n');
fprintf(1,'.  Computational complexity of LURV_A                      \n');
fprintf(1,'............................................................\n');
fprintf(1,'.    floating point operations           = %6.4e \n',noflops);
fprintf(1,'.    CPU time                            = %6.4e \n',time);
fprintf(1,'............................................................\n');
fprintf(1,'.                                                           \n');

%-----------------------------------------------------------------------
% End of function lurvdemo
%-----------------------------------------------------------------------