Sparse Matrix Computation with spdiags in MATLAB
100 likes | 134 Views
Learn how to manipulate sparse matrices efficiently using the spdiags function in MATLAB. Download code examples and explore converting sparse matrices to reduced forms.
Sparse Matrix Computation with spdiags in MATLAB
E N D
Presentation Transcript
Lecture 13 Sparse matrix computation Matlab function: spdiags Download spdiags01.m And spdiags02.m
n=5; e = ones(n,1); A = spdiags([e -2*e 3*e], -1:1, n, n) full(A) >> spdiags01 A = (1,1) -2 (2,1) 1 (1,2) 3 (2,2) -2 (3,2) 1 (2,3) 3 (3,3) -2 (4,3) 1 (3,4) 3 (4,4) -2 (5,4) 1 (4,5) 3 (5,5) -2 ans = -2 3 0 0 0 1 -2 3 0 0 0 1 -2 3 0 0 0 1 -2 3 0 0 0 1 -2
n=5; e1 = [1; 2; 3; 4; 5]; e2=[-1; -2; -3; -4; -5]; e3=[11; 12; 13; 14; 15]; A = spdiags([e1 e2 e3], -1:1, n, n); full(A) >> spdiags02 … ans = -1 12 0 0 0 1 -2 13 0 0 0 2 -3 14 0 0 0 3 -4 15 0 0 0 4 -5 full(A) – display matrix in usual format
Define diagonals –2, 0 and 1 A = spdiags([e1 e2 e3], [-2 0 1], n, n) >> full(A) ans = -1 12 0 0 0 0 -2 13 0 0 1 0 -3 14 0 0 2 0 -4 15 0 0 3 0 -5
Convert from sparse matrix to reduced form Taking into account nonzero diagonals only [B,d] = spdiags(A) extracts all nonzero diagonals from the m-by-n matrix A. B is a min(m,n)-by-p matrix whose columns are the p nonzero diagonals of A. d is a vector of length p whose integer components specify the diagonals in A. Example: >> [B d]=spdiags(A) B = 1 -1 0 2 -2 12 3 -3 13 0 -4 14 0 -5 15 d = -2 0 1
The same in standard representation of matrix >> A0=full(A) A0 = -1 12 0 0 0 0 -2 13 0 0 1 0 -3 14 0 0 2 0 -4 15 0 0 3 0 -5 >> [B1,d1]=spdiags(A0) B1 = 1 -1 0 2 -2 12 3 -3 13 0 -4 14 0 -5 15 d1 = -2 0 1
If a column of B is longer than the diagonal it's replacing, spdiags takes elements of super-diagonals from the lower part of the column of B, and elements of sub-diagonals from the upper part of the column of B. >> full(A) ans = -1 12 0 0 0 0 -2 13 0 0 1 0 -3 14 0 0 2 0 -4 15 0 0 3 0 -5 B = 1 -1 0 2 -2 12 3 -3 13 0 -4 14 0 -5 15
Inclass 1 Modify Use spdiags to create the following sparse matrix: 1 2 0 0 0 0 0 0 0 0 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1 2 0 0 0 0 0 0 0 -1 0 1
Newton’s method for system of two nonlinear equations Download newton2d01.m Function newton2d01 Set initial (x,y) point, maximum number of iterations, and convergence tolerance. for i= 1:maxiter Perform Newton update; check for conver-gence; print iterates end Define functions f1(x,y), and f2(x,y) and JacobianJacob(x,y).
function newton2d01 %solution of system of two nonlienar equations % f1(x,y) = exp(x)*y = 0 % f2(x,y) = x*cos(y) = 0 % using a Newton iteration. Calls function f1 and f2 and % Jacobian matrix Jacob. % % Setup x = 0.2; y = .5; maxiter = 25; tol = 1.e-6; err = inf; iter = 0; fprintf('iter = %3.0f, x = %6.4f, y = %6.4f, \n', 0, x, y); % Newton update. while err > tol & iter <= maxiter J = Jacob(x,y); rhs = -[f1(x,y);f2(x,y)]; dxy = J\rhs; x = x + dxy(1); y = y + dxy(2); err = norm([f1(x,y);f2(x,y)],Inf); iter = iter+1; disp(['Iteration number = ',num2str(iter),' x=',num2str(x),' y=',num2str(y),' Backward error=',num2str(err)]); end function z = f1(x,y) z = exp(x)*y; function z = f2(x,y) z = x*cos(y); function xx = Jacob(x,y) xx = [exp(x)*y, exp(x); cos(y), -x*sin(y)];