Matrix Computation
2022 Fall
Benchmark problems
Laplacian matrix of graphs
The Koblenz Network Collection, Stanford Large Network Dataset Collection
MATLAB (un-direct case), construct_Laplacian_matrix.m
E = load('butterfly.txt');
n = max(max(E));
m = size(E,2);
if (m == 2 )
mtx_L = sparse(E(:,1), E(:,2), 1, n, n);
else
mtx_L = sparse(E(:,1), E(:,2), E(:,3), n, n);
end
D_L = diag(mtx_L);
if ( ~isempty(find(D_L ~= 0, 1)) )
mtx_L = mtx_L - sparse(1:n,1:n,D_L);
end
mtx_ L = mtx_L + mtx_L';
D = diag(sum(mtx_L, 2));
mtx_L = D - mtx_L;
G = graph(E(:,1), E(:,2));
[bins,binsizes] = conncomp(G);
idx_first_graph = binsizes(bins) == max(binsizes);
mtx_L = mtx_L(idx_first_graph,idx_first_graph);
SG = subgraph(G, idx_first_graph);
plot(SG)
%
[V, D] = eig(full(mtx_L));
trimesh(E,V(:,2),V(:,3),[],V(:,2),'LineWidth',3);
hold on; pIDX = V(:,2)>0; nIDX = V(:,2)<0;
plot(V(pIDX,2),V(pIDX,3),'ko','MarkerFaceColor','r');
plot(V(nIDX,2),V(nIDX,3),'ko','MarkerFaceColor','b');
view([0,0,1]); axis off;
Matrix dimension:
Social/loc-brightkite_edges, dimension = 58,228
Social/petster-friendships-cat, dimension = 149,700
Social/loc-gowalla_edges, dimension = 196,591
coauthorship/DBLP co-authorship, dimension = 317,080
com-amazon, dimension = 334,863
Social/petster-friendships-dog, dimension = 426,820
Social/Catster_Dogster, dimension = 623,766
Infrastructure/Pennsylvania, dimension = 1,088,092
Social/Youtube friendship, dimension = 1,134,890
Infrastructure/Texas, dimension = 1,379,917
Social/hyves, dimension = 1,402,673
Computer/Skitter/out.as-skitter, dimension = 1,696,415
Social/flickr-links, dimension = 1,715,255
Infrastructure/California, dimension = 1,965,206
Social/flixster, dimension = 2,523,386
Social/Orkut, dimension = 3,072,441
Social/YouTube, dimension = 3,223,589
Authorship/dblp-author, dimension = 4,000,150
Social/LiveJournal links, dimension = 5,204,176
Power and inverse power methods
Power method (PowerMethod_Norm.m, PowerMethod_LinFun.m)
Inverse power method (InversePower.m, InversePower_var_shift.m, InversePower_RayQuo.m)
MATLAB (solving sparse linear system)
dim_n = size(mtx_L,1);
sigma = 1.0e-3;
mtx_shift_L = mtx_L - sigma * speye(dim_n);
mtx_factor.Perm_amd_vec = amd(mtx_shift_L);
W_reorder = mtx_shift_L(mtx_factor.Perm_amd_vec, mtx_factor.Perm_amd_vec);
mtx_factor.mtx_Perm_amd = sparse(mtx_factor.Perm_amd_vec,1:dim_n,ones(dim_n,1));
[mtx_factor.mtx_Low_L,mtx_factor.mtx_upper_U,mtx_factor.mtx_Perm_LU] = lu(W_reorder);
sol = mtx_factor.mtx_Perm_amd * (mtx_factor.mtx_upper_U \ (mtx_factor.mtx_Low_L \ (mtx_factor.mtx_Perm_LU * rhs(mtx_factor.Perm_amd_vec,:))));
Krylov method for solving eigenvalue problem
run_Lanczos_solver.m, Lanczos_basic.m, Lanczos_with_IR_Reorth.m, Lanczos_Steps_with_Reorth.m, Implicit_Restart_for_Lanczos.m
MATLAB ( eigs, Download )
ncv = 10;
[EV, EW] = eigs(@(x)Shift_A_inv_b(x,mtx_factor), dim_n, ncv, 'largestabs'); %,'SubspaceDimension',45);
EW = sigma + 1./diag(EW);
[~,idx_sort] = sort(abs(EW));
EW = EW(idx_sort,1);
EV = EV(:,idx_sort);
function [sol] = Shift_A_inv_b(rhs,mtx_factor)
sol = mtx_factor.mtx_Perm_amd * (mtx_factor.mtx_upper_U \ (mtx_factor.mtx_Low_L \ (mtx_factor.mtx_Perm_LU * rhs(mtx_factor.Perm_amd_vec,:))));
end
Conjugate gradient method
Lecture note (download)
Benchmark problem (Photonic crystal, 2D Laplacian PDE)
Video: (Prof. Wen-Wei Lin)
Preconditioning conjugate gradient method
Part I (page 1~9)
Part II (page 9~17, SSOR)
Part III (page 35~45, Chebychev Semi-Iteration, page 16~19, Perron-Frobenius)
Part IV (page 19~20, Perron-Frobenius)
Generalized Conjugate Gradient Methods for Non-symmetric Linear Systems
GMRES: Generalized Minimal Residual Algorithm for Solving Non-symmetric Linear Systems
Preconditioning
slide ( 2010, 2020 )
Video: (Prof. Wen-Wei Lin)