TD ``Calcul Parallèle''

Éric Goubault

TD 1 
Quelques exercices C

1.1  Calcul de la fonction de Fibonacci

On cherche à calculer et à afficher la suite de valeurs définie par, Ecrire successivement un programme C itératif, puis un programme C récursif, puis un programme C récursif utilisant un argument entier supplémentaire appelé ``accumulateur'' permettant de ne pas recalculer trop souvent des itérées anciennes (à vous de réfléchir!).

1.2  Racine carrée

Codez la fonction racine carrée à l'aide de la méthode de Newton, c'est-à-dire où Öa s'obtient comme la limite de la suite récurrente associée à la fonction qui à x associe ((x+a/x ) / 2 ).

TD 2 
Quelques exercices C/PVM

2.1  Tri fusion

Reprendre l'exercice en fin des notes de TDs de la dernière fois. Une autre option pour le tri fusion parallèle est la suivante: Pourrait-on proposer encore une autre méthode de tri fusion? (eu égard aux communications entre B et M en (1) et (2))

2.2  Calcul d'intégrale

Ecrire un programme C/PVM sur le modèle du calcul de p (maître/esclave) qui effectue l'intégration d'une fonction par la méthode des trapèzes,
 
 
ó
õ
x1

x0

f(x)dx 
@
x1-x0
n
æ
ç
è
f(x0)+f(x1)
2
+ n-1
å
i = 1
f æ
ç
è
x0+i x1-x0
n
ö
÷
ø
ö
÷
ø

2.3  Crible d'Erathostène

Tony Hoare a proposé il y a longtemps de calculer les nombres premiers inférieurs à n2 en utilisant n+2 processus, qui recoivent de leur prédecesseur une suite croissante de nombres naturels dont le premier, p , est un nombre premier qui est imprimé, et transmettent à leur successeur les nombres de cette suite, dans le même ordre, à l'exclusion de p et de ses multiples.

En utilisant la bibliothèque PVM, écrire un programme C qui affiche la suite des nombres premiers par la méthode du crible d'Erathostène: un processus est chargé de générer les entiers naturels, dont on élimine d'abord les multiples de 2, puis 3, 5 etc. au moyen de processus filtrants successifs.

On recommande chaudement de montrer le programme (qui compile néanmoins) à un chargé de TP avant de l'essayer...

2.4  Equation de la chaleur

On considère l'équation de la chaleur sur l'ouvert W = ]0,1[×]0,1[ ,
rc u
t
= div(k(x).Ñu)+S(t).dxs(x)
avec les conditions aux bords
u = 0
et la condition initiale
u(x,0) = 0
k est la conductivité thermique du matériau considéré, c est la capacité calorifique massique et r la densité. S(t) est une source de chaleur et dxs est une masse de Dirac au point xs . Ici, on va supposer k constant dans tout le domaine.

On considère un schéma explicite en temps et une discrétisation explicite d'ordre deux en espace. On approxime donc u par ses valeurs discrétisées en les sommets (i,j) d'une grille régulière de pas h (en abscisse et ordonnée), et en des temps n (pas de temps Dt ), soit ui,jn . L'équation discrétisée est maintenant:
rc un+1i-uni
Dt
k
h2
(uni+1,j+uni-1,j+uni,j+1+uni,j-1-4uni,j)
Le critère de stabilité (expliqué également dans le prochain cours) de l'équation discrétisée nous donne une condition pour le choix du pas de temps:
Dt <  rch2
4k
Pour les applications numériques, on prendra
4kDt
rch2
1
2

et on choisira h .

On souhaite écrire un programme maître/esclave permettant de calculer le schéma précédent,

TD 3 
Exercices C/PVM (suite)

3.1  Calcul d'une opération associative

Soit Å une operation associative quelconque sur des nombres entiers (par exemple, max, ou + etc.). On se propose d'implémenter un programme PVM parallèle calculant rapidement cette opération associative sur un tableau de nombres entiers A par la méthode ``diviser pour régner''. Pour cela on écrira un maître qui contient au début le tableau A (que l'on initialisera comme on voudra, et de taille 2k , avec k > 1 ), puis qui lance récursivement deux fils avec chacun une moitié du tableau, ¼, ainsi de suite jusqu'à tomber sur des fils qui n'ont que deux nombres entiers. Ceux-ci effectuent l'opération associative et renvoient le résultat à leurs parents, ¼, et ainsi de suite. On pourra s'inspirer du code de FFT en annexe de cet exercice qui utilise la même méthode. Peut-on optimiser en réduisant le nombre de fils?

Remarque:

Ceci est un exercice un peu académique dans le sens où certains de ces algorithmes dits de ``réduction'' sont déjà implémentés sous PVM.

Rappel:

int numt;
numt=pvm_spawn(char *task,char **argv,int flag,char
                          *where,int ntask,int *tids);
Exemple:

Calcul récursif d'une transformée de fourier rapide (les sources sont accessibles sur la page web du cours, dans la section ``corrigés'').

``fft-startup'' est à la racine de l'arbre des appels récursifs:

#include <stdlib.h>
#include <stdio.h>
#include "pvm3.h"
#define RESULT_TAG 1
#define MAX_LENGTH 128

main(argc, argv)
int argc;
char *argv[];
{
/* argv[i] is the unidimensional complex signal to be processed */
/* format is fft Re(x_1) Im(x_1) Re(x_2) ... Re(x_n) Im(x_n) */
int i,bufid,mtid,tid,numt;
float R[MAX_LENGTH];
        if ((mtid=pvm_mytid())<0)
           {
                printf("unexpected error occured!\n");
                exit(1);
           }
        setlinebuf(stdout);
       if (argc < 5)
          {
               printf("wrong number of parameters\n");
               exit(1);
          }
       numt = pvm_spawn("fft",&argv[1],PvmTaskDefault,"",1,
                                   &tid);
        if (numt != 1)
           {
                printf("could not spawn fft\n");
                exit(1);
           }
        bufid = pvm_recv(tid, RESULT_TAG);
        pvm_upkfloat(&R[0],argc-1,1);
       for (i=0;i<(argc-1)/2;i++)
          {
               printf(" %g %g",R[i],R[i+(argc-1)/2]); }
       printf("\n");
       pvm_exit();
       exit(0); }
Le programme récursif ``fft'' se charge du calcul effectif:
#include <math.h>
#include <stdio.h>
#include "pvm3.h"

#define RESULT_TAG 1
#define SIZE 64
#define MAX_LENGTH 128
#define PI 3.14159265358

main(argc, argv)
int argc;
char *argv[];
{
/* argv[i] is the unidimensional complex signal to be processed */
/* format is fft Re(x_1) Im(x_1) Re(x_2) ... Re(x_n) Im(x_n) */
/* sends to his parent an array of float with RESULT_TAG equal to 1 */
       int i,numt,bufid,mtid,ptid,lstid,rstid;
       float wr,wi;
       float Rr[MAX_LENGTH],Ri[MAX_LENGTH];
       float Xr[MAX_LENGTH];
       float Xi[MAX_LENGTH];
       float T1[MAX_LENGTH],T2[MAX_LENGTH];
       char *fargv1[MAX_LENGTH+1],*fargv2[MAX_LENGTH+1];
       if ((mtid=pvm_mytid())<0)
          {
               printf("unexpected error occured!\n");
               exit(1);
          }
       ptid=pvm_parent();
       setlinebuf(stdout);
       for (i=0;i<(argc-1)/2;i++)
          {
               sscanf(argv[2*i+1],"%g",&Xr[i]);
               sscanf(argv[2*(i+1)],"%g",&Xi[i]);
               if (i%2 == 0)
                  {
                       fargv1[i]=argv[2*i+1];
                       fargv1[i+1]=argv[2*(i+1)];
                  }
               else
                  {
                       fargv2[i-1]=argv[2*i+1];
                       fargv2[i]=argv[2*(i+1)];
                  }
          }
       fargv1[(argc-1)/2]=NULL;
        fargv2[(argc-1)/2]=NULL;
       if (argc==5)
          {
               Rr[1]=(Xr[0]-Xr[1])/2;
               Ri[1]=(Xi[0]-Xi[1])/2;
               Rr[0]=(Xr[0]+Xr[1])/2;
               Ri[0]=(Xi[0]+Xi[1])/2;
          }
       else
          {
               numt = pvm_spawn("fft",&fargv1[0],PvmTaskDefault,"",1,
                                             &lstid);
               if (numt != 1)
                  {
                       printf("could not spawn left son\n");
                       exit(1);
                  }
               numt = pvm_spawn("fft",&fargv2[0],PvmTaskDefault,"",1,
                                              &rstid);
               if (numt != 1)
                  {
                       printf("could not spawn right son\n");
                       exit(1);
                  }
               bufid = pvm_recv(lstid, RESULT_TAG);
               pvm_upkfloat(&T1[0],(argc-1)/2,1);

               bufid = pvm_recv(rstid, RESULT_TAG);
                pvm_upkfloat(&T2[0],(argc-1)/2,1);
               for (i=0;i<(argc-1)/4;i++)
                  {
                       wr=cos(-4*PI*i/(argc-1));
                       wi=sin(-4*PI*i/(argc-1));
                       Rr[i]=(T1[i]+wr*T2[i]-wi*T2[i+(argc-1)/4])/2;
                       Ri[i]=(T1[i+(argc-1)/4]+wr*T2[i+(argc-1)/4]+wi*T2[i])/2;
                       Rr[i+(argc-1)/4]=(T1[i]-wr*T2[i]+wi*T2[i+(argc-1)/4])/2;
                       Ri[i+(argc-1)/4]=(T1[i+(argc-1)/4]-wr*
                                   T2[i+(argc-1)/4]-wi*T2[i])/2;
                  }
          }
       pvm_initsend(PvmDataDefault);
       pvm_pkfloat(&Rr[0],(argc-1)/2,1);
       pvm_pkfloat(&Ri[0],(argc-1)/2,1);
       pvm_send(ptid,RESULT_TAG); pvm_exit(); exit(0); }

3.2  Client/Serveur

On considère une architecture client/serveur toute simple: un serveur de calculs S sait comment calculer certaines fonctions (on choisira ce que l'on veut) prenant un seul argument. Des clients, C1,¼,Cn , veulent avoir accès au calcul de ces fonctions à la demande. Clients et serveur seront chacun des processus PVM. Pour simplifier, on supposera que l'on spawn le serveur au début, et que c'est lui qui crée ses clients. Quel est le problème avec le recv bloquant de PVM? On pourra améliorer le programme avec un recv non-bloquant, décrit en annexe.

Remarque:

Cela peut servir entre autres à faire de l'équilibrage de charge ``centralisé'' où les clients demandent sitôt qu'il ont fini un calcul sur un sous-domaine d'autres calculs à effectuer au serveur.

Annexe:

Réception non-bloquante:

int bufid;
bufid=pvm_nrecv(int tid,int msgtag);
On peut vérifier si quelque chose est arrivé ou pas,
int info;
info=pvm_bufinfo(int bufid,int *bytes,int *msgtag,int *tid);
#include <stdio.h>
#include <pvm3.h>
#define MSG 1

main()
{
   int mtid,ctid,bufid,ok,i;
   if ((mtid=pvm_mytid())<0) { pvm_exit();
                               exit(1); }

   pvm_spawn("async2",(char **) 0,PvmTaskDefault,
                             (char *) 0,1,&ctid);
   fprintf(stdout,"I am spawning one child\n");
   fflush(stdout);

   ok = 0;
   i = -1;
   while (!ok)
      { i++;
        bufid = pvm_nrecv(ctid,MSG);
        if (bufid < 0) { pvm_perror("Pb with nrecv");
                         pvm_exit(); exit(1); }
        if (bufid != 0) ok=1;
        else if ((i % 60000) == 0)
                { fprintf(stdout,"I am computing...\n");
                  fflush(stdout);
                } }

   pvm_upkint(&i,1,1);
   fprintf(stdout,"I have received %d\n",i);
   fflush(stdout);
   pvm_exit();
   exit(1);
}
Esclave ``async2'':
#include <stdio.h>
#include <pvm3.h>
#define MSG 1

main()
{
   int i=42;
   sleep(4);
   fprintf(stdout,"I'm sending a message to my parent\n");
   fflush(stdout);
   pvm_initsend(PvmDataDefault);
   pvm_pkint(&i,1,1);
   pvm_send(pvm_parent(),MSG);
   pvm_exit();
   exit(1);
}
pvm> spawn -1 -> async1
[11]
1 successful
t40021
pvm> [11:t40021] I am spawning one child
[11:t40021] I am computing...
[11:t40021] I am computing...
[11:t40021] I am computing...
[11:t40021] I am computing...
[11:t40021] I am computing...
[11:t40022] I am sending a message to my parent
[11:t40022] EOF
[11:t40021] I have received 42
[11:t40021] EOF
[11] finished

TD 4 
Résolution de systèmes linéaires (C/PVM)

Au cours d'une résolution numérique de l'équation des ondes par un schéma implicite, on a à trouver la solution d'un système d'équations linéaires Ax = b . Après un algorithme de type pivot de Gauss que l'on n'étudiera pas ici, on se retrouve avec un système triangulaire supérieur Ux = b . On se propose ici de programmer la résolution parallèle de ce dernier système.

Soit donc à résoudre le système triangulaire supérieur
Ux = b
avec,
U =  æ
ç
ç
ç
ç
è
U1,1
U1,2
U1,3
¼
U1,n
U2,2
U2,3
¼
U2,n
¼
¼
Un,n
ö
÷
÷
÷
÷
ø

et Ui,i ¹ 0 pour tout 1 £ i £ n .

On procède par ``remontée'' c'est à dire que l'on calcule successivement,
 
 
xn
bn
Un,n
xi
bi-
n
å
j = i+1
Ui,jxj

Ui,i

pour i = n-1,n-2,¼,1 .

On écrit les équations précédentes en pseudo-code C de la façon suivante:

x[n]=b[n]/U[n,n];
for (i=n-1;i>=1;i--)
   {  
      x[i]=0;
      for (j=i+1;j<=n;j++)
         L: x[i]=x[i]+U[i,j]*x[j];
      x[i]=b[i]-x[i]/U[i,i];
   }
n est la taille de la matrice U représentée par le tableau U[i,j], où x[i] (i entre 1 et n) représente le vecteur ``inconnu'' x et b[i], le vecteur second membre du système d'équations linéaires.

Programmer le code de remontée correspondant en PVM. On supposera que Pn est le maître qui contient au début tout le tableau U et le tableau b et crée ses esclaves P1 , ¼, Pn-1 . On pourra comme dans le cas séquentiel commencer par calculer x[n]. Ensuite, l'idée sera de placer chacune des composantes des tableaux rentrant en compte pour le calcul de la ième équation linéaire ( 1 £ i £ n-1 )1 sur le ième processeur Pi . On pourra donc propager les valeurs déjà résolues des x[j], j ³ i+1 aux processeurs Pk , k £ i, qui effectueront une même étape de sommation élémentaire. A la fin il doit attendre le résultat de tous ses esclaves.


Footnotes:

1C'est celle qui sert à déterminer x[i] à partir des x[j] (j ³ i+1) lors de la remontée.


File translated from TEX by TTH, version 2.60.
On 11 Feb 2000, 15:48.