Matrix components

[arlk]

Thanks for the awesome package! It's been very helpful in putting together some complicated models. One of my simulations has a matrix of ODEs and I tried using ComponentArrays to handle it, but it seems to flatten it into a vector. MWE:

julia> c = ComponentArray(a = 1, b = rand(2,2))
ComponentArray{Float64}(a = 1.0, b = [0.24410936846253284, 0.7216454349529038, 0.7719341762466598, 0.5748953546713345])

julia> c[:b]
4-element view(::Array{Float64,1}, 2:5) with eltype Float64:
 0.24410936846253284
 0.7216454349529038
 0.7719341762466598
 0.5748953546713345

Instead I was expecting something like this:

julia> c[:b]
4-element view(::Array{Float64,2}, 2:5, 2:5) with eltype Float64:
 0.244109  0.771934
 0.721645  0.574895

[jonniedie]

Glad to hear it's been helpful! Yeah, I ran into this last week and forgot to make an issue for it, so thanks for making one. It isn't too hard to implement for the one-dimensional ComponentArray case, just hold onto the size information so we can reshape it after. The problem is what the behavior should be for higher dimensional ComponentArrays. So given the same c, we could extrude out a c2:

julia> c2 = c .* c'
5×5 ComponentArray{Tuple{Axis{(a = 1, b = 2:5)},Axis{(a = 1, b = 2:5)}},Float64,2,Array{Float64,2}}:
 1.0       0.244109   0.721645  0.771934  0.574895
 0.244109  0.0595894  0.17616   0.188436  0.140337
 0.721645  0.17616    0.520772  0.557063  0.414871
 0.771934  0.188436   0.557063  0.595882  0.443781
 0.574895  0.140337   0.414871  0.443781  0.330505

I struggle with what the shape of c2[:b, ;b] should be. For this case, since we have a square matrix, it's easy to say it should probably look something like this:

julia> reshape(c2[:b,:b], 2, 2, 2, 2)
2×2×2×2 reshape(view(::Array{Float64,2}, 2:5, 2:5), 2, 2, 2, 2) with eltype Float64:
[:, :, 1, 1] =
 0.0595894  0.188436
 0.17616    0.140337

[:, :, 2, 1] =
 0.17616   0.557063
 0.520772  0.414871

[:, :, 1, 2] =
 0.188436  0.595882
 0.557063  0.443781

[:, :, 2, 2] =
 0.140337  0.443781
 0.414871  0.330505

But for mixed cases like this one, I suppose we could just reshape them similarly:

julia> x = ComponentArray(a=rand(), b=rand(2), c=rand(3,3));

julia> x2 = x .* x';

julia> reshape(x2[:b,:c], 2, 3, 3)
2×3×3 reshape(view(::Array{Float64,2}, 2:3, 4:12), 2, 3, 3) with eltype Float64:
[:, :, 1] =
 0.211156  0.776949  0.0218616
 0.152617  0.561557  0.015801

[:, :, 2] =
 0.425994  0.854114  0.161747
 0.307897  0.61733   0.116906

[:, :, 3] =
 0.106867   0.759792  0.751145
 0.0772407  0.549157  0.542907

So maybe that's what I'll do. Just hold on to the input dimensions and reshape things when they come out. Because I agree, having matrix components would be useful. In fact, it's probably more useful than having broadcasted extrusion in the first place (since with DifferentialEquations.jl you can explicitly declare your jacobian prototype).

[jonniedie]

Closed by 8af621e. (Then reverted, then unreverted).