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

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

Closed
yewalenikhil65 opened this issue Jul 28, 2023 · 8 comments

Comments

@yewalenikhil65
Copy link

yewalenikhil65 commented Jul 28, 2023

I am trying following , and this works okay!

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

using Symbolics, ToeplitzMatrices

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

@shashi
Copy link
Member

shashi commented Aug 2, 2023

It seems the multiplication is switching to FFT algorithm for large convolutions. I'm not sure what we can do about that. What's your goal with this? We could potentially make fft a primitive in Symbolics. cc @ChrisRackauckas

@xtalax
Copy link
Contributor

xtalax commented Aug 2, 2023

I fully support this idea, there's a lot of simplification you can do in Fourier space

@ChrisRackauckas
Copy link
Member

A lot of pseudospectral PDE codes could be nicely represented with it. I think more generally the question is how to register array functions. If registration works well on that, we can start getting coverage.

@xtalax
Copy link
Contributor

xtalax commented Aug 2, 2023

Definitely agree on this point, I've run up against this a few times. Could just get the registration methods to wrap on Union{Num, Arr}

@yewalenikhil65
Copy link
Author

yewalenikhil65 commented Aug 2, 2023

It seems the multiplication is switching to FFT algorithm for large convolutions. I'm not sure what we can do about that. What's your goal with this? We could potentially make fft a primitive in Symbolics. cc @ChrisRackauckas

@shashi My aim was to create that F_expr vector in symbolic form as written in the code for arbitrary large N . It creates the set of quadratic equations I was interested in solving (they form system of equations related to deep gravity waves 🌊). Only thing is such matrix-vector multiplication wasn't possible due to error mentioned above for N >= 280 . I found a way around by helpful comment here by avoiding certain conditional statements

But it would be great if making fft a primitive in Symbolics could solve this issue. Certainly the said matrix-vector product is possible and user should be allowed to create them!

What looks overall from error trace is that plan_fft! method entertaining Vector{Num} should be included somewhere.

@shashi
Copy link
Member

shashi commented Aug 2, 2023

Array registration is a top priority for me right now.

@yewalenikhil65 even if we support fft node right now, you won't be able to access individual expressions after the convolution. It's best to request Toeplitz multiplication to add an option to multiplication where it does not take the path of calling fft.

@shashi shashi closed this as completed Aug 2, 2023
@ChrisRackauckas
Copy link
Member

Array registration is a top priority for me right now.

Yeah I think that's all you need. Once that exists, the community can register a whole standard library etc. and take it from there. But without that registration, these kinds of issues will just pop up. You don't need to worry about an fft node, just a registration process and I or anyone else will get to the fft node.

@xtalax
Copy link
Contributor

xtalax commented Aug 2, 2023

This probably needs a broad request to package maintainers to register important functions in package

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

No branches or pull requests

4 participants