4.2.2 Operationen mit dünnbesetzten Matrizen im Vektorrechner

Verschiedene Matrixspeicherungen
$\Longrightarrow$ verschiedene Implementierungen für die gleiche Operation.

Bsp.: $\underline{y} \,=\, A_{n\times n}\cdot\underline{x}$ in CRS 

DO $i \,=\, 1,\, n$

		      $y(i) \;:=\; 0$


DO $j \,=\, row\_ptr(i),\, row\_ptr(i+1)-1$

				     $y(i)\;=\;y(i) + val(j) \ast x(col\_ind(j))$


END DO 

END DO

Wie obige Doppelschleife vektorisieren ?

1.

Zeile $i$ in $A$ aufblasen und mit Nullen füllen
$\Longrightarrow$ innere Schleife $\widehat{=}$ $y_i \;=\; y_i + (A_{i,\ast},\underline{x})$
Widerspruch ist bei 3D-FEM Matrix kompletter Unsinn.

2.

GATHER/ SCATTER (sammeln/verteilen) auf der Matrixzeilen mittels des Indexvektor

     y = 0 

DO $i \,=\, 1,\, n$

		      $length \;:=\; row\_ptr(i+1)+1 - row\_ptr(i)\hspace{2em}$
[Register mit Vektorlänge]

		      $j1 \;:=\; row\_ptr(i)$

		      $y(i) \;:=\;y(i) + val(j1) \ast x(col\_ind(j1))$


END DO

Der Vektorrechner lädt das Indexfeld $col\_ind(j1)$ und adressiert damit den entsprechenden Eintrag in $\underline{x}$ ( GATHER). Analog läßt sich ein entsprechender Vektor speichern ( SCATTER).

• Ähnlich werden Matrix-Matrix-Operationen ausgeführt.
• Durch die indirekte Adressierung sind Speicherbankkonflikte und Page faults nicht zu vermeiden.
• Die Parallelisierung des Matrix-Vektor-Produkts für dünnbesetzte Matrizen sollte die (lokale) Besetztheitsstruktur der Matrix ausnutzen (Abschnitt 4.3).