JuliaLinearAlgebra / ToeplitzMatrices.jl

Fast matrix multiplication and division for Toeplitz matrices in Julia
MIT License
67 stars 24 forks source link

Symbolic (symmetric-toepltiz)matrix-vector multiplication doesn't work for higher dimensions of matrix #109

Open yewalenikhil65 opened 1 year ago

yewalenikhil65 commented 1 year ago

I am trying following , and this works okay! I am raising this issue both on Symbolics and ToeplitzMatrices

using Symbolics, ToeplitzMatrices

N =250

a = Symbolics.variables(:a, 0:N);

b= [1; (1:N).*a[2:end]]

A  = SymmetricToeplitz(a)
F_expr = A*b;

However, the code fails if i increase the code to N=280, even with FFTW


using Symbolics, ToeplitzMatrices,FFTW

N =280

a = Symbolics.variables(:a, 0:N);

b= [1; (1:N).*a[2:end]]

A  = SymmetricToeplitz(a)
F_expr = A*b;

ERROR: MethodError: no method matching plan_fft!(::Vector{Complex{Num}}, ::UnitRange{Int64})

Closest candidates are:
  plan_fft!(::StridedArray{T, N}, ::Any; flags, timelimit, num_threads) where {T<:Union{ComplexF64, ComplexF32}, N}
   @ FFTW ~/.julia/packages/FFTW/HfEjB/src/fft.jl:786
  plan_fft!(::AbstractArray; kws...)
   @ AbstractFFTs ~/.julia/packages/AbstractFFTs/hJ0Fz/src/definitions.jl:65

Stacktrace:
 [1] plan_fft!(x::Vector{Complex{Num}}; kws::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ AbstractFFTs ~/.julia/packages/AbstractFFTs/hJ0Fz/src/definitions.jl:65
 [2] plan_fft!(x::Vector{Complex{Num}})
   @ AbstractFFTs ~/.julia/packages/AbstractFFTs/hJ0Fz/src/definitions.jl:65
 [3] factorize(A::SymmetricToeplitz{Num, Vector{Num}})
   @ ToeplitzMatrices ~/.julia/packages/ToeplitzMatrices/5gW1c/src/linearalgebra.jl:150
 [4] mul!(y::Vector{Num}, A::SymmetricToeplitz{Num, Vector{Num}}, x::Vector{Num}, α::Bool, β::Bool)
   @ ToeplitzMatrices ~/.julia/packages/ToeplitzMatrices/5gW1c/src/linearalgebra.jl:37
 [5] mul!
   @ ~/software/julia-1.9.2/share/julia/stdlib/v1.9/LinearAlgebra/src/matmul.jl:276 [inlined]
 [6] *(A::SymmetricToeplitz{Num, Vector{Num}}, x::Vector{Num})
   @ LinearAlgebra ~/software/julia-1.9.2/share/julia/stdlib/v1.9/LinearAlgebra/src/matmul.jl:56
 [7] top-level scope
   @ REPL[36]:1
nsajko commented 1 year ago

This seems to be the culprit:

https://github.com/JuliaLinearAlgebra/ToeplitzMatrices.jl/blob/52775d6390525bce0ca09fba9fa10092adadca96/src/linearalgebra.jl#L19-L38

Modifying the code so the small case branch is taken unconditionally fixes the issue.

yewalenikhil65 commented 1 year ago

@nsajko shouldn't this issue be fixed by maintainers of ToeplitzMatrices.jl ? Or is it something that is expected to be done in isolation by users in their own compiled package of ToeplitzMatrices.jl ?

nsajko commented 1 year ago

shouldn't this issue be fixed by maintainers of ToeplitzMatrices.jl

I guess so. I'm not a maintainer, though, so don't ask me.

yewalenikhil65 commented 1 year ago

shouldn't this issue be fixed by maintainers of ToeplitzMatrices.jl

I guess so. I'm not a maintainer, though, so don't ask me.

Thanks for the help with that smaller case condition. It helped a lot