Solve Poisson Equation Using FFT I am trying to solve Poisson equation using FFT. The issue appears at wavenumber $k = 0$ when I want to get inverse Laplacian  which means division by zero.
We have 
${\nabla ^2}\phi  = f$ 
Taking FFT from both side we get:
$-k^2\hat\phi = \hat f $
or 
$\hat\phi = \frac{\hat f}{-k^2} $
Assuming that we want to solve this equation in periodic domain and using DFT using FFTW package, at $k=0$ we have a division by zero. Anybody knows how to deal with this singularity?

This is conceptually similar the integrating constants that show up when you are solving a differential equation by other methods. Usually these integrating constants are decided by your boundary conditions.
Instead of doing any division you can simply rewrite it as a least-norm problem: $$\|-k^2\hat \phi +\hat f\|_2^2+\text{any more terms you may want}$$
If your $\hat \phi,\hat f$ are stored in vectors the $k^2$ in the expression above will be  a diagonal weight matrix which multiplies the $\hat \phi$ vector.

Oops sorry I forgot where the FFTW comes into place, you can use it to calculate the matrix-vector product, $F$ below is the Fourier transform matrix.
$$\|-k^2 \hat \phi + F f\|_2^2$$
So $Ff$ means multiply $f$ vector by $F$ and it will be the output of calling your FFTW library with $f$ vector as the input.
And for other terms you have it might be $F^{-1}$ you want to multiply with but that is simply the corresponding IFFT routine in your library.
