+ All Categories
Home > Documents > Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1...

Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1...

Date post: 27-Dec-2019
Category:
Upload: others
View: 48 times
Download: 0 times
Share this document with a friend
77
Notiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele si distribuite ar trebui in primul rand sa ne raportam la diferentele lor fata de implementarile seriale: alocarea task-urilor; comunicarea intre procesoare a rezultatele intermediare obtinute prin calcul; sincronizarea calculelor efectuate de procesoare. De aici se obtine o clasificare a metodelor: metode sincrone si metode asincrone. Calculul paralel si distribuit reprezinta in prezent un domeniu larg de interes existand o activitate de cerc- etare intensa in acest sens motivata de numerosi factori. Intotdeauna a fost nevoie de o solutie de rezolvare a problemelor cu calcule foarte mari, dar, doar recent, tehnologiile avansate au luat in calcul folosirea calculului paralel in rezolvarea acestor probleme tinand cont si de calculatoarele paralele din ce in ce mai puternice. Dezvolatarea algoritmilor paraleli si distribuiti este determinata pe de o parte de legatura dintre nevoile de calcul noi si cele vechi, iar pe de alta parte de progresul tehnologic. Aplicatiile care privesc inteligenta artificiala, calculul simbolic si calculul numeric sunt cele care au jucat un rol important in dezvoltarea arhitecturilor paralele si distribuite. Nevoia de determinare rapida a unui numar mare de calcule numerice din cadrul aplicatiilor ce privesc rezolvarea ecuatiilor cu derivate partiale si a procesarii imaginilor a fost ceea ce a dus la dezvoltarea rapiditatii si a paralelizarii masinilor de calcul. Recent a crescut interesul pentru aplicarea calculului paralel si distribuit in aplicatii precum cele ce privesc analiza, simularea si optimizarea la scara larga a sistemelor interconectate. Alte exemple de aplicare a calcu- lului paralel sunt: rezolvarea sistemelor de ecuatii, programarea matematica si problemele de optimizare. O proprietate comuna a acestor probleme este aceea de a putea fi descompuse in subtask-uri. O distictie importanta este intre sistemele de calculul paralel si cele de calcul distribuit. Sistemele de calcul paralel consistau in mai multe procesoare care sunt localizate la o distanta mica unele de altele si care au fost concepute in scopul realizarii de taskuri de calcul in comun. Comunicarea dintre procesoare este de incredere si poate fi anticipata. Sistemele de calcul distribuit sunt diferite din mau multe puncte de vedere: procesoarele pot fi separate si comunicarea dintre ele este mult mai problematica; intarzierile de comunicare nu pot fi anticipate si legaturile de comunicatie pot sa nu fie de incredere; fiecare procesor poate fi angajat in propriile sale activitati in timp ce coopereaza cu alte procesoare pentru realizarea unui task comun. Sunt mai multi parametrii care pot fi folositi pentru a descrie sau pentru a clasifica sistemele de calcul paralel: tipul si numarul de procesoare. Sunt sisteme de calcul paralel care au multe procesoare numite masiv paralele si sunt sisteme cu un numar mic de procesoare numite coarse grained parallelism. prezenta sau absenta unei masini de control global. O clasificare privind acest aspect si referindu-ne la abilitatea diferitelor procesoare de a executa diferite instructiuni la un moment dat este in calculatoare paralele SIMD (Single Instruction Multiple Data) si MIMD (Multiple Instruction Multiple Data). 1
Transcript
Page 1: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Notiuni introductive

1 Arhitecturi paralele si distribuitePentru a realiza un studiu al metodelor numerice paralele si distribuite ar trebui in primul rand sa ne raportamla diferentele lor fata de implementarile seriale:

• alocarea task-urilor;

• comunicarea intre procesoare a rezultatele intermediare obtinute prin calcul;

• sincronizarea calculelor efectuate de procesoare. De aici se obtine o clasificare a metodelor: metodesincrone si metode asincrone.

Calculul paralel si distribuit reprezinta in prezent un domeniu larg de interes existand o activitate de cerc-etare intensa in acest sens motivata de numerosi factori. Intotdeauna a fost nevoie de o solutie de rezolvare aproblemelor cu calcule foarte mari, dar, doar recent, tehnologiile avansate au luat in calcul folosirea calcululuiparalel in rezolvarea acestor probleme tinand cont si de calculatoarele paralele din ce in ce mai puternice.

Dezvolatarea algoritmilor paraleli si distribuiti este determinata pe de o parte de legatura dintre nevoilede calcul noi si cele vechi, iar pe de alta parte de progresul tehnologic.

Aplicatiile care privesc inteligenta artificiala, calculul simbolic si calculul numeric sunt cele care au jucatun rol important in dezvoltarea arhitecturilor paralele si distribuite.

Nevoia de determinare rapida a unui numar mare de calcule numerice din cadrul aplicatiilor ce privescrezolvarea ecuatiilor cu derivate partiale si a procesarii imaginilor a fost ceea ce a dus la dezvoltarea rapiditatiisi a paralelizarii masinilor de calcul.

Recent a crescut interesul pentru aplicarea calculului paralel si distribuit in aplicatii precum cele ce privescanaliza, simularea si optimizarea la scara larga a sistemelor interconectate. Alte exemple de aplicare a calcu-lului paralel sunt: rezolvarea sistemelor de ecuatii, programarea matematica si problemele de optimizare. Oproprietate comuna a acestor probleme este aceea de a putea fi descompuse in subtask-uri.

O distictie importanta este intre sistemele de calculul paralel si cele de calcul distribuit.Sistemele de calcul paralel consistau in mai multe procesoare care sunt localizate la o distanta mica unele

de altele si care au fost concepute in scopul realizarii de taskuri de calcul in comun. Comunicarea dintreprocesoare este de incredere si poate fi anticipata.

Sistemele de calcul distribuit sunt diferite din mau multe puncte de vedere: procesoarele pot fi separate sicomunicarea dintre ele este mult mai problematica; intarzierile de comunicare nu pot fi anticipate si legaturilede comunicatie pot sa nu fie de incredere; fiecare procesor poate fi angajat in propriile sale activitati in timpce coopereaza cu alte procesoare pentru realizarea unui task comun.

Sunt mai multi parametrii care pot fi folositi pentru a descrie sau pentru a clasifica sistemele de calculparalel:

• tipul si numarul de procesoare. Sunt sisteme de calcul paralel care au multe procesoare numite masivparalele si sunt sisteme cu un numar mic de procesoare numite coarse grained parallelism.

• prezenta sau absenta unei masini de control global. O clasificare privind acest aspect si referindu-ne laabilitatea diferitelor procesoare de a executa diferite instructiuni la un moment dat este in calculatoareparalele SIMD (Single Instruction Multiple Data) si MIMD (Multiple Instruction Multiple Data).

1

Page 2: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

• operatii sincrone si asincrone. Distinctia in acest caz se refera la prezenta sau absenta unui ceasglobal folosit pentru a sincroniza operatiile diferitelor procesoare. Aceasta sincronizare este intot-deauna prezenta in SIMD. Operatiile sincrone sunt de preferat datorita catorva proprietati importante:comportamentul procesoarelor este mult mai usor de controlat si algoritmul este mult simplificat. Darin unele cazuri este destul de greu de sincronizat o retea de comunicare de date.

• legaturile dintre procesoare. Un aspect important al calculatoarelor paralele este mecanismul prin careprocesoarele fac schimb de informatii. Din acest punct de vedere exista doua cazuri: arhitecturi cumemorie comuna si arhitecturi care transmit mesaje. Plecand de la aceste doua tipuri de arhitecturiavem numerosi hibrizi. Cel mai des intalnit mod de a lucra cu memoria este switching systems (sistemede comutare) asa cum se poate observa in figura 1. Exista situatii in care nu avem memorie comuna, darfiecare procesor are propria sa memorie locala. In aceasta situatie procesoarele pot comunica printr-oretea de comunicare care consta in legaturi de conexiune care leaga anumite perechi de procesoare(2). Apare astfel problema de a stabili care procesoare sunt conectate. Ar fi de preferat daca toateprocesoarele ar fi direct conectate unele cu altele, sau daca pot fi legate printr-un numar de legaturiaditionale.

Figure 1: Sistem de conectare

Figure 2: Retea de comunicare

Structura retelei de conectare este foarte importanta atat in sistemele de calcul paralel cat si in cele decalcul distribuit dar cu o diferenta: in calculatoarele paralele reteaua de conectare este sub controlulunui designer si din aceasta cauza de obicei este regulata, iar in sistemele distribuite topologia reteleieste predeterminata si iregulata.

2

Page 3: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

2 Modele si exempleExista o varietate de modele de calcule paralele si distribuite. Un model simplu folosit in general in calcululparalel pentru algoritmi asincroni este acela al reprezentarii algoritmilor paraleli prin grafuri directe aciclice.Fiecare nod al grafului reprezinta operatiile pe care le face algoritmul si fiecare arc este folosit pentru areprezenta dependenta datelor. In particular, un arc (i, j) indica faptul ca operatia corespunzatoare nodului jfoloseste rezultatele operatiei corespunzatoare nodului i. Un astfel de exemplu este cel din figura 3.

Figure 3: Exemplu de graf direct orientat pentru (x1 + x2)(x2 + x3)

Un graf direct orientat realizeaza doar o reprezentare partiala a unui algoritm. In cadrul grafului suntspecificate operatiile care se efectueaza, operanzii si impune anumite conditii privind ordinea in care acesteoperatii se efectueaza. Pentru a determina un algoritm paralel trebuie sa specificam si ce procesor executa oanumita operatie si la ce moment o executa. Pentru aceasta trebuie impus ca un procesor sa execute cel multo operatie la un moment dat, si ca daca avem un nod (i, j) atunci operatia care ii corespunde nodului j poateincepe doar dupa ce operatia corespunzatoare nodului i a fost terminata.

Exista situatii cand pentru aceeasi problema de calcul exista mai multi algoritmi, deci sunt situatii incare avem mai multe grafuri direct orientate pentru rezolvarea unei probleme (figura 4). Apare astfel situatiadeterminarii unui graf ce minimizeaza timpul de calcul al solutiei problemei. Pentru cateva clase de problemeexista metode care construiesc un graf direct orientat cu un factor constant de optimizare.

In continuare consideram cateva exemple simple de taskuri de calcul numeric care pot fi reprezentate prinalgoritmi paraleli simpli. Toti algorimii considerati pot fi reprezentati printr-un graf direct orientat. Pentru ausura calculele presupunem ca fiecare operatie este executata intr-o unitate de timp si ca transmiterea de dateintre procesoare se realizeaza instantaneu.

2.1 Calcule cu vectori si matrici1. Adunarea scalarilor

Cea mai simpla sarcina de calcul este aceea de a aduna n scalari. Cel mai bun algoritm serial necesitan− 1 operatii, deci timpul de executie al algoritmului este n− 1.

In continuare prezentam un algoritm paralel construit sub presupunerea ca n este o putere a lui 2.Partitionam cei n scalari in n

2 perechi diferite si folosim n2 procesoare pentru a aduna cei doi scalari

din fiecare pereche. Deci, dupa o unitate de timp, avem de executat sarcina de a aduna n2 scalari.

Continuand in aceeasi maniera, dupa log n pasi vom mai avea de efectuat o singura adunare si calculul

3

Page 4: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 4: Graf direct orientat pentru (x1 + x2)(x2 + x3)

se va termina (figura 5). Acest algoritm poate fi generalizat la cazul in care n nu este putere a lui 2 siastfel timpul de executie devine ⌈log n⌉ folosind ⌊n

2 ⌋ procesoare.

Figure 5: Calcul paralel pentru adunarea a 16/15 scalari

Aceeasi problema a adunarii a n scalri poate fi rezolvata folosind graful direct orientat reprezentat infigura 6 si care aduce o imbunatatire a eficientei algoritmului.

2. Produsul scalar

Produsul scalarn∑

i=1

xiyi a doi vectori n-dimensionali poate fi calculata cu cu n procesoare in ⌈log n⌉+1

unitati de timp astfel: la primul pas fiecare procesor i calculeaza produsul xiyi si apoi avem nevoie de⌈log n⌉ unitati de timp pentru adunarea scalarilor folosind algoritmul anterior.

3. Adunarea si inmultirea matricilor

Suma a n matrici de dimensiune m×m poate fi calculta in ⌈log n⌉ unitati de timp si folosind m2⌊n2 ⌋

procesoare.

Inmultirea a doua matrici de dimensiune m × n si n × l consta in evaluarea a ml produse scalare devectori n-dimensionali si poate fi deci indeplinita in ⌈log n⌉+1 unitati de timp folosind nml procesoare.

4. Puterile unei matrici

4

Page 5: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 6: Algoritm alternativ pentru adunarea paralela a 16 scalari

Presupunem ca avem o matrice A de dimensiune n×n si dorim sa calculam Ak pentru k ∈ Z . Daca keste putere a lui 2, atunci task-ul poate fi indeplinit prin pasi repetitivi: calculam intai A2, apoi A2A2 =A4, etc. Dupa log k pasi se obtine Ak. Aceasta procedura presupune log k inmultiri consecutive dematrici si poate fi deci incheiata in log k(⌈log n⌉+ 1) unitati de timp folosind n3 procesoare.

O alta situatie de calcul a puterii unei matrici este reprezentata in figura 7.

Figure 7: Calcul paralel al puterilor unei matrici

2.2 Paralelizarea metodelor iterativePentru rezolvarea de sisteme de ecuatii, pentru optimizari si alte tipuri de probleme se folosesc in generalalgoritmi iterativi (metode de relaxare) de tipul:

x(t+ 1) = f(x(t)), t = 0, 1, . . . (1)

5

Page 6: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

unde fiecare x(t) este un vector n-dimensional si f : Rn → Rn.O notatie alternativa este x = f(x). Daca sirul {x(t)} generat de (1) este convergent la x∗ si daca functia

f este continua, atunci x∗ este punct fix al lui f , f(x∗) = x∗.Fie xi(t) componenta i a lui x(t) si fie fi componenta i a functiei f . Atunci putem scrie (1) ca:

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, n (2)

Algoritmul iterativ x = f(x) poate fi paralelizat prin luarea a n procesoare care sa calculeze o componentadiferita a lui x, asa cum se observa in (2). La fiecare pas, procesorul i cunoaste valorile tuturor componentelorlui x(t) de care depinde fi, calculeaza noua valoare a lui xi(t+1) si comunica aceasta valoare altor procesoarepentru a se putea trece la urmatoarea iteratie.

Comunicatiile necesare pentru executia unei iteratii de tipul (2) pot fi descrise cu ajutorul unui graf directca cel reprezentat in figura 8, numit si graf de dependenta. Nodurile reprezinta componentele lui x si pentrufiecare doua noduri i si j exista arcul (i, j), daca functia fj depinde de xi, asta inseamna ca procesorul itrebuie sa comunice valoarea lui xi(t) procesorului j.

Figure 8: Graf de dependenta

Graful din figura 8 este realizat pentru urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x2(t+ 1) = f2(x1(t), x2(t))x3(t+ 1) = f3(x2(t), x3(t), x4(t))x4(t+ 1) = f4(x2(t), x4(t))

Daca presupunem ca procesul iterativ reprezentat in ecuatia (1) se termina in T pasi, T ∈ N , structuraalgoritmului poate fi reprezentata printr-un graf direct orientat (9).

2.2.1 Iteratia Seidel-Gauss

Iteratia (2) in care toate componentele lui x sunt updatate simultan este numita iteratie de tip Jacobi. Intr-oforma alternativa, componentele lui x sunt updatate pe rand si valoarea cel mai recent calculata este folosita.Astfel ecuatia (2) devine:

xi(t+ 1) = fi(x1(t+ 1), . . . , xi−1(t+ 1), xi(t), . . . , xn(t)), i = 1, n (3)

Iteratia (3) se numeste algoritmul Seidel-Gauss bazat pe functia f .

Example 2.1 Consideram urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x2(t+ 1) = f2(x1(t+ 1), x2(t))x3(t+ 1) = f3(x2(t+ 1), x3(t), x4(t))x4(t+ 1) = f4(x2(t+ 1), x4(t))

6

Page 7: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 9: Graf direct orientat asociat grafului de dependenta

Figure 10: Paralelizarea iteratiei Seidel-Gauss

Exista mai multe alternative ale algoritmului Seidel-Gauss corespunzator functiei f , deoarece ordineain care sunt updatate componentele poate fi diferita. Atata timp cat viteza de convergenta nu se modificasemnificativ este natural sa schimbam ordinea de calcul a componentelor astfel incat paralelismul la fiecareiteratie este fie maximizat.

Example 2.2 Consideram urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x3(t+ 1) = f3(x2(t+ 1), x3(t), x4(t))x4(t+ 1) = f4(x2(t+ 1), x4(t))x2(t+ 1) = f2(x1(t+ 1), x2(t))

Proposition 2.1 Urmatoarele afirmatii sunt echivalente:

7

Page 8: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 11: Paralelizarea iteratiei Seidel-Gauss - alta ordine de updatare a componentelor

• Exista o ordine a variabilelor astfel incat o parte a algoritmului Seidel-Gauss sa poata fi executata inK pasi paraleli.

• Exista o colorare a grafului de dependenta care foloseste K culori si cu proprietatea ca nu exista unciclu pozitiv cu toate nodurile intr-un ciclu care au aceleasi culori.

• Exista o colorare a grafului de dependenta care foloseste cel mult K culori astfel incat nodurile adia-cente au culori diferite.

In cazul in care consideram exemplul din figura 8, doua culori sunt suficiente: h(1) = h(3) = h(4) = 1si h(2) = 2.

In cazul in care se folosesc schemele cu culori pentru o implementare paralela a algoritmului Seidel-Gauss, este usor de folosit cate un procesor diferit pentru fiecare componenta a lui x, deoarece fiecare procesorva fi idle atata timp cat sunt updatate variabile de culori diferite. Un remediu imediat ar fi acela de a folosimai putine procesoare, fiecare procesor avand asignate mai multe variabile de culori diferite.

Example 2.3 Intr-o varietate de algoritmi de forma x = f(x), folositi pentru rezolvarea de ecuatii cuderivate partiale sau in procesarea de imagini, fiecarei componente a vectorului x ii este asociat un punctparticular dintr-o anumita regiune dintr-un spatiu 2-dimensional.

Fie N multimea tuturor punctelor (i, j) ∈ R2 astfel incat i si j sunt numere intregi ce satisfac 0 ≤ i ≤ Msi 0 ≤ j ≤ M . Fie xij o componenta a vectorului x care corespunde punctului (i, j). Prin conectarea celormai apropiati vecini obtinem graful reprezentat in figura 12. Vedem acest graf ca si graf direct daca facemarcele bidirectionale si presupunem ca este un graf de dependenta asociat ecuatiei x = f(x). Asignam cateun procesor fiecarui punct (i, j). Acest procesor este responsabil de updatarea valorii lui xij si pentru aface acest lucru are nevoie sa cunoasca valorile componentelor lui x asociate punctelor vecine. Este astfelnatural de presupus ca procesoarele asociate punctelor vecine sunt unite printr-o conexiune directa.

In ceea ce priveste implementarea unei metode Seidel-Gauss asociate situatiei prezentate, observam cagraful din figura 12 poate fi colorat folosind doua culori asa cum se observa in figura 13. Daca asignamcate un procesor fiecarei componente xij , atunci acesta va fi idle jumatate din timp. Este deci rezonabil saasignam doua componente de culori diferite fiecarui procesor. In practica, numarul de puncte este foartemare si deci fiecare procesor va avea asignate mai mult de doua puncte.

3 Exercitii1. Se considera urmatoarea ecuatie:

8

Page 9: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 12: Graf de dependenta pentru un algoritm iterativ

Figure 13: Nodurile colorate ale grafului de dependenta

x(t+ 1) = a(t)x(t) + b(t)x(t− 1)

Se pune problema calculului lui x(t).

Datele de intrare sunt: x(−1), x(0), a(0), . . . , a(n− 1), b(0), . . . , b(n− 1). Gasiti un algoritm paralelcare sa rezolve aceasta problema presupunand ca avem un numar suficient de procesoare.

2. Se considera urmatoarea ecuatie:

x(t+ 1) = a(t)x(t) + b(t)x(t− 1)

Se pune problema calculului lui x(t).

Datele de intrare sunt: x(n−1), x(0), a(0), . . . , a(n−1), b(0), . . . , b(n−1). Gasiti un algoritm paralelcare sa rezolve aceasta problema presupunand ca avem un numar suficient de procesoare.

3. Propuneti un algoritm paralel pentru calculul valorii polinomului p(x) = anxn+ · · ·+a1x+a0 pentru

o valoare data x.

4. Se considera graful de dependenta din figura 14, unde toate arcele sunt interpretate ca fiind bidirec-tionale.

9

Page 10: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 14: Graf de dependenta - exercitiu

• Precizati cu cate culori poate fi colorat graful si care este numarul minim de culori.

• Verificati de cate procesoare aveti nevoie penru o implementare paralela.

10

Page 11: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Optimizari neconstranse si constranse

1 Optimizari necosntanseSe pune problema minimizarii unei functii de cost continue F : Rn → Rn. Pentru orice vector x∗ valoarea∇F (x∗) = 0 minimizeaza functia F. Deci pentru a afla vetorul pentru acre se obtine minimul functie F trebuierezolvat un sistem de forma ∇F (x) = 0.

Cle doua abordari de demonstrare a convergentei algoritmilor de rezolvare a unui astfel de sistem sunt:

• abordarea descendenta - valoarea functiei cost descreste catre valoarea minimala;

• abordarea alternativa - se introduce o norma si se construiesc iteratii cu proprietatea ca distanta intreiteratia curenta si punctul de minim descreste;

1.1 AlgoritmiPresupunem ca dorim sa rezolvam un sistem Ax = b, unde a este o matrice pozitiva si simetrica. Rezolvareasistemului este echivalenta cu minimizarea functie cost F definita prin

F (x) =1

2x′Ax− x′b

In acest context ∇F (x) = Ax− b si ∇2F (x) = A.Cei ami uzuali algoritmi de rezolvare a sistemului sunt:

• algoritmul Jacobi

x(t+ 1) = x(t)− γ[D(x(t))]−1∇F (x(t))

unde γ este un scalar pozitiv si D este o matrice diagonala (dii = ∇2F (x)).

• algoritmul Seidel-Gauss

xi(t+ 1) = xi(t)− γ∇iF (z(i, t))

∇2iF (z(i, t))

unde z(i, t) = (x1(t+ 1). . . . , xi−1(t+ 1), xi(t), . . . , xn(t)).

• algorimul gradientului

x(t+ 1) = x(t)− γ∇F (x(t))

• algorimul gradientului scalat

x(t+ 1) = x(t)− γ(D(t))−1∇F (x(t))

unde D(t) este o matrice de scalare.

1

Page 12: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 1: Directii descendente in algorimii gradientului si gradientului scalat

• metode de tip Newton

x(t+ 1) = x(t)− γ(∇2F (x(t)))−1∇F (x(t))

1.2 Analiza convergenteiPentru a putea folosi abordarea descendenta trebuie ca functia F sa indeplineasca urmatoarele conditii:

• F (x) ≥ 0, ∀x ∈ Rn;

• functia F este diferentiabil continua si exista o constanta K astfel incat:

||∇F (x)−∇F (y)||2 ≤ K||x− y||2

∀x, y ∈ Rn.

Lemma 1.1 Daca F este Lipschitz continua atunci:

F (x+ y) ≤ F (x) + y′∇F (x) +K

2||y||22

∀x, y ∈ Rn.

Proposition 1.1 Fie K1 si K2 doua constante pozitive. Consideram sirul {x(t)} genreta cu ajutorul algo-ritmului

x(t+ 1) = x(t) + γs(t)

unde s(t) satisface

||s(t)||2 ≥ K1||∇F (x(t))||2,∀tsi

s(t)′∇F (x(t)) ≤ −K2||s(t)||22,∀t.Daca 0 < γ < 2K2

K1, atunci

limt→∞

∇F (x(t)) = 0

2

Page 13: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

1.3 Algoritmi neliniariIdeea de baza a acestor algoritmi este de a fixa toate valorile vectorului x mai putin valoarea de pe pozitia i side a minimiza functia F relativ la xi. Exista doua alternative de implementare a acestor algoritmi:

• algoritmul Jacobi neliniar

Figure 2: Iteratia algoritmului Jacobi neliniar

xi(t+ 1) = argminxi

F (x1(t), ,̇xi−1(t), xi, xi+1(t), . . . , xn(t))

• algoritmul Seidel-Gauss neliniar

xi(t+ 1) = argminxi

F (x1(t+ 1), ,̇xi−1(t+ 1), xi, xi+1(t), . . . , xn(t))

2 Optimizari constranseSe pune problema minimizarii functiei cost F : Rn → Rn pentr o multime X ⊂ Rn. Presupunem ca F estediferentiabil continua si ca multimea X este nevida, inchisa si convexa.

Conditiile necesare si suficiente pentru ca vectorul x ∈ X sa fie optim sunt:

• Daca vectorul x ∈ X minimizeaza F peste X, atunci

(y − x)′∇F (x) ≥ 0,∀y ∈ X.

• Daca F este convexa pe multimea X, atunci conditia anterioara este suficienta pentru a minimiza functiaf pe multimea X.

In continuare folosim notatia [x]+ pentru a nota proiectia ortogonala a vectorului x pe multimea convexaX. In particular, [x]+ este definit prin relatia:

[x]+ = argminz∈X

||z − x||2

3

Page 14: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 3: Iteratia algoritmului Seidel-Gauss neliniar

Theorem 2.1 1. Pentru orice x ∈ Rn, exista un unic z ∈ X care minimizeaza ||z − x||2 si care va finotat cu [x]+;

2. Fiind dat un x ∈ Rn, un vector z ∈ X este egal cu [x]+ daca si numai daca

(y − x)′(x− z) ≤ 0,∀y ∈ X

3. Functia f : Rn → X definita prin f(x) = [x]+ este continua si neexpansiva, adica

||[x]+ − [y]+||2 ≤ ||x− y||2,∀x, y ∈ Rn

2.1 Algoritmul gradientuluiAlgoritmul gradientului bazat pe proiectii este dat de urmatoarea ecuatie:

x(t+ 1) = [x(t)− γ∇F (x(t))]+

unde γ este un scalar pozitiv.Urmatoarea propozitie arata ca pentru un scalar γ suficient de mic fiecare iteratie a algorimului descreste

valoarea functiei de cost, cu conditia sa nu fie un punct fix al functiei T (x) = [x− γ∇F (x)]+.

Proposition 2.1 • F (T (x)) ≤ F (x)− ( 1γ − K2 )||T (x)− x||22

• T (x) = x daca si numai daca (y − x)′∇F (x) ≥ 0, ∀y ∈ X . In particular, daca F este convexa pemultimea X, atunci T (x) = x daca si numai daca x minimizeaza pe F pe multimea X.

• Functia F este continua.

4

Page 15: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 4: Conditia satisfacuta de proiectia [x]+

Theorem 2.2 Daca 0 < γ < 2K si daca x∗ este limita sirului {x(t)} generat de x(t + 1) = [x(t) −

γ∇F (x(t))]+ atunci (y − x)′∇F (x) ≥ 0, ∀y ∈ X .In particular, daca F este convexa pe multimea X, atunci x∗ minimizeaza F pe multimea X.

2.2 Algoritmul gradientului scalatAlgoritmul gradientului scalat bazat pe proiectii este dat de urmatoarea ecuatie:

x(t+ 1) = [x(t)− γ(M(t))−1∇F (x(t))]+

unde γ este un scalar pozitiv si M este o matrice inversabila (de obicei se considera matrica Hessiana∇2F (x(t))).

5

Page 16: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 5: Iteratia algoritmului gradientului bazat pe proiectii

6

Page 17: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 6: Iteratia algoritmului gradientului scalat bazat pe proiectii

7

Page 18: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Problema drumului de cost minim

Problema celui mai scurt drum este des intalnita. Presupunem ca avem un graf si fiecarui arc (i, j) algrafului ii asociem un cost sau o lungime aij . Lungimea drumului (i, i1, . . . , ik, j) de la nudul i la nodul jcu arcele (i, i1), . . . , (ik, j) este definita de suma lungimii arcelor aii1 + · · ·+ aikj . Problema este de a gasidrumul de lungime minima de la nodul i la nodul j.

Aplicatii ale problemei drumului de lungime minima:

• probleme de cautare;

• probleme de control optim;

• subrutine in rezolvarea altor probleme mari.

Problema drumului minim este un caz particular al programarii dinamice care trateaza probleme de de-cizie optima. Pornind din nodul i, prima decizie este de a alege unarc (i, i1), apoi de a alege un arc (i1, i2),si asa mai departe pana la nodul 1. Pentru prima decizie trebuie sa facem o balanta intre a alege un arc delungime mica si dorinta de a evita alegerea unui nod j foarte departe de destinatie. Aceasta problema estecaptata in ecuatia:

x∗i = min

j(aij + x∗

j ), i = 2, . . . , n

x∗1 = 0

care va fi satisfacuta de drumul minim x∗i , i = 1, . . . , n.

Algoritmul de programare dinamica pentru problema drumului minim are forma:

xi = minj

(aij + xj), i = 2, . . . , n

x1 = 0

Acest algoritm este potrivit pentru o implementare paralela sau distribuita deoarece minimizarea peste jpoate fi implementata paralel pentru toate nodurile i ̸= 1.

1 Problema celui mai scurt drumPresupunem ca avem un graf ce contine n noduri numerotate de la 1 la n. Notam cu A(i) multimea nodurilorj pentru care exista un arc (i, j) de la nodul i. Nodul 1 este un nod special numit destinatie. Presupunmeca A(1) este vida. Pentru fiecare arc (i, j) avem un scalar aij care reprezinta lungimea arcului. Lungimeadrumului {(i, i1), . . . , (ik, j)} care incepe din nodul i si se termina in nodul j prin suma aii1 + · · ·+ aikj . Sepune problema gasirii delui mai scurt drum sau a unui drum de lungime minima pornind din fiecare nod i.

Presupunem ca exista un drum de la fiecare nod i = 2, . . . , n la nodul destinatie 1 si ca fiecare ciclu arelungime pozitiva.

Drumul de lungime minma reprezinta solutia sitemului (obtinand astfel ecuatia Bellman):

x∗i = min

j∈A(i)(aij + x∗

j ), i = 2, . . . , n

1

Page 19: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

x∗1 = 0

Iteratia

xi = minj∈A(i)

(aij + xj), i = 2, . . . , n

x1 = 0

cunoscuta sub denumirea de algoritmul Bellman-Ford converge la o solutie pentru un vector initial arbitrar xcu x1 = 0.

2 Algoritmul Bellman-FordIteratia k a argoritmului Bellman-Ford are forma:

xki = min

j∈A(i)(aij + xk−1

j ), i = 2, . . . , n

xk1 = 0

In ceea ce priveste conditiile initiale, presupunem ca x01 = 0 si x0

i ∈ R pentru i = 2, . . . , n. Spunem caalgoritmul se termina dupa iteratia k daca xk

i = xk−1i , i = 1, . . . , n.

Dandu-se nodurile i ̸= 1 si j ̸= 1 definim:

wkij =drumul de lungime dintre nodurile i si j si care are k arce

Lemma 2.1 Avem urmatoarea relatie:

xki = min

j=1,...,n(wk

ij + x0j ),∀i = 2, . . . , n, k ≤ 1

Theorem 2.1 • Exista un cel mai scurt drum de la fiecare nod i la 1 care are cel mult n-1 arce.

• Pentru orice set de conditii initiale algoritmul Bellman-Ford se termina dupa un numar finit k deiteratii.

xki este egal cu cea mai mica distanta xa

i st

• Daca x0i ≤ x∗

i , algoritmul Bellman-Ford se termina in cel mult m∗ iteratii, unde

m∗ = maxi=2,...,n

mi ≥ n− 1

si mi este cel mai mic numar de arce continut intr-un drum minim de la i la 1.

• Cea mia mica distanta x∗i reprezinta solutia unica a ecuatiei Bellman.

Figura 1 prezinta un exemplu in care algoritmul esueaza daca ciclurile au lungime nenegativa.Numarul de iteratii depinde foarte mult de conditiile initiale. Pentru anumite conditii initiale, acest numar

poate fi foarte mare si poate depinde si de dimensiunea arcelor (figura 2).

2

Page 20: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 1: Problema drumului de cost minim in cazul in care avem un ciclu de lungime 0

Theorem 2.2 Algoritmul Bellman-Ford se termina dupa cel mult k + 1 iteratii, unde

k =

m∗, daca β ≤ 0n− 1, daca β > 0si graful este aciclic

n− 2 + [ βL ], daca β > 0si graful are cicluri

m∗ este dat de

m∗ = maxi=2,...,n

mi ≥ n− 1

si mi este cel mai mic numar de arce continut intr-un drum minim de la i la 1.

In figura 3 este prezentat cazul in care graful are cicluri.

3 Alte metode paralele pentru determinarea drumului de cost minim• Algoritmul Floyd-Warshall

Incepe cu conditiile initiale:

x0ij =

{aij , dacaj ∈ A(i)∞, altfel

si genereaza secvential pentru k = 0, 1, . . . , n− 1 si pentru tpate nodurile i si j

xk+1ij =

{min{xk

ij , xkik+1 + xk

k+1j}, dacaj ̸= i

∞, altfel

• Algoritmul de dublare este dat de urmatoarele formule:

x1ij =

aij , dacaj ∈ A(i)0, dacai = j∞, altfel

xk+1ij =

{minm

{xkim + xk

jm}, dacaj ̸= i

0, dacai = j

3

Page 21: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 2: Relatia dintre numarul de iteratii si lungimile arcelor pentru algoritmul Bellman-Ford

4

Page 22: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 3: Numarul de iteratii - algoritmul Bellman-Ford

5

Page 23: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Notiuni introductive

1 Arhitecturi paralele si distribuitePentru a realiza un studiu al metodelor numerice paralele si distribuite ar trebui in primul rand sa ne raportamla diferentele lor fata de implementarile seriale:

• alocarea task-urilor;

• comunicarea intre procesoare a rezultatele intermediare obtinute prin calcul;

• sincronizarea calculelor efectuate de procesoare. De aici se obtine o clasificare a metodelor: metodesincrone si metode asincrone.

Calculul paralel si distribuit reprezinta in prezent un domeniu larg de interes existand o activitate de cerc-etare intensa in acest sens motivata de numerosi factori. Intotdeauna a fost nevoie de o solutie de rezolvare aproblemelor cu calcule foarte mari, dar, doar recent, tehnologiile avansate au luat in calcul folosirea calcululuiparalel in rezolvarea acestor probleme tinand cont si de calculatoarele paralele din ce in ce mai puternice.

Dezvolatarea algoritmilor paraleli si distribuiti este determinata pe de o parte de legatura dintre nevoilede calcul noi si cele vechi, iar pe de alta parte de progresul tehnologic.

Aplicatiile care privesc inteligenta artificiala, calculul simbolic si calculul numeric sunt cele care au jucatun rol important in dezvoltarea arhitecturilor paralele si distribuite.

Nevoia de determinare rapida a unui numar mare de calcule numerice din cadrul aplicatiilor ce privescrezolvarea ecuatiilor cu derivate partiale si a procesarii imaginilor a fost ceea ce a dus la dezvoltarea rapiditatiisi a paralelizarii masinilor de calcul.

Recent a crescut interesul pentru aplicarea calculului paralel si distribuit in aplicatii precum cele ce privescanaliza, simularea si optimizarea la scara larga a sistemelor interconectate. Alte exemple de aplicare a calcu-lului paralel sunt: rezolvarea sistemelor de ecuatii, programarea matematica si problemele de optimizare. Oproprietate comuna a acestor probleme este aceea de a putea fi descompuse in subtask-uri.

O distictie importanta este intre sistemele de calculul paralel si cele de calcul distribuit.Sistemele de calcul paralel consistau in mai multe procesoare care sunt localizate la o distanta mica unele

de altele si care au fost concepute in scopul realizarii de taskuri de calcul in comun. Comunicarea dintreprocesoare este de incredere si poate fi anticipata.

Sistemele de calcul distribuit sunt diferite din mau multe puncte de vedere: procesoarele pot fi separate sicomunicarea dintre ele este mult mai problematica; intarzierile de comunicare nu pot fi anticipate si legaturilede comunicatie pot sa nu fie de incredere; fiecare procesor poate fi angajat in propriile sale activitati in timpce coopereaza cu alte procesoare pentru realizarea unui task comun.

Sunt mai multi parametrii care pot fi folositi pentru a descrie sau pentru a clasifica sistemele de calculparalel:

• tipul si numarul de procesoare. Sunt sisteme de calcul paralel care au multe procesoare numite masivparalele si sunt sisteme cu un numar mic de procesoare numite coarse grained parallelism.

• prezenta sau absenta unei masini de control global. O clasificare privind acest aspect si referindu-ne laabilitatea diferitelor procesoare de a executa diferite instructiuni la un moment dat este in calculatoareparalele SIMD (Single Instruction Multiple Data) si MIMD (Multiple Instruction Multiple Data).

1

Page 24: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

• operatii sincrone si asincrone. Distinctia in acest caz se refera la prezenta sau absenta unui ceasglobal folosit pentru a sincroniza operatiile diferitelor procesoare. Aceasta sincronizare este intot-deauna prezenta in SIMD. Operatiile sincrone sunt de preferat datorita catorva proprietati importante:comportamentul procesoarelor este mult mai usor de controlat si algoritmul este mult simplificat. Darin unele cazuri este destul de greu de sincronizat o retea de comunicare de date.

• legaturile dintre procesoare. Un aspect important al calculatoarelor paralele este mecanismul prin careprocesoarele fac schimb de informatii. Din acest punct de vedere exista doua cazuri: arhitecturi cumemorie comuna si arhitecturi care transmit mesaje. Plecand de la aceste doua tipuri de arhitecturiavem numerosi hibrizi. Cel mai des intalnit mod de a lucra cu memoria este switching systems (sistemede comutare) asa cum se poate observa in figura 1. Exista situatii in care nu avem memorie comuna, darfiecare procesor are propria sa memorie locala. In aceasta situatie procesoarele pot comunica printr-oretea de comunicare care consta in legaturi de conexiune care leaga anumite perechi de procesoare(2). Apare astfel problema de a stabili care procesoare sunt conectate. Ar fi de preferat daca toateprocesoarele ar fi direct conectate unele cu altele, sau daca pot fi legate printr-un numar de legaturiaditionale.

Figure 1: Sistem de conectare

Figure 2: Retea de comunicare

Structura retelei de conectare este foarte importanta atat in sistemele de calcul paralel cat si in cele decalcul distribuit dar cu o diferenta: in calculatoarele paralele reteaua de conectare este sub controlulunui designer si din aceasta cauza de obicei este regulata, iar in sistemele distribuite topologia reteleieste predeterminata si iregulata.

2

Page 25: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

2 Modele si exempleExista o varietate de modele de calcule paralele si distribuite. Un model simplu folosit in general in calcululparalel pentru algoritmi asincroni este acela al reprezentarii algoritmilor paraleli prin grafuri directe aciclice.Fiecare nod al grafului reprezinta operatiile pe care le face algoritmul si fiecare arc este folosit pentru areprezenta dependenta datelor. In particular, un arc (i, j) indica faptul ca operatia corespunzatoare nodului jfoloseste rezultatele operatiei corespunzatoare nodului i. Un astfel de exemplu este cel din figura 3.

Figure 3: Exemplu de graf direct orientat pentru (x1 + x2)(x2 + x3)

Un graf direct orientat realizeaza doar o reprezentare partiala a unui algoritm. In cadrul grafului suntspecificate operatiile care se efectueaza, operanzii si impune anumite conditii privind ordinea in care acesteoperatii se efectueaza. Pentru a determina un algoritm paralel trebuie sa specificam si ce procesor executa oanumita operatie si la ce moment o executa. Pentru aceasta trebuie impus ca un procesor sa execute cel multo operatie la un moment dat, si ca daca avem un nod (i, j) atunci operatia care ii corespunde nodului j poateincepe doar dupa ce operatia corespunzatoare nodului i a fost terminata.

Exista situatii cand pentru aceeasi problema de calcul exista mai multi algoritmi, deci sunt situatii incare avem mai multe grafuri direct orientate pentru rezolvarea unei probleme (figura 4). Apare astfel situatiadeterminarii unui graf ce minimizeaza timpul de calcul al solutiei problemei. Pentru cateva clase de problemeexista metode care construiesc un graf direct orientat cu un factor constant de optimizare.

In continuare consideram cateva exemple simple de taskuri de calcul numeric care pot fi reprezentate prinalgoritmi paraleli simpli. Toti algorimii considerati pot fi reprezentati printr-un graf direct orientat. Pentru ausura calculele presupunem ca fiecare operatie este executata intr-o unitate de timp si ca transmiterea de dateintre procesoare se realizeaza instantaneu.

2.1 Calcule cu vectori si matrici1. Adunarea scalarilor

Cea mai simpla sarcina de calcul este aceea de a aduna n scalari. Cel mai bun algoritm serial necesitan− 1 operatii, deci timpul de executie al algoritmului este n− 1.

In continuare prezentam un algoritm paralel construit sub presupunerea ca n este o putere a lui 2.Partitionam cei n scalari in n

2 perechi diferite si folosim n2 procesoare pentru a aduna cei doi scalari

din fiecare pereche. Deci, dupa o unitate de timp, avem de executat sarcina de a aduna n2 scalari.

Continuand in aceeasi maniera, dupa log n pasi vom mai avea de efectuat o singura adunare si calculul

3

Page 26: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 4: Graf direct orientat pentru (x1 + x2)(x2 + x3)

se va termina (figura 5). Acest algoritm poate fi generalizat la cazul in care n nu este putere a lui 2 siastfel timpul de executie devine ⌈log n⌉ folosind ⌊n

2 ⌋ procesoare.

Figure 5: Calcul paralel pentru adunarea a 16/15 scalari

Aceeasi problema a adunarii a n scalri poate fi rezolvata folosind graful direct orientat reprezentat infigura 6 si care aduce o imbunatatire a eficientei algoritmului.

2. Produsul scalar

Produsul scalarn∑

i=1

xiyi a doi vectori n-dimensionali poate fi calculata cu cu n procesoare in ⌈log n⌉+1

unitati de timp astfel: la primul pas fiecare procesor i calculeaza produsul xiyi si apoi avem nevoie de⌈log n⌉ unitati de timp pentru adunarea scalarilor folosind algoritmul anterior.

3. Adunarea si inmultirea matricilor

Suma a n matrici de dimensiune m×m poate fi calculta in ⌈log n⌉ unitati de timp si folosind m2⌊n2 ⌋

procesoare.

Inmultirea a doua matrici de dimensiune m × n si n × l consta in evaluarea a ml produse scalare devectori n-dimensionali si poate fi deci indeplinita in ⌈log n⌉+1 unitati de timp folosind nml procesoare.

4. Puterile unei matrici

4

Page 27: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 6: Algoritm alternativ pentru adunarea paralela a 16 scalari

Presupunem ca avem o matrice A de dimensiune n×n si dorim sa calculam Ak pentru k ∈ Z . Daca keste putere a lui 2, atunci task-ul poate fi indeplinit prin pasi repetitivi: calculam intai A2, apoi A2A2 =A4, etc. Dupa log k pasi se obtine Ak. Aceasta procedura presupune log k inmultiri consecutive dematrici si poate fi deci incheiata in log k(⌈log n⌉+ 1) unitati de timp folosind n3 procesoare.

O alta situatie de calcul a puterii unei matrici este reprezentata in figura 7.

Figure 7: Calcul paralel al puterilor unei matrici

2.2 Paralelizarea metodelor iterativePentru rezolvarea de sisteme de ecuatii, pentru optimizari si alte tipuri de probleme se folosesc in generalalgoritmi iterativi (metode de relaxare) de tipul:

x(t+ 1) = f(x(t)), t = 0, 1, . . . (1)

5

Page 28: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

unde fiecare x(t) este un vector n-dimensional si f : Rn → Rn.O notatie alternativa este x = f(x). Daca sirul {x(t)} generat de (1) este convergent la x∗ si daca functia

f este continua, atunci x∗ este punct fix al lui f , f(x∗) = x∗.Fie xi(t) componenta i a lui x(t) si fie fi componenta i a functiei f . Atunci putem scrie (1) ca:

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, n (2)

Algoritmul iterativ x = f(x) poate fi paralelizat prin luarea a n procesoare care sa calculeze o componentadiferita a lui x, asa cum se observa in (2). La fiecare pas, procesorul i cunoaste valorile tuturor componentelorlui x(t) de care depinde fi, calculeaza noua valoare a lui xi(t+1) si comunica aceasta valoare altor procesoarepentru a se putea trece la urmatoarea iteratie.

Comunicatiile necesare pentru executia unei iteratii de tipul (2) pot fi descrise cu ajutorul unui graf directca cel reprezentat in figura 8, numit si graf de dependenta. Nodurile reprezinta componentele lui x si pentrufiecare doua noduri i si j exista arcul (i, j), daca functia fj depinde de xi, asta inseamna ca procesorul itrebuie sa comunice valoarea lui xi(t) procesorului j.

Figure 8: Graf de dependenta

Graful din figura 8 este realizat pentru urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x2(t+ 1) = f2(x1(t), x2(t))x3(t+ 1) = f3(x2(t), x3(t), x4(t))x4(t+ 1) = f4(x2(t), x4(t))

Daca presupunem ca procesul iterativ reprezentat in ecuatia (1) se termina in T pasi, T ∈ N , structuraalgoritmului poate fi reprezentata printr-un graf direct orientat (9).

2.2.1 Iteratia Seidel-Gauss

Iteratia (2) in care toate componentele lui x sunt updatate simultan este numita iteratie de tip Jacobi. Intr-oforma alternativa, componentele lui x sunt updatate pe rand si valoarea cel mai recent calculata este folosita.Astfel ecuatia (2) devine:

xi(t+ 1) = fi(x1(t+ 1), . . . , xi−1(t+ 1), xi(t), . . . , xn(t)), i = 1, n (3)

Iteratia (3) se numeste algoritmul Seidel-Gauss bazat pe functia f .

Example 2.1 Consideram urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x2(t+ 1) = f2(x1(t+ 1), x2(t))x3(t+ 1) = f3(x2(t+ 1), x3(t), x4(t))x4(t+ 1) = f4(x2(t+ 1), x4(t))

6

Page 29: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 9: Graf direct orientat asociat grafului de dependenta

Figure 10: Paralelizarea iteratiei Seidel-Gauss

Exista mai multe alternative ale algoritmului Seidel-Gauss corespunzator functiei f , deoarece ordineain care sunt updatate componentele poate fi diferita. Atata timp cat viteza de convergenta nu se modificasemnificativ este natural sa schimbam ordinea de calcul a componentelor astfel incat paralelismul la fiecareiteratie este fie maximizat.

Example 2.2 Consideram urmatoarea situatie:x1(t+ 1) = f1(x1(t), x3(t))x3(t+ 1) = f3(x2(t+ 1), x3(t), x4(t))x4(t+ 1) = f4(x2(t+ 1), x4(t))x2(t+ 1) = f2(x1(t+ 1), x2(t))

Proposition 2.1 Urmatoarele afirmatii sunt echivalente:

7

Page 30: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 11: Paralelizarea iteratiei Seidel-Gauss - alta ordine de updatare a componentelor

• Exista o ordine a variabilelor astfel incat o parte a algoritmului Seidel-Gauss sa poata fi executata inK pasi paraleli.

• Exista o colorare a grafului de dependenta care foloseste K culori si cu proprietatea ca nu exista unciclu pozitiv cu toate nodurile intr-un ciclu care au aceleasi culori.

• Exista o colorare a grafului de dependenta care foloseste cel mult K culori astfel incat nodurile adia-cente au culori diferite.

In cazul in care consideram exemplul din figura 8, doua culori sunt suficiente: h(1) = h(3) = h(4) = 1si h(2) = 2.

In cazul in care se folosesc schemele cu culori pentru o implementare paralela a algoritmului Seidel-Gauss, este usor de folosit cate un procesor diferit pentru fiecare componenta a lui x, deoarece fiecare procesorva fi idle atata timp cat sunt updatate variabile de culori diferite. Un remediu imediat ar fi acela de a folosimai putine procesoare, fiecare procesor avand asignate mai multe variabile de culori diferite.

Example 2.3 Intr-o varietate de algoritmi de forma x = f(x), folositi pentru rezolvarea de ecuatii cuderivate partiale sau in procesarea de imagini, fiecarei componente a vectorului x ii este asociat un punctparticular dintr-o anumita regiune dintr-un spatiu 2-dimensional.

Fie N multimea tuturor punctelor (i, j) ∈ R2 astfel incat i si j sunt numere intregi ce satisfac 0 ≤ i ≤ Msi 0 ≤ j ≤ M . Fie xij o componenta a vectorului x care corespunde punctului (i, j). Prin conectarea celormai apropiati vecini obtinem graful reprezentat in figura 12. Vedem acest graf ca si graf direct daca facemarcele bidirectionale si presupunem ca este un graf de dependenta asociat ecuatiei x = f(x). Asignam cateun procesor fiecarui punct (i, j). Acest procesor este responsabil de updatarea valorii lui xij si pentru aface acest lucru are nevoie sa cunoasca valorile componentelor lui x asociate punctelor vecine. Este astfelnatural de presupus ca procesoarele asociate punctelor vecine sunt unite printr-o conexiune directa.

In ceea ce priveste implementarea unei metode Seidel-Gauss asociate situatiei prezentate, observam cagraful din figura 12 poate fi colorat folosind doua culori asa cum se observa in figura 13. Daca asignamcate un procesor fiecarei componente xij , atunci acesta va fi idle jumatate din timp. Este deci rezonabil saasignam doua componente de culori diferite fiecarui procesor. In practica, numarul de puncte este foartemare si deci fiecare procesor va avea asignate mai mult de doua puncte.

3 Exercitii1. Se considera urmatoarea ecuatie:

8

Page 31: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 12: Graf de dependenta pentru un algoritm iterativ

Figure 13: Nodurile colorate ale grafului de dependenta

x(t+ 1) = a(t)x(t) + b(t)x(t− 1)

Se pune problema calculului lui x(t).

Datele de intrare sunt: x(−1), x(0), a(0), . . . , a(n− 1), b(0), . . . , b(n− 1). Gasiti un algoritm paralelcare sa rezolve aceasta problema presupunand ca avem un numar suficient de procesoare.

2. Se considera urmatoarea ecuatie:

x(t+ 1) = a(t)x(t) + b(t)x(t− 1)

Se pune problema calculului lui x(t).

Datele de intrare sunt: x(n−1), x(0), a(0), . . . , a(n−1), b(0), . . . , b(n−1). Gasiti un algoritm paralelcare sa rezolve aceasta problema presupunand ca avem un numar suficient de procesoare.

3. Propuneti un algoritm paralel pentru calculul valorii polinomului p(x) = anxn+ · · ·+a1x+a0 pentru

o valoare data x.

4. Se considera graful de dependenta din figura 14, unde toate arcele sunt interpretate ca fiind bidirec-tionale.

9

Page 32: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 14: Graf de dependenta - exercitiu

• Precizati cu cate culori poate fi colorat graful si care este numarul minim de culori.

• Verificati de cate procesoare aveti nevoie penru o implementare paralela.

10

Page 33: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Probleme de sincronizare in algoritmii paraleli si distribuiti

1 IntroducereIn orice algoritm paralel sau distribuit este necesar sa se coordoneze, intr-o anumita masura, activitatilediferitelor procesoare. Aceasta coordonare este adesea pusa in aplicare prin impartirea algoritmului in etape.Pe durata fiecarei etape, fiecare procesor trebuie sa execute o serie de calcule care depind de rezultatelecalculelor altor procesoare din fazele anterioare, cu toate acestea, calendarul de calcule al unui procesor intimpul unei faze poate fi independent de calendarul de calcule al altor procesoare in aceeasi faza. Toateinteractiunile au loc la sfarsitul fiecarei faze. Numim astfel de algoritmi algoritmi sincroni, si in continuare,ii vom compara cu alti algoritmi, numiti asincroni, pentru care nu exista notiunea de faze si coordonareacalculelor diferitelor procesoare este mai putin stricta.

2 Algoritmi sincroniPresupunem ca un sir de calcule este impartit in segmente consecutive, numite faze si numerotate 1, 2, . . . .Calculele in fiecare faza sunt impartite intre n procesoare ale unui sistem de calcul. In timpul unei faze t,fiecare procesor i face anumite calcule folosind datele problemei alaturi de informatiile pe care le primeste dela alte procesoare in timpul fazelor anterioare 1, 2, . . . , t − 1. Procesorul i trimite apoi informatii sub formade mesaje catre toate celelalte procesoare P (i, t) ⊂ {1, 2, . . . , n} si apoi trece la faza t+ 1.

O presupunere implicita aici este faptul ca toate calcule diferitelor procesoare pot fi efectuate independent,intr-o faza, si ca momentul lor de executie este irelevant. Acest lucru este necesar pentru algoritmii distribuiticare reproduc rezultatele algoritmului serial original la sfarsitul fiecarei faze. In unele cazuri, un procesorpoate sa nu stie de la ce procesoare sa astepte un mesaj in timpul unei faze, asta inseamna ca un procesor jpoate sa nu stie multimea {i|j ∈ P (i, t)}.

Se considera urmatorul proces iterativ:

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, . . . n (1)

unde variabila xi este actualizata de procesorul i.Un algoritm distribuit corespunzator poate fi pus in aplicare prin asocierea fazelor cu momentele de timp

t. Cerinta aici este ca fiecare procesor i sa actualizeze variabila xi folosind procesul iterativ (1) si apoi satrimita un mesaj cu valoarea actualizata tuturor procesoarelor j pentru care xi apare in mod explicit in functiafj . Astfel, multimea de procesoare P (i, t) care a primit un mesaj de la i in faza t este o multime a tuturorprocesoarelor j pentru care (i, j) este un arc in graful de dependenta.

Consideram cazul in care procesul iterativ (1) este implementat pe o masina cu memorie partajata. Unprocesor trimite un mesaj la toate procesoarele simultan scriind-ul in memoria partajata. Sa presupunem caun procesor va incepe un calcul a unei noi faze numai dupa ce toate mesajele fazei precedente au fost scrisein memoria partajata. Avem astfel un caz special al modelului precedent, unde P (i, t) = {1, 2, . . . , n} sireceptarea mesajului se realizeaza simultan pentru toate procesoarele.

Un astfel de algoritm distribuit se numeste sincron. Algoritmul sincron este matematic echivalent cu unalgoritm guvernat de un ceas global, unul pentru care inceputul fiecarei faze si sfarsitul receptionarii mesajelorse realizeaza simultan pentru toate procesoarele. In scopul de a pune in aplicare un algoritm sincron intr-unsistem distribuit, avem nevoie de un mecanism de sincronizare. Un astfel de algoritm este numit synchronizer.

1

Page 34: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

2.1 Sincronizarea globalaIdeea este de a lasa fiecare procesor detectat atunci cand toate mesajele trimise in timpul unei faze au fostprimite, si doar apoi sa inceapa calculul in cadrul etapei urmatoare. Cel mai simplu mod de a efectuareasincronizarea globala este prin utilizarea de timeout.

Presupunem ca nu exista nici un ceas global accesibil tuturor procesoarelor, dar in schimb, fiecare pro-cesor poate masura cu precizie lungimea oricarui interval de timp folosind un ceas local. Sa presupunem caexista o limita superioara cunoscuta Tp pentru timpului necesar fiecarui procesor i pentru a executa calculeledintr-o faza, si pentru mesajele aferente trimise de catre i. Sa presupunem, de asemenea, ca exista un intervalde timp (necunoscut) de lungime cunoscuta Tf in care toate procesoarele au inceput prima faza. Apoi, sin-cronizarea va fi efectuata in cazul in care fiecare procesor i incepe faza k in k(Tp + Tj) unitati de timp dupace acesta a inceput prima faza. Figura 1 ilustreaza acest proces. Dificultatea aceastei metode este faptul calimitele Tp si Tf pot fi prea conservatoare sau indisponibile.

Figure 1: Implementarea sincronizarii globale

O alta abordare este aceea prin care fiecare procesor trimite un mesaj de terminare a fazei curente catretoate procesoarele dupa ce stie ca toate mesajele sale pentru o anumita faza au fost primite. Presupunem cafiecare mesaj trimis este recunoscut printr-un mesaj de raspuns trimis de receptor catre transmitator, astfelincat fiecare procesor stie ca in cele din urma toate mesajele pe care le trimite in timpul unei faze au fost prim-ite, moment in care se poate emite un mesaj de terminare faza. Odata ce un procesor a trimis propriul mesajde terminare de faza si a primit mesajele corespunzatoare de la toate celelalte procesoare, se pot continuacalculele din etapa urmatoare.

Intr-un sistem cu memorie partajata, aceasta metoda este conceptual simpla prin utilizarea de variabilespeciale, care sunt accesibile tuturor procesoarelor (ca exemplu avem semafoarele, monitoarele). Ideea debaza este aceea ca aceste variabile iau la un moment dat valori care indica faptul ca un procesor a termint fazacurenta si poate trece la faza urmatoare.

Intr-un sistem bazat pe transmitere de mesaje, sincronizarea globala poate fi implementata cu ajutorulunui arbore ce are un nod special, desemnat ca radacina (root). Mesajele de terminare a fazelor procesoarelorsunt colectate de root, incepand de la frunze. Odata ce radacina a primit mesaje de incetare a fazelor dela toate procesoarele, poate trimite un mesaj de initiere a fazei urmatoare tuturor procesoarelor, si fiecareprocesor poate incepe o noua faza la primirea acestui mesaj (figura 2).

Se poate observa ca aceasta metoda necesita, in esenta, o acumulare intr-un singur nod, urmata de oemisiune a unui singur nod.

2

Page 35: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 2: Sincronizarea globala folosind un arbore

2.2 Sincronizarea localaIdeea principala a acestei metode este ca, daca un procesor stie mesajele pe care sa le astepte in fiecare faza,atunci poate incepe o noua faza dupa ce a primit toate aceste mesaje. In particular, procesorul j poate incepecalculul unei noi faze t + 1 dupa ce a primit mesaje de fazei precedente t de la toate procesoarele i dinmultimea {i|j ∈ P (i, t)}.

Nu este necesar pentru un procesor sa stie daca oricare alte mesaje trimise in timpul fazei t (inclusivproprii mesaje) au fost primite, deci nu este nevoie sa se piarda timpul de asteptare pentru ca aceste receptiisa fie finalizate si sa fie confirmate. Atunci cand procesorul j nu stie multimea {i|j ∈ P (i, t)}, dar in schimbcunoaste o multime Sj ⊃ {i|j ∈ P (i, t)}, schema poate fi modificata astfel incat toate nodurile i ∈ Sj cuj /∈ P (i, j) trimit un mesaj catre j in cadrul fazei t.

3 Algoritmi asincroniTimpul total de executie a multor algoritmi poate fi adesea substantial redus printr-o implemenatre asincrona.Analiza algoritmilor asincroni distribuiti ofera un contrast preliminar si informal cu algoritmii sincroni.

Avand in vedere un algoritm distribuit, pentru fiecare procesor, exista o multime de timpi in care proce-sorul executa calculele, si alti timpi la care procesorul primeste mesaje de la celelalte procesore. Algoritmuleste numit sincron daca el este matematic echivalent cu unul pentru care timpul de calcul, transmiterea sireceptionarea mesajelor sunt fixe.

Spunem ca algoritmul este asincron in cazul in care timpul (prin urmare si ordinea calculelor si receptareamesajelor de la procesoare) poate varia pe scara larga pentru doua executii diferite ale algoritmului.

Pentru a analiza in contrast, ne putem gandi la un algoritm distribuit ca la o colectie de algoritmi locali.Fiecare algoritm local este executat de un procesor diferit si foloseste ocazional date generate de alti algoritmilocali.

3

Page 36: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

In cazul cel mai simplu de algoritm sincron, calendarul operatiunilor desfasurate de fiecare procesor estecomplet determinat si este aplicat cu ajutorul unui ceas global.

Un tip mai complex de algoritmi sincroni este unul in care momentul exact al operatiunilor desfasuratede fiecare algoritm local nu este prestabilit, dar algoritmii locali inca mai trebuie sa astepte la momenteleprestabilite ca datele sa devina disponibile.

Un exemplu de algoritm asincron este atunci cand algoritmii locali nu asteapta date predeterminate pentrua deveni disponibile; ei continua sa calculeze, incearca sa rezolve problema data cu orice date se intampla safie disponibile la momentul dat.

Cel mai extrem tip de algoritm asincron este cel ce poate tolera modificari in datele problemei sau in sis-temul de calcul distribuit, fara a se reporni de la un punct initial. Aceasta situatie apare in principal in retelelede date, in cazul in care nodurile si legaturile de comunicare nu pot fi reparate sau diferiti algoritmi distribuiticare controleaza reteaua sunt executati. In unele retele, cum ar fi retelele mobile de pachete radio, schimbarilein topologia retelei pot fi frecvente. In alte retele, astfel de schimbari pot fi rare, dar timpul de executie alalgoritmului poate fi atat de mare astfel ca probabilitatea de aparitie a unei modificari topologice in timp cealgoritmul se executa este foarte mica. De exemplu exista situatii in care algoritmul trebuie sa functionezechiar daca au loc erori in ceea ce priveste nodurile sau legaturile retelei. Exista o serie de dificultati aici. Inprimul rand, trebuie sa se stie mereu care sunt nodurile si legaturile care au dat erori; aceste informatii potavea un impact asupra algoritmului distribuit. In al doilea rnd, un nod sau o legatura cu erori va afecta, deobicei, algoritmul. Pentru a face fata situatiei, algoritmul ar trebui sa fie capabil sa se adapteze la esec sau artrebui sa fie anulata eroarea si sa fie repornit algoritmul (figura 3).

Figure 3: Exemplu de dificultate de comunicare

4

Page 37: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

3.1 Metode de relaxare asincroneConsideram un sistem cu n procesoare si o problema de punct fix n-dimensionala, unde dorim sa aflam unvector x = (x1, . . . , xn) ce satisface:

xi = fi(x1, . . . , xn), i = 1, . . . n (2)

unde fi sunt functii date de n variabile. Presupunem ca fiecare procesor i updateaza variabila xi conformrelatiei (2) pornind de la un set de date initiale.

O implementare sincrona a algoritmului necesita ca procesorul i sa nu realizeze updatarea k pana nua primit toate rezultatele de care are nevoie in cadrul functiei fi din etapa k − 1. In aceasta situatie aparprobleme de sincronizare:

• Procesorul i dupa ce updateaza variabila xi trebuie sa ramana idle in timp ce asteapta mesaje de la alteprocesoare (figura 4). In particular, o comunicare lenta incetineste intregul calcul.

Figure 4: Exemplu de comunicare cu intarziere

Figure 5: Exemplu de comunicare lenta (canal de comunicare lent/ procesor lent)

• Procesoarele care sunt rapide, fie pentru ca au o putere de calcul mare, fie pentru ca au putine calculede efectuat, trebuie sa astepte ca procesoarele mai lente sa-si termine iteratiile (figura 5).

Intr-o versiune asincrona a algoritmului precedent exista o flexibilitate mai mare in ceea ce privestefolosirea informatiilor primite de la alte procesoare. In cadrul etapei k, pentru updatarea variabilei procesoru-lui i, pot fi folosite date de la etapele trecute, nu neaparat cele mai recente rezultate. De exemplu, procesorul

5

Page 38: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

i poate executa uptatarea din cadrul etapei k folosind rezultate de la procesorul j de la etapa k + 10, si inacelasi timp procesorul j poate sa execute faza k + 100 folosind rezultate din faza k − 10 a procesorului i.

Viteza de calcul si de comunicare poate fi diferita pentru procesoare diferite si de aceea pot exista multipleintarzieri de comunicare. Aceste viteze si intarzieri pot varia foarte mult pe masura ce algorimul progreseaza.Exista deasemena flexibilitate in ceea ce priveste conditiile initiale ale fiecarui procesor, variabila xi disponi-bila initial poate diferi de la un procesor la altul.

Procesoarele pot deasemenea executa mai multe iteratii daca nu sunt fortate sa astepte rezultate mairecente pentru diferite variabile. Exista pericolul ca astfel de rezultate sa nu fie eficiente deoarece folosescinformatii care nu mai sunt valabile. Totusi, atunci cand functioneza, algoritmii asincroni reduc penalizarilede sincronizare si obtin mai rapid o solutie a problemei (figura 6).

Figure 6: Diagrama de timp pentru algoritmi sincroni/asincroni

Un dezavantaj important al algoritmilor asincroni este acela ca ei pot distruge proprietatea de convergentape care o are un algoritm atunci cand este executat sincron. Astfel sunt cazuri in care trebuiesc impuse limitariin ceea ce privesti intarzierile care pot exista pentru a se garanta convergenta algoritmului. In toate cazurile,analiza algorimilor asincroni este considerata a fi mai dificila decat cea a algoritmilor sincroni.

Algoritmii asincroni pot oferi un avantaj important in ceea ce priveste metodele de relaxare si acest lucrupoate fi observat analizand metodele Jacobi si Seidel-Gauss pentru rezolvarea problemelor de punct fix.

Metoda Jacobi genereaza un sir x(t) = (x1(t), . . . , xn()t) cu ajutorul relatiei

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, . . . , n (3)

plecand de la un set de date initiale.Metoda Seidel-Gauss foloseste o formula ce utilizeaza cea mai recenta updatare a unei variabile:

x1(t+ 1) = fi(x1(t), . . . , xn(t))

x2(t+ 1) = fi(x1(t+ 1), x2(t), . . . , xn(t))

. . . . . . . . .

xn(t+ 1) = fi(x1(t+ 1), . . . , xn−1(t+ 1), xn(t))

6

Page 39: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 7: Metoda Jacobi - algoritm sincron si asincron

Metoda Jacobi este mai potrivita pentru implementari distribuite si este in esenta echivalenta cu algoritmidistribuiti sincroni. Metoda Seidel-Gauss nu este foarte potrivita pentru implementari paralele si distribuite.Este totusi posibil ca algoritmii asincroni sa mareasca viteza metodei Seidel-Gauss in comparatie cu metodaJacobi fara sa sacrifice paralelismul metodei Jacobi (figura 7).

7

Page 40: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Probleme de sincronizare in algoritmii paraleli si distribuiti

1 IntroducereIn orice algoritm paralel sau distribuit este necesar sa se coordoneze, intr-o anumita masura, activitatilediferitelor procesoare. Aceasta coordonare este adesea pusa in aplicare prin impartirea algoritmului in etape.Pe durata fiecarei etape, fiecare procesor trebuie sa execute o serie de calcule care depind de rezultatelecalculelor altor procesoare din fazele anterioare, cu toate acestea, calendarul de calcule al unui procesor intimpul unei faze poate fi independent de calendarul de calcule al altor procesoare in aceeasi faza. Toateinteractiunile au loc la sfarsitul fiecarei faze. Numim astfel de algoritmi algoritmi sincroni, si in continuare,ii vom compara cu alti algoritmi, numiti asincroni, pentru care nu exista notiunea de faze si coordonareacalculelor diferitelor procesoare este mai putin stricta.

2 Algoritmi sincroniPresupunem ca un sir de calcule este impartit in segmente consecutive, numite faze si numerotate 1, 2, . . . .Calculele in fiecare faza sunt impartite intre n procesoare ale unui sistem de calcul. In timpul unei faze t,fiecare procesor i face anumite calcule folosind datele problemei alaturi de informatiile pe care le primeste dela alte procesoare in timpul fazelor anterioare 1, 2, . . . , t − 1. Procesorul i trimite apoi informatii sub formade mesaje catre toate celelalte procesoare P (i, t) ⊂ {1, 2, . . . , n} si apoi trece la faza t+ 1.

O presupunere implicita aici este faptul ca toate calcule diferitelor procesoare pot fi efectuate independent,intr-o faza, si ca momentul lor de executie este irelevant. Acest lucru este necesar pentru algoritmii distribuiticare reproduc rezultatele algoritmului serial original la sfarsitul fiecarei faze. In unele cazuri, un procesorpoate sa nu stie de la ce procesoare sa astepte un mesaj in timpul unei faze, asta inseamna ca un procesor jpoate sa nu stie multimea {i|j ∈ P (i, t)}.

Se considera urmatorul proces iterativ:

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, . . . n (1)

unde variabila xi este actualizata de procesorul i.Un algoritm distribuit corespunzator poate fi pus in aplicare prin asocierea fazelor cu momentele de timp

t. Cerinta aici este ca fiecare procesor i sa actualizeze variabila xi folosind procesul iterativ (1) si apoi satrimita un mesaj cu valoarea actualizata tuturor procesoarelor j pentru care xi apare in mod explicit in functiafj . Astfel, multimea de procesoare P (i, t) care a primit un mesaj de la i in faza t este o multime a tuturorprocesoarelor j pentru care (i, j) este un arc in graful de dependenta.

Consideram cazul in care procesul iterativ (1) este implementat pe o masina cu memorie partajata. Unprocesor trimite un mesaj la toate procesoarele simultan scriind-ul in memoria partajata. Sa presupunem caun procesor va incepe un calcul a unei noi faze numai dupa ce toate mesajele fazei precedente au fost scrisein memoria partajata. Avem astfel un caz special al modelului precedent, unde P (i, t) = {1, 2, . . . , n} sireceptarea mesajului se realizeaza simultan pentru toate procesoarele.

Un astfel de algoritm distribuit se numeste sincron. Algoritmul sincron este matematic echivalent cu unalgoritm guvernat de un ceas global, unul pentru care inceputul fiecarei faze si sfarsitul receptionarii mesajelorse realizeaza simultan pentru toate procesoarele. In scopul de a pune in aplicare un algoritm sincron intr-unsistem distribuit, avem nevoie de un mecanism de sincronizare. Un astfel de algoritm este numit synchronizer.

1

Page 41: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

2.1 Sincronizarea globalaIdeea este de a lasa fiecare procesor detectat atunci cand toate mesajele trimise in timpul unei faze au fostprimite, si doar apoi sa inceapa calculul in cadrul etapei urmatoare. Cel mai simplu mod de a efectuareasincronizarea globala este prin utilizarea de timeout.

Presupunem ca nu exista nici un ceas global accesibil tuturor procesoarelor, dar in schimb, fiecare pro-cesor poate masura cu precizie lungimea oricarui interval de timp folosind un ceas local. Sa presupunem caexista o limita superioara cunoscuta Tp pentru timpului necesar fiecarui procesor i pentru a executa calculeledintr-o faza, si pentru mesajele aferente trimise de catre i. Sa presupunem, de asemenea, ca exista un intervalde timp (necunoscut) de lungime cunoscuta Tf in care toate procesoarele au inceput prima faza. Apoi, sin-cronizarea va fi efectuata in cazul in care fiecare procesor i incepe faza k in k(Tp + Tj) unitati de timp dupace acesta a inceput prima faza. Figura 1 ilustreaza acest proces. Dificultatea aceastei metode este faptul calimitele Tp si Tf pot fi prea conservatoare sau indisponibile.

Figure 1: Implementarea sincronizarii globale

O alta abordare este aceea prin care fiecare procesor trimite un mesaj de terminare a fazei curente catretoate procesoarele dupa ce stie ca toate mesajele sale pentru o anumita faza au fost primite. Presupunem cafiecare mesaj trimis este recunoscut printr-un mesaj de raspuns trimis de receptor catre transmitator, astfelincat fiecare procesor stie ca in cele din urma toate mesajele pe care le trimite in timpul unei faze au fost prim-ite, moment in care se poate emite un mesaj de terminare faza. Odata ce un procesor a trimis propriul mesajde terminare de faza si a primit mesajele corespunzatoare de la toate celelalte procesoare, se pot continuacalculele din etapa urmatoare.

Intr-un sistem cu memorie partajata, aceasta metoda este conceptual simpla prin utilizarea de variabilespeciale, care sunt accesibile tuturor procesoarelor (ca exemplu avem semafoarele, monitoarele). Ideea debaza este aceea ca aceste variabile iau la un moment dat valori care indica faptul ca un procesor a termint fazacurenta si poate trece la faza urmatoare.

Intr-un sistem bazat pe transmitere de mesaje, sincronizarea globala poate fi implementata cu ajutorulunui arbore ce are un nod special, desemnat ca radacina (root). Mesajele de terminare a fazelor procesoarelorsunt colectate de root, incepand de la frunze. Odata ce radacina a primit mesaje de incetare a fazelor dela toate procesoarele, poate trimite un mesaj de initiere a fazei urmatoare tuturor procesoarelor, si fiecareprocesor poate incepe o noua faza la primirea acestui mesaj (figura 2).

Se poate observa ca aceasta metoda necesita, in esenta, o acumulare intr-un singur nod, urmata de oemisiune a unui singur nod.

2

Page 42: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 2: Sincronizarea globala folosind un arbore

2.2 Sincronizarea localaIdeea principala a acestei metode este ca, daca un procesor stie mesajele pe care sa le astepte in fiecare faza,atunci poate incepe o noua faza dupa ce a primit toate aceste mesaje. In particular, procesorul j poate incepecalculul unei noi faze t + 1 dupa ce a primit mesaje de fazei precedente t de la toate procesoarele i dinmultimea {i|j ∈ P (i, t)}.

Nu este necesar pentru un procesor sa stie daca oricare alte mesaje trimise in timpul fazei t (inclusivproprii mesaje) au fost primite, deci nu este nevoie sa se piarda timpul de asteptare pentru ca aceste receptiisa fie finalizate si sa fie confirmate. Atunci cand procesorul j nu stie multimea {i|j ∈ P (i, t)}, dar in schimbcunoaste o multime Sj ⊃ {i|j ∈ P (i, t)}, schema poate fi modificata astfel incat toate nodurile i ∈ Sj cuj /∈ P (i, j) trimit un mesaj catre j in cadrul fazei t.

3 Algoritmi asincroniTimpul total de executie a multor algoritmi poate fi adesea substantial redus printr-o implemenatre asincrona.Analiza algoritmilor asincroni distribuiti ofera un contrast preliminar si informal cu algoritmii sincroni.

Avand in vedere un algoritm distribuit, pentru fiecare procesor, exista o multime de timpi in care proce-sorul executa calculele, si alti timpi la care procesorul primeste mesaje de la celelalte procesore. Algoritmuleste numit sincron daca el este matematic echivalent cu unul pentru care timpul de calcul, transmiterea sireceptionarea mesajelor sunt fixe.

Spunem ca algoritmul este asincron in cazul in care timpul (prin urmare si ordinea calculelor si receptareamesajelor de la procesoare) poate varia pe scara larga pentru doua executii diferite ale algoritmului.

Pentru a analiza in contrast, ne putem gandi la un algoritm distribuit ca la o colectie de algoritmi locali.Fiecare algoritm local este executat de un procesor diferit si foloseste ocazional date generate de alti algoritmilocali.

3

Page 43: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

In cazul cel mai simplu de algoritm sincron, calendarul operatiunilor desfasurate de fiecare procesor estecomplet determinat si este aplicat cu ajutorul unui ceas global.

Un tip mai complex de algoritmi sincroni este unul in care momentul exact al operatiunilor desfasuratede fiecare algoritm local nu este prestabilit, dar algoritmii locali inca mai trebuie sa astepte la momenteleprestabilite ca datele sa devina disponibile.

Un exemplu de algoritm asincron este atunci cand algoritmii locali nu asteapta date predeterminate pentrua deveni disponibile; ei continua sa calculeze, incearca sa rezolve problema data cu orice date se intampla safie disponibile la momentul dat.

Cel mai extrem tip de algoritm asincron este cel ce poate tolera modificari in datele problemei sau in sis-temul de calcul distribuit, fara a se reporni de la un punct initial. Aceasta situatie apare in principal in retelelede date, in cazul in care nodurile si legaturile de comunicare nu pot fi reparate sau diferiti algoritmi distribuiticare controleaza reteaua sunt executati. In unele retele, cum ar fi retelele mobile de pachete radio, schimbarilein topologia retelei pot fi frecvente. In alte retele, astfel de schimbari pot fi rare, dar timpul de executie alalgoritmului poate fi atat de mare astfel ca probabilitatea de aparitie a unei modificari topologice in timp cealgoritmul se executa este foarte mica. De exemplu exista situatii in care algoritmul trebuie sa functionezechiar daca au loc erori in ceea ce priveste nodurile sau legaturile retelei. Exista o serie de dificultati aici. Inprimul rand, trebuie sa se stie mereu care sunt nodurile si legaturile care au dat erori; aceste informatii potavea un impact asupra algoritmului distribuit. In al doilea rnd, un nod sau o legatura cu erori va afecta, deobicei, algoritmul. Pentru a face fata situatiei, algoritmul ar trebui sa fie capabil sa se adapteze la esec sau artrebui sa fie anulata eroarea si sa fie repornit algoritmul (figura 3).

Figure 3: Exemplu de dificultate de comunicare

4

Page 44: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

3.1 Metode de relaxare asincroneConsideram un sistem cu n procesoare si o problema de punct fix n-dimensionala, unde dorim sa aflam unvector x = (x1, . . . , xn) ce satisface:

xi = fi(x1, . . . , xn), i = 1, . . . n (2)

unde fi sunt functii date de n variabile. Presupunem ca fiecare procesor i updateaza variabila xi conformrelatiei (2) pornind de la un set de date initiale.

O implementare sincrona a algoritmului necesita ca procesorul i sa nu realizeze updatarea k pana nua primit toate rezultatele de care are nevoie in cadrul functiei fi din etapa k − 1. In aceasta situatie aparprobleme de sincronizare:

• Procesorul i dupa ce updateaza variabila xi trebuie sa ramana idle in timp ce asteapta mesaje de la alteprocesoare (figura 4). In particular, o comunicare lenta incetineste intregul calcul.

Figure 4: Exemplu de comunicare cu intarziere

Figure 5: Exemplu de comunicare lenta (canal de comunicare lent/ procesor lent)

• Procesoarele care sunt rapide, fie pentru ca au o putere de calcul mare, fie pentru ca au putine calculede efectuat, trebuie sa astepte ca procesoarele mai lente sa-si termine iteratiile (figura 5).

Intr-o versiune asincrona a algoritmului precedent exista o flexibilitate mai mare in ceea ce privestefolosirea informatiilor primite de la alte procesoare. In cadrul etapei k, pentru updatarea variabilei procesoru-lui i, pot fi folosite date de la etapele trecute, nu neaparat cele mai recente rezultate. De exemplu, procesorul

5

Page 45: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

i poate executa uptatarea din cadrul etapei k folosind rezultate de la procesorul j de la etapa k + 10, si inacelasi timp procesorul j poate sa execute faza k + 100 folosind rezultate din faza k − 10 a procesorului i.

Viteza de calcul si de comunicare poate fi diferita pentru procesoare diferite si de aceea pot exista multipleintarzieri de comunicare. Aceste viteze si intarzieri pot varia foarte mult pe masura ce algorimul progreseaza.Exista deasemena flexibilitate in ceea ce priveste conditiile initiale ale fiecarui procesor, variabila xi disponi-bila initial poate diferi de la un procesor la altul.

Procesoarele pot deasemenea executa mai multe iteratii daca nu sunt fortate sa astepte rezultate mairecente pentru diferite variabile. Exista pericolul ca astfel de rezultate sa nu fie eficiente deoarece folosescinformatii care nu mai sunt valabile. Totusi, atunci cand functioneza, algoritmii asincroni reduc penalizarilede sincronizare si obtin mai rapid o solutie a problemei (figura 6).

Figure 6: Diagrama de timp pentru algoritmi sincroni/asincroni

Un dezavantaj important al algoritmilor asincroni este acela ca ei pot distruge proprietatea de convergentape care o are un algoritm atunci cand este executat sincron. Astfel sunt cazuri in care trebuiesc impuse limitariin ceea ce privesti intarzierile care pot exista pentru a se garanta convergenta algoritmului. In toate cazurile,analiza algorimilor asincroni este considerata a fi mai dificila decat cea a algoritmilor sincroni.

Algoritmii asincroni pot oferi un avantaj important in ceea ce priveste metodele de relaxare si acest lucrupoate fi observat analizand metodele Jacobi si Seidel-Gauss pentru rezolvarea problemelor de punct fix.

Metoda Jacobi genereaza un sir x(t) = (x1(t), . . . , xn()t) cu ajutorul relatiei

xi(t+ 1) = fi(x1(t), . . . , xn(t)), i = 1, . . . , n (3)

plecand de la un set de date initiale.Metoda Seidel-Gauss foloseste o formula ce utilizeaza cea mai recenta updatare a unei variabile:

x1(t+ 1) = fi(x1(t), . . . , xn(t))

x2(t+ 1) = fi(x1(t+ 1), x2(t), . . . , xn(t))

. . . . . . . . .

xn(t+ 1) = fi(x1(t+ 1), . . . , xn−1(t+ 1), xn(t))

6

Page 46: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 7: Metoda Jacobi - algoritm sincron si asincron

Metoda Jacobi este mai potrivita pentru implementari distribuite si este in esenta echivalenta cu algoritmidistribuiti sincroni. Metoda Seidel-Gauss nu este foarte potrivita pentru implementari paralele si distribuite.Este totusi posibil ca algoritmii asincroni sa mareasca viteza metodei Seidel-Gauss in comparatie cu metodaJacobi fara sa sacrifice paralelismul metodei Jacobi (figura 7).

7

Page 47: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Algoritmi sincroni pentru sistemele de ecuatii liniare siinversarea matricilor

Fie A o matrice patratica reala (A ∈ Rn×n) si b un vector din Rn. Se considera sistemul:

Ax = b (1)

unde x este un vector necunoscut ce trebuie determinat.Exista o varietate de metode de rezolvare a acestui tip de sistem de obicei clasificate in metode directe

si metode iterative. Metodele directe determina solutia exacta a sistemului intr-un numar finit de operatii side obicei sunt de ordinul n3. Metodele iterative nu obtin solutia exacta a sistemului intr-un numar finit depasi, dar converg la o solutie asimptotica. Metodele iterative sunt uneori preferate metodelor directe deoarececonduc la o solutie cu o eroare acceptabila dupa un numar relativ mic de pasi.

O problema mai generala este aceea de calcul a inversei unei matrici (A−1), problema ce poate fi rezolvataatat prin metode directe cat si prin metrode iterative.

In continuare vom studia algoritmi paraleli pentru rezolvarea sistemelor (1) si pentru determinarea inverseimatricilor. O parte din algoritmii studiati vor fi doar implementari paralele ale unor algoritmi seriali, iarrestul vor fi algoritmi construiti special pentru o paralelizare buna. Vom analiza complexitatea si eficientametodelor directe, notiuni ce nu pot fi analizate in cazul metodelor iterative deoarece acestea nu se terminaintr-un timp finit, dar in legatura cu aceste metode vom discuta despre convergenta geometrica sau rataprogresiei geometrice. Aceasta inseamna ca vectorul {x(t)} generat de algoritmul iterativ are proprietatea||x(t)−x∗|| ≤ cρt, unde x∗ este solutia sistemului Ax = b, c este o constanta pozitiva, iar ρ este o constantapozitiva mai mica decat 1. Cu cat valoarea lui ρ este mai mica, cu atat algoritmul converge mai repede.

1 Algoritmi paraleli pentru rezolvarea sistemelor liniare cu structurispeciale

1.1 Matrici triunghiulare si metoda substitutiei inverseFie A o matrice reala triungiular inferioara (A ∈ Rn×n, aij = 0 pentru i < j). Ne propunem sa calculamA−1 in conditiile in care presupunem ca matricea este inversabila si ca aii ̸= 0 pentru orice i.

Consideram cazul in care aii = 1 pentru orice i. Notam A = I − L, unde L este o matrice strict inferiortriunghiulara (lij = 0 pentru i ≤ j).

Lemma 1.1 Daca A = I − L, unde L este o matrice strict inferior triunghiulara, atunci:

A−1 = I + L+ L2 + · · ·+ Ln−1.

Proof. Avem urmatoarea relatie

(I + L+ L2 + · · ·+ Ln−1)(I − L) = I − Ln = I

1

Page 48: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Un algoritm paralel pentru calculul inversei unei matrici este: se calculeaza paralel puterile matricii L(L2, . . . , Ln−1) si apoi se aduna rezultatele.

Un astfel de algoritm foloseste un numar foarte mare de procesoare, n4. Un algoritm mai eficient seobtine folosind urmatoarea lema:

Lemma 1.2 Daca A = I − L, unde L este o matrice strict inferior triunghiulara, atunci:

A−1 = (I + L)(I + L2)(I + L4) . . . (I + L2⌈log n⌉−2

)(I + L2⌈log n⌉−1

).

Obtinem astfel urmatorul algoritm: se calculeaza matricile L2, L4, . . . , L2⌈log n⌉−2

, L2⌈log n⌉−1

, adaugamidentitatea la fiecare din aceste matrici, apoi facem inmultirile dintre matricile obtinute. Timpul de executieeste la fel ca in cazul precedent, dar se folosesc n3 procesoare.

Presupunem ca nu avem relatia aii = 1 pentru orice i. Definim matricea diagonala D, astfel incat dii =aii, deci matricea D−1A este inferior triunghiulara si are elementele de pe diagonala principala egale cu 1.

Algoritmul in acest caz este: transformam A in D−1A, apoi inversam aceasta matrice si obtinem A−1D,apoi inmultim cu D−1 si obtinem A−1. Timpul de executie este acelasi ca in cazul precedent si se folosescn3 procesoare.

Un alt algoritm de calcul paralel al inversei unei matrici este imparte si cucereste. Impartim matricea Ain blocuri:

A =

[A1 0A2 A3

]unde A1 este de dimensiune ⌈n

2 ⌉ × ⌈n2 ⌉, A1 si A3 sunt inferior triunghiulare, iar inversa se calculeaza cu

formula:

A−1 =

[A−1

1 0−A−1

3 A2A−11 A−1

3

]Avem astfel urmatorul algoritm:

Algorithm 1 Algoritm pentru calculul inversei unei matriciif n = 1 then

A−1 este evidentend ifif n ̸= 1 then

partitioneaza matricea Ainverseaza A1 si A3

inmulteste A−13 cu A2 pentru a obtine A−1

3 A2

inmulteste A−13 A2 cu A−1

1 pentru a obtine A−13 A2A

−11

end if

Pentru executia algoritmului intr-un timp bun avem nevoie de n3 procesoare.O metoda practica de rezolvare a sistemului Ax = b este metoda substitutiei inverse, care se obtine prin

paralelizarea metodei secventiale de rezolvare a sistemelor liniare.Presupunem ca matricea A este inferior triunghiulara si ca avem a i-a ecuatie a sistemului Ax = b:

ai1x1 + ai2x2 + · · ·+ ainxn = bi.

Algoritmul foloseste n procesoare (figurile 1 si 2). Presupunem ca la pasul i avem calculate valorilepentru x1, . . . , xi−1 si aj1x1 + aj2x2 + · · · + aji−1xi−1 pentru j ≥ i. Procesorul i calculeaza valoarea luixi folosind formula:

2

Page 49: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

xi =1

aii(bi − ai1x1 − · · · − aii−1xi−1).

Apoi fiecare procesor j cu j ≥ i + 1 calculeaza aj1x1 + aj2x2 + · · · + ajixi prin adaugarea lui ajixi

la expresia deja cunoscuta aj1x1 + aj2x2 + · · · + aji−1xi−1. Algoritmul se termina la pasul n cand toatevalorile x1, . . . , xn sunt calculate.

Figure 1: Implementarea metodei substitutiei inverse

Figure 2: Implementarea alternativa pentru substitutia inversa

3

Page 50: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Metoda substitutiei inverse poate fi folosita si pentru calculul inversei unei matrici dupa cum urmeaza.Deoarece AA−1 = I , atunci coloana i, xi, a lui A−1 satisface Axi = ei, unde ei este coloana i a matricei I.Astfel inversa este obtinuta prin rezolvarea a s sisteme de n ecuatii.

1.2 Sisteme tridiagonaleSe considera sistemul Ax = b, unde A este matrice tridiagonala (aij = 0 pentru |i − j| > 1). Un astfel desistem are forma:

g1x1 + h1x2 = b1 (2)

fixi−1 + gixi + hixi+1 = bi, i = 2, 3, . . . , n− 1 (3)

fnxn−1 + gnxn = bn (4)

unde gi sunt elementele de pe diagonala principala, fi si hi sunt sub si deasupara diagonalei principale.Exista mai multe metode de rezolvare a acestui sistem, una din ele este reductia para-impara (figura 3).Folosim notatia x0 = xn+1 = 0 si din ecuatia (3) obtinem:

xi =1

gi(bi − fixi−1 − hixi+1) (5)

Folosim ecuatia (5) in care inlocuim i cu i− 1 si cu i+ 1 si obtinem:

figi−1

(bi−1 − fi−1xi−2 − hi−1xi) + gixi +hi

gi+1(bi+1 − fi+1xi − hi+1xi+2) = bi

care simplificata devine:

−fifi−1

gi−1xi−2 + (gi −

hi−1figi−1

− hifi+1

gi+1)xi −

hihi+1

gi+1xi+2 = bi −

figi−1

bi−1 −hi

gi+1bi+1 (6)

Se obtine astfel un sistem tridiagonal mai mic pentru a carui rezolvare se poate aplica aceeasi metoda.

2 Metode directe de rezolvare a ecuatiilor liniareFie A o matrice patratica n× n. Cele mai multe metode de rezolvarea a unui sistem Ax = b intai transformasistemul pana cand matricea A devine triunghiulara si apoi rezolva noul sistem cu metoda substitutiei inverse.Aceste transformari constau in inmultiri ale sitemului cu matricile M1, . . . ,Mk. Se obtine astfel un nousistem ce poate fi rezolvat cu metoda substitutiei inverse:

Mk . . .M1x = Mk . . .M1b.

Daca se doreste calculul inversei matricii A, aceasta poate fi calculata folosind relatia:

A−1 = (Mk . . .M1A)−1(Mk . . .M1).

2.1 Eliminarea GaussianaFie C0 = A si Ci = M i . . .M1A. Presupunem ca Ci−1 a fost deja calculata pentru 1 ≤ i ≤ n − 1 si areforma din figura 4.

Calculam matricea M i astfel incat Ci = M iCi−1 are aceeasi proprietate. Pentru aceasta consideram

M i = I−N i, unde N i are forma din figura 5, iar valorile nevide se calculeaza cu relatia N iji =

Ci−1ji

Ci−1ii

, j > i.

Se poate verifica faptul ca matricea Ci are forma din figura 4 si se poate calcula folosind direct relatia:

4

Page 51: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 3: Reductia para-impara

Cijk = Ci−1

jk −n∑

l=1

N ijlC

i−1lk = Ci−1

jk −N ijiC

i−1ik = Ci−1

jk −Ci−1

ji

Ci−1ii

Ci−1ik

In cazul in care Ci−1ii = 0, se poate folosi pivotrea partiala sau totala.

2.2 Triangularizarea matricilor folosind metoda rotatiilor GivensO diferenta notabila fata de triangularizarea prin metoda Gauss este aceea ca matricile M i sunt ortogonale(||M ix|| = ||x||, ∀x ∈ Rn). Aceasta relatie este echivalenta cu (M i)tM i = I .

Fie C o matrice data astfel incat Cil = Cjl = 0, pentru l < k si Cjk ̸= 0, deci are forma din figura 6.Se considera o matrice M, numita rotatie Givens, de forma 7 astfel incat:

c =Cik√

C2ik + C2

jk

s =Cjk√

C2ik + C2

jk

Acesta matrice M are urmatoarele proprietati:

• Toate liniile din MC sunt egale cu liniile din C mai putin liniile i si j;

• MCil = MCjl = 0, pentru l < k si MCjk = 0

3 Metode directe de calcul a inversei unei matriciFie A o matrice patratica nesingulara (n× n). Consideram polinomul caracteristic al matricei A:

5

Page 52: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 4: Matricea C

Figure 5: Matricea N

P (λ) = det(λI −A) =n∏

i=1

(λ− λi)

unde λi, i = 1, n, sunt valorile proprii ale matricei A. Fie ci, i = 1, n, coeficientii polinomului caracter-istic P (λ).

P (λ) = λn + c1λn−1 + · · ·+ cn−1λ+ cn

Se poate observa ca cn = (−1)nn∏

i=1

λi si deci matricea A este inversabila daca cn este nenul.

Folosind proprietatea lui Cayley-Hamilton (P (A) = An + c1An−1 + · · ·+ cn−1A+ cnI = 0), obtinem

ca inversa A−1 se calculeaza cu formula:

A−1 = − 1

cn(An−1 + c1A

n−2 + · · ·+ cn−1I)

Metoda directa de calcul a inversei unei matrici este data de urmatorul algoritm:

6

Page 53: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 6: Matricea C

Figure 7: Matricea M

7

Page 54: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Algorithm 2 Calculul inversei unei matriciCalculeaza Ak pentru k = 2, . . . , n

Calculeaza sk pentru k = 1, . . . , n folosind relatia sk =n∑

n=1λki = tr(Ak)

Rezolva sistemul urmator determinand c1, . . . , cn:

1 0 . . . . . . . . . 0

s1 2 0...

......

...... . . . . . . . . .

...sk−1 . . . s1 k 0 . . .

......

......

... 0sn−1 . . . sk−1 . . . s1 n

c1............cn

= −

s1............sn

Determina A−1 folosind relatia:

A−1 = − 1

cn(An−1 + c1A

n−2 + · · ·+ cn−1I)

8

Page 55: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Algoritmi sincroni pentru sistemele de ecuatii liniare siinversarea matricilor

Fie A o matrice patratica reala (A ∈ Rn×n) si b un vector din Rn. Se considera sistemul:

Ax = b (1)

unde x este un vector necunoscut ce trebuie determinat.Exista o varietate de metode de rezolvare a acestui tip de sistem de obicei clasificate in metode directe

si metode iterative. Metodele directe determina solutia exacta a sistemului intr-un numar finit de operatii side obicei sunt de ordinul n3. Metodele iterative nu obtin solutia exacta a sistemului intr-un numar finit depasi, dar converg la o solutie asimptotica. Metodele iterative sunt uneori preferate metodelor directe deoarececonduc la o solutie cu o eroare acceptabila dupa un numar relativ mic de pasi.

O problema mai generala este aceea de calcul a inversei unei matrici (A−1), problema ce poate fi rezolvataatat prin metode directe cat si prin metrode iterative.

In continuare vom studia algoritmi paraleli pentru rezolvarea sistemelor (1) si pentru determinarea inverseimatricilor. O parte din algoritmii studiati vor fi doar implementari paralele ale unor algoritmi seriali, iarrestul vor fi algoritmi construiti special pentru o paralelizare buna. Vom analiza complexitatea si eficientametodelor directe, notiuni ce nu pot fi analizate in cazul metodelor iterative deoarece acestea nu se terminaintr-un timp finit, dar in legatura cu aceste metode vom discuta despre convergenta geometrica sau rataprogresiei geometrice. Aceasta inseamna ca vectorul {x(t)} generat de algoritmul iterativ are proprietatea||x(t)−x∗|| ≤ cρt, unde x∗ este solutia sistemului Ax = b, c este o constanta pozitiva, iar ρ este o constantapozitiva mai mica decat 1. Cu cat valoarea lui ρ este mai mica, cu atat algoritmul converge mai repede.

1 Algoritmi paraleli pentru rezolvarea sistemelor liniare cu structurispeciale

1.1 Matrici triunghiulare si metoda substitutiei inverseFie A o matrice reala triungiular inferioara (A ∈ Rn×n, aij = 0 pentru i < j). Ne propunem sa calculamA−1 in conditiile in care presupunem ca matricea este inversabila si ca aii ̸= 0 pentru orice i.

Consideram cazul in care aii = 1 pentru orice i. Notam A = I − L, unde L este o matrice strict inferiortriunghiulara (lij = 0 pentru i ≤ j).

Lemma 1.1 Daca A = I − L, unde L este o matrice strict inferior triunghiulara, atunci:

A−1 = I + L+ L2 + · · ·+ Ln−1.

Proof. Avem urmatoarea relatie

(I + L+ L2 + · · ·+ Ln−1)(I − L) = I − Ln = I

1

Page 56: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Un algoritm paralel pentru calculul inversei unei matrici este: se calculeaza paralel puterile matricii L(L2, . . . , Ln−1) si apoi se aduna rezultatele.

Un astfel de algoritm foloseste un numar foarte mare de procesoare, n4. Un algoritm mai eficient seobtine folosind urmatoarea lema:

Lemma 1.2 Daca A = I − L, unde L este o matrice strict inferior triunghiulara, atunci:

A−1 = (I + L)(I + L2)(I + L4) . . . (I + L2⌈log n⌉−2

)(I + L2⌈log n⌉−1

).

Obtinem astfel urmatorul algoritm: se calculeaza matricile L2, L4, . . . , L2⌈log n⌉−2

, L2⌈log n⌉−1

, adaugamidentitatea la fiecare din aceste matrici, apoi facem inmultirile dintre matricile obtinute. Timpul de executieeste la fel ca in cazul precedent, dar se folosesc n3 procesoare.

Presupunem ca nu avem relatia aii = 1 pentru orice i. Definim matricea diagonala D, astfel incat dii =aii, deci matricea D−1A este inferior triunghiulara si are elementele de pe diagonala principala egale cu 1.

Algoritmul in acest caz este: transformam A in D−1A, apoi inversam aceasta matrice si obtinem A−1D,apoi inmultim cu D−1 si obtinem A−1. Timpul de executie este acelasi ca in cazul precedent si se folosescn3 procesoare.

Un alt algoritm de calcul paralel al inversei unei matrici este imparte si cucereste. Impartim matricea Ain blocuri:

A =

[A1 0A2 A3

]unde A1 este de dimensiune ⌈n

2 ⌉ × ⌈n2 ⌉, A1 si A3 sunt inferior triunghiulare, iar inversa se calculeaza cu

formula:

A−1 =

[A−1

1 0−A−1

3 A2A−11 A−1

3

]Avem astfel urmatorul algoritm:

Algorithm 1 Algoritm pentru calculul inversei unei matriciif n = 1 then

A−1 este evidentend ifif n ̸= 1 then

partitioneaza matricea Ainverseaza A1 si A3

inmulteste A−13 cu A2 pentru a obtine A−1

3 A2

inmulteste A−13 A2 cu A−1

1 pentru a obtine A−13 A2A

−11

end if

Pentru executia algoritmului intr-un timp bun avem nevoie de n3 procesoare.O metoda practica de rezolvare a sistemului Ax = b este metoda substitutiei inverse, care se obtine prin

paralelizarea metodei secventiale de rezolvare a sistemelor liniare.Presupunem ca matricea A este inferior triunghiulara si ca avem a i-a ecuatie a sistemului Ax = b:

ai1x1 + ai2x2 + · · ·+ ainxn = bi.

Algoritmul foloseste n procesoare (figurile 1 si 2). Presupunem ca la pasul i avem calculate valorilepentru x1, . . . , xi−1 si aj1x1 + aj2x2 + · · · + aji−1xi−1 pentru j ≥ i. Procesorul i calculeaza valoarea luixi folosind formula:

2

Page 57: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

xi =1

aii(bi − ai1x1 − · · · − aii−1xi−1).

Apoi fiecare procesor j cu j ≥ i + 1 calculeaza aj1x1 + aj2x2 + · · · + ajixi prin adaugarea lui ajixi

la expresia deja cunoscuta aj1x1 + aj2x2 + · · · + aji−1xi−1. Algoritmul se termina la pasul n cand toatevalorile x1, . . . , xn sunt calculate.

Figure 1: Implementarea metodei substitutiei inverse

Figure 2: Implementarea alternativa pentru substitutia inversa

3

Page 58: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Metoda substitutiei inverse poate fi folosita si pentru calculul inversei unei matrici dupa cum urmeaza.Deoarece AA−1 = I , atunci coloana i, xi, a lui A−1 satisface Axi = ei, unde ei este coloana i a matricei I.Astfel inversa este obtinuta prin rezolvarea a s sisteme de n ecuatii.

1.2 Sisteme tridiagonaleSe considera sistemul Ax = b, unde A este matrice tridiagonala (aij = 0 pentru |i − j| > 1). Un astfel desistem are forma:

g1x1 + h1x2 = b1 (2)

fixi−1 + gixi + hixi+1 = bi, i = 2, 3, . . . , n− 1 (3)

fnxn−1 + gnxn = bn (4)

unde gi sunt elementele de pe diagonala principala, fi si hi sunt sub si deasupara diagonalei principale.Exista mai multe metode de rezolvare a acestui sistem, una din ele este reductia para-impara (figura 3).Folosim notatia x0 = xn+1 = 0 si din ecuatia (3) obtinem:

xi =1

gi(bi − fixi−1 − hixi+1) (5)

Folosim ecuatia (5) in care inlocuim i cu i− 1 si cu i+ 1 si obtinem:

figi−1

(bi−1 − fi−1xi−2 − hi−1xi) + gixi +hi

gi+1(bi+1 − fi+1xi − hi+1xi+2) = bi

care simplificata devine:

−fifi−1

gi−1xi−2 + (gi −

hi−1figi−1

− hifi+1

gi+1)xi −

hihi+1

gi+1xi+2 = bi −

figi−1

bi−1 −hi

gi+1bi+1 (6)

Se obtine astfel un sistem tridiagonal mai mic pentru a carui rezolvare se poate aplica aceeasi metoda.

2 Metode directe de rezolvare a ecuatiilor liniareFie A o matrice patratica n× n. Cele mai multe metode de rezolvarea a unui sistem Ax = b intai transformasistemul pana cand matricea A devine triunghiulara si apoi rezolva noul sistem cu metoda substitutiei inverse.Aceste transformari constau in inmultiri ale sitemului cu matricile M1, . . . ,Mk. Se obtine astfel un nousistem ce poate fi rezolvat cu metoda substitutiei inverse:

Mk . . .M1x = Mk . . .M1b.

Daca se doreste calculul inversei matricii A, aceasta poate fi calculata folosind relatia:

A−1 = (Mk . . .M1A)−1(Mk . . .M1).

2.1 Eliminarea GaussianaFie C0 = A si Ci = M i . . .M1A. Presupunem ca Ci−1 a fost deja calculata pentru 1 ≤ i ≤ n − 1 si areforma din figura 4.

Calculam matricea M i astfel incat Ci = M iCi−1 are aceeasi proprietate. Pentru aceasta consideram

M i = I−N i, unde N i are forma din figura 5, iar valorile nevide se calculeaza cu relatia N iji =

Ci−1ji

Ci−1ii

, j > i.

Se poate verifica faptul ca matricea Ci are forma din figura 4 si se poate calcula folosind direct relatia:

4

Page 59: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 3: Reductia para-impara

Cijk = Ci−1

jk −n∑

l=1

N ijlC

i−1lk = Ci−1

jk −N ijiC

i−1ik = Ci−1

jk −Ci−1

ji

Ci−1ii

Ci−1ik

In cazul in care Ci−1ii = 0, se poate folosi pivotrea partiala sau totala.

2.2 Triangularizarea matricilor folosind metoda rotatiilor GivensO diferenta notabila fata de triangularizarea prin metoda Gauss este aceea ca matricile M i sunt ortogonale(||M ix|| = ||x||, ∀x ∈ Rn). Aceasta relatie este echivalenta cu (M i)tM i = I .

Fie C o matrice data astfel incat Cil = Cjl = 0, pentru l < k si Cjk ̸= 0, deci are forma din figura 6.Se considera o matrice M, numita rotatie Givens, de forma 7 astfel incat:

c =Cik√

C2ik + C2

jk

s =Cjk√

C2ik + C2

jk

Acesta matrice M are urmatoarele proprietati:

• Toate liniile din MC sunt egale cu liniile din C mai putin liniile i si j;

• MCil = MCjl = 0, pentru l < k si MCjk = 0

3 Metode directe de calcul a inversei unei matriciFie A o matrice patratica nesingulara (n× n). Consideram polinomul caracteristic al matricei A:

5

Page 60: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 4: Matricea C

Figure 5: Matricea N

P (λ) = det(λI −A) =n∏

i=1

(λ− λi)

unde λi, i = 1, n, sunt valorile proprii ale matricei A. Fie ci, i = 1, n, coeficientii polinomului caracter-istic P (λ).

P (λ) = λn + c1λn−1 + · · ·+ cn−1λ+ cn

Se poate observa ca cn = (−1)nn∏

i=1

λi si deci matricea A este inversabila daca cn este nenul.

Folosind proprietatea lui Cayley-Hamilton (P (A) = An + c1An−1 + · · ·+ cn−1A+ cnI = 0), obtinem

ca inversa A−1 se calculeaza cu formula:

A−1 = − 1

cn(An−1 + c1A

n−2 + · · ·+ cn−1I)

Metoda directa de calcul a inversei unei matrici este data de urmatorul algoritm:

6

Page 61: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 6: Matricea C

Figure 7: Matricea M

7

Page 62: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Algorithm 2 Calculul inversei unei matriciCalculeaza Ak pentru k = 2, . . . , n

Calculeaza sk pentru k = 1, . . . , n folosind relatia sk =n∑

n=1λki = tr(Ak)

Rezolva sistemul urmator determinand c1, . . . , cn:

1 0 . . . . . . . . . 0

s1 2 0...

......

...... . . . . . . . . .

...sk−1 . . . s1 k 0 . . .

......

......

... 0sn−1 . . . sk−1 . . . s1 n

c1............cn

= −

s1............sn

Determina A−1 folosind relatia:

A−1 = − 1

cn(An−1 + c1A

n−2 + · · ·+ cn−1I)

8

Page 63: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Metode iterative pentru rezolvarea sistemelor de ecuatii liniare

1 Metode iterative clasiceMetodele iterative sunt intens folosite, in special pentru rezolvarea de probleme mari, cum sunt cele dediscretizare a ecuatiilor cu derivate partiale.

Fie A o matrice patratica reala (A ∈ Rn×n) si b un vector din Rn. Se considera sistemul:

Ax = b (1)

unde x este un vector necunoscut ce trebuie determinat. Presupunem ca matricea A este inversabila si decisistemul are solutie unica. Presupunem ca a i-a ecuatie a sistemului este:

n∑j=1

aijxj = bi (2)

Metoda lui Jacobi de rezolvare a sistemului (1) este data de urmatorul algoritm: pornind cu x(0) ∈ Rn,se calculeaza x(t) pentru t = 1, 2, . . . cu ajutorul relatiei:

xi(t+ 1) = − 1

aii(∑j ̸=i

aijxj(t)− bi) (3)

Acest algoritm construieste un sir care converge la solutia sistemului. Exista situatii in care sirul construitcu metoda lui Jacobi nu este convergent (figura 1).

1

Page 64: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 1: Convergenta metodei Jacobi

2

Page 65: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Metoda Seidel-Gauss este urmatoarea: pornind cu x(0) ∈ Rn, se calculeaza x(t) pentru t = 1, 2, . . . cuajutorul relatiei:

xi(t+ 1) = − 1

aii(∑j<i

aijxj(t+ 1) +∑j>i

aijxj(t)− bi) (4)

In figura 2 este ilustrata o situatie de convergenta a algoritmului si una de divergenta.Variatii ale metodelor Jacobi si Seidel-Gauss se obtin folosind un scalar diferit de zero (numit parametru

de relaxare):

xi(t+ 1) = (1− γ)xi(t)−γ

aii(∑j ̸=i

aijxj(t)− bi) (5)

xi(t+ 1) = (1− γ)xi(t)−γ

aii(∑j<i

aijxj(t+ 1) +∑j>i

aijxj(t)− bi) (6)

Metoda lui Richardson se obtine prin rescrierea sistemului (1) sub forma x = x−γ(Ax−b) si construindsirul iterativ folosind formula:

x(t+ 1) = x(t)− γ(Ax(t)− b) (7)

O metoda mai generala se obtine folosind o matrice inversabila B, transformand sistemul intr-o formaechivalenta:

x = x−B(Ax− b)

si aplicand formula:

x(t+ 1) = x(t)−B(Ax(t)− b)

Theorem 1.1 Daca sirul {x(t)} generat de o metoda iterativa este convergent, atunci el converge la o solutiea sistemului Ax = b.

3

Page 66: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 2: Convergenta metodei Seidel-Gauss

4

Page 67: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

1.1 Implementari paralele ale metodelor iterative clasicePresupunem ca avem n procesoare si ca fiecare procesor i calculeaza valoarea lui xi(t) la fiecare iteratie t.Fiecare procesor i cunoaste intrarile de pe linia i a matricei A.

Pentru a calcula valoarea lui xi(t+1) procesorul i trebuie sa cunoasca valorile lui xj(t). Astfel se transmittoate valorile xj(t) de la iteratia t catre procesorul i chiar daca multe valori aij sunt nule.

In practica, p numarul de procesoare este adesea semnificativ mai mic decat n numarul de variabile. Inacest caz, mai multe variabile pot fi atribuite aceluiasi procesor.

Un punct final de interes se refera la oprirea algoritmului iterativ. De obicei conditia de terminare estelegata de evaluarea unei expresii, cum ar fi ||Ax(t)− b||, si de a opri algoritmul daca valoarea sa este destulde mica.

De exemplu, sa presupunem ca || · || este norma maximului. La fiecare iteratie, fiecare procesor cal-culeaza valoarea max |Axi−bi|. Aceste valori pot fi apoi comparate cu ajutorul unui arbore, fiecare procesorpropagand spre radacina arborelui propria valoare si valorile pe care le-a primit.

2 Analiza convergentei metodelor iterative clasiceFie D o matrice diagonala astfel incat dii = aii, si fie B o matrice cu zero pe diagonala principala B = A−D.Presupunand ca elementele de pe diagonala principala a matricei A sunt nenule, algoritmul Jacobi poate fiscris matriceal sub forma:

x(t+ 1) = D−lBx(t) +D−lb. (8)

Similar, algoritmul Jacobi modificat poate fi scris sub forma:

x(t+ 1) = [(1− γ)IγD−lB]x(t) + γD−lb. (9)

Pentru cazul in care se lucreaza cu algoritmul Seidel-Gauss sau Seidel-Gauss modificat, matricea A sepoate scrie sub forma A = L+D + U , unde L este strict inferior triunghiulara, D este diagonala, iar U estestrict superior triunghiulara. Astfel se obtin formulele echivalente:

x(t+ 1) = (1− γ)x(t)− γD−1[Lx(t+ 1) + Ux(t)− b] (10)

x(t+ 1) = (I + γD−1L)[(1− γ)I − γD−1U ]x(t) + γ(I + γD−1L)−1D−1 − b (11)

Analog, metoda Richardson poate fi descrisa cu ajutorul formulei:

x(t+ 1) = (I − γA)x(t) + γb (12)

Ecuatiile anterioare pot fi scrise unitar sub forma echivalenta:

x(t+ 1) = Mx(t) +Gb (13)

Deci pentru a stabili convergenta metodelor iterative ar trebui doar sa analizam convergenta sirului generatcu ajutorul relatiei precedente.

Theorem 2.1 Presupunand ca matricea I −M este inversabila, fie x∗ care verifica x∗ = Mx∗ + Gb si fiesirul generat cu ajutorul iteratiei x(t + 1) = Mx(t) + Gb. Atunci sirul converge la x∗ pentru orice x(0)daca si numai daca ρ(M) < 1.

Definition 2.1 O matrice patratica A este dominanta pe linii daca:∑j ̸=i

|aij | < |aii|,∀i. (14)

5

Page 68: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Remark 2.1 Daca A este dominanta pe linii, atunci matricea M asociata iteratiei Jacobi are raza spectralaρ(M) < 1.

Theorem 2.2 Daca A este dominanta pe linii, atunci metoda Jacobi converge la o solutie a sistemului Ax =b.

3 Metoda gradientului conjugatConsideram un sistem liniar de ecuatii: Ax = b. Presupunem ca A este o matrice patratica, simetrica sipozitiv definita.

Metodele directiilor conjugate sunt motivate de dorinta de a accelera viteza de convergenta a metodelorclasice iterative pentru aceasta clasa de probleme. In timp ce ele sunt garantate pentru a gasi solutia dupacel mult n iteratii, ele sunt privite ca fiind cele mai bune metode iterative, deoarece, de obicei, mai putin den iteratii sunt executate, in special pentru probleme mari. Aceste metode sunt, de fapt, aplicabile si pentruprobleme de optimizare noncuadratice. Pentru astfel de probleme, in general, se termina procesul dupa unnumar finit de iteratii, dar, atunci cand sunt aplicate in mod corespunzator, au convergenta atractiva si rata deconvergenta convenabila.

Presupunem ca b = 0. Functia cost a metodei este calculatat cu ajutorul relatiei:

F (x) =1

2x′Ax

Aceasta functie este strict convexa (deoarece A este pozitiv definita) si ia valoarea minima pentru x = 0,care este deasemenea solutia unica a sistemului Ax = 0.

Forma generala a metodei este data de relatia:

x(t+ 1) = x(t) + γ(t)s(t), t = 0, 1, . . .

unde s(t) este directia, iar γ(t) este un parametru definit de relatia F (x(t) + γ(t)s(t)) = minγ∈R

F (x(t) +

γs(t)).Vectorul directie este ales astfel incat:

s′(t)As(t) = 0, daca t ̸= r.

Facem urmatoarea notatie: g(t) = Ax(t) = ∇F (x(t)).

Proposition 3.1 Presupunem ca s(0), s(1), . . . , s(t) sunt diferite de zero si sunt conjugate. Atunci:

• Vectorii s(0), s(1), . . . , s(t) sunt liniar independenti.

• Daca 0 ≤ i ≤ k ≤ t atunci avem relatia:

g′(k + 1)s(i) = 0

• Vectorul x(k) satisface relatia F (x(k + 1)) ≤ F (x(k)), pentru orice k ≤ t.

3.1 Descrierea algoritmuluiSe construieste s(t) folosind formula:

s(t) = −g(t) +

t−1∑i=0

cis(i)

6

Page 69: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

unde coeficientii ci sunt alesi astfel incat s(t) este conjugat cu s(0), s(1), . . . , s(t − 1). Presupunand cas(0), s(1), . . . , s(t) sunt deja conjugate, pentru j = 0, . . . , t− 1 avem:

0 = s′(t)As(j) = −g′(t)As(j) +t−1∑i=0

cis′(i)As(j) = −g′(t)As(j) + cjs

′(j)As(j)

de unde se obtine:

cj =g′(t)As(j)

s′(j)As(j).

Deoarece g(j + 1)− g(j) = A(x(j + 1)− x(j)) = γ(j)As(j), obtinem:

cj =g′(t)(g(j + 1)− g(j))

s′(j)(g(j + 1)− g(j)).

Se observa astfel ca g(j) este o combinatie liniara intre s(0), s(1), . . . , s(j), si se obtine g′(t)g(j) = 0.Daca j < t− 1 atunci cj = 0, si obtinem:

s(t) = −g(t) + β(t)s(t− 1)

unde

β(t) =g′(t)g(t)

s′(t− 1)(g(t)− g(t− 1))=

g′(t)g(t)

g′(t− 1)g(t− 1)

Avem relatia:

s′(t)A(x(t) + γ(t)s(t)) = 0 =⇒ γ(t) = − s′(t)g(t)

s′(t)As(t).

Iteratia urmatoare in cadrul algoritmului gradientului conjugat x(t+ 1) este calculata prin relatia:

x(t+ 1) = x(t) + γ(t)s(t).

Proposition 3.2 Algorimul gradientului conjugat se termina in cel mult n pasi, asta inseamna ca exista unt ≤ n astfel incat g(t) = 0 si x(t) = 0.

3.2 Implementarea paralelaPresupunand ca la inceputul iteratie t, x(t), s(t− 1) si g′(t− 1)g(t− 1) au fost deja calculate, avem nevoiepentru sa evaluam vectorii g(t) = Ax(t), produsul vectorial g′(t)g(t), sa determinam β(t) si s(t), apoi saevaluam As(t), si in cele din urma sa calculam produsele vectoriale s′(t)(As(t)), s′(t)g(t).

Daca n procesoare sunt disponibile, este natural sa lasam procesorul i sa controleze componenta i avectorilor x(t), s(t) si g(t). Produsele scalare sunt calculate prin lasarea procesorului i sa calculeze produsulcomponentelor i si apoi acumuleze sumele partiale de-a lungul un arbore de procesoare. Aceasta este oacumulare intr-un singur nod si are nevoie de timp proportional cu diametrul retelei de interconectare. Apoi,valorile calculate ale produselor scalare sunt difuzate la toate procesoarele. Noi presupunem acum ca fiecareprocesor cunoaste intrarile dintr-o linie diferita a matricei A. Apoi, o matrice-vector Ax(t) poate fi calculataprin difuzarea vectorului x(t) si avand procesorul i care calculeaza produsul scalar al lui x cu linia i a matriceiA. Alternativ, procesor i ar putea calcula [A]jixi pentru fiecare j, si aceste valori ar putea fi reproduse la fiecareprocesor j, cu sume partiale formate de-a lungul drumului.

In cazul in care mai putin de n procesoare sunt disponibile, problemele implicate sunt aceleasi cu exceptiafaptului ca exista mai multe componente, si mai multe linii ale matricei A, alocate pentru fiecare procesor.

7

Page 70: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Algoritmi sincroni

1 Inversarea matricilorVom considera o metoda care se bazeaza pe algoritmul clasic al lui Newton pentru imbunatatiri iterative aleaproximarii inversei unei matrici patratice A. Aceasta metoda a fost obtinuta din dorinta de a calcula inversaunei matrici intr-un timp cat mai mic. Timpul de executie redus se bazeaza pe o proprietate dinstinctiva ametodelor Newton, convergenta patratica, convergenta ce se realizeaza in acelasi timp cu sirul {ρ2t}, unde ρeste o constanta pozitiva mai mica decat 1. Aceasta este mult mai rapida decat convergenta geometrica ce serealizeaza in aceasi timp cu sirul {ρt}.

Se da o matrice B de aceeasi dimensiune cu matricea A si definim o matrice reziduala R(B) = I −BA. Putem spune ca, in linii mari, R(B) masoara cat de departe este B de inversa matricei A. Consideramurmatorul algoritm pentru calculul inversei unei matrici:

Algorithm 1 Calculul inversei unei matriciFie B0 astfel incat ||I −B0A||2 < 1Calculeaza iterativ Bk utilizand formula

Bk+1 = (I +R(Bk))Bk = 2Bk −BkABk

Aceasta iteratie poate fi interpretata ca o metoda Newton de rezolvare a ecuatiei

X−1 −A = 0.

Avem urmatoarele relatii echivalente:

R(Bk+1) = I −Bk+1A = I − (I +R(Bk))BkA

R(Bk+1) = (I −BkA)−R(Bk)BkA = R(Bk)(I −BkA) = (R(Bk))2

Iterativ obtinem astfel ca R(Bk) = (R(B0))2k . Avem astfel inegalitatea:

||R(Bk)||2 = ||(R(B0))2k ||2 ≤ ||R(B0)||2

k

2

Deoarece am presupus ca ||R(B0)|| < 1, obtinem ca R(Bk) converge foarte repede la 0, echivalent cufaptul ca Bk converge foarte repede la A−1.

Deci algoritmul converge cu succes pentru o alegere favorabila a lui B0. Putem alege B0 folosind formula:

B0 =A′

tr(A′A)

unde tr(A′A) este urma matricei A′APentru orice matrice nesingulara A facem urmatoarea notatie:

κ(A) = ||A||2||A−1||2.

1

Page 71: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Aceasta valoare este numita numarul conditionat al lui A si joaca un rol important in analiza numericapentru a determina dificulatea de calcul al inversei matriei A.

Fie λ1 ≤ λ2 ≤ · · · ≤ λn valorile proprii ale matriceai A′A. Aceste valori proprii sunt reale si nenegativedeoarece matricea A′A este simetrica si nenegativa. Folosind valorile proprii obtinem: κ(A) = (λn

λ1)1/2

Proposition 1.1 Daca B0 = A′

tr(A′A) atunci ||I −B0A||2 ≤ 1− 1nκ2(A) .

Proof.Avem tr(A′A) = λ1 + · · ·+ λn.Rezulta ca a i-a valoare proprie ρi a matricei I −B0A = I − A′A

tr(A′A) este egala cu 1− λi

λ1+···+λn.

Deoarece λi ≤ λ1 + · · ·+ λn obtinem ca ρi ≥ 0.Deoarece λ1 + · · ·+ λn ≤ nλn obtinem:

ρi ≤ 1− λi

nλn≤ 1− λ1

nλn= 1− 1

nκ2(A)

Rezulta ca ρ(I −B0A) = maxi

|ρi| ≤ 1− 1nκ2(A) .

Deoarece I − B0A este matrice simetrica si ||I − B0A||2 = ρ(I − B0A) obtinem ||I − B0A||2 ≤1− 1

nκ2(A) .Algoritmul poate fi folosit si pentru rezolvarea sistemelor liniare Ax = b. Se calculeaza intai matricea

A−1 si apoi se obtine x = A−1b.

2 Lantul lui MarkovLanturile Markov sunt larg folosite ca modele probabilistice. In aplicatii, lanturile Markov sunt adesea intal-nite in spatii de dimensiuni mari. Adesea se doreste sa se calculeze distributia probabilista a acestor lanturi,si acest calcul necesita o implementare paralela.

Fie P matricea de tranzitie probabilista a unui lant Markov cu n stari. Matricea P are urmatoarele propri-etati:

• P ≥ 0

•n∑

j=1

pij = 1

Definition 2.1 Orice matrice care indeplineste cele doua proprietati se numeste matrice stocastica.

Definition 2.2 Un vector linie nenegativ π∗ care are suma componentelor egala cu 1 si are proprietateaπ∗ = π∗P se numeste distributie invarianta a lantului Markov asociata lui P.

Consideram urmatorul algoritm pentru determinarea distributiei invariante: consideram π(0) ≤ 0 cusuma componentelor egala cu 1, si calculam iterativ

π(t+ 1) = π(t)P

echivalent cu

π(t) = π(0)P t, t ≥ 0.

Dandu-se o matrice stocastica P, avem un graf direct orientat G = (N,A), unde N este multimea de starisi A = {(i, j), pij > 0, i ̸= j} multimea de arce (multimea tuturor tranzitiilor de stari care au o probabilitatepozitiva).

2

Page 72: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Definition 2.3 Spunem ca matricea P este ireductibila daca pentru orice i, j ∈ N exista un drum pozitiv ingraf de la i la j.

Definition 2.4 Matricea P se numeste periodica daca exista un k > 1 si niste multimi disjuncte N0, . . . , Nk−1

astfel incat daca i ∈ Nl si pij > 0 atunci j ∈ Nl+1(modk).

Definition 2.5 Spunem despre multimea P ca este aperiodica daca ea nu este periodica.

Definition 2.6 Spunem despre multimea P ca este primitiva daca exista un intreg pozitiv t astfel incat P t > 0.

Cateva exemple pot fi observate in figura 1.

Figure 1: Graf direct orientat asociat unei matrici stocastice

Proposition 2.1 O matrice stocastica P este primitiva daca si numai daca este ireductibila si aperiodica.

Proposition 2.2 Fie P o matrice stocastica.

• raza spectrala a lui P este egala cu 1 (ρ(P ) = 1);

• daca π este un vector cu suma elementelor egala cu 1, atunci vectorul πP are aceeasi proprietate.

Convergenta sirului iterativ π(t+ 1) = π(t)P se obtine folosind urmatoarea teorema:

Theorem 2.1 Fie P o matrice stocastica primitiva. Atunci:

• Exista un unic vector π∗ astfel incat π∗ = π∗P sin∑

i=1

π∗i = 1.

• Limita lui P t cand t tinde la infinit exista si este o matrice cu toate randurile egale cu π∗.

• Dacan∑

i=1

πi(0) = 1, atunci iteratia π(t+ 1) = π(t)P converge la π∗.

Procesul iterativ π = πP admite o implementare paralela in care fiecare procesor i calculeaza compo-nenta i a vectorului π si la fiecare pas comunica valoarea nou calculata procesorelor j pentru care pij ̸= 0.Se presupune ca fiecare procesor j cunoaste elementele coloanei j din matricea P.

O alta implementare se obtine daca fiecare procesor j cunoaste elementele liniei j din matricea P.

3

Page 73: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Metode iterative pentru probleme neliniare - contractii

Problemele neliniare sunt in general rezolvate prin metode iterative si analiza convergentei acestor metodeeste o problema importanta.

1 ContractiiMetodele iterative pot fi scrise in general sub forma:

x(t+ 1) = T (x(t)), t = 0, 1, . . .

unde T este o functie definita pe o submultime X a spatiului Rn cu valori in X si care are proprietatea:

||T (x)− T (y)|| ≤ α||x− y||, ∀x, y ∈ X (1)

unde α ∈ [0, 1).

Definition 1.1 O astfel de functie ce are proprietatea 1 se numeste contractie, iar procesul iterativ asociat eise numeste iteratie de tip contractie.

Fie o functie T : X → X . Orice vector x∗ ∈ X care verifica T (x∗) = x∗ se numeste punct fix al lui T siiteratia x = T (x) poate fi privita ca un algoritm de determinare a unui astfel de punct fix.

Se observa ca functiile contractie sunt continui: daca avem un sir {x(t)} care converge la un punct x∗ ∈ Xsi functia T este continua in punctul x∗, atunci x∗ este un punct fix al lui T.

Ca o alternativa la functiile de tip contractie se poate presupune ca o functie T : X → X are un punct fixx∗ ∈ X si are proprietatea

||T (x)− x∗|| ≤ α||x− x∗||, ∀x ∈ X (2)

unde α ∈ [0, 1).Se poate observa ca inegaliatea 2 este mai slaba decat conditia de contractie 1.

Definition 1.2 Orice functie ce are proprietatea 2 se numeste pseudocontractie.

In figura 1 se observa un exemplu de contractie (T : R → R, T (x) = x2 , α = 1

2 , T(0)=0) si un exemplude pseudocontractie (T : [0, 2] → [0, 2], T (x) = max{0, x− 1}.

Urmatoarea proprietate ne arata ca o contractie are un unic punct fix si ca iteratia corespunzatoare con-tractiei converge la acest unic punct fix.

Proposition 1.1 Presupunem ca T : X → X este o contractie de constanta α ∈ [0, 1) si ca X este osubmultime inchisa a lui Rn. Atunci:

• Functia T are un unic punct fix x∗ ∈ X;

• Pentru orice vector initial x(0) ∈ X , sirul {x(t)} generat de x(t+ 1) = T (x(t)) converge geometricla x∗. In particular:

||x(t)− x∗|| ≤ αt||x(0)− x∗||, ∀t ≥ 0

1

Page 74: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 1: Exemplu de contractie si pseudocontractie

Proposition 1.2 Presupunem ca X ⊂ Rn si ca functia T : X → X este o pseudocontractie de constantaα ∈ [0, 1) si cu un punct fix x∗ ∈ X . Atunci T nu mai are alt punct fix si sirul {x(t)} generat de x(t+ 1) =T (x(t)) satisface:

||x(t)− x∗|| ≤ αt||x(0)− x∗||, ∀t ≥ 0

pentru orice alegere a vectorului initial x(0) ∈ X . In particular sirul {x(t)} converge la x∗.

Proposition 1.3 Daca X ⊂ Rn este nevida, convexa si compacta si daca T : X → X este o functiecontinua, atunci exista x∗ ∈ X astfel incat T (x∗) = x∗ (un exemplu se poate vedea in figura 2).

Figure 2: Teorema de punct fix

2

Page 75: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

2 Contractii peste produse carteziene de multimi

In continuare presupunem ca X =m∏i=1

Xi, unde Xi ∈ Rni si n1 + · · ·+ nm = n. Orice vector x ∈ X poate

fi descompus in x = (x1, . . . , xm), cu xi ∈ Xi.Fie T : X → X o contractie si fie Ti : X → Xi componenta i a functie T:

T (x) = (T1(x), . . . , Tm(x))

Se poate observa ca:

||Ti(x)− Ti(y)||i ≤ maxj

||Tj(x)− Tj(y)||j = ||T (x)− T (y)|| ≤ α||x− y||

In figura 3 se poate observa o simulare a iteratiei Seidel-Gauss pentu contractii bloc.

Figure 3: Convergenta iteratiei Seidel-Gauss pentru contractii bloc

Ca o alternativa la gasirea unui punct fix al functiei T putem cauta o solutie a sistemului x = T (x). Acestsistem poate fi descompus in m sisteme mai mici:

xi = Ti(x1, . . . , xm), i = 1, . . . ,m

care pot fi rezolvate simultan.Fie Ri(x) multimea tuturor solutiilor ecuatiei i din sistem definita prin:

Ri(x) = {yi ∈ Xi, yi = Ti(x1, . . . , xi−1, yi, xi+1, . . . , xm)}

Se da un vector x(t) ∈ X , componenta i a vectorului urmator xi+1 este aleasa ca solutie a ecuatiei i asistemului:

xi(t+ 1) ∈ Ri(x(t))

Proposition 2.1 Presupunem ca X este inchisa si ca T : X → X este o contractie bloc. Atunci multimeaRi(x) are un singur element pentru fiecare i si pentru fiecare x ∈ X .

3

Page 76: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Definim functia Qi : X → Xi prin egalarea lui Qi(x) cu elementul unic al multimii Ri(x). Si unimaceste functii in cadrul functiei Q : X → X

Q(x) = (Q1(x), . . . , Qm(x))

Metoda de rezolvare a sistemului poate fi descrisa cu ajutorul algoritmului (figurile 4 si 5):

x(t+ 1) = Q(x(t)), t = 0, 1, . . .

Figure 4: Functia Q

Proposition 2.2 Daca T : X → X este o contractie bloc, atunci Q este deasemenea o contractie bloc. DacaX este inchisa, atunci metoda x(t+ 1) = Q(x(t)) converge la un punct fix al lui T.

Acelasi rezultat se obtine si in cazul in care T este o pseudocontractie si multimea X este inchisa siconvexa.

Rezultate similare cu cele de mai sus se obtin si pentru functii monotone.

4

Page 77: Notiuni introductive - popirlan.ropopirlan.ro/cris/ecnpd/cursuri.pdfNotiuni introductive 1 Arhitecturi paralele si distribuite Pentru a realiza un studiu al metodelor numerice paralele

Figure 5: Metoda Seidel-Gauss bazata pe functia Q pentru rezolvarea sistemului x = T (x)

5


Recommended