FFT with @view

General Usage
, ,
1

I want to do fft and want to write the result on a predefined array. However, there is a problem with this when the source and target arrays are both views:

using FFTW, LinearAlgebra
M = randn(100,100,100)
C = @view M[:,:,1]
B = plan_fft(C)
A = @view M[:,:,2]

mul!(A,B,C)

ERROR: MethodError: no method matching mul!(::SubArray{…}, ::FFTW.cFFTWPlan{…}, ::SubArray{…}, ::Bool, ::Bool)
The function `mul!` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  mul!(::AbstractVecOrMat, ::UniformScaling, ::AbstractVecOrMat, ::Number, ::Number)
   @ LinearAlgebra ~/.julia/juliaup/julia-1.11.1+0.x64.apple.darwin14/share/julia/stdlib/v1.11/LinearAlgebra/src/uniformscaling.jl:285
  mul!(::AbstractMatrix, ::AbstractVecOrMat, ::AbstractVecOrMat, ::Number, ::Number)
   @ LinearAlgebra ~/.julia/juliaup/julia-1.11.1+0.x64.apple.darwin14/share/julia/stdlib/v1.11/LinearAlgebra/src/matmul.jl:285
  mul!(::AbstractArray, ::Number, ::AbstractArray, ::Number, ::Number)
   @ LinearAlgebra ~/.julia/juliaup/julia-1.11.1+0.x64.apple.darwin14/share/julia/stdlib/v1.11/LinearAlgebra/src/generic.jl:132
  ...

Stacktrace:
 [1] mul!(C::SubArray{…}, A::FFTW.cFFTWPlan{…}, B::SubArray{…})
   @ LinearAlgebra ~/.julia/juliaup/julia-1.11.1+0.x64.apple.darwin14/share/julia/stdlib/v1.11/LinearAlgebra/src/matmul.jl:253
 [2] top-level scope
   @ REPL[6]:1
Some type information was truncated. Use `show(err)` to see complete types.

However, I can do mul!(A,B,C) if A and C are regular arrays.

Is there a way to do this using mul! ? I know there is an in place fft plan one can make, but that doesn’t suit my use case. I need this mul! approach to work otherwise I will redesign the code to avoid @view.

2

Fundamentally, you need a complex destination for your FFT, and then the docs for the AbstractFFT interface make it clear that the mul! interface requires both input and output to be complex arrays:

julia> using FFTW, LinearAlgebra

julia> M = randn(100, 100, 100);

julia> Mc = complex(M);

julia> C = @view Mc[:, :, 1];

julia> A = @view Mc[:, :, 2];

julia> B = plan_fft(C);

julia> mul!(A, B, C)
100×100 view(::Array{ComplexF64, 3}, :, :, 2) with eltype ComplexF64:
...
3 Likes
3

True, a hastily distilled example. Thank you. But in the problem code the target is complex. I will have to revisit … I must have misdiagnosed the issue

4

maybe you are using plan_rfft and you see something similar to

1 Like
5

thank you!