# Does there exists 1D FFT/IFFT implementation in pure Julia (i.e. no FFTW)

**URL:** <https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407>\
**Category:** Numerics\
**Created:** [August 16, 2017, 10:31am UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407 "2017-08-16T10:31:02Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![vgdev](https://avatars.discourse-cdn.com/v4/letter/v/47e85d/32.png) [@vgdev](https://discourse.julialang.org/u/vgdev)\
**Post date:** [August 16, 2017, 10:31am UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/1 "2017-08-16T10:31:02Z")

</div>

Hi,  
I need a pure julia FFT implementation as I am planning to use it in correspondence with the ForwardDiff package. Does this exist?

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [August 16, 2017, 10:37am UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/2 "2017-08-16T10:37:19Z")

</div>

You could adapt the `BigFloat` fft implementation in FastTransforms.jl: [https://github.com/MikaelSlevinsky/FastTransforms.jl/blob/master/src/fftBigFloat.jl](https://github.com/MikaelSlevinsky/FastTransforms.jl/blob/master/src/fftBigFloat.jl)

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [August 16, 2017, 11:59am UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/3 "2017-08-16T11:59:53Z")

</div>

I implemented an experimental one in: [https://github.com/JuliaLang/julia/pull/6193](https://github.com/JuliaLang/julia/pull/6193)

---

<div class="post-metadata">

**Author:** ![vgdev](https://avatars.discourse-cdn.com/v4/letter/v/47e85d/32.png) [@vgdev](https://discourse.julialang.org/u/vgdev)\
**Post date:** [August 16, 2017, 12:18pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/4 "2017-08-16T12:18:41Z")

</div>

Hmmm, it seems I will be running into problems with forwarddiffs dual numbers and the resulting complex numbers from the FFT.

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [August 16, 2017, 2:03pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/5 "2017-08-16T14:03:02Z")

</div>

Complex numbers are indeed not supported by ForwardDiff. If they were, a better way to differentiate through FFT would be to teach ForwardDiff the derivatives of FFTs: [https://github.com/JuliaDiff/ForwardDiff.jl/pull/165](https://github.com/JuliaDiff/ForwardDiff.jl/pull/165) and refs

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [August 16, 2017, 2:24pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/6 "2017-08-16T14:24:19Z")

</div>

Yes, since it is a linear transformation, the derivative of an FFT with respect to its inputs is simply another FFT.

---

<div class="post-metadata">

**Author:** ![RoyiAvital](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/royiavital/32/571_2.png) [@RoyiAvital](https://discourse.julialang.org/u/RoyiAvital)\
**Post date:** [August 16, 2017, 3:44pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/7 "2017-08-16T15:44:40Z")

</div>

For a vector x ou can always describe its Fourier Transform as y = D x where D is the [DFT Matrix](https://en.wikipedia.org/wiki/DFT_matrix).

You can even generalize it for cases which are “Zero Padded”.  
Since the Derivative is just {D}^{T} all you need is a function which generates those matrices.

Is there such built in function for that in Julia?

---

<div class="post-metadata">

**Author:** ![vgdev](https://avatars.discourse-cdn.com/v4/letter/v/47e85d/32.png) [@vgdev](https://discourse.julialang.org/u/vgdev)\
**Post date:** [August 16, 2017, 4:37pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/8 "2017-08-16T16:37:04Z")

</div>

Okey, so what I am really after is a way to find the derivative of the Volterra process \int\_0^tg\_\alpha(t-s) dW\_s wrt. \alpha \in (-0.5,0) where g\_\alpha(t-s) = (t-s)^\alpha and W a brownian motion. I am approximating this integral by looking at at the Riemann sum approximation \tilde{Y}(\frac{i}{n})= \sum\_{k =1}^i g(\frac{b\_k}{n})W\_{i-k} on an n sized grid, i = 0,1,2,\dots. This can be interpreted as a discrete convolution \tilde{Y}(\frac{i}{n}) = (g\*W)\_i, which can be computed efficiently by FFT (Ihave gotten this to work with ForwardDiff using a direct naive convolution, but its too slow). Furthermore, if I need to do several runs of this numerical scheme to find an estimate of \mathbb{E}[\tilde{Y}(\frac{i}{n})] \sim \frac{1}{N}\sum\_j^N\tilde{Y}^{(j)}(\frac{i}{n}). I have to perform this convolution each time I simulate a path of this process. So to speed up things I would like to compute the FFT of g(b\_k/n) once (since its deterministic) and multiply it with the FFT of a new set of random variates W, and then take the IFFT to arrive at the convolution.

The code to my current implementation (with some slight changes to the actual approximation scheme):

```julia
function gen_volterra(a,s,n,G,rand_1,rand_2)
  Z_1 = zeros(eltype(a),s); #first normal
  Z_2 = zeros(eltype(a),s); # 2nd normal
  #G = zeros(eltype(a),s+1); #Kernel function
  @fastmath @inbounds @simd for i = 1:s
    #G[i+1] = b(a,i+1)/(n^a); #fill kernel vector <<----this is where the error occurs
    Z_1[i] = get_sigma_Z1(n)*rand_1[i] #set correct variance of random variate
    Z_2[i] = get_sigma_Z2(a,n)*(get_rho(a)*rand_1[i]+sqrt(1-get_rho(a)*get_rho(a))*rand_2[i]); #set correct variance of random variate
  end
  return [0;(direct_conv(Z_1,G)[1:s]+Z_2)]; #Convolute
  #return G
end

function direct_conv{T}(a::Array{T}, b::Array{T})
    m = length(a)
    n = length(b)
    c = zeros(T,m+n-1)
    @fastmath @inbounds @simd for j=1:m
        @fastmath @inbounds @simd for k=1:n
            c[j+k-1] += a[j]*b[k]
        end
    end
    return c
end

```

and addtional helper functions:

```julia
function b(a, k)
  return (k^(a+1.0)-(k-1.0)^(a+1.0))/(a+1.0);
end

function get_sigma_Z1(n)
  return 1.0/sqrt(n);
end

function get_sigma_Z2(a,n)
  return 1.0/sqrt(((2.0*a+1.0)*n^(2.0*a+1.0)))
end

function get_rho(a)
  return sqrt(2.0*a+1.0)/(a+1.0)
end

function gen_G(a,n,s)
  i = range(1,s)
  return [0;map(i->b(a,i+1)/(n^a),i)]
end

```

Example of generating 1 path:

```julia
julia> a= -0.43
julia> n= 500
julia> s= 500
julia> G = gen_G(a,n,s)
julia> rand_1 = randn(s)
julia> rand_2 = randn(s)
julia> path = gen_volterra(a,s,n,G,rand_1,rand_2)

```

But the operation below will then be quite slow due to the number of floating point operations in the direct\_conv method.

```julia
m = ForwardDiff.derivative(a ->gen_volterra(a,s,n,G,rand_1,rand_2),-0.43)

```

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [August 16, 2017, 7:04pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/9 "2017-08-16T19:04:07Z")

</div>

> [@RoyiAvital](#):
>
> Is there such built in function for that in Julia?

You can just FFT the columns of the identity matrix. But usually it is preferable to work with the DFT matrix implicitly, since multiplying by it (or its transpose — it is symmetric) is computed efficiently by the `fft` function (and the inverse or conjugate by `ifft`).

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [August 17, 2017, 8:51am UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/10 "2017-08-17T08:51:31Z")

</div>

Why do you need AD for this? Just do the derivative wrt alpha, which is again a convolution that you can do with an FFT.

---

<div class="post-metadata">

**Author:** ![vgdev](https://avatars.discourse-cdn.com/v4/letter/v/47e85d/32.png) [@vgdev](https://discourse.julialang.org/u/vgdev)\
**Post date:** [August 21, 2017, 8:02pm UTC](https://discourse.julialang.org/t/does-there-exists-1d-fft-ifft-implementation-in-pure-julia-i-e-no-fftw/5407/11 "2017-08-21T20:02:32Z")

</div>

> [@antoine-levitt](#):
>
> Why do you need AD for this? Just do the derivative wrt alpha, which is again a convolution that you can do with an FFT.

Because this process is used as an input in a SDE for another process, which I didn’t include here for clairity. It is actually wrt to this other process I want the gradient, so the stuff outlined above is a subproblem as this is where the performance bottleneck is. But a solution might be to use [this](https://github.com/JuliaDiff/ForwardDiff.jl/pull/165)