mpi mnozenje matrica

mpi mnozenje matrica 22.01.2014. u 14:38 · pre 155 meseci #1

Radim neki mali projekat iz mpi-a i imam nekih problema. Pa ako mozete pomoci.
Kod je dolje ispod. Radi se o mnozenju matrica.
Ovo je tekst projekta

Pojašnjenje za MPI: Samo jedan čvor (master, npr. broj 0) ima matrice A i B (jer se one unose samo na jednom mjestu, to što se popunjavaju konstantama je stavljeno da se pojednostavi kod za vas, u stvarnom slučaju se uzimaju iz fajlova koji ne moraju biti dostupni svim čvorovima, itd) koje treba da dostavi (cijele ili dijelove) svim čvorovima koji množe matrice, a nakon množenja treba da na masteru rekonstruiše matricu C na osnovu parcijalnih rješenja sa svih čvorova. Kad se matrica C kompletira, saberu se svi elementi i dobije suma (ovo sabiranje je samo za provjeru, nas interesuje cijela matrica C).

Kad kompajliram kod kaze da je greska u scatterv. I ako mozete da pojasnite kako radi scatterv? Ja sam to pokusao nekako ali ne ide.

Code:

/* File taken from http://www.hp-see.eu/hp-see-computing-challenge */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <sys/time.h>
#include "mpi.h"

void ispis(long N, double *x)
{
   long i, j;

   for(i = 0; i < N; i ++)
   {
      for(j = 0; j < N; j ++){
         printf(" %lf",x[i * N + j]);

      }
         printf("\n");
   }
}

double gettime(void) {
   struct timeval tv;
   gettimeofday(&tv, NULL);
   return tv.tv_sec + 1e-6 * tv.tv_usec;
}

void matfill(long N, double *mat, double val) {
   long i, j;

   for(i = 0; i < N; i ++)
      for(j = 0; j < N; j ++)
         mat[i * N + j] = val;
}

void matmul(long N, double *a, double *b, double *c) {
   long i, j, k;
   int brojcvorova, suma=0;
   MPI_Comm_size(MPI_COMM_WORLD,&brojcvorova);
   int *temp=(int *) malloc(N * sizeof(int));
   int *temp2=(int *) malloc(N * sizeof(int));
   for(i=0;i<N;i++)
   {

     temp2[i]=suma+1;
     temp[i]=(N*N)/brojcvorova;
     suma+=temp[i];

  }
  temp2[0]=0;
   temp[brojcvorova-1]+=(N*N)%brojcvorova;


  MPI_Scatterv(a,temp2,temp,MPI_DOUBLE,c,N,MPI_DOUBLE,0,MPI_COMM_WORLD);
   for (i = 0; i < N; i ++)
      for (j = 0; j < N; j ++)
         for (k = 0; k < N; k ++)

            c[i * N + j] += a[i * N + k] * b[k * N + j];
}

int main(int argc, char **argv) {
   long N;
   double *A, *B, *C, t;
    int brojProcesa, ukupno,cvor;
   int  provjera=MPI_Init(&argc, &argv);

    if (provjera!=MPI_SUCCESS) {            //provjera da li je doslo do greske pri pokretanju mpi programa
        printf("Greska pri pokretanju MPI programa.\n");
        MPI_Abort(MPI_COMM_WORLD, provjera);
    }
    MPI_Comm_size(MPI_COMM_WORLD, &brojProcesa);
    MPI_Comm_rank(MPI_COMM_WORLD, &cvor);

    if(argc!=2)
    {
      if(cvor==0) printf("Netacan broj parametara!!\n");
      MPI_Finalize();
      return 1;
    }

    printf("Ukupan broj procesa %d Moj broj %d\n", brojProcesa,cvor);
    MPI_Barrier(MPI_COMM_WORLD);     //cekamo da svi dodju na isto mjesto


        N = atoi(argv[1]);
    double *bufer=(double *) malloc(N * N * sizeof(double));
        A = (double *) malloc(N * N * sizeof(double));
        B = (double *) malloc(N * N * sizeof(double));
        C = (double *) malloc(N * N * sizeof(double));
    if(cvor==0)
    {
        matfill(N, A, 1.0);
        matfill(N, B, 2.0);
        matfill(N, C, 0.0);
    }
    else { sleep(5); }
    MPI_Barrier(MPI_COMM_WORLD);

    MPI_Bcast(B,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
    MPI_Bcast(bufer,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);

   t = gettime();
   matmul(N, A, B, bufer);
   t = gettime() - t;
    MPI_Gather(bufer,N,MPI_DOUBLE,C,N,MPI_DOUBLE,0,MPI_COMM_WORLD);
    ispis(N,C);
   //fprintf(stdout, "%ld\t%le\t%le\n", N, t, (2 * N - 1) * N * N / t);
   //fflush(stdout);

   free(A);
   free(B);
   free(C);
  MPI_Finalize();
   return EXIT_SUCCESS;
}

Re: mpi mnozenje matrica 22.01.2014. u 18:44 · pre 155 meseci #2

MPI_Scatterv ima za zadatak da rasporedi elemente iz datog buffer-a po raspolozivim procesima, ali to pretpostavljam da vec znas. Parametre koje ova funkcija prima imas u dokumentaciji, a evo kako taj opis izgleda:

Citat:

sendbuf - Address of send buffer (choice, significant only at root).
sendcounts - Integer array (of length group size) specifying the number of elements to send to each processor.
displs - Integer array (of length group size). Entry i specifies the displacement (relative to sendbuf) from which to take the outgoing data to process i.
sendtype - Datatype of send buffer elements (handle).
recvbuf - Address of receive buffer (choice).
recvcount - Number of elements in receive buffer (integer).
recvtype - Datatype of receive buffer elements (handle).
root - Rank of sending process (integer).
comm - Communicator (handle).

Bacila sam pogled na tvoj kod i prvo sto sam uocila je sledece: ti za parametar sendcounts saljes niz temp2 koji ima duzinu N gde je N velicina matrice koju dobijas iz argumenata. Prva greska je to sto njegova duzina mora biti jednaka broju procesa u komunikatoru kom se obracas (naveden u parametru comm, u tvom slucaju MPI_COMM_WORLD). Trebao si zapravo napraviti niz koji je duzine brojcvorova. Isto vazi i za drugi parametar displs koji takodje prima niz koji bi trebao biti duzine brojcvorova.

Dalje, ne razumem zasto tako punis pomenuta dva niza. Verovatno si mislio da nesto drugo predstavljaju, tako da je najbolje da se sam pojasnis.

Sto se tice toga sto ti se kod ni ne kompajlira, mozes li prepisati ovde gresku koju ti izbacuje?

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 23.01.2014. u 01:10 · pre 155 meseci #3

E to me i bunilo. Znam sta radi scatterv ali nisam znao sta su ta dva druga niza sendcounts i displs. Koji su to nizovi? Ja sam shvatio da su to neki nizovi na koje treba podijeliti matricu.
Mozes li ti tacno pojasniti sta su ta dva parametra? Pa da ja pokusam dovrsiti kod. A gresku koju izbaci samo napise da nesto nije u redu kod scatter.
Hvala.

Re: mpi mnozenje matrica 23.01.2014. u 17:09 · pre 155 meseci #4

Citat:

Ja sam shvatio da su to neki nizovi na koje treba podijeliti matricu.

Evo sta konkretno predstavljaju pomenuti nizovi:
sendcounts - za svaki proces definise koliko elemenata treba da dobije. Dakle njegova duzina je jednaka broju procesa (cvorova), a i-ti element niza govori koliko elemenata ce proces ranka i dobiti.
displs - niz duzine koja je takodje jednaka broju procesa. i-ti element ovog niza govori od kog indeksa ce proces ranka i uzeti elemente.

Mislim da ces najlakse shvatiti iz primera. Recimo da imas 4 procesa, i niz [5, 2, 7, 9, 8, 3, 4, 6, 1] kog treba da podelis. Posto imas 4 procesa (cvora) nizovi sendcounts i displs treba da budu duzine 4. Definisimo ih npr na sledeci nacin:

sendcounts = [2, 3, 3, 1]
displs = [0, 2, 5, 8]

Sad da vidimo sta ce koji proces dobiti. Pocnimo od procesa ranka 0. Za njega je u sendcounts definisan broj 2, a u displs mu je definisan broj 0. To znaci da ce on dobiti dva elementa koja pocinju od nultog indeksa, tj dobice niz [5, 2]. Proces ranka 1 za sendcounts ima definisan broj 3, a za displs broj 2. To znaci da ce dobiti tri elementa, koja ce "uzeti" pocev od indeksa 2, te ce ovaj proces dobiti [7, 9, 8]. Analogno, proces ranka 2 ce dobiti [3, 4, 6], a proces ranka 3 dobija [1].

Iskrena da budem, ovu verziju scatter funkcije nisam koristila, ali sam tako shvatila iz dokumentacije koju sam ti navela u prethodnom postu.

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 23.01.2014. u 18:30 · pre 155 meseci #5

Tako sam nesto i pretpostavljao da su to sendcounts i displs samo nikako da napravim kako treba.
Sad sam prepravio kod. Ispise mi matricu ali ima jos nekih problema. Matrica nije dobro sabrana.
Kad je broj procesa (kod mene cvoroci kako sam pisao) veci od dimenzije matrice nece da zavrsi program. Moram ga prekinuti.
A kada je broj procesa veci od dimenzije matrice onda zavrsi sa greskom.

Code:

/* File taken from http://www.hp-see.eu/hp-see-computing-challenge */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <sys/time.h>
#include "mpi.h"

void ispis(long N, double *x)
{
   long i, j;

   for(i = 0; i < N; i ++)
   {
      for(j = 0; j < N; j ++){
         printf(" %lf",x[i * N + j]);

      }
         printf("\n");
   }
}

double gettime(void) {
   struct timeval tv;
   gettimeofday(&tv, NULL);
   return tv.tv_sec + 1e-6 * tv.tv_usec;
}

void matfill(long N, double *mat, double val) {
   long i, j;

   for(i = 0; i < N; i ++)
      for(j = 0; j < N; j ++)
         mat[i * N + j] = val;
}

void matmul(long N, double *a, double *b, double *c) {
   long i, j, k;
   int brojcvorova, suma=0;
   MPI_Comm_size(MPI_COMM_WORLD,&brojcvorova);
   int *sendcounts=(int *) malloc(brojcvorova * sizeof(int));
   int *displs=(int *) malloc(brojcvorova * sizeof(int));

   for(i=0;i<brojcvorova-1;i++)
   {

     sendcounts[i]=(N*N)/brojcvorova;
     displs[i]=suma;
     suma+=(N*N)/brojcvorova;

  }
  sendcounts[brojcvorova-1]=(N*N)/brojcvorova+(N*N)%brojcvorova;
  displs[brojcvorova-1]=suma;

  MPI_Scatterv(a,sendcounts,displs,MPI_DOUBLE,c,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
   for (i = 0; i < N; i ++)
      for (j = 0; j < N; j ++)
         for (k = 0; k < N; k ++)
            c[i * N + j] += a[i * N + k] * b[k * N + j];
}

int main(int argc, char **argv) {
   long N;
   double *A, *B, *C, t;
    int brojProcesa, ukupno,cvor;
   int  provjera=MPI_Init(&argc, &argv);

    if (provjera!=MPI_SUCCESS) {            //provjera da li je doslo do greske pri pokretanju mpi programa
        printf("Greska pri pokretanju MPI programa.\n");
        MPI_Abort(MPI_COMM_WORLD, provjera);
    }
    MPI_Comm_size(MPI_COMM_WORLD, &brojProcesa);
    MPI_Comm_rank(MPI_COMM_WORLD, &cvor);

    if(argc!=2)
    {
      if(cvor==0) printf("Netacan broj parametara!!\n");
      MPI_Finalize();
      return 1;
    }

    printf("Ukupan broj procesa %d, Moj broj %d\n", brojProcesa,cvor);
    MPI_Barrier(MPI_COMM_WORLD);     //cekamo da svi dodju na isto mjesto

        N = atoi(argv[1]);
    double *bufer=(double *) malloc(N * N * sizeof(double));
        A = (double *) malloc(N * N * sizeof(double));
        B = (double *) malloc(N * N * sizeof(double));
        C = (double *) malloc(N * N * sizeof(double));
    if(cvor==0)
    {
        matfill(N, A, 1.0);
        matfill(N, B, 2.0);
       // matfill(N, C, 0.0);
    matfill(N,bufer,0.0);
    }
    else { sleep(5); }
    MPI_Barrier(MPI_COMM_WORLD);

    MPI_Bcast(B,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
    MPI_Bcast(bufer,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);

    MPI_Barrier(MPI_COMM_WORLD);

   t = gettime();
   matmul(N, A, B, bufer);
   t = gettime() - t;

    MPI_Gather(bufer,N*N,MPI_DOUBLE,C,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
    MPI_Barrier(MPI_COMM_WORLD);
    if (cvor==0)
    {
    ispis(N,C); }
    else { sleep(5);}
   //fprintf(stdout, "%ld\t%le\t%le\n", N, t, (2 * N - 1) * N * N / t);
   //fflush(stdout);

   free(A);
   free(B);
   free(C);
  MPI_Finalize();
   return EXIT_SUCCESS;
}

Re: mpi mnozenje matrica 24.01.2014. u 17:31 · pre 155 meseci #6

Malo sam sad bolje pogledala kod, i primetila da najverovatnije nisi bas najbolje shvatio celu sustinu MPI komunikacije, kao ni kako paralelizovati mnozenje matrica :) Evo par gresaka koje sam na brzinu primetila:

- Ovim:

Code:
MPI_Bcast(B,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);

si svima poslao celu matricu B, sto je recimo OK, medjutim ovim:

Code:
MPI_Bcast(bufer,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);

nisi uradio nista. Svima si poslao prazan buffer u kog posle smestas rezultat.. Zasto? U tom trenutku u buffer-u nemas nijedan koristan podatak koji bi podelio drugim procesima.

- Ne vidim zasto si u funkciji matmul tako krenuo da delis matricu, a posle toga ne shvatam zasto si je posle onako mnozio i gde je tu paralelizacija. Ali hajde opet da pokazem na primeru. Recimo da je ovo matrica a:

i da imas 3 procesa. Nakon ove linije koda:

Code:
MPI_Scatterv(a,sendcounts,displs,MPI_DOUBLE,c,N*N,MPI_DOUBLE,0,MPI_COMM_WORLD);

promenljiva c ce imati sledece vrednosti u datim procesima:
proces 0 - 1 5 8 3 3
proces 1 - 4 1 7 1 6
proces 2 - 2 4 3 9 7 7

Nakon toga ti navodis bukvalno serijski algoritam mnozenja matrica, pri cemu tu istu promenljivu c uopste ne koristis u proracunu, vec je koristis za skladistenje rezultata, cime prebries sve sto si dobio iz MPI_Scatterv. Ali cak i da koristis promenljivu c u proracunu, ja zaista ne znam sta bi radio svaki proces sa tako polovicnim podacima koji nisu na smislen nacin podeljeni. Dalje, za proracun koristis matricu a, koja uopste nije popunjena, i matricu b, koja je OK jer si je poslao svima onim broadcast-om. Uvidjas da je greska sto je matrica a nepopunjena, medjutim, cak i da je popunjena, te tri ugnjezdene petlje ne predstavljaju nikakvu paralelizaciju buduci da svi procesi rade dupli posao mnozenja matrica.

Dakle, vidim da nisi shvatio osnovne koncepte, i zato bih ti toplo preporucila da procitas poglavlje "Distributed-Memory Programming with MPI" knjige "An Introduction to Parallel Programming" (Peter S. Pacheco). Ako to ne zelis, onda ti savetujem da uzmes nekog kolegu koji ce ti to objasniti. Ovako, preko foruma, vec postaje dosta neprakticno s obzirom da se ne radi o sitnoj greski, vec o celom konceptu paralelnog programiranja. :)

_{[Ovu poruku je menjao tuolarips dana 24.01.2014. u 19:24 GMT+1]}

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 25.01.2014. u 01:56 · pre 155 meseci #7

Ja sam dobio gotov kod koji trebam prepraviti tako da koristeci MPI uradim paralelizaciju mnozenja matrica.
Isao sam na ovaj nacin.
Kreiram matricu A i popunim je 1. Kreiram matricu B i popunim je 2. Takodje i matricu C sa 0.
Matrica se mnozi tako sto vrstu matrice A mnozim kolonom matrice B.
Zato sam mislio ovako. Svakom procesu dati jednu vrstu matrice A npr. prvu i citavu matricu B. Kada to pomnozim dobicu prvu vrstu matrice C. I tako redom. Kad to sve pomnozim svaki proces ce imati neku od vrsta matrice C.
I na kraju sa gather te sve vrste sastaviti u matricu C.
Nadam se da si razumjela na koji nacin sam htio.
Sto se tice onog bufera to sam stavio samo kao neki pomocni niz. Mada sam mogao odmah svima podijeliti i matricu C.
Problem je sad sto se ja prvi put susrecem sa MPI pa dok ne shvatim na koji nacin to sve ide.
Stavicu ovdje pocetni kod koji sam dobio. Sigurno da ti imas neki bolji nacin kako ovo uraditi. Pa ako imas neki laksi i prakticniji nacin reci, pa da pokusam.

Code:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <sys/time.h>

double gettime(void) {
   struct timeval tv;
   gettimeofday(&tv, NULL);
   return tv.tv_sec + 1e-6 * tv.tv_usec;
}

void matfill(long N, double *mat, double val) {
   long i, j;

   for(i = 0; i < N; i ++)
      for(j = 0; j < N; j ++)
         mat[i * N + j] = val;
}

void matmul(long N, double *a, double *b, double *c) {
   long i, j, k;

   for (i = 0; i < N; i ++)
      for (j = 0; j < N; j ++)
         for (k = 0; k < N; k ++)
            c[i * N + j] += a[i * N + k] * b[k * N + j];
}

int main(int argc, char **argv) {
   long N;
   double *A, *B, *C, t;

   N = atoi(argv[1]);
   A = (double *) malloc(N * N * sizeof(double));
   B = (double *) malloc(N * N * sizeof(double));
   C = (double *) malloc(N * N * sizeof(double));
   matfill(N, A, 1.0);
   matfill(N, B, 2.0);
   matfill(N, C, 0.0);

   t = gettime();
   matmul(N, A, B, C);
   t = gettime() - t;

   fprintf(stdout, "%ld\t%le\t%le\n", N, t, (2 * N - 1) * N * N / t);
   fflush(stdout);

   free(A);
   free(B);
   free(C);

   return EXIT_SUCCESS;
}

Re: mpi mnozenje matrica 25.01.2014. u 09:47 · pre 155 meseci #8

Imas OK ideju, samo sto dosta stvari nisi uzeo u obzir.. Taj problem je u paralelnom programiranju vec dosta obradjen, pa je meni sad, u nedostatku vremena, najlakse da te uputim na resenje. Malo naprednija od one sto sam ti vec dala je kniga slicnog naziva "Introduction to Parallel Computing, Second Edition" (Ananth Grama, Anshul Gupta, George Karypis, Vipin Kumar). U njoj imas deo "8.2 Matrix-Matrix Multiplication" koji ti opisuje resenje problema. Nadam se da ce ti pomoci.

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 25.01.2014. u 15:33 · pre 155 meseci #9

Nasao sam ove knjige i citam pomalo. Ali ovdje nigdje ne pise kako da uradim to. Objasnjeno je kako mogu da to uradim, na koje nacine. Da ima negdje neki primjer slican da mogu iz tog primjera da shvatim kako da uradim.
Naletio sam na neke primjere jos mnozenja ali su to sve slucajevi kada je matrica kreirana kao dvodimenzionalno polje i gdje je broj procesa i dimenzija "odgovaraju" tako da se moze to sve lako izracunati.

Da li bi mogla samo da napises kako da podijelim matricu A svim procesima tako da ne bude problema kad je broj procesa veci ili manji od dimenzije matrice??
Po meni ovo mnozenje matrica kako je napisano u matmul funkciji nije nekako prakticno. A mozda i grijesim jer nemam pojma.

Ako ti nije problem da samo dodas u ovaj gotov kod kako izdijeliti matricu A svim procesima.

Meni je rok do 30. januara da to uradim.

I trebao bih uraditi takodje preko open MP. Koji je laksi nacin?

Re: mpi mnozenje matrica 29.01.2014. u 22:58 · pre 155 meseci #10

Kao sto rekoh, imam jako malo slobodnog vremena, ali evo danas uspeh da ti "sklopim" nesto. Ispod ti je kod, moguce da ti neki momenti nece biti najjasniji, pa pitaj ako to bude slucaj, a ja cu ti odgovoriti kad ugrabim vremena. Ideja ovog algoritma je da svaki proces dobije odredjeni deo elemenata da proracuna. Uzmimo za primer mnozenje kvadratnih matrica dimenzije 4 u kom ucestvuje 5 procesa. Tada bi procesi dobili da racunaju sledece elemente:

proces 0: (0, 0) (0, 1) (0, 2)
proces 1: (0, 3) (1, 0) (1, 1)
proces 2: (1, 2) (1, 3) (2, 0)
proces 3: (2, 1) (2, 2) (2, 3)
proces 4: (3, 0) (3, 1) (3, 2) (3, 3)

Sto se tice podele matrica, ideja je da svaki proces dobije samo one elemente redova prve matrice i kolona druge matrice koji su mu neophodni za racun elemenata za koje je zaduzen. Zbog toga ce raspodela redova i kolona medju procesima biti ovakva:

proces 0: redovi-0, kolone-0,1,2
proces 1: redovi-0,1 kolone-0,1,2,3 (napomena: kolona 2 mu zapravo ne treba, ali je saljemo zbog ustede na komunikaciji. Naime, MPI ne podrzava da "preskocis" nesto u buffer-u kao sto bi nam u ovom slucaju bilo zgodno jer bismo trebali da saljemo elemente nulte i prve kolone, pa onda da preskocimo drugu i posaljemo trecu. Posto to nije podrzano, za ovako nesto bismo morali da razmenimo dve poruke sto generalno treba izbegavati jer je komunikacija medju procesima ono sto najvise usporava paralelne programe)
proces 2: redovi-1,2 kolone-0,1,2,3 (napomena: kolona 1 je nepotrebna, ali se salje iz istog razloga kog sam navela iznad)
proces 3: redovi-2 kolone-1,2,3
proces 4: redovi-3 kolone-0,1,2,3

Sto se tice poredjenja MPI-a i OpenMP-a, ne bih nijedan od ova dva okarakterisala laksim ili tezim. Imaju drugaciju svrhu, i u skladu s potrebama treba se odluciti za jedan od njih.

Code:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "mpi.h"

void matprint(char *matname, double *mat, long N) {
    long i, j;
    fprintf(stdout, "\n------------%s------------", matname);

    for (i = 0; i < N; i++) {
        fprintf(stdout, "\n");
        for (j = 0; j < N; j++)
            fprintf(stdout, "%f ", mat[i * N + j]);
    }

    fprintf(stdout, "\n--------------------------\n");
}

void matfill(long N, double *mat, double val) {
    long i, j;

    for(i = 0; i < N; i ++)
        for(j = 0; j < N; j ++) {
            mat[i * N + j] = val;
        }
}

int matmul(int commSize, int myRank, long N, double *rows, double *columns, double *result) {
    long i, j, count, index, start, end, startRow, startColumn, row, column;

    if (myRank >= N*N) {
        result[0] = 0;
        return 1;
    }

    count = (N*N)/commSize;
    if (count == 0) count = 1;

    start = myRank * count;
    end = start + count;

    if (myRank==commSize-1)
        end += (N*N)%commSize;

    startRow = start / N;

    if (start / N == (end - 1) / N) {
        startColumn = start % N;
    } else {
        startColumn = 0;
    }

    index = 0;
    for (i = start; i < end; i ++) {
        result[index] = 0;
        row = i / N - startRow;
        column = i % N - startColumn;
        for (j = 0; j < N; j ++) {
            result[index] += rows[row*N+j] * columns[column*N+j];
        }
        index++;
    }

    return end - start;
}

int dividerows(int commSize, int myRank, long N, double *mat, double *rows) {
    int i, start, end, startRow, endRow, count, rowsForMyRank;
    int *scounts = (int *)malloc(commSize*sizeof(int));
    int *displs = (int *)malloc(commSize*sizeof(int));

    count = (N*N)/commSize;
    if (count == 0) count = 1;

    for (i=0; i<commSize; i++) {
        if (i < N*N)
        {
            start = i * count;
            end = start + count - 1;

            if (i==commSize-1)
                end += (N*N)%commSize;

            startRow = start / N;
            endRow = end / N;

            scounts[i] = (endRow - startRow + 1) * N;
            displs[i] = startRow * N;
        } else {
            scounts[i] = 1;
            displs[i] = 0;
        }
    }

    MPI_Scatterv(mat, scounts, displs, MPI_DOUBLE, rows, N*N, MPI_DOUBLE, 0, MPI_COMM_WORLD);

    rowsForMyRank = scounts[myRank];

    free(scounts);
    free(displs);

    return rowsForMyRank;
}

int dividecolumns(int commSize, int myRank, long N, double *mat, double *columns) {
    int i, start, end, startColumn, endColumn, count, columnsForMyRank;
    int *scounts = (int *)malloc(commSize*sizeof(int));
    int *displs = (int *)malloc(commSize*sizeof(int));
    MPI_Datatype col, coltype;

    MPI_Type_vector(N, 1, N, MPI_DOUBLE, &col);
    MPI_Type_commit(&col);
    MPI_Type_create_resized(col, 0, 1*sizeof(double), &coltype);
    MPI_Type_commit(&coltype);

    count = (N*N)/commSize;
    if (count == 0) count = 1;
    for (i=0; i<commSize; i++) {
        if (i < N*N) {
            start = i * count;
            end = start + count - 1;

            if (i==commSize-1)
                end += (N*N)%commSize;

            if (start / N == end / N) {
                startColumn = start % N;
                endColumn = end % N;
            } else {
                startColumn = 0;
                endColumn = N-1;
            }

            scounts[i] = endColumn - startColumn + 1;
            displs[i] = startColumn;

            //fprintf(stdout, "\n\n----------- Process %d. start %d. end %d. startRow %d. endRow %d. scounts %d. displs %d. \n", myRank, start, end, startColumn, endColumn, scounts[i], displs[i]);
        } else {
            scounts[i] = 1;
            displs[i] = 0;
        }
    }

    MPI_Scatterv(mat, scounts, displs, coltype, columns, N*N, MPI_DOUBLE, 0, MPI_COMM_WORLD);

    columnsForMyRank = scounts[myRank] * N;

    free(scounts);
    free(displs);

    return columnsForMyRank;
}

int main(int argc, char **argv) {
    long N;
    int i, commSize, myRank, rowsSize, columnsSize, resultCount, currentDispls, count;
    int *rcounts, *displs;
    double *A, *B, *C, t, *rows, *columns, *localResult;

    int mpiInitSucceeded = MPI_Init(&argc, &argv);

    if (mpiInitSucceeded != MPI_SUCCESS) {
        printf("MPI failed to initialize.\n");
        MPI_Abort(MPI_COMM_WORLD, mpiInitSucceeded);
    }

    MPI_Comm_size(MPI_COMM_WORLD, &commSize);
    MPI_Comm_rank(MPI_COMM_WORLD, &myRank);

    N = 0;
    if (myRank == 0)
        N = atoi(argv[1]);

    MPI_Bcast(&N, 1, MPI_LONG, 0, MPI_COMM_WORLD);

    A = (double *) malloc(N * N * sizeof(double));
    B = (double *) malloc(N * N * sizeof(double));
    C = (double *) malloc((N * N < commSize ? commSize : N * N) * sizeof(double));

    if (myRank == 0) {
        matfill(N, A, 1.0);
        matfill(N, B, 2.0);
        matfill(N, C, 0.0);

        //t = MPI_Wtime();
        //matmul(N, A, B, C);
        //t = MPI_Wtime() - t;

        //fprintf(stdout, "%ld\t%le\t%le\n", N, t, (2 * N - 1) * N * N / t);
        //fflush(stdout);
    }

    rows = (double *)malloc(N*N*sizeof(double));
    rowsSize = dividerows(commSize, myRank, N, A, rows);

    columns = (double *)malloc(N*N*sizeof(double));
    columnsSize = dividecolumns(commSize, myRank, N, B, columns);

    /*if (myRank == 1) {
        fprintf(stdout, "\n\n----------- Process 1 rows(%d): \n", rowsSize);
        for (i=0; i<rowsSize; i++) {
            fprintf(stdout, " %f ", rows[i]);
        }
        fprintf(stdout, "\n-----------\n");
    }*/

    localResult = (double *)malloc(N*N*sizeof(double));
    resultCount = matmul(commSize, myRank, N, rows, columns, localResult);

    rcounts = (int *)malloc(commSize*sizeof(int));
    displs = (int *)malloc(commSize*sizeof(int));

    count = (N*N)/commSize;
    if (count == 0) count = 1;

    currentDispls = 0;
    for (i=0; i<commSize; i++) {
        displs[i] = currentDispls;
        if (i < N*N - 1) {
            rcounts[i] = count;
            currentDispls += count;
        } else if (i == N*N - 1){
            rcounts[i] = count + (N*N)%commSize;
            currentDispls += count + (N*N)%commSize;
        } else {
            rcounts[i] = 1;
            currentDispls++;
        }
    }

    MPI_Gatherv(localResult, resultCount, MPI_DOUBLE, C, rcounts, displs, MPI_DOUBLE, 0, MPI_COMM_WORLD);

    if (myRank == 0) matprint("Resulting matrix", C, N);

    free(A);
    free(B);
    free(C);
    free(rows);
    free(columns);
    free(localResult);

    MPI_Finalize();

    return EXIT_SUCCESS;
}

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 31.01.2014. u 14:24 · pre 154 meseci #11

Ti si ovo bas momacki odradila. Svaka ti cast na ovome i hvala puno.
Pregledao sam kod. Uspio sam skontati kako si radila i sta. Trebalo mi je malo vise vremena ali sam se izborio na kraju.
Hvala jos jednom.

Re: mpi mnozenje matrica 01.02.2014. u 11:30 · pre 154 meseci #12

Ako koristis gcc openmp je trivijalno.
Samo stavis

Code:

#pragma omp parallel for
for(...)....

I iskompajliras sa gcc -fopenmp ....

Re: mpi mnozenje matrica 01.02.2014. u 11:53 · pre 154 meseci #13

S tim sto bih napomenula da se promenljive j i k moraju proglasiti privatnim, inace ce doci do pogresnog (nedeterministickog) rezultata.

I da, u ovom slucaju se ispostavilo da je OpenMP resenje dosta jednostavnije. Nisam ni razmislila sta bi bilo OpenMP resenje ovog problema, pa sam na brzinu rekla ono sto generalno mislim, a to je da se ne moze reci da li je laksi OpenMP ili MPI. Deljena memorija moze da bude dosta problematicna jer se cesto neke stvari previde.

Potpis

"Time is a drug. Too much of it kills you." Terry Pratchett

Re: mpi mnozenje matrica 02.02.2014. u 10:14 · pre 154 meseci #14

Uradio sam Openmp. Bas je bilo mnogo lakse nego mpi. Sad ne znam da li je to uvijek tako ili ne.
Ali je bilo dosta lakse za shvatiti i uraditi.
I koristim gcc. Tako da nije bilo tesko.
Hvala na pomoci.

Srodne teme