Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

In electrostatics, the laplacian of the electrostatic potential $\Delta V(\mathbf{r})$ arising from a charge distribution $\rho(\mathbf{r})$ is

$$ \Delta V(\mathbf{r})=-\frac{\rho(\mathbf{r})}{\epsilon_0},$$

where $\epsilon_0$ is the vacuum permitivity.

I would like to use Fast Fourier Transforms (FFT) to solve this PDE using the FFTW library. FFTW computes $Y$, the Discrete Fourier Transform (DFT) of a 1D complex array $X$ of size $N$ with this definition for the forward transform: \begin{equation} Y_{k}\triangleq\sum_{n=0}^{N-1}X_{n}e^{-2i\pi nk/N}, \end{equation} and the following definition for the backward (or inverse) transform: \begin{equation} Y_{k}=\sum_{n=0}^{N-1}X_{n}e^{2i\pi nk/N}. \end{equation}

As far as I understand, the sampling interval in these definitions is equal to 1. Also, $FFT^{-1}(FFT(X))=N\cdot X$, which is not usual. The manual says this definition is not unitary and that normalization should be taken into consideration by the user.

My application is in 3 dimensions, with regular sampling in each direction, but the sampling frequency is not the same in all directions ($L_1$, $L_2$, $L_3 \neq 1$).

My question : could you help me get an expression of $\Delta V$ in k-space ? With a simplest definition of FFTs, it would be for instance $-\mathbf{k}^2 V(\mathbf{k})$. Here ... I can't get it.

Sorry if the question is a bit messy. I'm not a mathematician and this is quite new for me.

share|improve this question
1  
I think $\Delta V$ should be $-k^2 V$ in $k$-space. –  Fabian Jan 11 '12 at 16:35
    
@Fabian corrected. Thanks. –  max Jan 11 '12 at 17:36

1 Answer 1

up vote 3 down vote accepted

You define $L_1, L_2, L_3$ to be the sampling interval in the 3 directions. Thus a point $\mathbf{r}$ in real space is given by $$\mathbf{r} = (L_1 n_1, L_2 n_2, L_3 n_3) $$ where $n_j \in \{0,\dots, N-1\}$.

Let us understand how $\partial V(\mathbf r)/\partial r_1 $ is represented in Fourier space. As you sample your space, the derivative has to be replaced by a finite difference. So we have $$\frac{\partial V(\mathbf r)}{\partial r_1} \approx \frac{V(\mathbf{r} + L_1 \mathbf{e}_1) - V(\mathbf{r})}{L_1}= \frac{V_{n_1 +1, n_2, n_3} - V_{n_1-1,n_2,n_3}}{2L_1}$$ where I introduced the notation $V_{n_1,n_2,n_3} = V(\mathbf{r}=(L_1 n_1, L_2 n_2, L_3 n_3))$.

Given the definition of the Fourier transform in the question, you have $$\hat{V}_{k_1,k_2,k_3} = \sum_{n_1,n_2,n_3} V_{n_1,n_2,n_3} \exp\left(-2\pi i \sum_{j=1}^3 n_j k_j/N_j\right). \qquad\qquad(1) $$ Then $$V_{n_1,n_2,n_3} = \frac{1}{N_1 N_2 N_3} \sum_{k_1,k_2,k_3} \hat V_{k_1,k_2,k_3} \exp\left(2\pi i \sum_{j=1}^3 n_j k_j/N_j\right); \qquad\qquad(2)$$ Eq. (1) is the definition of $\hat{V}$ and $\hat{V}$ will be what the library will hopefully calculate for you. Eq. (2) then just follows from (1) --- and sadly does not coincide with the backtransform which the library offers.

Plugging (2) into the expression for $\partial V(\mathbf r)/\partial r_1$ yields $$ \frac{\partial V(\mathbf r)}{\partial r_1} = \frac{1}{N_1 N_2 N_3} \sum_{k_1,k_2,k_3} \frac{i\sin(2\pi k_1)}{L_1} \hat V_{k_1,k_2,k_3} \exp\left(2\pi i \sum_{j=1}^3 n_j k_j/N_j\right).$$ Thus $\partial/\partial k_1$ is represented in the Fourier domain as $i\sin(2\pi k_1)/L_1$. Analogue, you can easily show that $$ \frac\partial{\partial k_j} \mapsto \frac{i\sin(2\pi k_j)}{L_j}.$$

Using these results, we obtain $$ \Delta \mapsto - \sum_{j=1}^3\frac{\sin(2\pi k_j)^2}{L_j^2}. $$

share|improve this answer
    
That's very nice I'll now be able to demonstrate the Froueri transform (whatever the definition) of $\Delta$. Thank you @Fabian. What is unclear to me is considering my practical problem and the actual definition of backward FFT in FFTW, what should I divide the forward FFT of $-\rho(\mathbf{r})/\epsilon_0$ by in order to calculate $V(\mathbf{r})$? –  max Jan 11 '12 at 17:48
    
@MaximilienLevesque: keep the forward transform as defined in the package, so just solve $[\sum_j \sin^2(2\pi k_j)/L_j^2] \hat{V} = \hat{\rho}/\epsilon_0$. Just make sure that you include the factor $1/N_1 N_2 N_3$ when (at the end) transforming back to $V$. –  Fabian Jan 11 '12 at 18:07

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.