diff --git a/notes/ewald.tex b/notes/ewald.tex
deleted file mode 100644
index d595371..0000000
--- a/notes/ewald.tex
+++ /dev/null
@@ -1,598 +0,0 @@
-%% LyX 2.1.2 created this file.  For more info, see http://www.lyx.org/.
-%% Do not edit unless you really know what you are doing.
-\documentclass[10pt,english]{article}
-\usepackage{amsmath}
-\usepackage{amssymb}
-\usepackage{fontspec}
-\usepackage{unicode-math}
-\setmainfont[Mapping=tex-text,Numbers=OldStyle]{TeX Gyre Pagella}
-\usepackage[a5paper]{geometry}
-\geometry{verbose,tmargin=2cm,bmargin=2cm,lmargin=1cm,rmargin=1cm}
-\usepackage{array}
-\usepackage{multirow}
-\usepackage{esint}
-\usepackage[unicode=true,
- bookmarks=true,bookmarksnumbered=false,bookmarksopen=false,
- breaklinks=false,pdfborder={0 0 1},backref=false,colorlinks=false]
- {hyperref}
-\hypersetup{pdftitle={Accelerating lattice mode calculations with T-matrix method},
- pdfauthor={Marek Nečada}}
-
-\makeatletter
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% LyX specific LaTeX commands.
-\newcommand{\lyxmathsym}[1]{\ifmmode\begingroup\def\b@ld{bold}
-  \text{\ifx\math@version\b@ld\bfseries\fi#1}\endgroup\else#1\fi}
-
-%% Because html converters don't know tabularnewline
-\providecommand{\tabularnewline}{\\}
-
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% User specified LaTeX commands.
-\usepackage{unicode-math}
-
-% Toto je trik, jimž se z fontspec získá familyname pro následující
-\ExplSyntaxOn
-\DeclareExpandableDocumentCommand{\getfamilyname}{m}
- {
-  \use:c { g__fontspec_ \cs_to_str:N #1 _family }
- }
-\ExplSyntaxOff
-
-% definujeme novou rodinu, jež se volá pomocí \MyCyr pro běžné použití, avšak pro účely \DeclareSymbolFont je nutno získat název pomocí getfamilyname definovaného výše
-\newfontfamily\MyCyr{CMU Serif}
-
-\DeclareSymbolFont{cyritletters}{EU1}{\getfamilyname\MyCyr}{m}{it}
-\newcommand{\makecyrmathletter}[1]{%
-  \begingroup\lccode`a=#1\lowercase{\endgroup
-  \Umathcode`a}="0 \csname symcyritletters\endcsname\space #1
-}
-\count255="409
-\loop\ifnum\count255<"44F
-  \advance\count255 by 1
-  \makecyrmathletter{\count255}
-\repeat
-
-\renewcommand{\lyxmathsym}[1]{#1}
-
-\usepackage{polyglossia}
-\setmainlanguage{english}
-\setotherlanguage{russian}
-\newfontfamily\russianfont[Script=Cyrillic]{URW Palladio L}
-%\newfontfamily\russianfont[Script=Cyrillic]{DejaVu Sans}
-
-\makeatother
-
-\usepackage{xunicode}
-\usepackage{polyglossia}
-\setdefaultlanguage{english}
-\begin{document}
-\global\long\def\uoft#1{\mathfrak{F}#1}
-\global\long\def\uaft#1{\mathfrak{\mathbb{F}}#1}
-\global\long\def\usht#1#2{\mathbb{S}_{#1}#2}
-\global\long\def\bsht#1#2{\mathrm{S}_{#1}#2}
-\global\long\def\pht#1#2{\mathfrak{\mathbb{H}}_{#1}#2}
-\global\long\def\vect#1{\mathbf{#1}}
-\global\long\def\ud{\mathrm{d}}
-\global\long\def\basis#1{\mathfrak{#1}}
-\global\long\def\dc#1{\lyxmathsym{Ш}_{#1}}
-\global\long\def\rec#1{#1^{-1}}
-\global\long\def\recb#1{#1^{\widehat{-1}}}
-\global\long\def\ints{\mathbb{Z}}
-\global\long\def\nats{\mathbb{N}}
-\global\long\def\reals{\mathbb{R}}
-\global\long\def\ush#1#2{Y_{#1,#2}}
-\global\long\def\hgfr{\mathbf{F}}
-\global\long\def\hgf{F}
-\global\long\def\ph{\mathrm{ph}}
-\global\long\def\kor#1{\underline{#1}}
-\global\long\def\koru#1{\utilde{#1}}
-
-
-
-\title{Accelerating lattice mode calculations with $T$-matrix method}
-
-
-\author{Marek Nečada}
-\maketitle
-\begin{abstract}
-The $T$-matrix approach is the method of choice for simulating optical
-response of a reasonably small system of compact linear scatterers
-on isotropic background. However, its direct utilisation for problems
-with infinite lattices is problematic due to slowly converging sums
-over the lattice. Here I develop a way to compute the problematic
-sums in the reciprocal space, making the $T$-matrix method very suitable
-for infinite periodic systems as well.
-\end{abstract}
-
-\section{Formulation of the problem}
-
-Assume a system of compact EM scatterers in otherwise homogeneous
-and isotropic medium, and assume that the system, i.e. both the medium
-and the scatterers, have linear response. A scattering problem in
-such system can be written as
-\[
-A_{\alpha}=T_{\alpha}P_{\alpha}=T_{\alpha}(\sum_{\beta}S_{\alpha\leftarrow\beta}A_{\beta}+P_{0\alpha})
-\]
-where $T_{\alpha}$ is the $T$-matrix for scatterer α, $A_{\alpha}$
-is its vector of the scattered wave expansion coefficient (the multipole
-indices are not explicitely indicated here) and $P_{\alpha}$ is the
-local expansion of the incoming sources. $S_{\alpha\leftarrow\beta}$
-is ... and ... is ...
-
-...
-
-\[
-\sum_{\beta}(\delta_{\alpha\beta}-T_{\alpha}S_{\alpha\leftarrow\beta})A_{\beta}=T_{\alpha}P_{0\alpha}.
-\]
-
-
-Now suppose that the scatterers constitute an infinite lattice
-
-\[
-\sum_{\vect b\beta}(\delta_{\vect{ab}}\delta_{\alpha\beta}-T_{\vect a\alpha}S_{\vect a\alpha\leftarrow\vect b\beta})A_{\vect b\beta}=T_{\vect a\alpha}P_{0\vect a\alpha}.
-\]
-Due to the periodicity, we can write $S_{\vect a\alpha\leftarrow\vect b\beta}=S_{\alpha\leftarrow\beta}(\vect b-\vect a)$
-and $T_{\vect a\alpha}=T_{\alpha}$. In order to find lattice modes,
-we search for solutions with zero RHS
-\[
-\sum_{\vect b\beta}(\delta_{\vect{ab}}\delta_{\alpha\beta}-T_{\alpha}S_{\vect a\alpha\leftarrow\vect b\beta})A_{\vect b\beta}=0
-\]
-and we assume periodic solution $A_{\vect b\beta}(\vect k)=A_{\vect a\beta}e^{i\vect k\cdot\vect r_{\vect b-\vect a}}$,
-yielding
-\begin{eqnarray*}
-\sum_{\vect b\beta}(\delta_{\vect{ab}}\delta_{\alpha\beta}-T_{\alpha}S_{\vect a\alpha\leftarrow\vect b\beta})A_{\vect a\beta}\left(\vect k\right)e^{i\vect k\cdot\vect r_{\vect b-\vect a}} & = & 0,\\
-\sum_{\vect b\beta}(\delta_{\vect{0b}}\delta_{\alpha\beta}-T_{\alpha}S_{\vect 0\alpha\leftarrow\vect b\beta})A_{\vect 0\beta}\left(\vect k\right)e^{i\vect k\cdot\vect r_{\vect b}} & = & 0,\\
-\sum_{\beta}(\delta_{\alpha\beta}-T_{\alpha}\underbrace{\sum_{\vect b}S_{\vect 0\alpha\leftarrow\vect b\beta}e^{i\vect k\cdot\vect r_{\vect b}}}_{W_{\alpha\beta}(\vect k)})A_{\vect 0\beta}\left(\vect k\right) & = & 0,\\
-A_{\vect 0\alpha}\left(\vect k\right)-T_{\alpha}\sum_{\beta}W_{\alpha\beta}\left(\vect k\right)A_{\vect 0\beta}\left(\vect k\right) & = & 0.
-\end{eqnarray*}
-Therefore, in order to solve the modes, we need to compute the ``lattice
-Fourier transform'' of the translation operator,
-\begin{equation}
-W_{\alpha\beta}(\vect k)\equiv\sum_{\vect b}S_{\vect 0\alpha\leftarrow\vect b\beta}e^{i\vect k\cdot\vect r_{\vect b}}.\label{eq:W definition}
-\end{equation}
-
-
-
-\section{Computing the Fourier sum of the translation operator}
-
-The problem evaluating (\ref{eq:W definition}) is the asymptotic
-behaviour of the translation operator, $S_{\vect 0\alpha\leftarrow\vect b\beta}\sim\left|\vect r_{\vect b}\right|^{-1}e^{ik_{0}\left|\vect r_{\vect b}\right|}$
-that makes the convergence of the sum quite problematic for any $d>1$-dimensional
-lattice.%
-\footnote{Note that $d$ here is dimensionality of the lattice, not the space
-it lies in, which I for certain reasons assume to be three. (TODO
-few notes on integration and reciprocal lattices in some appendix)%
-} In electrostatics, one can solve this problem with Ewald summation.
-Its basic idea is that if what asymptoticaly decays poorly in the
-direct space, will perhaps decay fast in the Fourier space. I use
-the same idea here, but everything will be somehow harder than in
-electrostatics.
-
-Let us re-express the sum in (\ref{eq:W definition}) in terms of
-integral with a delta comb
-
-\begin{equation}
-W_{\alpha\beta}(\vect k)=\int\ud^{d}\vect r\dc{\basis u}(\vect r)S(\vect r_{\alpha}\leftarrow\vect r+\vect r_{\beta})e^{i\vect k\cdot\vect r}.\label{eq:W integral}
-\end{equation}
-The translation operator $S$ is now a function defined in the whole
-3d space; $\vect r_{\alpha},\vect r_{\beta}$ are the displacements
-of scatterers $\alpha$ and $\beta$ in a unit cell. The arrow notation
-$S(\vect r_{\alpha}\leftarrow\vect r+\vect r_{\beta})$ means ``translation
-operator for spherical waves originating in $\vect r+\vect r_{\beta}$
-evaluated in $\vect r_{\alpha}$'' and obviously $S$ is in fact
-a function of a single 3d argument, $S(\vect r_{\alpha}\leftarrow\vect r+\vect r_{\beta})=S(\vect 0\leftarrow\vect r+\vect r_{\beta}-\vect r_{\alpha})=S(-\vect r-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)=S(-\vect r-\vect r_{\beta}+\vect r_{\alpha})$.
-Expression (\ref{eq:W integral}) can be rewritten as
-\[
-W_{\alpha\beta}(\vect k)=\left(2\pi\right)^{\frac{d}{2}}\uaft{(\dc{\basis u}S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0))\left(\vect k\right)}
-\]
-where changed the sign of $\vect r/\vect{\bullet}$ has been swapped
-under integration, utilising evenness of $\dc{\basis u}$. Fourier
-transform of product is convolution of Fourier transforms, so (using
-formula (\ref{eq:Dirac comb uaFt}) for the Fourier transform of Dirac
-comb)
-\begin{eqnarray}
-W_{\alpha\beta}(\vect k) & = & \left(\left(\uaft{\dc{\basis u}}\right)\ast\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\right)(\vect k)\nonumber \\
- & = & \frac{\left|\det\recb{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\left(\dc{\recb{\basis u}}^{(d)}\ast\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\right)\left(\vect k\right)\nonumber \\
- & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\sum_{\vect K\in\recb{\basis u}\ints^{d}}\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\left(\vect k-\vect K\right).\label{eq:W sum in reciprocal space}
-\end{eqnarray}
-As such, this is not extremely helpful because the the \emph{whole}
-translation operator $S$ has singularities in origin, hence its Fourier
-transform $\uaft S$ will decay poorly. 
-
-However, Fourier transform is linear, so we can in principle separate
-$S$ in two parts, $S=S^{\textup{L}}+S^{\textup{S}}$. $S^{\textup{S}}$
-is a short-range part that decays sufficiently fast with distance
-so that its direct-space lattice sum converges well; $S^{\textup{S}}$
-must as well contain all the singularities of $S$ in the origin.
-The other part, $S^{\textup{L}}$, will retain all the slowly decaying
-terms of $S$ but it also has to be smooth enough in the origin, so
-that its Fourier transform $\uaft{S^{\textup{L}}}$ decays fast enough.
-(The same idea lies behind the Ewald summation in electrostatics.)
-Using the linearity of Fourier transform and formulae (\ref{eq:W definition})
-and (\ref{eq:W sum in reciprocal space}), the operator $W_{\alpha\beta}$
-can then be re-expressed as
-\begin{eqnarray}
-W_{\alpha\beta}\left(\vect k\right) & = & W_{\alpha\beta}^{\textup{S}}\left(\vect k\right)+W_{\alpha\beta}^{\textup{L}}\left(\vect k\right)\nonumber \\
-W_{\alpha\beta}^{\textup{S}}\left(\vect k\right) & = & \sum_{\vect R\in\basis u\ints^{d}}S^{\textup{S}}(\vect 0\leftarrow\vect R+\vect r_{\beta}-\vect r_{\alpha})e^{i\vect k\cdot\vect R}\label{eq:W Short definition}\\
-W_{\alpha\beta}^{\textup{L}}\left(\vect k\right) & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\sum_{\vect K\in\recb{\basis u}\ints^{d}}\left(\uaft{S^{\textup{L}}(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\left(\vect k-\vect K\right)\label{eq:W Long definition}
-\end{eqnarray}
-where both sums should converge nicely.
-
-
-\section{Finding a good decomposition}
-
-The remaining challenge is therefore finding a suitable decomposition
-$S^{\textup{L}}+S^{\textup{S}}$ such that both $S^{\textup{S}}$
-and $\uaft{S^{\textup{L}}}$ decay fast enough with distance and are
-expressable analytically. With these requirements, I do not expect
-to find gaussian asymptotics as in the electrostatic Ewald formula—having
-$\sim x^{-t}$, $t>d$ asymptotics would be nice, making the sums
-in (\ref{eq:W Short definition}), (\ref{eq:W Long definition}) absolutely
-convergent.
-
-The translation operator $S$ for compact scatterers in 3d can be
-expressed as
-\[
-S_{l',m',t'\leftarrow l,m,t}\left(\vect r\leftarrow\vect 0\right)=\sum_{p}c_{p}^{l',m',t'\leftarrow l,m,t}\ush p{m'-m}\left(\theta_{\vect r},\phi_{\vect r}\right)z_{p}^{(J)}\left(k_{0}\left|\vect r\right|\right)
-\]
-where $Y_{l,m}\left(\theta,\phi\right)$ are the spherical harmonics,
-$z_{p}^{(J)}\left(r\right)$ some of the Bessel or Hankel functions
-(probably $h_{p}^{(1)}$ in all meaningful cases; TODO) and $c_{p}^{l,m,t\leftarrow l',m',t'}$
-are some ugly but known coefficients (REF Xu 1996, eqs. 76,77). 
-
-The spherical Hankel functions can be expressed analytically as (REF
-DLMF 10.49.6, 10.49.1) 
-\begin{equation}
-h_{n}^{(1)}(r)=e^{ir}\sum_{k=0}^{n}\frac{i^{k-n-1}}{r^{k+1}}\frac{\left(n+k\right)!}{2^{k}k!\left(n-k\right)!},\label{eq:spherical Hankel function series}
-\end{equation}
- so if we find a way to deal with the radial functions $s_{k_{0},q}(r)=e^{ik_{0}r}\left(k_{0}r\right)^{-q}$,
-$q=1,2$ in 2d case or $q=1,2,3$ in 3d case, we get absolutely convergent
-summations in the direct space.
-
-
-\subsection{2d}
-
-Assume that all scatterers are placed in the plane $\vect z=0$, so
-that the 2d Fourier transform of the long-range part of the translation
-operator in terms of Hankel transforms, according to (\ref{eq:Fourier v. Hankel tf 2d}),
-reads
-
-\begin{multline*}
-\uaft{S_{l',m',t'\leftarrow l,m,t}^{\textup{L}}\left(\vect{\bullet}\leftarrow\vect 0\right)}(\vect k)=\\
-\sum_{p}c_{p}^{l',m',t'\leftarrow l,m,t}\ush p{m'-m}\left(\frac{\pi}{2},0\right)e^{i(m'-m)\phi}i^{m'-m}\pht{m'-m}{h_{p}^{(1)\textup{L}}\left(k_{0}\vect{\bullet}\right)}\left(\left|\vect k\right|\right)
-\end{multline*}
-Here $h_{p}^{(1)\textup{L}}=h_{p}^{(1)}-h_{p}^{(1)\textup{S}}$ is
-a long range part of a given spherical Hankel function which has to
-be found and which contains all the terms with far-field ($r\to\infty$)
-asymptotics proportional to$\sim e^{ik_{0}r}\left(k_{0}r\right)^{-q}$,
-$q\le Q$ where $Q$ is at least two in order to achieve absolute
-convergence of the direct-space sum, but might be higher in order
-to speed the convergence up.
-
-Obviously, all the terms $\propto s_{k_{0},q}(r)=e^{ik_{0}r}\left(k_{0}r\right)^{-q}$,
-$q>Q$ of the spherical Hankel function (\ref{eq:spherical Hankel function series})
-can be kept untouched as part of $h_{p}^{(1)\textup{S}}$, as they
-decay fast enough.
-
-The remaining task is therefore to find a suitable decomposition of
-$s_{k_{0},q}(r)=e^{ik_{0}r}\left(k_{0}r\right)^{-q}$, $q\le Q$ into
-short-range and long-range parts, $s_{k_{0},q}(r)=s_{k_{0},q}^{\textup{S}}(r)+s_{k_{0},q}^{\textup{L}}(r)$,
-such that $s_{k_{0},q}^{\textup{L}}(r)$ contains all the slowly decaying
-asymptotics and its Hankel transforms decay desirably fast as well,
-$\pht n{s_{k_{0},q}^{\textup{L}}}\left(k\right)=o(z^{-Q})$, $z\to\infty$.
-The latter requirement calls for suitable regularisation functions—$s_{q}^{\textup{L}}$
-must be sufficiently smooth in the origin, so that 
-\begin{equation}
-\pht n{s_{k_{0},q}^{\textup{L}}}\left(k\right)=\int_{0}^{\infty}s_{k_{0},q}^{\textup{L}}\left(r\right)rJ_{n}\left(kr\right)\ud r=\int_{0}^{\infty}s_{k_{0},q}\left(r\right)\rho\left(r\right)rJ_{n}\left(kr\right)\ud r\label{eq:2d long range regularisation problem statement}
-\end{equation}
- exists and decays fast enough. $J_{\nu}(r)\sim\left(r/2\right)^{\nu}/\Gamma\left(\nu+1\right)$
-(REF DLMF 10.7.3) near the origin, so the regularisation function
-should be $\rho(r)=o(r^{q-n-1})$ only to make $\pht n{s_{q}^{\textup{L}}}$
-converge. The additional decay speed requirement calls for at least
-$\rho(r)=o(r^{q-n+Q-1})$, I guess. At the same time, $\rho(r)$ must
-converge fast enough to one for $r\to\infty$.
-
-The electrostatic Ewald summation uses regularisation with $1-e^{-cr^{2}}$.
-However, such choice does not seem to lead to an analytical solution
-(really? could not something be dug out of DLMF 10.22.54?) for the
-current problem (\ref{eq:2d long range regularisation problem statement}).
-But it turns out that the family of functions
-\begin{equation}
-\rho_{\kappa,c}(r)\equiv\left(1-e^{-cr}\right)^{\text{\ensuremath{\kappa}}},\quad c>0,\kappa\in\nats\label{eq:binom regularisation function}
-\end{equation}
-might lead to satisfactory results; see below. 
-
-
-\subsubsection{Hankel transforms of the long-range parts, „binomial“ regularisation\label{sub:Hankel-transforms-binom-reg}}
-
-Let 
-
-\begin{eqnarray}
-\pht n{s_{q,k_{0}}^{\textup{L}\kappa,c}}\left(k\right) & \equiv & \int_{0}^{\infty}\frac{e^{ik_{0}r}}{\left(k_{0}r\right)^{q}}J_{n}\left(kr\right)\left(1-e^{-cr}\right)^{\kappa}r\,\ud r\nonumber \\
- & = & k_{0}^{-q}\int_{0}^{\infty}r^{1-q}J_{n}\left(kr\right)\sum_{\sigma=0}^{\kappa}\left(-1\right)^{\sigma}\binom{\kappa}{\sigma}e^{-(\sigma c-ik_{0})r}\ud r\nonumber \\
- & \underset{\equiv}{\textup{form.}} & \sum_{\sigma=0}^{\kappa}\left(-1\right)^{\sigma}\binom{\kappa}{\sigma}\pht n{s_{q,k_{0}}^{\textup{L}1,\sigma c}}\left(k\right).\label{eq:2D Hankel transform of regularized outgoing wave, decomposition}
-\end{eqnarray}
-From {[}REF DLMF 10.22.49{]} one digs 
-\begin{multline}
-\pht n{s_{q,k_{0}}^{\textup{L}1,c}}\left(k\right)=\frac{k^{n}\Gamma\left(2-q+n\right)}{2^{n}k_{0}^{q}\left(c-ik_{0}\right)^{2-q+n}}\hgfr\left(\frac{2-q+n}{2},\frac{3-q+n}{2};1+n;\frac{-k^{2}}{\left(c-ik_{0}\right)^{2}}\right),\\
-\Re\left(2-q+n\right)>0,\Re(c-ik_{0}\pm k)\ge0,\label{eq:2D Hankel transform of exponentially suppressed outgoing wave as 2F1}
-\end{multline}
-and using {[}REF DLMF 15.9.17{]} and  {[}REF DLMF 14.9.5{]}
-
-{\footnotesize{}
-\begin{multline}
-\pht n{s_{q,k_{0}}^{\textup{L}1,c}}\left(k\right)=\frac{k^{n}\Gamma\left(2-q+n\right)}{k_{0}^{q}\left(c-ik_{0}\right)^{2-q+n}}\left(\frac{-k^{2}}{\left(c-ik_{0}\right)^{2}}\right)^{-\frac{n}{2}}\left(1+\frac{k^{2}}{\left(c-ik_{0}\right)^{2}}\right)^{\frac{q}{2}-1}P_{q}^{-n}\left(\frac{1}{\sqrt{1+\frac{k^{2}}{\left(c-ik_{0}\right)^{2}}}}\right),\\
-k>0\wedge k_{0}>0\wedge c\ge0\wedge\lnot\left(c=0\wedge k_{0}=k\right)\label{eq:2D Hankel transform of exponentially suppressed outgoing wave expanded}
-\end{multline}
-}with principal branches of the hypergeometric functions, associated
-Legendre functions, and fractional powers. The conditions from (\ref{eq:2D Hankel transform of exponentially suppressed outgoing wave as 2F1})
-should hold, but we will use (\ref{eq:2D Hankel transform of exponentially suppressed outgoing wave expanded})
-formally even if they are violated, with the hope that the divergences
-eventually cancel in (\ref{eq:2D Hankel transform of regularized outgoing wave, decomposition}).
-
-One problematic element here is the gamma function $\text{Γ}\left(2-q+n\right)$
-which is singular if the arguments are negative integers, i.e. if
-$q-n\ge3$; but at least the necessary minimum of $q=1,2$ would be
-covered this way. The associated Legendre function can be expressed
-as a finite ``polynomial'' if $q\ge n$. In other cases, different
-expressions can be obtained from \ref{eq:2D Hankel transform of exponentially suppressed outgoing wave as 2F1}
-using various transformation formulae from either DLMF or \begin{russian}Прудников\end{russian}. 
-
-In fact, Mathematica is usually able to calculate the transforms for
-specific values of $\kappa,q,n$, but it did not find any general
-formula for me. The resulting expressions are finite sums of algebraic
-functions, Table \ref{tab:Asymptotical-behaviour-Mathematica} shows
-how fast they decay with growing $k$ for some parameters. The only
-case where Mathematica did not help at all is $q=2,n=0$, which is
-unfortunately important. But if I have not made some mistake, the
-expression (\ref{eq:2D Hankel transform of exponentially suppressed outgoing wave expanded})
-is applicable for this case.
-
-\begin{table}
-\begin{centering}
-{\footnotesize{}}%
-\begin{tabular}{cc|ccc}
-\multicolumn{2}{c|}{{\footnotesize{}$\kappa=0$}} &  & {\footnotesize{}$n$} & \tabularnewline
-\multicolumn{1}{c}{} &  & {\footnotesize{}0} & {\footnotesize{}1} & {\footnotesize{}2}\tabularnewline
-\hline 
-\multirow{2}{*}{{\footnotesize{}$q$}} & {\footnotesize{}1} & {\footnotesize{}2} & {\footnotesize{}1} & {\footnotesize{}1}\tabularnewline
- & {\footnotesize{}2} & {\footnotesize{}x} & {\footnotesize{}w} & {\footnotesize{}0}\tabularnewline
-\end{tabular}{\footnotesize{} \hspace*{\fill}}%
-\begin{tabular}{cc|ccccc}
-\multicolumn{2}{c|}{{\footnotesize{}$\kappa=1$}} & \multicolumn{5}{c}{{\footnotesize{}$n$}}\tabularnewline
- &  & {\footnotesize{}0} & {\footnotesize{}1} & {\footnotesize{}2} & {\footnotesize{}3} & {\footnotesize{}4}\tabularnewline
-\hline 
-\multirow{2}{*}{{\footnotesize{}$q$}} & {\footnotesize{}1} & {\footnotesize{}w} & {\footnotesize{}3} & {\footnotesize{}2} & {\footnotesize{}2} & {\footnotesize{}2}\tabularnewline
- & {\footnotesize{}2} & {\footnotesize{}x} & {\footnotesize{}1} & {\footnotesize{}w} & {\footnotesize{}1} & {\footnotesize{}1}\tabularnewline
-\end{tabular}{\footnotesize{} \hspace*{\fill}}%
-\begin{tabular}{cc|ccccc}
-\multicolumn{2}{c|}{{\footnotesize{}$\kappa=2$}} & \multicolumn{5}{c}{{\footnotesize{}$n$}}\tabularnewline
- &  & {\footnotesize{}0} & {\footnotesize{}1} & {\footnotesize{}2} & {\footnotesize{}3} & {\footnotesize{}4}\tabularnewline
-\hline 
-\multirow{2}{*}{{\footnotesize{}$q$}} & {\footnotesize{}1} & {\footnotesize{}0/w} & {\footnotesize{}3} & {\footnotesize{}4} & {\footnotesize{}3} & {\footnotesize{}3}\tabularnewline
- & {\footnotesize{}2} & {\footnotesize{}x} & {\footnotesize{}3} & {\footnotesize{}2} & {\footnotesize{}2} & {\footnotesize{}1}\tabularnewline
-\end{tabular}
-\par\end{centering}{\footnotesize \par}
-
-\protect\caption{Asymptotical behaviour of some (\ref{eq:2D Hankel transform of regularized outgoing wave, decomposition})
-obtained by Mathematica for $k\to\infty$. The table entries are the
-$N$ of $\protect\pht n{s_{q,k_{0}}^{\textup{L}\kappa,c}}\left(k\right)=o\left(1/k^{N}\right)$.
-The special entry ``x'' means that Mathematica was not able to calculate
-the integral, and ``w'' denotes that the first returned term was
-not simply of the kind $(\ldots)k^{-N-1}$.\label{tab:Asymptotical-behaviour-Mathematica}}
-\end{table}
-
-
-
-
-
-\subsection{3d (TODO)}
-
-\begin{multline*}
-\uaft{S_{l',m',t'\leftarrow l,m,t}\left(\vect{\bullet}\leftarrow\vect 0\right)}(\vect k)=\\
-\sum_{p}c_{p}^{l',m',t'\leftarrow l,m,t}\ush p{m'-m}\left(\theta_{\vect k},\phi_{\vect k}\right)\left(-i\right)^{p}\usht p{z_{p}^{(J)}}\left(\left|\vect k\right|\right)
-\end{multline*}
-
-
-
-\section{Major TODOs and open questions}
-\begin{itemize}
-\item Check if (\ref{eq:2D Hankel transform of exponentially suppressed outgoing wave expanded})
-gives a satisfactory result for the case $q=2,n=0$.
-\item Analyse the behaviour $k\to k_{0}$.
-\item Find a general algorithm for generating the expressions of the Hankel
-transforms.
-\item Three-dimensional case.
-\end{itemize}
-
-\section{(Appendix) Fourier vs. Hankel transform}
-
-
-\subsection{Three dimensions}
-
-Given a nice enough function $f$ of a real 3d variable, assume its
-factorisation into radial and angular parts 
-\[
-f(\vect r)=\sum_{l,m}f_{l,m}(\left|\vect r\right|)\ush lm\left(\theta_{\vect r},\phi_{\vect r}\right).
-\]
-Acording to (REF Baddour 2010, eqs. 13, 16), its Fourier transform
-can then be expressed in terms of Hankel transforms (CHECK normalisation
-of $j_{n}$, REF Baddour (1)) 
-\[
-\uaft f(\vect k)=\frac{4\pi}{\left(2\pi\right)^{\frac{3}{2}}}\sum_{l,m}\left(-i\right)^{l}\left(\bsht{f_{l,m}}{}\right)\left(\left|\vect k\right|\right)\ush lm\left(\theta_{\vect k},\phi_{\vect k}\right)
-\]
-where the spherical Hankel transform $\bsht l{}$ of degree $l$ is
-defined as (REF Baddour eq. 2)
-\[
-\bsht lg(k)\equiv\int_{0}^{\infty}\ud r\, r^{2}g(r)j_{l}\left(kr\right).
-\]
-Using this convention, the inverse spherical Hankel transform is given
-by (REF Baddour eq. 3)
-\[
-g(r)=\frac{2}{\pi}\int_{0}^{\infty}\ud k\, k^{2}\bsht lg(k)j_{l}(k),
-\]
-so it is not unitary. 
-
-An unitary convention would look like this:
-\begin{equation}
-\usht lg(k)\equiv\sqrt{\frac{2}{\pi}}\int_{0}^{\infty}\ud r\, r^{2}g(r)j_{l}\left(kr\right).\label{eq:unitary 3d Hankel tf definition}
-\end{equation}
-Then $\usht l{}^{-1}=\usht l{}$ and the unitary, angular-momentum
-Fourier transform reads
-\begin{eqnarray}
-\uaft f(\vect k) & = & \frac{4\pi}{\left(2\pi\right)^{\frac{3}{2}}}\sqrt{\frac{\pi}{2}}\sum_{l,m}\left(-i\right)^{l}\left(\usht l{f_{l,m}}\right)\left(\left|\vect k\right|\right)\ush lm\left(\theta_{\vect k},\phi_{\vect k}\right)\nonumber \\
- & = & \sum_{l,m}\left(-i\right)^{l}\left(\usht l{f_{l,m}}\right)\left(\left|\vect k\right|\right)\ush lm\left(\theta_{\vect k},\phi_{\vect k}\right).\label{eq:Fourier v. Hankel tf 3d}
-\end{eqnarray}
-Cool.
-
-
-\subsection{Two dimensions}
-
-Similarly in 2d, let the expansion of $f$ be 
-\[
-f\left(\vect r\right)=\sum_{m}f_{m}\left(\left|\vect r\right|\right)e^{im\phi_{\vect r}},
-\]
-its Fourier transform is then (CHECK this, it is taken from the Wikipedia
-article on Hankel transform) 
-\begin{equation}
-\uaft f\left(\vect k\right)=\sum_{m}i^{m}e^{im\phi_{\vect k}}\pht mf_{m}\left(\left|\vect k\right|\right)\label{eq:Fourier v. Hankel tf 2d}
-\end{equation}
-where the Hankel transform of order $m$ is defined as
-\begin{equation}
-\pht mg\left(k\right)=\int_{0}^{\infty}\ud r\, g(r)J_{m}(kr)r\label{eq:unitary 2d Hankel tf definition}
-\end{equation}
-which is already self-inverse, $\pht m{}^{-1}=\pht m{}$ (hence also
-unitary).
-
-
-\section{(Appendix) Multidimensional Dirac comb}
-
-
-\subsection{1D}
-
-This is all from Wikipedia
-
-
-\subsubsection{Definitions}
-
-\begin{eqnarray*}
-\lyxmathsym{Ш}(t) & \equiv & \sum_{k=-\infty}^{\infty}\delta(t-k)\\
-\lyxmathsym{Ш}_{T}(t) & \equiv & \sum_{k=-\infty}^{\infty}\delta(t-kT)=\frac{1}{T}\lyxmathsym{Ш}\left(\frac{t}{T}\right)
-\end{eqnarray*}
-
-
-
-\subsubsection{Fourier series representation}
-
-\begin{equation}
-\lyxmathsym{Ш}_{T}(t)=\sum_{n=-\infty}^{\infty}e^{2\pi int/T}\label{eq:1D Dirac comb Fourier series}
-\end{equation}
-
-
-
-\subsubsection{Fourier transform}
-
-With unitary ordinary frequency Ft., i.e.
-
-\[
-\uoft f(\vect{\xi})\equiv\int_{\mathbb{R}^{n}}f(\vect x)e^{-2\pi i\vect x\cdot\vect{\xi}}\ud^{n}\vect x
-\]
-we have 
-\begin{equation}
-\uoft{\lyxmathsym{Ш}_{T}}(f)=\frac{1}{T}\lyxmathsym{Ш}_{\frac{1}{T}}(f)=\sum_{n=-\infty}^{\infty}e^{-i2\pi fnT}\label{eq:1D Dirac comb Ft ordinary freq}
-\end{equation}
- and with unitary angular frequency Ft., i.e.
-\begin{equation}
-\uaft f(\vect k)\equiv\frac{1}{\left(2\pi\right)^{n/2}}\int_{\mathbb{R}^{n}}f(\vect x)e^{-i\vect x\cdot\vect k}\ud^{n}\vect x\label{eq:Ft unitary angular frequency}
-\end{equation}
-we have (CHECK)
-\[
-\uaft{\lyxmathsym{Ш}_{T}}(\omega)=\frac{\sqrt{2\pi}}{T}\lyxmathsym{Ш}_{\frac{2\pi}{T}}(\omega)=\frac{1}{\sqrt{2\pi}}\sum_{n=-\infty}^{\infty}e^{-i\omega nT}
-\]
-
-
-
-\subsection{Dirac comb for multidimensional lattices}
-
-
-\subsubsection{Definitions}
-
-Let $d$ be the dimensionality of the real vector space in question,
-and let $\basis u\equiv\left\{ \vect u_{i}\right\} _{i=1}^{d}$ denote
-a basis for some lattice in that space. Let the corresponding lattice
-delta comb be
-\[
-\dc{\basis u}\left(\vect x\right)\equiv\sum_{n_{1}=-\infty}^{\infty}\ldots\sum_{n_{d}=-\infty}^{\infty}\delta\left(\vect x-\sum_{i=1}^{d}n_{i}\vect u_{i}\right).
-\]
-
-
-Furthemore, let $\rec{\basis u}\equiv\left\{ \rec{\vect u}_{i}\right\} _{i=1}^{d}$
-be the reciprocal lattice basis, that is the basis satisfying $\vect u_{i}\cdot\rec{\vect u_{j}}=\delta_{ij}$.
-This slightly differs from the usual definition of a reciprocal basis,
-here denoted $\recb{\basis u}\equiv\left\{ \recb{\vect u_{i}}\right\} _{i=1}^{d}$,
-which satisfies $\vect u_{i}\cdot\recb{\vect u_{j}}=2\pi\delta_{ij}$
-instead.
-
-
-\subsubsection{Factorisation of a multidimensional lattice delta comb}
-
-By simple drawing, it can be seen that 
-\[
-\dc{\basis u}(\vect x)=c_{\basis u}\prod_{i=1}^{d}\dc{}\left(\vect x\cdot\rec{\vect u_{i}}\right)
-\]
-where $c_{\basis u}$ is some numerical volume factor. In order to
-determine $c_{\basis u}$, let us consider only the ``zero tooth''
-of the comb, leading to
-\[
-\delta^{d}(\vect x)=c_{\basis u}\prod_{i=1}^{d}\delta\left(\vect x\cdot\rec{\vect u_{i}}\right).
-\]
-From the scaling property of delta function, $\delta(ax)=\left|a\right|^{-1}\delta(x)$,
-we get
-\[
-\delta^{d}(\vect x)=c_{\basis u}\prod_{i=1}^{d}\left\Vert \rec{\vect u_{i}}\right\Vert ^{-1}\delta\left(\vect x\cdot\frac{\rec{\vect u_{i}}}{\left\Vert \rec{\vect u_{i}}\right\Vert }\right).
-\]
-
-
-From the Osgood's book (p. 375):
-
-\[
-\dc A(\vect x)=\frac{1}{\left|\det A\right|}\dc{}^{(d)}\left(A^{-1}\vect x\right)
-\]
-
-
-
-\subsubsection{Fourier series representation}
-
-
-
-
-\subsubsection{Fourier transform (OK)}
-
-From the Osgood's book https://see.stanford.edu/materials/lsoftaee261/chap8.pdf,
-p. 379
-
-\[
-\uoft{\dc{\basis u}}\left(\vect{\xi}\right)=\left|\det\rec{\basis u}\right|\dc{\rec{\basis u}}^{(d)}\left(\vect{\xi}\right).
-\]
-And consequently, for unitary/angular frequency it is
-
-\begin{eqnarray}
-\uaft{\dc{\basis u}}\left(\vect k\right) & = & \frac{1}{\left(2\pi\right)^{\frac{d}{2}}}\uoft{\dc{\basis u}}\left(\frac{\vect k}{2\pi}\right)\nonumber \\
- & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\dc{\rec{\basis u}}^{(d)}\left(\frac{\vect k}{2\pi}\right)\nonumber \\
- & = & \left(2\pi\right)^{\frac{d}{2}}\left|\det\rec{\basis u}\right|\dc{\recb{\basis u}}\left(\vect k\right)\nonumber \\
- & = & \frac{\left|\det\recb{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\dc{\recb{\basis u}}\left(\vect k\right).\label{eq:Dirac comb uaFt}
-\end{eqnarray}
-
-
-
-
-
-\subsubsection{Convolution}
-
-\[
-\left(f\ast\dc{\basis u}\right)(\vect x)=\sum_{\vect t\in\basis u\ints^{d}}f(\vect x-\vect t)
-\]
-
-
-
-\end{document}