Skip to content

Conversation

@dlfivefifty
Copy link
Contributor

The restriction to one-based indexing seems artificial for Diagonal. This PR removes the restriction.

@jishnub
Copy link
Member

jishnub commented Dec 15, 2025

IIRC, the main issue was that matrix multiplication wasn't well defined for offset arrays. This isn't insurmountable (discussion at JuliaArrays/OffsetArrays.jl#91 (comment)), but I imagine there might be more code that needs to be updated to properly support offset Diagonal matrices beyond just constructing one.

The previous attempt at getting offset matrix multiplication working more generally is at JuliaArrays/OffsetArrays.jl#270.

@dlfivefifty
Copy link
Contributor Author

Diagonal multiplication is easier because it reduces to broadcasting...

@mcabbott
Copy link
Collaborator

mcabbott commented Dec 16, 2025

Apart from *, how confident are you that no LinearAlgebra methods accepting Diagonal silently assume one-based indexing?

function Matrix{T}(D::Diagonal) where {T}
    B = Matrix{T}(undef, size(D))
    ...
    copyto!(B, D)
...    
function triu!(D::Diagonal{T}, k::Integer=0) where T
    n = size(D,1)
    if !(-n + 1 <= k <= n + 1)
...
function rmul!(A::AbstractMatrix, D::Diagonal)
    _muldiag_size_check(A, D)
    for I in CartesianIndices(A)
        row, col = Tuple(I)
        @inbounds A[row, col] *= D.diag[col]
...

In fact, that last one is already unsafe -- now #1525

julia> using OffsetArrays, LinearAlgebra

julia> rmul!(OffsetArray([1 2; 3 4], 5, 6), Diagonal([7, 8]))
2×2 OffsetArray(::Matrix{Int64}, 6:7, 7:8) with eltype Int64 with indices 6:7×7:8:
 7021710134324764715   6945210754644233948
 2618386329264742529  -4556322564421083720

@dkarrasch
Copy link
Member

Indeed, there are many methods that rely on the one-based indexing check at construction, like (additionally to the ones already mentioned) rmul!(::Diagonal, ::Tridiagonal) for instance. I wouldn't even call it "silently", because both Diagonal and Tridiagonal have the corresponding checks at construction, so why would we even check it in other places?

@dlfivefifty dlfivefifty changed the title Allow diagonal offset arrays WIP: Allow diagonal offset arrays Dec 16, 2025
@dlfivefifty
Copy link
Contributor Author

OK I've renamed this WIP. I think if we add tests for all the functions which specialise for Diagonal it should be OK.

In theory Tridiagonal, Bidiagonal and SymTridiagonal can be generalised as well but they are less clear what the conditions are.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants