朱莉娅功能表现混乱

问题描述 投票:2回答:2

我是朱莉娅的新手,因为我第一次尝试了一个温和的非平凡程序,下面写了下面的,它计算了从1到N范围内所有数字的除数sigma函数 - 除数sigma函数sigma(k,n)计算总和n的除数的第k次幂(如果k为0,则只是除数的数量)。由于sigma(k,n)是n的乘法函数,我们将其写为如下(原则上,这应该比分解我们范围内的每个整数更有效):

using Primes

thefunc(i) = (a::Int, b::Int)-> sum(BigInt(a)^(k*i) for k in 0:b)

function powsfunc(x, n, func)
    bound = convert(Int, floor(log(n)/log(x)))
    return Dict{Int, Number}(x^i => func(x, i) for i in 0:bound)
end

dictprod2(d1, d2, bd)=Dict(i*j=> d1[i]*d2[j] for i in keys(d1), j in keys(d2) if i*j<bd)

function dosigma(k, n)
    primefunc = thefunc(k)
    combfunc = (a, b) -> dictprod2(a, b, n)
    allprimes = [powsfunc(p, n, primefunc) for p in primes(n)]
    trivdict = Dict{Int, Number}(1=>1)
    theresult = reduce(combfunc, trivdict, allprimes)
    return theresult
end

好消息是以上的作品。坏消息是它的速度非常慢,dosigma(0, 100000)占用了10分钟的CPU时间,耗资150GB(!)。问题是:为什么?

performance julia
2个回答
1
投票

这是我的目标。第一个powsfunc

function powsfunc(x, n, func)
    bound = Int(floor(log(n)/log(x)))
    final_dict = Dict{Int, Int}()

    for i in 0:bound
        final_dict[x^i] = func(x,i)
    end

    return final_dict
end

我将字典改为{Int, Int}类型。这比{Int, Number}更精确,所以它应该稍微高效 - 正如克里斯在评论中指出的那样。但实际测试并没有表现出任何性能差异。

另一个变化是,我不是使用理解来定义Dict,而是逐个分配字典元素。这降低了内存估计 - 由@benchmark powsfunc(10, 1000, thefunc(0))测量,从1.68KiB到1.58KiB。不是100%肯定这是为什么。我想我在某处看到Dict的定义尚未完全抛光,但我找不到那个参考。也许更有知识的人可以阐明为什么会这样?

dictprod2可以获得巨大的收益。再一次,感谢Chris指出这一点。用Dict类型定义Int会产生巨大的差异。

function dictprod2(d1, d2, bd)
    res = Dict{Int, Int}()

    for (key1, val1) in d1
        for (key2, val2) in d2
            key1*key2>bd ? continue : true
            res[key1*key2] = val1 * val2
        end
    end 

    return res
end

我对此进行了基准测试:

primefunc = thefunc(0)
b1 = powsfunc(2, 1000, primefunc)
b2 = powsfunc(5, 1000, primefunc)

@benchmark dictprod2(b1, b2, 1000)

您的代码的结果是9.20KiB分配,中位运行时间为5.513μs。新代码仅使用1.97KiB,中位运行时间为1.68μs。

最后的功能,我几乎保持不变,因为我没有找到任何改进它的方法。

function dosigma(k, n)
    primefunc = thefunc(k)
    combfunc = (a, b) -> dictprod2(a, b, n)
    allprimes = [powsfunc(p, n, primefunc) for p in primes(n)]
    trivdict = Dict{Int, Int}(1=>1)
    theresult = reduce(combfunc, trivdict, allprimes)
    return theresult
end

检查@benchmark所有这些需要dosigma(0, 1000)从81ms和28MB到只有16ms和13MB。

我还运行了dosigma(0, 100000)并获得了85s和52GiB的分配。

我将把它留给专家来补充这个答案。我相信更有知识的人可以让这更快。


0
投票

我做了一些更多的思考/编码/分析(主要是在非常好的JuliaBox平台上),一行摘要是:

朱莉娅不是LISP

这意味着功能风格的编程只会让你到目前为止。

现在,让我们来看看代码(无论如何我们都想看到)。首先,这里是参考“愚蠢”的实现:

using Primes

function sigma(k, x)
   fac = factor( x)
   return prod(sum(i^(j*k) for j in 0:fac[i]) for i in keys(fac))
   end

allsigma(k, n) = [sigma(k, x) for x in 2::n]

在JuliaBox(也在我的笔记本电脑上),k = 0,n = 10000000需要大约40秒。精明的读者会注意到,由于溢出,这将打破较大的k。我解决这个问题的方法是将函数替换为:

function sigma(k, x)
   fac = factor( x)
   return prod(sum(BigInt(i)^(j*k) for j in 0:fac[i]) for i in keys(fac))
   end

对于相同的计算(在我的桌面上95秒),这需要131秒,因此稍微慢了3倍。相比之下,Mathematica中的相同计算如下:

foo = DivisorSigma[0, Range[10000000]] // Timing

需要27秒(在桌面上),这表明,Mathematica最有可能首先检查计算是否可以在fixnums中完成,然后以明显的方式进行。

现在我们继续进行“智能”实施。这里的工作假设是,Julia不是用于操作数据的函数语言级别,因此,考虑到这一点,我避免了以下代码的字典的创建和销毁:

using Primes

thefunc(i) = (a::Int, b::Int)-> sum(BigInt(a)^(k*i) for k in 0:b)

function biglist(n, func, thedict)
bot = Int(ceil(sqrt(n)))
theprimes = primes(bot, n)
    for i in theprimes
        thedict[i] = func(i, 1)
    end
    return theprimes
end

function powsfunc(x, n, func, motherdict)
    bound = convert(Int, floor(log(n)/log(x)))
    res = Int[]
    for i in 1:bound
        tmp = x^i
        push!(res, tmp)
        motherdict[tmp] = func(x, i)
    end
    return res
end

function makeprod(l1, l2, nn, dd)
    res = []
    for i in l1
        for j in l2 
            if i*j <= nn
                dd[i*j] = dd[i] * dd[j]
                push!(res, i*j)
            end
        end
    end
    return vcat(l1, l2, res)
end

function dosigma3(n, k)
    basedict = Dict{Int, BigInt}(1=>1)
    ff = thefunc(k)
    f2(a, b) = makeprod(a, b, n, basedict)
    smallprimes = reverse(primes(Int(ceil(sqrt(n))) -1))
    bl = biglist(n, ff, basedict)
    for i in smallprimes
        tmp = powsfunc(i, n, ff, basedict)
        bl = makeprod(bl, tmp, n, basedict)
    end
    return basedict
end

现在,dosigma3(100000, 0)需要0.5秒,比我的原始代码加速1500倍,比其他答案加速150倍。

dosigma3(10000000, 0也需要很长时间才能在JuliaBox上运行,但在前面提到的桌面上需要130秒,因此在愚蠢的实现中只有两倍。

对代码的检查表明,makeprod例程没有排序各种输入,这可能导致更快的终止。为了使它们更快,我们可以引入合并步骤,因此:

function domerge(l1, l2)
   newl = Int[]
   while true
        if isempty(l1)
            return vcat(l2, reverse(newl))
        elseif isempty(l2)
            return vcat(l1, reverse(newl))
        elseif l1[end]>l2[end]
                tmp = pop!(l1)
                push!(newl, tmp)
        else
                tmp = pop!(l2)
                push!(newl, tmp)
        end
    end
end

function makeprod2(l1, l2, nn, dd)
    res = Int[]
    for i in l1
         restmp = Int[]
         for j in l2 
            if i*j > nn
                break
            end
            dd[i*j] = dd[i] * dd[j]
            push!(restmp, i*j)
        end
        res = domerge(res, restmp)
    end
    return domerge(l1, domerge(l2, res))
end

function dosigma4(n, k)
    basedict = Dict{Int, BigInt}(1=>1)
    ff = thefunc(k)
    f2(a, b) = makeprod(a, b, n, basedict)
    smallprimes = reverse(primes(Int(ceil(sqrt(n))) -1))
    bl = biglist(n, ff, basedict)
    for i in smallprimes
        tmp = powsfunc(i, n, ff, basedict)
        bl = makeprod2(bl, tmp, n, basedict)
    end
    return basedict
end

然而,这已经花费了15秒已经为100000,所以我没有尝试10000000,并且表明我自己在朱莉娅的无能(毕竟我是新手)或朱莉娅在管理记忆方面的无能。我非常期待有见地的评论。

更新一个simpler sieving implementation将代码加速另一个因子5,并引入一些更多的缓存购买另一个因子10(!),因此计算第一个10000000(与上面相同的机器,没有使用的并行性)的sigma_k需要28(当第一次完成时为2秒(第二次完成时。所以,我猜智力并非完全没用。

© www.soinside.com 2019 - 2024. All rights reserved.