Лицензия MIT
© Алексей Александрович Щербаков, 2024
В данной лекции будем рассматривать открытые оптические системы - резонаторы и волноводы. Для краткости перепишем уравнения Максвелла в частотной области в операторной форме \begin{equation}\tag{1} \hat{\mathcal{L}}(\boldsymbol{r},\omega) \boldsymbol{F}(\boldsymbol{r},\omega) = \boldsymbol{J}(\boldsymbol{r},\omega) \end{equation} где $\hat{\mathcal{L}} = \omega\hat{\mathcal{P}} - \hat{\mathcal{D}}$, так что \begin{equation*} \hat{\mathcal{P}}(\boldsymbol{r},\omega) = \left( \begin{array}{cc} \varepsilon & -i\xi \\ -i\xi^T & \mu \end{array} \right), \;\; \hat{\mathcal{D}} = \left( \begin{array}{cc} 0 & \nabla\times \\ \nabla\times & 0 \end{array} \right) \end{equation*} Вектор поля включает электрическую и магнитную компоненты $\boldsymbol{F}(\boldsymbol{r},\omega) = (\boldsymbol{E},\;i\boldsymbol{H})^T$ и составной вектор тока $\boldsymbol{J}(\boldsymbol{r},\omega) = (-i\boldsymbol{J}_e,\;\boldsymbol{J}_m)^T$. Тогда можно обобщить запись уравнения на функцию Грина \begin{equation}\tag{2} \hat{\mathcal{L}}(\boldsymbol{r},\omega) \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega) = \hat{\mathbb{I}} \delta(\boldsymbol{r},\boldsymbol{r}') \end{equation} и интегральное решение уравнений Максвелла \begin{equation}\tag{3} \boldsymbol{F}(\boldsymbol{r},\omega) = \int dV' \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega) \boldsymbol{J}(\boldsymbol{r}',\omega) \end{equation}
Если возбуждение открытой системы происходит в течение короткого времени в момент $t=0$, то временной отклик системы запишется как \begin{equation*} \boldsymbol{F}(\boldsymbol{r},t) \sim \int d\omega \boldsymbol{F}(\boldsymbol{r},\omega) e^{-i\omega t},\;t>0 \end{equation*} где $\boldsymbol{F}(\boldsymbol{r},\omega)$ - решение уравнения в частотной области для заданного пространственного распределения источников в нулевой момент времени. Предполагая, что анлитическое продолжение функции Грина, и, соответственно, вектора полей имеет счетное количество полюсов $\omega_n = \Omega_n - i\Gamma_n$, $\Gamma_n > 0$, основная теорема о вычетах позволяет переписать решение как \begin{equation*} \boldsymbol{F}(\boldsymbol{r},t) \sim \sum_n \mathrm{Res}[\boldsymbol{F}(\boldsymbol{r},\omega),\omega_n] e^{-i\omega_n t},\;t>0 \end{equation*} так что затухание каждой компоненты определяется множителем $\exp(-\Gamma_n t)$.
Поскольку правая часть уравнения на функцию Грина $(2)$ не зависит от частоты, то полюса функции Грина должны быть собственными числами уравнения \begin{equation}\tag{4} \hat{\mathcal{L}}(\boldsymbol{r},\omega_n) \boldsymbol{F}_n(\boldsymbol{r}) = 0 \end{equation} то есть, уравнения на собственные значения (моды) уравнений Максвелла в заданном открытом резонаторе или волноводе. Эти числа, вообще говоря, комплексные. Вне рассматриваемого резонатора/волновода поле представляется в виде уходящих волн, которые в трехмерном случае имеют множитель $\exp(-i\omega_n t + ik_nr)$ где волновое число свободного пространства, окружающего резонатор $k_n=\omega_n\sqrt{\varepsilon\mu}$. Отсюда видно, что затухающие во времени состояния открытых систем описываются экспоненциально расходящимися пространственными функциями, поэтому при разложении по таким состояниям нужно руководствоваться правилами, учитывающими данную специфику. Такие состояния в литературе называются квазинормальными модами или резонансными состояниями.
Уравнение $(4)$ должно быть дополнено необходимыми граничными условиями. Решая это уравнение для комплексных частот имеют ввиду, что любые источники в заданной конфигурации рассматриваемых резонаторов локализованы, так что поле вне некоторого конечного объема имеет вид уходящих волн. Математически можно рассматривать уравнение на собственные значения как предельный случай общего уравнения $(1)$ с источником $\boldsymbol{J}(\boldsymbol{r},\omega) \sim \omega-\omega_n$ при $\omega\rightarrow\omega_n$. При этом $\boldsymbol{F}(\boldsymbol{r},\omega)\rightarrow \boldsymbol{F}_n(\boldsymbol{r})$.
Разложение функции Грина по резнансным состояниям основано на теореме Миттаг-Леффлера: если функция $f(z)$ является аналитической всюду, кроме счетного числа полюсов $z_n$ с вычетами $R_n$, и $\lim_{z\rightarrow\infty}f(z)/z^p=0$, то эта функция может быть представлена в виде \begin{equation*} f(z) = f_p(z) + \sum_n\frac{R_n}{z-z_n} \end{equation*} где $f_p$ - многочлен степени $p-1$. Применим эту теорему для записи функции Грина, учитывая, что $\lim_{\omega\rightarrow\infty}\mathcal{G}=0$ \begin{equation}\tag{5} \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega) = \sum_n \frac{\mathcal{Q}_n(\boldsymbol{r},\boldsymbol{r}')}{\omega-\omega_n} \end{equation} Подставим это выражение в уравнение на функцию Грина $(2)$, умножим на некотрое поле $\boldsymbol{F}_b(\boldsymbol{r})$, отличное от нуля в конечном объеме, и проинтегрируем обе части уравнения, что даст \begin{equation*} \sum_n \frac{\hat{\mathcal{L}}(\boldsymbol{r},\omega)\boldsymbol{A}_n(\boldsymbol{r})}{\omega-\omega_n} = \boldsymbol{F}_b(\boldsymbol{r}) \end{equation*} где \begin{equation*} \boldsymbol{A}_n(\boldsymbol{r}) = \intop_{V'} dV' \mathcal{Q}_n(\boldsymbol{r},\boldsymbol{r}') \boldsymbol{F}_b(\boldsymbol{r}') \end{equation*} В пределе $\omega\rightarrow\omega_m$ получаем систему независимых уравнений $\hat{\mathcal{L}}(\boldsymbol{r},\omega_n)\boldsymbol{A}_n(\boldsymbol{r}) = 0$, откуда следует, что получившиеся векторные амплитуды пропорциональны резонансным состояниям, задваемым уравнением $(4)$. А, поскольку выбор поля $\boldsymbol{F}_b(\boldsymbol{r})$ произволен в рамках оговоренных ограничений, то необходимо потребовать, чтобы применение оператора Грина проецировало на $n$-е резонансное состояние. \begin{equation*} \mathcal{Q}_n(\boldsymbol{r},\boldsymbol{r}') \sim \boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{X}(\boldsymbol{r'}) \end{equation*} Чтобы найти вектор $\boldsymbol{X}$, воспользуемся соотношением взаимности $\mathcal{G}(\boldsymbol{r}, \boldsymbol{r}',\omega) = \mathcal{G}^T (\boldsymbol{r}',\boldsymbol{r},\omega)$. Отсюда сразу получается, что $\boldsymbol{X}\sim\boldsymbol{F}_n(\boldsymbol{r}')$. Таким образом, выбирая соответствующую нормировку резонансных состояний, можно записать \begin{equation}\tag{6} \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega) = \sum_n \frac{\boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r'})}{\omega-\omega_n} \end{equation}
В ряде задач представление $(6)$ оказывается недостаточным, так что необходимо учесть разрезы на комплексной плоскости и выбор подходящих римновых поверхностей, например, за счет возникновения квадратного корня в дисперсионном уравнении для плоских волн. В более общем виде разложение тензора Грина записывается как \begin{equation}\tag{6a} \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega) = \sum_n \frac{\boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r'})}{\omega-\omega_n} + \sum_p \mathcal{G}_p(\boldsymbol{r},\boldsymbol{r}',\omega) \end{equation} где индекс $p$ соответствует вкладу отдельных разрезов: \begin{equation*} \mathcal{G}_p(\boldsymbol{r},\boldsymbol{r}',\omega) = \frac{1}{2\pi i} \int_{C_p} d\omega' \frac{\Delta \mathcal{G}(\boldsymbol{r},\boldsymbol{r}',\omega')}{\omega-\omega'} \end{equation*} Здесь $C_m$ - путь вдоль соответствущего разреза, и $\Delta \mathcal{G}$ - изменение значения функции Грина на различных римановых поверхностях. Необходимо отметить, что в численных решениях вклад разрезов может быть также записан в форме $(6)$.
Подставим полученное разложение функции Грина в объемное интегральное уравнение для вектора поля: \begin{equation}\tag{7} \boldsymbol{F}(\boldsymbol{r},\omega) = \sum_n \frac{\boldsymbol{F}_n(\boldsymbol{r})}{\omega-\omega_n} \int dV' \boldsymbol{F}_n(\boldsymbol{r'}) \cdot \boldsymbol{J}(\boldsymbol{r}',\omega) \end{equation} Как обсуждалось выше $\boldsymbol{F}(\boldsymbol{r},\omega) \underset{\omega\rightarrow\omega_{n}}{\rightarrow}\boldsymbol{F}_{n}(\boldsymbol{r})$ для источников вида $\boldsymbol{J}(\boldsymbol{r},\omega) = (\omega-\omega_n) \tilde{\boldsymbol{J}}(\boldsymbol{r})$. Тогда переходя к пределу $\omega\rightarrow\omega_n$ в уравнении $(7)$, получаем \begin{equation}\tag{8} \int dV' \boldsymbol{F}_n(\boldsymbol{r'}) \cdot \tilde{\boldsymbol{J}}(\boldsymbol{r}') = 1 \end{equation} Теперь возьмем уравнения $(1)$ и $(4)$ для поля некоторых локализованных источников и резонансных состояний, домножим первое на состояние $\boldsymbol{F}^{\ddagger}_{n}(\boldsymbol{r})$, а второе на решение для данных источников $\boldsymbol{F}(\omega,\boldsymbol{r})$, учтем явное представление оператора, и сложим: \begin{align*} \omega \boldsymbol{F}_{n}(\boldsymbol{r}) \cdot \hat{\mathcal{P}}(\boldsymbol{r},\omega) \boldsymbol{F}(\boldsymbol{r},\omega) - \omega_n \boldsymbol{F}(\boldsymbol{r},\omega) \cdot \hat{\mathcal{P}}(\boldsymbol{r},\omega_n) \boldsymbol{F}_{n}(\boldsymbol{r}) + \\ + \boldsymbol{F}(\boldsymbol{r},\omega) \cdot \hat{\mathcal{D}}(\boldsymbol{r},\omega_n) \boldsymbol{F}_{n}(\boldsymbol{r}) - \boldsymbol{F}_{n}(\boldsymbol{r}) \cdot \hat{\mathcal{D}}(\boldsymbol{r},\omega_n) \boldsymbol{F}(\boldsymbol{r},\omega) = \\ = (\omega-\omega_n) \boldsymbol{F}_{n}(\boldsymbol{r}) \cdot \tilde{\boldsymbol{J}}(\boldsymbol{r}) \end{align*} Подставляя во вторую строчку явный вид векторов поля и оператора $\hat{\mathcal{D}}$, и используя разложение дивергенции векторного произведения, получаем, \begin{equation*} \boldsymbol{F}(\boldsymbol{r},\omega) \cdot \hat{\mathcal{D}}(\boldsymbol{r},\omega_n) \boldsymbol{F}_{n}(\boldsymbol{r}) - \boldsymbol{F}_{n}(\boldsymbol{r}) \cdot \hat{\mathcal{D}}(\boldsymbol{r},\omega_n) \boldsymbol{F}(\boldsymbol{r},\omega) = i\nabla \cdot (\boldsymbol{E}_{n} \times \boldsymbol{H} - \boldsymbol{E} \times \boldsymbol{H}_n) \end{equation*} Проинтегрируем общее равенство по некотрому конечному объёму, включающему локализованные источники и перейдём к пределу $\omega\rightarrow\omega_n$. При этом воспользуемся разложением аналитических векторов поля и произведения $\omega \hat{\mathcal{P}}(\boldsymbol{r},\omega)$ продолжений вблизи точки $\omega=\omega_n$: \begin{equation*} \omega \hat{\mathcal{P}}(\boldsymbol{r},\omega) \approx \omega_n \hat{\mathcal{P}}(\boldsymbol{r},\omega_n) + (\omega-\omega_n) {[\omega \hat{\mathcal{P}}(\boldsymbol{r},\omega)]'}_{\omega=\omega_n} \end{equation*}
\begin{equation*} \boldsymbol{E}(\boldsymbol{r},\omega) \approx \boldsymbol{E}_n(\boldsymbol{r}) + (\omega-\omega_n) {[\boldsymbol{E}(\boldsymbol{r},\omega)]'}_{\omega=\omega_n} = \boldsymbol{E}_n(\boldsymbol{r}) + (\omega-\omega_n) \frac{1}{\omega_n} (\boldsymbol{r}\nabla) \boldsymbol{E}_n(\boldsymbol{r}) \end{equation*}Здесь штрих обозначает производную по частоте, а последнее равенство получается, если учесть, что на границе указанного объема поле представимо в виде разложения по уходящим на бесконечность волнам и имеет место соотношение $\boldsymbol{E}(\boldsymbol{r},\omega) = \boldsymbol{E}(\omega\boldsymbol{r})$. В итоге учитывая соотношение $(8)$ получаем \begin{equation}\tag{9} V_n + S_n = 1 \end{equation} где \begin{equation}\tag{10} V_n = \intop_V dV' \boldsymbol{F}_n(\boldsymbol{r}') \cdot {[\omega \hat{\mathcal{P}}(\boldsymbol{r}',\omega)]'}_{\omega=\omega_n} \boldsymbol{F}_n(\boldsymbol{r}') \end{equation}
\begin{equation}\tag{11} S_n = i \oint_{\partial V} d\boldsymbol{\Sigma} (\boldsymbol{E}_{n} \times \boldsymbol{H}'_n - \boldsymbol{E}'_n \times \boldsymbol{H}_n) \end{equation}На практике производная по частоте поля в резонансном состоянии может быть вычислена следующим образом. Вернемся к уравнению $(7)$ и рассмотрим его в окрестности состояния $m$, подставив соответствующее представление источника: \begin{equation*} \boldsymbol{F}(\boldsymbol{r},\omega) = \sum_n \frac{\omega-\omega_m}{\omega-\omega_n} \boldsymbol{F}_n(\boldsymbol{r}) \int dV' \boldsymbol{F}_n(\boldsymbol{r'}) \cdot \tilde{\boldsymbol{J}}_m(\boldsymbol{r}',\omega) = \sum_n \frac{\omega-\omega_m}{\omega-\omega_n} \boldsymbol{F}_n(\boldsymbol{r}) I_{nm} \end{equation*} Здесь $I_{mm} = 1$ согласно полученному выше соотношению. Поле вне резонансной системы раскладывается по уходящим векторным волнам. Запишем такое разложение для каждого резонансного сосотояния $\boldsymbol{F}_n(\boldsymbol{r}) = \sum_L \alpha_{Ln} \boldsymbol{O}(\boldsymbol{r},\omega)$ и продифференцируем по частоте в точке $\omega=\omega_m$: \begin{equation}\tag{12} \boldsymbol{F}_m'(\boldsymbol{r}) = \sum_L [\alpha_{Lm} \boldsymbol{O}'(\boldsymbol{r},\omega_m) + \alpha_{Lm}' \boldsymbol{O}(\boldsymbol{r},\omega_m)] \end{equation} При подстановке в соотношение ортогональности второе слагаемое не вносит вклада по следующей причине. Рассмотрим поля $\boldsymbol{F}_{1,2}(\boldsymbol{r},\omega)$, возбуждаемые локализованными в области резонатора источниками $\boldsymbol{J}_{1,2}(\boldsymbol{r},\omega)$. Тогда для Лоренц-взаимной среды \begin{equation*} \boldsymbol{F}_1 \cdot \hat{\mathcal{D}} \boldsymbol{F}_2 - \boldsymbol{F}_2 \cdot \hat{\mathcal{D}} \boldsymbol{F}_1 = \boldsymbol{F}_2 \cdot \boldsymbol{J}_1 - \boldsymbol{F}_1 \cdot \boldsymbol{J}_2 \end{equation*} Интегрирование по объему и применение теоремы Гаусса даёт \begin{equation*} i \oint_{\partial V} d\boldsymbol{\Sigma} (\boldsymbol{E}_2\times\boldsymbol{H}_1-\boldsymbol{E}_1\times\boldsymbol{H}_2) = \int_V dV (\boldsymbol{F}_2 \cdot \boldsymbol{J}_1 - \boldsymbol{F}_1 \cdot \boldsymbol{J}_2) \end{equation*} Подстановка в правую часть решений для поля через функцию Грина и применение теоремы взаимности даёт \begin{equation}\tag{13} i \oint_{\partial V} d\boldsymbol{\Sigma} (\boldsymbol{E}_2\times\boldsymbol{H}_1-\boldsymbol{E}_1\times\boldsymbol{H}_2) = 0 \end{equation} Это уравнение показывает, что при подстановке $(12)$ в $(9)$ второе слагаемое обнуляется, и достаточно вычислять производные уходящих векторных волн по частоте, что в случае плоских, сферических или цилиндрических волн можно сделать аналитически.
Производя вывод, такой же, как был использован для равенства $(13)$ для двух резонансных состояний $\boldsymbol{F}_{m,n}(\boldsymbol{r})$, получаем второе уравнение, выражающее ортогональность состояний: \begin{equation}\tag{14} \int_V dV \boldsymbol{F}_m \cdot [\omega_n \hat{\mathcal{P}}(\omega_n) - \omega_m \hat{\mathcal{P}}(\omega_m)] \boldsymbol{F}_n + i \oint_{\partial V} d\boldsymbol{\Sigma} (\boldsymbol{E}_{m} \times \boldsymbol{H}_n - \boldsymbol{E}_n \times \boldsymbol{H}_m) = 0 \end{equation}
Переходя к вопросу о полноте системы резонансных состояний, подставим разоложение $(6)$ в уравнение для функции Грина $(2)$ и используем явный вид оператора и уравнение на состояния $(4)$. После несложных преобразований получаем \begin{equation}\tag{15} \sum_n \frac{\omega \hat{\mathcal{P}}(\boldsymbol{r},\omega) - \omega_n \hat{\mathcal{P}}(\boldsymbol{r},\omega_n)}{\omega-\omega_n} \boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r}') = \hat{\mathbb{I}} \delta(\boldsymbol{r}-\boldsymbol{r}') \end{equation} Полученное равенство выражает замыкание системы резонансных состояний и является необходимым условием полноты. В случае отсутствия дисперсии $\hat{\mathcal{P}}(\boldsymbol{r},\omega)=\hat{\mathcal{P}}(\boldsymbol{r})$ у материалов оно упрощается \begin{equation}\tag{15a} \hat{\mathcal{P}}(\boldsymbol{r}) \sum_n \boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r}') = \hat{\mathbb{I}} \delta(\boldsymbol{r}-\boldsymbol{r}') \end{equation} В более реалистичном случае дисперсии Друде-Лоренца, когда \begin{equation*} \hat{\mathcal{P}}(\boldsymbol{r},\omega) = \hat{\mathcal{P}}_{\infty} (\boldsymbol{r}) + \sum_{i} \frac{\hat{\mathcal{P}}_{i}(\boldsymbol{r})}{\omega-\Omega_{i}} \end{equation*} линейная независимость лоренцианов приводит следующему выражению замыкания и дополнительным правилам сумм: \begin{equation}\tag{15b} \hat{\mathcal{P}}_{\infty}(\boldsymbol{r}) \sum_n \boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r}') = \hat{\mathbb{I}} \delta(\boldsymbol{r}-\boldsymbol{r}') \end{equation}
\begin{equation*} \hat{\mathcal{P}}_{i}(\boldsymbol{r}) \sum_n \frac{\boldsymbol{F}_n(\boldsymbol{r}) \otimes \boldsymbol{F}_n(\boldsymbol{r}')}{\omega_n-\Omega_i} = \hat{\mathbb{I}} \delta(\boldsymbol{r}-\boldsymbol{r}') \end{equation*}Наличие такого набора соотношений может свидетельствовать об избыточности базиса, тем не менее, на практике, когда разложение производится по небольшому числу состояний, удается получить достаточную точностью расчетов. Также открытым остаётся вопрос о полноте при наличии в системе особых точек, в которых происходит одновременное вырождение и собственных частот, и собственных векторов.
Рассмотрим некоторый объём, поверхность которого охватыват рассматриваемый резонатор. На этой поверхности любое поле может быть разложено по полному набору приходящих и уходящих векторным волнам, $\boldsymbol{I}_{L}(\boldsymbol{r},\omega)$ и $\boldsymbol{O}_{L}(\boldsymbol{r},\omega)$: \begin{equation}\tag{16} \boldsymbol{F}(\boldsymbol{r},\omega) = \sum_L [\alpha^{inc}_L(\omega) \boldsymbol{I}_{L}(\boldsymbol{r},\omega) + \alpha^{out}_L(\omega) \boldsymbol{O}_{L}(\boldsymbol{r},\omega)] \end{equation} Набор волн может быть выбран так, чтобы удовлетворять соотношениям \begin{equation*} \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{inc}_{L}\times\boldsymbol{H}^{out}_{L'} - \boldsymbol{E}^{out}_{L'}\times\boldsymbol{H}^{inc}_{L}) = \delta_{LL'} \end{equation*}
\begin{equation*} \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{\sigma}_{L}\times\boldsymbol{H}^{\sigma}_{L'} - \boldsymbol{E}^{\sigma}_{L'}\times\boldsymbol{H}^{\sigma}_{L}) = 0,\;\sigma=inc,out \end{equation*}тогда коэффциенты \begin{equation*} \alpha^{inc}_L = -i \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{out}_{L}\times\boldsymbol{H} - \boldsymbol{E}\times\boldsymbol{H}^{out}_{L}) \end{equation*}
\begin{equation*} \alpha^{sca}_L = i \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{inc}_{L}\times\boldsymbol{H} - \boldsymbol{E}\times\boldsymbol{H}^{inc}_{L}) \end{equation*}Взаимодействие резонатора с падающим полем может быть описано в терминах матрицы рассеяния: \begin{equation}\tag{17} \alpha^{out}_L= \sum_{L'} S_{LL'} \alpha^{inc}_{L'} \end{equation}
Выразим в явном виде компоненты матрицы рассеяния, выделив в физической системе некоторый фоновый оператор $\hat{\mathcal{P}}_b$ и добавку $\Delta\hat{\mathcal{P}} = \hat{\mathcal{P}} - \hat{\mathcal{P}}_b$, описывающую локальное измерение материальных параметров за счет присутствия резонатора. Запишем полное поле в системе как сумму фонового и рассеянного $\boldsymbol{F}_{tot} = \boldsymbol{F}_b + \boldsymbol{F}_{sca}$, где фоновое поле состоит из приходящих и уходящих волн и удовлетворяет уравнению $\hat{\mathcal{L}}_b \boldsymbol{F}_b = 0$. Отсюда можно записать уравнение на рассеянное поле как \begin{equation}\tag{18} \hat{\mathcal{L}}(\boldsymbol{r},\omega) \boldsymbol{F}_{sca}(\boldsymbol{r},\omega) = -\omega \Delta\hat{\mathcal{P}}(\boldsymbol{r},\omega) \boldsymbol{F}_b \end{equation} Решение этого уравнения через объемное интегральное уравнение с функцией Грина, разложенной по резнансным состояниям, приводит к выражению \begin{equation}\tag{19} \boldsymbol{F}_{tot}(\boldsymbol{r},\omega) = \boldsymbol{F}_b(\boldsymbol{r},\omega) - \sum_n \frac{I_n(\omega)}{\omega-\omega_n}\boldsymbol{F}_n(\boldsymbol{r}) \end{equation} где интеграл перекрытия \begin{equation}\tag{20} I_n(\omega) = \omega \int dV \boldsymbol{F}_n(\boldsymbol{r}) \cdot \Delta\hat{\mathcal{P}}(\boldsymbol{r},\omega) \boldsymbol{F}_b(\boldsymbol{r},\omega) \end{equation} Это представление используется для резонансного разложения поля в ближней зоне. Для представления поля в дальней зоне перейдём к разложениям по векторным волнам: $\boldsymbol{F}_{tot,L} = \boldsymbol{F}_{b,L} + \boldsymbol{F}_{sca,L}$. Соотствествующим образом можно предствить и матрицу рассеяния: \begin{equation}\tag{21} S_{LL'}(\omega) = S^b_{LL'}(\omega) + S^{sca}_{LL'}(\omega) \end{equation} Используя выражения для коэффициентов приходящих и уходящих волн, можно получить \begin{equation}\tag{22} S^b_{LL'}(\omega) = i \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{inc}_{L}\times\boldsymbol{H}_{b,L'} - \boldsymbol{E}_{b,L'}\times\boldsymbol{H}^{inc}_{L}) \end{equation}
\begin{equation}\tag{23} S^{sca}_{LL'}(\omega) = i \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{inc}_{L}\times\boldsymbol{H}_{sca,L'} - \boldsymbol{E}_{sca,L'}\times\boldsymbol{H}^{inc}_{L}) \end{equation}Подставляя выражение через разность материальных операторов, приходим к выражению \begin{equation}\tag{24} S_{LL'}(\omega) = S^b_{LL'}(\omega) - \sum_n \frac{\alpha_{L}^{(n)}(\omega)I^b_{n,L'}(\omega)}{\omega-\omega_n^p} \end{equation} где \begin{equation*} \alpha^{(n)}_L = -i \oint_{\partial V} d\boldsymbol{\Sigma} \cdot (\boldsymbol{E}^{inc}_{L}\times\boldsymbol{H}_n - \boldsymbol{E}_n\times\boldsymbol{H}^{inc}_{L}) \end{equation*}
\begin{equation*} I^b_{n,L}(\omega) = \omega \int dV \boldsymbol{F}_n(\boldsymbol{r}) \cdot \Delta\hat{\mathcal{P}}(\boldsymbol{r},\omega) \boldsymbol{F}_{b,L}(\boldsymbol{r},\omega) \end{equation*}Одним из способов нахождения разложения по резонансным состояниям сложных систем является использование известных или легко вычислимых параметров резонансных состояний в более простых системах. По аналогии с разложением задачи на фоновую и резонансную, предположим, что оператор $\hat{\mathcal{L}}_0$ (и его материальная часть $\hat{\mathcal{P}}_0$) соответствует базовой задаче с известным разложение, а её отличие от рассматриваемой задается оператором $\delta\hat{\mathcal{P}}_0$. Тогда задача на собственные состояния перепишется как \begin{equation*} \hat{\mathcal{L}}_0(\boldsymbol{r},\omega_{\nu}) \boldsymbol{F}_{\nu}(\boldsymbol{r}) = -\omega_{\nu} \delta\hat{\mathcal{P}}_0(\boldsymbol{r},\omega_{\nu}) \boldsymbol{F}_{\nu}(\boldsymbol{r}) \end{equation*} Здесь моды сложной системы нумеруются индексом $\nu$, чтобы различать их с модами базисной простой системы. Решением уравнения будет \begin{equation*} \boldsymbol{F}_{\nu}(\boldsymbol{r}) = -\sum_n \frac{\omega_{\nu}\boldsymbol{F}_{n}(\boldsymbol{r})}{\omega_{\nu}-\omega_n} \int_V dV' \boldsymbol{F}_{n}(\boldsymbol{r}') \delta\hat{\mathcal{P}}_0(\boldsymbol{r}',\omega_{\nu}) \boldsymbol{F}_{\nu}(\boldsymbol{r}') \end{equation*} Будем искать неизвестные состояния в виде разложения $\boldsymbol{F}_{\nu} = \sum_n c_n \boldsymbol{F}_{n}$. Подставляя его в обе части последнего уравнения, приходим к системе алгебраических уравнений, выражающих задачу на собственные значения: \begin{equation}\tag{25} (\omega_{\nu}-\omega_n) c_n = -\omega_{\nu} \sum_m V_{nm}(\omega_{\nu}) c_m \end{equation} где \begin{equation*} V_{nm} = \int_V dV' \boldsymbol{F}_{n}(\boldsymbol{r}') \delta\hat{\mathcal{P}}_0(\boldsymbol{r}',\omega_{\nu}) \boldsymbol{F}_{m}(\boldsymbol{r}') \end{equation*}
nanophotonics](https://iopscience.iop.org/article/10.1088/1361-6641/ac3290), Semicond. Sci. Technol. 37 013002 (2022) 2. Weiss, T. and Muljarov, E. A. How to calculate the pole expansion of the optical scattering matrix from the resonant states, Phys. Rev. B 98 085433 (2018) 3. P. T. Kristensen, K. Herrmann, F. Intravaia, and K. Busch, Modeling electromagnetic resonators using quasinormal modes, Adv. Opt. Photon. 12, 612-708 (2020)