Fourier operator

Fourier operator is linear integral operator with the exponential kernel.

Let some complex–valued function $$A$$ be defined for real values of the argument, id est, $$A : \mathbb R \mapsto \mathbb C$$. Then the Fourier operator, acting on $$A$$, gives function  $$B : \mathbb R \mapsto \mathbb C$$, or, simply, $$B=\mathrm{Fourier}~A$$, such that


 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(1) ~ ~ ~ \displaystyle

B(x)= (\mathrm{Fourier} A)(x)= \frac{1}{\sqrt{2\pi}} \int_{- \infty}^{\infty} \exp( - \mathrm{i} x y) A(y) ~ \mathrm{d} y $$

If the integral converges, then, function $$B$$ is called Fourier transform of function $$A$$.

Other notations
Slightly different notations are used in the article Fourier transform; operator $$\mathcal F$$ is defined such that


 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(2) ~ ~ ~ \displaystyle \mathcal F A (x)= \sqrt{2 \pi} (\mathrm{Fourier}~ A)(2 \pi x)$$

In such a way,


 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(3) ~ ~ ~ \displaystyle \mathcal F A (x)= \int_{- \infty}^{\infty} \exp( - 2 \pi \mathrm{i} x y) A(y) ~ \mathrm{d} y $$

Inverse operator
The inverse of the Fourier operator is its conjugate, id est,


 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(4) ~ ~ ~ \displaystyle {\mathrm {Fourier}}^{-1} = \mathrm{Fourier}^*$$

where $$\mathrm{Fourier}^*$$ acts on a function $$A$$ in such a way, that


 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(5) ~ ~ ~ \displaystyle \mathrm{Fourier}^*=\frac{1}{\sqrt{2\pi}} \int_{- \infty}^{\infty} \exp(\mathrm{i} x y) A(y) ~ \mathrm{d} y $$

Then, for the operator $$\mathcal F$$, the inverse acts on a fiction $$A$$ as follows:
 * $$ \!\!\!\!\!\!\!\!\!\!\!\!\!(6) ~ ~ ~ \displaystyle {\mathcal F}^{-1} ~A~ (x) = \sqrt{2 \pi} \mathrm{Fourier}^* ~ A~ (2 \pi x)

$$ $$= \int_{- \infty}^{\infty} \exp( 2 \pi \mathrm{i} x y) A(y) ~ \mathrm{d} y $$

Application of the Fourier operator
The Fourier transform is used in various sciences and, especially, in physics. In Quantum mechanics, function $$A$$ may have sense of the wave function in the coordinate representation, then $$B=\mathrm{Fourier}~ A$$ corresponds to the momentum representation. Many equations can be solved using the Fourier Transform. Especially efficient is the Fourier transform in the consideration of the paraxial approximation of the wave equations.

General properties of the Fourier operator
The Fourier operator is linear: for complex constants $$\alpha$$ and $$\beta$$,
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!(7) ~ ~ ~ \displaystyle \mathrm{Fourier} (\alpha A+\beta B)= \alpha \mathrm{Fourier}(A) + \beta \mathrm{Fourier} (B)$$

Norm of the Fourier operator
The norm $$||A||$$ of a function $$A$$ can be defined as
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (9) ~ ~ ~ ||A||^2 = \int_{-\infty}^{\infty} |A(x)|^2 ~\mathrm{d} x$$

The norm of any function is assumed to be non–negative; $$||A||\ge 0$$.

The Fourier operator is unitary, so, it preserves the norm of a function; if $$A=\mathrm{Fourier}(B)$$ then $$||A||=||B||$$

Iteration of the Fourier operator
The iterations of the Fourier operator do not provide a wide variety, because
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (10) ~ ~ ~ \mathrm{Fourier}^4= I$$

Any integrable continuous function $$A$$ can be used as a an origin to construct a self–Fourier function as a polynomial
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (11) ~ ~ ~ S(A)= A+ \mathcal{F} A+ \mathcal{F}^2 A+ \mathcal{F}^3 A$$

which is self-Fourier function; id set, $~S~$ is eigenfunction of the Fourier operator with eigenvalue unity.

Self–Fourier functions
In general function $A$ is called self-Fourier function, if
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (12) ~ ~ ~ \mathcal{F} A = A$$

All the self-Fourier functions are eigenfunction of the Fourier operator with eigenvalue unity. Any self–Fourier function can be represented as series,
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (13) ~ ~ ~ A(x) ~=~ \sum_{n=0}^{\infty} ~c_n~ \mathrm{HermiteH}_{4n}\!(x)~ \exp(-x^2/2)$$

where $$\mathrm{HermiteH}$$ are the Hermit polynomials , and $c$ are arbitrary complex coefficients; the sequence $c$ should decay quickly enough to provide the convergence of the series .

In particular, the first three self-Fourier functions $A_0, A_4,A_8$ can be defined with
 * $\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (13) ~ ~ ~ A_0(x)=\exp(-x^2/2)$
 * $\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (14) ~ ~ ~ A_4(x)=\exp(-x^2/2) (x^4-3x^2)$
 * $\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (15) ~ ~ ~ A_8(x)=\exp(-x^2/2) (x^8-14x^6+35x^4)$

The self-Fourier functions are good for testing of the numerical implementations of the Fourier operator $\mathcal{F}$.

Numerical implementation of the Fourier operator


In order to approximate the Fourier operator, the Discrete Fourier operator is used. In this section, the numerical implementation is discussed. The C++ implementation fafo.cin of operator $\mathrm{Fourier}$ is used to plot figures. This implementation is described below. The codes used to plot figures are supplied in the descriptions of the figures.

There is no need to describe functions $$f(x)\!=\!\sin(2x)$$, $$f(x)\!=\!\sin(kx)$$ , $$f(x)\!=\!\sin(2 \pi x/T)$$ and so on separately, it is sufficient to have function sin. In the similar way, it is sufficient to describe the discrete implementation of the Fourier operator at some grid. Then, other implementations can be built-up with minimal modifications of the algorithm, shifting and scaling the argument of the function and applying the corresponding modifications to the Fourier transform.

For the number $$N=2^m$$ of grid points, the step $$d$$ of the grid is chosen as follows:
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (16) ~ ~ ~ d = \sqrt{2 \pi /N}$$

The grid points $$x_j$$, $$j=0..N$$ are chosen for natural number $$~j~$$; $$~0\!\le ! j \!<\! N$$ in the following way:
 * \!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (17) ~ ~ ~ x_j = (j-N/2) d ~

At such a discretization, the node with number $$N/2$$ always corresponds to the argument zero. The 0th node is assumed to be at the left hand side of the sampled range. However, it can be placed also at the right hand side as well, and interpreted as $x_{N}$. The algorithm makes no difference between these points.

The same grid is used both for the function and for its Fourier–transform. Functions $$A$$ and $$B=\mathrm{Fourier}(A)$$ are represented with the arrays;
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (18) ~ ~ ~ A(x_j)=A_j~ ~; ~ ~ ~ B(x_j)=B_j~$$

The approximation of $$B$$ corresponds to the replacement of the integral to the discrete sum:
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (19) ~ ~ ~

B_j=\frac{d}{\sqrt{2 \pi}} \sum_{n=0}^{N-1} A_n \exp( -\mathrm{i} ~ x_j ~  x_n) $$

It is convenient to chose $$N=2^m$$ for some natural number $$m$$; then the fast implementation of such a summation is especially simple.

For the testing of the algorithm, the transform of discrete analogy of the Gaussian exponential $$A(x)\!=\!\exp(-x^2/2)$$  is suggested at the figure above.

Such a Gaussian is self-Fourier function, and it coincides with its Fourier-transform. However, at small number of the grid points, $$N\!=\!4$$, the deviation of the blue segmented line (that represents the array $$A$$) deviates from the red line, that represents the array $$B$$ as the numeric transform of $$A$$. At larger number $$N$$ of grid points, the strong zooming-in is necessary to see the deviation.

A little bit more complicated example with self-Fourier function
 * $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\displaystyle (20) ~ ~ ~ A(x)\!=\!\exp(-x^2/2)(−3x^2+x^4)$$

is shown in figure at right in the same notations, as in previous figure, but $$N\!=\!16$$. In this case, the segmented lines for arrays $$A$$ and $$B$$ overlap so perfectly, that the deviation is not seen. This deviation $$A\!-\!B$$ scaled with factor 100 is shown with the saw-like line along the abscissa axis. The deviation is of order of $$10^{-3}$$, and only in the 0th node of the grid it is a little bit larger, than $$10^{-3}$$.

The examples show, that the numerical approximation of the Fourier operator is much better, than the 1st order approximation of the function and its Fourier transform with the arrays $$A$$ and $$B$$, the deviation of the segmented line from the smooth Gaussian is much larger than the deviation of the initial array from its discrete transform.

Everyone may load the code and play with the discrete implementation of the Fourier operator for various input arrays, approximating various functions.

Keywords
Fourier, Integral transform, Linear algebra