The fft-based representation of a wave height field expresses the wave height $h(x, t)$ at the horizontal position $x = (x, z)$ as the sum of sinusoids with complex, time-dependent amplitudes:


where $t$ is the time and $k$ is a two-dimensional vector with components $k = (k_x, k_z)$, $k_x = 2πn/L_x$, $k_z = 2πm/Lz$, and $n$ and $m$ are integers with bounds $−N/2 ≤ n < N/2$ and $−M/2 ≤ m < M/2$. The fft process generates the height field at discrete points $x = (nLx/N,mLz/M)$. The value at other points can also be obtained by switching to a discrete fourier transform, but under many circumstances this is unnecessary and is not applied here. The height amplitude Fourier components, $\tilde{h}(k, t)$, determine the structure of the surface.

Another version of the same formula:


Here $X$ is a horizontal position of a point whose height we are evaluating. The wave vector $K$ is a vector pointing in the direction of travel of the given wave, with a magnitude k dependent on the length of the wave (l):

$$k = 2p /\lambda$$

And the value of $h (K,t)$ is a complex number representing both amplitude and phase of wave $K$ at time $t$. Because we are using discrete Fourier transformation, there are only a finite number of waves and positions that enters our equations. If $s$ is dimension of the heightfield we animate, and $r$ is the resolution of the grid, then we can write:

$$K = (kx,kz) = (2pn / s,2pm /s)$$

where $n$ and $m$ are integers in bounds $–r/2 < n,m < r/2$. Note that for FFT, $r$ must be power of two.


  • How many $k$ do I need? Does it mean, that $k$ takes all possible values (on the grid I mean)?

  • And what do I get as a result - complex number that shows height and (what?) in concrete point so I need make computations for all points on the grid or the whole grid itself?

  • How to apply fft to this formula?

  • 2
    $\begingroup$ I hope you mean that you're "trying hard", not "hardly trying". ;) $\endgroup$ Sep 4 '10 at 18:35
  • 9
    $\begingroup$ It appears that, rather than posting your question at MO, what you really need is a consultant: someone who is an expert in applied mathematics and, for appropriate remuneration, will analyze your model and provide you advice tailored to your specific problem. $\endgroup$ Sep 4 '10 at 21:32
  • 3
    $\begingroup$ It seems to me that this is not a research question for mathematicians; this is an implementation question for those who routinely do computational fluid dynamics as part of their job (but NOT applied mathematicians, who study the theoretical basis for methods) - I would suggest some kind of engineer or experimental physicist could help you, rather than a mathematician. All these numerical calculations/methods will be well-documented and implemented in standard hydrodynamics computer packages, and very familiar to those sorts of people (but not me). $\endgroup$
    – Zen Harper
    Sep 10 '10 at 2:58

It isn't totally clear to me what you are trying to do, mostly because I am not sure what problem you are trying to solve. I suspect you are looking for solutions to some wave / Laplacian equation, but without more specific details I can't comment very precisely. Nonetheless, I think I can perhaps help explain some of your questions. But first, let me recommend the following text book on spectral methods for engineers:


Now I will try to answer your questions in order:

1 It really depends on what you want to do, but if you have a real function on a grid of n x m numbers, then you need only n/2 x m/2 complex coefficients to represent by conjugate symmetry. This follows from the observation that for real f,

$ \hat{f}(\omega) = \overline{ \hat{f}(-\omega) }$

I think MATLAB has a special function called rfft that will compute this for you automatically. For a more specific answer you will probably need to analyze whatever PDE you are trying to solve and then figure out the number of coefficients you will need to get some desired accuracy.

2 In the Fourier domain, the magnitude of the kth coefficient determines the amplitude of the wave, and the arg value determines the phase. You can compute these using the square roots and arctangents respectively.

3 How you apply the FFT once again depends on what you are trying to do. Assuming you are trying to solve a PDE, what you basically will do is try to convert your PDE into a Fourier multiplier, then use it to solve a system of ODEs in order to get a solution to your final problem. Constructing the multiplier is usually straightforward and makes use of the nice fact that:

$ \mathcal{F}( \partial f )(\omega) = - i \omega \hat{f}(\omega) $

Though in the discrete case you do need to be a bit careful how you do interpolation. A common strategy is to use Shannon/Whittaker interpolation (aka sinc functions). You should read Boyd's book since he talks about them at some length. Once this is done, solving the system typically reduces to inverting a banded/diagonal matrix which is easy to do numerically.

Hope this helps, but I can't really be much more specific until you give me an actual problem to work with.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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