diff --git a/lepaper/T-matrix paper.bib b/lepaper/T-matrix paper.bib
deleted file mode 100644
index 2fc6909..0000000
--- a/lepaper/T-matrix paper.bib	
+++ /dev/null
@@ -1,60 +0,0 @@
-
-@article{xu_efficient_1998,
-  title = {Efficient {{Evaluation}} of {{Vector Translation Coefficients}} in {{Multiparticle Light}}-{{Scattering Theories}}},
-  volume = {139},
-  issn = {0021-9991},
-  abstract = {Vector addition theorems are a necessary ingredient in the analytical solution of electromagnetic multiparticle-scattering problems. These theorems include a large number of vector addition coefficients. There exist three basic types of analytical expressions for vector translation coefficients: Stein's (Quart. Appl. Math.19, 15 (1961)), Cruzan's (Quart. Appl. Math.20, 33 (1962)), and Xu's (J. Comput. Phys.127, 285 (1996)). Stein's formulation relates vector translation coefficients with scalar translation coefficients. Cruzan's formulas use the Wigner 3jm symbol. Xu's expressions are based on the Gaunt coefficient. Since the scalar translation coefficient can also be expressed in terms of the Gaunt coefficient, the key to the expeditious and reliable calculation of vector translation coefficients is the fast and accurate evaluation of the Wigner 3jm symbol or the Gaunt coefficient. We present highly efficient recursive approaches to accurately evaluating Wigner 3jm symbols and Gaunt coefficients. Armed with these recursive approaches, we discuss several schemes for the calculation of the vector translation coefficients, using the three general types of formulation, respectively. Our systematic test calculations show that the three types of formulas produce generally the same numerical results except that the algorithm of Stein's type is less accurate in some particular cases. These extensive test calculations also show that the scheme using the formulation based on the Gaunt coefficient is the most efficient in practical computations.},
-  number = {1},
-  journal = {Journal of Computational Physics},
-  doi = {10.1006/jcph.1997.5867},
-  author = {Xu, Yu-lin},
-  month = jan,
-  year = {1998},
-  pages = {137-165},
-  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/STV5263F/Xu - 1998 - Efficient Evaluation of Vector Translation Coeffic.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/VMGZRSAA/S0021999197958678.html}
-}
-
-@book{jackson_classical_1998,
-  address = {{New York}},
-  edition = {3 edition},
-  title = {Classical {{Electrodynamics Third Edition}}},
-  isbn = {978-0-471-30932-1},
-  abstract = {A revision of the defining book covering the physics and classical mathematics necessary to understand electromagnetic fields in materials and at surfaces and interfaces. The third edition has been revised to address the changes in emphasis and applications that have occurred in the past twenty years.},
-  language = {English},
-  publisher = {{Wiley}},
-  author = {Jackson, John David},
-  month = aug,
-  year = {1998},
-  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/3BWPD4BK/John David Jackson-Classical Electrodynamics-Wiley (1999).djvu}
-}
-
-@article{mie_beitrage_1908,
-  title = {Beitr{\"a}ge Zur {{Optik}} Tr{\"u}ber {{Medien}}, Speziell Kolloidaler {{Metall{\"o}sungen}}},
-  volume = {330},
-  copyright = {Copyright \textcopyright{} 1908 WILEY-VCH Verlag GmbH \& Co. KGaA, Weinheim},
-  issn = {1521-3889},
-  language = {en},
-  number = {3},
-  journal = {Ann. Phys.},
-  doi = {10.1002/andp.19083300302},
-  author = {Mie, Gustav},
-  month = jan,
-  year = {1908},
-  pages = {377-445},
-  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/RM9J9RYH/Mie - 1908 - Beiträge zur Optik trüber Medien, speziell kolloid.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/F5A7EX6R/abstract.html}
-}
-
-@book{kristensson_scattering_2016,
-  address = {{Edison, NJ}},
-  title = {Scattering of {{Electromagnetic Waves}} by {{Obstacles}}},
-  isbn = {978-1-61353-221-8},
-  abstract = {This book is an introduction to some of the most important properties of electromagnetic waves and their interaction with passive materials and scatterers. The main purpose of the book is to give a theoretical treatment of these scattering phenomena, and to illustrate numerical computations of some canonical scattering problems for different geometries and materials. The scattering theory is also important in the theory of passive antennas, and this book gives several examples on this topic. Topics covered include an introduction to the basic equations used in scattering; the Green functions and dyadics; integral representation of fields; introductory scattering theory; scattering in the time domain; approximations and applications; spherical vector waves; scattering by spherical objects; the null-field approach; and propagation in stratified media. The book is organised along two tracks, which can be studied separately or together. Track 1 material is appropriate for a first reading of the textbook, while Track 2 contains more advanced material suited for the second reading and for reference. Exercises are included for each chapter.},
-  language = {English},
-  publisher = {{Scitech Publishing}},
-  author = {Kristensson, Gerhard},
-  month = jul,
-  year = {2016},
-  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/3R7VYZUK/Kristensson - 2016 - Scattering of Electromagnetic Waves by Obstacles.pdf}
-}
-
-
diff --git a/lepaper/arrayscat.lyx b/lepaper/arrayscat.lyx
index ab54aa2..437e971 100644
--- a/lepaper/arrayscat.lyx
+++ b/lepaper/arrayscat.lyx
@@ -1,29 +1,40 @@
-#LyX 2.1 created this file. For more info see http://www.lyx.org/
-\lyxformat 474
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 584
 \begin_document
 \begin_header
+\save_transient_properties true
+\origin unavailable
 \textclass article
+\begin_preamble
+\DeclareUnicodeCharacter{0428}{Ш }
+\end_preamble
 \use_default_options true
 \maintain_unincluded_children false
 \language english
 \language_package default
-\inputencoding auto
-\fontencoding global
-\font_roman TeX Gyre Pagella
-\font_sans default
-\font_typewriter default
-\font_math auto
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
 \font_default_family default
-\use_non_tex_fonts true
+\use_non_tex_fonts false
 \font_sc false
-\font_osf true
-\font_sf_scale 100
-\font_tt_scale 100
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
 \graphics default
-\default_output_format pdf4
+\default_output_format default
 \output_sync 0
 \bibtex_command default
 \index_command default
+\float_placement class
+\float_alignment class
 \paperfontsize default
 \spacing single
 \use_hyperref true
@@ -38,7 +49,7 @@
 \pdf_colorlinks false
 \pdf_backref false
 \pdf_pdfusetitle true
-\papersize default
+\papersize a4paper
 \use_geometry false
 \use_package amsmath 2
 \use_package amssymb 1
@@ -59,6 +70,8 @@
 \suppress_date false
 \justification true
 \use_refstyle 1
+\use_minted 0
+\use_lineno 0
 \index Index
 \shortcut idx
 \color #008000
@@ -67,10 +80,14 @@
 \tocdepth 3
 \paragraph_separation indent
 \paragraph_indentation default
-\quotes_language english
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
 \papercolumns 1
 \papersides 1
 \paperpagestyle default
+\tablestyle default
 \tracking_changes false
 \output_changes false
 \html_math_output 0
@@ -117,6 +134,11 @@
 \end_inset
 
 
+\begin_inset FormulaMacro
+\newcommand{\uvec}[1]{\mathbf{\hat{#1}}}
+\end_inset
+
+
 \begin_inset FormulaMacro
 \newcommand{\ud}{\mathrm{d}}
 \end_inset
@@ -127,11 +149,25 @@
 \end_inset
 
 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
 \begin_inset FormulaMacro
 \newcommand{\dc}[1]{Ш_{#1}}
 \end_inset
 
 
+\end_layout
+
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\dc}[1]{\mathrm{III}_{#1}}
+\end_inset
+
+
 \begin_inset FormulaMacro
 \newcommand{\rec}[1]{#1^{-1}}
 \end_inset
@@ -162,6 +198,16 @@
 \end_inset
 
 
+\begin_inset FormulaMacro
+\newcommand{\vsh}[3]{\vect A_{#1,#2,#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\vshD}[3]{\vect A'_{#1,#2,#3}}
+\end_inset
+
+
 \begin_inset FormulaMacro
 \newcommand{\hgfr}{\mathbf{F}}
 \end_inset
@@ -213,7 +259,7 @@
 
 
 \begin_inset FormulaMacro
-\newcommand{\particle}{\mathrm{\Omega}}
+\newcommand{\particle}{\mathrm{\Theta}}
 \end_inset
 
 
@@ -232,15 +278,193 @@
 \end_inset
 
 
+\begin_inset FormulaMacro
+\newcommand{\rcoeffp}[1]{a_{#1}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoeffincp}[1]{a_{#1}^{\mathrm{inc.}}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoeff}{a}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoeffinc}{a^{\mathrm{inc.}}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoeffptlm}[4]{\rcoeffp{#1,#2#3#4}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoefftlm}[3]{\rcoeffp{#1#2#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rcoeffincptlm}[4]{\rcoeffincp{#1,#2#3#4}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\vswfrtlm}[3]{\vect v_{#1#2#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\outcoeff}{f}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\outcoeffp}[1]{f_{#1}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\outcoeffptlm}[4]{\outcoeffp{#1,#2#3#4}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\outcoefftlm}[3]{\outcoeffp{#1#2#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\vswfouttlm}[3]{\vect u_{#1#2#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\Tp}[1]{T_{#1}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\openball}[2]{B_{#1}\left(#2\right)}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\closedball}[2]{B_{#1}#2}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\tropr}{\mathcal{R}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\troprp}[2]{\mathcal{\tropr}_{#1\leftarrow#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\trops}{\mathcal{S}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\tropsp}[2]{\mathcal{\trops}_{#1\leftarrow#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\truncated}[2]{\left[#1\right]_{l\le#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\truncate}[2]{\left[#1\right]_{#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\dlmfFer}[2]{\mathsf{P}_{#1}^{#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\antidelta}{\gamma}
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+These are compatibility macros for the (...)-old files:
+\end_layout
+
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\vswfr}[3]{\vswfrtlm{#3}{#1}{#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\vswfs}[3]{\vswfouttlm{#3}{#1}{#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\svwfs}[3]{\vswfouttlm{#3}{#1}{#2}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\coeffrip}[4]{\rcoeffptlm{#1}{#4}{#2}{#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\coeffsip}[4]{\outcoeffptlm{#1}{#4}{#2}{#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\coeffr}{\rcoeffp{}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\coeffs}{\outcoeffp{}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\transop}{\trops}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\coeffripext}[4]{\rcoeffincptlm{#1}{#4}{#2}{#3}}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\Kp}{K}
+\end_inset
+
+
 \end_layout
 
 \begin_layout Title
-Some nice title about multiple scattering approach to photonic nanoparticle
- arrays (outline)
-\end_layout
-
-\begin_layout Author
-Marek Nečada
+QPMS Photonic Multiple Scattering suite (TODO better title)
 \end_layout
 
 \begin_layout Standard
@@ -254,6 +478,10 @@ The purpose of SIAM Journal on Scientific Computing (SISC) is to advance
 
 \begin_layout Quotation
 SISC papers are classified into three categories:
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -307,12 +535,31 @@ Due to space limitations, articles are normally limited to 20 journal pages.
 Category: Methods and Algorithms for Scientific Computing?
 \end_layout
 
+\begin_layout Abstract
+The (somewhat underrated) T-matrix multiple scattering method (TMMSM) can
+ be used to solve the electromagnetic response of systems consisting of
+ many compact scatterers.
+ It largely surpasses other methods in the number of scatterers it can deal
+ with, while retaining very good accuracy.
+\end_layout
+
+\begin_layout Abstract
+TODO REWRITE: We release a modern implementation of the method under GNU
+ General Public Licence, with several theoretical advancements presented
+ here, such as exploiting the system symmetries to further improve the efficienc
+y of the method, or extending it on infinite periodic systems.
+\end_layout
+
 \begin_layout Section
 Outline
 \end_layout
 
 \begin_layout Itemize
 Intro: 
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -343,6 +590,10 @@ my implementation
 \end_deeper
 \begin_layout Itemize
 Finite systems:
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -353,11 +604,19 @@ motivation (classes of problems that this can solve: response to external
 
 \begin_layout Itemize
 theory
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
 \begin_layout Itemize
 T-matrix definition, basics
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -385,6 +644,10 @@ Example results (or maybe rather in the end)
 \end_deeper
 \begin_layout Itemize
 Infinite lattices:
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -394,6 +657,10 @@ motivation (dispersion relations / modes, ...?)
 
 \begin_layout Itemize
 theory
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
@@ -432,10 +699,83 @@ My implementation.
 Maybe put the numerical results separately in the end.
 \end_layout
 
+\begin_layout Section
+TODOs
+\end_layout
+
+\begin_layout Itemize
+Consistent notation of balls.
+ How is the difference between two cocentric balls called?
+\end_layout
+
+\begin_layout Itemize
+Abstract.
+\end_layout
+
+\begin_layout Itemize
+Translation operators: rewrite in sph.
+ harm.
+ convention independent form.
+\end_layout
+
+\begin_layout Itemize
+Truncation notation.
+\end_layout
+
+\begin_layout Itemize
+Example results!
+\end_layout
+
+\begin_layout Itemize
+Figures.
+\end_layout
+
+\begin_layout Itemize
+Concrete comparison with other methods.
+\end_layout
+
+\begin_layout Itemize
+Fix and unify notation (mainly indices) in infinite lattices section.
+\end_layout
+
+\begin_layout Itemize
+Carefully check the transformation directions in sec.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Symmetries"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Itemize
+The text about symmetries is pretty dense.
+ Make it more explanatory and human-readable.
+\end_layout
+
+\begin_layout Itemize
+Check whether everything written is correct also for non-symmetric space
+ groups.
+\end_layout
+
 \begin_layout Standard
 \begin_inset CommandInset include
 LatexCommand include
 filename "intro.lyx"
+literal "true"
+
+\end_inset
+
+
+\begin_inset CommandInset include
+LatexCommand include
+filename "finite.lyx"
+literal "true"
 
 \end_inset
 
@@ -443,9 +783,19 @@ filename "intro.lyx"
 \end_layout
 
 \begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
 \begin_inset CommandInset include
 LatexCommand include
-filename "finite.lyx"
+filename "finite-old.lyx"
+literal "true"
+
+\end_inset
+
+
+\end_layout
 
 \end_inset
 
@@ -456,6 +806,38 @@ filename "finite.lyx"
 \begin_inset CommandInset include
 LatexCommand include
 filename "infinite.lyx"
+literal "true"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+\begin_inset CommandInset include
+LatexCommand include
+filename "infinite-old.lyx"
+literal "true"
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset CommandInset include
+LatexCommand include
+filename "symmetries.lyx"
+literal "true"
 
 \end_inset
 
@@ -466,6 +848,7 @@ filename "infinite.lyx"
 \begin_inset CommandInset include
 LatexCommand include
 filename "examples.lyx"
+literal "true"
 
 \end_inset
 
@@ -475,8 +858,10 @@ filename "examples.lyx"
 \begin_layout Standard
 \begin_inset CommandInset bibtex
 LatexCommand bibtex
-bibfiles "T-matrix paper"
+btprint "btPrintCited"
+bibfiles "tmpaper"
 options "plain"
+encoding "default"
 
 \end_inset
 
diff --git a/lepaper/comparison.lyx b/lepaper/comparison.lyx
new file mode 100644
index 0000000..6257c19
--- /dev/null
+++ b/lepaper/comparison.lyx
@@ -0,0 +1,109 @@
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
+\begin_document
+\begin_header
+\save_transient_properties true
+\origin unavailable
+\textclass article
+\use_default_options true
+\maintain_unincluded_children false
+\language english
+\language_package default
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
+\font_default_family default
+\use_non_tex_fonts false
+\font_sc false
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
+\graphics default
+\default_output_format default
+\output_sync 0
+\bibtex_command default
+\index_command default
+\float_placement class
+\float_alignment class
+\paperfontsize default
+\spacing single
+\use_hyperref true
+\pdf_title "Sähköpajan päiväkirja"
+\pdf_author "Marek Nečada"
+\pdf_bookmarks true
+\pdf_bookmarksnumbered false
+\pdf_bookmarksopen false
+\pdf_bookmarksopenlevel 1
+\pdf_breaklinks false
+\pdf_pdfborder false
+\pdf_colorlinks false
+\pdf_backref false
+\pdf_pdfusetitle true
+\papersize default
+\use_geometry false
+\use_package amsmath 1
+\use_package amssymb 1
+\use_package cancel 1
+\use_package esint 1
+\use_package mathdots 1
+\use_package mathtools 1
+\use_package mhchem 1
+\use_package stackrel 1
+\use_package stmaryrd 1
+\use_package undertilde 1
+\cite_engine basic
+\cite_engine_type default
+\biblio_style plain
+\use_bibtopic false
+\use_indices false
+\paperorientation portrait
+\suppress_date false
+\justification true
+\use_refstyle 1
+\use_minted 0
+\use_lineno 0
+\index Index
+\shortcut idx
+\color #008000
+\end_index
+\secnumdepth 3
+\tocdepth 3
+\paragraph_separation indent
+\paragraph_indentation default
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
+\papercolumns 1
+\papersides 1
+\paperpagestyle default
+\tablestyle default
+\tracking_changes false
+\output_changes false
+\html_math_output 0
+\html_css_as_file 0
+\html_be_strict false
+\end_header
+
+\begin_body
+
+\begin_layout Section
+Comparison to other methods
+\begin_inset CommandInset label
+LatexCommand label
+name "sec:Comparison"
+
+\end_inset
+
+
+\end_layout
+
+\end_body
+\end_document
diff --git a/lepaper/examples.lyx b/lepaper/examples.lyx
index 1a4f7b2..dee469a 100644
--- a/lepaper/examples.lyx
+++ b/lepaper/examples.lyx
@@ -1,29 +1,37 @@
-#LyX 2.1 created this file. For more info see http://www.lyx.org/
-\lyxformat 474
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
 \begin_document
 \begin_header
+\save_transient_properties true
+\origin unavailable
 \textclass article
 \use_default_options true
 \maintain_unincluded_children false
-\language finnish
+\language english
 \language_package default
-\inputencoding auto
-\fontencoding global
-\font_roman TeX Gyre Pagella
-\font_sans default
-\font_typewriter default
-\font_math auto
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
 \font_default_family default
-\use_non_tex_fonts true
+\use_non_tex_fonts false
 \font_sc false
-\font_osf true
-\font_sf_scale 100
-\font_tt_scale 100
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
 \graphics default
-\default_output_format pdf4
+\default_output_format default
 \output_sync 0
 \bibtex_command default
 \index_command default
+\float_placement class
+\float_alignment class
 \paperfontsize default
 \spacing single
 \use_hyperref true
@@ -58,6 +66,8 @@
 \suppress_date false
 \justification true
 \use_refstyle 1
+\use_minted 0
+\use_lineno 0
 \index Index
 \shortcut idx
 \color #008000
@@ -66,10 +76,14 @@
 \tocdepth 3
 \paragraph_separation indent
 \paragraph_indentation default
-\quotes_language swedish
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
 \papercolumns 1
 \papersides 1
 \paperpagestyle default
+\tablestyle default
 \tracking_changes false
 \output_changes false
 \html_math_output 0
diff --git a/lepaper/finite-old.lyx b/lepaper/finite-old.lyx
new file mode 100644
index 0000000..c6721bb
--- /dev/null
+++ b/lepaper/finite-old.lyx
@@ -0,0 +1,385 @@
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
+\begin_document
+\begin_header
+\save_transient_properties true
+\origin unavailable
+\textclass article
+\use_default_options true
+\maintain_unincluded_children false
+\language english
+\language_package default
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
+\font_default_family default
+\use_non_tex_fonts false
+\font_sc false
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures false
+\graphics default
+\default_output_format default
+\output_sync 0
+\bibtex_command default
+\index_command default
+\float_placement class
+\float_alignment class
+\paperfontsize default
+\spacing single
+\use_hyperref true
+\pdf_title "Sähköpajan päiväkirja"
+\pdf_author "Marek Nečada"
+\pdf_bookmarks true
+\pdf_bookmarksnumbered false
+\pdf_bookmarksopen false
+\pdf_bookmarksopenlevel 1
+\pdf_breaklinks false
+\pdf_pdfborder false
+\pdf_colorlinks false
+\pdf_backref false
+\pdf_pdfusetitle true
+\papersize default
+\use_geometry false
+\use_package amsmath 1
+\use_package amssymb 1
+\use_package cancel 1
+\use_package esint 1
+\use_package mathdots 1
+\use_package mathtools 1
+\use_package mhchem 1
+\use_package stackrel 1
+\use_package stmaryrd 1
+\use_package undertilde 1
+\cite_engine basic
+\cite_engine_type default
+\biblio_style plain
+\use_bibtopic false
+\use_indices false
+\paperorientation portrait
+\suppress_date false
+\justification true
+\use_refstyle 1
+\use_minted 0
+\use_lineno 0
+\index Index
+\shortcut idx
+\color #008000
+\end_index
+\secnumdepth 3
+\tocdepth 3
+\paragraph_separation indent
+\paragraph_indentation default
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
+\papercolumns 1
+\papersides 1
+\paperpagestyle default
+\tablestyle default
+\tracking_changes false
+\output_changes false
+\html_math_output 0
+\html_css_as_file 0
+\html_be_strict false
+\end_header
+
+\begin_body
+
+\begin_layout Subsection
+The multiple-scattering problem
+\begin_inset CommandInset label
+LatexCommand label
+name "subsec:The-multiple-scattering-problem"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+In the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix approach, scattering properties of single nanoparticles in a homogeneous
+ medium are first computed in terms of vector sperical wavefunctions (VSWFs)—the
+ field incident onto the 
+\begin_inset Formula $n$
+\end_inset
+
+-th nanoparticle from external sources can be expanded as 
+\begin_inset Formula 
+\begin{equation}
+\vect E_{n}^{\mathrm{inc}}(\vect r)=\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\sum_{t=\mathrm{E},\mathrm{M}}\coeffrip nlmt\vswfr lmt\left(\vect r_{n}\right)\label{eq:E_inc}
+\end{equation}
+
+\end_inset
+
+where 
+\begin_inset Formula $\vect r_{n}=\vect r-\vect R_{n}$
+\end_inset
+
+, 
+\begin_inset Formula $\vect R_{n}$
+\end_inset
+
+ being the position of the centre of 
+\begin_inset Formula $n$
+\end_inset
+
+-th nanoparticle and 
+\begin_inset Formula $\vswfr lmt$
+\end_inset
+
+ are the regular VSWFs which can be expressed in terms of regular spherical
+ Bessel functions of 
+\begin_inset Formula $j_{k}\left(\left|\vect r_{n}\right|\right)$
+\end_inset
+
+ and spherical harmonics 
+\begin_inset Formula $\ush km\left(\hat{\vect r}_{n}\right)$
+\end_inset
+
+; the expressions, together with a proof that the VSWFs span all the solutions
+ of vector Helmholtz equation around the particle, justifying the expansion,
+ can be found e.g.
+ in 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "chapter 7"
+key "kristensson_scattering_2016"
+literal "true"
+
+\end_inset
+
+ (care must be taken because of varying normalisation and phase conventions).
+ On the other hand, the field scattered by the particle can be (outside
+ the particle's circumscribing sphere) expanded in terms of singular VSWFs
+ 
+\begin_inset Formula $\vswfs lmt$
+\end_inset
+
+ which differ from the regular ones by regular spherical Bessel functions
+ being replaced with spherical Hankel functions 
+\begin_inset Formula $h_{k}^{(1)}\left(\left|\vect r_{n}\right|\right)$
+\end_inset
+
+, 
+\begin_inset Formula 
+\begin{equation}
+\vect E_{n}^{\mathrm{scat}}\left(\vect r\right)=\sum_{l,m,t}\coeffsip nlmt\vswfs lmt\left(\vect r_{n}\right).\label{eq:E_scat}
+\end{equation}
+
+\end_inset
+
+The expansion coefficients 
+\begin_inset Formula $\coeffsip nlmt$
+\end_inset
+
+, 
+\begin_inset Formula $t=\mathrm{E},\mathrm{M}$
+\end_inset
+
+ are related to the electric and magnetic multipole polarization amplitudes
+ of the nanoparticle.
+\end_layout
+
+\begin_layout Standard
+At a given frequency, assuming the system is linear, the relation between
+ the expansion coefficients in the VSWF bases is given by the so-called
+ 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix, 
+\begin_inset Formula 
+\begin{equation}
+\coeffsip nlmt=\sum_{l',m',t'}T_{n}^{lmt;l'm't'}\coeffrip n{l'}{m'}{t'}.\label{eq:Tmatrix definition}
+\end{equation}
+
+\end_inset
+
+The 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix is given by the shape and composition of the particle and fully
+ describes its scattering properties.
+ In theory it is infinite-dimensional, but in practice (at least for subwaveleng
+th nanoparticles) its elements drop very quickly to negligible values with
+ growing degree indices 
+\begin_inset Formula $l,l'$
+\end_inset
+
+, enabling to take into account only the elements up to some finite degree,
+ 
+\begin_inset Formula $l,l'\le l_{\mathrm{max}}$
+\end_inset
+
+.
+ The 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix can be calculated numerically using various methods; here we used
+ the scuff-tmatrix tool from the SCUFF-EM suite 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "SCUFF2,reid_efficient_2015"
+literal "true"
+
+\end_inset
+
+, which implements the boundary element method (BEM).
+\end_layout
+
+\begin_layout Standard
+The singular VSWFs originating at 
+\begin_inset Formula $\vect R_{n}$
+\end_inset
+
+ can be then re-expanded around another origin (nanoparticle location) 
+\begin_inset Formula $\vect R_{n'}$
+\end_inset
+
+ in terms of regular VSWFs, 
+\begin_inset Formula 
+\begin{equation}
+\begin{split}\svwfs lmt\left(\vect r_{n}\right)=\sum_{l',m',t'}\transop^{l'm't';lmt}\left(\vect R_{n'}-\vect R_{n}\right)\vswfr{l'}{m'}{t'}\left(\vect r_{n'}\right),\\
+\left|\vect r_{n'}\right|<\left|\vect R_{n'}-\vect R_{n}\right|.
+\end{split}
+\label{eq:translation op def}
+\end{equation}
+
+\end_inset
+
+Analytical expressions for the translation operator 
+\begin_inset Formula $\transop^{lmt;l'm't'}\left(\vect R_{n'}-\vect R_{n}\right)$
+\end_inset
+
+ can be found in 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "xu_efficient_1998"
+literal "true"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+If we write the field incident onto the 
+\begin_inset Formula $n$
+\end_inset
+
+-th nanoparticle as the sum of fields scattered from all the other nanoparticles
+ and an external field 
+\begin_inset Formula $\vect E_{0}$
+\end_inset
+
+ (which we also expand around each nanoparticle, 
+\begin_inset Formula $\vect E_{0}\left(\vect r\right)=\sum_{l,m,t}\coeffripext nlmt\vswfr lmt\left(\vect r_{n}\right)$
+\end_inset
+
+), 
+\begin_inset Formula 
+\[
+\vect E_{n}^{\mathrm{inc}}\left(\vect r\right)=\vect E_{0}\left(\vect r\right)+\sum_{n'\ne n}\vect E_{n'}^{\mathrm{scat}}\left(\vect r\right)
+\]
+
+\end_inset
+
+and use eqs.
+ (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:E_inc"
+
+\end_inset
+
+)–(
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:translation op def"
+
+\end_inset
+
+), we obtain a set of linear equations for the electromagnetic response
+ (multiple scattering) of the whole set of nanoparticles, 
+\begin_inset Formula 
+\begin{equation}
+\begin{split}\coeffrip nlmt=\coeffripext nlmt+\sum_{n'\ne n}\sum_{l',m',t'}\transop^{lmt;l'm't'}\left(\vect R_{n}-\vect R_{n'}\right)\\
+\times\sum_{l'',m'',t''}T_{n'}^{l'm't';l''m''t''}\coeffrip{n'}{l''}{m''}{t''}.
+\end{split}
+\label{eq:multiplescattering element-wise}
+\end{equation}
+
+\end_inset
+
+It is practical to get rid of the VSWF indices, rewriting (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:multiplescattering element-wise"
+
+\end_inset
+
+) in a per-particle matrix form 
+\begin_inset Formula 
+\begin{equation}
+\coeffr_{n}=\coeffr_{\mathrm{ext}(n)}+\sum_{n'\ne n}S_{n,n'}T_{n'}p_{n'}\label{eq:multiple scattering per particle p}
+\end{equation}
+
+\end_inset
+
+and to reformulate the problem using (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:Tmatrix definition"
+
+\end_inset
+
+) in terms of the 
+\begin_inset Formula $\coeffs$
+\end_inset
+
+-coefficients which describe the multipole excitations of the particles
+ 
+\begin_inset Formula 
+\begin{equation}
+\coeffs_{n}-T_{n}\sum_{n'\ne n}S_{n,n'}\coeffs_{n'}=T_{n}\coeffr_{\mathrm{ext}(n)}.\label{eq:multiple scattering per particle a}
+\end{equation}
+
+\end_inset
+
+Knowing 
+\begin_inset Formula $T_{n},S_{n,n'},\coeffr_{\mathrm{ext}(n)}$
+\end_inset
+
+, the nanoparticle excitations 
+\begin_inset Formula $a_{n}$
+\end_inset
+
+ can be solved by standard linear algebra methods.
+ The total scattered field anywhere outside the particles' circumscribing
+ spheres is then obtained by summing the contributions (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:E_scat"
+
+\end_inset
+
+) from all particles.
+\end_layout
+
+\end_body
+\end_document
diff --git a/lepaper/finite.lyx b/lepaper/finite.lyx
index 4084ce9..1f2d76c 100644
--- a/lepaper/finite.lyx
+++ b/lepaper/finite.lyx
@@ -1,29 +1,37 @@
-#LyX 2.1 created this file. For more info see http://www.lyx.org/
-\lyxformat 474
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 584
 \begin_document
 \begin_header
+\save_transient_properties true
+\origin unavailable
 \textclass article
 \use_default_options true
 \maintain_unincluded_children false
-\language finnish
+\language english
 \language_package default
-\inputencoding auto
-\fontencoding global
-\font_roman TeX Gyre Pagella
-\font_sans default
-\font_typewriter default
-\font_math auto
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
 \font_default_family default
-\use_non_tex_fonts true
+\use_non_tex_fonts false
 \font_sc false
-\font_osf true
-\font_sf_scale 100
-\font_tt_scale 100
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
 \graphics default
-\default_output_format pdf4
+\default_output_format default
 \output_sync 0
 \bibtex_command default
 \index_command default
+\float_placement class
+\float_alignment class
 \paperfontsize default
 \spacing single
 \use_hyperref true
@@ -58,6 +66,8 @@
 \suppress_date false
 \justification true
 \use_refstyle 1
+\use_minted 0
+\use_lineno 0
 \index Index
 \shortcut idx
 \color #008000
@@ -66,10 +76,14 @@
 \tocdepth 3
 \paragraph_separation indent
 \paragraph_indentation default
-\quotes_language swedish
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
 \papercolumns 1
 \papersides 1
 \paperpagestyle default
+\tablestyle default
 \tracking_changes false
 \output_changes false
 \html_math_output 0
@@ -84,45 +98,45 @@ Finite systems
 \end_layout
 
 \begin_layout Itemize
-
-\lang english
 motivation (classes of problems that this can solve: response to external
  radiation, resonances, ...)
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
 \begin_layout Itemize
-
-\lang english
 theory
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
 \begin_layout Itemize
-
-\lang english
 T-matrix definition, basics
+\begin_inset Separator latexpar
+\end_inset
+
+
 \end_layout
 
 \begin_deeper
 \begin_layout Itemize
-
-\lang english
 How to get it?
 \end_layout
 
 \end_deeper
 \begin_layout Itemize
-
-\lang english
 translation operators (TODO think about how explicit this should be, but
  I guess it might be useful to write them to write them explicitly (but
  in the shortest possible form) in the normalisation used in my program)
 \end_layout
 
 \begin_layout Itemize
-
-\lang english
 employing point group symmetries and decomposing the problem to decrease
  the computational complexity (maybe separately)
 \end_layout
@@ -130,14 +144,33 @@ employing point group symmetries and decomposing the problem to decrease
 \end_deeper
 \end_deeper
 \begin_layout Subsection
+Motivation/intro
+\end_layout
 
-\lang english
-Motivation
+\begin_layout Standard
+The basic idea of MSTMM is quite simple: the driving electromagnetic field
+ incident onto a scatterer is expanded into a vector spherical wavefunction
+ (VSWF) basis in which the single scattering problem is solved, and the
+ scattered field is then re-expanded into VSWFs centered at the other scatterers.
+ Repeating the same procedure with all (pairs of) scatterers yields a set
+ of linear equations, solution of which gives the coefficients of the scattered
+ field in the VSWF bases.
+ Once these coefficients have been found, one can evaluate various quantities
+ related to the scattering (such as cross sections or the scattered fields)
+ quite easily.
+ 
+\end_layout
+
+\begin_layout Standard
+The expressions appearing in the re-expansions are fairly complicated, and
+ the implementation of MSTMM is extremely error-prone also due to the various
+ conventions used in the literature.
+ Therefore although we do not re-derive from scratch the expressions that
+ can be found elsewhere in literature, we always state them explicitly in
+ our convention.
 \end_layout
 
 \begin_layout Subsection
-
-\lang english
 Single-particle scattering
 \end_layout
 
@@ -164,12 +197,40 @@ ity
 , and that the whole system is linear, i.e.
  the material properties of neither the medium nor the scatterer depend
  on field intensities.
- Under these assumptions, the EM fields in 
+ Under these assumptions, the EM fields 
+\begin_inset Formula $\vect{\Psi}=\vect E,\vect H$
+\end_inset
+
+ in 
 \begin_inset Formula $\medium$
 \end_inset
 
- must satisfy the homogeneous vector Helmholtz equation 
-\begin_inset Formula $\left(\nabla^{2}+k^{2}\right)\Psi\left(\vect r,\vect{\omega}\right)=0$
+ must satisfy the homogeneous vector Helmholtz equation together with the
+ transversality condition 
+\begin_inset Formula 
+\begin{equation}
+\left(\nabla^{2}+k^{2}\right)\Psi\left(\vect r,\vect{\omega}\right)=0,\quad\nabla\cdot\vect{\Psi}\left(\vect r,\vect{\omega}\right)=0\label{eq:Helmholtz eq}
+\end{equation}
+
+\end_inset
+
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+frequency-space Maxwell's equations
+\begin_inset Formula 
+\begin{align*}
+\nabla\times\vect E\left(\vect r,\omega\right)-ik\eta_{0}\eta\vect H\left(\vect r,\omega\right) & =0,\\
+\eta_{0}\eta\nabla\times\vect H\left(\vect r,\omega\right)+ik\vect E\left(\vect r,\omega\right) & =0.
+\end{align*}
+
+\end_inset
+
+
+\end_layout
+
 \end_inset
 
  
@@ -187,43 +248,197 @@ todo define
 \end_inset
 
  with 
-\begin_inset Formula $k=TODO$
+\begin_inset Formula $k=k\left(\omega\right)=\omega\sqrt{\mubg(\omega)\epsbg(\omega)}/c_{0}$
 \end_inset
 
- [TODO REF Jackson?].
- Its solutions (TODO under which conditions? What vector space do the SVWFs
- actually span? Check Comment 9.2 and Appendix f.9.1 in Kristensson)
+, as can be derived from the Maxwell's equations 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "???"
+key "jackson_classical_1998"
+literal "false"
+
+\end_inset
+
+.
+ 
+\end_layout
+
+\begin_layout Subsubsection
+Spherical waves
 \end_layout
 
 \begin_layout Standard
+Equation 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:Helmholtz eq"
+plural "false"
+caps "false"
+noprefix "false"
 
-\lang english
-Throughout this text, we will use the same normalisation conventions as
- in 
+\end_inset
+
+ can be solved by separation of variables in spherical coordinates to give
+ the solutions – the 
+\emph on
+regular
+\emph default
+ and 
+\emph on
+outgoing
+\emph default
+ vector spherical wavefunctions (VSWFs) 
+\begin_inset Formula $\vswfrtlm{\tau}lm\left(k\vect r\right)$
+\end_inset
+
+ and 
+\begin_inset Formula $\vswfouttlm{\tau}lm\left(k\vect r\right)$
+\end_inset
+
+, respectively, defined as follows:
+\begin_inset Formula 
+\begin{align}
+\vswfrtlm 1lm\left(k\vect r\right) & =j_{l}\left(kr\right)\vsh 1lm\left(\uvec r\right),\nonumber \\
+\vswfrtlm 2lm\left(k\vect r\right) & =\frac{1}{kr}\frac{\ud\left(krj_{l}\left(kr\right)\right)}{\ud\left(kr\right)}\vsh 2lm\left(\uvec r\right)+\sqrt{l\left(l+1\right)}\frac{j_{l}\left(kr\right)}{kr}\vsh 3lm\left(\uvec r\right),\label{eq:VSWF regular}
+\end{align}
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{align}
+\vswfouttlm 1lm\left(k\vect r\right) & =h_{l}^{\left(1\right)}\left(kr\right)\vsh 1lm\left(\uvec r\right),\nonumber \\
+\vswfouttlm 2lm\left(k\vect r\right) & =\frac{1}{kr}\frac{\ud\left(krh_{l}^{\left(1\right)}\left(kr\right)\right)}{\ud\left(kr\right)}\vsh 2lm\left(\uvec r\right)+\sqrt{l\left(l+1\right)}\frac{h_{l}^{\left(1\right)}\left(kr\right)}{kr}\vsh 3lm\left(\uvec r\right),\label{eq:VSWF outgoing}\\
+ & \tau=1,2;\quad l=1,2,3,\dots;\quad m=-l,-l+1,\dots,+l,\nonumber 
+\end{align}
+
+\end_inset
+
+where 
+\begin_inset Formula $\vect r=r\uvec r$
+\end_inset
+
+, 
+\begin_inset Formula $j_{l}\left(x\right),h_{l}^{\left(1\right)}\left(x\right)=j_{l}\left(x\right)+iy_{l}\left(x\right)$
+\end_inset
+
+ are the regular spherical Bessel function and spherical Hankel function
+ of the first kind, respectively, as in [DLMF §10.47], and 
+\begin_inset Formula $\vsh{\tau}lm$
+\end_inset
+
+ are the 
+\emph on
+vector spherical harmonics
+\emph default
+
+\begin_inset Formula 
+\begin{align}
+\vsh 1lm\left(\uvec r\right) & =\frac{1}{\sqrt{l\left(l+1\right)}}\nabla\times\left(\vect r\ush lm\left(\uvec r\right)\right)=\frac{1}{\sqrt{l\left(l+1\right)}}\nabla\ush lm\left(\uvec r\right)\times\vect r,\nonumber \\
+\vsh 2lm\left(\uvec r\right) & =\frac{1}{\sqrt{l\left(l+1\right)}}r\nabla\ush lm\left(\uvec r\right),\nonumber \\
+\vsh 3lm\left(\uvec r\right) & =\uvec r\ush lm\left(\uvec r\right).\label{eq:vector spherical harmonics definition}
+\end{align}
+
+\end_inset
+
+In our convention, the (scalar) spherical harmonics 
+\begin_inset Formula $\ush lm$
+\end_inset
+
+ are identical to those in [DLMF 14.30.1], i.e.
+\begin_inset Formula 
+\[
+\ush lm=\left(\frac{\left(l-m\right)!\left(2l+1\right)}{4\pi\left(l+m\right)!}\right)^{\frac{1}{2}}e^{im\phi}\dlmfFer lm\left(\cos\theta\right)
+\]
+
+\end_inset
+
+where importantly, the Ferrers functions 
+\begin_inset Formula $\dlmfFer lm$
+\end_inset
+
+ defined as in [DLMF §14.3(i)] do already contain the Condon-Shortley phase
+ 
+\begin_inset Formula $\left(-1\right)^{m}$
+\end_inset
+
+.
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO názornější definice.
+\end_layout
+
+\end_inset
+
+ 
+\end_layout
+
+\begin_layout Standard
+The convention for VSWFs used here is the same as in 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "kristensson_spherical_2014"
+literal "false"
+
+\end_inset
+
+; over other conventions used elsewhere in literature, it has several fundamenta
+l advantages – most importantly, the translation operators introduced later
+ in eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:translation op def"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ are unitary, and it gives the simplest possible expressions for power transport
+ and cross sections without additional 
+\begin_inset Formula $l,m$
+\end_inset
+
+-dependent factors (for that reason, we also call our VSWFs as 
+\emph on
+power-normalised
+\emph default
+).
+ Power-normalisation and unitary translation operators are possible to achieve
+ also with real spherical harmonics – such a convention is used in 
 \begin_inset CommandInset citation
 LatexCommand cite
 key "kristensson_scattering_2016"
+literal "false"
 
 \end_inset
 
 .
 \end_layout
 
-\begin_layout Subsubsection
-
-\lang english
-Spherical waves
-\end_layout
-
 \begin_layout Standard
-
-\lang english
 \begin_inset Note Note
 status open
 
 \begin_layout Plain Layout
+Its solutions (TODO under which conditions? What vector space do the SVWFs
+ actually span? Check Comment 9.2 and Appendix f.9.1 in Kristensson)
+\end_layout
 
-\lang english
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
 TODO small note about cartesian multipoles, anapoles etc.
  (There should be some comparing paper that the Russians at META 2018 mentioned.)
 \end_layout
@@ -234,37 +449,1596 @@ TODO small note about cartesian multipoles, anapoles etc.
 \end_layout
 
 \begin_layout Subsubsection
-
-\lang english
 T-matrix definition
 \end_layout
 
-\begin_layout Subsubsection
-Absorbed power
+\begin_layout Standard
+The regular VSWFs 
+\begin_inset Formula $\vswfrtlm{\tau}lm\left(k\vect r\right)$
+\end_inset
+
+ constitute a basis for solutions of the Helmholtz equation 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Helmholtz eq"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ inside a ball 
+\begin_inset Formula $\openball 0{R^{>}}$
+\end_inset
+
+ with radius 
+\begin_inset Formula $R^{>}$
+\end_inset
+
+ and center in the origin; however, if the equation is not guaranteed to
+ hold inside a smaller ball 
+\begin_inset Formula $B_{0}\left(R\right)$
+\end_inset
+
+ around the origin (typically due to presence of a scatterer), one has to
+ add the outgoing VSWFs 
+\begin_inset Formula $\vswfrtlm{\tau}lm\left(k\vect r\right)$
+\end_inset
+
+ to have a complete basis of the solutions in the volume 
+\begin_inset Formula $\openball 0{R^{>}}\backslash B_{0}\left(R\right)$
+\end_inset
+
+.
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+Vnitřní koule uzavřená? Jak se řekne mezikulí anglicky?
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+The single-particle scattering problem at frequency 
+\begin_inset Formula $\omega$
+\end_inset
+
+ can be posed as follows: Let a scatterer be enclosed inside the ball 
+\begin_inset Formula $B_{0}\left(R\right)$
+\end_inset
+
+ and let the whole volume 
+\begin_inset Formula $\openball 0{R^{>}}\backslash B_{0}\left(R\right)$
+\end_inset
+
+ be filled with a homogeneous isotropic medium with wave number 
+\begin_inset Formula $k\left(\omega\right)$
+\end_inset
+
+.
+ Inside this volume, the electric field can be expanded as
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+doplnit frekvence a polohy
+\end_layout
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{equation}
+\vect E\left(\omega,\vect r\right)=\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\left(\rcoefftlm{\tau}lm\vswfrtlm{\tau}lm\left(k\vect r\right)+\outcoefftlm{\tau}lm\vswfouttlm{\tau}lm\left(k\vect r\right)\right).\label{eq:E field expansion}
+\end{equation}
+
+\end_inset
+
+If there was no scatterer and 
+\begin_inset Formula $B_{0}\left(R_{<}\right)$
+\end_inset
+
+ was filled with the same homogeneous medium, the part with the outgoing
+ VSWFs would vanish and only the part 
+\begin_inset Formula $\vect E_{\mathrm{inc}}=\sum_{\tau lm}\rcoefftlm{\tau}lm\vswfrtlm{\tau}lm$
+\end_inset
+
+ due to sources outside 
+\begin_inset Formula $\openball 0R$
+\end_inset
+
+ would remain.
+ Let us assume that the 
+\begin_inset Quotes eld
+\end_inset
+
+driving field
+\begin_inset Quotes erd
+\end_inset
+
+ is given, so that presence of the scatterer does not affect 
+\begin_inset Formula $\vect E_{\mathrm{inc}}$
+\end_inset
+
+ and is fully manifested in the latter part, 
+\begin_inset Formula $\vect E_{\mathrm{scat}}=\sum_{\tau lm}\outcoefftlm{\tau}lm\vswfouttlm{\tau}lm$
+\end_inset
+
+.
+ We also assume that the scatterer is made of optically linear materials,
+ and hence reacts on the incident field in a linear manner.
+ This gives a linearity constraint between the expansion coefficients
+\begin_inset Formula 
+\begin{equation}
+\outcoefftlm{\tau}lm=\sum_{\tau'l'm'}T_{\tau lm}^{\tau'l'm'}\rcoefftlm{\tau'}{l'}{m'}\label{eq:T-matrix definition}
+\end{equation}
+
+\end_inset
+
+where the 
+\begin_inset Formula $T_{\tau lm}^{\tau'l'm'}=T_{\tau lm}^{\tau'l'm'}\left(\omega\right)$
+\end_inset
+
+ are the elements of the 
+\emph on
+transition matrix,
+\emph default
+ a.k.a.
+ 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix.
+ It completely describes the scattering properties of a linear scatterer,
+ so with the knowledge of the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix, we can solve the single-patricle scatering prroblem simply by substitut
+ing appropriate expansion coefficients 
+\begin_inset Formula $\rcoefftlm{\tau'}{l'}{m'}$
+\end_inset
+
+ of the driving field into 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:T-matrix definition"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+ The outgoing VSWF expansion coefficients 
+\begin_inset Formula $\outcoefftlm{\tau}lm$
+\end_inset
+
+ are the effective induced electric (
+\begin_inset Formula $\tau=2$
+\end_inset
+
+) and magnetic (
+\begin_inset Formula $\tau=1$
+\end_inset
+
+) multipole polarisation amplitudes of the scatterer, and this is why we
+ sometimes refer to the corresponding VSWFs as the electric and magnetic
+ VSWFs.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO mention the pseudovector character of magnetic VSWFs.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TOOD H-field expansion here?
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula $T$
+\end_inset
+
+-matrices of particles with certain simple geometries (most famously spherical)
+ can be obtained analytically [Kristensson 2016, Mie], but in general one
+ can find them numerically by simulating scattering of a regular spherical
+ wave 
+\begin_inset Formula $\vswfouttlm{\tau}lm$
+\end_inset
+
+ and projecting the scattered fields (or induced currents, depending on
+ the method) onto the outgoing VSWFs 
+\begin_inset Formula $\vswfrtlm{\tau}{'l'}{m'}$
+\end_inset
+
+.
+ In practice, one can compute only a finite number of elements with a cut-off
+ value 
+\begin_inset Formula $L$
+\end_inset
+
+ on the multipole degree, 
+\begin_inset Formula $l,l'\le L$
+\end_inset
+
+, see below.
+ We typically use the scuff-tmatrix tool from the free software SCUFF-EM
+ suite 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "reid_efficient_2015,SCUFF2"
+literal "false"
+
+\end_inset
+
+.
+ Note that older versions of SCUFF-EM contained a bug that rendered almost
+ all 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix results wrong; we found and fixed the bug and from upstream version
+ xxx onwards, it should behave correctly.
+ 
 \end_layout
 
 \begin_layout Subsubsection
-
-\lang english
 T-matrix compactness, cutoff validity
 \end_layout
 
+\begin_layout Standard
+The magnitude of the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix elements depends heavily on the scatterer's size compared to the
+ wavelength.
+ Fortunately, the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix of a bounded scatterer is a compact operator [REF???], so from certain
+ multipole degree onwards, 
+\begin_inset Formula $l,l'>L$
+\end_inset
+
+, the elements of the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix are negligible, so truncating the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix at finite multipole degree 
+\begin_inset Formula $L$
+\end_inset
+
+ gives a good approximation of the actual infinite-dimensional itself.
+ If the incident field is well-behaved, i.e.
+ the expansion coefficients 
+\begin_inset Formula $\rcoefftlm{\tau'}{l'}{m'}$
+\end_inset
+
+ do not take excessive values for 
+\begin_inset Formula $l'>L$
+\end_inset
+
+, the scattered field expansion coefficients 
+\begin_inset Formula $\outcoefftlm{\tau}lm$
+\end_inset
+
+ with 
+\begin_inset Formula $l>L$
+\end_inset
+
+ will also be negligible.
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO when it will not be negligible
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+A rule of thumb to choose the 
+\begin_inset Formula $L$
+\end_inset
+
+ with desired 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix element accuracy 
+\begin_inset Formula $\delta$
+\end_inset
+
+ can be obtained from the spherical Bessel function expansion around zero,
+ TODO.
+ 
+\end_layout
+
+\begin_layout Subsubsection
+Power transport
+\end_layout
+
+\begin_layout Standard
+For convenience, let us introduce a short-hand matrix notation for the expansion
+ coefficients and related quantities, so that we do not need to write the
+ indices explicitly; so for example, eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:T-matrix definition"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ would be written as 
+\begin_inset Formula $\outcoeffp{}=T\rcoeffp{}$
+\end_inset
+
+, where 
+\begin_inset Formula $\rcoeffp{},\outcoeffp{}$
+\end_inset
+
+ are column vectors with the expansion coefficients.
+ Transposed and complex-conjugated matrices are labeled with the 
+\begin_inset Formula $\dagger$
+\end_inset
+
+ superscript.
+\end_layout
+
+\begin_layout Standard
+With this notation, we state an important result about power transport,
+ derivation of which can be found in 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "sect. 7.3"
+key "kristensson_scattering_2016"
+literal "true"
+
+\end_inset
+
+.
+ Let the field in 
+\begin_inset Formula $\openball 0{R^{>}}\backslash B_{0}\left(R\right)$
+\end_inset
+
+ have expansion as in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:E field expansion"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+ Then the net power transported from 
+\begin_inset Formula $B_{0}\left(R\right)$
+\end_inset
+
+ to 
+\begin_inset Formula $\openball 0{R^{>}}\backslash B_{0}\left(R\right)$
+\end_inset
+
+ via by electromagnetic radiation is
+\begin_inset Formula 
+\begin{equation}
+P=\frac{1}{2k^{2}\eta_{0}\eta}\left(\Re\left(\rcoeffp{}^{\dagger}\outcoeffp{}\right)+\left\Vert \outcoeffp{}\right\Vert ^{2}\right)=\frac{1}{2k^{2}\eta_{0}\eta}\rcoeffp{}^{\dagger}\left(\Tp{}^{\dagger}\Tp{}+\frac{\Tp{}^{\dagger}+\Tp{}}{2}\right)\rcoeffp{}.\label{eq:Power transport}
+\end{equation}
+
+\end_inset
+
+In realistic scattering setups, power is transferred by radiation into 
+\begin_inset Formula $B_{0}\left(R\right)$
+\end_inset
+
+ and absorbed by the enclosed scatterer, so 
+\begin_inset Formula $P$
+\end_inset
+
+ is negative and its magnitude equals to power absorbed by the scatterer.
+\end_layout
+
+\begin_layout Subsubsection
+Plane wave expansion
+\end_layout
+
+\begin_layout Standard
+In many scattering problems considered in practice, the driving field is
+ a plane wave.
+ A transversal (
+\begin_inset Formula $\vect k\cdot\vect E_{0}=0$
+\end_inset
+
+) plane wave propagating in direction 
+\begin_inset Formula $\uvec k$
+\end_inset
+
+ with (complex) amplitude 
+\begin_inset Formula $\vect E_{0}$
+\end_inset
+
+ can be expanded into regular VSWFs 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "7.???"
+key "kristensson_scattering_2016"
+literal "false"
+
+\end_inset
+
+ as
+\begin_inset Formula 
+\[
+\vect E_{\mathrm{PW}}\left(\vect r,\omega\right)=\vect E_{0}e^{ik\uvec k\cdot\vect r}=\sum_{\tau,l,m}\rcoeffptlm{}{\tau}lm\left(\vect k,\vect E_{0}\right)\vswfrtlm{\tau}lm\left(k\vect r\right),
+\]
+
+\end_inset
+
+with expansion coefficients
+\begin_inset Formula 
+\begin{eqnarray}
+\rcoeffptlm{}1lm\left(\vect k,\vect E_{0}\right) & = & 4\pi i^{l}\vshD 1lm\left(\uvec k\right),\nonumber \\
+\rcoeffptlm{}2lm\left(\vect k,\vect E_{0}\right) & = & -4\pi i^{l+1}\vshD 2lm\left(\uvec k\right).\label{eq:plane wave expansion}
+\end{eqnarray}
+
+\end_inset
+
+where 
+\begin_inset Formula $\vshD{\tau}lm$
+\end_inset
+
+ are the 
+\begin_inset Quotes eld
+\end_inset
+
+dual
+\begin_inset Quotes erd
+\end_inset
+
+ vector spherical harmonics defined by duality relation with the 
+\begin_inset Quotes eld
+\end_inset
+
+usual
+\begin_inset Quotes erd
+\end_inset
+
+ vector spherical harmonics 
+\begin_inset Formula 
+\begin{equation}
+\iint\vshD{\tau'}{l'}{m'}\left(\uvec r\right)\cdot\vsh{\tau}lm\left(\uvec r\right)\,\ud\Omega=\delta_{\tau'\tau}\delta_{l'l}\delta_{m'm}\label{eq:dual vsh}
+\end{equation}
+
+\end_inset
+
+(complex conjugation not implied in the dot product here).
+ In our convention, we have
+\begin_inset Formula 
+\[
+\vshD{\tau}lm\left(\uvec r\right)=\left(\vsh{\tau}lm\left(\uvec r\right)\right)^{*}=\left(-1\right)^{m}\vsh{\tau}{l-}m\left(\uvec r\right).
+\]
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Cross-sections (single-particle)
+\end_layout
+
+\begin_layout Standard
+With the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix and expansion coefficients of plane waves in hand, we can state
+ the expressions for cross-sections of a single scatterer.
+ Assuming a non-lossy background medium, extinction, scattering and absorption
+ cross sections of a single scatterer irradiated by a plane wave propagating
+ in direction 
+\begin_inset Formula $\uvec k$
+\end_inset
+
+ and (complex) amplitude 
+\begin_inset Formula $\vect E_{0}$
+\end_inset
+
+ are 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "sect. 7.8.2"
+key "kristensson_scattering_2016"
+literal "true"
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{eqnarray}
+\sigma_{\mathrm{ext}}\left(\uvec k\right) & = & -\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\Re\left(\rcoeffp{}^{\dagger}\outcoeffp{}\right)=-\frac{1}{2k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\rcoeffp{}^{\dagger}\left(\Tp{}+\Tp{}^{\dagger}\right)\rcoeffp{},\label{eq:extincion CS single}\\
+\sigma_{\mathrm{scat}}\left(\uvec k\right) & = & \frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\left\Vert \outcoeffp{}\right\Vert ^{2}=\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\rcoeffp{}^{\dagger}\left(\Tp{}^{\dagger}\Tp{}\right)\rcoeffp{},\label{eq:scattering CS single}\\
+\sigma_{\mathrm{abs}}\left(\uvec k\right) & = & \sigma_{\mathrm{ext}}\left(\uvec k\right)-\sigma_{\mathrm{scat}}\left(\uvec k\right)=-\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\left(\Re\left(\rcoeffp{}^{\dagger}\outcoeffp{}\right)+\left\Vert \outcoeffp{}\right\Vert ^{2}\right)\nonumber \\
+ &  & =\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\rcoeffp{}^{\dagger}\left(\Tp{}^{\dagger}\Tp{}+\frac{\Tp{}^{\dagger}+\Tp{}}{2}\right)\rcoeffp{},\label{eq:absorption CS single}
+\end{eqnarray}
+
+\end_inset
+
+where 
+\begin_inset Formula $\rcoeffp{}=\rcoeffp{}\left(\vect k,\vect E_{0}\right)$
+\end_inset
+
+ is the vector of plane wave expansion coefficients as in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:plane wave expansion"
+
+\end_inset
+
+.
+ 
+\end_layout
+
 \begin_layout Subsection
-
-\lang english
 Multiple scattering
+\begin_inset CommandInset label
+LatexCommand label
+name "subsec:Multiple-scattering"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+If the system consists of multiple scatterers, the EM fields around each
+ one can be expanded in analogous way.
+ Let 
+\begin_inset Formula $\mathcal{P}$
+\end_inset
+
+ be an index set labeling the scatterers.
+ We enclose each scatterer in a ball 
+\begin_inset Formula $B_{\vect r_{p}}\left(R_{p}\right)$
+\end_inset
+
+ such that the balls do not touch, 
+\begin_inset Formula $B_{\vect r_{p}}\left(R_{p}\right)\cap B_{\vect r_{q}}\left(R_{q}\right)=\emptyset;p,q\in\mathcal{P}$
+\end_inset
+
+, 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO bacha, musejí být uzavřené!
+\end_layout
+
+\end_inset
+
+so there is a non-empty volume
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+jaksetometuje?
+\end_layout
+
+\end_inset
+
+ 
+\begin_inset Formula $\openball{\vect r_{p}}{R_{p}^{>}}\backslash B_{\vect r_{p}}\left(R_{p}\right)$
+\end_inset
+
+ around each one that contains only the background medium without any scatterers.
+ Then the EM field inside each such volume can be expanded in a way similar
+ to 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:E field expansion"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, using VSWFs with origins shifted to the centre of the volume:
+\begin_inset Formula 
+\begin{align}
+\vect E\left(\omega,\vect r\right) & =\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\left(\rcoeffptlm p{\tau}lm\vswfrtlm{\tau}lm\left(k\left(\vect r-\vect r_{p}\right)\right)+\outcoeffptlm p{\tau}lm\vswfouttlm{\tau}lm\left(k\left(\vect r-\vect r_{p}\right)\right)\right),\label{eq:E field expansion multiparticle}\\
+ & \vect r\in\openball{\vect r_{p}}{R_{p}^{>}}\backslash B_{\vect r_{p}}\left(R_{p}\right).\nonumber 
+\end{align}
+
+\end_inset
+
+Unlike the single scatterer case, the incident field coefficients 
+\begin_inset Formula $\rcoeffptlm p{\tau}lm$
+\end_inset
+
+ here are not only due to some external driving field that the particle
+ does not influence but they also contain the contributions of fields scattered
+ from 
+\emph on
+all other scatterers
+\emph default
+: 
+\begin_inset Formula 
+\begin{equation}
+\rcoeffp p=\rcoeffincp p+\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\tropsp pq\outcoeffp q\label{eq:particle total incident field coefficient a}
+\end{equation}
+
+\end_inset
+
+where 
+\begin_inset Formula $\rcoeffincp p$
+\end_inset
+
+ represents the part due to the external driving that the scatterers can
+ not influence, and 
+\begin_inset Formula $\tropsp pq$
+\end_inset
+
+ is a 
+\emph on
+translation operator
+\emph default
+ defined below in Sec.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "subsec:Translation-operator"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, that contains the re-expansion coefficients of the outgoing waves in origin
+ 
+\begin_inset Formula $\vect r_{q}$
+\end_inset
+
+ into regular waves in origin 
+\begin_inset Formula $\vect r_{p}$
+\end_inset
+
+.
+ For each scatterer, we also have its 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix relation as in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:T-matrix definition"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+,
+\begin_inset Formula 
+\[
+\outcoeffp q=T_{q}\rcoeffp q.
+\]
+
+\end_inset
+
+Together with 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:particle total incident field coefficient a"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, this gives rise to a set of linear equations
+\begin_inset Formula 
+\begin{equation}
+\outcoeffp p-T_{p}\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\tropsp pq\outcoeffp q=T_{p}\rcoeffincp p,\quad p\in\mathcal{P}\label{eq:Multiple-scattering problem}
+\end{equation}
+
+\end_inset
+
+which defines the multiple-scattering problem.
+ If all the 
+\begin_inset Formula $p,q$
+\end_inset
+
+-indexed vectors and matrices (note that without truncation, they are infinite-d
+imensional) are arranged into blocks of even larger vectors and matrices,
+ this can be written in a short-hand form
+\begin_inset Formula 
+\begin{equation}
+\left(I-T\trops\right)\outcoeff=T\rcoeffinc\label{eq:Multiple-scattering problem block form}
+\end{equation}
+
+\end_inset
+
+where 
+\begin_inset Formula $I$
+\end_inset
+
+ is the identity matrix, 
+\begin_inset Formula $T$
+\end_inset
+
+is a block-diagonal matrix containing all the individual 
+\begin_inset Formula $T$
+\end_inset
+
+-matrices, and 
+\begin_inset Formula $\trops$
+\end_inset
+
+ contains the individual 
+\begin_inset Formula $\tropsp pq$
+\end_inset
+
+matrices as the off-diagonal blocks, whereas the diagonal blocks are set
+ to zeros.
+\end_layout
+
+\begin_layout Standard
+In practice, the multiple-scattering problem is solved in its truncated
+ form, in which all the 
+\begin_inset Formula $l$
+\end_inset
+
+-indices related to a given scatterer 
+\begin_inset Formula $p$
+\end_inset
+
+ are truncated as 
+\begin_inset Formula $l\le L_{p}$
+\end_inset
+
+, laeving only 
+\begin_inset Formula $N_{p}=2L_{p}\left(L_{p}+2\right)$
+\end_inset
+
+ different 
+\begin_inset Formula $\tau lm$
+\end_inset
+
+-multiindices left.
+ The truncation degree can vary for different scatterers (e.g.
+ due to different physical sizes), so the truncated block 
+\begin_inset Formula $\tropsp pq$
+\end_inset
+
+ has shape 
+\begin_inset Formula $N_{p}\times N_{q}$
+\end_inset
+
+, not necessarily square.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+Such truncation of the translation operator 
+\begin_inset Formula $\tropsp pq$
+\end_inset
+
+ is justified by the fact on the left, TODO
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+If no other type of truncation is done, there remain 
+\begin_inset Formula $2L_{p}\left(L_{p}+2\right)$
+\end_inset
+
+ different 
+\begin_inset Formula $\tau lm$
+\end_inset
+
+-multiindices for 
+\begin_inset Formula $p$
+\end_inset
+
+-th scatterer, so that the truncated version of the matrix 
+\begin_inset Formula $\left(I-T\trops\right)$
+\end_inset
+
+ is a square matrix with 
+\begin_inset Formula $\left(\sum_{p\in\mathcal{P}}N_{p}\right)^{2}$
+\end_inset
+
+ elements in total.
+ The truncated problem 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can then be solved using standard numerical linear algebra methods (typically,
+ by LU factorisation of the 
+\begin_inset Formula $\left(I-T\trops\right)$
+\end_inset
+
+ matrix at a given frequency, and then solving with Gauss elimination for
+ as many different incident 
+\begin_inset Formula $\rcoeffinc$
+\end_inset
+
+ vectors as needed).
+\end_layout
+
+\begin_layout Standard
+Alternatively, the multiple scattering problem can be formulated in terms
+ of the regular field expansion coefficients,
+\begin_inset Formula 
+\begin{align*}
+\rcoeffp p-\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\tropsp pqT_{q}\rcoeffp q & =\rcoeffincp p,\quad p\in\mathcal{P},\\
+\left(I-\trops T\right)\rcoeff & =\rcoeffinc,
+\end{align*}
+
+\end_inset
+
+but this form is less suitable for numerical calculations due to the fact
+ that the regular VSWF expansion coefficients on both sides of the equation
+ are typically non-negligible even for large multipole degree 
+\begin_inset Formula $l$
+\end_inset
+
+, hence the truncation is not justified in this case.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO less bulshit.
+\end_layout
+
+\end_inset
+
+
 \end_layout
 
 \begin_layout Subsubsection
-
-\lang english
 Translation operator
+\begin_inset CommandInset label
+LatexCommand label
+name "subsec:Translation-operator"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+Let 
+\begin_inset Formula $\vect r_{1},\vect r_{2}$
+\end_inset
+
+ be two different origins; a regular VSWF with origin 
+\begin_inset Formula $\vect r_{1}$
+\end_inset
+
+ can be always expanded in terms of regular VSWFs with origin 
+\begin_inset Formula $\vect r_{2}$
+\end_inset
+
+ as follows:
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula 
+\begin{equation}
+\vswfrtlm{\tau}lm\left(k\left(\vect r-\vect r_{1}\right)\right)=\sum_{\tau'l'm'}\tropr_{\tau lm;\tau'l'm'}\left(k\left(\vect r_{2}-\vect r_{1}\right)\right)\vswfrtlm{\tau'}{l'}{m'}\left(\vect r-\vect r_{2}\right),\label{eq:regular vswf translation}
+\end{equation}
+
+\end_inset
+
+where an explicit formula for the (regular) 
+\emph on
+translation operator
+\emph default
+ 
+\begin_inset Formula $\tropr$
+\end_inset
+
+ reads in eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:translation operator"
+
+\end_inset
+
+ below.
+ For singular (outgoing) waves, the form of the expansion differs inside
+ and outside the ball 
+\begin_inset Formula $\openball{\left\Vert \vect r_{2}-\vect r_{1}\right\Vert }{\vect r_{1}}:$
+\end_inset
+
+ 
+\begin_inset Formula 
+\begin{eqnarray}
+\vswfouttlm{\tau}lm\left(k\left(\vect r-\vect r_{1}\right)\right) & = & \begin{cases}
+\sum_{\tau'l'm'}\trops_{\tau lm;\tau'l'm'}\left(k\left(\vect r_{2}-\vect r_{1}\right)\right)\vswfouttlm{\tau'}{l'}{m'}\left(\vect r-\vect r_{2}\right), & \vect r\in\openball{\left\Vert \vect r_{2}-\vect r_{1}\right\Vert }{\vect r_{1}}\\
+\sum_{\tau'l'm'}\tropr_{\tau lm;\tau'l'm'}\left(k\left(\vect r_{2}-\vect r_{1}\right)\right)\vswfrtlm{\tau'}{l'}{m'}\left(\vect r-\vect r_{2}\right), & \vect r\notin\closedball{\left\Vert \vect r_{2}-\vect r_{1}\right\Vert }{\left|\vect r_{1}\right|}
+\end{cases},\label{eq:singular vswf translation}
+\end{eqnarray}
+
+\end_inset
+
+where the singular translation operator 
+\begin_inset Formula $\trops$
+\end_inset
+
+ has the same form as 
+\begin_inset Formula $\tropr$
+\end_inset
+
+ in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:translation operator"
+
+\end_inset
+
+ except the regular spherical Bessel functions 
+\begin_inset Formula $j_{l}$
+\end_inset
+
+ are replaced with spherical Hankel functions 
+\begin_inset Formula $h_{l}^{(1)}$
+\end_inset
+
+.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO note about expansion exactly on the sphere.
+\end_layout
+
+\end_inset
+
+ 
+\end_layout
+
+\begin_layout Standard
+As MSTMM deals most of the time with the 
+\emph on
+expansion coefficients
+\emph default
+ of fields 
+\begin_inset Formula $\rcoeffptlm p{\tau}lm,\outcoeffptlm p{\tau}lm$
+\end_inset
+
+ in different origins 
+\begin_inset Formula $\vect r_{p}$
+\end_inset
+
+ rather than with the VSWFs directly, let us write down how 
+\emph on
+they
+\emph default
+ transform under translation.
+ Let us assume the field can be in terms of regular waves everywhere, and
+ expand it in two different origins 
+\begin_inset Formula $\vect r_{p},\vect r_{q}$
+\end_inset
+
+,
+\begin_inset Formula 
+\[
+\vect E\left(\vect r,\omega\right)=\sum_{\tau,l,m}\rcoeffptlm p{\tau}lm\vswfrtlm{\tau}lm\left(k\left(\vect r-\vect r_{p}\right)\right)=\sum_{\tau',l',m'}\rcoeffptlm q{\tau'}{l'}{m'}\vswfrtlm{\tau}{'l'}{m'}\left(k\left(\vect r-\vect r_{q}\right)\right).
+\]
+
+\end_inset
+
+Re-expanding the waves around 
+\begin_inset Formula $\vect r_{p}$
+\end_inset
+
+ in terms of waves around 
+\begin_inset Formula $\vect r_{q}$
+\end_inset
+
+ using 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular vswf translation"
+
+\end_inset
+
+, 
+\begin_inset Formula 
+\[
+\vect E\left(\vect r,\omega\right)=\sum_{\tau,l,m}\rcoeffptlm p{\tau}lm\sum_{\tau'l'm'}\tropr_{\tau lm;\tau'l'm'}\left(k\left(\vect r_{q}-\vect r_{p}\right)\right)\vswfrtlm{\tau'}{l'}{m'}\left(\vect r-\vect r_{q}\right)
+\]
+
+\end_inset
+
+and comparing to the original expansion around 
+\begin_inset Formula $\vect r_{q}$
+\end_inset
+
+, we obtain
+\begin_inset Formula 
+\begin{equation}
+\rcoeffptlm q{\tau'}{l'}{m'}=\sum_{\tau,l,m}\tropr_{\tau lm;\tau'l'm'}\left(k\left(\vect r_{q}-\vect r_{p}\right)\right)\rcoeffptlm p{\tau}lm.\label{eq:regular vswf coefficient translation}
+\end{equation}
+
+\end_inset
+
+For the sake of readability, we introduce a shorthand matrix form for 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular vswf coefficient translation"
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{equation}
+\rcoeffp q=\troprp qp\rcoeffp p\label{eq:reqular vswf coefficient vector translation}
+\end{equation}
+
+\end_inset
+
+(note the reversed indices; TODO redefine them in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular vswf translation"
+
+\end_inset
+
+, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:singular vswf translation"
+
+\end_inset
+
+? Similarly, if we had only outgoing waves in the original expansion around
+ 
+\begin_inset Formula $\vect r_{p}$
+\end_inset
+
+, we would get
+\begin_inset Formula 
+\begin{equation}
+\rcoeffp q=\tropsp qp\outcoeffp p\label{eq:singular to regular vswf coefficient vector translation}
+\end{equation}
+
+\end_inset
+
+for the expansion inside the ball 
+\begin_inset Formula $\openball{\left\Vert \vect r_{q}-\vect r_{p}\right\Vert }{\vect r_{p}}$
+\end_inset
+
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+CHECKME
+\end_layout
+
+\end_inset
+
+ and 
+\begin_inset Formula 
+\begin{equation}
+\outcoeffp q=\troprp qp\outcoeffp p\label{eq:singular to singular vswf coefficient vector translation-1}
+\end{equation}
+
+\end_inset
+
+outside.
+\end_layout
+
+\begin_layout Standard
+In our convention, the regular translation operator can be expressed explicitly
+ as (TODO CHECK CAREFULLY FOR POSSIBLE 
+\begin_inset Formula $(-1)^{m'}$
+\end_inset
+
+ AND SIMILAR FACTORS AND REWRITE IN TERMS OF SPHERICAL HARMONICS)
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+Teďka jsou tam zkopírovány výrazy pro C a D z Kristenssona, chybějí fase
+\end_layout
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{multline}
+\tropr_{\tau lm;\tau l'm'}\left(\vect d\right)=\frac{\left(-1\right)^{m+m'}}{2}\sum_{\lambda=\left|l-l'\right|}^{l+l'}\left(-1\right)^{\frac{l'-l+\lambda}{2}}\left(2\lambda+1\right)\sqrt{\frac{\left(2l+1\right)\left(2l'+1\right)\left(\lambda-\left(m-m'\right)\right)!}{l\left(l+1\right)l'\left(l'+1\right)\left(\lambda+\left(m-m'\right)\right)!}}\times\\
+\times\begin{pmatrix}l & l' & \lambda\\
+0 & 0 & 0
+\end{pmatrix}\begin{pmatrix}l & l' & \lambda\\
+m & -m' & m'-m
+\end{pmatrix}\left(l\left(l+1\right)+l'\left(l'+1\right)-\lambda\left(\lambda+1\right)\right)\times\\
+\times j_{\lambda}\left(d\right)P_{\lambda}^{m-m'}\left(\cos\theta_{\uvec d}\right)e^{i\left(m-m'\right)\phi_{\uvec d}},\\
+\tropr_{\tau lm;\tau'l'm'}\left(\vect d\right)=-i\frac{\left(-1\right)^{m+m'}}{2}\sum_{\lambda=\left|l-l'\right|+1}^{l+l'}\left(-1\right)^{\frac{l'-l+\lambda+1}{2}}\left(2\lambda+1\right)\sqrt{\frac{\left(2l+1\right)\left(2l'+1\right)\left(\lambda-\left(m-m'\right)\right)!}{l\left(l+1\right)l'\left(l'+1\right)\left(\lambda+\left(m-m'\right)\right)!}}\times\\
+\times\begin{pmatrix}l & l' & \lambda-1\\
+0 & 0 & 0
+\end{pmatrix}\begin{pmatrix}l & l' & \lambda\\
+m & -m' & m'-m
+\end{pmatrix}\sqrt{\lambda^{2}-\left(l-l'\right)^{2}}\sqrt{\left(l+l'+1\right)^{2}-\lambda^{2}}\times\\
+\times j_{\lambda}\left(d\right)P_{\lambda}^{m-m'}\left(\cos\theta_{\uvec d}\right)e^{i\left(m-m'\right)\phi_{\uvec d}},\qquad\tau\ne\tau'.\label{eq:translation operator}
+\end{multline}
+
+\end_inset
+
+The singular operator 
+\begin_inset Formula $\trops$
+\end_inset
+
+ for re-expanding outgoing waves into regular ones has the same form except
+ the regular spherical Bessel functions 
+\begin_inset Formula $j_{l}$
+\end_inset
+
+ in are replaced with spherical Hankel functions 
+\begin_inset Formula $h_{l}^{(1)}=j_{l}+iy_{l}$
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+In our convention, the regular translation operator is unitary, 
+\begin_inset Formula $\left(\tropr_{\tau lm;\tau'l'm'}\left(\vect d\right)\right)^{-1}=\tropr_{\tau lm;\tau'l'm'}\left(-\vect d\right)=\tropr_{\tau'l'm';\tau lm}^{*}\left(\vect d\right)$
+\end_inset
+
+,
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+todo different notation for the complex conjugation without transposition???
+\end_layout
+
+\end_inset
+
+ or in the per-particle matrix notation, 
+\begin_inset Formula 
+\begin{equation}
+\troprp qp^{-1}=\troprp pq=\troprp qp^{\dagger}\label{eq:regular translation unitarity}
+\end{equation}
+
+\end_inset
+
+.
+ Note that truncation at finite multipole degree breaks the unitarity, 
+\begin_inset Formula $\truncated{\troprp qp}L^{-1}\ne\truncated{\troprp pq}L=\truncated{\troprp qp^{\dagger}}L$
+\end_inset
+
+, which has to be taken into consideration when evaluating quantities such
+ as absorption or scattering cross sections.
+ Similarly, the full regular operators can be composed
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+better wording
+\end_layout
+
+\end_inset
+
+, 
+\begin_inset Formula 
+\begin{equation}
+\troprp ac=\troprp ab\troprp bc\label{eq:regular translation composition}
+\end{equation}
+
+\end_inset
+
+ but truncation breaks this, 
+\begin_inset Formula $\truncated{\troprp ac}L\ne\truncated{\troprp ab}L\truncated{\troprp bc}L.$
+\end_inset
+
+
 \end_layout
 
 \begin_layout Subsubsection
+Cross-sections (many scatterers)
+\end_layout
 
-\lang english
-Numerical complexity, comparison to other methods
+\begin_layout Standard
+For a system of many scatterers, Kristensson 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "sect. 9.2.2"
+key "kristensson_scattering_2016"
+literal "false"
+
+\end_inset
+
+ derives only the extinction cross section formula.
+ Let us re-derive it together with the many-particle scattering and absorption
+ cross sections.
+ First, let us take a ball circumscribing all the scatterers at once, 
+\begin_inset Formula $\openball R{\vect r_{\square}}\supset\particle$
+\end_inset
+
+.
+ Outside 
+\begin_inset Formula $\openball R{\vect r_{\square}}$
+\end_inset
+
+, we can describe the EM fields as if there was only a single scatterer,
+\begin_inset Formula 
+\[
+\vect E\left(\vect r\right)=\sum_{\tau,l,m}\left(\rcoeffptlm{\square}{\tau}lm\vswfrtlm{\tau}lm\left(\vect r-\vect r_{\square}\right)+\outcoeffptlm{\square}{\tau}lm\vswfouttlm{\tau}lm\left(\vect r-\vect r_{\square}\right)\right),
+\]
+
+\end_inset
+
+where 
+\begin_inset Formula $\rcoeffp{\square},\outcoeffp{\square}$
+\end_inset
+
+ are the vectors of VSWF expansion coefficients of the incident and total
+ scattered fields, respectively, at origin 
+\begin_inset Formula $\vect r_{\square}$
+\end_inset
+
+.
+ In principle, one could evaluate 
+\begin_inset Formula $\outcoeffp{\square}$
+\end_inset
+
+ using the translation operators (REF!!!) and use the single-scatterer formulae
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:extincion CS single"
+
+\end_inset
+
+–
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS single"
+
+\end_inset
+
+ with 
+\begin_inset Formula $\rcoeffp{}=\rcoeffp{\square},\outcoeffp{}=\outcoeffp{\square}$
+\end_inset
+
+ to obtain the cross sections.
+ However, this is not suitable for numerical evaluation with truncation
+ in multipole degree; hence we need to express them in terms of particle-wise
+ expansions 
+\begin_inset Formula $\rcoeffp p,\outcoeffp p$
+\end_inset
+
+.
+ The original incident field re-expanded around 
+\begin_inset Formula $p$
+\end_inset
+
+-th particle reads according to 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular vswf translation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{equation}
+\rcoeffincp p=\troprp p{\square}\rcoeffp{\square}\label{eq:a_inc local from global}
+\end{equation}
+
+\end_inset
+
+whereas the contributions of fields scattered from each particle expanded
+ around the global origin 
+\begin_inset Formula $\vect r_{\square}$
+\end_inset
+
+ is, according to 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:singular vswf translation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+,
+\begin_inset Formula 
+\begin{equation}
+\outcoeffp{\square}=\sum_{q\in\mathcal{P}}\troprp{\square}q\outcoeffp q.\label{eq:f global from local}
+\end{equation}
+
+\end_inset
+
+Using the unitarity 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular translation unitarity"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ and composition 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:regular translation composition"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ properties, one has
+\begin_inset Formula 
+\begin{align}
+\rcoeffp{\square}^{\dagger}\outcoeffp{\square} & =\rcoeffincp p^{\dagger}\troprp p{\square}\troprp{\square}q\outcoeffp q=\rcoeffincp p^{\dagger}\sum_{q\in\mathcal{P}}\troprp pqf_{q}\nonumber \\
+ & =\sum_{q\in\mathcal{P}}\left(\troprp qp\rcoeffincp p\right)^{\dagger}f_{q}=\sum_{q\in\mathcal{P}}\rcoeffincp q^{\dagger}f_{q},\label{eq:atf form multiparticle}
+\end{align}
+
+\end_inset
+
+where only the last expression is suitable for numerical evaluation with
+ truncated matrices, because the previous ones contain a translation operator
+ right next to an incident field coefficient vector (see Sec.
+ TODO).
+ Similarly, 
+\begin_inset Formula 
+\begin{align}
+\left\Vert \outcoeffp{\square}\right\Vert ^{2} & =\outcoeffp{\square}^{\dagger}\outcoeffp{\square}=\sum_{p\in\mathcal{P}}\left(\troprp{\square}p\outcoeffp p\right)^{\dagger}\sum_{q\in\mathcal{P}}\troprp{\square}q\outcoeffp q\nonumber \\
+ & =\sum_{p\in\mathcal{P}}\sum_{q\in\mathcal{P}}\outcoeffp p^{\dagger}\troprp pq\outcoeffp q.\label{eq:f squared form multiparticle}
+\end{align}
+
+\end_inset
+
+Substituting 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:atf form multiparticle"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:f squared form multiparticle"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ into 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:scattering CS single"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ and 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS single"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, we get the many-particle expressions for extinction, scattering and absorption
+ cross sections suitable for numerical evaluation: 
+\begin_inset Formula 
+\begin{eqnarray}
+\sigma_{\mathrm{ext}}\left(\uvec k\right) & = & -\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\Re\sum_{p\in\mathcal{P}}\rcoeffincp p^{\dagger}\outcoeffp p=-\frac{1}{2k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\Re\sum_{p\in\mathcal{P}}\rcoeffincp p^{\dagger}\left(\Tp p+\Tp p^{\dagger}\right)\rcoeffp p,\label{eq:extincion CS multi}\\
+\sigma_{\mathrm{scat}}\left(\uvec k\right) & = & \frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\sum_{q\in\mathcal{P}}\outcoeffp p^{\dagger}\troprp pq\outcoeffp q\nonumber \\
+ &  & =\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\sum_{q\in\mathcal{P}}\rcoeffp p^{\dagger}\Tp p^{\dagger}\troprp pq\Tp q\rcoeffp q,\label{eq:scattering CS multi}\\
+\sigma_{\mathrm{abs}}\left(\uvec k\right) & = & -\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\Re\left(\outcoeffp p^{\dagger}\left(\rcoeffincp p+\sum_{q\in\mathcal{P}}\troprp pq\outcoeffp q\right)\right).\nonumber \\
+\label{eq:absorption CS multi}
+\end{eqnarray}
+
+\end_inset
+
+
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+\begin_inset Formula $=\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\mbox{TODO}.$
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+An alternative approach to derive the absorption cross section is via a
+ power transport argument.
+ Note the direct proportionality between absorption cross section 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS single"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ and net radiated power for single scatterer 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Power transport"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset Formula $\sigma_{\mathrm{abs}}=-\eta_{0}\eta P/2\left\Vert \vect E_{0}\right\Vert ^{2}$
+\end_inset
+
+.
+ In the many-particle setup (with non-lossy background medium, so that only
+ the particles absorb), the total absorbed power is equal to the sum of
+ absorbed powers on each particle, 
+\begin_inset Formula $-P=\sum_{p\in\mathcal{P}}-P_{p}$
+\end_inset
+
+.
+ Using the power transport formula 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Power transport"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ particle-wise gives
+\begin_inset Formula 
+\begin{equation}
+\sigma_{\mathrm{abs}}\left(\uvec k\right)=-\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\left(\Re\left(\rcoeffp p^{\dagger}\outcoeffp p\right)+\left\Vert \outcoeffp p\right\Vert ^{2}\right)\label{eq:absorption CS multi alternative}
+\end{equation}
+
+\end_inset
+
+which seems different from 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS multi"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, but using 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:particle total incident field coefficient a"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, we can rewrite it as
+\begin_inset Formula 
+\begin{align*}
+\sigma_{\mathrm{abs}}\left(\uvec k\right) & =-\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\Re\left(\outcoeffp p^{\dagger}\left(\rcoeffp p+\outcoeffp p\right)\right)\\
+ & =-\frac{1}{k^{2}\left\Vert \vect E_{0}\right\Vert ^{2}}\sum_{p\in\mathcal{P}}\Re\left(\outcoeffp p^{\dagger}\left(\rcoeffincp p+\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\tropsp pq\outcoeffp q+\outcoeffp p\right)\right).
+\end{align*}
+
+\end_inset
+
+It is easy to show that all the terms of 
+\begin_inset Formula $\sum_{p\in\mathcal{P}}\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\outcoeffp p^{\dagger}\tropsp pq\outcoeffp q$
+\end_inset
+
+ containing the singular spherical Bessel functions 
+\begin_inset Formula $y_{l}$
+\end_inset
+
+ are imaginary, 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO better formulation
+\end_layout
+
+\end_inset
+
+ so that actually 
+\begin_inset Formula $\sum_{p\in\mathcal{P}}\Re\left(\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\outcoeffp p^{\dagger}\tropsp pq\outcoeffp q+\outcoeffp p^{\dagger}\outcoeffp p\right)=\sum_{p\in\mathcal{P}}\Re\left(\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\outcoeffp p^{\dagger}\troprp pq\outcoeffp q+\outcoeffp p^{\dagger}\outcoeffp p\right)=\sum_{p\in\mathcal{P}}\Re\left(\sum_{q\in\mathcal{P}\backslash\left\{ p\right\} }\outcoeffp p^{\dagger}\troprp pq\outcoeffp q+\outcoeffp p^{\dagger}\troprp pp\outcoeffp p\right)=\sum_{p\in\mathcal{P}}\sum_{q\in\mathcal{P}}\outcoeffp p^{\dagger}\troprp pq\outcoeffp q,$
+\end_inset
+
+ proving that the expressions in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS multi"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ and 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:absorption CS multi alternative"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ are equal.
 \end_layout
 
 \end_body
diff --git a/lepaper/infinite-old.lyx b/lepaper/infinite-old.lyx
new file mode 100644
index 0000000..8e8296f
--- /dev/null
+++ b/lepaper/infinite-old.lyx
@@ -0,0 +1,496 @@
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
+\begin_document
+\begin_header
+\save_transient_properties true
+\origin unavailable
+\textclass article
+\use_default_options false
+\maintain_unincluded_children false
+\language english
+\language_package none
+\inputencoding utf8
+\fontencoding default
+\font_roman "default" "default"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
+\font_default_family default
+\use_non_tex_fonts false
+\font_sc false
+\font_roman_osf false
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
+\graphics default
+\default_output_format default
+\output_sync 0
+\bibtex_command default
+\index_command default
+\float_placement class
+\float_alignment class
+\paperfontsize default
+\spacing single
+\use_hyperref false
+\papersize default
+\use_geometry false
+\use_package amsmath 1
+\use_package amssymb 0
+\use_package cancel 0
+\use_package esint 1
+\use_package mathdots 0
+\use_package mathtools 0
+\use_package mhchem 0
+\use_package stackrel 0
+\use_package stmaryrd 0
+\use_package undertilde 0
+\cite_engine basic
+\cite_engine_type default
+\biblio_style plain
+\use_bibtopic false
+\use_indices false
+\paperorientation portrait
+\suppress_date false
+\justification true
+\use_refstyle 0
+\use_minted 0
+\use_lineno 0
+\index Index
+\shortcut idx
+\color #008000
+\end_index
+\secnumdepth 3
+\tocdepth 3
+\paragraph_separation indent
+\paragraph_indentation default
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
+\papercolumns 1
+\papersides 1
+\paperpagestyle default
+\tablestyle default
+\tracking_changes false
+\output_changes false
+\html_math_output 0
+\html_css_as_file 0
+\html_be_strict false
+\end_header
+
+\begin_body
+
+\begin_layout Subsection
+Periodic systems and mode analysis
+\begin_inset CommandInset label
+LatexCommand label
+name "subsec:Periodic-systems"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+In an infinite periodic array of nanoparticles, the excitations of the nanoparti
+cles take the quasiperiodic Bloch-wave form 
+\begin_inset Formula 
+\[
+\coeffs_{i\nu}=e^{i\vect k\cdot\vect R_{i}}\coeffs_{\nu}
+\]
+
+\end_inset
+
+(assuming the incident external field has the same periodicity, 
+\begin_inset Formula $\coeffr_{\mathrm{ext}(i\nu)}=e^{i\vect k\cdot\vect R_{i}}p_{\mathrm{ext}\left(\nu\right)}$
+\end_inset
+
+) where 
+\begin_inset Formula $\nu$
+\end_inset
+
+ is the index of a particle inside one unit cell and 
+\begin_inset Formula $\vect R_{i},\vect R_{i'}\in\Lambda$
+\end_inset
+
+ are the lattice vectors corresponding to the sites (labeled by multiindices
+ 
+\begin_inset Formula $i,i'$
+\end_inset
+
+) of a Bravais lattice 
+\begin_inset Formula $\Lambda$
+\end_inset
+
+.
+ The multiple-scattering problem (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:multiple scattering per particle a"
+
+\end_inset
+
+) then takes the form
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula 
+\[
+\coeffs_{i\nu}-T_{\nu}\sum_{(i',\nu')\ne\left(i,\nu\right)}S_{i\nu,i'\nu'}e^{i\vect k\cdot\left(\vect R_{i'}-\vect R_{i}\right)}\coeffs_{i\nu'}=T_{\nu}\coeffr_{\mathrm{ext}(i\nu)}
+\]
+
+\end_inset
+
+or, labeling 
+\begin_inset Formula $W_{\nu\nu'}=\sum_{i';(i',\nu')\ne\left(i,\nu\right)}S_{i\nu,i'\nu'}e^{i\vect k\cdot\left(\vect R_{i'}-\vect R_{i}\right)}=\sum_{i';(i',\nu')\ne\left(0,\nu\right)}S_{0\nu,i'\nu'}e^{i\vect k\cdot\vect R_{i'}}$
+\end_inset
+
+ and using the quasiperiodicity, 
+\begin_inset Formula 
+\begin{equation}
+\sum_{\nu'}\left(\delta_{\nu\nu'}\mathbb{I}-T_{\nu}W_{\nu\nu'}\right)\coeffs_{\nu'}=T_{\nu}\coeffr_{\mathrm{ext}(\nu)},\label{eq:multiple scattering per particle a periodic}
+\end{equation}
+
+\end_inset
+
+which reduces the linear problem (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:multiple scattering per particle a"
+
+\end_inset
+
+) to interactions between particles inside single unit cell.
+ A problematic part is the evaluation of the translation operator lattice
+ sums 
+\begin_inset Formula $W_{\nu\nu'}$
+\end_inset
+
+; this is performed using exponentially convergent Ewald-type representations
+ 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "linton_lattice_2010"
+literal "true"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+In an infinite periodic system, a nonlossy mode supports itself without
+ external driving, i.e.
+ such mode is described by excitation coefficients 
+\begin_inset Formula $a_{\nu}$
+\end_inset
+
+ that satisfy eq.
+ (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:multiple scattering per particle a periodic"
+
+\end_inset
+
+) with zero right-hand side.
+ That can happen if the block matrix 
+\begin_inset Formula 
+\begin{equation}
+M\left(\omega,\vect k\right)=\left\{ \delta_{\nu\nu'}\mathbb{I}-T_{\nu}\left(\omega\right)W_{\nu\nu'}\left(\omega,\vect k\right)\right\} _{\nu\nu'}\label{eq:M matrix definition}
+\end{equation}
+
+\end_inset
+
+from the left hand side of (
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:multiple scattering per particle a periodic"
+
+\end_inset
+
+) is singular (here we explicitly note the 
+\begin_inset Formula $\omega,\vect k$
+\end_inset
+
+ depence).
+\end_layout
+
+\begin_layout Standard
+For lossy nanoparticles, however, perfect propagating modes will not exist
+ and 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ will never be perfectly singular.
+ Therefore in practice, we get the bands by scanning over 
+\begin_inset Formula $\omega,\vect k$
+\end_inset
+
+ to search for 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ which have an 
+\begin_inset Quotes erd
+\end_inset
+
+almost zero
+\begin_inset Quotes erd
+\end_inset
+
+ singular value.
+\end_layout
+
+\begin_layout Section
+\begin_inset ERT
+status collapsed
+
+\begin_layout Plain Layout
+
+{
+\end_layout
+
+\end_inset
+
+Symmetries
+\begin_inset ERT
+status collapsed
+
+\begin_layout Plain Layout
+
+}
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset CommandInset label
+LatexCommand label
+name "sm:symmetries"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+A general overview of utilizing group theory to find lattice modes at high-symme
+try points of the Brillouin zone can be found e.g.
+ in 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "chapters 10–11"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+; here we use the same notation.
+\end_layout
+
+\begin_layout Standard
+We analyse the symmetries of the system in the same VSWF representation
+ as used in the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix formalism introduced above.
+ We are interested in the modes at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point of the hexagonal lattice, which has the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ point symmetry.
+ The six irreducible representations (irreps) of the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ group are known and are available in the literature in their explicit forms.
+ In order to find and classify the modes, we need to find a decomposition
+ of the lattice mode representation 
+\begin_inset Formula $\Gamma_{\mathrm{lat.mod.}}=\Gamma^{\mathrm{equiv.}}\otimes\Gamma_{\mathrm{vec.}}$
+\end_inset
+
+ into the irreps of 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+.
+ The equivalence representation 
+\begin_inset Formula $\Gamma^{\mathrm{equiv.}}$
+\end_inset
+
+ is the 
+\begin_inset Formula $E'$
+\end_inset
+
+ representation as can be deduced from 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (11.19)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+, eq.
+ (11.19) and the character table for 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+.
+ 
+\begin_inset Formula $\Gamma_{\mathrm{vec.}}$
+\end_inset
+
+ operates on a space spanned by the VSWFs around each nanoparticle in the
+ unit cell (the effects of point group operations on VSWFs are described
+ in 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "schulz_point-group_1999"
+literal "true"
+
+\end_inset
+
+).
+ This space can be then decomposed into invariant subspaces of the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ using the projectors 
+\begin_inset Formula $\hat{P}_{ab}^{\left(\Gamma\right)}$
+\end_inset
+
+ defined by 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (4.28)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+.
+ This way, we obtain a symmetry adapted basis 
+\begin_inset Formula $\left\{ \vect b_{\Gamma,r,i}^{\mathrm{s.a.b.}}\right\} $
+\end_inset
+
+ as linear combinations of VSWFs 
+\begin_inset Formula $\vswfs lm{p,t}$
+\end_inset
+
+ around the constituting nanoparticles (labeled 
+\begin_inset Formula $p$
+\end_inset
+
+), 
+\begin_inset Formula 
+\[
+\vect b_{\Gamma,r,i}^{\mathrm{s.a.b.}}=\sum_{l,m,p,t}U_{\Gamma,r,i}^{p,t,l,m}\vswfs lm{p,t},
+\]
+
+\end_inset
+
+where 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+ stands for one of the six different irreps of 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+, 
+\begin_inset Formula $r$
+\end_inset
+
+ labels the different realisations of the same irrep, and the last index
+ 
+\begin_inset Formula $i$
+\end_inset
+
+ going from 1 to 
+\begin_inset Formula $d_{\Gamma}$
+\end_inset
+
+ (the dimensionality of 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+) labels the different partners of the same given irrep.
+ The number of how many times is each irrep contained in 
+\begin_inset Formula $\Gamma_{\mathrm{lat.mod.}}$
+\end_inset
+
+ (i.e.
+ the range of index 
+\begin_inset Formula $r$
+\end_inset
+
+ for given 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+) depends on the multipole degree cutoff 
+\begin_inset Formula $l_{\mathrm{max}}$
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+Each mode at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point shall lie in the irreducible spaces of only one of the six possible
+ irreps and it can be shown via 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (2.51)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+ that, at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point, the matrix 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ defined above takes a block-diagonal form in the symmetry-adapted basis,
+ 
+\begin_inset Formula 
+\[
+M\left(\omega,\vect K\right)_{\Gamma,r,i;\Gamma',r',j}^{\mathrm{s.a.b.}}=\frac{\delta_{\Gamma\Gamma'}\delta_{ij}}{d_{\Gamma}}\sum_{q}M\left(\omega,\vect K\right)_{\Gamma,r,q;\Gamma',r',q}^{\mathrm{s.a.b.}}.
+\]
+
+\end_inset
+
+This enables us to decompose the matrix according to the irreps and to solve
+ the singular value problem in each irrep separately, as done in Fig.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "smfig:dispersions"
+
+\end_inset
+
+(a).
+\end_layout
+
+\end_body
+\end_document
diff --git a/lepaper/infinite.lyx b/lepaper/infinite.lyx
index 26dfe37..9e74628 100644
--- a/lepaper/infinite.lyx
+++ b/lepaper/infinite.lyx
@@ -1,29 +1,37 @@
-#LyX 2.1 created this file. For more info see http://www.lyx.org/
-\lyxformat 474
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
 \begin_document
 \begin_header
+\save_transient_properties true
+\origin unavailable
 \textclass article
 \use_default_options true
 \maintain_unincluded_children false
-\language finnish
+\language english
 \language_package default
-\inputencoding auto
-\fontencoding global
-\font_roman TeX Gyre Pagella
-\font_sans default
-\font_typewriter default
-\font_math auto
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
 \font_default_family default
-\use_non_tex_fonts true
+\use_non_tex_fonts false
 \font_sc false
-\font_osf true
-\font_sf_scale 100
-\font_tt_scale 100
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
 \graphics default
-\default_output_format pdf4
+\default_output_format default
 \output_sync 0
 \bibtex_command default
 \index_command default
+\float_placement class
+\float_alignment class
 \paperfontsize default
 \spacing single
 \use_hyperref true
@@ -58,6 +66,8 @@
 \suppress_date false
 \justification true
 \use_refstyle 1
+\use_minted 0
+\use_lineno 0
 \index Index
 \shortcut idx
 \color #008000
@@ -66,10 +76,14 @@
 \tocdepth 3
 \paragraph_separation indent
 \paragraph_indentation default
-\quotes_language swedish
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
 \papercolumns 1
 \papersides 1
 \paperpagestyle default
+\tablestyle default
 \tracking_changes false
 \output_changes false
 \html_math_output 0
@@ -81,123 +95,201 @@
 
 \begin_layout Section
 Infinite periodic systems
+\begin_inset FormulaMacro
+\newcommand{\dlv}{\vect a}
+\end_inset
+
+
+\begin_inset FormulaMacro
+\newcommand{\rlv}{\vect b}
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+Although large finite systems are where MSTMM excels the most, there are
+ several reasons that makes its extension to infinite lattices (where periodic
+ boundary conditions might be applied) desirable as well.
+ Other methods might be already fast enough, but MSTMM will be faster in
+ most cases in which there is enough spacing between the neighboring particles.
+ MSTMM works well with any space group symmetry the system might have (as
+ opposed to, for example, FDTD with cubic mesh applied to a honeycomb lattice),
+ which makes e.g.
+ application of group theory in mode analysis quite easy.
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+Topology anoyne?
+\end_layout
+
+\end_inset
+
+ And finally, having a method that handles well both infinite and large
+ finite system gives a possibility to study finite-size effects in periodic
+ scatterer arrays.
 \end_layout
 
 \begin_layout Subsection
+Notation
+\end_layout
 
-\lang english
+\begin_layout Standard
+TODO Fourier transforms, Delta comb, lattice bases, reciprocal lattices
+ etc.
+\end_layout
+
+\begin_layout Subsection
 Formulation of the problem
 \end_layout
 
 \begin_layout Standard
-
-\lang english
-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
-\begin_inset Formula 
-\[
-A_{α}=T_{α}P_{α}=T_{α}(\sum_{β}S_{α\leftarrowβ}A_{β}+P_{0α})
-\]
+Let us have a linear system of compact EM scatterers on a homogeneous background
+ as in Section 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "subsec:Multiple-scattering"
+plural "false"
+caps "false"
+noprefix "false"
 
 \end_inset
 
-where 
-\begin_inset Formula $T_{α}$
+, but this time, the system shall be periodic: let there be a 
+\begin_inset Formula $d$
 \end_inset
 
- is the 
-\begin_inset Formula $T$
+-dimensional (
+\begin_inset Formula $d$
 \end_inset
 
--matrix for scatterer α, 
-\begin_inset Formula $A_{α}$
+ can be 1, 2 or 3) lattice embedded into the three-dimensional real space,
+ with lattice vectors 
+\begin_inset Formula $\left\{ \dlv_{i}\right\} _{i=1}^{d}$
 \end_inset
 
- is its vector of the scattered wave expansion coefficient (the multipole
- indices are not explicitely indicated here) and 
-\begin_inset Formula $P_{α}$
+, and let the lattice points be labeled with an 
+\begin_inset Formula $d$
 \end_inset
 
- is the local expansion of the incoming sources.
- 
-\begin_inset Formula $S_{α\leftarrowβ}$
+-dimensional integar multiindex 
+\begin_inset Formula $\vect n\in\ints^{d}$
 \end_inset
 
- is ...
- and ...
- is ...
-\end_layout
-
-\begin_layout Standard
-
-\lang english
-...
-\end_layout
-
-\begin_layout Standard
-
-\lang english
-\begin_inset Formula 
-\[
-\sum_{β}(\delta_{αβ}-T_{α}S_{α\leftarrowβ})A_{β}=T_{α}P_{0α}.
-\]
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Standard
-
-\lang english
-Now suppose that the scatterers constitute an infinite lattice
-\end_layout
-
-\begin_layout Standard
-
-\lang english
-\begin_inset Formula 
-\[
-\sum_{\vect bβ}(\delta_{\vect{ab}}\delta_{αβ}-T_{\vect aα}S_{\vect aα\leftarrow\vect bβ})A_{\vect bβ}=T_{\vect aα}P_{0\vect aα}.
-\]
-
-\end_inset
-
-Due to the periodicity, we can write 
-\begin_inset Formula $S_{\vect aα\leftarrow\vect bβ}=S_{α\leftarrowβ}(\vect b-\vect a)$
-\end_inset
-
- and 
-\begin_inset Formula $T_{\vect aα}=T_{\alpha}$
+, so the lattice points have cartesian coordinates 
+\begin_inset Formula $\vect R_{\vect n}=\sum_{i=1}^{d}n_{i}\vect a_{i}$
 \end_inset
 
 .
- In order to find lattice modes, we search for solutions with zero RHS
+ There can be several scatterers per unit cell with indices 
+\begin_inset Formula $\alpha$
+\end_inset
+
+ from set 
+\begin_inset Formula $\mathcal{P}_{1}$
+\end_inset
+
+ and (relative) positions inside the unit cell 
+\begin_inset Formula $\vect r_{\alpha}$
+\end_inset
+
+; any particle of the periodic system can thus be labeled by a multiindex
+ from 
+\begin_inset Formula $\mathcal{P}=\ints^{d}\times\mathcal{P}_{1}$
+\end_inset
+
+.
+ The scatterers are located at positions 
+\begin_inset Formula $\vect r_{\vect n,\alpha}=\vect R_{\vect n}+\vect r_{\alpha}$
+\end_inset
+
+ and their 
+\begin_inset Formula $T$
+\end_inset
+
+-matrices are periodic, 
+\begin_inset Formula $T_{\vect n,\alpha}=T_{\alpha}$
+\end_inset
+
+.
+ In such system, the multiple-scattering problem 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be rewritten as
+\end_layout
+
+\begin_layout Standard
 \begin_inset Formula 
-\[
-\sum_{\vect bβ}(\delta_{\vect{ab}}\delta_{αβ}-T_{α}S_{\vect aα\leftarrow\vect bβ})A_{\vect bβ}=0
-\]
+\begin{equation}
+\outcoeffp{\vect n,\alpha}-T_{\alpha}\sum_{\left(\vect m,\beta\right)\in\mathcal{P}\backslash\left\{ \left(\vect n,\alpha\right)\right\} }\tropsp{\vect n,\alpha}{\vect m,\beta}\outcoeffp{\vect m,\beta}=T_{\alpha}\rcoeffincp{\vect n,\alpha}.\quad\left(\vect n,\alpha\right)\in\mathcal{P}\label{eq:Multiple-scattering problem periodic}
+\end{equation}
 
 \end_inset
 
-and we assume periodic solution 
-\begin_inset Formula $A_{\vect b\beta}(\vect k)=A_{\vect a\beta}e^{i\vect k\cdot\vect r_{\vect b-\vect a}}$
+
+\end_layout
+
+\begin_layout Standard
+Due to periodicity, we can also write 
+\begin_inset Formula $\tropsp{\vect n,\alpha}{\vect m,\beta}=\tropsp{\alpha}{\beta}\left(\vect R_{\vect m}-\vect R_{\vect n}\right)=\tropsp{\alpha}{\beta}\left(\vect R_{\vect m-\vect n}\right)=\tropsp{\vect 0,\alpha}{\vect m-\vect n,\beta}$
 \end_inset
 
-, yielding
+.
+ Assuming quasi-periodic right-hand side with quasi-momentum 
+\begin_inset Formula $\vect k$
+\end_inset
+
+, 
+\begin_inset Formula $\rcoeffincp{\vect n,\alpha}=\rcoeffincp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}}$
+\end_inset
+
+, the solutions of 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem periodic"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ will be also quasi-periodic according to Bloch theorem, 
+\begin_inset Formula $\outcoeffp{\vect n,\alpha}=\outcoeffp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}}$
+\end_inset
+
+, and eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem periodic"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be rewritten as follows
 \begin_inset Formula 
-\begin{eqnarray*}
-\sum_{\vect bβ}(\delta_{\vect{ab}}\delta_{αβ}-T_{α}S_{\vect aα\leftarrow\vect bβ})A_{\vect a\beta}\left(\vect k\right)e^{i\vect k\cdot\vect r_{\vect b-\vect a}} & = & 0,\\
-\sum_{\vect bβ}(\delta_{\vect{0b}}\delta_{αβ}-T_{α}S_{\vect 0α\leftarrow\vect bβ})A_{\vect 0\beta}\left(\vect k\right)e^{i\vect k\cdot\vect r_{\vect b}} & = & 0,\\
-\sum_{β}(\delta_{αβ}-T_{α}\underbrace{\sum_{\vect b}S_{\vect 0α\leftarrow\vect bβ}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_{α}\sum_{\beta}W_{\alpha\beta}\left(\vect k\right)A_{\vect 0\beta}\left(\vect k\right) & = & 0.
-\end{eqnarray*}
+\begin{align}
+\outcoeffp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}}-T_{\alpha}\sum_{\left(\vect m,\beta\right)\in\mathcal{P}\backslash\left\{ \left(\vect n,\alpha\right)\right\} }\tropsp{\vect n,\alpha}{\vect m,\beta}\outcoeffp{\vect 0,\beta}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect m}} & =T_{\alpha}\rcoeffincp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}},\nonumber \\
+\outcoeffp{\vect 0,\alpha}\left(\vect k\right)-T_{\alpha}\sum_{\left(\vect m,\beta\right)\in\mathcal{P}\backslash\left\{ \left(\vect n,\alpha\right)\right\} }\tropsp{\vect 0,\alpha}{\vect m-\vect n,\beta}\outcoeffp{\vect 0,\beta}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect m-\vect n}} & =T_{\alpha}\rcoeffincp{\vect 0,\alpha}\left(\vect k\right),\nonumber \\
+\outcoeffp{\vect 0,\alpha}\left(\vect k\right)-T_{\alpha}\sum_{\left(\vect m,\beta\right)\in\mathcal{P}\backslash\left\{ \left(\vect 0,\alpha\right)\right\} }\tropsp{\vect 0,\alpha}{\vect m,\beta}\outcoeffp{\vect 0,\beta}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect m}} & =T_{\alpha}\rcoeffincp{\vect 0,\alpha}\left(\vect k\right),\nonumber \\
+\outcoeffp{\vect 0,\alpha}\left(\vect k\right)-T_{\alpha}\sum_{\beta\in\mathcal{P}}W_{\alpha\beta}\left(\vect k\right)\outcoeffp{\vect 0,\beta}\left(\vect k\right) & =T_{\alpha}\rcoeffincp{\vect 0,\alpha}\left(\vect k\right),\label{eq:Multiple-scattering problem unit cell}
+\end{align}
 
 \end_inset
 
-Therefore, in order to solve the modes, we need to compute the 
+so we reduced the initial scattering problem to one involving only the field
+ expansion coefficients from a single unit cell, but we need to compute
+ the 
 \begin_inset Quotes eld
 \end_inset
 
@@ -208,23 +300,272 @@ lattice Fourier transform
  of the translation operator,
 \begin_inset Formula 
 \begin{equation}
-W_{\alpha\beta}(\vect k)\equiv\sum_{\vect b}S_{\vect 0α\leftarrow\vect bβ}e^{i\vect k\cdot\vect r_{\vect b}}.\label{eq:W definition}
+W_{\alpha\beta}(\vect k)\equiv\sum_{\vect m\in\ints^{d}}\left(1-\delta_{\alpha\beta}\right)\tropsp{\vect 0,\alpha}{\vect m,\beta}e^{i\vect k\cdot\vect R_{\vect m}},\label{eq:W definition}
 \end{equation}
 
 \end_inset
 
+evaluation of which is possible but quite non-trivial due to the infinite
+ lattice sum, so we explain it separately in Sect.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "subsec:W operator evaluation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+As in the case of a finite system, eq.
+ can be written in a shorter block-matrix form,
+\begin_inset Formula 
+\begin{equation}
+\left(I-WT\right)\outcoeffp{\vect 0}\left(\vect k\right)=\rcoeffincp{\vect 0}\left(\vect k\right)\label{eq:Multiple-scattering problem unit cell block form}
+\end{equation}
+
+\end_inset
+
+ Eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem unit cell"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be used to calculate electromagnetic response of the structure to external
+ quasiperiodic driving field – most notably a plane wave.
+ However, the non-trivial solutions of the equation with right hand side
+ (i.e.
+ the external driving) set to zero, 
+\begin_inset Formula 
+\begin{equation}
+\left(I-WT\right)\outcoeffp{\vect 0}\left(\vect k\right)=0,\label{eq:lattice mode equation}
+\end{equation}
+
+\end_inset
+
+describes the 
+\emph on
+lattice modes.
+
+\emph default
+ Non-trivial solutions to 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:lattice mode equation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ exist if the matrix on the left-hand side 
+\begin_inset Formula $M\left(\omega,\vect k\right)=\left(I-W\left(\omega,\vect k\right)T\left(\omega\right)\right)$
+\end_inset
+
+ is singular – this condition gives the 
+\emph on
+dispersion relation
+\emph default
+ for the periodic structure.
+ Note that in realistic (lossy) systems, at least one of the pair 
+\begin_inset Formula $\omega,\vect k$
+\end_inset
+
+ will acquire complex values.
+ The solution 
+\begin_inset Formula $\outcoeffp{\vect 0}\left(\vect k\right)$
+\end_inset
+
+ is then obtained as the right 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+CHECK!
+\end_layout
+
+\end_inset
+
+ singular vector of 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ corresponding to the zero singular value.
+\end_layout
+
+\begin_layout Subsection
+Numerical solution
+\end_layout
+
+\begin_layout Standard
+In practice, equation 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem unit cell block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ is solved in the same way as eq.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ in the multipole degree truncated form.
+\end_layout
+
+\begin_layout Standard
+The lattice mode problem 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:lattice mode equation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ is (after multipole degree truncation) solved by finding 
+\begin_inset Formula $\omega,\vect k$
+\end_inset
+
+ for which the matrix 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ has a zero singular value.
+ A naïve approach to do that is to sample a volume with a grid in the 
+\begin_inset Formula $\left(\omega,\vect k\right)$
+\end_inset
+
+ space, performing a singular value decomposition of 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ at each point and finding where the lowest singular value of 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ is close enough to zero.
+ However, this approach is quite expensive, for 
+\begin_inset Formula $W\left(\omega,\vect k\right)$
+\end_inset
+
+ has to be evaluated for each 
+\begin_inset Formula $\omega,\vect k$
+\end_inset
+
+ pair separately (unlike the original finite case 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ translation operator 
+\begin_inset Formula $\trops$
+\end_inset
+
+, which, for a given geometry, depends only on frequency).
+ Therefore, a much more efficient approach to determine the photonic bands
+ is to sample the 
+\begin_inset Formula $\vect k$
+\end_inset
+
+-space (a whole Brillouin zone or its part) and for each fixed 
+\begin_inset Formula $\vect k$
+\end_inset
+
+ to find a corresponding frequency 
+\begin_inset Formula $\omega$
+\end_inset
+
+ with zero singular value of 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ using a minimisation algorithm (two- or one-dimensional, depending on whether
+ one needs the exact complex-valued 
+\begin_inset Formula $\omega$
+\end_inset
+
+ or whether the its real-valued approximation is satisfactory).
+ Typically, a good initial guess for 
+\begin_inset Formula $\omega\left(\vect k\right)$
+\end_inset
+
+ is obtained from the empty lattice approximation, 
+\begin_inset Formula $\left|\vect k\right|=\sqrt{\epsilon\mu}\omega/c_{0}$
+\end_inset
+
+ (modulo lattice points; TODO write this a clean way).
+ A somehow challenging step is to distinguish the different bands that can
+ all be very close to the empty lattice approximation, especially if the
+ particles in the systems are small.
+ In high-symmetry points of the Brilloin zone, this can be solved by factorising
+ 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ into irreducible representations 
+\begin_inset Formula $\Gamma_{i}$
+\end_inset
+
+ and performing the minimisation in each irrep separately, cf.
+ Section 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Symmetries"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, and using the different 
+\begin_inset Formula $\omega_{\Gamma_{i}}\left(\vect k\right)$
+\end_inset
+
+ to obtain the initial guesses for the nearby points 
+\begin_inset Formula $\vect k+\delta\vect k$
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Subsection
+Computing the Fourier sum of the translation operator
+\begin_inset CommandInset label
+LatexCommand label
+name "subsec:W operator evaluation"
+
+\end_inset
+
 
 \end_layout
 
-\begin_layout Subsection
-
-\lang english
-Computing the Fourier sum of the translation operator
-\end_layout
-
 \begin_layout Standard
-
-\lang english
 The problem evaluating 
 \begin_inset CommandInset ref
 LatexCommand eqref
@@ -233,20 +574,22 @@ reference "eq:W definition"
 \end_inset
 
  is the asymptotic behaviour of the translation operator, 
-\begin_inset Formula $S_{\vect 0α\leftarrow\vect bβ}\sim\left|\vect r_{\vect b}\right|^{-1}e^{ik_{0}\left|\vect r_{\vect b}\right|}$
+\begin_inset Formula $\tropsp{\vect 0,\alpha}{\vect m,\beta}\sim\left|\vect R_{\vect b}\right|^{-1}e^{ik_{0}\left|\vect R_{\vect b}\right|}$
 \end_inset
 
- that makes the convergence of the sum quite problematic for any 
+ that does not in the strict sense converge for any 
 \begin_inset Formula $d>1$
 \end_inset
 
 -dimensional lattice.
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
 \begin_inset Foot
 status open
 
 \begin_layout Plain Layout
-
-\lang english
 Note that 
 \begin_inset Formula $d$
 \end_inset
@@ -258,16 +601,20 @@ Note that
 
 \end_inset
 
- In electrostatics, one can solve this problem with Ewald summation.
+
+\end_layout
+
+\end_inset
+
+ In electrostatics, this problem can be solved with Ewald summation [TODO
+ REF].
  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.
+ We use the same idea here, but the technical details are more complicated
+ than in electrostatics.
 \end_layout
 
 \begin_layout Standard
-
-\lang english
 Let us re-express the sum in 
 \begin_inset CommandInset ref
 LatexCommand eqref
@@ -276,11 +623,14 @@ reference "eq:W definition"
 \end_inset
 
  in terms of integral with a delta comb
+\begin_inset FormulaMacro
+\renewcommand{\basis}[1]{\mathfrak{#1}}
+\end_inset
+
+
 \end_layout
 
 \begin_layout Standard
-
-\lang english
 \begin_inset Formula 
 \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}
@@ -382,8 +732,6 @@ W_{\alpha\beta}(\vect k) & = & \left(\left(\uaft{\dc{\basis u}}\right)\ast\left(
 status open
 
 \begin_layout Plain Layout
-
-\lang english
 Factor 
 \begin_inset Formula $\left(2\pi\right)^{\frac{d}{2}}$
 \end_inset
@@ -416,8 +764,6 @@ whole
 \end_layout
 
 \begin_layout Standard
-
-\lang english
 However, Fourier transform is linear, so we can in principle separate 
 \begin_inset Formula $S$
 \end_inset
@@ -463,7 +809,7 @@ reference "eq:W definition"
 
 \end_inset
 
- and legendre
+ and 
 \begin_inset CommandInset ref
 LatexCommand eqref
 reference "eq:W sum in reciprocal space"
@@ -484,7 +830,178 @@ W_{\alpha\beta}^{\textup{L}}\left(\vect k\right) & = & \frac{\left|\det\rec{\bas
 
 \end_inset
 
-where both sums should converge nicely.
+where both sums expected to converge nicely.
+ We note that the elements of the translation operators 
+\begin_inset Formula $\tropr,\trops$
+\end_inset
+
+ in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:translation operator"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be rewritten as linear combinations of expressions 
+\begin_inset Formula $\ush{\nu}{\mu}\left(\uvec d\right)j_{n}\left(d\right),\ush{\nu}{\mu}\left(\uvec d\right)h_{n}^{(1)}\left(d\right)$
+\end_inset
+
+ (TODO WRITE THEM EXPLICITLY IN THIS FORM), respectively, hence if we are
+ able evaluate the lattice sums sums
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+CHECK THE FOLLOWING EXPRESSION FOR CORRECT FUNCTION ARGUMENTS
+\end_layout
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{equation}
+\sigma_{\nu}^{\mu}\left(\vect k\right)=\sum_{\vect n\in\ints^{d}\backslash\left\{ \vect 0\right\} }e^{i\vect{\vect k}\cdot\vect R_{\vect n}}\ush{\nu}{\mu}\left(\uvec{R_{n}}\right)h_{n}^{(1)}\left(R_{n}\right),\label{eq:sigma lattice sums}
+\end{equation}
+
+\end_inset
+
+then by linearity, we can get the 
+\begin_inset Formula $W_{\alpha\beta}\left(\vect k\right)$
+\end_inset
+
+ operator as well.
+\end_layout
+
+\begin_layout Standard
+TODO ADD MOROZ AND OTHER REFS HERE.
+ 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "linton_one-_2009"
+literal "true"
+
+\end_inset
+
+ offers an exponentially convergent Ewald-type summation method for 
+\begin_inset Formula $\sigma_{\nu}^{\mu}\left(\vect k\right)=\sigma_{\nu}^{\mu(\mathrm{S})}\left(\vect k\right)+\sigma_{\nu}^{\mu(\mathrm{L})}\left(\vect k\right)$
+\end_inset
+
+.
+ Here we rewrite them in a form independent on the convention used for spherical
+ harmonics (as long as they are complex
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+lepší formulace
+\end_layout
+
+\end_inset
+
+).
+ The short-range part reads (UNIFY INDEX NOTATION)
+\begin_inset Formula 
+\begin{multline}
+\sigma_{n}^{m(\mathrm{S})}\left(\vect{\beta}\right)=-\frac{2^{n+1}i}{k^{n+1}\sqrt{\pi}}\sum_{\vect R\in\Lambda}^{'}\left|\vect R\right|^{n}e^{i\vect{\beta}\cdot\vect R}Y_{n}^{m}\left(\vect R\right)\int_{\eta}^{\infty}e^{-\left|\vect R\right|^{2}\xi^{2}}e^{-k/4\xi^{2}}\xi^{2n}\ud\xi\\
++\frac{\delta_{n0}\delta_{m0}}{\sqrt{4\pi}}\Gamma\left(-\frac{1}{2},-\frac{k}{4\eta^{2}}\right)Y_{n}^{m},\label{eq:Ewald in 3D short-range part}
+\end{multline}
+
+\end_inset
+
+
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+NEPATŘÍ TAM NĚJAKÁ DELTA FUNKCE K PŮVODNÍMU 
+\begin_inset Formula $\sigma_{n}^{m(0)}$
+\end_inset
+
+?
+\end_layout
+
+\end_inset
+
+and the long-range part (FIXME, this is the 2D version; include the 1D and
+ 3D lattice expressions as well)
+\begin_inset Formula 
+\begin{multline}
+\sigma_{n}^{m(\mathrm{L})}\left(\vect{\beta}\right)=-\frac{i^{n+1}}{k^{2}\mathscr{A}}\sqrt{\pi}2\left(\left(n-m\right)/2\right)!\left(\left(n+m\right)/2\right)!\times\\
+\times\sum_{\vect K_{pq}\in\Lambda^{*}}^{'}Y_{n}^{m}\left(\frac{\pi}{2},\phi_{\vect{\beta}_{pq}}\right)\sum_{j=0}^{\left[\left(n-\left|m\right|/2\right)\right]}\frac{\left(-1\right)^{j}\left(\beta_{pq}/k\right)^{n-2j}\Gamma_{j,pq}}{j!\left(\frac{1}{2}\left(n-m\right)-j\right)!\left(\frac{1}{2}\left(n+m\right)-j\right)!}\left(\gamma_{pq}\right)^{2j-1}\label{eq:Ewald in 3D long-range part}
+\end{multline}
+
+\end_inset
+
+where 
+\begin_inset Formula $\xi$
+\end_inset
+
+ is TODO, 
+\begin_inset Formula $\beta_{pq}$
+\end_inset
+
+ is TODO, 
+\begin_inset Formula $\Gamma_{j,pq}$
+\end_inset
+
+ is TODO and 
+\begin_inset Formula $\eta$
+\end_inset
+
+ is a real parameter that determines the pace of convergence of both parts.
+ The larger 
+\begin_inset Formula $\eta$
+\end_inset
+
+ is, the faster 
+\begin_inset Formula $\sigma_{n}^{m(\mathrm{S})}\left(\vect{\beta}\right)$
+\end_inset
+
+ converges but the slower 
+\begin_inset Formula $\sigma_{n}^{m(\mathrm{L})}\left(\vect{\beta}\right)$
+\end_inset
+
+ converges.
+ Therefore (based on the lattice geometry) it has to be adjusted in a way
+ that a reasonable amount of terms needs to be evaluated numerically from
+ both 
+\begin_inset Formula $\sigma_{n}^{m(\mathrm{S})}\left(\vect{\beta}\right)$
+\end_inset
+
+ and 
+\begin_inset Formula $\sigma_{n}^{m(\mathrm{L})}\left(\vect{\beta}\right)$
+\end_inset
+
+ .
+ Generally, a good choice for 
+\begin_inset Formula $\eta$
+\end_inset
+
+ is TODO; in order to achieve accuracy TODO, one has to evaluate the terms
+ on TODO lattice points.
+ (I HAVE SOME DERIVATIONS OF THE ESTIMATES IN MY NOTES; SHOULD I INCLUDE
+ THEM?)
+\end_layout
+
+\begin_layout Standard
+In practice, the integrals in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Ewald in 3D short-range part"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be easily evaluated by numerical quadrature and the incomplete 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+-functions using the series TODO and TODO from DLMF.
 \end_layout
 
 \end_body
diff --git a/lepaper/intro.lyx b/lepaper/intro.lyx
index f814a7c..e9f98f0 100644
--- a/lepaper/intro.lyx
+++ b/lepaper/intro.lyx
@@ -1,29 +1,37 @@
-#LyX 2.1 created this file. For more info see http://www.lyx.org/
-\lyxformat 474
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 583
 \begin_document
 \begin_header
+\save_transient_properties true
+\origin unavailable
 \textclass article
 \use_default_options true
 \maintain_unincluded_children false
-\language finnish
+\language english
 \language_package default
-\inputencoding auto
-\fontencoding global
-\font_roman TeX Gyre Pagella
-\font_sans default
-\font_typewriter default
-\font_math auto
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "TeX Gyre Pagella"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
 \font_default_family default
-\use_non_tex_fonts true
+\use_non_tex_fonts false
 \font_sc false
-\font_osf true
-\font_sf_scale 100
-\font_tt_scale 100
+\font_roman_osf true
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures false
 \graphics default
-\default_output_format pdf4
+\default_output_format default
 \output_sync 0
 \bibtex_command default
 \index_command default
+\float_placement class
+\float_alignment class
 \paperfontsize default
 \spacing single
 \use_hyperref true
@@ -58,6 +66,8 @@
 \suppress_date false
 \justification true
 \use_refstyle 1
+\use_minted 0
+\use_lineno 0
 \index Index
 \shortcut idx
 \color #008000
@@ -66,10 +76,14 @@
 \tocdepth 3
 \paragraph_separation indent
 \paragraph_indentation default
-\quotes_language swedish
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
 \papercolumns 1
 \papersides 1
 \paperpagestyle default
+\tablestyle default
 \tracking_changes false
 \output_changes false
 \html_math_output 0
@@ -81,6 +95,159 @@
 
 \begin_layout Section
 Introduction
+\begin_inset CommandInset label
+LatexCommand label
+name "sec:Introduction"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+The problem of electromagnetic response of a system consisting of many compact
+ scatterers in various geometries, and its numerical solution, is relevant
+ to many branches of nanophotonics (TODO refs).
+ The most commonly used general approaches used in computational electrodynamics
+, such as the finite difference time domain (FDTD) method or the finite
+ element method (FEM), are very often unsuitable for simulating systems
+ with larger number of scatterers due to their computational complexity.
+ Therefore, a common (frequency-domain) approach to get an approximate solution
+ of the scattering problem for many small particles has been the coupled
+ dipole approximation (CDA) where individual scatterers are reduced to electric
+ dipoles (characterised by a polarisability tensor) and coupled to each
+ other through Green's functions.
+ 
+\end_layout
+
+\begin_layout Standard
+CDA is easy to implement and has favorable computational complexity but
+ suffers from at least two fundamental drawbacks.
+ The obvious one is that the dipole approximation is too rough for particles
+ with diameter larger than a small fraction of the wavelength.
+ The other one, more subtle, manifests itself in photonic crystal-like structure
+s used in nanophotonics: there are modes in which the particles' electric
+ dipole moments completely vanish due to symmetry, regardless of how small
+ the particles are, and the excitations have quadrupolar or higher-degree
+ multipolar character.
+ These modes typically appear at the band edges where interesting phenomena
+ such as lasing or Bose-Einstein condensation have been observed – and CDA
+ by definition fails to capture such modes.
+\end_layout
+
+\begin_layout Standard
+The natural way to overcome both limitations of CDA mentioned above is to
+ include higher multipoles into account.
+ Instead of polarisability tensor, the scattering properties of an individual
+ particle are then described a more general 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix, and different particles' multipole excitations are coupled together
+ via translation operators, a generalisation of the Green's functions in
+ CDA.
+ This is the idea behind the 
+\emph on
+multiple-scattering 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix method 
+\emph default
+(MSTMM) (TODO a.k.a something??), and it has been implemented previously for
+ a limited subset of problems (TODO refs and list the limitations of the
+ available).
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO přestože blablaba, moc se to nepoužívalo, protože je težké udělat to
+ správně.
+\end_layout
+
+\end_inset
+
+ Due to the limitations of the existing available codes, we have been developing
+ our own implementation of MSTMM, which we have used in several previous
+ works studying various physical phenomena in plasmonic nanoarrays (TODO
+ examples with refs).
+ 
+\end_layout
+
+\begin_layout Standard
+Hereby we release our MSTMM implementation, the 
+\emph on
+QPMS Photonic Multiple Scattering
+\emph default
+ suite, as an open source software under the GNU General Public License
+ version 3.
+ (TODO refs to the code repositories.) QPMS allows for linear optics simulations
+ of arbitrary sets of compact scatterers in isotropic media.
+ The features include computations of electromagnetic response to external
+ driving, the related cross sections, and finding resonances of finite structure
+s.
+ Moreover, in QPMS we extensively employ group theory to exploit the physical
+ symmetries of the system to further reduce the demands on computational
+ resources, enabling to simulate even larger systems.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+(TODO put a specific example here of how large system we are able to simulate?)
+\end_layout
+
+\end_inset
+
+ Although systems of large 
+\emph on
+finite
+\emph default
+ number of scatterers are the area where MSTMM excels the most—simply because
+ other methods fail due to their computational complexity—we also extended
+ the method onto infinite periodic systems (photonic crystals); this can
+ be used for quickly evaluating dispersions of such structures and also
+ their topological invariants (TODO).
+ The QPMS suite contains a core C library, Python bindings and several utilities
+ for routine computations, such as TODO.
+ It includes extensive Doxygen documentation, together with description
+ of the API, making extending and customising the code easy.
+\end_layout
+
+\begin_layout Standard
+The current paper is organised as follows: Section 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Finite"
+
+\end_inset
+
+ is devoted to MSTMM theory for finite systems, in Section 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Infinite"
+
+\end_inset
+
+ we develop the theory for infinite periodic structures.
+ Section 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Applications"
+
+\end_inset
+
+ demonstrates some basic practical results that can be obtained using QPMS.
+ Finally, in Section 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sec:Comparison"
+
+\end_inset
+
+ we comment on the computational complexity of MSTMM in comparison to other
+ methods.
 \end_layout
 
 \end_body
diff --git a/lepaper/symmetries.lyx b/lepaper/symmetries.lyx
new file mode 100644
index 0000000..9016d58
--- /dev/null
+++ b/lepaper/symmetries.lyx
@@ -0,0 +1,1440 @@
+#LyX 2.4 created this file. For more info see https://www.lyx.org/
+\lyxformat 584
+\begin_document
+\begin_header
+\save_transient_properties true
+\origin unavailable
+\textclass article
+\use_default_options true
+\maintain_unincluded_children false
+\language finnish
+\language_package default
+\inputencoding utf8
+\fontencoding auto
+\font_roman "default" "default"
+\font_sans "default" "default"
+\font_typewriter "default" "default"
+\font_math "auto" "auto"
+\font_default_family default
+\use_non_tex_fonts false
+\font_sc false
+\font_roman_osf false
+\font_sans_osf false
+\font_typewriter_osf false
+\font_sf_scale 100 100
+\font_tt_scale 100 100
+\use_microtype false
+\use_dash_ligatures true
+\graphics default
+\default_output_format default
+\output_sync 0
+\bibtex_command default
+\index_command default
+\paperfontsize default
+\use_hyperref false
+\papersize default
+\use_geometry false
+\use_package amsmath 1
+\use_package amssymb 1
+\use_package cancel 1
+\use_package esint 1
+\use_package mathdots 1
+\use_package mathtools 1
+\use_package mhchem 1
+\use_package stackrel 1
+\use_package stmaryrd 1
+\use_package undertilde 1
+\cite_engine basic
+\cite_engine_type default
+\use_bibtopic false
+\use_indices false
+\paperorientation portrait
+\suppress_date false
+\justification true
+\use_refstyle 1
+\use_minted 0
+\use_lineno 0
+\index Index
+\shortcut idx
+\color #008000
+\end_index
+\secnumdepth 3
+\tocdepth 3
+\paragraph_separation indent
+\paragraph_indentation default
+\is_math_indent 0
+\math_numbering_side default
+\quotes_style english
+\dynamic_quotes 0
+\papercolumns 1
+\papersides 1
+\paperpagestyle default
+\tablestyle default
+\tracking_changes false
+\output_changes false
+\html_math_output 0
+\html_css_as_file 0
+\html_be_strict false
+\end_header
+
+\begin_body
+
+\begin_layout Section
+Symmetries
+\begin_inset CommandInset label
+LatexCommand label
+name "sec:Symmetries"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+If the system has nontrivial point group symmetries, group theory gives
+ additional understanding of the system properties, and can be used to reduce
+ the computational costs.
+ 
+\end_layout
+
+\begin_layout Standard
+As an example, if our system has a 
+\begin_inset Formula $D_{2h}$
+\end_inset
+
+ symmetry and our truncated 
+\begin_inset Formula $\left(I-T\trops\right)$
+\end_inset
+
+ matrix has size 
+\begin_inset Formula $N\times N$
+\end_inset
+
+,
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+nepoužívám 
+\begin_inset Formula $N$
+\end_inset
+
+ už v jiném kontextu?
+\end_layout
+
+\end_inset
+
+ it can be block-diagonalized into eight blocks of size about 
+\begin_inset Formula $N/8\times N/8$
+\end_inset
+
+, each of which can be LU-factorised separately (this is due to the fact
+ that 
+\begin_inset Formula $D_{2h}$
+\end_inset
+
+ has eight different one-dimensional irreducible representations).
+ This can reduce both memory and time requirements to solve the scattering
+ problem 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ by a factor of 64.
+\end_layout
+
+\begin_layout Standard
+In periodic systems (problems 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem unit cell block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:lattice mode equation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+) due to small number of particles per unit cell, the costliest part is
+ usually the evaluation of the lattice sums in the 
+\begin_inset Formula $W\left(\omega,\vect k\right)$
+\end_inset
+
+ matrix, not the linear algebra.
+ However, the lattice modes can be searched for in each irrep separately,
+ and the irrep dimension gives a priori information about mode degeneracy.
+\end_layout
+
+\begin_layout Subsection
+Excitation coefficients under point group operations
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO Zkontrolovat všechny vzorečky zde!!!
+\end_layout
+
+\end_inset
+
+In order to use the point group symmetries, we first need to know how they
+ affect our basis functions, i.e.
+ the VSWFs.
+\end_layout
+
+\begin_layout Standard
+Let 
+\begin_inset Formula $g$
+\end_inset
+
+ be a member of orthogonal group 
+\begin_inset Formula $O(3)$
+\end_inset
+
+, i.e.
+ a 3D point rotation or reflection operation that transforms vectors in
+ 
+\begin_inset Formula $\reals^{3}$
+\end_inset
+
+ with an orthogonal matrix 
+\begin_inset Formula $R_{g}$
+\end_inset
+
+:
+\begin_inset Formula 
+\[
+\vect r\mapsto R_{g}\vect r.
+\]
+
+\end_inset
+
+Spherical harmonics 
+\begin_inset Formula $\ush lm$
+\end_inset
+
+, being a basis the 
+\begin_inset Formula $l$
+\end_inset
+
+-dimensional representation of 
+\begin_inset Formula $O(3)$
+\end_inset
+
+, transform as 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "???"
+key "dresselhaus_group_2008"
+literal "false"
+
+\end_inset
+
+
+\begin_inset Formula 
+\[
+\ush lm\left(R_{g}\uvec r\right)=\sum_{m'=-l}^{l}D_{m,m'}^{l}\left(g\right)\ush l{m'}\left(\uvec r\right)
+\]
+
+\end_inset
+
+where 
+\begin_inset Formula $D_{m,m'}^{l}\left(g\right)$
+\end_inset
+
+ denotes the elements of the 
+\emph on
+Wigner matrix
+\emph default
+ representing the operation 
+\begin_inset Formula $g$
+\end_inset
+
+.
+ By their definition, vector spherical harmonics 
+\begin_inset Formula $\vsh 2lm,\vsh 3lm$
+\end_inset
+
+ transform in the same way,
+\begin_inset Formula 
+\begin{align*}
+\vsh 2lm\left(R_{g}\uvec r\right) & =\sum_{m'=-l}^{l}D_{m,m'}^{l}\left(g\right)\vsh 2l{m'}\left(\uvec r\right),\\
+\vsh 3lm\left(R_{g}\uvec r\right) & =\sum_{m'=-l}^{l}D_{m,m'}^{l}\left(g\right)\vsh 3l{m'}\left(\uvec r\right),
+\end{align*}
+
+\end_inset
+
+but the remaining set 
+\begin_inset Formula $\vsh 1lm$
+\end_inset
+
+ transforms differently due to their pseudovector nature stemming from the
+ cross product in their definition:
+\begin_inset Formula 
+\[
+\vsh 3lm\left(R_{g}\uvec r\right)=\sum_{m'=-l}^{l}\widetilde{D_{m,m'}^{l}}\left(g\right)\vsh 3l{m'}\left(\uvec r\right),
+\]
+
+\end_inset
+
+where 
+\begin_inset Formula $\widetilde{D_{m,m'}^{l}}\left(g\right)=D_{m,m'}^{l}\left(g\right)$
+\end_inset
+
+ if 
+\begin_inset Formula $g$
+\end_inset
+
+ is a proper rotation, but for spatial inversion operation 
+\begin_inset Formula $i:\vect r\mapsto-\vect r$
+\end_inset
+
+ we have 
+\begin_inset Formula $\widetilde{D_{m,m'}^{l}}\left(i\right)=\left(-1\right)^{l+m}D_{m,m'}^{l}\left(i\right)$
+\end_inset
+
+.
+ The transformation behaviour of vector spherical harmonics directly propagates
+ to the spherical vector waves, cf.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:VSWF regular"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:VSWF outgoing"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+:
+\begin_inset Formula 
+\begin{align*}
+\vswfouttlm 1lm\left(R_{g}\vect r\right) & =\sum_{m'=-l}^{l}\widetilde{D_{m,m'}^{l}}\left(g\right)\vswfouttlm 1l{m'}\left(\vect r\right),\\
+\vswfouttlm 2lm\left(R_{g}\vect r\right) & =\sum_{m'=-l}^{l}D_{m,m'}^{l}\left(g\right)\vswfouttlm 2l{m'}\left(\vect r\right),
+\end{align*}
+
+\end_inset
+
+(and analogously for the regular waves 
+\begin_inset Formula $\vswfrtlm{\tau}lm$
+\end_inset
+
+).
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO víc obdivu.
+\end_layout
+
+\end_inset
+
+ For convenience, we introduce the symbol 
+\begin_inset Formula $D_{m,m'}^{\tau l}$
+\end_inset
+
+ that describes the transformation of both types (
+\begin_inset Quotes eld
+\end_inset
+
+magnetic
+\begin_inset Quotes erd
+\end_inset
+
+ and 
+\begin_inset Quotes eld
+\end_inset
+
+electric
+\begin_inset Quotes erd
+\end_inset
+
+) of waves at once:
+\begin_inset Formula 
+\[
+\vswfouttlm{\tau}lm\left(R_{g}\vect r\right)=\sum_{m'=-l}^{l}D_{m,m'}^{\tau l}\left(g\right)\vswfouttlm{\tau}l{m'}\left(\vect r\right).
+\]
+
+\end_inset
+
+Using these, we can express the VSWF expansion 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:E field expansion"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ of the electric field around origin in a rotated/reflected system,
+\begin_inset Formula 
+\[
+\vect E\left(\omega,R_{g}\vect r\right)=\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\sum_{m'=-l}^{l}\left(\rcoefftlm{\tau}lmD_{m,m'}^{\tau l}\left(g\right)\vswfrtlm{\tau}l{m'}\left(k\vect r\right)+\outcoefftlm{\tau}lmD_{m,m'}^{\tau l}\left(g\right)\vswfouttlm{\tau}l{m'}\left(k\vect r\right)\right),
+\]
+
+\end_inset
+
+which, together with the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix definition, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:T-matrix definition"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be used to obtain a 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix of a rotated or mirror-reflected particle.
+ Let 
+\begin_inset Formula $T$
+\end_inset
+
+ be the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix of an original particle; the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix of a particle physically transformed by operation 
+\begin_inset Formula $g\in O(3)$
+\end_inset
+
+ is then 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+check sides
+\end_layout
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{equation}
+T'_{\tau lm;\tau'l'm'}=\sum_{\mu=-l}^{l}\sum_{\mu'=-l'}^{l'}\left(D_{\mu,m}^{\tau l}\left(g\right)\right)^{*}T_{\tau l\mu;\tau'l'm'}D_{m',\mu'}^{\tau l}\left(g\right).\label{eq:T-matrix of a transformed particle}
+\end{equation}
+
+\end_inset
+
+If the particle is symmetric (so that 
+\begin_inset Formula $g$
+\end_inset
+
+ produces a particle indistinguishable from the original one), the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix must remain invariant under the transformation 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:T-matrix of a transformed particle"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset Formula $T'_{\tau lm;\tau'l'm'}=T{}_{\tau lm;\tau'l'm'}$
+\end_inset
+
+.
+ Explicit forms of these invariance properties for the most imporant point
+ group symmetries can be found in 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "schulz_point-group_1999"
+literal "false"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+If the field expansion is done around a point 
+\begin_inset Formula $\vect r_{p}$
+\end_inset
+
+ different from the global origin, as in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:E field expansion multiparticle"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, we have
+\lang english
+
+\begin_inset Formula 
+\begin{align}
+\vect E\left(\omega,R_{g}\vect r\right) & =\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\sum_{m'=-l}^{l}\left(\rcoeffptlm p{\tau}lmD_{m',\mu'}^{\tau l}\left(g\right)\vswfrtlm{\tau}l{m'}\left(k\left(\vect r-R_{g}\vect r_{p}\right)\right)+\outcoeffptlm p{\tau}lmD_{m',\mu'}^{\tau l}\left(g\right)\vswfouttlm{\tau}l{m'}\left(k\left(\vect r-R_{g}\vect r_{p}\right)\right)\right).\label{eq:rotated E field expansion around outside origin}
+\end{align}
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Float figure
+placement document
+alignment document
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Caption Standard
+
+\begin_layout Plain Layout
+Scatterer orbits under 
+\begin_inset Formula $D_{2}$
+\end_inset
+
+ symmetry.
+ Particles 
+\begin_inset Formula $A,B,C,D$
+\end_inset
+
+ lie outside of origin or any mirror planes, and together constitute an
+ orbit of the size equal to the order of the group, 
+\begin_inset Formula $\left|D_{2}\right|=4$
+\end_inset
+
+.
+ Particles 
+\begin_inset Formula $E,F$
+\end_inset
+
+ lie on the 
+\begin_inset Formula $xz$
+\end_inset
+
+ plane, hence the corresponding reflection maps each of them to itself,
+ but the 
+\begin_inset Formula $yz$
+\end_inset
+
+ reflection (or the 
+\begin_inset Formula $\pi$
+\end_inset
+
+ rotation around the 
+\begin_inset Formula $z$
+\end_inset
+
+ axis) maps them to each other, forming a particle orbit of size 2
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+=???
+\end_layout
+
+\end_inset
+
+.
+ The particle 
+\begin_inset Formula $O$
+\end_inset
+
+ in the very origin is always mapped to itself, constituting its own orbit.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:D2-symmetric structure particle orbits"
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+TODO restructure this
+\end_layout
+
+\end_inset
+
+With these transformation properties in hand, we can proceed to the effects
+ of point symmetries on the whole many-particle system.
+ Let us have a many-particle system symmetric with respect to a point group
+ 
+\begin_inset Formula $G$
+\end_inset
+
+.
+ A symmetry operation 
+\begin_inset Formula $g\in G$
+\end_inset
+
+ determines a permutation of the particles: 
+\begin_inset Formula $p\mapsto\pi_{g}(p)$
+\end_inset
+
+, 
+\begin_inset Formula $p\in\mathcal{P}$
+\end_inset
+
+.
+ For a given particle 
+\begin_inset Formula $p$
+\end_inset
+
+, we will call the set of particles onto which any of the symmetries maps
+ the particle 
+\begin_inset Formula $p$
+\end_inset
+
+, i.e.
+ the set 
+\begin_inset Formula $\left\{ \pi_{g}\left(p\right);g\in G\right\} $
+\end_inset
+
+, as the 
+\emph on
+orbit
+\emph default
+ of particle 
+\begin_inset Formula $p$
+\end_inset
+
+.
+ The whole set 
+\begin_inset Formula $\mathcal{P}$
+\end_inset
+
+ can therefore be divided into the different particle orbits; an example
+ is in Fig.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:D2-symmetric structure particle orbits"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+ The importance of the particle orbits stems from the following: in the
+ multiple-scattering problem, outside of the scatterers 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+< FIXME
+\end_layout
+
+\end_inset
+
+ one has 
+\lang english
+
+\begin_inset Formula 
+\begin{align}
+\vect E\left(\omega,R_{g}\vect r\right) & =\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\sum_{m'=-l}^{l}\left(\rcoeffptlm p{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right)\vswfrtlm{\tau}l{m'}\left(k\left(\vect r-\vect r_{\pi_{g}(p)}\right)\right)+\outcoeffptlm p{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right)\vswfouttlm{\tau}l{m'}\left(k\left(\vect r-\vect r_{p}\right)\right)\right)\label{eq:rotated E field expansion around outside origin-1}\\
+ & =\sum_{\tau=1,2}\sum_{l=1}^{\infty}\sum_{m=-l}^{+l}\sum_{m'=-l}^{l}\left(\rcoeffptlm{\pi_{g}^{-1}(p)}{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right)\vswfrtlm{\tau}l{m'}\left(k\left(\vect r-\vect r_{p}\right)\right)+\outcoeffptlm{\pi_{g}^{-1}(p)}{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right)\vswfouttlm{\tau}l{m'}\left(k\left(\vect r-\vect r_{p}\right)\right)\right).
+\end{align}
+
+\end_inset
+
+This means that the field expansion coefficients 
+\begin_inset Formula $\rcoeffp p,\outcoeffp p$
+\end_inset
+
+ transform as 
+\begin_inset Formula 
+\begin{align}
+\rcoeffptlm p{\tau}lm & \mapsto\rcoeffptlm{\pi_{g}^{-1}(p)}{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right),\nonumber \\
+\outcoeffptlm p{\tau}lm & \mapsto\outcoeffptlm{\pi_{g}^{-1}(p)}{\tau}lmD_{m,\mu'}^{\tau l}\left(g\right).\label{eq:excitation coefficient under symmetry operation}
+\end{align}
+
+\end_inset
+
+Obviously, the expansion coefficients belonging to particles in different
+ orbits do not mix together.
+ As before, we introduce a short-hand block-matrix notation for 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "eq:excitation coefficient under symmetry operation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ (TODO avoid notation clash here in a more consistent and readable way!)
+\end_layout
+
+\begin_layout Standard
+
+\lang english
+\begin_inset Formula 
+\begin{align}
+\rcoeff & \mapsto J\left(g\right)a,\nonumber \\
+\outcoeff & \mapsto J\left(g\right)\outcoeff.\label{eq:excitation coefficient under symmetry operation block form}
+\end{align}
+
+\end_inset
+
+
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+
+\lang english
+The matrices 
+\begin_inset Formula $D\left(g\right)$
+\end_inset
+
+, 
+\begin_inset Formula $g\in G$
+\end_inset
+
+ will play a crucial role blablabla
+\end_layout
+
+\end_inset
+
+If the particle indices are ordered in a way that the particles belonging
+ to the same orbit are grouped together, 
+\begin_inset Formula $J\left(g\right)$
+\end_inset
+
+ will be a block-diagonal unitary matrix, each block (also unitary) representing
+ the action of 
+\begin_inset Formula $g$
+\end_inset
+
+ on one particle orbit.
+ All the 
+\begin_inset Formula $J\left(g\right)$
+\end_inset
+
+s make together a (reducible) linear representation of 
+\begin_inset Formula $G$
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Subsection
+Irrep decomposition
+\end_layout
+
+\begin_layout Standard
+Knowledge of symmetry group actions 
+\begin_inset Formula $J\left(g\right)$
+\end_inset
+
+ on the field expansion coefficients give us the possibility to construct
+ a symmetry adapted basis in which we can block-diagonalise the multiple-scatter
+ing problem matrix 
+\begin_inset Formula $\left(I-TS\right)$
+\end_inset
+
+.
+ Let 
+\begin_inset Formula $\Gamma_{n}$
+\end_inset
+
+ be the 
+\begin_inset Formula $d_{n}$
+\end_inset
+
+-dimensional irreducible matrix representations of 
+\begin_inset Formula $G$
+\end_inset
+
+consisting of matrices 
+\begin_inset Formula $D^{\Gamma_{n}}\left(g\right)$
+\end_inset
+
+.
+ Then the projection operators
+\begin_inset Formula 
+\[
+P_{kl}^{\left(\Gamma_{n}\right)}\equiv\frac{d_{n}}{\left|G\right|}\sum_{g\in G}\left(D^{\Gamma_{n}}\left(g\right)\right)_{kl}^{*}J\left(g\right),\quad k,l=1,\dots,d_{n}
+\]
+
+\end_inset
+
+project the full scattering system field expansion coefficient vectors 
+\begin_inset Formula $\rcoeff,\outcoeff$
+\end_inset
+
+ onto a subspace corresponding to the irreducible representation 
+\begin_inset Formula $\Gamma_{n}$
+\end_inset
+
+.
+ The projectors can be used to construct a unitary transformation 
+\begin_inset Formula $U$
+\end_inset
+
+ with components
+\begin_inset Formula 
+\begin{equation}
+U_{nri;p\tau lm}=\frac{d_{n}}{\left|G\right|}\sum_{g\in G}\left(D^{\Gamma_{n}}\left(g\right)\right)_{rr}^{*}J\left(g\right)_{p'\tau'l'm'(nri);p\tau lm}\label{eq:SAB unitary transformation operator}
+\end{equation}
+
+\end_inset
+
+where 
+\begin_inset Formula $r$
+\end_inset
+
+ goes from 
+\begin_inset Formula $1$
+\end_inset
+
+ through 
+\begin_inset Formula $d_{n}$
+\end_inset
+
+ and 
+\begin_inset Formula $i$
+\end_inset
+
+ goes from 1 through the multiplicity of irreducible representation 
+\begin_inset Formula $\Gamma_{n}$
+\end_inset
+
+ in the (reducible) representation of 
+\begin_inset Formula $G$
+\end_inset
+
+ spanned by the field expansion coefficients 
+\begin_inset Formula $\rcoeff$
+\end_inset
+
+ or 
+\begin_inset Formula $\outcoeff$
+\end_inset
+
+.
+ The indices 
+\begin_inset Formula $p',\tau',l',m'$
+\end_inset
+
+ are given by an arbitrary bijective mapping 
+\begin_inset Formula $\left(n,r,i\right)\mapsto\left(p',\tau',l',m'\right)$
+\end_inset
+
+ with the constraint that for given 
+\begin_inset Formula $n,r,i$
+\end_inset
+
+ there are at least some non-zero elements 
+\begin_inset Formula $U_{nri;p\tau lm}$
+\end_inset
+
+.
+ For details, we refer the reader to textbooks about group representation
+ theory
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+or linear representations?
+\end_layout
+
+\end_inset
+
+, e.g.
+ 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "Chapter 4"
+key "dresselhaus_group_2008"
+literal "false"
+
+\end_inset
+
+ or 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "???"
+key "bradley_mathematical_1972"
+literal "false"
+
+\end_inset
+
+.
+ The transformation given by 
+\begin_inset Formula $U$
+\end_inset
+
+ transforms the excitation coefficient vectors 
+\begin_inset Formula $\rcoeff,\outcoeff$
+\end_inset
+
+ into a new, 
+\emph on
+symmetry-adapted basis
+\emph default
+.
+ 
+\end_layout
+
+\begin_layout Standard
+One can show that if an operator 
+\begin_inset Formula $M$
+\end_inset
+
+ acting on the excitation coefficient vectors is invariant under the operations
+ of group 
+\begin_inset Formula $G$
+\end_inset
+
+, meaning that
+\begin_inset Formula 
+\[
+\forall g\in G:J\left(g\right)MJ\left(g\right)^{\dagger}=M,
+\]
+
+\end_inset
+
+then in the symmetry-adapted basis, 
+\begin_inset Formula $M$
+\end_inset
+
+ is block diagonal, or more specifically
+\begin_inset Formula 
+\[
+M_{\Gamma,r,i;\Gamma',r',j}^{\mathrm{s.a.b.}}=\frac{\delta_{\Gamma\Gamma'}\delta_{ij}}{d_{\Gamma}}\sum_{q}M{}_{\Gamma,r,q;\Gamma',r',q}^{\mathrm{s.a.b.}}.
+\]
+
+\end_inset
+
+Both the 
+\begin_inset Formula $T$
+\end_inset
+
+ and 
+\begin_inset Formula $\trops$
+\end_inset
+
+ operators (and trivially also the identity 
+\begin_inset Formula $I$
+\end_inset
+
+) in 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ are invariant under the actions of whole system symmetry group, so 
+\begin_inset Formula $\left(I-T\trops\right)$
+\end_inset
+
+ is also invariant, hence 
+\begin_inset Formula $U\left(I-T\trops\right)U^{\dagger}$
+\end_inset
+
+ is a block-diagonal matrix, and the problem 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+ can be solved for each block separately.
+\end_layout
+
+\begin_layout Standard
+From the computational perspective, it is important to note that 
+\begin_inset Formula $U$
+\end_inset
+
+ is at least as sparse as 
+\begin_inset Formula $J\left(g\right)$
+\end_inset
+
+ (which is 
+\begin_inset Quotes eld
+\end_inset
+
+orbit-block
+\begin_inset Quotes erd
+\end_inset
+
+ diagonal), hence the block-diagonalisation can be performed fast.
+ 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+Kvantifikovat!
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsection
+Periodic systems
+\end_layout
+
+\begin_layout Standard
+For periodic systems, we can in similar manner also block-diagonalise the
+ 
+\begin_inset Formula $M\left(\omega,\vect k\right)=\left(I-W\left(\omega,\vect k\right)T\left(\omega\right)\right)$
+\end_inset
+
+ from the left hand side of eqs.
+ 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:Multiple-scattering problem unit cell block form"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+, 
+\begin_inset CommandInset ref
+LatexCommand eqref
+reference "eq:lattice mode equation"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+ Hovewer, in this case, 
+\begin_inset Formula $W\left(\omega,\vect k\right)$
+\end_inset
+
+ is in general not invariant under the whole point group symmetry subgroup
+ of the system geometry due to the 
+\begin_inset Formula $\vect k$
+\end_inset
+
+ dependence.
+ In other words, only those point symmetries that the 
+\begin_inset Formula $e^{i\vect k\cdot\vect r}$
+\end_inset
+
+ modulation does not break are preserved, and no preservation of point symmetrie
+s happens unless 
+\begin_inset Formula $\vect k$
+\end_inset
+
+ lies somewhere in the high-symmetry parts of the Brillouin zone.
+ However, the high-symmetry points are usually the ones of the highest physical
+ interest, for it is where the band edges 
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+or 
+\begin_inset Quotes eld
+\end_inset
+
+dirac points
+\begin_inset Quotes erd
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+ are typically located.
+\end_layout
+
+\begin_layout Standard
+The transformation to the symmetry adapted basis 
+\begin_inset Formula $U$
+\end_inset
+
+ is constructed in a similar way as in the finite case, but because we do
+ not work with all the (infinite number of) scatterers but only with one
+ unit cell, additional phase factors 
+\begin_inset Formula $e^{i\vect k\cdot\vect r_{p}}$
+\end_inset
+
+ appear in the per-unit-cell group action 
+\begin_inset Formula $J(g)$
+\end_inset
+
+.
+ This is illustrated in Fig.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "Phase factor illustration"
+plural "false"
+caps "false"
+noprefix "false"
+
+\end_inset
+
+.
+\begin_inset Float figure
+placement document
+alignment document
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+
+\end_layout
+
+\begin_layout Plain Layout
+\begin_inset Caption Standard
+
+\begin_layout Plain Layout
+\begin_inset CommandInset label
+LatexCommand label
+name "Phase factor illustration"
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Plain Layout
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+More rigorous analysis can be found e.g.
+ in 
+\lang english
+
+\begin_inset CommandInset citation
+LatexCommand cite
+after "chapters 10–11"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+In the group-theoretical terminology, blablabla little groups blabla bla...
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+
+\lang english
+\begin_inset Note Note
+status open
+
+\begin_layout Plain Layout
+
+\lang english
+A general overview of utilizing group theory to find lattice modes at high-symme
+try points of the Brillouin zone can be found e.g.
+ in 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "chapters 10–11"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+; here we use the same notation.
+\end_layout
+
+\begin_layout Plain Layout
+
+\lang english
+We analyse the symmetries of the system in the same VSWF representation
+ as used in the 
+\begin_inset Formula $T$
+\end_inset
+
+-matrix formalism introduced above.
+ We are interested in the modes at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point of the hexagonal lattice, which has the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ point symmetry.
+ The six irreducible representations (irreps) of the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ group are known and are available in the literature in their explicit forms.
+ In order to find and classify the modes, we need to find a decomposition
+ of the lattice mode representation 
+\begin_inset Formula $\Gamma_{\mathrm{lat.mod.}}=\Gamma^{\mathrm{equiv.}}\otimes\Gamma_{\mathrm{vec.}}$
+\end_inset
+
+ into the irreps of 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+.
+ The equivalence representation 
+\begin_inset Formula $\Gamma^{\mathrm{equiv.}}$
+\end_inset
+
+ is the 
+\begin_inset Formula $E'$
+\end_inset
+
+ representation as can be deduced from 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (11.19)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+, eq.
+ (11.19) and the character table for 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+.
+ 
+\begin_inset Formula $\Gamma_{\mathrm{vec.}}$
+\end_inset
+
+ operates on a space spanned by the VSWFs around each nanoparticle in the
+ unit cell (the effects of point group operations on VSWFs are described
+ in 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "schulz_point-group_1999"
+literal "true"
+
+\end_inset
+
+).
+ This space can be then decomposed into invariant subspaces of the 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+ using the projectors 
+\begin_inset Formula $\hat{P}_{ab}^{\left(\Gamma\right)}$
+\end_inset
+
+ defined by 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (4.28)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+.
+ This way, we obtain a symmetry adapted basis 
+\begin_inset Formula $\left\{ \vect b_{\Gamma,r,i}^{\mathrm{s.a.b.}}\right\} $
+\end_inset
+
+ as linear combinations of VSWFs 
+\begin_inset Formula $\vswfs lm{p,t}$
+\end_inset
+
+ around the constituting nanoparticles (labeled 
+\begin_inset Formula $p$
+\end_inset
+
+), 
+\begin_inset Formula 
+\[
+\vect b_{\Gamma,r,i}^{\mathrm{s.a.b.}}=\sum_{l,m,p,t}U_{\Gamma,r,i}^{p,t,l,m}\vswfs lm{p,t},
+\]
+
+\end_inset
+
+where 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+ stands for one of the six different irreps of 
+\begin_inset Formula $D_{3h}$
+\end_inset
+
+, 
+\begin_inset Formula $r$
+\end_inset
+
+ labels the different realisations of the same irrep, and the last index
+ 
+\begin_inset Formula $i$
+\end_inset
+
+ going from 1 to 
+\begin_inset Formula $d_{\Gamma}$
+\end_inset
+
+ (the dimensionality of 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+) labels the different partners of the same given irrep.
+ The number of how many times is each irrep contained in 
+\begin_inset Formula $\Gamma_{\mathrm{lat.mod.}}$
+\end_inset
+
+ (i.e.
+ the range of index 
+\begin_inset Formula $r$
+\end_inset
+
+ for given 
+\begin_inset Formula $\Gamma$
+\end_inset
+
+) depends on the multipole degree cutoff 
+\begin_inset Formula $l_{\mathrm{max}}$
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Plain Layout
+
+\lang english
+Each mode at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point shall lie in the irreducible spaces of only one of the six possible
+ irreps and it can be shown via 
+\begin_inset CommandInset citation
+LatexCommand cite
+after "eq. (2.51)"
+key "dresselhaus_group_2008"
+literal "true"
+
+\end_inset
+
+ that, at the 
+\begin_inset Formula $\Kp$
+\end_inset
+
+-point, the matrix 
+\begin_inset Formula $M\left(\omega,\vect k\right)$
+\end_inset
+
+ defined above takes a block-diagonal form in the symmetry-adapted basis,
+ 
+\begin_inset Formula 
+\[
+M\left(\omega,\vect K\right)_{\Gamma,r,i;\Gamma',r',j}^{\mathrm{s.a.b.}}=\frac{\delta_{\Gamma\Gamma'}\delta_{ij}}{d_{\Gamma}}\sum_{q}M\left(\omega,\vect K\right)_{\Gamma,r,q;\Gamma',r',q}^{\mathrm{s.a.b.}}.
+\]
+
+\end_inset
+
+This enables us to decompose the matrix according to the irreps and to solve
+ the singular value problem in each irrep separately, as done in Fig.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "smfig:dispersions"
+
+\end_inset
+
+(a).
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_body
+\end_document
diff --git a/lepaper/tmpaper.bib b/lepaper/tmpaper.bib
new file mode 100644
index 0000000..a5c7b6f
--- /dev/null
+++ b/lepaper/tmpaper.bib
@@ -0,0 +1,426 @@
+
+@book{bohren_absorption_1983,
+  title = {Absorption and Scattering of Light by Small Particles},
+  abstract = {Not Available},
+  urldate = {2014-05-09},
+  url = {http://adsabs.harvard.edu/abs/1983asls.book.....B},
+  author = {Bohren, Craig F. and Huffman, Donald R.},
+  year = {1983},
+  keywords = {Particles,LIGHT SCATTERING,ABSORPTION,DUST,THEORY},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/HES6WJTP/(Wiley science paperback series) Craig F. Bohren, Donald R. Huffman-Absorption and scattering of light by small particles-Wiley-VCH (1998).djvu}
+}
+
+@misc{SCUFF2,
+  title = {{{SCUFF}}-{{EM}}},
+  url = {http://homerreid.dyndns.org/scuff-EM/},
+  author = {Reid, Homer},
+  year = {2018},
+  note = {http://github.com/homerreid/scuff-EM}
+}
+
+@article{xu_calculation_1996,
+  title = {Calculation of the {{Addition Coefficients}} in {{Electromagnetic Multisphere}}-{{Scattering Theory}}},
+  volume = {127},
+  issn = {0021-9991},
+  abstract = {One of the most intractable problems in electromagnetic multisphere-scattering theory is the formulation and evaluation of vector addition coefficients introduced by the addition theorems for vector spherical harmonics. This paper presents an efficient approach for the calculation of both scalar and vector translational addition coefficients, which is based on fast evaluation of the Gaunt coefficients. The paper also rederives the analytical expressions for the vector translational addition coefficients and discusses the strengths and limitations of other formulations and numerical techniques found in the literature. Numerical results from the formulation derived in this paper agree with those of a previously published recursion scheme that completely avoids the use of the Gaunt coefficients, but the method of direct calculation proposed here reduces the computing time by a factor of 4\textendash{}6.},
+  number = {2},
+  urldate = {2015-11-22},
+  journal = {Journal of Computational Physics},
+  doi = {10.1006/jcph.1996.0175},
+  url = {http://www.sciencedirect.com/science/article/pii/S0021999196901758},
+  author = {Xu, Yu-lin},
+  month = sep,
+  year = {1996},
+  pages = {285-298},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/8B2TWTJ2/1-s2.0-S0021999197956874-main (2).pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/NCD6BBNZ/Xu - 1996 - Calculation of the Addition Coefficients in Electr.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/NDSF7KI2/S0021999196901758.html}
+}
+
+@article{xu_efficient_1998,
+  title = {Efficient {{Evaluation}} of {{Vector Translation Coefficients}} in {{Multiparticle Light}}-{{Scattering Theories}}},
+  volume = {139},
+  issn = {0021-9991},
+  abstract = {Vector addition theorems are a necessary ingredient in the analytical solution of electromagnetic multiparticle-scattering problems. These theorems include a large number of vector addition coefficients. There exist three basic types of analytical expressions for vector translation coefficients: Stein's (Quart. Appl. Math.19, 15 (1961)), Cruzan's (Quart. Appl. Math.20, 33 (1962)), and Xu's (J. Comput. Phys.127, 285 (1996)). Stein's formulation relates vector translation coefficients with scalar translation coefficients. Cruzan's formulas use the Wigner 3jm symbol. Xu's expressions are based on the Gaunt coefficient. Since the scalar translation coefficient can also be expressed in terms of the Gaunt coefficient, the key to the expeditious and reliable calculation of vector translation coefficients is the fast and accurate evaluation of the Wigner 3jm symbol or the Gaunt coefficient. We present highly efficient recursive approaches to accurately evaluating Wigner 3jm symbols and Gaunt coefficients. Armed with these recursive approaches, we discuss several schemes for the calculation of the vector translation coefficients, using the three general types of formulation, respectively. Our systematic test calculations show that the three types of formulas produce generally the same numerical results except that the algorithm of Stein's type is less accurate in some particular cases. These extensive test calculations also show that the scheme using the formulation based on the Gaunt coefficient is the most efficient in practical computations.},
+  number = {1},
+  urldate = {2015-11-18},
+  journal = {Journal of Computational Physics},
+  doi = {10.1006/jcph.1997.5867},
+  url = {http://www.sciencedirect.com/science/article/pii/S0021999197958678},
+  author = {Xu, Yu-lin},
+  month = jan,
+  year = {1998},
+  pages = {137-165},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/STV5263F/Xu - 1998 - Efficient Evaluation of Vector Translation Coeffic.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/VMGZRSAA/S0021999197958678.html}
+}
+
+@book{mishchenko_light_1999-1,
+  title = {Light {{Scattering}} by {{Nonspherical Particles}}: {{Theory}}, {{Measurements}}, and {{Applications}}},
+  isbn = {978-0-08-051020-0},
+  shorttitle = {Light {{Scattering}} by {{Nonspherical Particles}}},
+  abstract = {There is hardly a field of science or engineering that does not have some interest in light scattering by small particles. For example, this subject is important to climatology because the energy budget for the Earth's atmosphere is strongly affected by scattering of solar radiation by cloud and aerosol particles, and the whole discipline of remote sensing relies largely on analyzing the parameters of radiation scattered by aerosols, clouds, and precipitation. The scattering of light by spherical particles can be easily computed using the conventional Mie theory. However, most small solid particles encountered in natural and laboratory conditions have nonspherical shapes. Examples are soot and mineral aerosols, cirrus cloud particles, snow and frost crystals, ocean hydrosols, interplanetary and cometary dust grains, and microorganisms. It is now well known that scattering properties of nonspherical particles can differ dramatically from those of "equivalent" (e.g., equal-volume or equal-surface-area) spheres. Therefore, the ability to accurately compute or measure light scattering by nonspherical particles in order to clearly understand the effects of particle nonsphericity on light scattering is very important.The rapid improvement of computers and experimental techniques over the past 20 years and the development of efficient numerical approaches have resulted in major advances in this field which have not been systematically summarized. Because of the universal importance of electromagnetic scattering by nonspherical particles, papers on different aspects of this subject are scattered over dozens of diverse research and engineering journals. Often experts in one discipline (e.g., biology) are unaware of potentially useful results obtained in another discipline (e.g., antennas and propagation). This leads to an inefficient use of the accumulated knowledge and unnecessary redundancy in research activities.This book offers the first systematic and unified discussion of light scattering by nonspherical particles and its practical applications and represents the state-of-the-art of this importantresearch field. Individual chapters are written by leading experts in respective areas and cover three major disciplines: theoretical and numerical techniques, laboratory measurements, and practical applications. An overview chapter provides a concise general introduction to the subject of nonspherical scattering and should be especially useful to beginners and those interested in fast practical applications. The audience for this book will include graduate students, scientists, and engineers working on specific aspects of electromagnetic scattering by small particles and its applications in remote sensing, geophysics, astrophysics, biomedical optics, and optical engineering.* The first systematic and comprehensive treatment of electromagnetic scattering by nonspherical particles and its applications* Individual chapters are written by leading experts in respective areas* Includes a survey of all the relevant literature scattered over dozens of basic and applied research journals* Consistent use of unified definitions and notation makes the book a coherent volume* An overview chapter provides a concise general introduction to the subject of light scattering by nonspherical particles* Theoretical chapters describe specific easy-to-use computer codes publicly available on the World Wide Web* Extensively illustrated with over 200 figures, 4 in color},
+  language = {en},
+  publisher = {{Academic Press}},
+  author = {Mishchenko, Michael I. and Hovenier, Joachim W. and Travis, Larry D.},
+  month = sep,
+  year = {1999},
+  keywords = {Science / Physics / General,Science / Applied Sciences,Science / Earth Sciences / Meteorology \& Climatology,Science / Physics / Geophysics,Science / Earth Sciences / Oceanography,Science / Earth Sciences / General}
+}
+
+@article{mackowski_calculation_1996,
+  title = {Calculation of the {{T}} Matrix and the Scattering Matrix for Ensembles of Spheres},
+  volume = {13},
+  issn = {1084-7529, 1520-8532},
+  language = {en},
+  number = {11},
+  urldate = {2015-11-09},
+  journal = {Journal of the Optical Society of America A},
+  doi = {10.1364/JOSAA.13.002266},
+  url = {https://www.osapublishing.org/josaa/abstract.cfm?uri=josaa-13-11-2266},
+  author = {Mackowski, Daniel W. and Mishchenko, Michael I.},
+  month = nov,
+  year = {1996},
+  pages = {2266},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/V59UV9H9/josaa-13-11-2266.pdf}
+}
+
+@misc{mackowski_mstm_2013,
+  title = {{{MSTM}} 3.0: {{A}} Multiple Sphere {{T}} -Matrix {{FORTRAN}} Code for Use on Parallel Computer Clusters},
+  url = {http://www.eng.auburn.edu/~dmckwski/scatcodes/},
+  author = {Mackowski, Daniel W.},
+  year = {2013},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/RQMQMC7H/mstm-manual-2013-v3.0.pdf}
+}
+
+@book{jackson_classical_1998,
+  address = {{New York}},
+  edition = {3 edition},
+  title = {Classical {{Electrodynamics Third Edition}}},
+  isbn = {978-0-471-30932-1},
+  abstract = {A revision of the defining book covering the physics and classical mathematics necessary to understand electromagnetic fields in materials and at surfaces and interfaces. The third edition has been revised to address the changes in emphasis and applications that have occurred in the past twenty years.},
+  language = {English},
+  publisher = {{Wiley}},
+  author = {Jackson, John David},
+  month = aug,
+  year = {1998},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/3BWPD4BK/John David Jackson-Classical Electrodynamics-Wiley (1999).djvu}
+}
+
+@misc{kristensson_spherical_2014,
+  title = {Spherical {{Vector Waves}}},
+  urldate = {2014-05-20},
+  url = {http://www.eit.lth.se/fileadmin/eit/courses/eit080f/Literature/book.pdf},
+  author = {Kristensson, Gerhard},
+  month = jan,
+  year = {2014},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/7MVDRPF2/Kristensson - 2014 - Spherical Vector Waves.pdf}
+}
+
+@incollection{mishchenko_t-matrix_1999,
+  title = {T-Matrix {{Method}} and {{Its Applications}}},
+  isbn = {978-0-08-051020-0},
+  abstract = {There is hardly a field of science or engineering that does not have some interest in light scattering by small particles. For example, this subject is important to climatology because the energy budget for the Earth's atmosphere is strongly affected by scattering of solar radiation by cloud and aerosol particles, and the whole discipline of remote sensing relies largely on analyzing the parameters of radiation scattered by aerosols, clouds, and precipitation. The scattering of light by spherical particles can be easily computed using the conventional Mie theory. However, most small solid particles encountered in natural and laboratory conditions have nonspherical shapes. Examples are soot and mineral aerosols, cirrus cloud particles, snow and frost crystals, ocean hydrosols, interplanetary and cometary dust grains, and microorganisms. It is now well known that scattering properties of nonspherical particles can differ dramatically from those of "equivalent" (e.g., equal-volume or equal-surface-area) spheres. Therefore, the ability to accurately compute or measure light scattering by nonspherical particles in order to clearly understand the effects of particle nonsphericity on light scattering is very important.The rapid improvement of computers and experimental techniques over the past 20 years and the development of efficient numerical approaches have resulted in major advances in this field which have not been systematically summarized. Because of the universal importance of electromagnetic scattering by nonspherical particles, papers on different aspects of this subject are scattered over dozens of diverse research and engineering journals. Often experts in one discipline (e.g., biology) are unaware of potentially useful results obtained in another discipline (e.g., antennas and propagation). This leads to an inefficient use of the accumulated knowledge and unnecessary redundancy in research activities.This book offers the first systematic and unified discussion of light scattering by nonspherical particles and its practical applications and represents the state-of-the-art of this importantresearch field. Individual chapters are written by leading experts in respective areas and cover three major disciplines: theoretical and numerical techniques, laboratory measurements, and practical applications. An overview chapter provides a concise general introduction to the subject of nonspherical scattering and should be especially useful to beginners and those interested in fast practical applications. The audience for this book will include graduate students, scientists, and engineers working on specific aspects of electromagnetic scattering by small particles and its applications in remote sensing, geophysics, astrophysics, biomedical optics, and optical engineering.* The first systematic and comprehensive treatment of electromagnetic scattering by nonspherical particles and its applications* Individual chapters are written by leading experts in respective areas* Includes a survey of all the relevant literature scattered over dozens of basic and applied research journals* Consistent use of unified definitions and notation makes the book a coherent volume* An overview chapter provides a concise general introduction to the subject of light scattering by nonspherical particles* Theoretical chapters describe specific easy-to-use computer codes publicly available on the World Wide Web* Extensively illustrated with over 200 figures, 4 in color},
+  language = {en},
+  booktitle = {Light {{Scattering}} by {{Nonspherical Particles}}: {{Theory}}, {{Measurements}}, and {{Applications}}},
+  publisher = {{Academic Press}},
+  author = {Mishchenko, Michael I. and {Travis, Larry D.} and Macke, Andreas},
+  editor = {Mishchenko, Michael I. and Hovenier, Joachim W. and Travis, Larry D.},
+  month = sep,
+  year = {1999},
+  keywords = {Science / Physics / General,Science / Applied Sciences,Science / Earth Sciences / Meteorology \& Climatology,Science / Physics / Geophysics,Science / Earth Sciences / Oceanography,Science / Earth Sciences / General},
+  pages = {147-172}
+}
+
+@article{schulz_point-group_1999,
+  title = {Point-Group Symmetries in Electromagnetic Scattering},
+  volume = {16},
+  issn = {1084-7529, 1520-8532},
+  language = {en},
+  number = {4},
+  urldate = {2016-08-04},
+  journal = {Journal of the Optical Society of America A},
+  doi = {10.1364/JOSAA.16.000853},
+  url = {https://www.osapublishing.org/abstract.cfm?URI=josaa-16-4-853},
+  author = {Schulz, F. Michael and Stamnes, Knut and Stamnes, J. J.},
+  month = apr,
+  year = {1999},
+  pages = {853},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/X9X48A6G/josaa-16-4-853.pdf}
+}
+
+@article{mishchenko_t-matrix_1994,
+  title = {T-Matrix Computations of Light Scattering by Large Spheroidal Particles},
+  volume = {109},
+  issn = {0030-4018},
+  abstract = {It is well known that T-matrix computations of light scattering by nonspherical particles may suffer from the ill-conditionality of the process of matrix inversion, which has precluded calculations for particle size parameters larger than about 25. It is demonstrated that calculating the T-matrix using extended-precision instead of double-precision floating-point variables is an effective approach for suppressing the numerical instability in computations for spheroids and allows one to increase the maximum particle size parameter for which T-matrix computations converge by as significant a factor as 2\textendash{}2.7. Yet this approach requires only a negligibly small extra memory, an affordable increase in CPU time consumption, and practically no additional programming effort. As a result, the range of particle size parameters, for which rigorous T-matrix computations of spheroidal scattering can be performed, now covers a substantial fraction of the gap between the domains of applicability of the Rayleigh and geometrical optics approximations.},
+  number = {1\textendash{}2},
+  urldate = {2017-01-18},
+  journal = {Optics Communications},
+  doi = {10.1016/0030-4018(94)90731-5},
+  url = {http://www.sciencedirect.com/science/article/pii/0030401894907315},
+  author = {Mishchenko, Michael I. and Travis, Larry D.},
+  month = jun,
+  year = {1994},
+  pages = {16-21},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/FT8KN354/mishchenko1994.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/TB425HGN/0030401894907315.html}
+}
+
+@article{mishchenko_t-matrix_1996,
+  series = {Light {{Scattering}} by {{Non}}-{{Spherical Particles}}},
+  title = {T-Matrix Computations of Light Scattering by Nonspherical Particles: {{A}} Review},
+  volume = {55},
+  issn = {0022-4073},
+  shorttitle = {T-Matrix Computations of Light Scattering by Nonspherical Particles},
+  abstract = {We review the current status of Waterman's T-matrix approach which is one of the most powerful and widely used tools for accurately computing light scattering by nonspherical particles, both single and composite, based on directly solving Maxwell's equations. Specifically, we discuss the analytical method for computing orientationally-averaged light-scattering characteristics for ensembles of nonspherical particles, the methods for overcoming the numerical instability in calculating the T matrix for single nonspherical particles with large size parameters and/or extreme geometries, and the superposition approach for computing light scattering by composite/aggregated particles. Our discussion is accompanied by multiple numerical examples demonstrating the capabilities of the T-matrix approach and showing effects of nonsphericity of simple convex particles (spheroids) on light scattering.},
+  number = {5},
+  urldate = {2017-01-18},
+  journal = {Journal of Quantitative Spectroscopy and Radiative Transfer},
+  doi = {10.1016/0022-4073(96)00002-7},
+  url = {http://www.sciencedirect.com/science/article/pii/0022407396000027},
+  author = {Mishchenko, Michael I. and Travis, Larry D. and Mackowski, Daniel W.},
+  month = may,
+  year = {1996},
+  pages = {535-575},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/8EA7QMDG/Mishchenko et al. - 1996 - T-matrix computations of light scattering by nonsp.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/HNWF8F6R/0022407396000027.html}
+}
+
+@article{hakala_lasing_2017,
+  title = {Lasing in Dark and Bright Modes of a Finite-Sized Plasmonic Lattice},
+  volume = {8},
+  copyright = {\textcopyright{} 2017 Macmillan Publishers Limited, part of Springer Nature. All rights reserved.},
+  issn = {2041-1723},
+  abstract = {Plasmonic dark modes are promising candidates for lasing applications. Here, Hakalaet al. show lasing at visible wavelengths in dark and bright modes of an array of silver nanoparticles combined with optically pumped dye molecules, opening up a route to utilization of all modes of plasmonic lattices.},
+  language = {en},
+  urldate = {2017-03-28},
+  journal = {Nature Communications},
+  doi = {10.1038/ncomms13687},
+  url = {http://www.nature.com/ncomms/2017/170103/ncomms13687/full/ncomms13687.html},
+  author = {Hakala, T. K. and Rekola, H. T. and V{\"a}kev{\"a}inen, A. I. and Martikainen, J.-P. and Ne{\v c}ada, M. and Moilanen, A. J. and T{\"o}rm{\"a}, P.},
+  month = jan,
+  year = {2017},
+  pages = {13687},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/73KCXGAP/ncomms13687.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/J6R8MHBH/ncomms13687-s1.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/X4PNEUQN/ncomms13687.html}
+}
+
+@article{mackowski_effective_2001,
+  series = {Light {{Scattering}} by {{Non}}-{{Spherical Particles}}},
+  title = {An Effective Medium Method for Calculation of the {{T}} Matrix of Aggregated Spheres},
+  volume = {70},
+  issn = {0022-4073},
+  abstract = {An effective medium approach is developed for describing the radiative scattering characteristics of large-scale clusters of spheres. The formulation assumes that the waves exciting each sphere in the cluster can be described by a regular vector harmonic expansion, centered about a common origin of the cluster, and characterized by an effective propagation constant mek. By combining this description with the multiple sphere interaction equations a `homogeneous' T matrix of the cluster is derived, which is analogous to using the effective propagation constant models of the Varadans in conjunction with Waterman's EBCM. However, it is shown that the homogeneous T matrix will not automatically satisfy energy conservation because it cannot account for dependent scattering effects among the spheres. A `discrete' formulation of the T matrix is then developed which retains the effective medium description of the exciting field yet provides for energy conservation. Illustrative calculations show that the effective medium T matrix can provide accurate predictions of the cross sections and scattering matrices of clusters containing a large number of uniformly packed spheres, yet this approximation uses a fraction of the computational time required for an exact solution.},
+  number = {4\textendash{}6},
+  urldate = {2017-06-05},
+  journal = {Journal of Quantitative Spectroscopy and Radiative Transfer},
+  doi = {10.1016/S0022-4073(01)00022-X},
+  url = {http://www.sciencedirect.com/science/article/pii/S002240730100022X},
+  author = {Mackowski, Daniel W.},
+  month = aug,
+  year = {2001},
+  pages = {441-464},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/9E7R7IRX/Mackowski - 2001 - An effective medium method for calculation of the .pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/D75CJ78C/S002240730100022X.html}
+}
+
+@book{dresselhaus_group_2008,
+  title = {Group {{Theory}}: {{Application}} to the {{Physics}} of {{Condensed Matter}}},
+  isbn = {978-3-540-32899-5},
+  abstract = {Every process in physics is governed by selection rules that are the consequence of symmetry requirements. The beauty and strength of group theory resides...},
+  urldate = {2017-10-31},
+  publisher = {{Springer, Berlin, Heidelberg}},
+  url = {//www.springer.com/us/book/9783540328971},
+  author = {Dresselhaus, Mildred S. and Dresselhaus, Gene and Jorio, Ado},
+  year = {2008},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/GFGPVB4A/Mildred_S._Dresselhaus,_Gene_Dresselhaus,_Ado_Jorio_Group_theory_application_to_the_physics_of_condensed_matter.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/E78682CJ/9783540328971.html}
+}
+
+@article{linton_lattice_2010,
+  title = {Lattice {{Sums}} for the {{Helmholtz Equation}}},
+  volume = {52},
+  issn = {0036-1445},
+  abstract = {A survey of different representations for lattice sums for the Helmholtz equation is made. These sums arise naturally when dealing with wave scattering by periodic structures. One of the main objectives is to show how the various forms depend on the dimension d of the underlying space and the lattice dimension \$d\_\textbackslash{}Lambda\$. Lattice sums are related to, and can be calculated from, the quasi-periodic Green's function and this object serves as the starting point of the analysis.},
+  number = {4},
+  journal = {SIAM Rev.},
+  doi = {10.1137/09075130X},
+  url = {http://epubs.siam.org/doi/10.1137/09075130X},
+  author = {Linton, C.},
+  month = jan,
+  year = {2010},
+  pages = {630-674},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/T86ATKYB/09075130x.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/ETB8X4S9/09075130X.html}
+}
+
+@book{bradley_mathematical_1972,
+  title = {The Mathematical Theory of Symmetry in Solids; Representation Theory for Point Groups and Space Groups},
+  isbn = {978-0-19-851920-1},
+  urldate = {2018-07-25},
+  publisher = {{Clarendon Press, Oxford}},
+  url = {http://gen.lib.rus.ec/book/index.php?md5=8539E3400CF65B6CC4FAC71B9DF286C5},
+  author = {Bradley, C. J. and Cracknell, A. P.},
+  year = {1972},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/SB5ZN5WH/C.J. Bradley, A.P. Cracknell - The mathematical theory of symmetry in solids_ representation theory for point groups and space groups (1972, Clarendon Press).djvu}
+}
+
+@article{moroz_quasi-periodic_2006,
+  title = {Quasi-Periodic {{Green}}'s Functions of the {{Helmholtz}} and {{Laplace}} Equations},
+  volume = {39},
+  issn = {0305-4470},
+  abstract = {A classical problem of free-space Green's function G 0{$\Lambda$} representations of the Helmholtz equation is studied in various quasi-periodic cases, i.e., when an underlying periodicity is imposed in less dimensions than is the dimension of an embedding space. Exponentially convergent series for the free-space quasi-periodic G 0{$\Lambda$} and for the expansion coefficients D L of G 0{$\Lambda$} in the basis of regular (cylindrical in two dimensions and spherical in three dimension (3D)) waves, or lattice sums, are reviewed and new results for the case of a one-dimensional (1D) periodicity in 3D are derived. From a mathematical point of view, a derivation of exponentially convergent representations for Schl{\"o}milch series of cylindrical and spherical Hankel functions of any integer order is accomplished. Exponentially convergent series for G 0{$\Lambda$} and lattice sums D L hold for any value of the Bloch momentum and allow G 0{$\Lambda$} to be efficiently evaluated also in the periodicity plane. The quasi-periodic Green's functions of the Laplace equation are obtained from the corresponding representations of G 0{$\Lambda$} of the Helmholtz equation by taking the limit of the wave vector magnitude going to zero. The derivation of relevant results in the case of a 1D periodicity in 3D highlights the common part which is universally applicable to any of remaining quasi-periodic cases. The results obtained can be useful for the numerical solution of boundary integral equations for potential flows in fluid mechanics, remote sensing of periodic surfaces, periodic gratings, and infinite arrays of resonators coupled to a waveguide, in many contexts of simulating systems of charged particles, in molecular dynamics, for the description of quasi-periodic arrays of point interactions in quantum mechanics, and in various ab initio first-principle multiple-scattering theories for the analysis of diffraction of classical and quantum waves.},
+  language = {en},
+  number = {36},
+  urldate = {2018-08-14},
+  journal = {J. Phys. A: Math. Gen.},
+  doi = {10.1088/0305-4470/39/36/009},
+  url = {http://stacks.iop.org/0305-4470/39/i=36/a=009},
+  author = {Moroz, Alexander},
+  year = {2006},
+  pages = {11247},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/268RXLJ4/Moroz - 2006 - Quasi-periodic Green's functions of the Helmholtz .pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/MGA5XR44/dlserr.pdf}
+}
+
+@article{linton_one-_2009,
+  title = {One- and Two-Dimensional Lattice Sums for the Three-Dimensional {{Helmholtz}} Equation},
+  volume = {228},
+  issn = {0021-9991},
+  abstract = {The accurate and efficient computation of lattice sums for the three-dimensional Helmholtz equation is considered for the cases where the underlying lattice is one- or two-dimensional. We demonstrate, using careful numerical computations, that the reduction method, in which the sums for a two-dimensional lattice are expressed as a sum of one-dimensional lattice sums leads to an order-of-magnitude improvement in performance over the well-known Ewald method. In the process we clarify and improve on a number of results originally formulated by Twersky in the 1970s.},
+  number = {6},
+  urldate = {2018-08-14},
+  journal = {Journal of Computational Physics},
+  doi = {10.1016/j.jcp.2008.11.013},
+  url = {http://www.sciencedirect.com/science/article/pii/S0021999108005962},
+  author = {Linton, C. M. and Thompson, I.},
+  month = apr,
+  year = {2009},
+  keywords = {Helmholtz equation,Ewald summation,Clausen function,Lattice reduction,Lattice sum,Schlömilch series},
+  pages = {1815-1829},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/YMRZHBY4/Linton ja Thompson - 2009 - One- and two-dimensional lattice sums for the thre.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/Z8CFQ6S9/S0021999108005962.html}
+}
+
+@book{olver_nist_2010,
+  edition = {1 Pap/Cdr},
+  title = {{{NIST}} Handbook of Mathematical Functions},
+  isbn = {978-0-521-14063-8},
+  urldate = {2018-08-20},
+  publisher = {{Cambridge University Press}},
+  url = {http://gen.lib.rus.ec/book/index.php?md5=7750A842DAAE07EBE30D597EB1352408},
+  author = {Olver, Frank W. J. and Lozier, Daniel W. and Boisvert, Ronald F. and Clark, Charles W.},
+  year = {2010},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/ZJ5LBQ8W/Olver ym. - 2010 - NIST handbook of mathematical functions.pdf}
+}
+
+@article{NIST:DLMF,
+  title = {{{NIST Digital Library}} of {{Mathematical Functions}}},
+  url = {http://dlmf.nist.gov/},
+  key = {DLMF},
+  note = {F.~W.~J. Olver, A.~B. Olde Daalhuis, D.~W. Lozier, B.~I. Schneider, R.~F. Boisvert, C.~W. Clark, B.~R. Miller and B.~V. Saunders, eds.}
+}
+
+@article{reid_efficient_2015,
+  archivePrefix = {arXiv},
+  eprinttype = {arxiv},
+  eprint = {1307.2966},
+  title = {Efficient {{Computation}} of {{Power}}, {{Force}}, and {{Torque}} in {{BEM Scattering Calculations}}},
+  volume = {63},
+  issn = {0018-926X, 1558-2221},
+  abstract = {We present concise, computationally efficient formulas for several quantities of interest -- including absorbed and scattered power, optical force (radiation pressure), and torque -- in scattering calculations performed using the boundary-element method (BEM) [also known as the method of moments (MOM)]. Our formulas compute the quantities of interest \textbackslash{}textit\{directly\} from the BEM surface currents with no need ever to compute the scattered electromagnetic fields. We derive our new formulas and demonstrate their effectiveness by computing power, force, and torque in a number of example geometries. Free, open-source software implementations of our formulas are available for download online.},
+  number = {8},
+  urldate = {2018-09-23},
+  journal = {IEEE Transactions on Antennas and Propagation},
+  doi = {10.1109/TAP.2015.2438393},
+  url = {http://arxiv.org/abs/1307.2966},
+  author = {Reid, M. T. Homer and Johnson, Steven G.},
+  month = aug,
+  year = {2015},
+  keywords = {Physics - Classical Physics,Physics - Computational Physics},
+  pages = {3588-3598},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/I2DXTKUF/Reid ja Johnson - 2015 - Efficient Computation of Power, Force, and Torque .pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/LG7AVZDH/1307.html}
+}
+
+@article{guo_lasing_2019,
+  title = {Lasing at \${{K}}\$ {{Points}} of a {{Honeycomb Plasmonic Lattice}}},
+  volume = {122},
+  abstract = {We study lasing at the high-symmetry points of the Brillouin zone in a honeycomb plasmonic lattice. We use symmetry arguments to define singlet and doublet modes at the K points of the reciprocal space. We experimentally demonstrate lasing at the K points that is based on plasmonic lattice modes and two-dimensional feedback. By comparing polarization properties to T-matrix simulations, we identify the lasing mode as one of the singlets with an energy minimum at the K point enabling feedback. Our results offer prospects for studies of topological lasing in radiatively coupled systems.},
+  number = {1},
+  urldate = {2019-01-10},
+  journal = {Phys. Rev. Lett.},
+  doi = {10.1103/PhysRevLett.122.013901},
+  url = {https://link.aps.org/doi/10.1103/PhysRevLett.122.013901},
+  author = {Guo, R. and Ne{\v c}ada, M. and Hakala, T. K. and V{\"a}kev{\"a}inen, A. I. and T{\"o}rm{\"a}, P.},
+  month = jan,
+  year = {2019},
+  pages = {013901},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/TDGW4CZ5/Guo ym. - 2019 - Lasing at $K$ Points of a Honeycomb Plasmonic Latt.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/8BW4R9F6/PhysRevLett.122.html}
+}
+
+@article{mie_beitrage_1908,
+  title = {Beitr{\"a}ge Zur {{Optik}} Tr{\"u}ber {{Medien}}, Speziell Kolloidaler {{Metall{\"o}sungen}}},
+  volume = {330},
+  copyright = {Copyright \textcopyright{} 1908 WILEY-VCH Verlag GmbH \& Co. KGaA, Weinheim},
+  issn = {1521-3889},
+  language = {en},
+  number = {3},
+  urldate = {2014-11-30},
+  journal = {Ann. Phys.},
+  doi = {10.1002/andp.19083300302},
+  url = {http://onlinelibrary.wiley.com/doi/10.1002/andp.19083300302/abstract},
+  author = {Mie, Gustav},
+  month = jan,
+  year = {1908},
+  pages = {377-445},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/RM9J9RYH/Mie - 1908 - Beiträge zur Optik trüber Medien, speziell kolloid.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/F5A7EX6R/abstract.html}
+}
+
+@book{kristensson_scattering_2016,
+  address = {{Edison, NJ}},
+  title = {Scattering of {{Electromagnetic Waves}} by {{Obstacles}}},
+  isbn = {978-1-61353-221-8},
+  abstract = {This book is an introduction to some of the most important properties of electromagnetic waves and their interaction with passive materials and scatterers. The main purpose of the book is to give a theoretical treatment of these scattering phenomena, and to illustrate numerical computations of some canonical scattering problems for different geometries and materials. The scattering theory is also important in the theory of passive antennas, and this book gives several examples on this topic. Topics covered include an introduction to the basic equations used in scattering; the Green functions and dyadics; integral representation of fields; introductory scattering theory; scattering in the time domain; approximations and applications; spherical vector waves; scattering by spherical objects; the null-field approach; and propagation in stratified media. The book is organised along two tracks, which can be studied separately or together. Track 1 material is appropriate for a first reading of the textbook, while Track 2 contains more advanced material suited for the second reading and for reference. Exercises are included for each chapter.},
+  language = {English},
+  publisher = {{Scitech Publishing}},
+  url = {http://gen.lib.rus.ec/book/index.php?md5=00CCB3E221E741ADDB2E236FD4A9F002},
+  author = {Kristensson, Gerhard},
+  month = jul,
+  year = {2016},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/ZRYZ4KLK/Kristensson - 2016 - Scattering of Electromagnetic Waves by Obstacles.pdf}
+}
+
+@article{ganesh_convergence_2012,
+  title = {Convergence Analysis with Parameter Estimates for a Reduced Basis Acoustic Scattering {{T}}-Matrix Method},
+  volume = {32},
+  issn = {0272-4979},
+  abstract = {Abstract.  The celebrated truncated T-matrix method for wave propagation models belongs to a class of the reduced basis methods (RBMs), with the parameters bein},
+  language = {en},
+  number = {4},
+  urldate = {2019-07-03},
+  journal = {IMA J Numer Anal},
+  doi = {10.1093/imanum/drr041},
+  url = {https://academic.oup.com/imajna/article/32/4/1348/654510},
+  author = {Ganesh, M. and Hawkins, S. C. and Hiptmair, R.},
+  month = oct,
+  year = {2012},
+  pages = {1348-1374},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/2CRM9IEU/ganesh2012.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/KLKJBTZU/Ganesh ym. - 2012 - Convergence analysis with parameter estimates for .pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/N5H8B7SF/654510.html}
+}
+
+@book{chew_fast_2000,
+  series = {Artech {{House Antennas}} and {{Propagation Library}}},
+  title = {Fast and {{Efficient Algorithms}} in {{Computational Electromagnetics}}},
+  isbn = {978-1-58053-152-8},
+  urldate = {2019-07-31},
+  publisher = {{Artech House Publishers}},
+  url = {http://gen.lib.rus.ec/book/index.php?md5=2A7D2CE03DB8CFC14E7189E9A441F759},
+  author = {Chew, Weng Cho and Jin, Jian-Ming and Michielssen, Eric and Song, Jiming},
+  year = {2000}
+}
+
+@article{pourjamal_lasing_2019,
+  title = {Lasing in {{Ni Nanodisk Arrays}}},
+  abstract = {Lasing in Ni Nanodisk Arrays},
+  language = {en},
+  urldate = {2019-07-31},
+  journal = {ACS Nano},
+  doi = {10.1021/acsnano.9b01006},
+  url = {https://pubs.acs.org/doi/suppl/10.1021/acsnano.9b01006},
+  author = {Pourjamal, Sara and Hakala, Tommi K. and Ne{\v c}ada, Marek and {Freire-Fern{\'a}ndez}, Francisco and Kataja, Mikko and Rekola, Heikki and Martikainen, Jani-Petri and T{\"o}rm{\"a}, P{\"a}ivi and van Dijken, Sebastiaan},
+  month = apr,
+  year = {2019},
+  file = {/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/C4SN68I6/Pourjamal ym. - 2019 - Lasing in Ni Nanodisk Arrays.pdf;/u/46/necadam1/unix/.mozilla/firefox/6m8fw48s.default/zotero/storage/S6AU6FV9/acsnano.html}
+}
+
+