Isend-Irecv non bloquant

De Wiki de Calcul Québec
Aller à : Navigation, rechercher
Autres langues :anglais 100% • ‎français 100%

Les exemples 5 et 6 illustrent la différence entre les fonctions MPI_Send/MPI_Recv et MPI_Isend/MPI_Irecv.

Chaque processus a pour mission de transférer un vecteur de nombres aléatoires (sendbuff) de dimension buffsize au processus suivant qui le reçoit dans recvbuff. Le dernier processus envoie son vecteur sendbuff au processus 0. La dimension buffsize est fournie en argument du programme à l'exécution. Cet exemple montre que les fonctions MPI_ISend et MPI_IRecv sont beaucoup plus appropriées que MPI_Send et MPI_Recv pour effectuer ce travail. Les fonctions MPI_Isend et MPI_Irecv sont non bloquantes, c'est-à-dire que la fonction retourne avant que la communication ne soit terminée. On peut donc placer MPI_Isend avant MPI_IRecv sans avoir peur d'un blocage (c'est-à-dire qu'un processus attende indéfiniment un message qui ne viendra pas, parce que le processus duquel vient le message attend aussi un message avec un MPI_Recv précédant le MPI_Send!!).

Puisque les fonctions MPI_Isend et MPI_Irecv n'attendent pas que la communication soit complétée, il faut placer en aval un MPI_Wait pour chacune des requêtes d'envoi et de réception du processus. Cette fonction, comme son nom l'indique, attend que les communications soient complétées avant de poursuivre. Dans le cas présent, chaque processus poste une requête d'envoi (MPI_Isend) et une requête de réception (MPI_Irecv). On voit donc clairement ici que les communications s'effectueront en même temps entre les différents processus.

Avant l'étape de communication, chaque processus calcule la somme du vecteur à envoyer. Les résultats de ces sommes sont ensuite rapatriés sur le processus 0 (avec un MPI_Gather) pour affichage à l'écran. Après l'étape de communication, chaque processus calcule la somme du vecteur reçu, et de la même façon ces sommes sont rapatriées sur le processus 0 pour affichage à l'écran. L'exemple 5 montre l'utilisation de MPI_Send et MPI_Recv pour effectuer le même travail. Essayez les deux pour comparer!

En Fortran

Fichier : isend-irecv.f
!--------------------------------------------------------
!  
!              sendbuff   recvbuff       sendbuff   recvbuff
!                                         
!              ########   ########       ########   ########
!              #      #   #      #       #      #   #      #
!          0   #  AA  #   #      #       #  AA  #   #  EE  #
!              #      #   #      #       #      #   #      #
!              ########   ########       ########   ########
!     T        #      #   #      #       #      #   #      #
!          1   #  BB  #   #      #       #  BB  #   #  AA  #
!     a        #      #   #      #       #      #   #      #
!              ########   ########       ########   ########
!     c        #      #   #      #       #      #   #      #
!          2   #  CC  #   #      #       #  CC  #   #  BB  #
!     h        #      #   #      #       #      #   #      #
!              ########   ########       ########   ########
!     e        #      #   #      #       #      #   #      #
!          3   #  DD  #   #      #       #  DD  #   #  CC  #
!     s        #      #   #      #       #      #   #      #
!              ########   ########       ########   ########
!              #      #   #      #       #      #   #      #
!          4   #  EE  #   #      #       #  EE  #   #  DD  #
!              #      #   #      #       #      #   #      #
!              ########   ########       ########   ########
!              
!                     AVANT                     APRÈS
!
! Auteur : Carol Gauthier
!          Centre de calcul scientifique
!          Université de Sherbrooke
!
! Dernière révision : août 2004
!--------------------------------------------------------
Program Exemple_MPI
 
  include 'mpif.h'
  integer ierr,ntasks,taskid,itask,status(MPI_STATUS_SIZE)
  integer i,j,k,buffsize,sentid,send_request,recv_request
  character argtmp*12
  real(8) inittime,recvtime,totaltime,rand,buffsum
 
  real(8),allocatable,dimension(:) :: recvtimes
  real(8),allocatable,dimension(:) :: buffsums
 
  real(8),allocatable,dimension(:) :: sendbuff
  real(8),allocatable,dimension(:) :: recvbuff
 
  !---------------------------------------------------------------
  ! Initialisation de MPI. Il est important de placer cet appel   
  ! au début du programme, immédiatement après les déclarations.
  call MPI_INIT( ierr )
 
  !---------------------------------------------------------------
  ! Obtenir le nombre de processus MPI et le numéro d'identification
  ! du processus courant.
  call MPI_COMM_SIZE(MPI_COMM_WORLD,ntasks,ierr)
  call MPI_COMM_RANK(MPI_COMM_WORLD,taskid,ierr)
 
  !---------------------------------------------------------------
 ! Dimension  du vecteur à partir des arguments
  call getarg(1,argtmp)
  read(argtmp,'(I12)')buffsize
 
  !---------------------------------------------------------------
  ! Affichage de la description de l'exemple.
  if ( taskid.eq.0 )then
    write(6,'(A)')
    write(6,'(A)')"##########################################################"
    write(6,'(A)')
    write(6,'(A)')" Exemple 6"
    write(6,'(A)')
    write(6,'(A)')" Communication point à point : MPI_Isend MPI_Irecv"
    write(6,'(A)')
    write(6,'(A,I12)')" Dimension du vecteur :",buffsize
    write(6,'(A,I5)')" Nombre total de processus :",ntasks
    write(6,'(A)')
    write(6,'(A)')"##########################################################"
    write(6,'(A)')
    write(6,'(A)')"                --> AVANT COMMUNICATION <--"
    write(6,'(A)')
  endif
 
   !---------------------------------------------------------------
   ! Allocation de la mémoire 
    allocate( sendbuff(0:buffsize-1) )
    allocate( recvbuff(0:buffsize-1) )
    allocate( recvtimes(0:ntasks-1) )
    allocate( buffsums(0:ntasks-1) )
 
   !-----------------------------------------------------------------
   ! Initialisation des vecteurs ou tableaux.
   call srand(taskid*10)
   do i=0,buffsize-1
    sendbuff(i)=rand();
   end do
 
   !-----------------------------------------------------------------
   ! Affichage avant communication.
 
   buffsum=0.0
   do i=0,buffsize-1
      buffsum=buffsum+sendbuff(i)
   end do
 
   call MPI_Gather(buffsum,1,MPI_REAL8,buffsums,1, MPI_REAL8,0,MPI_COMM_WORLD,ierr)
 
   if(taskid.eq.0)then
    do itask=0,ntasks-1
      sentid=itask+1
      if(itask.eq.(ntasks-1))sentid=0
      write(6,'(A,I3,A,I3,A,E14.8)')  &
      &       "Processus",itask," : somme du vecteur envoyé à ",sentid," =", buffsums(itask)
    end do
   end if
 
   !-----------------------------------------------------------------
   ! Communication
 
   inittime = MPI_Wtime();
 
   if ( taskid.eq.0 )then
     call MPI_Irecv(recvbuff,buffsize,MPI_REAL8,ntasks-1,MPI_ANY_TAG, &
     &             MPI_COMM_WORLD,recv_request,ierr)
     recvtime = MPI_Wtime()
     call MPI_Isend(sendbuff,buffsize,MPI_REAL8,taskid+1,0,MPI_COMM_WORLD,send_request,ierr)  
   else if( taskid.eq.(ntasks-1) )then
     call MPI_Isend(sendbuff,buffsize,MPI_REAL8,0,0,MPI_COMM_WORLD,send_request,ierr)
     call MPI_Irecv(recvbuff,buffsize,MPI_REAL8,taskid-1,MPI_ANY_TAG, &
     &             MPI_COMM_WORLD,recv_request,ierr)
     recvtime = MPI_Wtime()
   else
     call MPI_Irecv(recvbuff,buffsize,MPI_REAL8,taskid-1,MPI_ANY_TAG, &
     &             MPI_COMM_WORLD,recv_request,ierr)
     recvtime = MPI_Wtime()
     call MPI_Isend(sendbuff,buffsize,MPI_REAL8,taskid+1,0,MPI_COMM_WORLD,send_request,ierr)
   end if
 
   call MPI_Wait(send_request,status,ierr);
   call MPI_Wait(recv_request,status,ierr);
 
   totaltime=MPI_Wtime()
 
   !-----------------------------------------------------------------
   ! Affichage après communication.
 
   buffsum=0.0
   do i=0,buffsize-1
      buffsum=buffsum+recvbuff(i)
   end do
   if ( taskid.eq.0 )then
    do itask=0,ntasks-1
      buffsums(itask)=0.0
    end do     
   end if
 
   call MPI_Gather(buffsum,1,MPI_REAL8,buffsums,1, MPI_REAL8,0,MPI_COMM_WORLD,ierr)
 
   call MPI_Gather(recvtime,1,MPI_REAL8,recvtimes,1, MPI_REAL8,0,MPI_COMM_WORLD,ierr)
 
  if(taskid.eq.0)then
    write(6,*)
    write(6,'(A)')"##########################################################"
    write(6,*)
    write(6,'(A)')"                --> APRÈS COMMUNICATION <--"
    write(6,*)
    do itask=0,ntasks-1
      write(6,'(A,I3,A,F5.2,A,e14.8)')"Processus ",itask," : vecteur reçu à",recvtimes(itask), &
      &                             " secondes : somme = ",buffsums(itask)       
    end do
    write(6,*)
    write(6,'(A)')"##########################################################"
    write(6,'(A,F5.2,A)')" Temps total de communication : ",totaltime," secondes"
    write(6,'(A)')"##########################################################"
    write(6,*)
  end if
 
  !-----------------------------------------------------------------
  ! Libération de la mémoire
  deallocate(sendbuff)
  deallocate(recvbuff)
  deallocate(recvtimes)
  deallocate(buffsums)
 
  !-----------------------------------------------------------------
  ! Finalisation de MPI
 
  call MPI_FINALIZE( ierr )
 
end


En C

Fichier : isend_irecv.c
/*--------------------------------------------------------
 
              sendbuff   recvbuff       sendbuff   recvbuff
 
              ########   ########       ########   ########
              #      #   #      #       #      #   #      #
          0   #  AA  #   #      #       #  AA  #   #  EE  #
              #      #   #      #       #      #   #      #
              ########   ########       ########   ########
     T        #      #   #      #       #      #   #      #
          1   #  BB  #   #      #       #  BB  #   #  AA  #
     a        #      #   #      #       #      #   #      #
              ########   ########       ########   ########
     c        #      #   #      #       #      #   #      #
          2   #  CC  #   #      #       #  CC  #   #  BB  #
     h        #      #   #      #       #      #   #      #
              ########   ########       ########   ########
     e        #      #   #      #       #      #   #      #
          3   #  DD  #   #      #       #  DD  #   #  CC  #
     s        #      #   #      #       #      #   #      #
              ########   ########       ########   ########
              #      #   #      #       #      #   #      #
          4   #  EE  #   #      #       #  EE  #   #  DD  #
              #      #   #      #       #      #   #      #
              ########   ########       ########   ########
 
                     AVANT                     APRÈS
 
 Auteur : Carol Gauthier
          Centre de calcul scientifique
          Université de Sherbrooke
 
 Dernière révision : septembre 2005
-------------------------------------------------------- */
 
#include <malloc.h>
#include <stdlib.h>
#include <stdio.h>
#include <time.h>
#include "math.h"
#include "mpi.h"
 
int main(int argc,char** argv)
{
   /*===============================================================*/
   /* Déclaration des variables                                     */
   int          taskid, ntasks;
   MPI_Status   status;
   MPI_Request	send_request,recv_request;
   int          ierr,i,j,itask,recvtaskid;
   int	        buffsize;
   double       *sendbuff,*recvbuff;
   double       sendbuffsum,recvbuffsum;
   double       sendbuffsums[1024],recvbuffsums[1024];
   double       inittime,totaltime,recvtime,recvtimes[1024];
 
   /*===============================================================*/
   /* Initialisation de MPI. Il est important de placer cet appel   */
   /* au début du programme, immédiatement après les déclarations.   */
   MPI_Init(&argc, &argv);
 
   /*===============================================================*/
   /* Obtenir le nombre de processus MPI et le numéro d'identification */
   /* du processus courant taskid                                   */
   MPI_Comm_rank(MPI_COMM_WORLD,&taskid);
   MPI_Comm_size(MPI_COMM_WORLD,&ntasks);
 
   /*===============================================================*/
   /* Obtenir buffsize à partir des arguments.                      */
   buffsize=atoi(argv[1]);
 
   /*===============================================================*/
   /* Affichage de la description de l'exemple.                     */
   if ( taskid == 0 ){
     printf("\n\n\n");
     printf("##########################################################\n\n");
     printf(" Exemple 6 \n\n");
     printf(" Communication point à point : MPI_Isend MPI_Irecv \n\n");
     printf(" Dimension de chaque vecteur : %d\n",buffsize);
     printf(" Nombre total de processus : %d\n\n",ntasks);
     printf("##########################################################\n\n");
     printf("                --> AVANT COMMUNICATION <--\n\n");
   }
 
   /*=============================================================*/
   /* Allocation de la mémoire.                                   */ 
   sendbuff=(double *)malloc(sizeof(double)*buffsize);
   recvbuff=(double *)malloc(sizeof(double)*buffsize);
 
   /*=============================================================*/
   /* Initialisation des vecteurs ou tableaux                     */
   srand((unsigned)time( NULL ) + taskid);
   for(i=0;i<buffsize;i++){
     sendbuff[i]=(double)rand()/RAND_MAX;
   }
 
   /*==============================================================*/
   /* Affichage avant communication.                               */
 
   sendbuffsum=0.0;
   for(i=0;i<buffsize;i++){
     sendbuffsum += sendbuff[i];
   }   
   ierr=MPI_Gather(&sendbuffsum,1,MPI_DOUBLE,
                   sendbuffsums,1, MPI_DOUBLE,
                   0,MPI_COMM_WORLD);                  
   if(taskid==0){
     for(itask=0;itask<ntasks;itask++){
       recvtaskid=itask+1;
       if(itask==(ntasks-1))recvtaskid=0;
       printf("Processus %d : somme du vecteur envoyé à %d = %e\n",
               itask,recvtaskid,sendbuffsums[itask]);
     }  
   }
 
   /*===============================================================*/
   /* Communication.                                                */
 
   inittime = MPI_Wtime();
 
   if ( taskid == 0 ){
     ierr=MPI_Isend(sendbuff,buffsize,MPI_DOUBLE,
	           taskid+1,0,MPI_COMM_WORLD,&send_request);   
     ierr=MPI_Irecv(recvbuff,buffsize,MPI_DOUBLE,
	           ntasks-1,MPI_ANY_TAG,MPI_COMM_WORLD,&recv_request);
     recvtime = MPI_Wtime();
   }
   else if( taskid == ntasks-1 ){
     ierr=MPI_Isend(sendbuff,buffsize,MPI_DOUBLE,
	           0,0,MPI_COMM_WORLD,&send_request);   
     ierr=MPI_Irecv(recvbuff,buffsize,MPI_DOUBLE,
	           taskid-1,MPI_ANY_TAG,MPI_COMM_WORLD,&recv_request);
     recvtime = MPI_Wtime();
   }
   else{
     ierr=MPI_Isend(sendbuff,buffsize,MPI_DOUBLE,
	           taskid+1,0,MPI_COMM_WORLD,&send_request);
     ierr=MPI_Irecv(recvbuff,buffsize,MPI_DOUBLE,
	           taskid-1,MPI_ANY_TAG,MPI_COMM_WORLD,&recv_request);
     recvtime = MPI_Wtime();
   }
   ierr=MPI_Wait(&send_request,&status);
   ierr=MPI_Wait(&recv_request,&status);
 
   totaltime = MPI_Wtime() - inittime;
 
   /*===============================================================*/
   /* Affichage après communication.                                */
 
   recvbuffsum=0.0;
   for(i=0;i<buffsize;i++){
     recvbuffsum += recvbuff[i];
   }   
 
   ierr=MPI_Gather(&recvbuffsum,1,MPI_DOUBLE,
                   recvbuffsums,1, MPI_DOUBLE,
                   0,MPI_COMM_WORLD);
 
   ierr=MPI_Gather(&recvtime,1,MPI_DOUBLE,
                   recvtimes,1, MPI_DOUBLE,
                   0,MPI_COMM_WORLD);
 
   if(taskid==0){
     printf("##########################################################\n\n");
     printf("                --> APRÈS COMMUNICATION <-- \n\n");
     for(itask=0;itask<ntasks;itask++){
       printf("Processus %d : somme du vecteur reçu= %e : Temps=%f secondes\n",
               itask,recvbuffsums[itask],recvtimes[itask]);
     }  
     printf("\n");
     printf("##########################################################\n\n");
     printf(" Temps total de communication : %f secondes\n\n",totaltime);  
     printf("##########################################################\n\n");
   }
 
   /*===============================================================*/
   /* Libération de la mémoire                                      */
   free(recvbuff);
   free(sendbuff);
 
   /*===============================================================*/
   /* Finalisation de MPI                                           */
   MPI_Finalize();
}


Outils personnels
Espaces de noms

Variantes
Actions
Navigation
Ressources de Calcul Québec
Outils
Partager