Broadcasting functions with side effects

[willtebbutt]

using Distributions, FillArrays
rand.(Normal.(Zeros(3), Ones(3)))

yields something along the lines of

3-element Fill{Float64,1,Tuple{Int64}}:
 -1.0204222556659934
 -1.0204222556659934
 -1.0204222556659934

Unless I'm missing something obvious, it looks to me like the broadcasting implementation is making overly strong assumptions regarding the purity of whatever function is being broadcasted.

[dlfivefifty]

This has been discussed in JuliaLang issues (I don't remember which ones). I believe the consensus is that the behaviour of broadcasting is "defined" by the BroadcastStyle: there is not necessarily a guarantee that it is equivalent to [f(a[k]) for k in eachindex(a)].

I propose that the BroadcastStyle returns a FillStyle, and the definition of FillStyle is that it does the above behaviour.

[cstjean]

So the recommended solution is to use map(rand, Normal.(Zeros(3), Ones(3))?

[willtebbutt]

I propose that the BroadcastStyle returns a FillStyle, and the definition of FillStyle is that it does the above behaviour.

Okay. I take your point, but I really can't see how anyone would consider this behaviour intuitive. It's just completely different to what anyone expects to happen when you broadcast some function over an AbstractArray.

When I use a FillArray I expect it to have exactly the same behaviour as any other AbstractArray in all cases, just some operations are performed more efficiently, and their results may be represented more efficiently. The same goes for all other array types that exploit known-structure in some way (Diagonal, Toeplitz, etc) - their semantics are the same as a dense array, stuff just happens more quickly where possible.

For example were I a new user, not well-versed in the implementation-details of broadcast, what possible reason could I have for expecting anything other than this behaviour?

[willtebbutt]

At the risk of labouring the point,

using LinearAlgebra, Distributions, FillArrays
rand.(Normal.(Diagonal(Ones(10))))

produces a dense matrix as expected, and it would surprise me greatly if anyone would advocate for instead producing some structured output with one value on the diagonal and another everywhere else.

[willtebbutt]

Sorry, didn't mean to close the issue there.

[dlfivefifty]

Your Diagonal example is misleading as the resulting matrix is not diagonal, but perhaps you are right.

Is it possible to tell if a function is pure at compile time? Otherwise we can just have a list of pure functions like +, exp, etc. that we override broadcasted for

[willtebbutt]

Is it possible to tell if a function is pure at compile time?

I have no idea I'm afraid :(

Otherwise we can just have a list of pure functions like +, exp, etc. that we override broadcasted for

Sounds reasonable to me.

[dlfivefifty]

The current behaviour is consistent with sparse arrays in StdLib and so will remain unless StdLib changes:

julia> A = sparse(I, 3,3);  (_ -> randn()).(A)
3×3 SparseMatrixCSC{Float64,Int64} with 9 stored entries:
  [1, 1]  =  0.679325
  [2, 1]  =  1.54721
  [3, 1]  =  1.54721
  [1, 2]  =  1.54721
  [2, 2]  =  -0.963752
  [3, 2]  =  1.54721
  [1, 3]  =  1.54721
  [2, 3]  =  1.54721
  [3, 3]  =  -0.52672