JuliaMath / JuliaMath/AbstractFFTs.jl

Inplace irfft calculations with `mul!` mutate the fourier-space array

Open
#101 0 comments 1 reaction 0 assignees View on GitHub
Dominant language
Julia
Stars
134
Forks
40
PR merge metrics
No merged PRs in 30d

Description

I ran into some unexpected behavior with both `FFTW.jl` and `FastTransforms.jl` but *not* `CUFFT.jl`. I am thinking this repo is the appropriate location to raise this issue, but its possible it is an issue with both `FFTW.jl` and `FastTransforms.jl`

## Summary

When applying a irfft plan to a Fourier space matrix `uhat` with `mul!` and a preallocated output array `u`, `uhat` is unexpectedly mutated.

Given a `uhat` matrix of elements `ComplexF64` and shape `(3,4)`, these work as expected:

```julia
u = irfft(uhat, 4, 1:2)
```

```julia
plan = plan_irfft(uhat, 4, 1:2)
u = plan * uhat
```

this generates the correct `u`, but mutates `uhat`:

```julia
plan = plan_irfft(uhat, 4, 1:2)
u = zeros(4,4)
mul!(u, plan, uhat)
```

## Reproduction

```julia
using FFTW
using LinearAlgebra: mul!

# (3,4) array
uhat = ComplexF64.([
1 2 3 4;
1 2 3 4;
1 2 3 4;
])

#
# first, perform an irfft without a plan
#
u = irfft(uhat, 4, 1:2)

println("\n\n########\n########")
println("unplanned irfft:")
display(u)
println("û: unchanged:")
display(uhat)

# generate an irfft plan
plan = plan_irfft(uhat, 4, 1:2)

#
# calculate the irfft with the plan (allocates)
#
u_plan = plan * uhat

println("\n\n########\n########")
println("planned irfft with allocation:")
display(u_plan)
println("uhat: unchanged:")
display(uhat)

#
# calculate the irfft with the plan (reuse allocation)
#

u_plan_preallocation = zeros(4,4)
mul!(u_plan_preallocation, plan, uhat)

println("\n\n########\n########")
println("planned irfft with preallocation:")
display(u_plan)
println("uhat: HAS changed:")
display(uhat)
```

This outputs:

```
########
########
unplanned irfft:
4×4 Matrix{Float64}:
2.5 -0.5 -0.5 -0.5
0.0 0.25 0.0 -0.25
0.0 0.0 0.0 0.0
0.0 -0.25 0.0 0.25
û: unchanged:
3×4 Matrix{ComplexF64}:
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im

########
########
planned irfft with allocation:
4×4 Matrix{Float64}:
2.5 -0.5 -0.5 -0.5
0.0 0.25 0.0 -0.25
0.0 0.0 0.0 0.0
0.0 -0.25 0.0 0.25
uhat: unchanged:
3×4 Matrix{ComplexF64}:
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im
1.0+0.0im 2.0+0.0im 3.0+0.0im 4.0+0.0im

########
########
planned irfft with preallocation:
4×4 Matrix{Float64}:
2.5 -0.5 -0.5 -0.5
0.0 0.25 0.0 -0.25
0.0 0.0 0.0 0.0
0.0 -0.25 0.0 0.25
uhat: HAS changed:
3×4 Matrix{ComplexF64}:
10.0+0.0im -2.0-2.0im -2.0+0.0im -2.0+2.0im
10.0+0.0im -2.0-2.0im -2.0+0.0im -2.0+2.0im
10.0+0.0im -2.0-2.0im -2.0+0.0im -2.0+2.0im
```

Contributor guide

No contributing guide indexed for this repository

Assessment

This issue has not been assessed yet.

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.