1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
|
function test9
%TEST9 test cs_qr
%
% Example:
% test9
% See also: testall
% Copyright 2006-2007, Timothy A. Davis.
% http://www.cise.ufl.edu/research/sparse
rand ('state', 0) ;
index = UFget ;
[ignore f] = sort (max (index.nrows, index.ncols)) ;
f = f (1:100) ;
clf
% f = 185 ;
% f = 449 ;
% f = 186
for i = f
Prob = UFget (i) ;
disp (Prob) ;
A = Prob.A ;
[m n] = size (A) ;
if (m < n)
A = A' ;
end
[m n] = size (A) ;
sp = sprank (A) ;
% if (sprank (A) < min (m,n))
% continue
% end
Aorig = A ;
A = A (:, colamd (A)) ;
s1 = svd (full (A)) ;
tic ;
R = qr (A) ;
t1 = toc ; %#ok
% tic ;
% [Q,R] = qr (A) ;
% t1 = toc ;
[c,h,parent] = symbfact (A, 'col') ;
rnz = sum (c) ; %#ok
tic ;
[V2,Beta2,p,R2] = cs_qr (sparse(A)) ;
t2 = toc ; %#ok
v2 = full (V2) ;
if (any (spones (v2) ~= spones (V2)))
error ('got zeros!') ;
end
C = A ;
m2 = size (V2,1) ;
if (m2 > m)
C = [A ; sparse(m2-m, n)] ;
end
C = C (p,:) ;
% [H1,R1] = myqr (C) ;
% err1 = norm (R1-R2,1) / norm (R1)
% % [svd(A) svd(R1) svd(full(R2))]
% s2 = svd (full (R2)) ;
% err2 = norm (s1 - s2) / s1 (1)
% fprintf ('%10.6f %10.6f cs speedup %8.3f sprank %d n %d\n', ...
% t1, t2, t1/t2, sp, n) ;
% err2
% left-looking:
[V,Beta3,R3] = qr_left (C) ; %#ok
s3 = svd (full (R2)) ;
err3 = norm (s1 - s3) / s1 (1) ;
disp ('err3 = ') ; disp (err3) ;
if (err3 > 1e-12)
error ('!') ;
end
% right-looking:
[V,Beta4,R4] = qr_right (C) ; %#ok
s4 = svd (full (R2)) ;
err4 = norm (s1 - s4) / s1 (1) ;
disp ('err4 = ') ; disp (err4) ;
if (err4 > 1e-12)
error ('!') ;
end
% H2 = full (H2)
% R2 = full (R2)
subplot (2,4,1) ; spy (A) ; title ('A colamd') ;
subplot (2,4,2) ; spy (C) ; title ('A rperm') ;
subplot (2,4,3) ; treeplot (parent) ;
subplot (2,4,4) ; spy (Aorig) ; title ('Aorig') ;
subplot (2,4,5) ; spy (abs(R2)>0) ; title ('spqr R, no zeros') ;
subplot (2,4,6) ; spy (R) ; title ('matlab R') ;
subplot (2,4,7) ; spy (R2) ; title ('spqr R') ;
subplot (2,4,8) ; spy (V2) ; title ('spqr V') ;
drawnow
% if (err2 > 1e-9)
% error ('!') ;
% end
if (m2 > m)
fprintf ('added %d rows, sprank %d n %d\n', m2-m, sp, n) ;
end
% pause
end
|
