Introduction
ardea provides access to GPU devices and frameworks
natively in R. It currently provides rudimentary access to
OpenCL, CUDA, and Metal; letting
users compile, test, and implement their own device specific code.
ardea currently only requires existing OpenCL to
install appropriately, but will selectively build out access to
frameworks and devices that it can detect during package
installation.
This package currently does not provide an interface for vendor supplied libraries. Functionality is currently focused on ingesting a user supplied framework specific function, compiling that function and preparing it for dispatch, executing function dispatch and buffer management, and returning that function’s result to R.
## [1] TRUE
## [1] FALSE
## [1] TRUEUsing matrix multiplication as an example, a user can bring a naive
mm kernel function conforming to the a given framework’s rules and
requirements. In this case we can just write out our functions to temp
files, and walk them forward through the frameworks and devices
available on the current system. One of ardea’s current
quirks is that it expects the first argument of a function to be
the output. There is no keyword or tooling to enforce this, but
the underlying buffer management currently only returns that first input
from the device to R.
# write out a character vector to a tempfile for opencl, metal, and cuda
opencl_mm_naive <- '
__kernel void opencl_mm_naive(__global float* output,
const long M,
const long K,
const long N,
__global const float* A,
__global const float* B)
{
long row = get_global_id(0);
long col = get_global_id(1);
if (row < M && col < N) {
float sum = 0.0f;
for (long k = 0; k < K; k++) {
sum += A[row + k * M] * B[k + col * K];
}
output[row + col * M] = sum;
}
}
'
tmp01 <- tempfile(fileext = ".cl")
writeLines(text = opencl_mm_naive,
con = tmp01)
cuda_mm_naive <- '
extern "C" __global__ void cuda_mm_naive(float* output,
const long long M,
const long long K,
const long long N,
const float* A,
const float* B)
{
long long row = blockIdx.x * blockDim.x + threadIdx.x;
long long col = blockIdx.y * blockDim.y + threadIdx.y;
if (row < M && col < N) {
float sum = 0.0f;
for (long long k = 0; k < K; k++) {
sum += A[row + k * M] * B[k + col * K];
}
output[row + col * M] = sum;
}
}
'
tmp02 <- tempfile(fileext = ".cu")
tmp03 <- tempfile(fileext = ".ptx")
writeLines(text = cuda_mm_naive,
con = tmp02)
metal_mm_naive <- '
kernel void metal_mm_naive(device float* output [[buffer(0)]],
constant uint& M [[buffer(1)]],
constant uint& K [[buffer(2)]],
constant uint& N [[buffer(3)]],
device const float* A [[buffer(4)]],
device const float* B [[buffer(5)]],
uint2 id [[thread_position_in_grid]])
{
uint row = id.x;
uint col = id.y;
if (row < M && col < N) {
float sum = 0.0f;
for (uint k = 0; k < K; k++) {
sum += A[row + k * M] * B[k + col * K];
}
output[row + col * M] = sum;
}
}
'
tmp04 <- tempfile(fileext = ".metal")
tmp05 <- tempfile(fileext = ".metallib")
writeLines(text = metal_mm_naive,
con = tmp04)Matrix multiplication doesn’t need a lot of setup when using the
builtin %*% function, which is calling on your R
installation’s BLAS library and is highly optimized. However, because
we’re reaching out to a device almost entirely ab initio with
no helpers, we need to perform the work that those helpers would
perform.
# generic inputs, regardless of framework
var00 <- vector(mode = "numeric",
length = 250000)
var01 <- matrix(rnorm(250000),
nrow = 500,
ncol = 500)
var02 <- matrix(rnorm(250000),
nrow = 500,
ncol = 500)
dim01 <- nrow(var01)
dim02 <- ncol(var01)
# dim03 <- nrow(var02)
dim04 <- ncol(var02)
arg_list <- list(var00, # our expected output vector, we've just filled it with zeros
dim01,
dim02,
dim04,
as.vector(var01),
as.vector(var02))Every framework has its own unique verbiage and intricacies, though they’re all kind of doing similar things. OpenCL is designed to be be deployable across a variety of device vendors and device types, not just GPUs. Metal is Apple’s proprietary GPU / accelerator framework, though it shares some conceptual choices with OpenCL as until recently Apple’s hardware stack was somewhat diverse. CUDA is NVIDIA’s well known framework for accessing their devices, and includes an intimidating array of libraries and functionalities.
This means that frameworks share concepts, but not necessarily the same keywords, type names, specific paths, or user interactions.
opencl_arg_types <- c("float",
"long",
"long",
"long",
"float",
"float")
metal_arg_types = c("float",
"uint",
"uint",
"uint",
"float",
"float")
cuda_arg_types <- c("float",
"long",
"long",
"long",
"float",
"float")Though each framework handles program / module / library preparation
somewhat uniquely, an effort has been made to unify how users
interact with these steps. In this package OpenCL supports API-side
compilation and program generation, while Metal supports both API-side,
and system-side compilation and program management, and CUDA
only supports system-side compilation. Meaning that under the
hood Metal and CUDA are calling on system/vendor-specific compilers to
build their libraries/modules (i.e. clang and
nvcc respectively).
if (opencl_is_available() &
opencl_devices_exist()) {
print("OpenCL is available on this system!")
cl_dvcs <- opencl_device_information()
cl_ctx <- opencl_make_context(device = cl_dvcs[[1]])
cl_program <- opencl_make_program(cl_file = tmp01,
context = cl_ctx)
cl_knl <- opencl_make_kernelptr(program = cl_program,
kernel_names = "opencl_mm_naive")
}
if (cuda_is_available() &
cuda_devices_exist()) {
print("CUDA is available on this system!")
cu_dvcs <- cuda_device_information()
cu_ctx <- cuda_make_context(device = cu_dvcs[[1]])
cu_program <- cuda_make_program(cuda_file = tmp02,
context = cu_ctx,
ptx_file = tmp03)
cu_knl <- cuda_make_kernelptr(program = cu_program,
context = cu_ctx,
kernel_names = "cuda_mm_naive")
}
if (metal_is_available() &
metal_devices_exist()) {
print("Metal is available on this system!")
mtl_dvcs <- metal_device_information()
mtl_ctx <- metal_make_context(device = mtl_dvcs[[1]])
mtl_program <- metal_make_program(metal_file = tmp04,
metallib_file = tmp05,
context = mtl_ctx)
mtl_knl <- metal_make_kernelptr(program = mtl_program,
context = mtl_ctx,
kernel_names = "metal_mm_naive")
}## [1] "OpenCL is available on this system!"
## [1] "Metal is available on this system!"A simple wrapper function is supplied with ardea, though
it is mostly a thin wrapper for .Call() runner functions.
It does perform some error checking, but for the most part users are
being trusted with their judgements.
# our builtin optimized CPU implementation
system.time(res01 <- var01 %*% var02)
if (opencl_is_available() &
opencl_devices_exist()) {
print("OpenCL implementation:")
system.time(res02 <- simple_wrapper(framework = "opencl",
context_ptr = cl_ctx,
kernel_ptr = cl_knl$opencl_mm_naive,
arg_types = opencl_arg_types,
arg_list = arg_list,
problem_dims = as.integer(c(dim01,
dim04,
1L)),
group_dims = NULL,
workers_per = NULL))
plot(as.vector(res01),
res02,
pch = 46)
}
if (cuda_is_available() &
cuda_devices_exist()) {
print("CUDA implementation:")
system.time(res03 <- simple_wrapper(framework = "cuda",
context_ptr = cu_ctx,
kernel_ptr = cu_knl$cuda_mm_naive,
arg_types = cuda_arg_types,
arg_list = arg_list,
problem_dims = as.integer(c(dim01,
dim04,
1L)),
group_dims = NULL,
workers_per = NULL))
plot(as.vector(res01),
res03,
pch = 46)
}
if (metal_is_available() &
metal_devices_exist()) {
print("Metal implementation:")
system.time(res03 <- simple_wrapper(framework = "metal",
context_ptr = mtl_ctx,
kernel_ptr = mtl_knl$metal_mm_naive,
arg_types = metal_arg_types,
arg_list = arg_list,
problem_dims = as.integer(c(dim01,
dim04,
1L)),
group_dims = NULL,
workers_per = NULL))
plot(as.vector(res01),
res03,
pch = 46)
}
## user system elapsed
## 0.030 0.000 0.031
## [1] "OpenCL implementation:"
## [1] "Metal implementation:"This is currently the limit of ardea’s functionality.
This package began mostly as a curiousity project with Metal, and
turning it into a CRAN acceptable package became a bit of a labor of
love. If folks find this package useful I will consider expanding it to
include other vendors, or more complicated dispatch routines. As it
stands now, my primary concerns are conforming to CRAN submission
guidelines, ensuring user interactions are clear and predictable, and
working through a non-trivial TODO list.
## R version 4.5.1 (2025-06-13)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sequoia 15.7.9
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
##
## locale:
## [1] C/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## time zone: America/New_York
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ardea_0.0.9
##
## loaded via a namespace (and not attached):
## [1] digest_0.6.39 R6_2.6.1 fastmap_1.2.0 xfun_0.56
## [5] cachem_1.1.0 knitr_1.50 htmltools_0.5.9 rmarkdown_2.30
## [9] lifecycle_1.0.5 cli_3.6.6 sass_0.4.10 jquerylib_0.1.4
## [13] compiler_4.5.1 tools_4.5.1 evaluate_1.0.5 bslib_0.9.0
## [17] yaml_2.3.11 rlang_1.2.0 jsonlite_2.0.0