用显卡加速c语言程序,教程 | 如何在Julia编程中实现GPU加速

原标题：教程 | 如何在Julia编程中实现GPU加速

选自nextjournal

作者：Simon Danisch

参与：高璇、刘晓坤

GPU 的并行线程可以大幅提升速度，但也使得代码编写变得更复杂。而 Julia 作为一种高级脚本语言，允许在其中编写内核和环境代码，并可在大多数 GPU 硬件上运行。本文旨在介绍 GPU 的工作原理，详细说明当前的 Julia GPU 环境，以及展示如何轻松运行简单 GPU 程序。

为了简化操作，可以在 nextjournal 上注册账户，点击「edit」即可直接运行文章中的简单代码了。

注册地址：https://nextjournal.com/signup

首先，什么是 GPU？

GPU 是一种大型并行处理器，有几千个并行处理单元。例如，本文使用的 Tesla k80，能提供 4992 个并行 CUDA 核。GPU 在频率、延迟和硬件性能方面与 CPU 有很大的不同，但实际上 Tesla k80 有点类似于具有 4992 核的慢速 CPU。

能够启动的并行线程可以大幅提升速度，但也令使用 GPU 变得更困难。当使用这种未加处理的能量时，会出现以下缺点：

GPU 是一种有专属内存空间和不同架构的独立硬件。因此，从 RAM 到 GPU 内存(VRAM，显存)的传输时间很长。甚至在 GPU 上启动内核(调用调度函数)也会带来很大的延迟，对于 GPU 而言是 10us 左右，而对于 CPU 只有几纳秒。

在没有高级封装的情况下，建立内核会变得复杂。

低精度是默认值，高精度的计算可以很容易地消除所有性能增益。

GPU 函数(内核)本质上是并行的，所以编写 GPU 内核不比编写并行 CPU 代码容易，而且硬件上的差异增加了一定的复杂性。

与上述情况相关的很多算法都不能很好地迁移到 GPU 上。想要了解更多的细节，请看这篇博文：https://streamhpc.com/blog/2013-06-03/the application-areas-opencl- cuda-can- used/。

内核通常是用 C/ C++语言编写的，但这并不是写算法的最好语言。

CUDA 和 OpenCL 之间有差异，OpenCL 是编写底层 GPU 代码的主要框架。虽然 CUDA 只支持英伟达硬件，OpenCL 支持所有硬件，但并不精细。要看个人需求进行选择。

而 Julia 作为一种高级脚本语言，允许在其中编写内核和环境代码，同时可在大多数 GPU 硬件上运行！

GPUArrays

大多数高度并行的算法都需要同时处理大量数据，以克服所有的多线程和延迟损耗。因此，大多数算法都需要数组来管理所有数据，这就需要一个好的 GPU 数组库作为关键的基础。

GPUArrays.jl 是 Julia 为此提供的基础。它实现了一个专门用于高度并行硬件的抽象数组。它包含了设置 GPU、启动 Julia GPU 函数、提供一些基本数组算法等所有必要功能。

抽象意味着它需要以 CuArrays 和 CLArrays 的形式实现。由于继承了 GPUArrays 的所有功能，它们提供的接口完全相同。唯一的区别出现在分配数组时，这会强制用户决定这一数组是存在于 CUDA 还是 OpenCL 设备上。关于这一点的更多信息，请参阅「内存」部分。

GPUArrays 有助于减少代码重复，因为它允许编写独立于硬件的 GPU 内核，这些内核可以通过 CuArrays 或 CLArrays 编译到本地的 GPU 代码。因此，大多通用内核可以在从 GPUArrays 继承的所有包之间共享。

选择小贴士：CuArrays 只支持 Nvidia GPU，而 CLArrays 支持大多数可用的 GPU。CuArrays 比 CLArrays 更稳定，可以在 Julia 0.7 上使用。速度上两者大同小异。我建议都试一试，看看哪种最有效。

本文中，我将选择 CuArrays，因为本文是在 Julia 0.7 / 1.0 上编写的，CLArrays 暂不支持。

性能

用一个简单的交互式代码示例来快速说明：为了计算 julia 集合(曼德勃罗集合)，我们必须要将计算转移到 GPU 上。

using CuArrays, FileIO, Colors, GPUArrays, BenchmarkTools

using CuArrays:CuArray

"""

The function calculating the Julia set

"""

function juliaset(z 0, maxiter)

c = ComplexF32(- 0. 5, 0. 75)

z = z 0

fori in1:maxiter

abs2(z) > 4f 0&& return(i - 1) % UInt8

z = z * z + c

end

returnmaxiter % UInt8 # % is used to convert without overflow check

end

range = 100:50:2^ 12

cutimes, jltimes = Float64[], Float64[]

function run_bench( in, out)

# use dot syntax to apply `juliaset` to each elemt of q_converted

# and write the output to result

out .= juliaset.( in, 16)

# all calls to the GPU are scheduled asynchronous,

# so we need to synchronize

GPUArrays.synchronize(out)

end

# store a reference to the last results for plotting

last_jl, last_cu = nothing, nothing

forN inrange

w, h = N, N

q = [ComplexF32(r, i) fori= 1:-( 2.0/w) :-1, r=- 1.5:( 3.0/h) :1.5]

for(times, Typ) in((cutimes, CuArray), (jltimes, Array))

# convert to Array or CuArray - moving