        PROGRAM  PRODUIT_MPI
!
!===============================================================================
!===============================================================================
! CRIHAN     11/2000
! Guy Moebs   Guy.Moebs@crihan.fr
!
! Utilisation de la bibliotheque d'echange de messages MPI (www.mpi-forum.org)
!
! Programme SPMD
!
! Calcul du vecteur R = ( U*(tU) + LAMBDA*V*(tV) )  *  ( U + V )
!
!              i.e. RESULT = MAT * VECT
!
! 1.  Initialisation du vecteur par tous les processus
! 2.  Initialisation distribuee des ligne de la matrice par chaque processus
! 3.  Calcul du produit partiel par chaque processus :
!          lignes matrice * vecteur
! 4.  Tous les processus se redistribuent les produits partiels pour avoir
!     le resultat complet
!
!
!===============================================================================
!===============================================================================
!
!     DECLARATIONS
!
!===============================================================================
!===============================================================================
!
! MODULE(S) NECESSAIRE(S)
!
!     USE MPI
!
      use modscal
!
      IMPLICIT NONE
!
! COMMON(S) NECESSAIRE(S)
!
      include 'mpif.h'
!
! DECLARATION DES CONSTANTES
!
      INTEGER    ITYPE_REAL
      PARAMETER (ITYPE_REAL = MPI_DOUBLE_PRECISION)
!
! DECLARATION DES VARIABLES
!
      REAL(rp),  ALLOCATABLE,  DIMENSION(:,:)  ::  MAT
      REAL(rp),  ALLOCATABLE,  DIMENSION(:)    ::  VECT
      REAL(rp),  ALLOCATABLE,  DIMENSION(:)    ::  RCV_R
      REAL(rp),  ALLOCATABLE,  DIMENSION(:)    ::  RESULT
!
      INTEGER,   ALLOCATABLE,  DIMENSION(:)    ::  IPOS
      INTEGER,   ALLOCATABLE,  DIMENSION(:)    ::  IVAL
!
      REAL(rp)   ::  LAMBDA
!
      INTEGER    ::  N
      INTEGER    ::  NPROC, ME, NUM, I
      INTEGER    ::  IERR = 0
      INTEGER    ::  INDMIN, INDMAX
      INTEGER    ::  INTERV, NVAL, NVALMAX
!
      LOGICAL    ::  fichex
!
!===============================================================================
!===============================================================================
!
!     DEBUT DU PROGRAMME
!
!===============================================================================
!===============================================================================
!
!
!----Initialise MPI
!
      CALL MPI_INIT (IERR)
!
!
!----Donne le nombre de processus dans l'application
!
      CALL MPI_COMM_SIZE (MPI_COMM_WORLD, NPROC, IERR)
!
!
!----Mon rang dans l'application
!
      CALL MPI_COMM_RANK (MPI_COMM_WORLD, ME, IERR)
!
!
!----Barriere de synchronisation
!
      CALL MPI_BARRIER (MPI_COMM_WORLD, ierr)
!
      IF ( ME == 0 ) CALL system ("ps -u gm")
!
!
!----Lecture des dimensions du probleme
!
      IF ( ME == 0 ) THEN
!
         INQUIRE (file='mpi.dat', exist=fichex)
         IF ( .NOT. fichex ) THEN
            WRITE (6,1020) 'STOP STOP STOP STOP STOP'
            WRITE (6,1020) 'Le fichier d entree mpi.dat n est pas present'
            WRITE (6,1020) 'Le programme s arrete'
 1020 FORMAT (A)
            STOP
         END IF
!
         OPEN ( 19, FILE = 'mpi.dat', STATUS = 'OLD', FORM='FORMATTED' )
         READ (19, *) N
         CLOSE (19)
!
      END IF
!
!
!---Diffusion de la dimension
      CALL MPI_BCAST  ( N, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, IERR )
!
!
!----Determine les indices de calcul
!
      INTERV = N / NPROC
      INDMIN = ME * INTERV + 1
      INDMAX = INDMIN + INTERV - 1
      IF (ME == (NPROC-1) ) INDMAX = N
!
      NVAL    = INDMAX - INDMIN + 1
      NVALMAX = INTERV + MOD (N, NPROC)
!
      WRITE (6,1010) 'ME = ', ME, ' : ', &
                      INDMIN, ' => ', INDMAX, ' NVAL = ', NVAL
 1010 FORMAT ( 4(A, I4) )
!
!
!----Allocation dynamique
      ALLOCATE ( MAT(INDMIN:INDMAX,N), RESULT(INDMIN:INDMAX), VECT(N), &
                 RCV_R (N), IPOS (0:NPROC), IVAL (0:NPROC), STAT = IERR )
      IF (ierr /= 0) THEN
         WRITE (6,'(A)') ' Probleme allocation MAT, RESULT, VECT'
         WRITE (6,'(A)') ' IERR = ', IERR
         CALL MPI_ABORT (MPI_COMM_WORLD, 1, ierr)
      END IF
!
!
!----Initialisation des matrices
!
      CALL INIT_MAT_VECT (MAT, VECT, INDMIN, INDMAX, N)
!
!
!----Calcul du produit
!
      CALL PROD_LIGNE (MAT, VECT, RESULT, INDMIN, INDMAX, N)
!
!
!-----Envoi & reception du resultat
!
      IPOS (0) = 0
      DO NUM = 1, NPROC - 1
        IVAL (NUM-1) = INTERV
        IPOS (NUM)   = IPOS (NUM-1) + INTERV
      ENDDO
      IVAL (NPROC-1) = NVALMAX
!
      WRITE (6,1040) 'RESULTAT : PROC = ', ME, &
                      RESULT (INDMIN), RESULT (INDMAX)
 1040 FORMAT ( A, I4, 4E12.4 )
!
      CALL MPI_ALLGATHERV (RESULT (INDMIN), NVAL, ITYPE_REAL, RCV_R, &
                           IVAL, IPOS, ITYPE_REAL, MPI_COMM_WORLD, IERR)
!
      WRITE (6,1040) 'RESULTAT : PROC = ', ME, RCV_R (1), &
                      RCV_R (INDMIN),  RCV_R (INDMAX),  RCV_R (N)
!
!-----Liberation de la memoire
!
      DEALLOCATE ( MAT, VECT, RESULT, RCV_R, IPOS, IVAL )
!
!----Synchronisation finale
!
      CALL MPI_BARRIER (MPI_COMM_WORLD, IERR)
!
!
!----Quitter MPI
!
      CALL MPI_FINALIZE (IERR)
!
!===============================================================================
!===============================================================================
!
!     FIN DU PROGRAMME
!
!===============================================================================
!===============================================================================
!
      STOP
      END PROGRAM PRODUIT_MPI
!
!
!
!
!===============================================================
!
!     SUBROUTINE PROD_LIGNE : produit matrice vecteur distribue
!
!     MAT       (R ,E ) Matrice  du produit matriciel
!     VECT      (R ,E ) Vecteur  du produit matriciel
!     RESULT    (R, S ) Resultat du produit matriciel
!     INDMIN    (I ,E ) Indice initial des lignes distribuees
!     INDMAX    (I ,E ) Indice final   des lignes distribuees
!
!===============================================================
!
      SUBROUTINE PROD_LIGNE (MAT, VECT, RESULT, INDMIN, INDMAX, N)
!
! MODULE(S) NECESSAIRE(S)
!
      use modscal
!
      IMPLICIT NONE
!
! DECLARATION DES ARGUMENTS
!
      INTEGER, INTENT (IN)  ::  N
      INTEGER, INTENT (IN)  ::  INDMIN
      INTEGER, INTENT (IN)  ::  INDMAX
!
      REAL(rp),  DIMENSION(INDMIN:INDMAX,N), INTENT (IN)   ::  MAT
      REAL(rp),  DIMENSION(N),               INTENT (IN)   ::  VECT
      REAL(rp),  DIMENSION(INDMIN:INDMAX),   INTENT (OUT)  ::  RESULT
!
! DECLARATION DES VARIABLES LOCALES
!
      INTEGER  ::  I
      INTEGER  ::  J
!
! DEBUT DE LA SUBROUTINE
!
!
      RESULT (INDMIN:INDMAX) = zero
!
!------Produit en ligne : RESULT = MAT * VECT
!
      DO I = INDMIN, INDMAX
         DO J = 1, N
            RESULT (I) = RESULT (I) + MAT (I,J) * VECT(J)
         END DO
      END DO
!
! FIN DE LA SUBROUTINE
!
      RETURN
      END SUBROUTINE PROD_LIGNE
!
!
!
!
!===============================================================
!
!     SUBROUTINE INIT_MAT_VECT : initialisation matrice vecteur distribues
!
!     MAT       (R ,S ) Matrice  du produit matriciel
!     VECT      (R ,S ) Vecteur  du produit matriciel
!     INDMIN    (I ,E ) Indice initial des lignes distribuees
!     INDMAX    (I ,E ) Indice final   des lignes distribuees
!
!===============================================================
!
        SUBROUTINE INIT_MAT_VECT (MAT, VECT, INDMIN, INDMAX, N)
!
! MODULE(S) NECESSAIRE(S)
!
      use modscal
!
      IMPLICIT NONE
!
! DECLARATION DES ARGUMENTS
!
      INTEGER, INTENT (IN)  ::  N
      INTEGER, INTENT (IN)  ::  INDMIN
      INTEGER, INTENT (IN)  ::  INDMAX
!
      REAL(rp),  DIMENSION(INDMIN:INDMAX,N), INTENT (OUT)  ::  MAT
      REAL(rp),  DIMENSION(N),               INTENT (OUT)  ::  VECT
!
! DECLARATION DES VARIABLES LOCALES
!
      REAL(rp),  DIMENSION(N)    ::  U
      REAL(rp),  DIMENSION(N)    ::  V
!
      REAL(rp)   ::  LAMBDA
!
      INTEGER    ::  I, J
!
      REAL(rp)   ::  scal
!
! DEBUT DE LA SUBROUTINE
!
!
!----Initialisation des vecteurs U et V
!
      scal = 5.0_rp / 12.0_rp
      U (1:N) = scal
!
      DO I = 1, N
         V (I) = 0.2_rp + REAL (I, rp)
      END DO
!
!
!-----Calcul de LAMBDA
!
      LAMBDA = zero
      DO I = 1, N
         LAMBDA = LAMBDA + U (I) * V (I)
      END DO
!
!
!------Calcul de la matrice MAT = U*(tU) + LAMDA*V*(tV)
!
      DO J = 1, N
         DO I = INDMIN, INDMAX
             MAT (I,J) = U (I) * U (J) + LAMBDA * V (I) * V (J)
         END DO
      END DO
!
!
!------Calcul du vecteur VECT = V + U
      VECT (1:N) = U (1:N) + V (1:N)
!
! FIN DE LA SUBROUTINE
!
      RETURN
      END SUBROUTINE INIT_MAT_VECT

