E.

SCHEIBER

Laborator de ANALIZĂ NUMERICĂ

MATLAB

Braşov
Cuprins

I Matlab 4

1 Elemente de programare în MATLAB 5

1.1 Aplicaţii Matlab . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.2 Obiecte Matlab . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3 Elemente de programare in Matlab . . . . . . . . . . . . . . . . . . . 16
1.4 Funcţii(subprograme) în Matlab . . . . . . . . . . . . . . . . . . . . . 18
1.5 Grafică în Matlab . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
1.5.1 Grafică 2D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
1.5.2 Grafică 3D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
1.5.3 Animaţie în Matlab . . . . . . . . . . . . . . . . . . . . . . . . . 31

2 Algebră liniară numerică 34

2.1 Factorizarea unei matrice . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.2 Rezolvarea sistemelor algebrice de ecuaţii liniare . . . . . . . . . . . 36

3 Rezolvarea sistemelor şi ecuaţiilor algebrice 39

3.1 Rezolvarea sistemelor algebrice de ecuaţii neliniare . . . . . . . . . . 39
3.2 Rezolvarea ecuaţiilor algebrice . . . . . . . . . . . . . . . . . . . . . . 41
3.3 Rezolvarea ecuaţiilor polinomiale . . . . . . . . . . . . . . . . . . . . 42

4 Rezolvarea problemelor de interpolare 44

4.1 Interpolare polinomială . . . . . . . . . . . . . . . . . . . . . . . . . . 44
4.2 Interpolare cu funcţii spline cubice . . . . . . . . . . . . . . . . . . . 47

5 Derivare numerică 49
5.1 Derivarea funcţiilor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

6 Metoda celor mai mici pătrate 54

2
CUPRINS 3

7 Integrare numerică 59
7.1 Integrarea funcţiilor de o variabilă reală . . . . . . . . . . . . . . . . . 59
7.2 Calculul numeric al integralelor duble . . . . . . . . . . . . . . . . . . 62

8 Integrarea numerică a E.D.O. 63

8.1 Rezolvarea problemei cu valori iniţiale . . . . . . . . . . . . . . . . . 63
8.2 Rezolvarea problemelor bilocale . . . . . . . . . . . . . . . . . . . . . 65

9 Tranformata Fourier discretă 69

9.1 Calculul coeficienţilor Fourier . . . . . . . . . . . . . . . . . . . . . . 70

Bibliografie 72
Part I

Matlab

4
Capitolul 1

Elemente de programare în MATLAB

Cap. 2

Prezentăm pe scurt elemente ale limbajului de programare utilizat în pro-

dusul Matlab. Elementele prezentate în secţiunile următoare reprezintă doar
ghid de utilizare rapidă. Se presupune că utilizatorul posedă cunoştinţe de pro-
gramare, având experienţă de lucru în cel puţin un limbaj de programare proce-
durală.
Matlab este un produs comercial de calcul numeric disponibil în mediile Win-
dows, Linux şi Mac OS.
Matlab este disponibil în nor

[Link]

dar este necesară înregisrarea prealabilă.

Ca produs informatic Matlab are trăsături comune cu

• Scilab

• Octave

produse distribuite gratuit.

5
6 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

1.1 Aplicaţii Matlab

O aplicaţie Matlab este alcătuită din una sau mai multe funcţii Matlab. În
vederea executării aplicaţiei este nevoie în prealabil de declararea locaţiei funcţi-
ilor în variabila matlabpath. Această acţiune se face prin comanda Set Path din
meniul Matlab. Apelarea / pornirea aplicaţiei se va face dintr-un script în Matlab.
În varianta din nor, funcţiile client trebuie încărcate prin comanda Upload.

1.2 Obiecte Matlab

• Numere reale Matlab reprezintă în memorie orice număr în virgulă mo-
bilă.

• Constante
Tip Mnemonic Valoare
real pi π
real eps ≈ 2.220410−16
inf ∞
p
complex i ,j i= −1
NaN Not A Number

0/0
ans =

NaN

• Literali

– Identificatori - denumiri date de utilizator diverselor entităţi Matlab

(variabile, funcţii, etc).

* Primul caracter al unui identificator este literă.

* Următoarele caractere pot fi litere, cifre sau caracterul de sublin-
iere, .
* Se face distincţie între literele mari şi mici.
* Numărul maxim de caractere este 19.
1.2. OBIECTE MATLAB 7

– Literali numerici floor(10,4 Exemple:

1. Numărul real 10,23 se poate scrie:
10.23=0.1023e+2=0.1023E+2=1023e-2=1023E-2.
2. Numărul complex 1 + 5i se scrie 1 + 5 ∗ %i .
– Literali nenumerici
Un caracter sau un şir se caractere – string – se defineşte (scrie) cuprins
între apostroafe.
Concatenarea a două stringuri s1, s2 se programează prin [s1s2].

• Funcţii matematice uzuale:

Mnemonic Semnificaţie Mnemonic Semnificaţie
abs(x) |x| exp(x) ex
log(x) ln(x)
sin(x) sin(x) asin(x) arcsin(x)
sinh(x) sh(x) asinh(x) arcsh(x)
cos(x) cos(x) acos(x) arccos(x)
cosh(x) ch(x) acosh(x) arcch(x)
tan(x) tg(x) atan(x) arctg(x)
tanh(x) th(x) atanh(x) arcth(x)
cotg(x) ctg(x) coth(x) cth(x)
p sin x
sqrt(x) x sinc(x) x
floor(x) [x] ceil(x) dxe
R x −t 2 R∞
erf(x) p2 e dt gamma(x) Γ(x) = t x−1 e −t dt
π 0 0
gammaln(x) ln Γ(x)
real(z) ℜ(z) imag(z) ℑ(z)

• Şiruri de numere
Progresia aritmetică introduce definind primul termen (a), raţia (r ) şi o
marginea superioară sau inferioară (M ), după cum raţia este pozitivă sau
negativă, prin sintaxa
a :r :M
Dacă parametrul r lipseşte atunci raţia este 1.

a=0.2:0.3:1
0.2 0.5. 0.8
b=1:-0.3:0
1. 0.7 0.4 0.1
8 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

Funcţia Matlab linspace(a, b, n) defineşte progresia a+i b−a

n−1
, i = 0, 1, . . . , n−
1. Dacă se utilizează formatul linspace(a, b) atunci n = 100.
Un şir de numere (a i )1≤i ≤n se poate defini prin

1. i=1:n
a=formula termenului general, functie de i
2. a = [a 1 , a 2 , . . . , a n ],
unde a 1 , a 2 , . . . , a n sunt termenii şirului.

Indicele primului termen al unui şir este 1.

i=1:4
1 2 3 4
b=^2
1 4 9 16
b=[1,4,9,16]
a==b
1 1 1 1

Funcţii de variabilă vector

Mnemonic Semnificaţia
length(v) Numărul elementelor lui v
max(v) Cel mai mare element al lui v
min(v) Cel mai mic element al lui v
fliplr(v) Deplasarea elementelor spre dreapta cu o poziţie

Exemplul 1.1 Inversarea unui vector.

c=[1,2,3];
for i=1:length(c) c=fliplr(c); end
c

c =

3 2 1

• Matrice
1.2. OBIECTE MATLAB 9

– Definirea unei matrice.

Exemplu: Matricea µ ¶
1 2 3
a=
4 5 6
se defineşte în Matlab prin

a=[1,2,3;4,5,6]
Elementele unei linii se pot separa prin virgulă sau spaţiu. Elementele
unei matrice pot fi: constante/variabile numerice, constante booleene,
şiruri de caractere, polinoame, funcţii raţionale.
O progresie aritmetică pr og = a : r : M este interpretată ca un vector
linie (adică o matrice cu o singură linie).
Exemplu. Progresiile aritmetice a = 0.2, 0.6, 1.0, 1.4, 1.8 şi b = 0, −1, −2, −3, −4
se obţin prin

a=0.2 : 0.4 : 2
a=
! 0.2 0.6 1.0 1.4 1.8 !
b=-1 : -1 : -4
b=
! -1 -2 -3 -4 !

[ ] reprezintă matricea vidă.

Funcţia Matlab size(variabilaMatrice) returnează numărul liniilor şi
numărul coloanelor variabileiMatrice.
– Definirea unei matrice rară se face cu funcţia Matlab
a=sparse(i,j,v)
unde i , j , v sunt vectori de aceaşi dimensiune. Dacă

i = [i 1 , i 2 , . . . , i m ]
j = [ j1, j2, . . . , jm ]
v = [v 1 , v 2 , . . . , v m ]

atunci a i k , j k = v k , ∀k.
Funcţia full(a ) afişează o matrice rară a în formatul obişnuit.
Funcţia nnz(a) are ca valoare numărul elementelor nenule din a.
– Matrice speciale.
10 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

* zeros(m, n) defineşte o matrice nulă cu m linii şi n coloane;

* ones(m, n) defineşte o matrice cu toate elementele egale cu 1
având m linii şi n coloane;
* eye(m,n) defineşte o matrice unitate cu m linii şi n coloane;
* rand(m,n) defineşte o matrice cu m linii şi n coloane având ca
elemente numere (semi)aleatoare cuprinse între 0 şi 1.
– Selectarea unui element se obţine prin sintaxa

variabilaMatrice(numărLinie, numărColoană)

Prima linie şi prima coloană are numărul de ordine 1.

Elementul unei matrice rare se selectează prin intermediul funcţiei
full:
full(variabilaMatrice(numărLinie, numărColoană)).
– Operatori matriceali.
Simbol Semnificaţie
+ adunare
- scădere
∗ înmulţire de matrice
.∗ înmulţire pe componente (Produsul Hadamard)
ˆ ridicare la putere prin produs matriceal
.ˆ ridicarea la putere a componentelor
\ a\b = a −1 · b
bi , j
.\ a.\b = ( a )i , j
i,j
/ b/a = b · a −1
bi , j
./ b/a = ( a )i , j
i,j
’ transpunere / conjugare şi transpunere

– Funcţii matriceale.
Mnemonic Semnificaţia
triu(A) Matricea superior triunghiulară a lui A
tril(A) Matricea inferior triunghiulară a lui A
diag(A) Vectorul cu elementele diagonale ale lui A
diag(v) Matrice pe diagonală cu elementele vectorului v
size(A) Şir cu dimensiunile lui A

– Extinderea unei matrice. Date fiind o matrice a şi un vector linie v

1.2. OBIECTE MATLAB 11

* adăugarea liniei v ca prima linie a matricei a se obţine prin

[v; a]

* adăugarea liniei v ca ultimă linie a matricei a se obţine prin

[a; v]

* adăugarea vectorului v ca prima coloană a matricei a se obţine

prin
[v 0 , a] sau [v; a 0 ]0

* adăugarea vectorului v ca ultimă coloană a matricei a se obţine

prin
[a, v 0 ] sau [a 0 ; v]0

Trebuie observată diferenţa dintre cazul în care extinderea se face prin

linie faţă de cel în care se extinde prin coloană prin utilizarea “;” şi
respectiv “,”.
– Extragerea unei submatrice. O zonă compactă fixată prin liniile l 1 : l 2
şi coloanele c 1 : c 2 , inclusiv, se extrage din matricea a prin

a(l 1 : l 2 , c 1 : c 2 )

Pentru extragerea unei zone necompacte – aflată la intersecţia liniilor

l 1 , l 2 , . . . , l p cu coloanele c 1 , c 2 , . . . , c q – dintr-o matrice a având m linii
şi n coloane se definesc vectori linie
½
t r ue dacă i ∈ {l 1 , . . . , l p }
l i n = (l i n i )1≤i ≤m , l i ni =
f al se în caz contrar

½
t r ue dacă j ∈ {c 1 , . . . , c q }
col = (col j )1≤ j ≤n , col j =
f al se în caz contrar
Extragerea rezultă din
a(l i n, col )

• Elemente de calcul simbolic

– Simboluri pentru numere

12 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

a=sym(1/2)

a =

1/2

b=sym(1/3)

b =

1/3

(a+b)*a

ans =

5/12

– Variabile simbolice

syms a b x
f=a*x+b

f =

b + a*x

subs(f,x,2)

ans =

2*a + b

– Expresii simbolice
1. x=sym(sqrt(1+sym(sqrt(2))))

x =

(2^(1/2) + 1)^(1/2)
1.2. OBIECTE MATLAB 13

>> (x^2-1)^2

ans =

2
2. syms x y
expand((1+x)^3)
ans =

x^3 + 3*x^2 + 3*x + 1

factor((x+y)^7-x^7-y^7)

ans =

[7, x, y, x + y, x^2 + x*y + y^2, x^2 + x*y + y^2]

3. solve(x^3-1==0,x)

ans =
1
- (3^(1/2)*1i)/2 - 1/2
(3^(1/2)*1i)/2 - 1/2
– Calcul diferenţial şi integral:
Derivare simbolică, Calculul primitivei, Integrare simbolică

syms x
f=sin(x);
diff(f,x)

ans =

cos(x)

int(f,x)

ans =

-cos(x)
14 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

int(f,0,pi)

ans =

– Polinoame

c=[1,-3,2];
p=poly2sym(c)

p =

x^2 - 3*x + 2

a=sym(1/2)

a =

1/2

subs(p,a)

ans =

3/4

cp=sym2poly(p)

cp =

1 -3 2

q=poly2sym([1,-1])

q =

x - 1

simplifyFraction(p/q+q)
1.2. OBIECTE MATLAB 15

ans =

2*x - 3

gcd(p,q)

ans =

x - 1

lcm(p,q)

ans =

x^2 - 3*x + 2

[d,u,v]=gcd(p,q)

d =

x - 1

u =

v =

d-u*p-v*q

ans =

0
16 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

• Listă. O listă se defineşte prin

{e 1 , e 2 , . . . , e n }

cu elementele e 1 , e 2 , . . . , e n unde e i sunt obiecte Matlab.

Dacă l este o listă atunci

– l (i )
returnează al i -lea element al listei;

Exemplu.

l={1,true,’abc’,[1,2;3,4]}

l =

[1] [1] ’abc’ [2x2 double]

1.3 Elemente de programare in Matlab

Programarea în Matlab se face în cadrul funcţiilor definite în exterior. Pentru
editarea acestor funcţii Matlab posedă un editor propriu, dar se poate fi utilizat
orice editor de fişiere.
O linie de comentariu este % text - comentariu.
Fiecare instrucţiune introdusă în dreptul promptului Matlab este executată
şi rezultatul este afişat. Pe o linie pot fi puse mai multe instrucţiuni separate prin
, sau ;. În plus, încheierea unei comenzi cu caracterul ; are ca efect inhibarea
afişării rezultatului.
Pentru exprimarea condiţiilor se utilizează

Operatori relaţionali Operatori logici

Simbol Semnificaţie Simbol Semnificaţie
== = & şi
∼= 6 = | sau
<= ≤ ∼ negaţia
< <
>= ≥
> >
1.3. ELEMENTE DE PROGRAMARE IN MATLAB 17

• Instrucţiunea de atribuire
variabila=expresie

• Instrucţiuni condiţionate
1.
if condiţie
instrucţiuni
elseif condiţie
instrucţiuni
else
instrucţiuni
end

Exemplul 1.2 Rezolvarea ecuatiei x 2 + ax + b = 0 în linie de comandă.

a=3;b=2;delta=a^2-4*b;
if delta>=0 x=(-a+sqrt(delta))/2;y=(-a-sqrt(delta))/2;else x=-a/2;y=sqrt(-delta)/2;end

2.
switch expresie
case expr esi e 1
instrucţiuni
case expr esi e 2
instrucţiuni
..
.
otherwise
instrucţiuni
end

• Instrucţiuni de ciclare

1. for n = n1 : pas : n2 instrucţiuni end

Dacă pas = 1 atunci acest parametru este opţional.

Exemplul 1.3 5i =1 i
P

s=0;
for i=1:5 s=s+i; end

Exemplul 1.4 Împărţirea elementelor unei mulţimi în clase de echivalenţă

modulo 3.
18 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

a=[1,2,3,4,5,6,7,8,9];
c0=zeros(3,1);c1=zeros(3,1);c2=zeros(3,1);
i0=0;i1=0;i2=0;
for i=1:length(a)
switch a(i)-floor(a(i)/3)*3
case 0
i0=i0+1;c0(i0)=a(i);
case 1
i1=i1+1;c1(i1)=a(i);
case 2
i2=i2+1;c2(i2)=a(i);
end
end

[c0,c1,c2]

ans =

3 1 2
6 4 5
9 7 8

Exemplul 1.5 Şirul lui Fibonacci generat în linie de comandă.

F=[1,1];
n=15;
for i=3:n F(i)=F(i-1)+F(i-2); end

2. while condiţie instrucţiuni end

P5
Exemplul 1.6 i =1 i
s=0;
i=0;
while i<5 i=i+1;s=s+i; end

Instrucţiunea break realizează un salt necondiţionat la prima instrucţiune

aflată după instrucţiunea de ciclare.
Instrucţiunea continue realizează un salt necondiţionat la următorul pas
al ciclului.

1.4 Funcţii(subprograme) în Matlab

O funcţie Matlab este corect definită dacă produce rezultate pentru diferite
tipuri de variabile: numere, vectori, matrice. În acest scop, uzual operaţiile de
înmulţire, împărţire şi ridicare la putere se folosesc în varianta cu punct, adică
operaţiile se fac pe componente.
Funcţiile Matlab predefinite satisfac această cerinţă.
Funcţiile se pot defini în
1.4. FUNCŢII(SUBPROGRAME) ÎN MATLAB 19

• Linie de comandă.

– Funcţie anonimă

f nt = @(x)expr esi e

– Funcţie definită prin intermediul unui string

st r = 0 expr esi e 0
f nt = str2func([0 @(x)0 , st r ])

Funcţia f : R2 → R2 definită prin

µ ¶
x+y
f (x, y) =
x2 + y 2

se defineşte prin

– f=@(x,y) [x+y,x.^2+y.^2]

f =

@(x,y)[x+y,x.^2+y.^2]

>> f(2,3)

ans =

5 13
– str=’[x+y,x.^2+y.^2]’

str =

[x+y,x.^2+y.^2]

>> g=str2func([’@(x,y)’,str])

g =
20 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

@(x,y)[x+y,x.^2+y.^2]

>> g(2,3)

ans =

5 13

Variabilele funcţiei pot fi numere, vectori, matrice. În cazul în care vari-

abilele sunt tablouri, operaţiile algebrice se fac între componentele core-
spunzătoare.

• Exterior. Cu un editor de fişiere se crează subprogramul

function [y 1 , . . . , y n ]= f (x 1 , . . . , x m )
instructiuni Matlab

care se salvează într-un fişier obligatoriu cu numele f .m.

Ieşirea dintr-o funcţie se poate forţa prin instrucţiunea return.
Pentru a putea fi utilizată funcţia exterioară se va salva într-un catalog declarat
în fişerul matlabpath.m. Gestionarea acestui fişier se realizează apelând
meniul Set path.
Pentru funcţia din exemplul anterior, codul exemplului anterior este

function [u,v]=f(x,y)
u=x+y;
v=x^2+y^2;

cu utilizarea

[u,v]=f(2,3)

u =

v =
1.4. FUNCŢII(SUBPROGRAME) ÎN MATLAB 21

13

Exemplul 1.7 Funcţie pentru rezolvarea ecuaţiei de gradul al doilea.

function [u,v]=eq2(a,b)
delta=a*a-4*b;
if delta>=0
u=0.5*(-a+sqrt(delta));
v=0.5*(-a-sqrt(delta));
disp(’Real roots’)
else
u=-a/2;
v=sqrt(-delta);
disp(’Complex conjugate roots’)
end
Funcţia disp(obiectMatlab) afişează valoarea argumentului.

Exemplul 1.8 Funcţie pentru calculul valorii unui polinom.

function r=polyval1(x,c)
[l,n]=size(c)
r=c(1)
for i=2:n
r=r*x+c(i)
end
endfunction
function r=polyval2(x,c)
[l,n]=size(c)
i=0:n-1
u=x^i
r=u*c’
endfunction
Exemplul 1.9 Calculul coeficienţilor binomiali.

Variantă bazată pe formula de recurenţă (n ≥ k)


µ ¶  1 n = 0 sau k = 0
n µ ¶ µ ¶
= n − 1 n − 1
k + altfel
k k −1

22 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

function c=binomialcoeff(n,k)
c=0;
if k>n
return
end
if n==0 & k>0
c=1;
else
if k==0
c=1;
else
c=binomialcoeff(n-1,k-1)+binomialcoeff(n-1,k);
end
end
µ ¶
n n!
Varianta bazată pe formula = k!(n−k)! = exp(gammaln(n+1)−gammaln(k+
k
1) − gammaln(n − k + 1)).

function c=binomialcoeff1(n,k)
c=0;
if k>n then
return
end
c=exp(gammaln(n+1)-gammaln(k+1)-gammaln(n-k+1));

Exemplul 1.10 Funcţie pentru calculul primelor n termeni ai şirului Fibonacci.

function t=fib(n)
t=[1,1];
for i=3:n
t(i)=t(i-1)+t(i-2);
end

Exemplul
p 1.11 Funcţie pentru calculul limitei şirului (a n )n∈N definit prin a 0 =
p
2, a n+1 = a n + 2.

function [l,error]=lim(tol,nmi)
v=sqrt(2);
cont=true;
1.4. FUNCŢII(SUBPROGRAME) ÎN MATLAB 23

ni=1;
while(cont)
ni=ni+1;
l=sqrt(v+2);
d=abs(v-l);
v=l;
if(d<tol)|(ni>nmi)
cont=false;
end
end
if(d<tol)
error=0;
else
error=1;
end

O functie Matlab poate avea un număr variabil de parametri formali de in-

trare. Aceşti parametri se indică prin variabila varargin. Numărul argumentelor
suplimentare este dat de nargin (nar g i n = length(v ar ar g i n)).

Exemplul 1.12 Procedură pentru generarea nodurilor lui Cebîşev de speţa a doua
dintr-un interval [l , r ]. În lipsa parametrilor suplimentari intervalul implicit este
[−1, 1].

În intervalul [−1, 1], pentru n ∈ N∗ , nodurilor lui Cebîşev de speţa a doua sunt
definite ca punctele de extrem ale polinomului lui Cebîşev Tn (x) : cos kπ n
, k =0:
n.
Bijecţia între intervalele [−1, 1] şi [l , r ] este x 7→ l + r −l
2 (x + 1).

function x=chebpoints(n,varargin)
l=-1;
r=1;
if length(varargin)>0
if length(varargin)==2
l=varargin{1};
r=varargin{2};
disp([l,r])
else
error(’Wrong number of input parameters’)
end
end
k=0:n;
x=l+0.5*(r-l)*(cos(k*pi/n)+1);
24 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

1.5 Grafică în Matlab

1.5.1 Grafică 2D
Din mulţimea funcţiilor grafice ale produsului Matlab amintim

• plot. Cea mai simplă formă de utilizare este

pl ot (x, y)

unde x = (x i )1≤i ≤n , y = (y i )1≤i ≤n sunt doi vectori. Componentele lui x sunt

ordonate crescător. Se reprezintă funcţia dată de (x i , y i )i .
Alte formate de apelare sunt

pl ot (x, y, specLi ni e)pl ot (x 1 , y 1 , specLi ni e 1 , . . . , x n , y n , specLi ni e n )

Variabila specLinie este un string care precizează

– stilul
‘-’ linie continuă
‘–’ linie întreruptă
‘:’ curba reprezentată prin puncte
‘-.’ curbă reprezentată prin combinaţia segment punct
– marcaj (selecţe)
‘o’ cerc
‘+’ +
‘*’ asterisc
‘x’ x
‘s’ patrat
‘d’ diamant
– Culoare
‘r’ roşu
‘g’ verge
‘b’ albastru
‘c.’ cyan
‘m’ violet
‘y’ galben
‘k’ negru
‘w’ alb
1.5. GRAFICĂ ÎN MATLAB 25

Exemplul 1.13 Reprezentarea grafică a funcţiilor sin x şi cos x în [−π, π]

Reprezentarea grafică se obţine prin

t=-pi:0.1:pi;
x1=sin(t);
x2=cos(t);
plot(t’,x1’,t’,x2’,’r’)

Va rezulta imaginea din Fig. 1.1

Fig. 1.1: Grafică 2D: plot.

• fplot: Reprezentarea grafică într-un interval a unei funcţii. Utilizarea

funcţiei este
fplot( f )
fplot( f , i nt er v al )
26 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

În cazul neprecizării intervalului acesta este [−5, 5].

2
Exemplul 1.14 Reprezentarea funcţiei f (x) = xe −x în intervalul [−3, 3].

Scriptul petru obţinerea reprezentării grafice este

f=@(x)x.*exp(-x.^2);
interval=[-3,3];
fplot(f,interval)

Rezultă imaginea din Fig. 1.2.

Fig. 1.2: Grafică 2D: fplot.

Reprezentarea funcţiilor parametrice.

Exemplul 1.15 Reprezentarea funcţiei γ(t ) = (si n(2t ), cos(3t )).

1.5. GRAFICĂ ÎN MATLAB 27

syms t
x=sin(2*t);
y=cos(3*t);
fplot(x,y)

Rezultă imaginea din Fig. 1.3.

Fig. 1.3: Grafică 2D: fplot.

• quiver: Reprezentarea unui câmp vectorial U = u(x, y),V = v(x, y) într-

un interval bidimensional.
Sintaxa de utilizare

quiver(X , Y ,U ,V )

Şablonul de utilizare este indicat în exemplul următor.

28 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

Exemplul 1.16 Să se obţină imaginea câmpului vectorial (|y|, x − 2y) în in-
tervalul [−2, 2]2 .

Scriptul de rezolvare este

u=@(x,y)abs(y);
v=@(x,y)x-2*y;
x=linspace(-2,2,21);
y=linspace(-2,2,21);
[X,Y]=meshgrid(x,y);
U=u(X,Y);
V=v(X,Y);
quiver(X,Y,U,V)

Rezultă imaginea din Fig. 1.4.

Funcţiile meshgrid şi ndgrid:

Dacă

x=[1,2,3];
y=[10,20];

atunci
[X,Y]=ndgrid(x,y)
[U,V]=meshgrid(x,y)
X =
U =
1 1
2 2
1 2 3
3 3
1 2 3

V =
Y =
10 10 10
10 20
20 20 20
10 20
10 20
Aşadar meshgrid(x, y) = ndgrid(x 0 , y 0 ).
1.5. GRAFICĂ ÎN MATLAB 29

Fig. 1.4: Grafică 2D: quiver.

1.5.2 Grafică 3D
• mesh(X,Y,Z): Reprezentarea funcţiilor f : D ⊂ R2 → R.
X = (X i )1≤i ≤m şi Y = (Y1≤ j ≤n sunt vectori daţi în ordine crescătoare, iar
Z = (Zi , j ), i ∈ {1, . . . , m}, j ∈ {1, . . . , n} reprezintă valoarea funcţiei calculată
în (X i , Y j ).
Dacă suprafaţa este definită prin z = f (x, y), x ∈ [a, b], y ∈ [c, d ] atunci se
defineşte funcţia Matlab

function [x,y,z]=suprafata(m,n)
x=linspace(a,b,m);
y=linspace(c,d,n);
z=zeros(m,n)
for i=1:m
30 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

for j=1:n
z(i,j)=f(x(i),y(j));
end
end

[x,y,z]=suprafata(m,n);
mesh(x,y,z)

Exemplul 1.17 Reprezentarea semisferei x 2 + y 2 + z 2 = 1, z ≥ 0.

p
Semisfera superioară este dată de z = 1 − x 2 − y 2 , pentru care construim
funcţia Matlab

function[x,y,z]=sfera(r,m,n)
x=linspace(-r,r,m);
y=linspace(-r,r,n);
for i=1:m
for j=1:n
if x(i)^2+y(j)^2<=r^2
z(i,j)=sqrt(r^2-x(i)^2-y(j)^2);
else
z(i,j)=NaN;
end
end
end

Comenzile de execuţie sunt

[x,y,z]=sfera(1,20,20);
mesh(x,y,z)

Imaginea obţinută este dată în Fig. 1.5.

• fplot3: Reprezentarea grafică a curbelor în R3 .

Presupunem că x, y, z sunt funcţii Matlab definite într-un interval.

f pl ot 3(x, y, z)
f pl ot 3(x, y, z, i nt er v al )

unde x = X (t ), y = Y (t ), z = Z (t ), t ∈ i nt er v al . Intervalul implicit este [−5, 5].

Exemplul 1.18 Reprezentarea curbei γ(t ) = (t 2 sin(t ), t 2 cos(t ), t ), t ∈ [−5.5].

1.5. GRAFICĂ ÎN MATLAB 31

Fig. 1.5: Grafică 3D: mesh

x=@(t)t.^2.*sin(t);
y=@(t)t.^2.*cos(t);
z=@(t)t;
fplot3(x,y,z)

Se obţne imaginea din Fig. 1.6.

Calculând valorile funcţiilor x, y, z pe o reţea de puncte, stocate în tablourile
X , Y , Z , funcţia plot3(X , Y , Z ) generează imaginea graficului curbei.

1.5.3 Animaţie în Matlab

• comet. Reprezentare grafică animată.

Exemplul 1.19 Generarea dinamică a graficului funcţiei f (t ) = t 2 .

32 CAPITOLUL 1. ELEMENTE DE PROGRAMARE ÎN MATLAB

Fig. 1.6: Grafică 3D: mesh

t=linspace(0,1,10000);
comet(t,t.^2)

O valoare mai mare pentru al treilea argument al funcţiei linspace implică

o viteză mai mică a animaţiei.

• Funcţie care crează un fişier gif animat pentru desenarea graficul unei
funcţii f : [a, b] → R.

1 function t o g i f ( f , a , b , path , filename )

2 x = a:0.01:b;
3 figure ( 1 )
4 f i l e =[ path filename ] ;
5 t = linspace ( a , b) ;
6 x = f(t) ;
7 for n = 1 : length ( t )
1.5. GRAFICĂ ÎN MATLAB 33

8 plot ( t ( 1 : n ) , x ( 1 : n ) )
9 axis ( [ min( t ) max( t ) min( x ) max( x ) ] ) ;
10 drawnow
11 frame = getframe ( 1 ) ;
12 im = frame2im ( frame ) ;
13 [ imind ,cm] = rgb2ind (im , 2 5 6 ) ;
14 i f n == 1 ;
15 imwrite ( imind ,cm, f i l e , ’ g i f ’ , ’ Loopcount ’ , i n f ) ;
16 el se
17 imwrite ( imind ,cm, f i l e , ’ g i f ’ , ’ WriteMode ’ , ’append ’ ) ;
18 end
19 end
Capitolul 2

Algebră liniară numerică

Cap. 1 Cap. 3

2.1 Factorizarea unei matrice

Funcţiile

• [L,U , P ] = lu(A) calculează factorizarea LU a matricei A : P A = LU .

• [Q, R] = qr(A) calculează factorizarea QR a matricei A : A = QR

Exemplul 2.1 Să se calculeze factorizarea LU a matricei

 
1 2 −1 3 2

 2 4 −2 5 1 

A= −1 −2 1 −3 −4  .
 
3 6 2 10 7 
 

1 2 4 0 4

A=[1,2,-1,3,2;2,4,-2,5,1;-1,-2,1,-3,-4;3,6,2,10,7;1,2,4,0,4];
[L,U,P]=lu(A)
L =

34
2.1. FACTORIZAREA UNEI MATRICE 35

1.0000 0 0 0 0
0.6667 1.0000 0 0 0
0.3333 0 1.0000 0 0
0.3333 0 -0.5000 1.0000 0
-0.3333 0 0.5000 -1.0000 1.0000

U =

3.0000 6.0000 2.0000 10.0000 7.0000

0 0 -3.3333 -1.6667 -3.6667
0 0 3.3333 -3.3333 1.6667
0 0 0 -2.0000 0.5000
0 0 0 0 -2.0000

P =

0 0 0 1 0
0 1 0 0 0
0 0 0 0 1
1 0 0 0 0
0 0 1 0 0

Verificarea rezultatului se face calculând P ∗ A − L ∗U , rezultat care trebuie să fie

matricea nulă.

Exemplul 2.2 Să se calculeze factorizarea QR a matricei

 
6 6 1
X = 3 6 1 
2 1 1

X=[6,6,1;3,6,1;2,1,1];

[Q,R]=qr(X)
Q =

-0.8571 0.2857 -0.4286

-0.4286 -0.8571 0.2857
-0.2857 0.4286 0.8571
36 CAPITOLUL 2. ALGEBRĂ LINIARĂ NUMERICĂ

R =

-7.0000 -8.0000 -1.5714

0 -3.0000 -0.1429
0 0 0.7143

Verificarea
Q*R-X

ans =

1.0e-14 *

0 0.1776 0
0 0.0888 0.0222
0.0444 0.0444 0

2.2 Rezolvarea sistemelor algebrice de ecuaţii liniare

Pentru rezolvarea unui sistem algebric de ecuaţii liniare Ax = b, în cazul în
care

• numărul ecuaţiilor coincide cu numărul necunoscutelor;

• determinantul sistemului este diferit de zero;

se procedează după cum urmează:

1. se fixează matricele A şi b;

2. se calculează folosind una din variantele următoare

• x = A \b
• x = A −1 b
Expresia A −1 se poate introduce tastând Aˆ(-1).
• x = inv(A) ∗ b
• x = linsolve(A, b)

Funcţia Matlab det( A ) calculează determinantul matricei A.

2.2. REZOLVAREA SISTEMELOR ALGEBRICE DE ECUAŢII LINIARE 37

Exemplul 2.3 Să se rezolve sistemul algebric de ecuaţii liniare

x1 + 2x 2 + 3x 3 + 4x 4 11

 =

2x 1 + 3x 2 + 4x 3 + x4 = 12

 3x 1
 + 4x 2 + x3 + 2x 4 = 13

4x 1 + x2 + 2x 3 + 3x 4 = 14

Rezolvarea constă din

A = [1,2,3,4;2,3,4,1;3,4,1,2;4,1,2,3];
b = [11;12;13;14];
det(A)
ans=
160
x=inv(A)*b
x =

2.0000
1.0000
1.0000
1.0000

Cazul sistemelor incompatibile

În cazul unui sistem algebric de ecuaţii liniare incompatibil se poate calcula
elementul care minimizează funcţionala J (x) = kb − Axk22 . Acest element este
soluţia sistemului algebric de ecuaţii liniare compatibil A T Ax = A T b. Matricea
A † = (A T A)−1 A T este pseudoinversa matricei A.
Funcţia Matlab pinv(X ) calculează matricea X † - pseudoinversa matricei X ∈
M n,k (R).

Exemplul 2.4 Să se rezolve sistemul algebric de ecuaţii liniare


 2x − y + 3z = 7
x + y + z = 4
3x − 3y + 5z = 8


în sensul celor mai mici pătrate.

Rezolvarea este
38 CAPITOLUL 2. ALGEBRĂ LINIARĂ NUMERICĂ

A=[2,-1,3;1,1,1;3,-3,5];
b=[7;4;8];
x=pinv(A)*b
x =

1.4872
1.2949
1.5513
Capitolul 3

Rezolvarea sistemelor şi ecuaţiilor

algebrice

Cap. 2 Cap. 4

3.1 Rezolvarea sistemelor algebrice de ecuaţii neliniare

Pentru rezolvarea sistemului algebric de ecuaţii neliniare

 f 1 (x 1 , . . . , x n ) = 0

..
 .
f n (x 1 , . . . , x n ) = 0


sau scris sub formă concentrată f (x) = 0 cu

   
x1 f 1 (x 1 , . . . , x n )
 .  ..
x =  ..  f (x) =  ,
 
.
xn f n (x 1 , . . . , x n )

se utilizează funcţia Matlab

[x [,fval [,exitflag ] ] ]=fsolve(f,x0)

39
40 CAPITOLUL 3. REZOLVAREA SISTEMELOR ŞI ECUAŢIILOR ALGEBRICE

unde semnificaţia parametrilor este:

• x0 reprezintă o aproximaţie iniţială a soluţiei sistemului.

• f identificatorul funcţiei Matlab care defineşte sistemul neliniar

f (x) = 0.

• x aproximaţia calculată a soluţiei sistemului algebric.

• f v al reprezintă valoarea funcţiei f calculată în x; f v al = f (x).

• exi t f l ag indicatorul de răspuns al programului fsolve. În cazul rezolvării

cu succes, valoarea indicatorului este cuprins între 1 şi 4.

Astfel, rezolvarea constă din

1. definirea funcţei f ;

2. apelarea funcţiei Matlab fsolve.

Exemplul 3.1 Să se rezolve sistemul algebric de ecuaţii neliniare

 2
 10x 1 + x 1 − 2x 2 x 3 − 0.1 = 0
10x 2 − x 22 + 3x 1 x 3 + 0.2 = 0
10x 3 + x 32 + 2x 1 x 2 − 0.3 = 0


Rezolvarea constă din:

1. Definirea funcţiei (externe) care defineşte sistemul algebric

function q=fct(p)
x=p(1);
y=p(2);
z=p(3);
q(1)=10*x+x.^2-2*y.*z-0.1;
q(2)=10*y-y.^2+3*x.*z+0.2;
q(3)=10*z+z.^2+2*x.*y-0.3;

2. Rezolvarea propriu-zisă
3.2. REZOLVAREA ECUAŢIILOR ALGEBRICE 41

p0=[0,0,0];
[p,fval,exitflag]=fsolve(@fct,p0)

Equation solved.

fsolve completed because the vector of function values is near zero

as measured by the default value of the function tolerance, and
the problem appears regular as measured by the gradient.

<stopping criteria details>

p =

0.0099 -0.0200 0.0299

fval =

1.0e-07 *

0.1200 0.1717 0.1510

3.2 Rezolvarea ecuaţiilor algebrice

Dacă n = 1, adică dimensiunea spaţiului este 1, atunci f (x) = 0 cu f : I ⊂ R →
R reprezintă o ecuaţie algebrică.
Pentru rezolvarea ecuaţiilor algebrice se utilizează una din funcţiile Matlab
fsolve sau fzero.

Exemplul 3.2 Să se rezolve ecuaţia 2x = x 2 .

Rezolvarea ecuaţiei f (x) = 2x − x 2 = 0 este

f=@(x)2.^x-x.^2

f =

@(x)2.^x-x.^2

[x,fval,exitflag]=fsolve(f,1)
42 CAPITOLUL 3. REZOLVAREA SISTEMELOR ŞI ECUAŢIILOR ALGEBRICE

Equation solved.

fsolve completed because the vector of function values is near zero

as measured by the default value of the function tolerance, and
the problem appears regular as measured by the gradient.

<stopping criteria details>

x =

fval =

exitflag =

1
Dacă x0 = 5 atunci x = 4, iar pentru x0 = −1 găsim x = −0.7667.
Graficul funcţiei f este reprezentat în Fig. 3.1.
Utilizarea funcţiei fzero este identică cu cea a funcţiei fsolve.

3.3 Rezolvarea ecuaţiilor polinomiale

În cazul particular al ecuaţiilor polinomiale se cere determinarea tuturor rădăcinilor
reale sau complexe.
Pentru aflarea rădăcinilor unui polinom se procedează astfel:

1. Se defineşte un vector p cu coeficienţii polinomului, daţi în ordine de-

screscătoare a gradului.

2. Rădăcinile polinomului se calculează cu ajutorul funcţiei Matlab

x=roots(p)
3.3. REZOLVAREA ECUAŢIILOR POLINOMIALE 43

Fig. 3.1: y = 2x − x 2

unde x este un tablou cu rădăcinile polinomului p.

Exemplul 3.3 Să se determine rădăcinile polinomului p = x 4 + 2x 3 + 3x 2 + 2x + 1.

Rezolvarea constă din:

c=[1,2,3,2,1];
x=roots(p)
x =

1.0156 + 1.4152i
1.0156 - 1.4152i
 Capitolul 4

Rezolvarea problemelor de
interpolare

Cap. 3 Cap. 5

Fie F o familie interpolatoare de ordin n pe axa reală. Dându-se nodurile

(x i )1≤i ≤n şi numerele (y i )1≤i ≤n , dacă ϕ ∈ F este funcţia de interpolare care satis-
face condiţiile ϕ(x i ) = y i , 1 ≤ i ≤ n, se cere să se calculeze ϕ(z), unde z este un
punct dat.

4.1 Interpolare polinomială

Pentru F = P n−1 soluţia problemei de interpolare Lagrange este polinomul
n n z −x
X Y j
L(P n−1 ; x 1 , . . . , x n ; y 1 , . . . , y n )(z) = yi
i =1 j =1 xi − x j
j 6=i

Valoarea acestui polinom este calculat de funcţia Matlab prin formula baricen-
trică

1 function f =lagrange ( xd , x , y )

44
 4.1. INTERPOLARE POLINOMIALĂ 45

2 [mx, nx]= s i z e ( x ) ;
3 [my, ny]= s i z e ( y ) ;
4 i e r r o r =0;
5 i f ( nx~=ny ) | ( mx~ = 1 ) | (my~=1) ,
6 i e r r o r =1;
7 disp ( ’ data dimension e r r o r ’ )
8 abort
9 end
10 xx=sort ( x ) ;
11 for k =1: nx−1,
12 i f xx ( k)== xx ( k +1) ,
13 i e r r o r =1;
14 break ,
15 end
16 end
17 i f i e r r o r ~=0 ,
18 disp ( ’ data e r r o r ’ )
19 abort
20 end
21 [m, n]= s i z e ( xd ) ;
22 f =zeros (m, n ) ;
23 p=zeros (m, n ) ;
24 q=zeros (m, n ) ;
25 w=ones ( 1 , nx ) ;
26 for i =1:nx ,
27 for j =1:nx ,
28 i f i ~=j ,
29 w( i )=w( i ) * ( x ( i )−x ( j ) ) ;
30 end
31 end
32 end
33 for i =1:m,
34 for j =1:n ,
35 u=find ( x==xd ( i , j ) ) ;
36 i f ~isempty (u ) ,
37 f ( i , j )= y (u ) ;
38 el se
39 for k =1: nx ,
40 p( i , j )=p( i , j )+ y ( k ) / ( xd ( i , j )−x ( k ) ) /w( k ) ;
41 q ( i , j )=q ( i , j ) + 1 / ( xd ( i , j )−x ( k ) ) /w( k ) ;
42 end
43 f ( i , j )=p( i , j ) / q ( i , j ) ;
44 end
45 end
46 end

Semnificaţiile parametrilor formali şi a rezultatului sunt:

• xd este o matrice de numere. În fiecare element al matricei xd se cal-

culează valoarea polinomului de interpolare Lagrange.

• x = (x i )1≤i ≤n , şirul nodurilor de interpolare;

• y = (x i )1≤i ≤n , şirul valorilor interpolate;

• f este o matrice de aceleaşi dimensiuni ca xd şi are ca elemente valorile

 46 CAPITOLUL 4. REZOLVAREA PROBLEMELOR DE INTERPOLARE

polinomului de interpolare Lagrange calculate în elementele corespunză-

toare ale matricei xd .
Programul de mai sus, pentru calculul valorii polinomului de interpolare La-
grange, foloseşte formula baricentrică.
Exemplul 4.1 Să se calculeze L(P n ; x 0 , . . . , x n ; f )(z) pentru
f (x) = x 2
xi = i , i ∈ {0, 1, . . . , 5}
z = 0.5, 3, 5, 10
Rezolvarea este:

fct=@(x)x.^2
x=0:5;
y=fct(x);
z=[0.5,3.5,10];
f=lagrange(z,x,y)
f =

0.2500 12.2500 100.0000

Calculul simbolic al polinomului de interpolare Lagrange se calculează cu
funcţia

1 function lag =polyLagrange ( x , y )

2 syms s
3 z=poly2sym ( [ 1 , 0 ] , s ) ;
4 [ k , n]= s i z e ( x ) ;
5 [ k ,m]= s i z e ( y ) ;
6 i f m~=n
7 disp ( ’ Data e r r o r ’ )
8 return
9 end
10 v=sym( ’ v ’ , [ 1 , n ] ) ;
11 w=sym( ’w ’ , [ 1 , n ] ) ;
12 for i =1:n
13 v ( i )=sym( y ( i ) ) ;
14 end
15 for k =1:n−1
16 for i =1:n−k
17 w( i ) = ( ( z−x ( i ) ) * v ( i +1)−(z−x ( i +k ) ) * v ( i ) ) / ( x ( i +k)−x ( i ) ) ;
18 end
19 for i =1:n−k
20 v ( i )=w( i ) ;
21 end
22 end
23 lag =expand ( v ( 1 ) ) ;
4.2. INTERPOLARE CU FUNCŢII SPLINE CUBICE 47

Exemple de utilizare

x=[-2,-1,0,1,2];
y=[4,1,0,1,4];
polyLagrange(x,y)

ans =

s^2

z=[2,1,0,1,2];
polyLagrange(x,z)

ans =

- s^4/6 + (7*s^2)/6

syms s x
subs(z,s,x)

ans =

- x^4/6 + (7*x^2)/6

4.2 Interpolare cu funcţii spline cubice

Dacă F = S 3 , mulţimea funcţiilor spline cubice atunci pentru rezolvarea prob-
lemei de interpolare utilizăm funcţia Matlab

s = spline(x, y, z)

Semnificaţia parametrilor este

• x = (x i )1≤i ≤n , şirul nodurilor de interpolare;

• z punclele în care se calculează aproximaţia dată de funcţia spline;

• s aproximaţiile obţnute.
48 CAPITOLUL 4. REZOLVAREA PROBLEMELOR DE INTERPOLARE

Exemplul 4.2 Să se calculeze valorile funcţiei spline cubice de interpolare şi ale
derivatelor sale de ordin 1, 2, 3 pentru datele de interpolare

f (x) = x 3
xi = i , i ∈ {0, 1, . . . , 5}
z = 0.5, 3, 5, 10

Rezolvarea constă din

f=@(x)x.^3

f =

@(x)x.^3

>> x=0:5;
>> y=f(x);
>> z=[0.5,3.5,10];
>> s=spline(x,y,z)

s =

0.1250 42.8750 1000.0000

Capitolul 5

Derivare numerică

Cap. 4 Cap. 6

5.1 Derivarea funcţiilor

Matlab dispune de funcţia gradient care permite

• calculul derivatelor unei funcţii reale de o variabila reală;

Pe o reţea echidistantă de puncte x i = x 0 + i h, i ∈ {0, 1, . . . , n} se calculează
vectorul y = (y i )i , y i = f (x i ) şi se apelează

d y = gradient(y, h)

• calculul gradientului în cazul funcţiilor reale de mai multe variabile.

Considerăm cazul unei funcţii cu două variabile. Fie reţeaua bidimension-
ală definită prin şirurile echidistante x = (x i )0≤i ≤n şi y = (y j )0≤ j ≤n cu ratele
hx, h y. Se execută comenzile

[xx,yy]=meshgrid(x,y)
z=f(xx,yy)
[dzx,dzy]=gradient(f,hx,hy)

49
50 CAPITOLUL 5. DERIVARE NUMERICĂ

Elementul z i , j este f (x j , y i ).

Exemplul 5.1 Să se calculeze derivatele funcţiei f (x) = x 3 pentru x = −1 : 0.2 : 1.

f=@(x)x.^3;
x=-1:0.2:1;
y=f(x);
dy=gradient(y,0.2)
dy =

Columns 1 through 5

2.4400 1.9600 1.1200 0.5200 0.1600

Columns 6 through 10

0.0400 0.1600 0.5200 1.1200 1.9600

Column 11

2.4400
Observaţie 5.1

1. Pentru i = 0 se utilizează formula

4h f (x 0 ) f (x 1 ) − f (x 0 ) −0.512 − (−1)
d y0 = = = = 2.44.
h h 0.2
2. Pentru i ∈ {1, 2, . . . , n − 1} se utilizează formula
δ2h f (x i ) f (x i +1 ) − f (x i −1 )
d yi = =
2h 2h
Verificăm pentru i = 5
δ2h f (x 5 ) f (x 6 ) − f (x 4 )) 0 − (−0, 064
d y5 = = = = 0.16.
2h 2h 0.4
3. Pentru i = n se utilizează formula
∇h f (x n ) f (x n ) − f (x n−1 )
d yn = = =
h h
f (x 10 ) − f (x 9 ) 1 − 0.512
= = = 2.44.
h 0.2
5.1. DERIVAREA FUNCŢIILOR 51

Exemplul 5.2 Să se calculeze gradienţii funcţiei f (x) = x 3 y pe reţeaua de puncte

x = −0.4 : 0.2 : 0.4, y = 0 : 0.1 : 0.4.
f=@(x,y)x.^3.*y;
>> x=-0.4:0.2:0.4

x =

-0.4000 -0.2000 0 0.2000 0.4000

y=0:0.1:0.4

y =

0 0.1000 0.2000 0.3000 0.4000

[xx,yy]=meshgrid(x,y)

xx =

-0.4000 -0.2000 0 0.2000 0.4000

-0.4000 -0.2000 0 0.2000 0.4000
-0.4000 -0.2000 0 0.2000 0.4000
-0.4000 -0.2000 0 0.2000 0.4000
-0.4000 -0.2000 0 0.2000 0.4000

yy =

0 0 0 0 0
0.1000 0.1000 0.1000 0.1000 0.1000
0.2000 0.2000 0.2000 0.2000 0.2000
0.3000 0.3000 0.3000 0.3000 0.3000
0.4000 0.4000 0.4000 0.4000 0.4000

>> z=f(xx,yy)

z =

0 0 0 0 0
52 CAPITOLUL 5. DERIVARE NUMERICĂ

-0.0064 -0.0008 0 0.0008 0.0064

-0.0128 -0.0016 0 0.0016 0.0128
-0.0192 -0.0024 0 0.0024 0.0192
-0.0256 -0.0032 0 0.0032 0.0256

[dzx,dzy]=gradient(z,0.2,0.1)

dzx =

0 0 0 0 0
0.0280 0.0160 0.0040 0.0160 0.0280
0.0560 0.0320 0.0080 0.0320 0.0560
0.0840 0.0480 0.0120 0.0480 0.0840
0.1120 0.0640 0.0160 0.0640 0.1120

dxy =

-0.0640 -0.0080 0 0.0080 0.0640

-0.0640 -0.0080 0 0.0080 0.0640
-0.0640 -0.0080 0 0.0080 0.0640
-0.0640 -0.0080 0 0.0080 0.0640
-0.0640 -0.0080 0 0.0080 0.0640
Au loc aproximaţiile
∂ f (x j , y i ) ∂ f (x j , y i )
d zx i , j ≈ , d z yi , j ≈
∂x ∂y
şi formulele
i =0 j =0
f (x 1 ,y j )− f (x 0 ,y j )) f (x i ,y 1 )− f (x i ,y 0 ))
d zx(0, j ) = hx , d z y(i , 0) = hy ,

i ∈ {1, 2, . . . , n − 1} j ∈ {1, 2, . . . , n − 1}
f (x i +1 ,y j )− f (x i −1 ,y j )) f (x i ,y j +1 )− f (x i ,y j −1 ))
d zx(i , j ) = 2hx , d z y(i , j ) = 2h y ,

i =n j =n
f (x n ,y j )− f (x n−1 ,y j )) f (x i ,y n )− f (x i ,y n−1 ))
d zx(n, j ) = hx
, d z y(i , n) = hy
,

j ∈ {0, 1, . . . , n}, i ∈ {0, 1, . . . , n}.

5.1. DERIVAREA FUNCŢIILOR 53

Aceste formule nu ţin seama de faptul că în Matlab indicele iniţial este 1.
Capitolul 6

Construirea unei funcţii

de aproximare prin metoda
celor mai mici pătrate

Cap. 5 Cap. 7

Cazul liniar. Determinarea unui polinom de aproximare construit prin metoda

celor mai mici pătrate de grad m pentru datele (x i , y i )1≤i ≤n (m << n) se obţine
cu ajutorul funcţiei

coef =lq(m,x,y)

unde

• coe f – sunt coeficienţii polinomului în ordinea crescătoare a gradului;

• m – gradul polinomului de aproximare;

• x, y – doi vectori linie cu absisele şi respectiv, cu ordonatele datelor proble-

mei de aproximare.

Textul sursă al funcţiei lq este

54
 55

1 function coef=lq (m, x , y )

2 [mx, nx]= s i z e ( x ) ;
3 [my, ny]= s i z e ( y ) ;
4 i f ( ( nx~=ny ) | ( mx~ = 1 ) | (my~ =1 ) ) ,
5 disp ( ’ data dimension e r r o r ’ )
6 abort
7 end
8 u=zeros (m+1 ,nx ) ;
9 for i =1:m+1
10 for j =1: nx
11 u( i , j )= x ( j ) ^ ( i −1);
12 end
13 end
14 coef =(u * u’ ) ^ ( − 1 ) * u * y ’ ;

Exemplul 6.1 Să se calculeze polinoamele de aproximare de grad unu şi doi, con-
stituite prin metoda celor mai mici pătrate, pentru datele (x i , y i )0≤i ≤20 unde x i =
−2 + 0.2i , y i = f (x i ), f (x) = |x|. Să se reprezinte grafic funcţiile astfel obţinute.

x=-2:0.2:2;
y=abs(x);
c1=lq(1,x,y)

c1 =

1.0476
0.0000

>> c2=lq(2,x,y)

c2 =

0.3884
0.0000
0.4495

>> c3=lq(3,x,y)

c3 =

0.3884
56 CAPITOLUL 6. METODA CELOR MAI MICI PĂTRATE

0.0000
0.4495
0.0000
Reprezentarea grafică a celor trei aproximaţii
p1=@(x)c1(1)+x*c1(2);
p2=@(x)c2(1)+x*c2(2)+c2(3)*x.^2;
p3=@(x)c3(1)+x*c3(2)+c3(3)*x.^2+c3(4)*x.^3;
t=linspace(-2,2,100);
z=abs(t);
y1=p1(t);
y2=p2(t);
y3=p3(t);
plot(t,z,t,y1,t,y2,t,y3)
este dată în Fig. 6.1.

Fig. 6.1: Aproximarea funcţiei |x| prin polinome de grad 1,2 şi 3.

Pentru rezolvarea acestei probleme Matlab posedă funcţia

r ez f i t = fit(xD at a, yD at a, f ,0 St ar t Poi nt 0 , v)
 57

unde

• r ez f i t cu informaţii despre problemă şi cu rezultatele obţinute;

• xD at a, yD at a vectori coloană cu datele problemei de aproximare;

• f funcţie cu şablonul funcţiei de aproximare. Se obţine cu

f=fittype(‘sablon’)

• v vector linie cu valorile iniţiale ale parametrilor din şablon.

Rezolvarea exemplului anterior este

f=fittype(’a+b*x+c*x.^2’)

f =

General model:
f(a,b,c,x) = a+b*x+c*x.^2

rezfit=fit(x,y,f,’StartPoint’,[1,1,1])

rezfit =

General model:
fit1(x) = a+b*x+c*x.^2
Coefficients (with 95% confidence bounds):
a = 0.3884 (0.2735, 0.5033)
b = 7.282e-09 (-0.06313, 0.06313)
c = 0.4495 (0.391, 0.508)

Cazul neliniar.
Pentru determinarea funcţiei de aproximare y = ϕ(t , x 1 , . . . , x n ) care minimizează
expresia
m
Φ(x 1 , . . . , x n ) = [ϕ(t i , x 1 , . . . , x n ) − y i ]2
X
(6.1)
i =1

cu ajutorul funcţiei fit. Funcţia (x 1 , . . . , x n ) 7→ ϕ(t , x 1 , . . . , x n ) defineşte şablonul.

58 CAPITOLUL 6. METODA CELOR MAI MICI PĂTRATE

Exemplul 6.2 Să se calculeze funcţia de aproximare de forma y = ae bt în cazul

datelor

t=0:5;
y=2*exp(-t);

Rezolvarea este:

f=fittype(’a*exp(b*x)’)

f =

General model:
f(a,b,x) = a*exp(b*x)

x=(0:5)’;
y=2*exp(-x);
rezfit=fit(x,y,f,’StartPoint’,[1,1])

rezfit =

General model:
rezfit(x) = a*exp(b*x)
Coefficients (with 95% confidence bounds):
a = 2 (2, 2)
b = -1 (-1, -1)
Capitolul 7

Integrare numerică

Cap. 6 Cap. 8

7.1 Integrarea numerică a funcţiilor

de o variabilă reală
Pentru calculul integralei
Z b
I= f (x)dx
a

unde f : [a, b] → R este o funcţie continuă, Matlabf=@ oferă mai multe programe:

I = quad( f , a, b)
I = integral( f , a, b)

unde:

– I valuarea calculată a integralei;

– f funcţia Matlab de integrat;

59
60 CAPITOLUL 7. INTEGRARE NUMERICĂ

– a extremitatea stângă a intervalului de integrare;

– b extremitatea dreaptă a intervalului de integrare;

Exemplul 7.1 Să se calculeze integralele:

Z 1
1. 16x 15 dx
0
π
sin(x)
Z
2
2. dx
0 x

Rezolvările sunt:
1.

f=@(x)16*x.^15;
quad(f,0,1)

ans =

1.0000

2.

Definim funcţia Matlab exterioară

function y=fct(x)
if x==0
y=1;
else
y=sin(x)./x;
end

iar calculul este

quad(@fct,0,pi/2)

ans =

1.3708
7.1. INTEGRAREA FUNCŢIILOR DE O VARIABILĂ REALĂ 61

Utilizând integral, integralele exemplului anterior se calculează astfel:

1.

f=@(x)16*x.^15;
integral(f,0,1)

ans =

2.

integral(@fct,0,pi/2)

ans =

1.3708

• Dacă funcţia de integrat este cunoscută discret prin perechile (x i , y i )i , y i =

f (x i ) atunci integrarea se calculează cu

I = trapz(x, y),

unde x = (x i )i , y = (y i )i .

Exemplul 7.2 Simularea primei probleme din exemplul anterior.

f=@(x)16*x.^15;
x=0:10^(-6):1;
y=f(x);
trapz(x,y)

ans =

1.0000
62 CAPITOLUL 7. INTEGRARE NUMERICĂ

7.2 Calculul numeric al integralelor duble

Fie domeniul D definit prin

D = {(x, y) : a ≤ x ≤ b, yi n f ≤ y ≤ y sup(x)}. (7.1)

Integrala ZZ
f c t (x, y)dxdy (7.2)
D
se calculează în Matlab cu funcţia

I = integral2( f c t , a, b, yi n f , y sup)

unde parametrii corespund variabilelor introduse în (7.1) şi (7.2). yi n f , y sup pot
fi funcţii Matlab.
RR
Exemplul 7.3 Să se calculeze D x ydxdy unde domeniul D este delimitat de curbele
p
y = x 2 şi y = x.

Rezolvarea este

fct=@(x,y)x.*y;
yinf=@(x)x.^2;
ysup=@(x)sqrt(x);
integral2(fct,0,1,yinf,ysup)

ans =

0.0833
Capitolul 8

Integrarea numerică a ecuaţiilor

diferenţiale ordinare

Cap. 7 Cap. 9

8.1 Rezolvarea problemei cu valori iniţiale

Considerăm problema Cauchy

ẏ(x) = f (x, y(x)) x ∈ [x 0 , x f ], (8.1)

y(x 0 ) = y 0 , (8.2)

unde f : [x 0 , x f ] × R n → R n , despre care presupunem că admite o soluţie unică

y : [x 0 , x f ] → R n .
Matlab conţine mai multe funcţii pentru rezolvarea problemei cu valoare in-
iţială (8.1)-(8.2).
Funcţia ode45:
[x, y] = ode45(od e f un, t span, y 0 )
unde

• [x, y] este o matrice cu soluţia problemei;

63
64 CAPITOLUL 8. INTEGRAREA NUMERICĂ A E.D.O.

• od e f un este funcţia Matlab corespunzătoare membrului drept al ecuaţiei

(8.1);

• t span = [x 0 , x f ] intervalul de integrare;

• y 0 condiţia iniţială;

• y este
 matricea cu soluţiile calculate  în punctele vectorului x,
y 1 (x 1 ) y 1 (x 2 ) . . . y 1 (x m )
 y 2 (x 1 ) y 2 (x 2 ) . . . y 2 (x m ) 
 
y =
 ... ... ... ...


y n (x 1 ) y n (x 2 ) . . . y n (x m )

Exemplul 8.1 Să se rezolve problema Cauchy

ẏ(x) = y(x) x ∈ [0, 1]

y(0) = 1

Rezolvarea este

odefun=@(x,y)y;
[x,y]=ode45(odefun,[0,1],1);
plot(x,y)

Graficul soluţiei numerice este dat in Fig. 8.1.

Exemplul 8.2 Să se rezolve problema Cauchy

y˙1 = y 1 y 2 y 1 (1) = 1 x ∈ [1, 2]

y˙2 = − x1y 1 y 2 (1) = 1

Soluţia problemei Cauchy este y 1 (x) = x, y 2 (x) = x1 .

Rezolvarea numerică este
După definirea funcţiei

function deri=odefun(x,y)
deri=zeros(2,1);
deri(1)=y(1).*y(2);
deri(2)=-1/(x.*y(1));

rezolvarea este
8.2. REZOLVAREA PROBLEMELOR BILOCALE 65

Fig. 8.1: Graficul soluţiei problemei 8.1.

tspan=[1,2];
y0=[1;1];
[x,y]=ode45(@odefun,tspan,y0);
plot(x,y(:,1),x,y(:,2))
Graficul soluţiei numerice este dat in Fig. 8.2.

8.2 Rezolvarea problemelor bilocale

Fie problema bilocală

y 0 (x) = f (x, y(x)) x ∈ [a, b];

(8.3)
g (y(a), y(b)) = 0,

unde y : [a, b] → Rn , f : R × Rn → Rn şi g : Rn × Rn → R2n .

Se pot utiliza funcţiile Matlab:

sol = bv p4v(bv p f un, bc f un, sol i ni t )

sol = bv p4v(bv p f un, bc f un, sol i ni t )
66 CAPITOLUL 8. INTEGRAREA NUMERICĂ A E.D.O.

Fig. 8.2: Graficul soluţiei problemei 8.2.

unde

• bvpfun este funcţie Matlab definind funcţia f (x, y(x));

• bcfun este funcţie Matlab definind funcţia g (y(a), y(b));

• solinit este definit prin intermediul funcţiei bvpinit pe modelul

g=@(x)[expr_1,expr_2,\ldots]
x_retea=linspace(a,b,n)
solinit=bvpinit(x_retea,g)

expr i reprezintă aproximaţia iniţială pentru componenta y i (x).

Exemplul 8.3 Să se rezolve problema bilocală liniară

y 00 + y = 0,
y(0) = 0, y( π2 ) = 2.
8.2. REZOLVAREA PROBLEMELOR BILOCALE 67

Ecuaţia diferenţială este echivalentă cu sistemul diferenţial de ordinul 2

y 10 = y 2
y 20 = −y 1

Condiţiile la limită devin

y 1 (0) = 0
y 1 ( π2 ) = 2
Pentru această versiune funcţiile bvpfun şi bcfun sunt

function deri=bvpfun(x,y)
deri=zeros(2,1);
deri(1)=y(2);
deri(2)=-y(1);

function res=bcfun(ya,yb)
res=zeros(2,1);
res(1)=ya(1);
res(2)=yb(1)-2;
Rezolvarea este

g=@(x)[sin(x);cos(x)];
n=5;
x_retea=linspace(0,pi/2,n);
solinit=bvpinit(x_retea,g)

solinit =

solver: ’bvpinit’
x: [0 0.3927 0.7854 1.1781 1.5708]
y: [2x5 double]
yinit: @(x)[sin(x);cos(x)]

sol=bvp4c(@bvpfun,@bcfun,solinit)

solver: ’bvp4c’
x: [1x7 double]
y: [2x7 double]
yp: [2x7 double]
stats: [1x1 struct]
68 CAPITOLUL 8. INTEGRAREA NUMERICĂ A E.D.O.

[sol.x’,sol.y’]

ans =

0 0 2.0000
0.1963 0.3902 1.9616
0.3927 0.7654 1.8478
0.7854 1.4142 1.4142
1.1781 1.8477 0.7654
1.3744 1.9616 0.3902
1.5708 2.0000 0.0001

Pentru n = 30, graficul soluţiei numerice este în Fig. 8.3.

Fig. 8.3: Graficul soluţiei problemei 8.3.

Capitolul 9

Tranformata Fourier discretă

Cap. 8 ??

În Matlab transformata Fourier discretă şi transformata Fourier discretă in-

versă sunt calculate conform formulelor
n−1
y = Fn (x) x j w −k j
X
yk = k ∈ {0, . . . , n − 1}
j =0

şi respectiv
1 n−1
x = Fn−1 (y) yj wk j
X
xk = k ∈ {0, . . . , n − 1}
n j =0
2i π
unde x, y ∈ Cn iar w = e n . Prin Cn s-a notat mulţimea şirurilor de numere com-
plexe, periodice cu perioada n, Cn = {x = (x k )k∈Z : x k+n = x k , x k ∈ C, ∀k ∈ Z}.
În acest scop există funcţiile Matlab
• y=fft(x)
Calculul transformării discrete Fourier.

• x=ifft(y)
Calculul transformării discrete Fourier inverse.
unde semnificaţia parametrilor este x, y sunt şiruri de numere complexe reprezen-
tând originalul şi imaginea transformării.

69
70 CAPITOLUL 9. TRANFORMATA FOURIER DISCRETĂ

9.1 Calculul coeficienţilor Fourier

Exemplul 9.1 Să se calculeze coeficienţii Fourier ale funcţiei f (x) = sin x + cos3x.

Dacă integrala din expresia coeficientului Fourier complex

ak − i bk 1
Z 2π
ck = = f (x)e −i kx dx
2 2π 0

se calculează cu regula trapezelor, cu parametrul de discretizare n ∈ N∗ , atunci

se obţine aproximarea c = n1 Fn (ϕ), unde c = (c k )0≤k≤n−1 şi ϕ = (ϕ j )0≤ j ≤n−1 este
2π j
vectorul ale cărei componente sunt ϕ j = f ( n
).
Comenzile Scilab sunt

f=@(x)sin(x)+cos(3*x);
n=16;
h=2*pi/n;
x=0:h:2*pi-h;
ff=f(x);
c=fft(ff)/n;
a=c+conj(c);
b=(conj(c)-c)/i;
[a’,b’]

ans =

-0.0000 0
-0.0000 1.0000
-0.0000 0.0000
1.0000 0.0000
0.0000 0.0000
0.0000 0.0000
0.0000 0.0000
0.0000 -0.0000
0.0000 0
0.0000 0.0000
0.0000 -0.0000
0.0000 -0.0000
0.0000 -0.0000
1.0000 -0.0000
9.1. CALCULUL COEFICIENŢILOR FOURIER 71

-0.0000 -0.0000
-0.0000 -1.0000
Bibliografie

[1] BLAGA P., COMAN GH., POP S., TRÂMBIŢAŞ R., VASARU D., 1994, Anal-
iză numerică. Îndrumar de laborator. Univ. "Babeş–Bolyai" Cluj-Napoca
(litografiat).

[2] DINU M., LINCĂ Gh., 1999, Algoritmi şi teme speciale de analiză numerică.
Ed. Matrix rom, Bucureşti.

[3] GHINEA M., FIREŢEANU V., 1995, MATLAB calcul numeric, grafică, apli-
caţii. Rd. Teora, Bucureşti.

[4] MARINESCU GH., BADEA L., GRIGORE G., JAMBOR C., MAZILU P., RIZ-
ZOLI I., ŞTEFAN C., 1978, Probleme de analiză numerică. E.D.P., Bucureşti.

[5] MARINESCU GH., RIZZOLI I., POPESCU I., ŞTEFAN C., 1987, Probleme de
analiză numerică rezolvate cu calculatorul. Ed. Acad. R.S.R., Bucureşti.

[6] PLIS A. I., SLIVINA N.A., 1983, Laboratornyi praktikum po bysxei

matematike. Vysxa Xkola, Moskva.

[7] SCHEIBER E., LIXĂNDROIU D., CISMAŞIU C., 1982, Analiză numerică. În-
drumar de laborator. Univ. Braşov (litografiat).

[8] SCHEIBER E., SÂNGEORZAN L., GROVU M., 1993, Analiză numerică. În-
drumar de laborator. Univ. "Transilvania" Braşov (litografiat).

[9] SCHEIBER E., Laborator

de Analiză Numerică. Scilab. https:
//[Link], Tutorials.

72