Refactored Direct Methods

This commit is contained in:
Sergiusz Warga 2021-03-21 19:45:18 +01:00
parent 37feb7c484
commit 8cfa2c9c18
7 changed files with 37 additions and 53 deletions

View File

@ -1,4 +1,4 @@
function A = Alg1_outer_product_gaussian_elimination(A) function A = Alg1(A)
% Algorithm 1: Outer Product Gaussian Elimination % Algorithm 1: Outer Product Gaussian Elimination
% Performs a gaussian eliminaion on a square matrix A. % Performs a gaussian eliminaion on a square matrix A.
@ -23,3 +23,4 @@ for k = 1 : m-1
end end
end end

View File

@ -1,5 +1,5 @@
function [Q R] = Alg11(A) function [Q R] = Alg11(A)
% Algorithm 11: QR factorization via Householder algorithm. % Algorithm 11: QR factorization via Householder transformation.
[m, n] = size(A); [m, n] = size(A);

View File

@ -1,4 +1,4 @@
function [P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A) function [P, Q, L, U] = Alg2(A)
% Algorithm 2: Gaussian Elimination with Complete Pivoting. % Algorithm 2: Gaussian Elimination with Complete Pivoting.
% [P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A) % [P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A)
% computes the complete pivoting factorization PAQ = LU. % computes the complete pivoting factorization PAQ = LU.

View File

@ -1,4 +1,4 @@
function b = Alg3_forward_substitution(L, b) function b = Al3(L, b)
% Algorithm 3: Forward Substitution % Algorithm 3: Forward Substitution
% b = Alg3_forward_substitution(L, b) overwrites b with the solution to % b = Alg3_forward_substitution(L, b) overwrites b with the solution to
% Lx = b. % Lx = b.

View File

@ -1,4 +1,4 @@
function b = Alg4_back_substitution(U,b) function b = Alg4(U,b)
% Argorithm 4: Back Substitution % Argorithm 4: Back Substitution
% b = Alg4_back_substitution(U,b) returns vetor b with solution to the % b = Alg4_back_substitution(U,b) returns vetor b with solution to the
% Ux = b. % Ux = b.
@ -17,9 +17,9 @@ if length(b) ~= m
error('Vector b has wrong length!') error('Vector b has wrong length!')
end end
if det(U) < eps % if det(U) < eps
error('Matrix is not nonsingular!') % error('Matrix is not nonsingular!')
end % end
% The following algorithm is based on the Algrotihm 3.1.2 from [2]. % The following algorithm is based on the Algrotihm 3.1.2 from [2].

View File

@ -1,9 +1,9 @@
function A = Alg5_gauss_jordan_elimination(A) function A = Alg5(A)
% Algorithm 5: Gauss-Jordan Elimination % Algorithm 5: Gauss-Jordan Elimination
% A = Alg5_gauss_jordan_elimination(A) performs Gauss-Jordan elimination % A = Alg5_gauss_jordan_elimination(A) performs Gauss-Jordan elimination
% on an augmented matrix A. % on an augmented matrix A.
[m, n] = size(A); [m, ~] = size(A);
for k = 1 : m for k = 1 : m
@ -12,7 +12,7 @@ for k = 1 : m
A(k, :) = row; A(k, :) = row;
for i = 1 : m for i = 1 : m
if i ~= k if i ~= k && A(i, k) ~= 0
A(i, :) = A(i, :)-(A(i, k))*row; A(i, :) = A(i, :)-(A(i, k))*row;
end end
end end

View File

@ -1,54 +1,28 @@
clear all;
clc;
A = [2, -1, 0, 0;
-1, 2, -1, 0;
0, -1, 2, -1;
0, 0, -1, 2];
b = [0;0;0;5];
B = Alg1_outer_product_gaussian_elimination(A);
U = triu(B);
x = Alg4_back_substitution(U, b);
%% Problem 1 %% Problem 1
clear all; clear all;
clc;
A = [2, -1, 0, 0; A = [2, -1, 0, 0;
-1, 2, -1, 0; -1, 2, -1, 0;
0, -1, 2, -1; 0, -1, 2, -1;
0, 0, -1, 2]; 0, 0, -1, 2];
b = [0;0;0;5]; b = [0;0;0;5];
B = Alg5([A b])
B = gauss_jordan_elimination([A b])
[P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A);
b = P*b;
% Ly = b and Ux = y
y = forward_substitution(L, b);
x = Q*back_substitution(U, y);
% L*U
%% Problem 2 %% Problem 2
clear all;
clc;
A = [1, 1, 1; A = [1, 1, 1;
1, 1, 2; 1, 1, 2;
1, 2, 2]; 1, 2, 2];
b = [1;2;1]; b = [1;2;1];
[P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A); [P, Q, L, U] = Alg2(A);
b = P*b; b = P*b;
% Ly = b and Ux = y % Ly = b and Ux = y
y = Alg3_forward_substitution(L, b); y = Alg3(L, b); % Forward substitution
x = Q*Alg4(U, y) % Backward substitution
x = Q*Alg4_back_substitution(U, y)
% L*U
%% Problem 4 %% Problem 4
@ -60,25 +34,34 @@ bp = [0.168; 0.066];
kappa = cond(A) kappa = cond(A)
B = Alg5_gauss_jordan_elimination([A b]) B = Alg5([A b])
Bp = Alg5_gauss_jordan_elimination([A bp]) Bp = Alg5([A bp])
%% Problem 5 %% Problem 5
% AX = I3 clear all;
clc;
A = [2, 1, 2; A = [2, 1, 2;
1, 2, 3; 1, 2, 3;
4, 1, 2]; 4, 1, 2];
[P, Q, L, U] = Alg2_gaussian_elimination_with_complete_pivoting(A); [P, Q, L, U] = Alg2(A)
I = P*eye(3); I = P*eye(3);
% Ly = b and Ux = y % Ly = b and Ux = y
y = forward_substitution(L, I); y = Alg3(L, I); % Forward substitution
x = Q*Alg4(U, y) % Backward substitution
X = Q*back_substitution(U, y) inv(A) - x
inv(A)
%% Problem 6 %% Problem 6
clc;
A = [1 2 3 4;
-1 1 2 1;
0 2 1 3;
0 0 1 1];
[L, U, P] = Alg8(A)
det(A) - prod(diag(U))
%% Problem 10 %% Problem 10