{ "cells": [ { "cell_type": "markdown", "id": "6c42728a", "metadata": {}, "source": [ "# 02 - Troullier-Martins 伪化实践\n", "\n", "上一讲我们依靠 `solve_ae_atom` 得到 Al 的全电子 (AE) 轨道、能级以及径向网格,为赝势构造提供了可靠的参考。本笔记将以这些 AE 结果为基准,演示如何应用 Troullier-Martins (TM) 伪化方法生成范数守恒的赝波函数,并为后续的势反演与 KB 转换打下基础。\n" ] }, { "cell_type": "markdown", "id": "8d7e4d7b", "metadata": {}, "source": [ "## TM 理论速览\n", "\n", "TM 方法将内区波函数写成\n", "\\[\n", "u_{\\text{PS}}(r) = r^{l+1}\\exp\\big(P(r)\\big), \\qquad P(r) = \\sum_{i=0}^{N} a_{2i} r^{2i}\n", "\\]\n", "通过匹配 AE 轨道在伪化半径 $r_c$ 处的函数值与若干阶导数,并约束内区范数,使得赝波函数在内区光滑、外区与 AE 完全一致。\n", "\n", "核心条件:\n", "1. **连续性**:$u_{\\text{PS}}^{(k)}(r_c) = u_{\\text{AE}}^{(k)}(r_c)$,其中 $k=0,1,2$(若 `continuity_orders=2`)或 $k=0\\ldots4$(若 `continuity_orders=4`)。\n", "2. **范数守恒**:$\\int_0^{r_c} |u_{\\text{PS}}(r)|^2\\,dr = \\int_0^{r_c} |u_{\\text{AE}}(r)|^2\\,dr$,保证散射相位与 AE 一致。\n", "3. **光滑势**:指数型多项式消除了内区节点,使最终的赝势可微并利于平面波展开。\n", "\n", "在实际使用中可以对 `continuity_orders` 与 $r_c$ 做折中选择:\n", "- **continuity_orders=2**:未知系数更少,数值更稳定,适用于深束缚态或粗网格;\n", "- **continuity_orders=4**:多匹配两阶导数,可降低力量常数不连续导致的能量漂移;\n", "- **较小 $r_c$**:赝势更“硬”,平面波截断能需求高,但更接近 AE;\n", "- **较大 $r_c$**:赝势更“软”,收敛快,需注意是否引入幽灵态。\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "2860f055", "metadata": { "execution": { "iopub.execute_input": "2025-12-15T08:09:49.078659Z", "iopub.status.busy": "2025-12-15T08:09:49.078541Z", "iopub.status.idle": "2025-12-15T08:09:49.763194Z", "shell.execute_reply": "2025-12-15T08:09:49.762735Z" } }, "outputs": [], "source": [ "# 依赖导入与环境配置\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import platform\n", "from atomppgen.ae_atom import solve_ae_atom, AEAtomResult\n", "from atomppgen.tm import tm_pseudize, eval_derivatives_at, TMResult\n", "\n", "# 中文字体配置(兼容多平台)\n", "if platform.system() == 'Darwin': # macOS\n", " plt.rcParams['font.sans-serif'] = ['Arial Unicode MS', 'Heiti TC', 'STHeiti']\n", "elif platform.system() == 'Windows':\n", " plt.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei']\n", "else: # Linux / Colab\n", " plt.rcParams['font.sans-serif'] = ['DejaVu Sans', 'WenQuanYi Micro Hei']\n", "plt.rcParams['axes.unicode_minus'] = False\n", "\n", "# Colab 环境检测与本地路径加载(增强健壮性)\n", "try:\n", " import google.colab\n", " IN_COLAB = True\n", "except ImportError:\n", " IN_COLAB = False\n", "\n", "if IN_COLAB:\n", " !pip install -q git+https://github.com/bud-primordium/AtomSCF.git@main\n", " !pip install -q git+https://github.com/bud-primordium/AtomPPGen.git@main\n", "else:\n", " import sys\n", " import os\n", " # 假设 Notebook 位于 docs/source/tutorials\n", " project_root = os.path.abspath(\"../../../src\")\n", " if os.path.exists(project_root) and project_root not in sys.path:\n", " sys.path.insert(0, project_root)\n", " print(f\"已添加本地源码路径: {project_root}\")\n" ] }, { "cell_type": "markdown", "id": "e40895ca", "metadata": {}, "source": [ "## 1. 准备全电子参考态\n", "\n", "复用上一讲的参数(Al, Z=13, PZ81, n=1200, span=6.8)求解 AE 问题。注意这里我们使用 `lmax=2` 以覆盖 s, p, d 通道。\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "50fc5247", "metadata": { "execution": { "iopub.execute_input": "2025-12-15T08:09:49.764790Z", "iopub.status.busy": "2025-12-15T08:09:49.764675Z", "iopub.status.idle": "2025-12-15T08:10:11.154212Z", "shell.execute_reply": "2025-12-15T08:10:11.153548Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "AE 求解完成: True, 网格点数: 1200\n" ] } ], "source": [ "Z = 13\n", "al_ae = solve_ae_atom(\n", " Z=Z,\n", " xc=\"PZ81\",\n", " lmax=2,\n", " grid_type=\"exp_transformed\",\n", " grid_params={\"n\": 1200, \"total_span\": 6.8},\n", " spin_mode=\"LDA\",\n", ")\n", "print(f\"AE 求解完成: {al_ae.converged}, 网格点数: {len(al_ae.r)}\")\n" ] }, { "cell_type": "markdown", "id": "b9a3d891", "metadata": {}, "source": [ "## 2. 执行 TM 伪化\n", "\n", "我们将对 s/p/d 三个角动量通道分别进行伪化。注意:\n", "- **s (l=0)**:对应束缚价态 3s,选取 $r_c=2.0$ Bohr\n", "- **p (l=1)**:此处使用 AE 结果中能量最高的参考态(列表末尾,通常接近 0 的散射态)以强调散射性质匹配,选取 $r_c=2.2$ Bohr\n", "- **d (l=2)**:作为散射投影通道,同样使用能量最高的参考态(通常为近零散射态),选取 $r_c=2.2$ Bohr\n", "\n", "这里统一使用 `continuity_orders=4`(匹配至四阶导数)以获得更平滑的势。\n" ] }, { "cell_type": "code", "execution_count": 3, "id": "dad28cea", "metadata": { "execution": { "iopub.execute_input": "2025-12-15T08:10:11.157115Z", "iopub.status.busy": "2025-12-15T08:10:11.156989Z", "iopub.status.idle": "2025-12-15T08:10:11.172089Z", "shell.execute_reply": "2025-12-15T08:10:11.171636Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "l=0 rc=2.0: Norm Error=1.80e-16\n", "l=1 rc=2.2: Norm Error=2.30e-15\n", "l=2 rc=2.2: Norm Error=4.23e-13\n" ] } ], "source": [ "rc_map = {0: 2.0, 1: 2.2, 2: 2.2}\n", "tm_results = {}\n", "\n", "for l, rc in rc_map.items():\n", " # 选取该通道能量最高的参考态(列表末尾;p/d 通道通常为近零散射态)\n", " u_ae = al_ae.u_by_l[l][-1]\n", " eps = al_ae.eps_by_l[l][-1]\n", " \n", " # 调用 TM 伪化\n", " tm_res = tm_pseudize(\n", " r=al_ae.r,\n", " w=al_ae.w,\n", " u_ae=u_ae,\n", " eps=eps,\n", " l=l,\n", " rc=rc,\n", " continuity_orders=4,\n", " )\n", " tm_results[l] = {'tm': tm_res, 'u_ae': u_ae}\n", " \n", " print(f\"l={l} rc={rc}: Norm Error={tm_res.norm_error:.2e}\")\n" ] }, { "cell_type": "markdown", "id": "3d4ca370", "metadata": {}, "source": [ "## 3. 验证连续性\n", "\n", "为了确认伪波函数在 $r_c$ 处的平滑拼接,我们可以利用 `eval_derivatives_at` 工具检查 AE 与 PS 在接缝处的各阶导数是否吻合。\n" ] }, { "cell_type": "code", "execution_count": 4, "id": "aa8611bb", "metadata": { "execution": { "iopub.execute_input": "2025-12-15T08:10:11.173419Z", "iopub.status.busy": "2025-12-15T08:10:11.173333Z", "iopub.status.idle": "2025-12-15T08:10:11.176248Z", "shell.execute_reply": "2025-12-15T08:10:11.175750Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "--- l=1 (p-channel) at rc=2.2 ---\n", "u_AE: -0.010658, u_PS: -0.010658\n", "u'_AE: -0.001403, u'_PS: -0.001403\n", "u''_AE: 0.007498, u''_PS: 0.007498\n" ] } ], "source": [ "# 以 p 通道为例,检查 rc 处的导数匹配情况\n", "l_check = 1\n", "tm_p = tm_results[l_check]['tm']\n", "u_ae_p = tm_results[l_check]['u_ae']\n", "rc_p = tm_p.rc\n", "\n", "# 计算 AE 在 rc 处的导数(利用局部样条)\n", "derivs_ae = eval_derivatives_at(al_ae.r, u_ae_p, x0=rc_p, max_order=4)\n", "# TM 结果中已包含 PS 在 rc 处的导数信息(用于构建方程组)\n", "derivs_ps_check = tm_p.continuity_check\n", "\n", "print(f\"--- l={l_check} (p-channel) at rc={rc_p} ---\")\n", "print(f\"u_AE: {derivs_ae[0]:.6f}, u_PS: {derivs_ps_check['u']['ps']:.6f}\")\n", "print(f\"u'_AE: {derivs_ae[1]:.6f}, u'_PS: {derivs_ps_check['du']['ps']:.6f}\")\n", "print(f\"u''_AE: {derivs_ae[2]:.6f}, u''_PS: {derivs_ps_check['d2u']['ps']:.6f}\")\n" ] }, { "cell_type": "markdown", "id": "4f44effc", "metadata": {}, "source": [ "## 4. 可视化对比\n", "\n", "最后,我们将全电子波函数 $u_{\\text{AE}}(r)$ 与生成的 TM 伪波函数 $u_{\\text{PS}}(r)$ 绘制在同一张图上,直观展示内区无节点特性与外区重合度。\n" ] }, { "cell_type": "code", "execution_count": 5, "id": "cf056ed0", "metadata": { "execution": { "iopub.execute_input": "2025-12-15T08:10:11.177330Z", "iopub.status.busy": "2025-12-15T08:10:11.177257Z", "iopub.status.idle": "2025-12-15T08:10:11.376207Z", "shell.execute_reply": "2025-12-15T08:10:11.375714Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABdoAAAGOCAYAAACAIrHMAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjcsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvTLEjVAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAsjxJREFUeJzs3Qd8U1X7B/BfZtO994Cy9xIBAcWBgopb3Ih7b8W999/x6uvrFhRcqLhxD1BRUEFEpuxN6d4z4/4/z0mTpouW0jRN+vv6ufbk5t6bk5OQkzz33OfoNE3TQEREREREREREREREbaJv225ERERERERERERERMRAOxERERERERERERHRAeKIdiIiIiIiIiIiIiKiA8BAOxERERERERERERHRAWCgnYiIiIiIiIiIiIjoABgPZGciIiIiIiIiIiIi8i273Q6r1cqXwUtMJhMMBsM+t2GgnYiIiIiIiIiIiMgPaZqGvXv3oqioyNdVCXhRUVFISkqCTqdr8n4G2omIiIiIiIiIiIj8kCvInpCQgJCQkGaDwHRgJzMqKiqQk5OjbicnJze5HQPtRERERERERERERH6YLsYVZI+NjfV1dQJacHCw+ivBdmnvptLIcDJUIiIiIiKiVoxkIiIiIupMXDnZZSQ7eZ+rnZvLhc9AOxEREVEnM336dHXJ59NPP93ovsMPP1wtbXXfffepY1977bVN3j979mx1/76Wb775Bh3lp59+arE+smzbtq3etgsWLGjyeP/++2+9fZojbdzwMcxmMzIyMnDFFVegsLAQncUFF1yA7t27I9C9++676NatG4KCgnDZZZd16GN/9tlnmDZtmvu2670mf4mIiIh8jeliOkc7M3UMERERUSdSWlqKDz/8ECNGjMCrr76Km266qd2+ODscDsyZM0cd+6233sL//d//NTv65YUXXlDbNaV///7oKFKHJUuWuG9/+eWXePjhhxvVT/IkugLn0l4ffPABjjzyyEbHe++991r92DLR0SeffOK+XV1djVWrVuG2227DP//8g8WLF/NHTQeSfwtyQuG1115Djx49OvKh8cwzzzT5vhwwYECH1oOIiIiIOi+OaCciIiLqRN5//331V4LgGzZsaHZkdlvIsbZv347nnnsOJSUlmDt3brPbSgBxzJgxTS6RkZHoKBEREfUeu2fPnk3WT0Y5u4wbNw4ff/yxylnZVPvK9q0hx/R8jAkTJuCaa67BnXfeid9//x1//PFHOz5Takl2djYmTZqEY445Br169fJpg7nel/KXiIiIiNrOdUXtK6+80uiqzeauZr3//vs7ZZMz0E5ERETUibz++utqJPZRRx2FzMzMRl84G3J9AW3tsQcNGqQC0RI0bunY++vRRx+FyWRCfn5+vfXPPvusSrsi62VU/QMPPKCem0woJIFz+aIs69vLGWecgdzc3EZpPWQU+qZNm3D66acf0PFdI+m3bt2q/m7ZsgWnnHKKmoBKAq/Svj/++GO9ff7++28VJA4LC0N0dDTOPvts7N69u9EPjIbpbGQEt7zGLnl5eTjvvPMQFRWllquuuqrJHJEzZ87E4MGD1ckCGe3fmnQ3qampuPzyy9235USFnFQ5+OCD620nAeZzzz1XleWYN9xwgxphbrFYEBcXh5NPPlmd0Gnte8LVhqeddpp6TuHh4Tj++OOxbt26emlaxIMPPuhup6ZS5sh6uV/a09WuUq8///wThxxyiCpL+p/HH3+83n4yiZicRElJSUFoaChGjRqFL774wp1G6Oeff1aLK11MU6lj5GqHE044Qb0H5Bjy7/jXX39tVDdJQzN16lT1PGNiYnDRRRepK1k83ysTJ05UbS/3y4mF5cuX7/O1IyIiIvJXc+fORb9+/fDmm282uu+ss85S3+sbLrfeeis6IwbaiYiIiDoJyR8u6ShcOdolkPjpp5+qkbzNueeee+qlVmmOBBIlDYocW1x44YVYunRpswE8CbLabLZGS1OjxF0k+Cr3f/TRR42+PB977LEqEP3kk0/iv//9L+6++26V610CwI888gieeOIJtBcJMEt6G0kf0zBtjAQwpR4HwhUAlsC0kKBwcXGx+nEgI+klODplyhR30FyuTDjssMOQk5OjUve88cYbKkh/6KGHqv32x6mnnopvv/1WpTKR9D/ynmn4PO+9915ceuml6mTK559/rm5LOiJ5vIqKimaPfdxxx9U7QSDvD7nyQQK/rkCwBMZlvQSUxUknnYSvv/5avYbff/+9CqxL8Nl1cqA17wk5eSB1k/fi//73P3XVgUw8OnbsWBWw90wfdPHFF6uynDxoLTkRISdf5ASFpB6S495xxx2qTYTUT15DyQEv70tpM/mxJycMJLj+4osvYvjw4WqRx24qpZIE2SWQv2fPHsyaNUsdSxxxxBGNTrpccskl6kSTPI6kw5H3g9RHlJWVqcC6vIekfm+//bZ7XXl5eaufMxEREZE/yMnJUd+V5IpRScu4efPmevfLoBEZyNFw6bSTv2pERERE1CnceuutWkxMjFZdXa1ub9++XdPr9dojjzzi3mbChAlq2V8vvviiZjQataysLHW7vLxci4iI0C699NJ6273xxhuafEVsbhk4cOA+H0fqdsQRR7hvb968We03b948dXvKlCna5MmT6+3zv//9T3vrrbda9Txc9Vu4cGGj+2Sd6777779fi4+P12w2m/v+zMxMbfbs2e5jbN26dZ/Po1u3bprVanUvOTk52ty5c7Xo6Ght0KBB6th5eXnqWO+9955736KiIu3GG2/UVq9erW6fe+656nWVbT23SUhI0B566KF6z6thnaQO06dPV+Vvv/1WbTN//nz3/ZWVlep5ynZCHsNsNjd6XX/++We173PPPdfsc/7000/VNtu2bVO35X3Xr18/te7LL79U69599131PiosLNT27t2rHXroodqvv/5a7zjXXXedFhoa2ur3xF133aWZTCZt/fr17m2kveW9dvHFF7vXyT733Xef+7a0i+t5u0j7yXbSnp7tOnPmTPc2paWl6jlcc8016ra0p2wjz9/FbrdrY8aM0R588MEm/915vtfEaaedpsXGxmrFxcX1nkOfPn20ESNG1KvbeeedV6/O8jjyfhLLli1T2/z+++/u+zdt2qTNmDFD27lzZzOvHBEREXVV8l1w7dq16q8/ev7559V3KPnulZaWpt177731vuu5vgf7S3tzRDsRERFRJyCjamWEsqTPqKmpUaNYZVSrjEqWyR8PNLWKKyWNpC6RY8vxJN2JjCz2TFvh8vLLL6uRyw0XVw755kybNk2NAt67d697FLmkwHCNgJYR5TKSXUZ4yyh2SeciKTtktHF7cqWPWbhwobot+dRltLE859aS0dSS9sS1JCQkqBHakjrlq6++gsFgUCOyZaSzjFKW5y4jkGXS1P/85z8YOHCgOo6M0pk8ebJqB9eVAa70JPuTg1/ywku6FRl57iKpUFxt69pG3j/nnHNOvX2lvbt166ZeG9HwSgWJY8trI6OGfvjhB3e95f3Yp08f937y2snoc0nxkpiYiF9++UWlypG2lbZ+6aWX8Ntvv9Ubfd3Se0IeR9pC0s+46uNK2dJecxTIaHMX+TcgdXfVUeorr69nO+r1ejV6Xa4YaQ15fnIVg2fOdqPRqFIEyUh9z39jnnURkj7JVRe5EkOulJAR9pLGR64EkFFb8m8lLS3tAFqAiIiIqPN577331Hco+e4lV0rKd2nn+Ar/ZPR1BYiIiIgIKnCblZWlguqyNCTpQiTVRlusXr0ay5YtU2XJC92QfKG98sor663r27cvRo4cud+PJfnPJXAuaS/kr3x5lnzUrslKr7vuOhXolHQct99+O2677Tb07t1bTdAqwej2IgFLSSEjaVUkgCwnCKT99mfySklPIuk9hKTykecgOcGl/p4k+CypU+QxpC0lAC8/GOR1jI+PV6lRJJWIK51Iw3ZuLcnpLkFX+SHiSfKKu7jysCclJTXaX9ZJqhpJaSOpSzxJ+hJJ9yLBbQl8y4kPuXxXLuOV+ksgWX70fPfdd/VyYsoJnIcfflilwpHc8wcddJAK/u/Pe0KOL7nzJdjdkOcktwei4eXFEgR3nbySx5cTJg3bdX9IuzfX5kJS8LSmLnKfnLx46KGH1Hv31VdfVW0gJ06ef/75znuZNBEREXUaLyzchJKqxnP4eFuExYSrj2j9hPU7d+5UAx4klZ6QtH0vvPCCmuNGBnYI+f4sqTQbWrt2bb3vwJ0FA+1EREREnYAEOiWIK6PaPUlw88QTT1QjzNsaaJdgqIyglvzUEgT2JCOxZVLUhoH2tpKRylJfCTrLCHrJXS15t10kYC15tmWRnIwyelpyiLtGoLdXYFWceeaZatJN+cI+b948PPXUU/u1v4web83JBhnpLnnn5bHkpIaMQpbgs+T7lraV0d8S7L/55psb7et6vq7XpWEOfM+R0BLIloCujPiW4KyLBIo9txEyelzyjHuSUecy+lx+lMjVCZ5cgXcZSS0nDeRHjwR/JZ95QUGByjsuP3rkuK6R3xKIl/eP/Di69tpr1Yh5cfXVV6v9W/uekPaRKzeaen32NdGvtNm+2qu1pM1cE/V6Bttd8xc0lZO9qWO4Ruw3bHM5ptwv7/fWkJH98nkgE9r+9ddf6jNBguxyQsqVy52IiIioORJkL6l0XiHYmb333ntqgMakSZPUbRnwId+ZZE4jV6BdvkM29R1Rrk7sjBhoJyIiIvIxCTB/8cUXuOWWWzB+/PhG90sQWgKdMqJ5f8lEkDLKWr6kSjCzIZkcVUYtS8oRSYnSHiRViFz6KQFuCb5K2hIXmRxSRprLCHYJUMtIXRnJL8+9srKyXQPt0m4S7H7sscdUgNozNUh7kWCsnACZP3++Sn8iz00WCba7JrGVdpdRN8OGDXMHyOUEikxIK+llZHtXgFz2kVQiQvaRILeLBOv/7//+Tz2WKwWOpImRyUhdRo8erU4QyOgf+bHismjRIjVqSH607OsEgozEl6sOJBAu74fg4GD1mklAW06ISPC+Vy/nSCVJrSLP47777nNfKSGTrcqkqMIzcL2v94S0zzvvvKNS1HhecSAj5+UYMkq+KdJm0j6eJx7kee4veXyZpFeuKpHnL+R5XXTRRSq4LSdp5HnsayJgeT7yb1hGrrueg2wvJxek/q0diS6vrbwv5ISNjIaX95Qs8oNzX5MiExEREXmOLPeHx507dy6qqqpUukwX+W4rV0HKIAMhV5LKYCR/wUA7ERERkY/JiFUJiEs+56ZIMFxSSDSVUmbz5s0qUN9ckFwCd3J/c8eWAKgEo2XEvOcxJMjbMAWIZ0oV1+jlpkgKGEnFIXWeMWNGvVHJEvyVEdMyqloeb8eOHSqfuQSRZWRze5IgqeRPf/TRR1Vg2htpNyRILiO2Je2KBJzl5IGk+ZFAqaTGEbJeAuAyWlyuHJB6yGjljz/+WKWdEXK/rL/hhhvwwAMPqADygw8+qI7nctRRR6kRP3I1gIxil5zdcsJCRmNLewtJLSMBahlRL6lYJLgtqWLkNR4wYIB6L+2LjGyXtDtyia7U23XMIUOG4KefflKvp4srAC6j2eW4EgiWUf2u9DUSdHel2dnXe0JG+su/AXkPyAkX2U5OVEi+d7kaozkyOv7pp59WAXFpfznpIYH8plLQ7IucKJG86XIMea/IiHI5UbFu3TpVByHvTRnBLznj5YRJQ9JWMvpf3t8y6lxOZkh6pC1bttQ7EdISqYcE9U899VSVVknaT05CSA53WUdERETUkv1J3+IrGzZswN9//62+L8mgDhf5Di0pBptKF+MXOnx6ViIiIiKqZ9CgQdrAgQP32Sq9e/fWUlNTtfHjx2sTJkxwr58+fbrMFtTsflOmTNGio6O16urqZreZOHGiFhwcrBUUFGhvvPGGOt6+luuvv77FV/Daa69V265Zs6beervdrj366KNanz59NIvFoqWlpWlXXnmleuzWcNVv4cKFje6TdQ3ve/zxx9W6Tz/9tNExtm7d2uzjSBt369atVXWS45x11llaYmKiasfBgwdrs2fPrrfNsmXLtGOOOUYLCwvTwsPD1ev47bff1tvmu+++04YPH64FBQWp9nn//fe1448/Xr3GLmVlZaq9YmJi1HEuueQS7ZZbbmlU15kzZ6p6mEwmLSkpSbvmmmta3cYzZsxo1I433XSTWvfzzz/X2/aFF17Qevbsqercr18/7ZlnnnG/Dj/88EOr3hNi48aN2qmnnqpFRUVpISEh2ogRI7R333233jay73333VdvnTyePHd5fGnTFStWaKGhoeo13tdrLft4tqu0zaWXXqrFx8erxx83bpy2YMEC9/0//fST2sdsNmtz585t8r22cuVK7cQTT1Svi7wP5N/Vb7/95r5f6iD7uOrmIvXwfP3kOcjrHhcXp57LqFGjtPnz5zfxShEREVFXV1lZqa1du1b99Sf333+/Fhsb2+RvFPlddOyxx6rvSPIdOzc3t9FSXFzcKdtbJ//zdbCfiIiIiIiIiIiIiFpPUq/IpPRyVWJzV6N2Rv3798cxxxyj5jlqSNLGyFWecrWjXCnaFLnK03VlaGdqbwbaiYiIiIiIiIiIiPyMvwbaA7W9nbMTERERERERERERERFRmzDQTkRERERERERERER0ABhoJyIiIiIiIiIiIiI6AAy0ExEREREREREREREdAAbaiYiIiIiIiIiIiIgOAAPtREREREREREREREQHgIF2IiIiIiIiIiIiIqIDwEA7EREREREREREREdEBYKCdiIiIiIiIiIiIiLzu8MMPh06na3bZtm0bLrjgAlV+9NFHG+1fUlKCoKAgdO/evcnjy/4NjxkcHIyhQ4fi888/r7ftzJkzMXjwYISGhqrjXX311cjNzW3zc2OgnYiIiIiIiIiIiIi87rPPPlPBbFlmzJiBsWPHum/LkpGRobYzmUyNAuPiyy+/hN1ub/Fxli9f7j7munXrMGnSJJx22mnYsGGDuv/VV1/Fgw8+iMcffxzr16/He++9hxUrVuCoo45q1fGbwkA7EREREREREREREXldZGQk4uLi1BISEqIC6q7bsuj1znD1oYceiqVLlyI7O7ve/p9++qm6ryXR0dHuY8podRkdLyPbJVAv3njjDdx00004/vjjkZaWhjFjxmDOnDlYtWqVCtK3BQPtRERERERERERERNRppKenY9iwYZg/f757XXV1Nb7//nuccMIJ+308g8EATdNUsF04HA789NNPqKqqcm/Tq1cvrF27FkOGDGlTnRloJyIiIiIiIiIiIqJO5aSTTqqXPubHH39UwfeYmJj9Ok5FRQXuuOMOlRJGUsiIK664QqWxkdHs06dPx9tvv43CwkL0799f5YBvCwbaiYiIiIiIiIiIiKjTBdp/+OEHVFZWutPGnHzyya3aV0alR0VFqVQ1YWFhmDdvHubOnYvMzEx1/4UXXoivvvoKw4cPx7vvvotp06YhOTkZDz30UJvra2zznkRERERERERERETUuSx+HljyQsvbJQ8Fznmv/rp3zwKy/ml530OuBsZeA28aOnQoEhISVLqYKVOmqDQyd999NxYsWNDivhJEl9HqOp0O4eHhTY6CP/bYY9VSXFysAvovvvgi7r33XgwYMEBNnLq/GGgnIiIiIiIiIiIiChTVpUDpnpa3i0xtvK4ir3X7ymN0YPqYuLg4pKSkICMjo1X7SZBdJkFtyo4dO/DEE0/gySefVDnbZdS7BNZPPfVUjB49GgsXLmSgnYiIiIiIiIiIiKhLCwoHwlNa3i4krul1rdlXHqODAu3nnXeeGpF+yimntMsxjUYjXnjhBRx99NHq+C4y+l3SzURHR7ftuO1SOyIiIiIiIiIiIiLyPUnp0ta0Lg1TyfjYYYcdhqqqKrz22mtYtGhRuxxTRsbLBKgXX3wxcnNzceSRR8Jms6kJUZcuXYpXX321TcflZKhERERERERERERE1OkYjUYcd9xxKnXMoEGD2u24Eri/8cYb8fTTT6N///4YP348Vq1ahd9//73ZlDMt0WmaprVbDYmIiIiIiIiIiIjI62Sk99atW5GZmQmLxcIW93F7c0Q7EREREREREREREdEBYKCdiIiIiIiIiIiIiOgAMNBORERERERERERERHQAGGgnIkUmerjgggv8rjUOP/xwtbSkoKAA3bp1w7///ntAz3fv3r0488wzERkZqZbTTz8du3fvdt9/991344orrtjv4xIRUdfDvnf/Sb/bsP9m30tEFPhmz54NnU6Hbdu2HdBx2Pe2zj///IMTTjgBCQkJagLKSZMmqXUu7HuJmsZAOxF1Cdddd536otCvX782H8Nms6kvGCtXrsScOXPw+uuvqxmpjz76aFitVrXN7bffjs8//xzffPNNO9aeiIioa/a9Lg6HAzNmzMBHH33U6D72vURERO3X927YsAHjx49HWVkZZs6cqX73VlZWYuzYsdi0aRP7XqJ9MO7rTiKiQLB8+XJ88MEH2LFjxwEdZ968eSrILmfyhwwZotb16dNHleX45557LsLCwnDVVVfhjjvuwOTJk9vpGRAREXXNvlfISe3rr78ef/zxB4KDgxvdz76XiIio/fre//znP6pv/fLLLxESEqLWHXXUUepqgGeffRbPP/88+16iZnBEO5Gf+fvvvzFx4kSVtiQmJgbHHHOM6lBbUlRUhGuuuQYpKSkIDQ3FqFGj8MUXX9TbRkZl33rrrUhKSlIdqozUdqVacXnzzTdx8MEHq8eX4wwfPlwFoD0v6bNYLPjzzz9xyCGHqHJGRgYef/xx9zZyuZ9c9vfZZ59h6tSpCA8PV8/loosuQmlpab3He+GFF9C/f38EBQWhZ8+eeOqpp9Sotv3xyCOPYMqUKep57esyxOaW+++/X2339ddfuwPrLoMHD1brPEewX3zxxSoowFHtRET+T35UPvDAA7jyyitV3xcVFYXzzjsPeXl5Le7LvvfA+14xbdo0VFRUqO8Wcgl7U9j3EhEFDrvdjvvuuw+pqanqBKsMYPJM17kv7HsPvO/t27evOsHtCrIL+e0vr8f27dvd69j3EjVBIyK/UVpaqsXFxWlTp07VvvvuO+3LL7/UDjnkEC02NlYrKytrdj+bzaaNHTtWi46O1l544QXthx9+0KZNm6YZDAbtp59+Utt069ZN0+v12sSJE9Vx58yZo7Y/6KCD3MeZNWuW2uehhx7SFi5cqM2fP18bP3682m/Lli1qmzfeeEPdluM9//zz6rHOPPNMTT5u5s2bp7bZunWrui3PZcaMGdqCBQvUMWXd1Vdf7X68Rx99VK278sortW+//VZ75plntMjISO3GG290bzNhwgS1NEfaxWKxaK+//nq99VK/6dOnq3JOTo62ZMmSZpedO3eq7Q4++GDthBNOaPQYU6ZM0UaNGlVv3YgRI7SLLrqo2XoREZF/kP4iKipK9befffaZ9tprr2kxMTGq/9sX9r3t0/eKFStWNHmMhtj3EhEFhrvvvlv97nzwwQe1b775Rrviiis0o9GofhvKb8nmsO9tv763oc2bN6vX4M4776y3nn2v71VWVmpr165Vf8n37c1AO5EfWbZsmfpy8fvvv7vXbdq0SQWr99UpSkBc9vv000/d6+x2uzZmzBj15cXVAaelpWk1NTXube666y61X0lJibp98803a7fddlu9Yy9fvrxeEF0C7XJ75syZ9U4QSKd8zTXX1Au0n3feefWOJfUZNGiQKhcXF2shISHa+eefX2+bd999VwXyt23b1qpA+1dffaUea9WqVfXW7+uHenN69+7dqM7inHPO0fr06VNv3XXXXaelpqbu1/GJiKjzkf5CTgxLX+byySefqL5Ffvw3h31v+/S9Tb0ezR2DfS8Rkf8rKChQvwPlt2fDwU0tBdrZ93qn7y0vL1cDDkJDQ7WsrKx697Hv9T0G2jtXezNHO5EfkRQqcrnW8ccfj9NOO02ljZE0Mk888YS6X06eyWV2ngwGA3777TeYTCY1KYqLXq/HkiVL6m07ZswYtZ1Ljx491N/CwkKV3kXStojy8nKVUmbjxo1YsGCBe50nSRvjIvndEhMT97mNkNQwixcvVmWpm1wmLnnPZRJSF8kNJ6ljfv75Z5x//vkttpkrP520W3OaajdP0layyHZNpa2Ry+zkPk9paWnq8saamhqYzeYW60lERJ2X9LnSl7mceOKJqr+UnOHSF7Pv9V7fuz/Y9xIR+T9JwSm/A08++eR668844wx36lP+7u24vrekpASnnHKKSt/24YcfNkpLw76XqD7maCfyI5Ij7ZdfflEBc5nk5PTTT1cBbMltLl9GJPgsP/w9F1kneWRjY2Nb/MEqedc8ubZ3BZcluO7KDy9B8scee0zNPi4aBpo987kJo9HYKEi9r21cuW8nTZpU7/nI8xV79uxpVZvl5+erv3KioDlz5sxp1G6ey4MPPqi2k7y8DXPIi+LiYnWfJ9dt1+MTEZH/Sk9Pb9Q/Sr8qJ1TZ93q3790f7HuJiPyfKxe763efi8w15sK+t2P6XgneH3roofj1118xd+7cRic/BPteovo4op3Iz8go8zfeeAMzZ87EX3/9hbfeekvN+t27d2812enSpUvrbS8TmXz11Veq45Ugtmew3TWJ6ogRI1r12CeddJIKxktHe9BBB6nOWB7v7bffbudnWddhv/POO2qy0Yb2dabekwRCXJPixMXFNbmNnLho2G5NfamTeixbtqzR/Rs2bFCTy3pyBeRlklciIvJvDSc+ldFgsk76FekP2fd6r+/dH+x7iYj8X3R0tPqbk5OjfuM21Rez7/V+37ty5Uo16M1qteL777/HYYcd1uQ+7HupvdlsNjUI01/5b82JuqD58+fjwgsvxOrVq9UlWxLclUXOTGdnZ6uz1yNHjmy034QJE/Dkk0+qgPuUKVPcI9BlJLx8eZk3b16Ljy2BegkoP/fccyrFjMunn36q/jaVUuVAyGMEBQWpEQ3nnHOOe/3atWtx00034ZFHHkFycnKLx3Glv9m1a1ezXzjkS4nri8m+yBeNd999F//88w+GDh3qvrRR2uWee+6pt63UW0ZhyHMgIiL/JpeqS/o0SccmPvnkE/UjQNKZse/1bt+7P9j3EhH5P7lyWq58liu4x40b517/2Wefucvse73b927ZskV9x5FBdnL1QFMD31zY99KB2LZtGzIzM/H+++/jf//7n0pRJFdPnHrqqX7bsAy0E/nZlw4ZkS4fOrfddpvKFysjviX3+b4+iI499li17wUXXIBHH31UdcISMF63bh1eeumlVj22dMjdu3dXH36Sh02+/EjgQb4ACUld057k8WbMmKEC2JIX7vDDD0dWVhbuvfdelfN84MCBrTrO+PHjERwcrPLUDxs27IDqdNZZZ6l0OZIfUP5KwOX2229XdTnzzDPrbSt5eyUwT0RE/m/z5s0qP+mVV16J7du344477sBxxx2HI488stl92Pe2T9+7P9j3EhH5P0lTevfdd6tFAuryO/DLL79052ffF/a97dP3XnfddWqg3f/93/+hoKAAv//+u/u+iIgIDBgwwH2bfW/n5aiN0eiCg9W8ckKrqYEmc+AZjdB7zCXn3tZiga42C4JmtaoFBgP0HgMIm9tW5zHfX2vJIEYhA0MlViVB9/j4ePgz5mgn8iNyZlou25J0JJdccolK5bJmzRp1dr+5S7mEBOfly4kE4+ULi+wno7C/+eabRhOS7ouMXpfLyaZNm4bzzjsPZWVlKn2NBO6lg21vDz30kPrAlUlXJKAhgXf5orVw4UJYLJZWHUOC7DJ5rIzmP1AS4P/uu+/UF4vp06erdhg0aJBa5zmJrFzmKGdiGwbfiYjIP02dOlX1vTI3ipwAlj5A+qZ9Yd/bPn1va7HvJSIKHHJC++mnn1YpSuWKbElj8uyzz7a4H/veA+97ZRCfHEOugL/44otVvMBzueqqq9zbsu/t3NaPOEgt9sJC97r8119X67IfeqjethvGjVfrrXuy3OsK331Xrcu66+562246aqJaX7N5s3td0SeftKmOK1asUFdOSJaFo48+Gr169VIn2/yZTms4gyERUYCRvOqSimbr1q2NJrTzBjk58N5776nHdZ05JiIi/yRXc8lJ3tmzZ/u6Kn6FfS8RERH7XvK+qqoqFeuQ0eCeAxLX9euv/vZe/BuMtXPH5b38MnKf/S+ipp6OZI9g+7/DR0CrrETPH36AOc05H16BpCh+7HFETJmC1KeedG+74ZCxKnjfY/7nCKqdR6Hwgw8QfcYZ+113GQwqdZaMC/7e3i4MtBNRlyCjyyWv/X//+1+vPo6M8pe89zL6QvLaERGRf2Ogve3Y9xIREXUs9r1dT3OBX39IHdOjRw+Vjveyyy5zr5M0jddcc42ab0Am45UrLDIyMuAvgXamjiGiLuHFF19UqW/Wr1/v1cd5/PHHcfLJJzPITkREXR77XiIioo7Fvpdc9CEhavG8yl5nNjvXewTZ621bGzhX25pMzvUeQfZ9bbu/SkpK1GSow4cPd6+rqalRaYNvvvlm/P3331i0aBGSk5P96kXliHYiIiIiIiIiIiIiP9PSCOvOatGiRTjyyCNRWlrqrvf777+v5sCbNWsWOiuOaCciIiIiIiIiIiKiTuGff/5B3759650cWLVqFUaNGgV/xtQxRERERERERERERNQhrrnmGqxevbreusTERPc6u92OgoICv3s1GGgnIiIiIiIiIiIiIp+54IILsHnzZgwaNAgjR47Ehg0b/O7VMKKLcjgc2LNnD8LDw+tNDEBERNRVaZqmcuSlpKRA7zG5TXth30tERMS+l4iI2o9MICq/s2QEuCz+LCQkBPPnz6+3zhvPSX73SpuZTKZ2jwl32UC7BNnT09N9XQ0iIqJOZ+fOnUhLS2v347LvJSIiYt9LRETtp1u3bnj55ZdRWVnJZt1PQ4YMgdlsRnvqsoF2Gckutm/fjqioKF9XJ6DJWaLc3FzEx8d7ZYQksa19ge9rtnUgKioqUl/UXH1ke2Pf23H4GcW2DkR8X7OtAxH73sDBzyi2dSDi+7rzt7WMaM/Ozkb37t3rTSxKzbPZbGriVW/EKLtsoN11aUBERIRayLsfFlVVVaqdGWj3LrZ1x2Fbs60D9X0tvJVSjX1vx+FnFNs6EPF9zbYOROx7Awc/o9jWgYjv687f1rKPBOgNBoNaqHWpY7z1u5fDi4mIiIiIiIiIiIiIDgAD7UQBwlFTg7333ovyJ59SZSIiIvIP0m9n3XOPWtiHExERse8lIv/EQDtRoLDZUPzhR6j58ktVJiIiIj9hs6Fo3odqYR9ORETEvpeoranIyLft3GVztBMFGp3RiLjrrkNZebkqExERkX+Qfjv+huvdZSIiImLfS9QaZrNZ5XTfs2ePmkhVbntrzq1AmgzVld++tXntJa+7TDwr+fClvaWdm8Jv8kQBQmc2I/aKy2HPyVFlIiIi8g/Sb8ddcYWvq0FERNRlsO+lQCFB38zMTGRlZalgO7VuVHpeXh62bdu2XxPPipCQEGRkZDS7HwPtRERERERERERERH5IRldL8FdGatvtdl9Xp9MrKSnB8ccfj2XLliEsLKzV+8nod6PRuM8rBhhoJwoQchmLraAAjqIiaPHxvq4OERER7Ucfbi8sVGVDdDQv9yUiIvIy9r0UaCT4azKZ1EL7Jiljtm/frk5QWCwWtCcG2okChFZZic3jD1XlhGVLgf04K0dERES+7cM3jh2nyn2X/wVdSAhfDiIiIva9RORnGGgnClCaw4HSigqU2Qwoq7ahvNqGGpsD8Wteh2arAhx26OAAjMGA0QK9ORiGkGgYwuMRFJGA0IRusISE+/ppEBERERERERERdXoMtBMFQkB97yYUbfkLpgdOhn7vShS/MA4RVXuwLP5sfJ90ab3tH1z1BExaTYvHnZd2J9YlnoDIYBNiw8xICapGr4KfEZwxDNHdBsNobt/La4iIiLoqfUgI+v+7ztfVICIi6jLY9xKRNzDQTuSHivZsQtE/X0G3YzFi8pYiwpqHiCa2i67JarTOoTMAWsuPUW6MQqXVrpa9JVUoK/8HR26+Rd1n05mQHdobFbEDoU8bgYg+hyEqfQB0+zlbMxERERERERERUSBgoJ3ITyZqyS6pxurdxVi9pxiT/rke/UsXN7t9jS4IpZYUmKOSMSozGqFmI8KCjAgy6ZGV8F8YdIDOYIRD0wG2KmjWKjhqyuGoLIK+PBe6ynyEJPZGjN6EkkobbA4NsdV73Mc3alYklq0FZNk+D/gNKDXFIj/2YDi6H4qocRchJpz5ZYmIiIiIiIiIqGtgoJ2oE6sozELW4vcx3zwZ2aVW9/p1EePdgfYafTByo4ehMnYotN+yoZnD0f3hxxEbGopYAAMbHrTbma167O61fx0ODYUVNSjcYcOmDTaYclYjvHA1oip3QO8xND7cmo/wvd+gOG8ZHi8/DHHhQeidGI4+iWHoFR8Go4Gj3YmIiJriqKlBzlNPqXLCLbdAbzazoYiIiLyIfS8ReQMD7USdUO7GP1H28/NI3/0VempWhPaIBsJGuO8vyZiILTFWhPQ/Bgl9Dkaq0QRHRQXWP3iQul//UCtyw7SSXq9DbFgQYgcMB2SpVV1ehNxNy1G9ZTGCdi1GfOHfCHJUYGvYcECnQ15ZDfLK8rFkcz6m7bgbYRYzdANOROLIk2HmJKtERER1bDYUvvmWKibccAPAQDsREZF3se8lIi9goJ2oE9m94jvofv4/pBQuQ7zH+uGF38CWMR5D0yMxMCVSTVAKjK63r85oRMxll6GiolyVvS0oNAppQ48EZJERATYrsjf8iZpiG7rXhGBHQQUcGhBkL0PfokUwwA7s/Q41P92CbYlHQDf4NKSOPIGTqhIRUZcn/Xbs5Ze7+3MiIiLyLva9ROQN/CZP1Alkr1oA+4+PILVoWb31lYZw7OkxFZnjLsbI7j33eQyd2Yz4G65HTk6OKnc0vdGExAHjkAhgFIAqqx2bcsqwd+2vqDGEINheqrYzO6rQPetrIOtrVP4Yjl1pkxA8ajoS+o/nZKpERNQlSb+dcOMNvq4GERFRl8G+l4i8gYF2Ih8q2rUeJZ/fjoycBfXWF1jSUTT0MqRNuBA9/TTNisVkwKDUSAxKPR62w7dg58oFsK38EIm7vkWIvURtI8H37ts/BLZ/iLzgTKw/bh6G9slUE7cSERERERERERH5C0aziHzAanfg5/W5yFr6M6Z5BNkLLBkoHX0j0g+dhhijpIdpPU3TVJ52rbJSlTsTo8mM9IMmAwdNhsNajd1/fw3rinlI3vsjghyVaptqzYQvNlbh683rMCglEqMyY5AZFwqdTufr6hMREXmV9NvSfwtdcDD7PiIiIi9j30tE3sBAO1EH25Jbho+W70JBuRUIG4ONYQcjuXoLcg6egW5HXrLfAXYX+YG+ceTBqhy/bCkQFobOSG8KQuqok4FRJ6OqvARbl7yP4FXvYFnoEep+uwP4Z1exWs7Pehjm7ocg9fALYQmN9HXViYiIvEL68PUjnBOa913+F3QhIWxpIiIiL2LfS0TewEA7UQepLs3Huh/m4H3taPc6vV6HXYc+ifRBPdCjCwaSLaERyJx4KTDxUowvrULQ9kIs21aI8ho70irWon/uN0DuN6j66wlszTwNMUdchci0/r6uNhERERERERERUT16dDKVlZW46KKLEB4ejsTERDz66KPNbrtjxw4cd9xxCA0NRa9evfDBBx90aF2JWitn9Y+o+d8YDPvnAfQv/lWt6x4bguuP6o0jRg9vl9Hacql572VLEfX1V6rsb2LDLZg8KBm3H9sPZx2cjlH2Fe77LI5yZG5+E+EzD8HOl05F9lpnGxIREQUC6bdlJLsaze6HfTgREZG/Yd9LRF1iRPuMGTOwcuVK/Pnnn9i9ezemTp2KzMxMnH322fW2s9vtmDJlCnr27Km2//XXX3Heeeehf//+GDx4sM/qT+RJs9uw87OHkLbyOejhUOsmZ7+MXuNPwyG94ts1B6scSx8S4ve5XY0GPYamRwHTH0P+ptNQ+suLSN05HyatBnpoSM/+EfjgR+yJHIGaMdciY9RJ0BsMvq42ERFRm0m/zXQxREREHYd9LxEFfKBdRrPPmjULX3/9tQqYy3LjjTfipZdeahRo/+qrr7Bz504sXrwYYWFhKuD+8ccfY8mSJQy0U6dQXbgbhW9NR0bBUve6XZEHIeTMmRibkuDTuvmL2F4jEdvrdZQV5WL7gpeRuG42wq156r6U4uXAtxfinxWnwn78fzAsLUql4iEiIiIiIiIiIurSgfYVK1bAarVi7Nix7nXjxo3DY489pmaE9hyl+9NPP2HixIkqyO7y2WefdXidiZpSsul36D84F0k1zqCwA3psGnANepxyL4ymtk122hKtpga5z7+AyopyaLfcAlgsAfPihEXFo9ep98BafQu2/jIHEctfQmzlNnXfn6ETsGXZLvy4LhsT+iRgREaUGhVPRETkL1Qf/sKLqhx/9VXQmc2+rhIREVFAY99LRN7QqaJRu3btQlxcHMwePy5SUlJQVVWF/Pz8ettu3rwZ6enpuPXWW9Xfgw46CPPnz/dBrYnqy/ntLYS8cwLCaoPsJaZ47DzxffQ54yGvBdmFZrOh4NVXUfX2O6ociExBwcg8+grEzFiO3ZNnYV3KKdgSOkLdV1BuxSd/78b7H76LLV88jZqqCl9Xl4iIqFWk385/5RW1BGofTkRE1Jmw7yWigB/RLgH1oKCgeutctyWtjKeysjLMnj1b5XD/8ssv1Qj3U089Fb///rsKujdUXV2tFpeSkhL11+FwqIW8R9pXrkjoCu2869v/IuOP+923d4YPQ9A5byE9Mc3rz9+h1yPqvPNQUVmhyoHd3jokjzpVLZfklePnDbnYmF0GaBrGb38J3SpWo+SfF7Fj8BVIO+oKmIND270GXel97Wts645t6/bEvtd3+O/Gv9pa+u3oaee5y2Df4rW2pla+J9nWHYZ9b+Dgvxv/amv2vR3X1sS27my8+X7uVIF2i8Wigu2eXMHx4ODgeuuNRiOio6NV/na9Xo8hQ4bgxx9/xKuvvopXXnml0bEl/cwDDzzQaH1ubi5qamra/blQ/TdwcXGx+nCW1ypQLd1RghV7M3GjPgJhjhKsjp+CsOMfgUNnRk5OTofUQbvoQtQUFyOvqCig29qThNCP6xWCPfEGbFizVAXZRYQ1DxHLH0bpyhexuc8FiBxzPozm+p8jB6KrvK87A7Z1x5H3dHti3+s7/Hfjh2198cXqj/Th5OW2phaxrTsO+97AwX83ftjW7Hs7rq2JbR3AfW+nDbSnpqaqFDGSp91Um2Jj9+7dKgAfGxtbb9ukpKRG/9D79euH1audQbaG7rjjDtx00031RrRLypn4+HhERUV57TmR84NZ8utLWwfiB7O8D39Yl4PFO6uAiB54t9vDGBOahQEn3QJ9B+cKD/S23peEBGBY73TkDEpH9Y//h/TsH9X6cFsBBq79D0o3voW8g65D2hGXwWA68Ny3XbmtOxrbuuN4pm5rD+x7fYf/btjWgYjva7Z1IGLfGzj4GcW2DkR8X7OtA5HZi/MhdapA+7Bhw9RI9cWLF2PChAlq3aJFizBq1Kh6E6EKWffII4/UC8qvW7cOPXr0aPLYkoKmYVoaIQEyBsm8T16/QGxrraIAX/xbgiXbSuVJqnU9Dp6Mwf0TGr1nO0qgtnVrJfQZDfT5GPmb/0L5d48iI/sHtT7cmovw3+9B4YpXUTRmBrodOu2AT4R09bbuSGzrjtHe72X2vb7Ffzds60DE9zXbOtCw7w0s/IxiWwcivq/Z1oFG78UYTqeKDoWEhGD69Om4+eabsWbNGixcuBDPPfccrr76anV/Xl6eO5XM2WefDbvdjiuvvBKbNm3Ca6+9hm+++QaXX365j58FdRVa6V6UvngU0n+ZAZ3mzO90wpBkTByQ6JMgu6OiAusHDETh4UeoclcX2/MgZFz5EbLPXYAdCUe410dX7YR92Rz898eNWLunRF2RQERE5EvSb6/r118t7MOJiIjY9xKRf+pUgXbx9NNPY8CAARg9erQKpt95550444wz1H2SomHu3LmqLPnZFyxYgM2bN2PgwIEqD+xbb72FQYMG+fgZUFcJspe9ciwiyrZgeNF3mJjzBs4YmYaxveJ8XTVqILH3Qci46lPsnfoF9kSNVOu+S7oMOWU1eOv37Xjlly3Ynl/OdiMiIiIiIiIiosBIHSPCwsLw5ptvqqWhhiNPJSAvo96JOjrIXv7qsQgv26JuF5qSkHTYRRiQEe3TF0IXHIyevy5SV35ImepLGngoMOAH7Fz3B/S58UC+c9T/9vwKLJz/NiaWf4GQYx9ETI/hbDoiIupQ0m/3Xvybu0xERETse4nI/3S6QDtRp1aarYLsYaV1QfZdJ83D4IFDfF0zla7GGBMDvc3ms/zwnZ5Oh/QBY3C5pmFtVgm+Xb0XeaVVmLT3FSRXbYbjzSOxLeMUxJ5wP8LjM3xdWyIi6iJcfTgRERGx7yUi/9XpUscQdVqVhaiYNaVekH3HifMweJDvg+y0/wGNgSmRuGFiH5zZW0Ooo1St18OB7js+QtBLB2PzB3ehqryYTUtERERERERERC1ioJ2oNWrKUTH7NIQUbVA3C02J2HL8+xg6uPME2bWaGuS//Aoq33pblaller0OQ4cMQ9CNf2Pz0Bmo0oeq9WZHFXqufR7WZ0dg87cvwmGzsTmJiMhrpN/Oe/lltbAPJyIi8j72vUTkDQy0E7XEVoPKt89FSPZf6mapMRorj5iDg4YN61Rtp9lsyHvuOVTNmqXK1HpBwWHoecrdsF37N7b0OA92GNT6cGseei65AwX/GYXdS+ezSYmIyCuk38599r9qYR9ORETkfex7icgbmKOdqAWFBTmw52yBTE0mI55/G/MqJh0ypvO1m9GIyNNPQ2VllSrT/guLTkTY+S+gcMfVKPvyLqRnL1Dr4yo2499fXsAc62AcOzgJCeEWNi8REbUfoxFRU093l4mIiMjL2PcSkRfwmzzRPlTb7JizsgKlPV7Audvvxpq+1+LYIyZ2yslG9WYzkh58EDk5OapMbRedMQDRV36CvSt/hO77uxFf+i++Tr4SOXtLsSG7FKN7xOLIvnFsYiIiahfSbyc/9BBbk4iIqIOw7yUib2CgnagZmqbhw792IbukGjBG4dMhr+DKI3rDaGDGpa4iachR0AYdjg1//4LqnGSg0gqHBizZnI/SVV9hmGkXYk++Dfogud6BiIiIiIiIiIi6KgbaiZqStRI/F0Rh9e4SdTPIqMe0sZkINjtzd1PXodMb0PegI3CTzYFFG3Pxy4Zc2K3VmLzjP4it2Y2CrR+j7IiHkD7qxE55pQMREREREREREXkfh+YSNZTzL+xvHI/uX56DEFsRJHZ61qh0xIcHdeq2clRUYMNBI1E4+VhVpvZlNupxVP9E3HRMX0wKWY/omj1qfUzVdmR8fT52PT8FedtWs9mJiGi/Sb/97/ARamEfTkRE5H3se4nIGxhoJ/JUlgv7O1NhqClB94pVmJg9C0cPSES/pAi/aCetshKoqvJ1NQJaZLAJhx53DnLO+Aq7Qwe516fn/4ro2Ydh81vXobKkwKd1JCIi/yN9uOrHiYiIiH0vEfklpo4hcrFboc2bDkPxDnVzt6UPtg+/DWf2ifeLNtJZLOjx/XfIy8tXZfKuhH5jkB01D9s3fI/oxQ8joiYHBtjRc/MclD/3GbYePAPdjroCeiM/ZomIaN+k3+75ww/uMhEREXkX+14i8gaOaCdy+f5e6Lb/poolxlh8PuBpnDy6j9/k3dbp9TClpsKQnKTK1DFtnn7YNATf9De2DroWNTpneqFQWxEyl9yFf964Dtvzy/lSEBFRi/2JOS1VLezDiYiIvI99LxF5A6NxROKf94DfX1RFm86IuZmP4IRDR8Ji4uSn1DKTJQyZpz+Mqsv/wPbkY2vfRyZ8H3YiXv55Cz5YuhPFlVY2JRERERERERFRgGJOA6I9K6DNvx6uceufp9yIQaOPRlp0iF+1jWa1ouCdd1BVWgbtskuBoM49eWsgikjKRMTl72HPygX4d+WfKAxKVev/3lmEtVklOD6xAMOHDofREubrqhIRUSfrwwvffVeVo885BzqTyddVIiIiCmjse4nIGxhop66tPB/a++dCZ3NOIPpnzAkoH3QexvaMhT9+Uch9/P+c5QsvYKDdh1KGHImkQUcgdFsBvlubjYoaOxw1Fej9w0Wo+FGHksPuQ+rYs5kegIiI3H149mOPq3LU1KkMtBMREXkZ+14i8gYG2qlrW/Q0dMW7VHFHyEAs7DED1x6U5jd52esxGBB+/PGorq5SZfItvV6H0T1iMTgtEt+vzUbYkicRZc1R90X8cBX2/jUL5hOeQkyPEXypiIi6OoMBEVOmuMtERETEvpeI/A8D7dSl5Y25DZs3ZWFA0UK83e0RnDayB0LM/vnPQh8UhJQnn0BOTo4qU+cg76eThqUiJ/Ii7Jy/HukFi9X6pMK/4HjzKGzvcSYST3oIlsh4X1eViIh8RPrt1KeeZPsTERGx7yUiP8bJUKnLcjg0zFuRh09TbsKzfd7CgD590Dcp3NfVogCVkDkEadd8ia3HvI7CoDS1Tg8Hum2ZC+25Edjx7XPQ7JwwlYiIiIiIiIjIHzHQTl3Wr5vysKOgQpWDIxNw7OAkX1eJApxOr0fm2NMQcuNSbBx8C6r1wWp9sL0EGUvuQcF/xiB7wzJfV5OIiIiIiIiIiPYTA+3UtTgcwDd3oGjbSvywLlutknTsU0emIcjo3zlRHRUV2DRuPIpOOlmVqfMKsoSg92n3oOLyP7ElpTYnL4Dwil144+9ifPTXLpRV23xaRyIi6jjSb284ZKxa2IcTERGx7yUi/8RAO3UtS/4H/P4iQt88Gv3zv1erDukRi26xoQgE9sJCaMXFvq4GtVJ0YgZ6XPYOdp36GbJD+2NhwjQUmxKwbHshnv5uPX7blAe73cH2JCLqAqQPl4WIiIjY9xKRf/LPWR+J2mLXX8CPD6qiyVGFSkM4okNMOGZgYkC0p85iQffPPkNBQb4qk/9IG3I47AN/Q/jmHAT9W4BqmwNVVge+/3sTen56E7Sx1yL5kLOcl18QEVHAkX67x/zP3WUiIiJi30tE/oeBduoaqsuAjy8BHM50HAvjp2Fj+GhcODzV71PGeOb/DurdC4acCFUm/2IwGDC2TzIGZ8ThuzXZalT7kdmzkVS2BvjuCuxdNgshJzyJiMzhvq4qERF5pQ/vzXYlIiLqIOx7icgbGI2jruHbO4GCLaq4M7g/fki6GMMzotAnMdzXNSOqJ9xiwmkHpeGqCZlId+xyr08qWIqwOUdi19tXwlqax1YjIiIiIiIiIupEGGinwPfvl8DyOapYrQ/G+xn3whIUhOMHJyOQaFYriubNQ/UXX6gy+bf02DB0v+ZzbDrqNRQGpap1ejiQtuld2P47HHu+fx6anROmEhEFAum3Cz/4QC3sw4mIiNj3EpF/YqCdAltpNvD5te6bXyRfh/ygdEwamITQoMDKnCQ/zLPvux8VTz3NH+kBdDljr0PPQPANS/HvwBvViSIRbCtBym93oeCZMShYs8DX1SQionbow/fee59aGGgnIiLyPva9ROQNgRVpJPKkacBnVwEV+ermmohDsSxmCtKig3Fw9+jAayuDAWFHHonq6mpVpsBhCQ5Fv6n3I2/3dOycfyd67f1KrY8t24jKj6bhm6qFmDC4J4LNfN2JiPy2Dz/qKHeZiIiI2PcSkf9hoJ0C166lwKYfVbHUGIuP026DTq/DiUNToNPpEGj0QUFIff5/yMnJUWUKPHGpmYi9/F1s+ftHhPx4J5LK1+PHxIvw2/Zq/LV3vbpS46Bu0QH5/iYiCmTSb6e/8Lyvq0FERNRlsO8lIm9g6hgKXOmjgAu+RHlwCual34kKY5QayZ4eE+LrmhG1mQTRe4yYiNgbfsOaMU9jWeJpan1ZtR0fLd+NWd//hbw/5zmv6CAiIiIiIiIiog7BEe0U0LaHD8NrPd+GXW9GiNmAYwYk+bpKRO3CZDJh4ORLkFpRg69W7cWq3cVq/cD1zyMu/xPk/PkaQk96CqHpQ9jiRERERERERERexhHtFLA0TcMXK7NUkF0cMyAx4CZA9eSorMTmiUej+MyzVJm6hqgQM84ZnYFLDs1EH1MORud/ptYn5P2B4FmHY8+718BW5pyngIiIOifptzcdeZRa2IcTERGx7yUi/8RAOwWW/M3A0pmAw6FG+O4qdAacEyOCcHD3GAQ0TYNtzx44srOZNqQL6hkfhunHHYENE15AoTlFrdPDjpQNb8H67HDkLHgRcNh9XU0iImqKpsG6Z49amPqLiIioA7DvJSIvCNzhvdT1OBzAZ9cAOxZDW/URfo27U8b7qruOG5wMvT6wJ4jUBQUh4/33UFhYqMrU9egNevQ74hyUjToRa754Ar3XvwKzowrBtmIE/3IHCv6eDePxTyCi3+G+rioREXmQfrv7vA/cZSIiIvIu9r1E5A0c0U6BY9ksFWQXVQU7sbfG+UO1V0IYeieEIdDpDAYEDx4MY79+qkxdV1hoGAae+SDyL/gNGxKOda+PKV2PiPdOwo63rkSNzeHTOhIRUeM+XBb24URERN7HvpeIvIGBdgoMRTuAH+5335yXchus+mDodDKaPQk6KRB1MckZvdD7yrnYcNw87A3p417/V1UK/vP9BqzcVaTmMiAiIiIiIiIiogPDQDv5PwkUzr8eqClTN7d3PwPrLMNUeURGNJIjg9EVaDYbSuZ/gervv1dlIiEnmfqMOgZR1/+GVSMexKbwg7E05gQUV1ox98+dmLloK/bmFzEnMBGRD0m/XTx/vlrYhxMREbHvJSL/xBzt5P9WvANsXqCK9vAUvBl+kSqbDDocPSARXYVWU4Os225zlk85BTCbfV0l6kQsQWYMPvF65JZegd4r92BDtvPE1Ja8cux++y4YjAUIO/lpBKcO8nVViYi6HOnD98y4VZXDjzoKOiO/ohMREbHvJSJ/w2/x5N9KsoBvZdJTp8X97kZFtTMf+/hecYgMNqHL0OsRcsghqKmpUWWipsSHB+GCsd3x795SfLFyD0JzV+Cgwq/UfY7XDsOefuch8YQHYAiNZgMSEXUUvR6hYw9xl4mIiIh9LxH5Hwbayb9Txnx5M1BVrG5WDpiKr2uGqHKwyYDD+sSjK9FbLEifNRM5OTmqTLSvdDL9kyPURMGrf9uAwh3JiK7Jgh52pPw7BxWbPkPx2NuRNOFS6AzsJoiIOqIPz3j9dTY0ERFRB2HfS0TewCEz5L/WfQ6s/9JZDo3HVynXqdi7OKxPHCwmg0+rR9TZmQx6DD/sROiv+ROr+l6LGp3zBE2IrQjJv9yO/GfHoWDdz76uJhERERERERFRp8dAO/mvHocDI6arYuHhj2JZjk6Vw4IMOKRnrI8rR+Q/IiMiMPjsh5F3wa/YmDDJvT6u9F/EvH8ids08B+W5231aRyIiIiIiIiKizqzTBdorKytx0UUXITw8HImJiXj00Udb3KeoqAjJycmYPXt2h9SROglLJHDic8AVv2K+bZR79eF9ExBk7Hqj2R2Vldh6wokovuACVSbaXyndeqPXle9j05R5yA7p416ftutLfPflPCzamAub3cGGJSJqZ9Jvb54yRS3sw4mIiLyPfS8ReUOnS747Y8YMrFy5En/++Sd2796NqVOnIjMzE2efffY+99m7d2+H1pM6j53mnliXtVmVZfLTUZkx6JI0DTWbne3gzqFD1Ib87b1GHgPr0COw/rsXkb78aeQGpWNp+NHQVu3Fn1sLMHlQEgYkR6htiYionfrwTezDiYiIOgz7XiIK9EC7jGafNWsWvv76a/Tv318tN954I1566aVmA+2//PKLWmREO3UBlUVATRkQmeZe9f3abHf5yH4JKu90V6QLCkL67DdQWFikykQHwmQyoe/x16Nk3DlYt3IjkK8HNCCvrAZv/74Dp5S/j8xDTkZ8r5FsaCKiAyT9dsacOe4yEREReRf7XiLyhk4VkVyxYgWsVivGjh3rXjdu3DgsXboUWhMjdKurq3H55ZfjhRdegNls7uDakk98dxfwwhhg2euAw4GdBRXYmFOm7ooJNeGgbtFd9oXRGQwIGTUKpuHDVJmoPURExWPyYWNxzRG90CMuVK3LLFuOUZv/h9i3J2LbrOkozWH+diKiAyH9dujoUWphH05EROR97HuJKOBHtO/atQtxcXH1guYpKSmoqqpCfn6+us+T5G8fMWIEJk6c2OKxJSgvi0tJSYn663A41ELeI+0rJ0oOuJ03L4T+77dVUfvuXmi9J+PHddXuNCmH9Y6DDvI4XTdtSru1NbGtG0iKCMJF47phbVYJLB/dodbpoaH7zk9R89LX2NLvIiQfdyuCQqPa/d3D93XHae/PDva9vsN/N2zrQMT3Nds6ELHvDRz8jGJbByK+r9nWgcjhxZhZpwq0S0A9qMHlsq7bklbG07p16/Dyyy+rfO6t8dhjj+GBBx5otD43Nxc1NTUHVG9q+Q1cXFysAsB6fdsuotBZyxH32bXu2yVjbsG2HBtWbM1Rt8MtBqRarMjJcd7uijSbHTWLf0NFRQUcRx0Fg8nk6yoFtPZ4X/ujeCNgO/UFrP5tFnpsmIUQRxnMWjV6rHsJZRvmYu2AKxAz5jzoDe33/uuqbe0L0s7tiX2v7/DfjX+1tfTh1iVLVNl0yCHQdcFJ3VuD72u2dSBi3xs4+BnlX23Nvrfj2prY1oHe93rSaU3lZPGRefPm4ZprrkF2dna9gPqAAQNUQNw1ol2qfNhhh+Hcc8/FFVdcodZ1794d999/Py644IJWj6pLT09XI+Wjotp/BCbV/2CW1y8+Pr7tgfavb4Vu6WuqrHU/FNq0T/Hun7uwZo/zyoQThiZjTI/YLt3sjooKbBx5sCr3/PMPGMPCfF2lgNYe72t/V16Ug5wvHka3rXNh1Gzu9fnBmSg79B6kjToJunZoG7Z1xykqKkJsbKz64hEREXHAx2Pf6zv8d+Nfbe3Zh/dethT6kJB2rmVg4PuabR2I2PcGDn5G+Vdbs+/tuLYmtnWg972ddkR7amqqCnxLnnaZiE/s3r0bFotFNYDLjh078Ouvv6rR7Lfffrs7cH7VVVfhww8/xBdffNHo2DIyvuFoeSEfFPyw8D6dTtf2tt6+GKgNssMYDN2JzyG33Ia1e0vlwAi3GHFwZixfR6MRluHD1L8fvdHI9ujs7+sAEB6ThPDzn0fBrmtQ+sU96Lb3O7U+tnIrYr+7CJ/nvY7h449FesyBB4y6elt3lPZuX/a9vsV/N37U1kYjgocPV0X24V5ua2o1tnXHYN8bWPjvxo/amn1vx7U1sa07GW++lztVoH3YsGEwGo1YvHgxJkyYoNYtWrQIo0aNUv+wPQPyW7durbfv+PHjccMNN+C8887r8HqTF1krgc+uqbt91L1ATA/8tHSnKzU7Du0dB5OBH/h6iwXd3nlHpc+RMlFHiUnrh5gr5mHv6p+B7+5GUslKbAo7CEusvbHkp80YmhaJYwYmISaUk1YTETVF+u3uc99l4xAREXUQ9r1EFPCB9pCQEEyfPh0333wz5syZowKGzz33HF555RV1f15eHsLDw9UIOUkV40kC9JJaJikpyUe1J6/46TGgYLOznDYKGH058sqqsWJXkVoVYjZgVGYMG5+oE0gaNAHagJ+xY/H7WJwbBmjOE6T/7CrGmt3FOEv7CplHXoCQqERfV5WIiIiIiIiIqF11umHATz/9tMrJPnr0aJx99tm48847ccYZZ6j7JCfU3LlzfV1F6ii7/wIW/89ZNpiBk54H9Ab8ujHPPZp9fK84BHHCMKJOQ3KyZ4w/G+eeNEXNnRBqdk7o17N4MQaufBT654Zj68f3o6ai1NdVJSIiIiIiIiIKzBHtIiwsDG+++aZaGtrXvK3btm3zcs2ow1WXAaEJQNleYMJtQHxflFXbsHxHobo7yKjv8hOgenJUVWH7eefBarUhbu67nEiNfMqg12FszziMyIjGzxtyMeDL2Wq9xVGOzJXPoGzdHOwefj0yjr4SBlPj+TOIiLpeHz5Nlbu9/RZTwBEREbHvJSI/1OlGtBO59ZgAXP0HMOF2YNz1atXvm/NhtTtPuBzcPQbBtaNlSU0HjqrVa2Bfv16ViToDi8mASQOTEDn9PWzJOB12OP/NhlkLkPnnfSh9aji2/zQHmsPu66oSEfm4D1+tFvbhRERE7HuJyD91uhHtRPUERwFH3KGKNTYHlmzJV2W9DhjXK5aN5UFnNiP1pRdRXFSkykSdSURiN0RcNAv5225E2df3o1v292p9VPVuRP10HXKXPo+qCfcgbeQUlX6GiKgrkX477eWX3GUiIiJi30tE/oeBdupcZCS2TudcGpCUMRU1zlGvQ9OiEBXCH6KedEYjwiZMQEVOjioTdUax3Qch9soPkb1uMWzf3YfUwj/V+vjyDcBX0/Bh3ns4ePRh6BYb6uuqEhF1GOm3ww8/nC1ORETEvpeI/BiHDVLnsvg5YO7ZQMmeeqsdDk1NguoyvnecDypHRO0lsf9YpFz3HXaf8C5ywvqpdWsjxuOvqjS8/PMWvLlkG7KKK9ngREREREREROQXOOyVOo+8jcDCRwF7NfDS78B1fwPB0equtVklyC+vUeVeCWFIiQr2cWU7H81uR/mSJbAWFUGbNAlg+g3q5HQ6HVIPOh7a8MnY8dt7+DM/BqidXmBdVinWZxVjWtlsJEy4GAjiyTUiCvA+/PffVTl0zBjoDJyDhoiIiH0vEfkbBtqp86SM+fxaZ5BdDDvXHWTXNA2LPEazH8rR7E3Sqqux65JLnWW5/Nxk6oAXjujA6fQGZBx6Ls53aCpF1A/rclBcacWgogXot+N1OLbMxubEyQg64X5Ep/VlkxNRQPbhOy++RJX7Lv8LupAQX1eJiIgooLHvJSJvYKCdOoelM4EdS5zl6EzgiLvcd+0oqFCLSI60oHdCmK9q2bnp9Qjq2xc2m42j2ckv6fU6jOweg6HpUfhjSwEyPvvEuR4O9M7+CvaZ32Jrt1MRe+zdiEjq7uvqEhG1bx/ez5lGi1ekERERdQD2vUTkBQy0k+8Vbgd+uL/u9on/A8x1I7l+3VQ/N7ukm6DG9BYLun/yMXJyclSZyF+ZDHr1b73qis+x+dtnkbJ2JoLtpTDAjszt82B75VNsyTwTccfejoj4dF9Xl4jogEm/3eNT58lFIiIi8j72vUTkDZwMlXxL04D51wHWcuftkRcBmYe67y6qqMHaPSWqHG4xYkhqpK9qSkQdzBIWhZ6n3Q/HtSuwpsfFqNI7T8AZNSt6bHkblhcPwpZ3bkBZwV6+NkRERERERETkUwy0k2/9/Raw5SdnOSINmPhAvbt/31IAh+Ysj86MgdHAtyxRVxMcEYPYY26F7dq/saXvpajRO6/YMGvV6L5xNt747k98vSoL5dU2X1eViIiIiIiIiLooRi3Jd0r2AN/eXXf7hGcBS4T7ptXuwNJtBaos8fVRmTG+qKXfcFRVYcf0C1B6/Q2qTBRoQiIT0OPsp2C95m9s7TUdVp0ZK6OOwh5zJn7ZmIcnv12Pb1ZnoayS738i8i/Sb2+fdr5a2IcTERGx7yUi/8Qc7eQ7qz8Cqoud5aFnA72Prnf3PzuLUFFjV+UhqVEIt5h8UUv/4XCgculSd5koUIXGpCDzvOdQlnsL8jbnwrhXB5tDQ7XNgUXr92LEF5ORm3ooYifdioiEDF9Xl4ioZQ4HKtiHExERdRz2vUTkBQy0k++MvRaI6gb88gQw6dF6d2mahsWb8923D+kZ64MK+hed2YyU//wHxSXFqkwU6MLiMzAxPgMHV1rx0/ocLNtWiKF53yKhahuweRtsL72Hrd1PR/QxtyIquYevq0tE1Czpt1OffcZdJiIiIu9i30tE3sBAO/nWgBOB/icAOl291dvyK5BV7Ez/kBYdjPQY5ySI1Dyd0YjwyZNQmZOjykRdRWSwCScNS8XhfROw+6vPULPHArOjSk2amrl1LmyvzsOW9FMQPek2RKf29nV1iYgakX47YvJktgwREVEHYd9LRN7AHO3kew2C7GLx5jx3eSxHsxNRKwPuA067C9ZrVtROmhqs1hs1G3rsmIeI10Zj66zpyN++mu1JRERERERERO2KgXbqWGs/d+Zm17RmNymusGLtnhJVDrcYMTg1sgMr6L80ux0Vy5fDtmqVKhN1VaExyWrSVNt1/2BL/ytRrXdeEWOAHZk7P0X0G+Pxx0fPIKu40tdVJSKq14fLwj6ciIjI+9j3EpE3ML8EdZyyHGD+dUBlIbDmE+C0WYAxqNFmv2/Nh6M2Dj+qewyMBp4Pag2tuho7z5umyknLlgImTh5LXVtIVCJ6nPk4KotnYOt3zyD53zdgsZdBgx4/WQei6MdN6JMYhsP6xKNHXCh0TVxdQ0TUUX349nPOVeW+y/+CLoQp84iIiNj3EpG/YaCdOoamQffVDGeQXeiNTQbZ7Q4Nf213bqPXAaN6xPAVai2dDqaMDNhlNDsDhkRuwZGxyJz6MKpKb8aWH19ETtZOFJmT1X0bssvUcoT1Z/RLjUPaIVOhNxjYekTU8X14twx3mYiIiNj3EpH/YaCdOoRl89fQ/TvfeSMkFjjuqSa3W5dVgtIqmyoPSIlAhIWjsltLHxyMHt98jZycHFUmogafQ+HR6HHyXUizOaBtL8SijbkorLDC6KjGmI3/QcS6fBT89iiKR1yN9AnTYTRb2IRE1CGk3+717bdsbSIiog7CvpeIvIGBdvK+8jxE/Ppg3W0JsofGNbnpH1sL3OXRmRzNTkTtz2zU45CeseozZtXuYuz57R1E2PLVfTGV2xHz260o+fNJ5Ay+FGlHXQFLKOeJICIiIiIiIqJ9Y/Jr8n7KmK9nQF9VmzKm/4nAwFOa3DSvrBqbcspUOSbUhJ7xYXx1iMhr9HodhqZHYfIZV2LX8W8hK+og930R1lz0Wv4o8J8B2PLuzSjL2c5XgoiIiIiIiIiaxUA7edfqj6Bb+6kqasExwPFPN5t7dNm2utHsozJjOTHhfnJUV2PXFVei9PbbVZmIWken1yPt4BORfMMCZE+dj50Jh7vvk8lTe2yYieAXh2PTG5cjq7iSzUpE7U767R2XX64W9uFERETex76XiLyBqWPIe4p3A1/e5L6pHfckdGEJTW5qszvck6Aa9MBB3aL5yuwvux3lv/ziLhPR/ksceBgw8DAUbP0HJQueQdquL2DUrDDAjj3lOnz94yb0jA/FuF5x6JcUzhOCRNR+ffjP7MOJiIg6DPteIvICBtrJe76/F6gqVsXKXlMQNPDUZjddm1WCsmpncHhgSiTCgvjW3F86kwlJjzyMkpJSVSaitovJHIqYi2ejJHcXcha8gMSN72Nx3Gnqvs255WpJD6rAZP0fSD38QgQFh7O5iajNpN9OfvRRd5mIiIi8i30vEXkDo5nkPZMfB6yV0PYsR8n4exC/j03/2OKZNoaToLb1i0LkKaegOieHP9KJ2klEfBoiznwM1TX3Y8LOMizenIe8shp1X68dH6BH9kxULHsSm3ufg/gjr0VEQhrbnoja1IdHndr0HDZERETU/tj3EpE3MNBO3hMWD5z1DrSSPdAqDc1ulltajS155aocH2ZGj7hQvipE1KkEmYNwSM8gjOkRg3/3lmLxhiyMWfuJui/EXoKe/74M2/qZ2JY8CcHjrkTigPHNzkdBRERERERERIGHk6GSd0mgKTx5n5v8uZWToLYHzW5H1bp1sG3cpMpE1P50Oh36J0fg4gl9UX3GXGxPnQI7nCcSjZoN3fd8icR5U5Dzn7HY/uNrsFVX8GUgolb34bKwDyciIvI+9r1E5A0MtFP7Wvs5UJ7X6s1lEtTlO5yToBr1OozoFsVXpI206mpsP+10lF56qSoTkXfF9xmNbpe+g8qrV2Brv0tRYYhw35dQuhbdFt2Cmif74fdfvkVRhTPdDBFRU6Tf3nrKqWphH05EROR97HuJyBuYOobaz54VwIcXAsHRwIn/A/oe2+IuMglqRY1z9PWg1AiEmPmWbDOdDsaEBNgdDqasIOpAYfEZCDvrKVirHsDWRW8jdMUsJJSvd/6zdNjxdXYUPv92vRoJf0iPGPSMD1Mj44mIGvbhrjIRERF5GfteIvICRjWpfdSUAx9dAjhsQHkusPuvVgXa/9ruHM0uDurGSVAPhD44GD1/WoicnBxVJqKOZbKEIvPoy6EddSn2rFmEmiUvY4c1CjWGEEAD1u4pUcsJhW8iKSkVKYdfCEtoJF8mIlL9du9ffmZLEBERdRD2vUTkDQy0U/v4+jYgf6OznDwUOOzWFncprrBiY06ZKkeHmNAznpOgEpH/0+n1SBk8ARg8ATGVNbBuK8Sf2wpQUmlDiK0YB++aDdPOGlT99QS2ZJyIiPGXIq7XQb6uNhEREREREREdAAba6cCt+QT4+y1n2RQKnPY6YDS3uNvynYXQNGf5oG7RTKVARAEnItiMo/on4vC+CViXVYLcX+fApDnztVsc5eixbS6wbS6yIwahcsj5SBl3LszBYb6uNhERERERERHtJ06GSgemaAfw+fV1t497Eojr1eJumqbhr22F7lSkIzKi+UocIEd1NXbfcCPK7rtflYmo8zDodRiUGokjzrwO+ef9iG3dTkeN3uK+P7FkNbr/eiscT/XFljlXInfzcp/Wl4g6lvTbu66/QS3sw4mIiNj3EpF/YqCd2s5uc+Zlry523h50GjDsnFbtui2/AvnlzlGdPeJCER3a8gh4aun1sKPsu+9g/flnVSaizim210h0v3AWHDf+i82jHkROaG/3fRZ7GXpsfReODy7ECws2Yum2AlTb+O+ZKODZ7Sj99lu1sA8nIiJi30tE/ompY6jtfnkC2PmHsxyVAUx5xjk8vRU8J0Ed2Z2ToLYHncmEhLvvQmlpqSoTUedmCY9Gz+Ouhzb5WmSt+w1Vv89C6q6vYNaq8WfsSdhVVIVdy3fjy5VZGJoeibGhWUjsPbLVn7NE5D+k30685253mYiIiNj3EpH/YaCd2iZ/M/DLk86yzgCcNguwRLZq1yqrHat2FamyxaTHwJQIvgrtQH6YR59zDqw5OfyRTuRnk6cmDzwUGHgoqkoLsHXRHGTpDwMqnPdX2xzYufYPJG68CPkhPVDS70wkjD8foTEpvq46EbVjHx5z7rlsTyIiog7CvpeIvIGBdmqb2J7AaTOBz68Dxt0ApI9q9a6rdxejxu6cBXVoWhRMBmYwIiISlvAYZB53Iy7VNOwuqlSpY/7ZWYyRBV86P3ortiB2+WOwL38COxMOgzbsPKQefCIMJqbfIiIiIiIiIvIlBtqp7SQne8oIZ9qY/bDMI23MQd04CWp70RwO1GzbDntBPrS4OEDPExhE/kqn0yEtOkQtxw5Kxu5fD0fWis1ILl6h7jfAjvSchcB3C1G2IAbZmScj4pDpiO8xzNdVJ6I29uHWHTtU2ZSRoa50ISIiIu9h30tE3sBAOx2YmMz92jy3tBrb8535EBIjgpAWHcxXoJ1oVVXYetxxzrZdthQIC2PbEgUAi8mAnkecDxxxPvK3r0Hx4jcQv+UThFvz1P1htgKEbXwd2Pg6lqVfiJrD71ZXC4UGsYsn8qc+fPPkY1W57/K/oAsJ8XWViIiIAhr7XiLyBv4Kp9bb9ANQngcMPavNrbZ8R/3R7DJqk9qPPjxcnZknosAU220gYrs9BYftMez660s4/n4bKdkLYdRs6v5VxoHY8E+WmkC1T2I4hqcEo19iKMzBPPFG5A99OBEREbHvJSL/1emuS62srMRFF12E8PBwJCYm4tFHH21223nz5mHQoEEICwvD6NGjsWTJkg6ta5dSvAv46FLgk8uB+dcDtpr9PoSmaVix0zkJqsTXh6VHeaGiXZc+JAS9//gdUV9+ocpEFLj0RhPSRp+MjCs+RM31a7F15D3YET0Gm8JGqvsdGvDv3lKs+/FtaE/2xraZ07Drr6/gsDkD8kTUuUi/3Xfpn2phH05ERMS+l4j8U6cb0T5jxgysXLkSf/75J3bv3o2pU6ciMzMTZ599dr3tFi1ahPPPPx8zZ87EuHHj8Prrr2Py5MlYv349kpKSfFb/gGSrBuZdAFQWOG+XZkuUZ78PszW/AkUVVlXunRCGcIupvWtKRNTlhEQlInPKLQBuwXUlVVi+vRArdhWhpNKG4YXfIMhRge67Pgd2fY7Sb+KR0/0EhI08Ewm9RzEPNBEREREREZEvA+3btm3D999/j6VLl2Lv3r1qXUJCAkaMGIFjjz1WBcbbOpp91qxZ+Prrr9G/f3+13HjjjXjppZcaBdpnz56NM844A+eee666/eCDD+L999/HF198gUsuuaRNj09N0DTgqxnArqXO2zLx6SkvtWmizRU7nKPZxfAMToJKRNTeEiMsOHZwMiYPSsLW3DKgOBNVlWthsZep+8OtuQivzedeaElHUc8TETnyTMRkDuWLQURERERERHQA9ita+tNPP2HSpEno168fXnvtNdjtdhUMl/QtRqMR7733HoYOHaq2+eWXX/a7MitWrIDVasXYsWPd62S0ugT0Je2Ip2uvvRa33npro2OUl5fv9+PSPiybBSyf4ywbLcAZbwLB+x8kt9k1rN5drMpBRj0GJEew2duZo6YGWXfeifLHHldlIuq6ZP6LHgnh6HHhqzDcuhE7jnoJOxMOhx0G9zbRVTuRueYFxMw5DAvffgw/rc9BYTk/O4h8QfrtPbffoRb24UREROx7iSjAR7Sfc8452LFjBy699FJ8/PHHCA0NbXK7qqoqfP7557j77ruRkZGBt99+u9WV2bVrF+Li4mA2m93rUlJS1DHz8/PVfS7Dhg2rt++PP/6IDRs24Igjjmj141ELti8Gvr6t7vYJzwEpw9vUbJvyKlBtc6jk7ANSImA2drrpAfyfzYaSTz9zlh952Ne1IaJOwhQUgoxDzwEOPQeVRTnI/uMDmNd9iqSiZdDDeRJ7qWE4Ctdko7q81NfVJeqabDYUf/qpKibdew/g8V2YiIiI2PcSUYAF2i+44AIcc8wxqlxRUdHsdhaLRaV0keW7777br8pIQD0oKKjeOtdtSSvTHMnLLicCzjzzTAwZMqTJbaqrq9XiUlJSov46HA61UAPFu6D74HzoHM6J87QxV0MbPFUabL+bStp3VVa5uipBJydJ0iLZ5l6g6fWIu+kmlJWXQTMY2MZeJu9reU/z88P72NbtJygiDhlHXwUcfRVKcnci788PUL17DQrNKSpVWMOrxw4U+17f4b8b/2pr6cPjb77ZXWbf4r22ptZhW3ec9n4/s+/1Hf678a+2Zt/bcW1NbOvOxpvv51YH2l1BdjFgwADMnz8fgwcPbvU+rSFBegm2e3IFx4ODg5vcZ+vWrZg4caIa+S7pbJrz2GOP4YEHHmi0Pjc3FzVMs1GfrQqxn50DU3mu8zVIHYvCIVcBOTloi9IqK9btKYLJbEZ4kBFhWjlycpo/WUNt55hyPKqKi5FbWAh9G/Lo0360tcOB4uJi9aWDbe1dbGtvCULYwdMQdjAwvdKGf7PL8deWurk02gP7Xt/hvxs/bOsTpqg/uUXt++8wkPB9zbYORPL50Z7Y9/oOP6P8sK3Z93ZcWxPbOoD73gOeDFX+oXkjF3pqaqpKESN52k0mk1q3e/duFYCPjY1ttP3GjRtVqpj4+Hg1OWt4eHizx77jjjtw00031RvRnp6ervaNiopq9+fi1wq3QVfjfNNpUd1gOvstJITEtPlwGzflwmQOQrDFgkP6xCMpMbEdK0sN/21KbmZ5X7MT9C62dcdhW3tfAoA+3YDD+sSjPRNPse/1Hf67YVsHIr6v2daByDNtantg3+s7/IxiWwcivq/Z1oHI7MU0jW0KtJ9wwgk47rjj1F/Jw+4Kirvce++9baqM5F2XSVUXL16MCRMmqHWLFi3CqFGjVPDQU15eHo4++mg1kl2C7JGRkfs8tqSgaZiWRkgwkgHJBmJ7AJf8CHx0CXSTH4curC43flv8s6tEpYyR1/CgbjFsby/RHA7Yc3Oh5eVBx0B7h5D3ND9D2NaBpL37Q/a9vsXPKP9pa+nDbbnOKwmN8fHQccSY19qaWo9t3THY9wYW/rvxn7Zm39txbU1s687Gm+/lNgXa165di6FDh6rJUWVp+A+wrYH2kJAQTJ8+HTfffDPmzJmDnJwcPPfcc3jllVfcwXUZtS4/3O+55x6UlZWpPPAyAl7ucx1DFjpA4UnA9Plq8tIDkVtajd2Fzvz6SZEWtZB3aFVV2HLEkaqcsGwpEBbGpiYiIvKTPnzThMNVue/yv6Djd1kiIiL2vUTkd9oUaF+4cCG85emnn8ZVV12F0aNHIywsDHfeeaeaWFVIOow33nhDTcz64YcfqjQzffv2rbf/fffdh/vvv99r9QtY2WuAuD6AwePqhAMMsou/dxS6y8PSmaLH64xGNaEhERER+WEfTkREROx7ichvtfobvYw0l5HqPXv2bNX2mzdvxkMPPYTZs2fvV4UkuP7mm2+qpSGZfMFzElNqJznrgNcnA6kjgDPeBCz7TsPTWvJ6rdhZ5I7ZD01rn+NS0/QhIei78h91JYiUiYiIyD9Iv91/9SpfV4OIiKjLYN9LRN7Q6qQ055xzDo4//niVl/3tt9/Gtm3bGm2zfv16zJo1CxMnTlTbnnvuue1dX2pvJXuAd6YC1SXAlp+AhY+126F3FFSgsMKqyt2iLYgIrp/Ln4iIiIiIiIiIiKhLjWifNGkSVq9ejblz5+LFF19U6Vtkltbo6Gg1C3FhYaH6O3bsWFx22WUq3YtMbEqdWEUB8NYpQPFO5+3kocCRd7fb4f/ZVewuD0wKbbfjEhEREREREREREXUm+xUJl8D5tGnT1CKB9VWrVqk0FSIpKQnDhg1TqV/ID1SXOUey5/7rvB3VDTjnAyCofV4/h0PDql3OtDEGvQ69E5jKxNscNTXIfuwxVFZWwnH//dBbOPEsERGRv/ThOY8/rsoJt98Ovdns6yoREREFNPa9ROQNbRpyfuSRR+7z/gULFrS1PtQRbDXAB9OA3cuct0MTgPM/BcKT2u0htuSVoazarsp9k8JhMbY6SxG1lc2GornvOct3t9+VCURERORlNhsK352rigm33AIw0E5ERMS+l4i6RqB9woQJ9W7bbDY1+ek333yDu+66C/7ManfAqNdBJ7N3BiKHHfjkcmBz7cmQoEhg2idATI92fZh/dtaljRmiJkGtbtfjU2M6oxGxV12F8vJyVSYiIiL/IP123NVXu8tERETEvpeI/E+bvsnfd999Ta5/88038eGHH+Lmm2+Gv9E0Dd+vzcaijXkwG/U4dUQqBqZIgDiAOBzA/OuANR87bxuDgXPeB5IGtevD2OwOrNlTospBRj36JoajqICBdm/Tmc2Iu+ZqOHJyVJmIiIj8g/Tb8dde4+tqEBERdRnse4nIG9o1n8fhhx+OH3/8Ef5oXVYpFq7Phc2hoaLGjvf+3In8sgALDjusQOleZ1lvBM6YA3Q7pN0fZkN2GSqtzrQxA5Ij1IkLIiIiIiIiIiIiokDVphHtO3bsaLRO0lU8++yzyMzMhD/6dVNuvdsScF/wbw6mjkxHwDAGAWe+A3x4ITDkTKDPJK88zMraSVDFkPQAuyqgk1+VYS8pgaO0DFp8vK+rQ0RERPvRhztKS1VZHx4euCkMiYiIOgn2vUTUaQLt3bt3b/QDQD6kkpKS8O6778LfFFdYsTWvQpUjgo2w2jQ1InvV7mJMGZKCYLMBAcNkAc56F/DSD7hqmx3rspxpY0LMBvSKD/PK41BjWmUlNo1xXqGQsGwpEMa2JyIi8pc+fMOo0arcd/lf0IWE+LpKREREAY19LxF1mkD71q1b692WoHtwcDDi/XQU7ea8Mnd5ZLcYFWRfsjkfVruG1XuKcXD3GPjtxKcLHgYOmg5Ed69b78VRUv9mlaLGrqnyoNQIGA16OCQ3PBEREREREREREVGAalOgvVu3bggk2/PL3eXMuFA1gacE2sW/WSX+GWi3W4FPrgBWf+ic/PTCb4CIZK8/bL20MWlRXn88qqMLDkaff1YgJzdXlYmIiMg/SL/db9VK5w1jm76eExEREfteIvIxfpOXnPP5laox9DogPSYYZoMeYUEGlFXbsSmnDFa7AyaDH03oaa0EPpgObPzWebt4F5C1wuuB9soau5oIVURYjMiMDfXq4xEaXVmiM5mgMxqZ25WIiMiPqJSMJpOvq0FERNRlsO8lIm/wo+ixd9gdGnLLqlQ5PjwIQUaD+sDtkxiu1kkalO35zvztfqG6FHhnal2Q3VA7AWrfY73+0GuzitUksmJwWiT0cuaCiIiIiIiIiIiIKMB1+UB7QVk17LUpxBMjLO6G6ZVQN5Hktry61DKdWuleYPbxwLZFztvmMOC8j4C+kzvk4VfsLHaXhzJtTIfTamqQ8+RTqHjpZVUmIiIi/yD9dvYTT6qFfTgRERH7XiLyT10+0J5bVu1ujITwIHe5u0fak20eOdw7rZx1wMyJQNY/ztvB0cD0z4HMQzvk4UurrNic60wbExNqQlo0c4R3NM1mQ+Ebb6D6/fdVmYiIiPyD9NsFr7+uFvbhRERE7HuJyD91+RztuWV1I389R7RHhZgQGWxCcaUVOwsqVIoZQ2dNhbLlZ+D9aUB17YjyyAzg3HlAQr8Oq8Lq3SXQtLpJUFW+M+pQkps9+sILUVFRocpERETkH6TfjrnoIneZiIiI2PcSkf/p8t/kiyus7saICTW7yxIo7h4bgn92Fas87XuKKpEeE4JOSSY6dQXZk4cB53wAhCd2aBVW7ipyl5k2xjd0ZjMSZtyCnJwcVSYiIiL/IP124q0zfF0NIiKiLoN9LxF5Q5cPtBdVyoh2g2qM6JD6wclusaEq0C52FlZ03kD72OuAgq1AaRZw2iwgqC6/fEeQUf/bCyrc6XeSIuuuDCAiIiIiIiIiIiIKdF0+0F5SKbmsDQgy6mEx1U9Z75lnfHdhJToNWzVgrMsnD0nTctxTzrKh41/SNXuK3WljBqdGdvjjk5OmadCsVpXbVcpERETkH1S/7ZpfxWhkCj4iIiL2vUTkh7r8ZKgyGtuVk71hXnEZme1Ky767qJME2veuAl4YBaz/uv56CbD7IMguVu+uTVsjgfY0Btp9RausxIahw1A08WhVJiIiIv8g/fa/g4eohX04ERER+14i8k9dPtAuk5yKqGBTo8YxGfTuCVJzSqtRY3PAp/55D5h5NFC4Dfj4MiB3g2/rA6Ckyopt+c60MfHhQSp1DBEREREREREREVFX0uVTx7hENcjP7pIaFYys4iqVGmVvcRUyYn2Qp726DPhqBvDPu3XrYnsB5lD42prdJfXSxjS8KoA6ji44GL1+X4Lc3DxVJiIiIv8g/XafP/9wl4mIiIh9LxH5Hwbaa0WGNB7RLlKigoHthaq8q6ii4wPte1cD8y4A8jfWrRt+HnDc04DJ95OO1ksbw/zsPiUnOQwREdBXVfGEBxERkR/24URERMS+l4j8FwPttSIsTTeFzyZElWHiS2cC394F2Kud68xhwJRngCFnoDMorbJia365KseHmZEYwbQxRERERERERERE1PUw0F4rLKjpEe2uCVEllbukkOkQpXuBT64Atiz0qMhg4PTZQFwvdBarPdLGDGLaGJ/TamqQ9/IrqCwvh3bjDYDF91c8EBERUSv78FdeVeW4yy+Dztx0SkMiIiJqH+x7icgbuvxkqC6hQYYmG0gmRI0Lc47Uzi2tdk+e6lV6I5C9uu72qMuAi3/oVEF2sWaPR9qYtEif1oUAzWZD/osvomrOHFUmIiIi/yD9dt4LL6iFfTgRERH7XiLyTxzRXiu8mRHtrlHtOaXVsDk05JdVIyHCyyOFQ+OAKc86J0A98X9A74nobMqqbdiS50wbExdmRpK324RaZjQi6uyzUFlZqcpERETkJ4xGRJ9ztrtMRERE7HuJyP/wm3wLI9qFBJFXwjl6e29JVfsG2h0OYPlsoPckIDK1bn3/KUDPIwFzB0++2kprdhczbUwnozebkXjPPcjJyVFlIiIi8g/Sbyfde6+vq0FERNRlsO8lIm9g6hgAwSYDjIbmmyLRI7Dernnas9cArx8DfHEj8M1tje/vpEF2sWq3R9qYVKaNISIiIiIiIiIioq6LI9rVRKjNj2Z3pY5xyS5ph0B7dRnwyxPA4ucBze5ct24+kLUSSB6Czs4zbUxsqBnJHu1DRERERERERERE1NUw0C6Bdsu+myE6xIQgox7VNgf2HsiIdkkTs/J94McHgNKsuvWxvYEpz/hFkF0wbUzn5KiowPpRoyE5feL++B36sDBfV4mIiIj2pw8H0PfPP6AP6bxXNRIREQUC9r1E5A0MtKv87PtuBp1Op9LH7CioQGGFFVVWOyymfY+Cb2THH8A3twN7ltetMwQBh90CjLseMAbBX6zeU+IuD05j2phOxWbzdQ2IiIioLdiHExERdSz2vUTUzhhol0C7ueVmSIoMUoF2V/qYbrGhrW/l7+8Dfnu2/rq+xwHHPAzE9oQ/KZe0MbllqhwTakIK08Z0GjqLBT0WLkB+Xp4qExERkX+QfrvXzz+5y0RERMS+l4j8DwPtMhmqueXR6Z4Tokr6mP0KtGeMAX6rLScMACY9CvQ8Av5ozZ4SOLS6SVBltD91Djq9HqbEROh1OlUmIiIi/+rDiYiIiH0vEfkvBtol0N6KNDBJnoH2fU2IWlkEVBYCMZl16/pMBgadDnQbC4yYDhj8t9lX7S52lwelMm0MERERERERERERkf9GfDt4RHuSR4qUnJLqxhvUlAN/vOJMEZM0BJg+X5K7O++Tv6fPgr/zTBsjE8SmRgX7ukrkQaupQcGcN1FVXgbtiisAXnpORETkP334W2+pcsy0adCZzb6uEhERUUBj30tE3sBAeytHtIeYjYiwGFFSZVMj2jVNc6ZNqSoBls4ElrwAVOQ5N962CNiyEOh5JALJuiymjenMNJsNuU8/7Sxfcomvq0NERET70YfnPPmUKkeffTYD7URERF7GvpeIvIGBdgCWVgTaRUKEBSVVZaiosaOsKBfhK2YBf7wEVNWlU4FODww9G4jthUDDtDGdnNGIiJNPQlVllSoTERGRnzAaEXnyye4yERERse8lIv/Db/JqtHrrAu2Sp33v7h0Yn/ceQl/8DLCW1w+wDzwVOGwGkNAPgaaixoZNOXVpY9KimTams9GbzUh+9FHk5OSoMhEREfkH6bdTHn/M19UgIiLqMtj3EpE3MNDeytQxIjEiCKfuehz9SxfXrdQbgSFnAeNvBOICbxS7y9o9dWljZBJUlTaHiIiIiIiIiIiIiBho3+dkqFptZLk2qJwYYcFX8eeoQLtdb4JhxPnAuOuB6G4B/1byTBszODXSp3UhIiIiIiIiIiIi6kz06GQqKytx0UUXITw8HImJiXj00Ueb3Xb58uUYOXIkLBYLhg8fjj/++GP/H1AHBBkbNIPdBqz+CJh5FLDmE/fqhIggbAsbhvkp1+PNUfOBKf/pEkF2z7QxUUwb02k5KiqwcfQYFB0/RZWJiIjIP0i/vf7gUWphH05ERMS+l4j8U6dLHTNjxgysXLkSf/75J3bv3o2pU6ciMzMTZ599dr3tysrKcNxxx+GSSy7BRx99hJdffhlTpkzBli1bVJC+tYJMhro0KFUlwPI3gT9eBop3OtcteR4YeIoa1R5kNCAm1ITFmIqgGj00TesSKVTWZXmkjUlh2pjOzFFa6usqEBERURuwDyciIupY7HuJvEdzOFTcVNMc0Bx2Z1mts0PTmwCDWSUS0eQ/hwOoyFfbOqQsdzgccGgO6GR/td4uG8MengLNECRFZyKSqkLoi3fJv2j1OM47pCzHcT6ubOiADlUpo9Q+JcXOwcQBH2iX0eyzZs3C119/jf79+6vlxhtvxEsvvdQo0D5v3jwEBwfjoYceUsFuGfn+wQcfqPUyIr61go0GoGinM7guQfbqkvob2GqAykIgJMadPqag3IpqmwOFFVbEhAb+pJOrdjFtjD/QWSzI/Oor5BfkqzIRUVchX55q7PLlC7DLlyj1Jc5Ztjuct2UROvlPpy5oU/9z3RZ6ndwCDAYdjHodDLLonH+7wol18h3pt3t+87W7TEREROx724MKNtYGMlWgU27JXIMSlqwNQGpVpbXbSFDU7t5e/tptNlQW5qMUlUBIHDRjkDvLslZTDl3ZXmcgU44PR+1fTQVHnceXxwWsCYNrg67OffVFO2Ao31t7v2s7h8dtDTo4YLPEoCp+SF1AVgNCdi2CvqbMvb3zTjmw57E0VMQNQVVUr9rH1KCrKUPMls9r6+t8nur5q2Cucz/1KJoDe3ucDmtQtPs3RFj+GsTsWaiel2sf1/7uOsABqzEcm/pdXltf5++UjG0fIrJonbpfPS/3vh63NQeyY0bi78iJCAurgk6vU8cYt+I2GBxV9eonf3X1/mr4q/ulyIo6yB18jilbjwnrH3HfL23pqrtqC3Xb+fetIW+hRl/3uo7c8zYO2vuB2l72d9XXVVZ/oWFXyAC81esZ92PK63PZhsuRVrEOehXWrv3N1YTvEy/BgsQL3LeDbSW4d+1xrXpPv9DrVfXYLkOKfsDZO+5vcb9KQxgeHPiNKleXl3aNQPuKFStgtVoxduxY97px48bhscceazR6/Pfff8f48ePd6+Sv7CfpY/Yn0D5l68PAf39RZznq6X0McMjVQOYEd452V6B9XZbzBckuqQr4QHtljR2bcp1neiKDTUiPCfZ1lagZOr0e5u7dYAgJVmUios7K4dBQY3OoILgEyKus9trFo2zzXG9X28tJbtleyrJY7c51Vnvtt0Ivcgfe9R5B+NrFZNDDpILzepiMephc65oqG/Vqf1lnbqqstnMejwH+rtaHd/d1NYiIWsUZgFOlunUNstJ6BsBUUMu1bW2gz31bp3ePanSu1QBblUfgyhmocu1b73FlRKMhyLNi0GTgXIPtnA/v8bhStkSqYKN70+py6KpLVKBQgo2eAThX/dVfnQH2yPS6pyMBw9I90FWVwO6wo6SwALpi11XgdftLQMtuiYEtIqPueQIw71mmAps6j8CdK8joDuJBQ1VMf9hD4uumkasqQPDevxq0rzPY6bmf/C3qPgnQm9yBP0veGlgK//XYpi745vm3xhKHooxj3O0v/4/Z+iVMlbm1AUNnmzfcT/4riR+J4viD3K+fzl6NjHWvuiKPdft7PL5qE03Dzl7noDI42f3ahRetRcr2z2ofyxXRs6OmuhoVZlNtIE+DTR+ElQNvd78rJUCaueMjxBUs9wgwegQaPQKXuVFDsbb7+c79ai/nH7/6LlhqihoENhuUoeHvjAuwPXa8+7WJrNyJo9fe4RGgdLaTK7DpCljKuvcGvopyc4x738HZn2Ds7jdqg5kSpKwNxNYGRl3lPEs3zOr7qrt9ZfdzN9+CnqXL3M9NtnMFO9V7pvbvL3Fn4euUa9z/bHSaHY+umoB9iaj9OyvzGWwKP9i9vm/JElywbcY+91VtCh3uGrKo3ropu5/FuPwPW9x3ffhozM58ut66G9ffjoTqHS3u+0Xydfgt/gz37ciavbj937vQGp+WD0aupe672cH5v+Hg3f9tcb9CUyLmGE+vt67nth/Ro+TnFvfNqQT+rRqL4BL5rHG+YlNzfkKQo7LFfX8KPwHb7HUphHXlhUgqXY3WKCirQo2h7jNcTrxE1mS3uJ/RXomaBr/DdJpDve9aVn8/bT8GNbn/XbmP1Lp9nf/+vK9TBdp37dqFuLg4mM11weuUlBRUVVUhPz9f3ee57eDBg+vtL9uuWbOmyWNXV1erxaWkxDlyvVf+AiCo9qNZOukhZ0IbcyUQ38+5obuTc0oIM7tvZxVVom9iGALZmj1FsNf+wxmYHF57BnL/3pxy2Yfsoy7/IK9iW3cctjXbOhC19+d0c32v9mQvwGKAATpYdDpY3F+QnD8jxFfJV2FZzAnufaNqsnDV5itqv5K5ttO5v5R57v9aj/+iyOz8gSYGFy3AEblvQr72aTqD+uvQ6aHBUPvXebvEFI+P0u6o9xzG5b2P1MoNtfsY3Ns23Hd7yGCsjTys3r6H5H0Im7SDPKb8MNcZ4YARdnVb/hrV+j2WPigzOa+cEyZHJWJqstR9eqMJOoMsZhiNRuiNZrXOYDSrRQXuawP9ZkNd2ajToay0BPEVBpiNRrXOfUJAAvtq29p1tcF+vZ6j9tuC/UHHYVt3jbZ2jXhUVyjZquGoroRmt8KhFlttWf7aoDlscNis6nO4OnaAOoGr6q0BpuwV0JfnAHar2g4O5z71/jrsKI/sjYKEMWofeUzZv8eq/6rAoK52P1kkGKX+OiQoKuvsWN19OnIjBrmvnIou/hdjN/yf2lYWvebcT+9ocFuz49nBn8CmM6KyrLRD+l77//WEzeL6nK/tTRv8rpuT+QQ2ho9y3+5T+gemb50hvWe97ZrqLe4cXD+Idvye5zA+f16L9d0QNgpvNAii3bDhPCRWb29x3y+Sr8FvcWe6b0tg6Pb19QNczXmm95vIsWS6bx9c8DlO3f1ki/sVmRLwVL+P6q07Z/vdGFwbREvcx75Lo4/Hx2m311t3/5ozWxVEeyfjQayOPMJ9u1v5Slyx5Wq0xqwB36LGEOK+PTH7fRyVM6fF/baHDMIneXWjRsVVm15EeqWMzN23HxIuwo+JKe7bFnsp7lv7XKvq+43tIOwMqQtVDS5aifE7W65vlT4Ur4deVm9dj11/ILPwqxb3La5yYJ355HrrTsn9A5G2vBb3XVJwDHaiLsBprSpGYlnLbaQet7wCJTV1MSVHZesCnGZbGaqt9QeL6h02GDVri/vWnSSotR/hnboTFq5dW/fdUX2GNPi82a/gaBv3rTux4t7R68+1Ufvu576uUfC6/W6nBvWtLTpPC8lvIOdfdZ/6DVN7W6dDsEkG+xicx5GHCwpDqSneYx/nby7XcZy39agKSUZiuFmdWFRrdUBpRC9kG12/y2R752O49nH9NcV2Q++EUPezMzoM2Jkzwfl8pX6126qTsfLYHre7p6cjPizSfWVydGQ/bNGd4d5WTUfqOoZ6aINa79AHYUKfOLWqrNSMZ9AFAu0SUA8K8jgjLa9v7W1JK9OabRtu5yKj4h944IEm73NYolEx8FxUDDoHjuBY5z++nJwmtzXU1LgfY9PuXAyI7pgzIr6y5N8c9/NNttiQ00y77It8QS8uLlYfFnqOtPYazWZD1eefo6qyCo6pp8PgccKK2h/f1x2Hbd1x5LO6PTXX95pghVmCFKKZbtRRU1GvTw+tKUe4raBVj5sSoiEyxJnuxaADMiuKkVy1ucX9ii2p6Bahg92VgsahoX/lcvQsWtzivr9GlOEvc90IH3Hs3hdhasUPnllJ92FV2Li6+leuxLW7b0FrzOjxOaz6ulQjEwvmYkzx5ypIb4NBnViQgL4zyO8sV+mM2BTUC5/FXV6/vgVvIsae6xzhpzM6cyfKpcWyGFx/TciOHIaCqMHukf1mWJGe+xN0tffLX53BCJ3eVHuiQG6boDeYYI3IgN4UDBk0o64G0KwwOqrVfQZ1AsHkd1dltcdnlPTh1fPnq3LQCSdAZ+xUX9E7jUDtD1TeUofdGdyWoLDKKyrBYAnsOmA1R8ChN9cGcwGdtRzG0p3O3KVy6bsKBDv3keCvuuy/9niFyYeqfVRKLQ0IK1iFkJItzkvma7dTf2X0pKqHMxBcaknBmvDxCA4udI6j1DQM2fIaLDV5tYFm534S2JGgswos165fkXgKNkWNV5+l8lqFVOfgzM23qeCyK8isFtjqyuq2HU9kvIpcU6q7bcYVz8fU3P+12IZ5xiT8X/c36627OOtRDC5f0uK+v0Ucj08Tkuqte2Lz6zBrdcHq5vygOwSrwhLct3tUZiG55B+0hryX5fO7qrIuQOftvteo7TtgU1NVjUpDXd9bXVXdKMjenIa/w2222n6+BXaHo9G+zhHlLZOr4T33DbJWobWqqqtRqdXtW1PTuvqiifra5c3eCioFR8Pn2sJr4mKV+nrsW1Xd8vcL97ZVFaj2OJlubRCgbY58BjV6rq18bWzWuriJ81gt/3tyqamqQqXO87Vp5XPVtCbeh80/V4lLFm52noBwDKlf39rDtYq1NkbknvqvugbyieYKYrrGo9cFK+tuy4lAnb5G7Sf7W3VmlBhjnduowCLqBzdryyXmRITobfXSH5ZYUrDX3qv2cdEoQOqqS01IEuIskjGiNgapAdsjRnoEQWvr5g7GyndjB/QGI8IjotEtvC6rRJgxHhsrJnk8P1cwtOFtoF+Mc35EV51tGIm1ISbntrV1U4FSdfC645SFpGFEsrk2kOsM6O5ynIccazF0HoFbdyPq9LVXlOgRETMU4yOC3c/VaE/EquD73Md27i+B59q5G2v3l78HJ2TCYQpz19dSeQTWdk/32MYZwFU18thPrtA5PU6CwHVVMmbeig01VzjrJt9fauutHl+lqHSuizSE4FybCRERETDonc9j74Cvax/Ledu5j3N7eWz1fUinx2EGEw6X7+u1J0N1OBJ7detb9R6uS+Dici3KcW2L+0UBOKvhyoFP1es1dM2UZWhz7fDmWiFAH+dVGi0Z0XBF6jBg4LBW7evaqtgEr+lU3+ItFosKoHtynY2XfOyt2bbhdi533HEHbrrppnpn9tPT07F68B0Yc9p1CDEFy8vaophYB0JWFakP3UqYkZBQ98Uq0Mil+nsrc1WbhgcbMbx3Wpty1MqXfNkvPj4+oH4UdTaOigpsfM75QyT+0ktgDAvsqy18je9rtnUg8ryirD001/fmh/SE1SJfDuW7pdbgr/NnxOiBvTGsX39YTAYEmfSwlEXD8X56vUt16y45rn+587mH9gOinJdmK0sXQ9tucQeUGl5u6BIRYsFlRw2st05XEArIVcMtOKRPEoZPHAGrBMVsDtgcGoz/bd0P2YHdEhAZn6L2tdk1xOeEArtbtSvMIeEw6uq+zoXrqxBlz29xP73e0Og705CKJa06IfFN4uXY4Kgb4RZmLcDp/97Tqvr+r9dM7Anu6749vPAbnLHrkXrb2NXVAq7FeZKg3BSLOUPfUSdO9LW588fseAUZhX9A08sJBSMcclm8nCRQi6n2ZIER+bEjsKP76fXy7ffa8BqMjhpALyNc9HV/a0e8qJMG8qM15VBURWSqx5T3p6m6EJF7flHb6+RHqzy2/ACtqIQ1IsJ5kkF+GOkNqEweA53J7H5fG8v3wliWVfdDzvlKOMsS3Pqvc6Rf0GEHwZY80DNzIQxF26C3VdV9D/O80+OHoiM4BpoMGnFtojlgKKodFSorPNMfNMgvaotIh2YMcaczQFUxjMU7ai/l97zU3yNNgQR+dXrUJA6v20R+YBRuhrEiu25bz5QCrjrI+z0kDpWxg9x5PSV+E779e+itFc66qaBxbfqJ2mC0raIcFUFmFCePRXmEM++qHNNYVYDkze975Fx1BqmdKQlcizM9w7q+V8FqCnfvm5DzG1L3fO0MNtemJHCPDFPrnAHpcksiFve9wx20lr9jNz6BxJJVzsv8JdAs4ZXagLJaVxtEXh5/Mn5Ovki1k5zEk4D2ff8cqVIBtBTInJn5DDaHjax/qf72W1v89yZHvmtw/UvVp+yZj0GtuFT/37DRWJl0CIItVvf7btLeH1p1qf7K0LHIs3hcgm7TI6FqK1oj2GxEsKXus8lQ2bo5C4w6R+PfgYbW9WlmQ+Pfmipo04ogm0Xq67GvyRHS7GeZuirKfduAtLgw2IxhqC53jiL0dt+bG9IL1RbnZ5ZS74owtQKxiUnoGRHjXhNXkoy9Za5+sUECigb7902Nca9TAS30QpY7pNHgsdyfYTrUhA/AgAznleuutcVlQ2Gvdo0Nb7yP61gRyb0wJDHe/VwtNQbsKnOlodU1kZKg7m/vbmlIDU5wv7/j8vthh+5Id4DQFWSs97g6HarNURjbN9m9WnoHh2UCthbHqWCm1WaFyWR2fzart1FtUM0UMwxHpafVezo7DdOgd9QGkt2P3fhvn/TRSI3McO9oqbBgc/T1ddvIon5re+7rfOwpfXqrq/ddgcrwvFOwOb+nKrsCoa6+yR0A1uthC47FOWm96gUMrbtvx5bqwtp+xxVsdI0grZvTJjOqFy6O7uVuIzkht33XG/XaVz1+7WM7H98ZODwxfiAclrpApbE6GXuKRtUe27mtfAYXl5QgMipKxRjUfXoTbovv736eQl/2BPKt99Qe36DyXqv/9Ho4qmqQfZzzCoju//c8HoxNdT8d1etX/YFKeQU9gqHuQKdOB31tcHSqwYAz5HuEmwy8OL1V45BvazKEeG+L+8l1kI0ToMxuxSMCciqz7tqIWsd+v8/fvbm5uSqe06NRPEe+D57Yqsft02hN/TDrvtQfygJgdOsGpDRtUBv3kzQyY9q4b3qrtvJs67rYWd1VutR5f/d22kB7amqqShEjZ6ZNJufphd27d6ugemxsbKNts7Ky6q2TbWV9U2S0e8MR8GLAcVdCHxTa6jqa9XrEh1uQU1qNvLIa1THIj7ZA9G92sfrBI53I4NQoGAxt/xIoHZJ8UDDQ7kUmE8KOOUadcNKbTGzrDsD3dcdhW3eM9v6Mbq7vjb52IaKiZAxE8zzC5E4RfYAbW5dnsFGvPPoy5+LiysnZYESnc7BIgzY46Xmgprzx6E/1t+4YhrBEhFoaDI04/XUZvgbID2gZaWqv/esuy2LHmMGHALHp9RM0Yprz2LJNg/00j2XGsQNhdUDlq5cl7I80WCuTa1MjWJ0Bv9oRp+4f8nKFQIgFQ9Kj1D41tScGgja1ckRebRDaRUahtpZdZ2pxX4OEp9Ro2foBQ5mE3pOleAuSSle1+JgFFVYsxcR6645c9yqC7c45aPbl3YwHsCqqLtCXXrEGV226Ea3xwICvUGV0ZTUFjsx+A0dnz2pyW3mp7enRqlzz2XV4ud9r9e6/YtMV6FbR8vv/x4QL8UPSxe7bZnsFHlhzTKvq+1LPl7AjtC4t46CiBTh3R8s/9qv0IXhg0Hf11p268wkcXPhli/uujDwCc7s9VG/d7WvvaNWl+h+l3Y5lMeHu2wlVW3DjhtZdBPy+6WSUmOsGy4zPXYXxWZ+2uF92UHdsiimvt+6o4i1ILFvb4r6GmlJUyj9WF02v3uutoXKdevy7cdQL6LSwn/DYVwK9rSFJttzBqtr95aqY1j2uvS7gKrEvoxHV+uB6QWY5luYKPuudgWi5PyEqDJbQEDU6U06OhRgzsatqnPuEmvxF7V/nVTfO9bagSIzrHaf2kRNjsr8Wfjq2lB1Ud2WO7Kt3XmUD95U3RkRHdMe5id3Ubzr1uHogr9c7td89jHVX5Rid28ttg0rpZcKJIdE42eysr3QfegyAA1PVyTYZjbivFru+9m9RUVETAbP273tjr13QYt87tdGa3gBOa9XjXtlozd21y75JcpGhDVce9narHtMZsvbUAzjSObH0/oe9ujc1NrN1AcORM9xBMrkCXAbjNfedypm13MPQJ1r1mD0brUkCRg5t2769pU+s3y+2WvdT0Gbpp7Zxx3AgPaPxoKcW2lqJbH7+E0d1NconTVLliPgU6IMaBN4sdamLuzr+FmNbBxq9FwcBd6pA+7Bhw1T+0cWLF2PCBOdkDIsWLcKoUa6zl3XGjBmDhx9+2D1Jqvz97bffcO+9Lf8gOFAyIaoE2uVHaX5ZNRIiWjfSwt+s3l2XQmBwaqRP60It0wcFIfXZZ9SXOykTEVEzakc1qaBLS8LrpxPYLwPb+GM0ob8zwN8M19g10Shkctzdamnyx74agusM3idpGs42N7iWb8inkli0XnBfs9eo/Md2W426DFwuex8T0xfDIzPV5ePyXchRFYfsmPtV/uS6EwHOQL8szuM5TzgM7NUdmUExauS+7B9jzsDuyrHOkwEOa+0JASlLSonaEwSaDeXGGISaDe4UGLKvmhSsFWREfENyAqI1GubF3J9JlBoGNJu7kkLVxwCkjStU5R3GxiOX1JUcbbA/+7XvY+xHHtMG3Jd/t7Rvg9ew4USM+9Jwkq7W7tvU+0ZG9Ksp7txzOEjAWF83r0NtYFkXFI7YULPzqp3aqyuyw/rX5h01NFrUaFqVxsmA5ORUBEdFuveNLOuFrY7T1eeY2lYFnw3uYLK6BF4FlQ04ZmACDHqDMxCsA8KTT8fW4sHOqzL0dfs4r8Zw3pZAstkSi7NMiYiPi1XzOsj+umFvIlf+Tepd6aCc26oAtKR/kiC20YhjTBZMlvkkauvrtLdVbXxeozU9ZIq/tgU/BzgnNmyThKPavi8R+QX5zZz232d9XQ0iCjCdKtAeEhKC6dOn4+abb8acOXPUD8TnnnsOr7zyiro/Ly8P4eHh6gz96aefjttuuw333HMPLrvsMrz66quoqKjA1KmNz8O3t8SIIKyqvaQ7uyQwA+2SNmZjtnOkV0SwEd1iW5NYh4iIiDotlTeyNt96U6Ibj/qSEJmhdmn+AstQIO3GNo4gPLdVQbSoJsdFfqbSj8hJALtMjmi1wSEnBmpPDmgyQaLNikxzGK4NSaydJNGZfz+v+xvQbDX105Oo1CS1fyUQ67BjSOIo9A1Ndaf7MFYYsSXmXhU0V9vW5tKurKhAcJDJGTiuPd7YvinOdDbq5AAQHjYGW8JtHklf6ybLqgs4aygPTsXoTOfkuK40LqXWo7G13OMSa49juPfVNFgSh2J4ct3pF4PNgm1FU2pPEGgNcq/W/ZWjZGZ0R3xYtCtVKiKL+mOr/kz3Zf11KR9cE00531OSO3x8r9q0D7U5UA0hR2NLYXK946vMNa58pK7jRPbCMRmJtY/pTM+TG3wjCuwVtY/jTMPjTDHgnHi4sqISoeHh6JcwAr2iJJ2U89gGWyx2ZL5em/tUAsfOfVQQufavM4+pAefED4LOFOQMHsuI68orkV95ujv4rOriSgNUm05I5g8INZpxX3CUcx8VuJYEQt+o7VoK1cssDHUzMdQ66ne0RkZT/5IOP7xV+0qYup5ecoVDy1c5yMm6IDlZFxdad7IuZkirHpOIiIioq+pUgXbx9NNP46qrrsLo0aMRFhaGO++8E2eccYa6T/IUvfHGG7jgggvU5ABffPEFrrjiCjz55JMYNGgQvvzyS7VPR4xod8kuqcJgBN5o73VZJWqUmhiYEtmm3OxERERE3iQBTqPZAiMsQNPT9DQt7tg2PmIsMKh/vTXNpQpodKn+kLNanZagUfbQ4Y0nNmxKkxfIj3mnVftmNpmUoXVpZ5xZeD0MbjytVqv1uarZu/adliES6Na6FBeNRAbXJrBoC35HJiIiIqJOGmiXQPmbb76ploackyfVkWD833//jY7mGWjfW9L6mc39CdPG+B9HZSU2T5qsZgSP+/Yb6ENbP/cAERER+b4PFz2lD284qSMRERGx7yWiTq/TBdr9geRXNBl0sNo15ARgoF3SxmxwpY2xGNGdaWP8g6bBlpPjLhMREZGfYB9ORETEvpeI/B4D7W0gM9knhAdhd1EV8sprYLU7YJLp6QMxbUwq08b4C11QELp99CEKCgpVmYiIiPyD9NuZn3zsLhMRERH7XiLyPwy0t5FMgCqBdhk4nFtajZSo4IBMGzMoJcKndaHW0xkMsPTvD2NOjioTERGRf/XhRERExL6XiPxX4AzD7mBJAZqn3TNtTLhKG8M830RERERERERERET7wkB7O0yIGkh52tfvLa1LG5MSodLkkH/QrFYUf/IJqr/+RpWJiIjIP0i/XfTxJ2phH05ERMS+l4j8E1PHtMeI9uLACbSv8kgbMzg10qd1of0jP8z33nW3szz1dIA5XomIiPymD8+6805Vjpg8CTqTyddVIiIiCmjse4nIGxhob6OIYCOCjHpU2xzILq1GIKi2SdqYUlVm2hg/ZDAg9LDDUF1TrcpERETkR334hMPcZSIiImLfS0T+h4H2NtLpdEiKtGB7fgWKKqwqt7nF5N8/jP7NKoXVzrQx/kofFIS0l19CTk6OKhMREZF/kH4745VXfF0NIiKiLoN9LxF5A3O0H4DEiLpgZk5JdUCljRnEtDFERERERERERERErcJAeztNiLrXzydE9UwbExZkQGZsqK+rREREREREREREROQXGGhvp0B7tp8H2tfv9UwbEwm9XufrKtF+clRWYsvkY1F87nmqTERERP5B+u1NkyaphX04ERER+14i8k/M0X4AAinQzrQxAUDTYN2xw10mIiIiP+rDt7MPJyIiYt9LRP6MgfYDEBZkVGlWyqrtfh1ol7QxMqJdyPPpEce0Mf5IFxSE9LffQlFhoSoTERGRf5B+u9u777jLRERExL6XiPwPA+3tMKq9LLdcBdvLqm0q+O5vmDYmMOgMBoSMGIGynBxVJiIiIv/qw4mIiIh9LxH5L+ZoP0CBkD6GaWOIiIiIiIiIiIiI2o6B9i4eaPdMGxNqZtoYf6bZbCj95lvU/PSTKhMREZF/kH675Jtv1MI+nIiIiH0vEfkn/8tz0skk+XmgfcPeMljtzokzB6ZGQK/X+bpK1EZaTQ323HSTszxlCmA2sy2JiIj8pA/ffcONqtx3+V/QGfkVnYiIiH0vEfkbfos/QAkRdRNWZZdUw9+s3F3kLg9OjfRpXegA6fUIPvhgWGtqVJmIiIj8hF6PkIMPdpeJiIiIfS8R+R8G2g+QxWRAVIgJRRVW7C2ugqZp0On8Y1R4lbUubUxYkKSNCfN1legA6C0WZMyZjZycHFUmIiIi/yD9dre33vR1NYiIiLoM9r1E5A0cMtMOEsOdo9qrbQ6UVPpPbux1WSXutDGDUiOZNoaIiIiIiIiIiIioDRhobwdJkXWjh/f6UZ72f3bWpY0Zmhbl07oQERERERERERER+SsG2ttBgh9OiFpRY8PGnDJVjgw2oVtsiK+rRAfIUVWFbaecipKLL1FlIiIi8g/Sb285+RS1sA8nIiJi30tE/ok52ttBUoT/jWhfvbsEDmfWGAxNi/SbvPK0Dw4Hqtevd5eJiIjIj/rwf/91l4mIiIh9LxH5Hwba20F8eBAkTq1pQI6fBNpX7qpLGzMknWljAoEuKAhpM19DUVGRKhMREZF/kH47fdZMd5mIiIjY9xKR/2GgvR2YDHrEhZqRW1aDnNJqOBxap55YtLjSii155aocF2ZGikeOefJfOoMBoWPHojwnR5WJiIjIP0i/HTZunK+rQURE1GWw7yUib2CO9naSHBWs/lrtGrJLO/eo9tW7i9Xoe9ckqEwbQ0RERERERERERNR2DLS3k9TaQLvYU1SJzuwfz7QxaZE+rQu1H81mQ9nPP8O6ZIkqExERkX+Qfrv0p5/+v707AY6yvv84/t3NnRBykANCDi4J0FgBCwG0AoLV/tV/dRAVjyKoFfyTKqNg1dazg39H284o1tJapVrtTMXK/EUL44GAchgEBDksRyAkHDkIuc/d/c/vF7MQBQV2n+fZ59n3a+aZ/LIJm4dfNvvZ/e5vvz99kOEAAJC9AOyJ1jFB0jflRKG9rKZZLsyTkHSssU0OHut8IaBPUqxknLSRK+zN19Ym5bPv7hxfdplIdLTVpwQAAM4ww8tmzdbj/E2fiyuSh+gAABiJ7AVgBB7FG7CiXRXaQxWr2R3M7ZbYgh9Ie3uHHgMAADtleIF/DAAAyF4A9kOhPUhioyL0xqJVDW1ypLZFPF6fRITghqhbD9b6x6o/O5zDHRsref/8p1RUVOgxAACwB5Xb/Ze8afVpAAAQNsheAEZgyYwBq9o7vD45Whd6G6Ierm2WI1+fV25qvKQk0FoEAAAAAAAAAAJFoT2IslPi/ePyENwQddOBE5ugjshlNTsAAAAAAAAABAOF9iDKSj7RrqM8xPq0e70+f3/2CLfID7OTrD4lBJm3pUUO3Hyz1M2Zo8cAAMAeVG7vn3aTPshwAADIXgD2RI/2IMpKjhOXS8TnC70V7XsqG6S+pUOPh/TuKfHR/Oodx+uVls1b/GMAAGATXq80b97sHwMAALIXgP1QbQ36hqgxUlnfqjdE7fB4JVItHw8Bmw7U+MfDc2gb40Su6GjJeu45qa2t1WMAAGAPKrezFz7vHwMAALIXgP1QaA+y7JQ4XWhXG6Ierm2RnNQTfdut0tLukR2H6/Q4PjpChvROtPqUYABXZKQkTp4kzRUVegwAAOyU4ZOtPg0AAMIG2QvACKGx3NpB8k4qrB+obpJQsP1QrbR7fHqserOHyip7AAAAAAAAAHACKq5BltcrwT8+cKxRQsHm0s5NUJUROSmWnguM4/N4pOmzz6R98xY9BgAA9qByu3HDZ/ogwwEAIHsB2BP9JYIsIzFGYiLd0trh1SvafT6fuNQOqRapaWyTvZWdBf+0HtGSkxpn2bnAWL7WVjl424zO8cZikagophwAAJtkeOn06Xqcv+lzccVb33oQAAAnI3sBGIFCe5C53S7JTY2X3RUNUt/SITVN7ZKaYN2mVptKT2yCOiI32dKiPwzmckn0wIHS4enQYwAAYKMMHzTQPwYAAGQvAPsJqdYx+/btk4kTJ0psbKzk5+fLsmXLTvu9HR0d8uijj0pubq6kpKTItddeK4cOHZJQkNfr5D7t1rWP8Xp9Ury/xv+c7cLcVMvOBcZzx8VJ/3f+T5IWL9ZjAABgDyq3By5bpg8yHAAAsheAPYVMod3r9co111wj/fv3l6+++kqKiopk6tSpUlJScsrvf/rpp+XVV1+Vv//977JhwwZpb2+XKVOm6FYtoVRoLz1m3YaoeyobpLa5XY8HZ/SQpHhaiQAAAAAAAACAYwvtq1atkr1798rzzz8veXl5MmfOHBk7dqy88sorp/z+l19+WR555BG55JJLZPDgwfLnP/9Z1q9fr6/Datkp8f53/ao+7VYp3n/MPx7Vn9XsAAAAAAAAAODoQrsqko8cOVISEhL8l1100UV6tfqpLFq0SK688kr/5129xxsbrWvV0iU2KkKykmL1+HBtizS2dph+Dg2tHbLjUJ0eJ8ZGypDePU0/B5jL29IiB2+/Q+rvu1+PAQCAPajcLp05Ux9kOAAAZC8AewqZzVDLysqkb9++3S7LysqS8vLyU37/5MmTu33+xz/+UXr37i1Dhw6VUDAwvYeUH+8sdpZUNUpB3yRTf/6mAzXi/bqLzsjcZIlws7GW43m90rRunX8MAABswuuVxrVkOAAAZC8AOwuZQntLS4vExMR0u0x93tzc/L3/9q233pKnnnpKFi5cKNHR0af8ntbWVn10qaur8/eGV0ew9esVL6v/01np3n20Xob1SRSzqD71xSXVauAvtBvxfzxT6merc7LyHMKBLzJSMv/3Kamvq9dj5ttY3K7Nw1ybO9fBZHb24gT+buw11yq3+zz9tH/M34dxc40zw1ybh+x1Dv5u7DXXZK95cw3mOtQYeXu2rNC+YMECfXQpLCyUtLS0bt+jnpzHxcV95/W89957ctNNN8ltt90ms2bNOu33qUL8448//q3LKysrpa2tTYItweuVlpZmXevedqBCxmaZN9Ulx5rlYGWtHuekxIq3qVYqmqy9AdfW1uo7Z7c7ZLoVOZJ39Ghpqq2VymPHmGuj55rbtWmYa/Oo++pgMjt7cQJ/Nzac68LR+kPLsRN77MCgucb3Yq7NQ/Y6B383Npxrste8uQZz7eDsPZnLp/5aLHD8+HF9dHnjjTdk+fLlsnr1av9lv/71r6W4uFhWrFhxyutYunSp3HDDDXLjjTfqTVO/64/+VKvqcnJypLq6WpKTk8UIi1bvk9KvN0Odf0W+JMVFiRleW39Adh2u1+MbR+fI+Sa3rTnVHbMqqqSnp3PHzFw7Brdr5tqJVC736tVLP/Do2TPwvT2syF504j7KPMw1c+1E3K7NQ/Y6B383zLUTcbtmrp3oeJCf94bEinb1BPvkJ9ljxozRK9zVZqZdG6KuWbNGJkyYcMp//+mnn+oi+4wZM+TFF1/0b4Z6OqoNzTdb0yiqOG/Uq3LnZSRK6bHO1jcl1U0yMjdFjFbT2CZfHW1Qu8NKz7hIKeibLO4Q6M+ufj9GzjVEfB6PtG7fLp6aGnGlpTHXJuB2bR7m2hzBvo+2IntxAn839plrleEtO3boceywYeKKiAjyGToHt2vm2mnIXmfhPso+c032mjfXYK5DjZG35ZD5K1EF9dzcXCkqKpIDBw7o4rlazT5z5kz99fb2dqmqqvL3h7rzzjt1u5knn3xSr4xTX1OH+r5QMSC98wWDrj7tZli/r7qrNbuM6d+LTVDDiK+1VUpvuFHqZ83WYwAAYA8qt/dPvV4fZDgAAGQvAHuKDKVXE1QrmNtvv13y8/OlX79+8q9//Uvy8vL8K9gnTpwoJSUletX7zp079eUZGRndrmflypWnXQVvttzUeImJdEtrh1e+OtIgXq/P0NXlbR1eKd5fo8eRbpeM6p9q2M9CCHK5JDIrS7wejx4DAACbcLkkKivLPwYAAGQvAPsJmUK7MnjwYN0u5lRU8fzkdvIWtZY/K5ERbsnvnShby2qlud0jB441Sf+0E6vcg+2LsuP65yjnZydJj5iQ+vXCYO64OBn4wftSUVGhxwAAwB5Ubg/66EOrTwMAgLBB9gJwdOsYp1KF9i67DtcZ9nPUCw9rdlf5Px87oJdhPwsAAAAAAAAAcAKFdoPlZyb63wG884hxfdp3HK6TyvrOvtz90+IlJzXesJ8FAAAAAAAAADiBQrvBEmIida92RRXCK+pbDFnN/vFXlf7Pxw/u3rce4cHb2irlc4qk4eFf6zEAALAHldsH/2eOPshwAADIXgD2RBNvExRkJcmB6iY93nqwViYPiw3q9ZdUNUpZTbMe90mKlcGZPYJ6/bAJj0caPvrIPwYAADbK8A+/7tFOhgMAQPYCsCUK7SZQG5O+9+VhUfu3qg1LJw3NEFdXP5kgWNltNXt6UK8b9uGKipLMxx+T+vp6PQYAAPagcrv3E4/7xwAAgOwFYD8U2k2QFBclA9ISZG9lo1Q1tEn58WbJTglOD/W9lQ2yp6JBj1MTouT8vklBuV7Yj3pinjx1qrRVVPAkHQAAm2V4yvXXW30aAACEDbIXgBHo0W6SC3KS/eMtB48HrTf7iu1H/J9PHpopbjer2QEAAAAAAADATBTaTezTHvl1EXxz6XFp93gDvs6dh+vl4LHO3uyZPWPkguwTxXyEH5/XK62794inpESPAQCAnTJ8tz7IcAAAyF4A9kSh3SRx0RH+ti5NbR7ZWhbYqnaPt/tq9p8M681q9jDna2mR/T/7mdTNmKnHAADAHlRu77v6v/VBhgMAQPYCsCcK7SYaO7CXf7xub7Vu/XKu1L+vqG/V49zUeBnaJzEo5wh7i0hJEVcSffoBALBjhqsDAACQvQDsic1QTZSdEqePsppmKT/eojdHHZTR46yvp7a5XT7YeVSPXS6Rq37YR1xqgLDmjo+XQZ9+IhUVFXoMAADsQeX24HVrrT4NAADCBtkLwAisaDeRKoZfPCjN//mHO4+e9ap29f3vfHFIWjs6e3CP6pciOakUVQEAAAAAAADAKhTaTab6tKcnxujx/uom+c/RhrP695tKa2T7oTo9jo+OkMt/0NuQ8wQAAAAAAAAAnBkK7SZzu10yaUiG/3O1Or3d07k6/ftU1rfKO18c9n9+7Yi+Eh9N9x908ra2yqF586Xxt7/VYwAAYA8qt8vvn6cPMhwAALIXgD1RaLfAD7OTpH9aZ7uX6sY2eX9HZ7/179Lc5pHX1u33t4y5MC9FCvqy6SVO4vFI/bvvStsHH+oxAACwCY9H6pYt0wcZDgAA2QvAnlgObVGv9muG95XnPtotajH7mt1VkpMSL+dnJ522yP7ypyVS2dCmP8/sGaM3QAW63a6ioiT9Vw9IQ32DHgMAAHtQuZ354K/8YwAAQPYCsB8K7RbJ6Bkr/1XQR97Z2tkK5p8bD4rLJd9apa7axby+4YAcretsBdIjJkKmj+0nsVERlpw3Qpd6Yp76859LR0UFT9IBALBbhk+fbvVpAAAQNsheAEag0G6hsQN7SdnxZtlcelw6vD55fUOp3ixVtZaJcLvkP0frZeP+Gv21riL7HT8eICkJ0VaeNgAAAAAAAADgJBTaLW4hc93IbD1WxXZlW3mtPr4pIzFGbhmTJ+mJMaafJ+zB5/VKe3m5eKqqxZeWpnbetfqUAADAmWb4oc53OUZl9REXGQ4AgKHIXgBGoNBuMbfbJVMvzJb+aQmy/Msj0tTWfRPLqAiXjBuYJhOHpEtMJO1icHq+lhbZd9lP9DhzY7FIjx5MFwAANsnwvZMn63H+ps/FFR9v9SkBAOBoZC8AI1BoD5GV7aP6pcrwnGT56ki9HK1rEY/XJ5k9Y2VwZqLERVNgxxneluLixOfrbDUEAADsQ2U4AAAgewHYF4X2EBIV4daboX5zQ1TgTLjj42Xw5xuloqJCjwEAgD2o3B6yeZPVpwEAQNggewEYgSbOAAAAAAAAAAAEgEI7AAAAAAAAAAABoNAOOIS3rU2OPPKIND7zrB4DAAB7ULl9+De/0QcZDgAA2QvAnii0A07R0SG1S96Stnff1WMAAGATHR1y/M0l+iDDAQAgewHYE5uhAg7hioyUtF/+UhoaG/UYAADYg8rt9Hvv8Y8BAADZC8B+eCQPOIQrOlp6zbpLPBUVegwAAOxB5XbarFlWnwYAAGGD7AVgBFrHAAAAAAAAAAAQAArtgEP4fD7pOHZMvMeP6zEAALBXhquDDAcAgOwFYE+0jgEcwtfcLHsv/rEeZ2wsFunRw+pTAgAAZ5jhu8ddpMf5mz4XV3w88wYAgIHIXgBGCNtCe9dqobq6OnG7WdhvJK/XK/X19RIbG8tcGznPTU3S4PH4b9eRXq+RPy7scbs2D3NtHnXfoRi1opbsNQ9/N/aa629muLujI8hn6QzcrplrJyJ7nYP7KHvNNdlr3lyDuQ6n7A3bQnt1dbX+mJeXZ/WpAMHXty+zCiCgjExKSgr6DJK9wBno04dpAsIQ2QtYiOwFwlK1Ac97w7bQnpqaqj+WlpYaUkxA91eKcnJy5ODBg9KzZ0+mxkDMtXmYa+baiWprayU3N9efkcFG9pqH+yjm2om4XTPXTkT2Ogf3Ucy1E3G7Zq6dqNbA571hW2jvesuLKrJT/DWHmmfmmrl2Gm7XzLUTGfW2ULLXfNxHMddOxO2auXYistc5uI9irp2I2zVz7URuA5730mAJAAAAAAAAAIAAUGgHAAAAAAAAACAAYVtoj4mJkUcffVR/BHPtFNyumWsn4nbtnLnmd2ke5pq5diJu18y1E5G9zsF9FHPtRNyumWsnijHwea/L5/P5gn6tAAAAAAAAAACEibBd0Q4AAAAAAAAAQDBQaAcAAAAAAAAAIABhWWhvbm6WmTNnSmJiomRmZsqCBQusPiXHa2trk2HDhsnHH39s9ak4VklJiVx99dWSlJQk/fv3l6eeekq8Xq/Vp+VIu3fvlssvv1zfhwwePFhef/11q08pLNx1110yYcIEq0/Dsd555x1xuVzdjuuuuy5o10/2mo/sNR7Zax6y1xpkr7HIXuche41H9pqH7LUG2Wvv7I2UMDRv3jzZunWrfPbZZ1JeXi5Tp07Vhclp06ZZfWqO1NraKrfccovs3LnT6lNx9AM6VWQfMWKEbNq0Sfbs2SO33nqrJCcny+zZs60+PUdRL16ouS4sLNT3I7t27ZKbbrpJsrOzZfz48VafnmN98skn8pe//EUuueQSq0/FsXbs2CFXXnmlLF682H9ZMDeHIXvNRfYaj+w1D9lrDbLXeGSvs5C9xiN7zUP2WoPstX/2hl2hXa2o++tf/yr//ve/ZejQofqYO3euvPjiixTaDboB33zzzUZcNU6yYcMGXVwvLi6WuLg4GThwoNx77716pTWF9uA6evSoFBQUyMKFC/WKdvUinVplvXTpUgrtBj6gVq/qX3zxxUb9CIjoF0PVO4/S0tKCPh9kr7nIXnOQveYhe81H9pqD7HUOstccZK95yF7zkb32z96wbB2zZcsWaW9vl3Hjxvkvu+iii3SB0ufzWXpuTrR69WqZOHGiflUOxhkyZIh++4sqsndRb39pbGxk2oOsT58+smTJEl1kV9Q7Y1atWqVb9sAYqg2SerfGpZdeyhQb/ATxvPPOM+S6yV5zkb3mIHvNQ/aaj+w1B9nrHGSvOche85C95iN77Z+9YVloLysr069aREdH+y/LysqSlpYWqa6utvTcnGjWrFny+9//XhISEqw+FUdLT0+Xyy67rNvbFl9++WUKkwZTd86qhcyAAQOkqKjI6B8XllRrnj/96U/6fgTGz/XatWvl/PPP13sPPPLII/qF6WAge81F9pqD7LUG2Ws8stc8ZK9zkL3mIHutQfYaj+x1RvaGZaFdFdS/2Xun63P11nbA7jwej27Xo144uv/++60+HUd78803ZcWKFRIfHy/btm2z+nQcR73LSLWMeeyxxyQjI8Pq03E0tV9JfX29uN1u3atOvbChPqq+6sFA9sLpyF7zkL3GInvNQ/YCgSF7zUP2GovsdU72hmWP9tjYWP2E/2Rq9a9yctsNwK4blsycOVO3kXnvvff0271gnOHDh+uPDQ0NMmPGDNm3b59u2YPgeOmll/Qry7/4xS+YUoP17dtX6urq/C2RunoEqhftfve730lERERA10/2wsnIXnORvcYie81D9gLnjuw1F9lrLLLXOdkbloV2Nalqpa8q3kRFRflf0VBFgF69ell9ekBAr+jfeuut8vbbb8tbb70lkyZNYjYN2hRm3bp1cs0113TrFbh//36pqqrSb2dEcLzxxhuyefNmSUlJ0Z+rF0k7OjokOTlZtm7dKrm5uUx1EJ38YKPrdq3mvKamJuCNYsheOBXZaw6y1zxkr7nIXuDskb3mIHvNQ/Y6J3vDsnWMeiUuMjJS9+PpsmbNGhk9ejQrUWFr99xzjyxdulSWLVsmV111ldWn41iqoD5lyhQ5cuSI/7KNGzfq9jG8WBdc//jHP/SO4GojTXWo3pc/+tGP9FjtrYHgWblypd7QV707o4t6kUM90AjGgw2yF05F9pqD7DUP2Wseshc4N2SvOche85C9zsnesFzRroph06dPl/vuu0/+9re/SUVFhTz33HOyaNEiq08NOGcbNmyQF154QRYuXCgXXHCBXlmtqLe9dK0GRnCMGjVKH3feead+a1FJSYk88MAD+gGf6vOF4Ondu3e3z9VKdvXuo379+jHNQaZebFYPOFSbnieffNJ/u37wwQeDcv1kL5yI7DUP2Wsestc8ZC9w9she85C95iF7nZO9SlhWhVRxbNiwYVJYWCjTpk2Thx56SK6//nqrTws4Z0uWLNEf58yZo1uXdB0jRoxgVoNMFdNVe57o6Gh9Jz179mwpKiqSJ554grmGbSUkJMjy5ct1azW1+ly1oVIPPubOnRu0n0H2wmnIXvOQvXAishc4e2SvecheOFGCCc97XT61vS0AAAAAAAAAADgnYbmiHQAAAAAAAACAYKHQDgAAAAAAAABAACi0AwAAAAAAAAAQAArtAAAAAAAAAAAEgEI7AAAAAAAAAAABoNAOAAAAAAAAAEAAKLQDAAAAAAAAABAACu0AAAAAAAAAAASAQjsAQ4wfP17Wrl37rcsnTJggLpfLf8THx+vLtm3bdsbX3a9fP1m8eHGQzxgAAHsjewEAIHsBWIdCO4Cga2pqku3bt8vo0aNP+fV58+ZJZWWlVFRUyBdffCF5eXkyZcoUaW9v57cBAADZCwBAyON5L4BvivzWJQAQoDVr1siYMWMkMvLUdzFqFXtaWpoep6eny7PPPisZGRm6OD98+HDmHwAAshcAgJDG814A38SKdgBnbP/+/brdy2OPPSYJCQlSVFR0yu/78MMPZdKkSWd8vXFxcfpjcnKy/7LXXntNCgoK9NdGjRqlH8Sc7Msvv5SxY8dKbGysXHjhhbJ582Z9+ccffyzZ2dly11136YL+M888w28YAGBbZC8AAGQvz3sBe6DQDuCsFRcXy9atW+Xee+8NuNDe0dEhf/jDH3RhXPVeV1555RW5++675eGHH9YF9cmTJ8tPf/pTKSsr8/+7RYsWyfz58/XXk5KS9Pd3KS8vF4/Ho/u+T5s2jd8wAMD2yF4AAMhenvcCIc4HAGeopKTEp+42Pvjgg9N+T3V1tS8zM9Pn9XpP+fXx48f7YmJifElJSfpwu936ePjhh/3/Zvjw4b6HHnqo278rLCz0zZ8/X4/z8vJ8c+fO9X/t7bff9sXFxenxypUr9Tnu2bOH3ysAwPbIXgAAyF6F571A6GNFO4Czlpube9qvrVy5UiZMmKBbzJzOrFmzZMuWLfpQK/RUK5oFCxbICy+8oL++a9cuKSws7PZvxo0bpy/vMmTIEP9YtZxpbm4+43MEAMBuyF4AAMhenvcCoY3NUAGctaioqNN+7UzaxqjCeFebGPVx5MiRugf78uXLZc6cOfpyn08tTD/B6/XqNjNdoqOjz/kcAQCwG7IXAACy92weHwAwHyvaAQTV2W6EevIDhJ49e+pxfn6+rF+/vtvX161bJ0OHDg3aeQIA4BRkLwAAZC8A67GiHUDQqM1K29raZMCAAd/5fU1NTVJVVaXHdXV18v7778uKFStk2bJl+jK1ql1ttFpQUCBjx46VV199VW++unjxYn5bAACQvQAAWIbnvQBOh0I7gKCuqLv00ku/9/ueeeYZfSixsbEyaNAgeemll+SKK67Ql91xxx16hbvq215aWqpby6xatYoV7QAAkL0AAFiK570ATseldkQ97VcBAAAAAAAAAMB3okc7AAAAAAAAAAABoNAOAAAAAAAAAEAAKLQDAAAAAAAAABAACu0AAAAAAAAAAASAQjsAAAAAAAAAAAGg0A4AAAAAAAAAQAAotAMAAAAAAAAAEAAK7QAAAAAAAAAABIBCOwAAAAAAAAAAAaDQDgAAAAAAAABAACi0AwAAAAAAAAAQAArtAAAAAAAAAADIuft/9DoFJqFKMhEAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "r = al_ae.r\n", "fig, axes = plt.subplots(1, 3, figsize=(15, 4), sharey=True)\n", "for ax, (l, data) in zip(axes, tm_results.items()):\n", " tm = data['tm']\n", " label_map = {0: 's', 1: 'p', 2: 'd'}\n", " \n", " ax.plot(r, data['u_ae'], label='AE', linewidth=2, alpha=0.6)\n", " ax.plot(r, tm.u_ps, '--', label='TM PS', linewidth=2)\n", " ax.axvline(tm.rc, color='tab:red', linestyle=':', label='$r_c$')\n", " \n", " ax.set_xlim(0, 5)\n", " ax.set_xlabel('r / Bohr')\n", " ax.set_title(f\"{label_map[l]}-channel (l={l})\")\n", " ax.grid(alpha=0.3)\n", "\n", "axes[0].set_ylabel('u(r)')\n", "handles, labels_plot = axes[0].get_legend_handles_labels()\n", "fig.legend(handles, labels_plot, loc='upper right')\n", "fig.suptitle('Al: AE vs TM Pseudo-wavefunctions')\n", "plt.tight_layout()\n", "plt.show()\n" ] }, { "cell_type": "markdown", "id": "46cee7ae", "metadata": {}, "source": [ "## 小结与展望\n", "\n", "- 复用 AE 求解得到的径向网格与能量,针对 s/p/d 通道构建了 `continuity_orders=4` 的 TM 伪波函数。\n", "- 通过 `TMResult.norm_error` 与 `eval_derivatives_at` 的显式比较,确认范数误差低至 $10^{-6}$ 量级、各阶导数在 $r_c$ 处几乎重合。\n", "- 可视化结果显示内外区无缝拼接,验证了所选 $r_c$ 与网格的合理性,为后续的势反演提供了鲁棒初值。\n", "\n", "下一讲(03-potential-inversion)将基于这些赝轨道使用 `invert_semilocal_potential` 提取赝势,从而完成 TM 方法的第二步:由波函数反演出局域势并做 KB 投影。欢迎继续探索。\n" ] } ], "metadata": { "kernelspec": { "display_name": ".venv", "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.12.2" }, "nbformat": 4, "nbformat_minor": 5 }, "nbformat": 4, "nbformat_minor": 5 }