The math is sound. The method is not.
We recast the azimuthal mode-coupled surface equations as a block-Toeplitz convolution in the mode index and diagonalize this coupling via a discrete Fourier transform in \(n\), then briefly review the basic definition and core properties of the DFT used in this step.
Diagonalize the \(m\)-coupling (DFT over \(m\))
By “diagonalize the \(m\)-coupling” we mean: take the infinite tridiagonal system that couples the modal amplitudes \(\{\hat{\Phi}_{n-1},\hat{\Phi}_n,\hat{\Phi}_{n+1}\}\) (and similarly for \(\hat{\eta}_n\)) and perform a change of variables so that each new mode evolves independently. The tool that achieves this decoupling is a discrete Fourier transform in the azimuthal index \(n\) (playing the role of \(m\)).
We start from the coupled equations in \(n\). For each fixed \(k>0\), the surface equations are \[ \begin{aligned} & -i \omega \hat{\Phi}_n+\frac{U k}{2}\left(\hat{\Phi}_{n-1}-\hat{\Phi}_{n+1}\right) +\left(g+\frac{T}{\rho} k^2\right) \hat{\eta}_n+\hat{p}_n=0, \\ & -i \omega \hat{\eta}_n+\frac{U k}{2}\left(\hat{\eta}_{n+1}-\hat{\eta}_{n-1}\right) -k \hat{\Phi}_n=0. \end{aligned}\tag{1}\label{ref1} \]
This is an infinite block-tridiagonal (Toeplitz) system along the index \(n\): each equation for \(n\) only involves \(n-1\), \(n\), and \(n+1\).
We package the unknowns into a vector sequence and view the system as a discrete convolution in \(n\): \[ Q_n=\begin{bmatrix} \hat{\Phi}_n \\ \hat{\eta}_n \end{bmatrix}, \qquad F_n=\begin{bmatrix} -\hat{p}_n \\ 0 \end{bmatrix}. \tag{2} \]
Then \(\eqref{ref1}\) can be written as \[ A Q_n+B Q_{n-1}+C Q_{n+1}=F_n, \tag{3}\label{ref3} \] with constant \(2\times 2\) blocks \[ \begin{aligned} A &= \begin{bmatrix} -i \omega & g+\dfrac{T}{\rho} k^2 \\ -k & -i \omega \end{bmatrix},\\[0.5em] B &= \begin{bmatrix} \dfrac{U k}{2} & 0 \\ 0 & -\dfrac{U k}{2} \end{bmatrix}, \qquad C = \begin{bmatrix} -\dfrac{U k}{2} & 0 \\ 0 & \dfrac{U k}{2} \end{bmatrix}. \end{aligned}\tag{4} \]
Because the coefficients are independent of \(n\), \(\eqref{ref3}\) is a (matrix-valued) convolution in the discrete index \(n\). Discrete Fourier transforms diagonalize such convolutions.
Discrete Fourier Transform in \(n\) (Diagonalization Step)
We define, for each \(\beta\in[0,2\pi]\), \[ \begin{aligned} \widetilde{Q}(k,\beta) &= \sum_{n\in\mathbb{Z}} Q_n(k)\, e^{i n \beta},\\ \widetilde{F}(k,\beta) &= \sum_{n\in\mathbb{Z}} F_n(k)\, e^{i n \beta}. \end{aligned}\tag{5} \]
Applying \(\sum_n (\cdot)\,e^{i n \beta}\) to \(\eqref{ref3}\) and using the standard shift identities \[ \sum_n Q_{n-1} e^{i n \beta} = e^{-i\beta}\widetilde{Q}(k,\beta),\qquad \sum_n Q_{n+1} e^{i n \beta} = e^{i\beta}\widetilde{Q}(k,\beta), \] we obtain \[ \bigl(A + B e^{-i\beta} + C e^{i\beta}\bigr)\, \widetilde{Q}(k,\beta) = \widetilde{F}(k,\beta). \tag{6}\label{ref6} \]
Thus the original infinite coupled system in \(n\) is reduced to independent \(2\times 2\) systems, one for each \(\beta\).
We now compute the \(\beta\)-dependent symbol: \[ \begin{aligned} B e^{-i \beta}+C e^{i \beta} &= \begin{bmatrix} \dfrac{U k}{2}\bigl(e^{-i \beta}-e^{i \beta}\bigr) & 0 \\ 0 & -\dfrac{U k}{2}\bigl(e^{-i \beta}-e^{i \beta}\bigr) \end{bmatrix} \\ &= \begin{bmatrix} -i U k \sin \beta & 0 \\ 0 & i U k \sin \beta \end{bmatrix}. \end{aligned}\tag{7} \]
Therefore \(\eqref{ref6}\) reads \[ \begin{bmatrix} -i \omega - i U k \sin \beta & g+\dfrac{T}{\rho} k^2 \\ -k & -i \omega + i U k \sin \beta \end{bmatrix} \begin{bmatrix} \widetilde{\Phi}(k,\beta) \\ \widetilde{\eta}(k,\beta) \end{bmatrix} = \begin{bmatrix} \widetilde{p}(k,\beta) \\ 0 \end{bmatrix}, \] where \(\widetilde{\Phi},\widetilde{\eta},\widetilde{p}\) are the corresponding components of \(\widetilde{Q}\) and \(\widetilde{F}\). Each \(\beta\) now yields a closed \(2\times 2\) linear system, i.e. the azimuthal coupling has been diagonalized in Fourier space.
Intro to the Discrete Fourier Transform (DFT)
The discrete Fourier transform (DFT) maps a finite sequence of \(N\) equally spaced samples into \(N\) complex “frequency bins” describing how much of each discrete complex sinusoid is present. Because the input is finite and indexed modulo \(N\), these sinusoids and the resulting spectra are naturally periodic, and convolutions become circular.
Let \[ x[n],\quad n=0,1,\ldots,N-1. \]
DFT: \[ X[k] = \sum_{n=0}^{N-1} x[n]\, e^{-j 2\pi k n / N},\qquad k=0,\ldots,N-1. \]
Inverse DFT: \[ x[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k]\, e^{j 2\pi k n / N}. \]
Introducing the twiddle factor \(W_N = e^{-j 2\pi / N}\), we can write \[ X[k] = \sum_{n=0}^{N-1} x[n]\; W_N^{k n}. \]
Matrix form.
Define the DFT matrix \(F_N\) by \(\bigl[F_N\bigr]_{k,n} = e^{-j 2\pi k n / N}\). Then \[
X = F_N x,\qquad F_N^{-1} = \frac{1}{N} F_N^*,
\] where \((\cdot)^*\) denotes conjugate transpose.
Core Properties
Let \(\mathcal{F}\{\cdot\}\) denote the DFT operator and \(\mathcal{F}^{-1}\{\cdot\}\) its inverse.
Linearity
\[ \mathcal{F}\{a x + b y\} = a X + b Y. \]Circular shift
If \(y[n] = x\bigl((n - n_0)\bmod N\bigr)\), then \[ Y[k] = X[k]\, e^{-j 2\pi k n_0 / N}. \]Modulation (frequency shift)
If \(y[n] = x[n]\, e^{j 2\pi k_0 n / N}\), then \[ Y[k] = X\bigl((k - k_0)\bmod N\bigr). \]Conjugate symmetry (real sequences)
If \(x[n]\) is real-valued, then \[ X[k] = \overline{X[-k \bmod N]}, \] i.e. the spectrum is Hermitian.Circular convolution theorem
For circular convolution \((x \otimes h)[n]\), \[ \mathcal{F}\{x \otimes h\} = X \cdot H,\qquad \mathcal{F}^{-1}\{X \cdot H\} = x \otimes h, \] where \(\cdot\) denotes pointwise multiplication in the frequency index.