（（Julia）。+运算符似乎未使用我的自定义广播功能，为什么？

Question

我正在实现一个只有一个非零值的自定义矩阵，无论您执行什么操作，这都是矩阵中唯一一个可能为非零的单元格。我称它为SVMatrix（单值矩阵）。我到目前为止的代码是

struct SVMatrix{T} <: Base.AbstractMatrix{T}
    value::T
    index::Tuple{Int,Int}
    size::Tuple{Int,Int}
end

function Base.broadcast(+, A::SVMatrix, B::AbstractArray)
    SVMatrix(A.value+B[A.index...], A.index, A.size) 
end

function Base.getindex(A::SVMatrix{T}, i::Int) where {T}
    if i == A.index[1] + A.index[2]*A.size[1]
        A.value
    else
        0
    end
end


function Base.getindex(A::SVMatrix{T}, i::Vararg{Int,2}) where {T}
    if i == A.index
        return A.value
    else
        0
    end
end

function Base.size(A::SVMatrix)
    A.size
end

然后我以下列方式与。+运算符一起计时广播功能

function time(n::Int)
    A = SVMatrix(1.0, (3,4), (n, n))
    B = rand(n,n)
    @time broadcast(+, A, B)
    @time A .+ B
end

time(1000)
println()
time(1000)

并获得结果

 0.000000 seconds
 0.008207 seconds (2 allocations: 7.629 MiB, 47.51% gc time)

 0.000000 seconds
 0.008258 seconds (2 allocations: 7.629 MiB)

因此，似乎。+并未使用我的自定义广播功能，即使它说in the documentation表示

实际上，f。（args ...）等效于broadcast（f，args ...），提供了方便的语法来广播任何函数（点语法）。

为什么会得到这些结果？

Answer 1

这实际上是您不应扩展广播的一个很好的例子。

julia> struct SVMatrix{T} <: Base.AbstractMatrix{T}
           value::T
           index::Tuple{Int,Int}
           size::Tuple{Int,Int}
       end

julia> @inline function Base.getindex(A::SVMatrix{T}, i::Vararg{Int,2}) where {T}
           @boundscheck checkbounds(A, i...)
           if i == A.index
               return A.value
           else
               return zero(T)
           end
       end

julia> Base.size(A::SVMatrix) = A.size

julia> SVMatrix(1.0, (1,1), (2, 2)) .+ ones(2, 2)
2×2 Array{Float64,2}:
 2.0  1.0
 1.0  1.0

.+的结果应不是为[2 0; 0 0]！如果我们使用广播的实现（已更正以::typeof(+)的形式在DNF noted上调度），则当其他人使用它并期望它的行为与所有其他AbstractArray相同时，您的数组将被意外破坏。

现在，您可以返回巧妙重新计算的SVMatrix的操作是.*：

julia> SVMatrix(2.5, (1,1), (2, 2)) .* ones(2, 2)
2×2 Array{Float64,2}:
 2.5  0.0
 0.0  0.0

我们可以在O（1）的空间和时间中执行此操作，但是默认实现是循环所有值并返回密集的Array。这是朱莉娅的多重派遣大放异彩的地方：

julia> Base.broadcasted(::typeof(*), A::SVMatrix, B::AbstractArray) = SVMatrix(A.value*B[A.index...], A.index, A.size)

julia> SVMatrix(2.5, (1,1), (2, 2)) .* ones(2, 2)
2×2 SVMatrix{Float64}:
 2.5  0.0
 0.0  0.0

因为这是O（1）运算并且是巨大胜利，所以我们可以选择退出广播融合并立即重新计算新的SVMatrix-即使在“融合”表达式中也是如此。不过，您还没有完成！

• 需要对兼容形状实施错误检查。
• 需要允许广播SVMatrix(2.5, (1,1), (2, 2)) .* rand(2)之类的东西。
• [理想情况下，您将实现一个BroadcastStyle以允许在“参数列表中至少有一个BroadcastStyle”上进行分派。然后，您将实现SVMatrix，这将允许Base.broadcasted(::ArrayStyle{SVMatrix}, ::typeof(*), args...)出现在SVMatrix的任一侧，但这是高级主题。