Category Archives: Fortran

Il linguaggio di programmazione Fortran

Traduzioni da un linguaggio a un altro

btf3

Un discorso che devo prendere un po’ alla lunga, da lontano, roba di prima di voi giovani, ma poi arrivo al dunque, cioè all’oggi, promesso.

In Fortran π, pi, il rapporto tra le misure della circonferenza e il raggio non è definito. Ma tutti sanno che l’arcotangente di 1 rad vale π/4.
E allora diventa tutto semplice, gli antichi (non gli antichissimi ma non cambia molto) usavano il Fortran 77 (1978).
Vediamo se ricordo bene perché possono esserci diversi modi di scrivere, tenendo presente che i tipi sono impliciti e allora

      pi1 = 4.0 * atan(1.0)
      pi2 = 4.0 * atan(1)
      pi3 = 4 * atan(1.0)
      pi4 = 4.0 * datan(1.0)
      
      print *, pi1, pi2, pi3, pi4
      end

ap0

Partito molto male: 2 errori –no uno solo in 2 varietà– correggo.

La riga pi2 = 4.0 * atan(1) è irrecuperabile, a correggerla diventa uguale alla precedente.
La riga pi4 = 4.0 * datan(1.0) è più intrigante. Intanto la D iniziale sta per double precision, quello che in C si chiama double.

      pi1 = 4.0 * atan(1.0)
      pi3 = 4 * atan(1.0)
      pi4 = 4.0 * datan(1.d0)
      pi5 = 4.0d0 * datan(1.d0)
      pi6 = 4.0d0 * datan(dble(1.0))
      
      print "(f12.10)" , pi1, pi3, pi4, pi5, pi6
      end

ap1

OK, ottengo gli stessi valori. Notare la promozione a double nella moltiplicazione per il calcolo di pi4. Notare anche la funzione dble, cast a double direbbe il C-ista.

Ma non va come previsto: otteniamo sempre lo stesso numero di cifre corrette, sia usando i real*4 che i real*8 corrispondenti a float e double del C.
Devo dichiarare il tipo quando non rientra in quelli di default. Infatti questo viene automaticamente assegnato a integer (2 bytes nel Fortran IV, 4 nel 77) se la prima lettera della variabile è nell’intervallo [I-N], real*4 per tutte le altre. Si ricorda inoltre che al di fuori delle stringhe quotate si ha una conversione automatica al maiuscolo.
Ecco:

      double precision pi4
      real * 8 pi6

      pi1 = 4.0 * atan(1.0)
      pi3 = 4 * atan(1.0)
      pi4 = 4.0 * datan(1.d0)
      pi5 = 4.0d0 * datan(1.d0)
      pi6 = 4.0d0 * datan(dble(1.0))
      
      print "(f12.10)", pi1, pi3, pi4, pi5, pi6
      end

ap2

OK😀 Ah! dimenticavo: gli spazi non sono significativi, real*8 e real * 8 sono sinonimi, lo sarebbe pure r e a l * 8.

Lo so che nessuno usa più il Fortran, neanch’io, ma devo raccontare questo altrimenti poi non si capisce.

Neanche e, quella di Euler, dei logaritmi naturali, è predefinita. Ma facciamo come per π:

      real * 8 e2

      e1 = exp(1.0)
      e2 = dexp(1.d0)      

      print "(f12.10)", e1, e2
      end

ae

OK😀 notare che al posto di 1.0 posso scrivere 1. come posso scrivere .1 per 0.1; i vecchi lo facevano sempre.

Finora niente di sconvolgente ma ecco un passo ancora

      e = exp(1.0)
      pi = 4.0 * atan(1.0)
      ie2 = 2 * e
      ipi2 = 2 * pi
        
      print *, e, pi
      print *, ie2, ipi2
      end

app

OK, i cast a integer*4 tronca i real*4. Ma attenzione, questo avviene –ovviamente– dopo la moltiplicazione. Ovviamente è diverso da questo caso:

      e = exp(1.0)
      pi = 4.0 * atan(1.0)
      ie = e
      ipi = pi
      ie2 = 2 * ie
      ipi2 = 2 * ipi
        
      print *, e, pi
      print *, ie2, ipi2
      end

iapp

Se si vuole l’arrotondamento all’intero più vicino ci sono due modi:

      e = exp(1.0)
      ie1 = e + 0.5
      ie2 = int(e + 0.5)
      ie3 = nint(e)

      print *, ie1, ie2, ie3
      end

round

Siccome la conversione da real a integer avviene automaticamente per troncamento i vecchi usavano sommare 1/2 al valore da convertire. La funzione int() in questo caso è implicita per ie1 e esplicita per ie2. Esiste la funzione di arrotondamento nint() ma è usata raramente (mai).

Finora ho barato; no, non barato ma omesso di dire tutta la verità. Alcune funzioni intrinseche del Fortran sono da sempre (OK, quasi) polimorfe.
Già vista datan() ma risulta più chiaro con la funzione max():

      m1 = max(3, 1, 4, 1, 5)
      m2 = max0(3, 1, 4, 1, 5)
      m3 = max1(3.0, 1.0, 4.0, 1.0, 5.0)
      rm2 = amax0(3, 1, 4, 1, 5)
      rm3 = amax1(3.0, 1.0, 4.0, 1.0, 5.0)
      print *, m1, m2, m3
      print *, rm2, rm3
      end

fsp

Sì, conta la lettera iniziale (A per i real, D per i double) e l’eventuale numero finale (0 per integer, 1 per real). Come se non fosse abbastanza ci sono state variazioni (tante) tra le versioni 77 e 90 (e successive), RTFM (cit.).

Ci sarebbe da raccontare molto su common ma salto.
Come pure salto tutto il capitolo sulle matrici, c’è Octave visto recentemente.

Invece una cosa sempre in agguato: il passaggio di variabili ai sottoprogrammi. In Fortran ci sono le function, funzioni come quelle di C e Python e quasi tutti i linguaggi. Ma ci sono anche le subroutines –anzi sono nate prima– sottoprogrammi chiamate con la keyword call che non ritornano un valore come le funzioni ma modificano i dati passati. Tutti, perché vengono passati per indirizzo (questo è vero solo per le versioni non nuovissime ma il codice è spesso piuttosto attempato). Un esempio minimo:

*     main
      n = 5
      call foo(n, nq)
      print *, n, nq
      call foo(n, nq)
      print *, n, nq
      call foo(n, nq)
      print *, n, nq
      end
      
      subroutine foo(k, l)
      k = k + 1
      l = k ** 2
      end

pass

Sembra facile, traduco in Python:

### versione iniziale, non funziona

def foo(k, l):
      k = k + 1
      l = k ** 2

#  main
n = 5
foo(n, nq)
print(n, nq)
foo(n, nq)
print(n, nq)
foo(n, nq)
print(n, nq)

py-pass

Ehmmm… no; non basta aggiustare la sintassi delle singole righe. Il valore di l | nq dev’essere ritornato espressamente. Inoltre è inutile passarlo tra gli argomenti:

### ancora non funziona

def foo(k):
      k = k + 1
      l = k ** 2
      return l

#  main
n = 5
nq = foo(n)
print(n, nq)
nq = foo(n)
print(n, nq)
nq = foo(n)
print(n, nq)

py-pass1

No! non ci sono ancora: ovvio, gli passo sempre lo stesso parametro n = 5; la modifica all’interno di foo() è locale e in main non viene vista, correggo:

def foo(k):
      k = k + 1
      l = k ** 2
      return l

#  main
n = 5
nq = foo(n)
n += 1
print(n, nq)
n += 1
nq = foo(n)
print(n, nq)
n += 1
nq = foo(n)
print(n, nq)

py-pass2

OK! ma il codice è orrendo! Non voglio vedere niente di simile😦

Come visto anche per un caso molto banale la traduzione automatica, riga per riga non è possibile. Ci sarebbe da dare la colpa al fortrainer che se avesse fatto le cose per bene… scrivendo la subroutine foo() così:

      subroutine foo(k, l)
      j = k + 1
      l = j ** 2
      end

la variabile k sarebbe stata solo di input –o anche intent(in) per le versioni nuove del Fortran– non sarebbe ritornata modificata e il main avrebbe dovuto essere scritto meglio.

Adesso provate –Gedankenexperiment– ad applicare questo a codice reale👿 Ma secondo me va fatto. Prima o poi. Non necessariamente da me.
:mrgreen:

Il Fortran iniziale, I/O, FORMAT

wo-the-wizard

Il titolo è fuorviante, riassumo una conversazione con un collega dei tempi antichi riguardo al post Lo sviluppo iniziale dei linguaggi di programmazione.
Naturalmente riguarda il Fortran l’unico linguaggio che c’era sul ‘puter che usavamo. Il programma postato da Knuth e Trabb Pardo (K&T) è per la versione 0 (zero), forse solo virtuale. La prima versione con cui ho avuto a che fare è la IV, anzi un’evoluzione compatibile con la 66 ma c’era codice scritto per la versione II (1958), trasferito da schede a nastri o dischi o semplicemente ricopiato dal listato.
In ogni caso il programma postato da K&T deve essere aggiornato (ehi! 50 anni!): non solo le funzioni non hanno (dalla preistoria) il nome che deve terminare con F ma non ho nemmeno un lettore di schede😦
E non saprei nemmeno trovare le schede😦 e la perforatrice😦

forig

Però il programma prevede la lettura di dati formattati in modo preciso, cosa che da terminale non è agevole (bisogna essere proprio pistini) ma si può fare un formattatore, come i millemila come si usava una volta.

Ah! un altro problema: la gestione dei file: allora si aprivano da OS, il programma lo dava per scontato e alla fine si chiudevano, sempre da OS. I comandi Prime erano rispettivamente OPEN seguito dai parametri che non ricordo ma roba tipo 1 1 (unità 5 in lettura [aggiornamento: il secondo parametro è da verificare, probamilente non come pensavo]) e CLOSE 1 (unità 5); siccome di solito si dovevano chiudere più unità c’era il comodo CLOSE ALL. Su Prime (come con il DEC) i comandi potevano essere abbreviati, nel nostro caso al posto di OPEN bastava O, era il comando più importante iniziante con quella lettera, e C stava per CLOSE.
Quindi, se ricordo bene, il mio command-file per l’esecuzione del programma FMT sarebbe stato così:

O DATI-FREE 1 1
O DATI 2 2
SEG FMT
C ALL

L’istruzione di lancio SEG FMT valeva nell’ipotesi che il programma fosse stato compilato per la versione segmentata (prevista per i programmi grossi) usata (a sproposito) abitualmente. Però se uno era razionale vedeva che bastava la versione normale, anche più veloce quindi la riga sarebbe stata RUN FMT, RUN abbreviabile in R. OK, mode nostalgia OFF.

La versione corrente usa le normali istruzioni OPEN e CLOSE, non disponibili allora. Per evidenziare il nuovo uso il minuscolo, adesso equivalente al maiuscolo, allora vietato.

Ecco fmt.f

C FMT legge in free-format e scrive formattato

      DIMENSION A(11)
      
      open(1, file='dati-free')
      READ(1, *) A
      close(1)
      
      open(2, file='dati')
      WRITE(2, 2) A
      close(2)
2     FORMAT(6F12.4)

      END

Aggiornamento: per il Prime le prime 4 unità di I/O erano predefinite e riservate per cui avrei dovuto usare numeri più alti, tipicamente 5 e 6 al posto di 1 in OPEN e 2 in WRITE.

E questo è il file di dati dati-free

1 10 100 1000 10000 98765 987654
987654.3 2.7818 3.1415 -42

Compilo ed eseguo

f0

ottenendo il file dati

      1.0000     10.0000    100.0000   1000.0000  10000.0000  98765.0000
 987654.0000 987654.3125      2.7818      3.1415    -42.0000

Uh! visto il secondo dato della seconda riga? è 987654.3125 invece del previsto 987654.3000😳 ma il formato REAL*4 (quello usato per le variabili inizianti con una lettera fuori dell’intervallo [I-N] riservato agli interi) garantisce 7-8 cifre significative.

Per essere proprio sicuri-sicuri si può fare la verifica, con il file verifica.f

C VERIFICA

      DIMENSION A(11)
      
      open(1, file='dati')
      READ(1, 2) A
      close(1)
2     FORMAT(6F12.4)

      DO 3 I = 1, 11
        PRINT *, A(I)
3     CONTINUE     

      END

f1

OK, ci siamo. A questo punto torno al programma di K&P, opportunamente aggiornato (tpk.f).

C     THE TPK ALGORHTM, FORTRAN STYLE
*     aggiornato, modifiche in minuscolo

      FUN(T) = SQRT(ABS(T)) + 5.0**3
      DIMENSION A(11)
1     FORMAT(6F12.4)
      open(1, file = 'dati') 
      READ(1, 1) A
      close(1)
      DO 10 J =1, 11
      I = 11 - J
      Y = FUN(A(I+1))
      IF(400.0 - Y) 4, 8, 8
4     write(*, 5) I
5     FORMAT(I10, 10H TOO LARGE)
      GO TO 10
8     write(*, 9) I, Y
9     FORMAT(I10, F12.7)
10    CONTINUE
      STOP 52525
      end

f2

OK😀

C’è ancora una cosa emersa ieri: le funzioni, quando non ci sono si fanno. A dire tutta la verità i fortrainers preferivano –di gran lunga– le SUBROUTINEs ma questo sarebbe un altro rant, lunghissimissimo😉

Era tipico di IBM fornire funzioni che non si avevano (no, niente nomi). Per gioco supponiamo non fossero definite ABS e SQRT.
La prima è banale (e lasciata come esercizio) mentre la seconda è leggermente più impegnativa. Siccome il mondo è vario da adesso in poi continuo con Python 3.x, il Basic dei nostri giorni.

Invece di limitarci alla radice quadrata conviene probabilmente essere più generali; che ne dite di usare i logaritmi? OK, la Wiki (come faremmo senza?) ci dice tutto. Illo tempore probabilmente avrei cercato sul Dwight –hey! c’è online–, c’è il capitolo Logarithmic Functions, p.130
log-x
Per l’esponenziale c’è Exponential Functions a p.125

exp-x

Quindi ecco myfunc.py:

# modulo per Python 3.x; 
# per 2.x importare division from tuture

my_e = 2.71828182845904523536 

def my_log(x):
    r = (x - 1) / (x + 1)
    d = r
    v = d
    c = 1
    while d > 1e-8:
        c += 2
        t = 1
        for i in range(1, c + 1):
            t *= r
        d = 1/c * t
        v += d

    return round(2 * v, 8)

def my_exp(x):
    fact = 1
    t = x/fact
    e = 1 + t
    num = x
    c = 1
    while t > 1.e-8:
        c += 1
        fact *= c
        num *= x
        t = num / fact
        e += t
    st = str(e)
    ppoint = st.index(".")
    if ppoint < 8:
        e = round(e, 8 - ppoint)
    else:
        e = round(e, 1) 

    return e
    
def my_sqrt(x):
    return my_exp(1/2 * my_log(x))

py0

:mrgreen:

Variabili: uso e abuso

CiEUn tweet di un programmatore funzionale che riporta un tweet che –OK, è qui.

v0

OK, pythonisticamente è corretto. [2] è una lista contenente l’intero 2.

v1

In Python, come in taaanti linguaggi ‘normali‘ le variabili sono –come dice il nome– variabili e possono essere ridefinite.
Ma forse Python esagera, consideriamo il ciclo seguente:

v2

Sembra tutto OK, vero?
Ma forse… approfondiamo:

v3

Uhmmm… qui l’indice del ciclo viene modificato all’interno del ciclo stesso ma ciò non turba l’interprete che al giro successivo sa di dover usare l’elemento successivo della lista. La stessa cosa capita se usiamo range(), che infatti produce la lista usata dal ciclo:

v4

OK😀

Fin qui è tutto abbastanza scontato, tranne per i nuovi che usano solo linguaggi strettamente funzionali (un po’ li invidio) ma a me che sono vecchio è tornata in mente una cosa di tanto tempo fa. Adesso ve la racconto.
Illo tempore ogni computer aveva il suo sistema operativo, proprietario e di solito diverso da tutti gli altri. Anche i compilatori, di sicuro il Fortran e pochi altri. Quando è nato il PC IBM con il DOS di Microsoft è successo che (numero dei sistemi operativi)++, nient’altro; e peraltro tutt’altro che buono. Ma sto divagando.

Riprendo l’esempio che appalla i funzionali, in Fortran 77, quello corrente di allora si ha (uso i al posto di x per le convenzioni implicite del linguaggio):

* c0.f
      i = 1
      do 10, i = 2, 4
          print *, i
10    end do
      print *, i      
      
      end

v5

OK, a parte s/for/do/ e * per dire ‘default’ in print tutto come previsto.

Occorre anche considerare una cosa che tutti i fortrainers davano per scontato: all’uscita del ciclo il contatore è incrementato. Ma è inaffidabile, tra l’altro varia rispetto alla versione precedente del compilatore (Fortran IV).
Vediamo se si può, come in Python, …

* c1.f
      i = 1
      do 10, i = 2, 4
          i = 2 * i
          print *, i
10    end do
      print *, i
      
      end

v6No, nope! come previsto il compilatore GNU segnala la violazione.

Ci può essere un caso più sottile, quello capitatomi (io ero il debugger):

* c2.f
      i = 1
      do 10, i = 2, 4
          call doppio(i)
10    end do
      print *, i
      
      end
      
      subroutine doppio(n)
      n = 2 * n
      print *, n
      end

v7

OOPS!😳 cos’è successo? I fortrainers lo sanno: le variabili vengono sempre passate per indirizzo a funzioni e subroutines. In questo caso l’indice del ciclo viene modificato all’interno del ciclo stesso ma il compilatore non può saperlo. È cura del programmatore –che queste cose le dovrebbe sapere– diminuire il carico di lavoro al debugger, sempre che il bug venga rilevato, possibilmente presto😀

Al programmatore era stato senz’altro stato detto anche “cambia il nome alle variabili” che è una frase ce dev’essere interpretata e tradotta in linguaggio normale. Vediamola applicata:

* c3.f
      i = 1
      do 10, i = 2, 4
          call doppio(i)
10    end do
      print *, i
      
      end
      
      subroutine doppio(k)
      n = k
      n = 2 * n
      print *, n
      end

v8OK😀 C’era un altro motivo per cambiare le variabili, cosa che oggi probabilmente ha meno senso: le variabili locali sono near pointers mentre le globali sono far pointers. Non credo di dire cose nuove, ma forse dimenticate al momento opportuno😉

:mrgreen:

Breve paragone tra l’“overloading” in C++, C(11) e Fortran(90)

fvc

In C++ possiamo scrivere codice come questo:

#include <iostream>
 
int miafunz(int a, int b)
{
    std::cerr << "int, int\n";
    return a + b;
}
 
int miafunz(double a, int b)
{
    std::cerr << "double, int\n";
    return 2.0*a + 2.5*b;
}
 
int main()
{
    double d0 = 1.0;
    int di = 2;
    std::cout << miafunz(d0, di) << "\n";
    std::cout << miafunz(5, di) << "\n";
}

Lo compiliamo e lo eseguiamo:

fvc-compilesegui

Quello che accade dietro le quinte è che il compilatore in realtà genera due funzioni diverse e “capisce” qual è quella giusta da chiamare basandosi sul tipo degli argomenti.

Per il programmatore è tutto trasparente, mentre compilatore e linker hanno una gatta da pelare in più. Il problema è questo: le due funzioni sono in realtà diverse, pezzi di codice diverso, e ciascuna deve avere associato un simbolo diverso1, anche se il programmatore usa solo il nome che ha deciso (miafunz).

Se non è definita e imposta un’ABI specifica per tutte le implementazioni del C++, ci possono essere “incomprensioni” perché ogni compilatore “decora” i nomi delle funzioni come meglio crede. (Cfr. name mangling)

Detto altrimenti: il compilatore C++ si deve inventare due nomi per le due funzioni, laddove noi usiamo un unico nome. È abbastanza ragionevole ipotizzare che i nuovi nomi siano costruiti a partire dal nome da noi scelto, secondo un certo schema2 — ma non è affatto necessario che sia così. Nel caso di g++ (4.7.2),

fvc-readelfgcc

ed è lo stesso con clang++ (3.0), per la cronaca. È buono che ci sia una certa convergenza, ma purtroppo in realtà non possiamo farci affidamento. Invece di spiare l’ELF potete dire al g++ di fermarsi alla generazione del codice assembly (opzione -S) e dare uno sguardo ai simboli usati.

fvc-asmsym

Il C non ha di questi problemi

Il C non ha di questo problemi perché non permette l’overloading. Se scrivete

#include <stdio.h>

int miafunz(int a, int b)
{
    fprintf(stderr, "int, int\n");
    return a + b;
}

int miafunz(double a, int b)
{
    fprintf(stderr, "double, int\n");
    return 2.0*a + 2.5*b;
}


int main(void)
{
    double d0 = 1.0;
    int di = 2;
    printf("%d\n", miafunz(d0, di));
    printf("%d\n", miafunz(5, di));
    return 0;
}

il compilatore si ribella,

fvc-gccribelle

perché stiamo ridefinendo una funzione (ho usato gcc 5.2 perché evidenzia meglio gli errori, in stile clang). In C i simboli delle funzioni devono essere unici (all’interno del loro campo di visibilità3).

Potremmo scrivere

int miafunz_ii(int a, int b)
{
    fprintf(stderr, "int, int\n");
    return a + b;
}

int miafunz_di(double a, int b)
{
    fprintf(stderr, "double, int\n");
    return 2.0*a + 2.5*b;
}

ma poi starebbe a noi scegliere la funzione “giusta” per gli argomenti giusti. Fattibile, ma potenziale fonte di errori e comunque appesantisce la lettura (e la scrittura…). Un esempio più concreto è quello di certe funzioni matematiche per le quali in effetti nella libreria esistono diverse versioni a seconda del tipo di argomenti usati e restituiti; p.es. sin (double), sinf (float), sinl (long double).

Il C11 è venuto incontro a questa esigenza tramite una macro, _Generic, il cui uso però risulta limitato, cioè non dà la stessa elasticità dell’overloading del C++: in pratica è utile solo per tgmath.h (type generic math).

Riciclo un esempio da un mio post su un altro blog.

Il Fortran è diverso

Quando le persone sentono Fortran credono di essere tornati indietro nel tempo: è un linguaggio vecchio, ormai in disuso, è brutto, ecc… Non è così. Di sicuro è “di nicchia”. Ma comunque, a differenza di altri linguaggi di quell’epoca, è andato avanti e si è ammodernato. L’ultima versione dello standard è il Fortran 2008, ma è pianificato il Fortran 2015.

Una delle caratteristiche del Fortran moderno sono le interfacce. In pratica il Fortran ha lo stesso problema che il C(11) ha in parte risolto con _Generic; però la soluzione è molto più… generica!… elastica, elegante e (a mio immodesto avviso) intelligente.

Dettagli sintattici a parte, il programmatore non lascia al compilatore l’onere di generare i nomi delle funzioni: egli scrive le diverse funzioni come deve (p.es. miafunz_di, miafunz_ii) e specifica che miafunz è un’interfaccia a quelle funzioni (il compilatore di solito è nelle condizioni di poter dedurre da sé la funzione giusta in base al tipo degli argomenti, proprio come avviene per il C++; nei casi in cui non dovesse essere così, nell’interfaccia dovremmo specificare per intero gli argomenti e i loro tipi).

Riciclo di nuovo un esempio da un altro blog; nell’esempio my_sum è definita come interfaccia per le funzioni my_sum_r, my_sum_i, ecc…


  1. A codice oggetto generato, non è sempre necessario. Per esempio nell’esempio specifico in teoria non c’è nessun bisogno che i simboli delle funzioni rimangano, perché di fatto non ci interessa esportarli. Tant’è che tali simboli possono essere rimossi, per esempio con strip. Il discorso sarebbe diverso se stessimo scrivendo del codice che deve essere linkato (una libreria statica o dinamica, per esempio).

  2. Una soluzione ragionevole (perché intellegibile) è quella di codificare nel nome della funzione i tipi degli argomenti; vediamo infatti che miafunz(double,int) diventa _Z7miafunzdi e miafunz(int,int) diventa _Z7miafunzii.

  3. O «confine entro cui sono visibili». È la mia traduzione approssimativa di scope. Non voglio affrontare questioni linguistiche sottili, ma la traduzione “visibilità” non mi sembra del tutto corretta. La visibilità (visibility) è la proprietà di essere visibili (la proprietà di un oggetto di poter essere visto); il termine scope fa riferimento a quanto è “estesa” questa proprietà, cioè fin dove vale o non vale per un determinato oggetto; dunque si riferisce all’“area di pertinenza” di una proprietà (in generale, e della visibilità nello specifico). La traduzione “visibilità” di per sé non veicola l’idea di “estensione”. Prendiamo per esempio la frase «the name must be unique in its scope». La traduzione «il nome deve essere unico nella sua visibilità» non corrisponde all’originale; in italiano siamo costretti ad usare una perifrasi: «il nome deve essere unico» «all’interno dell’ambito in cui è visibile», o «all’interno dell’estensione della sua visibilità», o «entro la portata della sua visibilità», o «nel suo campo di visibilità» (forse l’opzione migliore), o altre espressioni simili. Per quanto mi riguarda, in casi come questi non ho nessun problema ad usare termini inglesi nell’italiano settoriale.

Argomenti passati a funzione – 2

sicNon so se solo noi… (inteso come quelli con cui lavoravo illo tempore) ma in seguito al post Argomenti passati a funzione ho avuto una sessione di ricordi con il mio boss di allora. Riassumo.

Intanto la funzione corretta (cquad.f):

1
2
3
4
5
6
7
* quad restituisce il quadrato del parametro - corretta
       function quad(x)
       y = x
       y = y * y
       quad = y
       
       end

Verifico, con questo main (c1.f):

1
2
3
4
5
6
7
8
* versione corretta #1
      x = 5
      qx = quad(x)
      print *, ' il quadrato di', x, ' vale ', qx
      qx = quad(8.0)
      print *, ' il quadrato di 8.0 vale ', qx
      
      end 

arg-c1

Nella funzione si copia il parametro in una variabile locale e si modifica solo questa. È così usabile come previsto nel main; qui a raccontarla tutta inizialmente avevo scritto 8 invece di 8.0

Però, mi fa osservare il Signor Capo, la prassi era quella di usare le subroutines, ecco squad.f:

1
2
3
4
       subroutine quad(x)
       x = x * x
       
       end

che verifico con il main c2.f:

1
2
3
4
5
6
* versione corretta #2
      x = 5
      call quad(x)
      print *, ' il quadrato di 5.0 vale ', x
      
      end 

arg-c2

OK, quasi. Ma questo era il modo standard di operare (almeno da noi). Notare che call quad(8.0) non avrebbe senso, non ottengo il risultato (sì un errore, invece).

Una googlata è stata molto fruttuosa: ho trovato il manuale di Prime F77 che riporta in due punti un warning a proposito di questo problema. Ma si sa che nessuno legge i manuali👿

pagina 5-1:

5-1

e pagina 5-4

5-4

Abbiamo pensato anche a intent ma ci risulta posteriore, Fortran 90. Può darsi che qualche implementazione l’avesse già allora ma dalla ricerca con Google non ho trovato niente. Trovato invece cose interessanti:

Il commento di Piero Patteri, esaminato con più cura di quanto ho risposto, si riferisce a un altro aspetto ancora, interessa a qualcuno? Devo parlarne? Secondo il Boss il Fortran si usa ancora, anche se molto meno, ci sono concorrenti, tanti. Anche Matlab.


L’immagine lassù è del manuale del Fortran. Trovata su Ebay (in realtà KiJiJi), con questa nota: Nave scuola per migliaia di programmatori, questo volume del 1981 (riedizione di quella del 1971) ha costituito un punto di riferimento per la nascente attività di programmazione per docenti e ricercatori.
E pare sia ancora nel catalogo di Zanichelli. No, solo nel catalogo storico.

Argomenti passati a funzione

pr1meL’altro giorno ero alle prese con la filosofia della programmazione funzionale, cosa di moda ultimamente ma –scopro– vecchia ormai. E mi è tornato in mente un caso di tento tempo fa; quasi-quasi lo racconto. Anzi vado, così verifico se nel frattempo i compilatori sono migliorati. E pensa te che ho dovuto reinstallare il Fortran. Dai c’è di peggio. Per esempio la politica. Il caso è molto semplice. Pochissime righe di codice, antico; qui uso il 77 ma andrebbe bene anche il IV (1966). L’unica cosa da tener presente è che gli argomenti a una funzione sono sempre passati per indirizzo (reference). Quindi quanto si torna dalla chiamata se le variabili usate come parametri sono state modificate risulteranno modificate. Mi prendo una piccola libertà per semplificare il problema: i programmatori allora prediligevano di gran lunga le subroutine, sottoprogrammi che non ritornano esplicitamente un valore ma si affidano ai side effects; invece qui uso una function, quelle del C, c’erano da sempre anche in Fortran. Creo una piccolissima (e ridondante per visualizzare il caso in oggetto) funzione quad() che ritorna il quadrato dell’argomento che le viene passato (quad.f):

1
2
3
4
5
6
7
* quad restituisce il quadrato del parametro
       function quad(x)
       
       x = x * x
       quad = x
       
       end

Iessainou! pessima. Se uno programma così lo licenzio a calcinculo. CMQ, la prima riga è un commento, la quinta specifica il valore da ritornare, calcolata nella riga precedente, inutile, anzi dannosa. Inutile in quanto il calcolo poteva farlo sulla stessa riga che definisce il valore della funzione. Questa funzione può essere chiamata da un main, così (bug1.f):

1
2
3
4
5
6
      x = 5
      print *, 'calcolo il quadrato di X, con X =', x
      qx = quad(x)
      print *, ' il quadrato di', x, ' vale ', qx
      
      end 

Compilo ed eseguo: bug1 Oops! anzi atz!👿 Ma lo sapevo, x è stata modificata nella funzione. C’era (e forse c’è ancora) di peggio, consideriamo questo caso (bug2.f):

1
2
3
4
5
6
      x = 5
      print *, 'calcolo il quadrato di X, con X =', x
      qx = quad(5)
      print *, ' il quadrato di', x, ' vale ', qx
      
      end 

Compilo, quad è già compilata, devo solo linkarla, risparmio tempo, i vecchi a queste cose erano attenti😉 ed eseguo: bug2 OK, ecco. Come previsto. A questo punto via con il debug. Cioè no, anzi tutto come previsto. Ma allora…👿

L’avvenire del Fortran e i suoi successori

fortran

Forse mi sono lasciato prendere la mano, non so se sia così necessario questo post. Ma siccome l’ho scritto lo pubblico😕
Tutto è partito da qui: Why Scientists Are Still Using FORTRAN in 2014 che rimanda a Ars technica, qui: Scientific computing’s future: Can any coding language top a 1950s behemoth? dove troviamo

But almost universally, the language in which these simulation codes are written is Fortran, a relic from the 1950s.

e

Wherever you see giant simulations of the type that run for days on the world’s most massive supercomputers, you are likely to see Fortran code.

Poi dice che subito dopo viene inventato il Lisp. Che non ha successo per la sua stranezza:

Weirdness, because the prefix notation used by Lisp made expressions in that language look a lot less like normal mathematics than did math rendered in Fortran. (Fortran stood, after all, for FORmula TRANslator.) A normal chemist or engineer is far more comfortable with y = (a + b)/c than with (setf y ((/ (+ a b) c))).

Ok, l’ho sentita tante volte che mi verrebbe voglia di ignorarla. Invece no! Perché è (quasi) tutto lì. Davvero, ve ne faccio il riassunto.

(setf y (/ (+ a b) c))

Parto dalla parentesi più interna: (+ a b); questa è una lista il cui primo elemento è una funzione. Questa funzione opera sui successivi elementi della lista e li somma. Questa somma (chiamiamola A+B viene restituita: Lisp restituisce sempre l’ultima espressione calcolata.
Quindi la nostra espressione diventa:

(setf y (/ A+B c))

dove (/ A+B c) è una lista in cui il primo elemento è la funzione / che restituisce il rapporto del secondo elemento della lista per tutti i successivi.

Il valore restituito viene assegnato alla variabile y, tramite setf (setf sta per set (quote e sì, si può scrivere anche setq o set ' –non sono esattamente la stessa cosa ma il risultato è lo stesso).

Uh! visto? il Lisp è uniforme. Si comporta sempre allo stesso modo: il primo elemento della lista (chiamato CAR o first) indica come trattare il resto della lista (CDR o rest) e restituisce il risultato.
Per cui il Lisp è fatto da liste dentro liste dentro liste ancora.  Tutto lì (ok, sto semplificando). Non so voi ma a me sembra molto semplice.

Per contro il Fortran riproduce la sintassi matematica che si usa a scuola:

y = (a + b) / c

Sì, lo ammetto, siamo abituati così e per un calcolo semplice come questo meglio il Fortran. Ma per un caso non banale si vede la differenza: con il Lisp ho una concatenazione di liste e con il Fortran una serie di assegnazioni (attenzione il segno = si deve leggere “diventa“, dice il Wirth (che scrive :=)).

11c

Poi una cosa ancora: io sono vecchio e non so nemmeno se si usano ancora ma quando frequentavo il Poli, roba anni ’70, i migliori usavano le calcolatrici HP, con RPNreverse polish notation— quasi come il Lisp (cioè no, esattamente il contrario), nel nostro caso si sarebbe scritto

a
b
+
c
/

Con risparmio di digitalizzazione di tasti –non servono le parentesi–, si consumano meno i diti e –sopratutto– i risultati intermedi sono disponibili nello stack (mi dicono che dovrei usare il termine “pila”).

Ma tornando alla programmazione allora quando è uscito il Fortran era una cosa nuova di pakka, un compilatore per linguaggio facile per gli umani, alternativo alle criptiche istruzioni assembly. E sui calcolatori computer c’era solo il Fortran. Molto diverso da quello di oggi, istruzioni molto elementari, quasi da assembly. E tanti GOTO, spesso mascherati come IF.

Ok, fine dello sfogo, torno al post di Ars Technica secondo il quale ci sono tre candidati a sostituire il Fortran:

Haskell, bello, sintassi da matematici, funzionale (cioè odia le variabili), chissà;
Clojure, Lisp & JVM, sì c’è Java sotto; ho scritto qualcosa in passato;
Julia, nuovissimo, ne so molto poco.

Ma ce ne sono altri, secondo me: il C++, Java, Python e altri ancora. Tutti ampiamente utilizzati.
E i nuovissimi Go e Rust. Di Go sono a conoscenza di un programma funzionante, per il Web, Rust promette bene (ne parla spesso robitex su questo blog).
Inoltre proprio ieri vosto questo post: #golang: The Next Great Teaching Language. Il blog è nuovo ma promette bene. E il post è illuminante.

Invece il Fortran non è amato dai giovani (quelli che conosco), lo trovano vecchio, inusuale. E se possono lo sostituiscono con qualcosa di più pratico, per esempio Matlab (a volte basta la versione free Octave). Ma sto parlando di chi usa il computer come un mezzo complementare. E lo usa per AutoCAD, fogli di calcolo Excel e come macchina da scrivere.

Mi mancano i programmatori di una volta, quando il mondo era nuovo😉

Precisione iii

Corto circuito

Questo post è un bellissimo corto circuito: conosco Orlando perché è un amico del GuIT che tra l’altro ho conosciuto l’anno scorso di persona al GuIT meeting di Napoli. E proprio Orlando ha commentato il post di Juhan Precisione su questo stesso blog, di qualche giorno fa su alcuni calcoli.
Questo dimostra che gli appassionati non vedono l’ora di scoprire nuove cose sul proprio argomento preferito.

Chiudo il circuito

Eccomi dunque a chiudere il circuito informatico presentrando un programma in Go che esegue le stesse operazioni eseguite in awk ed in Fortran:

package main

import (
    "fmt"
    "math/big"
)

func main() {
    s, u := big.NewInt(2), big.NewInt(1)
    for i := 1; i <= 10; i++ {
        // a = s - 1
        a := big.NewInt(0).Set(s).Sub(s, u)
        // s*a + 1 = s*(s - 1) + 1
        s.Mul(a, s).Add(s, u)
        fmt.Println(i, s)
    }
}

Tecnicamente, vorrei farvi notare che è possibile in Go concatenare le chiamate dei metodi come per esempio si legge alla riga 12 del codice sorgente, se il metodo stesso ritorna lo stesso oggetto.

Ho utilizzato il pacchetto dei grandi numeri big necessario perché dopo la sesta iterazione i numeri superano il limite massimo del tipo int64. Certo avrei potuto utilizzare il tipo uint64 ma potevo chiudere il circuito soltanto esagerando alla grande portando il ciclo fino a 10 iterazioni… ecco di seguito il risultato:

1 3
2 7
3 43
4 1807
5 3263443
6 10650056950807
7 113423713055421844361000443
8 12864938683278671740537145998360961546653259485195807
9 165506647324519964198468195444439180017513152706377497841851388766535868639572406808911988131737645185443
10 27392450308603031423410234291674686281194364367580914627947367941608692026226993634332118404582438634929548737283992369758487974306317730580753883429460344956410077034761330476016739454649828385541500213920807

A questo punto, ci sta propio bene un
Alla Prossima!
R.

Precisione – 2

precisione
Il post di ieri, Precisione è stato visitatissimo! E Orlando mi chiede una cosa che è difficile far stare in un commento. E allora invece di un commento ecco un intero post.

Orlando pubblica un programma in Fortran questo:

program test
  real :: s=2.0
  do i=1,7
    s=s*(s-1)+1
    write(*,10) s
  end do
10 format(F30.1)
end program test

che produce questo risultato:

                           3.0
                           7.0
                          43.0
                        1807.0
                     3263443.0
              10650057179136.0
 113423716646946804535918592.0

simile ma non uguale a quelli degli script da me pubblicati.
Allora, vediamo. Comincio con modificare lo script di bc, facendogli scrivere i risultai dei singoli passi del ciclo, ecco:

s = 2.0
for (i = 1; i <= 7; i++) {
    s = s * (s - 1) + 1
    print s, "\n"
}
quit

3.0
7.0
43.0
1807.0
3263443.0
10650056950807.0
113423713055421844361000443.0

Dove si vede che le difformità sono negli ultimi due passi, questi:

3.0
7.0
43.0
1807.0
3263443.0
* 10650056950807.0
F 10650057179136.0
         ^ 
* 113423713055421844361000443.0
F 113423716646946804535918592.0
          ^

dove * rappresenta l’output di bc, F quello del Fortran e ^ la prima cifra diversa tra i due. Si vede subito che il numero di cifre per i linguaggi normali (anche il C, dice infatti Orlando) sono 7-8. Più che sufficienti in condizioni normali come concludevo nel post di ieri.
Naturalmente si può fare meglio. In particolare per il programma Fortran ditato sopra basta modificare la seconda riga, così:

real(16) :: s=2.0

e ottengo lo stesso risultato per Fortran e bc (e gli altri programmi, presumibilmente, anzi certamente).

In passato (l’anno scorso) ho raccontato un po’ di cose sul Fortran (per via che si stava facendo manutenzione a cose preistoriche), si trova tutto nell’indice di OK, panico, cercare per Fortran.

In particolare, siccome non ricordo certo cosa ho scritto a giugno del 2012 consiglio una scorsa qui e qui.

Resto disponibile per modifiche e chiarimenti, anzi questo post possiamo considerarlo un prepost😉
E, un’ultima cosa (Steve J. finiva sempre così), lo sapete che Orlando ha un sito tutto suo? Questo: Orlando Iovino’s Personal Website. Da finire, certo😉

f2py – usare codice Fortran in Python

jakeTutto è partito da qui: Numba vs. Cython: Take 2.

Jake VanderPlas è un tipo tosto, fareste bene a seguirlo. Si presenta così: “I’m an NSF post-doctoral fellow at University of Washington, where I research and teach in Astronomy, Astrophysics, and Machine Learning“.

Jake accenna a f2py. Ehi! potrebbe essere utile a chi ha vecchia roba scritta in Fortran, magari si riesce a metterle un’interfaccia grafica o fare altre mescite di quelle che, fa … beh! Chissà, vediamo.

Intanto salta fuori che f2py è un programma a se stante che compila il codice Fortran (usando il compilatore che hai, dev’essercene uno) creando una libreria (.so), quelle che Python chiama moduli.
Ma è anche un modulo di NumPy, cioè Python. Ed è tutto scritto in Python.
Allora la prima cosa da fare (è da tanto che mi dicevo di farlo) è installare NumPy: f2py è lì. Installare NumPy richiede 952MB di spazio su disco! OK, c’è tempo per un caffè lungo.

Fatto, viene installato anche SciPy, ne ho già parlato in passato. Per vedere comincio a riprodurre l’esempio preso dalla documentazione, questa è la subroutine in Fortan 77:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
C FILE: FIB1.F
      SUBROUTINE FIB(A,N)
C
C     CALCULATE FIRST N FIBONACCI NUMBERS
C
      INTEGER N
      REAL*8 A(N)
      DO I=1,N
         IF (I.EQ.1) THEN
            A(I) = 0.0D0
         ELSEIF (I.EQ.2) THEN
            A(I) = 1.0D0
         ELSE 
            A(I) = A(I-1) + A(I-2)
         ENDIF
      ENDDO
      END
C END FILE FIB1.F

Si crea il modulo Python con il comando:

f2py -c fib1.f -m fib1

e questo è il file di Python:

#!/usr/bin/python
# -*- coding: utf-8 -*-

import numpy as np
import fib1 #quello in Fortran

a = np.zeros(10, 'd')
fib1.fib(a)
print a

f1

OK, vediamo se si può migliorare: Fibonacci restituisce interi, anzi ai tempi di Leonardo Pisano (Fibonacci per i matematti) i float non erano conosciuti. E meno che mai in doppia precisione (come si dice in Fortran).

Le modifiche sono minime:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
C FILE: FIB2.F
      SUBROUTINE FIB(A,N)
C
C     CALCULATE FIRST N FIBONACCI NUMBERS
C versione modificata, uso interi
      INTEGER N
      INTEGER A(N)
      DO I=1,N
         IF (I.EQ.1) THEN
            A(I) = 0.0D0
         ELSEIF (I.EQ.2) THEN
            A(I) = 1.0D0
         ELSE 
            A(I) = A(I-1) + A(I-2)
         ENDIF
      ENDDO
      END
C END FILE FIB2.F

che compilo con

f2-0

f2py riempie il terminale di scritte, roba da preoccuparsi se uno non sapesse che è normale. Ecco lo script Python

#!/usr/bin/python
# -*- coding: utf-8 -*-

#modificato, uso interi

import numpy as np
import fib2 #quello in Fortran

a = np.zeros(10, 'int')
fib2.fib(a)
print a

f2-1

OK. Vediamo qualche caso tipico, non quelli presenti nei tutorial, ma quelli che so che possono capitare, per esempio un “Hello World!” poco funzionale, chissà se funziona.

1
2
3
4
5
      subroutine hw
c scrive il classico saluto

      print *, "Hello World!"
      end

Ho dovuto usare la sintassi vecchia, quella fixed-form, in attesa di leggere il manuale. Ma un passo per volta.

#!/usr/bin/python
# -*- coding: utf-8 -*-

import hw
print "Adesso chiamo il Fortran"
hw.hw()
print "Fatto, OK?"

hw

OK. La pessima abitudine dei programmatori di una volta era di scrivere su file quando capita (visto che mi sto funzionalizzando?), vediamo:

1
2
3
4
5
6
7
8
      subroutine hw2
c scrive il classico saluto

      print *, "Hello World!"
      open(1, file = "prova-f2py.txt")
      write(1, '("Hello World!")')
      close(1)
      end

Per sicurezza ho controllato che con un main in Fortran funzionasse, non si sa mai… È OK.

#!/usr/bin/python
# -*- coding: utf-8 -*-

import hw2
print "Adesso chiamo il Fortran"
hw2.hw2()
print "Fatto, OK?"

hw2

OK, continua a funzionare. E se uno ha una subroutine che ne chiama un’altra in un altro file cosa succede?

1
2
3
4
      subroutine hw3

      call hw2  
      end

La hw2 è la stessa del caso precedente. Il comando per f2py è:

f2py -c hw3.f hw2.f -m hw3

e poi:

#!/usr/bin/python
# -*- coding: utf-8 -*-

import hw3
print "Adesso chiamo il Fortran"
hw2.hw2()
print "Fatto, OK?"

hw3

Chissà se posso passare una stringa al Fortran?

1
2
3
4
5
6
      integer function fs(st)
      character*(*) st

      print *, st
      fs = len(st)
      end
#!/usr/bin/python
# -*- coding: utf-8 -*-

import fs

msg = "f2py rockz!"
print "Adesso chiamo il Fortran"
res = fs.fs(msg)
print res
print "Fatto, OK?"

fs

OK! Un’ultima cosa, poi mi leggo tutto il manuale ma per intanto il free-form, non so se non si possa semplificare, questo è quello che sono riuscito a ottenere finora:

1
2
3
4
5
6
7
8
9
module fhw
  integer i
contains
  subroutine hw
    ! scrive il classico saluto

    print *, "Hello World!"
  end subroutine hw
end module fhw
#!/usr/bin/python
# -*- coding: utf-8 -*-

import hw
print "Adesso chiamo il Fortran"
hw.hw()
print "Fatto, OK?"

hw

No! è solo questione di mettere l’estensione .f90 e funziona tutto:

1
2
3
4
5
subroutine hw
! scrive il classico saluto

print *, "Hello World!"
end subroutine hw
#!/usr/bin/python
# -*- coding: utf-8 -*-

import hw
print "Adesso chiamo il Fortran"
hw.hw()
print "Fatto, OK?"

L’output è uguale al precedente.

Quasi dimenticavo:
NumPy è opera di Travis E. Oliphant, con il cotributo di numerosi collaboratori, vedi Origins of NumPy in Guide to NumPy.
f2py è stato scritto da Pearu Peterson.

Bellissima (anche se l’ho solo scorsa per adesso) la Guide to Numpy, se si ha fretta c’è anche An introduction to Numpy and Scipy. E poi c’è il Cookbook, il sito di SciPy dove ci sono ovviamente NumPy e f2py.

Quindi, conclusione provvisoria: si può usare? Forse, dipende. Leggo il manuale e vi dico.
Prossimamente… forse…😉