Move MIN,MAX,SQ macros to a common header
[qpms.git] / lepaper / infinite.lyx
blob9e74628092b8d0c5a4d346b4b15885fae623cdec
1 #LyX 2.4 created this file. For more info see https://www.lyx.org/
2 \lyxformat 583
3 \begin_document
4 \begin_header
5 \save_transient_properties true
6 \origin unavailable
7 \textclass article
8 \use_default_options true
9 \maintain_unincluded_children false
10 \language english
11 \language_package default
12 \inputencoding utf8
13 \fontencoding auto
14 \font_roman "default" "TeX Gyre Pagella"
15 \font_sans "default" "default"
16 \font_typewriter "default" "default"
17 \font_math "auto" "auto"
18 \font_default_family default
19 \use_non_tex_fonts false
20 \font_sc false
21 \font_roman_osf true
22 \font_sans_osf false
23 \font_typewriter_osf false
24 \font_sf_scale 100 100
25 \font_tt_scale 100 100
26 \use_microtype false
27 \use_dash_ligatures true
28 \graphics default
29 \default_output_format default
30 \output_sync 0
31 \bibtex_command default
32 \index_command default
33 \float_placement class
34 \float_alignment class
35 \paperfontsize default
36 \spacing single
37 \use_hyperref true
38 \pdf_author "Marek Nečada"
39 \pdf_bookmarks true
40 \pdf_bookmarksnumbered false
41 \pdf_bookmarksopen false
42 \pdf_bookmarksopenlevel 1
43 \pdf_breaklinks false
44 \pdf_pdfborder false
45 \pdf_colorlinks false
46 \pdf_backref false
47 \pdf_pdfusetitle true
48 \papersize default
49 \use_geometry false
50 \use_package amsmath 1
51 \use_package amssymb 1
52 \use_package cancel 1
53 \use_package esint 1
54 \use_package mathdots 1
55 \use_package mathtools 1
56 \use_package mhchem 1
57 \use_package stackrel 1
58 \use_package stmaryrd 1
59 \use_package undertilde 1
60 \cite_engine basic
61 \cite_engine_type default
62 \biblio_style plain
63 \use_bibtopic false
64 \use_indices false
65 \paperorientation portrait
66 \suppress_date false
67 \justification true
68 \use_refstyle 1
69 \use_minted 0
70 \use_lineno 0
71 \index Index
72 \shortcut idx
73 \color #008000
74 \end_index
75 \secnumdepth 3
76 \tocdepth 3
77 \paragraph_separation indent
78 \paragraph_indentation default
79 \is_math_indent 0
80 \math_numbering_side default
81 \quotes_style english
82 \dynamic_quotes 0
83 \papercolumns 1
84 \papersides 1
85 \paperpagestyle default
86 \tablestyle default
87 \tracking_changes false
88 \output_changes false
89 \html_math_output 0
90 \html_css_as_file 0
91 \html_be_strict false
92 \end_header
94 \begin_body
96 \begin_layout Section
97 Infinite periodic systems
98 \begin_inset FormulaMacro
99 \newcommand{\dlv}{\vect a}
100 \end_inset
103 \begin_inset FormulaMacro
104 \newcommand{\rlv}{\vect b}
105 \end_inset
108 \end_layout
110 \begin_layout Standard
111 Although large finite systems are where MSTMM excels the most, there are
112  several reasons that makes its extension to infinite lattices (where periodic
113  boundary conditions might be applied) desirable as well.
114  Other methods might be already fast enough, but MSTMM will be faster in
115  most cases in which there is enough spacing between the neighboring particles.
116  MSTMM works well with any space group symmetry the system might have (as
117  opposed to, for example, FDTD with cubic mesh applied to a honeycomb lattice),
118  which makes e.g.
119  application of group theory in mode analysis quite easy.
120 \begin_inset Note Note
121 status open
123 \begin_layout Plain Layout
124 Topology anoyne?
125 \end_layout
127 \end_inset
129  And finally, having a method that handles well both infinite and large
130  finite system gives a possibility to study finite-size effects in periodic
131  scatterer arrays.
132 \end_layout
134 \begin_layout Subsection
135 Notation
136 \end_layout
138 \begin_layout Standard
139 TODO Fourier transforms, Delta comb, lattice bases, reciprocal lattices
140  etc.
141 \end_layout
143 \begin_layout Subsection
144 Formulation of the problem
145 \end_layout
147 \begin_layout Standard
148 Let us have a linear system of compact EM scatterers on a homogeneous background
149  as in Section 
150 \begin_inset CommandInset ref
151 LatexCommand eqref
152 reference "subsec:Multiple-scattering"
153 plural "false"
154 caps "false"
155 noprefix "false"
157 \end_inset
159 , but this time, the system shall be periodic: let there be a 
160 \begin_inset Formula $d$
161 \end_inset
163 -dimensional (
164 \begin_inset Formula $d$
165 \end_inset
167  can be 1, 2 or 3) lattice embedded into the three-dimensional real space,
168  with lattice vectors 
169 \begin_inset Formula $\left\{ \dlv_{i}\right\} _{i=1}^{d}$
170 \end_inset
172 , and let the lattice points be labeled with an 
173 \begin_inset Formula $d$
174 \end_inset
176 -dimensional integar multiindex 
177 \begin_inset Formula $\vect n\in\ints^{d}$
178 \end_inset
180 , so the lattice points have cartesian coordinates 
181 \begin_inset Formula $\vect R_{\vect n}=\sum_{i=1}^{d}n_{i}\vect a_{i}$
182 \end_inset
185  There can be several scatterers per unit cell with indices 
186 \begin_inset Formula $\alpha$
187 \end_inset
189  from set 
190 \begin_inset Formula $\mathcal{P}_{1}$
191 \end_inset
193  and (relative) positions inside the unit cell 
194 \begin_inset Formula $\vect r_{\alpha}$
195 \end_inset
197 ; any particle of the periodic system can thus be labeled by a multiindex
198  from 
199 \begin_inset Formula $\mathcal{P}=\ints^{d}\times\mathcal{P}_{1}$
200 \end_inset
203  The scatterers are located at positions 
204 \begin_inset Formula $\vect r_{\vect n,\alpha}=\vect R_{\vect n}+\vect r_{\alpha}$
205 \end_inset
207  and their 
208 \begin_inset Formula $T$
209 \end_inset
211 -matrices are periodic, 
212 \begin_inset Formula $T_{\vect n,\alpha}=T_{\alpha}$
213 \end_inset
216  In such system, the multiple-scattering problem 
217 \begin_inset CommandInset ref
218 LatexCommand eqref
219 reference "eq:Multiple-scattering problem"
220 plural "false"
221 caps "false"
222 noprefix "false"
224 \end_inset
226  can be rewritten as
227 \end_layout
229 \begin_layout Standard
230 \begin_inset Formula 
231 \begin{equation}
232 \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}
233 \end{equation}
235 \end_inset
238 \end_layout
240 \begin_layout Standard
241 Due to periodicity, we can also write 
242 \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}$
243 \end_inset
246  Assuming quasi-periodic right-hand side with quasi-momentum 
247 \begin_inset Formula $\vect k$
248 \end_inset
251 \begin_inset Formula $\rcoeffincp{\vect n,\alpha}=\rcoeffincp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}}$
252 \end_inset
254 , the solutions of 
255 \begin_inset CommandInset ref
256 LatexCommand eqref
257 reference "eq:Multiple-scattering problem periodic"
258 plural "false"
259 caps "false"
260 noprefix "false"
262 \end_inset
264  will be also quasi-periodic according to Bloch theorem, 
265 \begin_inset Formula $\outcoeffp{\vect n,\alpha}=\outcoeffp{\vect 0,\alpha}\left(\vect k\right)e^{i\vect k\cdot\vect R_{\vect n}}$
266 \end_inset
268 , and eq.
270 \begin_inset CommandInset ref
271 LatexCommand eqref
272 reference "eq:Multiple-scattering problem periodic"
273 plural "false"
274 caps "false"
275 noprefix "false"
277 \end_inset
279  can be rewritten as follows
280 \begin_inset Formula 
281 \begin{align}
282 \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 \\
283 \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 \\
284 \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 \\
285 \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}
286 \end{align}
288 \end_inset
290 so we reduced the initial scattering problem to one involving only the field
291  expansion coefficients from a single unit cell, but we need to compute
292  the 
293 \begin_inset Quotes eld
294 \end_inset
296 lattice Fourier transform
297 \begin_inset Quotes erd
298 \end_inset
300  of the translation operator,
301 \begin_inset Formula 
302 \begin{equation}
303 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}
304 \end{equation}
306 \end_inset
308 evaluation of which is possible but quite non-trivial due to the infinite
309  lattice sum, so we explain it separately in Sect.
311 \begin_inset CommandInset ref
312 LatexCommand eqref
313 reference "subsec:W operator evaluation"
314 plural "false"
315 caps "false"
316 noprefix "false"
318 \end_inset
321 \end_layout
323 \begin_layout Standard
324 As in the case of a finite system, eq.
325  can be written in a shorter block-matrix form,
326 \begin_inset Formula 
327 \begin{equation}
328 \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}
329 \end{equation}
331 \end_inset
333  Eq.
335 \begin_inset CommandInset ref
336 LatexCommand eqref
337 reference "eq:Multiple-scattering problem unit cell"
338 plural "false"
339 caps "false"
340 noprefix "false"
342 \end_inset
344  can be used to calculate electromagnetic response of the structure to external
345  quasiperiodic driving field – most notably a plane wave.
346  However, the non-trivial solutions of the equation with right hand side
347  (i.e.
348  the external driving) set to zero, 
349 \begin_inset Formula 
350 \begin{equation}
351 \left(I-WT\right)\outcoeffp{\vect 0}\left(\vect k\right)=0,\label{eq:lattice mode equation}
352 \end{equation}
354 \end_inset
356 describes the 
357 \emph on
358 lattice modes.
360 \emph default
361  Non-trivial solutions to 
362 \begin_inset CommandInset ref
363 LatexCommand eqref
364 reference "eq:lattice mode equation"
365 plural "false"
366 caps "false"
367 noprefix "false"
369 \end_inset
371  exist if the matrix on the left-hand side 
372 \begin_inset Formula $M\left(\omega,\vect k\right)=\left(I-W\left(\omega,\vect k\right)T\left(\omega\right)\right)$
373 \end_inset
375  is singular – this condition gives the 
376 \emph on
377 dispersion relation
378 \emph default
379  for the periodic structure.
380  Note that in realistic (lossy) systems, at least one of the pair 
381 \begin_inset Formula $\omega,\vect k$
382 \end_inset
384  will acquire complex values.
385  The solution 
386 \begin_inset Formula $\outcoeffp{\vect 0}\left(\vect k\right)$
387 \end_inset
389  is then obtained as the right 
390 \begin_inset Note Note
391 status open
393 \begin_layout Plain Layout
394 CHECK!
395 \end_layout
397 \end_inset
399  singular vector of 
400 \begin_inset Formula $M\left(\omega,\vect k\right)$
401 \end_inset
403  corresponding to the zero singular value.
404 \end_layout
406 \begin_layout Subsection
407 Numerical solution
408 \end_layout
410 \begin_layout Standard
411 In practice, equation 
412 \begin_inset CommandInset ref
413 LatexCommand eqref
414 reference "eq:Multiple-scattering problem unit cell block form"
415 plural "false"
416 caps "false"
417 noprefix "false"
419 \end_inset
421  is solved in the same way as eq.
423 \begin_inset CommandInset ref
424 LatexCommand eqref
425 reference "eq:Multiple-scattering problem block form"
426 plural "false"
427 caps "false"
428 noprefix "false"
430 \end_inset
432  in the multipole degree truncated form.
433 \end_layout
435 \begin_layout Standard
436 The lattice mode problem 
437 \begin_inset CommandInset ref
438 LatexCommand eqref
439 reference "eq:lattice mode equation"
440 plural "false"
441 caps "false"
442 noprefix "false"
444 \end_inset
446  is (after multipole degree truncation) solved by finding 
447 \begin_inset Formula $\omega,\vect k$
448 \end_inset
450  for which the matrix 
451 \begin_inset Formula $M\left(\omega,\vect k\right)$
452 \end_inset
454  has a zero singular value.
455  A naïve approach to do that is to sample a volume with a grid in the 
456 \begin_inset Formula $\left(\omega,\vect k\right)$
457 \end_inset
459  space, performing a singular value decomposition of 
460 \begin_inset Formula $M\left(\omega,\vect k\right)$
461 \end_inset
463  at each point and finding where the lowest singular value of 
464 \begin_inset Formula $M\left(\omega,\vect k\right)$
465 \end_inset
467  is close enough to zero.
468  However, this approach is quite expensive, for 
469 \begin_inset Formula $W\left(\omega,\vect k\right)$
470 \end_inset
472  has to be evaluated for each 
473 \begin_inset Formula $\omega,\vect k$
474 \end_inset
476  pair separately (unlike the original finite case 
477 \begin_inset CommandInset ref
478 LatexCommand eqref
479 reference "eq:Multiple-scattering problem block form"
480 plural "false"
481 caps "false"
482 noprefix "false"
484 \end_inset
486  translation operator 
487 \begin_inset Formula $\trops$
488 \end_inset
490 , which, for a given geometry, depends only on frequency).
491  Therefore, a much more efficient approach to determine the photonic bands
492  is to sample the 
493 \begin_inset Formula $\vect k$
494 \end_inset
496 -space (a whole Brillouin zone or its part) and for each fixed 
497 \begin_inset Formula $\vect k$
498 \end_inset
500  to find a corresponding frequency 
501 \begin_inset Formula $\omega$
502 \end_inset
504  with zero singular value of 
505 \begin_inset Formula $M\left(\omega,\vect k\right)$
506 \end_inset
508  using a minimisation algorithm (two- or one-dimensional, depending on whether
509  one needs the exact complex-valued 
510 \begin_inset Formula $\omega$
511 \end_inset
513  or whether the its real-valued approximation is satisfactory).
514  Typically, a good initial guess for 
515 \begin_inset Formula $\omega\left(\vect k\right)$
516 \end_inset
518  is obtained from the empty lattice approximation, 
519 \begin_inset Formula $\left|\vect k\right|=\sqrt{\epsilon\mu}\omega/c_{0}$
520 \end_inset
522  (modulo lattice points; TODO write this a clean way).
523  A somehow challenging step is to distinguish the different bands that can
524  all be very close to the empty lattice approximation, especially if the
525  particles in the systems are small.
526  In high-symmetry points of the Brilloin zone, this can be solved by factorising
528 \begin_inset Formula $M\left(\omega,\vect k\right)$
529 \end_inset
531  into irreducible representations 
532 \begin_inset Formula $\Gamma_{i}$
533 \end_inset
535  and performing the minimisation in each irrep separately, cf.
536  Section 
537 \begin_inset CommandInset ref
538 LatexCommand ref
539 reference "sec:Symmetries"
540 plural "false"
541 caps "false"
542 noprefix "false"
544 \end_inset
546 , and using the different 
547 \begin_inset Formula $\omega_{\Gamma_{i}}\left(\vect k\right)$
548 \end_inset
550  to obtain the initial guesses for the nearby points 
551 \begin_inset Formula $\vect k+\delta\vect k$
552 \end_inset
555 \end_layout
557 \begin_layout Subsection
558 Computing the Fourier sum of the translation operator
559 \begin_inset CommandInset label
560 LatexCommand label
561 name "subsec:W operator evaluation"
563 \end_inset
566 \end_layout
568 \begin_layout Standard
569 The problem evaluating 
570 \begin_inset CommandInset ref
571 LatexCommand eqref
572 reference "eq:W definition"
574 \end_inset
576  is the asymptotic behaviour of the translation operator, 
577 \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|}$
578 \end_inset
580  that does not in the strict sense converge for any 
581 \begin_inset Formula $d>1$
582 \end_inset
584 -dimensional lattice.
585 \begin_inset Note Note
586 status open
588 \begin_layout Plain Layout
589 \begin_inset Foot
590 status open
592 \begin_layout Plain Layout
593 Note that 
594 \begin_inset Formula $d$
595 \end_inset
597  here is dimensionality of the lattice, not the space it lies in, which
598  I for certain reasons assume to be three.
599  (TODO few notes on integration and reciprocal lattices in some appendix)
600 \end_layout
602 \end_inset
605 \end_layout
607 \end_inset
609  In electrostatics, this problem can be solved with Ewald summation [TODO
610  REF].
611  Its basic idea is that if what asymptoticaly decays poorly in the direct
612  space, will perhaps decay fast in the Fourier space.
613  We use the same idea here, but the technical details are more complicated
614  than in electrostatics.
615 \end_layout
617 \begin_layout Standard
618 Let us re-express the sum in 
619 \begin_inset CommandInset ref
620 LatexCommand eqref
621 reference "eq:W definition"
623 \end_inset
625  in terms of integral with a delta comb
626 \begin_inset FormulaMacro
627 \renewcommand{\basis}[1]{\mathfrak{#1}}
628 \end_inset
631 \end_layout
633 \begin_layout Standard
634 \begin_inset Formula 
635 \begin{equation}
636 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}
637 \end{equation}
639 \end_inset
641 The translation operator 
642 \begin_inset Formula $S$
643 \end_inset
645  is now a function defined in the whole 3d space; 
646 \begin_inset Formula $\vect r_{\alpha},\vect r_{\beta}$
647 \end_inset
649  are the displacements of scatterers 
650 \begin_inset Formula $\alpha$
651 \end_inset
653  and 
654 \begin_inset Formula $\beta$
655 \end_inset
657  in a unit cell.
658  The arrow notation 
659 \begin_inset Formula $S(\vect r_{\alpha}\leftarrow\vect r+\vect r_{\beta})$
660 \end_inset
662  means 
663 \begin_inset Quotes eld
664 \end_inset
666 translation operator for spherical waves originating in 
667 \begin_inset Formula $\vect r+\vect r_{\beta}$
668 \end_inset
670  evaluated in 
671 \begin_inset Formula $\vect r_{\alpha}$
672 \end_inset
675 \begin_inset Quotes erd
676 \end_inset
678  and obviously 
679 \begin_inset Formula $S$
680 \end_inset
682  is in fact a function of a single 3d argument, 
683 \begin_inset Formula $S(\vect r_{\alpha}\leftarrow\vect r+\vect r_{\beta})=S(\vect 0\leftarrow\vect r+\vect r_{\beta}-\vect r_{\alpha})=S(-\vect r-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)=S(-\vect r-\vect r_{\beta}+\vect r_{\alpha})$
684 \end_inset
687  Expression 
688 \begin_inset CommandInset ref
689 LatexCommand eqref
690 reference "eq:W integral"
692 \end_inset
694  can be rewritten as
695 \begin_inset Formula 
697 W_{\alpha\beta}(\vect k)=\left(2\pi\right)^{\frac{d}{2}}\uaft{(\dc{\basis u}S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0))\left(\vect k\right)}
700 \end_inset
702 where changed the sign of 
703 \begin_inset Formula $\vect r/\vect{\bullet}$
704 \end_inset
706  has been swapped under integration, utilising evenness of 
707 \begin_inset Formula $\dc{\basis u}$
708 \end_inset
711  Fourier transform of product is convolution of Fourier transforms, so (using
712  formula 
713 \begin_inset CommandInset ref
714 LatexCommand eqref
715 reference "eq:Dirac comb uaFt"
717 \end_inset
719  for the Fourier transform of Dirac comb)
720 \begin_inset Formula 
721 \begin{eqnarray}
722 W_{\alpha\beta}(\vect k) & = & \left(\left(\uaft{\dc{\basis u}}\right)\ast\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\right)(\vect k)\nonumber \\
723  & = & \frac{\left|\det\recb{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\left(\dc{\recb{\basis u}}^{(d)}\ast\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\right)\left(\vect k\right)\nonumber \\
724  & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\sum_{\vect K\in\recb{\basis u}\ints^{d}}\left(\uaft{S(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\left(\vect k-\vect K\right)\label{eq:W sum in reciprocal space}\\
725  & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\sum_{\vect K\in\recb{\basis u}\ints^{d}}e^{i\left(\vect k-\vect K\right)\cdot\left(-\vect r_{\beta}+\vect r_{\alpha}\right)}\left(\uaft{S(\vect{\bullet}\leftarrow\vect 0)}\right)\left(\vect k-\vect K\right)\nonumber 
726 \end{eqnarray}
728 \end_inset
731 \begin_inset Note Note
732 status open
734 \begin_layout Plain Layout
735 Factor 
736 \begin_inset Formula $\left(2\pi\right)^{\frac{d}{2}}$
737 \end_inset
739  cancels out with the 
740 \begin_inset Formula $\left(2\pi\right)^{-\frac{d}{2}}$
741 \end_inset
743  factor appearing in the convolution/product formula in the unitary angular
744  momentum convention.
746 \end_layout
748 \end_inset
750 As such, this is not extremely helpful because the the 
751 \emph on
752 whole
753 \emph default
754  translation operator 
755 \begin_inset Formula $S$
756 \end_inset
758  has singularities in origin, hence its Fourier transform 
759 \begin_inset Formula $\uaft S$
760 \end_inset
762  will decay poorly.
764 \end_layout
766 \begin_layout Standard
767 However, Fourier transform is linear, so we can in principle separate 
768 \begin_inset Formula $S$
769 \end_inset
771  in two parts, 
772 \begin_inset Formula $S=S^{\textup{L}}+S^{\textup{S}}$
773 \end_inset
777 \begin_inset Formula $S^{\textup{S}}$
778 \end_inset
780  is a short-range part that decays sufficiently fast with distance so that
781  its direct-space lattice sum converges well; 
782 \begin_inset Formula $S^{\textup{S}}$
783 \end_inset
785  must as well contain all the singularities of 
786 \begin_inset Formula $S$
787 \end_inset
789  in the origin.
790  The other part, 
791 \begin_inset Formula $S^{\textup{L}}$
792 \end_inset
794 , will retain all the slowly decaying terms of 
795 \begin_inset Formula $S$
796 \end_inset
798  but it also has to be smooth enough in the origin, so that its Fourier
799  transform 
800 \begin_inset Formula $\uaft{S^{\textup{L}}}$
801 \end_inset
803  decays fast enough.
804  (The same idea lies behind the Ewald summation in electrostatics.) Using
805  the linearity of Fourier transform and formulae 
806 \begin_inset CommandInset ref
807 LatexCommand eqref
808 reference "eq:W definition"
810 \end_inset
812  and 
813 \begin_inset CommandInset ref
814 LatexCommand eqref
815 reference "eq:W sum in reciprocal space"
817 \end_inset
819 , the operator 
820 \begin_inset Formula $W_{\alpha\beta}$
821 \end_inset
823  can then be re-expressed as
824 \begin_inset Formula 
825 \begin{eqnarray}
826 W_{\alpha\beta}\left(\vect k\right) & = & W_{\alpha\beta}^{\textup{S}}\left(\vect k\right)+W_{\alpha\beta}^{\textup{L}}\left(\vect k\right)\nonumber \\
827 W_{\alpha\beta}^{\textup{S}}\left(\vect k\right) & = & \sum_{\vect R\in\basis u\ints^{d}}S^{\textup{S}}(\vect 0\leftarrow\vect R+\vect r_{\beta}-\vect r_{\alpha})e^{i\vect k\cdot\vect R}\label{eq:W Short definition}\\
828 W_{\alpha\beta}^{\textup{L}}\left(\vect k\right) & = & \frac{\left|\det\rec{\basis u}\right|}{\left(2\pi\right)^{\frac{d}{2}}}\sum_{\vect K\in\recb{\basis u}\ints^{d}}\left(\uaft{S^{\textup{L}}(\vect{\bullet}-\vect r_{\beta}+\vect r_{\alpha}\leftarrow\vect 0)}\right)\left(\vect k-\vect K\right)\label{eq:W Long definition}
829 \end{eqnarray}
831 \end_inset
833 where both sums expected to converge nicely.
834  We note that the elements of the translation operators 
835 \begin_inset Formula $\tropr,\trops$
836 \end_inset
838  in 
839 \begin_inset CommandInset ref
840 LatexCommand eqref
841 reference "eq:translation operator"
842 plural "false"
843 caps "false"
844 noprefix "false"
846 \end_inset
848  can be rewritten as linear combinations of expressions 
849 \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)$
850 \end_inset
852  (TODO WRITE THEM EXPLICITLY IN THIS FORM), respectively, hence if we are
853  able evaluate the lattice sums sums
854 \begin_inset Note Note
855 status open
857 \begin_layout Plain Layout
858 CHECK THE FOLLOWING EXPRESSION FOR CORRECT FUNCTION ARGUMENTS
859 \end_layout
861 \end_inset
864 \begin_inset Formula 
865 \begin{equation}
866 \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}
867 \end{equation}
869 \end_inset
871 then by linearity, we can get the 
872 \begin_inset Formula $W_{\alpha\beta}\left(\vect k\right)$
873 \end_inset
875  operator as well.
876 \end_layout
878 \begin_layout Standard
879 TODO ADD MOROZ AND OTHER REFS HERE.
881 \begin_inset CommandInset citation
882 LatexCommand cite
883 key "linton_one-_2009"
884 literal "true"
886 \end_inset
888  offers an exponentially convergent Ewald-type summation method for 
889 \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)$
890 \end_inset
893  Here we rewrite them in a form independent on the convention used for spherical
894  harmonics (as long as they are complex
895 \begin_inset Note Note
896 status open
898 \begin_layout Plain Layout
899 lepší formulace
900 \end_layout
902 \end_inset
905  The short-range part reads (UNIFY INDEX NOTATION)
906 \begin_inset Formula 
907 \begin{multline}
908 \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\\
909 +\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}
910 \end{multline}
912 \end_inset
915 \begin_inset Note Note
916 status open
918 \begin_layout Plain Layout
919 NEPATŘÍ TAM NĚJAKÁ DELTA FUNKCE K PŮVODNÍMU 
920 \begin_inset Formula $\sigma_{n}^{m(0)}$
921 \end_inset
924 \end_layout
926 \end_inset
928 and the long-range part (FIXME, this is the 2D version; include the 1D and
929  3D lattice expressions as well)
930 \begin_inset Formula 
931 \begin{multline}
932 \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\\
933 \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}
934 \end{multline}
936 \end_inset
938 where 
939 \begin_inset Formula $\xi$
940 \end_inset
942  is TODO, 
943 \begin_inset Formula $\beta_{pq}$
944 \end_inset
946  is TODO, 
947 \begin_inset Formula $\Gamma_{j,pq}$
948 \end_inset
950  is TODO and 
951 \begin_inset Formula $\eta$
952 \end_inset
954  is a real parameter that determines the pace of convergence of both parts.
955  The larger 
956 \begin_inset Formula $\eta$
957 \end_inset
959  is, the faster 
960 \begin_inset Formula $\sigma_{n}^{m(\mathrm{S})}\left(\vect{\beta}\right)$
961 \end_inset
963  converges but the slower 
964 \begin_inset Formula $\sigma_{n}^{m(\mathrm{L})}\left(\vect{\beta}\right)$
965 \end_inset
967  converges.
968  Therefore (based on the lattice geometry) it has to be adjusted in a way
969  that a reasonable amount of terms needs to be evaluated numerically from
970  both 
971 \begin_inset Formula $\sigma_{n}^{m(\mathrm{S})}\left(\vect{\beta}\right)$
972 \end_inset
974  and 
975 \begin_inset Formula $\sigma_{n}^{m(\mathrm{L})}\left(\vect{\beta}\right)$
976 \end_inset
979  Generally, a good choice for 
980 \begin_inset Formula $\eta$
981 \end_inset
983  is TODO; in order to achieve accuracy TODO, one has to evaluate the terms
984  on TODO lattice points.
985  (I HAVE SOME DERIVATIONS OF THE ESTIMATES IN MY NOTES; SHOULD I INCLUDE
986  THEM?)
987 \end_layout
989 \begin_layout Standard
990 In practice, the integrals in 
991 \begin_inset CommandInset ref
992 LatexCommand eqref
993 reference "eq:Ewald in 3D short-range part"
994 plural "false"
995 caps "false"
996 noprefix "false"
998 \end_inset
1000  can be easily evaluated by numerical quadrature and the incomplete 
1001 \begin_inset Formula $\Gamma$
1002 \end_inset
1004 -functions using the series TODO and TODO from DLMF.
1005 \end_layout
1007 \end_body
1008 \end_document