{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2badadd6",
   "metadata": {},
   "source": [
    "Лицензия MIT\n",
    "\n",
    "© Алексей Александрович Щербаков, 2024"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a7bef001",
   "metadata": {},
   "source": [
    "# Лекция 2.2. Фотонные кристаллы. Полный модальный базис в одномерном случае"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "746a64ed",
   "metadata": {},
   "source": [
    "В данной лекции рассмотрим аналитическое решение для одномерного фотонного кристалла, ячейка которого состоит из двух слоёв различных материалов с параметрами $\\varepsilon_{1,2}$ и $\\mu_{1,2}$. Данная модель является прямым аналогом модели Кронига-Пенни в электронной теории твердых тел, но в отличие от квантовомеханического случая, где периодическая система прямоугольных ям является лишь приближением потенциала, позволяющим получить аналитическое решение, в оптике кусочно-постоянное изменение материальных параметров среды описывает реально существующие структуры (если пренебрегать их внешними границами) - брэгговские зеркала."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55fbf405",
   "metadata": {},
   "source": [
    "## Уравнение Гельмгольца\n",
    "\n",
    "Пусть границы разделов материалов в одномерном фотонном кристалле с периодом $\\Lambda$ параллельны плоскостям $XY$ декартовой системы координат, так что периодичность задается направлением оси $Z$. Без ограничения общности можно считать поля не зависящими от координаты $y$, и уравнения Максвелла разделятся для ТЕ и ТМ поляризации соответственно:\n",
    "\\begin{equation}\\tag{1}\n",
    "    \\begin{array}{c}\n",
    "    -\\dfrac{\\partial E_{y}}{\\partial z} = i\\omega\\mu H_{x} \\\\\n",
    "    \\dfrac{\\partial E_{y}}{\\partial x} = i\\omega\\mu H_{z} \\\\\n",
    "    \\dfrac{\\partial H_{x}}{\\partial z}-\\dfrac{\\partial H_{z}}{\\partial x} = -i\\omega\\varepsilon E_{y}\n",
    "    \\end{array}\n",
    "\\end{equation}\n",
    "и\n",
    "\\begin{equation}\\tag{2}\n",
    "    \\begin{array}{c}\n",
    "    \\dfrac{\\partial H_{y}}{\\partial z} = i\\omega\\varepsilon E_{x} \\\\\n",
    "    \\dfrac{\\partial H_{y}}{\\partial x} = -i\\omega\\varepsilon E_{z} \\\\\n",
    "    \\dfrac{\\partial E_{x}}{\\partial z}-\\dfrac{\\partial E_{z}}{\\partial x} = i\\omega\\mu H_{y}\n",
    "    \\end{array}\n",
    "\\end{equation}\n",
    "Уравнения Гельмгольца, следующие из этих двух, систем можно записать в едином не зависящем от поляризации виде\n",
    "\\begin{equation}\\tag{3}\n",
    "    \\dfrac{\\partial}{\\partial z}\\left(\\dfrac{1}{\\eta\\left(z\\right)}\\dfrac{\\partial F}{\\partial z}\\right)+\\dfrac{1}{\\eta\\left(z\\right)}\\dfrac{\\partial^{2}F}{\\partial x^{2}}+\\omega^{2}\\zeta\\left(z\\right)F=0\n",
    "\\end{equation}\n",
    "где $F^{s}\\left(x,z\\right) = E_{y}\\left(x,z\\right)$, $\\eta^s=\\mu$ для ТЕ поляризации и $F^p\\left(x,y\\right) = H_{y}\\left(x,y\\right)$, $\\eta^p=\\varepsilon$ для ТМ поляризации. Также обозначим вторую ненулевую компоненту поля в плоскости слоёв символом $G$, так что $G^s = H_x$ и $G^p = E_x$. Эту компоненту можно найти из уравнений (1),(2): $G^{x,p} = (\\pm)_{s,p}(i/\\omega\\eta^{s,p}) dF^{s,p}/dz$. Уравнение Гельмгольца дополняется условиями непрерывности полей $F$ и $G$ на границах разделов сред между слоями кристалла и условием квазипериодичности."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4bbb27af",
   "metadata": {},
   "source": [
    "## Дисперсионное уравнение\n",
    "\n",
    "Для решения уравнения воспользуемся методом разделения переменных. Записывая $F=u\\left(x\\right)v\\left(z\\right)$, мы приходим к двум уравнениям на новые функции:\n",
    "\\begin{split}\\tag{4}\n",
    "    \\dfrac{1}{u\\left(x\\right)}\\dfrac{d^{2}u\\left(x\\right)}{dx^{2}} &= -\\beta^{2} \\\\\n",
    "    \\dfrac{\\eta\\left(z\\right)}{v\\left(z\\right)}\\dfrac{d}{dz}\\left(\\dfrac{1}{\\eta\\left(z\\right)}\\dfrac{\\partial v\\left(z\\right)}{\\partial z}\\right) & + \\left[ \\omega^{2}\\varepsilon\\left(z\\right)\\mu\\left(z\\right)-\\beta^{2} \\right] = 0\n",
    "\\end{split}\n",
    "Решение первого соответствует трансляционной инвариантности полей вдоль оси $X$ (направления вдоль плоскостей кристалла) и задается экспонентой\n",
    "\\begin{equation}\\tag{5}\n",
    "    u\\left(x\\right)\\sim\\exp\\left(i\\beta x\\right)\n",
    "\\end{equation}\n",
    "с константой распространения $\\beta$. Второй оператор в случае, когда материальные параметры являются действительными функциями координат, является самосопряженным. Следовательно, его собственные значения $\\beta^2$ вещественны, а собственные векторы, соответствующие этим значениям, ортогональны. Соотношение ортогональности записывается как\n",
    "\\begin{equation}\\tag{6}\n",
    "    \\dfrac{1}{\\Lambda}\\intop_{0}^{\\Lambda}\\dfrac{v_{m}\\left(z\\right) v^*_{n}\\left(z\\right)}{\\omega\\eta\\left(z\\right)}dz=\\delta_{mn}\n",
    "\\end{equation}\n",
    "\n",
    "В рассматриваемом нами случае двух однородных слоев с толщинами $d_{1,2}$ на периоде ($d_1+d_2=\\Lambda$) уравнение на функцию $v(z)$ в каждом слое сводится к уравнению второго порядка с постоянными коэффициентами\n",
    "\\begin{equation}\\tag{7}\n",
    "    \\dfrac{d^{2}v\\left(z\\right)}{dz^{2}}+\\kappa_{1,2}^{2}v\\left(z\\right)=0\n",
    "\\end{equation}\n",
    "где $\\kappa_{1,2}=\\sqrt{\\omega^{2}\\varepsilon_{1,2}\\mu_{1,2}-\\beta^{2}}$. Решения этого уравнения выразим в форме\n",
    "\\begin{equation}\\tag{8}\n",
    "    v\\left(z\\right)=\\begin{cases}\n",
    "\\mathcal{A}_{1}^{+}e^{i\\kappa_{1}z}+\\mathcal{A}_{1}^{-}e^{-i\\kappa_{1}z} & 0\\leq z\\leq d_{1}\\\\\n",
    "\\mathcal{A}_{2}^{+}e^{i\\kappa_{2}\\left(z-d_{1}\\right)}+\\mathcal{A}_{2}^{-}e^{-i\\kappa_{2}\\left(z-d_{1}\\right)} & d_{1}\\leq z\\leq\\Lambda\n",
    "\\end{cases}\n",
    "\\end{equation}\n",
    "Постоянные коэффициенты $\\mathcal{A}_{1,2}^{\\pm}$ связаны друг с другом условиями сопряжения и квазипериодичности:\n",
    "\\begin{split}\\tag{9}\n",
    "    v\\left(d_{1}-0\\right) = v\\left(d_{1}+0\\right) \\\\\n",
    "    e^{ik_{B}\\Lambda}v\\left(0\\right) = v\\left(\\Lambda\\right)\n",
    "\\end{split}\n",
    "Те же условия, записанные для второй компоненты поля $G$, дают аналогичные условия для функции $\\left(1/\\eta\\right)dv/dz$. В итоге получается система из четырёх уравнений\n",
    "\\begin{split}\\tag{10}\n",
    "    \\mathcal{A}_{1}^{+}e^{i\\kappa_{1}d_{1}}+\\mathcal{A}_{1}^{-}e^{-i\\kappa_{1}d_{1}} = \\mathcal{A}_{2}^{+}+\\mathcal{A}_{2}^{-} \\\\\n",
    "    e^{ik_{B}\\Lambda}\\left(\\mathcal{A}_{1}^{+} + \\mathcal{A}_{1}^{-}\\right) = \\mathcal{A}_{2}^{+}e^{i\\kappa_{2}d_{2}}+\\mathcal{A}_{2}^{-}e^{-i\\kappa_{2}d_{2}} \\\\\n",
    "    \\dfrac{\\kappa_{1}}{\\omega\\eta_{1}}\\left(\\mathcal{A}_{1}^{+}e^{i\\kappa_{1}d_{1}}-\\mathcal{A}_{1}^{-}e^{-i\\kappa_{1}d_{1}}\\right) = \\dfrac{\\kappa_{2}}{\\omega\\eta_{2}}\\left(\\mathcal{A}_{2}^{+}- \\mathcal{A}_{2}^{-}\\right) \\\\\n",
    "    e^{ik_{B}\\Lambda}\\dfrac{\\kappa_{1}}{\\omega\\eta_{1}}\\left(\\mathcal{A}_{1}^{+}-\\mathcal{A}_{1}^{-}\\right) = \\dfrac{\\kappa_{2}}{\\omega\\eta_{2}}\\left(\\mathcal{A}_{2}^{+}e^{i\\kappa_{2}d_{2}}-\\mathcal{A}_{2}^{-}e^{-i\\kappa_{2}d_{2}}\\right)\n",
    "\\end{split}\n",
    "Приравнивая детерминант к нулю, после несложных алгебраических преобразований получаем дисперсионное уравнение:\n",
    "\\begin{equation}\\tag{11}\n",
    "    \\cos\\left(\\kappa_{1}d_{1}\\right)\\cos\\left(\\kappa_{2}d_{2}\\right)-\\dfrac{1}{2}\\left(\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}+\\dfrac{\\eta_{1}}{\\kappa_{1}}\\dfrac{\\kappa_{2}}{\\eta_{2}}\\right)\\sin\\left(\\kappa_{1}d_{1}\\right)\\sin\\left(\\kappa_{2}d_{2}\\right) = \\cos\\left(k_{B}\\Lambda\\right)\n",
    "\\end{equation}\n",
    "которое для действительных $\\varepsilon_{1,2}$, $\\mu_{1,2}$ имеет бесконечный дискретный набор решений $\\beta_m^2$. Для заданных параметров кристалла имеется конечное число $\\beta_m^2 > 0$, соответствующих распространяющимся вдоль оси $x$ модам и бесконечно число $\\beta_m^2 < 0$, соответствующих затухающим решениям. Это можно проиллюстрировать, если построить график зависимости левой и правой частей уравнения от $\\beta^2$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "36d4023d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/xnp5ZAAAACXBIWXMAAAsTAAALEwEAmpwYAAA8ZklEQVR4nO3dd3ycZ5Xo8d+ZGc2oS5bVrOIuySUusRWnOomTkE4CJGETOmxuYJdclgV2FxYW9rLLvbDcyy67cMmadmGBhISQxCHFOCROIc2yLffe1Cyr96557h8zIwtHbWbeKa/mfD8ffyTNvJr38dg+fnSe85xHjDEopZSa/RyxHoBSSqno0ICvlFIJQgO+UkolCA34SimVIDTgK6VUgtCAr5RSCSLsgC8ipSLykogcFJEDIvJXE1wjIvLvInJcRPaKyLpw76uUUio4LgteYwT4vDFml4hkADtFZJsx5uC4a24Byvy/LgV+4P+olFIqSsKe4Rtjzhpjdvk/7wYOAcUXXHYn8HPj8yaQLSLzwr23UkqpmbNihj9GRBYCFwNvXfBUMVA77us6/2NnJ3iNB4AHANLS0tYvW7bMyiEqpWahI0eOAFBRURG1e3qN4UBDF4WZyeRleKJ23+ns3LmzxRiTN9FzlgV8EUkHHgc+a4zpCvV1jDGbgc0AlZWVpqqqyqIRKqVmq2uvvRaA7du3R+2ePYMjXPS1rfz9rct44OolUbvvdETkzGTPWVKlIyJJ+IL9L40xv53gknqgdNzXJf7HlFLKlgJ9yBwiMR7JzFlRpSPAj4FDxpjvTHLZFuAj/mqdy4BOY8w70jlKKWUXXhv2nbQipXMl8GFgn4hU+x/7e2A+gDHmIeBZ4FbgONAHfNyC+yqlVOz4A77YaIYfdsA3xrwGTPk7Nr6ffT4d7r2UUmqmhke9fPaRas51DfDv911MUXaKpa9vCKR0LH3ZiNKdtkqpWenXO2p5Zt9Zqs60872Xjlv++oGUjo3ivQZ8pdTs9GhVLcvnZXLP+hKe2FXP4Miopa8fWLS1U0pHA75SatZp7h5kb10nd6wp4vrlBfQPj7K/vtPSewTWbDWlo5RSMfTWqVYArlgyl8qFcwDYcbrd0nt4zdiqraWvG0mW7rRVSql4UHW6nVS3k5VFmbicDkpzUthn8Qw/MMXXGb5SSsXQ4cYuKgozcDl9Ia4sP4MTTT2W3uP8oq19Ir4GfKXUrGKM4UhjNxUFGWOPLc1P52RLL6MW7pYKlGXaKKOjAV8pNbu09AzR3jdM+fiAn5fO0IiX2rY+y+5jNKWjlFKxdfRcNwDLCs8H/IW5aQDUWBjwA4u2mtJRSqkYOdLoC/jl4wJ+8RzfLtuGjn7L7hOY4dso3mvAV0rNLkfPdTM3zU1u+vke9QUZHhwC9RYG/ICE6paplFLx5FhTD2UF6X/ymMvpoDAz2dKAfz6lYx8a8JVSs0pNWx8LctLe8XjxnBTq261P6dhogq8BXyk1ewwMj9LcPUhpzjs7YxZmpdDYNWDZvc63VrBPxNeAr5SaNerafVU4pTmp73guL91DS/egZfcaS+nYJ95rwFdKzR61bb6UTcmcd87w8zI89A6N0js4Ysm9jA0PQNGAr5SaNcZm+HMmmOFn+Kp2WnqsmeUbXbRVSqnYqW3vx+NyjAX38fL9jzVblNYZK8O3UcS3JOCLyE9EpElE9k/y/LUi0iki1f5fX7XivkopNV5tWx8lc1ImTLPkWR3wx1or2CfiW9Ue+f8B3wN+PsU1rxpjbrfofkop9Q617X2UTJDOgXEB36KUTsLW4RtjXgHarHgtpZQKVV17/4QlmQDZKUkAtPUOWXIvrcOf2uUiskdEnhORlVG8r1IqAYx6DR19wxRnTzzDdzkdZKUk0W5VwMd+Z9pG68SrXcACY0yPiNwKPAmUTXShiDwAPAAwf/78KA1PKWV3Q6NeAOZlJU96TU6am7a+YUvuNzbDt+TVoiMqM3xjTJcxpsf/+bNAkojkTnLtZmNMpTGmMi8vLxrDU0rNAkMjvoBfkDl5wJ+TauEMX+vwJyYiheJ/V0Rkg/++rdG4t1IqMQQC/lQz/Dmpbtr7rE3p2OkAFEtSOiLyMHAtkCsidcDXgCQAY8xDwN3AX4jICNAP3GuMse6sMaVUwgukdAqnCvhpbg6d7bLkfl4bLtpaEvCNMfdN8/z38JVtKqVURAyNeJmTmkRyknPSa3w5fKtSOvZbtNWdtkqpWWFoxEvhFPl7gOzUJAaGvfQPjYZ9P68u2iqlVGwMjXqnTOcAZCb7avG7B6yo1NEZvlJKxcRMZviZ/s1XXRYE/POtFcJ+qajRgK+Usj1jYHgGM/yMZN+yZddA+C2Sz6d07BPxNeArpWxvrEJnuhm+P6XT1W/FDF8PQFFKqagb23Q1bQ7fN8PvtmCGn7DtkZVSKpaG/TP8/An64I+XMbZoa0VKJ9At0z4RXwO+Usr2AgE/L33qgJ+ZEsjhW1ClY8ONVxrwlVK2Nzzqi745ae4pr0tJcuJ0iCVlmV4bHoCiAV8pZXvDo16SnA5czqlDmoiQmeyiq9+KHL4u2iqlVNQN+QP+TGQkJ1kyw9c6fKWUioHhES9JzplF3swUl0V1+GN1OmG/VrRowFdK2d7wqJn5DN9j0Qzf/1FTOkopFSXGGIZHvbhdMwtnmSkuS8oyseGibbSOOFRKxbnhUS/76ztxuxwsL8zEYZPkdPfgCF4TxAw/OcmSnbbn6/DtQwO+UopdNe381SO7qW3rB2BNSRbf+8A6SnMmPhA8nrR0DwLMOIefkWzNDN9oHb5Sym6ONHbzwR++hSD8+30X8z/fu4pTLb3cu/lNWnsGYz28aTWPBfwZpnSSk+geHGHUG96he4Hv1pSOUhFy7Fw3T+89y9HGbtwuB6tLsrh9ddG0XRLVxHoHR/iLX+wkzePisU9dPnYA+EXFmdz9gzf4+u8O8t17L47xKKfW0uM7wWrmKR1f2OsZHCHL3y45FF4bntI6KwP+03saMEBOqptl8zLInWa7tYp/HX1DfG3LAZ6qbsAhsHBuGkOjXrbsaeDbW4/wqWuW8Jnry3DaJO8cLx56+QQnW3p5+L9dNhbsAVaXZPOpaxbz7y8e5/2VpVy5NDeGo5xac/cAMPOUzlhP/P7hsAK+SdRFWxH5CXA70GSMuWiC5wX4LnAr0Ad8zBizy4p7T+RvfrOHgWHv2NdrS7P54KXzec/FxTOeBaj4caqllw/96C3OdQ3w4KalfOKqRWNb6E+19PKdbUf57h+OUV3bwQ8+tI5U96ycx1iusXOAH756knevKeLyJXPf8fxfblrK47vq+c62o/Ed8HsGEZEgUjrWdMxM5PbI/w+4eYrnbwHK/L8eAH5g0X0ntPWzV/PC567mV/dfyhduLKdvaIS/+c1ebv63V3jzZGskb60sdqa1l3seeoP+4VF+8xdX8IWbKv6kX8qi3DT+w593fvVYM5/8r51jrXLV1Da/cpKRUcPf3Fgx4fPJSU7u37iInWfa2VXTHuXRzVxL99CMZ/dwvmNmuA3UErYO3xjzCtA2xSV3Aj83Pm8C2SIyz4p7T2TB3DSW5mdwxdJcHryujK2fvZoffaSSEa/hAz98k++/dBxvmAs2KvI6+4b5+E93MOL18ugnL2Ntafak137g0vl8832refVYC//j6QPRG6RNdQ8M82hVLbevnsf8uZNX4txTWUpGsosfv3YqiqMLTnPPYFA/uad7fDP8vqFwZ/i+j3ZK6UQrv1EM1I77us7/2DuIyAMiUiUiVc3NzZbcXES4YUUBz35mI7evLuLbW4/w2V9XMzKqM8F4ZYzh75/YR01bHz/8SCVL8zOm/Z73X1LKJ69ezC/fquF3exuiMEr7+s3OOnoGR/j4lYumvC7d4+Lu9SVsO3COTgtq1yOhuTu4gJ/mcQLQOzga1n3tWIcfdwltY8xmY0ylMaYyLy/P0tdO87j47r1r+ZubKtiyp4H//vBu/fE/Tj2xu55n9p3lr99VziULc2b8fV+4qYK1pdn8w5P7aesdiuAI7csYw8Nv17CmNJs1U/zUFPCetcUMjXrZur8x8oMLQUvPYFApncAaT9gzfP9HG03woxbw64HScV+X+B+LOhHh05uW8pXblvPc/ka+8uS+scUXFR/OdQ3wtacOcMnCOXzqmiVBfW+S08G37lpN98AI33jmUIRGaG8HGro4eq6Hu9eXzOj61SVZLJybylN7YvJPdkper/EH/JmHslS3NTP884u29on40Qr4W4CPiM9lQKcx5myU7j2h+zcu5jPXLeXRqjoeevlkLIeiLvAvzx9hcMTLt+9eE1KZZUVhBvdvXMzju+o40NAZgRHa2xO760lyCrevmtkymohwx9piXj/RGncbsTr7hxkeNbiDCvi+GX7/cLgB3/fRPuHeooAvIg8DbwAVIlInIn8uIp8SkU/5L3kWOAkcB34I/KUV9w3XX7+rnHevKeJbzx/m5aPWrBeo8FTXdvD4rjo+cdUiFuamhfw6f3HtErJSkvj21iMWjs7+vF7D03sa2FSRz5xpToca713LCzAGth+Jr38nLT3B7bIFcLscJDmF3sFwUzoJOsM3xtxnjJlnjEkyxpQYY35sjHnIGPOQ/3ljjPm0MWaJMWaVMabKivuGS0T49t2rKS9I5/OPVo9t0VaxYYzh608fIC/Dw4PXLQ3rtbJSkvjLa5ew/UizluKOU13XQVP3ILfOcHYfsLIok7wMDy8eaYrQyELTHGQfnYBUt4u+IWtm+Hba6xd3i7bRlpzk5HsfWEf3wAiff2yP5vNjaPvRZnbVdPDXN5SPlc6F46NXLKQg08O/bjtqwehmhxcOnsPpEDZV5Af1fQ6HsKkij1eONo8dGB4PmgMz/Bm2Rg5IczvDnuF7x1I69on4CR/wAcoLMvjKbct55Wgzj+2si/VwEpIxhu++cIzi7JQZLyZOJznJyf1XLeatU23sq9NcPsC2g+e4dFEOWanBtxS4blk+3QMj7K7psH5gIQq2cVpAqseKGX7i7rS1vQ9euoANC3P4xjOHaPL35lDR88qxFqprO/jLTUtmfJDFTPzZhlLSPS5++KouzJ9u6eVYUw/vWlEQ0vdfvjgXEXj9RIvFIwtdc88gbqcDV5B5lTS3k16LNl7ZKeBLPKcwKisrTVVVCOn+574IjfuC/rb+4VH21nWQk+ambAYbfZR1DjR0MjjiZW1ptuU7F8+09nK2a4CLS7PxuJyWvradNHT2U9PWx9rSbJJDfB/21nfgcggr5mVZPLrQHG/uoat/mM897Fuc3/53l83o+w6e7cIYw8qi0H8fTd0DnGzpjczfq8JVcMs3Q/pWEdlpjKmc6Dmd4Y+TkuSkKDuF1t6hsPtsqJnrGRyhe3CEeVnJEdmmHmid3NiV2D+5dfQNkZrkDDnYA2T5e8nHS2vg4VFvSA0RHQKjVv0ebDTDn51tBUP8nxFg7tAI9/zvl8kf9fDkR6+0zTFvdvbVR6t5vqGRNz55PYTRrnYyHuDHP69id007r3/4ektTRnbRNzTCnf/j93ziykWsvnV5yK9z/NA5/vxnVfzqmku5YknsO2h+7ruvUpydDIV/43vg48/M6Pu+//Bu9td38tLHrw353i+8VcPfP7GPN++53jbnMSTe3/xppLpd/O3NFeyt6+TJ6vjbWTjbNHcP8rs9Z7l7fUlYvcmnc9+GUlp6hnjx8LmI3SOevX2qjeFRw1Vl4QXpSxblIAJvnZyqV2L0tPQMhnTehRVVOufr8MN6majSgD+B96wtZnVJFv/y/JGw+22oqf3qrRqGRr189IqFEb3PNeX5FGYm8/DbtdNfPAu9dqwFt8sRVF+iiWQmJ1FRkEF1bYc1AwvDqNfQ2jNIXkbwAd/KOnwN+DbncAhfuW0FjV0D/NcbZ2I9nFlraMTLL946w7UVeSzJS4/ovZwO4f2VJbxyrJm69r6I3isevXa8hcoFc0hOCn9x8eL52VTXdsS8xXhb7xBeQ0gBP83jq9IJp2hlrCzTRkl8DfiT2LAoh41lufznKyfD/tFPTez3Bxtp7h7kYxGe3QfcU1mKMfBUdWK1Tm7uHuRwY3fY6ZyAtaXZdPYPc6q115LXC1WgrUIoKZ1UtwtjYDCMbrnaLXOW+et3ldPWO8TPdZYfEY9W1VGcncLGMmvbYE+mNCeVSxbO4cnd9Qm1ozpQN3+VRccUXjx/DgDVMd6AFdh0FeoMHwhrMqcHoMwy6+bP4ZryPDa/ckJn+Rar7+jn1WPN3LW+JKoHj9+xtphjTT0cOtsdtXvG2hsnWslMdoVVcz7e0rx0MjwudtfG9tjDwAw/L8QZPhBWHl8PQJmF/vpd5bT3DfOzN07HeiizyuM76zAG7rGojcJM3bZqHi6HxGVv90h5+1QblyzMsew/VodDWFOaHfMWC4EZfm4oM/xAT/wwijJ0hj8LrS3N5pryPH7y2ikGwuyfrXy8XsOjVbVcuXQupTmTn6caCTlpbq4uz+Pp6oaYLzpGQ3P3ICdberlkUXjVORe6eH42hxu7Y1rF1tw9SEqScyx4ByPFgkNQvDZsiK8BfwY+ec1iWnqG+O2uxJkVRtKbJ1upa+/n/ZWl018cAXeuLaKhc4CdNbFNSURD1WlfvXy45ZgXunh+NqNeE9OmdC3+ksxQ+tGnWXSQOeii7axz+eK5rCrO4kevnkyIWWGkPbazjsxkFzetLIzJ/a9blo/b6YjbM1qt9PbpNpKTHKwqtrb3zaribAD2N3RZ+rrBaO4ZJDd95oe4jGfFMYea0pmlRIRPXrOYky29bDuUmDs1rdI/NMrWA43ctrrIkprwUGQkJ3Hl0rlsPdg466t1dpxuY21ptuXtJPIyPORneDhQH7sZfnN3aJuuANIsOMhcF21nsZtXFlKak8LmV7TNbjheOHSOvqFR7lxbFNNx3HxRIbVt/Rw8G7sZaqR1DwxzsKGLDRancwJWFmVyIIYz/JaeoZBq8AFSA2WZYVTpJGwdvojcLCJHROS4iHxxguc/JiLNIlLt/3W/FfeNJpfTwf1XLWbnmXZ2J0DuN1Keqm6gMDM5YkFopm5YXoBDmNVpnV01HXgNli/YBqwsyuJ4c09MihmGR7209Q6FPcPXOvwgiYgT+D5wC7ACuE9EVkxw6a+NMWv9v34U7n1j4a71JaR7XLoRK0QdfUO8fLSJd6+ZF/MupHPTPVyyMIetB2Zvim7HqTacDmGdf6OU1VYWZTLqNRxpjP6ehrbeISC0TVfga4UOvhRjqOKlRXQwrJjhbwCOG2NOGmOGgEeAOy143biT7nFx9/oSntl7Vg88D8Fz+xsZHjXcubY41kMBfGmdI+e6OdncE+uhRMTbp9tYWZQ5VpFitcBGrlikdcZq8ENM6TgcgtvlYGAk/J9ObDTBtyTgFwPjWxDW+R+70F0isldEfiMik9bjicgDIlIlIlXNzc0WDM9aH7psAUOjXn69oybWQ7GdLdUNLM5LY2VRZqyHAjBWJTQbZ/lDI1721HZYXo45XmlOChnJLg40RH/hNpy2CgEpSU4Gwpnh+yv2EiqlM0NPAwuNMauBbcDPJrvQGLPZGFNpjKnMy4tOj5VgLM1PZ2NZLr94s4aR0dAbLyWaxs4B3jzVyh1rikKqm46EouwULirO5A+zsPLq0NkuBke8EUvngK96bcW82CzcNofRViEgJclJfxjrD2OLtiG/QvRZEfDrgfEz9hL/Y2OMMa3GmEAO5EfAegvuGzMfuXwhjV0DbDs4+wJFpPxubwPGwB1rYludc6HrKvLZVdNOR99QrIdiqUBhwcXzsyN6n5VFWRxu7Ir65MeKGX5ykoOB4TC6ZSbioi2wAygTkUUi4gbuBbaMv0BE5o378g7gkAX3jZnrluVTnJ2i/XWC8PTes6wqzmJxhPveB+vaZfl4Dbx8NP7Sh+HYXdtBQaaHeRE+em9lUSYDw15OR7lVcnP3IBkeV1h7OZLDnOGP1eHbJ96HH/CNMSPAg8BWfIH8UWPMARH5uojc4b/sMyJyQET2AJ8BPhbufWPJ6RA+fPkC3jzZxtFzidN1MVQNHf3sqe3gllWx2Vk7lTUl2eSkudl+ZJYF/JoOLi6dE/H0WUVhBgBHGqO78N3cPUheZuize/D10wmnpPR8Hb59Ir4lOXxjzLPGmHJjzBJjzDf8j33VGLPF//mXjDErjTFrjDGbjDGHrbhvLN2zvoQkp/BIgh6ZF4ytB3y17jfHqJXCVJwO4ZryPLYfaWJ0lrTNaO0ZpKatj7URTueAb03LIXCkMbp5/ObuwbDy9wDJrvACPsbYanYPutM2ZHPTPdy4opDf7q5j0ILSrtns+f2NVBRkxF06J2DTsnza+4bZU9cR66FYInDe7MWl2RG/V3KSk4W5aRyJ8k+6zT2D5GeGl65KcYeb0rHXgi1owA/LvRtK6egbnpVlfVZp6Rlkx+k2broo/mb3AVeX5eIQ2H64KdZDscTumg6cDmFVibUN0yZTUZDB0XPRTek0dQ2EPcNPSXKGtfHKYGyVzgEN+GG5ckkuJXNStCZ/CtsOnsNr4jOdE5Cd6mbd/Dm8eGR2BPzq2g6WFWaMneoUaeUFGZxu7Q0reAajd3CE3qHRsCp0ADwWVOnEeMN40DTgh8HhEP6sspQ/Hm/lTIwPdI5Xz+9vZMHcVJbPy4j1UKa0aVk+++u7aOoaiPVQwjLqNVTXdrA2CumcgGWFGRgDx5uiM8sPlGTmhxnwU5LCy+H7Ujr2ivga8MN0T2UpDoFHq3Tx9kKd/cO8fqKFm1cWxv2Pvpsq8gFsX61zormHnsGRsYPGo6HcX6lzOEoLt2ObrmIc8H0pnbCGEHUa8MNUmJXMpop8Hquq0523F3jx8DmGR01c5+8Dls/LoDAzmZdsntap9p8zG+kNV+MtyEnF7XJErUTZik1XcL4OP9QzEYyxVw0+aMC3xL0b5tPUPchLNp8dWu35/Y0UZiaztiQ71kOZloiwsSyXPx5vsXV55u7adrJSklg0Ny1q93Q5HZTlp3MkSgu3gbRb2CkdtxOvgaEQJ2rGGE3pJKJNFXnkZ3h45G1dvA3oGxrh5aPN3LSyIOatkGdqY3keXQMj7LVxeebumg7WlGZH/T2vKMiIWi1+c88gTocwJzW04w0DArt0Q1241Rl+gnI5Hdy1voTtR5tp6rb3op9VXj7SzMCw1xbpnICrluYiAq8da4n1UELSMzjCkXPdUam/v1BFYQbnugaj0pOoudt3lm24/6klJ/nCX6h5fK+xVx8d0IBvmbvWlTDqNTy1uyHWQ4kLzx9oJCfNHfOTrYKRk+ZmZVEmr9o04O+t68CY6ObvAwILt9Gox2/qHiQ/I/weQeEegqKLtglsaX46a0uzeXxX3aw/GHs6gyOjvHioiXctL8DltNdfsY1leeyqaacnjKPvYmW3f8E2miWZARUFgZ46kU/rhHN4+XiBgB/qIShGZ/iJ7e71JRxu7I7pwc7x4PXjrXQPjnCzjdI5ARuX5jLiNbx5ojXWQwna7poOFuemkR1mbjsU87KSyUh2RaXFgm+GH37ATw5zhu/VXjqJ7d2ri3C7HPxmZ12shxJTz+9vJMPj4oqlc2M9lKCtXziH5CQHrx23V1rHGP+Gqxikc8BX5VSWnx7xlM6o19DaY80Mfyzgh5jD1xl+gstKTeJdKwrYsqeBoZHErMkfGfXy+4ONXLc8H48r9F7lseJxObls8VxeOWavEtu69n5aegajuuHqQhWFGRw71x3RlGZb7xBeE34NPvjKMgEGQ6zS8RqjrRUS3d3rSmjrHbL9Bp5QvX26jfa+YW6xYTon4KqluZxs7qW+oz/WQ5mx3VHskDmZsvwM2vuGx3bCRkKgCs6alI4v/IU6w/dt17BXxNeAb7GNZbnkZXh4PEHTOlv3N5Kc5ODq8vg7j3imAmN/zUaz/OqaDpKTHCwrjF3PonL/wu2xCKZ1rNplC+FX6YDO8BOey+ngvRcX8+LhJlojONOJR16v4fkDjVxTnhe1To2RUJafTkGmh1dsVJ65u7ad1cXZMa2KKi/wnXcQyRYLYwE/3bqyzFCrdLxezeErfDX5I17Dlj2JVZNfXdfBua5Bbrlo3vQXxzER4aqlebZpszA4MsqB+q6YLdgG5GV4yEpJiujCbZOFM3yPVumERkRuFpEjInJcRL44wfMeEfm1//m3RGShFfeNVxWFGawqzkq4ap2t+xtJcgqbluXHeihhu7o8l46+YQ40dMZ6KNM62NDF0KiXdTEO+CJCeUE6xyI8w8/wuMYWXMMxNsMPtUqHBJzhi4gT+D5wC7ACuE9EVlxw2Z8D7caYpcC/At8K977x7q51xRxo6Ipay9hYM8bw3P5GrliSS1ZKUqyHE7Yrl+YC2GLX7e6xDpmxq9AJKCvI4GgEK3Ws2nQFkOQUnA4JuZeOHWf4Eu4fjIhcDvyjMeYm/9dfAjDG/K9x12z1X/OGiLiARiDPTHPzjIwMs379+rDGFysjo4adNe0UZiazYG5qrIcTcX1Do+yt62BxXrolFRTxYF9dJ06nsGJeZqyHMqVjTT30DAzHRcBv7BzgdGsv6xbMwR2B9YSD/k2NK4r+9M+kuroagLVr1wb1ejtOt5GfEdq/0eNNPXQPjsS0MmoiL7/88k5jTOVEz1nxJ1IMjD/9o87/2ITXGGNGgE5gwl05IvKAiFSJSNXw8LAFw4sNl1OYk5pES88gidBpoa13CAHmpNp/dh+QlZpE98BI3OfxewZGSE+Oj/c91R1u5cvUhka9JLms+4/EIYI3jH+gNpvgE3elFMaYzcBmgMrKSrN9+/bYDigM2w6e47/9vIqvfrSS65cXxHo4EXXjv77M+lQ3v/7k5bEeimVeO9bCh378Fl//2CVxuy7R1DXAhv/5B75y23Lu37g41sOhpWeQyn9+gQdvX8GfX7XI0tc2xnDR17byZ5fM56vv/tOs8bXXXgtAsPHiqm+9yIZFOXzn/WuDHs9nHt7NvvpOXvrCtUF/byRNdbqcFf9V1gOl474u8T824TX+lE4WYL9mJUG6tiKPuWluHt81uxdvjzf1cPRcj603W02kcuEcPC5HXO+6HdtwFQfpHIDcdA85ae6ILNx2+w8vn5cVfklmQHIYxxx6jbHdDN+KgL8DKBORRSLiBu4FtlxwzRbgo/7P7wZenC5/PxskOR285+JiXjjYRHtv5PuEx8rz+88CcLPNyzEvlJzk5NLFc+N64XZ3TQdJTmFlUfysM/hOv7I+4Dd2+nbZFlga8B1htEdOwANQ/Dn5B4GtwCHgUWPMARH5uojc4b/sx8BcETkOfA54R+nmbHXXuhKGRr08vXf21uQ/t7+RdfOzKbTwH2K8uLosl+NNPZztjM82C7tq2llRlDXWCCwelBdkcPxcj+WVOmf9Ad/SGb7LyWCIfa+MMYlXlglgjHnWGFNujFlijPmG/7GvGmO2+D8fMMbcY4xZaozZYIw5acV97WBFUSYr5mXO2pr8mtY+DjR02X6z1WSuKovf8syRUS976zrirkqkvCCd7sGRsQBtlXP+1yvMtC7ge5IcIQd8rzcBZ/hqenevL2FvXWdEt5zHyvMHAumc2ZW/D6goyCAvwxOXAf9wYzcDw96YnHA1lbKCwOlX1v59D/wHUmBlwHc5GQz1ABQSdIavpnbn2iJcDpmVDdWe29/IRcWZlObMzr0GIsLGslz+eLwFb5yVZwYWbNfFyYJtQKSaqDV29ZOb7sZtYVlmcpIjjPbIU1fExCMN+FEwN93Ddcvy+e3uekZGZ0+f/LOd/eyu6Zi16ZyAjWW5tPUOcfBsfO2a3l3TTm66h5I5KbEeyp/ISXOTm+6OyAzf6nUiT5g5fHuFew34UXPX+hKauwfjMjUQquf3NwKzN50TEGizEG/lmdU1HVw8PzsuZ5nlBRkcbbJ4ht85QGGmtf+5eVyO0HvpGHDYLILabLj2takin5w096xavH1ufyPlBeksyUuP9VAiKj8jmeXzMnn1aPz8Z93eO8TJlt64y98H+Cp1rO2pc7ZzwNIKHfAF/JAXbRO1SkdNz+1ycOfaIrYdPEdHn/1r8pu7B9lxum3W1d5P5uqyXKrOtNE3NBLroQBQPXbCVXzl7wPKCtLpHRq17NSw/qFROvuHLU/pJCeFvmjrNfZrraABP4ruXu+vyZ8FffJ/f7ARY+DWVbM7nRNwVVkuw6OGt061xXoogK/+3iGwuiQr1kOZULnFlTqNXdaXZML5GX4oP4n4Nl7ZK+RrwI+ilUVZLJ+XyW92Xdh5wn6e39/Iotw0Kgpid6ReNF2yMAePyxE3aZ0dp9tYUZRJmifu2mEBUJ4fCPjW5PEDG98sT+kkOTHG15QtWEYPMVfTuXt9CXtqOyJ6SESktfYM8vqJVm65qNB2M5xQJSc52bAoh1fjYOF2aMRLdW0HlyzMifVQJpWVmkR+hse6GX5g01UEcvhASHl8Xz98e/3914AfZYGa/N/YuKHas/sbGfUa3r2mKNZDiaqry/I4FgdtFvY3dDIw7I3rgA++tI5VtfhnIxXw/S0pQqnFNwad4aup5aZ7uLYinyd22bcm/+k9DZTlp7OsMDHSOQGBNguvxbi0dod/HSHeA35ZQTrHm3os2bB2rmuAzGQXqW5rU1jnZ/jBL9zqDF/NyN3rS2iyaU3+2c5+dpxu491rimz3lz1cywozyE2PfZuFHafbWZSbZtlRf5FSXpBB//Aode3h/0TkK8m0foNZeCkdrdJRM3Ddsnzmprl5ZEdNrIcStGf2nsUYEi6dA76KjKvLcnkthm0WvF5D1Zk2LlkYn+WY41lZqdPQ0R+RbqweVxgHmZsEPMRcBc/tcnB3ZQkvHGriXJe1HQUjbcueBlYVZ7EoNy3WQ4mJjeWxbbNwvLmHjr7huE/ngC+lA3C0KfyAX9feT2mO9TP85KTwFm11p62akfsumc+o1/BYVe30F8eJUy297K3r5I4EnN0HBNosvHw0NtU6b9skfw+QmZzEvKzksBduuwaG6ewfpnSO9Q36AjP8UBZtfSde6QxfzcDC3DSuXDqXh9+ujbsujJP5nX/D2O1rEmN37UTyM5JZVZzFi4ebYnL/qtNt5GV4WDDXHt1Jywoywk7p1LX51gBKIhHwk0JftE3IE69U6O7bMJ/6jv64a8o1EWMMW/Y0sGFhTkQWz+xk07J8dte00xaDYyt3nG5nw8Ic2yyYl+f7KnVGw5jU1LX3AUSkK2i4i7aaw1czduOKQuamuXn47fhfvN1f38Wxph7uWJu46ZyA65fl4zXw8tHozvJrWvuo7+hnw6L4T+cElBdkMDjipaatL+TXCFT5RCLgB46GDGXR1hiTWDN8EckRkW0icsz/ccLSAREZFZFq/68LDzhPWHZavH18Vx1ulyMhq3MutKo4i9x0Dy8eju5PZn884SsHvXLp3KjeNxxjC7dhpHVq2/tIdTvJSXNbNawx4czwTQLO8L8I/MEYUwb8gckPJ+83xqz1/7pjkmsSkh0WbwdHRnmyup4bVxSQlZIU6+HEnMMhbKrI4+UjTVHdPPfH4y3kZ3hs1Y66bOz0q9ADfl17PyVzUiKSxhpbtA21Ssde8T7sgH8n8DP/5z8D3hPm6yWc8Yu34eQ5I+mlw0109A1z9/qSWA8lbly3LJ+ugRF2nmmPyv2MMbxxopUrl+baJn8PkO5xUZydElYTNV/Aj8wi9diibQgpHd8/V/v8WUD4Ab/AGHPW/3kjUDDJdckiUiUib4rIe6Z6QRF5wH9tVXNz/C9mWuGDly6gvqM/ZpUf0/nNzjoKMj1sLMuL9VDixlVluSQ5JWp/ZkfOddPaO8TlS+yTzgkoL0gPOaVjjKGurS9ixzgmhzHDn5XdMkXkBRHZP8GvO8dfZ3wNpSeboi4wxlQCHwD+TUSWTHY/Y8xmY0ylMaYyLy8xAsyNKwooykrmJ6+divVQ3qG5e5CXjjTznouLcdrtb3cEZSQnsWFRDi8cOheV+/3xeCtwfh+AnZQXZHCyuTek9Fdb7xDdgyMsmBuZjX5JTkEktBn+rMzhG2NuMMZcNMGvp4BzIjIPwP9xwumOMabe//EksB242LLfwSzgcjr46BULeeNkKwcb4uug7N/uqmPUa7h7naZzLnTjikJONPdy3IKdpNN5/XgLC+emUpxtv5LYsoIMhka9nAmhUud0ay8Ai3Ijk9IRkZCPOfQmWpUOsAX4qP/zjwJPXXiBiMwREY//81zgSuBgmPedde69ZD4pSU5++sf4meV7vYZfvV3DJQvnjC2+qfMCh7c/t68xovcZHBnljZOttpzdA2NdVQ+F0I7iZHMg4EduodrjcoaW0mEWzvCn8U3gXSJyDLjB/zUiUikiP/JfsxyoEpE9wEvAN40xGvAvkJWaxF3ri3mquoGWnsFYDweAV4+3cKa1jw9dtiDWQ4lLBZnJrF8wh2f3Rzbgv32qjb6hUa5blh/R+0RKWUE6bqeDffWdQX/vqZZeXA6JWA4ffP10Qm+PHIEBRVBYAd8Y02qMud4YU+ZP/bT5H68yxtzv//x1Y8wqY8wa/8cfWzHw2ehjVyxiaNTLL9+Mj41Yv3jzDHPT3GMzWfVOt1xUyKGzXZzxpx4i4cXDTXhcDq5YYs8ZvsflpKIwg311wQf80629zM9JJckZuT2iHpeTgRAPQLFTxRToTtu4sjQ/nWsr8vivN8+ENOOwUn1HP384dI73X1I6Vqus3ummlf60TgRn+duPNHP5krmkuO3757CqJIt99Z1BHxZ+srk34p1ZfTn80Gb4dqtj0IAfZ+6/ajEtPYM8vjO2B50/8nYNBvjAhvkxHUe8K81JZXVJVsQC/qmWXk619LKpwp7pnIDVxVl0D4xwpnXmC7der+F0ay8LIx3wkxxhHHFor4ivAT/OXLl0LmtKs/m/248zHKMjEAeGR3n47Ro2VeRTmmOProyxdPNFheyp7aC+w/qzbl/y1/nbPeCvKskCYG8QefzGrgEGhr1RmOGHtmjra49sLxrw44yI8N83LaWuvZ8t1Q0xGcMTu+tp6Rni/qsWxeT+dnPrRb520U/vsf7Pa9vBcyzNT2e+TdohT6a8IAO3y8H+IAL+8Sbf7tzFeZEN+KEu2moOX1ni+uX5LJ+Xyfe3H496uwWv1/DDV09yUXGmLXd1xsLC3DTWzc/m8Z11Qeeop9LUPcBbp1q5dZX9zx9IcjpYPi+TvXUdM/6ew42+Ms7lhZkRGpVP6Iu2msNXFhARHty0lJPNvTy77+z032ChFw6d42RzLw9cvcR2s5dYumt9Cceaejhg4ca5rfsb8Rq4fbX9Az748vj767tmfODPobPdFGYmMycCXTLHC33RVnP4yiI3X1RIWX4639l2NGq5fGMM399+gpI5KdyqpZhBuX1VEW6ng9/srLPsNX+39yxl+eljh4Hb3eqSLHoGRzjRPLNGaofOdrFsXuR/77rTVsWc0yH83c3LONXSyyNROiDlpSNN7Knt4MFNS3FFsO55NspKTeKGFfls2dNgyX/QTV0DvH26jdtmyeweoNJ/Du+O09N3GB0a8XKiuYfl8yKbzgHfISghVemgOXxloeuX53Ppohz+7YVj9AyORPRexhi+s+0o83NSuUvbIIfknvWltPUOsfVA+CWaz+47izFw2yzI3wcsnJtKbrqHHafbpr32ZEsPw6NmrC1DJLldDgZCWrTVHL6ykIjwpVuX09o7xOaXT0T0XlsPnGN/fRefub4sorsaZ7Ory/OYn5PKz14/HdbrGGP4dVUdK4syZ1UPIxFhw6I5vH1q+oAf6LsTjRm+x+VgKMQzbW02wdeAH+/WlmZz++p5bH71JLVhnAs6lcGRUf7Xc4dYmp/Oe/TM2pA5HcJHLl/AjtPtHGgIvo1AwP76Lg6d7eLeWbjp7ZKFOdR39NMwzZ6F6poOUt1OFke4Bh98M/xQAr5vhm+viK8B3wa+dOtyHCJ89an9lpb9Bfz0j6c509rHP9y+QnP3YbpnfSkpSc6wZvkP76ghOcnBHbPw/ODAAexvnmyd8rodp9tZN39OVP4+up1ORrxmxtVDAVqloyKiODuFz99YwUtHmnnW4la8TV0D/McfjnHD8nyuKU+MA2ciKSs1ifeuK+bJ6gaaQjiYvr13iN/uquOONUWz8vzg5YWZ5KZ72H5k8tPsugaGOdzYReXCOVEZk9t/kPlQkIvt3ghMviJNA75NfPTyBawqzuKrT+2nqTv4QDIRYwx//8Q+hr2Gr9y2wpLXVPDJqxcz6jU89PLJoL/3F2+eYWDYy/0bF0dgZLHncAjXVuTx8tHmSTcV7q7pwGt86Z9oCAT8oEszdYavIsXldPCd96+hd2iEzz+6J+gfPyfyWFUdLxxq4u9uXhbxBlWJZMHcNN57cTG/fOsM54KY5fcOjvCzN05zTXnerKm9n8iminw6+4eprp24PLPqdBtOh7C2NDsq4xmb4QcZ8LVbpoqosoIMvnr7Sl491sIPwqzaOdncw9d/d5DLFufw8SsWWjNANeYz15VhDHzrucMz/p4fvXqKlp4hPntDWQRHFntXleXidAjbDk58APyrx1q4qDiLNI8rKuPxOENN6WiVjoqw+zaU8u41RXx76xF+tze0Zl2dfcPc/7Mq3C4H//ueNTjsNk2xgflzU7l/4yJ+u7t+RnXn57oG2PzKCW5eWcjF86OTu46VrJQkrinP44ndde842Ly5e5A9dR1cH8XTvcZSOkEeZG7QKh0VYSLCt+9eTeWCOXzu13t4/XhLUN8/MDzKX/5qJ7Xtffznh9dTMsfeXRjj2ac3LaVkTgqfe7Sa7oHhSa8zxvCl3+5j1Bj+7pZlURxh7Ly/soRzXYO8euxP//4GNpzduLIgamMJfdE2wXbaisg9InJARLwiUjnFdTeLyBEROS4iXwznnsq3FfyHH6lkwdxUPvbTHTyzd2YN1noGR7j/Z1W8fqKVb75vddQWxRJVmsfFv/3ZWurb+/nCY3veMZsN+PFrp3jxcBN/e9OyiPd+jxfXLSsgJ83Nr8a1DTHG8GhVLSuLMlkW4Q6Z47mdoeXwTQL20tkPvA94ZbILRMQJfB+4BVgB3CciWhISpjlpbh771OWsKsniwYd38a3nD0/Z8e9AQyd3fu81Xj/Rwr/ctVrbJ0RJ5cIcvnLbCrYeOMfnH9tD/9D5PyNjDL986wz//MwhbrmokI8l0FqK2+XgQ5fOZ9vBc+yp7QB8RzkeaOjiw5ctiOpYPEmhBnxst2gb1qqIMeYQTPtjzQbguDHmpP/aR4A7gYPh3FtBdqqbX95/KV976gA/2H6Cp3bX84mrFnH98gLm56QyODLKntpOHttZy5O768lJc/OL+y+17WHYdvWJqxbRPzzKt7ceYW9dJx/YMJ85aW6e3tPAy0ebubo8j3/9s7UJt5Zy/9WLeWRHLZ/9dTVfuW05//DkfhblpvG+ddGdjIQ6w/facKdtNJbBi4HacV/XAZdOdrGIPAA8ADB//uzbWm615CQn37p7NXesLeL//P4I//zMIf75mUN/ck2a28nHrljEX11fRlbq7NvMYwef3rSU1SVZ/MvzR/jGs74/n7wMD3938zIeuHoxzgQL9gCZyUl8/4Pr+PhPd/DnP6siKyWJhz68fiynHi1ji7ah5PAjMaAImjbgi8gLwETN0b9sjHnK6gEZYzYDmwEqKyvtt5UtRq5cmsuVS3M52dzD26faaOgcwONyUFGQwaWLc8hI1kAfaxvL8thYlkdDRz99QyMsnJuW8K0sLlmYw4tfuIbdNR2smz+HvAxP1McQSh1+oMWJ3RZtpw34xpgbwrxHPVA67usS/2MqAhbnpbM4Lz3Ww1BTKMpOifUQ4kp+RjI3rYzdgTuekAK+76PdUjrRmF7sAMpEZJGIuIF7gS1RuK9SSk3L7XQCwQV879gMPyJDiphwyzLfKyJ1wOXAMyKy1f94kYg8C2CMGQEeBLYCh4BHjTEHwhu2UkpZI5ReOoFcs92WXsKt0nkCeGKCxxuAW8d9/SzwbDj3UkqpSDifw5/5TluvTXP4ib1ipJRKeKHstA3k8G0W7zXgK6USWyh1+Lpoq5RSNpTkFERCW7S1Ww5fA75SKqGJCG6nI6iNV2M5fJttvdKAr5RKeMEeZB6o0rFZRkcDvlJKeYIN+P5LNYevlFI243YGF/ATcuOVUkrNBm6XI8SNV/aK+BrwlVIJL9gcvlbpKKWUTbldjqA2XnltuvNKA75SKuEFm8NnbONVZMYTKRrwlVIJL/iUju+j5vCVUspmPC5niBuv7EUDvlIq4YW68Upn+EopZTO+gB9Ee2Sv1uErpZQteZxB1uGPFenYK+JrwFdKJbzgUzpah6+UUrYUfB2+72NC5fBF5B4ROSAiXhGpnOK60yKyT0SqRaQqnHsqpZTVEqWXTlhn2gL7gfcB/zmDazcZY1rCvJ9SSlku6JSOTXP44R5ifgjs95tWSqnx3C4HI16D12twzCAxb7QOf0oG+L2I7BSRB6J0T6WUmpFgDzIP5PCdNlu1nXaGLyIvAIUTPPVlY8xTM7zPVcaYehHJB7aJyGFjzCuT3O8B4AGA+fPnz/DllVIqdIGDzAdHvCQnOae93q7dMqcN+MaYG8K9iTGm3v+xSUSeADYAEwZ8Y8xmYDNAZWWlmegapZSykscf5Geaxz+/aGuviB/xlI6IpIlIRuBz4EZ8i71KKRUXPGMz/JnttjUJWpb5XhGpAy4HnhGRrf7Hi0TkWf9lBcBrIrIHeBt4xhjzfDj3VUopK43l8IOc4c+6lM5UjDFPAE9M8HgDcKv/85PAmnDuo5RSkRTsou2oNxDw7RXxdaetUirhBRZtZz7D9320WbzXgK+UUsGmdAJ1+HYry9SAr5RKeMHn8H0fNaWjlFI2Ewj4Mz31yq69dDTgK6USXvA5fF20VUopW/IEncP3fdSAr5RSNuNx+XbazvTUK7vW4WvAV0olvFAXbbW1glJK2cz5gD+z1go6w1dKKZsKdqet0UVbpZSyp2CrdAL/L2jAV0opm0ly+gJ38O2RIzakiNCAr5RKeCKC2+WY8cYrba2glFI25nHO/CBzba2glFI25nYFE/C1SkcppWzL7XIEsfHK91Hr8JVSyoY8Qczwjc7wlVLKvkJL6dgr4mvAV0op/AF/pu2RE7EOX0S+LSKHRWSviDwhItmTXHeziBwRkeMi8sVw7qmUUpHgDqpKJzHr8LcBFxljVgNHgS9deIGIOIHvA7cAK4D7RGRFmPdVSilLBZPSGWuPbLMkflgB3xjze2PMiP/LN4GSCS7bABw3xpw0xgwBjwB3hnNfpZSymtvlDPrEK5vFe1wWvtYngF9P8HgxUDvu6zrg0sleREQeAB7wf9kjIkdCHE8u0BLi90aTjtN6dhmrjtNauSIS9jjlwZlfW/StkG8Tyfd0wWRPTBvwReQFoHCCp75sjHnKf82XgRHgl6GOMMAYsxnYHO7riEiVMaYy3NeJNB2n9ewyVh2ntewyTojdWKcN+MaYG6Z6XkQ+BtwOXG8Cxal/qh4oHfd1if8xpZRSURRulc7NwN8Cdxhj+ia5bAdQJiKLRMQN3AtsCee+Simlghdulc73gAxgm4hUi8hDACJSJCLPAvgXdR8EtgKHgEeNMQfCvO9MhJ0WihIdp/XsMlYdp7XsMk6I0Vhl4iyMUkqp2UZ32iqlVILQgK+UUgliVgV8Efknf5uHahH5vYgU+R8XEfl3f2uHvSKyLg7GOmFbChFZKCL9/t/D2LpIvI3T/9yX/O/pERG5KYbDRETuEZEDIuIVkcpxj8fb+znhOP3Pxc37eSER+UcRqR/3Pt4a6zGNZ5f2LSJyWkT2+d/DqqgPwBgza34BmeM+/wzwkP/zW4HnAAEuA96Kg7HeCLj8n38L+Jb/84XA/liPbwbjXAHsATzAIuAE4IzhOJcDFcB2oHLc4/H2fk42zrh6PycY9z8CX4j1OCYZm9P/fi0G3P73cUWsxzXJWE8DubG6/6ya4RtjusZ9mQYEVqTvBH5ufN4EskVkXtQHOI6ZWVuKmJtinHcCjxhjBo0xp4Dj+NpoxIQx5pAxJtRd2VEzxTjj6v20GW3fMkOzKuADiMg3RKQW+CDwVf/DE7V3KI722KbwCXw/gQQsEpHdIvKyiGyM1aAmMH6c8f6ejhev7+d4dng/H/Sn9n4iInNiPZhx7PDeBRjg9yKy099GJqqs7KUTFdO1ejDGfBn4soh8CV/9/9eiOsBxQmxLcRaYb4xpFZH1wJMisvKCn17iYZxRN5NxTiAu3894NNW4gR8A/4QvYP0T8H/wTQBUcK4yxtSLSD6+/UuHjTGvROvmtgv4ZppWD+P8EngWX8CPSXuH6cY6UVsKY8wgMOj/fKeInADKgYgt8IQyTmLwngbxZz/+e+Lu/ZxEzFuQzHTcIvJD4HcRHk4wYv7ezZQxpt7/sUlEnsCXjopawJ9VKR0RKRv35Z3AYf/nW4CP+Kt1LgM6jTFnoz7AcSZrSyEieeI7QwARWQyUASdjM8op22dsAe4VEY+ILMI3zrdjMcapxNv7OYW4fj8vWPN6L7A/VmOZgC3at4hImohkBD7HVxAR1ffRdjP8aXxTRCoAL3AG+JT/8WfxVeocB/qAj8dmeH/ie/gqMraJ79icN40xnwKuBr4uIsP4fh+fMsa0xW6YE4/TGHNARB4FDuJL9XzaGDMaq0GKyHuB/wDygGdEpNoYcxNx9n5ONs54ez8n8C8ishZfSuc08MmYjmYcY8yIiATatziBn5jotG8JVgHwhP/fkQv4lTHm+WgOQFsrKKVUgphVKR2llFKT04CvlFIJQgO+UkolCA34SimVIDTgK6VUgtCAr5RSCUIDvlJKJYj/D11Gb9llM5FdAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "def get_kappa(eps : float, beta2 : np.ndarray) -> np.ndarray :\n",
    "    return np.sqrt(eps - beta2)\n",
    "\n",
    "def dispersion_fun(d1 : float, d2 : float, eps1 : float, eps2 : float, beta2 : np.ndarray, polarization : str) -> None :\n",
    "    kappa_d_1 = (2 * np.pi) * d1 * get_kappa(eps1, beta2)\n",
    "    kappa_d_2 = (2 * np.pi) * d2 * get_kappa(eps2, beta2)\n",
    "    if polarization == 's' or polarization == 'TE':\n",
    "        C = 0.5 * (kappa_d_1 / kappa_d_2 + kappa_d_2 / kappa_d_1)\n",
    "    elif polarization == 'p' or polarization == 'TM':\n",
    "        C = 0.5 * ((eps2/eps1) * kappa_d_1 / kappa_d_2 + (eps1/eps2) * kappa_d_2 / kappa_d_1)\n",
    "    return np.cos(kappa_d_1) * np.cos(kappa_d_2) - C * np.sin(kappa_d_1) * np.sin(kappa_d_2)\n",
    "\n",
    "k_B = 0.1 # Bloch wavevector modulo in the units of 2*pi/\\Lambda\n",
    "d1, d2 = 0.3, 0.4 # layer thicknesses in the units of wavelength\n",
    "eps1, eps2 = 1, 4\n",
    "beta2 = np.linspace(-30, 1.5*max(eps1,eps2), 500, dtype=complex)\n",
    "lh_part = np.real(dispersion_fun(d1, d2, eps1, eps2, beta2, 's'))\n",
    "rh_part = np.full(beta2.shape, np.cos(2*np.pi*k_B))\n",
    "\n",
    "plt.plot(beta2, lh_part, label='left-hand part')\n",
    "plt.plot(beta2, rh_part, label='right-hand part')\n",
    "plt.ylim([-2, 2])\n",
    "plt.axvline(x=0, color='k')\n",
    "plt.axhline(y=0, color='k')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5566ad38",
   "metadata": {},
   "source": [
    "## Поля мод\n",
    "\n",
    "Уравнения на постоянные коэффициенты $\\mathcal{A}_{1,2m}^{\\pm}$ можно переписать, введя новую константу:\n",
    "\\begin{equation}\\tag{12}\n",
    "    \\begin{array}{c}\n",
    "\\mathcal{A}_{1m}^{+}=\\sqrt{\\mathcal{C}_{m}}\\left(1-e^{ik_{B}\\Lambda}e^{-i\\kappa_{1m}d_{1}}e^{i\\kappa_{2m}d_{2}}\\right)\\left(\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}-1\\right)\\\\\n",
    "\\mathcal{A}_{1m}^{-}=\\sqrt{\\mathcal{C}_{m}}\\left(1-e^{ik_{B}\\Lambda}e^{i\\kappa_{1m}d_{1}}e^{i\\kappa_{2m}d_{2}}\\right)\\left(\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}+1\\right)\\\\\n",
    "\\mathcal{A}_{2m}^{+}=\\dfrac{1}{2}\\sqrt{\\mathcal{C}_{m}}\\left[\\mathcal{A}_{1m}^{+}e^{i\\kappa_{1m}d_{1}}\\left(1+\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}\\right)+\\mathcal{A}_{1m}^{-}e^{-i\\kappa_{1m}d_{1}}\\left(1-\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}\\right)\\right]\\\\\n",
    "\\mathcal{A}_{2m}^{-}=\\dfrac{1}{2}\\sqrt{\\mathcal{C}_{m}}\\left[\\mathcal{A}_{1m}^{+}e^{i\\kappa_{1m}d_{1}}\\left(1-\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}\\right)+\\mathcal{A}_{1m}^{-}e^{-i\\kappa_{1m}d_{1}}\\left(1+\\dfrac{\\eta_{2}}{\\eta_{1}}\\dfrac{\\kappa_{1m}}{\\kappa_{2m}}\\right)\\right]\n",
    "\\end{array}\n",
    "\\end{equation}\n",
    "Эту константу можно найти из условия ортогональности:\n",
    "\\begin{split}\n",
    "    1 &= \\dfrac{1}{\\Lambda}\\int_{0}^{\\Lambda} \\dfrac{v_{m}\\left(z\\right)v_{m}^{*}\\left(z\\right)}{\\omega\\eta\\left(z\\right)}dz= \\\\\n",
    "    &=\\left\\{ \\dfrac{1}{\\omega\\eta_{1}}\\dfrac{d_{1}}{\\Lambda}\\left[\\left(\\left|\\mathcal{A}_{1m}^{+}\\right|^{2}e^{-\\Im\\kappa_{1m}d_{1}}+\\left|\\mathcal{A}_{1m}^{-}\\right|^{2}e^{\\Im\\kappa_{1m}d_{1}}\\right) \\dfrac{\\sinh\\left(\\Im\\kappa_{1m}d_{1}\\right)}{\\Im\\kappa_{1m}d_{1}} + 2\\Re\\left(\\mathcal{A}_{1m}^{+}\\mathcal{A}_{1m}^{-*}e^{i\\Re\\kappa_{1m}d_{1}}\\right)\\dfrac{\\sin\\left(\\Re\\kappa_{1m}d_{1}\\right)}{\\Re\\kappa_{1m}d_{1}}\\right]+\\right. \\\\\n",
    "    &\\left. + \\dfrac{1}{\\omega\\eta_{2}}\\dfrac{d_{2}}{\\Lambda}\\left[\\left(\\left|\\mathcal{A}_{2m}^{+}\\right|^{2}e^{-\\Im\\kappa_{2m}d_{2}}+\\left|\\mathcal{A}_{2m}^{-}\\right|^{2}e^{\\Im\\kappa_{2m}d_{2}}\\right) \\dfrac{\\sinh\\left(\\Im\\kappa_{2m}d_{2}\\right)}{\\Im\\kappa_{2m}d_{2}}+2\\Re\\left(\\mathcal{A}_{2m}^{+}\\mathcal{A}_{2m}^{-*}e^{i\\Re\\kappa_{2m}d_{2}}\\right)\\dfrac{\\sin\\left(\\Re\\kappa_{2m}d_{2}\\right)}{\\Re\\kappa_{2m}d_{2}}\\right] \\right\\} ^{-1}\n",
    "\\end{split}\n",
    "\n",
    "Разложение произвольного поля по найденным собственным решениям будет иметь вид:\n",
    "\\begin{split}\\tag{13}\n",
    "    F &= \\sum_{m=1}^{\\infty}v_{m}\\left(z\\right)\\left(a_{m}^{+}e^{i\\beta_{m}x}+a_{m}^{-}e^{-i\\beta_{m}x}\\right) \\\\\n",
    "    G=\\left(\\mp_{s,p}\\right)\\dfrac{1}{i\\omega\\eta}\\dfrac{\\partial F}{\\partial x} &= \\left(\\mp_{s,p}\\right) \\sum_{m=1}^{\\infty} v_{m}\\left(z\\right) \\dfrac{\\beta_{m}}{\\omega\\eta} \\left(a_{m}^{+}e^{i\\beta_{m}x}-a_{m}^{-}e^{-i\\beta_{m}x}\\right)\n",
    "\\end{split}\n",
    "где $a_m^{\\pm}$ - постоянные коэффициенты разложения."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "74fe4a47",
   "metadata": {},
   "source": [
    "## Предел эффективной среды\n",
    "\n",
    "Вернёмся к полученному дисперсионному уравнению и рассмотрим распространение волн в $\\Gamma$-точке. Применяя формулы для тригонометрических функций двойного угла уравнение можно факторизовать следующим образом:\n",
    "\\begin{split}\n",
    "    \\dfrac{\\eta_{1}}{\\kappa_{1}}\\dfrac{\\kappa_{2}}{\\eta_{2}}\\left[\\sin\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)\\cos\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)+\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\sin\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)\\cos\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)\\right]\\times \\\\\n",
    "    \\times\\left[\\cos\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)\\sin\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)+\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\sin\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)\\cos\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)\\right]=0\n",
    "\\end{split}\n",
    "Отсюда\n",
    "\\begin{equation*}\n",
    "    \\left[\\begin{array}{c}\n",
    "\\tan\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)=-\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\tan\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)\\\\\n",
    "\\tan\\left(\\dfrac{\\kappa_{2}d_{2}}{2}\\right)=-\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\tan\\left(\\dfrac{\\kappa_{1}d_{1}}{2}\\right)\n",
    "\\end{array}\\right.\n",
    "\\end{equation*}\n",
    "Можно показать, что эти два типа решений соответствуют симметричным и антисимметричным модам. Переходя к пределу $k\\Lambda\\rightarrow 0$ подставим $\\tan\\alpha\\approx\\alpha$, что даёт:\n",
    "\\begin{equation*}\n",
    "    \\left[\\begin{array}{c}\n",
    "\\kappa_{1}d_{1}=-\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\kappa_{2}d_{2}\\\\\n",
    "\\kappa_{2}d_{2}=-\\dfrac{\\kappa_{1}}{\\eta_{1}}\\dfrac{\\eta_{2}}{\\kappa_{2}}\\tan\\kappa_{1}d_{1}\n",
    "\\end{array}\\right.\\Rightarrow\\left[\\begin{array}{c}\n",
    "\\dfrac{d_{1}}{d_{2}}=-\\dfrac{\\eta_{2}}{\\eta_{1}}\\\\\n",
    "\\kappa_{2}^{2}d_{2}\\eta_{1}=-\\kappa_{1}^{2}d_{1}\\eta_{2}\n",
    "\\end{array}\\right.\n",
    "\\end{equation*}\n",
    "Первое уравнение решений не имеет. Подставляя во второе уравнение явный вид $\\kappa_{1,2}$, приходим к соотношению\n",
    "\\begin{equation}\\tag{14}\n",
    "    \\left(\\omega^{2}\\varepsilon_{2}\\mu_{2}-\\beta^{2}\\right)\\eta_{1}d_{2}=-\\left(\\omega^{2}\\varepsilon_{1}\\mu_{1}-\\beta^{2}\\right)\\eta_{2}d_{1}\n",
    "\\end{equation}\n",
    "В случае ТЕ поляризации для немагнитной среды получаем\n",
    "\\begin{equation}\\tag{15}\n",
    "    \\beta^{2}=\\omega^{2}\\mu_{0}\\dfrac{\\varepsilon_{1}d_{1}+\\varepsilon_{2}d_{2}}{d_{1}+d_{2}}\n",
    "\\end{equation}\n",
    "то есть, эффективная диэлектрическая проницаемость вычисляется как средняя проницаемость на периоде. Для ТМ поляризации\n",
    "\\begin{equation}\\tag{16}\n",
    "    \\beta^{2}  = \\omega^{2}\\mu_{0} \\left( \\dfrac{\\left(1/\\varepsilon_{1}\\right)d_{1} + \\left(1/\\varepsilon_{2}\\right) d_{2}}{d_{2}+d_{1}} \\right)^{-1}\n",
    "\\end{equation}\n",
    "и вычислять нужно обратную среднюю проницаемость.\n",
    "\n",
    "В случае $\\beta=0$, когда мы хотим рассмотреть распространение волн перпендикулярно слоям кристалла, раскладываем дисперсионное уравнение по малому параметру $k_B\\Lambda$, что дает не зависящее от поляризации решение\n",
    "\\begin{equation}\\tag{17}\n",
    "    k_{B}^{2}=\\omega^{2}\\mu_{0}\\dfrac{\\varepsilon_{1}d_{1}+\\varepsilon_{2}d_{2}}{d_{1}+d_{2}}\n",
    "\\end{equation}\n",
    "Таким образом, фотонный кристалл в приближении эффективной среды будет вести себя как одноосный однородный анизотропный кристалл, который в главных осях описывается тензором диэлектрической проницаемости\n",
    "\\begin{equation}\\tag{18}\n",
    "    \\hat{\\varepsilon}=\\mathrm{diag}\\left\\{ \\dfrac{\\varepsilon_{1}d_{1}+\\varepsilon_{2}d_{2}}{d_{1}+d_{2}};\\thinspace\\dfrac{\\varepsilon_{1}d_{1}+\\varepsilon_{2}d_{2}}{d_{1}+d_{2}};\\thinspace\\left(\\dfrac{\\left(1/\\varepsilon_{1}\\right)d_{1}+\\left(1/\\varepsilon_{2}\\right)d_{2}}{d_{2}+d_{1}}\\right)^{-1}\\right\\} \n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f72765d2",
   "metadata": {},
   "source": [
    "### Литература\n",
    "\n",
    "1. Botten, I. C., Craig, M. S., McPhedran, R. C., Adams, J. L., & Andrewartha, J. R. [The Dielectric Lamellar Diffraction Grating](https://www.tandfonline.com/doi/abs/10.1080/713820571). Optica Acta: International Journal of Optics, 28(3), 413–428 (1981)\n",
    "2. Tishchenko, A.V. [Phenomenological representation of deep and high contrast lamellar gratings by means of the modal method](https://link.springer.com/article/10.1007/s11082-005-1188-2#citeas). Opt Quant Electron 37, 309–330 (2005)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}