<- Fedora Core 3: un'installazione difficile - Indice Generale - Copertina - OpenVPN, VPN in pochi passi ->

Sistemi Liberi

Il g77 e un esempio di applicazione in fisica classica come strumento di studio scolastico nelle scienze

di Gabriele Zucchetta

Le scienze computazionali

L'avvento dei computer ha messo a disposizione dei ricercatori un nuovo strumento di studio e di analisi, la cui potenzialità è cresciuta parallelamente agli sviluppi tecnologici che il campo informatico ha avuto.

La Matematica, la Chimica, la Fisica, la Biologia sono le aree più importanti dove l'approccio computazionale ha aperto nuove strade nella ricerca. Con il termine scienze computazionali possiamo indicare quei campi di ricerca che fanno uso del computer. Oltre al campo della ricerca pura, l'approccio computazionale ha permesso di sviluppare strumenti di calcolo per la risoluzione di problemi pratici come, ad esempio, le progettazioni di strutture o le previsioni meteorologiche.

Nell'ambito scolastico, l'approccio computazionale, in affiancamento allo studio teorico e a quello sperimentale di laboratorio, è uno strumento con il quale lo studente può misurare le proprie capacità, approfondire le proprie conoscenze ed imparare ad analizzare problemi nuovi trovando soluzioni.

Infatti, l'approccio computazionale ad un problema non richiede solo conoscenze in merito al problema stesso, ma necessita di un approccio interdisciplinare ed è questo a rendere questo tipo di studio molto stimolate per lo studente.

Gli strumenti richiesti sono, a differenza delle esperienze in laboratorio, alla portata di chiunque: è sufficiente un computer, oramai in ogni casa, e di un compilatore. Tra i compilatori disponibili, il g77 risponde bene alle esigenze degli studenti, come verrà mostrato, appunto, nel corso dell'articolo.

Il g77

La programmazione dei primi computer avveniva attraverso istruzioni dette codice macchina, un sistema estremamente complesso e solo alla portata di persone molto specializzate. Pertanto, sono stati creati dei linguaggi alternativi, sufficientemente semplici da utilizzare, ai quali è stato dato il nome di "linguaggi di programmazione". Il compilatore si occupa di tradurre in codice macchina le istruzioni scritte in questi linguaggi: un file che contiene la serie di istruzioni per un programma è detto codice sorgente. Il g77 o GNU-FORTRAN [1] è un compilatore FORTRAN (versione 77) sviluppato dalla Free Software Foundation [2] sotto licenza GPL [3].

Ma la storia del FORTRAN affonda le sue origini ben più lontano nel tempo. Infatti, negli anni '50 è nata, in casa IBM, la prima versione del FORTRAN, acronimo di "FORmula TRANslator", cioè "traduttore di formule". Il fatto che la sintassi del linguaggio sia stata sviluppata sulla base del calcolo scientifico ha permesso al FORTRAN di diventare il linguaggio di programmazione delle scienze computazionali.

Nel corso degli anni si sono susseguite varie versioni dello standard FORTRAN (II, III, IV, FORTRAN 66, ecc.) fino alla versione 77, al cui standard risponde il g77. Nel 1990 è stato definito un nuovo standard, al quale sono state apportati miglioramenti nel 1990 e nel 1995. Questi due nuovi standard sono stati denominati, rispettivamente, FORTRAN 90 e FORTRAN 95. Recentissimamente, nel 2003, è stato definito un nuovo standard denominato FORTRAN 2003.

L'uso del compilatore è decisamente semplice. Supponendo di avere il compilatore g77 correttamente istallato nel proprio sistema (si faccia riferimento ai manuali della propria distribuzione per l'installazione del compilatore g77 che, generalmente, fa parte dei pacchetti di sviluppo), i comandi per compilare un programma e per eseguirlo sono, rispettivamente, i seguenti:

g77 programma.f -o programma
./programma

Dove programma.f è il file che contiene il codice sorgente del programma. Per scrivere un file sorgente si può usare un qualsiasi editor di testi, tra cui vim [4] è uno dei più potenti.

Negli ultimi anni, altri linguaggi di programmazione stanno velocemente sostituendo il FORTRAN: in particolare il C [7] e, più recentemente, Java [8]. Questi hanno il vantaggio di essere più versatili, nel senso che si prestano bene per la programmazione di qualsiasi problema, non solo quelli scientifici. Pertanto, la loro diffusione in ogni campo convince vecchi e nuovi programmatori scientifici a preferirli al FORTRAN.

Tuttavia, uno strumento appositamente creato per problemi matematici è più adatto allo studio, risultando più veloce sia da insegnare che da apprendere. In particolare, anche per la sua non eccessiva complessità, il FORTRAN 77 è particolarmente adatto a questo scopo.

Per completezza, è bene dire che lo sviluppo del compilatore g77 al fine di adeguarlo allo standard 95 ha generato due progetti: il g95 [5] ed il gfortran [6]. Per chi lavora a livello professionale, lo standard 95 è sicuramente più interessante, pertanto il g95 o il gfortran sono strumenti più potenti rispetto all'ormai datato g77.

Applicazione

Di seguito, viene proposto un programma su un problema scientifico per far "toccare con mano" il significato di approccio computazionale. Il problema che viene affrontato riguarda la legge della Gravitazione Universale ed è la derivazione e l'analisi del moto dei pianeti. Si farà ricorso alla Meccanica, che è un ramo della Fisica che studia il comportamento di sistemi sottoposti all'azione di forze.

La legge della gravitazione universale

Sin dall'antichità, si riteneva che la Terra fosse il centro dell'Universo e che attorno ad essa ruotassero gli altri pianeti. Questa teoria, detta geocentrica, fu sviluppata da Tolomeo nel II secolo dopo Cristo.

In netto contrasto con questa teoria, l'astronomo Copernico sviluppò, nel XVI secolo, la teoria eliocentrica che vede i pianeti, compresa la Terra, muoversi attorno al Sole. La teoria eliocentrica ebbe il suo definitivo sopravvento su quella geocentrica quando Keplero, all'inizio del XVII secolo, sulla base di dati scientifici, formulò le sue famose tre leggi:

1. i pianeti descrivono intorno al Sole orbite ellittiche di cui il Sole occupa uno dei due fuochi;
2. le aree descritte dal raggio vettore tracciato dal sole ai pianeti sono proporzionali ai tempi impiegati a descriverle;
3. i quadrati dei tempi impiegati dai pianeti a descrivere le proprie orbite sono proporzionali ai cubi dei semiassi maggiori delle ellissi.

Più tardi, alla vine del XVII secolo, il fisico Isaac Newton generalizzò il lavoro di Keplero arrivando alla Legge di Gravitazione Universale, enunciandola nel 1684:

Qualsiasi corpo dell'Universo attrae ogni altro corpo con una forza diretta lungo la linea che congiunge i centri dei due corpi, di intensità direttamente proporzionale al prodotto delle loro masse ed inversamente proporzionale al quadrato della loro distanza.

Matematicamente, si ha:

(1) F = G · ((m1 · m2) / r2)

Dove:

• F = forza gravitazionale tra i due corpi;
• m₁ = massa del primo corpo;
• m₂ = massa del secondo corpo;
• r = distanza tra i due corpi;
• G = costante di gravitazione universale.

Occorre tener presente che la 1 è il modulo di una forza centrale attrattiva.

L'approccio numerico

Conoscendo l'espressione della forza che caratterizza il sistema (equazione 1) e facendo ricorso alle leggi della Meccanica Classica (leggi sviluppate dallo stesso Newton), possiamo derivare un algoritmo che ci permette di calcolare il moto dei pianeti o, detto ancora più esplicitamente, di calcolarne le posizioni istante per istante. Infine, note le posizioni istantanee, possiamo elaborare i dati, calcolare quantità fisiche (distanza massima, distanza media, tempo di rivoluzione) fino ad estrapolare leggi fisiche.

Per fare ciò, oltre a conoscere le forze, abbiamo bisogno di uno schema per passare dalle condizioni ad un istante a quelle di un istante successivo e pensare di ripetere questo schema varie volte fino a studiare il sistema per una arbitraria quantità di tempo. Detto in termini matematici, uno schema siffatto permette di integrare le equazioni del moto.

Per semplicità, altrimenti dovremmo considerare anche l'interazione reciproca tra i vari pianeti orbitanti, ci poniamo nella condizione in cui vi sia un pianeta isolato che orbita intorno ad un altro pianeta. Tenendo presente che l'orbita avviene su un piano, possiamo fissare un opportuno sistema di riferimento cartesiano centrato sul pianeta attorno al quale l'altro orbita e scomporre la forza nelle sue componenti lungo gli assi. Facendo riferimento alla figura 1, si possono ricavare, sfruttando la similitudine tra i triangoli MOP e MON rispetto ai triangoli QOR e QOS:

(2) ax = -a · (x / r)
(3) ay = -a · (y / r)

A questo punto, uguagliando l'equazione 1 con quella della seconda legge di Newton della meccanica (F = m · a), si può ricavare l'accelerazione a cui il pianeta orbitante è sottoposto:

(4) a = G · (M / r2)

dove con M indichiamo la massa del pianeta attorno al quale avviene l'orbita. Si noti che la massa del pianta orbitante non entra nell'espressione dell'accelerazione. Sostituendo la 4 nella 2 e 3 si ottengono, rispettivamente, le equazioni 5 e 6:

(5) ax = -G · (M / r2) · (x / r) = -G · ((M · x) / r3)

Figura 1: scomposizione dei vettori accelerazione in componenti cartesiane.

(6) ay = -G · (M / r2) · (y / r) = -G · ((M · y) / r3)

Benché quello che stiamo studiando non sia, evidentemente, un moto uniformemente accelerato, possiamo pensare di scegliere un intervallo di tempo tra un istante e l'altro in modo tale che, applicando al moto le leggi del moto uniformemente accelerato, l'errore che si commette sia trascurabile. Sotto questa ipotesi, prendiamo in considerazione le leggi matematiche che descrivono il moto uniformemente accelerato:

(7) st = s0 + v0t + (1/2) · a · t2
(8) vt = v0 + at

dove:

• s_t = posizione al tempo t;
• v_t = velocità al tempo t;
• s₀ = posizione iniziale (al tempo 0);
• v₀ = velocità iniziale (al tempo 0);
• a = accelerazione (costante nel tempo).

Riscriviamo queste due formule applicandole al nostro caso, nell'intervallo di tempo tra t e t + Δt ( dove Δt è la differenza di tempo tra un istante e quello successivo):

(9) x(t + Δt) = x(t) + vx(t)Δt + (1/2) · ax(t)(Δt)2
(10) y(t + Δt) = y(t) + vy(t)Δt + (1/2) · ay(t)(Δt)2
(11) vx(t + Δt) = vx(t) + ax(t)(Δt)2
(12) vy(t + Δt) = vy(t) + ay(t)(Δt)2

Questo algoritmo, conosciuto come algoritmo di Eulero, può essere facilmente migliorato sostituendo l'accelerazione iniziale con quella media tra l'istante t e t + Δt. Sulla base di ciò, le equazioni precedenti diventano:

(13) x(t + Δt) = x(t) + vx(t)Δt + (1/2) · ax(t)(Δt)2
(14) y(t + Δt) = y(t) + vy(t)Δt + (1/2) · ay(t)(Δt)2
(15) vx(t + Δt) = vx(t) + (1/2) · (ax(t) + ax(t + Δt))(Δt)2
(16) vy(t + Δt) = vy(t) + (1/2) · (ay(t) + ay(t + Δt))(Δt)2

Questa versione è chiamata algoritmo di Verlet.
Preme far presente che l'algoritmo di Verlet è valido anche per leggi diverse da quella gravitazionale ed è la particolare legge di forza che determina l'espressione dell'accelerazione da usare (nel nostro caso, equazioni 5 e 6). Il sorgente del programma riportato di seguito implementa quanto derivato in questa sezione; il sistema preso in considerazione è quello Sole-Terra. I dati necessari allo studio di questo sistema sono riportati nella tabella 1 e sono utilizzati all'interno del programma usato per questo studio.

Il sorgente

 1       program grav
 2 C
 3 C Dichiarazioni variabili:
 4 C
 5       integer i
 6       real deltat,massa,g,posx,posy,velx,vely,accx,accy,r
 7 C
 8 C Definizione parametri:
 9 C  - massa del sole (massa) = 1.989e30 (Kg)     
10 C  - cost. gravitazione universale (g) = 
11 C          = 6.67259e-11 (m**3 Kg**(-1) s**(-2)) 
12 C
13       parameter(massa=1.989e30,g=6.67259e-11)
14       deltat= 10000
15       posx = 152.1e+9
16       posy = 0.0
17       velx = 0.0
18       vely = 29.29e+3
19 C
20 C Calcolo della distanza e delle accelerazioni 
21 C della posizione di partemza
22 C
23       r = sqrt(posx**2.0+posy**2.0)
24 C
25       accx = 0.0 - (massa*g*posx)/(r**3.0)
26       accy = 0.0 - (massa*g*posy)/(r**3.0)
27 C
28 C Inizio del ciclo del calcolo dalla traiettoria del pianeta
29 C 
30       do i=1,10000
31 C
32 C Calcolo delle accelerazioni e delle velocità al tempo t
33 C 
34       posx = posx+(deltat*velx)+(0.5*deltat**2.0*accx)
35       posy = posy+(deltat*vely)+(0.5*deltat**2.0*accy)
36       velx = velx+(0.5*deltat*accx)
37       vely = vely+(0.5*deltat*accy)
38 C
39 C Calcolo della distanza, delle accelerazioni e delle 
40 C velocità al tempo t+deltat
41 C
42       r = sqrt(posx**2.0+posy**2.0)
43 C
44       accx = 0.0 - (massa*g*posx)/(r**3.0)
45       accy = 0.0 - (massa*g*posy)/(r**3.0)
46 C
47       velx = velx+(0.5*deltat*accx)
48       vely = vely+(0.5*deltat*accy)
49 C
50 C Stampa delle posizioni del pianeta orbitante
51 C
52       write(*,*) posx,posy 
53       enddo
54       end

I numeri su ogni riga non fanno parte del sorgente, ma sono stati aggiunti per agevolarne il commento.

Compilazione, esecuzione e analisi dei dati

La compilazione avviene, come già detto, nel seguente modo:

g77 grav.f -o grav

Così facendo, si ottiene un file (grav) eseguibile, che può essere lanciato semplicemente così:

./grav

In questo modo, l'output del programma (due serie di lunghe colonne di numeri) viene indirizzato sullo schermo, che è lo standard output. Per poter successivamente elaborare questi dati, è bene reindirizzare l'output su un file (dati.dat) eseguendo il programma nel seguente modo:

./grav > dati.dat

Analisi del codice sorgente

Al fine di entrare nel dettaglio dell'implementazione dell'algoritmo, vedremo alcuni aspetti del sorgente, tenendo presente che ciò è solo un'introduzione al linguaggio, la cui conoscenza richiede ben più nozioni. Si faccia riferimento alla sezione "Conclusioni" per informazioni su guide e manuali.

In primo luogo, il codice può essere scritto a partire dalla settima colonna fino alla settantaduesima, eredità dovuta al supporto (le schede perforate) che era in uso negli anni in cui il linguaggio è nato. Per tale motivo, ci sono sei spazi vuoti prima della prima lettera di ogni riga di codice.

Le righe di commento possono essere scritte inserendo in prima colonna la lettera "C". È bene non lesinare nell'uso dei commenti, in quanto permettono agli altri, e anche a noi stessi, di comprendere ciò che viene fatto dal programma, rendendo meno complesse le fasi di correzione degli errori e di miglioramento e sviluppo del codice.

Il codice inizia e termina, rispettivamente, con la dichiarazione del nome (riga 1) e della fine del programma (riga 54). Queste dichiarazioni sono sempre presenti nel programma e rappresentano il programma più semplice che sia possibile scrivere. Infatti, un sorgente composto da queste due sole istruzioni può essere compilato ed eseguito.

Nelle righe 5 e 6 vengono dichiarate le variabili. Questa operazione permette di far assegnare al nostro programma la giusta quantità di memoria. Le variabili contengono numeri che vengono calcolati durante l'esecuzione del programma e i loro valori possono variare man mano che il programma viene eseguito. Nel nostro caso, la riga 5 specifica una variabile di tipo intero, mentre la riga 6 specifica dieci variabili di tipo reale. La descrizione del significato di ogni variabile è riportata nella tabella 2.

Le unità di misura riportate nella terza colonna della tabella non sono esplicite nel programma ma occorre tener presente questo aspetto ed usare numeri coerenti rispetto ad un determinato sistema di misura. Si possono simulare altri sistemi modificando le seguenti variabili o parametri: massa (riga 13), posx (riga 15), posy (riga 16), velx (riga 17), vely (riga 18). Si fa presente che, dopo la modifica del sorgente, il programma va ricompilato, in modo da creare un nuovo eseguibile con le modifiche apportate.

Nella riga 13, le variabili massa e g vengono definite come parametri: viene loro assegnato un valore che non può essere cambiato durante l'esecuzione del programma. L'uso dei parametri è frequente e permette di rendere semplici e indolori alcune modifiche al programma, sia da parte dell'autore che di altri programmatori.

Nelle righe successive, dalla 14 alla 18, vengono assegnati dei valori alle variabili. L'assegnazione viene effettuata tramite il simbolo "=" che ha un significato ben diverso dal simbolo che tutti studiano a scuola in algebra matematica. Questo simbolo indica al programma che alla variabile posta a sinistra del simbolo deve essere dato il valore di ciò che è a destra. A parole, la riga 14 significa: "assegna il valore 10000 alla variabile deltat". Anche la riga 23 assegna un valore ad una variabile, ma il valore è il risultato di un'espressione e a parole questa riga significa: "assegna il valore della radice quadrata (sqrt) della somma (+) dei quadrati (**2) delle variabili posx e posy". Si tratta, ovviamente, dell'applicazione del teorema di Pitagora per il calcolo dell'ipotenusa di un triangolo rettangolo. Il calcolo viene eseguito usando i valori che le variabili hanno in quel momento. Come si può notare, la stessa cosa viene rieseguita alla riga 41, dove però i valori di posx e posy sono diversi, pertanto diverso risulta anche il valore che in quel momento del programma viene assegnato a r.

Il cuore del programma sono le righe che vanno dalla 30 alla 53, nelle quali vi sono una serie di istruzioni che vengono ripetute un certo numero di volte. Questo è possibile in quanto queste istruzioni sono all'interno di un ciclo do, che viene controllato dalla variabile i. Senza entrare troppo nella spiegazione, possiamo dire che il ciclo viene ripetuto 10000 volte. Il tempo necessario ad un normale PC per eseguire tutte le operazioni decine di migliaia di volte è di pochi secondi, ben diverso sarebbe se queste fossero eseguite a mano, anche con l'ausilio di una calcolatrice. Questa è una misura delle potenzialità dei computer per la comprensione di problemi scientifici.

L'analisi dei comandi all'interno del ciclo si lascia al lettore come esercizio. Si tenga presente che non è altro che la messa in pratica dell'introduzione al problema fatta nel paragrafo "La legge della gravitazione universale".

L'ultima operazione del ciclo do (riga 52) consiste nello stampare (comando "write") le variabili posx e posy. In altre parole, vengono scritte le coordinate della posizione del pianeta calcolate ad ogni passo. Questi valori sono i numeri di cui abbiamo bisogno per la nostra analisi e vengono stampati per ogni esecuzione del ciclo do.

Se, diversamente, si volesse solo avere il valore della distanza tra i due pianeti, la riga 52 sarebbe dovuta essere:

      write(*,*) r

Appare molto intuitivo il motivo di tutto ciò. Il lettore può, come esercizio, modificare il programma e vedere il risultato del cambiamento.

Analisi dei dati

Come detto, i dati ottenuti dal programma sono stati scritti nel file dati.dat. Di seguito, si riportano le prime dieci righe di questo file:

 1  1.52099701E+11  292900000.
 2  1.52098832E+11  585798848.
 3  1.52097391E+11  878695552.
 4  1.52095375E+11  1.17158886E+09
 5  1.52092787E+11  1.46447782E+09
 6  1.52089625E+11  1.75736115E+09
 7  1.52085889E+11  2.05023795E+09
 8  1.5208158E+11  2.34310707E+09
 9  1.52076698E+11  2.63596749E+09
10  1.52071242E+11  2.92881792E+09

Come si vede, vi sono due colonne di dati, dove ogni riga riporta a sinistra il valore della variabile posx e a destra quella di posy. Con un opportuno programma, è possibile plottare i dati ed avere una rappresentazione grafica di ciò che si è calcolato. Si è utilizzato Grace [9], di cui si è già parlato in un precedente articolo del giornale [10]. Dalla figura 2 si può vedere quello che è quasi impossibile percepire sfogliando il file dei dati: il nostro Pianeta effettua, attorno al Sole, un'orbita ellittica. Sebbene questo risultato sia ben noto, conferma la validità delle nostre manipolazioni matematiche sul problema fisico e dei calcoli effettuati dal nostro programma.

Figura 2: Rappresentazione grafica dei dati prodotti dal programma.

Conclusioni

Prima di passare alle conclusioni, mi preme ringraziare il professor Giorgio Pastore per il suo prezioso contributo a questo articolo, segnalandone il progetto [11] "Fare scienza al computer: Sperimentazione di un laboratorio di simulazione numerica ad uso degli insegnanti e degli allievi dell'ultimo biennio delle scuole superiori", il cui titolo è autoesplicativo sui contenuti.

Detto ciò, possiamo passare alle conclusioni. Questa è stata una particolare introduzione al g77 in chiave non professionale, ma come strumento di studio. Uno strumento valido se si pensa al percorso mentale che ogni studente deve fare per affrontare i tanti problemi che i suoi studi gli propongono a livello teorico. Un percorso che richiede un'approfondita conoscenza del problema, collegamenti interdisciplinari e che, per tali motivi, è fortemente stimolante. Approfondimenti sul FORTRAN possono essere fatti in rete che, come al solito, è ricca di materiale. In particolare, sono reperibili varie guide e manuali, come ad esempio la Professional Programmer's Guide to Fortran77 di Clive G. Page [12]. Di sicuro interesse possono essere le FAQ del newsgroup comp.lang.fortran [13].

Questo strumento è anche il frutto del lavoro di tanti programmatori che, nell'ambito del progetto GNU, lavorano con uno spirito di collaborazione ed è a disposizione di chiunque nella totale libertà che solo i progetti sotto licenza GPL possono dare.

Riferimenti bibliografici

[1] http://www.gnu.org/software/fortran/fortran.html

[2] http://www.fsf.org/

[3] http://www.fsf.org/licensing/licenses/gpl.html

[4] http://www.vim.org/

[5] http://g95.sourceforge.net/

[6] http://gcc.gnu.org/fortran/

[7] http://gcc.gnu.org/

[8] http://gcc.gnu.org/java/

[9] http://plasma-gate.weizmann.ac.il/Grace/

[10] http://journal.pluto.it/pj0109/Grace-articolo.html

[11] http://www.democritos.it/education/iatfys.php

[12] http://www.star.le.ac.uk/cgp/prof77.html

[13] http://www.faqs.org/faqs/fortran-faq/

<- Fedora Core 3: un'installazione difficile - Indice Generale - Copertina - OpenVPN, VPN in pochi passi ->