Tehnici avansate

MATLAB furnizează un limbaj de nivel înalt care ne permite să gândim abstract despre matrice, vectori şi funcţii. Uneori, totuşi, trebuie să acordăm atenţie unor aspect concrete, cum ar fi viteza de execuţie sau complexitatea programelor. În acest context câteva trucuri pot fi utile.
Mai devreme sau mai tarziu, vom observa că MATLAB este uneori lent. El tinde să acorde prioritate timpului de dezvoltare a programelor în dauna timpului de execuţie. Dacă codul este prea lent ăşi nu vă satisfice, puteţi folosi profiler-ul pentru a detecta eventualele gâtuiri. Dacă se programează cu atenţie şi cu focalizare pe timpul de execuţie se pot obţine câştiguri importante în eficienţă. La limită, se pot scrie subrutinele consumatoare de timp în C sau Fortran şi să se lege codul compilat cu MATLAB. A se vedea în help la “external interfaces.” Acest capitol sugerează câteva moduri de a obţine îmbunătăţiri substanţiale.
La fel, desi vectorii şi matricele sunt o paradigm bună în multe situaţii, exista probleme pentru care ele sunt consumatoare de resurse. Vom vedea şi alte metode de gestionare a datelor.
Cuprins

Prealocarea memoriei

MATLAB ascunde procesul complicat şi costisitor de alocare a memoriei pentru variabile. Această generozitate convenabilă poate cauza irosirea a mult timp de execuţie. Considerăm codul următor, care implementează metoda lui Euler pentru ecuaţia diferenţială vectorială şi memorează valorile la fiecare moment de timp:
tic
A = rand(200);
y = ones(200,1);
dt = 0.001;
for n = 1:(1/dt)
y(:,n+1) = y(:,n) + dt*A*y(:,n);
end
toc
Elapsed time is 0.090000 seconds.
Aceasta necesită aproximativ 0.109679 secde pe un laptop fabricat în 2018. Aproape tot timpul se consumă pe o sarcină necomputaţională.
Când MATLAB întâlneşte instrucţiunea y = ones(200,1), reprezentând condiţia iniţială, el cere sistemului de operare să-I aloce un bloc de memorie pentru 200 de numere în virgulă flotantă. La prima execuţie a ciclului, devine clar că avem nevoie de spaţiu pentru 400 de numere, deci va fi nevoie de un nou bloc de aceasta dimensiune. La iteraţia următoare, dimensiunea nu mai corespunde şi este nevoie să fie alocată mai multă memorie. Acest program mic necesită 1001 de alocări de memorie individuale pentru a creşte dimensiunea!
Schimbând a doua linie în y = ones(200,1001); nu se schimbă nimic din matematica rezolvării problemei, ci se alocă toata memoria necesară o singură dată. Acest pas se numeşte prealocare. Cu prealocare programul necesită în jur de 0.092783 secunde. Pe acelaşi calculator ca mai sus viteza creşte de aproape 3 ori.

Vectorizarea

Vectorizarea se referă la evitarea ciclurilor for şi while. De exemplu, presupunem că x este un vector coloană şi că dorim să calculăm matricea D astfel încât . Implementarea standard necesită două cicluri imbricate:
x=(1:5)';
n = length(x);
D = zeros(n); % preallocation
for j = 1:n
for i = 1:n
D(i,j) = x(i) - x(j);
end
end
Ciclurile se pot scrie aici în orice ordine. Ciclul interior poate fi înlocuit cu o operaţie vectorială:
n = length(x);
D = zeros(n); % preallocation
for j = 1:n
D(:,j) = x - x(j);
end
Putem evita ciclul, trecând de la vectori la tablouri bidimensionale. Implementarea este un pic mai subtilă:
n = length(x);
X = x(:,ones(n,1)); % copy columns to make n by n
D = X - X.';
Cea de-a doua linie este un truc introdus în secţiunea referitoare la matrice (numit în cercurile MATLAB “Tony’s trick”). De notat şi utilizarea lui .’ în ultima linie pentru compatibilitate cu elemente complexe.
De ce este necesară vectorizarea? Sunt două motive: viteză şi stil. Niciunul nu este simplu. Până în jurul lui 2002, vectorizarea atentă ducea întotdeauna la îmbunătăţiri notabile de viteză. Dar, situaţia s-a schimbat într-o oarecare măsură datorată unei tehnici numite accelerare JIT. Accelerarea, care se aplică automat, poate înlătura penalizarea de viteză care apare traditional în ciclurile MATLAB. Nu orice ciclu poate fi optimizat, iar vectorizarea codului nu este critică pentru viteză în orice situaţie. De exemplu, pentru timpul în milisecunde pentru fiecare din cele trei metode de mai sus a fost 0.025, 0.010 şi respectiv 0.025. În acest caz, un nivel mediu de vectorizare s-a dovedit cel mai rapid.
Stilul este ceva subiectiv. Programatorilor vechi, dar care sunt noi în MATLAB, multlipicarea ciclurilor imbricate peste operaţii scalare li se pare naturală şi inevitabilă. Totuşi, MATLAB face uşor accesul la vectori şi matrice şi nu la elementele lor individuale, private ca obiecte atomice. Superioritatea lui A*B faţă de un ciclu imbricat triplu pentru a calcula un produs de matrice este evidentă; multe probleme par să fie situate între aceste două extreme. În codurile de mai sus, de exemplu, fiecare versiune este mai scurtă de cât cea de dinaintea ei. În timp ce a doua nu este mai puţin clară decât prima, versiunea a treia ne cere să ne gândim un pic cand o întâlnim prima dată. Spre deosebire de primele două versiuni, a treia nu este uşor ajustată pentru a ţine cont de antisimetria evidentă din rezultat.
O altă ilustrare clasică a competiţiei implicată de vectorizare provine din eliminarea gaussiană. Iată implementarea clasică, fără nici o vectorizare:
n = length(A);
for k = 1:n-1
for i = k+1:n
s = A(i,k)/A(k,k);
for j = k:n
A(i,j) = A(i,j) - s*A(k,j);
end
end
end
Din nou, putem începe vectorizarea cu ciclul cel mai interior, în j. Fiecare iteraţie a acestui ciclu este independentă de toate celelalte. Acest paralelism este o indicaţie serioasă că putem utiliza în loc o operaţie vectorială:
n = length(A);
for k = 1:n-1
for i = k+1:n
s = A(i,k)/A(k,k);
cols = k:n;
A(i,cols) = A(i,cols) - s*A(k,cols);
end
end
Noua versiune face ideea algoritmică de operaţie pe linii mai clară şi mai de preferat. Totuşi, partea cea mai interioară a ciclului rămas este de asemenea vectorizabilă:
n = length(A);
for k = 1:n-1
rows = k+1:n;
cols = k:n;
s = A(rows,k)/A(k,k);
A(rows,cols) = A(rows,cols) - s*A(k,cols);
end
Trebuie să vă încordaţi un pic muşchii de algebra liniară pentru a vedea că produsul exterior vectorial din următoarea până la ultima linie este adecvat. Această observaţie este clarificatoare şi lămuritoare, dar nu conduce la o îmbunătăţire incontestabilă a stilului.
Pentru a compara viteza acestor trei versiuni, fiecare a fost rulată de 20 de ori pentru diferite valori ale dimensiunii n a matricei pe un laptop în MATLAB 2019a. Rezultatele în milisecunde pe factorizare se dau în tabela 6.1. Rezultatele acestui experiment argumentează impotriva unei recomandări uniforme a vectorizării bazate pe viteza de execuţie! Experimentul indică neajunsurile modului classic de gândire despre performanţa algoritmilor doar în termeni de numar de flops: nici una din liniile de mai sus nu este reprezentativă pentru complexitatea de .
O aplicaţie necontroversată a vectorizării este utilizarea operaţiilor pe tablouri şi a funcţiilor utilitare la evaluarea unei expresii matematice. În comparaţiile următoare, este dificil să susţinem că versiunile cu ciluri sunt superioare corespondentelor vectorizate, odată ce înţelegeţi funcţiile de utilizat:
y = t.*sin(t.^2); y = zeros(size(t));
for i = 1:length(t)
y(i) = t(i)*sin(t(i)ˆ2);
end
%------------------------------------------------------------
dz = diff( z([1:end 1]) ); dz = zeros(size(z));
for j = 1:length(z)-1
dz(j) = z(j+1)-z(j);
end
dz(end) = z(end)-z(1);
%------------------------------------------------------------
e = 1+sum(1./cumprod(1:20)); e = 1; p = 1;
for j = 1:20,
p = p*j;
e = e + 1/p;
end
Concluzia este că ciclurile nu trebuie programate neglijent. Scrierea de cod care lucrează la nivel vectorial este uşoară şi natural în cele mai multe cazuri. Dar dacă profiler-ul indică că se consumă mult timp într-un loc unde nivelul de vectorizare este selectabil, experimentarea este singurul mod de a vedea care nivel este cel mai bun. Pentru introducere în tipuri mai sofisticate de vectorizare, a se vedea help-urile online pentru accumarray, arrayfun şi bsxfun.

Mascare (Masking)

Un tip special de vectorizare este mascarea. Să considerăm un vector x de valori în care vrem să evaluăm o functie definită pe porţiuni
Metoda standard cu cicluri este:
x=-1:0.1:1;
f = zeros(size(x));
for j = 1:length(x)
if abs(x(j)) <= 0.5
f(j) = 1 + cos(2*pi*x(j));
end
end
Modul mai scurt utilizează o mască:
f = zeros(size(x));
mask = (abs(x) < 0.5);
f(mask) = 1 + cos(2*pi*x(mask));
Masca este un index logic al lui x (vezi pagina ??). Vă puteţi referi, la nevoie, la punctele nemascate utilizând ˜mask .
Considerăm o nouă versiune a sumei de numere prime de la pagina??32. Iată cum putem număra numerele prime mai mici ca 100 şi cum le putem însuma:
isprm = isprime(1:100); % which are prime?
sum( isprm ) % how many primes
ans =
25
sum( find(isprm) ) % sum the primes
ans =
1060
Aici find converteşte un index logic într-unul absolut, i.e., furnizează un vector de numere prime.

Excepţii de domeniu

Este important să cunoștem regulile de domeniu și de vizibilitate când scriem cod format din mai multe funcții.
Cea mai puțin utilă și potențial cea mai vătămătoare violare a regulilor de domeniu provine de la variabilele globale. Orice variabilă poate fi declarată globală înainte de a i se atribui pentru prima dată o valoare în domeniul ei curent. După aceasta, orice alt spațiu de lucru poate declara variabila ca fiind globală și să-i acceseze sau să-i schimbe valoarea. Necazul cu variabilele globale este că nu se potrivesc cu proiectele mari și nici cu cele de mărime moderată; conflicte de nume și comportări neașteptate tind să apară rapid. În anumite momente valorile globale au fost mai mult sau mai puțin necesare și pe mai pot găsi încă exemple de coduri care le utilizează. La ora actuală sunt mecanisme mai dezirabile și mai stabile. În particular, variabilele globale nu se vor utiliza pentru a transmite parametrii între funcții.
O excepție interesantă la regulile de domeniu o constituie variabila persistentă. Când o funcție declară o variabilă persistentă, valoarea variabilei este păstrată între apelurile acelei funcții particulare. În timp ce acest mecanism necesită puțină atenție pentru a asigura o utilizare corectă, este mult mai puțin genral decât o variabilă globală, deoarece variabila persistentă este vizibilă doar în spațiul de lucru al funcției care o declară.
Variabilele persistente se pot utiliza pentru a înregistra informații despre starea internă a unei funcții, sau pentru a păstra rezultate preliminare costisitoare, care ar putea fi reutilizate ulterior. Considerăm exemplul care calculează numerele lui Fibonacci:
function y = fib(n)
persistent f
if length(f) < 2, f = [1 1]; end
for k = length(f)+1:n
f(k) = f(k-2) + f(k-1);
end
y = f(1:n);
La primul apel al funcției, f va fi vid. Spre deosebire de alte variabile, variabilelor persistente li se atribuie o valoare inițială—concret, matricea nulă. Astfel, lui f i se vor atribui primele două valori din șir. După aceasta, f se extinde după cum este nevoie pentru a calcula primele n valori, iar șirul este returnat apelantului. La viitoarele apeluri ale lui fib, orice valoare calculată anterior este mai degrabă accesată decât calculată. Dacă funcția ar fi fost apelată cu secvența de lungime , efortul de calcul este proporțional cu nu cu . Același efect se poate realiza punând f ca argument de intrare, dar care face apelul de funcție reponsabil cu tratarea de date relevante doar în interiorul lui fib. Abordarea de aici este autoconținută.

Siruri de caractere

Un șir de caractere MATLAB este delimitat între apostrofuri. Dacă apostroful este conținut în șir el trebuie dublat, ca în 'It''s Cleve''s fault'. Un șir de caractere este de fapt un vector linie de coduri de caractere și astfel șirurile de caractere pot fi concatenate utilizând sintaxa matricială:
str = 'Hello world';
str(1:5)
ans = 'Hello'
double(str)
ans = 1×11
72 101 108 108 111 32 119 111 114 108 100
char(ans)
ans = 'Hello world'
['Hello',' ','world']
ans = 'Hello world'
['Hello';'world']
ans = 2×5 char array
'Hello' 'world'
Concatenarea verticală nu este atât de directă, deoarece toate liniile unui tablou trebuie să aibă aceeași lungime. Se poate utiliza char în locul completării șirurilor cu spații înainte de concatenare. Dacă, totuși, doriți să creați colecți de șiruri, atunci un mecanism mai bun este tabloul de celule. Sunt disponibile multe funcții pentru tratarea șirurilor de caractere; a se vedea help-ul pe strfun. Iată câteva exemple:
upper(str)
ans = 'HELLO WORLD'
strcmp(str,'Hello world')
ans = logical
1
strfind('world',str)
ans = []
De un interes particular este conversia între numere și reprezentarea lor prin caractere. Se poate converti un șir, cum ar fi ’3.14’ în numărul corespunzător cu str2num sau str2double . De notat, diferența importantă între double, care convertește caracterele în codurile lor și str2double , care interpretează șirul ca un număr (îl convertește într-un număr)!
Pentru operația inversă de conversie a unui număr într-un șir de caractere, avem funcțiile num2str și sprintf. Ambele acceptă șiruri de conversie în stilul din C cu specificatori de conversie (%d, %f, etc.), dar din nefericire fără estensii pentru numere complexe. În num2str, se aplică un singur specificator de conversie tuturor elementelor unui tablou, în timp ce sprintf este mai flexibil, putând aplica specificatori de conversie multiplii, în mod ciclic, asupra datelor, procedând în ordinea uzuală pe linii la prelucrarea tablourilor. Pentru ieșiri pe ecran sau într-un fișier se utilizează fprintf .
O utilizare comună a funcției sprintf este de a crea șiruri, cum ar fi nume de fișiere, care conțin în ele întregi consecutivi, ca în exemplul:
for n = 1:12
fn = sprintf('mydata_%.2i',n);
load(fn), data(:,:,n) = A;
end
Acesta va încărca fișierele de date cu numele mydata_01.mat , …, mydata_12.mat, memorând o matrice numită A ca nivele ale unui tablou tridimensional. De notat că
load mydata_01
load('mydata_01')
|
sunt funcțional identice, dar al doilea mod este mai potrivit când numele de fișier este memorat într-o variabilă. Această ,,dualitate comandă/funcție” are loc pentru orice comandă MATLAB care admite argumente suplimentare separate prin spații. O altă utilizare comună este la formatarea unor tabele la ieșire. Script-ul următor afișează aproximațiile Taylor succesive ale lui :
clear
T=zeros(7,1);
x=0.25; n=1:6; c=1./cumprod([1 n]);
for n=1:7
T(n)=polyval(c(n:-1:1),x);
end
fprintf('\n T_n(x) |T_n(x)-exp(x)|\n');
T_n(x) |T_n(x)-exp(x)|
fprintf('----------------------------------\n');
----------------------------------
fprintf('%15.12f %8.3e\n', [T;abs(T-exp(x))] )
1.000000000000 2.840e-01 1.250000000000 3.403e-02 1.281250000000 2.775e-03 1.283854166667 1.713e-04 1.284016927083 8.490e-06 1.284025065104 3.516e-07 1.284025404188 1.250e-08

Tablouri de celule

Gruparea obiectelor de tipuri și dimensiuni diferite este o necesitate frecventă în programare. De exemplu, să presupunem că dorim să tabelăm polinoamele Cebîșev , , , , ș.a.m.d. În MATLAB, un polinom se poate reprezenta ca un vector de coeficienți, ordonați descrescător după puterile variabilei, deci lungimea depinde de gradul polinomului. Pentru o colecție de polinoame Cebîșev consecutive, începând cu , am putea utiliza un tablou triunghiular, dar acest lucru nu este nici convenabil, nici general.
Tablourile de celule sunt un mecanism de grupare a obiectelor eterogene într-o singură variabilă. Ele sunt indexate ca tablourile numerice obișnuite, dar elementele lor pot fi orice, inclusiv alte tablouri de celule. Un tablou de celule se creaza sau se indexează prin acolade {} în loc de paranteze. Tablourile de celule pot avea orice mărime și dimensiune, iar elementele lor nu trebuie să fie de aceeași mărime sau tip. Datorită generalității lor, tablourile de celule sunt doar niște containere; ele nu suportă nici un fel de operații aritmetice.
str = { 'Goodbye', 'cruel', 'world' }
str = 1×3 cell
'Goodbye' 'cruel' 'world'
str{2}
ans = 'cruel'
T = cell(1,9);
T(1:2) = { [1], [1, 0] };
for n=2:8
T{n+1}=[2*T{n} 0] - [0 0 T{n-1}];
end
T
T = 1×9 cell
 123456789
11[1,0][2,0,-1][4,0,-3,0][8,0,-8,0,1][16,0,-20,0,5,0][32,0,-48,0,18,0,-1][64,0,-112,0,56,0,-7,0][128,0,-256,0,160,0,-32,0,1]
T{4}
ans = 1×4
4 0 -3 0
Al doilea exemplu evidențiază o subtilitate la extragerea elementelor dintr-un tablou de celule. Referirea într-o atribuire la T(1:2) se referă la un subtablu extras din tabloul de celule T; astfel membrul drept al atribuirii este un tablou de celule compatibil. Construcția T{n+1} se referă la un singur element al tabloului, în acest caz un vector. În rezumat,
T{4}
ans = 1×4
4 0 -3 0
T(4)
ans = 1×1 cell array
{1×4 double}
Astfel, membrul drept al atribuirii din ciclu nu trebuie inclus între acolade. Daca facem aceasta nu se obține o eroare la primul pas al ciclului — elementul T{3} va fi el însuși un tablou de celule cu un element — dar va cauza o eroare la al doilea pas, când tabloul de celule este concatenat ilegal cu un vector. Cum se va comporta referința la T{4:5}, care se va referi la doi vectori? Răspunsul este că MATLAB acționează ca și cum cei doi vectori se introduc in linia de comanda ts as though those two vectors had just been entered in the command line at that spot:
T{4:5}
ans = 1×4
4 0 -3 0
ans = 1×5
8 0 -8 0 1
Această sintaxă se dovedește surprinzător de utilă, în particular idiomul T{:} , care este interpretat ca o listă de elemente ale lui T, separate prin virgulă. De exemplu, utilizând tabloul de celule str crea anterior, obținem
char(str{:}) % = char(’Goodbye’,’cruel’,’world’)
ans = 3×7 char array
'Goodbye' 'cruel ' 'world '
Combinată cu cat (comanda pentru concatenare de tablouri), această sintaxă se poate utiliza pentru a converti celule ce conțin date de dimensiune compatibilă în tablouri numerice:
c = { [3 4], [5 6] };
cat(1, [1 2], c{:} )
ans = 3×2
1 2 3 4 5 6
cat(2, [1 2], c{:} )
ans = 1×6
1 2 3 4 5 6
e = {}; cat(2, [1 2], e{:} )
ans = 1×2
1 2
Se observă că pentru tabloul de celule vid e, expresia e{:} nu este numai vidă, dar înghite virgula din fața ei, care altfel ar cauza o eroare de sintaxă. Operația inversă de conversie a unui tablou în unul de celule se poate realiza prin num2cell :
anum2cell(1:6)
ans = 1×6 cell
 123456
1123456
num2cell(1:6,2)
ans = 1×1 cell array
{1×6 double}
Al doilea argument al lui num2cell specifică care dimensiune să se ,,împacheteze" într-o celulă; valoarea implicită este 1.
Tabloul de celule special varargin se utilizează pentru a transmite argumentele opționale funcțiilor. În acest mod, se pot scrie funcții care au o comportare care depinde de numărul de argumente transmise și care poate fi ignorată. De exemplu, presupunem că dorim să scriem o funcție care returneză specificația de culoare pentru albastru, atât în modelul de culoare RGB (implicit) sau în modelul HSV:
function b = blue(varargin)
if nargin < 1
varargin = {'rgb'};
end
switch(varargin{1})
case 'rgb'
b = [0 0 1];
case 'hsv'
b = [2/3 1 1];
otherwise
error('Unrecognized color model.')
end
Atât blue cât și blue('rgb') returnează prima opțiune, în timp ce blue('hsv') returnează cea de-a doua. Un al doilea mod de a utiliza varargin este când se poate da un număr variabil de parametrii de intrare. De exemplu, considerăm:
function s = add(s,varargin)
for n = 1:nargin-1
s = s + varargin{n};
end
Aici, primul argument de intrare se atribuie numelui s, iar toate celelate tabloului de celule varargin . Cuvăntul cheie nargin returnează numărul total de intrări, atât cele cu nume cât și cele opționale. Pentru a asigura aceeai funcționalitate și la ieșire, se poate utiliza varargout și nargout .