unh2010:iam931:hw4soln
====== Differences ====== This shows you the differences between two versions of the page.
| Both sides previous revision Previous revision Next revision | Previous revision | ||
|
unh2010:iam931:hw4soln [2010/10/28 19:48] gibson |
unh2010:iam931:hw4soln [2010/10/29 08:14] (current) gibson |
||
|---|---|---|---|
| Line 1: | Line 1: | ||
| - | {{:unh2010:iam931:svderrs.png?400}} {{:unh2010:iam931:qrerrs.png?400}} | ||
| - | <code> | ||
| - | % A matlab program that constructs a matrix via SVD factors and then compares those to the computed SVD | ||
| - | m = 50; % dimension of matrices | ||
| - | S = 100; % number of tests | ||
| - | for s = 1:S | ||
| - | |||
| - | % construct A from random SVD factors | ||
| - | [U, tmp] = qr(randn(m,m)); | ||
| - | [V, tmp] = qr(randn(m,m)); | ||
| - | S = diag(sort(rand(m,1).^6, 'descend')); | ||
| - | A = U*S*V'; | ||
| - | | ||
| - | % compute SVD of A | ||
| - | [U2,S2,V2] = svd(A); | ||
| - | |||
| - | % change signs of cols of U2 if they don't match U, ditto V | ||
| - | for k=1:m | ||
| - | uk = U(:,k); | ||
| - | if uk'*U2(:,k) < 0 | ||
| - | U2(:,k) = -1*U2(:,k); | ||
| - | V2(:,k) = -1*V2(:,k); | ||
| - | end | ||
| - | end | ||
| - | kappa = cond(A); | ||
| - | loglog(kappa, norm(U-U2), 'r.') | ||
| - | hold on | ||
| - | loglog(kappa, norm(V-V2), 'b.') | ||
| - | loglog(kappa, norm(S-S2), 'k.') | ||
| - | loglog(kappa, norm(A-U2*S2*V2'), 'g.'); | ||
| - | end | ||
| - | |||
| - | legend('norm(U-U2)','norm(V-V2)','norm(S-S2)', 'norm(A-U2*S2*V2)''', ... | ||
| - | 'Location', 'NorthWest'); | ||
| - | legend boxoff | ||
| - | xlabel('cond(A)') | ||
| - | axis([1 10^20 10^-16 1]) | ||
| - | set(gcf, 'color', [1 1 1]) | ||
| - | set(gcf,'inverthardcopy','off') | ||
| - | |||
| - | title('SVD errors, A has singvals (0<random<1).^6') | ||
| - | </code> | ||
| - | |||
| - | To modify this program to do similar tests on QR decomposition, substitute in this code segment | ||
| - | <code> | ||
| - | R = triu(randn(m,m)); | ||
| - | [Q, tmp] = qr(randn(m,m)); | ||
| - | |||
| - | for k=1:m | ||
| - | if R(k,k) < 0 | ||
| - | R(k,k) = -1*R(k,k); | ||
| - | end | ||
| - | end | ||
| - | </code> | ||
| - | and appropriate changes in the plotting and labeling. | ||
unh2010/iam931/hw4soln.1288320490.txt.gz · Last modified: 2010/10/28 19:48 by gibson
Page Tools
Except where otherwise noted, content on this wiki is licensed under the following license: CC Attribution-Share Alike 3.0 Unported
