Laboraratorul 1. INTRODUCERE ÎN MATLAB

Bibliografie:

1. M. Ghinea, V. Fireţeanu, Matlab: Calcul numeric- Grafică-Aplicaţii, ed. Teora,

Bucureşti, 1998.

2. I. Iatan , Îndrumător de laborator în Matlab 7.0, Ed. Conspress, Bucureşti, 2009.

Scopuri:

1) Aplicarea Matlab-ului in calcule matematice fundamentale.

2) Calcul numeric în Matlab cu aplicaţii în Algebră.

3) Rezolvarea in Matlab a ecuaţiilor si sistemelor de ecuatii diferenţiale.

4) Notiunea fisier in Matlab: creare, scriere, adaugare de date, citire.

5) Grafica in Matlab: reprezentari grafice 2D, reprezentarea grafica a histogramelor,

reprezentari grafice 3D.

MATLAB-ul este un pachet de programe de o performanţă remarcabilă care are o vastă

aplicabilitate atât în domeniul ştiinţei cât al ingineriei.

Pentru lansarea în execuţie a programului se acţionează dublu click pe pictograma

Matlab de pe Desktop sau se selectează StartAll Programs Matlab.

APLICAŢII ÎN CALCULE MATEMATICE FUNDAMENTALE

1) Să se calculeze expresiile următoare:

a) 33 25172517

>> (17+5*sqrt(2))^(1/3)+(17-5*sqrt(2))^(1/3)

ans =

5.0367

b)

33

1

1

1

1

i

i

i

i

>> ((1+i)/(1-i))^3-((1-i)/(1+i))^3

ans =

0 - 2.0000i

c) 53log xB , pentru 7x ;

Folosim formula de schimbare a bazei logaritmice

2

a

xx

b

b

alog

loglog .

>> x=7;

>> log(x^(1/5))/log(3)

ans =

0.3542

d) xxxxC 7cos5cos3coscos 4444 , pentru 8

x ;

>> x=pi/8;

>> C=cos(x)^4+cos(3*x)^4+cos(5*x)^4+cos(7*x)^4

C =

1.5000

2) Să se sorteze în ordine descrescătoare elementele vectorului

7820176.0 x , cu precizarea indicelui fiecărui element.

>> x=[-0.76 -1 20 8 -7]

x =

-0.7600 -1.0000 20.0000 8.0000 -7.0000

>> [y,I]=sort(x,’descend’)

y =

20.0000 8.0000 -0.7600 -1.0000 -7.0000

I =

3 4 1 2 5

3) Se consideră matricea

34.71278

78547.5

79.007.4

39971

A .

Se cere:

a) Transformaţi matricea A într-un vector coloană b ;

>> A=[1 -7 99 3;4.7 0 0.9 -7;5.7 4 5 78;-78 12 -7.4 3];

>> b=A(:);

3

b) Să se extragă submatricea D de dimensiune 23x , ce constă din elementele situate pe

ultimele trei linii şi primele două coloane ale matricei A .

>> D=A(2:4,1:2)

D =

4.7000 0

5.7000 4.0000

-78.0000 12.0000

CALCUL NUMERIC ÎN MATLAB CU APLICAŢII ÎN ALGEBRĂ

4) Rezolvaţi ecuaţiile neliniare:

a) arctgx1 x

>> f=@ (x) 1+x-atan(x);

>> x=fzero(f,[-3 3])

x =

-2.1323

>> f(x)

ans =

2.2204e-016

b) 0cos2 xx

>> f=@ (x) x.^2-cos(pi*x);

>> s=fsolve(f,[-0.5,0.5])

s =

-0.4384 0.4384

>> f(s)

ans =

1.0e-011 *

0.962 0.1001

5) Rezolvaţi sistemele neliniare:

a)

026

036

33

33

yyx

xyx (considerând ca punct iniţial 5.0,5.0 )

Pasul 1. Definim funcţia vectorială f .

4

>>f=@(x) [x(1)^3+x(2)^3-6*x(1)+3;

x(1)^3-x(2)^3-6*x(2)+2];

Pasul 2. Rezolvăm sistemul neliniar.

>> s = fsolve(f,[0.5 0.5])

s =

0.53237037226762 0.35125744755245

Pasul 3. Efectuăm verificarea soluţiei pe care am obţinut-o.

>> f(s)

ans =

1.0e-009 *

0.29623858921468

0.19358559200100

b)

0

1

9

2

222

zyx

xyz

zyx

(considerând ca punct iniţial 6.1,2.0,5.2 ).

>> f=@(x) [x(1)^2+x(2)^2+x(3)^2-9;x(1)*x(2)*x(3)-1;x(1)+x(2)-x(3)^2]

f =

@(x) [x(1)^2+x(2)^2+x(3)^2-9;x(1)*x(2)*x(3)-1;x(1)+x(2)-x(3)^2]

>> s = fsolve(f,[2.5 0.2 1.6])

s =

2.49137569683072 0.24274587875742 1.65351793930053

>> f(s)

ans =

1.0e-011 *

0.11226575225010

0.13493650641294

-0.05244693568329

REZOLVAREA IN MATLAB A ECUAŢIILOR SI SISTEMELOR DE ECUATII DIFERENŢIALE

6) Să se integreze ecuaţia de tip Riccati

xxyxyyx 212 22 , 0x

>> y=dsolve('x*Dy=y^2-(2*x+1)*y+x^2+2*x','x')

5

y =

(-x-1+x^2*C1)/(-1+x*C1)

7) Să se integereze ecuaţia diferenţiala Euler:

xyyx 1

>> y=dsolve('x*D3y+D2y=1+x','x')

y =

1/12*x^3+x*log(x)*C1-C1*x+1/2*x^2+C2*x+C3

8) Să se rezolve următorul sistem de ecuaţii diferenţiale liniare neomogen:

33

23

213

2312

2321

xyyy

xyyy

xxyyy

>> [y1,y2,y3]=dsolve('Dy1=y2+y3-x-x^2','Dy2=3*y1+y3-2-x^2','Dy3=3*y1+y2+x-3','x')

y1 =

1+2/3*C2*exp(3*x)-C3*exp(-2*x)

y2 =

x+exp(-x)*C1+C2*exp(3*x)+C3*exp(-2*x)

y3 =

-exp(-x)*C1+C2*exp(3*x)+C3*exp(-2*x)+x^2

NOTIUNEA FISIER IN MATLAB

9) Construiti un program in Matlab, cu ajutorul caruia:

- sa se creeze un fisier in care vor fi salvate doua matrice, una de tipul ( ) si

cealalta ( ) ;

- sa se citeasca continutul fisierului in alte doua matrice.

Se va construi un script in Matlab cu numele fisier.m, avand continutul:

W=[1 6; 7 8]; w=[-2 3 5; 6 7 8];

fid=fopen('pond1fb','w');

fprintf(fid,'%f\t',W);

fprintf(fid,'%f\t',w);

fclose(fid);

fid=fopen('pond1fb','r');

6

a=fscanf(fid,'%f\t',[2,3]);

b=fscanf(fid,'%f\t',[2,2]);

fclose(fid);

In linia de comanda se scrie

>>fisier

GRAFICA IN MATLAB

10) Reprezentati grafic in Matlab printr-o histograma, notele obtinuite de cei 50 de studenti

inscrisi la cursul de Statistica. Notele se vor citi dintr-un fisier text.

Se construieste un script in Matlab, cu denumirea fis.m, ce contine:

fid=fopen('pond.txt','r');

u=fscanf(fid,'%f\t',[1,50]);

fclose(fid);

In linia de comanda scriem:

>> fis

>> for j=1:max(u)

v=find(u(1:length(u))==j);

n(j)=length(v);

end

>> sum(n)

ans =

50

>> v=1:10;

>> hist(u,v)

>> grid on

7

11) Reprezentaţi grafic următoarele funcţii în plan:

21

2arcsin)

x

xxfa

, 5,5x

Pasul 1. Se fixează intervalul pe care va fi reprezentată funcţia şi un anumit pas.

>>x=-5:0.1:5;

Pasul 2. Definim funcţia ce urmează să fie reprezentată.

>> f=@(x) asin(2*x./(1+x.^2));

Pasul 3. Realizăm reprezentarea grafică.

>> plot(f(x),'m','LineWidth',4)

0 20 40 60 80 100 120

-2

-1.5

-1

-0.5

0

0.5

1

1.5

2

8

xxfb arcsin) , xxg arccos , 1,1x

>>x=-1:0.01:1;

>> f=@(x) asin(x);

>> g=@(x) acos(x);

>> plot(x,f(x),'r',x,g(x),'b','LineWidth',3)

12) Reprezentati arcul de parabola 2: xyAB , care uneste punctele 1,1A si 4,2B .

Secventa Matlab urmatoare permite reprezentarea arcului AB .

>> x=-4:.1:4;

>> y=x.^2;

>> plot(x,y,1,1,'or',2,4,'or')

13) Scrieţi un fişier “function” în Matlab pentru a reprezenta grafic funcţia

-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1-2

-1

0

1

2

3

4

-4 -3 -2 -1 0 1 2 3 40

2

4

6

8

10

12

14

16

A(1,1)A(1,1)A(1,1)A(1,1)

B(2,4)

9

0,0

0,1

cos

x

xx

xxf , 5.0,5.0x , 01.0h

Se selectează succesiv File->New->M-file şi se scriu următoarele instructiuni

function r=f(x)

if x~=0

r=x*cos(1/x);

elseif x==0

r=0;

end

end

Se salvează fişierul cu f.m apoi în linia de comanda se scrie:

>> x=-0.5:0.01:0.5;

>> for k=1:length(x)

y(k)=f(x(k));

end

>> plot(x,y)

14) Reprezentati grafic in Matlab corpul solid, de forma tetraedrului din primul octant,

marginit de planele 1 zyx , 0x , 0y , 0z .

>> x=[1 0 0 0 1 0]; y=[0 0 1 0 0 1]; z=[0 0 0 1 0 0];

>> plot3(x,y,z,1,0,0,'ob',0,1,0,'ob',0,0,1,'ob',0,0,0,'ob')

-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5-0.4

-0.3

-0.2

-0.1

0

0.1

0.2

0.3

0.4

10

15) Realizati in Matlab graficul corpului, limitat de suprafetele: xzy 42 22 si 2x .

Secventa Matlab urmatoare permite reprezentarea corpului omogen.

>> [y,z]=meshgrid(-3:.03:3,-2.5:.03:2.5);

>> x=y.^2/4+z.^2/2;

>> [m,n]=size(x);

>> x1=2*ones(m,n);

>> plot3(x,y,z,x1,y,z)

00.2

0.40.6

0.81

0

0.2

0.4

0.6

0.8

1

0

0.2

0.4

0.6

0.8

1

xy

z

A(1,0,0)

O(0,0,0)

C(0,0,1)

B(0,1,0)

01

23

45

6

-5

0

5

-3

-2

-1

0

1

2

3

xy

z

Laboraratorul 2. Aplicatii ale testelor Chauvenet si Young

Bibliografie:

1. Văduva, I. Modele de simulare, Editura Universitatii din Bucureşti, 2004.

2. Trandafir, R. Note de curs, Facultatea de Hidrotehnica, an III, an AIA, 2011-2012.

3. M. Marusteri, Notiuni fundamentale de biostatisticã. Note de curs. University Press,

Târgu Mureş, 2006, http://www.umftgm.ro/info/Curs_Notiuni_fundamentale.pdf.

4. M. Ghinea, V. Fireţeanu, Matlab: Calcul numeric- Grafică-Aplicaţii, ed. Teora,

Bucureşti, 1998.

5. I. Tutunea, Algoritmi, logica si programare¸ Reprografia Universitatii Craiova, 1993.

6. I. Iatan , Îndrumător de laborator în Matlab 7.0, Ed. Conspress, Bucureşti, 2009.

Scopuri:

1) Aplicarea testului Chauvenet privind eliminarea datelor afectate de erori aberante;

2) Aplicarea testului Young pentru identificarea erorilor sistematice.

TESTUL CHAUVENET

Fiind date valorile, nxx ,......,1 rezultate din măsurători, se consideră că valoarea ix este

afectată de erori aberante dacă verifică relaţia:

zxxi || , (1)

unde:

x = n

xn

jj

1 este media aritmetică a valorilor observate,

n

ii xx

n 1

21 semnifica abaterea standard a valorilor observate,

z se ia din tabele (în functie de numarul n de valori din sir, ce reprezinta

dimensiunea sirului sau volumul esantionului) sau se aproximează conform relatiei

următoare: 2213,3604,31

862,0435,0

aa

az

, cu

n

na

4

12 .

Observatie. Este suficient ca verificarea relatiei (1) sa fie efectuata doar pentru

valorile extreme (minima si maxima) din cadrul esantionului.

Daca, în urma aplicarii testului, rezulta ca una dintre valorile testate este afectata de

erori aberante, valoarea respectiva este eliminata din cadrul esantionului, se recalculeaza

valorile mediei si abaterii standard pentru valorile ramase si se reia verificarea conditiei (1),

2

algoritmul aplicându-se pâna când conditia respectiva nu mai este verificata pentru nici una

dintre cele doua valori extreme ale esantionului.

Exemplul 1. În simularea utilizării unei imprimante conectate la o reţea se urmăreşte

repartiţia numărului maxim de fişiere care sunt in lista de aşteptare pentru listare pe o

perioada de 30 de zile, înregistrându-se valorile:

săpt L Ma Mi J V săpt L Ma Mi J V

I 15 19 13 18 20 IV 19 21 13 23 17

II 12 21 22 19 21 V 18 20 21 13 15

III 13 13 13 14 17 VI 20 19 12 18 12

Să se verifice dacă acest eşantion are valori aberante şi dacă există, să se elimine! Scrieţi un

program Matlab pentru testul Chauvenet. Datele de intrare se iau dintr-un fişier text.

Componentele

esantionului

(ordonate

crescator)

Valori

extreme

|| xxi z

1 15 12 12 5.0333 8.3158

2 19 12 12 5.0333

3 13 12 12 5.0333

4 18 13 23 5.9667

5 20 13

6 12 13

7 21 13

8 22 13

9 19 13

10 21 14

11 13 15

12 13 15

13 13 17

14 14 17

15 17 18

16 19 18

17 21 18

18 13 19

19 23 19

20 17 19

21 18 19

22 20 20

23 21 20

24 13 20

25 15 21

26 20 21

3

3

27 19 21

28 12 21

29 18 22

30 12 23

0333.17x , 5183.3

In concluzie, nici una dintre valorile extreme nu este afectata de erori aberante.

Exemplul 2. A fost măsurată greutatea a 15 indivizi adulţi. Rezultatele măsurătorilor

sunt cele din tabelul următor. Să se verifice dacă acest eşantion are valori aberante şi dacă

există, să se elimine! Scrieţi un program Matlab pentru testul Chauvenet. Datele de intrare se

iau dintr-un fişier text.

Prima aplicare a criteriului Chauvenet

Componentele

esantionului

(ordonate

crescator)

Valori

extreme

|| xxi z

1 58 35 35 38.9333 56.53065

2 60 50 160 86.06667

3 80 55

4 77 58

5 83 60

6 75 65

7 82 70

8 79 75

9 50 77

10 35 79

11 70 80

12 160 80

13 80 82

14 65 83

15 55 160

9333.73x , 652.26

După cum rezultă din tabel, valoarea 160 va trebui să fie eliminată din datele supuse

prelucrării.

A doua aplicare a criteriului Chauvenet

Componentele

esantionului

(ordonate

crescator)

Valori

extreme

|| xxi z

1 35 35 32.7857 29.19275

2 50 83 15.2143

3 55

4 58

4

5 60

6 65

7 70

8 75

9 77

10 79

11 80

12 80

13 82

14 83

7857.67x , 93443.13

De data aceasta trebuie eliminata din tabel valoarea 35.

A treia aplicare a criteriului Chauvenet

Componentele

esantionului

(ordonate

crescator)

Valori

extreme

|| xxi z

1 50 50 20.30769 22.64275

2 55 83 12.69231

3 58

4 60

5 65

6 70

7 75

8 77

9 79

10 80

11 80

12 82

13 83

30769.70x , 95715.10

Utilizand Matlab 7.9, vom construi programul corespunzator testului Chauvenet:

Etapa 1. Din fisierul dateb.txt se vor citi datele initiale.

Etapa 2. Se scrie functia sterge.m, ce permite stergerea unei componente a unui vector.

function v=sterge(u,k,n) for i=1:k-1 v(i)=u(i); end if n~=0 for j=k:n-1

5

5 v(j)=u(j+1); end end end

Etapa 3. Se scrie programul propriu- zis.

fid=fopen('dateb.txt','r'); u=fscanf(fid,'%f\t',[1,15]); fclose(fid); m=mean(u); sigma=std(u); n=length(u); a=(2*n-1)/(4*n) z=(0.435-0.862*a)/(1-3.604*a+3.213*a^2); it=1; u=sort(u); mi=min(u); ma=max(u); p1=length(find(u==mi)); p2=length(find(u==ma)); for i=1:p1 w(i)=mi; end for i=p1+1:p1+p2 w(i)=ma; end q=find(abs(w-m)>z*sigma); while (length(q)~=0) it=it+1; for k=1:length(q) h=find(u==w(q(k))); for l=1:length(h) v=sterge(u,h(l),n); n=n-1; u=v; end end m=mean(u); sigma=std(u); n=length(u); a=(2*n-1)/(4*n); z=(0.435-0.862*a)/(1-3.604*a+3.213*a^2); u=sort(u); mi=min(u); ma=max(u); p1=length(find(u==mi)); p2=length(find(u==ma)); for i=1:p1 w(i)=mi; end for i=p1+1:p1+p2 w(i)=ma; end q=find(abs(w-m)>z*sigma); end

6

TESTUL LUI YOUNG

Problema depistarii si eliminarii erorilor sistematice este mai dificila datorita

multitudinii de factori care se interconditioneaza. Vom prezenta testul lui Young, test care nu

ofera posibilitatea eliminarii erorilor sistematice, ci doar pe aceea a aprecierii influentei

cauzelor sistematice asupra datelor de sondaj.

Pasul 1. Intrare: nxx ,......,1 - şir de date experimentale şi probabilitatea de acceptare

(coeficient de încredere); se calculează mărimile:

.

1

1

2

2

21

11

2

M

xxn

n

iii

Pasul 2. Se compară M cu valoarea critică (VCI) inferioară şi valoarea critică superioară

(VCS), VCI<M<VCS, valori luate din tabele în funcţie de n şi sau determinate cu relaţiile:

a) Dacă 95,0 atunci

341,0919,8

2

057,1317,3

003,0081,0491,0

neVCS

nnVCI

b) Dacă 99,0 atunci

388,1574,33

33,2

33,2

882,0484,3

427,411

269,1883,192

neVCS

n

nVCI

Dacă inegalitatea este satisfăcută, atunci se consideră că şirul de date experimentale are un

caracter aleator (nu este afectat de erori sistematice) cu probabilitatea .

Pasul 3. Stop!

Observatie. În tabele se dau valori pentru corespunzătoare diferitelor valori ale lui n

)25( n

Exemplul 3. Testati daca esantionul de date experimentale, ce consta din valorile ix

ale timpului de latenţă a instalării efectului hipnotic în cazul amobarbitalului (vezi tabelul

urmator) are sau nu un caracter aleator cu probabilitatea 95.0 :

7

7

Componentele esantionului (secunde)

1 16.1

2 15.5

3 13.4

4 22.8

5 12.1

6 11.3

7 11.6

8 6.3

9 8.8

10 7.1 function [VCI, VCS]=valcrit(n,al) if al==0.95 VCI=0.491+0.081*n-0.003*n^2; VCS=3.317-1.057*exp(-8.919*n^(-0.341)); elseif al==0.99 VCI=(192.883+1.269*n^2.33)/(411.427+n^2.33); VCS=3.484-0.882*exp(-33.574*n^(-1.388)); end end

fid=fopen('data.txt','r'); x=fscanf(fid,'%f\t',[1,10]); fclose(fid) n=10;al=0.95 sigma=var(x) delta=0; for i=1:n-1 delta=delta+(x(i+1)-x(i))^2; end delta=delta/(n-1); M=delta/sigma [VCI,VCS]=valcrit(10,al); if VCI<M<VCS disp(['esantionul datelor experimentale are un caracter aleator cu

probabilitatea:', num2str(al)]); else disp('esantionul datelor experimentale este afectat de erori

sistematice') end