Übungsblatt 1

18.04.2019 / B. Leder

Wissenschaftliches Rechnen III / CP III

Aufgabe 1.1: Methode der Konjugierten Gradienten in Matlab

Implementieren Sie den Algorithmus “Methode der Konjugierten Gradienten” zur Lösung des Gleichungssystems \[Ax=b\] in Matlab. Die Iteration ist vollständig bestimmt durch \[\begin{align} r^{(0)} &= b - A x^{(0)} \nonumber \\ p^{(0)} &= r^{(0)} \nonumber \\ \alpha_k &= \frac{{r^{(k)}}^T r^{(k)}}{{p^{(k)}}^T A p^{(k)}} \nonumber \\ x^{(k+1)} &= x^{(k)} + \alpha_k p^{(k)} \nonumber \\ r^{(k+1)} &= r^{(k)} - \alpha_k A p^{(k)} \nonumber \\ \beta_k &= \frac{{r^{(k+1)}}^T r^{(k+1)}}{{r^{(k)}}^T r^{(k)}} \nonumber \\ p^{(k+1)} &= r^{(k+1)} + \beta_k p^{(k)} \nonumber \end{align}\]

Die Herleitung der einzelnen Gleichungen erfolgt in der Vorlesung. Die Iteration erzeugt eine Folge von Vektoren \(\{x^{(k)},\,r^{(k)},\,p^{(k)}\}\) (Approximation der Lösung, Residuum, konjugierte Suchrichtung) und Koeffizienten \(\{\alpha_k,\, \beta_k\}\) mit \(k=0,1,2,...\).

Die Anwendung der Matrix \(A\) auf einen Vektor \(A u\) soll explizit, also als Matrix-Vektor-Multiplikation erfolgen. Die Matrix \(A\) soll dem diskreten Laplace-Operator mit Dirichlet-Randbedingungen in zwei Dimensionen entsprechen (\(A \sim -\Delta\)). Sie können dafür die auf der Website bereitgestellte Matlab-Routine laplace.m benutzen. Die Matrix wird hier als sparse erzeugt (help spdiags), so dass der Speicherbedarf \(O(n)\) ist und Matrix-Vektor-Multiplikation \(O(n)\) Operationen benötigt (\(n\): Matrixgröße). Sie können die Struktur von \(A\) mit spy(A) inspizieren.

Benutzen Sie als Abbruchbedingungen eine maximale Iterationszahl (zur Sicherheit) und eine gewünschte Reduktion der anfänglichen Residuumsnorm um einen Faktor \(\epsilon\)

\[||r^{(k)}||<\epsilon ||r^{(0)}||\]

Achten Sie bei der Implementierung auf Effizienz, d.h. vermeiden Sie die Wiederholung von Berechnungen (benutzen Sie Hilfsvariablen/Zwischenspeicher). Validieren Sie Ihr Programm mit einem geeigneten Test.

10 Punkte

Aufgabe 1.2: Konjugierten Gradienten in Matlab - GPU Version

Ändern Sie Ihr Programm aus Aufgabe 1.1 dahingehend, dass die Iteration vollständig auf der GPU ausgeführt wird, siehe Aufgabe 0.1. Die Erstellung von \(A\) und \(b\) soll weiterhin in der CPU geschehen. Beachten Sie, dass alle Variablen für alle Zwischenschritte der Iteration im Device-Speicher liegen (siehe isa(A,'gpuArray'). Validieren Sie die GPU Version.

Benutzen Sie die CPU und die GPU Version zur Bestimmung von

  1. Latenz für das Kopieren zwischen Host und Device
  2. Laufzeiten von CPU und GPU Version
  3. Speedup mit und ohne Latenzen

Messen Sie diese Größen in Abhängigkeit der Matrixgröße \(n\leq 1024\). Verwenden Sie für CPU und GPU Version das gleiche \(\epsilon\). Plotten und interpretieren Sie die Ergebnisse.

10 Punkte

Aufgabe 1.3 (freiwillig): Eigenschaften der Methode der Konjugierten Gradienten

Beweisen Sie die folgenden Eigenschaften der Residuen \(r^{(k)}\) und Suchrichtungen \(p^{(k)}\), die durch die Iteration in Aufgabe 1.1 erzeugt werden: \[{r^{(k)}}^T p^{(i)} = 0 \quad i<k\] \[{p^{(k)}}^T A p^{(i)} = 0 \quad i\neq k\] \[{r^{(k)}}^T r^{(i)} = 0 \quad i\neq k\] \[{p^{(k)}}^T r^{(k)} = {r^{(k)}}^T r^{(k)} \] Hinweis: Beweisen Sie die ersten beiden Zeilen durch Induktion.

0 Punkte