Allow offset indices in some linear algebra #43552
Closed
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
This wants to make
A * B
work with OffsetArrays with minimal changes. The easy part is always to callsimilar(..., axes(A), ...)
and neversize
.The messier part is that it adapts the fallback
generic_matmatmul!
to allow offsets. The plan would be never to hit that fallback, by adding a methodgeneric_matmatmul!(C::OffsetArray, ...)
which removes offsets, and then callsmul!
. That will allow an offset adjoint strided array (etc) will dispatch to the right place, after the offset is removed, without duplication or ambiguities. So perhaps the fallback need not exist?Third, if we want to un-wrap the OffsetArrays and call
mul!
again, then we needα,β
notMulAdd(α,β)
. I thought that the extra stage (and instability?) of making that twice was causing an extra allocation, so I added a step to callgeneric_matmatmul!
without it, first. This does not seem to remove the allocaiton but might be the right thing to do anyway.Finally, it might also be worth making
generic_matmatmul!
accept functions likeadjoint
instead of characters'C'
. In many cases this information is in the type domain when it is called, butα,β
are not. (Last commit which does this unfortunately adds an ambiguity, right now.)