riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per...

25
1 ,0DWODE±1R]LRQLJHQHUDOL ±&RV¶q0DWODE Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio di programmazione ad alto livello, sulla stessa linea del Fortran, ma con particolari facilitazioni nelle elaborazioni di matrici. - Possibilità di costruire grafici in 2 e 3 dimensioni - Ampia disponibilità di programmi interni per la risoluzione dei problemi dell’Analisi Numerica - Disponibilità di pacchetti per svariati tipi di applicazioni (Toolbox). Tra queste applicazioni ricordiamo l’elaborazione numerica dei segnali, la simulazione di sistemi dinamici, il calcolo simbolico ecc. Matlab lavora in modo interattivo, cioè l’utente digita una istruzione ed ha immediatamente la risposta. Il prompt su cui si digita l’istruzione è la coppia di caratteri >> : >> (comando Matlab. Per eseguire, digitare Enter) Per uscire dalla sessione di lavoro interattiva: >> quit E’ anche possibile lavorare a programma. In questo caso l’utente scrive l’intero programma mediante un editor, quindi lo registra su un file , e infine lo manda in esecuzione. $ULWPHWLFDHIXQ]LRQLHOHPHQWDUL. Operatori aritmetici . + - / * ^ Funzioni elementari . sin cos asin acos tan atan, atan2 sinh asinh cosh ….. exp log log10 sqrt abs conj imag real fix (parte intera), floor, (intero inferiore dell’argomento), ceil (intero superiore), round (intero più vicino) rem (resto della divisione tra i 2 argomenti), sign, ecc. I dati numerici sono rappresentati in memoria nella forma a virgola mobile (floating point), su parole di 64 bit (16 decimali) per i PC. Tuttavia il formato esterno può essere deciso dall’utente con i seguenti comandi: >> format short % 5-6 cifre (opzione di default) - es: 3.1416 >> format short e % 5-6 cifre, forma esponenziale - es: 0.25 e-4 >>format long % 16 cifre >> format long e % 16 cifre, forma esponenz. ….. (altri formati) Costanti predefinite . pi = π, realmax, realmin: numeri massimo (2 1023 ) e minimo (2 -1022 ) rappresentabili in floating-point , i (oppure j): unità immaginaria (ad esempio 1+2i è un mumero complesso), Inf=Esempio: >> exp(4.27)+cos(pi/4) % questo è un comando Matlab 72.2287 % questo è l’output corrispondente Il valore di una espressione può essere assegnato ad una variabile: >> a= exp(4.27)+cos(pi/4) a= 72.2287 >> a-2 70.2287 Negli identificatori delle variabili lettere maiuscole e minuscole sono considerate diverse. Così c1 e C1 sono due variabili GLVWLQWH. Se un’istruzione termina con il carattere ; l’output viene soppresso. Sulla stessa riga si possono digitare più istruzioni separate da virgola (o da ; se si vuole eliminare l’output), ad esempio: >> ceil(1.27), 1.27; exp(1.27) 2 3.5609

Transcript of riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per...

Page 1: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

1

�������������� ����������������

������������������Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio di programmazione ad alto livello, sulla stessa linea del Fortran, ma con particolari facilitazioni nelle elaborazioni di matrici. - Possibilità di costruire grafici in 2 e 3 dimensioni - Ampia disponibilità di programmi interni per la risoluzione dei problemi dell’Analisi Numerica - Disponibilità di pacchetti per svariati tipi di applicazioni (Toolbox). Tra queste applicazioni

ricordiamo l’elaborazione numerica dei segnali, la simulazione di sistemi dinamici, il calcolo simbolico ecc.

Matlab lavora in modo interattivo, cioè l’utente digita una istruzione ed ha immediatamente la risposta. Il prompt su cui si digita l’istruzione è la coppia di caratteri >> : >> (comando Matlab. Per eseguire, digitare Enter) Per uscire dalla sessione di lavoro interattiva: >> quit E’ anche possibile lavorare a programma. In questo caso l’utente scrive l’intero programma mediante un editor, quindi lo registra su un file , e infine lo manda in esecuzione. ��������������������� ���������������. Operatori aritmetici. + - / * ^ Funzioni elementari. sin cos asin acos tan atan, atan2 sinh asinh cosh ….. exp log log10 sqrt abs conj imag real fix (parte intera), floor, (intero inferiore dell’argomento), ceil (intero superiore), round (intero più vicino) rem (resto della divisione tra i 2 argomenti), sign, ecc. I dati numerici sono rappresentati in memoria nella forma a virgola mobile (floating point), su parole di 64 bit (≅ 16 decimali) per i PC. Tuttavia il formato esterno può essere deciso dall’utente con i seguenti comandi: >> format short % 5-6 cifre (opzione di default) - es: 3.1416 >> format short e % 5-6 cifre, forma esponenziale - es: 0.25 e-4 >>format long % 16 cifre >> format long e % 16 cifre, forma esponenz. ….. (altri formati) Costanti predefinite. pi = π, realmax, realmin: numeri massimo (21023) e minimo (2-1022) rappresentabili in floating-point , i (oppure j): unità immaginaria (ad esempio 1+2i è un mumero complesso), Inf=∞ Esempio: >> exp(4.27)+cos(pi/4) % questo è un comando Matlab 72.2287 % questo è l’output corrispondente Il valore di una espressione può essere assegnato ad una variabile: >> a= exp(4.27)+cos(pi/4) a= 72.2287 >> a-2 70.2287 Negli identificatori delle variabili lettere maiuscole e minuscole sono considerate diverse. Così c1 e C1 sono due variabili ��������. Se un’istruzione termina con il carattere ; l’output viene soppresso. Sulla stessa riga si possono digitare più istruzioni separate da virgola (o da ; se si vuole eliminare l’output), ad esempio: >> ceil(1.27), 1.27; exp(1.27) 2 3.5609

Page 2: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

2

Se una istruzione è troppo lunga per essere contenuta in una riga, prima di andare a capo digitare 3 punti (…)

����������������������� ����!"

�����#����� �������$��$����%����������Matlab (= Matrix Laboratory) è un linguaggio particolarmente potente nella elaborazione di vettori e matrici, che si scrivono in conformità agli esempi seguenti: Vettore riga: [1 2 3 0.3] (oppure [1,2,3,0.3]) Vettore colonna: : [1; 2; 3; 0.3] Matrice: [1 2 5; 0 3 9; 0.5 4 0] Come separatore tra 2 elementi di un array si può anche usare la virgola (invece dello spazio). Se una matrice (o qualsiasi istruzione) è troppo lunga, prima di andare a capo 3 punti (…) Se un elemento di una matrice è una espressione, non lasciare spazi dentro l’elemento. Matrici e vettori non vengono dimensionati. Ciò permette il ridimensionamento in corso di lavoro, come è mostrato nell’esempio: >> A=[1 2; 3 4] A= 1 2

3 4 >> A(3,3)=5 A= 1 2 0 3 4 0 0 0 5 Come si vede, agli elementi non definiti è stato assegnato il valore 0. Indici: Debbono essere >0. A(i,j) , v(i) indicano un elemento della matrice A o del vettore v. Trasposta. Con A’ si indica la trasposta di una matrice A (o vettore) Gli operatori +, -, * indicano le corrispondenti operazioni tra array, se esse sono possibili. Esempio: se u=[1 2 3], v=[3 2 0] il prodotto u*v non è possibile, ma u*v’=7. Il prodotto di una matrice (o vettore) per uno scalare è sempre possibile, ad es. 3*v=[9 6 0]. L’operazioni a+c (a array, c numero) somma ad ogni elemento di a il numero c. Analogam. Per a-c. A^n (A matrice quadrata, n intero >0) moltiplica A n volte per se stessa. Operatori di divisione. Vi sono 2 operatori di divisione , / (slash) e \ (backslash) ������: Se A(nxn) è quadrata e regolare e b è un vettore, x=A\b calcola la soluzione del sistema lineare Ax=b. Se invece B è una matrice nxn avente come colonne b1,b2,…bn , allora X=A\B sarà a sua volta una matrice, avente come colonne x1,x2,…xn le soluzioni dei sistemi lineari Ax1=b1, Ax2=b2, … Axn=bn (si dice anche che X è soluzione dell’equazione matriciale AX=B). Se A(mxn) con m>n, x=A\b calcola la soluzione di Ax=b nel senso dei minimi quadrati. Se m<n vi saranno infinite soluzioni, ed x=A\B ne fornisce una base (in questo caso x è una matrice). ����: X=B/A indica la matrice tale che XA=B (poco usato). Operazioni elemento per elemento (gli array A e B debbono avere le stesse dimensioni): A .*B , A./B Prodotto e divisione elemento per elemento A .^B ciascun elemento di A viene elevato a potenza, con esponente il corrispond. di B

Esempio:

==

=

=

2215

107*2^,

169

412.^,

43

21����� ,

In ogni caso, i risultati sono disposti in un array che ha le stesse dimensioni di A (e di B).

Page 3: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

3

c . / B (c costante, B vettore o matrice) Risultato: [c/b1, c/b2,….] Tutte le funzioni definite in 1.2 possono essere applicate ad una matrice (o vettore), in tal caso la funzione si intende applicata ����������������. Operatore : A(m:k, i:j) separa la sottomatrice di A con le righe da m a k, e le colonne da i a j. A(:, i:j) sottomatrice delle colonne da i a j. A(m:k, :) sottomatrice delle righe da m a k. Ad esempio, nel caso della matrice A=[1 2 0; 3 4 1; 2 6 3]: >> b=A(2:3, :) b= 3 4 1 2 6 3 Aggiunta alla matrice A di una colonna: A=[A u] aggiunge ad A il vettore u come ultima colonna ( u deve essere vettore colonna, se è vettore riga basterà fare A=[A u’] ) Aggiunta di una riga: A=[A; v] aggiunge ad A il vettore riga v come ultima riga Esempio: >> B=[1 2 3] >> for k=2:3 , B=[B; k k+1 k+2]; end >> B

B=

543

432

321

Azzeramentodi righe o colonne: A(i,:) =[] azzera la i-ma riga di A. A(i:j, :)=0 azzera le righe da i a j A(: , m:k]=[] azzera tutte le colonne da m a k ������&�� ����������'���� �"�Se A e’ un vettore calcola la somma degli elementi. Se A e’ una matrice calcola il vettore delle somme di A per colonne. $��( �" Se A e’ un vettore calcola il prodotto degli elementi. Se A e’ una matrice calcola il vettore dei prodotti di A per colonne. )�*�+,�� � �" Assegna ad m ed n i valori delle dimensioni di A. #,(�� �"� Calcola il determinante -,��. �"� Calcola l’inversa di A �,��� �" Vettore degli autovalori di A )�*#+,��� �" V: matrice degli autovettori, D matrice diagonale con gli autovalori �, ���� �*�". Inizializzazione a 0 di una matrice A m x n: #,�!� �*�"��costruisce D con tutti 1 in posizione diagonale (se D è quadrata, D=eye(n) è la matrice

identica) #,(��� ." costruisce una matrice diagonale con il vettore v sulla diagonale principale. /,���( �"� matrice quadrata ad elementi numeri casuali compresi tra 0 ed 1 /,����( �*�" Matrice random di dimens. m x n Matrici speciali: 0�(����(, 0��, ��.0��, ���$�� ,…. Ecc. Esempio. L’equazione matriciale AX=I (matrice identica) ha evidentem. la soluzione X=A-1. Calcolare l’inversa di una matrice sia con l’operatore \ sia con la funzione inv. >> A=[1 2 –3;-2 2 0;4 –2 –1] A= 1 2 -3 -2 2 0 4 -2 -1 >> A\eye(3) -0.3333 1.3333 1.0000 -0.3333 1.8333 1.0000 -0.6667 1.6667 1.0000

Page 4: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

4

>> inv(A) -0.3333 1.3333 1.0000 -0.3333 1.8333 1.0000 -0.6667 1.6667 1.0000 ��1������� �����(��.���������.,����'�����'���2�����������oppure .,���'��2��Costruisce il vettore riga min, min+incr, min+2*incr,…max. Se l’incremento manca, si intende uguale ad 1. Esempio: a=0:10:100 produce il vettore riga a=[0 10 20 30 …..100] ����$��� ���*���2*��$����" crea un vettore di n punti reali equidistanziati tra min e max Esempi: a=linspace(0, 2, 41) produce il vettore riga a=[0 0.05 0.1,….1.95 2] sum(1:1:10) risultato:55 sum((1:1:10).^2) risultato: 385 (somma i quadrati dei primi 10 numeri interi) ��3��&�� �����(�����������Se A è matrice quadrata nxn ed f(x) una funzione analitica, nell’algebra della matrici si può definire la matrice quadrata nxn f(A) (questa teoria non fa parte del corso). In matlab: �4��� �" Radice di A, cioè una matrice B tale che B*B=A (Attenz: diverso da sqrt(A), i cui elementi sono le radici degli elementi di A) �2$� �" eA. ���� �" log(A) ….. ��5��6�����0��(������������Sono sequenze di caratteri racchiuse tra apici. Es: nome = ’Matteo’. In questo caso nome è un vettore di 6 caratteri. Si possono considerare anche matrici: nomi=[ ‘Mario’; ‘Maria’ ; ‘Luigi’] però tutti devono avere lo stesso numero di caratteri (se non ce l’hanno, riempire con blank). Funzioni: �0�� ��������*���������*7�" con stringa1, stringa2,… anche di lunghezze diverse tra loro. Forma la matrice, aggiungendo automaticamente (se necessario) i blank. �.�� �������" Se stringa è un valido comando Matlab, lo esegue. Es: >> eval(‘x=cos(pi/12)’) % Calcola l’espressione, e ne assegna il valore ad x. La funzione eval è utile quando in una sessione di lavoro si deve calcolare più volte una data espressione. Esempio: >> f=’cos(x)+cos(2*x)+cos(3*x)’ >> x=0.1; a=2+eval(f) >> x=0.2; b=3+eval(f) ��8��&�� �����(���������(������$���Vi sono numerose funzioni, che danno il tempo in diversi formati. Ad esempio: ���� (… calcolo..) t=��� La variabile t contiene il tempo di CPU (in secondi) dal tic, cioè il tempo impiegato per il calcolo. (���. Data del giorno, ecc.

Page 5: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

5

1��9�������� ���������������Oltre che in modo interattivo, Matlab consente di lavorare con un programma memorizzato, che viene scritto mediante un editor incorporato. I programmi Matlab sono detti Script Files, e consentono anche la definizione di funzioni. Sia gli Script Files che le funzioni debbono essere registrati con l’estensione .m, e vengono detti pertanto m-files. Premettiamo la definizione degli operatori relazionali e logici, e delle strutture di controllo. 1���:$������������ ����������������I 2 valori logici (costanti logiche) sono: 0 Falso, 1 Vero (ma ogni intero ≠0 è considerato vero) Operatori Relazionali: < minore, <= minore o uguale, > maggiore, >= maggiore o uguale, == uguale, ~= diverso A primo e secondo membro di una operazione relazionale debbono esserci 2 valori numerici. Esempio: >> t=2>3 % t=0 (falso) Operatori Logici: & AND , | OR , ~ NOT a primo e secondo membro di una operazione logica debbono esserci 2 valori logici. Più in generale gli operatori relazionali e logici si possono applicare alle matrici, ed il risultato sarà una matrice a valori logici avente le stesse dimensioni di quelle poste a confronto. Esempio: se x=[0 5 3 7] ed y=[0 2 8 7] >> m=(x>y)&(x>4) % risultato: m=[0 1 0 0] >> p= x|y % risultato: p=[0 1 1 1] (attenz: qualsiasi intero ≠0 è considerato vero) Altre funzioni logiche: ��! 2"���fornisce 1 (vero) se ogni elemento dell’array x è non-zero (…..altre funzioni) 1���6���������(���������������������$�������� ����'�Cicli For: for n=j:k for n=j:m:k Blocco di istruzioni blocco di istruzioni end end Nella prima forma l’iteratore n va da j a k con passo 1. Nella seconda da j a k con passo m. Cicli while: while espress. logica Blocco di istruzioni end Il ciclo viene eseguito finchè l’espress. logica rimane vera (=1). Selezioni if – elseif - else: if espress. logica 1° blocco di istruz. elseif espress. logica 2° blocco di istruz. else 3° blocco di istruz. end come in Fortran: la elseif e/o la else possono mancare. In ogni blocco di istruzioni possono essere annidati altri blocchi (cicli o else-if). Se l’espressione logica è complicata, racchiuderla tra parentesi Ad es. (x>0 & n=100) Le for, while, if possono essere usate anche in modo interattivo nella finestra di lavoro principale (command window). In tal caso l’intera struttura va in esecuzione solo quando è digitata la end. break provoca l’uscita da un ciclo for o da un ciclo while.

Page 6: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

6

Selezioni Switch-case-Otherwise: switch flag case valore1 1° blocco di istruz. case valore2 2° blocco di istruzioni ………. otherwise ultimo blocco di istruzioni end Flag è una qualsiasi variabile. Se il suo attuale valore è valore1 viene eseguito il 1° blocco. Se valore2 il 2° blocco, e così via. 1�1��;(�������(������ �����(���$��������� 6���$�������"�Uno script file è un file di comandi matlab , che può essere chiamato ed eseguito in qualsiasi punto del Command Window. Gli script files sono memorizzati nella Directory di Lavoro (Current Directory), che appare nel piccolo riquadro in alto a destra della finestra Matlab, e per default è c:\matlab6p1\work. Può però essere cambiata (ad es. in A:\) agendo sul pulsante (…) a destra. Per crearlo File → New ��-File (appare uno schermo editor)→scriverlo: % commento su cosa fa x=…. ; ( ; se non si vuole visualizzare il risultato nella esecuzione) z=…. (istruzioni) Per salvarlo: File → Save as (finestra di dialogo. Digitare il nome con estensione .m) ����� Per eseguirlo >>nome (in qualsiasi punto della finestra principale) – Viene subito eseguito. Può essere richiamato per correggerlo dal Command Window con <<��(��������� . Se il file nome.m non esiste, >> edit nome.m consente anche di crearlo. Comandi utili in uno Script-File: (��$ .�����" Visualizza valore sullo schermo. Valore puo’ essere una variabile (anche un array), od una ‘stringa di caratteri’. Se si vogliono visualizzare più elementi, porli in un array. �,��$�� =��������" L’elaborazione si arresta e viene visualizzata la stringa . Riparte dopo l’introduzione da tastiera di un valore, che viene assegnato ad a. Nel caso si voglia introdurre un array occorre digitare anche le parentesi quadre, ad es. [1 2; 4 1] $����. Arresto finche’ viene premuto un tasto. $���� �": arresto per n secondi. >�!���(�� passa il controllo al Command Window. Per mettere in rilievo la particolare situazione, il prompt si presenterà così: k>> ><<������� ritorna allo Script File, dove era uscito con keyboard. Le variabili di uno Script File sono globali, cioè condivise col Command Window (e viceversa) Dunque una sessione di lavoro può essere interattiva, a programma o mista. Una sessione di lavoro interattiva (o parte di essa) può essere conservata mediante il comando (���!�. Supponiamo che debba essere conservata nel file miofile.m. Allora: >>diary on >> diary miofile.m ( sessione di lavoro che si vuole conservare) >> diary off (parte che non si vuole conservare) >> diary on (si vuole conservare di nuovo, sempre in miofile.m) >>diary off

Page 7: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

7

Per rieseguirlo bisogna cancellare tutti gli output generati durante l’esecuzione precedente, e gli eventuali messaggi di errore. 1�3������������� �����Una funzione è una parte di codice che esegue un determinato calcolo (ad esempio l’inversione di una matrice), e che può essere chiamata tutte le volte che occorre sia dal Command Window sia da un programma. Le funzioni Matlab scritte dal programmatore vengono registrate in appositi m-files. La prima istruzione utile di una funzione (dopo eventuali commenti) deve essere:

���������).��������(�����$��+,�����?��� ���� .��������(����$��". Le variabili di output sono i dati che la funzione calcola e restituisce al Command Window o al programma chiamante. Se vi è una sola variabile la coppia di parentesi [] può essere omessa. Le variabili di input e di output possono essere variabili semplici o array, tra le variabili di input ci possono essere anche stringhe di caratteri. �0�������(��������� ���� (ad esempio dal Command Window): >> [var. output]=nome(var.input) oppure nome(var. input) (anche inclusa in una espressione) La chiamata delle funzioni interne Matlab segue ovviamente questa sintassi. Ad esempio: >> [m, n]=size([1 2 3; 7 6 3]) m=2 n=3 Le funzioni vengono scritte (o modificate) con l’editor di Matlab come script file. Sono salvate con l’estensione m. Esempio 1 % calcolo del fattoriale di n>0 function fatt=fattoriale(n) fatt=1; for i=1:n, fatt=fatt*i; end Esempi di chiamate: >> fattoriale(5) ans = 120 >> 2+2 % il valore dell’espressione è assegnato alla variabile ans, che può essere ans = 4 % riutilizzata. >> zz=7+fattoriale(ans) % il risultato è ora assegnato a zz, non ad ans zz = 31 Esempio 2 % questa funzione riceve r come variabile di input. Costruisce 2 vettori x ed y contenenti le % coordinate di 100 punti della circonferenza di centro 0 e raggio r, fa il grafico di tale % circonferenza e restituisce al punto di chiamata i vettori x ed y. function [x,y]=circlefn(r); % [x,y] sono le variabili di output (in questo caso 2 vettori) theta=linspace(0, 2*pi,100); x=r*cos(theta); y=r*sin(theta); plot(x,y); axis(‘equal’) Esempi di chiamate: >> R=5 >> [x,y]=circlefn(R ); % Fa il grafico, e passa al Command Window i vettori x,y >>[cx,cy]=circlefn(4); % cerchio di raggio 4 “ “ “ “ >> circlefn(1); % fa solo il grafico del cerchio di raggio 1 Si noti che nella chiamata le variabili di output se non servono possono essere omesse. L’omissione del carattere ; dopo le chiamate comporta comunque il display dei vettori x, y sul monitor.

Page 8: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

8

Le variabili interne ad una funzione sono locali, e vengono azzerate dopo l’uscita. Per dichiarare le variabili globali (cioè comuni alla funzione, e Command Window o script file): ������2*�!7, sia nella funzione che nella parte chiamante. Anche se chiamate da uno script file, le funzioni devono essere in un m-file ��������.. 1�5���&�� �����$����������������������

Matlab possiede molte funzioni interne al sistema per la risoluzione di problemi di Calcolo Numerico, ad esempio il calcolo di un integrale, la risoluzione di un’equazione differenziale, il calcolo delle radici di una equazione eccetera. L’uso di queste funzioni comporta che l’utente prepari una sua propria funzione descrivente il �����������problema (ad esempio: la funzione integranda nel calcolo di un integrale definito). La risoluzione completa del problema appare dunque costituita da 3 moduli:

Nel primo modulo vi è il programma (o semplicemente l’istruzione del Command Window) che chiama la funzione Matlab risolvente il problema (ad esempio l’integratore). Tra i parametri della chiamata deve essere presente e ���������������� l’identificatore della funzione particolare (ad esempio l’integrando), preparata dall’utente. La funzione Matlab utilizzerà questo identificatore per chiamare, quando è necessario, la funzione particolare. Se la funzione particolare dipende da parametri c1, c2..dello Script File che possono variare durante l’elaborazione, essi possono essere dichiarati globali (sia nello Script File che nella funzione). Esempio. Si vuole risolvere l’equazione differenziale y ‘= sen(x*y)+cos(x), con la condizione iniziale y(0)=1, nell’intervallo [0,2] (problema di Cauchy). Una importante funzione Matlab per il problema di Cauchy è la ode23, la cui sintassi è: [x, y]=ode23(‘funzione_utente’, intervallo, y0) La risoluzione della particolare equazione comporta dunque che l’utente abbia preparato il seguente m-file: % particolare equazione differenziale function y1=fxy(x,y) y1=sin(x*y)+cos(x) dopo di che il problema sarà risolto con la seguente chiamata della ode23: >> [xa,ya]=ode23(‘fxy’,0,2,1) Il risultato sarà una coppia di vettori xa,ya contenenti rispettivamente le ascisse e le ordinate della soluzione calcolata in un numero adeguato di nodi.

3���@�������3����@������������(����������� La funzione piu’ importante per la grafica e’ la plot, che richiede come argomenti principali due vettori (vettx e vetty) contenenti le coordinate (xi,yi) dei punti da graficare. La sintassi e’: $��� .���2*�.���!*�=�$ �����" ( ‘opzioni’ puo’ mancare) vettx e vetty devono avere lo stesso numero di componenti. Opzioni indica il colore, ed il tipo o marcatore di linea. Esempi colore: m magenta, r rosso, g verde, b blu, w bianco, k nero, y giallo Esempi tipo e marcatore di linea: - continua (default), -- tratteggiata, : punteggiata, -. Punto-linea, + segni più, o cerchietti, * asterischi, x lettere x,….

Script File, oppure Command Window

Funzione di calcolo numerico di Matlab

Funzione particolare (dell’utente)

Page 9: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

9

Il grafico appare in una finestra apposita (Graphics Window) sovrapposta al Command Window, e che puo’ essere inviata alla stampante mediante l’istruzione print. Vi sono regole di default che consentono di ottenere una visualizzazione ottimale del grafico, riguardo alle scale degli assi, la forma del riquadro ecc. Tuttavia vi sono istruzioni che consentono di personalizzare il grafico aggiungendo scritte, commenti, e variando le regole di default della presentazione. Tali istruzioni debbono seguire la plot, e la loro esecuzione modifichera’ via via nel modo voluto il grafico creato dalla plot. Ne ricordiamo alcune: xlabel(‘scritta’) etichetta l’ asse x, ylabel(‘scritta’) etichetta l’asse y, title(‘scritta’) pone un titolo in testa al grafico, text( coord.x, coord.y, ‘scritta’) inserisce ‘scritta’ nel grafico con inizio alle coord. specificate. grid: inserisce nel grafico una griglia. axis([xmin,xmax,ymin,ymax]) limita il grafico all’intervallo specificato axis(‘equal’) Stessa scala sulle ascisse e le ordinate axis(‘square’) Nel grafico i 2 assi hanno la stessa lunghezza print Il grafico (nella attuale fase di formazione) viene stampato (La stampa può essere effettuata anche dal menù file di Graphics Window) Esempio 1: Grafico di un cerchio di raggio 1 >> theta=linspace (0, 2*pi, 100); >> x=cos(theta); y=sin(teta); % sono creati i vettori x e y, lunghi 100 >> plot(x,y, ‘y’) % effettua il grafico sulla finestra Graphics Window. >> axis(‘equal’) % stessa scala su ascisse e ordinate (se no appare una ellisse) >> xlabel(‘x’) % Etichetta l’asse x >> ylabel(‘y’) % “ y >> title(‘cerchio’) % Titolo del grafico >> print % scrive il grafico sulla stampante di default Esempio 2: Grafico di più funzioni sulla stessa figura: >> x=linspace(0, pi/2, 50); >> plot(x,sin(x)) % grafico di sin(x) tra 0 e ������� ������è un vettore lungo 50) >> y=1-x.^2; % elevaz. al quadrato del vettore x punto per punto >> plot(x, sin(x),x, y,’g*’), grid % grafico di 2 funzioni con griglia. La lineadella 2°

% funzione è verde e formata da asterischi.

0 0.5 1 1.5 2

-1.5

-1

-0.5

0

0.5

1

3����@�������(����� �����(������������������������ fplot(‘funzione’,[xmin,xmax]) Con la fplot basta dare l’intervallo del grafico, il vettore x viene calcolato automatic. da Matlab. La funzione deve essere data come costante carattere. Usare operatori del tipo ‘punto per punto’.

Page 10: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

10

Esempio: >> fplot(‘exp(-0.1*x) .* sin(x)’,[0,3] ) % attenz. .* : anche la funzione è un vettore ! >> xlabel(‘x’), ylabel(‘y=e^(x/10) sin(x)’), title(‘grafico’) 3�1�����.��(����.�����@�������(�������� �������$�������Il grafico delle curve di livello di una superficie z=f(x,y) richiede una discretizzione mediante un adeguato numero di nodi (xi,yi) del rettangolo in cui si vuole il grafico. Le funzioni necessarie sono: )A*B+,���0���( �2'(2'2*��!'(!'!" Crea 2 matrici X ed Y descriventi la discretizzazione. Ogni coppia di elementi corrispondenti (xij,yij) fornisce le coordinate di un nodo della discretizzazione. )�*0+,������� A*B*C" Crea il grafico. X,Y, sono le 2 matrici precedenti. Z è la matrice che contiene i valori della z=f(x,y) nei nodi della discretizzazione. C è una matrice che contiene le coordinate dei punti che contour utilizza per costruire il grafico. h è un vettore associato alle curve di livello che appariranno nel grafico. Se C ed h non interessano possono essere omesse, tuttavia sono necessarie se si vuole etichettare le curve con le rispettive quote mediante la successiva istruzione ����� �*0"��Contour(X,Y,Z, [a,b]) limita le curve a quelle comprese tra le quote a e b. Esempio: curve di livello della funzione f(x,y)= 2x2-3xy+x-y2+1 >> [X, Y]=meshgrid(-2: 0.2: 2, -2: 0.2: 2); % Crea in un rettangolo la griglia di nodi xi, yj >> Z= 2*X .^2 – 3* X .*Y +X- Y .^2 +1; % In ogni nodo della griglia calcola la funzione Z ( attenz: operazioni ‘elemento per elemento!) >> [C,h]=contour(X,Y,Z); % Crea il grafico (C è la matrice dei nodi, ) >> clabel(C,h) % etichetta le curve di livello (……. Altre scritte se occorre) >> print % Stampa il grafico

Per la f(x,y)=0 (funzione implicita): [C,h]=contour(X,Y,Z, [0 0]) % serve solo la curva di livello 0 della z=f(x,y) 3�3��@�������1#�� �����"�Oltre che mediante le sue curve di livello, una funzione z=f(x,y) può essere visualizzata (ed in modo più naturale) mediante il grafico della sua superficie, il che si ottiene mediante la funzione mesh. Ad esempio vogliamo graficare nel rettangolo -2������-2�������������

Page 11: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

11

[ ]))(( 222 ������ +−−= Il programma sarà: >> [X,Y]=meshgrid(-2:0.1:2); % poiché è specificato solo l’intervallo sull’asse x, l’intervallo sull’asse y sarà preso uguale. >> Z= X.*exp(-((X-Y.^2).^2+Y.^2)); >> mesh(X,Y,Z), xlabel(‘x’), ylabel(‘y’), zlabel(‘z’)….

Ovviamente l’aspetto della superficie cambia a seconda del punto P(x,y,z) dello spazio da cui la guarda l’osservatore. Tale punto può essere definito dai 2 angoli ��������� � � � ����� ������primo dei quali specifica la rotazione orizzontale dall’asse y per raggiungere la proiezione di P sul piano (x,y), ed il secondo l’angolo tra la retta congiungente P con (0,0,0) ed il piano (x,y). Matlab pone per default questi angoli a -37.5° e 30°. Tuttavia questi valori possono essere cambiati mediante l’istruzione .��D �*� "����� Così se, ad esempio, l’ultima istruzione mesh fosse seguita da view(40, -30) la superfice assumerebbe un aspetto diverso (sarebbe vista ‘dal di sotto’ anziché ‘dal di sopra’) Grafico di curve dello spazio: $���1 2*!* *��$ �����"�� dove x, y, z sono i vettori delle coordinate (xiyizi) dei punti della curva, ed anche ‘opzioni’ sono analoghe a quelle viste per la plot. La plot3 può essere seguita, come la mesh, dalla funzione view. Esempio: >> z=0:0.2:20; >> plot3(sin(z),cos(z),z), title(‘Spirale’) Matlab possiede moltissime altre funzioni per la grafica e la creazione di figure in movimento. Per approfondire Matlab: W.J.Palm III – Matlab 6 per l’ingegneria e le scienze (Mc Graw Hill) Testo di consultazione per Analisi Numerica: G. Pesamosca – Metodi dell’Analisi Numerica (Ediz. Kappa)

Page 12: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

12

��EF�-�8��

&������������������(��('

Command Window: è lo spazio di lavoro Launch Pad: Elenca i Toolbox e altri programmi installati (Simile ad Esplora Risorse), Command History. Lista dei comandi immessi in precedenza. Un comando può essere traferito da questa finestra nel command Window (selezionarlo, poi KS e trascinarlo). ���G������'�&������;(��������D���H����H��(�D�Però può cambiare secondo la finestra attiva (cioè le 3 standard + WorkSpace+Editor+Graphic Window) ��(������(�����G�(������������������������(��(: Matlab �� � � ���op Layout � ��� ��������� ����!����à: default (standard), Command Vindow only, Simple (Command Vindow e Command History), ecc H��>�6$���. E’una finestra (normalmente nascosta) che contiene i nomi e le dimensioni di tutte le variabili usate nella sessione. Per vederla: Matlab �� � "�#��� �� La stessa cosa per vedere altre finestre, eventualmente non presenti nella finestra standard (es.: Matlab �� � $�����%�&� � ���������(�������������(��������������: Matlab '&�� (� �#�(����&�"�&���� �# �Command History, oppure Work Space…) Ma per azzerare il Command Window anche: >> clc (più semplice) ;������ �����(��������������: X ��������(�������!:

Matlab File Edit View Vindow

Command Window Launch Pad

Command History

Current Directory …

Page 13: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

13

è la directory dove vengono registrati o letti i programmi che vanno in esecuzione, cioè gli m-files. Normalmente è settata su c:\matlab6p1\work. Per cambiarla, selezionare sul pulsante (…) immediatamente a destra, la directory (o il dischetto) voluto (ad esempio a:\). ������(����(�D: Possibilità di ritornare su un comando precedente per modificarlo, rieseguirlo, ecc. spostandosi con le freccette ( � �� � ��)�) La variabile ���. Può essere riutilizzata: ad esempio se l’ultimo risultato vale 0.5, 4*ans varrà 2. ;������ �����(�����������.�������: >> clear ;������ �����(�����.��������.�*�.�, .. >> clear v1, v2… 4���: uscita da Matlab

Page 14: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

14

��������������$$���� �������������������������

�I�������$��� ��������$$������� ���������������4��(������

�I�������$��� �����$������������Siano y=[y0,y1,…yn] gli n+1 valori di una funzione negli n+1 nodi x=[x0,x1,…xn]. L’istruzione a=polyfit(x, y, n) assegna al vettore a i coefficienti del polinomio interpolatore Pn(x)=anx

n+ an-1xn-1+ …. +a1x+ a0 (1)

I coefficienti sono sistemati in a così: a=[an,an-1,…a0] (2) Se x è un vettore di elementi x1,x2,…xm ed a il vettore (2), l’istruzione y=polyval(a,x) costruisce il vettore y=[y1,y2,…ym] dei valori di Pn(x) nei punti x1,x2,…xm . ;����� ��. Costruire la tavola di f(x)= |cos(x)| in [0,2] e passo h=0.25. Graficare questi punti con cerchietti. Costruire quindi l’interpolatore di grado 8 P8(x) sui dati della tavola, e graficare (sulla stessa figura precedente, ma con passo 0.05) anche P8(x) ed f(x). Nota. La funzione 0��( tiene la figura precedente attiva, in modo che su di essa si possano inserire altri grafici. 0��(���� annulla l’effetto di hold. >> x=0:0.25:2; y=abs(c0s(x)); >> plot(x,y,’o’) >> hold >> a=polyfit(x,y,8) >> xplot = -0.1:0.1:2 >> plot(xplot, polyval(a,xplot)), title(‘punti tabulate e interpolatore di grado 8’) >> pause >> plot(xplot, abs(cos(xplot))) Si nota che f(x) è approssimata in modo molto modesto da P8(x). Si tenga però presente che f(x) è solo di classe C0, quindi l’errore sarà elevato. La formula del resto non è valida, l’errore non è maggiorabile. I����$$������� �����$���������������������4��(������Sia m il numero dei nodi e dei valori della funzione, ed x=[x0,x1,…xm], y=[y0,y1,…ym] i 2 vettori. Sia inoltre m il grado del polinomio Pn(x) che si vuol costruire. Se m>n+1, allora la a= polyfit (x,y,n) calcola Pn(x) nel senso dei minimi quadrati, cioè Pn(x) sarà il polinomio tale che

( ) min)(1

22 =−= ∑=

P

L

LLQ���σ

Il vettore a è ordinato sempre come in (2). Anche la y=polyval(a,x) ha lo stesso significato del punto precedente. ;����� ���(fatto a lezione): Determinare la retta che approssima nel senso dei minimi quadrati i seguenti dati (xi,yi): (-1 0), (0 2.5), (1 3), (2 4), (3 3.5). >> x=[-1 0 1 2 3]; y=[0 2.5 3 4 3.5]; >> a=polyfit(x,y,1) % non è necessario determinare la matrice �* a = 0.85 1.75

Page 15: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

15

I�1��$$������� ���������������4��(���������$������������Matlab possiede una istruzione diretta per risolvere il problema dell’approssimazione ai minimi quadrati solo nel caso polinomiale. Nel caso si voglia approssimare i dati sperimentali (x1,y1), (x2,y2),…. (xm,ym) secondo un modello non polinomiale: y = a1 ϕ1(x)+ a2 ϕ2(x)+ ….. + an ϕn(x) (n < m) occorrerà procedere così: 1) – costruire la matrice Φ(m x n) di elementi ϕi,j = ϕj (xi), 2) – a= Φ\y con y (m x 1) vettore ������ dei dati sperimentali. Infatti, se Φ è quadrata e regolare l’operazione Φ\y fornisce (come già visto) la soluzione del sistema Φx = y. Se invece Φ ha più righe che colonne Φ\y fornisce la soluzione del sistema delle equazioni normali, cioè Ha=D, con H=ΦTΦ, e D=ΦTy. (vedere esercizio) I�3� �����$��� �����$��������������������

L’interpolazione polinomiale a tratti è un’altra tecnica che consente di approssimare mediante polinomi di grado basso una funzione continua f(x), di cui sono noti i valori nei nodi (x1,y1), (x2,y2),…. (xm,ym). Vi sono due varianti principali: Interpolazione lineare a tratti. Se x* è un punto compreso tra xi ed xi+1, la f(xi) viene approssimata mediante interpolazione lineare tra xi ed xi+1. Interpolazione a tratti di terzo grado. Se x* è un punto compreso tra xi ed xi+1, la f(xi) viene

approssimata mediante l’interpolatore di terzo grado tra xi-1 , xi , xi+1 ed xi+2: xi-1 xi x* xi+1 xi+2

Ovviamente, se x* cade in [x0, x1] o in [xn-1 xn] l’interpolatore di 3° grado viene costruito in altro modo (nel primo caso, sui nodi x0, x1 , x2 , x3). La funzione che consente l’interpolazione a tratti è: yn = interp1(x, y, xn,’metodo’) in cui: x, y sono rispettivamente i vettori dei nodi xi, e dei valori yi, xn è il nuovo valore x* (ma può essere anche un array!) in cui si vuole valutare f(x), yn è l’approssimazione calcolata,e ‘metodo’ può essere ‘linear’ se si vuole l’interpolazione lineare o ‘cubic’ se si vuole l’interpolazione di terzo grado. Se ‘metodo’ è omesso l’interpolazione sarà per default lineare. Un altro importante ‘metodo’ di approssimazione a tratti, non visto a lezione, è ‘spline’. In ogni intervallo [xi,xi+1]è costruito un polinomio di 3° grado. L’insieme di questi polinomi garantisce non solo che sia yi= f(xi) nei nodi, ma anche che la derivata di f(x) sia continua nei nodi. Si può anche dire che la interp1 implementa alcuni criteri per ricostruire una funzione continua, a partire da alcuni suoi campioni discreti. ��������. a)- Costruire il vettore dei valori di f(x)=exp(-x2) nei nodi –2, -1.5, -1, …1.5, 2. b)- Graficare in [-2,2] l’approssimazione ottenuta mediante la interp1 con metodo ‘linear’ (per il grafico, si costruiscano 2 vettori x, y di 81 elementi, cioè con passo 0.05) c)- Graficare in [-2,2] l’errore exp(-x2)-yappr. (yappr è l’approssimaz. lineare a tratti). Ripetere i punti a), b), c) con il metodo ‘cubic’ e con il metodo ‘spline’. Per le approssimazioni ‘linear’ e ‘cubic’ (ma specialmente ‘linear’) si notano dei punti angolosi in corrispondenza dei nodi, molto evidenti nella funzione di errore exp(-x2)-yappr.. Il motivo è che quando x* passa da xi-ε ad xi+ ε l’approssimazione viene effettuata con un ������� polinomio interpolatore. Poiché i 2 polinomi relativi ad [xi-1,xi] ed [xi,xi+1] assumono xi lo stesso valore yi , è garantita la continuità della approssimazione, mentre non c’è nulla che possa garantire anche la continuità della derivata prima. Tale continuità è invece garantita dalle approssimazioni di tipo ‘spline’.

Page 16: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

16

J�#���.� �������K��(������

�J��� #������� �� �������(���.� �����diff(X)���Se X è un vettore (riga o colonna) di dimensione n , fornisce il vettore (riga o colonna) delle differenze, di dimensione n-1. Se X è una matrice n x m, fornisce la matrice delle differenze (colonna per colonna), che sarà di dimensione (n-1) x m. diff(X,k) come sopra, ora calcola la differenza k-ma. Dati 2 vettori X, Y tali che (x1,y1), (x2,y2), ….(xn,yn) siano i valori di una tavola di una funzione f(x), allora diff(Y) ./ diff(X)sarà un vettore di dimensione n-1, le cui componenti (yi+1-yi)/ (xi+1-xi) forniscono un’approssimazione della derivata f ‘(x) nei nodi. J���K��(������. integrale=quad (funz, a, b, tol, trace, p1, p2,…) Calcolo dell’integrale definito tra a e b di una funzione f(x), data analiticamente. Argomenti obbligatori: funz : è l’identificatore dell’integrando. Può essere dato come function-file, e in questo caso deve essere preceduto da @. Può essere dato anche come funzione definita entro la quad, o compresa tra apici oppure mediante la funzione inline. Esempio: quad(‘exp(-x.^2).*sin(x)’, ……. oppure quad(inline(‘exp(-x.^2).*sin(x)’),……... ������.: siccome x è considerata internamente come un vettore costruito sui nodi, usare gli operatori vettoriali .* ./ (ad es: x .*sin(x)+2) a, b : sono gli estremi di integrazione. Argomenti opzionali: tol : indica l’errore assoluto massimo tollerato. Il valore di default è 10-6. trace : se è un valore +�,��è visualizzata una tabella indicativa di come si è svolta l’integrazione (ad esempio: mostra il numero di valutazioni dell’integrando in vari subintervalli) . p1, p2 …. : valori numerici di eventuali parametri che compaiono nella funzione integranda. Se p1,p2… sono presenti, l’integrando deve essere del tipo funz(x,p1,p2…), cioè x 1° parametro. Se compare uno qualsiasi degli argomenti opzionali, debbono comparire anche i precedenti. Quad applica il metodo di Simpson adattivo, con un massimo di 10 suddivisioni. [int nvf]=quad1( ‘funzione’, a, b, tol, trace, p1, p2,…) Esattamente come quad, ma lavora con un algoritmo più sofisticato riuscendo ad effettuare il calcolo anche se l’integrando è poco regolare. . Il parametro di uscita nvf è il numero complessivo di valutazione dell’integrando. Se però vi sono singolarità evidenti, conviene spezzare l’intervallo di integrazione sulla su di esse integrale=trapz(y) l’argomento y è un vettore di dati. trapz calcola l’integrale dei dati con il metodo dei trapezi, supponendo il passo h=1. Se h ≠1, moltiplicare il risultato per h. integrale=dblquad(‘funzione’,c,d,a,b,tol,metodo)

Calcolo dell’integrale doppio, esteso ad un rettangolo: ∫ ∫=E

D

G

F

������� ),(int

���. esattamente come in quad. ����(�. È il metodo con cui sono calcolati internamente gli integrali in una dimensione (quad o quad1). Per default è quad

Page 17: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

17

Esercizio. E’ data f(x)= x sin(x). Effettuare il grafico e l’integrazione con la traccia in [-1,1] >> fplot(‘x .* sin(2*x)’,[-1,1]), grid, title(‘f(x)=x sin(2 x)’) >> int=quad(@fint, -1,1,1e-5,1) l’m-file della funzione integrando è: % funzione integrando fint function f=fint(x) f=x .*sin(2*x); Tuttavia si poteva anche evitare l’m-file e chiamare quad così: >> int=quad(‘x.*sin(2*x)’, -1,1,1e-5,1) % oppure quad(inline(‘x.*sin(2*x)’)…. La traccia mostra una notevole difficoltà del calcolo nell’intorno dell’origine, dove la |fiv(x)| è più alta. Il grafico della fiv(x) può essere visualizzato con le seguenti istruzioni Mathematica: f=x Sin[2 x] D[f, {x,4}] Plot[%, {x,-1,1}] L��/�(����(������4�� ����������������. L������������(��������(����(������4�� ������ 2",M. La funzione Matlab è fzero: sol=fzero(@funz, x0, opzioni, p1, p2,…) Argomenti obbligatori: ���� è l’identificatore della f(x). Può essere una Function File od una funzione matlab. ���� deve essere preceduto da @ ( ad esempio: fzero(@cos, 1.5 ) ). Può anche essere una funzione definita direttamente entro la fzero mediante la funzione inline, come dall’esempio seguente che calcola la radice di x3-2x-5=0 vicina a 2: >> x=fzero(inline(‘x.^3-2*x-5’),2) x=2.0946 Attenzione: x è inteso come un vettore. Definire le operazioni elemento per elemento. (per compatibilità con la versione matlab precedente la funzione può essere data semplicemente racchiusa tra apici, ad es: fzero(‘cos’,1.5) oppure fzero( ‘x.^3-2*x-5’,2) ) �� prima approssimazione della radice. Ma si può dare anche l’intervallo in cui cade: [x0 x1]. Argomenti opzionali: ������� sono parametri che consentono di definire l’errore max tollerato, oppure il tipo di display (ad es. il risultato di ogni iterata). Vedere Help per dettagli. ������,… eventuali parametri della f(x), che sarà del tipo: funzione(x, p1, p2…). Se ������� non è presente, indicarlo con [] (prima di p1, p2…) [sol, valf]=fzero(@funz, x0, opzioni, p1, p2,…) Come la precedente, in più fornisce il valore valf di f(x) sull’approssimazione sol della radice. L’algoritmo è basato sul metodo di Brent, alternato se occorre con iterate che possono essere a seconda dei casi bisezioni, oppure secanti. L’algoritmo ��� calcola radici di molteplicità 2. Esempio. Calcolo della radice di f(x)=sin(x) – ex+5=0, approssimaz. iniziale x0=1. % funzione frad function f=frad(x) f = sin(x)-exp(x)+5; % Chiamata da Command Window >> x=fzero(@frad, 1)

Page 18: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

18

L������������(�������������(����(�������4�� ����������������Sia anx

n+ an-1xn-1+…. a1x+ a0=0 l’equazione. Il calcolo delle n radici è fatto dalla roots:

v=roots([an an-1 …. a1 a0]) Il risultato è il vettore che contiene le n radici (reali e complesse). La funzione: a=poly(v) ricostruisce dal vettore v delle radici il vettore a dei coefficienti del polinomio. Esempio. Radici di x3-6x2-72x-27=0 >> v=[1 -6 -72 -27]; >> roots(v) 12.1229 -5.7345 -0.3884 % attenz. : esce un vettore colonna �M��;4�� �����(������� �������(������. �M����F����� �����$�����$���. Matlab possiede diverse funzioni per la risoluzione del problema di Cauchy. Tutte richiedono che i sistemi o le singole equazioni di ordine superiore al primo siano scritti nella forma di sistemi del primo ordine:

Le funzioni principali sono ode23 ed ode45, che implementano rispettivamente i metodi di Runge-Kutta del 2/3 ordine e del 4/5 ordine. Sintassi della ode23 (e dalla ode45) nella loro forma più semplice: [x,y]=ode23(@funzione_utente, intervallo, y0) dove: x è il vettore colonna dei valori della variabile indipendente in cui la soluzione è ottenuta. y è una matrice di n colonne (n=numero di equazioni). Ciascuna colonna contiene nell’ordine i valori di y1, y2,….yn corrispondenti ai valori della variabile situati nel vettore x. funzione_utente è l’m-file creato dall’utente, che dovrà calcolare i secondi membri del sistema (od equazione). Nella sua forma più semplice, sarà del tipo: function der = funzione_utente(x ,y) e dovrà porre nel vettore ������ der i secondi membri del sistema, calcolati sulla variabile indipendente x, e sugli elementi del vettore colonna y intervallo è l’intervallo della x in cui è richiesta la soluzione (vettore riga, ad es. [0 1]. y0 è il vettore delle condizioni iniziali (ad es [y10 y20… yn0] ). Si tenga presente che i nodi xi (componenti del vettore di output x) non sono prevedibili a priori, in quanto gli integratori lavorano �����������������. Se però si desidera avere la soluzione in m determinati punti xi, si deve porre il parametro di input intervallo così: [x0 x1 x2 ….xm] La ode23 e la ode45 lavorano con passo variabile, calcolando ad ogni passo una stima dell’errore ed accettando il risultato se l’errore relativo è < 10-3, e l’errore assoluto è (per ogni componente della soluzione) < 10-6. La sintassi completa degli integratori Matlab è la seguente: [x,y]=ode23(@funzione_utente, intervallo, y0,opzioni, p1, p2,…) ( analogo per ode45) dove:

==

====

0,021

0,2022122

0,1012111

)(),,,(

)(),,,(

)(),,,(

QQQQQ

Q

Q

������������

������������

������������

���

Page 19: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

19

opzioni è un argomento opzionale, creato con la odeset (vedere Help ). Mediante questo argomento si può ad esempio variare la tolleranza relativa ed assoluta (RelTol ed AbsTol) rispetto ai valori di default 10-3 e 10-6. Un’altra possibilità è quella di generare in output il grafico delle componenti della soluzione (ma questo si può fare anche utilizzando la plot sugli output x, y). p1, p2,… sono eventuali parametri della funzione_utente, che sarà in questo caso del tipo function der=funzione_utente(t, y, p1, p2,..) Altre funzioni Matlab per la soluzione delle equazioni differenziali: ode23 - equazioni non-stiff, metodo Runge-Kutta di ordine 4/5 ode45 - equazioni non-stiff, metodo Runge-Kutta di ordine 4/5 ode15s - equazioni stiff, metodo di ordine variabile. ode113 – equazioni non stiff, algoritmo Adams Bashfort-Moulton di ordine variabile ode23t - equaz. moderatamente stiff, metodo dei trapezi Le equazioni stiff sono quelle difficili da risolvere numericamente. Per un’equazione del 1° ordine

y’=f(x,y) ‘stiff’ significa che la ����

∂∂ ),( è in modulo elevata.

;����� ����. Calcolare in [0,2] la soluzione del problema di Cauchy y’=x sin(y)+1, y(0)=1. Visualizzare la soluzione in forma di tabella e di grafic0. >> [x,y]=ode23(@fode1, [0,2],1]; >> plot(x,y), grid, title(‘Soluzione equaz. differ. y’=x sin(y)+1’) >> soluzione=[x y] e la funzione_utente fode1 sarà: % m-file della funzione-utente function y1=fode1(x,y) y1=x*sin(y)+1 ; ;����� ����. Calcolare in [0,1] la soluzione del problema di Cauchy y’’=3 y’+2 y2+y-1, condiz. iniziali y(0)=0, y’(0)=1/2. Grafico della soluzione e della derivata prima, e visualizzazione della soluzione e della derivata prima in forma di tabella. Matlab richiede che l’equazione del 2° ordine sia scritta in forma di sistema del 1° ordine, cioè:

−++=′=′

123 12122

21

�����

��

==

2/1)0(

0)0(

2

1

>> [x,y]=ode23(@fzutente, [0,1],[0 1/2]); >> plot(x, y(:,1),’r’, x, y(:,2),’g’), grid, title(‘Soluzione equaz. differ. y’=x sin(y)+1’) >> [x y(:,1) y(:,2)] e la funzione_utente fzutente sarà: % m-file della funzione-utente function y1=fzutente(x,y) y1(1,1)=y(2) ; % Attenzione ! il vettore delle derivate deve essere un vettore colonna y1(2,1)=3*y(2)+2*x*y(1)^2+y(1)-1 ;

Page 20: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

20

pstat vettore a 6 componenti, rappresentanti una statistica sulle ‘performances’ del risolutore:

1- numero di passi corretti, 2- numero di passi non riusciti, 3- numero di valutazioni della funzione_utente, ……(vedere la struttura opzioni, creata con odeset). Questo parametro di uscita è opzionale.

funzione_utente è la funzione (M-file) che descrive il sistema. tspan - Intervallo di integrazione y0 - vettore colonna delle n condizioni iniziali I parametri di input successivi sono opzionali. opzioni – e’ una struttura creata prima della chiamata del risolutore, mediante la funzione odeset: opzioni=odeset(‘opzione1’, valore1, ‘opzione2’,valore2,….) Esempi di opzioni sono: reltol tolleranza relativa . Valori possibili: 1e-3 (default), 1e-4, … abstol tolleranza assoluta. Val. possibili: 1e-6 (default), 1e-4,… outputfcn Valori possibili per l’output, ad es: ‘odeplot’ (grafico), ‘odephas2’ (grafico nel piano delle fasi), ‘odeprint’ (visualizza la soluzione man mano che procede)... stat deve essere ‘on’ se è richiesta la statistica nell’ output del risolutore . Esempio: opzioni=odeset(‘rltol’,1e-4, ‘stats’,’on’,’outputfcn’,’odeplot’)

Page 21: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

21

p1,p2,…eventuali parametri che compaiono dentro il sistema di equazioni differenziali. Funzioni di utilità nell’integrazione delle equazioni differenziali: odefile – è un help sulla sintassi dei vari integratori. odeset – Importante. Già visto sopra. odeplot, odephas2, odeprint: varie funzioni, che possono essere chiamate dalla opzione outputfcn �M�1���F��$ �����=�.�����. Se si vuole arrestare il calcolo al verificarsi di una certa condizione (ad esempio: yn=h , e/o y’n=k) occore porre ad ‘on’ l’opzione events: opzioni=odeset(…‘events’, ‘on’,….) In tal caso, l’integratore passa alla funzione_utente un parametro (flag) che viene posto ad ‘events’ quando deve essere verificata la condizione. Quindi: ������ intervallo=…; y0=….; opzioni= odeset(…‘events’, ‘on’,…. ) [t,y,te,ye]=ode23(‘funz_utente’, intervallo, y0, opzioni) (te, ye sono il valori della var.indipendente e della y al verificarsi dell’evento) funz_utente: function [y1, isfinale, derfinale]=funz_utente(t,y,flag) if isempty(flag) % se ode23 non ha assegnato a flag il valore ‘events’, calcola y1 secondo l’equazione (o % il sistema) differenziale. else switch flag case ‘events’ y1=y – h % assegna ad y1 la condiz. d’arresto ���������������. Se l’equaz. è del 2°

% ordine, assegnare la condizione d’arresto se occorre anche alla derivata isfinale=1 % se la condizione è proprio una condizione di uscita dalla ode23. otherwise Error(‘errore non previsto’) end end derfinale deve figurare tra gli argomenti anche se non è utilizzato.

Page 22: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

22

�����������������������EF�-�

Le variabili simboliche debbono essere dichiarate mediante la �!��: >> syms x y z1 z2 % crea 4 variabili simboliche x y z1 z2 Sulle variabili simboliche si può operare con gli operatori e le usuali funzioni matematiche, formando espressioni simboliche. Non è necessario dichiarare simboliche le espressioni create mediante sole variabili simboliche. Ad es: >> f1=x*y^2-sqrt(y+x) >> f2=3*x-2*y^2 >> f3=f1+f2 E’ possibile manipolare le espressioni simboliche mediante funzioni analoghe a quelle usate in Matematica, cioè expand, factor, collect, … Es: >> syms x y >> expand(x+y)^2 x^2+2*x*y+y^2 >> expand(sin(x+y)) sin(x)*cos(y)+cos(x)*sin(y) >> factor(x^2-1) (x+1)*(x-1) In una espressione simbolica E è possibile sostituire una varabile (o una intera parte di E) con un’altra variabile, un numero od un’altra espressione. La funzione Matlab è: ��� ;*�.*��"��dove v (vecchio) è la parte che viene sostituita, n (nuovo) è la parte che sostituisce. Esempio: >> syms x y a >> esp=x^3-2*y+1; >> t=subs(esp,y,a) x^3-2*a+1 >> z=subs(t,x,2) 9-2*a La funzione poly2sym(v) trasforma il vettore v contenente i coefficienti di un polinomio in un polinomio simbolico. Ad es: >> p=[1 3 -2 5]; >> poly2sym(p) x^3+3*x^2-2*x+5 Nota: il calcolo simbolico di Matlab presuppone che la variabile indipendente sia sempre x. Tuttavia questa posizione di default può essere cambiata. Se E è una espressione simbolica dipendente da una sola variabile simbolica (x), la funzione � $���(E, [xmin xmax]) grafica E nell’intervallo specificato. Il grafico può essere migliorato utilizzando le solite funzioni axis, xlabel, ylabel, grid….. Esempio: syms x f=x^2-6*x+3 ezplot(f, [-2 2]), xlabel(‘x’), ylabel(‘f’) Calcolo simbolico dell’Analisi Matematica. Derivazione: (��� ;" Calcola la derivata dell’espressione simbolica E rispetto alla variabile di default x Esempio: >> syms x >> diff((1+x)/(1+x^2)) 1/(1+x^2)-2*(1+x)/(1+x^2)^2*x (��� ;*�"��deriva E n volte rispetto ad x (n=1,2,…). (��� ;* " deriva E rispetto alla variabile simbolica z (anche diversa da x) ��� ;" integrale indefinito di E rispetto alla variabile di default x. Esempio:

Page 23: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

23

>> syms x >> f=int(x*cos(x)) cos(x)+x*sin(x) ��� ;*�*" (con E espressione simbolica, a, b valori numerici) calcola l’integrale definito di E tra a e b, rispetto alla variabile di default x. Esempio: >> syms x >> i=int(x*sin(x),-1,1) 2*sin(1)-2*cos(1) >> double(i) % la funzione (���� ;"� converte l’espressione E dalla forma esatta alla numerica decimale (equivalente della funz. N in Matematica).

0.6023

Page 24: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

24

Esercizi 1– Di una f(x) conosciamo i seguenti valori (xi,yi): (0 0), (0.25 0.52), (1 1.56), (1.5 2.1).

a) stima della f ‘(x) in x=0.5

b) stima dell’ ∫5.1

0)( ����

MATLAB >>syms x ; >> t=[0 0.25 1 1.5]; y=[0 0.52 1.56 2.1]; >> a=polyfit(t,y,3) 0.2987 -1.0667 2.3280 -0.0000 >> p=poly2sym(a) p= 112/375*x^3-16/15*x^2+291/125*x-22034213020631/316912650057057350374175801344 >> der=diff(p) 112/125*x^2-32/15*x+291/125 >> d05=subs(der,x,0.5) 1.4853 >> i=double(int(p,0,1.5)) 1.7970 2 - Stabilire graficamente il numero e la localizzazione approssimata delle radici della equazione x(1-exp(-x/2))-1=0, che può essere posta nella forma x=x exp(-x/2)+1. Stabilire (sempre graficamente) la convergenza o meno verso di esse della iterazione xn+1=xnexp(-xn/2)+1. Applicare questa iterazione per calcolare una radice, arrestando il calcolo quando |xn+1-xn|<10-5. MATLAB (m-file) syms x fx=x*(1-exp(-x/2))-1 ezplot(fx,[-2,2]),grid pause fi=x*exp(-x/2)+1 df=diff(fi) ezplot(df,[-2,2]) pause x0=2; test=1; iter=0 while test, x1=x0*exp(-x0/2)+1; test=abs(x1-x0)>1.e-5; iter=iter+1; x0=x1; end [x0 iter] Stesso calcolo, fatto con la funzione Matlab fzero: >> fzero(inline(‘x*(1-exp(-x/2))-1’), 2) 1.72836 3- Stessa equazione dell’esercizio 2. Calcolo della 1° radice con il metodo di Newton, dopo aver verificato che la condizione di convergenza è soddisfatta, e con una scelta opportuna dall’approssimazione iniziale x0. Condizione di arresto: |xn+1-xn|<10-6 & |f(xn+1)|<10-6 MATLAB >> syms x >> f=x*(1-exp(-x/2))-1 >> f1=diff(f) % la derivate prima serve per il metodo di Newton >> f2=diff(f,2) >> ezplot(f, [-2,2]), grid % il grafico mostra la prima radice tra -2 ed 1 >> ezplot(f2, [-2,2]) % il grafico mostra che f ’’(x ) > 0 tra -2 e -1. Quindi il metodo di

Newton converge prendendo x0=-2 (estremo di Fourier) % m-file del metodo di Newton x0= -2; test=1; iter=0; eps=1e-6; while test,

Page 25: riccardo.broglia/.../CommMat/dispensamatlab.pdf · 1 ˘ Matlab è un sistema software integrato per il calcolo tecnico e scientifico. Le sue principali caratteristiche sono: - Linguaggio

25

f=x0*(1-exp(-x0/2))-1; f1=1-exp(-x0/2)+x0*exp(-x0/2)/2; x1=x0-f/f1; test= ~(abs(x1-x0) < eps & abs(f) < eps); x0=x1; iter=iter+1; end [x0 iter] -1.207 5

4 – E’ data la funzione f(x)=1/(1+x2). Trovare |)(|max22

�� ,9

[

,9

≤≤−= (massimo assoluto)

MATLAB >> syms x >> f=1/(1+x^2) >> d4=diff(f,4) >> d4=simplify(d4) (* espressione identica a quella di Mathematica *) >> ezplot(abs(d4), [-2,2]) (* il massimo per x=0 non è rilevabile perché fuori grafico *) >> ezplot(abs(d4,[-0.5, 0.5]), grid (* grafico più dettagliato, ora si vede che in x=0 d4 ~ = 24 *) [Matlab non possiede una funzione per il calcolo dei massimi e minimi relativi di una f(x)]