Mise en oeuvre d'OpenMP dans le logiciel MEFISTO

Alain PERRONNET Laboratoire Jacques-Louis Lions UPMC Paris


BUTS de cette mise en oeuvre:
1. Récupérer plus rapidement les résultats d'un calcul lourd de plusieurs heures de calcul;

2. Exploiter au mieux les potentialités de mon PC:

   . la mémoire partagée de 8 Giga octects;

   . les 8 potentiels threads (processus légers) éxecutables en parallèle (i-core 7 d'intel);

3. Ne pas avoir à modifier l'algorithmique des programmes éléments finis du logiciel,
   notamment, ne pas avoir à coder une méthode de sous-domaines pour laquelle
   toute la structure du programme est à revoir;

4. En cas de succès, utiliser la machine HPC du LJLL avec ses 160 possibilités de threads
   en parallèle pour obtenir encore plus vite les résultats;

5. Voir le côut d'implémentation d'OpenMP et ses difficultés de mise en oeuvre;

6. Tirer une conclusion sur la nécessité ou non d'utiliser un second langage MPI, ... , pour
   améliorer les temps de restitution d'un calcul "éléments finis".

POURQUOI OpenMP?
- Disponible et gratuit sur mon PC (option -openmp aux compilateurs GNU);

- Utilise la mémoire partagée du PC;

- Langage peu compliqué à mettre en oeuvre par ajoût de directives dans le fortran;

- L'application "surveillance du système" visualise en continu l'emploi des threads en parallèle;

- Permet de tout faire sur le PC: frappe, exécution, analyse du résultat.

Programme test pour l'ajoût d'OpenMP:

Un résultat numérique, a été obtenu par une exécution séquentielle avec 81178 secondes d'attente:

183 600 tétraèdres Taylor-Hood, 34 562 sommets; 260 484 noeuds; 816 014 degrés de liberté Vitesse Pression;
ρ=1 Re=300; dt=0.05 T=120; => 2400 pas de temps calculés;
Résolution par gradient conjugué sur des matrices condensées demandant 3 738 810 double precision et entiers;
Le module de la vitesse est ici représenté dans le plan Y=0 par des couleurs au temps 120.
Exécution séquentielle en 81178 secondes:


Exécution avec 8 threads en 26676 secondes:



L'algorithme programmé:
Navier-Stokes Option6: Méthode des caractéristiques et pas fractionnaires:

A L'ETAPE 0 IL FAUT RECUPERER U0 et CONSTRUIRE F0

A L'ETAPE n+1, au temps tn+1=t0 + (n+1)  δ, le problème consiste a trouver
 {u(tn+1), p(tn+1)} SOLUTION de


(8) ( ρ - dt μ Δ ) {ui*(tn+1,x)} =
      ρ {ui(tn,X(tn;tn+1,x))} -dt CoGrPr gradi p(tn,x) + dt {fgi(tn+1)})

      où X(tn;tn+1,x) est la position de la molécule du fluide a tn
      qui sera en x à tn+1 (méthode des caractéristiques)

      Condition aux limites
       Ui(tn+1,X) = Ui(tn+1,X)      avec X sur Γ 
      -dt μ dUi(tn+1,γ)/dn = 0   avec γ sur Frontiere-Γ 
                                                               pour  i=1,...,NDIM


(10)  -dt CoGrPr Δ {p(tn+1)-p(tn)}=-(ρ-dt μ Δ)Div{u*(tn+1)}
      Condition aux limites
       Pn+1,m+1(γ) = P(tn+1,γ) pour Γ sur la frontiere
      -dt CoGrPr dPn+1,m+1/dn =-n . (ρ-dt μ Δ){u*(tn+1)}


(12) ( ρ - dt μ Δ ) {u(tn+1)-u*} = - dt CoGrPr GRAD(p(tn+1)-p(tn))
     Condition aux limites
      Ui(tn+1,X)-u* = Ui(tn+1,X)-u*   avec X sur Γ 
     -dt μ dUi(tn+1,γ)/dn = 0      avec γ sur Frontiere-Γ 
                                                                pour  i=1,...,NDIM

     p(tn+1,Sommets) = p(tn,Sommets) + ( p(tn+1,Sommets)-p(tn,Sommets) )
 
     U(tn+1,Noeuds) = U(tn,Noeuds) + ( U(tn+1,Noeuds)-U*(tn,Noeuds) )

Ce qui conduit, préalablement aux itérations sur le temps,
à assembler une fois, à partir des matrices calculées sur chacun des nbe éléments finis du maillage:

- une matrice VG "éléments finis" en vitesse, d'ordre le nombre de noeuds vitesse nbnv
  (nombre de sommets et milieux des arêtes des tétraèdres),
  stockée sous forme condensée dans 3 tableaux ( i, j, VGij )

- une matrice PG "éléments finis" en pression, d'ordre le nombre de sommets pression nbsp,
  stockée sous forme condensée dans 3 tableaux ( i, j, PGij )

Ensuite, à chaque pas de temps,
à assembler les vecteurs seconds membres, à partir des vecteurs calculés
sur chacun des nbe éléments finis du maillage de:

- (8)   3 systèmes linéaires "éléments finis", de nombre de lignes nbnv;

- (10)  1 système linéaire "éléments finis", de nombre de lignes nbsp;

- (12)  3 systèmes linéaires "éléments finis", de nombre de lignes nbnv;

- (13)  1 combinaison linéaire de 2 vecteurs de nbnv composantes;

- (14)  1 combinaison linéaire de 2 vecteurs de nbsp composantes.

D'où la modification des sous-programmes les plus utilisés et coûteux:
- Résolution du gradient conjugué avec matrice condensée (sp gcaxb.f)
=> . fonction produit scalaire de 2 vecteurs;     (sp proscd.f => proscd.f95 lib util)

DOUBLE PRECISION FUNCTION PROSCD( V1, V2, N )
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! BUT : PRODUIT SCALAIRE DE 2 VECTEURS V1 V2 REELS DOUBLE PRECISION
! ----- DE N COMPOSANTES AVEC NBTHREAD THREADS
!       V1 et V2 SONT DANS LA MEMOIRE PARTAGEE A TOUS LES THREADS
! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! AUTEURS: Philippe PARNAUDEAU  LJLL UPMC Paris             Octobre 2012
!          Alain    PERRONNET   LJLL UPMC Paris & St Pierre du Perray
!.......................................................................
!$ USE OMP_LIB
   IMPLICIT  NONE
   INTEGER           N, K
   DOUBLE PRECISION  V1(N), V2(N), S
!$ INTEGER           NBTHREADS

!///////////////////////////////////////////////////////////////////////
!!$OMP PARALLEL 
!!$OMP MASTER
!      print *,'PROSCDomp: Nombre THREADS=',OMP_GET_NUM_THREADS()
!!$OMP END MASTER
!!$OMP END PARALLEL
!///////////////////////////////////////////////////////////////////////

      PROSCD = 0.d0
!///////////////////////////////////////////////////////////////////////
!$    NBTHREADS = OMP_GET_NUM_THREADS()

!$OMP PARALLEL PRIVATE(K,S) SHARED(N,V1,V2)
!$OMP DO SCHEDULE( STATIC, N/NBTHREADS ) REDUCTION(+:PROSCD)

      DO K = 1, N
         S = V1(K) * V2(K)
         PROSCD = PROSCD + S
      ENDDO

!$OMP END DO
!$OMP END PARALLEL
!///////////////////////////////////////////////////////////////////////

!     print *,'PROSCD: (V1,V2)=',PROSCD,' pour',N,' composantes'

      RETURN
END FUNCTION PROSCD


   . fonction combinaison linéaire de 2 vecteurs; (sp cl2ved.f => cl2ved.f95 lib util)

   . sous programme produit matrice vecteur;      (sp magcvx.f => magcvx.f95 lib reso)

SUBROUTINE MAGCVX( NTDL, NPDLFX, LPLIGN, LPCOLO, AG, X,  Y )
! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! BUT:     Y=AG*X avec Y,X VECTEURS (NTDL)  AG MATRICE MORSE SYMETRIQUE
! ----     VARIANTE DE magcve.f AVEC CHAQUE LIGNE D'UN DEGRE DE LIBERTE
!          CONSIDEREE COMME LA LIGNE DE LA MATRICE UNITE MAIS STOCKEE
!          SANS MODIFICATION POUR PRENDRE EN COMPTE LE STOCKAGE
!          SYMETRIQUE DE LA MATRICE
!
!
! ENTREES:
! --------
! NTDL   : NOMBRE DE LIGNES DE LA MATRICE ET DU VECTEUR
! NPDLFX : NPDLFX(N)=0 SI LE DL N EST LIBRE, >0 SI LE DL N EST FIXE
!          NO TEMOIN D'UN DEGRE DE LIBERTE FIXE'
! LPLIGN : LES POINTEURS SUR LES LIGNES DE LA MATRICE MORSE AG
! LPCOLO : LES NUMEROS DES COLONNES DES COEFFICIENTS DE LA MATRICE MORSE
! AG     : LES COEFFICIENTS DE LA MATRICE MORSE AG
! X      : VECTEUR DE NTDL COMPOSANTES
!
! SORTIE ou MODIFIE:
! ------------------
! Y      : VECTEUR DE NTDL COMPOSANTES
!          ATTENTION: Y DOIT ETRE DIFFERENT DE X A L'APPEL
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! AUTEUR : Alain PERRONNET LJLL UPMC & St Pierre du Perray  Octobre 2012
!23456---------------------------------------------------------------012
!$ USE OMP_LIB
   IMPLICIT  NONE
   INTEGER           NTDL
   DOUBLE PRECISION  AG(1:*), X(NTDL), Y(NTDL), AJ, XI, YI, P
   INTEGER           NPDLFX(NTDL), LPLIGN(0:NTDL), LPCOLO(1:*)
   INTEGER           I, J, NODIAG0, NODIAG, NOCOL
!$ INTEGER           NBTHREADS

!///////////////////////////////////////////////////////////////////////
!$OMP PARALLEL PRIVATE( I, J, NODIAG0, NODIAG, NOCOL, NBTHREADS, AJ, XI, YI, P ) &
!$OMP          SHARED ( NTDL, NPDLFX, LPLIGN, LPCOLO, AG, X,  Y )
!$    NBTHREADS = OMP_GET_NUM_THREADS()

!$OMP DO SCHEDULE( STATIC, NTDL/NBTHREADS )
!     INITIALISATION DE TOUS LES Y(I) OBLIGATOIRE AVANT SOMMATION
      DO I = 1, NTDL
         IF( NPDLFX(I) .EQ. 0 ) THEN
!           LIGNE D'UN DEGRE DE LIBERTE LIBRE
            Y(I) = 0D0
         ELSE
!           LA LIGNE I DU DL FIXE EST EN FAIT LA LIGNE IDENTITE
            Y(I) = X(I)
         ENDIF
      ENDDO
!$OMP END DO

!$OMP DO SCHEDULE( STATIC, NTDL/NBTHREADS )
      DO I = 1, NTDL

!        NUMERO +1 DU COEFFICIENT DIAGONAL I-1
         NODIAG0 = LPLIGN(I-1) + 1

!        NUMERO DU COEFFICIENT DIAGONAL I
         NODIAG = LPLIGN(I)

!        I-EME COEFFICIENT DE X
         XI = X(I)

!        I-EME COEFFICIENT DE Y
         YI = Y(I)

         IF( NPDLFX(I) .EQ. 0 ) THEN

!           LIGNE D'UN DEGRE DE LIBERTE LIBRE
!           PARCOURS DES COLONNES AVANT LA DIAGONALE DE LA LIGNE I DE AG
            DO J = NODIAG0, NODIAG-1

!              NO DE LA COLONNE DU COEFFICIENT J
               NOCOL = LPCOLO(J)

!              CONTRIBUTION DU COEFFICIENT AG(J) AVANT LA DIAGONALE
               AJ = AG(J)
               P  = AJ * X(NOCOL)
               YI = YI + P

!              CONTRIBUTION DU COEFFICIENT APRES LA DIAGONALE PAR SYMETRIE
               IF( NPDLFX(NOCOL) .EQ. 0 ) THEN
!                 COLONNE ET LIGNE D'UN DEGRE DE LIBERTE LIBRE
!$OMP ATOMIC
                  Y(NOCOL) = Y(NOCOL) + AJ * XI
               ENDIF
!
            ENDDO

!           CONTRIBUTION DU COEFFICIENT DIAGONAL SUPPOSE 1 SI DL FIXE
            P  = AG(NODIAG) * XI
            YI = YI + P

         ELSE

!           LIGNE D'UN DEGRE DE LIBERTE FIXE
            DO J = NODIAG0, NODIAG-1

!              NO DE LA COLONNE DU COEFFICIENT J
               NOCOL = LPCOLO(J)

               IF( NPDLFX(NOCOL) .EQ. 0 ) THEN
!                 COLONNE D'UN DEGRE DE LIBERTE LIBRE
!$OMP ATOMIC
                  Y(NOCOL) = Y(NOCOL) + AG(J) * XI
                ENDIF

            ENDDO

!           LA LIGNE I DU DL FIXE EST EN FAIT LA LIGNE IDENTITE
            YI = XI

         ENDIF

!        RENDEZ VOUS APRES SOMME DES COEFFICIENTS SOUS-DIAGONAUX DANS YI
         Y(I) = YI

      ENDDO
!$OMP END DO
!$OMP END PARALLEL
!///////////////////////////////////////////////////////////////////////

      RETURN
END SUBROUTINE MAGCVX


- Construction de la matrice PG  (sp p13daggc.f => p13daggc.f95 lib flui et asmegc.f lib reso)
- Construction de la matrice VG  (sp agvitth.f  => agvitth.f95  lib flui)

SUBROUTINE AGVITTH( Rho,    DtMhu,  XYZNOE, &
                    NUTYEL, NBNOEF, NBELEM, NUNOEF, &
                    NORESO, NBLVG,  NBCVG,  LPLIGN, LPCOLO, VG )
! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! BUT:     CONSTRUIRE ET ASSEMBLER LA MATRICE GLOBALE DU PROBLEME
! ----     VG=( Rho - dt Mhu Laplacien ) = Integrale Rho P2 P2
!                                        + deltat Mhu grad P2 grad P2 dX
!          avec INTEGRATION EXACTE SUR L'ELEMENT FINI
!
! ENTREES:
! --------
! Rho    : COEFFICIENT DE LA MATRICE DE MASSE      P2  P2
! DtMhu  : COEFFICIENT DE LA MATRICE DE VISCOSITE DP2 DP2
! XYZNOE : XYZNOE(3,NBLVG)  3 COORDONNEES DES NOEUDS DU MAILLAGE
! NUTYEL : NUMERO DU TYPE D'EF ( 15 TRIANGLE, 20 TETRAEDRE TAYLOR-HOOD )
! NBNOEF : NOMBRE DE NOEUDS D'UN EF DE CE TYPE
!          ( 6 POUR LE TRIANGLE, 10 POUR LE TETRAEDRE)
! NBELEM : NOMBRE D'EF DU MAILLAGE
! NUNOEF : NUNOEF(NBELEM,NBNOEF) NO DES 6 NOEUDS DES NBELEM EF

! NORESO : CODE RESOLUTION DES SYSTEMES LINEAIRES AVEC LA MATRICE VG
!          1 FACTORISATION COMPLETE DE CHOLESKY ET MATRICE PROFIL VG
!          2 GRADIENT CONJUGUE AVEC UN STOCKAGE MORSE DE VG
! NBLVG  : NOMBRE DE NOEUDS DU MAILLAGE (SOMMETS ET MILIEUX DES ARETES)
! NBCVG  : NOMBRE DE COEFFICIENTS DE VG
!          (SI NORESO=1   NBCVG=0
!          (SI NORESO=2   NBCVG=NOMBRE DE COLONNES STOCKEES)
! LPLIGN : POINTEURS SUR LES COEFFICIENTS DIAGONAUX DE VG
! LPCOLO : NUMERO DE COLONNE DES COEFFICIENTS DE LA MATRICE MORSE VG
!
! SORTIE :
! --------
! VG     : MATRICE GLOBALE ASSEMBLEE
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! AUTEUR : ALAIN PERRONNET LJLL UPMC & St Pierre du Perray  Octobre 2012
!23456---------------------------------------------------------------012
!$ USE OMP_LIB
      IMPLICIT NONE
      include"./incl/donflu.inc95"
      include"./incl/a___npef.inc95"

      INTEGER           NORESO, NBLVG, NBCVG, LPLIGN(0:NBLVG), &
                        LPCOLO(NBCVG), NUTYEL, NBNOEF, NBELEM, &
                        NUNOEF(NBELEM,NBNOEF)
      DOUBLE PRECISION  VG(NBCVG), Rho, DtMhu
      REAL              XYZNOE(3,NBLVG)

      INTEGER           I, J, M, NUELEM, NONOEF(10)
      DOUBLE PRECISION  AE(55), S
!$    INTEGER           NBTHREADS

!     MISE A ZERO DE LA MATRICE VG
      CALL AZEROD( NBCVG, VG )

!///////////////////////////////////////////////////////////////////////
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE( NUELEM, I, J, M, NONOEF, AE, S )
!$    NBTHREADS = OMP_GET_NUM_THREADS()

!     BOUCLE SUR LES EF de TAYLOR-HOOD
!$OMP DO SCHEDULE( STATIC, NBELEM/NBTHREADS )
      DO NUELEM = 1, NBELEM

!        NO DES NOEUDS DE L'ELEMENT FINI NUELEM
         DO M = 1, NBNOEF
            NONOEF(M) = NUNOEF(NUELEM,M)
         ENDDO

         IF( NUTYEL .EQ. 15 ) THEN

!           CALCUL DE LA MATRICE DE VISCOSITE Ae DE LA VITESSE P2
            CALL FAE2P2P2( Rho, DtMhu, NONOEF, NBLVG, XYZNOE, AE )

         ELSE

!           CALCUL DE LA MATRICE DE VISCOSITE Ae DE LA VITESSE P2
            CALL FAE3P2P2( Rho, DtMhu, NONOEF, NBLVG, XYZNOE, AE )

         ENDIF

!        ASSEMBLAGE DE LA MATRICE ELEMENTAIRE DANS LA MATRICE GLOBALE
         IF( NORESO .EQ. 1 ) THEN
!           MATRICE SYMETRIQUE PROFIL
            CALL ASMEPC( NBNOEF, NONOEF, 1, AE, AE, 1, &
                         LPLIGN, VG )
         ELSE
!           MATRICE SYMETRIQUE MORSE
!           LES MATRICES ELEMENTAIRE ET GLOBALE SONT ICI SYMETRIQUES
            CALL ASMEGC( NBNOEF, NONOEF, 1, AE, AE, 1, &
                         LPLIGN, LPCOLO, VG )
         ENDIF

      ENDDO     ! FIN DE LA BOUCLE DES EF POUR ASSEMBLER VG
!$OMP END DO
!$OMP END PARALLEL
!///////////////////////////////////////////////////////////////////////

!      call affvect( 'AGVITTH: MATRICE VITESSE VG', 20,    VG )
!      call afl1ve(  'AGVITTH: MATRICE VITESSE VG', NBCVG, VG )

      RETURN
END SUBROUTINE AGVITTH


      SUBROUTINE ASMEGC( NBDL,   NODL,   NCODAE, AES, AE,
     %                   NCODAG, LPLIGN, LPCOLO, AG )
C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C BUT :   ASSEMBLER LA MATRICE ELEMENTAIRE EN LA MATRICE GLOBALE
C -----   MORSE SYMETRIQUE OU NON TOUTE EN MEMOIRE CENTRALE
C         (DIFFERENCE AVEC LE SP ASABGC PAS DE VECTEUR ELEMENTAIRE)
C
C ENTREES:
C --------
C NBDL   : NOMBRE DE   DEGRES DE LIBERTE DE L ELEMENT FINI
C NODL   : NO DES NBDL DEGRES DE LIBERTE DE L ELEMENT FINI
C NCODAE : CODE DE STOCKAGE DE LA MATRICE ELEMENTAIRE
C           0 : MATRICE DIAGONALE
C          -1 : MATRICE NON SYMETRIQUE
C           1 : MATRICE SYMETRIQUE
C AES    : MATRICE ELEMENTAIRE SYMETRIQUE (NBDL * (NBDL+1) / 2 )
C AE     : MATRICE ELEMENTAIRE  (NBDL,NBDL)
C NCODAG : CODE DE STOCKAGE DE LA MATRICE GLOBALE MORSE
C           0 : MATRICE DIAGONALE
C          -1 : MATRICE NON SYMETRIQUE
C           1 : MATRICE SYMETRIQUE
C LPLIGN : TABLEAU POINTEUR SUR LES COEFFICIENTS DIAGONAUX
C LPCOLO : TABLEAU DES NUMEROS DE COLONNES
C
C SORTIES:
C --------
C AG     : MATRICE MORSE TOUTE EN MEMOIRE CENTRALE
C ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C AUTEUR :   A. PERRONNET LABORATOIRE ANALYSE NUMERIQUE PARIS   MAI 1989
C MODIF OMP: A. PERRONNET LJLL UPMC & St Pierre du Perray  Decembre 2012
C ......................................................................
C$    use OMP_LIB
      INTEGER           NODL(1:NBDL), LPLIGN(*), LPCOLO(1:*)
      DOUBLE PRECISION  AES(1:*), AE(NBDL,NBDL), AG(1:*)
C
      IF( NCODAE .EQ. 0 .AND. NCODAG .EQ. 0 ) THEN
C
C        LES MATRICES ELEMENTAIRE ET GLOBALE SONT DIAGONALES
C        ===================================================
         DO 10 I=1,NBDL
C           LE NO GLOBAL DU I-EME DEGRE DE LIBERTE LOCAL
            IEG = NODL( I )
C           ASSEMBLAGE DE AE( IEG )
C$OMP ATOMIC
            AG( IEG ) = AG( IEG ) + AES( I )
 10      CONTINUE
C
      ELSE IF( NCODAE .GT. 0 .AND. NCODAG .GT. 0 ) THEN
C
C        LES MATRICES ELEMENTAIRE ET GLOBALE SONT SYMETRIQUES
C        ====================================================
         L = 0
         DO 40 I=1,NBDL
C
C           LE NO GLOBAL DU I-EME DEGRE DE LIBERTE LOCAL
            IEG = NODL( I )
C
C           ASSEMBLAGE DE AES DANS AG
C           LA BOUCLE SUR LES COLONNES DE LA MATRICE ELEMENTAIRE
            DO 30 J=1,I
               JEG = NODL( J )
C              ADRESSES DANS LA MATRICE MORSE
               IF(JEG.LE.IEG) THEN
                  K1  = LPLIGN(IEG) + 1
                  K2  = LPLIGN(IEG+1)
                  KEG = JEG
               ELSE
                  K1  = LPLIGN(JEG) + 1
                  K2  = LPLIGN(JEG+1)
                  KEG = IEG
               END IF
C              RECHERCHE DU NUMERO DE COLONNE JEG PAR DICHOTOMIE
               KI=K1
               KS=K2
 21            KM=(KI+KS)/2
               IF(LPCOLO(KM)-KEG) 22,23,24
 22            KI=KM
               GO TO 25
 24            KS=KM
 25            IF(KI+1-KS)        21,26,30
 26            IF(LPCOLO(KI)-KEG) 27,28,30
 28            KM=KI
               GO TO 23
 27            IF(LPCOLO(KS)-KEG) 30,29,30
 29            KM=KS
C
 23            IAJEG = KM
               L     = L + 1
C$OMP ATOMIC
               AG( IAJEG ) = AG( IAJEG ) + AES( L )
 30         CONTINUE
 40      CONTINUE
         RETURN
C
      ELSE IF( NCODAE .LT. 0 .AND. NCODAG .LT. 0 ) THEN
C
C        LES MATRICES ELEMENTAIRE ET GLOBALE SONT NON SYMETRIQUES
C        ========================================================
         DO 140 I=1,NBDL
C           LE NO GLOBAL DU I-EME DEGRE DE LIBERTE LOCAL
            IEG = NODL( I )
C           ADRESSES DANS LA MATRICE MORSE
            K1 = LPLIGN(IEG) + 1
            K2 = LPLIGN(IEG+1)
C
C           LA BOUCLE SUR LES COLONNES DE LA MATRICE ELEMENTAIRE
            DO 130 J=1,NBDL
               JEG = NODL( J )
C              RECHERCHE DU NUMERO DE COLONNE JEG PAR DICHOTOMIE
               KI=K1
               KS=K2
 121           KM=(KI+KS)/2
               IF(LPCOLO(KM)-JEG) 122,123,124
 122           KI=KM
               GO TO 125
 124           KS=KM
 125           IF(KI+1-KS)        121,126,130
 126           IF(LPCOLO(KI)-JEG) 127,128,130
 128           KM=KI
               GO TO 123
 127           IF(LPCOLO(KS)-JEG) 130,129,130
 129           KM=KS
C
 123           IAJEG = KM
C              ASSEMBLAGE DE AE(I,J)
C$OMP ATOMIC
               AG( IAJEG ) = AG( IAJEG ) + AE( I , J )
 130        CONTINUE
 140     CONTINUE
C
      ELSE IF( NCODAE .EQ. 0 .AND. NCODAG .NE. 0 ) THEN
C
C        LA MATRICE ELEMENTAIRE EST DIAGONALE MAIS
C        LA MATRICE GLOBALE EST SYMETRIQUE OU NON SYMETRIQUE
C        ===================================================
         DO 170 I=1,NBDL
C           LE NO GLOBAL DU I-EME DEGRE DE LIBERTE LOCAL
            IEG = NODL( I )
C           SA POSITION SUR LA DIAGONALE DE LA MATRICE GLOBALE MORSE
            IEG = LPLIGN( IEG + 1 )
C           ASSEMBLAGE DE AES( IEG )
C$OMP ATOMIC
            AG( IEG ) = AG( IEG ) + AES( I )
 170     CONTINUE
C
      ENDIF
      RETURN
      END


- Construction des 3 seconds membres du système en vitesse  (sp bgvitth.f    => bgvitth.f95   )
- Construction du    second  membre  du système en pression (sp bgdivlath3.f => bgdivlath3.f95)
- Construction des 3 seconds membres du système en vitesse  (sp bggrath.f    => bggrath.f95   )

- Maximum de la vitesse aux noeuds du maillage              (sp maxvit.f    => maxvit.f95   )
- Flux normal de la vitesse sur les surfaces frontalières   (sp fluxvitsf.f => fluxvitsf.f95)
- Presentation de Vx,Vy,Vz,Pression en 1 tableau aux noeuds (sp vxyzpre.f   => vxyzpre.f95  )

Liste des problèmes rencontrés:
L'instruction EQUIVALENCE entre variables et un tableau semble poser problème.
=> Suppression de EQUIVALENCE

. Utiliser sur le PC la totalité des threads, fait qu'il n'est plus possible de
travailler simultanément... ou à retarder l'éxécution en cours.
Exécuter le programme avec 2 threads de moins que le total disponible semble le bon choix.

. Obtention de RÉSULTATS DIFFÉRENTS pour 2 passages différents ou un nombre différent de threads demandés
<= Accès d'au moins 2 threads à la même mémoire, modification de la valeur
   dans chaque thread et écriture sur la mémoire, la dernière écriture
   effaçant toutes les précédentes réécritures.
=> ajouter la directive ATOMIC  (cf asmegc.f);
   bien déclarer privées les variables prenant des valeurs différentes selon le thread d'exécution.


. Exemple de précision de l'assemblage d'un vecteur dans un vecteur global:
  (ou d'une matrice élémentaire dans une matrice globale)

At Time   0.500000E-01 : VELOCITY PRESSURE computation ===============================================
0 thread:
BGVITTH VECTEUR BG
  V(   1)=  0.70214E-03  V(   2)=  0.10883E-02  V(   3)=  0.21064E-02  V(   4)=  0.18134E-02  V(   5)=  0.14237E-02
  V(   6)=  0.14043E-02  V(   7)=  0.14237E-02  V(   8)=  0.18137E-02  V(   9)=  0.21064E-02  V(  10)=  0.21356E-02
  V(  11)=  0.18134E-02  V(  12)=  0.14043E-02  V(  13)=  0.14237E-02  V(  14)=  0.18137E-02  V(  15)=  0.21064E-02
  V(  16)=  0.21356E-02  V(  17)=  0.18134E-02  V(  18)=  0.14043E-02  V(  19)=  0.14237E-02  V(  20)=  0.18137E-02
 BGVITTH VECTEUR BG Min  Vect(n) = -4.96056672457384495E-002  ncMIN=      206126
 BGVITTH VECTEUR BG Mean|Vect(n)|=  3.75167242676738245E-003
 BGVITTH VECTEUR BG Max  Vect(n) =  4.96056672457384357E-002  ncMAX=      206038
 gcaxb: ||B0|| =   29864.296659327549     
 gcaxb: ||X0|| =   259249.00000000000     
 gcaxb: ||AV0||=   29864.297149371512     
 gcaxb: ||R0|| =  1.80881864330702718E-004
 CG Iteration   0  ||Ax-b||**2=  1.80881864330702718E-004  EpsZero=  9.99999975E-05  NTDL=      260484
 CG Iteration  12  from ||Ax-b||**2=  1.80881864330702718E-004  to   1.55254977087381506E-008  => CONVERGENCE with Eps=  9.99999975E-05


2 threads:
BGVITTH VECTEUR BG
  V(   1)=  0.7021434E-03  V(   2)=  0.1088266E-02  V(   3)=  0.2106430E-02  V(   4)=  0.1813437E-02  V(   5)=  0.1423708E-02
  V(   6)=  0.1404287E-02  V(   7)=  0.1423708E-02  V(   8)=  0.1813692E-02  V(   9)=  0.2106430E-02  V(  10)=  0.2135562E-02
  V(  11)=  0.1813437E-02  V(  12)=  0.1404287E-02  V(  13)=  0.1423708E-02  V(  14)=  0.1813692E-02  V(  15)=  0.2106430E-02
  V(  16)=  0.2135562E-02  V(  17)=  0.1813437E-02  V(  18)=  0.1404287E-02  V(  19)=  0.1423708E-02  V(  20)=  0.1813692E-02
 BGVITTH VECTEUR BG Min  Vect(n) = -4.96056672457384426E-002  ncMIN=      206134
 BGVITTH VECTEUR BG Mean|Vect(n)|=  3.75167 403274361162E-003
 BGVITTH VECTEUR BG Max  Vect(n) =  4.96056672457384357E-002  ncMAX=      206038
 gcaxb: ||B0|| =   29864.296659864296     
 gcaxb: ||X0|| =   259249.00000000000     
 gcaxb: ||AV0||=   29864.297149371512     
 gcaxb: ||R0|| =  1.79903504976940392E-004
 CG Iteration           0  ||Ax-b||**2=  1.79903504976940392E-004  EpsZero=  9.99999975E-05  NTDL=      260484
 CG Iteration          12  from ||Ax-b||**2=  1.79903504976940392E-004  to   1.53240086394572182E-008  => CONVERGENCE with Eps=  9.99999975E-05
 GCAXB: Solution Xl= Min  Vect(n) =   0.0000000000000000       ncMIN=        3789
 GCAXB: Solution Xl= Mean|Vect(n)|=  0.99354119774549166     
 GCAXB: Solution Xl= Max  Vect(n) =   1.1648122716996268       ncMAX=       15670
 

4 threads:
BGVITTH VECTEUR BG
  V(   1)=  0.70214E-03  V(   2)=  0.10883E-02  V(   3)=  0.21064E-02  V(   4)=  0.18134E-02  V(   5)=  0.14237E-02
  V(   6)=  0.14043E-02  V(   7)=  0.14237E-02  V(   8)=  0.18137E-02  V(   9)=  0.21064E-02  V(  10)=  0.21356E-02
  V(  11)=  0.18134E-02  V(  12)=  0.14043E-02  V(  13)=  0.14237E-02  V(  14)=  0.18137E-02  V(  15)=  0.21064E-02
  V(  16)=  0.21356E-02  V(  17)=  0.18134E-02  V(  18)=  0.14043E-02  V(  19)=  0.14237E-02  V(  20)=  0.18137E-02
 BGVITTH VECTEUR BG Min  Vect(n) = -4.960566724573844 26E-002  ncMIN=      206134
 BGVITTH VECTEUR BG Mean|Vect(n)|=  3.75167 854757440687E-003
 BGVITTH VECTEUR BG Max  Vect(n) =  4.96056672457384357 E-002  ncMAX=      206038
 gcaxb: ||B0|| =   29864.296660974403     
 gcaxb: ||X0|| =   259249.00000000000     
 gcaxb: ||AV0||=   29864.297149371512     
 gcaxb: ||R0|| =  1.75908139914939373E-004
 CG Iteration   0  ||Ax-b||**2=  1.75908139914939373E-004  EpsZero=  9.99999975E-05  NTDL=      260484
 CG Iteration  12  from ||Ax-b||**2=  1.75908139914939373E-004  to   1.5 0616849379181280E-008  => CONVERGENCE with Eps=  9.99999975E-05


8 threads:
 BGVITTH VECTEUR BG
  V(   1)=  0.70214E-03  V(   2)=  0.10883E-02  V(   3)=  0.21064E-02  V(   4)=  0.18134E-02  V(   5)=  0.14237E-02
  V(   6)=  0.14043E-02  V(   7)=  0.14237E-02  V(   8)=  0.18137E-02  V(   9)=  0.21064E-02  V(  10)=  0.21356E-02
  V(  11)=  0.18134E-02  V(  12)=  0.14043E-02  V(  13)=  0.14237E-02  V(  14)=  0.18137E-02  V(  15)=  0.21064E-02
  V(  16)=  0.21356E-02  V(  17)=  0.18134E-02  V(  18)=  0.14043E-02  V(  19)=  0.14237E-02  V(  20)=  0.18137E-02
 BGVITTH VECTEUR BG Min  Vect(n) = -4.96056672457384426E-002  ncMIN=      206134
 BGVITTH VECTEUR BG Mean|Vect(n)|=  3.7516 8264737622001E-003
 BGVITTH VECTEUR BG Max  Vect(n) =  4.96056672457384357E-002  ncMAX=      206038
 gcaxb: ||B0|| =   29864.296665405171     
 gcaxb: ||X0|| =   259249.00000000000     
 gcaxb: ||AV0||=   29864.297149371512     
 gcaxb: ||R0|| =  1.74607224328496849E-004
 CG Iteration   0  ||Ax-b||**2=  1.74607224328496849E-004  EpsZero=  9.99999975E-05  NTDL=      260484
 CG Iteration  12  from ||Ax-b||**2=  1.74607224328496849E-004  to   1.52030438117300345E-008  => CONVERGENCE with Eps=  9.99999975E-05



Calcul de la PRESSION  au temps 0.500000E-01: Un peu plus loin dans le calcul --------------------------
0 thread:
 Sd Membre BGDIVLATH
  V(   1)= -0.39522E-07  V(   2)= -0.10829E-05  V(   3)= -0.13978E-06  V(   4)= -0.28370E-09  V(   5)= -0.25759E-05
  V(   6)= -0.19790E-07  V(   7)= -0.52695E-06  V(   8)= -0.33372E-05  V(   9)= -0.35791E-07  V(  10)= -0.75587E-06
  V(  11)= -0.36139E-05  V(  12)= -0.47480E-07  V(  13)= -0.86679E-06  V(  14)= -0.36815E-05  V(  15)= -0.52893E-07
  V(  16)= -0.90329E-06  V(  17)= -0.36727E-05  V(  18)= -0.54083E-07  V(  19)= -0.90897E-06  V(  20)= -0.40244E-08
 Sd Membre BGDIVLATH Min  Vect(n) = -5.13728986339159338E-002  ncMIN=        8755
 Sd Membre BGDIVLATH Mean|Vect(n)|=  4.14183930871127131E-004
 Sd Membre BGDIVLATH Max  Vect(n) =  7.03476976908405521E-002  ncMAX=        1294
 gcaxb: ||B0|| =  0.52557651214980372     
 gcaxb: ||X0|| =   0.0000000000000000     
 gcaxb: ||AV0||=   0.0000000000000000     
 gcaxb: ||R0|| =  0.52557651214980372     
 CG Iteration    0  ||Ax-b||**2=  0.52557651214980372       EpsZero=  9.99999975E-05  NTDL=       34562
 CG Iteration  385  from ||Ax-b||**2=  0.52557651214980372       to   5.01251797321817188E-005  => CONVERGENCE with Eps=  9.99999975E-05


2 threads:
Sd Membre BGDIVLATH
  V(   1)= -0.3977577E-07  V(   2)= -0.1100148E-05  V(   3)= -0.1433028E-06  V(   4)= -0.2982607E-09  V(   5)= -0.2619714E-05
  V(   6)= -0.2128133E-07  V(   7)= -0.5440158E-06  V(   8)= -0.3385362E-05  V(   9)= -0.3923357E-07  V(  10)= -0.7801152E-06
  V(  11)= -0.3641960E-05  V(  12)= -0.5116731E-07  V(  13)= -0.8851787E-06  V(  14)= -0.3685920E-05  V(  15)= -0.5503834E-07
  V(  16)= -0.9101525E-06  V(  17)= -0.3672567E-05  V(  18)= -0.5466859E-07  V(  19)= -0.9087995E-06  V(  20)= -0.4013950E-08
 Sd Membre BGDIVLATH Min  Vect(n) = -5.13729538537972125E-002  ncMIN=        8755
 Sd Membre BGDIVLATH Mean|Vect(n)|=  4.14028744100109749E-004
 Sd Membre BGDIVLATH Max  Vect(n) =  7.03213023014770261E-002  ncMAX=        1294
 gcaxb: ||B0|| =  0.52478994868193551     
 gcaxb: ||X0|| =   0.0000000000000000     
 gcaxb: ||AV0||=   0.0000000000000000     
 gcaxb: ||R0|| =  0.52478994868193551     
 CG Iteration           0  ||Ax-b||**2=  0.52478994868193551       EpsZero=  9.99999975E-05  NTDL=       34562
 CG Iteration         386  from ||Ax-b||**2=  0.52478994868193551       to   4.85809208086402056E-005  => CONVERGENCE with Eps=  9.99999975E-05


4 threads:
Sd Membre BGDIVLATH
  V(   1)= -0.37228E-07  V(   2)= -0.97176E-06  V(   3)= -0.13132E-06  V(   4)= -0.25782E-09  V(   5)= -0.22712E-05
  V(   6)= -0.17670E-07  V(   7)= -0.48075E-06  V(   8)= -0.28745E-05  V(   9)= -0.29693E-07  V(  10)= -0.64752E-06
  V(  11)= -0.30471E-05  V(  12)= -0.36600E-07  V(  13)= -0.69278E-06  V(  14)= -0.30054E-05  V(  15)= -0.39864E-07
  V(  16)= -0.69819E-06  V(  17)= -0.28803E-05  V(  18)= -0.42557E-07  V(  19)= -0.70597E-06  V(  20)=  0.40211E-08
 Sd Membre BGDIVLATH Min  Vect(n) = -5.14745954505346839E-002  ncMIN=        8755
 Sd Membre BGDIVLATH Mean|Vect(n)|=  4.13627939534453562E-004
 Sd Membre BGDIVLATH Max  Vect(n) =  7.02700759595715069E-002  ncMAX=        1294
 gcaxb: ||B0|| =  0.52196590004469734     
 gcaxb: ||X0|| =   0.0000000000000000     
 gcaxb: ||AV0||=   0.0000000000000000     
 gcaxb: ||R0|| =  0.52196590004469734     
 CG Iteration    0  ||Ax-b||**2=  0.52196590004469734       EpsZero=  9.99999975E-05  NTDL=       34562
 CG Iteration  387  from ||Ax-b||**2=  0.52196590004469734       to   4.25113919730824033E-005  => CONVERGENCE with Eps=  9.99999975E-05


8 threads:
Sd Membre BGDIVLATH
  V(   1)= -0.39124E-07  V(   2)= -0.10691E-05  V(   3)= -0.14780E-06  V(   4)= -0.36750E-09  V(   5)= -0.25785E-05
  V(   6)= -0.25082E-07  V(   7)= -0.56258E-06  V(   8)= -0.33645E-05  V(   9)= -0.45016E-07  V(  10)= -0.79682E-06
  V(  11)= -0.36052E-05  V(  12)= -0.54816E-07  V(  13)= -0.89516E-06  V(  14)= -0.36018E-05  V(  15)= -0.56036E-07
  V(  16)= -0.91764E-06  V(  17)= -0.35529E-05  V(  18)= -0.53368E-07  V(  19)= -0.90093E-06  V(  20)= -0.43771E-08
 Sd Membre BGDIVLATH Min  Vect(n) = -5.13300870491007616E-002  ncMIN=        8755
 Sd Membre BGDIVLATH Mean|Vect(n)|=  4.13554653304912994E-004
 Sd Membre BGDIVLATH Max  Vect(n) =  7.03503944847874846E-002  ncMAX=        1294
 gcaxb: ||B0|| =  0.52250333576920138     
 gcaxb: ||X0|| =   0.0000000000000000     
 gcaxb: ||AV0||=   0.0000000000000000     
 gcaxb: ||R0|| =  0.52250333576920138     
 CG Iteration    0  ||Ax-b||**2=  0.52250333576920138       EpsZero=  9.99999975E-05  NTDL=       34562
 CG Iteration  384  from ||Ax-b||**2=  0.52250333576920138       to   5.06110138140815387E-005  => CONVERGENCE with Eps=  9.99999975E-05


Moralité, faire plusieurs éxécutions avec un nombre différent de threads et
vérifier que les résultats sont "identiques" à un nombre suffisant de chiffres près.


Modes de calcul des temps de l'exécution:
Au début du sous programme fluinth6.f ou fluinbf6.f, déclaration des variables:
      DOUBLE PRECISION  DATE00, DATE
      integer           t0,  t1, nbclockps
      DOUBLE PRECISION  t_cpu_0, t_cpu_1
C$    DOUBLE PRECISION  tt0, tt1


Au début du module fluinth6.f ou fluinbf6.f, initialisation des temps:
CALL SECONDES1970( DATE00 ) 
<=> C: *DATE00 = (double) time( NULL )             => la seconde du départ (0h 1 juin 1970)

call system_clock (count=t0, count_rate=nbclockps) => temps du départ

call cpu_time(t_cpu_0)                             => temps CPU de départ

C$OMP PARALLEL SHARED( tt0 )
C$    tt0 = OMP_GET_WTIME()                        => temps du départ
C$OMP END PARALLEL


A la fin du module fluinth6.f ou fluinbf6.f, les affichages:
CALL SECONDES1970( DATE )
print *,'fluinth6: Execution seconds=',DATE-DATE00     => fluinth6: Execution seconds=1449.

call system_clock( count=t1, count_rate=nbclockps )
PRINT *,'fluinth6 elapsed real execution seconds=',
   real( t1-t0, kind=8 ) / real( nbclockps, kind=8 )       => fluinth6 elapsed real execution seconds=   1448.625000 

      call cpu_time(t_cpu_1)
      PRINT *,'fluinth6: cpu_time=',t_cpu_1-t_cpu_0         => fluinth6: cpu_time= 11076.28 secondes (cumulées des threads)

C$OMP PARALLEL SHARED( tt0, tt1 )
C$    tt1 = OMP_GET_WTIME()
C$    print *,'OMP seconds=', tt1-tt0
C$OMP END PARALLEL                                         => OMP seconds=   1448.6192778999994     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005     
                                                              OMP seconds=   1448.6193055740005

La commande time $MEFISTO/pp/ppflui                        => real    25m4.324s   =  1504 secondes
                                                              user    184m20.899s = 11061 secondes cumulées des threads
                                                              sys     0m16.709s

Temps CPU: CALL TEMPSCPU(*tclock) = ((double) clock()) /
 ( (double) CLOCKS_PER_SEC )                               => TOTAL COMPUTATION TIME= 11076.28 CPU SECONDS cumulées des threads


Ci-dessous: L'impression à chaque pas de temps de Last Command Real Time=SECONDES est le resultat dans le programme principal des instructions:
 30   CALL  SECONDES1970( DATE )
      SECONDES = DATE - DATE0
      DATE0    = DATE
      WRITE(IMPRIM,20030) SECONDES
20030 FORMAT(/'Last Command Real Time=',G15.6,' seconds')

Résultats pour 0.05s le pas de temps, 183 600 tétraèdres TAYLOR-HOOD, 34 562 sommets, 260 484 noeuds et 816 014 d.l. Vitesse(P2)-Pression(P1):
- Exécution séquentielle   AUCUNE directive  OpenMP: 34s/Pas de Temps de 0.05s     Last Command Real Time=  81178 s      CONVERGENCE des GC au temps=120 s
  27/11/2012
- Exécution parallèle 1 thread  + directives OpenMP: 37.16s/Pas de Temps de 0.05s  Last Command Real Time=  89212 s      CONVERGENCE des GC au temps=120 s
  27/11/2012                                       => Déccélération 0.91 

- Exécution parallèle 2 threads + directives OpenMP: 25s/Pas de Temps de 0.05s     Last Command Real Time=  72596 s      CONVERGENCE des GC au temps=120 s
  27/11/2012                                       => Accélération  1.36                                    => Accélération  1.12

- Exécution parallèle 4 threads + directives OpenMP: 17s/Pas de Temps de 0.05s     Last Command Real Time=  35962 s      NON CONVERGENCE du  GC au temps= 114 s
  1/12/2012                                        => Accélération  2        35962 / 114 * 120 =  37855 s   => Accélération  2.14

- Exécution parallèle 4 threads + directives OpenMP: 15.35s/Pas de Temps de 0.05s  Last Command Real Time=  36882 s      CONVERGENCE des GC au temps= 120 s
  5/12/2012  avec integrale moyenne volumique=0    => Accélération  2.21                                    => Accélération  2.20

- Exécution parallèle 8 threads + directives OpenMP: 13.88s/Pas de Temps de 0.05s  Last Command Real Time=  23019 s      NON CONVERGENCE du  GC au temps= 82.8 s
  2/12/2012                                        => Accélération  2.45     23019 /  82 * 120 =  33686 s   => Accélération  2.41

- Exécution parallèle 8 threads + directives OpenMP: 11.10s/Pas de Temps de 0.05s  Last Command Real Time=  26676 s      CONVERGENCE du  GC au temps= 120 s
  11/12/2012 avec moyenne=0 et vitesse convectée   => Accélération  3.06                                    => Accélération  3.04
  Time=120.003  VelocityMean=0.999268 VelocityMax=1.51191   PressureMean=-34.0433 PressMax-Min=3.28687   Flux-=-84.0000 Flux+=84.5853 Flux-+ =0.585287

- Exécution parallèle 6 threads + directives OpenMP: 12.53s/Pas de Temps de 0.05s  Last Command Real Time= 15079s        CONVERGENCE du  GC au temps= 60 s
  13/12/2012 avec moyenne=0 et vitesse convectée   => Accélération  2.71                                    => Accélération 2.69
  Time=59.9994  VelocityMean=0.993302 VelocityMax=1.46847   PressureMean=-0.349340E-01 PressMax-Min=3.11558   Flux-=-84.0000 Flux+=84.5316 Flux-+=0.531648


- Exécution séquentielle   AUCUNE directive  OpenMP: 26.12s/Pas de Temps de 0.05s  Last Command Real Time=31375s        CONVERGENCE des GC au temps=60 s
  19/12/2012 avec moyenne=0 et vitesse convectée
  Time=59.9994  VelocityMean=0.997427 VelocityMax=1.51769   PressureMean=-0.283292E-01 PressMax-Min=3.19100   Flux-=-84.0000  Flux+=84.9504  Flux-+=0.950398

- Exécution parallèle 6 threads + directives OpenMP: 11.38s/Pas de Temps de 0.05s  Last Command Real Time=13697s        CONVERGENCE du  GC au temps= 60 s
  21/12/2012 avec moyenne=0 vitesse convectée projetée => Accélération  2.30                                => Accélération 2.29
  Time=59.9994  VelocityMean=0.997125  VelocityMax=1.63501  PressureMean=-0.265324E-01 PressMax-Min=3.13565   Flux-=-84.0000  Flux+=   84.8085 Flux-+ =  0.808514


- Exécution sequentielle 1 thread + directives OpenMP: 26.53s/Pas de Temps de 0.1s  Last Command Real Time=33481s
  02/02/2013 avec moyenne=0 vitesse convectée projetée CHOLESKY L tL
  Time=119.999  VelocityMean=0.994483  VelocityMax=1.42976  PressureMean=-0.339109E-01 PressMax-Min=3.23492   Flux-=-84.0000  Flux+=84.5586  Flux-+=0.558603

- Exécution parallèle 6 threads + directives OpenMP: 20.99s/Pas de Temps de 0.1s  Last Command Real Time=26857s
  01/02/2013 avec moyenne=0 vitesse convectée projetée CHOLESKY L tL => Accélération 1.26
  Time=119.999  VelocityMean=0.994428  VelocityMax=1.42610  PressureMean=-0.293486E-01 PressMax-Min=3.28135   Flux-=-84.0000  Flux+=84.2240  Flux-+=0.223994    

Résultats pour 0.05s le pas de temps, 260 040 tétraèdres BREZZI-FORTIN, 48 785 sommets, 48785+260484 noeuds et 975 260 d.l. Vitesse Pression:
- Exécution séquentielle   AUCUNE directive  OpenMP: 20.36s/Pas de Temps de 0.05s Last Command Real Time= 24444s        CONVERGENCE des GC au temps=60 s
  20/12/2012 option 6: moyenne=0 et vitesse convectée Gradient conjugué
  Time=59.9994  VelocityMean=0.988612  VelocityMax=2.11612    PressureMean= -0.440098E-01 PressMax-Min=11.8017   Flux-=-84.2143  Flux+=84.3352  Flux-+=0.120867    

- Exécution parallèle 6 threads + directives OpenMP: 9.17s/Pas de Temps de 0.05s  Last Command Real Time= 11027 s       CONVERGENCE du GC au dernier temps final= 60 s
  13/12/2012 option 6: moyenne=0 et vitesse convectée Gradient conjugué   => Accélération 2.22
  Time=59.9994  VelocityMean=0.989035  VelocityMax=2.10014    PressureMean=-0.451065E-01  PressMax-Min=12.0009  Flux-=-84.2143  Flux+=84.2649  Flux-+=0.506429E-01
 
- Exécution parallèle 8 threads + directives OpenMP: 8.15s/Pas de Temps de 0.05s  Last Command Real Time= 19558 s       CONVERGENCE du  GC au temps= 120 s
  12/12/2012 option 6: moyenne=0 et vitesse convectée Gradient conjugué   => Accélération 2.50
  Time=120.003  VelocityMean=1.00563   VelocityMax=2.11311    PressureMean=-0.173629E-01  PressMax-Min=11.7314  Flux-=-84.2143  Flux+=89.5580  Flux-+=5.34367

- Exécution parallèle 6 threads + directives OpenMP: 7.93s/Pas de Temps de 0.05s  Last Command Real Time= 9540 s        CONVERGENCE du GC au dernier temps final= 60 s
  21/12/2012 option 6: moyenne=0 vitesse convectée projeté Gradient conjuguée  => Accélération 2.57
  Time=59.9994  VelocityMean=0.988307  VelocityMax=2.11318    PressureMean=-0.417128E-01 PressMax-Min=11.7768   Flux-=-84.2143  Flux+=84.3717  Flux-+=0.157407

Résultats pour différents pas de temps, 183 600 tétraèdres TAYLOR-HOOD, 34 562 sommets, 260 484 noeuds et 816 014 d.l. Vitesse(P2)-Pression(P1):
Exécution séquentielle   AUCUNE directive OpenMP: 12.6s/Pas de Temps de 0.05s   Last Command Real Time=30285s
Time=120.003   VelocityMean=0.997868  VelocityMax=1.45806    PressureMean=-0.328888E-01  PressMax-Min=3.41286     Flux-=-84.0000  Flux+=84.1472  Flux-+=0.147205    
19/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué

Exécution parallèle 6 threads + directives OpenMP: 7.83s/Pas de Temps de 0.05s   Last Command Real Time=18846s  => Accélération 1.61 
Time=120.003   VelocityMean=0.999435  VelocityMax=1.51210    PressureMean=-0.344601E-01  PressMax-Min=3.35377     Flux-=-84.0000  Flux+=84.3044  Flux-+=0.304383    
16/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué

Exécution parallèle 8 threads + directives OpenMP: 6.91s/Pas de Temps de 0.05s   Last Command Real Time=16644s  => Accélération 1.82 
Time=120.003   VelocityMean=0.997824  VelocityMax=1.48398    PressureMean=-0.242325E-01  PressMax-Min=3.29233     Flux-=-84.0000  Flux+=84.2119  Flux-+=0.211924    
15/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué

Résultats pour différents pas de temps, 183 600 tétraèdres BREZZI-FORTIN, 34 562 sommets, 34562+183600 noeuds et 689048 d.l. Vitesse(P1B)-Pression(P1):
Exécution séquentielle   AUCUNE directive OpenMP: 7.32s/Pas de Temps de 0.05s   Last Command Real Time=17586s
Time=120.003   VelocityMean=0.993988  VelocityMax=1.68325    PressureMean=-0.259710E-01  PressMax-Min=6.29447     Flux-=-84.4598  Flux+=84.4894  Flux-+= 0.295768E-01
19/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué

Exécution parallèle 6 threads + directives OpenMP: 2.98s/Pas de Temps de 0.05s   Last Command Real Time=7161s  => Accélération 2.46
Time=120.003   VelocityMean=0.994055  VelocityMax=1.72503    PressureMean=-0.308770E-01  PressMax-Min=6.08664     Flux-=-84.4598  Flux+=84.0525  Flux-+=-0.407390    
15/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué

Exécution parallèle 8 threads + directives OpenMP: 3.03s/Pas de Temps de 0.05s   Last Command Real Time=7282s  => Accélération 2.42
Time=120.003   VelocityMean=0.994084  VelocityMax=1.65286    PressureMean=-0.256841E-01  PressMax-Min=6.35291     Flux-=-84.4598  Flux+=84.6701  Flux-+= 0.210231    
16/01/2013  option 6: moyenne=0 vitesse convectée projetée Gradient conjugué


Exécution séquentielle 1 thread + directives OpenMP: 8.38s/Pas de Temps de 0.05s  Last Command Real Time=20150s
Time=120.003   VelocityMean=0.995933  VelocityMax=1.72782    PressureMean=-0.283742E-01  PressMax-Min=6.20095     Flux-=-84.4598  Flux+=84.6238  Flux-+= 0.163970    
04/02/2013  option 6: moyenne=0 vitesse convectée projetée Facorisation de Cholesky

Exécution parallèle 6 threads + directives OpenMP: 3.44s/Pas de Temps de 0.05s   Last Command Real Time=8299s  => Accélération 2.43
Time=120.003   VelocityMean=0.994844  VelocityMax=1.69037    PressureMean=-0.295341E-01  PressMax-Min=6.33560     Flux-=-84.4598  Flux+=84.3652  Flux-+=-0.946078E-01
03/02/2013  option 6: moyenne=0 vitesse convectée projetée Facorisation de Cholesky

Résultats pour différents pas de temps, 12470 triangles TAYLOR-HOOD, 6359 sommets, 25188 noeuds et 56735 d.l. Vitesse(P2)+Pression(P1):
Exécution séquentielle   AUCUNE directive OpenMP: 0.273s/Pas de temps de 0.05s   Last Command Real Time=657s
Time=120.003   VelocityMean=0.997915  VelocityMax=1.49592    PressureMean=-0.512409E-01  PressMax-Min=1.49956     Flux-=-12.0000   Flux+=11.9962  Flux-+=-0.380065E-02
14/01/2013   option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée   Gradient conjugué

Exécution séquentielle 1 thread + directives OpenMP: 0.388/Pas de temps de 0.05s  Last Command Real Time=932s  => Décélération 0.705
Time=120.003   VelocityMean=1.00076   VelocityMax=1.48156    PressureMean=-0.651925E-01  PressMax-Min=1.38706     Flux-=-12.0000   Flux+=12.0086  Flux-+= 0.858030E-02
04/02/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée   Gradient conjugue

Exécution parallèle 6 threads + directives OpenMP: 0.221s/Pas de Temps de 0.05s   Last Command Real Time=532s  => Accélération 1.24
Time=120.003   VelocityMean=0.998086  VelocityMax=1.48343    PressureMean=-0.373372E-01  PressMax-Min=1.39900     Flux-=-12.0000  Flux+=11.9650  Flux-+=-0.349863E-01
Time=120.003   VelocityMean=0.991301  VelocityMax=1.50448    PressureMean=-0.452973E-01  PressMax-Min=1.49634     Flux-=-12.0000  Flux+=11.7144  Flux-+=-0.285577  second calcul
14/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée   Gradient conjugue

Exécution parallèle 8 threads + directives OpenMP: 0.203s/Pas de Temps de 0.05s   Last Command Real Time=491s  => Accélération 1.34
Time=120.003   VelocityMean=0.999784  VelocityMax=1.47131    PressureMean=-0.857504E-01  PressMax-Min=1.33769     Flux-=-12.0000  Flux+=11.9990  Flux-+=-0.102494E-02    
Time=120.003   VelocityMean=0.991367  VelocityMax=1.49490    PressureMean= 0.362189E-02  PressMax-Min=1.49861     Flux-=-12.0000  Flux+=11.7034  Flux-+=-0.296629  second calcul
19/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée   Gradient conjugue


Exécution séquentielle 1 thread + directives OpenMP: 0.278/Pas de temps de 0.05s  Last Command Real Time=670s
Time=120.003   VelocityMean=0.978581  VelocityMax=1.61972    PressureMean=-0.100306      PressMax-Min=1.12007     Flux-=-12.0000   Flux+=12.4576  Flux-+ =0.457566  
04/02/2013  option 6:  pas fractionnaire, moyenne=0, vitesse convectée projetée  Factorisation de Cholesky   Battement derriere le cylindre à partir de t=100s

Exécution parallèle 6 threads + directives OpenMP: 0.206/Pas de temps de 0.05s    Last Command Real Time=498s  => Accélération 1.35
Time=120.003   VelocityMean=1.00036   VelocityMax=1.50442    PressureMean=-0.656686E-01  PressMax-Min=1.43816     Flux-=-12.0000  Flux+=11.9323  Flux-+=-0.676643E-01  
31/01/2013  option 6:  pas fractionnaire, moyenne=0, vitesse convectée projetée  Factorisation de Cholesky   Battement derriere le cylindre à partir de t=37s


Exécution séquentielle 1 thread + directives OpenMP: 0.414s/Pas de temps de 0.05s Last Command Real Time=1019s
Time=119.953   VelocityMean=0.999817  VelocityMax=1.47891    PressureMean=0.756415  PressMax-Min=1.37185    Flux-=-12.0000  Flux+=11.8570  Flux-+=-0.142951    
04/02/2013  option 5:  schema implicite et vitesse convectée en sortie Factorisation de Crout  Battement derriere le cylindre à partir de t=52s

Exécution parallèle 6 threads + directives OpenMP: 0.354/Pas de temps de 0.05s    Last Command Real Time=874s  => Accélération 1.17
Time=119.953   VelocityMean=0.997329  VelocityMax=1.47905    PressureMean=0.748277  PressMax-Min=1.41327    Flux-=-12.0000  Flux+=11.7531  Flux-+=-0.246917    
04/02/2013  option 5:  schema implicite et vitesse convectée en sortie Factorisation de Crout  Battement derriere le cylindre à partir de t=40s



Exécution parallèle 6 threads + directives OpenMP: 0.324/Pas de temps de 0.1s     Last Command Real Time=412s  => Accélération ?
Time=119.999   VelocityMean=0.991056  VelocityMax=1.45998    PressureMean=0.809014       PressMax-Min=1.44206     Flux-=-12.0000  Flux+=11.8869  Flux-+=-0.113080    
28/01/2013  option 5:  schema implicite et vitesse convectée en sortie  Gradient conjugue  Battement derriere le cylindre à partir de t=30s

Résultats pour différents pas de temps, 12470 triangles BREZZI-FORTIN, 6359 sommets, 6359+12470 noeuds et 44017 d.l. Vitesse(P1+Bulle)+Pression(P1):
Exécution séquentielle   AUCUNE directive OpenMP: 0.185s/Pas de temps de 0.05s  Last Command Real Time=443s
Time=120.003    VelocityMean=1.01189  VelocityMax=1.49265    PressureMean=0.523721E-02   PressMax-Min=1.39074     Flux-=-12.0000  Flux+=12.4421  Flux-+=0.442092    
14/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée   Gradient conjugue

Exécution parallèle 4 threads + directives OpenMP: 0.125s/Pas de temps de 0.05s   Last Command Real Time=300s  => Accélération 1.48
Time=120.003    VelocityMean=1.01103  VelocityMax=1.56386    PressureMean=-0.478190E-01  PressMax-Min=1.60427     Flux-=-12.0000  Flux+=12.4787  Flux-+=0.478721    
18/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée  Gradient conjugue

Exécution parallèle 6 threads + directives OpenMP: 0.113s/Pas de temps de 0.05s   Last Command Real Time=272s  => Accélération 1.64
Time=120.003    VelocityMean=1.00852  VelocityMax=1.56016    PressureMean=-0.424279E-01  PressMax-Min=1.58574     Flux-=-12.0000  Flux+=12.3827  Flux-+=0.382734    
18/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée  Gradient conjugue

Exécution parallèle 8 threads + directives OpenMP: 0.119s/Pas de temps de 0.05s   Last Command Real Time=286s  => Accélération 1.55
Time=120.003    VelocityMean=1.01139  VelocityMax=1.56533    PressureMean=-0.435052E-01  PressMax-Min=1.59472     Flux-=-12.0000  Flux+=12.5037  Flux-+=0.503736    
18/01/2013  option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugue

Exécution parallèle 6 threads + directives OpenMP: 0.120/Pas de temps de 0.05s    Last Command Real Time=288s  => Accélération 1.54
Time=120.003    VelocityMean=1.01001  VelocityMax=1.52277    PressureMean=-0.743607E-02  PressMax-Min=1.46389     Flux-=-12.0000  Flux+=12.3771  Flux-+=0.377110    
31/01/2013  option 6:  schema implicite et vitesse convectée projetee  Factorisation Cholesky   Battement derriere le cylindre à partir de t=17s


Exécution séquentielle   AUCUNE directive OpenMP: 0.189s (0.356s)/Pas de temps de 0.1s   Last Command Real Time=456s  (863s avec l'option -g au lieu de -O)
Time=239.905    VelocityMean=1.00193  VelocityMax=1.46058    PressureMean=0.866541      PressMax-Min=1.49268   Flux-=-12.0000  Flux+=12.2245  Flux-+=0.224476     
25/01/2013  option 5:  schema implicite et vitesse convectée en sortie Factorisation de Crout  Battement derriere le cylindre à partir de t=130s

Exécution parallèle 6 threads + directives OpenMP: 0.156s (0.283s)/Pas de temps de 0.1s  Last Command Real Time=378s (689s avec l'option -g au lieu de -O)  => Accélération 1.21
Time=239.905    VelocityMean=1.00419  VelocityMax=1.44304    PressureMean=0.741196     PressMax-Min=1.36728    Flux-=-12.0000  Flux+=12.3214  Flux-+=0.321442   option -O 
Time=239.905    VelocityMean=1.00117  VelocityMax=1.47497    PressureMean=0.892127     PressMax-Min=1.52457    Flux-=-12.0000  Flux+=12.2182  Flux-+=0.218224   option -g  
26/01/2013  option 5:  schema implicite et vitesse convectée en sortie Factorisation de Crout  Battement derriere le cylindre à partir de t=30s

Résultats pour 568 008 tétraèdres BREZZI-FORTIN, 103170 sommets, 103170+568008 noeuds et 2 116 704 d.l. Vitesse(P1Bulle)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] et [PG] matrices morses de 791 120 double precision
12/2/2013 Sur PCap: Exécution parallèle 6 threads + directives OpenMP: 10.36s/Pas de temps de 0.1s   Last Command Real Time=12459s
Time=119.999    VelocityMean=0.981748  VelocityMax=1.67351    PressureMean=-0.355406E-01 PressMax-Min=10.6557   Flux-=-84.1674  Flux+=84.1960  Flux-+=0.286818E-01

Résultats pour 568 008 tétraèdres TAYLOR-HOOD, 103170 sommets, 791120 noeuds et 2 476 530 d.l. Vitesse(P2)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] matrice morse de 11 475 722 et [PG] de 791 120 double precision
 
12/2/2013 Sur PCap: Exécution parallèle 6 threads + directives OpenMP: 10.36s/Pas de temps de 0.1s   Last Command Real Time=52312s
Time=119.999    VelocityMean=0.990824  VelocityMax=1.44413    PressureMean=-0.368969E-01 PressMax-Min=4.35235   Flux-=-84.0000 Flux+=84.3260 Flux-+=0.325976    

Résultats pour 183 600 tétraèdres TAYLOR-HOOD, 34 562 sommets, 260 484 noeuds et 816 014 d.l. Vitesse(P2)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] matrice morse de 3 738 810 et [PG] de 260 484 double precision
 
15/4/2013 Sur PCap :Exécution parallèle  6 threads + directives OpenMP: 6.12s/Pas de temps de 0.1s  Last Command Real Time=7375s  fluinth6: OMP seconds=7342s
Time=119.999   VelocityMean=0.995118  VelocityMax=1.47832   PressureMean=-0.334461E-01 PressMax-Min=3.18330  Flux-=-84.0000 Flux+=84.6626 Flux-+ =0.662656


16/4/2013 Sur hpc1: Exécution parallèle  1 thread  + directives OpenMP: 21.90s/Pas de temps de 0.1s Last Command Real Time=26320s fluinth6: OMP seconds=26277s
Time=119.999   VelocityMean=0.994124  VelocityMax=1.47347   PressureMean=-0.328748E-01 PressMax-Min=3.14507  Flux-=-84.0000 Flux+=84.6764 Flux-+=0.676442

15/4/2013 Sur hpc1: Exécution parallèle  6 threads + directives OpenMP: 4.28s/Pas de temps de 0.1s  Last Command Real Time=5183s  fluinth6: OMP seconds=5141s
Time=119.999   VelocityMean=0.995839  VelocityMax=1.47238   PressureMean=-0.334718E-01 PressMax-Min=3.29742  Flux-=-84.0000 Flux+=84.5903 Flux-+=0.590335

16/4/2013 Sur hpc1: Exécution parallèle  6 threads + directives OpenMP: 4.10s/Pas de temps de 0.1s  Last Command Real Time=4967s  fluinth6: OMP seconds=4925s
Time=119.999   VelocityMean=0.994420  VelocityMax=1.45713   PressureMean=-0.359343E-01 PressMax-Min=3.14133  Flux-=-84.0000 Flux+=84.5976 Flux-+=0.597564

17/4/2013 Sur hpc1: Exécution parallèle  8 threads + directives OpenMP: 3.32s/Pas de temps de 0.1s  Last Command Real Time=4027s  fluinth6: OMP seconds=3986s
Time=119.999   VelocityMean=0.993810  VelocityMax=1.47010   PressureMean=-0.328934E-01 PressMax-Min=3.11924  Flux-=-84.0000 Flux+=84.6651 Flux-+=0.665069

18/4/2013 Sur hpc1: Exécution parallèle  9 threads + directives OpenMP:14.16s/Pas de temps de 0.1s Last Command Real Time=17039s  fluinth6: OMP seconds=16998s
Time=119.999   VelocityMean=0.994806  VelocityMax=1.47239   PressureMean=-0.327808E-01 PressMax-Min=3.14816  Flux-=-84.0000 Flux+=84.6542 Flux-+=0.654204

18/4/2013 Sur hpc1: Exécution parallèle 10 threads + directives OpenMP:15.64s/Pas de temps de 0.1s Last Command Real Time=18807s  fluinth6: OMP seconds=18765s
Time=119.999   VelocityMean=0.994110  VelocityMax=1.47184   PressureMean=-0.335721E-01 PressMax-Min=3.14987  Flux-=-84.0000 Flux+=84.6889 Flux-+=0.688952

15/4/2013 Sur hpc1: Exécution parallèle 16 threads + directives OpenMP:15.53s/Pas de temps de 0.1s  Last Command Real Time=18682s fluinth6: OMP seconds=18639s
Time=119.999   VelocityMean=0.994337  VelocityMax=1.47141   PressureMean=-0.369067E-01 PressMax-Min=3.13744  Flux-=-84.0000 Flux+=84.6554 Flux-+=0.655436    

16/4/2013 Sur hpc1: Exécution parallèle 16 threads + directives OpenMP:15.44s/Pas de temps de 0.1s  Last Command Real Time=18580s fluinth6: OMP seconds=18537s
Time=119.999   VelocityMean=0.993318  VelocityMax=1.46962   PressureMean=-0.342513E-01 PressMax-Min=3.13219  Flux-=-84.0000 Flux+=84.6503 Flux-+=0.650289

17/4/2013 Sur hpc1: Exécution parallèle 16 threads + directives OpenMP:15.49s/Pas de temps de 0.1s  Last Command Real Time=18629s fluinth6: OMP seconds=18586s
Time=119.999   VelocityMean=0.992589  VelocityMax=1.43660   PressureMean=-0.343703E-01 PressMax-Min=3.11290  Flux-=-84.0000 Flux+=84.6458 Flux-+=0.645764

Résultats pour 183 600 tétraèdres BREZZI-FORTIN, 34562 sommets, 34562+183600 noeuds et 689 048 d.l. Vitesse(P1Bulle)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] et [PG] matrices morses de 260 484 double precision
13/4/2013 Sur PCap: Exécution parallèle 6 threads + directives OpenMP: 3.16s/Pas de temps de 0.1s  Last Command Real Time=3797s  fluinbf6: OMP seconds=3787s
Time=119.999   VelocityMean=0.993474  VelocityMax=1.63596   PressureMean=-0.305358E-01 PressMax-Min=6.14201   Flux-=-84.4598  Flux+=84.5794  Flux-+=0.119550    


16/4/2013 Sur hpc1: Exécution parallèle  1 thread  + directives OpenMP:11.65s/Pas de temps de 0.1s  Last Command Real Time=13997s fluinbf6: OMP seconds=13982s
Time=119.999   VelocityMean=0.993470  VelocityMax=1.63522   PressureMean=-0.294923E-01 PressMax-Min=6.13848   Flux-=-84.4598  Flux+=84.5669  Flux-+=0.107066

14/4/2013 Sur hpc1: Exécution parallèle  6 threads + directives OpenMP: 2.31s/Pas de temps de 0.1s  Last Command Real Time=2789s  fluinbf6: OMP seconds=2775s
Time=119.999   VelocityMean=0.993576  VelocityMax=1.63611   PressureMean=-0.309648E-01 PressMax-Min=6.12647   Flux-=-84.4598  Flux+=84.5990  Flux-+=0.139128     

16/4/2013 Sur hpc1: Exécution parallèle  6 threads + directives OpenMP: 2.36s/Pas de temps de 0.1s  Last Command Real Time=2845s  fluinbf6: OMP seconds=2830s
Time=119.999   VelocityMean=0.992593  VelocityMax=1.64120   PressureMean=-0.306766E-01 PressMax-Min=6.14202   Flux-=-84.4598  Flux+=84.5300  Flux-+=0.0701232

17/4/2013 Sur hpc1: Exécution parallèle  8 threads + directives OpenMP: 2.04s/Pas de temps de 0.1s  Last Command Real Time=2456s  fluinbf6: OMP seconds=2442s
Time=119.999   VelocityMean=0.994047  VelocityMax=1.63202   PressureMean=-0.279622E-01 PressMax-Min=6.12356   Flux-=-84.4598  Flux+=84.6500  Flux-+=0.190123  

18/4/2013 Sur hpc1: Exécution parallèle  9 threads + directives OpenMP:10.72s/Pas de temps de 0.1s  Last Command Real Time=12874s fluinbf6: OMP seconds=12860s
Time=119.999   VelocityMean=0.992872  VelocityMax=1.62581   PressureMean=-0.302102E-01 PressMax-Min=6.08857   Flux-=-84.4598  Flux+=84.6379  Flux-+=0.178043

18/4/2013 Sur hpc1: Exécution parallèle 10 threads + directives OpenMP:12.71s/Pas de temps de 0.1s  Last Command Real Time=15251s fluinbf6: OMP seconds=15265s
Time=119.999   VelocityMean=0.995780  VelocityMax=1.62721   PressureMean=-0.259526E-01 PressMax-Min=6.13514   Flux-=-84.4598  Flux+=85.0126  Flux-+=0.552746

14/4/2013 Sur hpc1: Exécution parallèle 16 threads + directives OpenMP:13.36s/Pas de temps de 0.1s  Last Command Real Time=16045s fluinbf6: OMP seconds=16030s
Time=119.999   VelocityMean=0.997257  VelocityMax=1.62605   PressureMean=-0.307908E-01 PressMax-Min=6.08522   Flux-=-84.4598  Flux+=85.7643  Flux-+=1.30448   

16/4/2013 Sur hpc1: Exécution parallèle 16 threads + directives OpenMP:13.39s/Pas de temps de 0.1s  Last Command Real Time=16045s fluinbf6: OMP seconds=16087s
Time=119.999   VelocityMean=0.987033  VelocityMax=1.61448   PressureMean=-0.268118E-01 PressMax-Min=6.07616   Flux-=-84.4598  Flux+=84.7808  Flux-+=0.320921 

Résultats pour 374 100 tétraèdres TAYLOR-HOOD, 69 949 sommets, 528 948 noeuds et 1 656 793 d.l. Vitesse(P2)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] matrice morse de 7 607 145 et [PG] de 528 948 double precision
 
20/4/2013 sur hpc1: Exécution parallèle 1 thread  + directives OpenMP:53.17s/Pas de temps de 0.1s  Last Command Real Time=63809s fluinth6: OMP seconds=63898s
Temps=119.999  VitesseMoy=0.980524  VitesseMax=1.47844   PressionMoy=-0.319144E-01 PressionMax-Min=4.18973  Flux-=-84.0000 Flux+=84.0194 Flux-+=0.193857E-01

19/4/2013 sur hpc1: Exécution parallèle 8 threads + directives OpenMP: 7.97s/Pas de temps de 0.1s  Last Command Real Time= 9648s fluinth6: OMP seconds= 9563s
Temps=119.999  VitesseMoy=0.975916  VitesseMax=1.41335   PressionMoy=-0.295486E-01 PressionMax-Min=4.18274  Flux-=-84.0000 Flux+=84.0159 Flux-+=0.158892E-01

21/4/2013 sur hpc1: Exécution parallèle 9 threads + directives OpenMP:30.78s/Pas de temps de 0.1s  Last Command Real Time=37023s fluinth6: OMP seconds=36936s
Temps=119.999  VitesseMoy=0.979432  VitesseMax=1.41096   PressionMoy=-0.345733E-01 PressionMax-Min=4.21081  Flux-=-84.0000 Flux+=84.0305 Flux-+=0.305389E-01

24/4/2013 sur hpc1: Exécution parallèle 16threads + directives OpenMP:31.51s/Pas de temps de 0.1s  Last Command Real Time=37901s fluinth6: OMP seconds=37812s
Time=119.999   VelocityMean=0.977942 VelocityMax=1.45538 PressureMean=-0.391074E-01 PressMax-Min=4.21460    Flux-=-84.0000 Flux+=83.8474 Flux-+=-0.152548

26/4/2013 sur hpc1: Exécution parallèle 8 threads +PROSCD=sum(V1(:)*V2(:)): 8.28s/Pas de temps de 0.1s  Last Command Real Time=10017s fluinth6: OMP seconds=9932s mem=3421876kb walltime=02:47:40
Time=119.999  VelocityMean=0.976584 VelocityMax=1.40277 PressureMean=-0.319825E-01 PressMax-Min=4.12076    Flux-=-84.0000 Flux+=84.0570 Flux-+=0.570123E-01

Résultats pour 374 100 tétraèdres BREZZI-FORTIN, 69 949 sommets, 69949+374100=444 049 noeuds, 1 402 096 d.l. Vitesse(P1Bulle)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée Gradient conjugué sur [VG] et [PG] de 528 948 double precision
 
20/4/2013 sur hpc1: Exécution parallèle 1 thread  + directives OpenMP:25.98s/Pas de temps de 0.1s  Last Command Real Time=31213s fluinbf6: OMP seconds=31180s
Temps=119.999 VitesseMoy=0.991423  VitesseMax=1.82233   PressionMoy=-0.239710E-01 PressionMax-Min=14.9848  Flux-=-84.2264 Flux+=85.6134 Flux-+=1.38698

20/4/2013 sur hpc1: Exécution parallèle 8 threads + directives OpenMP: 4.17s/Pas de temps de 0.1s  Last Command Real Time= 5038s fluinbf6: OMP seconds=5009s
Temps=119.999 VitesseMoy=0.991187  VitesseMax=1.82292   PressionMoy=-0.237756E-01 PressionMax-Min=14.9877  Flux-=-84.2264 Flux+=85.6154 Flux-+=1.38894

21/4/2013 sur hpc1: Exécution parallèle 9 threads + directives OpenMP:21.76s/Pas de temps de 0.1s  Last Command Real Time=26145s fluinbf6: OMP seconds=26116s
Temps=119.999 VitesseMoy=0.990637  VitesseMax=1.83617   PressionMoy=-0.223753E-01 PressionMax-Min=15.0625  Flux-=-84.2264 Flux+=85.7173 Flux-+=1.49088

24/4/2013 sur hpc1: Exécution parallèle 9 threads + directives OpenMP:24.92s/Pas de temps de 0.1s  Last Command Real Time=29934s fluinbf6: OMP seconds=29905s
Time=119.999  VelocityMean=0.991666 VelocityMax=1.83452 PressureMean=-0.351805E-01 PressMax-Min=15.0508   Flux-=-84.2264 Flux+=86.2805 Flux-+=2.05407

26/4/2013 sur hpc1: Exécution parallèle 8 threads + PROSCD=sum(V1(:)*V2(:)):4.22s/Pas de temps de 0.1s  Last Command Real Time=5097s fluinbf6: OMP seconds=5066s mem=2823132kb walltime=01:25:36
Time=119.999  VelocityMean=0.991337 VelocityMax=1.82340  PressureMean=-0.233241E-01 PressMax-Min=14.9890   Flux-=-84.2264 Flux+=85.6566 Flux-+=1.43013
Résultats pour 374 100 tétraèdres BREZZI-FORTIN, 69 949 sommets, 69949+374100=444 049 noeuds, 1 402 096 d.l. Vitesse(P1Bulle)+Pression(P1):
option 6: pas fractionnaire, moyenne=0, vitesse convectée projetée, Factorisation de CHOLESKY sur [VG] et [PG] de 69 445 473 double precision
2/5/2013 sur hpc1: Exécution parallèle 1 thread  + directives OpenMP:27.50s/Pas de temps de 0.1s  Last Command Real Time=33175s fluinbf6: OMP seconds=33143s mem=3892180kb walltime=09:13:36
Time=119.999  VelocityMean=0.994639 VelocityMax=1.82467  PressureMean=-0.441357E-01 PressMax-Min=15.0359   Flux-=-84.2264 Flux+=86.0558 Flux-+=1.82937

2/5/2013 sur hpc1: Exécution parallèle 8 threads + directives OpenMP: 6.90s/Pas de temps de 0.1s  Last Command Real Time= 8450s fluinbf6: OMP seconds= 8423s mem=2555940kb walltime=02:21:12
Time=119.999  VelocityMean=0.994744 VelocityMax=1.82479  PressureMean=-0.447416E-01 PressMax-Min=15.0070   Flux-=-84.2264 Flux+=86.0814 Flux-+=1.85496

CONCLUSIONS:
- L'algorithme lui-même, initialement programmé, n'a pas été modifié.

- L'apprentissage d'OpenMP, par simple ajoût de directives peu nombreuses, N'EST PAS DIFFICILE.

- Par contre, les CALCULS, avec ou sans OpenMP, NE DONNENT PAS les MÊMES RÉSULTATS,
  . soit, par un conflit mémoire NON RÉSORBÉ en l'absence d'une directive adéquate (ATOMIC ou CRITICAL) avec des différences importantes;
  . soit, parce que les calculs sont effectués dans un ORDRE DIFFÉRENT avec des différences portant sur le 3-ème chiffre et 
                que les opérations + - * ne sont pas associatives sur les registres de calcul.
  L'ajoût de directives OpenMP AJOUTE UNE DIFFICULTÉ (s'assurer d'obtenir à peu près les mêmes résultats)
  par rapport à la programmation purement séquentielle.
  Mais cela reste un problème purement technique et non algorithmique.

  Si l'algorithme est instable, l'ajoût de directives OpenMP peut le révéler, le détecter,
  le mettre en valeur, lors des exécutions avec un nombre différent de threads.

- L'accélération de la restitution des résultats est d' autant plus grande que le nombre de degés de liberté
  est grand devant le nombre de threads.
  L'exécution avec des threads de problème 2d naturellement plus pauvre en nombre de degrés de liberté est
  moins favorable que celle des problèmes en 3d plus volumineux en nombre de degrés de liberté.

- Les temps de restitution du résultat DIMINUENT de 1 à 6 en fonction du nombre CROISSANT de threads utilisés de 1 à 8
  jusqu'à un certain SEUIL (8 threads), puis AUGMENTENT très fortement!

  Avec un seul thread, OpenMP ralentit l'exécution. Cela s'explique par le fait qu'OpenMP ajoute des instructions,
  fait des calculs intermédiaires, des accès mémoire différents que dans le cas séquentiel, des synchronisations,...

  Mais plus génant, un processeur a une mémoire associée proche et tout accès à la mémoire éloignée associée à un autre
  processeur peut être largement plus couteuse. Notamment sur hpc1, ce temps d'accès aux mémoires associées aux autres
  processeurs va de 1 à 5.5!

  Les 4 derniers exemples montrent que le NOMBRE OPTIMAL de threads est 8 ce qui est le nombre 8 de threads par processeur
  sur PCap ou hpc1 où la mémoire nécessaire pour les calculs n'a pas excèdé 8 Gigaoctets,
  ce qui est bien loin sur hpc1 de la trentaine de  Gigaoctets libres de la mémoire associée à chacun de ses 20 processeurs.
  Le passage à 9 threads est une vraie catastrophe! Les temps de restitution s'accroissent de 1 à 5!

Moralité:
Si la mémoire nécessaire aux calculs n'excède pas celle associée à un processeur, le nombre optimal
de threads à choisir pour OpenMP est égal au nombre de threads par processeur (ici 8);


Finalement et globalement, cette adjonction de directives OpenMP a permis d'obtenir PLUS RAPIDEMENT le résultat 
d'un long calcul en sollicitant le nombre maximal de threads par processeur, sur notre PC standard (ou sur hpc1),
sans aucun côut financier supplémentaire et cette restitution plus rapide des résultats,
RENTABILISE BIEN l'EFFORT d'ajouter des directives OpenMP simples au code Fortran ou C...

Ainsi, pour atteindre une accélération de plusieurs fois le nombre de threads par processeur,
il sera nécessaire de faire appel à d'autres outils, de type MPI ou équivalents
et adapter les algorithmes de calcul (Méthodes de sous-domaines avec coloriage optimisé, ...)
à ces outils conduisant à des coûts d'analyse et programmation importants.


Alain Perronnet, Dernière mise à jour: 29 juillet 2013