RcppArrayFire 0.0.2: Rcpp-Integration für ArrayFire

Das RcppArrayFire-Paket nutzt Rcpp um ein Interface zwischen R und der ArrayFire-Bibliothek bereitzustellen. ArrayFire ist eine quelloffene Bibliothek, die GPUs und andere Hardwarebeschleuniger über CUDA oder OpenCL nutzen kann.

Die offizielle R-Anbindung veröffentlicht die ArrayFire-Datenstrukturen als Objekte in R, was eine große Menge an Quellcode erfordern würde um alle Funktionen im C/C++ API von ArrayFire zu unterstützen. RcppArrayFire, das auf RcppFire von Kazuki Fukui aufbaut, verwendet den Ansatz von Paketen wie RcppArmadillo oder RcppEigen um eine direkte Kommunikation zwischen R und ArrayFire auf der C++-Ebene zu ermöglichen.

Installation

RcppArrayFire wird unter Linux entwickelt und wurde bisher auch nur unter Linux getestet.

Voraussetzungen

Um RcppArrayFire nutzen zu können müssen die Entwicklungswerkzeuge für R, Rcpp sowie die ArrayFire-Bibliothek mit zugehörigen Header-Dateien installiert sein. Auf einem aktuellen Debian basierten oder abgeleiteten System ist folgendes möglich:

sudo apt-get install r-base-dev r-cran-rcpp libarrayfire-dev libarrayfire-unified3 

Dadurch werden das Unified- und CPU-Backend installiert. Das CUDA-Backend wurde bisher nicht für Debian paketiert, und die Nutzung des OpenCL-Backends wird durch einen Fehler im clBLAS-Paket behindert. Daher ist es bei ernsthafter Anwendung derzeit besser entweder die Bibliothek selbst zu bauen oder den offiziellen Installer zu nutzen:

sudo apt-get install r-base-dev r-cran-rcpp libglfw3
# download installer from http://arrayfire.com/download/
sudo sh arrayfire-installer.sh --exclude-subdir --prefix=/usr/local

Im letzten Befehl muss noch der Name des Installers sowie (optional) der Prefix angepasst werden. Für die Nutzung von GPUs werden zusätzliche Treiber benötigt. Bei vielen integrierten GPUs von Intel kann Beignet zum Einsatz kommen:

sudo apt-get install beignet-opencl-icd ocl-icd-libopencl1

Die Installation von CUDA oder anderen OpenCL Treibern würde den Umfang sprengen. Die ArrayFire-Dokumentation enthält weitere Informationen.

Paket-Installation

RcppArrayFire ist nicht auf CRAN, kann aber über drat installiert werden:

#install.packages("drat")
drat::addRepo("daqana")
install.packages("RcppArrayFire")

Über das Argument --with-arrayfire kann bei der Installation gegebenenfalls mitgeteilt werden, wo die ArrayFire-Bibliothek installiert wurde.:

install.packages("RcppArrayFire", configure.args = "--with-arrayfire=/opt/arrayfire-3")

Nutzung

Berechnung von pi durch Simulation

Die Berechnung von pi ist ein einfaches Beispiel für eine numerische Simulation. Dabei werden viele Punkte zufällig innerhalb des Einheitsquadrats erzeugt. Eine Näherung für pi ergibt sich aus dem Verhältnis der Punkte innerhalb des Einheitskreises zur Gesamtzahl der Punkte. In R könnte man dieses Verfahren wie folgt umsetzen:

piR <- function(N) {
    x <- runif(N)
    y <- runif(N)
    4 * sum(sqrt(x^2 + y^2) < 1.0) / N
}

set.seed(42)
system.time(cat("pi ~= ", piR(10^7), "\n"))
## pi ~=  3.140899
##    user  system elapsed 
##   0.784   0.076   0.860

Mit Hilfe des inline-Pakets oder cppFunction() aus dem Rcpp-Paket, die beide auch von RcppArrayFire unterstützt werden, ist es einfach C++-Code mit R zu verbinden. Hier eine Umsetzung dieses Verfahrens in C++ mit ArrayFire:

src <- '
double piAF (const int N) {
    array x = randu(N, f32);
    array y = randu(N, f32);
    return 4.0 * sum<float>(sqrt(x*x + y*y) < 1.0) / N;
}'
Rcpp::cppFunction(code = src, depends = "RcppArrayFire", includes = "using namespace af;")

RcppArrayFire::arrayfire_set_seed(42)
system.time(cat("pi ~= ", piAF(10^7), "\n"))
## pi ~=  3.141066
##    user  system elapsed 
##   0.000   0.000   0.018

Auffällig ist:

(1) Die Syntax ist nahezu identisch. Außer der Notwendigkeit Datentypen anzugeben und einem abweichenden Funktionsnamen zur Erzeugung der Zufallszahlen, fallen das Argument f32 von randu und die explizite Angabe von floatauf. Diese weisen ArrayFire an mit einfacher Genauigkeit zu rechnen, nachdem nicht alle Geräte mit doppelter Genauigkeit umgehen können. Ist doppelte Genauigkeit unterstützt und soll verwendet werden, so ist f64 und double anzugeben.

(2) Die Ergebnisse sind verschieden, da ArrayFire einen anderen Zufallszahlengenerator nutzt.

(3) Die Lösung mit ArrayFire ist deutlich schneller. Es kommt aber gelegentlich vor, dass die erste Ausführung einer Funktion wegen der verwendeten just-in-time Kompilierung langsamer ist als spätere Aufrufe.

Arrays als Argumente

Bisher war die Argumente und Rückgabewerte der Funktion einfache Typen wie double oder int. Es ist jedoch auch möglich Arrays direkt zu verwenden. Das Matrixprodukt X‘ X einer zufälligen Matrix X lässt sich in R so berechnen:

set.seed(42)
N <- 40
X <- matrix(rnorm(N * N * 2), ncol = N)
tXXR <- t(X) %*% X

In ArrayFire lässt die die Matrixmultiplikation über die passende matmul Funktion umsetzen:

src <- '
af::array squareMatrix(const RcppArrayFire::typed_array<f32>& x) {
    return af::matmulTN(x ,x);
}'
Rcpp::cppFunction(code = src, depends = "RcppArrayFire")
tXXGPU <- squareMatrix(X)

all.equal(tXXR, tXXGPU)
## [1] "Mean relative difference: 1.372856e-07"

Nachdem ein Objekt vom Typ af::array unterschiedliche Datentypen enthalten kann, wird hier die Wrapperklasse RcppArrayFire::typed_array<> eingesetzt um den gewünschten Datentyp beim Schritt von R nach C++ festzulegen. Die Verwendung einfach genauer Fließkommazahlen erklärt die Abweichung zwischen den Ergebnissen. Verwenden wir das CPU-Backend, so können wir sicher mit doppelter Genauigkeit arbeiten und erhalten ein identisches Ergebnis:

src <- '
af::array squareMatrixF64(const RcppArrayFire::typed_array<f64>& x) {
    return af::matmulTN(x ,x);
}'
Rcpp::cppFunction(code = src, depends = "RcppArrayFire")

RcppArrayFire::arrayfire_set_backend("CPU")
tXXCPU <- squareMatrixF64(X)
RcppArrayFire::arrayfire_set_backend("DEFAULT")

all.equal(tXXR, tXXCPU)
## [1] TRUE

Nutzung in einem Paket

Bei ernsthafter Nutzung empfiehlt es sich die Funktionen permanent zu definieren. Um dies zu ermöglichen, enthält RcppArrayFire die Funktion RcppArraFire.package.skeleton(). Diese Funktion initialisiert ein R-Paket mit einem configure-Skript um mit ArrayFire und RcppArrayFire zu linken. Neue Funktionalität kann dann in Form von C++-Funktionen im src-Ordner hinzugefügt werden. Funktionen die aus R heraus aufgerufen werden sollen, können mit dem Attribut [[Rcpp::export]] gekennzeichnet werden. Die Rcpp-Vignetten über Attribute und Pakete enthalten weitere Details.

Ausblick

Neben Tests auf anderen Plattformen außer Linux sollten zukünftige Versionen auch zusätzliche Funktionalität bringen: