## 5 Non-Uniform Discrete Fourier Transforms  
The main computational cost of the Sinkhorn algorithm comes from the matrix-vector product with the kernel matrix (8). Let \(k,\ell\in[K]\) and \(\alpha\in\mathbb{R}^{n_{\ell}}\). We briefly describe a fast summation method for the computation of \(\beta=\mathbf{\mathcal{K}}^{(k,\ell)}\alpha\), i.e.,  
\[\beta_{i_{k}}=\sum_{i_{\ell}=1}^{n_{\ell}}\alpha_{i_{\ell}}\,\exp\left(-\tfrac{ 1}{\eta}\left\|\mathbf{x}_{i_{k}}^{k}-\mathbf{x}_{i_{\ell}}^{\ell}\right\|_{2}^{2} \right),\quad i_{k}\in[n_{k}]. \tag{11}\]We refer to [32] for a detailed derivation and error estimates. The main idea is to approximate the kernel function  
\[\kappa(x)\coloneqq\exp(-\tfrac{1}{\eta}x^{2}),\quad x\in\mathbb{R}, \tag{12}\]  
using a Fourier series. In order to ensure the fast convergence of the Fourier series, we extend \(\kappa\) to a periodic function of certain smoothness \(p\in\mathbb{N}\). Let \(\epsilon_{\mathrm{B}}>0\) and \(\tau>\epsilon_{\mathrm{B}}+\max\left\{\|x_{k}^{k}-x_{\ell}^{\ell}\|_{2}^{2}: i_{k}\in[n_{k}],i_{\ell}\in[n_{\ell}]\right\}\). For \(x\in\mathbb{R}\), we define the regularized kernel  
\[\kappa_{\mathrm{R}}(x)\coloneqq \begin{cases}\kappa(x),&|x|\leq\tau-\epsilon_{ \mathrm{B}},\\ \kappa_{\mathrm{B}}(x),&\tau-\epsilon_{\mathrm{B}}<|x|\leq\tau,\\ \kappa_{\mathrm{B}}(\tau),&|x|>\tau,\end{cases}  \tag{13}\]  
where \(\kappa_{\mathrm{B}}\) is a polynomial of degree \(p\) that fulfills the Herminte interpolation conditions \(\kappa_{\mathrm{B}}^{(r)}(\tau-\epsilon_{\mathrm{B}})=\kappa^{(r)}(\tau- \epsilon_{\mathrm{B}})\) for \(r=0,\ldots,p-1\), and \(\kappa_{\mathrm{B}}^{(r)}(\tau)=0\) for \(r=1,\ldots,p-1\); see Figure 1.  
Then, we define a \(2\tau\)-periodic function on \(\mathbb{R}^{d}\) by  
\[\tilde{\kappa}(x)\coloneqq\kappa_{\mathrm{R}}(\|x\|),\quad x\in[-\tau,\tau)^ {d}. \tag{14}\]  
By construction, \(\tilde{\kappa}\) is \(p\) times continuously differentiable and we have \(\tilde{\kappa}(x_{k}^{k}-x_{i_{\ell}}^{\ell})=\exp(-\tfrac{1}{\eta}\|x_{k}^{k }-x_{i_{\ell}}^{\ell}\|_{2}^{2})\) for all \(i_{k}\in[n_{k}]\), \(i_{\ell}\in[n_{\ell}]\). Let \(M\in\mathbb{N}\). We approximate \(\tilde{\kappa}\) using the \(2\tau\)-periodic Fourier expansion of degree \(2M\),  
\[\tilde{\kappa}(x)\approx\sum_{m\in\{-M,\ldots,M-1\}^{d}}\hat{\kappa}(m)\,\exp (\mathrm{i}\tfrac{\pi}{\tau}m^{\intercal}x),\quad x\in\mathbb{R}^{d},\]  
with the discrete Fourier coefficients \(\hat{\kappa}(m)\in\mathbb{C}\), which can be efficiently approximated by the fast Fourier transform (FFT)  
\[\hat{\kappa}(m)\coloneqq\frac{1}{(2M)^{d}}\sum_{x\in\tfrac{\pi}{\eta}\{-M, \ldots,M-1\}^{d}}\hat{\kappa}(x)\,\exp(-\mathrm{i}\tfrac{\pi}{\tau}m^{ \intercal}x),\quad m\in\{-M,\ldots,M-1\}^{d}. \tag{15}\]  
This yields an approximation of (11) by  
\[\beta_{i_{k}}=\sum_{j=1}^{n_{\ell}}\tilde{\kappa}(x_{i_{k}}^{k}-x _{i_{\ell}}^{\ell})\,\alpha_{i_{\ell}}\approx\sum_{i_{\ell}=1}^{n_{\ell}}\sum _{m\in\{-M,\ldots,M-1\}^{d}}\hat{\kappa}(m)\,\exp\Bigl{(}\mathrm{i}\tfrac{ \pi}{\tau}m^{\intercal}(x_{i_{k}}^{k}-x_{i_{\ell}}^{\ell})\Bigr{)}\,\alpha_{i_ {\ell}}\\ =\sum_{m\in\{-M,\ldots,M-1\}^{d}}\hat{\kappa}(m)\left(\sum_{i_{ \ell}=1}^{n_{\ell}}\alpha_{i_{\ell}}\,\exp\Bigl{(}-\mathrm{i}\tfrac{\pi}{\tau }m^{\intercal}x_{i_{\ell}}^{\ell}\Bigr{)}\right)\,\exp\Bigl{(}\mathrm{i} \tfrac{\pi}{\tau}m^{\intercal}x_{i_{\ell}}^{k}\Bigr{)}. \tag{16}\]  
Figure 1: Regularized kernel \(\kappa_{\mathrm{R}}\) for \(\eta=1/4\), periodicity length \(\tau=1\), boundary interval \(\epsilon_{\mathrm{B}}=0.2\), and smoothness \(p=1\).