Script for Fast Fourier Transformation
Conversion of available experimental data from time domain to frequency domain.
Clearing workspace and command window-
clear;
clc;
Importing data.csv-
data = readtable("data.csv");
Spliting into vectors-
time = data{:,1};
V = data{:,2};
Determining dt, length of the vectors and sampling frequency-
l1 = length(V);
dt = (time(end)-time(1))/(l1-1);
fs = 1/dt;
FFT-
fft_V = fft(V, l1)*(2/l1);
abs_V = abs(fft_V);
freq = 0:(1/time(end)):fs/2-(1/time(end));
Ploting raw data (Given in time domain)-
plot(time, V)
xlabel('Time')
xlim([0 125])
ylabel('Heat Flux')
ylim([-10 50])
grid on
Ploting FFT-
figure
plot(freq, abs_V(1:length(freq)))
xlabel('Frequency')
ylim([-10 110])
xlim([-0.025 0.1])
grid on
Script for Gauss Elimination method
This matlab script can solve a system of linear equations by Gauss elimination method with partial pivoting. Fun fact is, this script can show the calculation steps. So, A lot less effort is needed while preparing assignments.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Title: Gauss Elimination with partial pivoting %
% Author: Dhiman Roy, Dept. of ME (BUET) %
% Date: August 23, 2021 %
% Licensed under Creative Commons %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
clc;
global showSteps;
showSteps = true; %can be true or false. true shows calculation steps
A = input('Please Enter the Co-efficient Matrix, A:\n' );
b = input('Please Enter the Constants(Must be a column vector, b:\n');
solByGEPP = gaussElimPartPiv(A, b)
function [x] = gaussElimPartPiv(A, b)
global showSteps;
[m, n] = size(A);
if m~=n
error("'A' must be a square matrix.");
else
aug = [A b];
if showSteps == 1
fprintf("Gauss elimination with partial pivoting.\nAugmented Matrix:\n");
disp(aug);
fprintf("\n");
end
for k = 1: n-1
%partial pivoting
for p = k+1:n
if abs(A(k,k) < A(p,k))
aug([k p],:) = aug([p k],:);
end
end
%forward elimination
for i = k+1: n
if showSteps == 1
fprintf("Row%d = Row%d - (%0.4f)/(%0.4f) * Row%d\n", i, i, aug(i,k), aug(k,k), k);
end
factor = aug(i,k)/aug(k,k);
aug(i,k:n+1) = aug(i,k:n+1)-factor*aug(k,k:n+1);
end
%display augmented matrix after each step
if showSteps == 1
fprintf("\n");
disp(aug)
fprintf("\n");
end
end
%back substitution
x = zeros(n,1);
x(n) = aug(n,n+1)/aug(n,n);
for i = n-1:-1:1
x(i) = (aug(i,n+1)-aug(i,i+1:n)*x(i+1:n))/aug(i,i);
end
end
end
Script for Gauss Jordan method
This matlab script can solve a system of linear equations by Gauss Jordan method with partial pivoting. Co-efficeint matrix and contant vector are be prompted for user input. Fun fact is, this script can show the calculation steps.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Title: Gauss Jordan with partial pivoting %
% Author: Dhiman Roy, Dept. of ME (BUET) %
% Date: August 23, 2021 %
% Licensed under Creative Commons %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
clc;
global showSteps;
showSteps = true; %can be true or false. true shows calculation steps
A = input('Please Enter the Co-efficient Matrix, A:\n' );
b = input('Please Enter the Constants(Must be a column vector, b:\n');
solByGJPP = gaussJordanPartPiv(A, b)
function [x] = gaussJordanPartPiv(A, b)
global showSteps;
[m, n] = size(A);
if m~=n
error("'A' must be a square matrix.");
else
aug = [A b];
if showSteps == 1
fprintf("Gauss Jordan with partial pivoting.\nAugmented Matrix:\n");
disp(aug);
fprintf("\n");
end
for k = 1: n-1
%partial pivoting
for p = k+1:n
if abs(A(k,k) < A(p,k))
aug([k p],:) = aug([p k],:);
end
end
%forward elimination
for i = k+1: n
if showSteps == 1
fprintf("Row%d = Row%d - (%0.4f)/(%0.4f) * Row%d\n", i, i, aug(i,k), aug(k,k), k);
end
factor = aug(i,k)/aug(k,k);
aug(i,k:n+1) = aug(i,k:n+1) - factor*aug(k,k:n+1);
end
%display augmented matrix after each step
if showSteps == 1
fprintf("\n");
disp(aug);
end
end
%elimination from upper triangular matrix
if showSteps == 1
fprintf("Elimination from upper triangular matrix:\n\n");
end
for k = n:-1:2
for i = k-1:-1:1
if showSteps == 1
fprintf("Row%d = Row%d - (%0.4f)/(%0.4f) * Row%d\n", i, i, aug(i,k), aug(k,k), k);
end
factor = aug(i,k)/aug(k,k);
aug(i,n+1:-1:i) = aug(i,n+1:-1:i) - factor*aug(k,n+1:-1:i);
end
if showSteps == 1
fprintf("\n");
disp(aug)
end
end
%converting to Identity Matrix
for i = 1:n
aug(i,:) = aug(i,:)/aug(i,i);
end
if showSteps == 1
fprintf("\nFinal form of the augmented matrix:\n\n");
disp(aug);
fprintf("\n");
end
end
x = aug(:,n+1);
end
Script for Gauss Seidel method
This matlab script can solve a system of linear equations by Gauss Seidel method, an iterative method. Co-efficeint matrix, contant vector, initial guesses and tolerence are be prompted for user input.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Title: Gauss Seidel method %
% Author: Dhiman Roy, Dept. of ME (BUET) %
% Date: August 23, 2021 %
% Licensed under Creative Commons %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
clc;
global showSteps;
showSteps = false; %true shows calculation steps
A = input('Please Enter the Co-efficient Matrix, A:\n' );
b = input('Please Enter the Constants(Must be a column vector, b:\n');
solByGaussSiedel = gaussSeidel(A, b)
function [x] = gaussSeidel(A, b)
global showSteps;
xIni = input("Enter initial guesses (Must be a column vector:\n");
eLimit = input("Enter tolerance:\n");
[m,n] = size(A);
for k = 1: n-1
for p = k+1:n
if abs(A(k,k) < A(p,k))
b([k p]) = b([p k]);
A([k p],:) = A([p k],:);
end
end
end
check = ones(n,1);
for i = 1:n
if 2*abs(A(i,i)) < sum(abs(A(i,:)))
check(i) = 0;
end
end
if all(check) ~= 1
error("Partial Pivoting does not make the coefficient matrix a diagonally dominant matrix.")
else
disp("Partial Pivoting makes the coefficient matrix a diagonally dominant matrix.")
end
x = xIni;
xOld = x;
for k = 1:100
err = 0;
for i = 1:n
x(i) = (b(i)-A(i,[1:i-1,i+1:n])*x([1:i-1,i+1:n]))/A(i,i);
err = err + abs((x(i)-xOld(i))/x(i));
xOld(i) = x(i);
end
if showSteps == 1
fprintf("Iteration %d\n",k);
disp(x);
disp(err);
end
if err < eLimit
break
end
end
end
Script for Jacobi method
Another iterative method for solving a system of linear equations.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Title: Jacobi method %
% Author: Dhiman Roy, Dept. of ME (BUET) %
% Date: August 23, 2021 %
% Licensed under Creative Commons %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
clc;
global showSteps;
showSteps = false; %true shows calculation steps
A = input('Please Enter the Co-efficient Matrix, A:\n' );
b = input('Please Enter the Constants(Must be a column vector, b:\n');
solByJacobi = jacobi(A, b)
function [x] = jacobi(A, b)
global showSteps;
xIni = input("Enter initial guesses (Must be a column vector:\n");
eLimit = input("Enter tolerance:\n");
[m,n] = size(A);
for k = 1: n-1
for p = k+1:n
if abs(A(k,k) < A(p,k))
b([k p]) = b([p k]);
A([k p],:) = A([p k],:);
end
end
end
check = ones(n,1);
for i = 1:n
if 2*abs(A(i,i)) < sum(abs(A(i,:)))
check(i) = 0;
end
end
if all(check) ~= 1
error("Partial pivoting does not make the coefficient matrix a diagonally dominant matrix.")
else
disp("Partial pivoting makes the coefficient matrix a diagonally dominant matrix.")
end
xOld = xIni;
xNew = zeros(n,1);
for k = 1:100
err = 0;
for i = 1:n
xNew(i) = b(i)/A(i,i)-A(i,[1:i-1,i+1:n])*(xOld([1:i-1,i+1:n])/A(i,i));
err = err + abs((xNew(i)-xOld(i))/xNew(i));
end
if showSteps == 1
fprintf("Iteration %d\n",k);
disp(xNew);
abs((xNew-xOld)/xNew);
end
if err < eLimit
break;
end
xOld = xNew;
end
x = xNew;
end
Script for SOR method
Solves a system of linear equations by successive over relaxation, prompts for co-efficient matrix, contant vector, initial guesses and relaxation factor.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Title: SOR method %
% Author: Dhiman Roy, Dept. of ME (BUET) %
% Date: August 23, 2021 %
% Licensed under Creative Commons %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
clc;
global showSteps;
showSteps = false; %true shows calculation steps
A = input('Please Enter the Co-efficient Matrix, A:\n' );
b = input('Please Enter the Constants(Must be a column vector, b:\n');
solBySOR = sor(A, b)
function [x] = sor(A, b)
global showSteps;
xIni = input("Enter initial guesses (Must be a column vector:\n");
eLimit = input("Enter tolerance:\n");
w = input("Enter Relaxation Factor:\n");
[m,n] = size(A);
for k = 1: n-1
for p = k+1:n
if abs(A(k,k) < A(p,k))
b([k p]) = b([p k]);
A([k p],:) = A([p k],:);
end
end
end
check = ones(n,1);
for i = 1:n
if 2*abs(A(i,i)) < sum(abs(A(i,:)))
check(i) = 0;
end
end
if all(check) ~= 1
error("Partial Pivoting does not make the coefficient matrix a diagonally dominant matrix.")
else
disp("Partial Pivoting makes the coefficient matrix a diagonally dominant matrix.")
end
x = xIni;
xOld = x;
for k = 1:100
err = 0;
for i = 1:n
x(i) = (1-w)*x(i) + w*(b(i)-A(i,[1:i-1,i+1:n])*x([1:i-1,i+1:n]))/A(i,i);
abs((x(i)-xOld(i))/x(i))
err = err + abs((x(i)-xOld(i))/x(i));
xOld(i) = x(i);
end
if showSteps == 1
fprintf("Iteration %d\n",k);
disp(x);
disp(err);
end
if err < eLimit
break
end
end
end