Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Structured matrices don't extend missing dimensions correctly in broadcasting #54087

Open
jishnub opened this issue Apr 15, 2024 · 1 comment · May be fixed by #54190
Open

Structured matrices don't extend missing dimensions correctly in broadcasting #54087

jishnub opened this issue Apr 15, 2024 · 1 comment · May be fixed by #54190
Labels
domain:broadcast Applying a function over a collection domain:linear algebra Linear algebra kind:bug Indicates an unexpected problem or unintended behavior

Comments

@jishnub
Copy link
Contributor

jishnub commented Apr 15, 2024

julia> D1 = Diagonal([1]); D2 = Diagonal([1,2]);

julia> D1 .+ D2
2×2 Diagonal{Int64, Vector{Int64}}:
 2  
   3

julia> Matrix(D1) .+ D2
2×2 Matrix{Int64}:
 2  1
 1  3

Ideally, these should produce identical results.
In the structured case, the 1x1 matrix D1 isn't expanded to ones(size(D2)), but instead it's expanded to Diagonal(ones(Int, size(D2,1))). I'm unsure if there's a way to resolve this while preserving type-stability.

@jishnub jishnub added the domain:linear algebra Linear algebra label Apr 15, 2024
@N5N3 N5N3 added the kind:bug Indicates an unexpected problem or unintended behavior label Apr 16, 2024
@mbauman
Copy link
Sponsor Member

mbauman commented Apr 16, 2024

Oh, yikes. The problem is in similar, but I'm not sure what to do about it:

function Base.similar(bc::Broadcasted{StructuredMatrixStyle{T}}, ::Type{ElType}) where {T,ElType}
inds = axes(bc)
fzerobc = fzeropreserving(bc)
if isstructurepreserving(bc) || (fzerobc && !(T <: Union{SymTridiagonal,UnitLowerTriangular,UnitUpperTriangular}))
return structured_broadcast_alloc(bc, T, ElType, length(inds[1]))
elseif fzerobc && T <: UnitLowerTriangular
return similar(convert(Broadcasted{StructuredMatrixStyle{LowerTriangular}}, bc), ElType)
elseif fzerobc && T <: UnitUpperTriangular
return similar(convert(Broadcasted{StructuredMatrixStyle{UpperTriangular}}, bc), ElType)
end
return similar(convert(Broadcasted{DefaultArrayStyle{ndims(bc)}}, bc), ElType)
end

julia> bc = Broadcast.instantiate(Broadcast.broadcasted(+, D1, D2))
Base.Broadcast.Broadcasted{LinearAlgebra.StructuredMatrixStyle{Diagonal}}(+, ([1;;], [1 0; 0 2]))

julia> similar(bc, Int)
2×2 Diagonal{Int64, Vector{Int64}}:
 4333067184           ⋅
          ⋅  4333069200

I think you're right that this will need to break type stability. It could alternatively be an error, but that also seems quite broken. This probably happens with all the structured matrices.

@dkarrasch dkarrasch added the domain:broadcast Applying a function over a collection label Apr 17, 2024
N5N3 added a commit to N5N3/julia that referenced this issue Apr 22, 2024
1. size-1 StructuredMatrix should behave like scalar during broadcast. Thus their `fzero` should return the only element.
(fix JuliaLang#54087)

2. But for simple broadcast with only one StructuredMatrix, we can keep stability as the structure is "preserved" even for size-1 case. Thus `count_structedmatrix` is added to check that.

3. `count_structedmatrix` is fused to keep `Bidiagonal` stability.
(replace JuliaLang#54067)
@N5N3 N5N3 linked a pull request Apr 22, 2024 that will close this issue
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
domain:broadcast Applying a function over a collection domain:linear algebra Linear algebra kind:bug Indicates an unexpected problem or unintended behavior
Projects
None yet
Development

Successfully merging a pull request may close this issue.

4 participants