Skip to content

Commit c0e4dbe

Browse files
robolrobol
committed
'SVD' variant for spy in @hodlr and @hss.
1 parent d331301 commit c0e4dbe

File tree

2 files changed

+93
-12
lines changed

2 files changed

+93
-12
lines changed

‎@hodlr/spy.m

Lines changed: 48 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,15 @@
1-
function spy(H)
1+
function spy(H, variant)
22
%SPY Inspect the rank structure of H.
33
%
44
% SPY(H) draw a plot displaying the low-rank structure in the off-diagonal
55
% blocks of H.
6+
%
7+
% SPY(H, 'svd') also plots bars representing the singular values of the
8+
% low-rank blocks, to indicate the decay present in those blocks.
9+
10+
if ~exist('variant', 'var')
11+
variant = 'standard';
12+
end
613

714
m = size(H, 1);
815
n = size(H, 2);
@@ -16,15 +23,16 @@ function spy(H)
1623
axis([ 0 n 0 m ]);
1724
hold on;
1825

19-
spy_draw_block(H, 0, 0, hodlrrank(H));
26+
spy_draw_block(H, 0, 0, hodlrrank(H), variant);
2027

2128
hold off;
2229

2330
end
2431

25-
function spy_draw_block(H, xoffset, yoffset, rk)
32+
function spy_draw_block(H, xoffset, yoffset, rk, variant)
2633
diag_color = [ 0.1, 0.3, 0.95 ];
2734
offdiag_color = [ 0.95, 0.95, 1.0 ];
35+
svd_color = [1.0, 0.3, 0.4];
2836

2937
if is_leafnode(H)
3038
mm = size(H, 1);
@@ -35,10 +43,10 @@ function spy_draw_block(H, xoffset, yoffset, rk)
3543

3644
fill(xb, yb, diag_color);
3745
else
38-
spy_draw_block(H.A11, xoffset, yoffset, rk);
46+
spy_draw_block(H.A11, xoffset, yoffset, rk, variant);
3947
mm = size(H.A11, 1);
4048
nn = size(H.A11, 2);
41-
spy_draw_block(H.A22, xoffset + nn, yoffset + mm, rk);
49+
spy_draw_block(H.A22, xoffset + nn, yoffset + mm, rk, variant);
4250

4351
mm2 = size(H.A22, 1);
4452
nn2 = size(H.A22, 2);
@@ -48,21 +56,54 @@ function spy_draw_block(H, xoffset, yoffset, rk)
4856
xoffset + nn + nn2, xoffset + nn ];
4957
yb = [yoffset, yoffset, yoffset + mm, yoffset + mm ];
5058
thisrk = size(H.U12, 2);
51-
59+
5260
fill(xb, yb, offdiag_color);
61+
62+
if strcmp(variant, 'svd')
63+
% Draw the SVD plot of the block in the background
64+
[~, R] = qr(H.U12, 0); [~, S] = qr(H.V12, 0);
65+
s = log(svd(R * S'));
66+
nsvd = length(s);
67+
for j = 1 : nsvd
68+
sjm = (s(1) - s(j)) / (s(1) - s(end)) * mm;
69+
xbj = [xoffset + nn + nn2 * (j-1)/nsvd, xoffset + nn + nn2 * j / nsvd, ...
70+
xoffset + nn + nn2 * j / nsvd, xoffset + nn + nn2 * (j-1)/nsvd ];
71+
ybj = [yoffset + sjm, yoffset + sjm, yoffset + mm, yoffset + mm ];
72+
fill(xbj, ybj, svd_color, 'LineStyle','none');
73+
end
74+
end
75+
5376
text(xoffset + nn + nn2 / 2, yoffset + mm / 2, ...
5477
sprintf('%d', thisrk), ...
5578
'HorizontalAlignment', 'center');
56-
79+
5780
% Low rank block (2, 1)
5881
xb = [xoffset, xoffset + nn, xoffset + nn, xoffset ];
5982
yb = [yoffset + mm, yoffset + mm, ...
6083
yoffset + mm + mm2, yoffset + mm + mm2 ];
6184
thisrk = size(H.U21, 2);
6285

6386
fill(xb, yb, offdiag_color);
87+
88+
if strcmp(variant, 'svd')
89+
% Draw the SVD plot of the block in the background
90+
[~, R] = qr(H.U21, 0); [~, S] = qr(H.V21, 0);
91+
s = log(svd(R * S'));
92+
nsvd = length(s);
93+
for j = 1 : nsvd
94+
sjm = (s(1) - s(j)) / (s(1) - s(end)) * mm2;
95+
xbj = [xoffset + nn * (j-1)/nsvd, xoffset + nn * j / nsvd, ...
96+
xoffset + nn * j / nsvd, xoffset + nn * (j-1)/nsvd ];
97+
ybj = [yoffset + mm + sjm, yoffset + mm + sjm, yoffset + mm + mm2, yoffset + mm + mm2 ];
98+
fill(xbj, ybj, svd_color, 'LineStyle','none');
99+
end
100+
end
101+
64102
text(xoffset + nn / 2, yoffset + mm + mm2 / 2, ...
65103
sprintf('%d', thisrk), ...
66104
'HorizontalAlignment', 'center');
67105
end
106+
end
107+
108+
function spy_draw_svd_block()
68109
end

‎@hss/spy.m

Lines changed: 45 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,19 @@
1-
function spy(H)
1+
function spy(H, variant)
22
%SPY Inspect the rank structure of H.
33
%
44
% SPY(H) draw a plot displaying the low-rank structure in the off-diagonal
55
% blocks of H.
6+
%
7+
% SPY(H, 'svd') also plots bars representing the singular values of the
8+
% low-rank blocks, to indicate the decay present in those blocks.
9+
10+
if ~exist('variant', 'var')
11+
variant = 'standard';
12+
end
13+
14+
if strcmp('variant', 'svd')
15+
H = hss_proper(H);
16+
end
617

718
m = size(H, 1);
819
n = size(H, 2);
@@ -16,15 +27,16 @@ function spy(H)
1627
axis([ 0 n 0 m ]);
1728
hold on;
1829

19-
spy_draw_block(H, 0, 0, hodlrrank(H));
30+
spy_draw_block(H, 0, 0, hssrank(H), variant);
2031

2132
hold off;
2233

2334
end
2435

25-
function spy_draw_block(H, xoffset, yoffset, rk)
36+
function spy_draw_block(H, xoffset, yoffset, rk, variant)
2637
diag_color = [ 0.1, 0.3, 0.95 ];
2738
offdiag_color = [ 0.95, 0.95, 1.0 ];
39+
svd_color = [1.0, 0.3, 0.4];
2840

2941
if H.leafnode == 1
3042
mm = size(H, 1);
@@ -35,10 +47,10 @@ function spy_draw_block(H, xoffset, yoffset, rk)
3547

3648
fill(xb, yb, diag_color);
3749
else
38-
spy_draw_block(H.A11, xoffset, yoffset, rk);
50+
spy_draw_block(H.A11, xoffset, yoffset, rk, variant);
3951
mm = size(H.A11, 1);
4052
nn = size(H.A11, 2);
41-
spy_draw_block(H.A22, xoffset + nn, yoffset + mm, rk);
53+
spy_draw_block(H.A22, xoffset + nn, yoffset + mm, rk, variant);
4254

4355
mm2 = size(H.A22, 1);
4456
nn2 = size(H.A22, 2);
@@ -50,6 +62,20 @@ function spy_draw_block(H, xoffset, yoffset, rk)
5062
thisrk = rank(H.B12);
5163

5264
fill(xb, yb, offdiag_color);
65+
66+
if strcmp(variant, 'svd')
67+
% Draw the SVD plot of the block in the background
68+
s = log(svd(H.B12));
69+
nsvd = length(s);
70+
for j = 1 : nsvd
71+
sjm = (s(1) - s(j)) / (s(1) - s(end)) * mm;
72+
xbj = [xoffset + nn + nn2 * (j-1)/nsvd, xoffset + nn + nn2 * j / nsvd, ...
73+
xoffset + nn + nn2 * j / nsvd, xoffset + nn + nn2 * (j-1)/nsvd ];
74+
ybj = [yoffset + sjm, yoffset + sjm, yoffset + mm, yoffset + mm ];
75+
fill(xbj, ybj, svd_color, 'LineStyle','none');
76+
end
77+
end
78+
5379
text(xoffset + nn + nn2 / 2, yoffset + mm / 2, ...
5480
sprintf('%d', thisrk), ...
5581
'HorizontalAlignment', 'center');
@@ -61,6 +87,20 @@ function spy_draw_block(H, xoffset, yoffset, rk)
6187
thisrk = rank(H.B21);
6288

6389
fill(xb, yb, offdiag_color);
90+
91+
if strcmp(variant, 'svd')
92+
% Draw the SVD plot of the block in the background
93+
s = log(svd(H.B21));
94+
nsvd = length(s);
95+
for j = 1 : nsvd
96+
sjm = (s(1) - s(j)) / (s(1) - s(end)) * mm2;
97+
xbj = [xoffset + nn * (j-1)/nsvd, xoffset + nn * j / nsvd, ...
98+
xoffset + nn * j / nsvd, xoffset + nn * (j-1)/nsvd ];
99+
ybj = [yoffset + mm + sjm, yoffset + mm + sjm, yoffset + mm + mm2, yoffset + mm + mm2 ];
100+
fill(xbj, ybj, svd_color, 'LineStyle','none');
101+
end
102+
end
103+
64104
text(xoffset + nn / 2, yoffset + mm + mm2 / 2, ...
65105
sprintf('%d', thisrk), ...
66106
'HorizontalAlignment', 'center');

0 commit comments

Comments
 (0)