【问题标题】:Julia: Products of sequences in a vectorized wayJulia:向量化的序列乘积
【发布时间】:2017-04-06 00:40:42
【问题描述】:

学习从 Python 到 Julia,我正在尝试转换我拥有的旧代码,即计算此表达式序列的乘积:

我有两个版本的 Python 代码,一个使用 for 循环实现,另一个使用广播。 for 循环版本是:

import numpy as np
A = np.arange(1.,5.,1)
G = np.array([[1.,2.],[3.,4.]])

def calcF(G,A):
    N = A.size
    print A
    print N
    F = []
    for l in range(N):
        F.append(G/A[l])
        print F[l]
        for j in range(N):
            if j != l:
                F[l]*=((G - A[l])/(G + A[j]))*((A[l] - A[j])/(A[l] + A[j]))
    return F

F= calcF(G,A)
print F

而我从对我的问题here的回复中学到的矢量化版本是这个函数:

def calcF_vectorized(G,A):
    # Get size of A
    N = A.size

    # Perform "(G - A[l])/(G + A[j]))" in a vectorized manner
    p1 = (G - A[:,None,None,None])/(G + A[:,None,None])

    # Perform "((A[l] - A[j])/(A[l] + A[j]))" in a vectorized manner
    p2 = ((A[:,None] - A)/(A[:,None] + A))

    # Elementwise multiplications between the previously calculated parts
    p3 = p1*p2[...,None,None]

    # Set the escaped portion "j != l" output as "G/A[l]"
    p3[np.eye(N,dtype=bool)] = G/A[:,None,None]

    Fout = p3.prod(1)

    # If you need separate arrays just like in the question, split it
    return np.array_split(Fout,N)

我试图天真地将 Python for 循环代码转换为 Julia:

function JuliacalcF(G,A)
    F = Array{Float64}[]
    for l in eachindex(A)
        push!(F,G/A[l])
        println(A[i])
        for j in eachindex(A)
            if j!=l
                F[l]*=((G - A[l])/(G + A[j]))*((A[l] - A[j])/(A[l] + A[j]))
            end
        end
    end
    #println(alpha)
    return F
end
A = collect(1.0:1.0:5.0)
G = Vector{Float64}[[1.,2.],[3.,4.]]
println(JuliacalcF(G,A))

但是有没有办法像numpy 广播矢量化版本那样以智能的方式做到这一点?

【问题讨论】:

  • 您确定需要矢量化版本吗?如果您对性能(以及可能矢量化)感兴趣,我建议您首先profile 您的代码并查看performance tips

标签: python algorithm numpy julia array-broadcasting


【解决方案1】:

另外,请查看More-DotsLoop Fusion,其中通过示例描述了矢量化。

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 2021-08-05
    • 2016-01-09
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2013-07-19
    相关资源
    最近更新 更多