Tullio.jl icon indicating copy to clipboard operation
Tullio.jl copied to clipboard

Error when specifying the range of an index with a UnitRange

Open jgidi opened this issue 1 year ago • 4 comments

Hi!

If we perform an operation with @tullio that needs to specify the range of an index manually, as follows:

using Tullio

N = 10
A = rand(5)

@tullio B[i, j] := A[i] (j in 1:N)

It errors with a message: ERROR: MethodError: no method matching similar(::Vector{Float64}, ::Type{Float64}, ::Tuple{Base.OneTo{Int64}, UnitRange{Int64}})

I have found two workarounds so far.

  1. Write the value of N explicitly (which is very limiting):
@tullio B[i, j] := A[i] (j in 1:10)
  1. Use Base.OneTo instead of a UnitRange:
@tullio B[i, j] := A[i] (j in Base.OneTo(N))

This seems to be because similar does not accept UnitRanges as axes, and probably Tullio made some conversion in the first case when N was provided explicitly, but I don't know why the same conversion is not performed when the range is given as 1:N. Also, I can confirm this happens at least for the julia versions 1.6.7 and 1.9.0.

While the second workaround circumvents the problem completely, I think this behavior is odd and worth mentioning.

Thanks in advance!

jgidi avatar May 15 '23 04:05 jgidi

Yes this really ought to work. It can see the 1 so it should be smart enough to make a OneTo for you. (In fact I thought there was already code for this, but must be wrong.)

Note that if you write 2:N it has to let this though, and make an OffsetArray:

julia> using OffsetArrays

julia> @tullio B[i, j] := A[i] (j in 2:N)
5×9 OffsetArray(::Matrix{Float64}, 1:5, 2:10) with eltype Float64 with indices 1:5×2:10:
 0.580193  0.580193  0.580193  0.580193  0.580193  0.580193  0.580193  0.580193  0.580193
 ...

mcabbott avatar May 15 '23 04:05 mcabbott

I can effectively reproduce your example. It seems that the code for that functionality only gets available after loading OffsetArrays, even for a UnitRange starting at 1.

After using OffsetArrays my first example works, but an OffsetArray is returned, which shouldn't be the case:

julia> using Tullio                                                    
                                                                       
julia> N = 10;                                                         
                                                                       
julia> A = rand(5);

julia> @tullio B[i, j] := A[i] (j in 1:N)
ERROR: MethodError: no method matching similar(::Vector{Float64}, ::Type{Float64}, ::Tuple{Base.OneTo{Int64}, UnitRange{Int64}})
...

julia> using OffsetArrays

julia> @tullio B[i, j] := A[i] (j in 1:N)
5×10 OffsetArray(::Matrix{Float64}, 1:5, 1:10) with eltype Float64 with indices 1:5×1:10:
 0.897518   0.897518   0.897518   …  0.897518   0.897518   0.897518
 ...

jgidi avatar May 15 '23 04:05 jgidi

OK there is code to make OneTo, but it's too restrictive:

https://github.com/mcabbott/Tullio.jl/blob/5d3105552e87efe6de981fb4d78e23e7df587fbb/src/macro.jl#L633-L637

julia> r = :(1:N)
:(1:N)

julia> r.args[1] == :(:) && length(r.args) == 3
true

julia> r.args[2] == 1
true

julia> r.args[3] isa Integer
false

julia> @tullio B[i, j] := A[i] (j in 1:10)  # with literal 10, this path makes OneTo(10):
5×10 Matrix{Float64}:
 0.573777  0.573777  0.573777 ...

mcabbott avatar May 15 '23 12:05 mcabbott

This may seem like an odd question, but do we really need to check that r.args[3] isa Integer? I just tested removing that check on L635 and everything seems to work great. Obviously, if we pass a non-integer N it will fail,

julia> N = 5.5;

julia> @tullio B[i, j] := A[i] (j in 1:N)
ERROR: MethodError: no method matching Base.OneTo(::Float64)
...

but I think that error is pretty obvious to the user (at least more than it is now) and I don't think there would be side effects besides yielding the same error if a literal Float is passed as 1:5.5.

If you like the idea I can propose this (very small) PR.

jgidi avatar May 15 '23 13:05 jgidi