Massively parallel computations of the Bose-Hubbard model with time-dependent potentials

Abstract

The present thesis describes an efficient algorithm to determine the time evolution of the Bose-Hubbard model with time-dependent potentials. The Bose-Hubbard model was interpreted as a system of particles in a periodic structure at absolute zero. The system described here had M potential wells and N particles.

The main part of the algorithm presented here was the calculation of a matrix-vector product, which was carried out matrix-free, i.e. without storing the matrix explicitly. Due to the 'on-the-fly' nature of the matrix-vector product, it was straightforward to incorporate a time dependence of the potentials into the system. The order of Fock states was an important part of the derivation of this algorithm. Using combinatorics a state indexing routine was found. For neighbouring site transitions of particles between two states the access time to the amplitudes of these states is of order O(1). Therefore the time complexity of the matrix-vector product is O(MD) with the Hilbert space dimension D. It is also possible to use the neighbouring site transitions multiple times to perform an arbitrary transition, so more complicated operators can be implemented. The time evolution itself was carried out using two different Runge-Kutta methods, namely the classical 4th order Runge-Kutta method and an embedded Runge-Kutta method, which can provide 4th and 5th order time steps simultaneously.

A performance improvement was achieved by executing the presented algorithm in parallel utilizing a graphics card as a general purpose computational device. The GPU was chosen as a platform for these calculations, as it is capable to execute a large amount of tasks simultaneously. The implementation was done using CUDA, an API provided by Nvidia. The performance of this implementation was compared to cuSPARSE, a library for parallel sparse matrix computations also developed by Nvidia. The presented algorithm was found to outperform the well optimized cuSPARSE implementation by a factor of 165% on a consumer level graphics card.

To prove the relevancy of the GPU implementation, a CPU one of the very same algorithm was made. The execution time of the CPU implementation surpassed the GPU one by a factor of 200% with consumer level hardware and 2300% on a professional graphics card.

The presented algorithm was also compared to a different method of matrix-free matrix-vector multiplication involving hash functions and a Bisection search to retrieve the index of a given Fock state. Due to this required search algorithm, the time complexity of this method is O(MDlog D), which has a worse scaling behaviour than the presented algorithm.

Last but not least, the functionality of the presented method was demonstrated by application to different systems and comparison to another work, where the Bogoliubov-Backreaction method was implemented.

Die vorliegende Arbeit beschreibt einen effizienten Algorithmus zur Berechnung der Zeitentwicklung des Bose-Hubbard Modells mit zeitabhängigen Potentialen. Das Bose-Hubbard Modell wurde interpretiert als ein System von Teilchen in einer periodischen Gitterstruktur bestehend aus M Potentialtöpfen und N Teilchen am absoluten Temperaturnullpunkt.

Hauptbestandteil des vorgestellten Algorithmus ist eine Matrix-Vektor-Multiplikation, welche "Matrix-frei" ausgeführt wurde, d.h. ohne die Matrix explizit zu speichern. Wegen der dynamischen Natur dieses Matrix-Vektor-Produkts war es unkompliziert eine Zeitabhängigkeit der Potentiale in das System einzubauen. Die Reihenfolge der Fock-Zustände war ein wichtiger Teil der Herleitung des Algorithmus. Mithilfe von Kombinatorik wurde eine Zustandsnummerierungsmethode gefunden. Für Nächste-Nachbar-Übergänge von Teilchen zwischen zwei Zuständen ist die Zugriffszeit auf die Amplituden dieser Zustände von der Ordnung O(1). Deshalb ist die Zeitkomplexität des Matrix-Vektor-Produkts von der Ordnung O(MD), wobei D die Dimension des Hilbertraums ist. Es ist außerdem möglich durch hintereinander Ausführung von mehreren Nächste-Nachbar-Sprüngen beliebige Übergänge zu berechnen, sodass auch kompliziertere Operatoren damit umgesetzt werden können. Die Zeitentwicklung selbst wurde mithilfe von zwei verschiedenen Runge-Kutta-Verfahren durchgeführt, der klassischen Runge-Kutta-Methode 4ter Ordnung und einer eingebetteten Runge-Kutta-Methode. Die eingebettete Runge-Kutta-Methode bestimmt Zeitschritte 4ter und 5ter Ordnung gleichzeitig.

Durch parallele Ausführung des vorgestellten Algorithmus auf einer Grafikkarte wurde eine Verbesserung der Laufzeit erzielt. Die GPU wurde als Plattform für diese Berechnungen gewählt, da sie in der Lage ist eine große Anzahl an Rechnungen gleichzeitig auszuführen. Die Implementation wurde mittels CUDA realisiert, einer von Nvidia bereitgestellten API. Die Performanz dieser Implementation wurde mit cuSPARSE verglichen, einer Bibliothek für dünn besetzte Matrizen, welche ebenfalls von Nividia entwickelt wurde. Der vorgestellte Algorithmus übertraf die Leistung der wohl optimierten cuSPARSE Implementation um einem Faktor von 165% auf einer kommerziell erhältlichen Grafikkarte.

Um die Relevanz der GPU Implementation zu zeigen, wurde eine CPU Umsetzung desselben Algorithmus erstellt. Die Laufzeit der CPU Implementation übertraf die der GPU um einen Faktor von 200% auf einer gewöhnlichen und um 2300% auf einer professionellen Grafikkarte.

Der vorgestellte Algorithmus wurde außerdem mit einer anderen Matrix-freien Methode zur Matrix-Vektor-Multiplikation verglichen, welche zur Berechnung des Index eines Fock-Zustands eine Hashfunktion und eine Bisektion erfordert. Wegen des benötigten Suchalgorithmuses ist die Zeitkomplexität von der Ordnung O(MDln D). Somit hat diese Methode ein schlechteres Skalierungsverhalten als der vorgestellte Algorithmus.

Zum Schluss wurde die Funktionalität der beschriebenen Methode demonstriert, indem sie auf unterschiedliche Systeme angewandt und außerdem mit Ergebnissen einer anderen Arbeit verglichen wurde, in der die Bogoliubov-Backreaction-Methode verwendet wurde.