Calcolo Numerico - Matlab

Annamaria Mazzia

Matlab e Foglio Elettronico:una breve introduzione

Dipartimento di Ingegneria Civile Edile e AmbientaleUniversit degli Studi di Padova

Creative Commons Attribuzione-Non commerciale-Non opere derivate 3.0 Italy License

a.a. 2011/2012

INDICE

Indice ii

I Primi passi inMATLAB 1I.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1I.2 Avvio di MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2

I.2.1 Matrici e vettori in MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2I.3 Comandi utili . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5I.4 MATLABper scrivere ed eseguire programmi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

I.4.1 Strutture . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6I.5 Dati di input . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

I.5.1 Programma sul metodo delle bisezioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9I.6 Dati di output . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10I.7 Grafici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11I.8 Sulle potenzialit di MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12I.9 Applicazioni di MATLABnel Calcolo Numerico . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

I.9.1 Sullinstabilit numerica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13I.9.2 Sullinterpolazione e approssimazione di dati . . . . . . . . . . . . . . . . . . . . . . . . . 14I.9.3 Confronto tra schemi per equazioni differenziali ordinarie . . . . . . . . . . . . . . . . . 20

II Tabelle e grafici con il foglio elettronico 23II.1 Il foglio elettronico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23II.2 Inserire funzioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23II.3 Formule ricorsive . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28II.4 Fare grafici nel foglio elettronico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31

II.4.1 Grafici multicolonna . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36II.5 Inserire dati da files . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

ii

APPEN

DICE I

PRIMI PASSI IN MATLAB

Luniverso non potr essere lettofinch non avremo imparato illinguaggio ed avremo familiarizzatocon i caratteri con cui scritto. scritto in linguaggio matematico, ele lettere sono triangoli, cerchi edaltre figure geometriche, senza lequali umanamente impossibilecomprendere una singola parola.

Galileo Galilei

I.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1I.2 Avvio di MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2

I.2.1 Matrici e vettori in MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2I.3 Comandi utili . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5I.4 MATLABper scrivere ed eseguire programmi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

I.4.1 Strutture . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6I.5 Dati di input . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

I.5.1 Programma sul metodo delle bisezioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9I.6 Dati di output . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10I.7 Grafici . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11I.8 Sulle potenzialit di MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12I.9 Applicazioni di MATLABnel Calcolo Numerico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

I.9.1 Sullinstabilit numerica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13I.9.2 Sullinterpolazione e approssimazione di dati . . . . . . . . . . . . . . . . . . . . . . . . . . 14I.9.3 Confronto tra schemi per equazioni differenziali ordinarie . . . . . . . . . . . . . . . . . . 20

I.1 Introduzione

MATLABpu essere considerato un linguaggio di programmazione ad alto livello e, allo stesso tempo,un ambiente interattivo per lo sviluppo di algoritmi, per la visualizzazione e lanalisi di dati, e per il calcolonumerico.

1

I. PRIMI PASSI IN MATLAB

Per chi ha gi imparato a programmare in linguaggio FORTRAN, il passaggio a MATLAB semplice (ba-sta capire come cambiano le strutture di programmazione e il modo di programmare). In pi, MATLABhale sue functions, permette calcoli molto pi sofisticati di quelli di una calcolatrice, ed capace di fare graficie visualizzare i risultati in maniera diretta. Molto spesso la conoscenza di MATLAB un requisito richiestoper molte posizioni lavorative in ambito ingegneristico: a tal fine, utile imparare qualcosa di MATLAB .Ed altrettanto importante prendere atto del fatto che non esiste un linguaggio di programmazione che vadabene per risolvere tutti i problemi (quindi non basta conoscere e saper usare un solo linguaggio di program-mazione). Il FORTRAN si rivela la scelta migliore per risolvere algoritmi numerici complicati (un esempiotra tanti: risolvere problemi accoppiati di flusso e trasporto in mezzi porosi dove le incognite del problemasono dellordine delle migliaia e migliaia). MATLAB esso stesso un programma complesso (originaria-mente scritto in FORTRAN e successivamente riscritto in linguaggio C) che va bene per risolvere programmidelaborazione numerica, per lavorare conmatrici (MATLAB infatti sta per MATrix LABoratory laboratoriomatriciale) e per la grafica.

MATLAB un prodotto della The MathWorks : vedremo in particolare la versioneMATLAB7.10.0.499 (R2010a).

Un prodotto simile a MATLABma open source GNU Octave (si vada sul sito http://www.octave.org).

Nel seguito vedremo come si lavora inMATLABtenendo presente che quasi tutti i comandi che daremovalgono alla stessa maniera anche per Octave (eccezion fatta per linterfaccia grafico).

I.2 Avvio di MATLAB

Tralasciamo la parte relativa alla procedura di installazione diMATLAB(che dipende dal computer e dalsistema operativo in uso) e diamo per scontato che il programma sia presente sul proprio computer.

In ambiente Linux, per avviare MATLABbasta digitare il comando matlab da una finestra di termina-le. Una volta avviato MATLABcomparir il prompt dei comandi in una finestra come quella mostrata inFigura I.1. Per uscire dallambiente si digita exit dalla finestra dei comandi oppure dal menu File o ancoradallicona di chiusura (x in alto a destra della finestra).

Per imparare a usareMATLAB , inizialmente si prende familiarit con la finestra dei comandi eseguendocalcoli come se fosse una calcolatrice. Vediamo che il risultato viene assegnato ad una variabile detta ans eche nella Workspace si trovano informazioni su di essa (si veda Figura I.2).

I.2.1 Matrici e vettori inMATLAB

I vettori si possono scrivere come vettori riga o vettori colonna.Per esempio:

>>x=[1 2]>> x=[1 2]

x =

1 2

>> y=[10 20]

y =

1020

Un vettore riga x viene scritto usando le parentesi quadre e scrivendo le componenti del vettore una dopolaltra. Un vettore colonna pu essere creato facendo la trasposta di un vettore riga (mediante il simbolo ),oppure mettendo un punto e virgola dopo ciascuna componente:

2

I.2. Avvio di MATLAB

Figura I.1: Lambiente MATLAB consiste di una serie di finestre, alcune delle quali vengono aperte di default:la Command Window (la finestra dei comandi), la Current Folder (la directory corrente in cui si sta la-vorando), la Workspace (lo spazio di lavoro) e la Command History (la finestra con la storia dei comandidati).

Figura I.2: Primi calcoli in MATLAB : osservare cosa succede nella Command Window nella Workspace enella Command History.

3


>> y=[10;20]

y =

1020

In MATLAB i vettori altro non sono che un caso particolare di matrice a n righe e 1 colonna (vettorecolonna) o a 1 riga e n colonne (vettore riga). Quindi per scrivere unamatrice si scrivono i valori dellamatriceriga per riga andando a capo con il punto e virgola:

>> A=[1 2 3 4;5 6 7 8;9 10 11 12;13 14 15 16]

A =

1 2 3 45 6 7 89 10 11 12

13 14 15 16

Per matrici (e quindi per i vettori) si possono fare le operazioni di somma e prodotto in maniera moltosemplice. Date due matrici A e B si ha

C=A+B: matrice sommaC=A-B: matrice differenzaC=A*B: matrice prodotto (deve essere la matrice A di dimensione nm e la matrice B di dimensionem r altrimenti si ha unmessaggio di errore).C=A: matrice traspostaC=A.*B : matrice i cui elementi sonoC (i , j )= A(i , j )B(i , j )

Osserviamo che, per indicare un valore della matrice A basta specificare lindice di riga e di colonna: peresempio

>> A(2,2)

ans =

6

Per indicare gli elementi di tutta la colonna i si usa A(:,i), mentre A(i,:) indica gli elementi della riga i.>> A(:,2)

ans =

26

1014

>> A(2,:)

ans =

5 6 7 8

Loperatore due punti pu dunque essere usato per estrarre unintera riga o colonna da una matrice ounintera sottomatrice. Se vogliamo estrarre le ultime due righe e colonne dellamatriceA, si digita il comando

4

I.3. Comandi utili

>> M=A(3:4,3:4)

M =

11 1215 16

In questomodo abbiamo creato lamatriceM che ha come elementi quelli della sottomatrice diA con le ultimedue righe e colonne.

Per scrivere una matrice vuota A, invece, si usa listruzione A=[ ].Se si vuole risolvere un sistema lineare Ax = b in MATLABsi pu semplicemente usare una function

propria di MATLABnascosta nelloperatore \ (backslash): basta digitare il comando x= A \b, con bvettore colonna, per avere in x il vettore incognito. Ad esempio (si confronti con lesercizio ??):>> A=[0.2 1 0.2; 1 6.5 1.75; 0 2 2.25]

A =

0.2000 1.0000 0.20001.0000 6.5000 1.7500

0 2.0000 2.2500

>> b=[2.8 ; 19.25; 10.75]

b =

2.800019.250010.7500

>> x=A\b

x =

123

I.3 Comandi utili

Per lavorare meglio sia nella finestra dei comandi sia, successivamente, per scrivere ed eseguire dei veri epropri programmi, risultano utili le seguenti funzioni:

who fornisce lelenco di tutte le variabili presenti nella finestra dei comandi (lo si pu vedere anchenella Workspace;whos fornisce lelenco di tutte le variabili insieme allo spazio in memoria da esse occupato;help pu essere usato da solo o insieme al nome di una function di MATLABo creata dallutente, emostra tutte le informazioni utili per capire come usare MATLABo quella function;clear se usata da sola cancella tutte le variabili presenti nella finestra dei comandi, se invece seguita da un elenco di variabili (messe una di seguito allaltra senza virgole) cancella tutte le variabilidellelenco;il punto e virgola ; messo dopo unistruzione non fa vedere il risultato dellistruzione nella finestradei comandi;il simbolo % il simbolo per scrivere commenti: ci che viene scritto dopo % rappresenta uncommento;diary permette di salvare, su un file che viene chiamato diary, il contenuto di ci che viene scrittonella finestra dei comandi. Il file diary si trova nella directory corrente in cui si sta lavorando;

5


diary off chiude il file aperto mediante listruzione diary;diary filediary comando analogo a diary, ma il file non si chiamer diary bens filediary,dove filediary rappresenta il nome del file definito dallutente;save filesave salva tutte le variabili presenti nel Workspace, nel file filesave.mat confilesave nome scelto dallutente;load filesave ripristina lo stato del Workspace, per esempio dopo aver chiuso la sessione diMATLAB , caricando tutte le variabili che erano state salvate in filesave.mat.

Ce una notevole differenza di significato tra i files generati dallistruzione diary e quelli generati me-diante save. Ci che salviamo con il comando diary paragonabile a ci che pu essere scritto alla lavagnae di cui prendiamo appunti: scriviamo un file di testo che possiamo elaborare o da cui possiamo trarre deirisultati per i nostri problemi. Al contrario, quello che viene salvato mediante il comando save costituitoda una o pi variabili, che possono essere riutilizzate nuovamente inMATLABnella finestra dei comandi, insessioni successive a quella in cui sono state create, senza dover rieseguire tutti i comandi che hanno portatoalla creazione delle variabili stesse.

Il comando save utile anche per salvare singolematrici o una lista dimatrici nella directory corrente. Cisoffermiamo sulla possibilit di salvare i dati in formato ASCII, inmodo da poter utilizzare i dati che salviamonon solo in MATLABma anche in programmi scritti con altri linguaggi di programmazione. Supponiamodi avere una matrice A e di volerla salvare nel file matriceA.dat. Scriveremo il comando

save matriceA.dat A -asciiQuando vogliamo caricare in MATLAB la matrice dal file useremo il comandoload(matriceA.dat)e in questo modo avremo la matrice nella variabile matriceA (il nome del file senza lestensione dopo il

punto).

I.4 MATLABper scrivere ed eseguire programmi

A differenza del FORTRAN, MATLABnon ha subroutines e functions ma solo functions che possonoavere uno o pi variabili di ingresso e una o pi variabili in uscita.

I files in cui vengono scritti i programmi si chiamano M-files dal fatto che hanno lestensione .m:prova.m, fun.m, . . . I files vanno scritti con leditor di testo di MATLABo con un altro editor di testo informato ASCII (quello che facciamo per il FORTRAN).

Quello che in FORTRAN corrisponde al programma principale, in MATLABprende il nome di script,per il fatto che presenta una lista di comandi da eseguire in ambiente MATLAB , tra cui la chiamata allefunctions proprie di MATLABo definite dallutente (sempre tramite M-files). Inoltre non bisogna dare nes-suna istruzione di inizio e fine programma (cio non esiste lequivalente delle istruzioni program e stop,end con cui apriamo e chiudiamo un programma principale in FORTRAN). importante vedere subito unaimportante differenza tra FORTRAN e MATLAB : mentre in FORTRAN, siamo abituati a scrivere un pro-gramma in un unico file, includendo in esso anche le subroutines e le functions che ci servono, poich poidobbiamo compilare per avere il programma eseguibile, in MATLAB le functions e gli scripts vengonosempre scritti su files distinti.

Per eseguire uno script salvato in nomefile.m non c bisogno di dover compilare il programma (co-me invece si deve fare per il FORTRAN): basta lanciare il programma nella finestra dei comandi mediante ilnome del file (senza lestensione .m) o mediante il comando run(nomefile.m) o run(nomefile).Una function invece pu essere eseguita tramite uno script oppure direttamente nella finestra deicomandi, specificando le variabili di ingresso e di uscita, come vedremo con qualche esempio tra breve.

Per quanto riguarda il nome da dare alle variabili, MATLAB sensibile alle lettere minuscole emaiuscole: una variabile chiamata A diversa da a, una variabile XOLD diversa da xold.

importante osservare, inoltre, che nello scrivere M-files le variabili utilizzate non vanno dichiarate,come invece va fatto in FORTRAN.

I.4.1 Strutture

Per capire le differenze di programmazione tra FORTRAN eMATLAB , partiamo dagli operatori logici.

6

I.4. MATLABper scrivere ed eseguire programmi

Operatore MATLAB Operatore FORTRAN Significato< .LT minore> .GT. maggiore= .GE. maggiore o uguale== .EQ. uguale= .NE. non uguale&& .AND. congiunzione|| .OR. disgiunzione .NOT. negazione

Tabella I.1: Operatori relazionali e logici

Vediamo poi come si scrivono i vari cicli in FORTRAN e in MATLAB :

in FORTRAN

i f ( espr . logica ) then{ is truzione 1a }{ istruzione 2a }{ . . . . }

else{ is truzione 1b }{ istruzione 2b }{ . . . . }

end i f

in MATLAB

i f ( espr . logica ){ istruzione 1a }{ istruzione 2a }{ . . . . }

else{ is truzione 1b }{ istruzione 2b }{ . . . . }

end

in FORTRAN

i f ( espr . logica1 ) then{ is truzione 1a }{ istruzione 2a }{ . . . . }

else i f ( espr . logica2 ) then{ is truzione 1b }{ istruzione 2b }{ . . . . }

. . . .else

{ is truzione 1z }{ istruzione 2z }{ . . . . }

end i f

in MATLAB

i f ( espr . logica1 ){ istruzione 1a }{ istruzione 2a }{ . . . . }

e l se i f ( espr . logica2 ){ istruzione 1b }{ istruzione 2b }{ . . . . }

. . . .else

{ is truzione 1z }{ istruzione 2z }{ . . . . }

end

Osserviamo che in MATLABnon ci sono i then e che elseif va scritto tutto attaccato.

7


in FORTRAN

do while ( espressione logica ){ istruzione 1 }{ istruzione 2 }{ . . . }{ is truzione n }

end do

in MATLAB

while ( espressione logica ){ istruzione 1 }{ istruzione 2 }{ . . . }{ is truzione n }

end

Nel ciclo while si elimina il do.

in FORTRAN

do ind= val iniz , val f in , incr{ i s t ruz ion i }

end do

in MATLAB

for ind= va l in iz : incr : va l f in{ i s t ruz ion i }

end

Il ciclo do in FORTRAN diventa ciclo for in MATLAB . Lordine delle variabili relative al valore inizialedel ciclo, valore finale e incremento non pi come in FORTRAN e ci sono i due punti al posto della virgola:osservare bene le differenze.

C unaltra struttura che si trova in MATLABe non in FORTRAN 77: la costruzione switch-case.

switch ( espressione ) % ( scalare o str inga )case { valore1 } % eseguita se l espress ione e valutata al valore1 ){ i s t ruz ion i }{ . . . }case { valore2 } % ( eseguita se l espress ione e valutata al valore2 ){ i s t ruz ion i }{ . . . }

otherwise{ i s t ruz ion i }{ . . . }

end

Il ciclo con switch confronta i valori dati nellespressione di input (subito dopo switch) con ciascunvalore assegnato a case ed esegue le istruzioni relative al case in cui valore ed espressione coincidono.

Nellesempio che riportiamo, a seconda del valore assegnato alla variabile scelta cambiano i valori daassegnare alle variabili a e b:sce l ta=test1 ;switch sce l tacase {test1}x0= 0 . 1 ;x1= 0 . 2 ;case {test2}x0= 0 . 0 ;x1= 1 . 0 ;otherwisedisp (nessuncasotestscelto)end

Osserviamo che scelta una variabile di stringa di caratteri, il nome del caso test scritto tra apici; pervisualizzare unmessaggio sulla CommandWindow, abbiamo usato la function di MATLABchiamata disp.

8

I.5. Dati di input

I.5 Dati di input

Quando si lavora nella finestra dei comandi, per assegnare il valore alle variabili basta scrivere il nomedella variabile seguito dal simbolo di uguale e dal valore (o dai valori) da assegnare (a seconda che si tratti diuna variabile scalare, matrice, vettore...):

>> a=10.5;>> A=[1 2 3; 4 5 6; 7 8 9];

Nelleseguire uno script, alla stessa maniera, si pu assegnare il valore alle variabili direttamente allinternodello script. Tuttavia, se si vuole dare maggiore generalit al programma e si vogliono dare in input i valoridelle variabili, conviene usare la function input. Vediamo con un esempio:a=input(scriviilvaloredellavariabilea)Il messaggio contenuto tra apici viene visualizzato sulla finestra dei comandi e il prompt aspetter che lu-tente scriva il valore da assegnare ad a. La function input si pu dunque paragonare alle istruzioni di writee read in FORTRAN per assegnare il valore alle variabili.

Questa procedura pu essere utilizzata sia per assegnare il valore a variabili scalari, sia permatrici e vetto-ri. Tuttavia, se i dati di input sonomolto pesanti (ad esempiomatrici di dimensionimolto elevate), convieneusare inmaniera opportuna la funzione di input unitamente al comando load - scrivendo una volta per tut-te la matrice di input in un file da caricare ogni volta che si vuole eseguire il programma con quella matrice.Per esempio, abbiamo scritto nel file A.dat i valori della matrice e vogliamo dare in input questa matrice alnostro programma. Invece di scrivereA= input(matriceA);nello script scriveremoA= input(filediinputconlamatriceA, s);A= load(A);

Osserviamo che abbiamo usato input inserendo, oltre alla stringa tra apici file di input con lamatrice A, anche s: questa opzione serve perch noi scriveremo il nome del file su cui si trova memo-rizzata la matrice e questo file viene letto come una stringa di caratteri per cui inizialmente alla variabile Aviene associato il file con la matrice (in questo caso A.dat). Con il comando successivo, viene caricato il filee memorizzato di nuovo sulla variabile A, che quindi diventa la nostra matrice.

Si osservi la differenza che c nel lasciare o meno spazi bianchi prima dellapice che chiude la frase cheviene visualizzata sulla CommandWindow tramite input.

utile sapere anche che, se unistruzione troppo lunga e si vuole andare a capo, si utilizzano tre puntini... sulla riga che si vuole interrompere e si prosegue a scrivere sulla riga successiva.

I.5.1 Programma sulmetodo delle bisezioni

Consideriamo lalgoritmo del metodo delle bisezioni, per calcolare gli zeri di unopportuna funzione f .

Ad esempio vogliamo risolvere il problema f (x)= 0 con f (x)= (x2)2 sin(x) nellintervallo [1,3].

Dobbiamo scrivere il programma principale, in un file che chiamiamo bisez.m, e la function legata allafunzione f (x), in un file che chiamiamo fun.m.

Scriviamo il programma principale:

% programma per i l calcolo degl i z e r i di una funzione mediante% i l metodo di bisezionea=input (primoestremodellintervalloa ) ;b=input (secondoestremodellintervallob ) ;aux=fun (a ) * fun (b ) ;i f aux>=0

disp (estremidellostessosegno)break

enditmax=100;t o l l =1.e10;

9


i t e r =0;c=(a+b ) * 0 . 5 ;scarto=abs (ba ) * 0 . 5 ;while i t e r t o l l

i t e r = i t e r +1;aux=fun (a ) * fun ( c ) ;i f aux>0

a=c ;else

b=c ;end

c=(a+b ) * 0 . 5 ;scarto=abs (ba ) * 0 . 5 ;

endi f fun ( c)==0 | | scarto < t o l l

sprintf (%s%15.8e , soluzioneapprossimatac= , c )else

sprintf (%s , raggiuntonumeromaxdiiterazioni)end

Per poter essere eseguito, dobbiamo scrivere in un file chiamato fun.m la function fun. Si ha:function y=fun ( x )% funzione per lo schema de l l e bi sez ioni% input x% output y=(x /2)^2 sin ( x )y=(x/2).^2 sin ( x ) ;

Nel programma principale abbiamo semplicemente tradotto lalgoritmo del metodo delle bisezioni. Ab-biamo considerato dei comandi che non abbiamo ancora visto inMATLAB , per la stampa dei risultati (me-diante sprintf) e abbiamo usato il comando break per interrompere lesecuzione del programma se unacondizione non verificata.

Le righe di commento scritte dopo listruzione function y=fun(x) vengono visualizzate sullaCommandWindow se, una volta salvata la function, digitiamo listruzione help fun

I.6 Dati di output

Lambiente MATLABpermette di visualizzare tutte le variabili presenti al termine dellesecuzione di unprogramma o di una funzione. Tuttavia, alcune volte vogliamo visualizzare subito dei risultati sulla Com-mand Window o su un file separato rispetto allambiente di lavoro. Abbiamo gi visto, nellesempio prece-dente il comando sprintf: esso serve per stampare un messaggio sulla Command Window utilizzando uncerto tipo di formato.

Nellesempio avevamo usatosprintf(%s%15.8e, soluzioneapprossimatac=, c)

Allinterno delle parentesi tonde si possono riconoscere due parti: nella primaparte si scrive, tra apici e unodiseguito allaltro, il formato da utilizzare per la stampa delle variabili che sono scritte subito dopo. Nel nostrocaso, le variabili da stampare sono la stringa di caratteri soluzione approssimata c= e il valore dellavariabile c. Per la stringa di caratteri si usa il formato delle stringhe che dato da %smentre per la variabile c,volendo scriverla mediante un formato esponenziale con 8 cifre nella mantissa, usiamo il formato %15.8e.

Rispetto al FORTRAN, il formato per una variabile preceduto dal simbolo % ma per il resto non ci sonogrosse differenze. Per la variabili intere o reali possiamo decidere quante cifre utilizzare per rappresentare lavariabile (%15.8e nellesempio). Riassumiamo nella tabella seguente i principali tipi di formato e i simboliche servono per andare a capo o creare nuove linee: Se si vogliono scrivere i risultati su un file, occorre aprireil file e associarlo ad una variabile mediante la function fopen. Ad esempiofid= fopen(risultati.txt,w)Con questa istruzione aprimamo il file di risultati dal nome risultati.txt (w indica che il file di

10

I.7. Grafici

formato Significato%s stringhe di caratteri%d formato intero%f formato fisso%e formato esponenziale (del tipo 3.5e+00)%E formato esponenziale (del tipo 3.5E +00)\ n nuova linea\ r per andare a capo

Tabella I.2: Il formato

Figura I.3: Finestra con il grafico

scrittura) associandolo alla variabile fid. Per scrivere sul file, al posto della function sprintf si user lafunction fprintf che differisce dalla prima per il fatto che bisogna indicare la variabile associata al file discrittura dati. Il comando sprintf di prima diventa ora:fprintf(fid , %s%15.8e, soluzioneapprossimatac=, c)Per chiudere il file si usa listruzione fclose(fid). Quindi fopen e fclose sono equivalenti alle

istruzioni open e close che abbiamo visto in FORTRAN.

I.7 Grafici

Supponiamo di voler fare il grafico di una serie di dati (xi , yi ), i = 1, . . . ,n. Sulla Command Window (oallinterno di uno script) basta digitare il comandoplot(x,y)

Si aprir una nuova finestra con il grafico (vedi Figura I.3).Potremopoimodificare il colore, il tipo di linea, inserire titolo, legenda,...operandodirettamente sulmenu

della finestra del grafico, o inserendo i comandi opportuni tramite la CommandWindow1.Si possono sovrascrivere grafici luno sullaltro utilizzando il comando hold on. Oppure si possono

affiancare grafici mediante il comando subplot. Lasciamo gli approfondimenti allhelp on line.Per fare il grafico di una funzione, si possono seguire diverse strade.

1Ci sarebbe tanto da dire a riguardo ma lasciamo che il lettore curioso approfondisca largomento utilizzando lhelp on line diMATLAB . In Octave, invece, le modifiche ai grafici non possono essere fatte usando la finestra del grafico.

11


Se si ha a disposizione la function (propria di MATLABo scritta su un file .m), si pu costruire il vettorecon il valore della funzione in un numero determinato di punti equidistanti sullintervallo in cui si desideravisualizzarla. A tal proposito utile la function linspace che permette di discretizzare un intervallo chiuso[a,b] in un prefissato numero di punti. Useremo allora le seguenti istruzioni:

x=linspace ( 0 , 2 ) ;% discretizziamo l in te rva l lo [0 ,2 ] in 100 part i% uguali% x=linspace (50 ,0 ,2) d i s c r e t i zza l in te rva l lo in%50 part i uguali

y=myfun( x ) ; % valuto la function myfun nel ve t tore xplot ( x , y )

Questo approccio comodo quando la funzione da visualizzare ha una espressione complicata o quandostiamo lavorando allinterno di uno script. La function myfun deve essere scritta in modo che sia possibilevalutarla direttamente su un vettore (che quello che facciamo tramite listruzione y=myfun(x)). Le opera-zioni dimoltiplicazione, divisione ed elevamento a potenza devono essere vettorizzate, facendo precendere ilsimbolo dimoltiplicazione, divisione e elevamento a potenza dal simbolo del punto ( .*, ./, . ) permet-tendo, in tal modo, che le operazioni vengano fatte componente per componente del vettore. Le operazionidi somma e differenza sono vettorizzate per definizione.

Ad esempio: la function myfun definita tramite le istruzionifunction [ y ]= f ( x )y=log ( x ) * x

non vettorizzata. Se invece scriviamo

function [ y ]= f ( x )y=log ( x ) . * x

allora la function vettorizzata.Unaltra via di scrivere la function come una funzione inline nella Command Window e di utilizzare la

function ezplot per farne il grafico.Supponiamo di voler fare il grafico della funzione f (x)= ex10sin(x)1. Per avere questa funzione nella

Command Window durante la sessione di lavoro (quindi non come una function scritta su file), scriviamo ilcomandofun = inline(exp(x)-10*sin(x)-1)Sulla CommandWindow compaiono le seguenti righe:

>> fun=inline(exp(x)-10*sin(x)-1)

fun =

Inline function:fun(x) = exp(x)-10*sin(x)-1

Se vogliamo farne il grafico nellintervallo [0,2] scriveremoezplot(fun, 0, 2)e viene direttamente creato il grafico della funzione fun. La function ezplot pu essere utilizzata anchecon funzioni intrinseche di MATLABo definite dallutente in forma vettorizzata (per esempioezplot(sin,0,2)crea il grafico della funzione sin(x) nellintervallo [0,2]).

I.8 Sulle potenzialit di MATLAB

MATLABha un grande numero di functions predefinite che permettono di eseguire applicazioni in di-versi settori propriamente matematici, per approssimare zeri di funzioni, per lavorare su matrici sparse, ri-solvere equazioni differenziali, lavorare su dati statistici, fare grafici in due e tre dimensioni... Basta vedere ilDemo di MATLABper trovare la function di cui si ha bisogno.

12

I.9. Applicazioni di MATLABnel Calcolo Numerico

possibile inoltre installare pacchetti specifichi per specifiche aree di applicazioni, quali sistemi dicontrollo, biologia computazionale, finanza computazionale, meccatronica...

In questa breve presentazione di MATLAB , tuttavia, cercheremo di vedere come utilizzarlo per capiremeglio alcuni dei problemi trattati in CalcoloNumerico, rielaborando o ripresentando alcuni esempi gi visti.

I.9 Applicazioni di MATLABnel Calcolo Numerico

I.9.1 Sullinstabilit numerica

Riprendiamo lesempio sullinstabilit numerica visto nel Capitolo sulla rappresentazione dei numeri alcalcolatore, per cui vogliamo approssimare il valore degli integrali espressi mediante la formula

yn =Z10

xn

x+10 dx

per valori di n = 1,2, . . . ,30.Avevamo visto due formule ricorsive che ci permettevano di approssimare yn , una era instabile e laltra

era stabile.Proviamo a rivedere gli algoritmi e scriviamo un programmaMATLABda eseguire per verificare quanto

avevamo detto.

% Eserc iz io sul l in s tab i l i ta numerica% calcolo del l in tegra le y_n= int_0^1 x^n / ( x+10) dx%% yins t : ve t tore con i valor i del l algoritmo ins tab i l e% ys t : ve t tore con i valor i del l algoritmo s tab i l e%% algoritmo ins tab i l e%yinst (1)= log (11) log ( 10 ) ; %corrisponde al valore in i z ia l efor i =1:30

yinst ( i +1)= 1/ i 10*yinst ( i ) ;end%% algoritmo s tab i l e% s i r ichiede che i l valore del l in tegra le y_n1 sia approssimato% con una accuratezza data dal valore di input to ln1= input (indicen1 ) ;to l=input (tolleranzatol ) ;k= log10 ( to l ) + n1 ;k= f i x ( k+1) ; % f ix e una function che e f f e t tua l arrotondamento

% del numero in modo da avere un valore interoyst ( k )=0 ;for j =k1:1:1

yst ( j ) =1/10*(1/ j yst ( j +1 ) ) ;end

Uno volta eseuito lo script, nella CommandWindow si hanno i due vettori che possono essere confrontati traloro. Osserviamo che lequivalente programma FORTRAN sarebbe meglio scriverlo non utilizzando i vettorie facendosi stampare i risultati intermedi su un file. In questo caso, potendo avere i risultati sulla CommandWindow, ci conviene usare dei vettori. Il valore iniziale y0 si avr nella prima componente dei vettori chevengono creati. Perci si faccia attenzione agli indici utilizzati (per yst si usa j e j+1: perch?).

Volendo, si pu modificare lo script facendo uso della function di MATLABsingle che converte ilrisultato in singola precisione in modo da confrontare i due algoritmi con i calcoli in singola precisione.

% yinst s ing : ve t tore del l algoritmo ins tab i l e lavorando% in singola precis ione

13


Figura I.4: Algoritmo instabile: schermata del confronto tra luso o meno della function single

% yst s ing : ve t tore del l algoritmo s tab i l e lavorando% in singola precis ioneyinsts ing (1)= single ( log (11 ) ) s ingle ( log ( 1 0 ) ) ;for i =1:30

yinsts ing ( i +1)= single (1/ i ) s ingle (10* yinsts ing ( i ) ) ;endn1= input (indicen1 ) ;to l=input (tolleranzatol ) ;k= log10 ( to l ) + n1 ;k= f i x ( k+1) ;ysts ing (k )=0;for j =k1:1:1ysts ing ( j ) =s ingle (1/10)* s ingle (1/ j ysts ing ( j +1 ) ) ;end

Come si pu osservare dalle Figure I.4 e I.5, i risultati ottenuti dallalgoritmo instabile cambiano a secondache si usi o meno la function single mentre abbiamo gli stessi risultati (consideriamo le cifre corrette insingola precisione, usando il formato format short e) per lalgoritmo stabile.

I.9.2 Sullinterpolazione e approssimazione di dati

In MATLABesistono gi delle function che permettono di interpolare e approssimare delle serie di dati.C la function polyfit che, dati i vettori contenenti le ascisse e le ordinate da interpolare o approssima-

re, di dimensione n, e il gradom del polinomio che si vuole creare, fornisce in output il vettore contenente icoefficienti del polinomio interpolante (o approssimante) in ordine decrescente [am am1 . . . a0] da cui il po-linomio p(x)= amxm+am1xm1+ . . .a0: perm = n1 si ha il polinomio di interpolazione, perm < n1 si

14


Figura I.5: Algoritmo stabile: schermata del confronto tra luso o meno della function single

ha il polinomio di approssimazione. Lalgoritmo si basa sul processo di minimizzazione nel senso dei minimiquadrati.

Esempio:

>> x=[7 8 9 10];>> y=[3 1 1 9];>> p=polyfit(x,y,3)

p =

1.0000 -23.0000 174.0000 -431.0000Significa che il polinomio di interpolazione :

p(x)= x323x2+174x431Una volta ricavati i coefficienti, si pu fare un grafico del polinomio utilizzando la function polyval.

>> xx=linspace(x(1), x(4));>> yy=polyval(p,xx);>> plot(x,y,o, xx,yy)

Con polyval si valuta il polinomio, i cui coefficienti sono dati dal vettore p, nei punti di xx. Abbiamo usatola function plot per rappresentare sullo stesso grafico due curve, quella dei dati x,y (grafico che facciamoper punti utilizzando dei cerchietti) e quella del polinomio.

15


Scriviamo ora delle function che ci permettano di ottenere il polinomio di interpolazione sia usandolapproccio delle funzioni base monomiali che porta alla costruzione della matrice di Vandermonde, siacostruendo i polinomi di Lagrange o utilizzando le differenze divise di Newton.

Usando le funzioni base monomiali, scriviamo la seguente function, interpmonom:function p=interpmonom(x , y )% function p=interpmonom(x , y )% interpolazione monomiale% dati i valor i x e y da interpolare s i cos t ru i s ce i l ve t tore p% dei c o e f f i c i e n t i del polinomio di interpolazione% applicando i l metodo dei c o e f f i c i e n t i indeterminati%% se x e y non sono gia v e t t o r i colonna l i rendiamo ta l i% mediante l e due i s t ruz ioni succes s ivex=x ( : ) ;y=y ( : ) ;i f length ( x)~= length ( y )% length e una function che misura la lunghezza del ve t tore% ( s i confronti la di f ferenza tra length e s i z e )

error (MATLAB:interpmonom , . . .ivettorixeynonhannolastessalunghezza)

elsen=length ( x )1;% V matrice di Vandermonde cos t rui ta in maniera r i cor s ivaV( : , 1 )= ones (n+1 , 1 ) ;for i =2:n+1

V( : , i )= x . *V ( : , i 1);endp=V\y ;% i l ve t tore p contiene i c o e f f i c i e n t i del polinomio interpolatore% in ordine crescente p0 p1 p2 . . .% se vogliamo usare la function del MATLAB polyval per valutare% ta le polinomio in piu punti , dobbiamo s c r i v e r l i in ordine decrescentefor i =1:n+1aux ( i )=p(n+2 i ) ;endp=aux ;end

Osserviamo che abbiamo usato listruzione error per mostrare un messaggio di errore e far interromperelesecuzione della function, nel caso in cui i dati di input x e y non abbiano la stessa lunghezza. La stringaMATLAB:interpmonom una stringa di identificazione dellerrore (pu essere anche omessa), mentre lastringa i vettori x e y non hanno la stessa lunghezza quella che viene visualizzata durante lesecuzione delcodice.

La matrice V stata costruita in maniera ricorsiva. Una volta calcolato il vettore p possiamo valutare ilpolinomio di interpolazione mediante la polyval.

Riprendendo lesempio di prima, con x,y,xx,yy gi dati:

>> p=interpmonom(x,y);>> plot(x,y,o, xx,yy)

Per quanto riguarda linterpolazione di Lagrange, si considerino le due functions che chiamamolagrange e interplagrange rispettivamente. La prima valuta li -simo polinomio di Lagrange e laltravaluta il polinomio di interpolazione di Lagrange in un assegnato punto (o nelle componenti di un vettore).

function yval=lagrange ( xval , x , i )% function yval=lagrange ( xval , x , i )% function che calcola i l polinomio isimo di Lagrange

16


% valutandolo in xval% xval puo e s s e re uno scalare o un vet tore% x vet tore de l l e a s c i s s e da interpolarexval=xval ( : ) ;n=length ( x ) ;yval =ones ( length ( xval ) , 1 ) ; % s i crea un vet tore di t u t t i 1for j =1:n

i f j ~= iyval=yval . * ( xvalx ( j ) ) / ( x ( i ) x ( j ) ) ;

endend

function yval=interplagrange ( xval , x , y )% function yval=interplagrange ( xval , x , y )% dati i v e t t o r i x e y da interpolare% la function implementa l interpolazione di Lagrange valutandola% in xval% xval puo e s s e re uno scalare o un vet tore% questa function chiama la function lagrange ( xval , x , i )i f length ( x)~= length ( y )

error (MATLAB:interplagrange , . . .ivettorixeynonhannolastessalunghezza)

elsexval=xval ( : ) ;yval=zeros ( length ( xval ) , 1 ) ; % s i crea un vet tore di t u t t i zeron=length ( x )1; % n grado del polinomiofor i =1:n+1

yval = yval + lagrange ( xval , x , i ) * y ( i ) ;endend

Eseguiamo lesempio di prima per vedere cosa si ottiene:

>> yy=interplagrange(xx,x,y);>> plot(x,y,o, xx,yy)

Ora la function interplagrange sostituisce luso delle due function polyfit, polyval ointerpmonom, polyval.

Calcoliamo ora il polinomio di interpolazione mediate le differenze divise di Newton. Scriviamo duefunctions, la prima che scrive la tabella delle differenze divise, la seconda che valuta il polinomio di in-terpolazione implementando lalgoritmo di Horner2 in modo da minimizzare il numero delle operazioni daeseguire.

function table=divdi f ( x , y )% function table=divdi f (x , y )% x asc i s s e dei dati da interpolare% y ordinate dei dati da interpolare% table tabel la de l l e d i f f e renze div i sex=x ( : ) ;y=y ( : ) ;n=length ( x ) ;m=length ( x ) ;i f n~=m

error (MATLAB:differenze_divise ,erroresuidati)else

2William Horner (1786-1837) fu un matematico inglese, ricordato essenzialmente per il suo metodo sulle equazioni algebriche.Spieghiamo lalgoritmo solo per rendere comprensibile la function che scriviamo.

17


table=zeros (n ,n ) ; % iniziamo la tabel la come una matrice di ze r itable ( : , 1 )= y ;for j =2:n

for k=2: jtable ( j , k)= ( table ( j , k1) table ( j 1,k1) )/ . . .

( x ( j ) x ( jk+1) ) ;end

endend

La diagonale principale della matrice table ha i coefficienti a0, a1, . . .an del polinomio di interpolazione.Il polinomio dunque:

p(x)= a0+a1(xx0)+a2(xx0)(xx1)+ . . .+an(xx0)(xx1) (xxn1)

Per valutarlo in un punto x eseguiamo i seguenti passaggi, applicando, in tal modo, lalgoritmo di Horner:

p = anp = p(xxn1)+an1= an(xxn1)+an1

p = p(xxn2)+an2= an(xxn1)(xxn2)+an1(xxn2)+an2...

p = p(xx0)+a0= an(xxn1)(xxn2) . . . (xx1)(xx0)+ . . .+a1(xx1)(xx0)+a0= a0+a1(xx0)(xx1)+ . . .+an(xx0)(xx1) . . . (xxn2)(xxn1)

La function da scrivere dunque:

function yval= interpdivdi f ( xval , x , table )% function yval=interpdivdi f ( xval , x , table )% x asc i s s e dei dati da interpolare% table tabel la de l l e d i f f e renze divise ,% ottenuta dalla function divdi f ;% servono i valor i del la diagonale principale% xval scalare o ve t tore in cui calcolare i l polinomio% interpolatore% yval valore (o valor i ) del polinomio interpolatore% valutato in xval% nel calcolo del polinomio interpolatore s i applica% l algoritmo di Hornerxval=xval ( : ) ;x=x ( : ) ;n=length ( x )1;yval=table (n+1 ,n+1)*ones ( length ( xval ) , 1 ) ;for j =n:1:1

yval=yval . * ( xvalx ( j ) ) + table ( j , j ) ;end

Le functions appena descritte possono essere usate inmaniera del tutto equivalenti per risolvere il proble-ma dellinterpolazione. Ci sono per dei casi in cui i risultati ottenuti dalla polyfit e dalla interpmonomnon sono corretti in quanto la matrice di Vandermonde ad esse legate malcondizionata.

Abbiamo gi descritto gli effetti del malcondizionamento nellinterpolazione. Quando applichiamo lafunction polyfit ai dati di quel problema, si ha il seguente messaggio di avvertimento (un warning):

18


>> poli=polyfit(x,y,5)Warning: Polynomial is badly conditioned. Add points with distinct X

values, reduce the degree of the polynomial, or trycenteringand scaling as described in HELP POLYFIT.

Questo ci dice che i risultati che avremo non saranno buoni e, effettivamente, se facciamo il grafico dei datidel problema e del polinomio di interpolazione ottenuto con la polyfit, si nota subito che il polinomio completamente errato (vedi Figura I.8) Anche nellapplicare la function interpmonom si ha il messaggio di

Figura I.6: Uso della function polyfit per lesempio malcondizionato.

avvertimento

>> pmon=interpmonom(x,y)Warning: Matrix is close to singular or badly scaled.

Results may be inaccurate. RCOND = 5.537667e-31.

Se aggiungiamo al grafico precendente la curva corrispondente al polinomio ottenuto con la interpmonomosserviamo come le due curve siano tra loro distinte e non interpolano i dati.

Figura I.7: Uso della function interpmonom per lesempio malcondizionato.

Proviamo invece ad applicare lalgoritmo di Lagrange o delle differenze divise di Newton e aggiungiamole nuove curve sul precedente grafico:

>> ylagr=interplagrange(xx,x,y);>> table=divdif(x,y);

19


>> ynewt=interpdivdif(xx,x,table);>> plot(xx,ylagr,xx,ynewt)

Figura I.8: Uso delle functions corrispondenti agli algoritmi di Lagrange e delle differenze divise di Newtonnellesempio malcondizionato. Osserviamo che la Figura ?? relativa allo stesso problema stata ottenutaeseguendo le stesse functions (per semplicit abbiamo omesso i risultati ottenuti dalla interpmonom) inambiente Octave.

Si nota subito che le due curve sono tra loro coincidenti e interpolatorie!

I.9.3 Confronto tra schemi per equazioni differenziali ordinarie

Quando abbiamo introdotto alcuni metodi per la risoluzione di equazioni differenziali ordinarie, inparticolare i metodi di Eulero esplicito, di Eulero implicito e di Crank-Nicolson.

Li abbiamo anche confrontati per capire meglio il concetto di stabilit, considerando lequazione testy 0 =y con y(0)= 1 (si veda a proposito la Figura ??).

Scriviamo ora un programmaMATLABche, per questa equazione test,ci permetta di poter confrontare i tre metodi tra loro e con la soluzione esattacrei dei grafici simili a quelli di Figura ??salvi i valori numerici dei diversi metodi e della soluzione esatta in un file.

clear% h = passo di discre t izzazione temporale% Tfin = tempo f inale di osservazione% n = numero di vol te in cui verra applicato ciascuno dei metodi% proposti% f id = variabi le associata al f i l e di r i s u l t a t i% t = vet tore dei tempi% ye = vet tore dei valor i a s soc ia t i al metodo di Eulero e sp l i c i t o% ( di lunghezza n+1 perche la prima componente corrisponde% al valore del la soluzione al tempo zero )% yi = vet tore dei valor i a s soc ia t i al metodo di Eulero implici to% ycn = vet tore dei valor i a s soc ia t i al metodo di CrankNicolson% yex = vet tore dei valor i a s soc ia t i al la soluzione esat tah=input (passoh ) ;Tfin=input (tempofinalediosservazione ) ;n=Tfin /h ;f id=fopen (ode_a_confronto.txt ,w ) ;t (1 )=0 ;y (1)=1 ;ye (1)=y ( 1 ) ;

20


yi (1)=y ( 1 ) ;ycn(1)=y ( 1 ) ;yex (1)=y ( 1 ) ;for i =1:n

ye ( i +1)= ye ( i ) h*ye ( i ) ;y i ( i +1)= yi ( i )/(1+h) ;ycn ( i +1)= (2h)/(2+h)* ycn ( i ) ;t ( i +1)= t ( i )+h ;yex ( i +1)=exp( t ( i +1 ) ) ;

endplot ( t , yex , k , linewidth ,2 )hold on %comando per sovrascr ivere i g ra f i c i l uno sul l a l t roplot ( t , ye ,b , linewidth ,2 )plot ( t , yi ,r ,linewidth ,2 )plot ( t , ycn , g ,linewidth ,2 )legend (soluzioneesatta , Euleroesplicito , Euleroimplicito , Crank-Nicolson , )% la legenda , lo spessore de l l e l inee , o i l colore possono es se re% messi s ia direttamente dalla f ine s t ra del la figura% sia dal programma stesso , come in questo caso ( solo in Octave s i% deve operare direttamente dal programma per poter cambiare% le proprieta del la figura ) .hold of ffor i =1:n+1fpr int f ( f id , \n%5.2f%12.6e%12.6e%12.6e%12.6e , t ( i ) , yex ( i ) , . . .

ye ( i ) , y i ( i ) , ycn ( i ) ) ;endfclose ( f id ) ;

Questo programma specifico per lequazione test assegnata. La sua esecuzioneci permette ci comprenderemeglio il concetto di stabilit dei metodi studiati per la soluzione di equazioni differenziali ordinarie.

21

APPEN

DICE II

TABELLE E GRAFICI CON IL FOGLIO ELETTRONICO

Io dico di aver capito unequazionequando sono in grado di predire lepropriet delle sue soluzioni senzaeffettivamente risolverla.

Paul Adrien Maurice Dirac

II.1 Il foglio elettronico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23II.2 Inserire funzioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23II.3 Formule ricorsive . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28II.4 Fare grafici nel foglio elettronico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31

II.4.1 Grafici multicolonna . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36II.5 Inserire dati da files . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

II.1 Il foglio elettronico

Il foglio elettronico pu essere usato nellambito del Calcolo Numerico per implementare ed eseguiresemplici algoritmi ricorsivi. Inoltre pu essere utilizzato per importare dati e per fare grafici. Nel seguito,daremo alcune indicazioni su come usare il foglio elettronico, in particolare la versione in lingua inglese diOpenOffice 2.4. Osserviamo subito che nella versione inglese i numeri decimali vanno scritti facendo uso delpunto . ad esempio 2.1023 a differenza della versione in lingua italiana dove i numeri decimali, di default,vanno scritti facendo uso della virgola ad esempio 2,1023.1 Tenendo conto di questa differenza saremo ingrado di usare correttamente il foglio elettronico.

II.2 Inserire funzioni

Un foglio elettronico si presenta come una tabella di cellette, ognuna delle quale viene individuata dalnumero che contraddistingue la corrispondente riga e dalla lettera che indica la corrispondente colonna. Adesempio, la cella C7 si trova sulla colonna C e sulla riga 7.

1 Se vogliamo usare la versione italiana di OpenOffice, si pu scegliere di scrivere i numeri decimali usando il punto . anzich lavirgola , andando, mediante il Menu, su Strumenti > Opzioni > Impostazioni della Lingua > Lingue e scegliere come LinguaLocale (nella seconda riga a partire dallalto) quella Italiana (Svizzera) o unaltra lingua che permette di ottenere come separatore dellecifre decimali il punto.

23

II. TABELLE E GRAFICI CON IL FOGLIO ELETTRONICO

Figura II.1: Si vuole lavorare sulla cella C7. La lettera C e il numero 7 che individuano la colonna C e la riga7 sono evidenziate rispetto alle altre colonne e righe. La cella C7 ha il bordo pi scuro rispetto alle altrecelle. Sono state cerchiate in rosso le due righe chiamate Name Box (Casella del Nome) e Input Line (Riga didigitazione).

Quando, con il mouse, ci spostiamo su una cella e clicchiamo il tasto sinistro del mouse, allora il contor-no della cella diventa di colore nero grassetto: vuol dire che possiamo lavorare su quella cella (si veda laFiguraII.1). Inoltre, sopra la tabella delle celle, sulla sinistra, viene visualizzato nella Name Box (nella Caselladel Nome) il nome della cella su cui vogliamo lavorare.

Se in una cella scriviamo un numero o applichiamo una formula numerica che avr come risultato unnumero, la cella pu essere considerata come una variabile e sar chiamata con la coppia lettera-numero chela individua (per esempio la cella C7).

Per scrivere un numero o applicare una formula numerica, una strada da seguire la seguente:Una volta che siamo sulla cella (per esempio la cella C7 evidenziata rispetto alle altre) scriviamo =e poi il numero che vogliamo scrivere, ad esempio 0.1 (oppure 0,1 se la versione di Office quellaitaliana). Vediamo che quello che scriviamo sulla cella compare anche nella riga di Input Line (o rigaper immettere i comandi, Riga di digitazione) che si trova sopra la tabella delle celle. Nel momento incui scriviamo =, accanto a quella riga compare un simbolo a forma di x dal colore rosso (al posto delsimbolo di sommatoria che cera prima) e il simbolo di accettazione in verde (vedasi le immagini inFigura II.2 e II.3 e le si confronti con la FiguraII.1).Se invece di un numero vogliamo applicare una formula numerica, possiamo o scrivere direttamente ilvalore della funzione (per esempio cos(0.1) + 0.3 ) oppure possiamo avvalerci dellaiuto delle funzionipredefinite che possiamo utilizzare nel foglio elettronico. Dal menu del foglio elettronico, cliccandosu Insert (Inserisci) troviamo le due voci Function. . . (Funzione) e Function list (Lista di funzioni). Diqui possiamo vedere quali sono le funzioni che possiamo utilizzare. Ricordiamo che nella versioneitaliana, la funzione valore assoluto si chiama ASS (ABS in quella inglese), la funzione seno si chiamaSEN (SIN in quella inglese). Nel momento in cui iniziamo a digitare il nome di una funzione, moltoprobabilmente ci verr gi data indicazione sulla funzione che inizia in quel modo, se labbiamo giutilizzata precedentemente. Inoltre, nel Name Box, dove prima era individuata la cella su cui stavamolavorando, vengono visualizzate le funzioni pi recenti che abbiamo utilizzato (si vedano le Figure II.4e II.5).

24

II.2. Inserire funzioni

Figura II.2: Scriviamo = nella cella C7 per inserire un numero o per valutare una funzione in un numero.

Figura II.3: Inseriamo il valore 0.1 nella cella C7.

25


Figura II.4: Inseriamo il valore cos(0.1) nella cella C7 partendo dal menu Insert (Inserisci) e Function(Funzione).

26

II.2. Inserire funzioni

Figura II.5: Inseriamo il valore cos(0.1) nella cella C7. Nella figura in alto si vede lultimo passaggio di cosaaccade dopo aver inserito la funzione da menu. Nella figura in basso si vede cosa succede quando si scrivedirettamente la formula.

27


Figura II.6: Inseriamo in C2 il valore cos(B2).

II.3 Formule ricorsive

Supponiamo ora di voler implementare lalgoritmo dello schema di punto fisso, applicato alla funzioneg (x)= cos(x).

Nella colonna A scriveremo le iterazioni k che facciamo nellimplementare lo schema. Per implementarela formula xk = g (xk1) useremo le colonne B e C: in B metteremo i valori delle approssimazioni xk , in Cscriveremo il valore g (xk ).

Nella riga 1 delle colonne A, B e C, per ricordare ci che stiamo facendo, scriviamo semplicemente k, xk eg(xk), rispettivamente.

Al primo passo, nella cella A2 diamo il valore 0 per indicare literazione k = 0, a B2 assegniamo il valore ini-ziale per x0, per esempio 0.0 e in C2 assegniamo il valore g (x0), cio cos(B2), che ci servir al passo successivo.Osserviamo che abbiamo usato la cella B2 come una variabile avendo dato a C2 il valore cos(B2). Possiamoscrivere a mano cos(B2) oppure, inserendo la funzione cos() e cliccando con il mouse sulla variabile-cellain cui valutare la funzione (i cui contorni diventeranno di colore rosso), automaticamente viene inserito ilnome della cella allinterno delle parentesi (in questo caso la cella B2). Dopo aver schiacciato il tasto di inviocomparir il valore numerico della funzione calcolata (si veda la Figura II.6).

Al passo successivo, il passo k deve essere uguale a 1: applichiamo alla cella A3 la formula A2+1 in mododa ottenere 1 (vedremo perch conviene applicare questa semplice formula). La cella B3 ha come significatox1 = g (x0)= cos(B2)=C2, il valore che si trova nella cella C2. Nella cella C3 andremo a valutare cos(B3), inmodo da utilizzarlo nella cella B4 al passo successivo.

Siamo ora in grado di scrivere queste formule in maniera ricorsiva effettuando una sola operazione dicopia-incolla o di trascinamento. . .

Difatti, ora sulla riga 4 noi dobbiamo:in A4, incrementare il valore delliterazione che viene effettuata (la cella A4 deve essere uguale ad A3+1);in B4, aggiornare la variabile xk mettendo nella cella B4 il valore che si trova nella cella C3;in C4, aggiornare la variabile g (xk ) calcolando cos(B4).

28

II.3. Formule ricorsive

Figura II.7: Operazione di copia-incolla per applicare le formule ricorsive dalla riga 3 alla riga 4. Nella figurain alto, abbiamo annerito le tre cellette della riga 3. Nella figura in basso, dopo aver fatto Copy delle trecellette della riga 3, anneriamo le tre cellette della riga 4 su cui incolleremo le formule ricorsive.

29


Figura II.8: Risultato delloperazione di copia e incolla dalle cellette della riga 3 alle cellette della riga 4.

Questa operazione andr ripetuta fino a quandonon ci accorgiamodi essere arrivati a convergenza - vedremopoi come.

Per aggiornare la riga 4 basta semplicemente annerire con il mouse le cellette A3-B3-C3 (dove ci sonodelle formule che si riferiscono alla riga precedente), spostandoci con il mouse e cliccando contemporanea-mente con il tasto sinistro del mouse stesso. Clicchiamo quindi con il tasto destro del mouse e si aprir unafinestra di comandi: in particolare, clicchiamo alla voce Copy (Copia) e poi con il mouse anneriamo (con lastessa tecnica descritta prima) le tre cellette successive.

Clicchiamo di nuovo il tasto destro del mouse e facciamo Paste (Incolla). Vedremo che le cellette ora sonoaggiornate e con le stesse formule che avremmo scritto amano celletta per celletta (in A4 troviamo la formulaA3+1, in B4 troviamo C3, in C4 troviamo cos(B4) ). Si vedano le Figure II.7 e II.8. Se anzich annerire solole cellette della riga successiva, annerrissimo le cellette di molte pi righe, lalgoritmo verr applicato insuccessione (talmente rapida che ci sembra istantanea) per molti pi passi (tanti quante sono le righe su cuiabbiamo applicato la formula ricorsiva mediante loperazione di copia e incolla). Si veda la Figura II.9.

Loperazione di copia e incolla si pu fare anche in unaltra maniera mediante loperazione di trascina-mento del mouse. Una volta che abbiamo annerito le cellette in cui si trovano le formule ricorsive da ap-plicare sulle righe successive, osserviamo che langolo in basso a destra della cornice nera che si venuta acreare ha un bordino pi spesso. Ci avviciniamo con il tasto sinistro del mouse su quel bordino e poi tra-sciniamo il mouse verso il basso. Vediamo che, via via che andiamo verso il basso, i bordi delle cellette sucui vogliamo applicare le formule ricorsive diventano di colore rosso. Una volta che lasciamo il tasto sini-stro del mouse, le formule vengono applicate e troviamo i valori numerici nelle cellette. Notiamo che nellaversione 2.3 di OpenOffice, una volta che lasciamo il tasto sinistro del mouse, per completare loperazione ditrascinamento dobbiamo cliccare su Formule nella finestra che automaticamente viene aperta. Dopo avercompletato questa operazione, si conclude correttamente loperazione di trascinamento. Osserviamo che cisono al pi due cifre decimali nella rappresentazione dei numeri. Per visualizzare meglio i numeri, convienecambiare formato: anneriamo le cellette delle colonne B e C, clicchiamo il tasto destro del mouse e, sullafinestra che si apre, andiamo su Format Cells (Formatta Celle). Si apre una nuova finestra da cui clicchiamosu Numbers (Numeri) e scegliamo la categoria del formato scientifico. Qui possiamo definire manualmente

30

II.4. Fare grafici nel foglio elettronico

Figura II.9: Risultato delloperazione di trascinamento (o di copia e incolla) su pi righe.

il numero di zeri che vogliamo nella mantissa (vedasi Figura II.10). 2 Se vogliamo ora aggiungere delle itera-zioni allalgoritmo, mediante unoperazione di copia e incolla o di trascinamento, oltre alle formule anche ilformato viene incollato/trascinato nelle celle successive. Per capire quando arrestare lalgoritmo, utilizziamoil test di arresto sugli scarti. Sulla colonna D andremo a scrivere il valore assoluto della differenza tra dueapprossimazioni successive e applicheremo tale formula in maniera ricorsiva. Faremo quindi ABS(B3-B2)e trascineremo questa formula nelle celle successive con il formato che ci sembrer pi opportuno (si ve-da la Figura II.11). In base alla tolleranza prefissata, quando lo scarto diverr minore della tolleranza, alloraarresteremo lalgoritmo (e quindi ci fermeremo nelloperazione di copia e incolla o di trascinamento delleformule da una cella alla successiva).

II.4 Fare grafici nel foglio elettronico

Una volta terminata limplementazione del metodo di punto fisso, possiamo pensare di fare un graficosemilogaritmico degli scarti in funzione delle iterazioni.

Dobbiamo quindi vedere come fare nel foglio elettronico un grafico in cui sullasse delle ascissemettiamoi valori delle iterazioni e sullasse delle ordinate inseriamo i valori assoluti degli scarti. Inoltre, dobbiamorendere il grafico semilogaritmico per studiare il profilo di convergenza del metodo.

Incominciamo annerendo solo la colonna che contiene le iterazioni dal valore 1 fino alla fine delleiterazioni che abbiamo effettuato: questa operazione viene fatta tenendo premuto il tasto sinistro del mouse.

2Osserviamo che il formato scelto pu andare stretto nella cella/celle selezionate: in tal caso, al posto del numero saranno visua-lizzati i simboli ###. Per avere una corretta visualizzazione dei numeri, si dovr allargare la colonna corrispondente in cui si trovano lecelle per le quali non si vede il formato corretto. Cliccando con il tasto sinistro del mouse sulla cella che porta la lettera della colonna chedesideriamo allargare, vedremo che tutte le celle della colonna selezionata si anneriscono. Ci si sposta poi, con il mouse, verso il bordodestro della cella (per esempio ci troviamo sulla cella della colonna C al contorno con la cella della colonna D) e vedremo comparire unadoppia freccia, grazie alla quale si potr allargare tutta lampiezza della colonna verso destra. Allargheremo la colonna fino a quandotutte le celle saranno visualizzate correttamente nel loro formato.

31


Figura II.10: Finestra per modificare il formato numerico dei risultati. Viene scelto un formato scientificodefinito dallutente.

Figura II.11: Schema del punto fisso, in cui stata aggiunta, nella colonna D, il valore assoluto dello scartotra due approssimazioni successive. Si noti che abbiamo scelto un formato numerico diverso per gli scartirispetto a quello usato per xk e cos(xk).

32


Figura II.12: Selezione delle colonne per fare il grafico.

Dobbiamo poi passare alla colonna degli scarti, dal valore ottenuto alliterazione 1 fino al valore ottenutoallultima iterazione. Per non perdere traccia della colonna delle iterazioni, prima di iniziare a premere conil tasto sinistro del mouse sulla colonna degli scarti, dobbiamo schiacciare anche il tasto Ctrl della tastiera delcomputer. Solo dopo andiamo ad annerire la colonna degli scarti. In questo modo lasciamo nera anche lacolonna delle iterazioni (si veda la Figura II.12).

A questo punto, andiamo sul Menu principale del foglio elettronico alla voce Insert (Inserisci) e di qui suChart (Diagramma).

Si apre una finestra che ci permette la scelta del grafico. Noi sceglieremo il grafico XY (Scatter) (nellaversione italiana XY (Dispersione) ) con solo linee - Lines Only (vedasi Figura II.13). Procedendo (con Next -Avanti), si vede che lintervallo dei punti su cui fare il grafico stato gi scelto e compare nella riga relativaal Data Range (Area dati). Per la serie dei dati, Data Series, si vede quali sono i valori delle ascisse (X-values)e quali le ordinate (Y-values). Si vede anche la scritta Column D che indica il fatto che abbiamo preso per leordinate i valori presi dalla colonna D (vedremo che poi uscir il nome Column D nella legenda della figura evedremo come dare un nome pi significativo) (si veda la Figura II.14). Andando ancora avanti, nella ChartElements (Elementi del diagramma), si pu decidere cosa scrivere come titolo, sottotitolo, etichette sugli assix e y. Nella Figura II.15 abbiamo dato titolo, sottotitolo. . .A questo punto si clicca su Finish (Fine) e sul foglioelettronico compare il nostro grafico (che non ancora semilogaritmico), mostrato nella Figura II.16.

Per fare modifiche al grafico, bisogna cliccare una volta con il tasto sinistro del mouse posizionato sulgrafico, in modo da vedere in alto a sinistra unicona con il simbolo dellncora. Ci si posiziona con il mouseallinterno del grafico e si clicca una volta con il tasto destro. Compare una finestra e qui si clicca la voce Edit(Edita). Ora possibile cambiare il tipo di assi o le linee del grafico, i colori, la legenda. . . .

Ad esempio, andando con il mouse sullasse delle Y e cliccando due volte con il tasto sinistro del mouse,si apre una finestra da cui poter cambiare le linee, i colori, i caratteri, ma soprattutto la scala del grafico.

Da Scale (Scala) infatti possibile scegliere la scala logaritmica come si vede nella Figura II.17. Una voltacliccato su OK, si torna al grafico che non pi quello di prima ma in scala semilogaritmica (si veda la Fi-

33


Figura II.13: Finestra per la scelta del grafico

Figura II.14: Finestra del Data Series

Figura II.15: Titolo, sottotitolo, etichette.

34


Figura II.16: Grafico

Figura II.17: Scelta della scala

gura II.18). Possiamo fare altre modifiche (per esempio cliccando sulla curva del grafico (Data Series ColumnD), possiamo cambiare lo spessore della curva). Inoltre, possiamo cambiare i caratteri usati per il titolo, icolori dello sfondo. E cos via.

Quello che ci interessa ora come cambiare la legenda.Clicchiamo due volte con il tasto sinistro del mouse sul grafico. Poi clicchiamo con il tasto destro: si apre

unafinestra con varie opzioni di scelta tra cuiData Ranges. Clicchiamo suData Ranges e si aprir una finestra,quella che vediamo nella Figura II.19. Per cambiare la scritta Column D, dobbiamo inserire nella casellaRange for Name (Area per Nome) il nome di una cella del foglio elettronico in cui abbiamo scritto un nome dilegenda che riteniamo pi appropriato (per esempio, nella cella E1 abbiamo scritto punto fisso). Cliccandosullicona che si trova a destra della finestra del Range for Name, possiamo poi direttamente cliccare sullacella E1 in modo da inserire automaticamente il nome della cella sulla casella Range for Name (il risultato diquesta operazione mostrata nella Figura II.20. Vediamo che non compare pi il nome Column Dma puntofisso sia su questa finestra sia, dopo aver dato lOK, sul grafico (vedasi la Figura II.21).

35


Figura II.18: Grafico semilogaritmico

Figura II.19: Finestra del Data Ranges

II.4.1 Grafici multicolonna

Osserviamo che la procedura eseguita per fare il grafico semilogaritmico appena descritto pu essere ri-petuta in maniera analoga per fare grafici multicolonna (per esempio, sullasse delle x riportiamo le iterate,mentre sullasse delle ordinate riportiamo i valori assoluti degli scarti relativi ai metodi di Newton-Raphson,Regula Falsi, tangente fissa, secante fissa. . . ): una volta che abbiamo inserito i dati relativi alle ascisse e alleordinate della prima serie di dati (per esempio quelli delmetodo di Newton-Raphson), ementre stiamo anco-ra lavorando sulla finestra che ci permette di costruire il grafico, dopo la scelta del tipo di grafico (Chart Type)e Data Range, al passo del Data Series possiamo aggiungere (Add) altre serie di dati (per esempio iterazioni escarti del metodo della Regula Falsi) per fare un grafico multicolonna. Si vedano le Figure II.22 e II.23.

36


Figura II.20: Finestra del Data Ranges dopo aver operato la modifica

Figura II.21: Grafico finale semilogaritmico

37


Figura II.22: Dopo aver schiacciato sul tasto Add (Aggiungi), per aggiungere altri dati nel grafico, bisognaselezionare i nuovi dati.

Figura II.23: Dopo aver selezionato le cellette opportune dei valori delle ascisse da aggiungere al grafico, ivalori vengono visualizzati sulla finestra. Analogamente si procede per lasse y.

II.5 Inserire dati da files

Spesso il foglio elettronico viene utilizzato per fare grafici importando dati numerici da files esterni per esempio files di tipo testo generati dallesecuzione di programmi scritti in FORTRAN o altri linguaggi diprogrammazione.

Per inserire correttamente i dati su pi colonne, si pu seguire questo percorso: dal Menu principale delfoglio elettronico si va su Insert (Inserisci) > Sheet from file (Foglio da file). A questo punto si apre unafinestra da cui possibile selezionare il percorso in cui si trova il file da aprire. Una volta che importiamoil file, si apre unaltra finestra (la finestra di Text import (importa testo)) in cui dobbiamo selezionare alcuniparametri che ci permetteranno di importare i dati inserendoli su pi colonne (se i dati sono scritti su picolonne, che la cosa che ci importa per poter fare i successivi grafici).

Le opzioni che dobbiamo selezionare sono, oltre a Tab, che gi troviamo indicato, ancheMerge delimiterse Space (si veda la Figura II.24). Una volta dato lOK, si apre unaltra finestra (Insert Test) in cui digitiamo dinuovoOK. Troviamo i dati del filemessi in colonna come erano nel file di partenza. Il formato numerico non molto carino (solo due cifre decimali, anche se andando a vedere celletta per celletta possiamo osservare chei valori sono memorizzati cos come erano nel file originario), ma selezionando le celle possiamo cambiare

38

II.5. Inserire dati da files

Figura II.24: Importare dati da un file esterno.

il formato in modo da visualizzarlo correttamente. A questo punto, abbiamo i nostri dati in tabella sul foglioelettronico e possiamo fare i grafici o altre operazioni tramite il foglio elettronico.

39

IndicePrimi passi in MATLABIntroduzioneAvvio di MATLABMatrici e vettori in MATLAB

Comandi utiliMATLABper scrivere ed eseguire programmiStrutture

Dati di inputProgramma sul metodo delle bisezioni

Dati di outputGraficiSulle potenzialit di MATLABApplicazioni di MATLABnel Calcolo NumericoSull'instabilit numericaSull'interpolazione e approssimazione di datiConfronto tra schemi per equazioni differenziali ordinarie

Tabelle e grafici con il foglio elettronicoIl foglio elettronicoInserire funzioniFormule ricorsiveFare grafici nel foglio elettronicoGrafici multicolonna

Inserire dati da files

Calcolo Numerico - Matlab

Documents

Transcript of Calcolo Numerico - Matlab