function [A, b] = matrixs(n)
e = ones(n, 1);
A = sparse(n, n);
i = 1:n;
a = spdiags([i'.*e], 0, n, n);
a1 = spdiags([1/2*e], -1, n, n);
a2 = spdiags([1/2*e], 1, n, n);
a3 = spdiags([1/2*e], -2, n, n);
a4 = spdiags([1/2*e], 2, n, n);
A = a + a1 + a2 + a3 + a4;
b = ones(n, 1);
end
n_values = [4, 10, 100, 1000, 10000, 100000];
for n = n_values
[A, b] = matrixs(n);
figure;
spy(A);
title(['Nonzero Matrix A for n = ', num2str(n)]);
figure;
spy(inv(full(A)));
title(['Inverse of Matrix A for n = ', num2str(n)]);
eigenvalues = eigs(A, 10, 'lm');
figure;
plot(real(eigenvalues), imag(eigenvalues), 'o');
title(['Eigenvalues of Matrix A for n = ', num2str(n)]);
% LU decomposition
tic;
[L, U, P] = lu(A);
lu_time = toc;
% Cholesky factorization
tic;
R = chol(A);
chol_time = toc;
fprintf('LU time for n = %d: %f seconds\n', n, lu_time);
fprintf('Cholesky time for n = %d: %f seconds\n', n, chol_time);
% Jacobi and Gauss-Seidel methods
tol = 1e-8;
max_iter = 1000;
% Jacobi method
x_jacobi = zeros(n, 1);
1
for iter = 1:max_iter
x_old = x_jacobi;
for j = 1:n
x_jacobi(j) = (b(j) - A(j, :) * x_old + A(j, j) * x_old(j)) /
A(j, j);
end
if norm(x_jacobi - x_old, 2) < tol
break;
end
end
% Gauss-Seidel method
x_gs = zeros(n, 1);
for iter = 1:max_iter
x_old = x_gs;
for j = 1:n
x_gs(j) = (b(j) - A(j, :) * x_gs + A(j, j) * x_gs(j)) / A(j, j);
end
if norm(x_gs - x_old, 2) < tol
break;
end
end
fprintf('Jacobi method l2-norm for n = %d: %f, iterations: %d\n', n,
norm(x_jacobi, 2), iter);
fprintf('Gauss-Seidel method l2-norm for n = %d: %f, iterations: %d\n',
n, norm(x_gs, 2), iter);
end
2
3
LU time for n = 4: 0.000093 seconds
Cholesky time for n = 4: 0.000058 seconds
Jacobi method l2-norm for n = 4: 0.892679, iterations: 12
Gauss-Seidel method l2-norm for n = 4: 0.892679, iterations: 12
4
5
LU time for n = 10: 0.000092 seconds
Cholesky time for n = 10: 0.000057 seconds
Jacobi method l2-norm for n = 10: 0.932109, iterations: 11
Gauss-Seidel method l2-norm for n = 10: 0.932109, iterations: 11
6
7
LU time for n = 100: 0.000113 seconds
Cholesky time for n = 100: 0.000066 seconds
Jacobi method l2-norm for n = 100: 0.968049, iterations: 11
Gauss-Seidel method l2-norm for n = 100: 0.968049, iterations: 11
8
9
LU time for n = 1000: 0.000386 seconds
Cholesky time for n = 1000: 0.000252 seconds
Jacobi method l2-norm for n = 1000: 0.972561, iterations: 11
Gauss-Seidel method l2-norm for n = 1000: 0.972561, iterations: 11
10
Warning: 6 of the 10 requested eigenvalues converged. Eigenvalues that did not converge are NaN.
11
LU time for n = 10000: 0.003092 seconds
Cholesky time for n = 10000: 0.001903 seconds
Jacobi method l2-norm for n = 10000: 0.973022, iterations: 11
Gauss-Seidel method l2-norm for n = 10000: 0.973022, iterations: 11
12
Error using full
Requested 100000x100000 (74.5GB) array exceeds maximum array size preference (5.0GB). This might
cause MATLAB to become unresponsive.
Related documentation
13
