julia - 分阶段编程——Jake Bolewski 的演讲

Question

版本:从 Julia v0.4 开始(我使用 0.5.0-dev+433 (2015-09-29 15:39 UTC))

引用:Jake Bolewski: Staged programming in Julia

问题:在观看了 Jakes Bolewski 关于 StaticVec 的演讲后，我没有理解 length 函数示例背后的想法。

julia> type StaticVec{T,N}
        vals::Vector{T}
       end

julia> StaticVec(T,vals...) = StaticVec{T,length(vals)}([vals...])
StaticVec{T,N}

julia> v= StaticVec(Float64,1,2,3)
StaticVec{Float64,3}([1.0,2.0,3.0])

非分阶段长度:

julia> function Base.length{T,N}(v::StaticVec{T,N}) 
        N
       end
length (generic function with 58 methods)

julia> code_llvm(length, (StaticVec{Float64,3},))

define i64 @julia_length_21889(%jl_value_t*) {
top:
  ret i64 3
}

和分阶段长度版本

julia> @generated function Base.length{T,N}(v::StaticVec{T,N}) 
        :(N)
       end
length (generic function with 58 methods)

julia> code_llvm(length, (StaticVec{Float64,3},))

define i64 @julia_length_21888(%jl_value_t*) {
top:
  ret i64 3
}

给出相同的llvm代码。

我想我理解分阶段编程背后的想法，但在这个特定的例子中我不明白说话者的意图。谁能给我解释一下吗？

Answer 1

该示例可能不是最佳选择，因为正如您所指出的，它根本不需要生成的函数。陈嘉豪最近写了一篇blog post有一个使用生成函数进行有效数据平滑的出色示例。下面是一段与他类似的代码，它通过剥离前面和后面的 M 循环迭代来进一步提高效率，从而避免主循环体中的分支:

immutable SavitzkyGolayFilter{M,N} end

wrapL(i, n) = ifelse(1 ≤ i, i, i + n)
wrapR(i, n) = ifelse(i ≤ n, i, i - n)

@generated function smooth!{M,N}(
    ::Type{SavitzkyGolayFilter{M,N}},
    data::AbstractVector,
    smoothed::AbstractVector,
)
    # compute filter coefficients from the Jacobian
    J = Float64[(i-M-1)^(j-1) for i = 1:2M+1, j = 1:N+1]
    e₁ = [1; zeros(N)]
    C = J' \ e₁

    # generate code to evaluate filter on data matrix
    pre  = :(for i = 1:$M end)
    main = :(for i = $(M+1):n-$M end)
    post = :(for i = n-$(M-1):n end)
    for loop in (pre, main, post)
        body = loop.args[2].args
        push!(body, :(x = $(C[M+1]) * data[i]))
        for j = reverse(1:M)
            idx = loop !== pre ? :(i-$j) : :(wrapL(i-$j,n))
            push!(body, :(x += $(C[M+1-j]) * data[$idx]))
        end
        for j = 1:M
            idx = loop !== post ? :(i+$j) : :(wrapR(i+$j,n))
            push!(body, :(x += $(C[M+1+j]) * data[$idx]))
        end
        push!(body, :(smoothed[i] = x))
    end
    quote
        n = length(data)
        n == length(smoothed) || throw(DimensionMismatch())
        @inbounds $pre; @inbounds $main; @inbounds $post
        return smoothed
    end
end

smooth{S<:SavitzkyGolayFilter,T}(::Type{S}, data::AbstractVector{T}) =
    smooth!(S, data, Vector{typeof(1.0*one(T))}(length(data)))

例如，为 smooth(SavitzkyGolayFilter{3,4}, rand(1000)) 生成的代码如下:

n = length(data)
n == length(smoothed) || throw(DimensionMismatch())
@inbounds for i = 1:3
    x =   0.5670995670995674   * data[i]
    x +=  0.02164502164502159  * data[wrapL(i - 3, n)]
    x += -0.1298701298701297   * data[wrapL(i - 2, n)]
    x +=  0.32467532467532445  * data[wrapL(i - 1, n)]
    x +=  0.32467532467532473  * data[i + 1]
    x += -0.12987012987013022  * data[i + 2]
    x +=  0.021645021645021724 * data[i + 3]
    smoothed[i] = x
end
@inbounds for i = 4:n-3
    x =   0.5670995670995674   * data[i]
    x +=  0.02164502164502159  * data[i - 3]
    x += -0.1298701298701297   * data[i - 2]
    x +=  0.32467532467532445  * data[i - 1]
    x +=  0.32467532467532473  * data[i + 1]
    x += -0.12987012987013022  * data[i + 2]
    x +=  0.021645021645021724 * data[i + 3]
    smoothed[i] = x
end
@inbounds for i = n-2:n
    x =   0.5670995670995674   * data[i]
    x +=  0.02164502164502159  * data[i - 3]
    x += -0.1298701298701297   * data[i - 2]
    x +=  0.32467532467532445  * data[i - 1]
    x +=  0.32467532467532473  * data[wrapR(i + 1, n)]
    x += -0.12987012987013022  * data[wrapR(i + 2, n)]
    x +=  0.021645021645021724 * data[wrapR(i + 3, n)]
    smoothed[i] = x
end
return smoothed

正如您想象的那样，这会生成非常高效的机器代码。我希望这能在一定程度上澄清生成函数的概念。