80

将性能与 R 进行比较的 Julia 示例似乎特别令人费解https://github.com/JuliaLang/julia/blob/master/test/perf/perf.R

您可以从以下两种算法中获得的最快性能是多少(最好说明您所做的更改以使其更像 R)?

## mandel

mandel = function(z) {
    c = z
    maxiter = 80
    for (n in 1:maxiter) {
        if (Mod(z) > 2) return(n-1)
        z = z^2+c
    }
    return(maxiter)
}

mandelperf = function() {
    re = seq(-2,0.5,.1)
    im = seq(-1,1,.1)
    M = matrix(0.0,nrow=length(re),ncol=length(im))
    count = 1
    for (r in re) {
        for (i in im) {
            M[count] = mandel(complex(real=r,imag=i))
            count = count + 1
        }
    }
    return(M)
}

assert(sum(mandelperf()) == 14791)

## quicksort ##

qsort_kernel = function(a, lo, hi) {
    i = lo
    j = hi
    while (i < hi) {
        pivot = a[floor((lo+hi)/2)]
        while (i <= j) {
            while (a[i] < pivot) i = i + 1
            while (a[j] > pivot) j = j - 1
            if (i <= j) {
                t = a[i]
                a[i] = a[j]
                a[j] = t
            }
            i = i + 1;
            j = j - 1;
        }
        if (lo < j) qsort_kernel(a, lo, j)
        lo = i
        j = hi
    }
    return(a)
}

qsort = function(a) {
  return(qsort_kernel(a, 1, length(a)))
}

sortperf = function(n) {
    v = runif(n)
    return(qsort(v))
}

sortperf(5000)
4

2 回答 2

107

这个问题的关键词是“算法”:

您可以从以下两种算法中获得的最快性能是多少(最好说明您所做的更改以使其更像 R)?

如“你能在 R 中以多快的速度制作这些算法?” 这里讨论的算法是标准的 Mandelbrot 复杂循环迭代算法和标准的递归快速排序内核。

当然有更快的方法来计算这些基准中提出的问题的答案——但不是使用相同的算法。您可以避免递归,避免迭代,并避免其他 R 不擅长的事情。但是,您不再比较相同的算法。

如果您真的想在 R 中计算 Mandelbrot 集或对数字进行排序,是的,这不是您编写代码的方式。您要么尽可能地对其进行矢量化——从而将所有工作推入预定义的 C 内核中——要么只是编写一个自定义 C 扩展并在那里进行计算。无论哪种方式,结论是 R 的速度不够快,无法单独获得真正好的性能——您需要让 C 完成大部分工作才能获得良好的性能。

这正是这些基准测试的重点:在 Julia 中,您永远不必依赖 C 代码来获得良好的性能。你可以只用纯 Julia 编写你想做的事情,它会有很好的性能。如果迭代标量循环算法是做你想做的最自然的方法,那就去做吧。如果递归是解决问题最自然的方法,那也没关系。无论是通过非自然向量化还是编写自定义 C 扩展,您都不会被迫依赖 C 来获得性能。当然,您可以在自然的情况下编写矢量化代码,因为它通常是在线性代数中;如果您已经有一些库可以满足您的需求,则可以调用 C。但你不必这样做。

我们确实希望对跨语言的相同算法进行最公平的比较:

  1. 如果有人在 R 中有使用相同算法的更快版本,请提交补丁!
  2. 我相信julia 网站上的 R 基准已经是字节编译的,但是如果我做错了并且比较对 R 不公平,请告诉我,我会修复它并更新基准。
于 2012-05-23T01:00:27.583 回答
47

嗯,在 Mandelbrot 示例中,矩阵 M 的维度已转置

M = matrix(0.0,nrow=length(im), ncol=length(re))

因为它是通过count在内部循环中递增(的连续值im)来填充的。我的实现创建了一个复数向量mandelperf.1并对所有元素进行操作,使用索引和子集来跟踪向量的哪些元素尚未满足条件Mod(z) <= 2

mandel.1 = function(z, maxiter=80L) {
    c <- z
    result <- integer(length(z))
    i <- seq_along(z)
    n <- 0L
    while (n < maxiter && length(z)) {
        j <- Mod(z) <= 2
        if (!all(j)) {
            result[i[!j]] <- n
            i <- i[j]
            z <- z[j]
            c <- c[j]
        }
        z <- z^2 + c
        n <- n + 1L
    }
    result[i] <- maxiter
    result
}

mandelperf.1 = function() {
    re = seq(-2,0.5,.1)
    im = seq(-1,1,.1)
    mandel.1(complex(real=rep(re, each=length(im)),
                     imaginary=im))
}

对于 13 倍的加速(结果相等但不相同,因为原始返回数字而不是整数值)。

> library(rbenchmark)
> benchmark(mandelperf(), mandelperf.1(),
+           columns=c("test", "elapsed", "relative"),
+           order="relative")
            test elapsed relative
2 mandelperf.1()   0.412  1.00000
1   mandelperf()   5.705 13.84709

> all.equal(sum(mandelperf()), sum(mandelperf.1()))
[1] TRUE

快速排序示例实际上并没有排序

> set.seed(123L); qsort(sample(5))
[1] 2 4 1 3 5

但我的主要加速是围绕枢轴矢量化分区

qsort_kernel.1 = function(a) {
    if (length(a) < 2L)
        return(a)
    pivot <- a[floor(length(a) / 2)]
    c(qsort_kernel.1(a[a < pivot]), a[a == pivot], qsort_kernel.1(a[a > pivot]))
}

qsort.1 = function(a) {
    qsort_kernel.1(a)
}

sortperf.1 = function(n) {
    v = runif(n)
    return(qsort.1(v))
}

加速 7 倍(与未校正的原件相比)

> benchmark(sortperf(5000), sortperf.1(5000),
+           columns=c("test", "elapsed", "relative"),
+           order="relative")
              test elapsed relative
2 sortperf.1(5000)    6.60 1.000000
1   sortperf(5000)   47.73 7.231818

由于在最初的比较中,对于 mandel,Julia 比 R 快 30 倍,对于快速排序快 500 倍,因此上述实现仍然没有真正的竞争力。

于 2012-04-02T00:17:13.613 回答