Esercitazioni di C++ per Fisici Annotazioni dal Corsofontana/+.pdf · 2011. 7. 14. · dati, il...

Esercitazioni di C++ per Fisici

Annotazioni dal Corso

Andrea Fontana

Dipartimento di Fisica Nucleare e TeoricaUniversita di Pavia

A.A. 2004-2005

1

Indice

1 Numeri Complessi 41.1 Esercizio 1.1: struttura in stile C e suo uso . . . . . . . . . . . . 51.2 Esercizio 1.2: aggiunta di funzioni alla struttura . . . . . . . . . 81.3 Esercizio 1.3: trasformazione da C a C++ . . . . . . . . . . . . . 121.4 Esercizio 1.4: compilazione con Makefile . . . . . . . . . . . . . . 15

2 Frattali 172.1 Esercizio 2.1: classe Complesso e suo uso . . . . . . . . . . . . . . 182.2 Esercizio 2.2: overload di metodi e di operatori . . . . . . . . . . 222.3 Esercizio 2.3: classe FrattC per mappa logistica . . . . . . . . . . 292.4 Esercizio 2.4: classe FrattC per insiemi di Julia e Mandelbrot . . 392.5 Esercizio 2.5: libreria di ROOT e output grafico . . . . . . . . . . 47

3 Sistema Solare 613.1 Esercizio 3.1: classe Vettore e suo uso . . . . . . . . . . . . . . . 623.2 Esercizio 3.2: overload di metodi e di operatori . . . . . . . . . . 663.3 Esercizio 3.3: classi CorpoCeleste e SistemaSolare (con elementi

di design in UML) . . . . . . . . . . . . . . . . . . . . . . . . . . 703.4 Esercizio 3.4: simulazione del Sistema Solare con output numeri-

co e grafico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 823.5 Esercizio 3.5: aggiunta di una Sonda . . . . . . . . . . . . . . . . 91

4 Bibliografia 100

3

Introduzione alle Esercitazioni

Queste lezioni riguardano la programmazione in linguaggio C++ e consisto-no in un ciclo di esercitazioni a supporto del Corso di Informatica per la FisicaII tenuto dalla Prof. Adele Rimoldi presso l’ Universita di Pavia a partire dall’Anno Accademico 2001-2002. Il corso di Esercitazioni di Informatica II e uncorso di laboratorio e viene effettuato dagli studenti utilizzando il calcolatoreper scrivere programmi in C++ e per utilizzarli successivamente in ambienteUnix: l’ obiettivo principale e quello di fornire agli studenti uno strumento perrisolvere problemi di calcolo.

Le esercitazioni sono suddivise in tre parti: partendo dal concetto di strut-tura in C, viene introdotta la nozione di classe che gioca un ruolo centrale inC++ e nell’ ambito della programmazione ad oggetti e successivamente ven-gono illustrate le nozioni fondamentali del C++, quali l’ incapsulamento deidati, il sovraccarico degli operatori e infine l’ ereditarieta degli oggetti. Questiargomenti vengono presentati attraverso esercizi, intesi nel senso tradizionaledel termine, che riguarderanno argomenti di interesse fisico e di attualita: nellaprima parte verra mostrato come poter utilizzare numeri complessi sia in un pro-gramma C sia in un programma C++, sottolineando la possibilita di integrarei due linguaggi, pur utilizzando la semplice nozione di struttura.

Nella seconda parte, dopo aver definito una classe C++ per rappresentarenumeri complessi, gli oggetti con essa costruiti saranno utilizzati per definire unanuova classe che rappresenti un insieme frattale. Il concetto di frattale, breve-mente introdotto, sara’ illustrato iterando semplici algoritmi in C++: verrannoriprodotti l’ insieme di biforcazione della mappa logistica e gli insiemi di Juliae Mandelbrot, accennando brevemente anche all’ utilizzo di librerie grafiche perla visualizzazione dei risultati del programma.

La terza parte infine si propone di scrivere in C++ una simulazione delSistema Solare. Si tratta di un programma che utilizzera varie classi: una classevettore (analoga alla classe numero complesso) e le classi per descrivere il motodi corpi celesti, del sistema solare e anche di sonde interplanetarie in base allalegge di gravitazione di Newton. Un semplice utilizzo di una libreria graficaconcludera il corso presentando un modello animato del nostro sistema solare.

Buon lavoro e buon divertimento con il C++!

Andrea FontanaPavia, Marzo 2005

5

1 Numeri Complessi

La prima parte di queste esercitazioni si ricollega direttamente al Corso di Ceffettuato dall’ autore durante il I semestre. L’ obiettivo e di illustrare come siapossibile utilizzare un programma in C per lavorare con numeri complessi, cioecon un tipo di dato che non e intrinsecamente incluso nel linguaggio: questoporta alla scrittura di una struttura Complesso e ad una estensione rudimentaledel linguaggio C per l’ utilizzo di nuovi tipi di dati.

Si vuole inoltre sottolineare la necessita di definire delle funzioni in modoesplicito per effettuare delle operazioni algebriche sui numeri definiti e l’ impos-sibilita di ricorrere alla notazione intuitiva dei simboli algebrici usuali: questaimpostazione servira per introdurre successivamente il concetto di classe in mo-do naturale e permettera di apprezzare il vantaggio della definizione dei metodidi interfaccia della classe e del sovraccarico degli operatori.

Infine verra mostrato come sia possibile integrare il linguaggio C++ nel lin-guaggio C, utilizzando in modo perfettamente trasparente le nuove istruzioni delC++ in un programma che utilizza strutture C: in questo senso, come peraltrosuggerito dal nome del linguaggio, il C++ puo essere visto come una evoluzionee un potenziamento del C.

Il capitolo viene chiuso da un richiamo all’ utilizzo del comando make percompilare e linkare programmi suddivisi su piu files, come richiesto dalla pro-grammazione in C++.

4

1.1 Esercizio 1.1: struttura in stile C e suo uso

Scrivere una struttura in C che rappresenti un numero complesso e usarla in unprogramma C che permetta di effettuare le operazioni algebriche elementari connumeri complessi.

• Listato del programma

file Struttura1.c:

#include <stdio.h>#include <math.h>

main(){double den;

struct Complesso{double x;double y;

};

struct Complesso a;struct Complesso b = {3.,5.};struct Complesso c;

a.x = 0.;a.y = 1.;

printf("a: (%lf,%lf)\n",a.x,a.y);printf("b: (%lf,%lf)\n\n",b.x,b.y);

/* coniugio */c.x = b.x;c.y = -b.y;printf("cc: (%lf,%lf)\n",c.x,c.y);/* modulo */printf("cm: %lf\n",sqrt(c.x*c.x + c.y*c.y));/* fase */printf("cf: %lf\n\n",atan2(c.y,c.x));

...

5

.../* somma */c.x = a.x + b.x;c.y = a.y + b.y;printf("c+: (%lf,%lf)\n",c.x,c.y);/* sottrazione */c.x = a.x - b.x;c.y = a.y - b.y;printf("c-: (%lf,%lf)\n",c.x,c.y);/* moltiplicazione */c.x = a.x*b.x - a.y*b.y;c.y = a.x*b.y + a.y*b.x;printf("c*: (%lf,%lf)\n",c.x,c.y);/* divisione */den = b.x*b.x + b.y*b.y;c.x = (a.x*b.x - a.y*b.y)/den;c.y = (a.x*b.y + a.y*b.x)/den;printf("c/: (%lf,%lf)\n",c.x,c.y);

}

• Note esplicative

– Il programma richiede l’ uso delle funzioni di I/O e delle funzionimatematiche della libreria standard del C e quindi e’ necessario in-cludere all’ inizio i rispettivi header files chiamando la direttiva peril preprocessore #include.

– Una struttura in C consiste in un nuovo tipo di variabile che risul-ta composto da piu variabili di tipo elementare (int, float, doubleo char). In questo caso definiamo un tipo di variabile Complesso,costituito da due variabili double che rappresentano la parte realee la parte immaginaria del numero complesso. Si noti il ; (punto evirgola) finale alla chiusura del corpo della struttura!

– Definiamo poi 3 numeri complessi a, b e c e inizializziamo b conil numero complesso 3+5i attraverso l’ uso delle parentesi graffe.La variabile a viene invece inizializzata accedendo separatamente aidue membri x e y e i risultati vengono poi stampati con due funzio-ni printf, trattando sempre in modo esplicitamente indipendente leparti reale e immaginaria.

– La variabile c e’ usata per memorizzare il coniugato del numero b eper stamparne il valore, nonche’ il modulo e la fase: si noti che inogni caso dobbiamo accedere ai membri della struttura trattandoliseparatamente in termini dei tipi elementari del C. In altri termini,non esiste una funzione printf che stampi direttamente un numerocomplesso considerandolo come un’ entita unica. La funzione sqrt()

6

calcola la radice quadrata del suo argomento e atan2() l’ arcotangen-te di un angolo in un triangolo rettangolo di cateti y e x, tenendoconto del segno dei due argomenti per determinare il quadrante delrisultato: entrambe appartengono alla libreria matematica del C chedeve essere quindi linkata con l’ opzione -lm del compilatore:gcc -o Struttura1 Struttura1.c -lm

– Il calcolo delle quattro operazioni algebriche fondamentali su a e bcon memorizzazione dei risultati in c puo essere fatto esplicitamentecalcolando separatamente la parte reale e la parte immaginaria delrisultato secondo il seguente schema:

c = (Re a + Re b) + i(Im a + Im)¯

c = (Re a− Re b) + i(Im a− Im b)c = (Re a ∗ Re b− Im a ∗ Im b) + i(Re a ∗ Im b + Re b ∗ Im a)

c =(Re a ∗ Re b + Im a ∗ Imb) + i(Re b ∗ Im a− Re a ∗ Im b)

(Re b)2 + (Im b)2

Si noti l’ uso della variabile d’ appoggio den per il calcolo del deno-minatore e in generale la necessita di scrivere tutti i conti in dettaglioper fare una semplice moltiplicazione o divisione: non sarebbe megliopoter scrivere direttamente c=a*b o c=a/b e lasciare al programmalo svolgimento dei conti in dettaglio? In C questo non e’ possibilee l’ evoluzione del C in C++ ha tra le sue motivazioni la risposta aquesto problema.

7

1.2 Esercizio 1.2: aggiunta di funzioni alla struttura

Scrivere un programma che permetta di effettuare operazioni su numeri com-plessi attraverso la chiamata a funzioni.

• Listato dei programmi

file Complesso.h:

struct Complesso{double x;double y;

};

file Struttura2.c:

#include <stdio.h>#include <math.h>#include "Complesso.h"

/* Prototipi delle funzioni */struct Complesso Coniugio(struct Complesso c);double Modulo(struct Complesso c);double Fase(struct Complesso c);struct Complesso Addizione(struct Complesso c1,

struct Complesso c2);struct Complesso Sottrazione(struct Complesso c1,

struct Complesso c2);struct Complesso Moltiplicazione(struct Complesso c1,

struct Complesso c2);struct Complesso Divisione(struct Complesso c1,

struct Complesso c2);

main(){struct Complesso a;struct Complesso b = {3.,5.};struct Complesso c;

a.x = 0.;a.y = 1.;

printf("a: (%lf,%lf)\n",a.x,a.y);printf("b: (%lf,%lf)\n\n",b.x,b.y);

...

8

...

/* coniugio */c = Coniugio(b);printf("cc: (%lf,%lf)\n",c.x,c.y);/* modulo */printf("cm: %lf\n",Modulo(c));/* fase */printf("cf: %lf\n\n",Fase(c));/* addizione */c = Addizione(a,b);printf("c+: (%lf,%lf)\n",c.x,c.y);/* sottrazione */c = Sottrazione(a,b);printf("c-: (%lf,%lf)\n",c.x,c.y);/* moltiplicazione */c = Moltiplicazione(a,b);printf("c*: (%lf,%lf)\n",c.x,c.y);/* divisione */c = Divisione(a,b);printf("c/: (%lf,%lf)\n",c.x,c.y);

}

struct Complesso Coniugio(struct Complesso c){struct Complesso z;z.x = c.x;z.y = -c.y;

return z;}

double Modulo(struct Complesso c){return sqrt(c.x*c.x + c.y*c.y);

}

double Fase(struct Complesso c){return atan2(c.y,c.x);

}...

9

...struct Complesso Addizione(struct Complesso c1,

struct Complesso c2){struct Complesso z;z.x = c1.x + c2.x;z.y = c1.y + c2.y;

return z;}struct Complesso Sottrazione(struct Complesso c1,

struct Complesso c2){struct Complesso z;z.x = c1.x - c2.x;z.y = c1.y - c2.y;

return z;}

struct Complesso Moltiplicazione(struct Complesso c1,struct Complesso c2)

{struct Complesso z;z.x = c1.x*c2.x - c1.y*c2.y;z.y = c1.x*c2.y + c1.y*c2.x;

return z;}

struct Complesso Divisione(struct Complesso c1,struct Complesso c2)

{struct Complesso z;double den;

den = c2.x*c2.x + c2.y*c2.y;z.x = (c1.x*c2.x - c1.y*c2.y)/den;z.y = (c1.x*c2.y + c1.y*c2.x)/den;

return z;}

10


– Scriviamo la struttura in un file separato di nome Complesso.h: que-sto non e strettamente necessario, ma risulta comodo in vista di quan-to si fara abitualmente nel seguito definendo delle classi in C++ epermette inoltre di definire variabili globali di tipo Complesso.

– Il file che contiene la funzione main() include questo file .h all’ inizioe pertanto il file prende il nome di header file (file di intestazione).Si noti l’ uso delle parentesi acute <> e del doppio apice: il prepro-cessore cerca i files tra parentesi acute nelle directory di sistema diUnix (/include, /usr/include, /usr/local/include) e i files tra doppiapici nella directory corrente (files scritti dall’ utente).

– Definiamo prima di main() i prototipi delle funzioni che vogliamodefinire, specificando il numero e il tipo dei parametri e il tipo resti-tuito dalla funzione: poiche Complesso e a tutti gli effetti un nuovotipo di variabile lo possiamo passare alle funzioni come ogni altrotipo elementare. In particolare le funzioni di coniugio, addizione,sottrazione, moltiplicazione e divisione accettano in ingresso uno odue numeri complessi e restituiscono un numero complesso, mentreil calcolo di modulo e fase restituisce naturalmente solo un numeroreale dato un numero complesso in ingresso. In tutti i casi di que-sto esempio, la chiamata alla funzione viene fatta per semplicita pervalore, ma sarebbe possibile effettuarla anche per indirizzo passandoalla funzione l’ indirizzo della variabile.

– Il programma main() e analogo al precedente, con la differenza chele operazioni vengono svolte chiamando le rispettive funzioni primadella stampa del risultato: inoltre un commento, racchiuso tra */ e*/, illustra l’ operazione in corso.

– La funzione main() e seguita dalla dichiarazione di tutte le funzio-ni prima definite attraverso i prototipi. Le funzioni di coniugio e leoperazioni algebriche definiscono al loro interno una variabile localez in cui viene memorizzato il numero complesso calcolato e poi re-stituito a main(). Le funzioni di modulo e fase restituiscono invecesubito il risultato che viene calcolato sulla riga dell’ istruzione re-turn. Le operazioni fatte all’ interno delle varie funzioni accedonocomunque sempre ai singoli tipi elementari che sono membri dellastruttura Complesso.

11

1.3 Esercizio 1.3: trasformazione da C a C++

Utilizzare la struttura Complesso in un programma C++ che utilizzi numericomplessi in modo sia statico sia dinamico facendo uso delle funzioni di I/O delC++.

• Listato del programma

file Struttura3.C:

#include <iostream>#include <stdio.h>#include "Complesso.h"

using namespace std;

main(){int i;struct Complesso a = {2.,7.};struct Complesso *z;struct Complesso w[100];

struct Complesso *b;

z = &a;printf("(%lf,%lf)\n",a.x,a.y);printf("(%lf,%lf)\n",(*z).x,(*z).y);printf("(%lf,%lf)\n\n",z->x,z->y);

for(i=0;i<100;i++){w[i].x = i;w[i].y = -i;printf("(%lf,%lf)\n",w[i].x,w[i].y);

}

b = new Complesso();b->x = 1.;b->y = 2.;printf("\n(%lf,%lf)\n",b->x,b->y);delete b;b = NULL;printf("(%lf,%lf)\n",b->x,b->y);

std::cout << b->x << "\t" << b->y << std::endl;}

12


– In questo programma integriamo delle istruzioni C++ in un program-ma che utilizza le strutture Complesso definite finora: questo mostracome i linguaggi C e C++ siano collegati tra loro e come sia possibileil loro utilizzo simultaneo. Il C++ visto cosı potrebbe sembrare soloun’ estensione del C (come suggerirebbe il nome ad un osservatorecasuale) e a questo punto del corso questo punto di vista va bene:la potenza del C++ rispetto ad altri linguaggi risiede nei concetti diclasse e di oggetto che verranno ampiamente trattati nella secondae terza parte delle dispense e che riguardano la programmazione adoggetti ed un approccio alla programmazione completamente diversoda quello visto finora.

– Prima di introdurre le novita del C++, vediamo come sia possibiledefinire dei puntatori a strutture e degli array di strutture con leusuali notazioni: dato un puntatore a Complesso z e’ sempre possibileaccedere ai membri della struttura con la notazione (*z).x, in cui laparentesi e necessaria in quanto l’ operatore . ha precedenza maggioredell’ operatore *. Questa notazione puo essere abbreviata con lanotazione z → x:

(∗z).x = z → x (1)

in cui si usa l’ operatore freccia → per accedere alla memoria defe-renziata dal puntatore.

– Veniamo ora a due aspetti nuovi del C++: le funzioni di I/O e lagestione dinamica della memoria.

– L’ header file iostream.h contiene le definizioni necessarie per utilizza-re le funzioni di I/O del C++ di cui vedremo solamente l’ equivalentedi printf(), cioe cout. Come si vede e possibile includere anche ilfile stdio.h del C e usare anche printf() nello stesso programma. Ilcomando std::cout appartiene al namespace std e puo essere interpre-tato a questo livello come un comando analogo a printf che accettadelle variabili, delle stringhe e anche dei caratteri di controllo dell’output (come \n per andare a capo o \t per il tabulatore) e che liconcatena usando l’ operatore <<: in realta cout e un oggetto C++basato su una classe e verra reinterpretato correttamente piu avanti.Il comando endl e un modo alternativo per andare a capo equivalen-te a \n. Il namespace std, definito all’ inizio del programma con l’istruzione using puo essere sottinteso (std sta per standard, cioe ilnamespace di default) in quanto per i programmi che vedremo nelseguito l’ unico oggetto cout e quello di default e non ci sono rischidi sovrapposizioni con altri oggetti di altri namespace.

– La gestione dinamica della memoria consente di ottimizzare la capa-cita del programmatore di accedere alle variabili e di usare la me-moria. Tutte le variabili definite finora in C e in C++ sono dette

13

variabili statiche e il loro scope viene deciso dal compilatore: unavariabile locale interna ad una funzione viene rimossa dalla memoriaquando l’ esecuzione della funzione termina e anche una variabile glo-bale viene cancellata al termine di un programma in cui era definita.La zona di memoria in cui le variabili statiche sono definite prende ilnome di stack. Esiste anche un altra zona di memoria, detta heap,in cui lo scope di una variabile, cioe la durata della sua esistenzain memoria, viene decisa dal programmatore con le istruzioni newe delete: new permette di allocare memoria e delete di liberare lamemoria, in modo totalmente indipendente dallo scope automaticodeciso dal compilatore. Una variabile dinamica, definita in una fun-zione, continua ad esistere anche dopo l’ esecuzione della funzione:in realta cio esisteva gia in C con due funzioni di libreria, malloc()e free(), ma non era implementato in modo molto efficiente ed e orastato sostituito dal metodo del C++.

– Per ogni istruzione new, deve essere presente la corrispondente istru-zione delete. Dopo un delete inoltre, e buona norma inizializzare ilpuntatore a NULL (cioe rimuovere ogni indirizzo di memoria): solocosı la zona di memoria allocata con new non e piu accessibile e lavariabile b viene realmente rimossa. Le ultime due istruzioni del pro-gramma provocano un crash in fase di esecuzione, in quanto si stacercando di accedere ad una zone di memoria che non e piu allocata!

14

1.4 Esercizio 1.4: compilazione con Makefile

Scrivere un Makefile per compilare e linkare il programma dell’ esercizio prece-dente con il comando di Unix make.

• Listato del file Makefile

file Makefile:

# istruzioni per creare il file eseguibile

Struttura3: Struttura3.og++ -o Struttura3 Struttura3.o

# istruzioni per creare il file oggetto

Struttura3.o: Struttura3.C Complesso.hg++ -c Struttura3.C


– Il comando make permette di semplificare la compilazione e il linkingdi un programma C o C++: quando un eseguibile dipende da piufiles sorgente, spesso la modifica di uno solo di questi files richiedecomunque la compilazione di tutti gli altri se si utilizza direttamente ilcompilatore dal prompt di Unix. Questo risulta scomodo, in quantoanche files non modificati vengono ricompilati inutilmente con unnotevole spreco di risorse e di tempo.

– Ad esempio, nel caso degli esercizi precedenti, potremmo pensaredi scrivere le funzioni per agire sulla struttura Complesso in un fileseparato chiamato Funzioni.C e potremmo compilarlo e linkarlo conl’ istruzione:g++ -o Struttura3 Funzioni.C Struttura3.COra, se facciamo una modifica nel programma main (che sta nel fileStruttura3.C) ma non nelle funzioni (file Funzioni.C), per creare ilnuovo eseguibile dobbiamo eseguire ancora il comando precedente equindi dobbiamo necessariamente ricompilare sia il file Struttura3.C(modificato) sia il file Funzioni.C (non modificato). Questa operazio-ne, specialmente in un progetto con molti files, risulta inefficiente edinutile e quindi e meglio non farla.

– Il comando make e la scrittura di un file detto Makefile permettonodi ovviare a questo problema controllando la data e l’ ora dei filessorgente, oggetto ed eseguibile e ricompilando solo i files sorgentecon data o ora successive a quelle dei rispettivi files oggetto (cioe piurecenti e quindi modificati dopo l’ ultima compilazione) e linkandopoi alla fine tutti i files oggetto con i files di libreria per produrre l’eseguibile.

15

– Le istruzioni per il comando make sono contenute in un file che deveessere chiamato Makefile e che deve risiedere nella stessa directorydei files dai compilare e seguono una sintassi particolare: per ognifile da creare si scrive una riga con il nome del file stesso seguito da :(due punti) e dai nomi di tutti i files da cui il file da creare dipende.Ad esempio, il nostro eseguibile Struttura3 dipende dal file oggettoStruttura3.o e quindi scriveremo:Struttura3: Struttura3.oIl concetto di dipendenza tra files e cruciale per comprendere il fun-zionamento di make: un file eseguibile dipende da uno o piu filesoggetto e un file oggetto dipende a sua volta da un file sorgente.Infatti il file Struttura3.o dipende dai files Struttura3.C e anche daComplesso.h, nel senso che ogni modifica a questi due files comportala rigenerazione del file Struttura3.o per essere inclusa nel programmaeseguibile.Struttura3.o: Struttura3.C Complesso.h

– Dopo ognuna delle righe che descrivono le dipendenze tra files dob-biamo specificare le istruzioni con cui i files da creare (targets) devonoessere creati a partire dai files sorgenti scrivendo esattamente le stes-se istruzioni che scriveremmo dal prompt di Unix per compilare e perlinkare dopo un carattere di tabulazione. Tutti i files oggetto creativengono linkati alle librerie producendo l’ eseguibile:g++ -o Struttura3 Struttura3.o(-o e l’ optione del compilatore/linker per rinominare il file eseguibileche di default viene chiamato a.out) e tutti i files sorgenti vengonocompilati nei rispettivi files oggetto:g++ -c Struttura3.C(-c e l’ opzione per compilare senza linkare). Il carattere tab e es-senziale per il funzionamento del comando make in quanto viene uti-lizzato per distinguere le istruzioni da eseguire dalle dipendenze trafiles.

– Come si vede, in un Makefile risultano nettamente separate la fasedi compilazione dei sorgenti e la fase di link dei vari oggetti allelibrerie: queste due operazioni sono sempre presenti, ma invocandoil compilatore dalla riga di comando vengono spesso incorporate inun’ unica operazione.

16

2 Frattali

La seconda parte di queste esercitazioni si propone di scrivere una classe C++che rappresenti numeri complessi e di utilizzarla per lo studio di particolarisuccessioni numeriche complesse generate da semplici algoritmi: questo portaalla generazione di insiemi frattali e fornisce lo spunto per scrivere un’ altra classeche rappresenti un oggetto frattale utilizzando la classe dei numeri complessi.

La classe Complesso e costituita da dati e da metodi (funzioni) membri diclasse e permette di effettuare le operazioni viste nella prima parte in modo piuefficiente e sintetico: in particolare verra visto l’ utilizzo di metodi inline e ilsovraccarico (overload) di metodi e operatori e quindi la possibilita di utilizzarei tradizionali operatori algebrici (+, -, * e /) nel calcolo di espressioni C++che coinvolgono numeri complessi. Sara inoltre fatto un cenno all’ overload deimetodi di I/O appartenenti alla classe iostream della libreria standard C++ perla stampa di numeri complessi.

Il concetto di frattale gioca un ruolo molto importante nella Fisica modernae verra introdotto in modo elementare attraverso le nozioni di invarianza discala e di autosomiglianza: saranno esaminati due esempi, la mappa logisticae il suo diagramma di biforcazione e gli insiemi di Julia/Mandelbrot, ottenutidall’ iterazione di espressioni complesse. Cio comportera la scrittura di unanuova classe FrattC che utilizza oggetti della classe Complesso e la progressivaaggiunta di nuovi metodi necessari nei vari casi: la trattazione sara limitataai concetti fondamentali e sara incentrata sugli aspetti di programmazione inC++ del problema.

Lo studio dei frattali proposti verra effettuato inizialmente in modo nume-rico e successivamente in modo grafico. Si mostrera infatti come utilizzare unalibreria grafica per rappresentare grafici ottenuti da valori numerici e distribu-zioni di punti nel piano complesso: la libreria utilizzata fa parte di un insiemedi tools C++ denominato ROOT e sviluppato ai laboratori di Fisica del Cerndi Ginevra e per un suo utilizzo piu completo si rimanda a testi piu avanzatiche esulano dagli obiettivi di questo corso.

17

2.1 Esercizio 2.1: classe Complesso e suo uso

Scrivere una classe in C++ che rappresenti un numero complesso e utilizzarlain un programma in C++.


file Complesso.h:

class Complesso{

public:Complesso();~Complesso();

void setRe(double a) {x=a;};void setIm(double b) {y=b;};double Re() {return x;};double Im() {return y;};

double Mod();double Phase();Complesso Conjugate();

private:double x;double y;

};

18

file Complesso.C:

#include "Complesso.h"#include <math>

Complesso::Complesso(){

x=0;y=0;

}

Complesso::~Complesso(){}

double Complesso::Mod(){

return sqrt(x*x + y*y);}

double Complesso::Phase(){

return atan2(y,x);}

Complesso Complesso::Conjugate(){

return Complesso(x,-y);}

file UsoComplesso1.C:

#include <iostream>#include "Complesso.h"

main(){

Complesso z;z.SetRe(7);z.SetIm(5);cout << z.Mod() << "\t" << z.Phase() << endl;

return 0;}

19


– La scrittura di una classe permette di definire un nuovo tipo di va-riabile (in questo caso Complesso) del tutto analoga a int, float odouble. Si tratta di un costrutto molto potente del C++ che sta alcuore della programmazione ad oggetti: con oggetto si intende unavariabile di un tipo definito in base a una certa classe. La nozionedi classe estende quella di struttura fornendo la possibilita di rag-gruppare non solo dei dati, ma anche delle funzioni che agiscano suquesti dati e introduce quindi il concetto di oggetto: un oggetto e unnuovo tipo di variabile definito dall’ utente e costituito da variabili eda funzioni che agiscono su queste variabili.

– Una classe e pertanto a tutti gli effetti un miniprogramma, con varia-bili e funzioni proprie, che puo essere utilizzato in altri programmi,estendendo la libreria di default del C++ non solo per quanto ri-guarda le funzioni (come gia succede in altri linguaggi come il C),ma piuttosto per quanto riguarda i tipi di variabili: cio e gia possi-bile con le strutture che, come abbiamo visto, non permettono perodi raggruppare funzioni, ma solo dati. Le variabili di una classe so-no locali e possono essere viste dall’ esterno solamente utilizzando lefunzioni della classe, che prendono il nome di metodi. Questo fattoe evidenziato dalle parole chiave private: e public: che separano levariabili e le funzioni, rispettivamente: le variabili sono private e sonoaccessibili attraverso le funzioni, che sono public, e possono quindiessere chiamate con l’ operatore .(punto). Ad esempio l’ istruzione:

Complesso c;c.Re();

restituisce la parte reale del numero, mentre scrivere

c.x;

darebbe un messaggio di errore perche si sta cercando di accederead una variabile protetta. Questa caratteristica del C++ prende ilnome di incapsulamento dei dati e serve a garantire la separazionetra dati e metodi all’ interno delle classi.

– La parola chiave class e seguita da un nome che identifichera poi ilnuovo tipo di variabile (in questo caso Complesso) e dalle parentesigraffe che identificano il corpo della classe (si noti il punto e virgolafinale).

– Due funzioni sono molto importanti e portano lo stesso nome dellaclasse: il costruttore Complesso() e il distruttore ˜Complesso(). Ilcostruttore viene chiamato quando definiamo una variabile Comples-so e serve a inizializzare le variabili membro di una classe, ponendolead esempio uguali a zero, come in questo caso. Il distruttore viene

20

chiamato alla fine dello scope dell’ oggetto dichiarato e di defaultcancella le variabili dalla memoria: nel caso di un oggetto statico cioavviene quando la variabile esce di scope, mentre per un oggetto di-namico il distruttore viene invocato dall’ istruzione delete applicataad un oggetto della classe.

– Tutti i metodi del file di implementazione sono definiti specificandoil nome della classe prima del nome del metodo, separandoli con ilsimbolo :(due punti) ripetuto due volte. Le funzioni per agire sui datimembro e lavorare con numeri complessi sono normali funzioni cheservono a impostare la parte reale e immaginaria, a leggerne il con-tenuto e a calcolare modulo, fase e coniugato del numero complessodato.

– I metodi per leggere e scrivere le parti reali e immaginarie (che pren-dono anche il nome di getters e setters) sono funzioni standard peraccedere in lettura e in scrittura ai membri privati della classe e sonotutti definiti per semplicita su un’ unica riga nel file di definizionedella classe (funzioni inline): una funzione definita in linea vieneeseguita piu rapidamente di una funzione che viene definita nel fi-le di implementazione della classe perche il compilatore sostituisceogni chiamata a questa funzione con il codice binario in linguaggiomacchina della funzione nel punto esatto del programma in cui lafunzione viene chiamata. Invece se si definisce la funzione nel filedi implementazione, il compilatore non sostituisce ogni chiamata conil codice binario, bensı con un’ istruzione di salto alla locazione dimemoria a partire dalla quale e definita la funzione: questo saltocomporta alcune operazioni supplementari di salvataggio e ripristinodi registri di memoria (variabili di sistema) che vengono effettuate adogni chiamata e che si traducono in un calo di prestazioni del pro-gramma. Con i metodi inline, il tempo morto dovuto alle chiamatedi sistema si riduce, ma la dimensione del codice eseguibile aumen-ta perche per ogni chiamata viene scritto nel programma eseguibile,ripetendolo, lo stesso codice: e dunque vantaggioso definire come me-todi inline funzioni molto brevi poiche esse non occuperanno moltospazio e saranno eseguite rapidamente.

– Per chiamare le funzioni dei vari oggetti, si usa la sintassi indicata:cioe il nome dell’ oggetto seguito da punto e dalla funzione, analo-gamente a quanto si faceva con i membri di una struttura. Oggettidiversi definiti in base alla stessa classe hanno le stesse funzioni, mai valori da queste calcolate sono diversi a seconda dei valori delle va-riabili associate agli oggetti. Si ribadisce il fatto che, a differenza diquanto succedeva con le strutture, e necessario utilizzare i metodi diclasse per accedere sia in lettura sia in scrittura alle variabili (private)membro di classe attraverso l’ uso dell’ operatore .(punto).

21

2.2 Esercizio 2.2: overload di metodi e di operatori

Aggiungere all’ esercizio precedente la possibilita di utilizzare gli operatori al-gebrici e di I/O attraverso l’ overload dei metodi e degli operatori della classeComplesso.


file Complesso.h:

#include <iostream>

class Complesso{

public:Complesso();Complesso(double a, double b);~Complesso();

void setRe(double a) {x=a;};void setIm(double b) {y=b;};double Re() {return x;};double Im() {return y;};

double Mod();double Phase();Complesso Conjugate();

Complesso operator+(const Complesso &c);Complesso operator-(const Complesso &c);Complesso operator*(const Complesso &c);Complesso operator/(const Complesso &c);Complesso& operator+=(const Complesso &c);Complesso& operator-=(const Complesso &c);

friend ostream& operator<<(ostream& o,const Complesso& c);

private:double x;double y;

};

22

file Complesso.C:

#include "Complesso.h"#include <math>

Complesso::Complesso(){

x=0;y=0;

}

Complesso::~Complesso(){}

Complesso::Complesso(double a, double b){

x=a;y=b;

}

double Complesso::Mod(){

return sqrt(x*x + y*y);}

double Complesso::Phase(){

return atan2(y,x);}

Complesso Complesso::Conjugate(){

return Complesso(x,-y);}

Complesso Complesso::operator+(const Complesso &c){

Complesso z;z.x = x + c.x;z.y = y + c.y;

return z;}...

23

...Complesso Complesso::operator-(const Complesso &c){

Complesso z;z.x = x - c.x;z.y = y - c.y;

return z;}Complesso Complesso::operator*(const Complesso &c){

Complesso z;z.x = x*c.x - y*c.y;z.y = x*c.y + y*c.x;

return z;}

Complesso Complesso::operator/(const Complesso &c){

Complesso z;z.x = (x*c.x + y*c.y)/(c.x*c.x + c.y*c.y);z.y = (y*c.x - x*c.y)/(c.x*c.x + c.y*c.y);

return z;}

Complesso& Complesso::operator+=(const Complesso &c){

x += c.x;y += c.y;

return *this;}

Complesso& Complesso::operator-=(const Complesso &c){

x -= c.x;y -= c.y;

return *this;}...

24

...ostream& operator<<(ostream& o, const Complesso& c) {o << "(" << c.x << "," << c.y << ")";return o;

}

file UsoComplesso2.C:

#include <iostream>#include "Complesso.h"

main(){

Complesso y(7,5);Complesso w(3,1);Complesso z;

z = y + w;cout << z << endl;

return 0;}


– In C++ e possibile dare piu definizioni di una funzione con lo stessonome, ma con un diverso numero e/o tipo di parametri in ingresso:il compilatore e in grado di chiamare la funzione corretta semplice-mente controllando la lista dei parametri (ma, si badi bene, non iltipo restituito) e questo fornisce al programmatore uno strumento inpiu per il controllo degli oggetti in un programma. La funzione defi-nita in questi diversi modi prende il nome di funzione sovraccaricata(overloaded) e come vedremo subito questa tecnica e ampiamenteusata per i metodi di classe.

– Il costruttore di default viene invocato al momento della dichiarazio-ne di una variabile Complesso che viene poi inizializzata con la partereale e immaginaria a zero. Per inizializzare il numero complesso adun altro valore si possono poi invocare le funzioni SetRe() e SetIm()passando loro i valori da utilizzare, ma e anche possibile inizializza-re un numero complesso ad un qualunque valore scrivendo un altrocostruttore, con lo stesso nome e con due parametri in ingresso checorrispondono alla parte reale e immaginaria:

Complesso(double a, double b);

ove la funzione, definita nel file Complesso.C, semplicemente assegnai valori in ingresso a e b ai membri di classe x e y. Poiche x e y sono

25

membri di classe e la funzione e un metodo di classe implementatocon l’ operatore ::, i membri di classe sono variabili locali di tutti imetodi di classe e non e necessario ridefinirli entro la funzione. E’possibile quindi inizializzare un numero complesso al momento delladichiarazione con un’ istruzione del tipo:

Complesso z(5.,12.);

– In C++ per ogni tipo elementare come int e float e definita nellalibreria di default la classe corrispondente con dati membro e metodidi classe. Quando si sommano due numeri tra loro con la notazione:

a+b;

in realta si chiama un metodo denominato operator+() che appartie-ne all’ oggetto corrispondente al primo addendo:

a.operator+(b);

e che restituisce come risultato la somma di b ad a. Lo stesso av-viene per gli altri operatori algebrici attraverso la chiamata ai meto-di: operator-(), operator*(), operator/(), operator+=() e operator-=(). Per il modo in cui il metodo viene chiamato, si capisce che lenotazioni:

a + b e b + a

sono profondamente diverse in quanto nel primo caso si chiama unmetodo dell’ oggetto a mentre nel secondo si chiama un metodo dell’oggetto b. Non e quindi immediato o banale trattare la proprietacommutativa tra oggetti di un programma C++.

– In base a quanto visto e dunque possibile sovraccaricare i metodi ope-rator+() e simili per la classe complesso e specificare al programmache si intende operare sia sulla parte reale sia sulla parte immaginariadegli operandi: sara poi finalmente possibile utilizzare il simbolo +(e simili) per effettuare la somma fra due numeri complessi!

– L’ argomento dei metodi sovraccaricati della classe Complesso presen-ta una novita rispetto a quanto visto finora nella chiamata a funzioni.In C e in C++ una funzione puo essere chiamata per valore o perindirizzo con le differenze che abbiamo visto nel primo semestre: inparticolare nella chiamata per indirizzo, che viene effettuata usandodei puntatori, si hanno i due vantaggi di poter modificare dalla fun-zione chiamata il valore di variabili nel programma chiamante e diottimizzare l’ esecuzione del programma nel caso in cui l’ argomentodella funzione non sia costituito da un singolo valore scalare (comecon array, strutture o classi). Oltre ai puntatori, il C++ introducele referenze per accedere al contenuto di variabili: una referenza e

26

una variabile sinonimo (o alias) di un’ altra variabile, viene dichia-rata interponendo il simbolo & tra il tipo e il nome della variabile epuo essere pensata come un nome alternativo della variabile stessa.Nella chiamata a funzione, e possibile passare come argomento del-la funzione una referenza alla variabile: la funzione lavorera alloradirettamente sulla locazione di memoria della variabile senza farneuna copia, analogamente a quanto succede con i puntatori ma conuna notazione piu semplice. Per evitare che la funzione a cui vie-ne passata una referenza modifichi il programma chiamante, si faprecedere la dichiarazione della referenza dalla parola chiave const:ogni tentativo della funzione di modificare la locazione di memoriaa cui punta la referenza genera cosı un errore di compilazione ed epiu facile intercettare errori di programmazione. L’ utilizzo di constapplicato a referenze e molto usato in C++ nella chiamata a funzioniperche presenta i vantaggi della chiamata per indirizzo e minimizzail rischio di errori difficili da individuare che si avrebbero con l’ usodi un semplice puntatore.

– Il sovraccarico degli operatori impliciti += e -= richiede l’ utilizzo diuna referenza anche nel tipo restituito dalla funzione che deve esserequindi Complesso&: questa operazione e necessaria se si vuole ren-dere l’ operatore implicito cumulabile, cioe se lo si vuole riutilizzaredue o piu volte nella medesima istruzione. Consideriamo ad esempiole seguenti istruzioni:

int a=5,b=2,c=1;a-=b-=c;

in cui, lavorando su interi, l’ operatore -= viene concatenato. Inquesto caso, dopo l’ operazione indicata, i valori delle variabili sonorispettivamente a=4, b=1 e c=1. Scrivendo invece:

(a-=b)-=c;

si avrebbero i valori a=2, b=2 e c=1. Queste istruzioni funzionanocorrettamente perche l’ operatore -= definito per tipi int nella libre-ria standard restituisce una referenza ad intero e non un intero: senon fosse cosı, a=a-b sarebbe calcolato correttamente, ma il valorerestituito sarebbe salvato in un altro oggetto ad un altro indirizzo. Aquesto valore (e non al vero a) verrebbe sottratto c, con il risultatosbagliato di avere a=3, b=2 e c=1.

– Ogni oggetto definito secondo una classe possiede un puntatore auto-matico, nascosto all’ utente, denominato this che serve per puntareall’ oggetto stesso: questo e comodo quando l’ oggetto contiene unmetodo che ne modifica lo stato e che restituisce l’ oggetto stesso mo-dificato. Nel caso di operator+=(), si vuole aggiungere un numeroall’ oggetto dato e il risultato deve essere memorizzato nell’ ogget-to stesso, cioe nella variabile puntata da this dopo l’ operazione che

27

per l’ appunto e *this (si ricorda che l’ operatore * applicato ad unpuntatore restituisce il contenuto della cella di memoria puntata).

– A causa dell’ incapsulamento, gli attributi privati di un oggetto sonoinaccessibili a oggetti appartenenti a classi diverse: tuttavia in alcunicasi questo e permesso consentendo di rendere disponibili questi at-tributi ad alcune classi fidate, definendole come amiche (friend). Leclassi definite come amiche avranno accesso diretto ai dati membridella classe entro cui questo legame e definito: inoltre si puo conce-dere l’ accesso ai membri privati di classe anche ad alcune funzioniche vanno dichiarate facendo precedere il nome del metodo dalla pa-rola chiave friend e questo permette di ottimizzare l’ I/O dei nostrioggetti Complesso. Infatti per stampare su schermo il valore di unnumero dobbiamo ancora scrivere:

cout << a.Re() << " " << a.Im() << endl;

mentre sarebbe molto piu comodo e desiderabile scrivere un’ istru-zione del tipo:

cout << a << endl;

in cui il numero complesso a e realmente trattato come un oggetto ase stante. Questo e possibile a patto di definire come friend il meto-do operator<<() (<< e uno degli operatori ridefinibili dall’ utentein C++), un metodo che viene utilizzato nella classe ostream perdefinire il simbolo cout nella stampa di un numero complesso: questedefinizioni sono contenute nell’ header file iostream.h che e necessarioincludere all’ inizio del file. La funzione, definita friend della classeComplesso, accetta due parametri: il primo, una referenza ad un og-getto di tipo ostream, rappresenta il canale di output su cui e direttal’ uscita del programma (in genere lo standard output), mentre ilsecondo e una referenza all’ oggetto da ridirigere sul canale. Nel filedi implementazione, la funzione accede direttamente agli attributi diComplesso indirizzando il numero complesso sul canale di output o:cio rende il programma molto generale in quanto l’ output puo esserereindirizzato anche su un file o su una stampante senza modificare ilcodice.

28

2.3 Esercizio 2.3: classe FrattC per mappa logistica

Utilizzando la classe Complesso, scrivere una classe che rappresenti una mappalogistica e scrivere un programma per calcolarne le prime iterazioni stampandoil risultato sullo schermo.

• Teoria

– Una mappa logistica permette di generare una successione di numeriiterando una funzione secondo lo schema:

xn+1 = fp(xn) (2)

ove il parametro p e il valore iniziale x0 sono dati e il problemaconsiste nello studiare il comportamento della successione per mezzodi un programma in C++. Nel nostro caso, siamo interessati allafunzione:

xn+1 = pxn(1− xn) (3)

che prende storicamente il nome di mappa (o funzione) logistica 1 eche genera i seguenti numeri:

x1 = px0(1− x0)x2 = px1(1− x1)x3 = px2(1− x2)...

– Dato un valore iniziale compreso nell’ intervallo [0,1], si puo mostrareche i valori della successione cadono sempre nel medesimo intervalloe che per uno stesso valore iniziale, ad esempio x0 = 0.2, la successio-ne presenta un comportamento estremamente diverso a seconda deivalori del parametro p, come riassunto dallo schema seguente:∗ p = 0.90: decadimento a zero;

Questo valore di p produce la seguente successione di numeri:0.2 0.144 0.1109 0.08877 0.0728 0.06075 0.05135 0.04384 0.037730.03268 0.02845 0.02487 0.02183 0.01922 0.01696 0.01501 0.01330.01182 0.01051 0.009358 0.008343...

∗ p = 2.60: convergenza a un valore costante;Questo valore di p produce la seguente successione di numeri: 0.20.416 0.6317 0.6049 0.6214 0.6117 0.6176 0.6141 0.6162 0.61490.6157 0.6152 0.6155 0.6153 0.6154 0.6154 0.6154 0.6154 0.61540.6154 0.6154...

∗ p = 3.20: oscillazione con periodo 2;Questo valore di p produce la seguente successione di numeri:0.2 0.512 0.7995 0.5129 0.7995 0.513 0.7995 0.513 0.7995 0.5130.7995 0.513 0.7995 0.513 0.7995 0.513 0.7995 0.513 0.7995 0.5130.7995...

1Questa funzione rappresenta bene un modello matematico di crescita di una popolazionecon un rate di crescita variabile in funzione delle dimensioni della popolazione ed e statabattezzata in questo modo dai matematici.

29

∗ p = 3.52: oscillazione con periodo 4;Questo valore di p produce la seguente successione di numeri: 0.20.5632 0.8659 0.4086 0.8506 0.4473 0.8702 0.3975 0.843 0.46580.8759 0.3827 0.8316 0.493 0.8798 0.3722 0.8225 0.5139 0.87930.3735 0.8237...

∗ p = 4.00: variazione caotica;Questo valore di p produce la seguente successione di numeri: 0.20.64 0.9216 0.289 0.8219 0.5854 0.9708 0.1133 0.402 0.9616 0.14780.5039 0.9999 0.0002463 0.000985 0.003936 0.01568 0.06174 0.23170.7121 0.82...

I vari comportamenti sono riassunti in Fig. 1 seguenti:

Figura 1: Comportamento della mappa logistica per diversi valori del parametrop: decadimento a zero, tendenza ad un valore costante, oscillazione tra due valorie oscillazione caotica.

– La sequenza di valori assunti dalla mappa prende il nome di orbi-ta, mentre i valori a cui la successione di stabilizza vengono dettiattrattori: nel caso di convergenza ad un singolo valore si parla diattrattore di periodo 1, nel caso di oscillazione tra due valori di at-trattore di periodo 2, e cosı via fino ad arrivare ad un attrattorecaotico in cui la successione non si stabilizza mai su una sequenza dinumeri ripetibile.

– E interessante rappresentare il comportamento di questa mappa inun grafico con il parametro p sull’ asse delle ascisse e i valori as-

30

sunti dalla mappa (dopo un certo numero di iterazioni iniziale dettotransiente, necessario perche il comportamento della successione sistabilizzi): si ottiene un grafico che prende il nome di diagrammadi biforcazione e che si presenta come in Fig. 2 per quanto riguardavalori del parametro nell’ intervallo [3,4].

Figura 2: Diagramma di biforcazione della funzione logistica per valori delparametro nell’ intervallo [3,4].

Questo grafico consente di esaminare l’ effetto globale del parametrodi controllo p sul comportamento della mappa logistica. Per valori dip compresi in [0,1], l’ orbita decade a zero, mentre per valori nell’ in-tervallo [1,3] converge ad un singolo valore costante che costituisce unattrattore di periodo 1. A questo punto, nell’ intervallo [3,3.449490...]si incontra un attrattore di periodo 2 e poi in [3.449490...,3.544090...]un attrattore di periodo 4 seguiti in sequenza da altri attrattori di pe-riodo 8, 16, 32 e cosı via: questo comportamento prosegue finche nonsi arriva ad un attrattore di periodo infinito per p=3.569945 in cuiil comportamento e caotico. Tra questo valore e p=4 si trova infineuna ricca varieta di attrattori, sia periodici sia caotici: ad esempio,in corrispondenza del valore p=3.828435 si ha un attrattore di peri-odo 3, come si puo vedere dalla finestra che si apre nel diagrammadi biforcazione. Inoltre il diagramma e un frattale, cioe presenta unastruttura autosimile a livelli di zoom successivi, cioe all’ aumentaredella precisione con cui si effettua il calcolo. Questa e una carat-teristica di questa funzione scoperta solamente una ventina di annifa da un ricercatore dell’ IBM che si chiamava Benoit Mandelbrot:la caratteristica principale dei frattali e di essere simili a se stessi alvariare della scala a cui li si analizza e noi tutti siamo circondati daoggetti frattali anche nella vita di tutti i giorni (si pensi al profilo di

31

una nube, a un cavolfiore, al cervello o semplicemente ad un viaggioin treno attraverso le grandi pianure centrali degli Stati Uniti (co-me in Fig. 3)) e la Fisica sta rivelando le molteplici applicazioni diquesta intuizione in vari campi. La situazione che si presenta e estre-mamente varia ed interessante, nonstante la semplicita della formulautilizzata nel nostro caso: un vasto campo di ricerca, che va oltre gliscopi di questo esercizio (si veda la bibliografia per alcuni testi in cuil’ argomento e approfondito), si apre ora davanti a noi e l’ obiettivodel programma che scriveremo in questo esercizio e di generare questesequenze di numeri, nonche di calcolare il diagramma di biforcazione.

Figura 3: L’ invarianza di scala, tipica dei frattali, non e poi cosı lontana dallanostra vita di tutti i giorni!

32


file FrattC1.h:

#include "Complesso.h"

class FrattC{

public:FrattC();~FrattC();FrattC(int N);

int d() {return durata;};void inizializza(Complesso X, Complesso P);double evolvi();void stampa();void biforca();

private:Complesso *x;Complesso p;int passo;int durata;double epsilon;

};

33

file FrattC1.C:

#include "FrattC.h"#include <iostream>#include <iomanip>

FrattC::FrattC(){}

FrattC::~FrattC(){

if(x!=NULL){delete[] x;x=NULL;}

}

FrattC::FrattC(int N){

durata=N;x=new Complesso[N];epsilon=1e-6;

}

void FrattC::inizializza(Complesso X, Complesso P){

passo=0;x[passo]=X;p=P;

}

double FrattC::evolvi(){

x[passo+1]=p*x[passo]*(Complesso(1.,0.)-x[passo]);

passo++;return( x[passo].Re() );

}...

34

...void FrattC::biforca(){double y;

for(double par=2.9;par<=4.;par+=0.001){Complesso p(par,0.);

inizializza(Complesso(0.2,0.),p);for(int i=1;i<durata;++i){

y=evolvi();if(i>50)

cout << p.Re() << "\t" << y << endl;}

}}

void FrattC::stampa(){

for(int i=0;i<durata;++i)cout << setprecision(4) <<setw(8) << x[i] << "\n";

}

file UsoFrattC1.C:

#include "Complesso.h"#include "FrattC.h"#include <iostream>

main(){

Complesso Z(0.2,0.), R(3.00,0.);FrattC a(100);a.inizializza(Z, R);for(int i=0;i<a.d()-1;++i)

a.evolvi();

cout << "Premi RETURN per finire..." << endl;cin.get();return 0;

}

35


– La classe presentata, FrattC, permette di studiare numericamente lesuccessioni generate dalla mappa logistica calcolandone un numeroarbitrario di valori (limitato dalla memoria del computer e dai tempidi calcolo) e di stampare il risultato sullo schermo o in un file peruna successiva rappresentazione grafica. Poiche nel seguito estende-remo la classe anche all’ uso di numeri complessi con lo studio dellerappresentazioni complesse di altre funzioni iterate (insiemi di Juliae di Mandelbrot), la classe utilizza oggetti di tipo Complesso la cuiparte immaginaria per il momento viene posta a zero.

– I membri privati di classe sono costituiti da:

*xun puntatore a Complesso *x che viene allocato nel costrutto-re dell’ oggetto e che costituisce un array dinamico di numericomplessi, usato per memorizzare i valori della successione;

pun numero Complesso p che rappresenta il parametro (per orareale, vedi esercizio seguente);

passodurata

due interi per memorizzare il passo di iterazione corrente e ilnumero totale di iterazioni (durata);

εun numero reale ε per memorizzare la precisione a cui si vuo-le arrivare con il calcolo: si puo decidere di fermarsi nel calco-lo o dando un numero totale di passi (durata) o decidendo peresempio in base ad un criterio del tipo |xn+1 − xn| < ε.

– Per quanto riguarda i metodi, si hanno due costruttori (quello di de-fault non viene praticamente usato), il distruttore, una funzione diinizializzazione e una di evoluzione, noche una funzione per la stampadei valori della successione. Il metodo biforca viene usato per studia-re il comportamento della successione al crescere del parametro p eserve per calcolare i valori di p in corrispondenza dei quali il compor-tamento cambia (ad esempio da un valore limite a una oscillazionetra 2 o 4 o piu valori) e per produrre il diagramma di biforcazione.

– Il secondo costruttore accetta un intero N in ingresso che specifica ilnumero di passi che l’utente vuole eseguire e che viene memorizzatonella variabile durata, alloca in modo dinamico memoria per un arraydi N variabili di tipo Complesso e inizializza anche il parametro dellaprecisione di calcolo.

– La funzione inizializza() ha il compito di memorizzare il valore delparametro p e di calcolare il primo passo dell’ iterazione: si noti cheabbiamo usato lettere maiuscole per le variabili in ingresso (X e P)

36

e variabili minuscole per i membri di classe (x e p). Si sarebberopotute scrivere queste istruzioni all’ interno del costruttore in cuiviene allocata la memoria e questa e una scelta del programmatore,ma per il modo in cui la classe e utilizzata si vedra che e piu conve-niente procedere cosı. Ovviamente non si tratta dell’ unico modo dirisolvere il problema e la scelta deve essere fatta in fase di design delprogramma: non esistono metodi da seguire in questo caso, poichesolamente l’ esperienza e la pratica, unite ad un’ idea, possono pro-durre un codice efficiente e facile da programmare e da usare oppureno. La progettazione di una classe e di un programma ad oggetti euna fase molto importante e delicata e alcuni suggerimenti verrannodati nella terza parte del corso accennando alla metodologia UML(Unified Modeling Language).

– La funzione evolvi() esegue un singolo passo di iterazione calcolandoil valore xn+1 a partire dal valore di xn, incrementando di uno ilpasso corrente e restituendo il nuovo valore calcolato (come partereale di numero complesso). La funzione cosı scritta e pensata peressere chiamata dall’ interno di un ciclo in un programma main()per un numero di volte pari alla durata voluta dall’ utente: a questoproposito, la funzione inline d() serve proprio a restituire il valoredella durata di un oggetto di tipo FrattC dopo che questo sia statoinizializzato.

– La funzione biforca() e costituita da un ciclo for sul parametro pche viene utilizzato per calcolare il diagramma di biforcazione: all’interno di questo ciclo l’ oggetto di tipo FrattC viene inizializzato efatto evolvere un certo numero di volte per ogni valore del parametrop e, dopo un transiente di 50 iterazioni, vengono stampati il valore dip e i valori della successione corrispondenti. Si noti che il parametro ereale, ma viene inserito (embedded) in un oggetto di tipo complessoe trattato dal codice come un numero complesso attraverso la nostraclasse per motivi di generalita che risulteranno piu chiari nel prossimoesercizio. Inoltre l’ oggetto usato per memorizzare la successione esempre lo stesso e la medesima memoria viene riutilizzata ad ogniciclo per scrivere i valori della successione.

– Infine la funzione stampa() e costituita da un ciclo su tutti gli ele-menti dell’ array x membro di FrattC e dalla successiva stampa ditutte le componenti. Viene fatto uso di alcuni manipolatori dell’input/output, definiti nel file iomanip.h incluso all’ inizio del pro-gramma: questi servono a formattare l’ output, ossia a scegliere ilformato dei numeri stampati su schermo. In particolare, setpreci-sion(4) serve a scegliere di visualizzare i valori con 4 cifre decimali,mentre setw(8) serve a fissare a 8 caratteri la larghezza (width) delcampo entro cui ogni numero verra scritto.

– Una nota di stile: il parametro ε che avevamo previsto di utilizzarein fase di progettazione non viene in realta utilizzato. Questo e stato

37

fatto intenzionalmente per mostrare come un programma possa evol-vere con il nascere di nuove idee e l’ utilizzo di ε per decidere quandofermarsi nel calcolo della successione (invece di fissare sempre 100passi come abbiamo fatto) e lasciato come esercizio.

38

2.4 Esercizio 2.4: classe FrattC per insiemi di Julia eMandelbrot

Includere nella classe precedente, preservandone il piu possibile la struttura, lapossibilita di iterare una mappa complessa che generi degli insiemi frattali nelpiano complesso noti come insiemi di Julia e di Mandelbrot.

• Teoria

– Consideriamo un’ altra mappa da iterare:

zn+1 = z2n + c

dove z e un generico numero complesso z = a + ib e c e anch’ essocomplesso. Come nel caso della mappa logistica, la soluzione di que-sta mappa complessa esibisce una grande varieta di comportamenti aseconda dei valore del parametro complesso c e del punto iniziale z0 eun buon punto di partenza consiste nell’ investigare il comportamentodella successione di numeri per alcune condizioni iniziali notevoli.

– Fissiamo c al valore (-1,0) e studiamo le due successioni generate apartire da z0 = (0, 0) e z0 = (0, 1): la successione che parte da (0,0)produce un’ orbita periodica che oscilla tra due valori (0,0) e (-1,0),mentre la successione che parte da (0,1) continua ad aumentare diampiezza e diverge verso l’ infinito. Facendo iniziare la successioneda altri punti, si trova che in generale per alcuni valori essa oscillatra 2 valori e che per altri valori essa tende all’ infinito: e possibiledisegnare nel piano complesso i punti per i quali la successione nonfugge all’ infinito e l’ insieme che si ottiene (il cosiddetto insiemeprigioniero della mappa) e rappresentato in Fig. 4 . L’ insieme deipunti sul bordo di questa regione del piano prende il nome di insiemedi Julia e si puo mostrare che per alcuni valori del parametro c essoe connesso e per altri e disconnesso, come si vede dalla figura per idue casi esaminati).

– L’ iterazione della mappa complessa per determinare l’ insieme diJulia per un dato parametro c richiede un notevole sforzo computa-zionale (come vedremo eseguendo il programma) e spesso e sufficientedeterminare se un dato valore di c produrra un insieme connesso op-pure un insieme disconnesso: e sufficiente allora iterare la mappa perogni valore di c partendo sempre dal punto z0 = (0, 0). Se il puntoscappa all’ infinito, l’ insieme di Julia corrispondente a quel valoredi c e disconnesso, mentre insiemi di Julia connessi corrispondonoa orbite finite. Disegnando nel piano complesso c tutti i punti cheportano ad insiemi di Julia connessi, si ottiene il cosiddetto insie-me di Mandelbrot: il bordo di questo insieme di punti presenta unastruttura infinitamente complessa, di natura frattale, come illustratoin Fig. 5. Successivi zoom nella regione del bordo rivelano sempre

39

ulteriori dettagli e assegnando un diverso colore al punto a secondadella periodicita dell’ orbita permette di generare delle figure moltobelle che sono anche state definite artistiche 2: sono disponibili inrete molti programmi per la generazione di frattali sia in ambienteLinux sia in ambiente Windows e fra questi si segnala il program-ma Fractint (molto veloce tra l’ altro per il fatto di effettuare tutti icalcoli in aritmetica intera!).

2Con l’ avvento dei primi personal computers negli anni ottanta, l’ applicazione di questialgoritmi porto ad un boom della computer graphics: un ottimo libro, ricco anche di fotografiee disegni, e il Peitgen-Richter, La bellezza dei Frattali, Bollati-Boringhieri

40

Figura 4: Due casi estremi di insiemi di Julia connesso (in alto) e disconnesso(in basso), corrispondenti rispettivamente ai valori del parametro c=(-1,0) ec=(0.5,0.1).

41

Figura 5: L’ insieme di Mandelbrot nel piano complesso c e zoom successivi chene dimostrano la natura frattale.

42


file FrattC2.h:


class FrattC{


int d() {return durata;};void inizializza(Complesso X, Complesso P);double evolvi(int flag);void stampa();void biforca();


};

43

file FrattC2.C:

#include "FrattC.h"#include <iostream>#include <iomanip>

FrattC::FrattC(){}FrattC::~FrattC(){

if(x!=NULL){delete[] x; x=NULL;}

}

FrattC::FrattC(int N){


}

void FrattC::inizializza(Complesso X, Complesso P){


}

double FrattC::evolvi(int flag){

if(flag==1){

x[passo+1]=p*x[passo]*(Complesso(1.,0.)-x[passo]);passo++;return( x[passo].Re() );

} else if(flag==2){

x[passo+1] = x[passo]*x[passo] + p;passo++;return( x[passo].Mod() );

}}...

44

...void FrattC::biforca(){double y;



y=evolvi(1);if(i>50)

cout << p.Re() << "\t" << y << endl;}

}}


for(int i=0;i<durata;++i)cout << setprecision(4) <<setw(8) << x[i] << "\n";

}

file UsoFrattC2.C:


main(){

//Complesso R(-0.5,0.); //JuliaComplesso R(0.,0.); //MandelbrotFrattC b(100);for(int i=0;i<B.d()-1;++i)b.evolvi(2);


}

45


– Sviluppiamo la classe FrattC scritta precedentemente per generare lesuccessioni di numeri degli insiemi di Julia e di Mandelbrot, utiliz-zando a pieno gli oggetti Complesso e mantenendo il piu possibile lastruttura della classe utilizzata: questo e un esempio, seppur moltosemplice, di reuse cioe di riutilizzo di codice gia scritto ed e una tec-nica di programazione molto usata anche in C++. Si tratta quindidi modificare la classe aggiungendo l’ algoritmo della nuova mappada iterare: questo comporta solamente una modifica alla funzioneevolvi per poter scegliere quale mappa iterare, la mappa logistica ola mappa di Julia/Mandelbrot.

– La nuova funzione evolvi accetta in ingresso un parametro intero, det-to flag, che viene utilizzato per decidere quale delle due mappe iterareattraverso un’ istruzione if: se flag=1 viene iterata la mappa logistica,mentre se flag=2 viene invece iterata la mappa di Julia/Mandelbrot.Il valore restituito dal metodo e il nuovo valore iterato (numero reale)per la mappa logistica oppure il modulo del nuovo valore (complesso)per la mappa di Julia: questo e utile perche, come abbiamo visto, l’insieme di Julia e definito dai punti che non divergono all’ infinito chesono identificabili controllando il valore del loro modulo. Inoltre lafunzione restituisce sempre un valore double in entrambi i casi. Tut-te le altre funzioni della classe rimangono inalterate: in particolarela funzione biforca() viene usata solo nel caso della mappa logistica,mentre l’ inizializzazione e la stampa sono comuni ai due casi.

– Il risultato viene ancora stampato in forma numerica come successio-ne di numeri complessi: nel prossimo esercizio verra mostrato comeaggiungere a questo programma chiamate a funzioni grafiche in mo-do da visualizzare gli insiemi di Julia e di Mandelbrot, nonche ildiagramma di biforcazione della mappa logistica.

46

2.5 Esercizio 2.5: libreria di ROOT e output grafico

Utilizzando gli oggetti della libreria grafica di ROOT, migliorare l’ output delprogramma precedente visualizzando in modo grafico i seguenti punti: compor-tamento della mappa logistica, diagramma di biforcazione, insieme di Julia einsieme di Mandelbrot.


file FrattC3.h:


class FrattC{


int d() {return durata;};void inizializza(Complesso X, Complesso P);double evolvi(int flag);double calcola();void calcola(Complesso X);void stampa();void biforca();


};

47

file FrattC3.C:

#include "FrattC.h"#include <iostream>#include <iomanip>#include <TCanvas.h>#include <TF1.h>#include <TH2.h>#include <TLine.h>FrattC::FrattC(){}FrattC::~FrattC(){

if(x!=NULL){delete[] x; x=NULL;}

}FrattC::FrattC(int N){


}void FrattC::inizializza(Complesso X, Complesso P){


}double FrattC::evolvi(int flag){

if(flag==1){

x[passo+1]=p*x[passo]*(Complesso(1.,0.)-x[passo]);passo++;return( x[passo].Re() );

} else if(flag==2){

x[passo+1] = x[passo]*x[passo] + p;passo++;return( x[passo].Mod() );

}}...

48

...TCanvas *c1;TF1 *l, *r;TLine *l1,*l2;

double FrattC::calcola(){

double y;char logistica[20];

if(c1==NULL) c1=new TCanvas("c1","Frattali",600,600);

c1->Draw();c1->cd();

x[passo+1]=p*x[passo]*(Complesso(1.,0.)-x[passo]);

sprintf(logistica,"%lf*x*(1.-x)",p.Re());l=new TF1("logistica",logistica,0.,1.);l->SetLineColor(2);r=new TF1("linea","x",0.,1.);r->SetLineColor(4);

if(passo==0) y=0; else y=x[passo].Re();l1=new TLine(x[passo].Re(),y,

x[passo].Re(),x[passo+1].Re());l1->SetLineWidth(2);l2=new TLine(x[passo].Re(),x[passo+1].Re(),

x[passo+1].Re(),x[passo+1].Re());l2->SetLineWidth(2);

r->Draw("same");l->Draw("same");l1->Draw("same");l2->Draw("same");

c1->Update();

passo++;return( x[passo].Re() );

}...

49

...TH2F *h2;

void FrattC::calcola(Complesso X){

double y,z,w;

if(c1==NULL) c1=new TCanvas("c1","Frattali",600,600);

c1->Draw();c1->cd();if(h2==NULL){//h2=new TH2F("h2","Julia",

1000,-2.,2.,1000,-1.,1.);h2=new TH2F("h2","Julia",

1000,-2.,0.5,1000,-1.,1.);h2->SetStats(0);

}

for(y=-2.;y<=0.5;y+=.01){for(z=-1.;z<=1.;z+=.01){//Julia//inizializza(Complesso(y,z),X);//Mandelbrotinizializza(X,Complesso(y,z));for(int i=0;i<d()-1;++i)w=evolvi(2);

if(w<100.) h2->Fill(y,z);}

}h2->Draw();c1->Update();

}...

50

...TCanvas *c2;TH2F *h3;

void FrattC::biforca(){double y;

if(c2==NULL) c2=new TCanvas("c2","Biforcazione",800,600);

if(h3==NULL){h3=new TH2F("h3","Biforcazione",

1000,2.9,4.1,1000,0.,1.);h3->SetStats(0);

}



y=evolvi(1);c2->cd();if(i>50)

//cout << p.Re() << "\t" << y << endl;h3->Fill(p.Re(),y);

}h3->Draw();

}c2->Update();

}


for(int i=0;i<durata;++i)cout << setprecision(4) << setw(8)<< x[i] << "\n";

}

51

file UsoFrattC3.C:


main(){

//Complesso R(-0.5,0.); //JuliaComplesso R(0.,0.); //MandelbrotFrattC b(100);for(int i=0;i<b.d()-1;++i)b.calcola(R);


}


– Con questo esercizio si vuole mostrare come sia possibile, utilizzan-do un’ opportuna libreria di funzioni esterna, aggiungere un outputgrafico ad un programma preesistente: si introdurranno altri oggettiC++ presenti in questa libreria che verranno utilizzati per rappre-sentare in modo grafico il comportamento della mappa logistica, ildiagramma di biforcazione e infine per disegnare (anche se solo inbianco e nero per semplicita) gli insiemi di Julia e di Mandelbrot.

– La libreria utilizzata e quella di ROOT, un tool per l’ analisi da-ti sviluppato presso il laboratorio CERN di Ginevra e attualmenteutilizzato in vari progetti di Fisica. Altre librerie sono disponibili fa-cilmente in internet, ma non verranno presentate per brevita: quellaqui proposta vuole essere solo un esempio di come procedere per dareai propri programmi un output grafico in ambiente Unix.

– Seguendo la filosofia del reuse, estendiamo la classe FrattC aggiun-gendo le nuove funzionalita grafiche: molte strade possono essereseguite a questo punto e la via presentata e solo una di queste, mail lettore e invitato a pensarne e sperimentarne altre utilizzando glistrumenti C++ oramai noti.

– La modifica della classe consiste nell’ aggiunta di due metodi sovrac-caricati tra di loro e denominati calcola(), uno per la mappa logisticae l’ altro per la mappa di Julia/Mandelbrot, e nella modifica del me-todo biforca per la visualizzazione grafica dei risultati. Analizziamoprima il metodo calcola() per la funzione logistica considerando 4

52

nuovi oggetti grafici della libreria di ROOT: TCanvas, TF1, TLi-ne e TH2. Per poter utilizzare questi oggetti dobbiamo ricordareche ognuno di essi e definito da una classe e che quindi sara necessarioincludere il corrispondente header file all’ inizio del programma: sitratta di header files di ROOT che possono essere trovati usando leparentesi acute come per gli header di sistema (a patto che ROOTsia stato installato correttamente!). L’ oggetto TCanvas serve adisegnare il foglio (canvas in inglese) su cui gli altri oggetti graficiverranno disegnati e viene definito sullo heap attraverso un punta-tore definito in modo globale e denominato c1: il costruttore dell’oggetto TCanvas accetta 4 parametri in ingresso secondo il prototipo

TCanvas("c1","Frattali",600,600);

ove le prime due stringhe definiscono il nome e il titolo del canvas,mentre i due numeri danno la dimensione in x e in y del canvas in pi-xel. Nome e titolo dell’ oggetto servono perche ROOT e organizzatosecondo una struttura gerarchica di classi e ogni oggetto della libreriaha tra i suoi metodi una stringa che ne contiene il nome e una strin-ga che ne contiene il titolo: questi attributi, apparentemente di pocaimportanza, vengono invece utilizzati in modo avanzato per otteneredelle informazioni sugli oggetti contenuti in memoria al runtime (conla tecnica del RTTI, Run Time Type Information, che esula dai nostriobiettivi). L’ oggetto, che viene definito solo una volta controllandoche il puntatore c1 non sia nullo, viene poi disegnato con il metodoDraw() e selezionato con il metodo cd(): un’ oggetto di tipo TCanvasha moltissimi metodi di classe e diversi dati membro che a noi non ser-vono per questo esercizio, ma che sono descritti sul manuale online diROOT accessibile al sito http://root.cern.ch/root/Categories.html.Un altro oggetto molto importante e denominato TF1 e rappresen-ta il grafico di una funzione unidimensionale. Il costruttore di unoggetto TF1 e definito come segue:

TF1("linea","x",0.,1.);

ove la prima stringa rappresenta il nome dell’ oggetto e la secondaserve a passare all’ oggetto la definizione della funzione da disegnare(in questo caso la retta y=x). Nel caso della mappa logistica e neces-sario stampare all’ interno della stringa di definizione della funzioneanche il valore numerico del parametro p e a questo scopo viene fattouso della funzione di libreria C sprintf() che funziona come una printfo una fprintf, ma che stampa il risultato nella stringa che viene pas-sata come primo argomento: si ribadisce con questo esempio che Ce C++ possono essere integrati benissimo tra loro dal momento chealcune cose vengono risolte ancora in modo efficiente dal pur supera-to linguaggio C! Tra i vari metodi della classe TF1, viene fatto usodella funzione SetLineColor(n), con n intero, per cambiare il colore

53

della funzione disegnata: in questo caso 2 per rosso e 4 per blu. Sinoti infine che le due funzioni, denominate l e r, sono ancora definitesullo heap attraverso un puntatore globale.Veniamo ora alla parte animata della nostra grafica, cioe l’ utilizzo dilinee per seguire in modo grafico l’ evoluzione della mappa logisticaattraverso il plot della retta y=x inclinata di 45 gradi. Comincia-mo dall’ asse delle x con un valore x0 e proiettiamolo in su fino adintersecare la funzione logistica, calcolando in questo modo il valoresuccessivo della successione x1. Per trovare il terzo valore, cioe perporre x1 sull’ asse delle x, ne proiettiamo il valore sulla retta y=x ericalcoliamo il valore dalla funzione logistica: la procedura viene ri-petuta per tutti gli altri punti della successione ed ha come risultatoi grafici mostrati in Fig. 6 per valori del parametro p=0.90, 2.60,3.20 e 4.00.Il costruttore dell’ oggetto TLine deve essere invocato con la seguentesintassi:

TLine(x1,y1,x2,y2);

in cui bisogna specificare le coordinate del punto inziale (x1, y1) e delpunto finale (x2, y2) tra cui tracciare la linea. Tra i vari metodi di-sponibili in questa classe utilizziamo SetLineWidth(n), con n intero,per disegnare le linee con uno spessore diverso da quello usato per lefunzioni per maggiore chiarezza del grafico. Si noti inoltre l’ utilizzodi una if per tracciare la prima iterazione a partire dall’ asse dellex (cioe con y=0). Le due funzioni e le due linee vengono successi-vamente disegnate sul foglio e dal momento che vogliamo disegnarlesullo stesso foglio dobbiamo utilizzare il metodo Draw() con l’ opzio-ne same): tutte le opzioni di Draw() possono essere studiate sempresul manuale di ROOT online. Infine per aggiornare il display ad ognipasso facciamo uso del metodo di TCanvas Update() che ridisegnasul foglio tutti gli oggetti: la funzione infine termina esattamente co-me la corrispondente funzione evolvi incrementando di un’ unita ilpasso e restituendo il valore della successione calcolato.

– La rappresentazione grafica del diagramma di biforcazione viene trat-tata in una versione modificata del metodo biforca in cui viene de-finito un nuovo TCanvas (di dimensioni diverse dal precedente) e incui viene introdotto l’ ultimo degli oggetti grafici della libreria diROOT che esamineremo, l’ oggetto TH2F, che rappresenta un isto-gramma bidimensionale: un istogramma bidimensionale puo esserepensato come un foglio di carta millimetrata sovrapposto al canvas eserve sostanzialmente per definire un sistema di riferimento cartesia-no per riferirsi ai vari punti del piano attraverso delle cellette di unagriglia di cui possiamo scegliere il passo. Nel nostro caso vogliamorappresentare il diagramma di biforcazione nella sua regione piu inte-ressante, per valori di p compresi tra 2.9 e 4.1 ad esempio: definiamo

54

Figura 6: Evoluzione grafica della mappa logistica nei casi di convergenza a zero,convergenza ad un valore fissato, oscillazione bistabile e oscillazione caotica.

dunque il costruttore specificando nome e titolo dell’ oggetto seguitidal numero di divisioni della griglia sull’ asse delle ascisse e dai limitiinferiore e superiore del grafico e dalle stesse informazioni per l’ assedelle ordinate

TH2F("h3",Biforcazione",1000,2.9,4.1,1000,0.,1.);

In questo caso abbiamo quindi una griglia con 1000 divisioni sull’ assedelle x nell’ intervallo [2.9,4.1] e con 1000 divisioni sull’ asse delle ynell’ intervallo [0,1]. Il metodo SetStats(0) serve per non visualizzarealcuni parametri statistici relativi all’ istogramma che per ora nonutilizzeremo quali ad esempio medie e varianze lungo gli assi (si vedail manuale di ROOT). La struttura della funzione e del tutto iden-

55

tica all’ esercizio precedente, ma invece di stampare il risultato sulloschermo, accendiamo un pixel sull’ istogramma, cioe carichiamo unvalore nella celletta corrispondente alle coordinate del punto calcola-to colorandola cosı in nero: questo e il ruolo del metodo Fill(x,y), cioecaricare un valore nell’ istogramma in corrispondenza delle coordina-te (x,y). L’ istogramma viene infine disegnato con Draw() e il foglioaggiornato con Update(), esattamente come prima, con il risultatofinale riportato in Fig. 7.

– La rappresentazione grafica degli insiemi di Julia e Mandelbrot vienegestita dal metodo calcola(Complesso X) che sovrascrive il metodocalcola(): ancora definiamo un nuovo canvas e un nuovo istogrammabidimensionale per visualizzare i due insiemi, denominati rispettiva-mente c1 e h2. Vogliamo una funzione generale che possa visualizzarei due insiemi e ci garantiamo questa generalita mantenendo del codicecommentato da utilizzare per l’ uno o per l’ altro caso: questa e un’abitudine dei programmatori molto comune, lasciare delle istruzioninel codice commentandole per un uso futuro. Consideriamo anzitut-to l’ algoritmo per calcolare l’ insieme di Julia: la funzione ricevein input il numero complesso X che rappresenta il parametro c (daconsiderarsi fissato) e il calcolo si basa su due cicli sulle variabili y ez annidati uno nell’ altro, in cui y e z rappresentano la parte reale eimmaginaria del punto iniziale della successione. L’ oggetto FrattCviene inizializzato con questi valori e fatto evolvere per tutta la suadurata fornita dal metodo d() (si noti che l’ indice parte da zero!: ilmetodo evolvi(n) viene invocato con n=2 e restituisce quindi il mo-dulo del nuovo valore calcolato per viene utilizzato per capire se lasuccessione diverge oppure no confrontandolo con 100. Se il valoresi mantiene finito la corrispondente cella nel grafico bidimensionaleviene riempita con il metodo Fill() e il calcolo prosegue in modo ana-logo per tutti i punti della griglia: il risultato tipico che si ottiene perl’ insieme di Julia con c=(-0.5,0) e riportato in Fig. 8.

– Per quanto riguarda l’ insieme di Mandelbrot la struttura logica delprogramma e esattamente identica tranne che ora il parametro Xin ingresso alla funzione va considerato come il punto iniziale dellasuccessione e i due cicli devono essere interpretati come sulla par-te reale e immaginaria del parametro c. Questo si riflette bene nelmodo in cui l’ oggetto viene inizializzato, chiamando la funzione congli argomenti invertiti rispetto al caso dell’ insieme di Julia: infatti,ricordando le definizioni, l’ insieme di Julia si ottiene iterando su za c fissato, mentre l’ insieme di Mandelbrot iterando su c a z fissa-to. La parte grafica resta anch’ essa invariata e quindi vediamo chelo stretto legame matematico tra l’ insieme di Julia e l’ insieme diMandelbrot si traduce in un numero minimo di modifiche da fare alprogramma per calcolare l’ uno o l’ altro insieme. Una particolaritasu cui vale la pena di riflettere e anche l’ esistenza di un legame tra

56

l’ insieme di Mandelbrot e la mappa logistica che si puo dimostra-re (magari anche usando questi programmi e lo sviluppo di quest’idea viene lasciato come esercizio): se percorriamo l’ asse delle ascis-se sull’ insieme di Mandelbrot (Fig. 9) incontriamo delle gemme incorrispondenza di punti che prendono il nome di punti di gemmazio-ne e questi stessi punti coincidono con le ascisse del diagramma dibiforcazione percorso con l’ asse delle ascisse cambiato di segno in cuicambia il ciclo dell’ orbita della mappa logistica! In piu i punti in cuicompare una piccola copia dell’ insieme di Mandelbrot totale cadonoin corrispondenza delle finestre del diagramma di biforcazione che neinterrompono il velo caotico.

– Se ROOT e stato installato correttamente, la libreria in cui si tro-vano gli oggetti grafici descritti e i relativi header files si trovanonella directory /usr/local/root: per creare il file eseguibile dobbiamocomunicare queste informazioni al compilatore e al linker e il modomigliore per farlo e di scrivere un Makefile e di compilare poi conil comando Unix make. Il file proposto da utilizzare per compilareil nostro programma eseguibile con output grafico potrebbe essere ilseguente:file Makefile:

ROOTLIBS = $(shell root-config --libs)ROOTGLIBS = $(shell root-config --glibs)#Fratt++: Fratt++.o FrattC.o Complesso.o

g++ -o Fratt++ Fratt++.o FrattC.o Complesso.o$(ROOTLIBS) $(ROOTGLIBS)

#Fratt++.o: Fratt++.C

g++ -c Fratt++.C#FrattC.o: FrattC.C

g++ -c FrattC.C#Complesso.o: Complesso.C

g++ -c Complesso.C

ove abbiamo chiamato Fratt++.C la versione finale e definitiva delnostro programma main(). Si tratta di definire due simboli logicientro il makefile che contengono il risultato di un’ utility di ROOTche fornisce il nome di tutte le librerie nel formato voluto dal linkercon un’ istruzione del tipo root-config –libs che viene eseguita almomento del link.

57

Figura 7: Il diagramma di biforcazione della mappa logistica nell’ intervallo [3,4]generato dal programma Fratt++ con la grafica di ROOT.

58

Figura 8: L’ insieme di Julia corrispondente a c=(-0.5,0) generato da Fratt++con ROOT.

59

Figura 9: L’ insieme di Mandelbrot generato da Fratt++ in bianco e nero conROOT.

60

3 Sistema Solare

La terza parte di queste esercitazioni fa riferimento ad un corso di teoria diprogrammazione Object Oriented tenuto presso l’ Universita La Sapienza diRoma dal Prof. G. Organtini e disponibile in rete. Il materiale ivi presentatoe stato riorganizzato per praticita in esercizi progressivi che si propongono l’ambizioso compito di scrivere in C++ una simulazione del nostro Sistema Solarepartendo dalla descrizione del moto bidimensionale dei pianeti sotto l’ azionedella forza gravitazionale di Newton.

Nei primi due esercizi si riprende l’ argomento della progettazione di unaclasse e dello sviluppo della sua implementazione, ripercorrendo la traccia pre-cedentemente fornita dalla classe Complesso, ma approfondendo alcuni aspettiavanzati. In particolare viene illustrato un utilizzo piu avanzato dei costruttorinella copia di oggetti in oggetti dello stesso tipo attraverso il costruttore copia.

Successivamente vengono definite due nuove classi, CorpoCeleste e Sistema-Solare: la prima classe ha tra i dati membri vettori bidimensionali per specificareposizione e velocita di un pianeta, nonche uno scalare per la massa. Tra i metodisi ha una funzione che calcola la forza agente sul pianeta a causa del sole e ditutti gli altri pianeti e una funzione di evoluzione in cui l’ equazione di motodel pianeta viene integrata in modo approssimato ed elementare: il risultato,ottenuto prima in modo numerico e poi in modo grafico attraverso la libreriaROOT a partire da condizioni iniziali realistiche ottenute da tavole astronomi-che, e il moto dei pianeti come lo conosciamo e il simulatore potrebbe essereutilizzato anche per verificare le leggi di Keplero. L’ implementazione utilizza ascopo illustrativo alcuni oggetti di tipo container della libreria STL (StandardTemplate Library) e ne illustra anche l’ utilizzo attraverso degli iteratori. Conquesti esercizi si vuole sottolineare l’ importanza di un buon design di una classeutilizzando come ausilio il linguaggio grafico UML (Unified Modeling Langua-ge), a cui viene fatto qualche breve cenno finalizzato al nostro scopo pratico dirappresentare un moto bidimensionale in un piano cartesiano.

L’ ultimo esercizio riguarda l’ ereditarieta tra classi e il polimorfismo in C++ed ha come obiettivo la scrittura di una classe Sonda derivata da CorpoCeleste:una sonda puo essere pensata come un piccolo pianeta dotato di razzi, cioesoggetto sia alla forza gravitazionale degli altri pianeti sia alla spinta eventuale dirazzi e il suo moto puo essere calcolato con un’ opportuna funzione di evoluzionetemporale che tenga conto di queste due forze. L’ applicazione della classepermette di lanciare un oggetto Voyager nel nostro sistema solare ed al limiteanche di simulare una delle missioni Apollo (la sfida e lanciata... ;-)).

Non vengono trattati argomenti squisitamente computazionali, come la ri-soluzione di equazioni differenziali e l’ integrazione delle equazioni di moto, peri quali si rimanda a corsi e a testi piu avanzati.

61

3.1 Esercizio 3.1: classe Vettore e suo uso

Scrivere una classe che rappresenti un vettore bidimensionale sulla base carte-siana.


file Vettore1.h:

#ifndef _VETTORE#define _VETTORE

#include <math.h>#include <iostream>

class Vettore {public:Vettore();Vettore(double a, double b);Vettore(const Vettore& v); // costruttore copia

void setX(double a) { x = a; };void setY(double b) { y = b; };double X() { return x; };double Y() { return y; };double mod() { return sqrt(x*x+y*y); };double phase() { return atan2(y,x); };Vettore unit();

protected:double x;double y;

};

#endif

62

file Vettore1.C:

#include "Vettore.h"

Vettore::Vettore() {x = 0;y = 0;

}

Vettore::Vettore(double a, double b) {x = a;y = b;

}

Vettore::Vettore(const Vettore& v) {x = v.x;y = v.y;

}

Vettore Vettore::unit() {Vettore s(x/mod(), y/mod());return s;

}


– La classe Vettore viene definita in modo analogo alla classe Com-plesso sviluppata nella seconda parte del corso e in queste note sivogliono approfondire alcuni aspetti non ancora incontrati, ma im-portanti nella scrittura di una classe: tra questi sono rilevanti lacompilazione selettiva attraverso particolari direttive al preprocesso-re, il costruttore copia (copy constructor) e la possibilita di una classedi autoreferenziarsi (ad esempio attraverso un metodo di classe chechiama il costruttore della stessa classe).

– La definizione di una classe scritta in un header file viene inclusa suc-cessivamente nel file di implementazione ed eventualmente anche inaltri files da compilare: e opportuno assicurarsi che la classe sia defi-nita una ed una sola volta ed evitare ripetizioni e messaggi di errori infase di compilazione dovuti a dichiarazioni multiple. Il preprocessorefornisce la direttiva #ifndef per effettuare la compilazione selettivadegli header files di classi in modo da garantire che il codice relati-vo alla definizione di una classe sia passato al compilatore una solavolta: #ifndef controlla se un simbolo logico, che prende il nome ditag e che viene in genere indicato con caratteri maiuscoli e precedutoda un carattere di underscore, sia stato definito e in caso negativo lodefinisce assieme alla dichiarazione della classe come qui di seguito:

63


class Vettore{

...}#endif

Una direttiva #endif e necessaria alla fine della definizione dellaclasse. Il tag risulta definito dopo la prima dichiarazione della classee anche se l’ header file viene incluso piu volte la classe non vieneridefinita piu volte.

– Tra i costruttori della classe Vettore compare un nuovo tipo di co-struttore che prende come argomento una referenza costante ad unoggetto di tipo Vettore e che prende il nome di costruttore copia. InC++ la copia di un oggetto viene effettuata dal compilatore copiandomembro a membro tutti i dati membro incapsulati nell’ oggetto e cionon comporta problemi per oggetti con membri statici come i nostritipi Complesso e Vettore. Ci sono tuttavia problemi a copiare ogget-ti di una classe i cui membri sono definiti come puntatori: infatti inquesto caso solo il puntatore viene copiato e non i dati da esso defe-renziati, a meno di definire un copy constructor che specifica che cosafare. In questo caso l’ oggetto contiene membri statici e il costrut-tore copia non sarebbe necessario, ma viene definito comunque perillustrare come utilizzarlo: il costruttore copia viene implicitamenteinvocato quando si assegna un oggetto della classe ad un altro ogget-to, passando questo secondo oggetto come argomento ed effettuandoesplicitamente la copia membro a membro. Come esempio consi-deriamo il caso piu interessante del costruttore copia per la classeFrattC in cui uno dei dati membro e un puntatore ad un array di-namico di numeri complessi: in questo caso il costruttore copia deveallocare dinamicamente la memoria e riferirla al puntatore membrodi classe nonche copiarvi i valori del vettore contenuto nell’ oggettoiniziale attraverso un ciclo sulle singole componenti come illustratonell’ esempio seguente:

FrattC::FrattC(const FrattC& t){

x = new Complesso[durata=t.durata];for(int i=0;i<durata;++i)

x[i] = t.x[i];}

– Con questa classe e possibile descrivere vettori bidimensionali, manon e difficile estenderla a rappresentare anche vettori tridimensio-nali aggiungento la terza componente ai dati e ai metodi membri di

64

classe. In questo caso l’ utente che utilizza oggetti Vettore non de-ve modificare il codice a causa dell’ estensione a tre dimensioni, inquanto l’ interfaccia della classe non viene modificata da questa ag-giunta e la rappresentazione del vettore (in due o tre dimensioni) eincapsulata entro la classe e risulta nascosta all’ utente. Questo e unodei vantaggi dell’ incapsulamento e della separazione tra interfaccia eimplementazione in C++ e per meglio apprezzarne l’ utilita vediamoun altro esempio: i due possibili modi di rappresentare un vettoretridimensionale, cioe la rappresentazione cartesiana e la rappresen-tazione polare. Nel nostro caso abbiamo scelto di usare la rappre-sentazione cartesiana e di rappresentare un vettore attraverso le suedue componenti x e y, ma potremmo anche scegliere di utilizzare larappresentazione polare caratterizzando il vettore con modulo e fase:si tratta di una scelta di implementazione che non deve comportaremodifiche nel codice dell’ utente che accede agli oggetti attraversoi metodi di interfaccia. Infatti supponendo di scegliere la rappre-sentazione polare per l’ implementazione del vettore e di avere comevariabili private il modulo m e la fase f del vettore, potremmo scriverei metodi getters in questo modo:

double Vettore::X(){

return m*cos(f);}

double Vettore::Y(){

return m*sin(f);}

e l’ utente continuerebbe ad usare la stessa interfaccia senza accor-gersi di cambiamenti nell’ implementazione.

– Tra i metodi nell’ implementazione si noti la funzione unit() che ha loscopo di normalizzare il vettore restituendo il versore corrispondentedi modulo unitario. Questa funzione utilizza al suo interno uno deicostruttori della classe e il metodo mod() precedentemente definitoinline: si tratta di un operazione consentita in C++ e non presentaproblemi di autoreferenza. Quindi in C++ all’ interno della definizio-ne di una classe si possono utilizzare oggetti appartenenti alla stessaclasse che si sta definendo!

65

3.2 Esercizio 3.2: overload di metodi e di operatori

Completare la classe Vettore aggiungendo gli operatori sovraccaricati per ef-fettuare operazioni algebriche vettoriali. Si consideri anche la possibilita dieffettuare il prodotto scalare.


file Vettore2.h:


#include <math.h>#include <iostream>

class Vettore {public:Vettore();Vettore(double a, double b);Vettore(const Vettore& v); // costruttore copiavoid setX(double a) { x = a; };void setY(double b) { y = b; };double X() { return x; };double Y() { return y; };double mod() { return sqrt(x*x+y*y); };Vettore unit();

Vettore operator+(const Vettore& v);Vettore& operator+=(const Vettore& v);Vettore operator-(const Vettore& v);Vettore& operator-=(const Vettore& v);Vettore operator-();double operator*(const Vettore &v);Vettore operator*(const double &a);Vettore operator/(const double &a);friend ostream& operator<<(ostream& o,

const Vettore& v);protected:double x;double y;

};#endif

66

file Vettore2.C:

#include "Vettore.h"

Vettore::Vettore() {x = 0;y = 0;

}

Vettore::Vettore(double a, double b) {x = a;y = b;

}

Vettore::Vettore(const Vettore& v) {x = v.x;y = v.y;

}

Vettore Vettore::unit() {Vettore s(x/mod(), y/mod());return s;

}

Vettore Vettore::operator+(const Vettore& v) {Vettore s(x+v.x, y+v.y);return s;

}

Vettore& Vettore::operator+=(const Vettore& v) {x += v.x;y += v.y;return *this;

}

Vettore Vettore::operator-(const Vettore& v) {Vettore s(x-v.x, y-v.y);return s;

}...

67

...

Vettore& Vettore::operator-=(const Vettore& v) {x -= v.x;y -= v.y;return *this;

}Vettore Vettore::operator-() {Vettore s(-x, -y);return s;

}

double Vettore::operator*(const Vettore &v) {return x*v.x + y*v.y;

}

Vettore Vettore::operator*(const double &c) {Vettore s(x*c, y*c);return s;

}Vettore Vettore::operator/(const double &c) {Vettore s(x/c, y/c);return s;

}

ostream& operator<<(ostream& o, const Vettore& v) {o << "(" << v.x << "," << v.y << ")";return o;

}


– In questo esercizio vengono sovraccaricati gli operatori algebrici edi output della classe vettore in modo analogo a quanto fatto inprecedenza per la classe Complesso.

– Una particolarita a cui si e gia fatto cenno e l’ utilizzo del costruttoredella classe Vettore all’ interno di uno dei metodi della classe stessa,come in operator+() ad esempio. Si noti inoltre che il metodo opera-tor*() e stato sovraccaricato due volte: una volta per gestire il casodel prodotto di uno scalare per un vettore e una volta per calcolareil prodotto scalare tra due vettori.

– Supponiamo ora di aver scritto una classe Vettore per rappresenta-re un vettore tridimensionale (come abbiamo visto si tratta di unamodifica semplice dell’ implementazione). Come potremmo fare per

68

includere anche il prodotto vettore tra le operazioni algebriche? Do-vremmo sovraccaricare il metodo operator*() passandogli un oggettodi tipo Vettore sia per calcolare il prodotto scalare sia per calcolareil prodotto vettore: l’ oggetto restituito nei due casi e diverso ed eun double o un Vettore, ma questa strada non e percorribile perchecon l’ overloading non si distinguono due funzioni dal tipo restituito,bensı solo e solamente dal numero e tipo di parametri. Quindi cio chesi puo fare e o utilizzare un flag passato alla funzione per sceglieretra i due tipi di prodotto:

operator*(const Vettore& v, int flag);

oppure definire in modo generale due nuove funzioni, non sovracca-ricando nessuno degli operator predefiniti:

double Vettore::ProdScal(const Vettore& v);Vettore Vettore::ProdVett(const Vettore& v);

La codifica di queste funzioni viene lasciata come esercizio.

69

3.3 Esercizio 3.3: classi CorpoCeleste e SistemaSolare (conelementi di design in UML)

Implementare le classi CorpoCeleste e Sistema Solare per la descrizione dipianeti in moto nel Sistema Solare secondo la legge di gravitazione di Newton.

• Design con UML del modello di sistema solare

– Come abbiamo visto negli esempi precedenti, nella scrittura di unprogramma in C++ e molto importante la fase di design e di proget-tazione delle classi e un notevole sforzo va sempre compiuto per ten-tare di anticipare tutti i possibili modi di utilizzo di un oggetto e perpresentare all’ utente un’ implementazione razionale ed un’ interfac-cia logica e semplice da utilizzare. Questa fase della programmazionee cosı importante che negli ultimi anni si e tentato di standardizzar-la creando un apposito linguaggio grafico chiamato UML (UnifiedModeling Language). In questo esercizio si vuole fare un breve cen-no a questa tecnica applicandola alla progettazione di due classi cherappresentino un corpo celeste del sistema solare e il sistema solarestesso rispettivamente.

– Nell’ approccio UML il design della classe segue due fasi: una primafase porta alla definizione delle relazioni tra l’ utente e gli oggettidella classe che si vuole progettare e tra gli oggetti delle varie classidel progetto e consiste nella scrittura di un diagramma in cui questerelazioni sono illustrate in modo grafico secondo una data simbologiae che prende il nome di Use Case. La seconda fase riguarda invecele singole classi e la definizione preliminare di tutti i dati membro edi tutti i metodi da definire per rappresentare gli oggetti in questionee le azioni che con essi l’ utente puo intraprendere e consiste in unarappresentazione grafica degli oggetti e delle relazioni che intercorro-no tra le classi corrispondenti: a questo schema viene dato il nomedi Modello Logico.

– Costruiamo anzitutto lo Use Case per il modello di Sistema Solare: siprocede disegnando una colonna per l’ utente rappresentato schema-ticamente da un omino e una colonna per ciascuna classe, nel nostrocaso CorpoCeleste e SistemaSolare: si tracciano poi delle linee trat-teggiate verticali in ogni colonna che rappresentano il fluire del tempoe servono per rappresentare le interazioni tra l’ utente e le classi infunzione del tempo. Queste interazioni vengono poi rappresentate dalinee orizzontali orientate da chi compie l’ azione verso chi la subiscee poste una dopo l’ altra in ordine cronologico dall’ alto in basso:queste frecce terminano sulla linea di evoluzione della classe oggettoin corrispondenza di un piccolo rettangolo che indica un’ azione che l’oggetto corrispondente deve compiere attraverso un metodo di classe.Nel caso del sistema solare, l’ utente deve istanziare alcuni oggetti ditipo CorpoCeleste e un oggetto di tipo SistemaSolare attraverso due

70

istruzioni new, nonche aggiungere i vari corpi al sistema solare at-traverso un metodo aggiungiCorpo (a questo livello di progettazionepossiamo pensare di memorizzare i pianeti in un array di oggetti ditipo CorpoCeleste e, per motivi che saranno chiari in seguito, deci-diamo di utilizzare un container di tipo vector di STL e di utilizzarloper memorizzare puntatori a oggetti CorpoCeleste). Pensiamo infinedi lasciar evolvere il sistema solare e che quest’ ultimo faccia muoveread ogni passo i singoli pianeti integrando l’ equazione di moto, notele masse e le condizioni di moto iniziali dei vari pianeti. Tutte questeidee e informazioni sono racchiuse in modo sintetico nello Use Casedi Fig. 10.

– Il modello logico schematizza a grandi linee la struttura delle due clas-si CorpoCeleste e SistemaSolare per quanto riguarda dati e funzionimembro: le classi sono rappresentate da riquadri che ne riportano ilnome, i dati membro e i metodi e sono unite da linee che ne indicanole dipendenze, come indicato in Fig. 11. Una freccia nel modellologico rappresenta una classe (ad esempio CorpoCeleste) che dipen-de da un’ altra classe (Vettore), mentre una linea terminata con unrombo indica una relazione di aggregazione, nel senso che un sistemasolare e composto da corpi celesti.

71

Figura 10: Use case UML per il modello di sistema solare.

72

Figura 11: Modello logico UML per il modello di sistema solare.

73


file CorpoCeleste.h:

#ifndef TAG_CORPOCELESTE#define TAG_CORPOCELESTE

#include "vettore.h"#include <vector>

const double G=6.672e-11;

class CorpoCeleste{public:CorpoCeleste();CorpoCeleste(double mass);CorpoCeleste(double mass, vettore x0);CorpoCeleste(double mass, vettore x0, vettore v0);void setMass(double mass) {M=mass;};void setX0(vettore x0) {x=x0;};void setV0(vettore v0) {v=v0;};vettore X() {return x;};vettore V() {return v;};double m() {return M;};vettore F(vector<CorpoCeleste*> lista);void evolve(double dT, vettore F);void evolve(double dT,

vector<CorpoCeleste*> lista);protected:vettore x;vettore v;double M;

};#endif

74

file CorpoCeleste.C:

CorpoCeleste::CorpoCeleste() {M = 0;vettore zero(0,0);x = zero;v = zero;

}

CorpoCeleste::CorpoCeleste(double mass) {M = mass;vettore zero(0,0);x = zero;v = zero;

}

CorpoCeleste::CorpoCeleste(double mass, vettore x0) {M = mass;x = x0;vettore zero(0,0);v = zero;

}

CorpoCeleste::CorpoCeleste(double mass, vettore x0,vettore v0) {

M = mass;x = x0;v = v0;

}

vettore CorpoCeleste::F(vector<CorpoCeleste*> list) {vettore R(0,0);for (int i=0; i<list.size(); i++) {vettore r = list[i]->X()-x;double d = r.mod();vettore dir = r.unit();if (d!=0) {R += dir*G*M*list[i]->m()/(d*d);

}}return R;

}...

75

...void CorpoCeleste::evolve(double dT, vettore F){vettore a=F/M;v+=a*dT;x+=v*dT;

}

void CorpoCeleste::evolve(double dT,vector<CorpoCeleste*> lista)

{evolve(dT,F(lista));

}


– La classe corpo celeste deve rappresentare un pianeta per quantoconcerne il moto nel sistema solare ed e quindi caratterizzata da duevettori costruiti con la nostra classe Vettore che ne rappresentanoposizione e velocita in un piano e da uno scalare che ne rappresenta lamassa. Queste grandezze servono a caratterizzare lo stato di un corpoceleste e costituiscono i dati private. Le funzioni public comprendonovari costruttori tra loro sovraccaricati, alcuni metodi di tipo gettere setter, una funzione per il calcolo della forza sul pianeta correntea causa della presenza di tutti gli altri corpi celesti e due metodiche si occupano dell’ evoluzione temporale e dell’ integrazione dell’equazione di moto secondo il semplice algoritmo iterativo:

ai = F/m

vi = vi−1 + ai ∗ dt

xi = xi−1 + vi ∗ dt

ove dt e lo step di evoluzione temporale che sulla scala da noi consi-derata corrispondera a 86400 s (cioe a un giorno terrestre).

– I costruttori permettono di inizializzare massa, posizione e velocita invari modi ponendo eventualmente a zero le quantita non specificatedirettamente nella chiamata attraverso la definizione di un vettorelocale nullo chiamato appunto zero. Il distruttore non e stato di-chiarata esplicitamente e verra aggiunto, nella sua forma di default,dal compilatore: la classe contiene comunque oggetti statici che sonorimossi dalla memoria appena l’ oggetto che li contiene esce di sco-pe. Infine getters e setters sono scritti in modo standard inline e nonnecessitano di spiegazioni.

76

– Vediamo ora il metodo F(), che ha il compito di calcolare la forzadi attrazione gravitazionale sul corpo celeste corrente a causa dellapresenza massiva degli altri corpi celesti presenti nella simulazione eche costituisce il cuore di tutta la simulazione. In questa funzioneutilizziamo per la prima volta oggetti della libreria STL quali con-tainers di tipo vector e iteratori: il loro utilizzo verra illustrato nelpunto successivo dopo aver illustrato l’ algoritmo. La funzione Friceve in ingresso un vector di puntatori a oggetti di classe Corpo-Celeste, cioe una lista di indirizzi di memoria a cui si trovano glioggetti che rappresentano i vari pianeti. La lista contiene puntatoria oggetti e non gli oggetti stessi per due motivi: anzitutto per motividi efficienza e rapidita di calcolo per non passare alla funzione troppidati ad ogni chiamata e poi anche per poter definire una classe Sondache ereditera da CorpoCeleste e per poter realizzare le proprieta dipolimorfismo come vedremo nell’ ultimo esercizio. Nella funzione, sidefinisce une vettore R, posto inizialmente a zero, per calcolare laforza risultante sul pianeta attraverso la formula della forza di New-ton tra il pianeta corrente e tutti gli altri corpi celesti: questa forza ecalcolata con un ciclo sulle varie componenti del vector list la cui di-mensione non e inizialmente nota alla funzione e viene calcolata conil metodo STL size(). Si vede qui un primo vantaggio dei containersSTL rispetto ai tradizionali arrays: in particolare con un vector nonsolo non e necessario conoscerne la dimensione, ma e anche possibi-le estenderne la lunghezza e inserire nuovi oggetti in una posizioneintermedia con i metodi push back() e insert() rispettivamente. Perogni pianeta del sistema solare diverso da quello corrente, si defini-sce il vettore distanza r e se ne valutano modulo e direzione: pereseguire queste operazioni si utilizzano i metodi e gli operatori dellaclasse Vettore definiti in precedenza e il risultato viene memorizzatonel vettore risultante R che viene restituito dalla funzione. Si osservicome il controllo di valutare la forza generata da un pianeta diversoda quello corrente venga effettuato attraverso la verifica di distanzanon nulla.

– Facciamo un breve cenno all’ utilizzo di oggetti container di STLcome vector. L’ oggetto viene definito attraverso un template incui il tipo di variabile da memorizzare nel vector denominato lista econsiderato un parametro, nel nostro caso un puntatore a CorpoCe-leste. L’ utilizzo di templates permette di utilizzare lo stesso codicegia scritto per un tipo anche per altri tipi: ad esempio potremmopensare di scrivere una classe Complesso con dati membri int, floate double e potremmo quindi definirla come template e utilizzare poii corrispondenti oggetti con chiamate del tipo:

Complesso<int> a;Complesso<float> b;Complesso<double> t;

77

Analogamente vector e pensato come container per tipi di qualun-que classe (di libreria o definita dall’ utente) ed ha quindi caratteredi estrema generalita e utilita. La Standard Template Library e co-stituita da containers, algoritmi ed iteratori e fornisce un numerodi metodi per effettuare operazioni su vector. Tra le funzioni piucomode e utilizzate si segnalano:

size()restituisce la lunghezza;

[]operator[ ] e usato per accedere alle singole componenti con lanotazione familiare degli arrays del C;

push back()aggiunge un elemento nell’ ultima posizione estendendo di un’unita la dimensione del vector;

clear()restituisce la lunghezza;

sort()riordina gli oggetti secondo un criterio che puo essere definitodall’ utente. Nel caso di un vector numerico viene effettuato ilriordinamento quicksort in ordine crescente;

insert()inserisce una componente in una posizione a scelta dell’ utente;

Per accedere alle componenti si utilizza un indice intero come nel casodegli arrays standard. Ma e anche possibile utilizzare un iteratoreche altro non e che un ulteriore oggetto C++ di STL che generalizzail semplice indice intero, come verra illustrato nel prossimo esercizio.

– L’ evoluzione temporale del moto procede attraverso il calcolo perpassi successivi della posizione di tutti i pianeti seconda l’ algoritmoillustrato in precedenza: questo e stato implementato in due meto-di sovraccaricati e denominati evolve() dei quali il secondo chiamail primo e che usano un passo temporale dT e una lista di pianeti.Questo approccio e stato seguito per evidenziare con maggior chia-rezza il calcolo della forza risultante con il metodo F() a cui vienepassato il vector di puntatori a CorpoCeleste che, come vedremo,definira lo stato dell’ oggetto SistemaSolare ad un certo istante. Laforza viene poi divisa per la massa e integrata due volte per calcolareil vettore posizione: ci si fermi un momento a riflettere quante ope-razioni queste semplici istruzioni comportano e quante istruzioni dicodice scritte finora sia nella classe Vettore che in CorpoCeleste e Si-stemaSolare vengono invocate per ottenere questo risultato. Il C++,attraverso le nozioni di incapsulamento dei dati e di sovraccarico deglioperatori, permette all’ utente di effettuare queste operazioni in mo-do apparentemente banale, ma adesso siamo in grado di apprezzaretutto il lavoro che sta dietro questo risultato!

78

file SistemaSolare.h:

#ifndef TAG_SISTEMASOLARE#define TAG_SISTEMASOLARE

#include "CorpoCeleste.h"#include <vector>

class SistemaSolare {public:SistemaSolare();int nCorpi() {return corpi.size();};void aggiungiCorpo(CorpoCeleste *corpo){corpi.push_back(corpo);

}void setMaxTime(double T) {TMax=T;};void setDeltaT(double T) {dT=T;};void evolvi(double dT);void evolvi(void (*f)(vector<CorpoCeleste *>));private:double TMax;double dT;vector<CorpoCeleste *> corpi;

};#endif

79

file SistemaSolare.C:

#include "SistemaSolare.h"#include <iomanip>

SistemaSolare::SistemaSolare(){corpi.clear();dT=0;TMax=0;

}

void SistemaSolare::evolvi(double dT){for(int i=0;i<corpi.size();i++)

corpi[i]->evolve(dT, corpi);}

void SistemaSolare::evolvi(void(*f)(vector<CorpoCeleste *>))

{double T=0;while(T<TMax){evolvi(dT);(*f)(corpi);T+=dT;

}}


– La classe SistemaSolare contiene il vector di puntatori a CorpoCelestedenominato corpi che caratterizza il sistema solare e due parametridouble che controllano l’ evoluzione temporale: il passo dT e la vitamassima del sistema TMax. Tra i metodi, oltre a costruttori e di-struttore di default e non esplicitamente codificati, si hanno alcunefunzioni ausiliarie definite inline che servono per definire i parame-tri temporali e per aggiungere elementi al vector dei pianeti o perconoscerne la lunghezza, cioe il numero di corpi presenti nel siste-ma (queste ultime due funzioni fanno uso dei metodi STL size() epush back di cui sopra). Il costruttore viene usato per inizializzarea zero il vector con il metodo clear() e per resettare il passo e ladurata della simulazione. Ancora, come nella classe CorpoCeleste,abbiamo due metodi per l’ evoluzione temporale denominati evolvi()e sovraccaricati e ancora uno dei due chiama l’ altro: ricordando il

80

modello logico che avevamo scritto progettando la classe in UML, sipuo capire come il metodo evolvi() del sistema solare invochi, uno aduno, i metodi evolve() dei vari pianeti passando ad ognuno di essi ilvector di puntatori agli altri corpi celesti in modo che la forza risul-tante possa essere valutata. Il secondo metodo evolvi(), piu esterno,contiene un ciclo while() che regola l’ evoluzione del sistema e checostituisce il centro motore di tutta la simulazione: questo metodoriceve in ingresso un puntatore ad una funzione generica f() che asua volta ha come parametro un vector di puntatori a CorpoCele-ste e restituisce un void. Questa ulteriore funzione viene utilizzatanel prossimo esercizio per generare l’ output sia numerico sia graficodella simulazione e viene chiamata attraverso il puntatore a scopodidattico per illustrare un nuovo modo di uso dei puntatori.

– Un puntatore a funzione e semplicemente un puntatore alla locazionedi memoria a cui si trova il codice eseguibile della funzione chiamatae puo essere dichiarato anche se, in senso stretto, una funzione non euna variabile. Per chiamare una funzione f() attraverso un puntatoresi deve anzitutto definire il puntatore con la sintassi:

void (*swap)(int a, int b)

in cui ad esempio viene definito un puntatore swap ad una funzioneche ha due argomenti di tipo int e che restituisce un void. Per chia-mare la funzione poi e sufficiente utilizzare la notazione standard peri puntatori con l’ operatore *: se swap e il puntatore alla funzione,allora *swap e la funzione e la chiamata sara ad esempio del tipo:

(*swap)(a, b)

Le parentesi sono necessarie affinche la chiamata alla funzione vengaeffettuata correttamente. Se scrivessimo:

void *swap(int a, int b) /* ERRORE! */

significherebbe che swap e una funzione che ritorna un puntatore avoid e non che e un puntatore ad una funzione che restituisce void:un errore di questo tipo potrebbe avere effetti disastrosi sul buonfunzionamento di un programma!

81

3.4 Esercizio 3.4: simulazione del Sistema Solare con out-put numerico e grafico

Scrivere un programma di simulazione del Sistema Solare utilizzando le classisviluppate negli esercizi precedenti con output numerico e grafico utilizzando lalibreria di ROOT.


file Pianeti.C:

#include "CorpoCeleste.h"#include "SistemaSolare.h"#include <iostream>#include <iomanip>

void output(vector<CorpoCeleste *> corpi) {vector<CorpoCeleste *>::const_iterator corpo =

corpi.begin();while (corpo != corpi.end()) {cout << setprecision(3) << setw(10)

<< (*corpo)->X().X()<< setprecision(3) << setw(10)<< (*corpo)->X().Y();

cout << endl;corpo++;

}}

int main() {SistemaSolare ss;CorpoCeleste sole(1.98e30, vettore(0., 0.),

vettore(0., 0.));CorpoCeleste mercurio(3.28e24, vettore(6.99e10, 0.),

vettore(0., 43479.17));CorpoCeleste venere(4.83e24, vettore(1.09e11, 0.),

vettore(0., 34840.23));CorpoCeleste terra(5.98e24, vettore(1.52e11, 0.),

vettore(0., 29476.35));CorpoCeleste marte(6.37e23, vettore(2.49e11, 0.),

vettore(0., 23025.48));CorpoCeleste giove(1.9e27, vettore(0., 8.16e11),

vettore(-12723.41, 0.));CorpoCeleste saturno(5.67e26, vettore(1.50e12, 0.),

vettore(0., 9370.69));...

82

...ss.aggiungiCorpo(&sole);ss.aggiungiCorpo(&mercurio);ss.aggiungiCorpo(&venere);ss.aggiungiCorpo(&terra);ss.aggiungiCorpo(&marte);ss.aggiungiCorpo(&giove);ss.aggiungiCorpo(&saturno);ss.setMaxTime(20*3.16e7);ss.setDeltaT(86400.);ss.evolvi(&output);return 0;

}

– Cominciamo con lo scrivere un programma main() in cui si definiscaun oggetto SistemaSolare composto da un certo numero di oggettiCorpoCeleste e con il farlo evolvere secondo i criteri illustrati soprae stampando le posizioni successive dei pianeti sullo schermo. Suc-cessivamente aggiungeremo a questo programma delle chiamate allalibreria grafica di ROOT con lo scopo di visualizzare il moto deipianeti durante l’ evoluzione.

– Scriviamo una funzione di output che accetta in ingresso un vectordi puntatori a oggetti CorpoCeleste e che utilizza un iteratore STLper eseguire un ciclo si tutti i pianeti stampandone su schermo lecoordinate di posizione x e y. E questa la funzione il cui punta-tore passeremo al metodo evolvi() della classe SistemaSolare, comeabbiamo visto nell’ esercizio precedente. Un iteratore costituisce l’estensione dell’ indice di un vettore e deve essere dichiarato con unaistruzione del tipo:

tipocontainer<tipooggetti>::const_iteratornomeiterator = nomecontainer.begin();

in cui nel nostro caso il container e un vector di puntatori a Corpo-Celeste e l’ iteratore, denominato corpo, viene posizionato all’ iniziodel vector dalla funzione begin(). All’ interno del ciclo di stampaottenuto con un while() che controlla il raggiungimento dell’ ultimoelemento con il metodo end() si utilizza l’ iteratore per accedere allesingole componenti: l’ iteratore si sostituisce all’ oggetto del vettorea cui esso fa riferimento come se fosse un puntatore ad esso con lanotazione *iteratore e poiche nel nostro caso l’ oggetto contenuto nelvector e un puntatore si utilizza l’ operatore -¿.Si osservi attentamente la notazione di un’ istruzione del tipo:

(*corpo)->X().X()

83

in cui si chiamano in cascata due metodi diversi, appartenenti a classidiverse ma chiamati allo stesso modo: il primo metodo X() appartienea CorpoCeleste e restituisce il vettore posizione, mentre il secondometodo X() appartiene a Vettore e restituisce la componente x delvettore. Si noti come l’ iteratore venga incrementato in maniera deltutto innocua usando l’ operatore postfisso ++ come su un normaleindice intero: in realta cio e possibile perche nella classe iterator ilmetodo operator++() e stato sovraccaricato, come ben sappiamo.

– Scritte le classi, il corpo del programma main() risulta facile da com-prendere in quanto le istruzioni si susseguono in ordine logico secondoquanto, potremmo dire, si farebbe con carta e penna per rappresenta-re il sistema solare: vengono chiamati il costruttore di SistemaSolaree dei vari pianeti inizializzando gli oggetti in modo opportuno. Ivari pianeti vengono poi aggiunti al sistema solare con il metodo ag-giungiCorpo (si ricordi l’ uso del puntatore). Il sistema viene quindiinizializzato e lasciato evolvere per 20 anni (un anno contiene circaπ 107 secondi) a passi di un giorno con il risultato di ottenere la stam-pa su schermo delle coordinate dei vari pianeti ad ogni passo: questivalori possono essere salvati in un file e analizzati con un programmadi tipo foglio elettronico oppure si puo visualizzare la simulazione inmodo grafico come faremo qui di seguito.

– Tabella di caratteristiche dei pianeti del Sistema Solare (dati reali):

Massa Perielio Afelio Eccentricita(Kg) (1e6 Km) (1e6 Km)

Mercurio 3.30e23 45.9 69.7 0.206Venere 4.87e24 107.4 109 0.007Terra 5.97e24 147.1 152.1 0.017Marte 6.42e23 206.7 249.1 0.093Giove 1.90e27 740.9 815.7 0.048Saturno 5.68e26 1347 1507 0.056Urano 8.68e25 2735 3004 0.047Nettuno 1.02e26 4456 4537 0.009Plutone 1.27e22 4425 7375 0.25

84

file PianetiAnimati.C:

#include "CorpoCeleste.h"#include "SistemaSolare.h"#include <iostream>#include <iomanip>#include <math.h>

#include <TROOT.h>#include <TView.h>#include <TCanvas.h>#include <TPolyMarker.h>

TCanvas *screen;TPolyMarker **pm;

void output(vector<CorpoCeleste *> corpi) {vector<CorpoCeleste *>::const_iterator corpo

= corpi.begin();int count = 0;while (corpo != corpi.end()) {double x = (*corpo)->X().X()*400/1e15;double y = (*corpo)->X().Y()*400/1e15;pm[count]->SetPoint(0, x, y);float size = 0.2*(log10((*corpo)->m())-20);if (size<0) {size = .5;

}pm[count++]->SetMarkerSize(size);corpo++;

}screen->Modified();screen->Update();

}

TROOT gRoot("PianetiAnimati", "Pianeti Animati");

int main() {cout << "Scegli fra i seguenti set inserendo un

numero" << endl;cout << "1) Sistema Solare normale" << endl;cout << "2) Sistema Solare con Giove molto piu’

massiccio del normale" << endl;...

85

...cout << "3) Sistema Solare con Giove pesante e

diretto verso il Sole " << endl;cout << "4) Sistema Solare con Giove e Saturno

pesanti " << endl;cout << "Input set :";int set = 0;cin >> set;

SistemaSolare ss;CorpoCeleste sole(1.98e30);CorpoCeleste mercurio(3.28e24);CorpoCeleste venere(4.83e24);CorpoCeleste terra(5.98e24);CorpoCeleste marte(6.37e23);CorpoCeleste giove(1.9e27);CorpoCeleste saturno(5.67e26);

sole.setX0(vettore(0., 0.));sole.setV0(vettore(0., 0.));mercurio.setX0(vettore(6.99e10, 0.));mercurio.setV0(vettore(0., 43479.17));venere.setX0(vettore(1.09e11, 0.));venere.setV0(vettore(0., 34840.23));terra.setX0(vettore(1.52e11, 0.));terra.setV0(vettore(0., 29476.35));marte.setX0(vettore(2.49e11, 0.));marte.setV0(vettore(0., 23025.48));giove.setX0(vettore(0., 8.16e11));giove.setV0(vettore(-12723.41, 0.));saturno.setX0(vettore(1.50e12, 0.));saturno.setV0(vettore(0., 9370.69));

if (set == 1) {} else if (set == 2) {cout << "Giove ha la massa vicina a quella del

Sole" << endl;giove.setMass(1.e30);

} else if (set == 3) {cout << "Giove ha la massa vicina a quella del

Sole e i due corpi si vanno incontro" << endl;giove.setMass(1.e30);giove.setV0(vettore(1000., -12723.41));sole.setV0(vettore(-1000., 12723.41));

...

86

...} else if (set == 4) {cout << "Giove e Saturno hanno una massa vicina a

quella del Sole" << endl;giove.setMass(1.e30);saturno.setMass(1.e30);

}

ss.aggiungiCorpo(&sole);ss.aggiungiCorpo(&mercurio);ss.aggiungiCorpo(&venere);ss.aggiungiCorpo(&terra);ss.aggiungiCorpo(&marte);ss.aggiungiCorpo(&giove);ss.aggiungiCorpo(&saturno);ss.setMaxTime(20*3.16e7);ss.setDeltaT(86400.);

screen = new TCanvas("c1", "Sistema Solare",0, 0, 800, 800);

screen->SetFillColor(1);TView *view = new TView(1);view->SetRange(-400,-400,-400,400,400,400);pm = new TPolyMarker*[ss.nCorpi()];for (int i=0; i<ss.nCorpi(); i++) {pm[i] = new TPolyMarker(1);pm[i]->SetMarkerStyle(20);pm[i]->SetMarkerColor(10-i);pm[i]->SetMarkerSize(20);pm[i]->Draw();

}screen->Update();ss.evolvi(&output);

cout << "Simulazione finita...";return 0;

}...


– Modifichiamo ora il programma aggiungendo un output grafico conla libreria di ROOT in modo da visualizzare su schermo il motodei pianeti: le modifiche sono molto simili a quelle introdotte nell’esempio dei frattali e consistono in aggiunte di oggetti grafici su cuivisualizzare le posizioni dei vari corpi celeste calcolate per mezzo del

87

nostro algoritmo. Gli oggetti grafici nuovi, oltre al gia noto TCan-vas, sono tre in questo esercizio: TROOT, TView e TPolyMarker esaranno illustrati nel seguito sia per la funzione output() che per lafunzione main(). Inoltre in questo programma viene implementatoun semplice menu (non grafico) per poter utilizzare il simulatore invarie configurazioni iniziali: con il sistema solare normale o in altresituazioni anomale in cui il programma e in grado di rispondere al-la domanda Che cosa succederebbe se?, ad esempio aumentando lamassa di un pianeta o dandogli una velocita iniziale perturbativa delmoto: la struttura di questo menu e elementare e si basa semplice-mente sull’ uso di cout e cin per I/O di valori dal programma percaricare una variabile di scelta set.

– La classe TROOT definisce un oggetto globale di ROOT denominatogROOT che serve per l’ utilizzo di alcuni menu grafici una volta che l’esecuzione del programma e terminata: definire un oggetto TROOTprima di un programma main() crea un’ area di memoria globale acui tutti gli oggetti grafici caricati successivamente fanno riferimentoe che rimane definita per tutta l’ esecuzione del programma. L’ uso diTROOT costituisce un argomento avanzato che esula dagli obiettividelle esercitazioni e per approfondimenti si puo consultare il manualedi ROOT disponibile online al sito http://root.cern.ch.

– Per definire un sistema di riferimento all’ interno di un TCanvasavevamo visto l’ oggetto TH2F che definisce un istogramma bidi-mensionale. Un altro modo di definire un sistema di assi cartesiani efornito dalla classe TView qui utilizzata: in questo caso il costruttore invocato con argomento 1 per indicare che gli assi da definire sonocartesiani (2 indicherebbe polari, 3 cilindrici e cosı via...) e succes-sivamente viene chiamato il metodo SetRange(x1,y1,z1,x2,y2,z2) chedelimita il punto di vista nel parallelepipedo compreso tra i punti(x1,y1,z1) e (x2,y2,z2).

– Il programma utilizza un TCanvas di ROOT come foglio su cui rap-presentare la posizione dei pianeti attraverso un nuovo oggetto graficodi tipo TPolyMarker: un TPolyMarker serve per rappresentare unalinea spezzata tra uno o piu punti in cui in corrispondenza di ognisingolo punto viene disegnato un simbolo come un cerchietto o altro.Nel programma si definisce un array di puntatori aTPolyMarker o, inaltri termini, un puntatore a puntatore di tipo TPolyMarker per me-morizzare gli indirizzi dei polymarkers corrispondenti ai vari pianeti:pm e un vettore di puntatori a TPolyMarker in cui ogni componen-te viene allocata dinamicamente come un TPolyMarker costituito daun singolo punto di cui si specificano le coordinate (x,y) e il tipo ecolore del simbolo grafico con gli opportuni metodi. Questa opera-zione viene fatta per ogni corpo celeste attraverso un ciclo for chetermina quando l’ ultimo pianeta del vector di SistemaSolare e stato

88

raggiunto. Ogni simbolo i pianeta viene quindi disegnato sul canvasche viene aggiornato con il metodo Update().

– La possibilita di definire puntatori a puntatori non deve stupire ospaventare, anche se puo essere fuorviante all’ inizio. Per capire beneche cosa si sta facendo e utile ricordare sempre la relazione esisten-te tra puntatori e arrays e cioe che il nome di un vettore e ilpuntatore al suo primo elemento: in questo caso e quindi equi-valente interpretare **pm come un puntatore a puntatore o come unarray di puntatori alla cui prima componente pm punta e la secondainterpretazione genera meno confusione ed e piu semplice da capi-re. Questo vettore di puntatori viene poi utilizzato nella funzioneoutput() per aggiornare le dimensioni del simbolo di ogni pianeta inmodo proporzionale alla massa (con un algoritmo del tutto arbitra-rio) e soprattutto per aggiornare la posizione del simbolo sul canvascon i nuovi valori di coordinate (opportunamente riscalati affinchesiano contenuti nel grafico). Il risultato grafico e riportato in Fig. 12durante l’ evoluzione dei pianeti.

89

Figura 12: Fotografia del nostro sistema solare durante l’ evoluzione temporaledel simulatore.

90

3.5 Esercizio 3.5: aggiunta di una Sonda

Aggiungere al programma precedente una sonda interplanetaria sfruttando il po-limorfismo del linguaggio C++ attraverso la definizione di una classe Sonda cheeredita da CorpoCeleste e che sia dotata di razzi e con un metodo evolvi() ditipo virtual che realizzi un late binding in esecuzione.


file Sonda.h:

#ifndef _SONDA#define _SONDA

#include "CorpoCeleste.h"

class Sonda : public CorpoCeleste {public:Sonda() : CorpoCeleste() { init(); };Sonda(double mass) : CorpoCeleste(mass) { init(); };Sonda(double mass, vettore x0) :

CorpoCeleste(mass, x0) { init(); };Sonda(double mass, vettore x0, vettore v0) :

CorpoCeleste(mass, x0, v0){ init(); };

void setOrigin(CorpoCeleste *origin) { o = origin; };void setStartTime(double start) { Tstart = start; };void accendiRazzi(double T, vettore v) {Taccensioni.push_back(T+Tstart);vRazzi.push_back(v); };

virtual void evolve(double dT, vettore F);protected:void init() { Tstart = -1; T = 0;

starting = 1; Nstep =0; };vector<double> Taccensioni;vector<vettore> vRazzi;CorpoCeleste *o;double Tstart;double T;char starting;int Nstep;

};

#endif

91

file Sonda.C:

#include "Sonda.h"

void Sonda::evolve(double dT, vettore F) {T += dT;if (Tstart < 0) return;if (T < Tstart) {x = o->X();

} else {if (starting) {x = o->X();v += o->V();starting = 0;

}if ((Taccensioni.size()>0)&&

(Taccensioni.size()>Nstep)&&(T >= Taccensioni[Nstep])) {

v += vRazzi[Nstep++];}v += F/M*dT;x += v*dT;

}}

file Pianeti.C:

#include "CorpoCeleste.h"#include "SistemaSolare.h"#include "Sonda.h"#include <iostream.h>#include <iomanip.h>

#include <TROOT.h>#include <TView.h>#include <TCanvas.h>#include <TPolyMarker.h>

#include <math.h>

TCanvas *screen;TPolyMarker **pm;

...

92

...void output(vector<CorpoCeleste *> corpi) {vector<CorpoCeleste *>::const_iterator corpo

= corpi.begin();int count = 0;while (corpo != corpi.end()) {double x = (*corpo)->X().X()*400/1e15;double y = (*corpo)->X().Y()*400/1e15;pm[count]->SetPoint(0, x, y);float size = 0.2*(log10((*corpo)->m())-20);if (size<0) size = .5;pm[count++]->SetMarkerSize(size);corpo++;

}screen->Modified();screen->Update();

}

TROOT gRoot("PianetiAnimati", "Pianeti Animati");

int main(int argc, char **argv) {SistemaSolare ss;cout << "Scegli fra i seguenti set inserendo un

numero" << endl;cout << "1) Sistema Solare con Voyager che orbita

tra i pianeti interni" << endl;cout << "2) Sistema Solare con Giove molto piu’

massiccio del normale" << endl;cout << "3) Sistema Solare con Giove pesante e

diretto verso il Sole " << endl;cout << "4) Sistema Solare con Giove e Saturno

pesanti " << endl;cout << "Gli altri numeri simulano un normale

sistema solare" << endl;cout << "Input set :";int set = 0;cin >> set;CorpoCeleste sole(1.98e30);CorpoCeleste mercurio(3.28e24);CorpoCeleste venere(4.83e24);CorpoCeleste terra(5.98e24);CorpoCeleste marte(6.37e23);CorpoCeleste giove(1.9e27);CorpoCeleste saturno(5.67e26);

...

93

...Sonda voyager(110000.);

sole.setX0(vettore(0., 0.));sole.setV0(vettore(0., 0.));mercurio.setX0(vettore(6.99e10, 0.));mercurio.setV0(vettore(0., 43479.17));venere.setX0(vettore(1.09e11, 0.));venere.setV0(vettore(0., 34840.23));terra.setX0(vettore(1.52e11, 0.));terra.setV0(vettore(0., 29476.35));marte.setX0(vettore(2.49e11, 0.));marte.setV0(vettore(0., 23025.48));giove.setX0(vettore(0., 8.16e11));giove.setV0(vettore(-12723.41, 0.));saturno.setX0(vettore(1.50e12, 0.));saturno.setV0(vettore(0., 9370.69));

voyager.setV0(vettore(24500., 64000.));if (set == 1) {cout << "Il voyager ha una velocita’ minore e

rimane tra i pianeti interni" << endl;voyager.setV0(vettore(-24500., 20000.));


Sole" << endl;giove.setMass(1.e30);


Sole e i due corpi si vanno incontro" << endl;giove.setMass(1.e30);giove.setV0(vettore(1000., -12723.41));sole.setV0(vettore(-1000., 12723.41));

} else if (set == 4) {cout << "Giove e Saturno hanno una massa vicina a

quella del Sole" << endl;giove.setMass(1.e30);saturno.setMass(1.e30);

}voyager.setOrigin(&terra);voyager.setStartTime(180*86400.);voyager.accendiRazzi(400*86400,vettore(-5000., 0.));

...

94

...

ss.aggiungiCorpo(&sole);ss.aggiungiCorpo(&mercurio);ss.aggiungiCorpo(&venere);ss.aggiungiCorpo(&terra);ss.aggiungiCorpo(&marte);ss.aggiungiCorpo(&giove);ss.aggiungiCorpo(&saturno);ss.aggiungiCorpo(&voyager);ss.setMaxTime(20*3.16e7);ss.setDeltaT(86400.);

screen = new TCanvas("c1", "Sistema Solare",0, 0, 800, 800);

screen->SetFillColor(1);TView *view = new TView(1);view->SetRange(-400,-400,-400,400,400,400);pm = new TPolyMarker*[ss.nCorpi()];for (int i=0; i<ss.nCorpi(); i++) {pm[i] = new TPolyMarker(1);pm[i]->SetMarkerStyle(20);pm[i]->SetMarkerColor(10-i);pm[i]->SetMarkerSize(20);pm[i]->Draw();

}

screen->Update();ss.evolvi(&output);cout << "Simulazione finita...";cin >> set;return 0;

}


– Affrontiamo in quest’ ultimo esercizio il tema dell’ ereditarieta e delpolimorfismo tra classi in C++ estendendo il simulatore del SistemaSolare con l’ aggiunta di un nuovo tipo di oggetto per rappresentareil moto di una sonda interplanetaria. Una sonda e un corpo celestesoggetto alla forza gravitazionale di Newton, ma dotato di razzi chepossono fornire un’ ulteriore spinta per un determinato intervallo ditempo modificandone quindi la legge di moto. Questi fatti possonoessere tradotti in C++ e implementati nel nostro simulatore attra-verso la scittura di una nuova classe, che chiameremo Sonda, che

95

eredita tutte le proprieta della classe CorpoCeleste, ma con in piudelle variabili per gestire l’ accensione e lo spegnimento dei razzi edun metodo evolvi() che tenga conto delle modifiche all’ equazione dimoto.

– L’ ereditarieta tra classi in C++ viene esplicitata con la definizionedella nuova classe Sonda in cui si specifica la relazione con la classeCorpoCeleste secondo la notazione:

\item{} class classe_figlia: public classe_madre {...}

Da questo momento in poi tutte le caratteristiche di Corpoceleste(dati e funzioni membro di classe) appartengono anche a Sonda. Laparola chiave public nella dichiarazione di ereditarieta serve garantirealla classe Sonda la possibilita di accedere ai dati protected dichiaratinella classe CorpoCeleste: si ricorda che in una classe dati e funzionipublic sono accessibili dall’ esterno, dati e funzioni private sono ac-cessibili sono entro metodi della stessa classe, mentre dati e funzioniprotected sono equivalenti a private ma garantiscono l’ accesso dafunzioni di una classe ereditaria.

– Oltre ai dati membro di CorpoCeleste, la classe Sonda contiene unvector STL di double per memorizzare i tempi di accensione dei razzie un vector STL di vettori per memorizzare i valori di spinta dei razziche modificheranno la velocita del corpo celeste. Inoltre aggiungiamoun puntatore a CorpoCeleste per memorizzare l’ indirizzo di memoriadel pianeta da cui la Sonda inizialmente parte e alcuni valori per ilcontrollo dell’ evoluzione del moto: Tstart per l’ istante di lancio, Tper il tempo della sonda, un flag starting per sapere se i razzi sonoaccesi o spenti e infine Nstep per contare il numero di accensionidei razzi. Aggiungiamo anche un metodo init() nella parte protecteddella classe per l’ inizializzazione di queste variabili di stato.

– Il costruttore della classe Sonda deve istanziare comunque un oggettodi tipo CorpoCeleste e tutti i costruttori di Sonda tra loro sovrac-caricati chiamano i corrispondenti costruttori di CorpoCeleste vistiprecedentemente ed eseguono anche la funzione init(), comune a tut-ti e dichiarata come protected(). Altri metodi della classe Sondasono di tipo setters per definire il pianeta di partenza e i parame-tri per l’ accensione dei razzi, con istruzioni che dovrebbero esserecomprensibili. Si noti l’ uso di push back per caricare dati nei duevectors dei tempi di accensione e delle velocita addizionali dovute airazzi. Molto importante e discutere il metodo evolvi() della classeSonda che, come si vede, e stato dichiarato virtual e che realizza inpratica il polimorfismo del C++: esso verra illustrato prima veden-done l’ implementazione e spiegandone le modalita di chiamata e diesecuzione.

96

– L’ unica funzione non definita inline di Sonda e il metodo evolvi() incui viene integrata l’ equazione di moto della sonda tenendo conto diun’ eventuale spinta dovuta ai razzi che possono essere accesi e spentia piacimento. In questa funzione vengono effettuati diversi controlliattraverso una serie di istruzioni if: se Tstart e negativo, cioe se irazzi non vengono accesi la sonda non parte e il metodo esce senzafare nulle, mentre per T inferiore a Tstart (positivo) la sonda e fermasul pianeta da cui viene lanciata e si muove in maniera solidale adesso: si capisce quindi l’ utilita di un puntatore a CorpoCeleste comemembro di classe di Sonda. Si effettua poi un controllo sul fatto cheil lancio sia stato effettuato oppure no attraverso il flag starting: illancio viene effettuato fornendo alla sonda una velocita iniziale pari aquella del pianeta iniziale e ponendo successivamente il flag a zero pernon rientrare piu in questo caso dell’ if (il lancio viene effettuato unasola volta). Infine si controlla l’ accensione dei razzi con un’ istruzio-ne if che effettua un triplice controllo attraverso due AND booleane:se il vector dei tempi di accensione contiene almeno un elemento, sela sua dimensione e superiore al numero di volte che si sono accesi irazzi e se il tempo T e maggiore dei tempi di accensione specificatiin tale vector, allora la velocita della sonda viene incrementata conla spinta fornita dai razzi e caricata nel corrispondente vector e ilvalore del numero di accensioni e incrementato di un’ unita. L’ equa-zione di moto viene infine integrata con la stessa tecnica usata perCorpoCeleste.

– Analizziamo ora le proprieta di polimorfismo in relazione all’ eredita-rieta tra classi. Nel nostro caso abbiamo oggetti di tipo CorpoCele-ste e oggetti di tipo Sonda e per entrambi risulta definito un metodoevolvi(): cerchiamo di capire che cosa succede eseguendo il seguenteframmento di codice:

Sonda *Voyager = new Sonda();CorpoCeleste *Terra = new CorpoCeleste();Voyager.evolvi();Terra.evolvi();

Nelle ultime due istruzioni vengono chiamati i metodi evolvi() di Son-da e di CorpoCeleste rispettivamente: ma a causa dell’ ereditarietauna Sonda e a tutti gli effetti un CorpoCeleste e, se non avessimo de-finito un metodo evolvi() in Sonda, sarebbe stato chiamato il metodoevolvi() di CorpoCeleste. Tutto avrebbe funzionato ancora, ma senzala gestione dei razzi e quindi con un risultato sbagliato: quindi unaclasse Sonda che eredita da CorpoCeleste deve avere esplicitamenteun metodo evolvi() che verra chiamato su oggetti di tipo Sonda()(questa proprieta prende il nome di early binding).

97

Ma consideriamo ora il seguente frammento di codice:

Sonda *Voyager = new Sonda();CorpoCeleste *pVoyager = Voyager;pVoyager.evolvi();

Che cosa succede adesso? pVoyager e un puntatore a CorpoCeleste e,sfruttando l’ ereditarieta, anche Voyager in quanto CorpoCeleste e unpuntatore ad un CorpoCeleste che e anche una Sonda: ma il metodoevolvi() chiamato dall’ ultima istruzione e quello di CorpoCeleste enon quello di Sonda perche pVoyager non ha modi di sapere di esse-re stato incluso con l’ ereditarieta dentro Voyager! Affinche questachiamata esegua il metodo evolvi() di Sonda e sufficiente dichiarareil metodo evolvi() di Sonda virtual. Un metodo virtuale controllase l’ oggetto a cui appartiene ha una relazione di ereditarieta conaltri oggetti e chiama la funzione della classe giusta: in questo ca-so il compilatore si accorge che pVoyager e incluso (embedded) inVoyager e chiama quindi il metodo evolvi() di Sonda come ci aspet-teremmo. Questa capacita di invocare a run-time la corretta funzionetra piu classi legate da una relazione di ereditarieta prende il nome dilate binding (contrapposto all’ early binding che si avrebbe senzadichiarare il metodo virtual).

– Il late binding e il polimorfismo funzionano con l’ utilizzo di punta-tori, come negli esempi precedenti, e questo spiega il motivo per cuiin questi esercizi si e sempre lavorato con un vector di puntatori aCorpoCeleste e non con un vector di oggetti CorpoCeleste sempli-cemente: in quest’ ultimo caso, oltre ad una perdita di prestazionidel programma, l’ ereditarieta con polimorfismo tra CorpoCeleste eSonda non avrebbe funzionato.

– Un’ ultima nota sul programma main: anche la Sonda Voyager vieneaggiunta al Sistema Solare con il metodo aggiungiCorpo() in quantoCorpoCeleste e deve essere inizializzata opportunamente usando imetodi di Sonda() visti sopra prima di far evolvere il simulatore. Ilrisultato finale, con il Voyager in orbita, e riportato in Fig. 13. Buonaesplorazione! ;-)

98

Figura 13: Il Voyager e decollato!

99

4 Bibliografia

1. Kerningham & Ritchie, Linguaggio C, Ed. Jackson 1989

2. G. Organtini, Elementi di programmazione object oriented, Dispensa

3. A. Rimoldi, Linguaggio C++, Dispensa

100

Esercitazioni di C++ per Fisici Annotazioni dal Corsofontana/+.pdf · 2011. 7. 14. · dati, il...

Documents

Transcript of Esercitazioni di C++ per Fisici Annotazioni dal Corsofontana/+.pdf · 2011. 7. 14. · dati, il...