====== Macierze c.d. ======
* tworzenie macierzy i indeksowanie elementów macierzy
* operacje macierzowe i funkcje działające na macierzach
* macierze logiczne i operacje na macierzach logicznych
* optymalizacja kodu
* wektoryzacja - zamiana pętli po elementach na mechanizmy indeksowania
* prealokacja - dynamiczne modyfikowanie rozmiarów macierzy jest kosztowne
===== Definiowanie macierzy =====
* macierze można tworzyć za pomocą operatora indeksowania ''%%:%%''
* uwaga: ilość elementów w wierszach musi być jednakowa
A = [1 2 3; 4 5 6]
B = [5:8 ; 10:13]
A =
1 2 3
4 5 6
B =
5 6 7 8
10 11 12 13
===== Funkcje tworzące macierze =====
* ''%%rand(n,m)%%'' losowe wartość z przedziału $[0, 1]$
* ''%%zeros(n,m)%%'' macierz zerowa
* ''%%ones(n,m)%%'' macierz zawierająca jedynki
* ''%%eye(n,m)%%'' macierz jednostkowa (jedynki na diagonali)
a = rand(3)
b = rand(3, 2)
a =
0.306403 0.834324 0.986851
0.031466 0.872381 0.695912
0.942068 0.691957 0.619879
b =
0.244435 0.975030
0.421339 0.069198
0.647612 0.217188
c = zeros(5)
d = ones(2,5)
c =
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
d =
1 1 1 1 1
1 1 1 1 1
e = eye(4, 6)
e =
Diagonal Matrix
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
0 0 0 1 0 0
**Ćwiczenie** utwórz macierz o wymiarach 100x100, która posiada wartości 1 w pierwszej i czwartej ćwiartce oraz wartości 0 w pozostałych miejscach, zgodnie z poniższym przykładem macierzy 4x4\\
1 1 0 0
1 1 0 0
0 0 1 1
0 0 1 1
===== Indeksowanie macierzy =====
* indeksowanie 2 wymiarów macierzy ''%%(N,M)%%'' od ''%%1%%'' do ''%%end%%''
* indeksowanie liniowe ''%%(N)%%'' od ''%%1%%'' do ''%%end%%'' po elementach kolejnych kolumn
X = [1, 2, 3; 4, 5, 6]
X =
1 2 3
4 5 6
% indeksowanie 2-wymiarowe
a = X(1,2) % 1-szy wiersz 2-ga kolumna
b = X(1:2, 1)
c = X(1:2, 1:2)
a = 2
b =
1
4
c =
1 2
4 5
e = X(:, 1) % cała kolumna
f = X(end, :) % cały wiersz
e =
1
4
f =
4 5 6
% indeksowanie liniowe
a = X(1)
b = X(end)
c = X(1:end) # wszystkie elementy w wektorze wierszowym
d = X(:) # wszystkie elementy w wektorze kolumnowym
a = 1
b = 6
c =
1 4 2 5 3 6
d =
1
4
2
5
3
6
===== Rozmiary wektorów i macierzy =====
* ''%%length(w)%%'' ilość elementów wektora lub najdłuższy wymiar macierzy
* ''%%size(m)%%'' zwraca ilość wierszy i kolumn
* ''%%numel(m)%%'' ilość elementów
w = 1:100;
a = length(w)
b = size(w)
c = numel(w)
a = 100
b =
1 100
c = 100
X = rand(13, 100);
a = length(X)
[w k] = size(X)
c = numel(X)
a = 100
w = 13
k = 100
c = 1300
===== Iteracja po elementach =====
* W wielu przypadkach iteracja po elementach macierzy za pomocą pętli nie jest konieczna
* Indeksowanie za pomocą ''%%:%%'', operatory macierzowe (''%%*%%'', ''%%+%%'', itd.) oraz funkcje dedykowane do obliczeń na macierzach (np. ''%%sum(X)%%'') będą z reguły wydajniejsze od pętli
% Iterowanie po elementach macierzy
w = rand(1, 10);
for i=1:length(w)
fprintf('Element %d wektora = %f\n', i, w(i))
end
Element 1 wektora = 0.468979
Element 2 wektora = 0.160596
Element 3 wektora = 0.254776
Element 4 wektora = 0.487126
Element 5 wektora = 0.721925
Element 6 wektora = 0.443289
Element 7 wektora = 0.345285
Element 8 wektora = 0.392169
Element 9 wektora = 0.011494
Element 10 wektora = 0.963213
% Iterowanie po elementach macierzy
m = rand(2 ,3);
[w k] = size(m);
for i=1:w
for j=1:k
fprintf('Element (%d,%d) macierzy = %f\n', i, j, m(i, j))
end
end
Element (1,1) macierzy = 0.192715
Element (1,2) macierzy = 0.097834
Element (1,3) macierzy = 0.704280
Element (2,1) macierzy = 0.938989
Element (2,2) macierzy = 0.119784
Element (2,3) macierzy = 0.608151
**Ćwiczene** napisz funkcję ''%%suma%%'', która wyznaczy sumę wszystkich elementów macierzy podanej w argumencie. Ćwiczenie wykonaj z użyciem instrukcji pętli, iterując po wszystkich elementach macierzy.
W jaki sposób osiągnąć to samo za pomocą funkcji matlabowej ''%%sum%%''?
**Ćwiczenie** napisz skrypt/funkcję o nazwie ''%%suma2.m%%'', która uruchamia operację sumowania z poprzedniego ćwiczenia wielokrotnie (np. n=10000) dla macierzy 100x100 o losowych wartościach i mierzy czas wykonania tej operacji. Tą samą operację powtórz za pomocą funkcji ''%%sum%%'' i porównaj czasy obliczeń z obu podejść.
Czas wykonywania operacji można wyznaczyć za pomocą poleceń ''%%tic%%'' i ''%%toc%%'', przykład:
tic
sum(x)
toc
Elementy macierzy można również zsumować za pomocą iloczynu z macierzami jedynkowymi:
$$
\left[\begin{array}{ccc}
1 & 1 & 1 \\
\end{array}\right] \cdot\left[\begin{array}{ll}
3 & 1 \\
2 & 1 \\
1 & 0
\end{array}\right]
\cdot\left[\begin{array}{ll}
1 \\
1 \\
\end{array}\right]
$$
Sprawdź również szybkość obliczeń sumy tą metodą.
===== Zmiana rozmiaru macierzy =====
* ''%%reshape(A, M, N)%%'' lub ''%%reshape(A, [M N])%%'' zmienia kształt macierzy na ''%%MxN%%'' (elementy układane kolumnami)
* ''%%fliplr(m)%%'' , ''%%flipup(m)%%'', ''%%rot90()%%'' odbicie lewo-prawo, góra-dół i obrót 90 stopni
w = 1:12;
a = reshape(w, 3, 4)
a1 = reshape(a, [2 6])
a =
1 4 7 10
2 5 8 11
3 6 9 12
a1 =
1 3 5 7 9 11
2 4 6 8 10 12
b = fliplr(a)
c = flipud(a)
d = rot90(a)
b =
10 7 4 1
11 8 5 2
12 9 6 3
c =
3 6 9 12
2 5 8 11
1 4 7 10
d =
10 11 12
7 8 9
4 5 6
1 2 3
**Ćwiczenie** korzystając z funkcji ''%%randperm%%'' wygeneruj macierz 10x10 zawijającą liczby od 1 do 100 umieszczone w losowej kolejności.
===== Puste macierze =====
m = []; # pusta macierz (wektor)
a = size(m)
b = length(m)
c = numel(m)
m = [m 1]
size(m)
a =
0 0
b = 0
c = 0
m = 1
ans =
1 1
===== Usuwanie elementów macierzy =====
* dla wektorów można usuwać dowolne elementy
* dla macierzy można usuwać tylko całe wiersze lub całe kolumny
w = 1:12
m = reshape(w, 3, 4)
w(1) = [] % usuwanie elementu
w([3 6]) = []
w(6:end) = []
% m(1,1) = []
m(1, :) = []
m(:, 2) = []
w =
1 2 3 4 5 6 7 8 9 10 11 12
m =
1 4 7 10
2 5 8 11
3 6 9 12
w =
2 3 4 5 6 7 8 9 10 11 12
w =
2 3 5 6 8 9 10 11 12
w =
2 3 5 6 8
m =
2 5 8 11
3 6 9 12
m =
2 8 11
3 9 12
===== Operacje na macierzach =====
* operacje skalarne na elementach macierzy: ''%%+%%'', ''%%-%%'', ''%%*%%'', ''%%/%%'', ''%%.^%%'', ''%%.\%%''
A = [1 2 3; 4 5 6]
A =
1 2 3
4 5 6
a = A + 3 # to samo co A .+ 3
b = A * 5 # to samo co A .* 5
c = A / 5 # to samo co A ./ 5
a =
4 5 6
7 8 9
b =
5 10 15
20 25 30
c =
0.20000 0.40000 0.60000
0.80000 1.00000 1.20000
d = A .^ 2
e = A .\ 5 # to samo co 5 ./ A
d =
1 4 9
16 25 36
e =
5.00000 2.50000 1.66667
1.25000 1.00000 0.83333
* operacje macierzowe: ''%%*%%'', ''%%/%%'', ''%%\%%''
* ''%%A * B%%'' iloczyn macierzy (iloczyn skalarny dla wektorów)
* ''%%A / B%%'' dzielenie macierzy, wynikiem jest ''%%X%%'' takie, że ''%%X * B = A%%''
* ''%%A \ B%%'' dzielenie macierzy, wynikiem jest ''%%X%%'' takie, że ''%%A * X = B%%''
* potęgowanie ''%%A^n%%'' tylko dla macierzy ''%%A%%'' kwadratowych
A = [1 2 3; 4 5 6];
B = [1 2 3; 4 5 6];
a = A * B' # ilość wierszy A == ilość kolumn B
b = A / B # wymiary A i B takie same
c = A \ B # wymiary A i B takie same
a =
14 32
32 77
b =
1.0000e+00 4.7673e-16
2.2979e-16 1.0000e+00
c =
0.83333 0.33333 -0.16667
0.33333 0.33333 0.33333
-0.16667 0.33333 0.83333
d = A(1:2,1:2)^2
d =
9 12
24 33
* opearcje macierzowe element po elemencie: ''%%+%%'', ''%%-%%'', ''%%.*%%'', ''%%./%%'', ''%%.^%%'', ''%%.\%%''
* wymiary macierzy muszą się zgadzać
A = [ 1 2 3; 4 5 6 ];
B = [ 3 4 6; 3 4 3 ];
a = A + B # to samo co A .+ B
b = A .* B
c = A ./ B
a =
4 6 9
7 9 9
b =
3 8 18
12 20 18
c =
0.33333 0.50000 0.50000
1.33333 1.25000 2.00000
d = A .^ B
e = A .\ B # to samo co B ./ A
d =
1 16 729
64 625 216
e =
3.00000 2.00000 2.00000
0.75000 0.80000 0.50000
* operacje macierzy i wektora element po elemencie powielane wzdłuż pasującego wymiaru
A = ones(2, 3) + 1
w1 = [1 2 3]
w2 = [1; 2]
A =
2 2 2
2 2 2
w1 =
1 2 3
w2 =
1
2
a = A + w1 # sumowanie elementów w wierszach macierzy
b = A + w2 # sumowanie elementów w kolumnach macierzy
a =
3 4 5
3 4 5
b =
3 3 3
4 4 4
a = A .* w1
b = A .* w2
a =
2 4 6
2 4 6
b =
2 2 2
4 4 4
===== Macierze i funkcje =====
Wiele funkcji Matlaba jest przystosowana do działania na macierzach o dowolnych kształtach
w = 1:6
X = reshape(w, 2, 3)
a = sin(w)
b = cos(X)
w =
1 2 3 4 5 6
X =
1 3 5
2 4 6
a =
0.84147 0.90930 0.14112 -0.75680 -0.95892 -0.27942
b =
0.54030 -0.98999 0.28366
-0.41615 -0.65364 0.96017
**Ćwiczenie** popraw program ''%%pole_kola.m%%'' tak aby funkcja zwracała poprawny wynik dla argumentu macierzowego. Wynikiem funkcji będzie macierz zawierająca w każdym z elementów $P_{ij}$ wartość pola dla promienia $R_{ij}$. Postaraj się nie korzystać z instrukcji pętli.
**Ćwiczenie** zmierz czas wykonywania poniższych instrukcji.
* porównaj szybkość działania przy zastosowaniu prealokowanego wektora ''%%x%%''
* przyśpiesz działanie stosując wektoryzację (pozbądź się instrukcji pętli)
n=100000
tic
for i=1:n
x(i) = sqrt(i);
end
toc
sum(x)
n = 100000
Elapsed time is 0.43374 seconds.
ans = 21082008.97392
===== Wektory i macierze logiczne =====
* macierze zawierające wartości logiczne (typ ''%%logical%%''), wartość ''%%0%%'' to fałsz, wartość ''%%1%%'' (lub niezerowa) to prawda
* indeksowanie macierzy za pomocą macierzy logicznych
* funkcje ''%%false()%%'' i ''%%true()%%'' pozwalają tworzyć macierze logiczne
* funkcja ''%%logical()%%'' zamienia macierz liczbową na logiczną
X = [1 , -5, 7 ; -1, 42, -12];
ind = X < 0
whos ind
ind =
0 1 0
1 0 1
Variables in the current scope:
Attr Name Size Bytes Class
==== ==== ==== ===== =====
ind 2x3 6 logical
Total is 6 elements using 6 bytes
# indeksowanie macierzą logiczną
a = X(ind)
b = X( X > 0)
X( X > 0) = 111
a =
-1
-5
-12
b =
1
42
7
X =
111 -5 111
-1 111 -12
A = true(2)
B = false(2, 3)
whos A B
A =
1 1
1 1
B =
0 0 0
0 0 0
Variables in the current scope:
Attr Name Size Bytes Class
==== ==== ==== ===== =====
A 2x2 4 logical
B 2x3 6 logical
Total is 10 elements using 10 bytes
# wartości logiczne mogą być traktowane jak liczby całkowite
C = A + 1
whos A C
C =
2 2
2 2
Variables in the current scope:
Attr Name Size Bytes Class
==== ==== ==== ===== =====
A 2x2 4 logical
C 2x2 32 double
Total is 8 elements using 36 bytes
**Ćwiczenie** dlaczego poniższe wyrażenie nie działa? Spróbuj to naprawić.
w = 1:5;
w([1 0 1 0 1]) = -1
error: w(0): subscripts must be either integers 1 to (2^63)-1 or logicals
===== Przydatne funkcje operujące na wartościach logicznych =====
* ''%%any()%%'' zwraca 1 gdy wektor posiada niezerowy element, \\ dla macierzy zwraca wektor z wartością 1 gdy w kolumnie jest niezerowa wartość\\
* ''%%all()%%'' zwraca 1 dla każdej kolumny której wszystkie elementy są niezerowe
* ''%%find()%%'' zwraca indeksy (indeksowanie liniowe) niezerowych elementów
* ''%%isequal()%%'' sprawdza czy macierze mają jednakowe elementy
w1 = [ 0 0 1 2 ];
w2 = zeros(1,10);
a = all(w1)
b = any(w1)
c = all(w2)
d = any(w2)
a = 0
b = 1
c = 0
d = 0
E = [ 1 2 ; 0 0 ; 3 * eye(2)]'
all(E)
any(E)
E =
1 0 3 0
2 0 0 3
ans =
1 0 0 0
ans =
1 0 1 1
# indeksy liniowe
k = find(E)
a = E(k)
k =
1
2
5
8
a =
1
2
3
3
# indeksy wierszy i kolumn
[w k] = find(E)
w =
1
2
1
2
k =
1
1
3
4
# szukanie zerowych elementów
k2 = find(~E)
k2 =
3
4
6
7
# szukanie elementów spełniających warunki
k3 = find(E > 2)
k3 =
5
8
X = E
isequal(E, X)
% to samo za pomocą funkcji all()
all(all(X == E))
X(1,1) = -X(1,1)
isequal(E, X)
all(all(X == E))
X =
1 0 3 0
2 0 0 3
ans = 1
ans = 1
X =
-1 0 3 0
2 0 0 3
ans = 0
ans = 0
===== Operatory logiczne =====
* operatory logiczne działające element po elemencie: ''%%&%%'' and, ''%%|%%'' or, ''%%~%%'' not
X = magic(5)
X(X > 5 & X < 10)
X =
17 24 1 8 15
23 5 7 14 16
4 6 13 20 22
10 12 19 21 3
11 18 25 2 9
ans =
6
7
8
9
**Ćwiczenie** utwórz macierz ''%%X%%'' o wymiarach 10x10 posiadającą losowe wartości z zakresu od -1 do 1. Następnie zastąp wszystkie elementy o wartości ujemnej w macierzy ''%%X%%'' liczbą ''%%0%%''.
===== Zadanie 5 Monte Carlo =====
Napisz funkcję o nazwie ''%%mc_pi%%'', która wyznaczy przybliżenie wartości liczby $\pi$ za pomocą całkowania metodą Monte Carlo. Zadanie sprowadza się do wyznaczenia pola koła o promieniu 1.
$$ P_{\circ} = \pi r^2 = \pi $$
Pole ćwierci koła można wyznaczyć licząc stosunek liczby $K$ losowych punktów wewnątrz ćwierci koła o promieniu 1 do liczby $N$ ilości wszystkich wylosowanych punktów w kwadracie o polu 1.
$$4 \cdot \frac{K}{N} \approx \pi $$
{{zajecia:pp2_2020_2:mc_pi.png| image.png}}
Algorytm:
- Wylosuj $N$ losowych punktów $(x, y)$, gdzie $x, y \in [0, 1]$
- Wyznacz liczbę punktów $K$, która trafiła w tarczę koła, tzn. spełniające warunek \\ $$\sqrt{x^2 + y^2} < 1$$
- Wyznacz przybliżenie liczby $\pi$ równe $4 \cdot \frac{K}{N}$
Spróbuj wykonać to zadanie bez wykorzystywania instrukcji pętli (wektoryzacja kodu) oraz zadbaj o wcześniejsze zaalokowanie wykorzystywanych na potrzeby obliczeń macierzy lub wektorów.
Funkcja ''%%mc_pi()%%'' przyjmuje jeden argument ''%%N%%'', który jest liczbą całkowitą określającą ilość losowań.\\
Gdy argument nie jest podany to używana jest domyślna wartość ''%%N=1000%%''