Tecnicamente's Weblog

Just another WordPress.com weblog

Archivi delle etichette: C

Fortran95, il ciclo “DO” su matrici di grandi dimensioni

Se state usando il Fortran95 i motivi sono semplici, o vi stanno obbligando oppure volete usare un linguaggio di programmazione che non ha eguali nelle applicazioni del calcolo numerico ed in particolar modo nei calcoli dove si fa uso pesante di matrici.

Questo breve articolo ha due scopi:

  1. mostrare che bisogna organizzare le matrici nel modo migliore per permettere ai cicli di lettura di essere molto efficienti;
  2. mostrare le ottimizzazioni intrinseche del linguaggio e del compilatore;

Supponiamo di avere due matrici di grande dimensione la [A] e la [B], il calcolo che si vuole eseguire è c(i,j) = a(i,j) + b(i,j). Voglio trovare la matrice somma [C]=[A]+[B].

Operazione estremamente semplice, supponiamo ancora che decidiamo di usare il linguaggio C e che i dati siano interi, vediamo cosa accade usando il seguente codice:

#include<stdio.h>
#include<time.h>

#define X 10000
#define Y 10000

int main(void){
 int i, j;
 float t0, t1;
 //int a[X][Y], b[X][Y], c[X][Y];
 int **a, **b, **c;

/* allocazione della memoria */
 a = (int **)calloc(X,sizeof(int *));
 b = (int **)calloc(X,sizeof(int *));
 c = (int **)calloc(X,sizeof(int *));

 for(i=0;i<X;i++){
    a[i] = (int *)calloc(Y,sizeof(int));
    b[i] = (int *)calloc(Y,sizeof(int));
    c[i] = (int *)calloc(Y,sizeof(int));
 }

/* inizializzazione delle matrici */
 for(i=0;i<X;i++){
    for(j=0;j<Y;j++){
       a[i][j] = 1;
       b[i][j] = 2;
    }
 }
 
/* somma scorrendo per le colonne */
 t0 = clock();
 for(i=0;i<X;i++){
   for(j=0;j<Y;j++){
     c[i][j] = a[i][j] + b[i][j];
   }
 }
 
 t1 = clock();
 printf("%d\n",c[X-1][Y-1]);
 printf("tempo richiesto somme su colonne = %f\n",(double)(t1 - t0) / CLOCKS_PER_SEC);

/* somma scorrendo le righe */
 t0 = clock();
 for(j=0;j<Y;j++){
   for(i=0;i<X;i++){
     c[i][j] = a[i][j] + b[i][j];
   }
 }
 t1 = clock();
 printf("%d\n",c[X-1][Y-1]);
 printf("tempo richiesto somme su colonne = %f\n",(double)(t1 - t0) / CLOCKS_PER_SEC);

 return 0;
}

Il cui risultato sul mio povero PC è rappresentato nella seguente immagine:

C

Ok ma perchè ho messo un codice in C ed un test del risultato?

Ovviamente per mostrare come in Fortran95 sia comodo descrivere lo stesso algoritmo e per vedere quali sono i risultati in termini di tempo di accesso, quindi qui il codice:

program prova
implicit none
  integer :: i,j
  integer, parameter :: x=10000, y=10000
  real :: t0,t1
  integer, dimension(x,y) :: a,b,c

  a=1 ; b = 2

!SOMMA SCORRENDO LE COLONNE
  call cpu_time(t0)
  do i=1,x
    do j=1,y
      c(i,j) = a(i,j) + b(i,j)
    enddo
  enddo
  call cpu_time(t1)
  print*, c(x,y)
  print*, "tempo richiesto somme su colonne = ",t1-t0
  print*,""

!SOMMA SCORRENDO LE RIGHE
  do j=1,y
    do i=1,x
      c(i,j) = a(i,j) + b(i,j)
    enddo
  enddo
  call cpu_time(t1)
  print*, c(x,y)
  print*, "tempo richiesto somme su righe = ", t1-t0
  print*,""

!SOMMA IN MODO NATIVO DEL FORTRAN
  call cpu_time(t0)
  c=a+b
  call cpu_time(t1)
  print*, c(x,y)
  print*, "tempo richiesto somme tra matrici = ", t1-t0

end program prova

Con il seguente risultato:F

Come vedete la principale differenza tra i due linguaggi è che in Fortran le MATRICI ed i VETTORI sono tipi nativi, quindi è possibile fare delle operazioni direttamente su di esse.

Occorre tenere presente che ho usato i seguenti compilatori gcc e gfortran, quindi il backend è lo stesso, il perchè della differenza di risultato lo lascio agli esperti, io non so spiegarlo.

E’ interessante notare come una implementazione non oculata degli accessi agli array in C può essere disastrosa.

Tutti i modi usati in Fortran sono migliori, in particolare l’uso delle operazioni a livello “nativo” sulle matrici permette di ridurre notevolmente le righe di codice oltre al non trascurabile fatto avere tempi di calcolo migliori che in C.

Una interessante possibilità dei compilatori fortran è quella di parallelizzare il codice, cosa che purtroppo non riesco a testare :(.

 

Buon Coding in Fortran.

Annunci

Importare file dati in Excel

E’ indubbio che Excel come tutti i fogli di calcolo sia ormai uno strumento indispensabile nell’analisi dei dati, e sicuramente tutti vi siete trovati ad affrontare il triste problema della “virgola”.

Infatti nelle impostazioni Italiane (ma forse anche straniere), i formati decimali sono separati da “,” e non da “.” .

Questo può creare dei problemi quando si devono importare file numerici realizzati da altri software o scaricati da data-logger, proprio come era successo a me anni fa.

Il problema di per se è facilmente risolvibile quando il file da importare contiene solo numeri decimali con separazione basata su “.”, in questi casi si apre il file con un qualsiasi editor e si esegue la sostituzione automatica di “.” con “,” in tal modo il gioco è fatto e si importa direttamente in Excel un file leggibile secondo il formato che esso si attende.

Tuttavia esistono casi in cui l’operazione precedente non è possibile, in particolare quando il file da importare è mostruosamente grande, ed in tal caso occore un editor di testi serio, o quando il file contiene dati promiscui, come ad esempio il file che dovevo importare io:

000000    25/10/2004    17.07.35    21.487    1.612    26.90    27.17
000001    25/10/2004    17.17.38    21.223    1.596    26.05    27.04
000002    25/10/2004    17.27.35    20.820    1.574    25.08    25.89
000003    25/10/2004    17.37.35    19.464    1.515    23.74    24.90
000004    25/10/2004    17.47.35    16.722    1.426    22.43    23.74
000005    25/10/2004    17.57.35    11.868    1.233    21.33    22.81
000006    25/10/2004    18.07.35    6.576    0.809    20.38    21.99
000007    25/10/2004    18.17.35    2.836    0.382    19.77    21.67
000008    25/10/2004    18.27.35    0.746    0.110    19.30    21.56
000009    25/10/2004    18.37.36    0.090    0.013    19.26    21.62
000010    25/10/2004    18.47.36    0.010    0.001    18.95    21.56
000011    25/10/2004    18.57.36    0.003    0.001    18.92    21.63
000012    25/10/2004    19.07.36    0.003    0.001    19.25    21.62
000013    25/10/2004    19.17.36    0.002    0.001    19.11    21.45

Nel mio caso oltre ad essere un file piuttosto grosso,  (mesi di acquisizione dati), mi serviva modificare i separatori decimali solo dalla quarta colonna in poi, per cui il metodo pocanzi menzionato non era applicabile.

Per risolvere il problema sono bastate poche righe di codice C, dove imponevo di mutare i “.” in “,” dal secondo punto letto in poi per ogni riga.

Quindi ho prodotto il seguente codice, dove con DOT_LIMIT pari a 2 richiedo la modifica a partire dal terzo “.” letto, e DOT_UP_LIMIT pari a 6 per azzerare il contatore a fine linea.

Riposto il codice, per quanto banale  spero possa essere utile a qualcuno.

#define DOT_LOW_LIMIT 2
#define DOT_UP_LIMIT 6
#include <stdio.h>
#include <stdlib.h>

int main(int argc, char *argv[])
{
  FILE *fpin;
  FILE *fpout;

  int dot_count;
  char carattere;

  if (argc !=3 )
  {
           fprintf(stderr,"Uso: %s <file_in> <file_out>\n",argv[0]);
           exit(0);
  }

  fpin = fopen(argv[1],"r");
  fpout = fopen(argv[2],"w");

  dot_count = 0;

  while(fscanf(fpin,"%c",&carattere) != EOF ){
           if (carattere == EOF) break;
           if (carattere == '.'){
                         dot_count++;
                         if(dot_count > DOT_LOW_LIMIT){
                                     fprintf(fpout,"%c",',');
                                     fprintf(stdout,"%c",',');
                         }
                         else{
                              fprintf(fpout,"%c",carattere);
                              fprintf(stdout,"%c",carattere);
                         }
                         if(dot_count == DOT_UP_LIMIT) dot_count = 0;
           }
           else{
                fprintf(fpout,"%c",carattere);
                fprintf(stdout,"%c",carattere);
           }
  }

  fclose(fpin);
  fclose(fpout);

  return 0;
}

Ottenendo questo risultato:

 

000000    25/10/2004    17.07.35    21,487    1,612    26,90    27,17
000001    25/10/2004    17.17.38    21,223    1,596    26,05    27,04
000002    25/10/2004    17.27.35    20,820    1,574    25,08    25,89
000003    25/10/2004    17.37.35    19,464    1,515    23,74    24,90
000004    25/10/2004    17.47.35    16,722    1,426    22,43    23,74
000005    25/10/2004    17.57.35    11,868    1,233    21,33    22,81
000006    25/10/2004    18.07.35    6,576    0,809    20,38    21,99
000007    25/10/2004    18.17.35    2,836    0,382    19,77    21,67
000008    25/10/2004    18.27.35    0,746    0,110    19,30    21,56
000009    25/10/2004    18.37.36    0,090    0,013    19,26    21,62
000010    25/10/2004    18.47.36    0,010    0,001    18,95    21,56
000011    25/10/2004    18.57.36    0,003    0,001    18,92    21,63
000012    25/10/2004    19.07.36    0,003    0,001    19,25    21,62
000013    25/10/2004    19.17.36    0,002    0,001    19,11    21,45

Fortran parte 3: Fortran Vs. C/C++ o compilatore Vs. compilatore?

Basta fare una ricerca dove si vogliono recuperare informazioni circa le performances offerte dal Fortran95 contro C/C++ e si trovano centianai di benchmark a favore dell’uno e dell’altro linguaggio.

Se poi si va su gruppi di discussione si assiste a vere e proprie guerre di religione per partito preso e che quasi sempre si concludono con insulti.

Non è mia intenzione schierarmi dall’una o dall’altra parte, credo semplicemente che i software, come qualsiasi altra cosa, sono dei prodotti che devono essere ingegnerizzati.

Quindi quando c’è da metter giù due righe di codice si può partire alla “garibaldina”, ma quando le cose sono serie ed impegnative occorre spendere del tempo per svolgere delle analisi e successivamente impostare l’architettura del software e la scelta degli strumenti da usare.

Tra gli strumenti c’è sicuramente il linguaggio di programmazione, difficilmente un solo linguaggio è adatto a coprire tutti gli scopi.

Attenzione, con un linguaggio come ad esempio il C o il C++ si può fare praticamente tutto, non a caso sono “general purpose“, e devo essere sincero, sono un grande ammiratore di questi linguaggi e delle loro potenzialità.

Ma in un ambiente di produzione il fattore tempo è fondamentale, il ricorso a strumenti che facilitano certe parti del progetto è, a mio parere, obbligatorio.

Il Fortran non si pone come un linguaggio per risolvere qualsiasi problema, il Fortran si pone come un linguaggio per la soluzione di problemi numerici.

In tale ottica offre una serie di strumenti che altri linguaggi non hanno, per capirci dovendo spiegare quale è il linguaggio più somigliante al fortran, io dire Matlab.

La sostanziale differenza è che Fortran necessita di compilatori, Matlab mi pare sia in pseudocodice.

Ma torniamo alla questione Fortran Vs. C/C++, voglio solo fare degli esempi per dare una idea di quello che mi piace del Fortran ma anche del C/C++.

Si prenda in considerazione l’esempio seguente, una matrice B contiene dei valori numerici, ed una seconda matrice A è nulla.

Si vuole che A assuma gli stessi elementi di B quando questi sono pari.

Il problema proposto è solo uno spunto per mostrare le capacità del linguaggio.

In C si può scrivere il seguente programma:

#define X 15000
#define Y 15000
#define ESECUZIONE 5

#include<stdio.h>
#include<stdlib.h>
#include<malloc.h>

int main(void)
{
  int i,j,k;
  int **a;
  int **b;

  for(k=0;k<ESECUZIONE;k++)
  {
    a = malloc(sizeof(int)*X);
    b = malloc(sizeof(int)*X);

    for(i=0;i<Y;i++)
    {
      a[i] = malloc(sizeof(int)*Y);
      b[i] = malloc(sizeof(int)*Y);

      for(j=0;j<Y;j++)
      {
        b[i][j] = j+1;
      }
    }

    for(i=0;i<X;i++)
    {
      for(j=0;j<Y;j++)
      {
        if((b[i][j]%2)==0)
        {
          a[i][j] = b[i][j];
        }
      }
    }
    printf("%d\n",a[10-1][10-1]);

    for(i=0;i<X;i++)
    {
      free(a[i]);
      free(b[i]);
    }
    free(a);
    free(b);
    return 0;
}

return 0;
}

Come si potrebbe scrivere lo stesso programma in Fortran95?
Ecco la risposta:

program esempio
implicit none

integer, parameter :: X = 15000
integer, parameter :: Y = 15000
integer, parameter :: ESECUZIONE = 5

integer :: i,j,k

integer, allocatable :: A(:,:), B(:,:)

do k=1,ESECUZIONE
 allocate(A(X,Y))
 allocate(B(X,Y))

 A = 0

 do j=1,Y
     B(:,j) = (/(i,i=1,X)/)
 enddo

 where(mod(B,2)==0) A=B

 write(*,*) B(10,10)

 deallocate (A)
 deallocate (B)

enddo
end program esempio

o più sinteticamente:

program esempio_2
implicit none

integer, parameter :: X = 15000
integer, parameter :: Y = 15000
integer, parameter :: ESECUZIONE = 5

integer :: i,j,k

integer, allocatable :: A(:,:), B(:,:)

do k=1,ESECUZIONE
	allocate(A(X,Y))
	allocate(B(X,Y))

	A = 0

	forall(j=1:Y)
		B(:,j)=j
	end forall

	where(mod(B,2)==0) A=B

	write(*,*) B(10,10)

        deallocate (A)
	deallocate (B)
enddo
end program esempio_2

La cosa evidente, a prescindere da fatto che abbiate capito qualcosa nel codice Fortran :), è che in Fortran le linee di codice richieste sono effettivamente poche, e vi assicuro che se si conosce il linguaggio sono anche estremamente leggibili.

Vi ho mostrato due versioni in Fortran95, la differenza tra i due sta nel fatto che nel primo codice uso un ciclo DO…ENDDO per riempire la matrice B lungo le colonne, ciò è necessario per velocizzare gli accessi in memoria visto che il Fortran organizza i dati in memoria seguendo l’ordine delle collone, a differenza del C/C++ che seguno le righe.

Nel secondo codice il ciclo DO…ENDDO precedende è sostituito dal un ciclo FORALL…END FORALL, tale ciclo è molto particolare, e dovrebbe dare il meglio di se su architetture multiprocessore.

Infatti FORALL implementa il calcolo parallelo nativamente e direttamente dal compilatore, cosa non da poco, purtroppo io dispongo di un sistema non proprio performante e per di più monoprocessore infatti eccovi la scheda del mio sistema:

Adesso passiamo alla copilazione ed al “benchmark” dei codici.

Per la precisione userò una versione per usi non commerciali del compilatore Fortran di Intel e l’immancabile gcc per il codice in C.

Compilo i primi due codici ed ottengo:

Quale è la conclusione?

La velocità dipende:

– dal tempo impiegato a scrivere tutte le righe di codice

– da come il programmatore implementa il codice

– dal tempo impiegato a debuggare

– dalle ottimizzazioni realizzate dal compilatore

Circa l’ultimo punto ho una domanda, come mai gfortran che poi è gcc impiega più tempo rispetto alla versione in C, mentre il compilatoreo ifort impiega praticamente lo stesso tempo di gcc?

Ovviamente questo test non è un test valido per qualificare le prestazioni dei compilatori.

Se qualcuno prova ed eseguire i test su una macchina dual core mi farebbe piacere avere i risultati.

In tal caso il codice fortran da usare è il secondo dove è presente l’istruzione FORALL.

Appena recupero i sorgenti della mia tesi, vedo di eseguire un test sui diversi compilatori, a suo tempo impiegava 60 ore, compilatore Compaq Visual Fortran e un bel PC da 3200 Mhz.