Automatische Generierung von effizienten, parallelen Implementierungsvarianten für ODE-Verfahren aus deren Datenflussgraphen mit Kernelfusion und Tiling für GPUs und CPUs

Titelangaben

Werner, Tim:
Automatische Generierung von effizienten, parallelen Implementierungsvarianten für ODE-Verfahren aus deren Datenflussgraphen mit Kernelfusion und Tiling für GPUs und CPUs.
Bayreuth , 2023 . - 343 S.
( Dissertation, 2023 , Universität Bayreuth, Fakultät für Mathematik, Physik und Informatik)
DOI: https://doi.org/10.15495/EPub_UBT_00007002

Volltext

Link zum Volltext (externe URL):

Angaben zu Projekten

Abstract

Das Lösen von Anfangswertproblemen mit Hilfe von ODE-Verfahren ist in vielen wissenschaftlichen Bereichen und für ingenieurtechnische Anwendungen zwingend erforderlich. Bedauerlicherweise gibt es eine große Anzahl von relevanten Anfangswertproblemen, ODE-Verfahren und Prozessorarchitekturen, wobei potentiell jede Kombination von Anfangswertproblem, ODE-Verfahren und Prozessorarchitektur von einer speziell dafür optimierten Implementierungsvariante profitiert. Wegen dieser großen Anzahl an Kombinationen ist es sehr zeitaufwändig, per Hand eine optimierte Implementierungsvariante für jede davon zu erstellen. Eine solche optimierte Implementierungsvariante ist jedoch erforderlich, sofern man die Rechenleistung eines Prozessors gut ausnutzen möchte.


In dieser Arbeit befassen wir uns damit, wie wir automatisch optimierte Implementierungsvarianten von ODE-Verfahren für das Lösen von Anfangswertproblemen generieren können. Dabei ist es unser Ziel ein breites Spektrum an ODE-Verfahren, allgemeine Anfangswertprobleme sowie solche mit einer ausnutzbaren Zugriffsstruktur und GPUs sowie CPUs als optimierte Zielplattformen zu unterstützen.
Hierfür zeigen wir zunächst, dass sich viele häufig verwendete ODE-Verfahren, wie die RK-Verfahren, die PIRK-Verfahren oder die Peer-Verfahren, zu einem Datenflussgraphen bestehend aus Auswertungen der problemspezifischen rechten Seite (RHS), Linearkombinationen (LC) und Reduktionen (RED) als Grundoperationen abstrahieren lassen. Bei vielen Verfahren besteht der Datenflussgraph zusätzlich aus einer Abhängigkeitskette aus RHS→LC-Gliedern, in welcher voneinander unabhängige Grundoperationen zum selben RHS→LC-Glied zusammengefasst werden.

Diesen Datenflussgraphen können wir nun leicht zur Generierung einer Basisimplementierung für das Lösen eines Anfangswertproblems verwenden, welche eine Grundoperation nach der anderen über je ein paralleles Kernel berechnet. Ein solches Kernel nutzt die Systemparallelität innerhalb einer Grundoperation aus, um deren Berechnungen auf sämtlichen Rechenkernen der Zielplattform auszuführen. Für diese Basisimplementierung stellen wir von uns durchgeführte experimentelle Untersuchungen auf einer Test-GPU und einer Test-CPU vor. Diese experimentellen Untersuchungen zeigen, dass eine solche Basisimplementierung für dichtbesetzte Anfangswertprobleme mit einer hohen arithmetischen Intensität, wie das N-Körperproblem, bereits eine gute Performance erzielen kann. Dahingegen ist sie für dünnbesetzte Anfangswertprobleme mit einer geringen arithmetischen Intensität, wie unser Testproblem BRUSS2D oder viele andere Stencil-Probleme, durch die DRAM-Bandbreite gebunden und kann die Rechenwerke der Test-GPU und Test-CPU nur schlecht auslasten.


Deswegen verwenden wir als Nächstes die Technik der Kernelfusion, um für allgemeine, dünnbesetzte Anfangswertprobleme mit einer geringen arithmetischen Intensität die Datenwiederverwendung und damit die Performance zu erhöhen. Mit Hilfe dieser Technik der Kernelfusion verschmelzen wir jeweils mehrere Grundoperationen im Datenflussgraphen zu einem Kernel, um innerhalb eines solchen fusionierten Kernels die geladenen Daten über den On-Chip-Speicher wiederzuverwenden, die erzeugten Daten über den On-Chip-Speicher zu transferieren und ein unnötiges Zurückschreiben von den erzeugten Daten in den DRAM zu vermeiden. Dabei können wir jeweils sämtliche Grundoperationen in einem RHS→LC-Glied zu einem Kernel verschmelzen. Jedoch darf eine RHS-Operation bei einem allgemeinen Anfangswertproblem für die Auswertung einer ihrer Komponenten beliebig auf ihren Argumentvektor zugreifen. Deshalb können wir keine LC→RHS-Abhängigkeiten und Grundoperationen aus unterschiedlichen RHS→LC-Gliedern miteinander verschmelzen. Zusätzlich wenden wir vor der Fusion zwei Enabling-Transformationen an, und zwar die Splitting-Transformation zum Aufspalten von LC-Operationen und die Cloning-Transformation zur redundanten Berechnung einer Grundoperation. Als Nächstes stellen wir vor, wie sich eine Variante mit Kernelfusion sowohl auf GPUs als auch auf CPUs mit architekturspezifischen Optimierungen automatisch generieren lässt. Schließlich führen wir eine experimentelle Untersuchung der per Kernelfusion und per Enabling-Transformationen generierten Implementierungsvarianten von mehreren expliziten ODE-Verfahren mit dem Testproblem BRUSS2D auf mehreren modernen GPUs und CPUs durch. Diese Untersuchungen zeigen, dass sich auf den getesteten GPUs durch die Kernelfusion und Enabling-Transformationen eine Beschleunigung zwischen 1.62 und 3.32 gegenüber der jeweiligen Basisimplementierung des ODE-Verfahrens erzielen lässt, während sich auf den getesteten CPUs eine Beschleunigung zwischen 1.98 und 4.14 gegenüber der jeweiligen Basisimplementierung des ODE-Verfahrens erzielen lässt. Diese Beschleunigung ist darauf zurückzuführen, dass die Kernelfusion und die Enabling-Transformationen zusammen die zwischen Prozessor und DRAM transferierte Datenmenge je nach ODE-Verfahren auf 50% bis 25% gegenüber der jeweiligen Basisimplementierung reduzieren können. Allerdings zeigen die Experimente auch, dass die per Kernelfusion und Enabling-Transformationen generierten Implementierungsvarianten auf den getesteten modernen GPUs und CPUs weiterhin durch die DRAM-Bandbreite gebunden sind.

Da das beliebige Zugriffsmuster der RHS-Operationen für ein allgemeines Anfangswertproblem weitere Lokalitätsoptimierungen verhindert, spezialisieren wir uns, um solche weiteren Lokalitätsoptimierungen durchzuführen, als Nächstes auf dünnbesetzte rechte Seiten mit einer beschränkten Zugriffsdistanz. Bei diesen darf eine RHS-Operation für die Auswertung einer gegebenen Komponente nur auf diejenigen Komponenten des Argumentvektors, die sich innerhalb der beschränkten Zugriffsdistanz um der gegebenen Komponente herum befinden, zugreifen. Dies erlaubt es uns LC→RHS-Abhängigkeiten über eine partitionierte Fusion zu verschmelzen. Diese partitionierte Fusion können wir nun entlang der LC→RHS-Abhängigkeitskette des ODE-Verfahrens anwenden, wodurch sich zwei Tiling-Schemata ergeben, und zwar das trapezoidale Tiling-Schema und das hexagonale Tiling-Schema. Innerhalb eines Kernels mit Tiling kann ein Prozessor nun wieder analog, wie bei einem per Kernelfusion generierten Kernel, Daten wiederverwenden. Wir ermöglichen dieses Tiling nicht nur über die Stufen eines ODE-Verfahrens, sondern auch über dessen Zeitschritte. Ebenso stellen wir vor, wie mehrere Rechenkerne eines Prozessors an einem Tile zusammenarbeiten können, um darüber die maximale Tile-Größe und die Datenwiederverwendung zu erhöhen. Zusätzlich entwickeln wir einen Autotuning-Ansatz, um für beide Tiling-Schemata die optimalen Tile-Größen zu bestimmen. Als Nächstes gehen wir darauf ein, wie sich eine Variante mit Tiling sowohl auf GPUs als auch auf CPUs mit architekturspezifischen Optimierungen implementieren lässt. Schließlich führen wir für mehrere explizite ODE-Verfahren eine experimentelle Untersuchung der Implementierungsvarianten mit Tiling und mit BRUSS2D als Testproblem auf mehreren modernen GPUs und CPUs durch. Diese Untersuchung zeigt, dass sich durch das Tiling sowohl die Datenwiederverwendung als auch die Performance sowohl auf den getesteten GPUs als auch auf den getesteten CPUs gegenüber den Implementierungsvarianten mit Kernelfusion stark erhöhen lassen. So erzielen für die getesteten ODE-Verfahren die Implementierungsvarianten mit Tiling gegenüber einer durch Kernelfusion generierten Referenzimplementierung bei einer kleinen Zugriffsdistanz von 32 auf den getesteten GPUs eine Beschleunigung von 1.8 bis 8.2 und auf den getesteten CPUs eine Beschleunigung von 7.5 bis 14.1. Ebenso zeigen die Untersuchungen, dass es bereits bei kleinen bis mittelgroßen Zugriffsdistanzen für die Performance einer Variante mit Tiling sehr vorteilhaft ist, wenn ein Tile nicht nur von einem einzelnen Kern abgearbeitet wird, sondern wenn jeweils mehrere Rechenkerne an einem Tile zusammenarbeiten. Durch dieses Zusammenarbeiten von mehreren Rechenkernen an einem Tile lässt sich zum Beispiel auf der getesteten Volta-GPU und Skylake-CPU bei einer beschränkten Zugriffsdistanz von 131 072 noch eine Beschleunigung von 1.4 beziehungsweise 1.6 gegenüber der per Kernelfusion generierten Referenzimplementierung erzielen.


Wir stellen zudem ein von uns entwickeltes automatisiertes Framework vor, welches sowohl die Generierung als auch die Ausführung der Varianten der Basisimplementierung, der Varianten mit Kernelfusion sowie der Varianten mit Tiling anhand des Datenflussgraphen eines ODE-Verfahrens durchführt. Dieses Framework unterstützt sowohl GPUs als auch CPUs als optimierte Zielplattformen und implementiert zusätzlich eine Reihe von hilfreichen Funktionalitäten, wie ein automatisches Speichermanagement, eine NUMA-Unterstützung für CPUs und eine Autotuning-Funktionalität für das Tiling. Dieses Framework beinhaltet zudem eine DSL, mit welcher ein Benutzer leicht ODE-Verfahren zum Framework hinzufügen kann.

Abstract in weiterer Sprache

Solving initial value problems using ODE methods is essential in many scientific fields and for engineering applications. Unfortunately, there are a large number of relevant initial value problems, ODE methods and processor architectures, with potentially every combination of initial value problem, ODE method and processor architecture benefiting from a specially optimised implementation variant. Because of this large number of combinations, it is very time-consuming to manually create an optimised implementation variant for each of them. However, such an optimised implementation variant is necessary if one wants to make good use of the computing power of the processor.


In this thesis, we address how we can automatically generate optimised implementation variants of ODE methods for solving initial value problems. Our goal is to support a wide range of ODE methods, general initial value problems with an arbitrary memory access pattern, initial value problems with limited access distance, and GPUs as well as CPUs as optimised target platforms. To this end, we first show that many commonly used ODE methods, such as RK methods, PIRK-Methods or Peer-Methods, can be abstracted by a dataflow graph consisting of evaluations of the problem-specific right-hand side (RHS), linear combinations (LC) and reductions (RED) as basic operations. For many methods, the data flow graph additionally consists of a dependency chain of RHS→LC-links, in which basic operations being independent of each other are inserted into the same RHS→LC-link.

We can now easily use this data flow graph to generate a basic implementation for solving an initial value problem, which computes one basic operation at a time via one parallel kernel each. Such a parallel kernel exploits the system parallelism within a basic operation to perform its computations on all computational cores of the target platform. For this basic implementation, we conduct an experimental evaluation on a test GPU and test CPU. This experimental evaluation shows that such a basic implementation can already achieve a good performance for densely populated initial value problems with a high arithmetic intensity, such as the N-body problem. On the other hand, for sparse initial value problems with a low arithmetic intensity, such as our test problem BRUSS2D or many other stencil problems, it is bound by the DRAM bandwidth and can only poorly utilise the execution units of the processor.


Therefore, we next use the technique of kernel fusion to increase the data reuse for general sparse initial value problems with a low arithmetic intensity. Using this technique, we can fuse groups of basic operations in the data flow graph into one kernel each. Such a fused kernel reuses its loaded data via the on-chip memory of the processor, transfers its generated data via the on-chip memory of the processor and avoids unnecessary write-backs of its generated data to DRAM. In doing so, we can fuse all the basic operation of a RHS→LC-link into one kernel each. However, for general initial value problems, an RHS operation may arbitrarily access its argument vector for evaluating one of its components. Therefore, we cannot fuse LC→RHS-dependencies and basic operations from different RHS→LC-links. In addition, we apply two enabling transformations before the fusion, namely the splitting transformation to split LC operations and the cloning transformation to redundantly compute a basic operation. Next, we present how a variant can be automatically generated with kernel fusion on both GPUs and CPUs with architecture-specific optimisations. Finally, for several explicit ODE methods, we perform an experimental evaluation of the implementation variants generated via kernel fusion and enabling transformations with the test problem BRUSS2D on several modern GPUs and CPUs. This evaluation shows that on GPUs, depending on the method and the GPU used, the kernel fusion and enabling transformations can achieve a speedup between 1.62 and 3.32 compared to the respective basic implementation of the method, while on CPUs, depending on the method and the CPU used, a speedup between 1.98 and 4.14 can be achieved compared to the respective basic implementation of the method. This speedup is due to the fact that, depending on the method, the kernel fusion and enabling transformations can reduce the amount of data transferred between the processor and DRAM to between 50% and 25% compared to the respective basic implementation of the method. However, the experiments also show that the implementation variants generated by kernel fusion and enabling transformations are still bound by the DRAM bandwidth on the tested modern GPUs and CPUs.

Since the arbitrary access pattern of RHS operations for general initial value problems prevents further locality optimisations, we next specialise in RHS with a limited access distance. Here, an RHS evaluation of a single component may only access those components of the argument vector that are within the limited access distance around this given component. This allows us to fuse LC→RHS-dependencies via a partitioned fusion. We can now apply this partitioned fusion along the dependency chain of the ODE method, resulting in two tiling schemes, namely the trapezoidal tiling scheme and the hexagonal tiling scheme. In a kernel generated via tiling, a processor can now reuse data again analogously to a kernel generated via kernel fusion. We realize this tiling not only across the stages of an ODE method, but also across its time steps. We also introduce how multiple cores of a processor can cooperate on a tile to increase the maximum tile size and thereby the data reuse. Next, we address how a variant with tiling can be implemented on both GPUs and CPUs with architecture-specific optimisations. Additionally, we develop an autotuning approach, which determines the optimal tile sizes for both tiling schemes. Finally, for several explicit ODE methods, we conduct an experimental evaluation of the implementation variants generated via tiling with BRUSS2D as a test problem on several modern GPUs and CPUs. This evaluation shows that the variants generated via tiling can greatly increase both performance and data reuse on both the test GPUs and test CPUs compared to the variants generated via kernel fusion. For example, for the tested ODE methods, for a small limited access distance of 32 and on the tested GPUs the variants with tiling achieve a speedup of 1.8 to 8.2 compared to a reference implementation generated by kernel fusion. Also on the tested CPUs for the same small limited access distance of 32 the variants with tiling achieve a speedup of 7.5 to 14.1 compared to a reference implementation generated by kernel fusion. The experimental evaluation also shows that even for small to medium access distances, it is very advantageous for the performance of the variants with tiling if a tile is not only processed by a single core, but if several cores cooperate on one tile. By having several cores of the Volta GPU or of the Skylake CPU cooperate on the same tile, for a large limited access distance of 131 072 a speedup of 1.4 and 1.6 respectively can be achieved compared to the reference implementation generated by kernel fusion.


For this thesis, we implemented the generation and the execution of the variants of the basic implementation, the variants with fusion and the variants with tiling via the data flow graph abstraction of an ODE method by an automatic framework. This framework supports both GPUs and CPUs as target platforms and additionally implements several of helpful functionalities, such as automatic memory management, NUMA support for CPUs and an autotuning functionality for the tiling. Additionally, this framework includes an DSL, which allows a user to easily add new ODE methods to the framework.