Autystyczne notatki

MATLAB + Simulink

11 min read

Octave or Scilab - FOSS alternatives to MATLAB

You can dock and undock each menu element (To undo layout changes: Home → Layout → Default).

MATLAB is case-sensitive (A ≠ a)

Warning!

MATLAB uses the dot (.) as decimal separator, not the comma (,). Decimals are written like this: 5.7, not like this: 5,7.

Confirm commands by pressing the Enter key. The results of commands ending with a semicolon (;) won’t show up in the command window. You can write multiple commands in the same line, provided that you separate them with a semicolon. The ans variable holds the answer/result of operations for which a variable wasn’t defined and changes after each command.

Basic arithmetic functions: +, -, *, /, ^. Also, sqrt, which by default calculates imaginary numbers (try doing sqrt(-1) in MATLAB).

Basic commands: who - list variables in workspace. whos - list variables in workspace, alongside their datatypes and size. clear - remove variables from workspace (can be used for singular variables or all of them). clc - clear the command window. plot - used for plotting graphs. help - try using the help command to see what this command does.

You’ll want to start most scripts with the clear and clc commands, in this order: clear all;clc;close all

You can make comments by preceding them with a % sign or multiple lines of comments between %{ and }%.

Table/Matrix/Array:

Syntax: {table_name}=[{row1_element1} {row1_element2}; {row2_element1} {row2_element2}] E.g. A=[1 2; 3 4]. You can also use the comma (,) to separate elements in a row, like this: A=[1, 2; 3, 4].

Difference between table and matrix operations:

The . is the table operator and is used right before an arithmetic operation, like this: C=A.*B, D=C.^B. There is not much difference for the + and - operations, but:

  • When doing multiplication (*), the following table operation: C=A.*B will multiply the A and B tables’ elements that are in the same row and column, while a matrix operation: C=A*B will follow the matrix multiplication rules. No difference when multiplying by a constant.
  • Similarly, when doing exponentiation (^), the following table operation: C=A.^2 will raise each element of the table A to the power of 2, while a matrix operation: C=A^2 will be equivalent to multiplying the A matrix by itself.

Other methods of creating arrays (tables/matrices):

  • a=zeros(3) - creates a 3 by 3 array filled with zeroes.
  • b=ones(3) - creates a 3 by 3 array filled with ones.
  • c=ones(3)*7 - creates a 3 by 3 array filled with sevens.
  • d=ones(3,4) - creates a 3 by 4 array filled with ones.
  • e=rand(6) - creates a 6 by 6 array with randomized elements between 0 and 1.
  • f=eye(3) - creates an identity matrix (all ones on the diagonal).
  • g=magic(5) - creates a magic matrix (sum of all rows, columns and diagonals is the same)
  • i=50*rand(3) - creates an array with randomized elements between 0 and 50.
  • j=150*rand(3)-50 - creates an array with randomized elements between -50 and 100 (alternatively: j=randi([-50,100], 3)).

Matrix transposition: ' (rows become columns), e.g. n=m’.

  • Matrix inversion: inv(m). If you multiply a matrix by its inverse, you obtain an identity matrix.

You can select an array element this way: {array_name}({row},{column}). E.g. element=A(2,3). Note that the first position is the row and the second is the column. You can also replace a specific element in the array, e.g. A(2,3)=50. By utilizing the “all set elements” operator (:), we can:

  • Choose all elements of a row, like this: A(2,:). Understand it like this: A(second row, all elements (columns)).
  • Choose all elements of a column, like this: A(:,4). Understand it like this: A(all elements (rows), fourth column).
  • Choose all elements of a range of columns/rows, like this: A(3:4,2:4). Indexation in MATLAB starts from 1.
Vectors:

The default step (h) is equal to 1. E.g. x=5:15 will create the vector: x=[5,6,7,8,9,10,11,12,13,14,15]. You can declare a different step value, e.g. x=5:3:15, where 3 is the step value h.

Vector arithmetic:

  • Choose a specific vector element, e.g. x(4) is the fourth element of the x vector.
  • Replace vector element/assign to variable, e.g. x(4)=12, var=x(4).
  • Remove vector element, e.g. x(4)=[].
  • Choose the last vector element, e.g. x(end).
  • Adding an element to the end of the vector, e.g. x(end+1)=17.
  • Calculate the length of the vector (number of elements), e.g. k=length(x).
  • Choose a range of vector elements, e.g. x(1:5). This way you’re de facto choosing a sub-vector from the vector x, i.e. x(4:8) <=> x([4 5 6 7 8]). You can use this to choose any vector elements, e.g. x([2 5 11]) selects elements indexed 2, 5, and 11.
  • You can make use of steps to create new vectors from pre-existing ones with a specific pattern, e.g. d=x(1:2:end) creates a vector d that contains every second element of vector x.
  • You can concatenate vectors, e.g. z=[x y] and also add values to them, like this:zz=[3 11 21 x(3:8) 4:8 12 y(1:5)].
Common errors:
  • Mistaking an element’s index with its value.
  • Chosen index exceeds the number of array elements.
  • Using the comma as decimal separator.
  • Indexing from zero or negative numbers/non-integers.
  • Unpaired parentheses/incorrect parentheses type.
  • Forgot to close a loop with end (The IF statement also requires an end).
First program:

M-files are specific to MATLAB.

To create a new script, go to: Home → New (Script), or press CTRL+N. You can then run the script by right-clicking it and choosing the “Run section” option.

File naming shenanigans:

  • File names can’t start with a number.
  • File names can’t contain spaces or dots, also avoid using #, %, $, etc.
  • Do not use MATLAB-specific commands (e.g. length, for, while, plot) to name your files.

Pro tip: Use camelCase for naming files and variables! Also use names instead of single-letter variables, e.g.: Alice, Bob, Charlie, and so on.

Loops and conditionals:
  • end <=> } in other programming languages.
  • Use x=x+1 instead of x++ (as in C).
  • There is no built-in do while loop in MATLAB.
FOR loop:

Usually used when you know how many times the loop will be executed, e.g.:


clear all
clc
close all
for k=1:10
	x(k)=50-k
end
WHILE loop:

Usually used when you don’t know how many times the loop will be executed, e.g.:


clear all
clc
close all
count=0
a=0
while a<=999
	a=1000*rand(1)
	count=count+1
end
IF condition:

clear all
clc
x=rand(1,20)
for k=1:20
	if x(k) <= 0.3
	x(k)=0;
	elseif x(k) > 0.8
	x(k)=1;
	else
	x(k)=x(k);
	end
end
disp(x)
Multiplication table:

clear all
clc
A=zeros(10);
for k=1:10
	for n=1:10
		A(k,n)=k*n;
	end
end
disp(A)

Understand it like this: The internal loop for n=1:10 is executed in its entirety for each time the external loop for k=1:10 gets executed once. In other words, we start with the first row (k=1) and for that row, we multiply it by columns 1 to 10, like this: row 1 * column 1, row 1 * column 2, row 1 * column 3, ... row 2 * column 1, row 2 * column 2, row 2 * column 3, ...

Ordered elements table:

clear all
clc
A=zeros(10);
a=0;
count=1;
for k=1:10
	for n=1:10
		A(k,n)=n+a;
	end
	a=n*count;
	count=count+1;
end
disp(A)
Chessboard:

clear all
clc
A=zeros(10);
range=1;
count=1;
for k=1:10
	%{we want to switch the range of every second row, which happens each time the count is divisible by 2:}%
	if mod(count,2) == 0
		range=2;
	else
		range=1;
	end
	for n=range:2:(range+9)
		if A(k,n) == 0
			A(k,n)=1;
		end
	end
	count=count+1;
end
disp(A)
Quadratic equations:

clear all;
close all;
clc;

%{requesting user input:}%
a=input("Podaj współczynnik przy najwyższej potędze: ")
b=input("Podaj współczynnik przy pierwszej potędze: ")
c=input("Podaj wyraz wolny: ")

delta=b^2-4*a*c;

if delta==0
	x0=(-1*b)/(2*a)
elseif delta>0
	x1=(-1*b+sqrt(delta))/(2*a)
	x2=(-1*b-sqrt(delta))/(2*a)
elseif delta<0
	x1=(-1*b+i*sqrt(-1*delta))/(2*a)
	x2=(-1*b-i*sqrt(-1*delta))/(2*a)
end
Approximation and interpolation:
Polynomial approximation:

clear all;
clc;
close all;

%generating points:
x=-50:20:70;
y=10*rand(1,length(x));

%{making the graph smoother:}%
xx=-50:70;

%creating the matrices:
Y=y';

%1st. deg.:
X=[ones(length(x),1), x'];
A=inv(X'*X)*X'*Y;
w1=A(1)+A(2)*xx;

%2nd. deg.:
X=[ones(length(x),1), x', x'.^2];
A=inv(X'*X)*X'*Y;
w2=A(1)+A(2)*xx+A(3)*(xx.^2);

%4th. deg.:
X=[ones(length(x),1), x', x'.^2, x'.^3, x'.^4];
A=inv(X'*X)*X'*Y;
w3=A(1)+A(2)*xx+A(3)*(xx.^2)+A(4)*(xx.^3)+A(5)*(xx.^4);

%6th. deg.:
X=[ones(length(x),1), x', x'.^2, x'.^3, x'.^4, x'.^5, x'.^6];
A=inv(X'*X)*X'*Y;
w4=A(1)+A(2)*xx+A(3)*(xx.^2)+A(4)*(xx.^3)+A(5)*(xx.^4)+A(6)*(xx.^5)+A(7)*(xx.^6);

%{plotting the functions, with descriptions:}%
plot(x,y,"*",xx,w1,"red",xx,w2,"green",xx,w3,"blue",xx,w4,"magenta")
legend("Węzły aproksymacji", "Aproksymacja funkcją liniową","Aproksymacja wielomianem drugiego stopnia","Aproksymacja wielomianem czwartego stopnia","Aproksymacja wielomianem szóstego stopnia")
grid minor
Approximation with MATLAB functions:

clear all;
clc;
close all;

%generating points:
x=-50:20:70;
y=10*rand(1,length(x));

%making the graph smoother:
xx=-50:70;

%MATLAB functions:
p=polyfit(x,y,2);
f=polyval(p,xx);

%creating the matrices:
Y=y';

%2nd. deg.:
X=[ones(length(x),1), x', x'.^2];
A=inv(X'*X)*X'*Y;
w2=A(1)+A(2)*xx+A(3)*(xx.^2);

plot(x,y,"*",xx,w2,"red",xx,f,"--")
legend("Węzły aproksymacji", "Aproksymacja średniokwadratowa","Aproksymacja funkcjami MATLAB'a")
grid minor
Approximation with other base functions:

clear all;
clc;
close all;

%generating points:
x=125:125:625;
y=[80, 36.56, 26.34, 24.33, 20.23];

%making the graph smoother:
xx=125:625;

%creating the matrices:
Y=y';

%Logarithmic, 1st. deg.:
X=[ones(length(x),1), (log(x))'];
A=inv(X'*X)*X'*Y;
w1=A(1)+A(2)*log(xx);

%Logarithmic, 2nd. deg.:
X=[ones(length(x),1), (log(x))', (log(x))'.^2];
A=inv(X'*X)*X'*Y;
w2=A(1)+A(2)*log(xx)+A(3)*(log(xx).^2);

%Rational, 1st. deg.:
X=[ones(length(x),1), (1./x)'];
A=inv(X'*X)*X'*Y;
w3=A(1)+A(2)*(1./xx);

%Rational, 2nd. deg.:
X=[ones(length(x),1), (1./x)', (1./x)'.^2];
A=inv(X'*X)*X'*Y;
w4=A(1)+A(2)*(1./xx)+A(3)*((1./xx).^2);

plot(x,y,"*",xx,w1,"red",xx,w2,"green",xx,w3,"blue",xx,w4,"magenta")
legend("Węzły aproksymacji", "Aproksymacja funkcją logarytmiczną pierwszego stopnia","Aproksymacja funkcją logarytmiczną drugiego stopnia","Aproksymacja funkcją wymierną pierwszego stopnia","Aproksymacja funkcją wymierną drugiego stopnia")
grid minor
Calculating the absolute error, Root Mean Square Error (RMSE) and the coefficient of determination (R^2):

clear all;
clc;
close all;

%generating points:
x=125:125:625;
y=[80, 36.56, 26.34, 24.33, 20.23];

%making the graph smoother:
xx=125:625;

%creating the matrices:
Y=y';
%Logarithmic, 1st. deg.:
X=[ones(length(x),1), (log(x))'];
A=inv(X'*X)*X'*Y;
yy=A(1)+A(2)*log(xx);

plot(x,y,"*",xx,yy,"red")
legend("Węzły aproksymacji", "Aproksymacja funkcją logarytmiczną pierwszego stopnia")
grid minor

%defining the vector step:
h=125;

%calculating errors:
D=yy(1:h:end)-y
RMSE=sqrt(sum(D.^2)/length(D))
R=1-(sum(D.^2)/sum((y-mean(y)).^2));
disp("R^2 = "); disp(R)
Non-linear equations:

Given equation:


xx=0:0.1:1.5;
f=xx-cos(xx);
y=xx;
g=cos(xx);

figure(1)
plot(xx,f,"red")
grid minor

figure(2)
plot(xx,y,"green",xx,g,"blue")
grid minor
Iterative method:

clear all;
close all;
clc;

%initial values:
x(1)=1;
k=1;

%{keep looping until condition is met:}%
while true
	x(k+1)=cos(x(k));
	if abs(x(k+1)-x(k)) < 0.001
		%display results:
		x
		k
		h=abs(x(k+1)-x(k))
		break
	else
		k=k+1;
	end
end
Relaxation method:

clear all;
close all;
clc;

format short%long

x(1)=1;
k=1;
a=sin(0.75);

while true
	x(k+1)=(a*x(k)+cos(x(k)))/(a+1);
	if abs(x(k+1)-x(k)) < 0.001
		%display results:
		x
		k
		h=abs(x(k+1)-x(k))
		break
	else
		k=k+1;
	end
end
Newton’s method:

clear all;
close all;
clc;

format long

x(1)=1;
k=1;

while true
	h(k)=-1*(x(k)-cos(x(k)))/(1+sin(x(k)));
	x(k+1)=x(k)+h(k);
	if abs(h(k)) < 0.001
		%display results:
		x
		k
		precision=h(k)
		break
	else
		k=k+1;
	end
end
Secant method:

clear all;
close all;
clc;

format long

x(1)=0;
x(2)=1;
k=2;

while true
	h(k)=-1*(x(k)-cos(x(k)))*(x(k)-x(k-1))/((x(k)-cos(x(k)))-(x(k-1)-cos(x(k-1))));
	x(k+1)=x(k)+h(k);
	if abs(h(k)) < 0.001
		%display results:
		x
		k
		precision=h(k)
		break
	else
		k=k+1;
	end
end
Numerical integration:
Left rectangle rule:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
y=0.5*x.^3-5*x.^2+80;
n=(b-a)/h;

%plotting the function:
xx=a:b;
yy=0.5*xx.^3-5*xx.^2+80;
plot(xx,yy)
grid on

%{left rectangle rule with step h:}%
Rl=0;
for k=1:n
    Rl=Rl+y(k);
end
Rl=Rl*h

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_Rl=100*(Rl-Q)/Q
Right rectangle rule:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
y=0.5*x.^3-5*x.^2+80;
n=(b-a)/h;

%{right rectangle rule with step h:}%
Rp=0;
for k=2:(n+1)
    Rp=Rp+y(k);
end
Rp=Rp*h

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_Rp=100*(Rp-Q)/Q
Midpoint rule:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
yc=0.5*(x+0.5*h).^3-5*(x+0.5*h).^2+80;
n=(b-a)/h;

%{midpoint rule with step h:}%
Rc=0;
for k=1:n
    Rc=Rc+yc(k);
end
Rc=Rc*h

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_Rc=100*(Rc-Q)/Q
Simpson’s rule:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
y=0.5*x.^3-5*x.^2+80;
n=(b-a)/h;

%{Simpson's rule with step h:}%
S=0;
for k=1:(n/2)
    S=S+y(2*k-1)+4*y(2*k)+y(2*k+1);
end
S=S*(h/3)

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_S=100*(S-Q)/Q
Trapezoidal rule:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
y=0.5*x.^3-5*x.^2+80;
n=(b-a)/h;

%{trapezoidal rule with step h:}%
T=0;
for k=2:n
    T=T+y(k);
end
T=h*((y(1)+y(n+1))/2+T)

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_T=100*(T-Q)/Q
Richardson extrapolation:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-4;
b=12;
h=2;
x=a:h:b;
y=0.5*x.^3-5*x.^2+80;
n=(b-a)/h;

%{trapezoidal rule with step h:}%
T=0;
for k=2:n
    T=T+y(k);
end
T=h*((y(1)+y(n+1))/2+T);

%{trapezoidal rule with step 2h:}%
T2=0;
m=0;
count=1;
for m=1:(n/2)
    T2=T2+y(2*m-count)+y(2*m+count);
end
T2=h*T2;

%{Richardson extrapolation:}%
E=T+(1/3)*(T-T2)

%{calculating the error:}%
f=@(x) 0.5*x.^3-5*x.^2+80;
Q=integral(f,a,b);
d_E=100*(E-Q)/Q
Exercise:

clear all;
clc;
close all;

%{defining the range, step, and function:}%
a=-pi/4;
b=pi/2;
n=6;
h=(b-a)/n;
x=a:h:b;
y=sin(cos(x));

%{left rectangle rule with step h:}%
Rl=0;
for k=1:n
    Rl=Rl+y(k);
end
Rl=Rl*h

%{right rectangle rule with step h:}%
Rp=0;
for k=2:(n+1)
    Rp=Rp+y(k);
end
Rp=Rp*h

%{midpoint rule with step h:}%
yc=sin(cos((x+0.5*h)));
Rc=0;
for k=1:n
    Rc=Rc+yc(k);
end
Rc=Rc*h

%{Simpson's rule with step h:}%
S=0;
for k=1:(n/2)
    S=S+y(2*k-1)+4*y(2*k)+y(2*k+1);
end
S=S*(h/3)

%{trapezoidal rule with step h:}%
T=0;
for k=2:n
    T=T+y(k);
end
T=h*((y(1)+y(n+1))/2+T)

%{trapezoidal rule with step 2h:}%
T2=0;
m=0;
count=1;
for m=1:(n/2)
    T2=T2+y(2*m-count)+y(2*m+count);
end
T2=h*T2;

%{Richardson extrapolation:}%
E=T+(1/3)*(T-T2)

%{calculating the error:}%
f=@(x) sin(cos(x));
Q=integral(f,a,b);
d_Rl=100*(Rl-Q)/Q
d_Rp=100*(Rp-Q)/Q
d_Rc=100*(Rc-Q)/Q
d_S=100*(S-Q)/Q
d_T=100*(T-Q)/Q
d_E=100*(E-Q)/Q
First-degree differential equations:

Work in progress.

Differential equations of higher degrees:

Work in progress.

Graph View

Backlinks