diff --git a/apps/elph.cpp b/apps/elph.cpp index 2894544..7ad2099 100644 --- a/apps/elph.cpp +++ b/apps/elph.cpp @@ -67,6 +67,14 @@ int main(int argc, char const *argv[]) ElectronPhonon.read_mesh_geometry_from_msh_file(mesh_band_input_file); ElectronPhonon.read_mesh_bands_from_msh_file(mesh_band_input_file); + unsigned int nb_bands = ElectronPhonon.get_number_bands(); + std::cout << "Number of bands: " << nb_bands << std::endl; + if (my_options.nrLevels > nb_bands) { + std::cout << "Number of bands requested is greater than the number of bands in the mesh file. Resetting to " << nb_bands << std::endl; + my_options.nrLevels = nb_bands; + } + + ElectronPhonon.compute_electron_phonon_rates_over_mesh(); ElectronPhonon.export_rate_values("rates_all.csv"); diff --git a/apps/mpi_DOS_MeshBZ.cpp b/apps/mpi_DOS_MeshBZ.cpp index f9872fb..8ef5e47 100644 --- a/apps/mpi_DOS_MeshBZ.cpp +++ b/apps/mpi_DOS_MeshBZ.cpp @@ -258,7 +258,7 @@ int main(int argc, char *argv[]) { std::cout << "Exporting DOS values to file: " << out_file_bands << std::endl; export_multiple_vector_to_csv(out_file_bands + ".csv", list_header, dos_at_bands); - const std::string python_plot_dos = std::string(CMAKE_SOURCE_DIR) + "/python/plot_density_of_states.py"; + const std::string python_plot_dos = std::string(CMAKE_SOURCE_DIR) + "/python/plots/plot_density_of_states.py"; bool call_python_plot = plot_with_python.isSet(); if (call_python_plot) { std::string python_call = "python3 " + python_plot_dos + " --file " + out_file_bands + ".csv"; diff --git a/python/plots/plot_phonon_rate.py b/python/plots/plot_phonon_rate.py index cf02808..31e1495 100644 --- a/python/plots/plot_phonon_rate.py +++ b/python/plots/plot_phonon_rate.py @@ -24,7 +24,7 @@ def scatter_plot_rates(filename): data = np.loadtxt(filename) energy = data[:,0] energy -= np.min(energy) - for i in range(1, data.shape[1]): + for i in range(2, data.shape[1]): ax.scatter(energy, data[:,i], label=f"Mode {i}", s=1) ax.set_xlabel("Energy (eV)") ax.set_ylabel("Rate (s$^-1$)") @@ -34,24 +34,43 @@ def scatter_plot_rates(filename): ax.legend() def plot_rates(filename): - fig, ax = plt.subplots() data = np.loadtxt(filename) energy = data[:,0] energy -= np.min(energy) dos = data[:,1] - ax.plot(energy, dos, label="DOS") + nb_modes = data.shape[1] - 2 + fig, ax = plt.subplots(nb_modes, figsize=(6, 4*nb_modes)) + for i in range(1, data.shape[1]-3): + ax[i-2].plot(energy, data[:,i], label=f"Mode {i}", color="black") + ax[i-2].set_xlabel("Energy (eV)") + ax[i-2].set_ylabel("Rate (s$^-1$)") + ax[i-2].set_title(f"Mode {i}") + + # max_rate = np.max(data[:,2:]) + # dos /= np.max(dos) + # ax.plot(energy, dos, label="DOS", color="black", linestyle="--") + + # ax.legend() + fig.tight_layout() + + +def plo_dos(filename): + data = np.loadtxt(filename) + energy = data[:,0] + energy -= np.min(energy) + dos = data[:,1] + fig, ax = plt.subplots(figsize=(6, 4)) + ax.plot(energy, dos, label="DOS", color="black", linestyle="--") ax.set_xlabel("Energy (eV)") ax.set_ylabel("DOS") fig.tight_layout() ax.legend() - - - if __name__ == "__main__": parser = ArgumentParser() parser.add_argument("-f", "--filename", type=str, required=True, help="Filename to plot") args = parser.parse_args() - # scatter_plot_rates(args.filename) + scatter_plot_rates(args.filename) plot_rates(args.filename) + plo_dos(args.filename) plt.show() \ No newline at end of file diff --git a/python/plots/test_gradient.ipynb b/python/plots/test_gradient.ipynb new file mode 100644 index 0000000..b864c22 --- /dev/null +++ b/python/plots/test_gradient.ipynb @@ -0,0 +1,110 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "# 3D plotting\n", + "from mpl_toolkits.mplot3d import Axes3D" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAGbCAYAAADX6qdpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADj40lEQVR4nOx9d3wb9f3+o2nJe894xTNOHGd5ZjCaElpGKStAIZCwS9pCoKVsWvYopRAgLQ0tO6GUL6uMQiAhw1neO957yJK8JFvzfn/k97lKssaddLJk557Xq6++cKS7j6S7z3Pv9TwCiqIo8ODBgwcPHi4g9PUCePDgwYPH/ABPGDx48ODBgxF4wuDBgwcPHozAEwYPHjx48GAEnjB48ODBgwcj8ITBgwcPHjwYgScMHjx48ODBCDxh8ODBgwcPRuAJgwcPHjx4MAJPGDx48ODBgxF4wuDBgwcPHozAEwYPHjx48GAEnjB48ODBgwcj8ITBgwcPHjwYgScMHjx48ODBCDxh8ODBgwcPRuAJgwcPHjx4MAJPGDx48ODBgxF4wuDBgwcPHozAEwYPHjx48GAEnjB48ODBgwcj8ITBgwcPHjwYgScMHjx48ODBCDxh8PApzGYz9Ho9TCYTKIry9XJ48ODhBGJfL4DHmQmKomAymaDX6zE9PQ2RSASRSASxWEz/v0AggEAg8PVSefDg8f/BEwaPOYfZbIbRaITJZAIACIVCCAQCGI1GGAwGmih4AuHBw7/AEwaPOQNFUTCbzTAYDKAoit78BQIBhEIhhEIh/TqKongC4cHDzyCg+MQxjzkARVEwGAx0VEE2e1LDIGTh6L2EbE6cOIG8vDwEBwfzBMKDxxyDjzB4eB0kqjCZTHT6iQ0IEQiFQhgMBj6FxYOHj8ATBg+vgRS2u7q6oNPpkJGR4fEGTt7Pp7B48Jh78ITBwyuwTEHNzMxgenqas83aNotqGYGQf+cJhAcP7sETBg/OYTKZYDAYYDabrSIBLsBkg+cJhAcP74AnDB6cgWzKRqMRwP/aZQUCAadDeWyPxRMIDx7cgCcMHpyAFLbNZjMAWEUVXBIGFxs4TyA8eLgHnjB4eARS2DYajXQKynZT9XWE4Qo8gfDgwQw8YfBwG7azFc5aZh1t8mQOQyaTMTrnXGzQzgikpqYGCQkJiIqK4gmExxkHnjB4uAWy0TuKKiwhFArtEsbU1BRqamowOTkJuVyOiIgIhIeHIyIiAgEBAQ6PN9ezprZzIGRKnY9AeJxp4AmDBytYFrYpimI8iGe5yVMUhf7+fjQ1NSE5ORn5+fmYmpqCWq1Gb28vGhsbERgYaEUgUqnUmx+LFXgpEx5nKnjC4MEYjkQDXcGyhmE0GtHQ0AClUokVK1YgKioKer0e0dHRiI6OBgAYDAaMjY1hbGwM3d3daGhoQFBQEMLDw+niuj+Br4HwOFPAEwYPl7AnGshmoyOvHR8fR01NDeRyOdauXYuAgAC6q8oSEokEMTExiImJAfA/AlGr1TAYDKivr0dQUBAiIiLoKEQikXDzYTkATyA8Fip4wuDhFI5EA9lCq9Xi+PHjyMjIQHp6OqtjWBKISqWi369Wq9He3g6tVovg4GArAhGLvXNpu1M/4QmEx0IBTxg8HMJT0UAA0Ol0tJZUUVERwsPDPVqTQCCAWCxGVFQUYmNj6XOQCKStrQ1arRYhISE0eXiTQNwBTyA85iv85y7i4TdgMlvBBKOjo6itrUVgYCBCQkI8JgvL9VkiICAAcXFxiIuLA3CaQNRqNdRqNVpbWzEzMzOLQEQiESdr4QI8gfCYL+AJg4cVSArqxIkTWLp0KWQyGetNyWw2o7W1FT09PViyZAlEIhG6uro4WR+TtQQEBCA+Ph7x8fEAgJmZGZpAWlpaoNPpEBoaSndghYWFzRsC0Wq1OHXqFJYtWwaxWMwTCI85BU8YPGhYigYqlUorVzym0Gq1qKmpgdlsRmlpKYKDgzE8POzTSW+ZTIaEhAQkJCQAAKanp6FWqzE2Noampibo9XqEhobSNZDQ0FCnBDLXm7IlgRiNRqhUKn4OhIdPwBMGD4eigfY6mJxhcHAQDQ0NSExMRE5ODr3p+puWlFwuh1wuR2JiIiiKwvT0NF0DGRgYgNFonEUgXCruegLyPVoSGp/C4jFX4AnjDIcj0UBH09n2YDKZ0NTUhOHhYeTn59O1BAJ/1pISCAQIDAxEYGAgTSBarZYmkL6+PphMJoSFhSEiIgImk2nOJ83trdn2v/kaCI+5AE8YZyhcFbaZbvKTk5Oorq6GVCpFWVkZ5HK5w/NxAW9vcAKBAEFBQQgKCkJSUhJNIKQGMj09jcbGRgwMDNA1kJCQkDmLQJh8jzyB8PAWeMI4A8FENNAVYVAUhd7eXrS0tCAtLQ0ZGRkON01/jjBcwZJAFi1ahOPHj9O1ELVajZ6eHlAURZNHeHg4QkJCvLbxulNX4gmEB1fgCeMMg60bnqMNQSgUOqxh6PV6NDQ0YGxsDKtXr0ZkZKTTc3K56fh6AxMIBJDJZIiJiUFycjIoiqJ1sMbGxuhuMNK+GxERgeDgYL/6DngC4eEueMI4Q8BWNNBRVKBSqVBbW4uQkBCsXbuWkSigqwI6ORfTzcjXNQRLCAQChISEICQkBCkpKaAoCpOTk3QNpLOzEwKBgCaPiIgIBAUFub3xeuOz8wTCgyl4wjgD4I5ooC1hUBSF9vZ2dHZ2Ijs7GykpKYw3C1evcyfN4q8QCAQIDQ1FaGgoUlJSYDabaQJRKpXo6OiAUCi0IpDAwEBWZDkXdRxXBKJQKBAVFYXAwECeQM4g8ISxgOGJaKBlSmpmZgY1NTXQ6/UoLi5GaGgoq3Vw3VbrywiD7bmFQiHCwsIQFhaG1NRUmkDUajUUCgXa2togEomsdLBcEYgv50CA/z08BAQE0P/GRyBnBnjCWKDwVDSQbMwjIyOoq6tDbGwsVq9e7ZYmk683eX+CJYGkpaXBbDZjYmICarUaw8PDaG1thVgstvICkcvl9G/nD98j+T3FYjEkEgmfwjqDwBPGAgQXooEA0NvbC7Vajby8PCQmJrq9noUUYXANkp4KDw9Heno6TCYTTSBDQ0M4deoUpFKplQqvP2y6pA4G8DWQMwk8YSwgcCUaODU1BY1GA6PRiLKyMgQGBnq0roW2yXtzk7NMTwGnu9rGx8ehVqsxODiI8fFxAEBjYyP9OqZ+6FzCbDY7/B54Alm44AljgYDJbAWTYxDrVIlEgoyMDI/JAnBOGO7MFCwk8nEFkUiEyMhIunVZqVSisbERAQEB6O/vR3NzM2QymVUR3ZkfOldgU3znCWThgCeMeQ4SVSgUCoSEhLidsrC1TiXtoFyA65TUmQyhUAiRSISMjAwAp383Ymdr6YduOUjINYGQDd/d34InkPkLnjDmMchNptPpcPz4cZxzzjlu3VD2rFO7urrmrO7gTiR0psJ2oxaLxbP80EkKi/ihBwYGWnVhMZmdcbUGAJzJofAEMn/AE8Y8hW1hG2C/kVIUha6uLrS1tc2yTnU26e0O+AiDOzj7DiQSySwCsRwi1Gg0Hvuhk9/SW78FTyD+C54w5hmcFbbZbMo6nQ51dXXQaDQoLCyc5YbnjTQSV0NnZ3qEwQaWfujAaVkXQiAdHR3QaDS0Hzrp1nJFIN4mDFs4IpCpqSlUVFSgrKyMTtXxBOJd8IQxj+CssM0mIiDWqZGRkSgrK7O7Qfhr3eFM3wA8JV2pVIrY2FjaD12v19M6WO3t7bQfumUNxHb2huuUFFtYEgiZB+EjkLkBTxjzBK5EA5kQhq11alJSktPWSK5SUnyEwS243PSkUqldP/SxsTErP3RLO1tyXfh687W8F/gU1tyAJww/B1PRQFcbvD3rVGdgY6DkCv40pewpfP0ZvH1+R37oY2NjtB86uXbUajXCw8N95odOCMMWPIF4Dzxh+DHsueE5uqCdbfCOrFOdwVs1DK6OdaZiroUa7fmhj4yMYHJyEs3Nzaz90LmEs+FBS/AEwh14wvBDuCMaaC8l5co61Rm8lZLiAr5+yvc1fLmJyeVyREdHo7OzE2VlZXQEYs8PPTw8HGFhYV6rdZjNZrfIiScQ98EThp/BXdFA2w2eqXWqI/hrSupMm/S2hT98dsuHGLlcDrlcTvuhT09P0wTS398Po9GIsLAwugYSGhrKGYEwjTBcgScQ5uAJw4/giWgg2eDZWKc6g7+mpPwBvtwk/ME7xNEaBAIBAgMDERgYOMsPfWxsDH19fTCZTAgLC6NTWJ74oTuqYXgKnkAcgycMPwCZrSApKHd0oIRCIfR6PaqrqxlbpzrDXKak9Ho9JBIJ40hqoRCPO/AXwmCyUQsE1n7oFEVBo9HQBGLph25pZ8uUBLxFGLZgQiAqlYp2XlzIBMITho/BhWggcPrmaW5uRnh4OGPrVGfwxsZsezyz2YyWlhZ0d3dDIpHQT522HhA8/AvupoIEAgGCg4MRHBxs5YdOBgm7u7tpAiHXgTM/9LkiDFvYI5Curi6kpKRYEcVCjEB4wvARPHHDsz1Oe3s7pqamkJiYiPz8fE4uSq5rGLYEND09jZqaGphMJhQXF8NgMNAmQsQDIjIycpYCKx9h+EeEwVXtgDyVWxIIqYG48kP3FWHYgkTjEolkwRtK8YThA1jOVgDs3fAILK1Tw8LCEBkZydkF6M0IQ6FQoLa2FnFxccjNzYXZbKafLImJ0Pj4OFQqlZUCa2RkJGZmZtxy/Vso8BfC8FbtgBAI8UMnBGLPD31mZsbn3wWBZccWkxSWWq3Gli1b8N1333mcDZhLnLl3no9gNpsxMjICkUjEKl9rC1vrVDKUxxW4rGFYHq+lpQU9PT1YunQp3Vljex5bDwgioKdSqWg/7LGxMToCCQsL89nwmC/g602Sq+4kVxAKhQgNDUVoaKhdP/SxsTEAQH19PU0irvzQvQWTyeTwGrRHIBMTEzh69Oi8u255wpgjWIoGdnV1ITw8HKGhoayPYzKZcOrUKfT391tZp3KtLstlSoqgrq6O8aS5JWwF9AAgNDQUarUaTU1NMBgMCA0NpQnEk84bf4c/pON8FeXY+qF3dHRgcnISQUFBGBkZQVtbG8RisVUKa65qYWzSYwKBAFqt1mfk5gl4wpgD2Ctsu7O5T01NoaamBgKBYJZ1KtcpJC6PNzo6CrPZDJlMhoKCAo+fqkQiET19THr/VSoV1Gq1VecNIRDLvLen8PWG7S8pKV+vgaxDJpMhPT19lh+6ZS3MkkBkMhnnaycPg2yua61Wy+l1OVfgCcPLsCcaKBKJWBEGRf3POjU5ORnZ2dmznma4jjC4SElRFIW2tjZ0dXVBJBIhMzPTY7KwvcEse/9J66Zl3ru9vZ32ySYEwnaI0d/g603GX4rNtutw5Ic+NjaGwcFBtLS0QCqVWnmBcHEtWEr3MAWJMOYbeMLwEpyJBrLZ3G2tUy3TMpbwt5SUTqdDTU0NdDodSkpKcOLEiTl5OrdXOJ2YmIBKpaI3jYCAACsCmU9FR19HOGQNviYt4PRG7awBwrYWRgiETKE3NzfT14InfujkvmPzMKTRaPiUFI/TcCUayHRzt2ed6gjeKFK7uzkplUrU1NQgKioKq1atolsIfaElRbpqiEGU0WikO7CIhWlQUBBNHvb8H/wJ/rBZ+8MaAPbFd1sCIdeCWq2mu/Hkcjl9HTAlEMtUM1OQlNR8g//eGfMQlrMVjnwrgNMXlsFgcHocR9apjsB1kdodAiIzIZ2dncjNzcWiRYvodXNFGJ5uVGKxGFFRUYiKigJg7UBn6f9AIhB76qtnujSIv6ak2ML2WjAajfS10NPTQ7dzWxKIvWjU2b3uCCTCmG/gCYMjsJnYFgqF9Ots4co61RGcHdMdsCUgvV6PmpoaTE9Po7i4eFYHmCvCYLMRckmMtg50RH1VpVKhoaGBFs8jBOIPKSFfwx9IC+CeuMRisUM/dMto1NKNUCqV0tpvbMDXMM5gsBUNdFT0Hh0dRV1dHSIiIhxapzqCq6iFLdhEBCqVCjU1NYiIiMDKlSvtpnS4jDC8uWlb+j9YiuepVCr09PTAaDSis7MTWq3WJ33//rBZ+8MaAO9HOrbt3ESNYGxsDJ2dnbQfOtn4DQYD43uWT0mdgXBXNNA23cPGOtXZMee6rZaiKHR2dqK9vR05OTlITk52uu755odhTzzvyJEjCAwMhEKhoPv+ScE0MjISMpnMq2vyh816oaSk2EIikczyQycdWCaTCQcPHkRwcLBVF5ajehhPGGcYPBENtIww2FqnOoI3uqScHU+v16Ourg5TU1MoKipCWFiY0+P5Sw3D03OLRCLExcUhMjISZrOZLpoODAygpaUFMpmMJg+Sslho8AfSAnxPXCSdKRAIoNPpUFBQYOWHPj09TdfDiCIBIRA+JXWGgAvRQLIZu2Od6uqYXMHZBq9Wq1FTU4OwsDDGqTMuNxh/qSMIhUJ6M1i8eLFV0dQyZUEIxHLDcBf+sFn7wxoA3xOG5TpEIpFdP3RyPRA/9JmZGXz11VeYnJzEokWL3D7nDz/8gOeeew4VFRUYHBzE//3f/+GSSy5x+p79+/djx44daGhoQHJyMh588EHccMMNrM7LEwYLcCUaCJye2m5sbGRtneoIc5GSsuzeysrKQmpqKuPPP19qGJ7Atmiq1+vp+gfZMCz9r92xL/WHzdpb4oPzdR2Oit4ymcyKQKanp1FbW4vh4WEcPXoU4+PjOHHiBM4++2xccMEFKC0tZXxOjUaDgoICbNu2DZdeeqnL13d2duKCCy7AbbfdhnfffRf79u3DTTfdhISEBGzatInxeXnCYAjLwralkBhbTE5OorW1FSaTCevWreNs6tjbKSmDwYC6ujpMTEyw6t4i8OeNng3YfAapVIq4uDj6gcDWvtRkMtEdN5GRkU69Hyzha8KYK/FBJuvwB8Jg6i0ul8tRXFyMPXv24KqrrkJBQQEyMzPx/fff4/PPP2dFGD/5yU/wk5/8hPHrd+3ahfT0dPzpT38CACxZsgSHDh3Cn//8Z54wuISlaKA7/daWxyHWqbGxsZiYmOBUosKbKanx8XFUV1cjODgYZWVlbuXlz4QIwxVs/a+J+5xarUZXVxcEAoFVAd2ecJ4/fHZ/iHIA/yEMd9tqk5OTsW3bNmzbts1LK/sfysvLsXHjRqu/bdq0CXfeeSer4/CE4QQURUGv11tFFe7cKHq9Hg0NDbR1KkVRaGho4HSt3pIj7+rqQmtrKzIzM5GWlub2RuFoo/eHjccXEAis3ecspbuJ8qqlC2FkZCQCAgL8YrP2l1SQvxAG0wjDEnPdJTU0NDQr9R0XF4eJiQlMT08zfnjlCcMBSLvsd999hzVr1rjsAnIEUiAOCQmhrVPVajWnQ3YA95PeJKrq7OzEmjVraEE3d+EsMmCzAfp6s/QWbKW7LXWP+vr60NTURM98yOVyVj3/XIO4y/ka/kIYbJVqgdPpSb6tdgHAtrDNVlnW8jhEJiM7OxspKSn0ZufuMZ2By5TU+Pg4amtrQVEUJ/7gALepJF+nZeaCtCx1jzIyMuip446ODoyPj+PgwYNWLZvh4eFzZsbjD1EO4D+EwXYdZCB0LgkjPj4ew8PDVn8bHh5GaGgoq9Q4TxgWsBUNJGkottHAzMwMamtrodPp7MpkcF1vALiTIyd1lkWLFqGnp4ezOYKFMIfhS5CpY6VSCYlEgkWLFtH1j+bmZtqml6SvvGki5Q+EQdrbfb0O4HSEwTbiIi3Xc4XS0lJ88cUXVn/75ptvWBXaAZ4wADgXDWQbDVhapxKlVlt4gzA8TUkZjUbU19dDrVZj1apVkMvl6Onp4XCF82/S2x9BNmvLnn9iIkUIpK+vD2az2co4iGkHFtM1+PrJnlwD/mBx6gstqampKbS1tdH/3dnZierqakRGRiIlJQX33Xcf+vv78dZbbwEAbrvtNuzcuRO/+93vsG3bNnz33Xf44IMP8J///IfVec94wnA1sc00wiB+1bbWqfZANncun9Q8IaHJyUlUVVVBLpejrKwMAQEBmJmZ4XSNXNVY/OGJ0t9gaSKVlJRkZSJFhggthww9tS71hyd7yyyAr+GLovfJkydxzjnn0P+9Y8cOAMD111+Pf/7znxgcHLR64EtPT8d//vMf3HXXXfjLX/6CRYsW4e9//zurllrgDCcMs9kMvV7vtF2WSYThzDrVHsjTiDsXmrNjuiNHTpz80tPTkZGRYSVHTl7jbXtTtsfnIwzn35cjEylL61JPjIP8ISXljsudt8A2wtDr9TAajR6lpM4++2yn98E///lPu++pqqpy+5zAGUoYbEQDRSKRwwiDiXWqPXiLMNhspEajEY2NjRgdHcXKlSvp6WQCS8LgAs5qGGzO4euNytdk5c5mbWkiRbyviWQFMQ4KCgqyIhBnEiY8YViD7X2s0WgAgO+Smg8wm80wGo2MRQMdpaSYWqc6OibgXrHMEdgUvaemplBdXQ2JRIKysjK7Cqtck9pC1JKarxCJRFbGQUS2W61Wo7293a5onuU14A/dSSQt5mviAthHGFqtFgBPGH4Nd0UD7aWk2Fin2gMhqblUlyXo7+9HY2MjUlNTkZmZ6fBCn8sIg+1xzmR44+neVrabmEip1Wo0NTXBYDDQHVgRERF+U8PwNWkRsH2o0mq1kMlkflGwZ4szgjA8EQ20jDDcsU51dlyuCcNZkdpkMqGpqQnDw8OMIiJ/JQzgzI4w5iId5MhEiliXGo1G9Pb2wmAwICIiAkFBQXNOIP5EGGwH94g9q69J1x0seMJg64ZnCxJhuGud6gje0H4C7G8oJAUlFouxdu1aRiY//koY8/Em4xpz+R3YM5E6duwYgoKCoFQq0d7ebmUiRTqwvA1/aO0lYJuSmq9+3sACJgyuRAOFQiGmpqZw5MgRt6xTnR2X6wgDmP3kNTAwgIaGBqSkpCArK4vxhU2iMK7WyEcY3MDXn51cF7GxsYiOjrYykRocHERLSwsCAgIQGRlJE4g3TKT8KcJwJyXFRxh+BE/c8CxhNpuhVqsxPj6OpUuXumWd6gjeJAzg9FNPc3MzhoaGUFBQQOen2YDLTZ6PMLiBP3QoWa7Bcr4DgJWJVHd3NxoaGhjblrKBvxAGqY2yLXrzEYafgLTLehJVAP+zTp2ZmUFUVJRH7lj24I7kiDOQz2k2m6HRaKzmQtxNEcwFYZhMJpw6dQpGoxFRUVGIiIhwGcH5w1P2mXx+ZxukIxMptVqN1tZWzMzMICQkhI5A3DGRcrWGuQR5QGNbw5iPHVLAAiIM28K2J2RhaZ0aHx8PlUrF5VIBcK8uSz7ryMgITp06haSkJOTk5Hh0U3EZBdkjDK1Wi6qqKgiFQoSEhKCzsxP19fX0hkJsTS0/g683S1/D12RJ1sD0d3BmIjUwMACj0YiwsDCaQEJCQhgd218IwzKLwRTzVakWWCCEQVJHEokEEonE7f5sy04iYp1KdHm4BtcRBtlIWlpasHz5cr+0fbU81sjICGpra5GUlISMjAx6A9DpdFCpVFCr1WhoaIDRaLTyhCCdYGcq/C0lxRbOTKS6u7sBwKqA7ijX7w+tvWQdADvC4IvePoLlbEVtbS3S09ORkJDg1rEmJydRXV0NqVRqlcbhemMn4FLinKTPALAeInQGb6SkKIpCW1sburq6sGzZMiQkJMBkMtHfRUBAgFVLp0ajoQmko6ODjhwHBwdpU6EzDb7eKLnqULJnIjU1NQWVSgWFQoG2tja6A4tEIKTDj0uVBE9AWmrZ/CZzLW3OJeYtYdgWtj3xrSCS3mlpacjIyLC6GbzhXQFwl+4ZHh5GXV0dEhMTodFoON1Aue6SMhgMOHnyJKanp1FSUoKQkBCX7yEbCtFE6uzsxNDQEG0qFBQURKev5sITwtfRja/PD3jv6V4oFCI0NBShoaGzTKT6+/vR3NwMmUyGyMhImEwmv/ku2JInX8OYY9gTDRSLxawjAVvr1MjIyFmv8VaE4SlhmM1mnDp1Cn19fVi6dCkSEhIwPDzslWFALqDX6zE8PIzo6GiUlZW51SkjFAoRHBwMmUyG1atX05IWKpUKLS0t0Ol0dD6ceEL4+mmca8z3lBQbWJpIAac7sEj6SqlUwmAw4Pjx43T0MZcmUgTuuO1ptdo59cLgEvOKMJyJBjoTCbQHe9ap9uCPEcb09DSqq6thNptRWlpKP614YxjQU8KgKAp9fX0YHBxEaGgoVqxY4fFmQ9ZkKWlBPCFI+opIO5N0RmRk5JwMlM0FzhTCsIVYLEZMTAxiYmIgkUgwNTWFmJgYqNVq+oEhNDSU/s1DQ0O9Xhh31wvDnTZ3f8C8IQxXooFMIwFn1qn2wJaImMLdzZ0YNMXFxWHJkiVWTzdcF6k9TUmZTCY0NjZCoVAgISHBYecam0K2o9/K0hNi0aJFMJvNmJychEqlwtDQEE6dOkWnM0hB1R98qdnC1xEG+a183aFEfMWJiRQAqweG/v5+mM1mqw4sLk2kLNcx114YvsS8IQyTyQSTyeSwA4pJSsqVdao9+EtKymw2o7W1FT09PVi6dKldgyZv6VO5A61Wi+rqagiFQpSVlaGvrw/T09OcrIvJmoRCIcLCwhAWFob09HR6oEylUqGjo4NWZHVnHsDXG7YvQc7v6yjHXu1ALpcjKSmJNpGybJggJlLh4eH0b+6JiRSBOxHG9PQ0n5LyNghRuONbATCzTnV0XF+31c7MzKCmpgYGgwGlpaUOLzZ/SUkpFArU1tYiISEBubm5dGThy0lv24EyosiqUqnop1GymURGRvq1dIM/EJavvxtXT/b2GiZIxElMpKRSqVUHljsNI+5EGHxb7RzA1WyFSCSCXq+f9Xc21qn2QJ6yuR4UYprqGh0dRU1NDWJjY5GXl+fyJuGaMNgcz7Jl1jYK8jctKVtFVmJpSgT1JBKJVf3DG3pI7sAfUlKA742LSEqKKWwjTtKBpVKpaBOpwMBA+jcPDw9ndHx3axh8hOFluLpJ7G3AbK1THR0X4H6yVCgUwmAwOPx3y813yZIljKRJuJ4eZ3M8vV6P2tpaaLVauy2zvo4wXB3T0tLU3mYSHByMyMhImM1mr0ScTOHrlJS/eGl7ej/admAZDAY6Zdne3g6tVkubSBHFAXsPa+52SfERho9hSRjuWqc6Oi5w+sLgQjSNwFn6SKfToaamBjqdjtG8ApNjugOmm/z4+DiqqqoQGhqK0tJSu09m/hZhOIPtZkL0kFQqFQwGA6qrq63SV94opjqDP0QYviYMrgvvEomE7sACYKU40NTUBL1eT5tIkZZtcr+xWQfxF+GL3j4GIQxPrFPtwVLUj0s42tyVSiVqamoQFRXFqtbi7JjugklKigzQMTGU8tcIwxUs9ZCUSiWysrJgMBigUqnQ1dVFK7YSAmHiN+Iu/CUl5WvC8LaWlK3iANHAIlEnRVEIDw+HyWSCVCpl9bvwhDEHYJKS0ul0OHLkiNvWqY7O643WWtvN3bLdNzc3F4sWLWJ9U3qjhuFokye6WyMjI1i1ahXtD+3Osdh+Tl+mZQQCAWQyGeLi4ma17xI/CLlcblVM5TIy9YeUlLtabVyvY67qKJYt26QDi9S8ent7MT4+jkOHDllpnjmb+eHFB30MiqIwOjoKjUaD7Oxsl0+6tbW12L17N3p7e5GXl4ebb74ZqampDl/vjdZaS8LQ6XSora3F9PQ043ZfR8ecixoGaZkVCAQoLS1lNBDH1Qbj643KFvbad8k0cnt7O6anpzkdJvOHCMMffgNfig9a1rympqbo+R7LmZ+AgAArAiFNE3xKyscg1qkTExMICAjA4sWLnb5+z5492Lp1K0QiEYxGI77++mu88sor+OKLL1BaWmr3Pd5orSUkpFKpUFNTg4iICKxcudKjp9G5qGHYa5l191juwpdP2a7ObTmNDJxu31WpVFbtu5bRhzvtuzxh+I/4oNlshlgsRnh4OG3bbDQa6aaJnp4eNDY2IigoCA0NDRCLxZDJZB53Sb3yyit47rnnaIO0l19+GUVFRQ5f/+KLL+K1115DT08PoqOjcfnll+Opp55inT6d14QxOjqKuro6REREoKCggFZsdYSJiQncfvvttHcGcDq1otPpcMstt6C2ttbuzeCtCEOr1aKiogI5OTlITk72+Eb0Zg3DMmWWl5eHpKQk1sfiIiXlD5sVG8hkMiQmJtJy3rZqrBKJhK59MLEz9fWG7Q9T3oD/yJvba6sVi8WIioqi07R6vR5jY2P4/PPP8a9//QszMzO46KKLsHHjRvzoRz9inT7fu3cvduzYgV27dqG4uBgvvvgiNm3ahJaWFruSI++99x5+//vf44033kBZWRlOnTqFG264AQKBAC+88AKrz+v7X54hLC8OMltRVVWFrKwsFBQUICAgwOWm/u2339qdNiZT1M3NzXbfx3WEodfr0dnZCb1ej6KiIpfyJEzBtTQISUkZDAZUVlZiYGAAJSUlrMmC67X5Oo/vLkgqIzU1FStXrsT69euRm5sLiUSC7u5uHDp0CMePH0dbWxtUKpXd69nXn91fNmp/MVBiEulIpVLExsbisccew5dffgmhUIhbb70VPT09uPbaa9He3s7qnC+88AJuvvlmbN26FXl5edi1axcCAwPxxhtv2H39kSNHsHbtWlxzzTVIS0vDeeedh6uvvhrHjx9ndV5gHkYYxPvBZDJZTT2TwrSzJzB7g31M/p3LCIOIHspkMsjlcoSFhXFyXMA7EcbMzAyOHDmC4OBghy2zTI/lbLNjugn5w2bF1RpEItGsJ1GSvmpqaoLBYLBS3yXXuq8jDH/4DfyFMNjOYWi1Wsjlctxwww3Ytm0b6wcAvV6PiooK3HffffTfhEIhNm7ciPLycrvvKSsrwzvvvIPjx4+jqKgIHR0d+OKLL3DdddexOjcwzwhjaGgI9fX1SExMRE5OjtUPZTlg5+gHPPvssx12PBEhM3vgIsKgKApdXV1oa2tDVlYWQkJCUF9f79ExbcF16kyj0WBsbAyZmZlYvHixRxsFH2G4hlQqpcX0SHGUEAhp3zUajRgdHYVUKvVq+64j+FNKaj6ug8iCkHuJ7T01OjoKk8k0y1EzLi7OYYbkmmuuwejoKNatW0en42+77Tbcf//9rM4NzKOUlMlkQmdnJ/Lz8+1KZFgO2DlCfHw8zcy2P5TBYMBZZ51ldxP3tK3WYDCgqqoK3d3dKCwsRFpamlcK6VxtymazGfX19RgfH0dsbCwyMjI8fqr050lvf4RAIEBQUBCSk5NRUFCA9evXIz8/H8DpxoPy8nIcPXoUp06dwujoKF2T8zb4lJQ12EYY09PTcz7lvX//fjz55JN49dVXUVlZiY8++gj/+c9/8Nhjj7E+1ryJMEQikcMuJuB/2jauNvYHHngAOTk5uO+++9DX14eYmBgoFAoAQFdXF8466yzs3r0bl1xyidWx3SWM8fFxVFdXIzg4GGVlZXRRk+v0EVfHJF4bAJCYmMhZJwofYXgGorQqEomwdOlSBAQE0FIWra2tmJmZQWhoqJV5lDc2VD4lZQ22WlK2EQZbREdHQyQSYXh42Orvw8PDDjMkDz30EK677jrcdNNNAID8/HxoNBrccssteOCBB1it3/ffOAs4+2BMB+wEAgGuuOIK/PKXvwQAnHfeeTTTCgQCaDQaXHXVVXjsscfozdedaICkoI4fP46UlBSsWrXKqgPGHwljdHQUR44cQWhoKIqLiyGRSDj39ObiOGcyyIZNpCxycnJQWlqKkpISxMfHQ6PRoKamBocOHUJtbS36+vqg1Wo5JWt/2Kj9ZR1s23s9ncGQSqVYvXo19u3bZ7WGffv2OXyg1mq1s74rsma218W8iTCYgE3qiISF09PTuOeee3D8+HF89tlnCA4OxtTUFJ544gnU1dVh9+7drCMMg8GA+vp6jI2NYc2aNYiIiJj1Gm+06ro76U1RFDo6OtDR0WEldOhKINGd8/jTceYjHH12Wy8IMn0+MjKC1tZWepCMRCDuNi/4Q0rKbDb7DWG4E2F4OrS3Y8cOXH/99VizZg2Kiorw4osvQqPRYOvWrQCALVu2ICkpCU899RQA4KKLLsILL7yAlStXori4GG1tbXjooYdw0UUXsc4gnLGEQaaTtVotBAIB/v73v6OsrAzt7e3Iy8tDW1sbPv30U5x11ll4+umnkZKSwui44+PjqKmpQWBgoFPrV3KRcRlauzPpbTAYUFtbi6mpqVlT5lymkfgIgzu4+g4EAgFCQ0MRGhqKtLQ0mEwmOn3V1dWFhoYGK/Oo8PBwxtegP6Sk/ElinaIo1hGGpzWMzZs3Q6FQ4OGHH8bQ0BBWrFiBr776ii6E9/T0WH03Dz74IAQCAR588EH09/cjJiYGF110EZ544gnW5z5jCYP8aFqtFgAQFhaGvXv3Yv369WhsbMSWLVvwzTffoLGxEddeey2efvpp5ObmOjweRVHo7e1FS0sLFi9e7LKryBuy6WxTUhMTE6iqqqLrK7ZPnf5IGAAfYbDdsG3bd3U6HS2k19jYCKPRaKW+GxQU5PAc/kAY5Br3B8Jguw6upM23b9+O7du32/23/fv3W/23WCzGI488gkceecTj884rwnDHE8MRLFNSBMuWLcOrr76KrVu34q233sLf//53/O1vf8Px48exfft2aDQabN++fdY6jEYj6uvroVarGQnxAdYRBldgQxj9/f1obGx0Sm5c1lkWipaUr8mKi/MHBARYte9qNBqaQDo7OyESiazSV5ZTyP5QbPY3wmATYXCRkvIl5hVhuALRh2ICkpKynfy++uqrcezYMezatQv33HMPvv/+e/zhD3/Axx9/jN/+9reoqanBzp076R74yclJVFVVQS6Xo6ysjPGIvzcIg8lTvNlsRlNTEx3KOpN/5yMM/4I3pMUtrUyTk5NhNptpHSTiKRMUFESnr0wmk89J218Uc8l3wWYdvmir5RILjjDcqWHY4tlnn0VVVRWOHTuG66+/Hm+++SZSU1Oxc+dOvPPOO2hpacGePXtgNpvR3NyM9PR01rMK5EKbywiDtMxSFIXS0lKXF64/EoavNwl/gDe/A+LtQRo1DAYDHX2cOnUKMzMzEIvF6Orqott35/o38YcoB/hfwZvN59doNLSv/HyE7791FmCSkmK6AdtLSRFIpVK8++67iImJQW1tLR577DFcdtll+PTTTxEREYETJ06guLgYn332GVauXInMzEy3bhquW2udHU+pVOLIkSMICQlBcXExo6ccfyQMgI8w5hISiQSxsbHIzc1FWVkZUlNTIZfLMTk5ierqahw8eBB1dXXo7++3ey95A/7QqUXW4Y4963z18wYWYITBNCVlW/S2xaJFi/DWW2/hggsuwEcffYS0tDQ8+eST+Prrr3HllVeiq6sLv//97xEcHIwtW7a4tV6uW2vtbcqOWmaZro/LGsZCiTB8tQZ/cLsTiUQIDAzE0qVLrdp3h4eHaR8IS/Vdd9t3ncHfIgw2mM9+3sA8izBcgauUFME555yDP/7xjwCAl19+GV988QWGh4fx/vvv4+KLL4Zer8ctt9yCe+65xy1pBq7lQWw3eCJJ0tfXh6KiIlZkAcxNhEFRFNra2nD8+HG0tLRAoVC4/C7P1AiDwF/EB0n7blpaGlatWoX169cjOzsbIpEInZ2dOHjwIE6ePIn29nao1WrOrnV/IQx3I4z5TBjzKsJgkpJiG2EYDAYYDAaHT0J33303Dh48iK+//hq33norvv32W+Tk5GDPnj148skn8fjjj2Pnzp1oaGjAO++8w6hDisCbKSlSjA8MDERpaalLnwV78DZhGAwG1NTUQKvVIjU1FRMTE2hra8PMzAyt0hoVFYXg4GC3xdoWEvyBKJ211YrFYkRHR9M5ep1OR4snNjQ0wGQy0e27ERERTtt3ncFfCMPdCIPvkvITEF9vJrBk+enpaYeEodFocOutt6Kurg4DAwO466678Nlnn0EkEuHBBx9Efn4+tm3bhu+//x7r1q3Dhx9+iKVLlzJag7cIY2BgAA0NDUhLS3O7vsL1+mwJY2pqCpWVlQgKCkJxcTEoikJCQgKA07+HUqmESqVCd3c3hEIhTR6BgYF+sXH6Av6QkmKzWQcEBCAhIQEJCQl0+65KpYJSqUR7ezvEYrFV+opph6G/EIa7EcZ8Jgzff+scgk1KKiAggL7xHKWlBgYGUF5ejkWLFuGBBx5AYGAgvvvuOzpNBQA/+9nPcODAAaSlpaGzsxMbNmzAJ598wmgNXBMGkS5uamrCihUrkJWV5TeS5GR9wGmhtPLyciQkJGDVqlWzyFoul2PRokVYvnw51q9fj2XLlkEmk6G3txcnT54EAM7THPMJ/pKSYgPSvpuSkoIVK1Zg/fr1yMvLQ0BAAHp7e3H48GEcO3YMra2tUCqVTu9jfyEMtkq1FEVhenp6XhPGgoswmBKGQCBAYGAgNBrNLMIwmUxobm6m/XLlcjmGh4fx6quv4oYbbsAzzzyDwsJCXHjhhQBOD/wdPnwY1157Lb7//nts3rwZDz74IO6//36nFzaXhDEzM4OGhgaYzWasW7eOkzwp1ykp4mzY1dWF/Px8Wl3T2Tks2zwzMjIwMTGBkydP0p/XZDJZDZnN5/ywK/hDZMXVpLdIJKJ/s4yMDBgMBjp91dLSAp1OZ2UeZdm+6y+E4c46NBoN3yU1V+By0hsATRiW7YAajQbV1dUQCoUoKyuDXC6HVquF2WzGVVddhePHj+PVV1/FjTfeiCNHjiAjIwMAEBUVhc8++wz33XcfXn75ZTz++OOoq6vD3//+d4SEhNg9P1eEoVQqUVNTg4iICExOTnK2aXI5J0J+l8HBQZSUlDj8TlyBpC3y8vIAgPbIJiJ7MpnMKs0hFs+rS9wp/CElxVY7iSkkEgni4uIQFxdHP4kTAunu7oZAIKB/U4PB4BeEwTbCAPiit1+BLWHYTnsTR7+kpCTk5OTQFyUR9TObzXj66adRWVmJo0eP4qqrrsKBAwfoC0AsFuO5555Dfn4+tm/fjk8++QRtbW3417/+hcWLF886v6dttRRFobOzE+3t7cjNzUV0dDSGh4c5ewp0R8zQHqamplBVVQUAKCkpcasAbwuiVhoSEkL7ZBuNRlpkz1Xx3JPz+hq+rmF4+/wk+g8MDMSiRYtgNpvp9t2hoSGMj49DKBSipaUFkZGRCA8P90r7riu4W/TmIww/gbuEodFo0NjYiIGBAeTn58+yP7QUCiRDfaWlpairq8P27duxe/duq5toy5YtyMnJwVVXXYWGhgasW7cO77zzDs4999xZx3X3Cd5oNKKurg7j4+MoKipCWFgYXfDnijC4SEmNjIygtrYWiYmJ6Onp8fjGdva5bLt0yFOqUqm0Kp4TAuGCuOYS/kBWvhAfFAqFCAsLQ1hYGNLT09HV1QWlUgmBQID29nZMT0/T6ruRkZEIDQ2dkwiEbdHbaDRCr9fzEcZcwdWFKhaLWaekAKCyshLFxcUoKyuz+2Na2r+KxWIkJSXh7bffxk9/+lO89957KCkpwS233GL1nuLiYhw+fBibN2/GiRMncNFFF+GZZ57BHXfcYdUi6g5h2OpXWbr4AdzleD1JSVEUhfb2dtpWNyIiAj09PZxtOEw2T0uPCEuNpL6+PjQ1NSE4OJgmj7CwML9IcziDv6Sk/OF7ksvlyM7OBnC6fqdSqaBWq1FXVwez2WylvuuJw50zmEwmVg9AU1NTAMBHGP4Ctikekt8Wi8UoKSlxeCNYFtsIzjrrLDz22GO4//77cffdd2PFihUoKiqyel9iYiK++eYbbN++He+88w7uuece1NTU4OWXX4ZMJnMrwnDWMsu1oKG7KSmj0Yja2lpMTk7S9Qq9Xg/A86dkT1qELYvner1+1oyAvxfP/UVa3N/WIJPJkJiYiMTERFAURde1RkdH0d7eDolEYvXbchVZuuO2B/CEMadwliZhmpIym804deoUPeQXGRnplv3rXXfdhePHj+Pjjz/G1VdfjaNHj85Sf5XJZHj99ddRUFCAe++9F2+//TYtXsim6G02m9HS0oL+/n4UFBQgNjbW7joB7lIX7qSkNBoNKisrIZPJrAYGXa2N7Sbk6WeUSqVWEt9kk1EoFLRDXVRUlN3iuS83TF9v1v5CWo42aoFAYFXXMplMdGTZ09ODxsZGOrKMjIxEWFiY20V8tjWM6elpSKXSed2IMX9XbgdisZh2wXJ0URPFVrPZjISEBFRVVTmVByGwF70IBAL87W9/Q0NDA1pbW7FlyxZ8/vnnsy5AgUCAX/3qV8jLy8O1116L48ePY+3atXj66adRUFDg8twzMzOorq6GyWRymDYj5wG4izDYEgapVyQnJyM7O9vqN3BFGEzTLd7YrOxtMkSh1bZ4TlGUz2oJ/lLD8HVKik2UY9m+CwB6vZ7+bZuammAwGKzad9k0RrCNMDQajdfSY3OFBUUYtrUGW4yMjKCurg5xcXFYsmQJ3drJRGXTUfooNDQUe/fuxbp16/D999/j0UcfxWOPPWb3GD/60Y9w6NAhXHHFFWhsbMRNN92EBx54AL///e8dnlelUqG6uhrR0dFYunSp0wtUIBBwPp3N5FiWAofLli2jJ7Ztj0VeywW8uXmKRCKHxXOTyYSqqio6+pjL4rm/PN37wxrcbZ6QSqVW7btarZYmkK6uLjp1SQiE+N7Ygzt+3v6Y6mSDM4IwyMBYT08Pli5disTERACuFWst4aw+kpeXh127dmHLli147rnnUFRUhIsuusjuazMyMnDgwAHceOON+PTTT/Hoo49idHQUTz/9tNWaKYpCV1cX2trakJOTg+TkZEY3KpfDdkxqGKRba2JiYpYnuO26AN/VMDyBZfH8wIEDyMrKglar9Unx3Nebtb+QFldNHUFBQQgKCqLbdycmJqBSqTA4OIiWlhbI5XI6LWmbmmQ7hzHfZUGAeUgYzjZEYkpkubHPzMygpqYGBoMBpaWlVgUnNoThqkB95ZVX4tixY3jllVfoob7MzEy7rw0JCcGePXtwzz334NVXX50lXmjZMltYWIjw8HCX6yPwhv6To01Co9GgqqoKAQEBjAUO50OE4QwkfZWQkOCweB4eHk5HIHK5nLMN1l9SUguFMGwhFAoRHh5O329Go5GOPkj7bmhoKB19uEMYXF4PvsC8IwxnsC1Oj46OoqamBrGxscjLy5v14zqyabUHJh1YTz31FCorK1FeXk4P9Tl6ohAKhfjVr36F5ORkPP7447R44VtvvYWZmRnIZDKrllmm8IYPt71NQqFQoKamBosWLUJ2drbLG3g+RxjOYFs812g0UCqVjIrnbLGQN2t/XINYLEZMTAzdyDI9PU0TSF9fHwwGA7q6ujA9Pc3o4WC++3kDC4wwgP9JnBPNImemQc5c9+wd1xVhWA711dfXY/v27XjjjTccXkQikQgbNmzAgQMHcPnll6OzsxM//vGP8dRTT+G2225zW+SNy7ZawHqTt6xXWKb3mKzL9liewB+etm1BBPaCg4NnFc/JE6ojfSQ25/AlzmTSksvlkMvldPvuwYMHERwcTD8cSKVSK1ka24e9hVDD8P0EDku4uliFQiEaGxsxNDSEkpISp6ZBTEyUCJjOTCQmJuLtt9+GSCTC+++/j7/97W9O12o2m5GXl4d//OMfWL58OWZmZnDXXXfhiSeecGvj50rOA5i9yRuNRlRXV6O3txfFxcWMycL2eEz/7ulrvQE23y0pnmdnZ6OkpAQlJSWIi4ujBy8PHTqEhoYGDA4OMpLl9weSPJMJwxIkXZuUlISVK1diw4YNyM3NhVgsRnd3Nw4dOoQTJ06gra0NKpUKMzMznCjVvvLKK0hLS4NMJkNxcTGOHz/u9PVjY2O44447kJCQgICAAGRnZ+OLL75w+/zzjjCcQalU0ib1paWlLgXu2EQYbIYCN2zYgCeeeAIAcM899+DYsWMOj2k0GnHixAmYzWZ6yA8AHn/8cVx99dWYnJxkdE7LY3KdkjKbzdBqtTh69CgMBgPKysocFrddHY+rOYz5CFI8z8/Px/r165Gfnw+5XI6+vj4cPnwYx48fpzcYe7+hP2zW/tJW6+s1kHWQNLdIJEJUVBSysrJQVFSEdevWITk5GXq9HseOHUNKSgree+899Pf30239bLF3717s2LEDjzzyCCorK1FQUIBNmzZhZGTE7uv1ej1+/OMfo6urCx9++CFaWlrw+uuvIykpye3P7PtvnQMQm08yMJacnMwoV+yNCIPgN7/5DS699FIYDAZcc801dn9UrVaLiYkJyOVylJSUICwsDM8//zz++te/QiqV4pNPPsHZZ5+Njo4OxuflkjDITTk6Oory8nJERUVhzZo1breRctHBxXVqy1cgBdbFixejsLAQ69atQ0pKCvR6PRobG3Hw4EHU1NSgr68PWq3WL2RBAP9oq/UX0nKm3EtqW3l5edi0aRO+/vprZGZmYmRkBOvXr0dCQoLTdnp7eOGFF3DzzTdj69atdGdmYGAg3njjDbuvf+ONN6BSqfDxxx9j7dq1SEtLw1lnncVo9ssR5j1h6HQ6nDx5EgMDAyguLkZgYCDjDZOrtlp7EAgE+Otf/4rs7Gz09/fjuuuuoyfLSctsR0cHpFIp8vPzrS6866+/Ht988w3i4+Np8cLvvvuO8Xm53kzr6+uRm5uLJUuWeHSjcrm2+U4YtrDcYNauXYvVq1cjPDwcCoUCx44dQ3l5OXp6emAymdzyj+cK/hDl+EOEQfYYJusQCoVYuXIl4uLicOmll0KpVOLDDz/EunXrGJ9Pr9ejoqICGzdutDruxo0bUV5ebvc9n376KUpLS3HHHXcgLi4Oy5Ytw5NPPumRQva8IwzLi1WlUuHIkSOQSCR0moSNYi3bojfbJ/eQkBDs3bsXQUFBOHDgAB599FEYjUbU1NSgq6sL2dnZEIlEdm9AIl5YWFgIlUqFiy66CDt37nS5UXIVYRA9KAAoKCjwKIwl4IowfL1hzYW8Nymck/w4EdozGo04ePAgKioq0NXVhYmJiTklT395uvf1Gsgew2YdpK1WKpVi/fr1tAEbE4yOjsJkMs1S0o6Li8PQ0JDd93R0dODDDz+EyWTCF198gYceegh/+tOf8PjjjzM+ry3mZZeUZaeO7VAbG8Jgk5Ii9Qa2WLJkCXbt2oXrrrsOzz//PIKDg3HuueeitLQUGo3G6eaelJSEb775BnfccQfefffdWeKFjtbpKWFotVpUVVVBLBbTmxcX4CMM90CK50KhEFqtFgUFBVbmQpay7ZGRkYy9sd2BP6Sk/IUwiLICU8z14J7ZbEZsbCz+9re/QSQSYfXq1ejv78dzzz2HRx55xK1jzrsIg4RmfX19KCoqQkpKitUF7E6EwbSG4W4od8UVV2Dbtm0AgOeffx7h4eEICAhgFLXIZDL8/e9/x3PPPQehUIi3334b5513HgYGBuy+3lPCIPWKiIgIFBYWct51tRAiDF+BpIM8LZ5zsQZfwh8Ig62OFOAZYURHR0MkEmF4eNjq78PDw7TVsS0SEhLoLAbBkiVLMDQ0RKtHs8W8IwyVSgWRSISysjKEhYXN+nd/SkkBpy+s5uZmXHzxxSgsLMTU1BSuvvpqaDQaxps7ES/89NNPERERQYsX2mupc3dTJu59VVVVyMnJQV5eHoRCIee+3nyE4T7sfWZ7xfPU1FS6eP7DDz+gpqYGvb29VsVzT9bgD4Th6zW467bnLmFIpVKsXr0a+/bto/9mNpuxb98+lJaW2n3P2rVr0dbWZrXHnDp1CgkJCW43rsw7wkhISMCKFSscio95K8Jwx05Vp9PhxIkTGB0dxfr16/HBBx/Qhexf/vKXrIfsNm7ciEOHDmHJkiUYHBzExo0b8fbbb89aJ1tiM5lMqK2tRVdXF4qKiqxmV7jsugK42eh9vVn4Eq4+OxHXI8XzNWvWICIiAqOjo3TxvLm5GQqFwq0UK1/D+N8a3IkwPEnv7tixA6+//jrefPNNNDU14fbbb4dGo8HWrVsBnHb6vO++++jX33777VCpVPjNb36DU6dO4T//+Q+efPJJ3HHHHW6vYV7WMJzdNCKRCDMzM4yOw0YahG2EoVarUV1djcjISKxevRpisRhBQUF45513sGnTJuzduxerVq1CTk4Oq6e2jIwM/PDDD9i2bRs+++wz3HzzzaitrcVTTz0FsVjMeoOfnp5GZWUlxGIxysrKZuW/+QjDf8D26d5y8jwlJQUmkwljY2NQKpWztJGioqJcTp470xWbK5A1+Jow3I0wPJn03rx5MxQKBR5++GEMDQ1hxYoV+Oqrr+hCeE9Pj9WakpOT8fXXX+Ouu+7C8uXLkZSUhN/85je499573V7DvCMMVxerO0Vvo9EIg8HgVDKZaYRBURS6u7vR2tqK7OzsWTWWdevW4cknn8S9996LBx54AI8//jjOO+88Vjch6b56/PHH8eSTT+Lll19GfX093nnnHVZRi1KpRHV1NRISEpCbm2v3BvBHwjhTIwxPvzsyXBYVFQXgf7LtKpUKvb29EAgETovn/jALwqad1dvrYBNhECl1TxtItm/fTg/32mL//v2z/lZaWoqjR496dE5LzLuUlCu4k5ICXKelmEQYpBW1s7MTa9asQWpqqt2b69e//jUuu+wyGI1GPPvssxgcHGS0XksIhUI8/PDDeP/99xEUFESLF3Z0dLjcWMgcSGVlpVW9wtF5fGXI5AxnYoQBcLtZWxbP161bRxfP+/v77RbPyXfuy83aXwjDFxGGP+CMJgypVEr/6K4Iw1WEodFocPToUeh0OpSVlSEiIsLhawUCAXbt2oWcnByoVCrccMMNbg9j/fznP8f+/fuRmpqKzs5O3HDDDfj2228dvt5kMqGurg6dnZ0oLCx0qrVF1sq1NpW/HMcd+JKovJkOsiyer1mzBuvXr0dqaioMBgNdPCdzOdPT0z77HvyJMOayS8pfcEYThkAgYFz4dnbc4eFhlJeXIyYmBmvWrGHUB088MWQyGQ4ePOh2XzQA5Ofn4/DhwzjrrLMwPT2NHTt22BUvnJ6exrFjx6DValFWVsbIZ2OuUlJnapqJDeZyk5ZIJLQzJSmek67EiooKj4vn7oJc076+XtxNSfGEMcfgsoYBMG+ttZeSMpvNaGlpQV1dHZYtW4acnBxWTz5LlizBr3/9awDAn/70J3z88ceM32uL6OhofP7557jmmmsAAI899hiuueYaTE1NATjdjlxeXo6wsDAUFRUxHu6aC8LQ6XRobW3FwMAAo4YFb8ifzAf4quBMiudk2n/9+vX0td7e3k5Pnnd2dnp98px0SPmaMNimpGZmZmA2mzkbgvUV5l3R2xXYEgbTaW/blJROp0NNTQ30ej1KSkrcvhA2bNiAyclJvPbaa7j55puRl5dHy0CwhUQiwf3334+MjAw888wz+Pjjj9Ha2oq//OUvmJ6eRm5uLpKTk1kd09s1jImJCVRWViIoKAgTExM4deoUAgMD6eKsty1P5xv8oeAsFoutiuczMzNQKpVWxfOIiAjaOIrLyXN/6JAC2EcYZH+Z7xHGvCQMZ0+Y7kYYbIrepGU2IiICq1at8shFTSgU4oEHHkBtbS0OHz6Mq666CgcPHnT7whIKhbjwwgvxox/9CJs3b0ZDQwN+/vOfY/fu3azJAvBuhDE0NIS6ujosXryYrqUQW0ylUklbnpLNJyoqipZE8WWE4atN29dRlaMIRyaT0Z7nZrMZk5OTUCqV6O/vR1NTE4KCgmjyCA8P92jD94cZDIB9hKHRaKxS4PMV85IwnIE47jEF05QUkcjo6upCa2srsrKyHHZBsYFQKIRIJMI777yD0tJSNDY24vbbb8ebb77ptuMeRVFYsWIFXn75ZTz88MNobm7G1VdfjWeffZYeGGR7PC5g6RHe3t6Ozs5OFBQUICYmBgaDARRFQSKRIDY2FrGxsaAoClNTU1AqlRgaGqKjD5PJhImJCYSEhPjF5jFX8IcZCFfft1AoRFhYGMLCwrB48WIYDAa6dbepqQkGgwERERF0625gYCCrz+QvhGE2m1lNS5MOKV+n0jzFgiQMiqIYX1hsBAiB0wqQZHqWC5CUT0JCAt59911s2rQJH3zwAYqLi92ayBQKhZiZmcGRI0ewePFiHDlyBNu3b8d7772Hu+++mxYvZJom4HrS22QyoaamBuPj4ygpKUFISIjTQnhISAhCQkKQlpYGg8EAtVqNxsZGWh7eXvSxkOFrwmB7flI8j4uLoz3PVSoVRkdH0d7ebmVrGhkZ6TJa9xfCcCfCcOX5PR/g+2/eDbia9AbAeJNjEmFoNBqcPHkSADglC8B6Q167di2eeuopAMC9996LI0eOsDoWRVEYGxvD2NgYMjMzsXTpUgQGBmL37t149tlnIRQK8dZbb+G8885jPPvB9exEe3s7dDodI0dEW5DoQyqVYunSpVi1ahVCQ0MxNDSE8vJyHDt2zGuie/4AX0cYnmo4keJ5SkoKVq5cSRfPRSIROjo6GBXP/YUw3KlhzPd0FLBAIwzgdC6cjeueI8IYHh5GXV0dFi1ahKmpKY/qFfZgW0zfvn07jh8/jn/961+49tprceTIEYdqlJYwm81obGyEQqFASEgIUlJS6H8TCAT49a9/jby8PFx33XU4duwY1q5di71796KwsNDpcbkiDLVajfHxcYSHh9MquJ7AUfShVCrR2NjosPYxn+GvNQx3YTl5npWVhZmZGahUKiiVSvT29gIAHXlERUUhICDAL4QHAfZzGKSl1h/W7gkWHGGQljtPPTHMZjNaW1vR09OD/Px8xMfHo7e3l/MnV9t2XYFAgNdeew11dXVobm7Gddddhy+//NIpUc3MzKCqqgoAkJmZOUsCmYCIF15++eVoamrCxo0b8corr+Daa691eGwuUlL9/f1obGxEYGAg4uPjOSEL283TUe1jeHiYrn2QjcfTwqsv4esIw5vfm0wmQ2JiIhITE62K5wMDA2hubkZQUBBNGu4MznEJtimphRJhzM+7xgU8VawlnhsjIyMoLS2ln/DdUax1BXsbcnBwMPbu3Yvg4GAcPHgQDz30kMP3q9VqlJeXIzg4GEVFRZDJZE43eCJeeNFFF0Gn0+Gmm27C7373O4eNAp5EGBRFobm5Gc3NzVi5ciVnRT9XxyDRR1paGlatWoV169YhPT0dRqOR9suura1Ff38/Y6FKf4CvU1JzeX5SPLecPE9LS4PJZIJWq6U9z3t7e6HRaOY8+mKbktJoNPO+pRaYpxEGl8N7toQxNjaG6upqhIeHY+XKlVZP9u56YjiDoyf4nJwc/O1vf8M111yDP//5zygqKsLPf/5zq9f09PSgpaXFSuSQSURAxAsfe+wxPPXUU3jppZdo8cLIyEir17pLGMSKVqvVoqSkBEFBQejp6XG6NjYbEps1cRl9+FoaxJfwJWGR4rnBYIBIJEJmZiadvmpvb4dEIqFbdyMiIpwKiXKBMzXCmJeE4Qru2rSSDTgzMxNpaWmzbo65ijAILr30Utx555148cUXccsttyAvLw85OTl0vWJkZASrV6+22uSZbvBCoRCPPPII8vPzcdNNN+G7777D2rVr8e9//xt5eXmsj2cJjUaDyspKyOVylJSU0DevP2hJ2dY+jEYjvfHY1j4iIyPp68Mf4A8Rhq9TeeTJ3p5su0qlQmdnJxoaGmjZ9sjISISGhnL+vfERxgKCWCxmHWH09/ejvb191gZsibmMMAgef/xxVFRU4ODBg9i8eTO+/fZb2kWrtLR01obGtuZw6aWXIisrC5dffjk6OzuxYcMGvPHGG7j44ovdOh6RTE9KSkJOTo7VjeqPkh5isdgq+tBoNBgdHaWjD7lcThdmmWhveRu+rmH4umhrr45iK9tOiudk8hyYXTz3FHyEMY/g6qJlEwmQlJNGo0FpaanTbpq5jjDI+t5++22UlpaiubkZ11xzDZ577jksW7bM7hOOO0VqIl74i1/8AgcOHMCVV16Jhx56CPfddx+rTZ5EaEuWLLGrgsulH4Y3iMfScMhR9EFRFO2jPNfRh6/J1tcRDsCs8O6seN7S0mKVggwLC3OreO5OWy0fYfgpmKakRkZG6HkEuVzusvWSrewIEzAhofj4ePz5z3/GL37xC/zwww/44YcfUFBQ4PB47mwsRLzw3nvvxauvvorHHnsMdXV1uOeeexgNUzU1NWF4eNjpnIo/RhjOYC/6OH78OJRKJbq6umZFH95O1/h6w/b1+QH2nVr2Js9J+zWZPA8PD6dTkEwaMyiKcivC4HJ+y1dYkIThKiVFURRaW1vR3d2NxYsXA/COTSsTuIoIyGYcEhKCRx55BA8//DDuu+8+rFq1CmvXrp31erY+4ZaQSCR44YUXkJ+fj1//+tf4+OOP0dDQgBdeeAE5OTl236PX61FdXQ2DwWA3RWa7Nn+OMFydMzg4GAKBAEuXLoVEIqGjj6amJhiNRlryIioqymvRh69TUr6uYXhaR7FtgNBqtVAqlVbFc/IbOiqek/uLTYQxPT3t0ntmPmBeEoYnXVJ6vR41NTWYmZlBaWkpJiYmADAjDG+lpBy1tOp0OlRVVcFsNqOsrAznnnsu6uvr8cEHH+Daa69FeXn5rKE+LuYmtm7ditzcXFx11VVobW3Ftddeiz179uCcc86xet3k5CQqKysRGhrKSIRxvkUYzmAv+lAqlRgZGUFra6tXog9fP+H7+vzA6c2aq+FZgUCAoKAgBAUFsSqeu0MYGo2Gr2H4Kxxt7OPj46iqqkJYWBhKS0shFosZiw8C3okwRCIR9Hr9rL+PjY2hqqoKUVFRWLp0KX1xvvrqq6irq0NTUxOuvfZafPnll1ZPQVxpP5WWluLw4cP42c9+hoaGBlx44YV47rnncPvtt0MgEGBkZAS1tbVITU1FZmamV9ph5+I4XMCy9pGammpV++Ay+vD1Z/YXwvBWlMO0eE4kbdimpPgahp/CNsKgKAp9fX1obm6e1TLLRnxwrorefX19aGpqsquIS4b61q5di0OHDuGhhx7C008/Tf+7pSKspzf3okWL8M9//hMPP/wwvvzyS+zYsQM1NTW466670NvbS0/AM4WjCEMgELBW0PVneDP68HWE4euU1FymxSyL5xRFYWJiAiqVilZSOH78OB19hIeHO404eMLwY4jFYuh0OgCn29+IxtKqVavopwcCpn4YgPdrGGazGc3NzRgcHLS7VoLs7Gy8/vrruOqqq/Diiy+isLAQl112GX08gLunwcDAQDz66KM4++yzcd999+HNN9/EsWPH8MEHH7AiC4B7IcP5AFfRh8FgsGr5dBZ9+PoJ31/baucCAoGALp5HR0ejsrIS6enpUCqVaG5upovn5He0LZ4vlDmMeSkNwrSGodVqcfToUWg0GpSVldndgNmkpLwZYeh0Opw4cQJqtRqlpaUOyYLgkksuwY4dOwAAt956K5qbm+njAczVel2BfNe33nornn76aQQHB6O5uRnnn38+Tpw4wfpY87XozRVI9LFkyRKUlZXRXtkjIyM4evQojh49itbWVruKu77+zL4mLMA/Cu+kjmL5OxYWFiIqKgoqlQonTpzAkSNH0NTUhFOnTkGhUHg8h/HKK68gLS0NMpkMxcXFOH78OKP37dmzBwKBAJdcconb57bEvCQMVxAKhdBqtSgvL0dkZCStsWQPbFJS3mqr1ev1KC8vR0BAAEpKShhfWH/84x+xYcMGTE1N4aqrrsLk5CR9Q3NJGDqdDuXl5Tj77LNx5MgR5ObmYmBgABs3bsQ777zD6ljzHVxu2iT6SE1NxapVq7B+/XpkZGTAZDKhqakJP/zwA2pqatDX10c/0PApKd9HObbCh6R4npycjBUrVmD9+vXIzc2FWCzGm2++iaysLExNTeHf//43jh49ynoP2bt3L3bs2IFHHnkElZWVKCgowKZNmzAyMuL0fV1dXbjnnnuwfv16tz6nPSw4wqAoCkqlEuPj41iyZAmWLFni9CK3TEm52gy8kZIaHx/H+Pg4UlNTUVBQwKrzggz1JSYmorm5Gbfddht9M3G1sU1NTUGtViMtLQ35+fnIzs7GDz/8gAsvvJCReKEl+AjDOcRiMWJiYpCbm0tHH+Hh4XT0MTQ0BLVa7TO/D3/YrP0lwnC2BlI8z8rKwhNPPIGKigqIRCL09fXhpz/9KaKjo7Fv3z7G53vhhRdw8803Y+vWrcjLy8OuXbsQGBiIN954w+F7TCYTfvGLX+APf/gDPTrABRYUYRCV2fHxcYSEhCAxMdHlewhhmM1mu91KluAyJUXmKwYGBhAYGIj09HS3bsa4uDi8++67EIvF+Pe//42dO3fSx/cEZFZlcHAQISEhVusLDQ3FBx98gPvuuw8A8NJLL+Hiiy+GSqVyekxnG72vNyJ/g73og7gTOoo+vA0+JXUabKXVFy9ejNHRUfzlL3+BQqHAV199hRUrVjB6L9nTNm7cSP9NKBRi48aNKC8vd/i+P/7xj4iNjcWNN97IeJ1MMC8Jw95FOz4+jiNHjkAkEmHJkiWMj2WZ/nGVluIqwtDr9Th58iSUSqXLCIgJSktL8eyzzwIA7r//fjQ2Nnq0TqPRiOrqagwODiI9Pd3u8BIRL3zvvfcQGBhIixc2NjY6PK4zwmATMfh60/IFxGIxZDIZoqOj7UYfpPahVCo5T5sS+EtKytdrYDvlbTKZMDMzg6CgIIhEIhQXF7usURKMjo7CZDIhLi7O6u9xcXEYGhqy+55Dhw5h9+7deP311xmvkSnmJWEA1ptGb28vjh8/jpSUFKxYsQJSqZTxTSORSOhBIFdPalxEGITYpFIpSkpKIJfLOSGh22+/HVdeeSVMJhOee+45xhastpiensaxY8dgMBjoeoqzzfzSSy/F/v37kZqaSosXfvrpp3Zfu1C6pHxJWKQF2VHto7m5mfaK4Dr64COM/62B7dAegDnpkpqcnMR1112H119/HdHR0Zwff1631ZLi4MjIiFUbKtvitFwux+TkpMuby9MIY2BgAA0NDcjIyKBTPFwN2hGnvoaGBjQ0NODGG2/Et99+y8oXQK1Wo6qqCvHx8cjNzaXdC11tzsuXL8fhw4dxzTXX4IcffrASL7S8ubmsYZyJcPTdkdpHTEyM07kPJvMCzsB2o/QG/IEw3NGRAtwjjOjoaIhEolkumkQA0xbt7e3o6urCRRddRP+N7C9isRgtLS3IyMhgvQ6CeRthaLVaHDt2DJOTk7NaZtkSBtNZDHcjDFKvaGxsxIoVK7B48WJ60+OKMIDTF+SePXsQGBiIY8eO4YEHHmD83r6+Ppw8eRKZmZnIy8ujbwimm3x0dDT+85//4PbbbwcAPPbYY7jmmmswNTVl9Tq+6O0+mDzhezP64COM/62BrVKtWCx2S1ZdKpVi9erVVkVys9mMffv2obS0dNbrc3NzUVdXh+rqavp/F198Mc455xxUV1cjOTmZ9RosMS8jDIqicPLkSURFRdFPwpZwJ8IAmNUw2BIGEefT6/UoLS2d9ZTBJWEAQFZWFu6++2489thjeOmll1BUVITLL7/c4evNZjNaWlowMDBgd1iQjZihRCLBn//8Z+Tn5+M3v/kNPv74Y7S1teGDDz7A4sWL3VbS5fE/sN2w7UUfKpUKCoWCjj4s3QadbYT+UMPwhzW464XhLtnu2LED119/PdasWYOioiK8+OKL0Gg02Lp1KwBgy5YtSEpKwlNPPQWZTIZly5ZZvZ/4uNj+3R3MS8IQCARWTm62EIlEoCiK8dMI0wiDbUpqYmIClZWVCAsLcyjO541hwPXr1+OOO+7AK6+8gltvvRVLly612whgMBishBjtzX+4s8lv27YNS5YswVVXXYX6+nqsW7cO7777LtLS0hZEDcNX8PQJ33LqPCUlBUajkZb6JtPKxG3Q3tQ531b7vzXMpdve5s2boVAo8PDDD2NoaAgrVqzAV199RRfCe3p65uw7mZeEAQABAQEON2/yYzJ9EmA67c1mcyf1isWLF1uloOytlSvtJ8t13nPPPaivr8eBAwdw1VVX4dChQ7RoGvA/G9XAwECUlJQ4VAB1N/1DxAs3b96MiooKXHjhhfj973+PCy64wO3PZbmmMxFckyTb6MPXKSk2D4HehMlkYqWYy4Xb3vbt27F9+3a7/7Z//36n7/3nP//p0bktMW9rGM5gSRhMwCYl5SrCIHpQjY2NKCgoQEZGhtObjGspD3JMoVCIt99+G0lJSWhpacGtt95Kbzijo6MoLy9HbGysS1lyT/w1Fi1ahG+//RZXX301TCYTnnjiCTz55JO0zpftedjgTIwwAO+RJYk+UlJSsHLlSrr2Qa7ngwcPYmxsDOPj44xUEbwB8pv7A2G4I22+EB505i1hOPvySfcRW19vJkVv8pRjD2TIRqFQoLS0FLGxsS7P7Q3CIFFBbGws3n33XUgkEnz00Ud48cUX0dXVhaqqKixZsmSW57aj9XmyOcvlcrzxxht4+umnIRQK8fnnn+O8886zavvV6XSorKxEZWUlurq6MDk56fScC+HGcwdz+YRvO3VeWFgIsVgMjUaDY8eOoby8HKdOnfLq3IctyD3ia8Jwp+i9ELwwgHmcknIFNgVqpikpcpHYC4uJmVBISAjttcEE3oowyPFKSkrw3HPP4c4778SDDz4IkUiELVu20IUwV+CiI0kgEODOO+9EbGwsfv3rX+PYsWNYu3Yt9u7di9zcXFRUVCAsLAwxMTFQqVTo7u6GWCymc+mRkZGzblBfRBi+jmp8dX6ilSSTyRAXF4fY2FiHtQ9ic+oN+AthsC16LxTzJIAnDADsUlLA7Bzm4OAg6uvrkZ6e7jIFZQuuxQKB2Z1XW7duxX/+8x988803eP7553H55ZfPKWEQbNiwAa+//jr+8Ic/oLm5GRs3bsQvf/lLbNu2DWlpaTCZTEhKSoLZbMbY2BhGR0fR1taGmZkZekOKjo72eYTh68E9X4FEOLa1D2JzSmofMpnMyu+Dq9kNck37+vd3J8JYCNLmAE8YAJhHGORCJcelKAqnTp1Cb28vCgoKGKWg7B2T69Zay7rD5OQkKioqcP/996O/vx+NjY34xS9+gf/+97+Mhvq4XJtAIEBSUhJ++OEHXHXVVfjuu+/w5z//GWazGX/84x+tzkk8IoDTN9zo6ChNIEKhEDqdDiEhIQgLC/P5E+dcwddFZ3uRta3NqaPOK0ufCE/Oz9ZsyxtgG2FMT0/zhOFreOLrbQumEYZAIKAL36QldXp6GiUlJQgODma2cDvgmjBI3WF4eBi1tbV0p9YHH3yAsrIylJeX47777sPzzz/v8lhcRhiEyPr6+rBjxw4UFBTgz3/+M/7yl7+grq4Ou3fvduhZkpKSQm9INTU1MJvNaGhogMlkQmRkJKKjoxEVFQWpVMrJWv0R/pASc3XfOYo+CNl7En34Q4cUWQcfYSwwuEMYTE2USL0iKCjI6TwIm7VyHWEoFAqo1WorG9XMzEzs3r0bV1xxBXbu3ImioiJceeWVLo/F1UZlNpsxNTUFo9GIsrIybNy4EYWFhbjpppvw3Xff4dxzz8WePXuwZMkSTExMQK/XIyQkxGpClojwBQUFITU1FZOTk1Aqlejv70dzczOCg4Np8ggJCfH50yiX8HWEwfb8jqIPlUqFlpYW6PV6VtGHvxCGu11SCwELmjCYeDQA7GxaAaCurg7p6enIzMzk5AYWCAScdZqYTCZMTEzAYDCguLgYoaGhVv9+0UUX4Z577sHzzz+P22+/HcuWLUNeXp7D43EV/UxPT6OtrQ0AUFxcTJPspZdeiszMTFxxxRXo7OzEueeei7vuugsxMTHQ6/UICwtDfn4+8vLy6O/a8v9DQ0MRGhqK9PR06PV6+mm2p6eH9iUgxVg2vfP+Cl8ThicbtqfRh78QBtt1aLVaJCQkeHFFc4d5ewcxSUkx3eiYEAbxhzAYDMjIyEBWVhbzxboAVxHGzMwMKisrYTKZsGjRollkQfDoo4+ioqIC33//PTZv3ozDhw87fC0XEcbY2BjdQabX62dFZMuXL8eBAwfwi1/8AocPH8bjjz+Oiy66CJdeeikmJiZw4MABSKVSq+/c3pqkUikSEhKQkJAAs9mM8fFxjI6OoqOjAw0NDQgPD6cL5/Pxic/XKSkuJ73diT78hTDYRhh8SmoegMuUlMFgQG1tLR1ahoWFcbZOgJun+LGxMVRVVSE6OhqhoaFObyyxWIy33noLJSUlaG1txS233IL333/f7mZg6eDnzmZBOsiysrIgl8vpKMMW0dHR2LVrF2699VYcOXIEn332Gb788ku88MILEIvFtMovUdB1BaFQiIiICERERCArKwvT09MYHR2FUqlEe3s7/TQbHR2N8PBwv9iIXGG+paTYgEn0QUiD7YbNJdyZNl9Icxj+f5e4Ca5SUlNTU7SzVWlpKSQSiVd8vT2VTT9x4gTS09OxbNkyRhFLTEwM3n//fUgkEnz88cd48cUXHa4NYP90S1EU2tvbUV9fj4KCAqSlpbmMVnQ6HS644ALaXUwgEODRRx9FSEgIXdOwPD4byOVy2nN5w4YNyMzMpFWEDx48iNraWgwMDNidQnf3nN6ArwljLoiVRB+WU+dZWVmgKAp6vR4HDx5EdXU1ent753zqnNxXc6kl5U+YtxEGl11SjgiDdBmlpaXR9Qp3FGtdwV3CsGzrXbFiBWJiYgAwTyMVFRXh+eefx29+8xs8+OCDWLVqFc466yyr11jOiTDdLMxmM+rr66FSqazqKK5+M7lcjiNHjmDfvn1ISEiAUCjE2NgYFAoFkpOT6Q4oTzdNkUhk9TQ7NTUFpVKJgYEBtLS0ICgoiC6ch4aG+k3h3NeE5SvxQbFYjOjoaBiNRhiNRixZsoSzziu2IPc+31a7wOBOhEFSUhRFoa2tDV1dXVZdRuS4XHY0Ae4RBmkv1Wg0s9p6hUIhDAYDo+PccsstOHbsGN577z1cd911KC8vR1JSEv3vlikpJtDr9aiqqoLJZEJJSQlkMpnVsewdh6IozMzM4JlnnsHXX3+NmJgYZGVlQaVSYXR0FIcOHcJjjz1mdZNy2eobEhKCkJAQpKWlQa/X0+ft6+uDQCCw2ox8iYWckmJ6fpFINKv2MTY2BqVSSdc+SK3K07kPe3Bn2pyvYcwDsCEMyzkMy3pFSUmJlcIr4B05craEodVqUVlZiYCAADpN5u7xBAIBdu7cibq6OtTV1dFDfeRpnk1KampqChUVFQgNDcXy5ctnPenZIwyKoqBSqXDdddfh2LFjyMnJQU5ODsLDw9Hb24vR0VEcPnzYirS9uWlJpVLEx8cjPj4eZrMZExMTGB0dRVdXFx2B9vX1IT4+3ieCcr7esH1Z67EX5ZLoIzo62mnnVWRkJCIiIjyOPkj9hOnvQNbkyZyWP2FBE4Y7KamjR49CLpfb3YjJcb0RYTBdq0qlQlVVFRITE5GTk2P3BmarMBsYGIg9e/agrKwMR48exX333Yc//elP9LEA19Ilo6OjqK6uRkpKCrKyshwW0C0Jg6IodHR04PLLL8fw8DDy8vLws5/9DNnZ2XTNgkQsf/3rX3HXXXfRQ4lzkZ4RCoUIDw9HeHg4MjMz6etjYmICvb29kEqlVoVzbxdifZ2S8nWE4Sotatt5ZTKZ6KnzU6dOcRJ9uNOpxdcw/ADemPSemJhAbGwssrOzHR7fGxEGUxLq7e1Fc3MzcnNznVotuqMwm5GRgd27d+Pyyy/HK6+8gqKiImzevJlRSoqsa8mSJVi0aJHD11kShtlsxokTJ3DFFVdAKBQiLy8PV155Jd2vTjrRzjrrLAwNDeGVV17BTTfdhODgYBiNRggEAhiNRlrKfS5A0mtLly6FSCSiNyPLNlBS+7BMxXEFf9iwfX1+Nr+1SCTiPPpgKwsCLKwuqXlLGK7AlDAoisLIyAiA0+2zOTk5Lo871zUM4kkwODiI1atX0xpL7h7PES688EL87ne/w7PPPksP9S1dutRp7aGlpQX9/f2M1kWOYzab8X//93+4+eabERMTg5ycHFx55ZV225WXLl2KRYsWYXBwEP/4xz/w85//HEqlEsuWLYPZbLbKKc8ledhuRhqNBkqlEkNDQzh16hSCgoLoJ9mwsDBO5xd8AX/wovCEsNhGH3K53O652Lb0EqLiIww/gLNuICaEYTQaUVtbi7GxMQCni96unuKEQiHj2ghTONvg9Xo9ampqoNPpHNqo2sIT06NHHnkEJ0+exHfffUc79dlbH/nupqamUFJSwviGMJlMeOGFF/Doo48iPT0dy5cvx89//nMr+Q9LCIVCrF27Fn19ffjzn/+MgoICrF27FsHBwXRPvtlsBkVR9O9CBB3J/3sbltanqampMBgMdOG8rq4OFEXRG1FUVJTbUjK+jDDIfebrCIOrtJ+96EOlUkGpVKKtrQ0BAQH072UZfbBdg8FggMlk4msY/g5XhEEsSmUyGdavXw/g9E2h0+mcphNEIpHTXn134IgwpqamUFlZieDgYKc2qvaO526+WyQS4c0330RpaSk91HfTTTdZHW9mZgYVFRWQSCQoKSlhJPhH6g47d+7E999/j+zsbBQWFuL88893+blyc3MRHx+PoaEhVFVVYdOmTfTntPQTsfyf5W8/19GHRCJBXFwc4uLiQFEUXTjv6elBU1MTQkND6dpHUFAQqwKqr+AP0uLeSolZRh/JyclOow+2a9BoNADARxj+DmeEMTIygtraWiQnJyM7O9sqYtBqtS4JwxspKds2WIVCgZqaGqdFZGfH82SNZKjvRz/6ET755BNERkaiuLgYADA+Po7KykpER0dj6dKljDZhe51Qa9euRUZGBnp6eiCXyxEYGIigoKBZ5GEwGDA4OIgVK1bQtYw77rjDbveaJXlQFAWTyUQTFfmNSeQxl9FHWFgYwsLCkJGRAZ1OR0+cWxpFRUdHM8qj8xGG9ye8nUUfarUaAoEALS0ts6IPe+AJw4/ANiVFunI6OjqwbNkyusAqkUggkUhgMBig1Wqd5uK91VZr6bHR1dWFtrY2LF26FImJiW4dz1NSKywsxJ/+9Cf86le/wptvvokLLrgAhYWFqK2tRUZGBtLT0xltHhRFobOzE5dffjkGBweRl5eHSy65BNnZ2aAoiv7OSQ1AIpEgKCiITr0NDQ0hODgYZWVlqK+vx9DQEHbt2oXf/va3Tj8/YO2QaJm6snxansvUFQAEBAQgKSmJNoqyfZIlRlEkj24Jf0hJ+bqG4akyNFvYRh/d3d1QKBQAMCv6IG6Dlr+RVquFXC6fF9IzTLAwPoUdEMIgF7rRaER1dTX6+vpQXFw8Sz2SjU0r1xEGOSaZkO7q6kJRUZFbZAF47sNNcNNNN+EXv/gFzGYzbr75Znz33XdYvnw5Fi9ezGjjMpvNOH78OM4++2yoVCosXboU1113HbKzswGcvhmlUinCw8ORmJiI1NRURERE0FFFf38/LWcOAOvWrUN8fDxeeuklTExMMP4cQqEQEokEAQEBkEqlkEqlVu6JRqMRer0eRqPR6W/LdUpIKBQiKioK2dnZKC0tRWFhISIiIjAyMoKjR4/i2LFjaGtrg1qtpsnOV/CXCMPXGy9FUZDL5cjJyUFpaSmKiooQFRUFpVKJEydOoLy8HC0tLRgdHYXRaKT15zz53l555RWkpaVBJpOhuLgYx48fd/ja119/HevXr6d11DZu3Oj09WyxYAmDpDbMZjM0Gg2OHj0Ko9GI0tJSu8qsTCXOvRFhkBbR48ePY2pqCqWlpR4JHHpS9LY9zl/+8hekp6dDrVZj586diIiIcPk+y06on/zkJ5DL5Vi+fDluuOEGpzLPZIo3ICAAFEUhMjIScrkcarUaXV1diIiIQGJiIkwmE/7617+69ZmEQiFEIhGkUilkMhkCAgIgFotpkiXkQYqVXD8cOIKlftKqVauwfv16WrK9vr4ehw4dgtFohFKptNLUmiv4Sw3D14RhuQaBQIDAwEBao2z9+vX0w9CpU6dw9tln484770RISAhaWlrcIvy9e/dix44deOSRR1BZWYmCggJs2rSJ7uy0xf79+3H11Vfj+++/R3l5OZKTk3Heeeehv7/f/Q9tgQVLGOQJcmRkBOXl5YiOjsbq1asdFmjZ+HpzvYkQOQq5XI6ioiKPe/i58rDQ6/VobGzEvffei5CQEFRUVODee+91+h6y6b7wwgu4/vrrkZSUhDVr1mDLli0uSZCiKIyOjmJsbAyJiYl0eiY5ORnJyckIDg6mpVr+9Kc/oaqqin76dhe20YdYLHYr+uAaYrEYsbGxyMvLw7p167BixQoAp6/nw4cP4+TJk+js7MTk5OScRB6+nvIG/IMwnLXVktoHiT5ee+01rFy5EmNjY1ixYgUyMjLwq1/9itXv9cILL+Dmm2/G1q1bkZeXh127diEwMBBvvPGG3de/++67+OUvf4kVK1YgNzcXf//732E2m7Fv3z63Pq8t5n0Nw9W/1dXVYdmyZS7TO0xTUlxHGENDQ+jq6kJAQACWL1/OyRMcF4Sh0WhQUVGB4OBgZGRk4KWXXsLWrVvx2muvoaioCFdfffWs95Ausx07duDdd9+lUy3nnXeeyxvdbDZjeHgYRqMRSUlJs3LVEokEYWFhWLt2LRoaGjA0NIR3330XF154IQAgKioKMTExHrWtOiqck3QhaUywjDzmqnBOJOuXLl0KqVQKpVIJpVJpZRRFCufeMIry9dAeWYOvCcNsNjP6fgUCAfLz83HOOeeguroahw4dwv79+9HQ0MD4e9Tr9aioqMB9991H/00oFGLjxo20grYrELkjVzNSTDGvCcMRjEYj6urqAJw257HUIXIEpjatXKnVEvnvzs5OpKSkYGxsjNPhLk+eOon8SFJSEnJycnDkyBH86Ec/wu9//3s8/fTT+OUvf4n8/HwsW7aMfg9FUVCr1bjuuutw9OhR5OTk4LzzzkNRUZHLz2U0GjE0NAShUIjExESnXSdCoRDr1q1DX18f9u7di4ceeoiOTDo7O1FfX4/w8HBER0cjJibG7e4U28L59PQ0Ghoa6K4YX858BAQEIDExEYmJiTCbzbT4Xnt7O6anp+nPz6X4nq+nzAH/IAyTycTKN57IggQGBuKnP/0pfvrTnzJ+7+joKEwmE+Li4qz+HhcXh+bmZkbHuPfee5GYmEhbBniKBUcYRJhPKpVCIpHM6jRxhLlMSRFCm5iYQElJCd22xxU8iTD6+/vR2NhoJT9CCOihhx7CyZMn8e2332Lz5s04cuQIwsLC6M6uyy67jO6E+tnPfuZyah447YExNDQEuVyOmJgYRptSXl4ekpOTMTQ0hL/+9a/43e9+R+s9EaOk0dFR2iiJkIe7RklkZiciIoK2s/XFzIe9TVsoFCIyMhKRkZHIysqaJX8hl8ut1HbdXZc/EIa/pMXYuu35Shbk6aefxp49e7B//37OpGrmNWHYXsAKhQK1tbW0MN/Bgwc99sSwhacpqenpaVRVVUEkEqG0tBRSqRQzMzOc5sfdIQxiQdvT04NVq1YhKiqK/jdSRBeJRPjnP/+J0tJStLe34+abb8b777+PiooKXHHFFRAIBFi6dKmVJpQzaLVaDA8P0wJ/TDckMv3d29uLl156CbfddhvdyECMksgAFtk86+rqYDab6dRNdHQ0oyfF8fFxOtqy9HC3TV0R8vDmzAeTqDEwMJAuxBLb09HRUTQ2NsJkMtGWp1FRUQ6n6x2d2x82a1+vga2WlCfCg9HR0RCJRBgeHrb6+/DwsMusyfPPP4+nn34a3377LZYvX+7W+e1hXhMGAen1b29vR15eHu3n4I5irTfbatVqNaqqquhiJrnwuCpSE7DtkjKZTKitraUjHlsZA8s23ejoaLz//vs499xz8emnn+LGG2/Ep59+6lITyhZk+jkmJmbWEB4TWEYZu3btwu9+97tZrxGJRIiNjUVsbCwoisLk5CQUCgV6e3vR2NiI0NBQ2kjJ3sS1QqFAXV0dMjMzkZKSYncdjmY+LEkE4C51xeYp39b2dGpqCqOjo7RRVHBwME2gISEhTo/N1zD+t4a58vOWSqVYvXo19u3bh0suuYQ+/759+7B9+3aH73v22WfxxBNP4Ouvv8aaNWvcOrcjzHvCMBqNqK+vx9jYGIqKiqw2K3cUa5lEGO74+pJUT3Z2NlJSUqxuPq4L6ZYS4K5u8pmZGVRVVUEoFNIRjy1sayKrV6/Gc889h9/85jf417/+hcWLF7vUhCIgU98TExNISEhgnDK09xkdRRn2QArHoaGh9MS1QqHA6OgoOjo6IJVKERMTQxeOBwcHcerUKSxdunRWDtnVuhwVzj1NXXmSFrI0iiLtuqRwXl1dTae2oqOjERkZOauw6w8pKX8gLbbig2QOw13s2LED119/PdasWYOioiK8+OKL0Gg02Lp1KwBgy5YtSEpKwlNPPQUAeOaZZ/Dwww/jvffeQ1paGoaGhgCA1jrzFPOaMKanp3H8+HGIxWKUlpbO2qy4sGm1heWTJFNZDKLounLlSkRHR9s9JtcpKXJuZzfYxMQEKisrERkZiWXLljn8PJYRC+mEqq6uRkBAABYvXozg4GCcf/75LsnCbDZDoVBAp9MhKSmJVfHQHphEGY4QEBCARYsWYdGiRbR2kEKhQGNjIz3nkJqa6pHLnquJc7aFcy7bZ6VSKRISEpCQkACz2Yzx8XEolUp0dnaioaEBYWFhVoVzf0lJzYU0iKs1sHXbc/Yg4wqbN2+GQqHAww8/jKGhIaxYsQJfffUV/RDT09NjtZ7XXnsNer0el19+udVxHnnkETz66KNur4NgXhMGkfFwZCTkrZQUcPpJw1V7HXHv02q1ThVdvZGSApxf3CMjI6ipqcHixYtdTm5bRiyWnVDZ2dkIDw/H1NQUPv74Y1x77bUOvxOTyUQ/7SQlJXFy45OOKaZRhiOQ/vnIyEia1OLj46FWq9Hd3Y2QkBC6cO4qdeNqve6KJXpz0looFNKTwaRxgEQfJPoKCgqCyWRi/YTNJfwlwmBDGNPT04zqec6wfft2hymo/fv3W/13V1eXR+dyhXlNGNHR0U7z5d5ISZEL1tVxSWeNXC5HSUmJ09kArgnDMsKwBUVR6O7uRmtrq5WeljOQCKOrqwuXX345BgYG6E6o2NhY/OMf/8DAwAC+/fZbnH/++bPeT6Q+pFIpYmNjOX1SXbJkCR1lvPbaay4HCx2B1HFmZmasvMj1ej1GR0ehUChosUBCHpGRkW5vnq5mPmx9PuZSmkMul8+Kvvr6+qDT6XDw4EG6cB4dHc2qcO4p/KGGwZYwF5J5EjDPCcMVvJGSEggELlNIJC9M5hhc3eRsag5MYLkRWcJsNqOpqQnDw8MoLCxknG4RCAQ4efIk7r77brudUBdffDE++OADVFZWIikpCfn5+fR7Z2ZmMDQ0hJCQEERGRnK+4VlGGS+//DJuu+021rIqer0e1dXVEAgEWLNmjRW5S6VSq5kH0nXU0tICnU5H5/1jYmLcbl1kkroivyUhk7k2iiJT77m5uRgdHZ1lFBUdHY3Q0FCvEZo7dUNvwJ2i90LxwgDmOWFwadPKNCUFOC9S9/T0oKWlxaVdqe3xAO5ytPZ8uA0GA6qrq2kjJibFZkJiJ0+exG9/+1vExcUhIyMDl1xyCWJjY+nXZWZmYu3atTh8+DC+/PJLxMbGIi4uDlNTU1AoFIiMjPRIG8sVLKOMXbt2sYoytFotqqqqEBISQluvOgIRCySCgVqtFgqFAkNDQ3TXESEPTzZP2+jDaDTi1KlTtAucZdvuXPl8kM2aFE/T0tJgMBjo1FVNTQ0EAoFV4ZxLZVl/UMsF2Kek+AhjHoENYZCnQyaEYS/CsHx6X7NmDSORPgJvEIZlZxMZZpTJZIyNmMiT7UsvvYRHHnkE6enpWLZsGc4991wYDAZ0d3dDJpPRU6zr16/HwMAAOjs78dFHH+HSSy+FVqtFbGys170A3I0yJiYmUFVVhbi4OEaRoCUsZa/J5kkGBisrKyEUCul5j6ioKLflOiiKQn19PfR6Pa2FZm/mw9sT5/aiX4lEgvj4eMTHx4OiKLpw3t3dTbctk8/PxijKHuZSisXZGiiK4iOMhQqRSMRY2ZNpSgqYHWGQlIbBYGD89G57POD00wtXT2WkLkJmP+Lj45Gbm8u4s0un0+Huu+/GO++8g6ysLJSWlmLTpk30++35WJx11lm0yczXX3+NK6+8krMJU1dYsmQJUlJSGEcZ5Kl48eLFSE1N9TiVIpFIZnUdKRQKtLe3o66uzip1xfT6MBgMqKqqgkAgwOrVq+lrw9HMB7kmvUEerrqkBAIBPYCZkZGBmZkZOvro7OyEVCq1a3nKFP5CGGzXoNFo3G4d90cseMLguoZBjksunsnJSVRWViI0NBSrVq1y60mS3NxcF76Hh4fR0dGB7OxspKamMnqfvU6oTZs2zdKEImKAYWFhMJlMmJ6ehkajQXFxMfbt24e+vj4cPXoUGzZsmJObnMxl9PT0uIwyBgYG0NTUhLy8PI87WBythXQdkdQVKZyfOnUKgYGB9MxHWFiY3e+HzMfI5XLk5+fb3WC9OfNhC7b1NZlMRhtFmUwmWu/K1igqOjqa0UOFP6SkyP3JlOyIWx8fYfgJfF3DIK2paWlpVrIR7oDLTinyxNnR0YEVK1YgJiaG8fssO6GWLFmCSy65xKUmlEgkgkwmg1qtRlJSEs455xzs27cP5eXlkEqlSE9Pp1NX3lBSJbCMMl577TX8/ve/t/v5urq6sGLFCiv5E28iMDAQKSkpSElJoT0tiAUvAKvUlUQigVarRUVFBa1dxWST5HrmwxaetLQSNd2oqCgrvauRkRG0trYiMDDQqnBub10mk4lOtfoK7qxhenp6wdizAvOcMAD2Nq2OwLStFjh9cw4NDWF0dJT2Z/AUXBGGyWRCfX09TCYT8vPzGZOF2WzGyZMn3dKEIgKCgYGBiI6ORkJCAkZHR1FTU4Njx44hOTkZk5OTGB0dpXv6AwMDIZVKOd0ALKe/X375Zdx+++10lEEGKEmNyR05Ei4gFosRFxeHuLg4Ou9vqbQbHBwMrVaLuLg45OXl+WTmwx646uCzrP0QAlWpVLTeFzHOIoVzMtzpD4ODpODN5nvwRBrEHzHvCcMZvJGSMplM0Gq1mJycnCVF4gm4IAydToeqqiraRpJpqG82m/HJJ5/gpptuYq0J5UhA8LzzzsPQ0BCGh4exb98+XHvttRAIBHTdY2xsDEKhkCYP0gHkKWznMn7/+9/TJDo1NYWioiK/ySlb5v0zMzMxNDSEhoYGyGQyDA0NYWxszGOlXcBx6spW78qZWKK3NmxiFEX0viYmJqBUKtHb24umpiZ6aJKr68MTsG1KMZvNPGHMJ3CdkiJ5ZbPZjNTUVE5bRT0ljKmpKVRUVCAsLAz5+fkoLy93eTySqnj55Zfx8MMPIz09nbEmFAC6KyY2NnZWnlYikeCyyy7DG2+8YTXUR/SMzGYzZmZmoNFooFAoYDabaaXVwMBAjwbiLKOMm266CZ2dnaAoCoWFhR7LkXgLIyMjaGxsRE5ODhYtWjTryZso7RKTKHc/h6vUlSOxxLmYshYIBHRdbPHixdDpdHThvKurC2azGc3NzYiKivJoaNJduNNSC4CvYcwXcJmSGh8fR2VlJa2tw/XTlieEMTo6iurqaqSmptK1FFcmSkw6oZy9V6VSYXJy0qmAYHh4uNVQX2JiIi21LBQKaXKgKAp6vR4ajYbuLgoICLBKXbGBZZTx8MMP45ZbbnFYOPYHkCL8smXLaI0ge0/eo6Oj6O7uttJ6cqS0yxRMC+dc+9gzgaVRlFKpRGNjI0QiEdra2qDT6RAeHk7XPuYianRnaA8AP4fhT+CqhmEZYdjmawcHB1FfX4/MzEykpaWhoaGB8xvIXcIgg4JLly61sqF1djzSCbVlyxaUl5c77ISyB7PZjJGREej1eiQmJrrczDMzM7Fu3TocOnSIFk2zVX8lLnIBAQGIjIyE0WiERqOhjaXEYjFNHjKZjNHkfHFxMXp7e/HRRx/h8ccf91uy6O7uRkdHB1auXOnQRtPyyZu0rJKZD1ul3cjISI9SV8D/ugBJxKHX6zE2Nobw8HDo9fo5dxkkkEqlyMrKcmoU5azzzFO444UhFArnrLV8LjDvCcMZ3CEM4HTqSS6Xg6IotLW1obu7GwUFBfR0M9fqsgB7iXOKotDc3IyBgQG7g4KOCMPdTijAfQFBMtTX0dGBjz76CFu3bnV6E4nFYnqDNJvNdMsuMZKxTF3Zu4Gnp6cREhKCxMRE2pXP0hfZH0Curf7+fqxatYpVelMmk1lpPZHUVVNTEwwGg5VJlLtaT5azQfX19ZBKpUhLS6OvK2/OfNiDbUrMkVEUeZgjhXNP0nf21uCOF4avay9c4owgDCYdHpYhrVarhUQiQV1dHSYnJ2eZCgmFQrpNkcu1MiUho9GImpoaaLValJaW2g15bQnDUubDnU4ovV6PoaEhBAQEICYmhtXmIBAIcPHFF+Mf//gH1Go1PvvsM1x++eWMbiRSGA8KCqLTaBqNBmq1GiMjI5DL5fTmIZFIMDU1hZGREURHR+Oss87CwMAAdu7cidtvv90jqXIuQVEUmpqaoFQqUVhY6FFRVCQSzTJJUigU6O/vp4vGJPpgq7Sr0+lQWVmJwMBA5Ofnz9nMhz0405GyNYqanJyEUqlEf38/mpubacmWqKgoj9SG3RUe5AnDj+Dsx7As7Ln6oUUiEaRSKfR6PVQqFZqbmyGRSFBSUjLrCUUkEkGn03m+eAswTUlNT0/TnuXOVHAtU3WedkJNT09jeHjYIwHBwMBAXHrppXjrrbfQ2tqKI0eOYO3atayOIRAIIJPJIJPJEBUVBYPBAI1GQ0+bkweEyMhIhISE0L7kpGPKH6IMs9mMuro6aDQaFBYWcpqusDRJWrx4sZXSbldXFyulXXKdhYWFzZoF8fbMhz0wFR60NMqyNYrq6emxmgmxZxTFxRoIFpqOFLAACMMZLL0rmDwZBAYGQq/Xo7y8HGvWrHEopeGtlJSrY5LCe0xMjMuBLnI8TzuhiK0pGaryBAkJCdi0aRO++OILHDhwAAkJCVi8eLHbx5NIJAgPD0dYWBhGR0cxNTVFDxCOjY0hKCgIhYWF6O3t9Ysog0SGRqMRa9as8XrHlj2lXYVCYaW0S6IPS+Iig4PR0dHIzc1lVDPicubDHtxVqnVkFNXR0YGGhgarwrmrzd2dGgYfYcwjkB/XaDQyujnJBko8t50dd66L3kNDQ7S/dFpaGqObmBQsd+zY4VYn1NjYGMbGxhAfH8/Zk9KKFSvQ39+PmpoafPLJJ9i2bZtH7ckURWFkZAQ6nQ6LFi2CRCIBRVGYnp6mDbbi4uIwNDSEZ555Bo8++uicejgQ6PV6VFVVQSwWY/Xq1V6deLcHS6VdiqKg0WgwOjqKwcFBq7RNUFAQWlpakJiYiKysLNabHRczH/bAhbS5pWSLpVEUaR4ICAigycPe3Is7KamFNIMBLHDCIN4VrjZ3s9mMlpYW+mJw9SQ9lxEGRVHo6OhAR0eHVeHdFQQCAcbGxvCrX/2KdScURVFQKBSYnp5GYmIi5xvspk2bMDw8jKGhIXz00Ue47rrr3NpATSYThoeHYTabkZiYSB9DIBDQdY2oqChs2LABw8PD2L17N4qLi5GYmEjnvOeiKEnSOyEhIU6tcOcKAoFglkw5IY/Ozk6IRCIYDAZamt5dcnN35sMevOGFYWsUpVKpoFQq0dTUBKPRiIiICLr2ERAQwLro7amftz9i3hMGEz0pZ5u7wWBATU0NZmZmEBERgYGBAZfT3nMVYZjNZtTX10OpVKK4uJhxSoi0zd56662YmJjAkiVLcPHFF2PJkiUu32u5CSclJXnlSVgsFuPSSy/FP/7xDwwODuKbb77BT37yE1bHMBqNGBwchFgsRmJiosPNRCAQID8/H0ePHsXQ0BBqa2tRWFgIhUJBP1WStExERATnm9LU1BQqKysRHR2NJUuW+GV6QiKRQCaTYXx8HJmZmXSKr7W1FdPT024p7dqDJ2KJ3jZPstc8oFQqMTg4SHudUBRF/z+T35GvYcxDiEQihx1NGo0GFRUVCAoKsuqEYuLr7e22WpLCMJlMKC0tZSzzQVEUKioqcP3110MgECAvLw/nnHMOAgIC0NfXR3ccSSSSWRe9wWDA0NCQy02YC5Chvr1796KqqgqJiYkoKChg9F69Xo/BwUHI5XLExMS4vHkFAgHtl7Fr1y78+te/plVUVSoVFAoF6uvrraapo6OjPZaaHx8fR1VVFRYtWoSMjAy/JAvgf1Lv2dnZtOlXZGQksrOz6dSVrdJuTEwMwsLCPBoYBKzvJWeFc7b1A09g2TxgaRTV0dGBkZERKJVKq8K5o+uET0nNQzhKSZHp6OTkZGRnZ0MgEDAWIPRGhGG5TvJUGhISgvz8fMaGR2azGZ9++iluuukmREdHW3VCmUwmehhOrVZDLBYjMDAQQUFBkMlk9CYcFBSE6OjoOdncMjIysH79ehw8eBBff/014uLiXAo5EsvX0NBQREREMF6n5fT3q6++ivvvv3/WU+XExATt3U0KouTf2T4pkk04IyODsbS8LzAyMoL6+nosWbLEbns1ecBITU2FwWCgCba6uhrAbKVdd0CIwDL6sPX5IAODRqPRK227zkCMohQKBUJDQxEWFkbLlTQ2NiIsLMyqcE6uSZ4w/BBsJc4pikJ3dzdaW1tnTUczFSBkMxDIFEKhkH6Sqa6uxqJFi2gicwXyRLZz50489NBDdjuhRCIR3W5oOwxHIpPg4GBERUXN6ZPwunXrMDAwgPb2dnqoz1HaQ6PRYGRkxC3LV4FAgPXr19MdU7/85S+tOqYsp6lJQZQ8WRMJbqZP1sPDw/QmbHl9+RsGBwdpSRImtTGJRDJLaVehUNBKu4RgSfHcXdimrgYGBjAyMoJly5bNKpzPJXmYzWaIxWK7RlFEcVgqlSIiIgJNTU2YmpryKCX1yiuv4LnnnsPQ0BAKCgrw8ssvo6ioyOHr//Wvf+Ghhx5CV1cXsrKy8Mwzz+CnP/2p2+e3h3lPGK5gubmbzWY0NjZCoVCgsLBwVosl+XFnZmZcHtMbKanJyUkMDg6y8gMnOkw7duzA22+/zagTynIYbmxsDCqVCnK5HDMzM+jq6oJcLqf/3duSGmSo74033sDY2Bg+++wzerDQEuPj41CpVB5ZvlrOZZAowxHkcjmSk5PpSWLiYUGerAl52BaF+/r6cOrUKSxfvpyxtLwvYLnO6Oho1u+3VNrNysqaRbByuZwTpd3BwUGcOnWK9i+Zi5kPR7CXFrM1ilKr1WhoaMDvfvc7jI2NITk5Gbt27cIFF1yA5ORkxufau3cvduzYgV27dqG4uBgvvvgiNm3ahJaWFrvkfuTIEVx99dV46qmncOGFF+K9997DJZdcgsrKSixbtszjz04goJwp1M0DmM1mGAwGh/9O+skTEhLomsCqVavs1gRuuukmvPPOO3jiiSdw9913Ozzm5OQkjh49ih//+MecfAaKonDixAmo1WqsWbOGsbGPrSbU4sWLcd5556G4uJhRJ5RSqcTU1BTi4+Mhk8lAURQ9DKfVaqHT6TwSAWSDwcFBvPXWWzCZTNiwYQPWrVtn9RknJibodXqCpqYmvPfeexgcHERjYyPruQxL+1WFQkE3S8TExGBmZgZ9fX1YsWIFK0/3uQbRr/LWOonSrkKhwOjoqNtKu4TUVqxYYVdny3bmw3Ir80b0ceLECaSmpjKKxsxmM66//nr6Pjp8+DAefPBBPProo4zOVVxcjMLCQuzcuZM+XnJyMn71q1/NMgYDgM2bN0Oj0eDzzz+n/1ZSUoIVK1Zg165dzD4gA5wREYZWq8WRI0cQHh7uVLWUTUqKqwjDaDTSEiTEtpIJSGrtsssuY60JZSkgmJSUROeeBQIBpFIpHVYbjUbav8JSBDAoKAgBAQGcpq4SEhJw/vnn4z//+Q9++OEHeqjPsr2XC8KyjDJeeeUVPPDAA6zeb9nLT4rCJC2j0+kQFBQElUoFkUjkkQyFN0BRFDo7O9HT08Nav4oNuFDa7e3tRWtrK1auXOmQ1JzNfFhGH2xnPhyBTVstqXP++Mc/xm9/+1uo1WpG5mzA6aaOiooKK2UCoVCIjRs3ory83O57ysvLsWPHDqu/bdq0CR9//DGjczLFvCcMVzekXq/H8PAwMjMzsXjxYqevZ2rTKhQK6YKcJxfgzMwMKisrIRKJsHjxYiiVSpfvseyEuuKKKwAAS5cuxRVXXMEoX240GjE0NASBQOBSQFAsFlvVPbRaLbRaLS1AaGl+xMWTXEFBAfr7+1FdXY1PPvkEF110ER3yc9Xea1nLeOWVV3DHHXd4NP0tl8sxNTUFoVCIoqIiaLVaunBuqXEUERHhU8VcInY4MDCA1atXz5njoCOlXcu2ZkIepK25p6cH7e3tWLVqFePfhsuZD0dwZ9KbpE/JQwYTjI6OwmQyzVJ1jouLQ3Nzs933DA0N2X09uVe5wrwnDEcgA29qtRrR0dHIyMhw+R6mXVKWF6W7G+XExASdLlu6dCkGBwddFtJddUK5AumEkslkiI2NZfX0KxQK6WEviqJo8yOlUgmj0Uh3XHlifgT8z6lvaGgI3333HW644QbOZ0Fyc3Np7293ogwCk8mE2tpazMzMoLCwEAEBAQgLC6NlKIgUh6WKLCkKz6WRE7GnVSgUWLNmjU87dxwp7TY0NMBoNEIul0Or1aKgoMAjIvdk5sMR3FWrXUhYkIRhMplQV1eHsbExJCQkMP6R2aSkyHnc2cxGRkbolsv09HR6It1ZmsteJ1R+fj4uvfRSRpPY09PTGBoaQlhYGKt2VHsgLchEAl6v10Or1c4yPyLzHmxAURRKSkrw5ZdfQqlUYt++fZx3epC5jJ6eHrejDIPBQBfA16xZM+tzWkpx5OTk0Cqyvb29dCsmFwZIrkAaPcbGxrBmzRq/sacFrIflcnNzcerUKfT19UEul6O6uppW2o2JiUFwcDAnMx+A+2KJbCOM6elptwgjOjoaIpGIlvMnGB4edth2Hh8fz+r17mLeE4btRWSZ5iktLUVfXx80Gg2jYzFNSZFzsm2tJV4UbW1tyM/Pt/oxXRkese2EsgSXAoK2sDQ/InUPS/MjiURCRx6u6h5kxiIuLg4///nPsWfPHlRXVyMpKYnxUB9TeBJlENlvmUyG5cuXu3wgsVWR1el0dNHcctrc044iW9gq4/pCQ4spurq6MDg4iMLCQoSGhkKn09EmUWyVdl3BNvqwnfmwfZ3la5mel6IotyMMqVSK1atXY9++fbjkkkvoc+/btw/bt2+3+57S0lLs27cPd955J/23b775BqWlpazP7wzznjAsMTY2hqqqKjrNQ35srmxaCZhEBLawbOktKiqalUJypiVl2QmVk5PDqhNKrVZjfHycUwFBZ7A0PzKZTPS8x+DgIIRCIZ26ksvlVusnMxYRERG0guiGDRvwww8/MB7qYwPLKIPMZTDJMWu1WlRWViI8PNylYrAjBAQEWKVlSMsu8e4mG6Mnw3AkXabT6eZEGdcTtLe3o7e316q2EhAQQLerWqb3mpubodfrHSrtsoWz1JXlzAcAxpIgBBqNxm0/7x07duD666/HmjVrUFRUhBdffBEajQZbt24FAGzZsgVJSUl46qmnAAC/+c1vcNZZZ+FPf/oTLrjgAuzZswcnT57E3/72N7fO7wgLgjAEAgH6+/vR0NCArKwspKam0j+sO657TLoZ2BARSV/o9XqUlJTYTQvYIwxPOqEsBQSTkpJ8smGIRCKrugdRkFUoFDCbzbRAoNlshlKpRGxsrNUNtnbtWgwMDKCtrc3lUJ87sI0yHnzwQaevn5ycRGVlJeLj4xkPVbqCSCSa1VFkOQxHWnbZ6DgZjUZUV1eDoiisXr3aY4kTb4GiKLS3t6Ovr89pId42vedIaTcmJgahoaFeSV2pVCp6zQaDgVHh3BMtqc2bN0OhUODhhx/G0NAQVqxYQVscA6etmS3PXVZWhvfeew8PPvgg7r//fmRlZeHjjz/mdAYDWABzGBRFob6+Hj09PSgoKJg1LDUwMIDe3l4UFxe7PNann36KK6+8EiUlJdi/f7/T1+7fvx8FBQUun0qJt0BgYCAKCgoc1jxIdHTOOefY7YTKzMxk3AllKSAYHx8/51LarkBSbBqNBhMTEzCZTAgICEBISAjtnEcwPT2Nf/zjHxgbG0NGRgauvPJKTvP9ZC5jYGAAjY2NDn9P8vukpqbSdSdvY3p6mk5dqdVqBAUF0U/VjqbNDQYDqqqqIBKJsGLFCr/1Midk0d/fj9WrV7v9JE4MkhQKBZRKJYRCoVXqiotrX6VSobq6GllZWYiLi2M080FRFKKjo1FbW4vc3FyP1+Av8K+dxA2QzcfWRpWATYRBnuBc1TAAZhGGSqWihfVcGdFYGh550gk1lwKC7oLMe0xMTEAgECA+Pp6ufSiVSkilUiudK+LU197ejkOHDmH9+vWcrYVJlEHSRVlZWaymdT2FXC5HSkoKUlJSaNkYhUKBqqoqq40xKioKIpGI7t+Xy+VO5418DYqi0NraisHBQY+7tmwNksbGxjhV2lWr1aiurkZOTg6SkpIAwOo+tS2ck8jDaDTCYDC4TYT+inlPGEKhEMuXL3dYT3CHMJikpFzVMPr7+9HY2IicnBykpKS4PB4hDHudUD//+c8Z5WlJ0dgXmlBsQAYHDQYDEhMT6YiC1D3IvMfAwAAtY3Luuefiv//9Lw4ePIjExERGbdJMQOYyLDumLKMMMhG+dOlSzjtO2IAI4MXHx9MbI1GQ1el0CA8Px9TUFEJDQ7F8+XK/fFAATpPFqVOnMDw8zHmLr1AoRGRkpEOlXSKsyVRpd2xsDNXV1cjOzqbJgpwHmJ26siQRMlPFpQWvP2DeE4YrzHUNgzw99fT0YOXKlYx1eoRCIfR6PbZv347333+fdSeUJ8J8cwmTyUQPEyUmJs56CiYT0iEhIXTdQ6PRIDo6GpmZmWhra8PHH3+MrVu32pWLcAeE1G2jjJ6eHrS1tdE6Rv4C241RpVKhtrYWQqEQSqUSJ06c4KQdlWvYzoN4uwnDVmmXiAQyUdol0vSZmZkudd1sC+eko9ET+Xd/BU8YFmDaVuvouGT+Y3x83GGKzB5IC97zzz+P2tpa5OTkYN26dSgtLWXUCTU+Pg61Wu2RMN9cwGAwYHBwEFKpFLGxsS6J0NI5jzwZvvfee1AoFNi7dy8uvPBChIaGIigoyKNctW2U8ctf/hIqlYouxvozAWu1WjQ0NCAxMRHZ2dm0e55CoUBXVxckEonVtLmvIg+KotDc3IzR0VGfzINYRmiWSrsdHR10cwEhEJPJhMrKSmRkZLBOQZrNZtx2223o6urCqVOn/OpBgwvwhGEBT1JSpDdfIBCgtLSUcVcS6YS6/PLL0d/fjyVLluD8889HQkICPYhD5hhsDeUtBQQTEhL8OvzV6XQYHBx0O10mEAgQFBSEK6+8Em+88QbUajUqKipQWFhI1z0sRRLZHt8yyvjDH/6Ayy67DGvWrPHrHPTk5CQqKiqsDJqkUikSExORmJhIq6cqFAp6ktoX0+YURaGpqQkqlcovhgftKe0SocTW1lZQFIXw8HAEBwezUnMwmUzYvn07qqurceDAAcZ2yvMJ875LCjjdRuiIFGZmZrB//35s2rTJ5SaiUCjoJwqtVuv0QqmqqkJ4eDjS09PpGzciIgLLli1jVGwknVAnT57ElVdeCWB2J5SlBIdGo6FbUUkxWKFQwGg0Ij4+3m9bJ4HT3+Xw8DB9k3oapnd0dGDPnj0AgJ/+9KfIz8+nRRK1Wu0scyim52tqasL777+P7u5u1NTU+LWXBenaSktLQ3p6usvXUxRFD3AqFApMTU0hLCzMytvcG6AoCo2NjbQSsz8/1ExOTuLkyZOIjo6GUChkpbRrNptx55134vvvv8f333/PqG45H3FGRBgAMxkPy5yqq7F+EmEoFArU1NQgLS2NsQ0n004oSwmOqKgouhVVrVZDr9fTyqn+nCclm1RMTAxngneLFy/GWWedhQMHDtBDfQkJCQgJCZllDgWAjs4CAwMdPgSYzWZaB2pwcBC7d+/GQw89xMl6uYZlmyfTlIlAIKCFJIkIIHmqbm9vh0wmszKI4iJ1RVEUGhoaMD4+7vdkMTU1hYqKCqSmpmLx4sUAYNeFkZAsMYkSCAQwm82499578c0332D//v0LliyAMyDCMJvN+O9//0v7WjsDeYIHTssrOzPAqa+vx/T0NMbGxrBs2TK79pb2QFrwSHGVbScUERAkchy2vhXu6Dd5A5a1lbi4OM4LnBRF4cMPP0RrayvCwsKwbdu2WakOiqKg0+noCI2I25HUFXmAIIV4gUCA8fFxvP/++/QgKFeFda5AWnxzc3M5i4As/SsUCgUAWLXsulMfMpvNaGhowOTkJFavXu3XsiQajQYnT56kU3uOYKm0q1Kp8MUXX0CtVkMkEuH777/HgQMHkJWVNYcrn3ssCMIwmUx0H7Q9fP3111i/fj2jTSs8PBwzMzNoaWlx6MVsNptRXl5Oa/QwlS0mMyN333033nrrLWRkZKC0tBTnnXceozQWSe3YCgha+lZMT0/T+k1BQUFu5fM9hW1txVubxczMDO3Ut3jxYmzevNnpZyUiiRqNhiZZmUwGjUZDF+IFAgF2796NY8eO4dZbb8XDDz/slbW7A2L9umzZsllS1lzBsiCsUCig1WqtZDiY1B/MZjPq6+sxNTU1b8giKSmJcYYAOL3nfPXVV9i5cyeOHz8OoVCI888/HxdeeCGuuuoqn9dpvIUFn5IC2HdKzczMOOyUIrILMzMziIuLY0UWY2Nj2LJlC44cOcJKEwoAbUJjL7Vjz7dCo9FYzTGwzee7C5Km0+l0VuZM3oBMJsNll12GN998Ex0dHTh48CA2bNjg8PXEHCo8PBxGoxGTk5NQq9X01C6ZpiYaU6+++iq2b9/uF1HGwMAAmpubvW79alsQJlIuw8PDaGlpoafNHclwEMFDrVbr9xpWRIWBzPWwlftvbGxEQ0MDjh49Coqi8Pnnn+Of//wnrrrqKi+u2rfgCcMGzqa9ifAc0fTX6XSMjmk2m9HT02PVCcVGE4oICCYkJLh8crH1rbCXzyfif1y3WFpKkrgyZ+IKcXFxOP/88/H555/j0KFDSExMRGZmJqO1jo+P09Ea+Z6GhoYQEBCAxMREDA4O4qWXXmJsq+ktEPc5R1al3kRgYCBSU1PpWQaSkqmsrIRQKLTyNhcIBLQ/yOrVq/2aLKanp1FRUYG4uDhkZmayIguKovDSSy/hL3/5C7755htaSXnFihUu9cjmO3jCsIGj1tqxsTFaeC43Nxfd3d0u5zXsdUKxccejKAojIyOYmZlxS0DQdo6B5PO5Nj0CTkdeg4ODPpEkWb58Ofr7+1FVVYVPP/0U27Ztc+pvQbxBiDIu8L8hL9KZtnr1agwMDGDnzp0499xzkZGRgZiYmDnfBDs7O9HV1cXKfc5bkEgks2Q4FAoFWlpaoNPpIBaLIRAI5gVZnDx5EjExMaxFJCmKwl//+lc888wz+Oqrr7B69WovrtT/4J/6ASzh6gf3dNp7cHAQJ06cQEZGBi1p7UpLyrIT6ic/+QkCAgKwfPly3HDDDYwFBAcGBmAwGDhRmxUIBJDJZIiKikJycjIWLVqEgIAAjI+Po7u7GwMDAxgfH3daC3IEvV6P/v5+yGQyxMfH+2Q47Mc//jESExMxMzODf//73zAYDHZfR6KI6Ohouxsw6Uxbs2YNUlJSIJPJ8H//93/o6+vDDz/8gBMnTqCrq4uxP7O7IJaq3d3dWLNmjc/JwhZk2jwnJwclJSW0f3lAQADKy8tx/PhxdHZ2YmpqCv5UJp2ZmaGdLnNycliTxZtvvolHH30Un332GUpKSry4Uv8EH2HYwHLamyhqdnZ2oqCgwGoQx5mWlG0nVFpaGpYvX864E4rtRDRbkAEvqVSKiIgIGAwGuu5hOQRHiubOwKWTnycQi8W49NJL8cYbb2B4eBj//e9/ccEFF1i9ZmJigpZRdzV3YDn9/f777+Ohhx5CYGAgXQxua2tDYGAgYmJiEBsb65Gsti0s9ZYKCwv9enrfZDKhpqYGwGkTH4lEQpsfkUlqqVTqF9PmMzMzOHnyJCIjI12KgdqCoii8//77uPfee/HJJ59wKoA5n8AThg1ISkqj0aCurg4qlYp+grKEMy0pvV6Pe+65B2+++SatCcW0E8oXAoISicTK9IiQh1qthlgspsnD1jGP2I5GRUVx7uTnDkJDQ/Gzn/0Me/bsQU1NDZKSkrBixQoAp1OKarUa8fHxjDtYsrOz6envl19+GY888giSk5ORnJwMo9EIpVKJkZERu/l8d1N8tlPRc2F65S5MJhOqq6thNpuxatUquv3W0vyI+HYrFArU19dbDcJFR0fPWQu4Tqejh2uXLFnC+r7697//jTvvvBMffvghzj33XC+t0v+xIAiDSUqKaaqF3KANDQ1IS0tDaWmp3bZAexGGp51QZAP2pYCgpfgf6bjSarWzHPP0er1f6lelp6djw4YN9FBfbGwsAgICMDk5icTERFYtnpZRxmuvvYbt27fT2kBisRhxcXG0P4KlI5zBYEBUVBRiY2NZbYqkHXVyctLvB90sTZpWrlzpcFbD0rebDMKNjo7Sg3Dh4eH0v3uLHInse1hYGPLy8liTxaefforbb78d77//Ps4//3yvrHG+YEEQhiuwsVMlN7fJZEJRUZHDJ0XbCMPTTih/FBB01nFFhhxJrcaf5LTLysowMDCA1tZWfPjhhzj//PORlpbmVh3IMsrYuXMnHnnkkVmvsXWEm5qawsjIiNWmGBsb69STgQhXTk9Po7Cw0K+LxkajkfbkWLlyJeNoSiAQ0JFsRkYGpqen6dRVa2srneJjKj/OBIQsQkJCsHTpUtbH/PLLL3HjjTfizTffxMUXX+zxeuY7zhjCYJKSGh0dxdTUFAAgIiLC6Y1AjmnPHY9tJ9To6Cg0Gg3rJ+C5BCkGT01N0ekXnU4HtVqNkZEReoI6KCjI58Y9AoEAF1xwAXbv3k3rfDFptXV0LEdRhqPXkwiNbIqk7kE8GUjdgxSKSWrHZDJhzZo1fjGp7whcOvrJ5fJZKT6FQmElP+7JtDkhi6CgILfI4rvvvsP111+Pv/3tb7j88stZn38h4owhDFcpqZ6eHrS0tNASH65aZknUYjab8dlnn+HGG29EdHQ0srOzsXnzZkYpJbPZjOHhYRiNRq8PuXkKs9mMoaEhesZCLBbTdRZLu9XR0VGfy5SQls+zzz4bX375JTo7O3Ho0CGnQ33OkJ2djdTUVKdRhiPYuuaRJ+qTJ09CLBYjKioK4+PjkEgkVnUAf4TBYEBlZSUkEgkKCgo4fTCwTPGR1C5pLqirq6OnzWNiYhg3jlRWViIwMBDLli1jHQEfPHgQV199NV5++WVcc8017n6sBQf/vTpZgEkNQ6/X2/03otM/MDCA1atX4/vvvwfgmjCIDeNLL72Ehx56iHUnlNFoxNDQEIRCoV0jIX+C7Vptbz7LjititarRaKBSqazsVudCpsRkMmFwcBAikQhLly4FAHz22WeshvpsIRAIsG7dOnR3dzOKMhzBdo5hZGQEzc3NMJlMmJmZQWNjI2JjY+0a+vgaBoMBFRUVCAgIQEFBgVdTkAKBABEREYiIiKCd8xQKBYaGhtDS0oLg4GCaPEiUZrvWyspKBAQEID8/n/Vajx49iiuuuALPPvssbrjhBr8W95xrLAjCAE5fZI76vR3VMIxGI2pqaqDRaFBSUkJPQAPOCYOiKIhEIuzevRtffvklsrKysGbNGvzkJz9hdKPrdDoMDQ1BLpcjJibGry9IInbIdK1isdhux9XY2JjbsuNMQdqRAwICaF2o/Px89Pf3o7KyEp9++im2bt3KWM7FEpZRxssvv+zx9LfBYEBHRwctiU8aHjo7O2lDH1L38HXxW6/X0woHvrB/JdFqWloa9Ho9nbrq7u6GWCy2atmlKApVVVWQSqVuEVtFRQUuvfRSPPbYY7jtttv8+t70BRYMYTiDvRoGGeCRSCR0/zjg2kTJXifU+vXrkZmZid7eXshkMvoCt5de4NobwpsgLb4hISG09AMb2HZc2cqUkClzLmRKCLEFBQXNakfeuHEjhoaGMDAwgI8++ghbtmxh/QRvG2X86le/cttNjchShIeH04OghGQzMzNn6TeRJ+rY2Ng5t1wldYDAwEC3nta5hlQqtYrSSHdaU1MTDAYDPTzozlpra2vxs5/9DPfffz9+/etf+/W96SucMYRhWcMYHx9HZWUlYmJi6BuWwJlNq9lsRm9vLy677DK6E+pnP/sZcnNzAZx+atRoNJiamoJSqZyVy3cmIOhv4Noj3FIE0dIYipjUWHpWsE3PEWILDQ21OzxoO9T39ddf44ILLmC9IXARZWg0GlRWViI6Otrh8JilfpNer6frHt3d3ZBIJHTkER4e7tUNnMwuBAcHu1UH8DYsu9MyMzNx8uRJGI1GiEQiHD58GKGhoVYGUc5+78bGRlx00UW488478dvf/pYnCwc4YwiDRBhDQ0Ooq6tDZmYm0tLSZl0Y9iIMe51QeXl5uPLKK606oSQSCR05WEqOq1QqCIVCmM3meUEWZCI6JibGKxaljoyhSKHTnmeFI5CIzRWxWQ711dbWIikpCStXrmS97vXr17sdZUxOTqKyspKupTDZlGwtV8kQXF1dHX09edJJ5AiELEg7qr+RhSVIl5lYLEZhYSFEIhF0Oh3dndbR0YGAgAD6u7Il2lOnTuHCCy/ELbfcgoceeognCydYMIThqoZhMpnQ0dGB9vZ2LF++3KGfgK2WFJkz+Pzzz7Ft2zbGnVBEcjw4OJgWECS2qmNjYz71q3AES2VcNhPRnoCkEAICAhAZGUlHaZOTk1YdV8Sr2xJk3oEpCaenp+Oss87C/v378d///hdxcXGsTYiysrLcijLGx8dRVVWFlJQUpKenu/Wb2w7BEd+K9vZ2upOIRB+etGeTdG1YWJhb7ahzCUIWZICQRKcBAQFYtGgRFi1aZJdop6amMDQ0hFWrVuHqq6/GL37xCzz22GN+/Vn9AQuGMJxBIBBAq9Wiu7sbRUVFTjd6Wy0po9GIV199FQ888ADrTiji5AYAycnJdPHdl34VjkBRFBQKBaanpzkRO3QXzqI0S2OomZkZqFQqxMfHs5oQLi0txcDAAE6dOoWPPvoI27ZtY/V+d6IMtVqN6upqLF682KEpF1vY+laQTqLBwUE0NzezSsdYgugtRUREuDUVPZcgOlYmk8lpS7K9afNPPvkEf/nLX9Dd3Y3ExEQkJyejq6uLtmflYR8LwnEPOF0/sNcJpdfrceLECUxNTeGss85yutFrNBo899xzePrppxEdHY2//vWv+OKLL/D222+zdsfT6/UYGhpyKiBIURS9IZKIhmyIcrl8zm5Wy3mQhIQEv5wFsCRajUYDiqIQFBSEsLAw1kQ7MzODf/zjH1Cr1UhPT8fmzZtZpVwoisIbb7yBY8eO4cYbb8Qf/vAHh68dHR1FbW0tsrOzsWjRIsbn8AR6vZ5Ox5BammXdw9F3RYrxkZGRbuktzSXMZjNqampgMBjcml8ZGBjApk2bUFxcjJKSEnz22Wc4cOAAPvzwQ36i2wkWDGHY8/XWaDSoqKiATCbD+Pg4fvzjHzt8/+TkJO655x58//336OnpgVgsRkBAAHQ6HXJzc+mLi8lNND09jeHhYVbdRZaFYI1GA7PZbEUe3sohk7kFoVCIuLg4v54HsbR+jfh/7Z1pWJNn+vbPsKPsQsK+qYAVZa24W6e2uAGJY+vSKlZnOh1rtXbmbWun1ta2ttaOR6u2tXZqsVVHyyKgIC4IYxVE2VQEVHZkCfuShJCQPO8Hj+f5JwhIQgjb/TuOfmhMyJ1AnjP3fV3XeVpaMrUPQPVgqLq6OkRERKCrqwtz5szBggULVFrLgwcPcPz4cVRWVuLevXuwtrbu8Tny8vIwZcqUfme+axqZTMa0odJ53Yp1D/r3TWdE9FWMHy7I5XLcuXMHnZ2d8Pf3V7njjc/nY/HixQgKCsIvv/zCvAetra0wMDAYtfGqmmD4fZXUEE1NTcjJyYGjoyMcHBxw/fr1Pu9/5swZpKeng81mo6KiAjKZjAmFmTFjRr+979V1cO1eCO4t7Gj8+PEaE4+e5haGK3SYVPfoV3WDodhsNpYuXYr4+Hhcv34d9vb2mDx5cr/X072W0X2XUVNTg4KCAnh7eyvZ4msbXV1dsNlssNlspQnqBw8eoLOzE1ZWVrCwsEBlZSXYbLbKGRHaho6ApUOuVBWLhoYGhISEwNfXF0ePHlX6+xgqw8+RxPBtfRgAjx49QlZWFjw8PODp6Qk9PT2meN0bqamp0NPTY4J36EyIcePGPXXqG/i/gnF9fT04HM6A7L67hx3RNYWWlhaUlZWhpqYGbW1t/bZs74nOzk5UVVVh3Lhxw14saFsSqVQKe3t7pYvEQIKhvL29mcS0+Ph4NDc393tNdC3D1tYWhw8fRkNDA/Nvjx49QkFBwRMZKkMNPUHt4eGBOXPmICgoCOPGjUNxcTE6OzvR1taG8vJyZtc23KDdfEUikVo7i+bmZoSFhWHy5Mk4fvz4oB+9Xr16FSEhIbC3tweLxUJsbOxTH5Oamgp/f38YGhpi0qRJiIiIGNQ1qsqoEgyKonD//n0UFhbC398fTk5OAMD8YfR2gVUUE1owpFIpxGIx2tvbkZaWxngS9ZaBUV9fj7a2Ntjb22vUppnuIrKysmIuiMbGxmhvb1c7KU8kEqG6uhoWFhawtrYe1mJBH5nJ5fKn1lfoYChLS0s4OjrCyckJ48ePh1AoREVFBR49eoTm5mZIJBKmo27RokVwcHBAZ2dnn0l9PUHvMvT19XHw4EEAQHl5OR4+fAg/Pz+1B/u0AYvFAovFQm1tLZydnTFv3jzY29ujqakJ6enpSEtLw8OHD9HS0jIsEvMoisK9e/cgEAjUioBtbW0Fl8uFvb09Tp06pRXrFaFQCB8fH3z33Xf9un9paSmWLVuGhQsXIjc3F2+//Tb+8pe/4MKFC4O80v4zamoYtH1Be3s7/P39leYH5HI5Ll68iOeee+6JojdFUZDJZIiIiMC3334LR0dHtLe3M/YdjY2N0NXVxYQJE2BhYQErKytMnjwZXl5ecHV1hY6ODmprayGTybReMFb0bRKLxf0y/Wtvb0d9ff2ImAehc8LpYbWBHMUp2pSIRCKlYCiJRIKjR49CJBJh2rRpWL58eb9FVLGWkZCQAIFAAH9//2ERKNUXAoEAWVlZPc6EKDrH1tfXaywcSl1osWhra0NAQIDKLcMCgQBcLhfjx4/H2bNnh8RqhcVi4cyZM+Byub3e57333kNCQgLy8vKY21avXo2WlhYkJSVpYZVPZ9TUMEpLSyGRSDBz5swnvn3o6OgwNtKK0DsLuVwOHo+H69evIysrCzo6OqAoChMmTMCaNWswf/58nDt3DufOnUNpaSmKi4tx48YNWFhYwMHBAa6urkPiNNrdt0nR9E+xBZV+P1paWtDS0qJyK+pQQNdXjIyMNOK39bRgqIULFyIxMRF3796Fg4MD/P39+/VzFWsZhw4dwqFDhwZl2FGTCAQCZGZmwtHRERMnTuxxMl4xHIqueyiGQ9GJeYPdfk1RFPLz89Ha2orAwECVxUIkEuGll16Cvr4+zpw5M+S+XH2Rnp6ORYsWKd0WHByMt99+e2gW1AOjZochlUrR1dXV67fQ5ORkPPvss8w3P3pnQYf/sFgstLW1ITExERkZGdDX18eCBQvwwgsvMB+Krq4uXL9+HXFxcYiNjUVzczMsLCxgbm4OCwsLTJo0CZ6enpg0adKQ5lr09G2adtcdzpkbNJ2dnaipqVHbw0oV6GAokUiEW7duMcFAL730Etzc3J66q6EdAOLj41FeXo6CgoIeO6aGC3Q+iJOTEyZOnKjSYymKYpo66urqIBAIBjUxr3tcraoXe7FYjFWrVkEoFCIpKWlId3392WF4eHjgtddew44dO5jbEhMTsWzZMohEomHRvTVqdhhPS9VTtAehdxWKYgE8to9YvXo1Vq9e3ePP0NPTw4IFC/DMM89g6dKlaG9vR0ZGBuLi4pCbm4vS0lLcunUL5ubmcHd3h6enJzw8PLT+i1b8Nk0PD9L27rW1tcNiULA3Ojo6UFtbqzVzRhaLxXhYLV68GO3t7SgqKsLZs2cRHBwMCwuLXoOh6NoVm81mdhkHDhzA7t27B3XN6tLW1obs7Gy4uLjAzc1N5ccrhkO5u7tDLBYz4vHw4UMmHMrGxgZmZmYD+t3RsQPqikVnZyfWrVuHlpYWXLp0adgfEY4URo1gPA0dHR10dXVBLpczwqEoFv2lsrISDx8+hLe3N2xtbcHj8fDFF18gNzcXsbGxiIuLw+3bt1FWVobs7GyYmZnBzc2NEQ9tHlcoTpo7OztDR0fnCcdYxRbUoRYP2vBQ1ZZkTcFisRAaGoqIiAg0NTUhOzsbS5YsecKmhHYi5vP5TOfW/PnzUV5ejsOHD2Pr1q3DbpdBG266ubnB1dVVIz/TyMiIScyTSqVobGxEXV0dsrOzlaarraysVB6MvH//PhoaGtQSC6lUitdeew1VVVVITk6GhYWFiq9saLC1tWU+lzR8Ph9mZmbDYncBjCHBoB1rabGgu0T6C0VRePjwIaqrq+Hv76/0R8hiseDn5wc/Pz98/PHHyM/PZ46t7ty5g/LycuTm5sLMzAwuLi7w9PSEp6fnoF4UpVIpamtrnygY09+mra2te3SMpcVD22ZzdDGetvAeKoyMjLBixQocO3YM5eXluHPnDp577rkngqFYLBZ0dHTAZrOhq6vLmFny+fxht8ugxUKT1iTd0dfXh62tLWxtbZVsx/Pz89HV1cXErVpbW/fZoURRFB48eID6+noEBASofKHs6urCX//6Vzx8+BApKSnDulOtO7NmzUJiYqLSbZcuXcKsWbOGaEVPMmpqGHK5vNeWSIqikJGRASMjI7i4uPSY0tUXMpkMeXl5EAgE8PPzU+mstqioiNl50KJhYWEBMzMzODk5wdPTE15eXmqF+vQG3eFFC8PTXitFUczUtEAg6Pfwm6ZoaWlBc3Oz1gwP+8O9e/cQFxcHAFi5ciU8PDwAPP47q6mpgUwmg4GBAVMjGjduHPh8Pk6fPo2Kigrcu3cPNjY2Q/kSADx+b3NycjBx4kQ4Oztr/fkpikJ7ezvq6upQX18PoVAIS0tLZveh+PumKApFRUWoqalBYGCgyjURmUyGN998Ezdu3MD//ve/IZuupxEIBCgqKgIA+Pn5Yf/+/Vi4cCGsrKzg7OyMHTt2oKqqCr/++iuAx4073t7eePPNN7Fx40ZcuXIFW7duRUJCAoKDg4fypTCMGsGgL3o93S6TydDc3IyKigo0NjbC2NgYbDYbHA7nqYE0EokEubm5YLFY8PHxGVBXSEVFBeLi4hAXF4eMjAyYmpoyRXN7e3tGPAZynKGJGgAtHkKhEBKJRCW7cVWghx3b2tpga2s77DpYLl68iMzMTBgaGuK1116Dubk5E//K4XAYy3r6mE8oFOLChQvIycnByy+/jK+++mpIrVaam5uRk5ODyZMnMzNJQ01HRwdT92hpaVGKW+Xz+aiurkZgYCDGjx+v0s+Vy+V4++23kZKSgpSUlCERx+6kpqZi4cKFT9weHh6OiIgIbNiwAWVlZUhNTVV6zPbt25Gfnw9HR0fs3LkTGzZs0N6in8KoFoyeOqFkMhkaGhrA5/PR0NAAAwMDcDgcsNnsJwp1QqEQOTk5MDMzw9SpUzX64a+pqcHZs2cRGxuL69evY9y4cUy3FYfDYcSDw+H0+6Kvqt13f6DtxoVCITo7O/s169EfKIpCQ0MDRCIR7Ozshswdty9kMhmOHz+Oqqoq2NjYYNGiRTA2Nu71d0LPC/z+++8oLS3Ff/7zH0ycOLFfRzGapqmpCbm5uVo1PVQVqVSKhoYGZvdBURRsbW1hb28PS0vLfh+LyuVyvPvuu0hISEBqaqpaBX1C/xi1gtFbJ5QitDEb/Qerp6fH7DxkMhnu3LnDtB8OZkG4oaEBCQkJiIuLw5UrV2BoaMjsEGxsbBjxoC0GeoI+1uFwOIM2Y9F9UJC2T1Gc9egPFEUxBePh6o5L097ejp9//hkikQiTJk3CypUr+7yQURSFX375BRkZGVi7di1ef/111NXVQSgUaiyv4mk0Njbi9u3b8PT0hIODw6A9j6YoKSlBeXk5Jk2axNSy5HK5Ut2jt78RuVyOnTt34vfff0dqaqpKfmAE1Rk1ggE8Prun0/FULW7L5XJGPPh8PmQyGaysrODu7q7V7O3W1lacP38esbGxuHTpEnR1dZljK2tra6Zg7uTkxAwY0g6u2jzW6T7r0X1QsLf3i/aFoq0+hrM7LvD4W3BOTg6Sk5NBURSCg4MZ/6neePjwIX777TelWoZIJGK+mLS2tsLMzIwRD1WPX/qCFgv6C8Zwp6ysDGVlZQgICGB2xXRmBf1+iUQiWFlZMUdX9N84RVH4/PPP8fPPPyM1NRVTpkwZypcyJhhVgiEWi5ldBaBeJ1RpaSnKysrg5uaGjo4O1NXVgcViMY6fqmyVB4pAIMDFixcRFxeHpKQkyGQyWFpawtzcHFZWVvDw8ACHw8GECROUHFy1jWJWhUgkYkKhTExMYGhoyPwOFK3UbW1th3XsJ/C4llNdXQ0TExOm60ZHRwfr1q3r85u74i5jw4YN+PTTT5X+XTE+tLGxkTGAZLPZKjdkKEJnbwylnboqlJeXo6SkBAEBAX12DHYX20uXLjGf7ePHjyMlJQXTp0/X4srHLqNGMJqampicBLrlURXkcjkzVerr68t826GtEfh8Purq6kBRFGxsbMDhcFTuLx8IYrEYV65cQWxsLBISEiAWi5ljK0tLS8bfys3NbUiPeBSLwIqhUEZGRmhubh4RVurA/02bm5mZMR1sMTExuH//PkxNTbFx48Y+dwY97TJ6oquriznHb2hogL6+PmxsbMBms5/Inu6L+vp63LlzB1OnToWtra3qL1jLVFRUoLi4GP7+/irZikskEpw4cQJHjx5Fbm4u7O3t8fLLLyMsLAxz5swZ9jvWkc6oEYyjR4/ib3/7G+bPnw8ul4vQ0NB+O7FKpVLcuXMHUqkUvr6+vR7r0HkCisdWiuKhrT/W9vZ2REREIC0tDWlpaWhvb2fEQ9GiZOLEiUNaTKZDodra2iAQCMBisWBiYjLooVADRSwWo6amBpaWlkrzNp2dnYiIiEBjYyNcXV2xevXqXl+D4i4jPDwcn3322VOfVy6Xo6mpSakIbG1tDTabrRR21J26ujrcvXsX3t7evWbVDycqKytRVFSkslgAj9/XH3/8Ebt370ZsbCwTt0r7gA23gcnRxqgRDPo4KTIyEjExMcjOzsacOXMY8eits6WjowM5OTkwNjbGtGnT+v3tnKIotLa2MuIhlUoZ8ejrwz1QBAIBcnJymBhNiqKQlpbGtOs2NDQw3VYWFhaYOHEi4281FG2rYrEYtbW1MDU1ZazGhUIhZDLZkA4K9oZIJAKfz4eVlVWPF7P6+npERERAKpVi1qxZPbZN0tC7jPLycuTn56s0l6H490UHR9FFYBsbG+b4kc/nIy8vD9OmTRtW2Ru98ejRIzx48OCJ4df+QFEUjh07hvfffx8JCQmYN2+e0r8N913raGDUCIYiFEWhvLwcUVFRiImJwa1btzBr1iyEhYUhLCwMdnZ2YLFYyMnJQUtLC+zs7ODh4aH2RYseTuLz+eDz+ejs7GSOFfrq8FCV5uZm5Obm9tq5JZfLkZmZyQwKPnr0iBEOc3NzJYsSbbjV0hff7t/UFQcFhUIhpFIpM+vRk2eTtqCtSaytrftsS+5tqK87FEUhIiICN27c6Pcuo7efQ6+NNv2ztLSEoaEhamtr4ePjMyyGBJ9GVVUV7t+/Dz8/P5UHVSmKwsmTJ/HOO+8gPj6+T6EmDB6jUjAUoSgKlZWViI6ORkxMDNLT0xEUFARvb2/89ttv+PHHH7FixQqNfTuhHT3pmkdHRwcmTJgADoczoF58Oh+6v0NYFEXh9u3bzM6juLiY2XmYmZnB1dWV6bgaDCsOVWZCug8KGhkZKXk2aQN6vf21Juk+1GdlZdXj/YqKivDbb7+hrKxM5V1Gb3R0dKC4uBg1NTUA8ETH1XD8pl1dXY3CwkL4+vr2+l71RVRUFDZv3oyoqCgsXrx4EFbYM9999x327dvHCPPBgwcxY8aMXu//zTff4IcffkBFRQWsra2xcuVKfPHFF8NuKFVdRr1gKEJRFKqrq7F161bExcVh4sSJsLS0RFhYGHg8HpycnDT+YaMvRHw+H0KhkBEPxWOFp6FoeKjusUNBQQHjb5Wfn8/kaJibm8PZ2ZkRD03kGre1taGxsRFsNlvlllF1QqE0tV5VZlhkMhlOnDiBR48ewcbGBuHh4T3WizS1y1CEvvj6+PjA1NSUKZo3NjbCyMiIEQ9zc/NhIR50vrmPj49a3k7x8fHYtGkT/vvf/yI0NHQQVtgzp0+fxvr163H48GEEBQXhm2++QWRkJO7fv9/j5/DkyZPYuHEjjh49itmzZ+PBgwfYsGEDVq9ejf3792tt3YPJmBIMAPjHP/6B3377DXFxcXB1dUVMTAyio6Pxxx9/wNfXF1wuF1wuF66urhr/sNHHCnw+HwKBgBnkYrPZvV5siouL8ejRI/j6+mrMdbO4uJjZedCOuvTOw9HRkRkUVPWbIN0UQIc0DdQXSjEUqqOjg5n1MDExgb6+vkZ+P62trWhqalJrve3t7Th69CiEQiG8vb0REhLS45oUdxn37t0bUK2BPtbp6eLbfRCVNkfUdju4IrW1tcjPz1dbLM6fP4/169fj2LFjWLly5SCssHeCgoLw7LPP4tChQwAeH/k6OTnhrbfewvvvv//E/bds2YKCggIkJyczt/3jH/9ARkYGrl27prV1DyZjTjBOnDiBWbNmwd3dnbmNoijU1dXhzJkziI6ORmpqKqZNm8aIx2BMend0dDDHVm1tbbC0tGQ+3IaGhpDL5cjPz0dzczP8/f01OtylyKNHjxjxSE9PV/K3srOzU/K36us9UBwgtLOz0/gkc18Rq4qzHqrQ3NzM1LDUPTIoLy/HyZMnQVEUXnzxRQQGBj5xH8Vdxvr16/H555+r9Vx0wbg/xzqKjrF1dXWQyWRKHVfaOOqjC/LTp09X6yjuypUrWL16NY4cOYK1a9cOwgp7RyKRYNy4cYiKilIKPQoPD0dLSwtTw1Lk5MmT2Lx5My5evIgZM2agpKQEy5Ytw7p16/DBBx9ocfWDx5gTjKdBX/ho8bhy5QqmTJkCLpcLHo+HyZMna1w8xGIxs/Ogp4Bp592AgACtnX/W1tbi3LlziI2NxR9//AFjY2OmaE77W3l6esLW1lbpPaCDhMRiMezs7AZ9gLC74R89KNjfUCiKotDU1IT29naNiNuNGzdw5cqVPof6BrrLoI8l1S0Yt7W1MeLR0dGhZFMyGK3XdKuvumJx9epVvPTSSzh48CDCw8O1frRWXV0NBwcHpKWlKdmLv/vuu/jf//6HjIyMHh934MAB/POf/wRFUejq6sIbb7yBH374QVvLHnSIYPQBfWGJi4tDdHQ0Ll++DA8PD2bn4eXlpfE/5La2NuTm5kImk0EmkzEFTQ6Ho1Xr78bGRiV/KwMDAyV/Kw8PD3h5ecHOzg51dXXo6uoaEl8oOmKVFg8AjHgYGxs/8fuhvxAIhUKNmR5SFIUzZ86gsLCw16G+gewy6CE3Pz8/jRxL0kej9fX1aGtrg7m5ObO71cTfGD1EqG6r740bN8DlcrF371688cYbQ1KHUUcwUlNTsXr1anz22WcICgpCUVERtm3bhr/+9a/YuXOnNpc/aBDB6Cf0+Xx8fDyio6Nx8eJFuLu7MzuPZ555ZsB/2LQ7rrm5OaZOnYquri6mlbKpqQmmpqaMeGijLZamtbUVSUlJiIuLw8WLF6Gjo8McW1laWsLZ2Rl+fn5wcXEZ0nkKelCQFo/uoVAsFmvQdkKKQ30uLi5Ys2bNE++FOrsM2j5DnSG3/qAYs9rc3Izx48cz4vE06/+eaGhowO3bt9UeIszKykJISAg++eQTbN26dciK9uocSc2bNw8zZ87Evn37mNuOHz+O119/HQKBYNjMGg0EIhhq0trairNnzyI6OhoXLlyAs7Mz023l7e2t8h9Ha2srcnJyYG9v3+Oxl0QiYT7YjY2NGD9+PGPLrs2EOqFQiEuXLiEmJgbnz58HRVHMzoP2t/L09ISrq+uQ2jRQFIXOzk5GPLq6upj1DJadekNDA3755RdIpVLMnDkTf/rTn55Ykyq7jLKyMpSWlg6aWHRH0W68sbERBgYGSjYlT7t408aHzzzzjFr2JLdv38ayZcvw/vvv4//9v/835B1eQUFBmDFjBg4ePAjg8VGos7MztmzZ0mPROyAgAIsWLcLevXuZ2/773/9i06ZNaG9vHxW2JUQwNEBbWxsSEhIQHR2N8+fPw97enhEPHx+fp4pHfX097t69i4kTJ/YrQpP+YPP5fCYQSlE8BvuD1tHRgaysLBgbG6O2thbx8fFISEhAR0fHE/5Wnp6ecHd3H3J/q5qaGkilUujq6g7qoGB+fj5iY2MBAH/+85/h6emp9O/93WWUlpaivLwc/v7+Q5JvLpPJlGxKWCyWUkZ39/eMzt9Q1/gwPz8fS5YswdatW/Hhhx8OuVgAj9tqw8PD8eOPP2LGjBn45ptv8Pvvv6OwsBAcDgfr16+Hg4MDvvjiCwDAxx9/jP379+PIkSPMkdTf//53BAQE4PTp00P8ajQDEQwNIxAIkJiYiOjoaCQmJsLGxoYRD39//yfEo6qqCoWFhWqbxtHmdXQgFN2Hz+FwBuR82hvt7e3Izs6Gra0tPDw8mJ8vlUpx9epVxMXFIT4+Hm1tbcPC34q2U6fDeWjB6B4KZWJignHjxmnkmOrSpUu4desWDA0NsWHDBqV2UtreIj09HevWrcOePXueeHxxcTEqKyuVLL+HErlcrmRTIpVKlbIqaLsadS3VHzx4gMWLF2PTpk347LPPhoVY0Bw6dIgZ3PP19cWBAwcQFBQEAHjuuefg6uqKiIgIAI8/i59//jl+++03JnQrJCQEn3/+ucZa4ocaIhiDiFAoRFJSEqKjo5GQkABLS0uEhoaCx+MhICAAH330Edzd3bFixQq1pl+7Q6cJ0t8KDQwMGPHoniaoDnQ+tIuLC9zc3Hr9eTKZDOnp6YiNjUV8fDzq6uqULEpof6vJkycPageYTCZDbW0tWCxWr3bqmgqF6v68J0+eRGVlZY9DffQuo7S0FPn5+cwuQ3HuJjAwUKtHjf2FdjKgxUMoFDJi7OHhoXLHWUlJCRYvXoxVq1Zh3759o+KcfzRDBENLdHR04MKFC4iOjkZ8fDxMTEwgFouxZ88erF27VuPnm72lCfb3PLo79LGZqvnQcrkcWVlZjL9VZWUlIx5mZmZwd3cfFH8rOntDMX+7P48RiUQQCARKg4JPC4XqCcWhvqlTpyI0NJR5fE+7DFosqqqqEBAQMCzFojstLS3Izs7GhAkTIJFIVA6GKi8vx+LFixESEoIDBw4QsRgBEMHQMiKRCC+//DLu3LmD2bNn4+LFizA2NkZoaCi4XC5mz56tcfGgbbP5fD5zHk3XPOj8kL6oqalBfn7+gLMWKIrCnTt3mEHBhw8fMrsOc3NzuLi4wMvLCx4eHgM6iunq6kJNTQ2zw1JnZ9U9FEpXV1flQcGKigqcOHGix6E+xV3GvXv30NraipqaGgQGBg7akKYmaW1tRXZ2NiZOnAhnZ2cAqgVDVVdXIzg4GM8//zwOHz5MxGKEQARDy7z11luMKaClpSU6OzuRnJyMqKgoxMXFwcDAACEhIeDxeJgzZ47Gi8X0BDA9ZU4XMzkcTo/2ERUVFSgqKlLb2qEvCgsLGX+re/fuKVmU0P5WXl5eKnUISaVS1NTUwMjICDY2Nho5D+8+KMhisfqc9VAkIyMDycnJ0NHRwauvvgpHR0cAyruMZcuWYe3atQgICBgRYtHW1oasrCy4u7v32qRB19bq6+vR0NAAPT09iEQiSCQS+Pr6IjQ0FEFBQfjll19GRffQWIEIhpZpbm6GsbFxj2f3UqkUKSkpiIyMRGxsLFgsFiMe8+bN0/gENUVRaG5uZqbM6TRBNpsNKysrlJWVobKyEn5+foPe1llaWsocW2VlZanlbyWRSFBTU4Px48djwoQJg1I87WlQkJ716CkUiqIoxMbGoqCgACYmJti4cSNz3ETvMoqLi5GdnQ1XV1eNr1fTtLe3IysrC66urv1eL73D/e9//4u9e/dCLBaDw+Hg66+/xpIlS7Q6U0QYGEQwhildXV1ITU1FVFQUzpw5A5lMhuXLl4PH42HBggUa7zSiA3voTA+JRAIWi4XJkyfDwcFBq98Cq6qqmG6rtLQ0jB8/XsnfysvLC56enko7iO6RqtrotOk+69FbKFRPQ30sFgsNDQ04c+YMsrKy8OqrrzLtmcMVgUCAzMxMpulBVZqbm7Fs2TLY2Nhg+vTpiI+PR1VVFU6cOAEejzcIKyZoGiIYI4Curi788ccfzM5DLBYz4vHcc89p1OhPLpfj7t27aG1thbW1NZqamiCRSGBtbc1kemhTPPh8PhISEhAbG4urV6/CyMiIEQ/a38rd3R1yuZyZ/xgKeguFott1m5ubERERAYlEgqCgIEyfPh0dHR0Qi8U4ceIEU8sYrhGrtFg4OzsrGXf2l9bWVoSGhoLNZiMmJgaGhoagKIrJCBkJaYEEIhgjDplMhmvXrjE7D4FAgGXLloHL5eL5558fUJuqTCbD7du3IZFI4O/vDwMDAyZNkD62EovFSuKhzYG8pqYmJCYmIi4uDsnJydDX12fEY8KECZgyZQq8vLzg4OAw5L38PYVC1dTU4MKFCwCA+fPnY+bMmdDV1WVqGcN1lyEUCpGZmQkHBwdMmjRJ5ccLBAJwuVyMHz8eZ8+eHTVhQmMRIhgjGHregRaPlpYWLFmyBDweD4sWLVLJSE4qlSInJwc6Ojrw9fXtUQjoqFD62IpOE6TbKAfbpVaRtrY2REVF4dSpU8jOzgaLxWLadSdMmMA46zo7Ow95B45UKmXadW/cuIH79+9DX18f69evB4fDQXFxMX799ddhucugxcLe3h6TJk1SWYhFIhH+/Oc/AwASExO1WtRXNS2vpaUF//rXvxATE4Ompia4uLjgm2++wdKlS7W25uEOEYxRglwuR0ZGBpNj3tDQgMWLF4PH4+HFF1/ss7AoFouRnZ2N8ePHw9vbu99HTrR40DnTiuIx2NPcdNbC1KlTYWZmhsuXLyM2Nhbnz5+HVCoddv5WdOZKR0cHUlJSUFNTA3NzcyxduhTm5uaIi4vDjRs38Morr+DLL78ckjV2RyQSITMzE7a2tmrZ+ovFYrz88svo6OhAUlKSVqfWVU3Lk0gkmDNnDthsNj744AM4ODigvLwcFhYW8PHx0dq6hztEMEYhcrkcmZmZjHjU1tbixRdfBI/HQ3BwsNJQGG3rYGVlhSlTpqj9bVwkEjHi0d7eDktLS2bWQ9PiQUeUTps27Ymshc7OTqSkpCAuLg4JCQkQCoVKzrqTJ0+Gl5cX3NzctLYjosVCIpHAzs4OYrEYR48ehUAggIeHB+bOnYvi4mJcvHgRJSUluHHjhlrf5jVJR0cHMjMzwWazlSxg+ktnZydeeeUV1NfX49KlS1q3xlA1Le/w4cPYt28fCgsLtbpTHmkQwRjlyOVy5OTkICoqCtHR0Xj06BFeeOEF8Hg8WFlZYfv27YiIiIC/v7/GLlAdHR1MzYP2lKJzzAd6fk0HCfVnLoRuFqA7rpqbm3v0t5o0adKg7YgoigKfz4dUKoW9vT2zw6msrMSJEycgl8vxwgsvICAgABEREcjIyMD8+fPxxhtvDGgyfyDQYmFjYwNPT0+Vn1sqlSI8PBxlZWVITk7W+PzO01DHmnzp0qWwsrLCuHHjEBcXBxsbG6xduxbvvfcemRNRgAjGGEIul+POnTuIiorCr7/+itraWsybNw9r167F0qVLNeI31Z3uaYJ0WA+Hw1FZPGi7b3WChGQyGTIyMphZDz6fr+RvRVuUTJ48WWNBVbRY0OFS3S88N2/exOXLl6Gjo4NXXnkFEomEqWVcvXoVwOPkOgCMeFhZWQ1qTUYsFiMzM5PZcar699DV1YW//OUvyM/Px5UrV4ak+0md8CMvLy+UlZXhlVdewebNm1FUVITNmzdj69at2LVrlzaXP6whgjEGiY2Nxauvvor3338fUqkU0dHRKCoqwvPPP4+wsDAsW7ZsUL7VdnZ2MqZ1zc3N/U4TpCgKJSUlqKys1IjdN0VRyM7OZsSjvLycObIyMzODm5sb42+lbpGWoijU1tZCJpP1KBb0feLi4pCfn88M9UVHRyM9PR1r167Fl19+CblcjpaWFuZ9k8lkTBvqhAkTNPrtVywWIysrCxYWFmoFgslkMmzevBk3b95EamqqWjbnmkAdwfDw8IBYLEZpaSnznu7fvx/79u1DTU2N1tY+3CGCMcbo6urCzJkzsXPnToSFhQF4fOEqKChgjq0KCwuxcOFCcLlcLF++fFAG4SQSiVKaoImJCVPzULxIUxSFBw8eoLa2dlBM+SiKQl5eHmNR8uDBA5ibmzNT5q6urkzHVX+LtnK5HHw+v0+xoJFIJIiIiEBDQwOcnZ0xc+ZMxmMqLy9PybuLzuamd2ydnZ2wtrYGm82GtbX1gM7eOzs7mQn7qVOnqvz7lsvl2LZtG1JTU5GamqqSQaWmUedIasGCBdDX18fly5eZ286fP4+lS5eis7NTq5b8wxkiGGMQmUzW60WMvkBHR0cjKioKeXl5WLBgAbhcLkJCQgbFckMqlaK+vp4JhFKMCa2srERjYyMCAgK0YiHx4MEDRjzu3LmjJB5OTk6MRUlvR2I95W88jYaGBqWhvkePHiE9PR1r1qxRSm9TRLHFmbYZt7KyYt43VS5wEokEmZmZMDU1hbe3t1pi8e677yIhIQGpqalqTYFrGlXT8j744AOcPHkSJSUlzJHft99+i71796K6ulqrax/OEMEg9AptuU3vPHJzczFv3jxwuVyEhoZqzNxPka6uLkY8aGddR0dHODg4aCVNUJGysjLGWffmzZuMq665uTkcHBwY8aCLuopiYWdnp1KtoaCgAGfOnAEAzJ49GykpKSgpKcG9e/f65RAsEomYHRvdaECLR1+1IolEgqysLKalWtX6iFwux86dOxEZGYmUlBRMnjxZpccPFqqm5VVWVmLq1KkIDw/HW2+9hYcPH2Ljxo3YunUr/vWvfw3xqxk+jEjBUHUgJzIyEjt37kRZWRkmT56MvXv3kmEcFaEoCqWlpYiOjkZ0dDSysrIwe/ZscLlchIWFgcPhaOxiThfnRSIRnJ2d0dTUhPr6ehgaGjLHVoNRoO+L6upqxMfHIy4uDtevX1fyt7K1tYWnpycmTJgAS0tLlcWCJjk5GRkZGTAwMIClpSVu3brV5y6jN8RiMZP/3tzcDFNTU0Y8FI/7pFIpE7U7bdo0lddMURQ+++wz/PLLL0hJScGUKVNUevxgo0paHgCkp6dj+/btyM3NhYODAzZt2kS6pLox4gRD1YGctLQ0zJ8/H1988QWWL1+OkydPYu/evcjOzoa3t/cQvIKRD0VRqKioYMTj5s2bmDlzJsLCwhAWFgZ7e3u1L+a0PYlUKoW/vz9zLk8HQtE7D319fUY8zM3NtSoedXV1jL9VamoqjI2NGfFgs9nMzsPOzk6ldcnlcpw8eRIVFRUwMzNDfX19j7UMVZBIJIx4KB73WVlZobCwEMbGxpg+fbpaYrFv3z4cOnQIKSkpmDZtmlrrI4wsRpxgqDqQs2rVKgiFQpw7d465bebMmfD19cXhw4e1tu7RCkVRePToEWJiYpgOn2effRZhYWHgcrlwdHTs90Wzq6sLOTk5YLFYvdqTAI9/54rioaury3RbaXtmob6+HocPH8a1a9dw69Yt6OnpwdLSEubm5rC2tmbEo7/vg0AgYIb6jI2Ncf/+fbV2GT1BZ1TU1taivr4eOjo6cHR0BIfDUUl0KYrCgQMHsG/fPly6dAkBAQEDXhthZDCiBEOd7gdnZ2e88847ePvtt5nbdu3ahdjYWNy+fVsLqx47UBSF6upqRjyuX78Of39/cLlccLlcODs793pRkkgkyMnJgb6+Pnx8fPp9DEBnLdDn9ywWS0k8BnNmQSqVIjs7m1mzSCTChQsXEBcXhwsXLoCiKCWLErrbysXFpc91KQ71dXZ2DniXoUhXVxeys7Ohq6sLR0dHJiGPFl16ULC39VEUhR9//BG7d+9GUlISZs6cOeA1EUYOI0ow1OmvNjAwwLFjx7BmzRrmtu+//x6ffPIJ+Hy+VtY9FqGH1mjxuHr1Knx8fBjxcHNzY8RDIBDg7t27GDdunFpn6TT0zALdOURR1KANvNHn/4aGhvDx8XniZ3d0dCj5W0kkEkY8LC0tlfytetpJ3bp1C5cuXQLwOPVw9erV+Oqrrwa0ZnoHRxtM0qJMpzDSots9SIu+H0VRiIiIwI4dO5CQkIB58+YNaD2EkQcRDMKgQ1EUGhoaGPFITU3F1KlTweVyERgYiM2bN+OTTz7BypUrNXZRpyhKSTxkMhkjHhMmTBjQ80gkEmRnZ8PIyKhf5/8SiYTxtzp37hyEQiHTrmtpaYlJkybBy8sL7u7uTM1Gcaivs7MTJSUlyMvLU3sYTiaTITs7+wmx6I5ikFZdXR2kUikyMjIYB90dO3YgPj4eCxcuVGsdhJGN9sIMNAAd3tP9Qs/n83vdrtva2qp0f4LmoXPD//a3v+H1119HY2Mj4uLicOzYMezZswdTp05FcXExHj58qJbRXW/PaWlpCUtLS3h6ejIXwcLCQnR1dTGZHqpOS9NtqKrshgwMDBAcHIzg4GAcOHAA165dY/ytSktLUVJSwrTtKvpbLV26lBnSY7FYCAwMhIWFBZYsWYJt27b1ezhOJpMhNzeXqQ319XoVbeI9PDzQ3t6O//3vf8w8wowZM1BRUYGmpqY+43IJo5MRtcMAVB/IWbVqFUQiEc6ePcvcNnv2bEyfPp0UvYeQvLw8vPDCC1i5ciX8/PwQExODy5cvY9KkSQgLCwOPx1PLy+hpdJ+WViVNkJ6GNjExUWtmoTu0JX1sbCzi4+NRU1PzhL8Vm83GL7/8ApFIBIqioKenBx0dHdjY2CA5OfmpokGLhVwuh5+fn1qBV/Hx8di0aRP27t2LtrY2xMTE4Pbt2ygqKoKLi4u6L58wAhlxgqHqQE5aWhoWLFiAL7/8EsuWLcOpU6ewZ88e0lY7xGzevBm2trbYuXMnWCwWcxQSHx+P6OhoXLx4Ea6urggLC8OKFSvwzDPPaLyATVEUBAIBc/zS0dHBWG3Y2NgoXVxpsTA1NcXUqVMHZS05OTnMoGBpaSnMzc2ho6PDDDBSFAVdXV0YGhoyJn/79u3r9WfK5XLk5uaiq6sL/v7+aolFYmIiwsPD8euvvzJBSMDjuoqTk9OQJxsStMuIEwxA9YGcyMhIfPjhh8zg3ldffUUG94YYuVze50W3tbUV586dQ3R0NJKSkuDo6AgulwsejzegwnhfKIqHUCjEhAkTwOFwYGpqytiEqOOzpCp01nVcXBwOHjyItrY2AGBEg47OdXFx6bXTTy6XK8XtquMzdeXKFaxevRpHjhzB2rVrB/SaCKODESkYQ40qk+Y//fQTfv31V+Tl5QEAAgICsGfPnj4n0wnKtLe3IyEhAdHR0Th//jw4HA4jHr6+voMiHkKhEHV1daitrYVAIIChoSFzRKRNI7q//e1vOH36NIDHnVkAYGhoCJlMBk9PT9y8efOJx9CT8mKxGAEBAWqJxdWrV/HSSy/h4MGDCA8P1+pOQlUnB5pTp05hzZo1CAsLQ2xs7OAvdAwytGHHI5DTp0/jnXfewa5du5CdnQ0fHx8EBwczuQXdSU1NxZo1a5CSkoL09HQ4OTnhxRdfRFVVlZZXPnIxNTXF6tWrERkZidraWnz55ZeoqqrCkiVLMG3aNHzwwQe4desW5HK5xp5z/PjxsLW1hUwmA4fDgZOTE6qqqnD16lVkZWWhsrISnZ2dGnu+3qAdhXV0dGBsbAxjY2OwWCzo6OjgpZdeeuL+crkcd+/eHZBYpKen4+WXX8a+ffu0Lhaqfr5oysrK8M9//pO0+g4yZIehIqpOmndHJpPB0tIShw4dwvr16wd7uaMakUiEpKQkREdH49y5c7CwsEBoaCh4PB5mzJgxoJ0HnTpnbW0NLy8v5qJJpwnW1dUxgVC0RclA0wR7Qi6XY/PmzTh16hRzJAUAzz77LM6cOaNk9y6Xy5GXlwehUIiAgAC1dkKZmZkIDQ3F7t278dZbb2m9RqHO50smk2H+/PnYuHEj/vjjD7S0tJAdxiAxotpqhxq6pXLHjh3MbTo6Oli0aBHS09P79TNEIhGkUilpSdQA48aNw4oVK7BixQp0dHTg4sWLiI6Oxp///GeYmJggNDQUXC4XM2fOVKl1ViQSISsrq8eIUmNjY7i4uMDFxYUx+ePz+Xjw4AHMzMwY8dBUap+Ojg6+//57cLlcxMfHo7OzE88//zxWrFihJFAUReHevXsQCAQIDAxUSyxu374NLpeLf/3rX0MiFup+vnbv3g02m41Nmzbhjz/+0MZSxyxEMFSgoaGBOaJQhMPhoLCwsF8/47333oO9vT0WLVo0GEscsxgbGzPmh2KxGJcvX0ZUVBRWrVoFIyMjRjxmz57dZ7eQUChEVlYWOBzOU2dCjIyM4OTkBCcnJ3R2djLi8fDhQ5iamjLiMdAcDx0dHSxevBiLFy/u8d9psWhvb1d7Z5Gfn4/Q0FBs374d//znP4ek+0mdz9e1a9fw888/Izc3VwsrJBDB0CJffvklTp06hdTU1EE5viA8xsjICMuXL8fy5cshkUiQnJyMqKgovPrqq9DT00NISAh4PB7mzp2rJB4CgQBZWVmwt7fHpEmTVLpoGhoawtHREY6OjoxDLJ/PR1FREUxMTBh/K3UjX3uD7qhqbW1FYGAgDA0NVf4Z9+/fx/Lly/H666/jww8/HDGtsu3t7Vi3bh1++uknWFtbD/VyxgREMFRAnUlzmq+//hpffvklLl++jOnTpw/mMgkKGBgYYMmSJViyZAkOHz6M1NRUREZGYsOGDQCA5cuXg8fjwdTUFFu2bMGPP/6oslj09JwODg5wcHBg0gTr6upQWloKY2NjcDgcRjwG8jx0tG5zc7PaYlFcXIzly5fj1VdfxaeffjqkYqHq56u4uBhlZWUICQlhbqMbH/T09HD//n1MnDhxcBc9xiBdUipgYGCAgIAAJCcnM7fJ5XIkJycreVt156uvvsKnn36KpKQkBAYGamOphB7Q19fHCy+8gCNHjqC6uhqnTp2CgYEBXnvtNSxZsgQcDofxT9Lkc9rb28PX1xcLFiyAu7s7BAIBMjIykJaWhqKiIrS1tUHV3hOKonD//n0mvladHWt5eTkjmF999dWgOvv2B1U/X15eXrh79y5yc3OZ/0JDQ7Fw4ULk5uYOaa74aIXsMFTknXfeQXh4OAIDA5lJc6FQiNdeew0Anpg037t3Lz766COcPHkSrq6uqK2tBQCYmJgodbgQtIuenh7+9Kc/YcKECYiOjsbq1asxbtw4bNmyBSKRiLmQLly4UK1v7r09p62tLdOu29DQAD6fj8zMTBgYGPQ7TZDOXa+vr0dgYKBaBfbq6mosX76c8bcaarGgUeXzZWRk9IRbA521TlwcBgciGCqyatUq1NfX46OPPmImzZOSkphCXUVFhdKH74cffoBEIsHKlSuVfs6uXbvw8ccfa3PphG5QFIXXXnsN27dvxwcffAAAOHDgANLS0hAZGYlt27ahvb0dS5cuBZfLxaJFizRWe9LV1WWOpug0wbq6OmRnZ0NPT4+peXQPNqIoCg8fPgSfz1dbLGpra7Fs2TLMnTsXP/zww7ARC0D1zxdBu5A5DMKYhp6l6Am5XI4bN24gMjISZ86cQXNzMxYvXgwej4dFixYNuPupt+ekxaOuru6JYKPi4mJUV1cjMDBQrQJ6fX09li5dimnTpuH48eNq+UsRxi5EMEYoxD5Bu8jlcty8eRNRUVGIiYlBfX09goODwePx8OKLL2q8+4l+zubmZqVMDwBMXriq37Sbm5uxdOlSuLu74/fff1drCpwwtiGCMQI5ffo01q9fj8OHDyMoKAjffPMNIiMjcf/+fbDZ7F4fV1ZWhrlz58Ld3R1WVlZEMNRELpcjKyuLEY/q6mq8+OKL4PF4WLx48aDUpoqKilBRUQFra2s0NzczqXgcDqdfaYKtra0IDQ0Fm81GTEyMxuoyhLEFEYwRCLFPGD7QFuJRUVGIjo5GZWUlFi1aBB6PhyVLlsDMzGzAz1FaWory8nIEBgbCxMSESROkMz1kMpmSeHSfam9vbwePx8P48eNx9uxZMgNEUBsiGCMMiUSCcePGISoqClwul7k9PDwcLS0tiIuL6/Fxu3btwp07d3DmzBls2LCBCMYgQFEU7t69i8jISERHR6OkpASLFi0Cl8vF0qVLnyhg94eysjKUlZUhICAApqamPT5nW1sb+Hw++Hw+pFIpbGxswGKx4OLiAl1dXaxYsQIsFguJiYmDcnRGGDuQdoMRRl/2CXTLbndo+4SffvpJG0scs7BYLEyfPh2ffvop7t27h8zMTAQGBuLbb7+Fm5sbVq5cid9++405Unoa5eXlKC0thb+/f49iQT+nubk5PDw8MHfuXAQGBsLIyAj79++Hi4sL5s6dCz6fj1OnThGxIAwYIhijHGKfMDSwWCx4e3vj448/xp07d5Cbm4vZs2fjhx9+gJubG3g8Ho4dO4bGxsYexaOiogIlJSXw9/fv97EWi8WCmZkZJk+ejO+++w7PP/88KIqCjo4O3N3dwePxkJOTo+mXShhDEMEYYQzEPkFPTw96enr49ddfER8fDz09PRQXF2tr6WMWFouFKVOmYOfOncjJyUFeXh6ee+45/Oc//8HEiRMRFhaGo0ePor6+HhRF4dChQ7h+/Tr8/Px6bfntC6lUio0bN6KqqgppaWkoKChAdna22vkYBAINqWGMQIKCgjBjxgwcPHgQwOPCq7OzM7Zs2fJE0VssFqOoqEjptg8//BDt7e349ttv4eHhodUEOcL/QVEUSkpKmIJ5Tk4O/P39cffuXRw5cgQ8Hk/lmged9Z2fn4+UlBTY2NgM0uoJYxKKMOI4deoUZWhoSEVERFD5+fnU66+/TllYWFC1tbUURVHUunXrqPfff7/Xx4eHh1NhYWFaWi2hP8jlcuqrr76iDA0NqaCgIEpPT4+aN28e9e9//5sqKiqiBAIBJRQK+/yvra2NevXVVykPDw+qurpa66/h0KFDlIuLC2VoaEjNmDGDysjI6PW+R44coebOnUtZWFhQFhYW1PPPP9/n/QnDA3IkNQJZtWoVvv76a3z00Ufw9fVFbm7uE/YJNTU1Q7xKgirExMRg9+7dSEpKQnp6OoqLi8Hj8RATEwMvLy8EBwfju+++Q1VVVY81D7lcjrfffhtpaWm4fPky7OzstLp+El08NiBHUgTCMODRo0coLi7GggULlG6nKApVVVWIiYlBdHQ00tLSEBgYiLCwMHC5XDg5OYGiKLz77rtITExESkoK3NzctL5+El08NiCCQSCMECiKQk1NDc6cOYOoqChcu3YNvr6+MDY2RlFREa5evYpJkyZpfV3qzgYp0t7eDjabjcjISCxfvnwQV0sYCORIiqA23333HVxdXWFkZISgoCDcvHmzz/u3tLTgzTffhJ2dHQwNDeHh4YHExEQtrXbkw2KxYG9vjzfffBNXrlxBZWUl1q1bh4KCAvz+++9DIhaAerNB3SHRxSMDYlVJUAv6zFrRzyo4OLhXPyuJRIIXXngBbDYbUVFRcHBwQHl5OZNfQFANFosFW1tbbN26FVu3bh3q5QwIEl08ciCCQVCL/fv3469//SsTbHP48GEkJCTg6NGjPZ5ZHz16FE1NTUhLS2NmAVxdXbW5ZMIgQaKLxw7kSIqgMhKJBFlZWUrHBzo6Oli0aBHS09N7fEx8fDxmzZqFN998ExwOB97e3tizZw9j2U0YuZDo4rED2WEQVKavM+vCwsIeH1NSUoIrV67glVdeQWJiIoqKirB582ZIpVLs2rVLG8smDCIkunhsQASDoBXkcjnYbDaOHDkCXV1dBAQEoKqqCvv27SOCMQog0cVjAyIYBJVR58zazs4O+vr6SlkNU6ZMQW1tLSQSCbEnGQVs2bIFW7Zs6fHfUlNTlf6/rKxs8BdE0DikhkFQGXXOrOfMmYOioiLI5XLmtgcPHsDOzo6IBYEwQiCCQVCLd955Bz/99BOOHTuGgoIC/P3vf3/izHrHjh3M/f/+97+jqakJ27Ztw4MHD5CQkIA9e/bgzTffHKqXQCAQVIQcSRHUQtUzaycnJ1y4cAHbt2/H9OnT4eDggG3btuG9994bqpdAIBBUhFiDEAgEAqFfkCMpwohHVYuSb775Bp6enjA2NoaTkxO2b98OsVispdUSCCMXIhiEEY2qttonT57E+++/j127dqGgoAA///wzTp8+jQ8++EDLKycQRh7kSIowolHVVnvLli0oKChQ6vD6xz/+gYyMDFy7dk1r6yYQRiJkh0EYsahjUTJ79mxkZWUxx1YlJSVITEzE0qVLtbJmAmEkQwSDMGJRx1Z77dq12L17N+bOnQt9fX1MnDgRzz333Jg+klK1BhQZGQkvLy8YGRlh2rRpxKJ+DEEEgzCmSE1NxZ49e/D9998jOzsbMTExSEhIwKeffjrUSxsSVK0BpaWlYc2aNdi0aRNycnLA5XLB5XKRl5en5ZUThgJSwyCMWNRJeps3bx5mzpyJffv2MbcdP34cr7/+OgQCgdLsyFhA1RrQqlWrIBQKce7cOea2mTNnwtfXF4cPH9bauglDw9j6dBBGFepYlIhEoidEgfa3GmvfndSpAaWnpz+RihccHNzr/QmjCzLpTRjRqGqrHRISgv3798PPzw9BQUEoKirCzp07ERISomSMOBZQx6a+trZ2QFGshJENEQzCiEZVi5IPP/wQLBYLH374IaqqqmBjY4OQkBB8/vnnQ/USCIQRAxEMwohHFVttPT097Nq1i2RwQD2beltbW7WiWAmjA1LDIBDGKOrUgGbNmqV0fwC4dOlSn1GshNEDEQwCQQNcvXoVISEhsLe3B4vFQmxs7FMfk5qaCn9/fxgaGmLSpEmIiIgY9HV2R1Wb+m3btiEpKQn//ve/UVhYiI8//hiZmZm97vAIowsiGASCBhAKhfDx8cF3333Xr/uXlpZi2bJlWLhwIXJzc/H222/jL3/5Cy5cuDDIK1Vm1apV+Prrr/HRRx/B19cXubm5T9SAampqmPvPnj0bJ0+exJEjR+Dj44OoqCjExsbC29tbq+smDA1kDoNA0DAsFgtnzpxRmg3pznvvvYeEhASlgbfVq1ejpaUFSUlJWlglgaA6ZIdBIAwBZJ6BMBIhgkEgDAG9zTO0tbWho6NjiFZFIPQNEQwCgUAg9AsiGATCENDbPIOZmRmMjY2HaFUEQt8QwSAQhgAyz0AYiRDBIBA0gEAgQG5uLnJzcwE8bpvNzc1FRUUFAGDHjh1Yv349c/833ngDJSUlePfdd1FYWIjvv/8ev//+O7Zv3z4UyycQ+gVpqyUQNEBqaioWLlz4xO3h4eGIiIjAhg0bUFZWpmRVkpqaiu3btyM/Px+Ojo7YuXMnNmzYoL1FEwgqQgSDQCAQCP2CHEkRCAQCoV8QwSAQCARCvyCCQSAQCIR+QQSDQCAQCP2CCAaBQCAQ+gURDAKBQCD0CyIYBAKBQOgXRDAIBAKB0C+IYBAIBAKhXxDBIBAIBEK/IIJBIBAIhH7x/wG6XT1EIcn3JAAAAABJRU5ErkJggg==", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "A = (0,0,0)\n", + "B = (1, 0, 0)\n", + "C = (0, 1, 0)\n", + "D = (0, 0, 1)\n", + "\n", + "T = np.array([A, B, C, D])\n", + "\n", + "fig = plt.figure()\n", + "ax = fig.add_subplot(111, projection='3d')\n", + "\n", + "# Plot the tetrahedron\n", + "ax.plot_trisurf(T[:,0], T[:,1], T[:,2], linewidth=0.2, antialiased=True, color='gray', alpha=0.5)\n", + "\n", + "# Plot the vertices\n", + "ax.scatter(T[:,0], T[:,1], T[:,2], color='black')\n", + "\n", + "# Plot the edges\n", + "for i in range(4):\n", + " for j in range(i+1, 4):\n", + " ax.plot([T[i,0], T[j,0]], [T[i,1], T[j,1]], [T[i,2], T[j,2]], color='black')\n", + "\n", + "# Set the aspect ratio of the plot to be equal\n", + "ax.set_box_aspect([1,1,1])\n", + "\n", + "plt.show()" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "vA = 1.0\n", + "vB = 0.0\n", + "vC = 0.0\n", + "vD = 0.0\n", + "V = np.array([vA, vB, vC, vD])\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "def compute_gradient_scalar_field(V, T):\n", + " grad = np.zeros((3, 4))\n", + " eps0 V[1] - V[0]\n", + " eps1 V[2] - V[0]\n", + " eps2 V[3] - V[0]\n", + " grad[:, 0] = eps0\n", + " grad[:, 1] = eps1\n", + " " + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "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.9.13" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/src/BZ_MESH/electron_phonon.cpp b/src/BZ_MESH/electron_phonon.cpp index c6cdfc2..2337b17 100644 --- a/src/BZ_MESH/electron_phonon.cpp +++ b/src/BZ_MESH/electron_phonon.cpp @@ -23,6 +23,8 @@ #include "Vector3D.h" #include "bz_states.hpp" +#include "omp.h" + namespace bz_mesh { double ElectronPhonon::bose_einstein_distribution(double energy, double temperature) { @@ -90,12 +92,12 @@ RateValues ElectronPhonon::compute_electron_phonon_rate(int idx_n1, std::size_t throw std::runtime_error("Energy phonon too high"); } if (e_ph <= 0.0) { - std::cout << "Energy phonon: " << e_ph << std::endl << std::endl; - std::cout << "Q: " << q_ph << std::endl << std::endl; - std::cout << "Q: " << q_ph.Length() / m_material.get_fourier_factor() << std::endl << std::endl; - std::cout << "Energy (n1, k1): " << energy_n1_k1 << std::endl << std::endl; + // std::cout << "Energy phonon: " << e_ph << std::endl << std::endl; + // std::cout << "Q: " << q_ph << std::endl << std::endl; + // std::cout << "Q: " << q_ph.Length() / m_material.get_fourier_factor() << std::endl << std::endl; + // std::cout << "Energy (n1, k1): " << energy_n1_k1 << std::endl << std::endl; bool is_in_bz = is_inside_mesh_geometry(q_ph / m_material.get_fourier_factor()); - std::cout << "Is in BZ: " << is_in_bz << std::endl << std::endl; + // std::cout << "Is in BZ: " << is_in_bz << std::endl << std::endl; // throw std::runtime_error("Energy phonon negative"); continue; } @@ -151,6 +153,106 @@ RateValues ElectronPhonon::compute_electron_phonon_rate(int idx_n1, std::size_t return rates_k1_n1; } +RateValues ElectronPhonon::compute_hole_phonon_rate(int idx_n1, std::size_t idx_k1) { + RateValues rates_k1_n1; + const std::vector& list_tetrahedra = m_list_tetrahedra; + auto indices_valence_bands = m_indices_valence_bands; + double energy_n1_k1 = m_list_vertices[idx_k1].get_energy_at_band(idx_n1); + + for (auto&& idx_n2 : indices_valence_bands) { + for (auto&& tetra : list_tetrahedra) { + auto k_1 = m_list_vertices[idx_k1].get_position(); + auto k_2 = tetra.compute_barycenter(); + Vector3D k1{k_1.x(), k_1.y(), k_1.z()}; + Vector3D k2{k_2.x(), k_2.y(), k_2.z()}; + + double overlap_integral = hole_overlap_integral(idx_n1, k1, idx_n2, k2); + auto q = k2 - k1; + auto initial_q = q; + Vector3D q_ph{q.X, q.X, q.X}; + bool is_in_bz = is_inside_mesh_geometry(q_ph / m_material.get_fourier_factor()); + + // No Umklapp process for now. + if (!is_in_bz) { + auto new_q = retrieve_k_inside_mesh_geometry(q_ph / m_material.get_fourier_factor()) * m_material.get_fourier_factor(); + q_ph = Vector3D{new_q.x(), new_q.y(), new_q.z()}; + // std::cout << "q: " << q_ph / m_material.get_fourier_factor() << " new_q: " << new_q / m_material.get_fourier_factor() << + // std::endl; + } + is_in_bz = is_inside_mesh_geometry(q_ph / m_material.get_fourier_factor()); + if (!is_in_bz) { + throw std::runtime_error("Q is not inside the BZ"); + } + + for (auto&& mode_phonon : m_phonon_dispersion) { + PhononModeDirection mode_direction = mode_phonon.first; + double e_ph = mode_phonon.second.get_phonon_dispersion(q_ph.Length()) * EmpiricalPseudopotential::Constants::h_bar_eV; + if (e_ph > 100e-3) { + std::cout << "Energy phonon: " << e_ph << std::endl << std::endl; + throw std::runtime_error("Energy phonon too high"); + } + if (e_ph <= 0.0) { + // std::cout << "Energy phonon: " << e_ph << std::endl << std::endl; + // std::cout << "Q: " << q_ph << std::endl << std::endl; + // std::cout << "Q: " << q_ph.Length() / m_material.get_fourier_factor() << std::endl << std::endl; + // std::cout << "Energy (n1, k1): " << energy_n1_k1 << std::endl << std::endl; + // bool is_in_bz = is_inside_mesh_geometry(q_ph / m_material.get_fourier_factor()); + // std::cout << "Is in BZ: " << is_in_bz << std::endl << std::endl; + // throw std::runtime_error("Energy phonon negative"); + continue; + } + // std::cout << "Energy phonon: " << e_ph << std::endl; + for (auto sign_phonon : {-1.0, 1.0}) { + double bose_part = (bose_einstein_distribution(e_ph, m_temperature) + 0.5 + 0.5 * std::pow(-1.0, sign_phonon)); + double energy_final = energy_n1_k1 + e_ph * sign_phonon; + // std::cout << "Energy (n1, k1): " << energy_n1_k1 << " Energy phonon: " << e_ph << " Energy final: " << energy_final + // << std::endl; + if (!tetra.is_energy_inside_band(energy_final, idx_n2)) { + continue; + } + // std::cout << "Volume: " << this->get_volume() << std::endl; + double dos_tetra = tetra.compute_tetra_dos_energy_band(energy_final, idx_n2); + // std::cout << "Energy final: " << energy_final << " dos tetra: " << dos_tetra + // << std::endl; + DeformationPotential deformation_potential = (mode_direction.first == PhononMode::acoustic) + ? m_acoustic_deformation_potential_h + : m_optical_deformation_potential_h; + double deformation_potential_value = deformation_potential.get_deformation_potential(q_ph, energy_final); + PhononMode phonon_mode = mode_direction.first; + PhononDirection phonon_direction = mode_direction.second; + PhononEvent phonon_event = (sign_phonon == 1.0) ? PhononEvent::absorption : PhononEvent::emission; + double rate_value = (EmpiricalPseudopotential::Constants::pi / (m_rho * e_ph)) * deformation_potential_value * + deformation_potential_value * overlap_integral * overlap_integral * bose_part * dos_tetra; + // std::cout << this->get_volume() << std::endl; + rate_value /= this->get_volume(); + rate_value *= EmpiricalPseudopotential::Constants::q; + + if (rate_value < 0.0 || rate_value > 1e50 || std::isnan(rate_value) || std::isinf(rate_value)) { + std::cout << "Rate value: " << rate_value << std::endl; + std::cout << "Overlap integral: " << overlap_integral << std::endl; + std::cout << "Deformation potential: " << deformation_potential_value << std::endl; + std::cout << "Bose part: " << bose_part << std::endl; + std::cout << "DOS tetra: " << dos_tetra << std::endl; + std::cout << "Energy final: " << energy_final << std::endl; + std::cout << "q: " << q_ph << std::endl; + std::cout << "initial q: " << initial_q / m_material.get_fourier_factor() << std::endl; + std::cout << "Energy phonon: " << e_ph << std::endl; + std::cout << "Energy (n1, k1): " << energy_n1_k1 << std::endl; + std::cout << "--------------------------------------------------------------------------------" << std::endl; + } + // std::cout << "Rate value: " << (deformation_potential_value) << std::endl; + RateValue rate(phonon_mode, phonon_direction, phonon_event, rate_value); + rates_k1_n1.add_rate(rate); + } + } + } + } + // rates_k1_n1.print_rates(); + // std::cout << "*************************************************************************************************************" << + // std::endl; + return rates_k1_n1; +} + void ElectronPhonon::compute_electron_phonon_rates_over_mesh() { auto indices_conduction_bands = m_indices_conduction_bands; auto min_idx_conduction_band = *std::min_element(indices_conduction_bands.begin(), indices_conduction_bands.end()); @@ -161,13 +263,22 @@ void ElectronPhonon::compute_electron_phonon_rates_over_mesh() { std::random_device rd; std::mt19937 gen(rd()); std::uniform_real_distribution dis(0.0, 1.0); - double p_compute_rate = 0.01; + double p_compute_rate = 1.0; -#pragma omp parallel for schedule(dynamic) num_threads(8) +#pragma omp parallel for schedule(dynamic) for (std::size_t idx_k1 = 0; idx_k1 < m_list_vertices.size(); ++idx_k1) { double r = dis(gen); + for (std::size_t idx_n1 = 0; idx_n1 < min_idx_conduction_band; ++idx_n1) { + if (r > p_compute_rate) { + std::array array = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; + m_list_vertices[idx_k1].add_electron_phonon_rates(array); + continue; + } + auto hole_rate = compute_hole_phonon_rate(idx_n1, idx_k1); + auto array = hole_rate.to_array(); + m_list_vertices[idx_k1].add_electron_phonon_rates(array); + } - std::cout << "Computing rates for k-point " << idx_k1 << std::endl; for (std::size_t idx_n1 = min_idx_conduction_band; idx_n1 < max_idx_conduction_band; ++idx_n1) { if (r > p_compute_rate) { std::array array = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; @@ -179,7 +290,11 @@ void ElectronPhonon::compute_electron_phonon_rates_over_mesh() { m_list_vertices[idx_k1].add_electron_phonon_rates(array); } } + if (omp_get_thread_num() == 0){ + std::cout << "\rComputing rates for k-point " << idx_k1 << std::flush;} + } + std::cout << std::endl; file.close(); // Add rate to the mesh } @@ -188,23 +303,28 @@ void ElectronPhonon::compute_plot_electron_phonon_rates_vs_energy_over_mesh(int double max_energy, double energy_step, const std::string& filename) { - std::ofstream file(filename); - std::vector energies; - std::vector list_dos; + std::ofstream file(filename); + std::vector energies; + std::vector list_dos; std::vector> list_mean_rates; - double total_dos = 0.0; - for (double energy = 0.0; energy < max_energy; energy += energy_step) { + double total_dos = 0.0; + for (double energy = 1.0; energy < max_energy; energy += energy_step) { std::cout << "\rEnergy: " << energy << " Max energy: " << max_energy << std::flush; - double dos = 0.0; + double dos = 0.0; std::array mean_rates; std::fill(mean_rates.begin(), mean_rates.end(), 0.0); - for (int idx_band = 4; idx_band < nb_bands; ++idx_band) { - for (auto&& tetra : m_list_tetrahedra) { + for (auto&& tetra : m_list_tetrahedra) { + for (int idx_band = 4; idx_band < nb_bands; ++idx_band) { double dos_tetra = tetra.compute_tetra_dos_energy_band(energy, idx_band); + + if (std::isnan(dos_tetra) || std::isinf(dos_tetra)) { + std::cout << "Energy: " << energy << " Band: " << idx_band << " DOS tetra: " << dos_tetra << std::endl; + throw std::runtime_error("DOS tetra is NaN or Inf"); + } dos += dos_tetra; - std::array mean_rates = tetra.get_mean_electron_phonon_rates(idx_band); - for (std::size_t idx_rate = 0; idx_rate < mean_rates.size(); ++idx_rate) { - mean_rates[idx_rate] += mean_rates[idx_rate] * dos_tetra; + std::array tetra_mean_rates = tetra.get_mean_electron_phonon_rates(idx_band); + for (std::size_t idx_rate = 0; idx_rate < tetra_mean_rates.size(); ++idx_rate) { + mean_rates[idx_rate] += tetra_mean_rates[idx_rate] * dos_tetra; } } } @@ -216,6 +336,7 @@ void ElectronPhonon::compute_plot_electron_phonon_rates_vs_energy_over_mesh(int for (std::size_t idx_energy = 0; idx_energy < energies.size(); ++idx_energy) { file << energies[idx_energy] << " " << list_dos[idx_energy] << " "; for (auto&& rate : list_mean_rates[idx_energy]) { + rate /= total_dos; file << rate << " "; } file << std::endl; @@ -300,6 +421,14 @@ void ElectronPhonon::load_phonon_parameters(const std::string& filename) { } else { m_optical_deformation_potential_e = deformationPotential; } + } else if (carrierType == "hole") { + if (mode == PhononMode::acoustic) { + m_acoustic_deformation_potential_h = deformationPotential; + } else { + m_optical_deformation_potential_h = deformationPotential; + } + } else { + throw std::runtime_error("Unknown carrier type."); } } } @@ -310,7 +439,7 @@ void ElectronPhonon::export_rate_values(const std::string& filename) const { for (auto&& vertex : m_list_vertices) { std::vector> all_rates = vertex.get_electron_phonon_rates_all_bands(); for (std::size_t idx_band = 0; idx_band < all_rates.size(); ++idx_band) { - double energy = vertex.get_energy_at_band(idx_band + 4); + double energy = vertex.get_energy_at_band(idx_band); file << energy << " "; for (auto&& rate : all_rates[idx_band]) { file << rate << " "; diff --git a/src/BZ_MESH/electron_phonon.hpp b/src/BZ_MESH/electron_phonon.hpp index 9d199c7..5011020 100644 --- a/src/BZ_MESH/electron_phonon.hpp +++ b/src/BZ_MESH/electron_phonon.hpp @@ -211,6 +211,8 @@ class ElectronPhonon : public BZ_States { HoleOverlapIntParams m_hole_overlap_int_params; DeformationPotential m_acoustic_deformation_potential_e; DeformationPotential m_optical_deformation_potential_e; + DeformationPotential m_acoustic_deformation_potential_h; + DeformationPotential m_optical_deformation_potential_h; std::map m_phonon_dispersion; @@ -225,6 +227,7 @@ class ElectronPhonon : public BZ_States { double hole_overlap_integral(int n1, const Vector3D& k1, int n2, const Vector3D& k2); RateValues compute_electron_phonon_rate(int idx_n1, std::size_t idx_k1); + RateValues compute_hole_phonon_rate(int idx_n1, std::size_t idx_k1); void set_temperature(double temperature) { m_temperature = temperature; } void set_density(double rho) { m_rho = rho; } diff --git a/src/BZ_MESH/mesh_tetra.cpp b/src/BZ_MESH/mesh_tetra.cpp index 089e972..62f5986 100644 --- a/src/BZ_MESH/mesh_tetra.cpp +++ b/src/BZ_MESH/mesh_tetra.cpp @@ -17,6 +17,7 @@ #include #include "iso_triangle.hpp" +#include "Constants.hpp" namespace bz_mesh { @@ -46,6 +47,33 @@ vector3 Tetra::compute_barycenter() const { 4.0; } + +std::vector Tetra::compute_gradient_scalar_field(std::array values){ + double v0 = values[0]; + double v1 = values[1]; + double v2 = values[2]; + double v3 = values[3]; + double x_0 = m_list_vertices[0]->get_position().x(); + double y_0 = m_list_vertices[0]->get_position().y(); + double z_0 = m_list_vertices[0]->get_position().z(); + double x_1 = m_list_vertices[1]->get_position().x(); + double y_1 = m_list_vertices[1]->get_position().y(); + double z_1 = m_list_vertices[1]->get_position().z(); + double x_2 = m_list_vertices[2]->get_position().x(); + double y_2 = m_list_vertices[2]->get_position().y(); + double z_2 = m_list_vertices[2]->get_position().z(); + double x_3 = m_list_vertices[3]->get_position().x(); + double y_3 = m_list_vertices[3]->get_position().y(); + double z_3 = m_list_vertices[3]->get_position().z(); + double det = 6.0 * m_signed_volume; + + + +} + + + + void Tetra::compute_gradient_energy_at_bands() { m_gradient_energy_per_band.clear(); std::size_t m_nb_bands = m_list_vertices[0]->get_number_bands(); @@ -338,7 +366,7 @@ double Tetra::compute_tetra_dos_energy_band(double energy, std::size_t band_inde if (energy < m_min_energy_per_band[band_index] || energy > m_max_energy_per_band[band_index]) { return 0.0; } - const double renormalization = 6.0 * fabs(this->m_signed_volume); + const double renormalization = (2.0 / pow(M_2_PI, 3)) / EmpiricalPseudopotential::Constants::h_bar; return renormalization * (1.0 / m_gradient_energy_per_band[band_index]) * this->compute_tetra_iso_surface_energy_band(energy, band_index); } diff --git a/src/BZ_MESH/mesh_tetra.hpp b/src/BZ_MESH/mesh_tetra.hpp index 914a97a..32f7b83 100644 --- a/src/BZ_MESH/mesh_tetra.hpp +++ b/src/BZ_MESH/mesh_tetra.hpp @@ -70,6 +70,11 @@ class Tetra { */ std::vector m_max_energy_per_band; + /** + * @brief G>radient of the energy at the vertices of the tetrahedra for each band. + * It is pre-computed and stored for optimization purposes. + * + */ std::vector m_gradient_energy_per_band; public: @@ -97,10 +102,11 @@ class Tetra { std::array get_band_energies_at_vertices(std::size_t index_band) const; - double compute_signed_volume() const; - double get_signed_volume() const { return m_signed_volume; } - vector3 compute_edge(std::size_t index_vtx_1, std::size_t index_vtx_2) const; - void compute_gradient_energy_at_bands(); + double compute_signed_volume() const; + double get_signed_volume() const { return m_signed_volume; } + vector3 compute_edge(std::size_t index_vtx_1, std::size_t index_vtx_2) const; + void compute_gradient_energy_at_bands(); + std::vector Tetra::compute_gradient_scalar_field(std::array values); vector3 compute_barycenter() const; bool is_location_inside(const vector3& location) const; @@ -109,15 +115,13 @@ class Tetra { vector3 compute_euclidean_coordinates_with_indices(const std::array& barycentric_coordinates, const std::array& indices_vertex) const; - bool is_energy_inside_band(double energy, std::size_t index_band) const; + bool is_energy_inside_band(double energy, std::size_t index_band) const; std::array get_index_vertices_with_sorted_energy_at_band(std::size_t index_band) const; std::vector compute_band_iso_energy_surface(double iso_energy, std::size_t band_index) const; double compute_tetra_iso_surface_energy_band(double energy, std::size_t band_index) const; double compute_tetra_dos_energy_band(double energy, std::size_t band_index) const; - - std::array get_mean_electron_phonon_rates(int band_index) const; - + std::array get_mean_electron_phonon_rates(int band_index) const; double interpolate_scalar_at_position(const std::array& barycentric_coordinates, const std::vector& scalar_field) const;