I have the following convolution as part of a numerical simulation.
$$T(r)=\int \mathrm{d}^3r_2\, p(r_2)f(r_2)\alpha(r-r_2)\, .$$
My problem is that the analytical expressions for $f$ and $p$ do exist but, I have the expression for $\alpha$ only in the Fourier domain in the form of $\alpha(k)$. I planned to evaluate in the following way:
- Construct a grid using the mesh grid of $100\times100\times100$ using meshgrid and linspace in numpy
ran = linspace(-1,1,N_r)
x,y,z = meshgrid(ran,ran,ran) #position space
- Construct the components xf, yf, zf in the Fourier domain from x, y, z
xf = fftn(x)
yf = fftn(y)
zf = fftn(z)
- Find the Fourier transform of $f(r)\times p(r)$ using fftn in numpy
- Multiply it with $\alpha(k)$
- Taking the inverse Fourier transform using ifftn in numpy.
I am not very sure that the above method works and I actually failed to verify it properly. I tried using scipy.ndimage.convolve to compare the results with the inverse Fourier transform of the product in the Fourier domain. Is it correct as to what I am doing with the code? And is there a way where I can verify that a method is working with a simpler example?
Trying to verify:
I have tried the following to test the theory. Seems like it does not work. I expect that the result res_1 and res_2 to be the same. I also used the function real to truncate the tiny imaginary part which results from the fftn and ifftn functions.
x = linspace(-1,1,10)
xf = fftn(x)
def f(x):
return x**2+x**3*sin(x)
def g(k):
return k**2+k**3/(3-k**2)
g_k = g(xf)
g_x = real(ifftn(g_k))
res_1 = img_con(g_x,f(x))
res_2 = real(ifftn(g(xf)*fftn(f(x))))
print(res_1)
print(res_2)
Am I doing something which is wrong?