{ "cells": [ { "attachments": {}, "cell_type": "markdown", "id": "349d0370", "metadata": {}, "source": [ "# The shortest path problem using QUBO\n", "This notebook implements a Perceval code for finding the shortest path on a _directed, weighted graph_ using the Quantum Unconstrained Binary\n", "Optimization (QUBO) model. It is mainly implementing the algorithm in https://arxiv.org/pdf/2112.09766.pdf and represents the work done during the 2022 LOQCathon. \n", "\n", "Authors: Beata Zjawin, Benjamin Pointard, Nathan Claudet, Noé Delorme, Rina Ismailati.\n", "\n", "The task is to find the shortest path from start to finish (such that the sum of the weights is minimized) on a simple 5-edge graph:\n", "\n", "![Graph_5_edges](../_static/img/graph_5edges.png)\n", "\n", "Start with importing the necessary packages.\n" ] }, { "cell_type": "code", "execution_count": 11, "id": "ab150578", "metadata": {}, "outputs": [], "source": [ "import math\n", "\n", "import perceval as pcvl\n", "from perceval.components import PS, BS, GenericInterferometer\n", "import numpy as np\n", "from scipy.optimize import minimize\n", "import matplotlib.pyplot as plt\n", "plt.rcdefaults()" ] }, { "attachments": {}, "cell_type": "markdown", "id": "7c25b79a-d88c-4bf2-9398-93c1434d06cd", "metadata": { "tags": [] }, "source": [ "## Implementation in Perceval\n", "To find the shortest path, first we need to compute the objective function of the weighted graph. The constraints specified by the shortest path problem can be written in the form of quadratic penalty functions, as required by the QUBO formalism [1]:\n", "\n", "\\begin{align*}\n", "H_s&=\\Big(\\sum_{j} x_{s,j} - \\sum_{k} x_{k,s} - 1\\Big)^2, \\\\\n", "H_f&=\\Big(\\sum_{j} x_{f,j} - \\sum_{k} x_{k,f} + 1\\Big)^2, \\\\\n", "H_i&=\\Big(\\sum_{j} x_{i,j} - \\sum_{k} x_{k,i}\\Big)^2. \\\\\n", "\\end{align*}\n", "\n", "Here, $s$ and $f$ correspond to the start and finish node, respectively. Together with the original objective function:\n", "\n", "\\begin{equation*}\n", "H_c= \\sum_{(i,j)\\in edges} c_{i,j} x^{2}_{i,j},\n", "\\end{equation*}\n", "\n", "the constraints form the final QUBO objective function:\n", "\n", "\\begin{equation*}\n", "H= \\alpha \\Big(H_s + H_f + \\sum_{i \\notin \\{s,f\\}} H_i\\Big) + H_c,\n", "\\end{equation*}\n", "\n", "where $\\alpha$ is a scaling factor satisfying $\\alpha > \\sum_{(i,j)\\in edges} c_{i,j}$.\n", "\n", "\n", "\n", "For our simple example, the objective function is given by (see the Appendix to generate Hamiltonians from arbitrary adjacency matrices):\n" ] }, { "cell_type": "code", "execution_count": 12, "id": "dcde3aa3-da33-4d7d-9846-7efa216550c6", "metadata": {}, "outputs": [], "source": [ "H1 = [[2., 32., -32., -32., 32., 0.],\n", " [0., 1., 32., 0., -32., -32.],\n", " [0., 0., 35., 32., -64., -32.],\n", " [0., 0., 0., 2., -32., 32.],\n", " [0., 0., 0., 0., 35., 32.],\n", " [0., 0., 0., 0., 0., 4.]]" ] }, { "attachments": {}, "cell_type": "markdown", "id": "6d0d9f59", "metadata": {}, "source": [ "Next, we implement a variational algorithm to find the optimal path using the QUBO model.\n", "\n", "Following the approach introduced in Ref. [2], the output of the simulation, i.e. the Fock states $|n> = |n_1,...,n_M>$, can be mapped into bit strings $b=(b_1,...b_M)$ encoding possible paths on the graph by applying the parity function:\n", "\n", "\\begin{equation*}\n", "\\rho_j: b_i^{(j)} = \\text{mod}[n_i,2] \\oplus j.\n", "\\end{equation*}\n", "\n", "For each evaluation, we need to test for $j=0$, $j=1$, $n=M$ and $n=M−1$, as it is not possible to know a priori which configuration is the optimal one. Here, the number of modes, $M$, corresponds to the number of possible paths and $n$ is the number of photons.\n", "\n", "In our example, $b=(b_{sa},b_{sb},b_{ab},b_{af},b_{ba},b_{bf})$, with $b_i \\in \\{0,1\\}$, where subscripts correspond to graph edges.\n", "\n", "In the cell below, we set up the functions necessary for the simulation, such as implementation of the parity function and sampling of the circuit." ] }, { "cell_type": "code", "execution_count": 13, "id": "3441d4b4", "metadata": {}, "outputs": [], "source": [ "def parify_samples(samples, j):\n", " \"\"\"apply the parity function to the samples\"\"\"\n", " def _parity(output, j):\n", " m = len(output)\n", " parity = [0]*m\n", " for i in range(0,m):\n", " parity[i] = (output[i] + j) % 2\n", " return pcvl.BasicState(parity)\n", " \n", " new_samples = pcvl.BSCount()\n", " for idx,sample in enumerate(samples):\n", " new_sample = _parity(sample,j)\n", " if new_sample in new_samples:\n", " new_samples[_parity(sample,j)] += samples[sample]\n", " else:\n", " new_samples[_parity(sample,j)] = samples[sample]\n", " return new_samples\n", " \n", "def set_parameters_circuit(parameters_circuit, values): \n", " \"\"\"set values of circuit parameters\"\"\"\n", " for idx, p in enumerate(parameters_circuit):\n", " parameters_circuit[idx].set_value(values[idx])\n", "\n", "def compute_samples(circuit, input_state, nb_samples, j):\n", " \"\"\"sample from the circuit\"\"\"\n", " p = pcvl.Processor(\"SLOS\", circuit)\n", " p.with_input(input_state)\n", " p.min_detected_photons_filter(0)\n", "\n", " sampler = pcvl.algorithm.Sampler(p) \n", " samples = sampler.sample_count(nb_samples)['results']\n", " \n", " return parify_samples(samples,j)" ] }, { "attachments": {}, "cell_type": "markdown", "id": "bc3125ec-1fe5-4df4-a5bf-0cb9a7de5df0", "metadata": {}, "source": [ "The central part of solving the shortest path problem with Perceval is initializing a universal circuit and optimizing its parameters. We will work with the following generic interferometer:\n", "\n", "![generic_6photons_interferometer](../_static/img/generic_6photons_interferometer.svg)" ] }, { "attachments": {}, "cell_type": "markdown", "id": "2bb3985b-3105-4d5b-893f-8d89f41d0051", "metadata": {}, "source": [ "We sample the circuit for all configurations of $j$ and $n$. In this example, we set all the initial parameters to $\\pi$, which is close to the optimal solution $-$ a good initial guess can dramatically improve simulation time. The loss function is given by\n", "\n", "\\begin{equation}\n", "E(j;\\theta,\\psi) = \\sum_{b^{(j)}} \\beta_{b^{(j)}} ,\n", "\\end{equation}\n", "\n", "where $\\theta$ and $\\psi$ stand for the parameters of the linear interferometer: beam-splitters and phase-shifters, respectively. For optimization, we use the Powell minimisation algorithm from the scipy.optimize package." ] }, { "cell_type": "code", "execution_count": 14, "id": "d0d7c887", "metadata": {}, "outputs": [], "source": [ "def test_configuration(circuit, nb_modes, j, n, H, nb_samples):\n", " \"\"\"output the samples for a given configuration of j and n (inludes minimisation of the loss function)\"\"\"\n", " parameters_circuit = circuit.get_parameters()\n", "\n", " input_state = pcvl.BasicState([1]*nb_modes)\n", " if n!=nb_modes:\n", " input_state = pcvl.BasicState([1]*(nb_modes-1)+[0])\n", "\n", " def loss(parameters):\n", " set_parameters_circuit(parameters_circuit, parameters)\n", " samples = compute_samples(circuit, input_state, nb_samples, j)\n", " E = 0\n", " for sample in samples:\n", " b = np.array([sample[i] for i in range(len(sample))])\n", " b_prime = np.dot(H, b)\n", " E += samples[sample]/nb_samples*np.dot(b.conjugate(), b_prime)\n", " \n", " return E.real\n", "\n", " # init_parameters = [2*(math.pi)*random.random() for _ in parameters_circuit] # initialize with random initial parameters\n", " init_parameters = [math.pi for _ in parameters_circuit] # initialize with a good guess\n", " best_parameters = minimize(loss, init_parameters, method='Powell', bounds=[(0,2*math.pi) for _ in init_parameters]).x \n", " set_parameters_circuit(parameters_circuit, best_parameters)\n", " samples = compute_samples(circuit, input_state, nb_samples, j)\n", " return samples\n", "\n", "def shortest_path(H, nb_samples):\n", " \"\"\"run the universal circuit and optimize the parameters\"\"\"\n", " nb_modes = len(H)\n", " circuit = GenericInterferometer(\n", " nb_modes,\n", " lambda i: BS(theta=pcvl.P(f\"theta{i}\"),phi_tr=pcvl.P(f\"phi_tr{i}\")),\n", " phase_shifter_fun_gen=lambda i: PS(phi=pcvl.P(f\"phi{i}\")))\n", " js = [0,1]\n", " ns = [nb_modes, nb_modes-1]\n", " configuration_samples = []\n", " for j in js:\n", " for n in ns:\n", " current_sample = test_configuration(circuit, nb_modes, j, n, H, nb_samples)\n", " print(f\"Configuration for (j,n)=({j},{n})\")\n", " print(current_sample)\n", " configuration_samples.append(([j,n],current_sample))\n", " \n", " return configuration_samples" ] }, { "attachments": {}, "cell_type": "markdown", "id": "75ab2dba-e47b-4d53-a92f-c13c077aeb69", "metadata": {}, "source": [ "Run the algorithm on 10000 samples (this may take a few minutes):" ] }, { "cell_type": "code", "execution_count": 15, "id": "7f440a3e", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Configuration for (j,n)=(0,6)\n", "{\n", " |1,1,0,0,0,0>: 1347\n", " |1,0,1,0,0,0>: 1186\n", " |0,0,0,0,0,0>: 33\n", " |0,1,0,1,0,0>: 898\n", " |0,1,0,0,1,0>: 1048\n", " |0,1,0,0,0,1>: 992\n", " |0,0,1,1,0,0>: 888\n", " |0,0,1,0,1,0>: 908\n", " |0,0,1,0,0,1>: 888\n", " |0,0,0,0,1,1>: 15\n", " |1,0,0,0,0,1>: 24\n", " |1,1,1,0,0,1>: 3\n", " |1,1,0,1,1,0>: 540\n", " |1,1,0,1,0,1>: 8\n", " |1,1,0,0,1,1>: 4\n", " |1,0,0,1,0,0>: 3\n", " |1,0,1,1,1,0>: 530\n", " |1,0,1,1,0,1>: 10\n", " |1,0,0,1,1,1>: 9\n", " |0,1,1,1,1,0>: 3\n", " |0,1,1,0,0,0>: 5\n", " |0,1,0,1,1,1>: 305\n", " |0,0,0,1,1,0>: 5\n", " |0,0,1,1,1,1>: 312\n", " |1,0,0,0,1,0>: 7\n", " |1,1,1,0,1,0>: 1\n", " |1,1,1,1,1,1>: 1\n", " |1,0,1,0,1,1>: 3\n", " |0,1,1,1,0,1>: 3\n", " |0,1,1,0,1,1>: 6\n", " |0,0,0,1,0,1>: 15\n", "}\n", "Configuration for (j,n)=(0,5)\n", "{\n", " |1,1,1,0,0,0>: 2\n", " |1,1,0,1,0,0>: 3\n", " |1,1,0,0,1,0>: 2\n", " |1,0,0,0,0,0>: 39\n", " |1,0,0,1,1,0>: 17\n", " |0,0,0,1,0,0>: 16\n", " |0,0,1,0,0,0>: 23\n", " |0,0,1,1,1,0>: 4\n", " |0,0,0,0,1,0>: 18\n", " |0,0,0,0,0,1>: 4\n", " |0,1,1,0,0,1>: 32\n", " |0,1,0,0,0,0>: 6554\n", " |0,1,0,1,1,0>: 3227\n", " |0,1,0,1,0,1>: 11\n", " |0,1,0,0,1,1>: 27\n", " |0,1,1,1,0,0>: 5\n", " |0,1,1,0,1,0>: 10\n", " |0,1,1,1,1,1>: 6\n", "}\n", "Configuration for (j,n)=(1,6)\n", "{\n", " |0,0,0,1,1,0>: 6\n", " |0,1,0,1,0,0>: 1\n", " |1,0,0,1,0,0>: 9979\n", " |1,1,0,0,0,0>: 2\n", " |1,1,1,1,0,0>: 1\n", " |1,0,1,0,0,0>: 3\n", " |1,0,0,0,0,1>: 7\n", " |1,0,1,1,0,1>: 1\n", "}\n", "Configuration for (j,n)=(1,5)\n", "{\n", " |1,1,0,0,1,0>: 2\n", " |0,1,0,0,1,1>: 3\n", " |0,1,0,1,1,0>: 9950\n", " |0,1,0,0,0,0>: 26\n", " |0,1,0,1,0,1>: 3\n", " |0,1,1,1,0,0>: 1\n", " |1,0,0,1,1,0>: 15\n", "}\n" ] } ], "source": [ "nb_samples = 10000\n", "data = shortest_path(H1, nb_samples)" ] }, { "cell_type": "markdown", "id": "308e2867", "metadata": {}, "source": [ "Plot the results for better visualization:" ] }, { "cell_type": "code", "execution_count": 16, "id": "f29aa6de", "metadata": {}, "outputs": [], "source": [ "def plot_samples(samples):\n", " \n", " values = list(samples.values())\n", "\n", " keys = samples.keys()\n", "\n", " key_list = []\n", " for x in keys:\n", " s = \"\".join(str(c) for c in list(x))\n", " key_list.append(s)\n", "\n", " y_pos = np.arange(len(key_list))\n", " barlist = plt.bar(y_pos, values, align='center', alpha=0.8)\n", " index = key_list.index('100100')\n", " barlist[index].set_color('m')\n", " plt.yscale('log')\n", " plt.xticks(y_pos, key_list)\n", " plt.xticks(rotation = 70)\n", " plt.ylabel('Sample values')\n", " plt.title('Sampling the Fock States')\n", "\n", " plt.show()" ] }, { "cell_type": "code", "execution_count": 17, "id": "f13d9a95", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjoAAAHcCAYAAADfvQDCAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABWh0lEQVR4nO3deVRU5eMG8GdAFgXEBUVRBHEHFRQR0RQ0lHApt7RSQzL0m7gkaWl9kyzLykwtxyUX0DTXzD0tMcMFv7KES64oEqKAiLIqyMz7+6PD/ERwGWbwMneezzmc49w7zjwvd8DHe997r0IIIUBEREQkQyZSByAiIiKqKiw6REREJFssOkRERCRbLDpEREQkWyw6REREJFssOkRERCRbLDpEREQkWyw6REREJFssOkRERCRbLDpERkyhUOCTTz7RPI6MjIRCocC1a9cky/Soa9euQaFQ4JtvvpE6it4pFApMmjRJ6hhEssaiQ6SjM2fOYPjw4XBycoKlpSWaNGmCvn374vvvv5c6mkHZt29fmdIlBYVCUeFXo0aNJM1VkeLiYixevBidOnVC7dq1UadOHbi5uWH8+PG4cOGC5nnHjx/HJ598grt371b6vZYuXYrIyEjdQxNJoIbUAYgM2fHjx9G7d280a9YMISEhaNSoEVJTU3HixAksXrwYkydPljqiVsaMGYPXXnsNFhYWz/299+3bB6VSKXnZ6du3L958880yy2rWrClRmscbNmwYfv31V7z++usICQnBgwcPcOHCBezZswfdu3dH27ZtAfz7GZ0zZw7Gjh2LOnXqVOq9li5dCjs7O4wdO1Z/AyB6Tlh0iHTw+eefw9bWFrGxseX+EcnMzJQmlA5MTU1hamoqdQxJtW7dGqNHj5Y6xhPFxsZiz549+Pzzz/Hhhx+WWbdkyRKd9t4QyQ0PXRHp4MqVK3Bzc6vwf8oNGzYs8zgiIgJ9+vRBw4YNYWFhAVdXVyxbtqzc33N2dsbAgQNx+PBhdOnSBTVr1kSHDh1w+PBhAMD27dvRoUMHWFpawtPTE3/99VeZvz927FhYW1vj6tWrCAgIgJWVFRwcHPDpp59CCPHE8VQ0R6c0z9GjR9G1a1dYWlrCxcUF69atK/f3T58+DV9fX9SsWRNNmzbF3LlzERER8dR5P2PHjoVSqQRQ9vDRo3744Qe0aNECFhYW8PLyQmxsbLnnXLhwAcOHD0e9evVgaWmJLl26YNeuXU8ctzYyMzMxbtw42Nvbw9LSEu7u7li7dm2556nVaixevFizrRo0aICXXnoJcXFxT3z9uXPnwsTE5ImHPq9cuQIA6NGjR7l1pqamqF+/PgDgk08+wYwZMwAAzZs313xfS7fFs3wmnZ2d8ffff+PPP//U/H0/Pz/N+rt37+Ldd9+Fo6MjLCws0LJlS3z11VdQq9VlXmfTpk3w9PSEjY0NateujQ4dOmDx4sVP/F4Q6QP36BDpwMnJCTExMTh79izat2//xOcuW7YMbm5uePnll1GjRg3s3r0bEydOhFqtRmhoaJnnJiUl4Y033sCECRMwevRofPPNNxg0aBCWL1+ODz/8EBMnTgQAzJs3DyNGjMDFixdhYvL//29RqVR46aWX0K1bN3z99dfYv38/wsPDUVJSgk8//VTrcSYlJWH48OEYN24cgoKCsGbNGowdOxaenp5wc3MDAKSlpaF3795QKBSYNWsWrKyssGrVqmc6DDZhwgTcuHEDv//+O3788ccKn/PTTz8hLy8PEyZMgEKhwNdff42hQ4fi6tWrMDMzAwD8/fff6NGjB5o0aYKZM2fCysoKW7ZsweDBg/Hzzz9jyJAhT81y//59ZGVllVlmY2MDCwsL3Lt3D35+fkhKSsKkSZPQvHlzbN26FWPHjsXdu3cxdepUzd8ZN24cIiMjERgYiLfffhslJSU4cuQITpw4gS5dulT43v/973/xxRdfYMWKFQgJCXlsRicnJwDAhg0b0KNHD9SoUfGv8qFDh+LSpUvYuHEjFi5cCDs7OwBAgwYNADzbZ3LRokWYPHkyrK2t8dFHHwEA7O3tAQCFhYXw9fVFWloaJkyYgGbNmuH48eOYNWsWbt68iUWLFgEAfv/9d7z++ut48cUX8dVXXwEAzp8/j2PHjpX5nhFVCUFElfbbb78JU1NTYWpqKnx8fMT7778vDhw4IIqLi8s9t7CwsNyygIAA4eLiUmaZk5OTACCOHz+uWXbgwAEBQNSsWVOkpKRolq9YsUIAEH/88YdmWVBQkAAgJk+erFmmVqvFgAEDhLm5ubh165ZmOQARHh6ueRwRESEAiOTk5HJ5oqOjNcsyMzOFhYWFeO+99zTLJk+eLBQKhfjrr780y27fvi3q1atX7jUrEhoaKir6lZScnCwAiPr164vs7GzN8p07dwoAYvfu3ZplL774oujQoYO4f/9+mbF3795dtGrV6onvL8S/34+KviIiIoQQQixatEgAEOvXr9f8neLiYuHj4yOsra1Fbm6uEEKIQ4cOCQBiypQp5d5DrVaXeb/Q0FAhhBDvvfeeMDExEZGRkU/NqVarha+vrwAg7O3txeuvvy6USmWZz0ap+fPnP/b7/6yfSTc3N+Hr61vuuZ999pmwsrISly5dKrN85syZwtTUVPzzzz9CCCGmTp0qateuLUpKSp46NiJ946ErIh307dsXMTExePnll3Hq1Cl8/fXXCAgIQJMmTcodLnl4QmtOTg6ysrLg6+uLq1evIicnp8xzXV1d4ePjo3ns7e0NAOjTpw+aNWtWbvnVq1fLZXv4tOXS05iLi4tx8OBBrcfp6uqKnj17ah43aNAAbdq0KfO++/fvh4+PDzw8PDTL6tWrh1GjRmn9fhUZOXIk6tatq3lcmqc0Q3Z2Ng4dOoQRI0YgLy8PWVlZyMrKwu3btxEQEIDLly8jLS3tqe/zyiuv4Pfffy/zFRAQAODfCdONGjXC66+/rnm+mZkZpkyZgvz8fPz5558AgJ9//hkKhQLh4eHlXv/RQ3JCCEyaNAmLFy/G+vXrERQU9NSMCoUCBw4cwNy5c1G3bl1s3LgRoaGhcHJywsiRI595jo42n8mKbN26FT179kTdunU13++srCz4+/tDpVIhOjoaAFCnTh0UFBTg999/f6ZcRPrEQ1dEOvLy8sL27dtRXFyMU6dO4ZdffsHChQsxfPhwJCYmwtXVFQBw7NgxhIeHIyYmBoWFhWVeIycnB7a2tprHD5cZAJp1jo6OFS6/c+dOmeUmJiZwcXEps6x169YAUKlr5DyaBwDq1q1b5n1TUlLKlLNSLVu21Pr9niVDaekpzZCUlAQhBD7++GN8/PHHFb5GZmYmmjRp8sT3adq0Kfz9/Stcl5KSglatWpU5TAgA7dq106wH/p1D4+DggHr16j1lVMC6deuQn5+PZcuWlSlQT2NhYYGPPvoIH330EW7evIk///wTixcvxpYtW2BmZob169c/9TW0+UxW5PLlyzh9+rTmUNijSifkT5w4EVu2bEFgYCCaNGmCfv36YcSIEXjppZeecbRElceiQ6Qn5ubm8PLygpeXF1q3bo3g4GBs3boV4eHhuHLlCl588UW0bdsW3377LRwdHWFubo59+/Zh4cKF5SZuPu7Mp8ctF0+ZZKwrqd5Xmwyl38Pp06dr9sA8Sl+lS5969OiBxMRELFmyBCNGjHimcvSoxo0b47XXXsOwYcPg5uaGLVu2IDIy8rFzdwBo/ZmsiFqtRt++ffH+++9XuL60XDds2BCJiYk4cOAAfv31V/z666+IiIjAm2++WeFEbiJ9YtEhqgKlk01v3rwJANi9ezeKioqwa9euMnsm/vjjjyp5f7VajatXr2r+oQGAS5cuAfj3LJqq4OTkhKSkpHLLK1pWkYrOstJG6R4sMzOzx+6R0ZWTkxNOnz4NtVpdZq9O6QX6SicJt2jRAgcOHEB2dvZTi0vLli3x9ddfw8/PDy+99BKioqJgY2NTqXxmZmbo2LEjLl++jKysLDRq1Oix31dtPpOPe40WLVogPz//mb7f5ubmGDRoEAYNGgS1Wo2JEydixYoV+Pjjj6tlASX54BwdIh388ccfFe7V2LdvHwCgTZs2AP5/b8TDz83JyUFERESVZVuyZInmz0IILFmyBGZmZnjxxRer5P0CAgIQExODxMREzbLs7Gxs2LDhmf6+lZUVAFT6GjANGzaEn58fVqxYoSmYD7t161alXvdh/fv3R3p6OjZv3qxZVlJSgu+//x7W1tbw9fUF8O/F/IQQmDNnTrnXqOjz0rFjR+zbtw/nz5/HoEGDcO/evSfmuHz5Mv75559yy+/evYuYmBjUrVtXczjpcd9XbT6TVlZWFW6XESNGICYmBgcOHKgwS0lJCQDg9u3bZdaZmJigY8eOAICioqLHDZNIL7hHh0gHkydPRmFhIYYMGYK2bduiuLgYx48fx+bNm+Hs7Izg4GAAQL9+/TT/o50wYQLy8/OxcuVKNGzYsMJ/lHVlaWmJ/fv3IygoCN7e3vj111+xd+9efPjhh4+dT6Gr999/H+vXr0ffvn0xefJkzenlzZo1Q3Z29lP32Hh6egIApkyZgoCAAJiamuK1117TKoNSqcQLL7yADh06ICQkBC4uLsjIyEBMTAyuX7+OU6dOVXp8ADB+/HisWLECY8eORXx8PJydnbFt2zYcO3YMixYt0uyJ6d27N8aMGYPvvvsOly9fxksvvQS1Wo0jR46gd+/eFd7fqlu3bti5cyf69++P4cOHY8eOHZrT5h916tQpvPHGGwgMDETPnj1Rr149pKWlYe3atbhx4wYWLVqkKTKl39ePPvoIr732GszMzDBo0CCtPpOenp5YtmwZ5s6di5YtW6Jhw4bo06cPZsyYgV27dmHgwIGayw0UFBTgzJkz2LZtG65duwY7Ozu8/fbbyM7ORp8+fdC0aVOkpKTg+++/h4eHh2Z+E1GVkehsLyJZ+PXXX8Vbb70l2rZtK6ytrYW5ublo2bKlmDx5ssjIyCjz3F27domOHTsKS0tL4ezsLL766iuxZs2aCk/nHjBgQLn3wkOnIpcqPfV6/vz5mmVBQUHCyspKXLlyRfTr10/UqlVL2Nvbi/DwcKFSqcq95rOcXl5RHl9f33KnHP/111+iZ8+ewsLCQjRt2lTMmzdPfPfddwKASE9Pf9y3UQghRElJiZg8ebJo0KCBUCgUmlPNKxrj4/ILIcSVK1fEm2++KRo1aiTMzMxEkyZNxMCBA8W2bdue+P6lr/fo9/hRGRkZIjg4WNjZ2Qlzc3PRoUMHzennj45n/vz5om3btsLc3Fw0aNBABAYGivj4+Ce+386dO0WNGjXEyJEjy22vhzN8+eWXwtfXVzRu3FjUqFFD1K1bV/Tp06fCcX722WeiSZMmwsTEpMz2fdbPZHp6uhgwYICwsbERAMps97y8PDFr1izRsmVLYW5uLuzs7ET37t3FN998o7nMwrZt20S/fv1Ew4YNhbm5uWjWrJmYMGGCuHnz5hO/10T6oBDiOc4mJKIqN3bsWGzbtg35+flSRwEAvPvuu1ixYgXy8/ON/vYSRPT8cY4OEenNo3NLbt++jR9//BEvvPACSw4RSYJzdIhIb3x8fODn54d27dohIyMDq1evRm5u7mOva0NEVNVYdIhIb/r3749t27bhhx9+gEKhQOfOnbF69Wr06tVL6mhEZKQ4R4eIiIhki3N0iIiISLZYdIiIiEi2jH6Ojlqtxo0bN2BjY6PzJeiJiIjo+RBCIC8vDw4ODuVutPswoy86N27cKHdHaCIiIjIMqampaNq06WPXG33RKb1ke2pqKmrXri1xGiIiInoWubm5cHR0fOpNcI2+6JQerqpduzaLDhERkYF52rQTTkYmIiIi2WLRISIiItli0SEiIiLZYtEhIiIi2WLRISIiItli0SEiIiLZYtEhIiIi2WLRISIiItmSTdEpLCyEk5MTpk+fLnUUIiIiqiZkU3Q+//xzdOvWTeoYREREVI3IouhcvnwZFy5cQGBgoNRRiIiIqBqRvOhER0dj0KBBcHBwgEKhwI4dO8o9R6lUwtnZGZaWlvD29sbJkyfLrJ8+fTrmzZv3nBITERGRoZC86BQUFMDd3R1KpbLC9Zs3b0ZYWBjCw8ORkJAAd3d3BAQEIDMzEwCwc+dOtG7dGq1bt36esYmIiMgAKIQQQuoQpRQKBX755RcMHjxYs8zb2xteXl5YsmQJAECtVsPR0RGTJ0/GzJkzMWvWLKxfvx6mpqbIz8/HgwcP8N5772H27NkVvkdRURGKioo0j0tv856Tk8O7l9MTFaUXoeRuidQxtFajTg1YNLKQOgYRkV7l5ubC1tb2qf9+13iOmbRWXFyM+Ph4zJo1S7PMxMQE/v7+iImJAQDMmzdPc9gqMjISZ8+efWzJKX3+nDlzqjY4yU5RehFOv3QaJXcMsOjUrYGO+zuy7BCRUarWRScrKwsqlQr29vZlltvb2+PChQuVes1Zs2YhLCxM87h0jw7Rk5TcLUHJnRIozBUwsZD8iO8zUxepUXKnBCV3S1h0iMgoVeuio62xY8c+9TkWFhawsOAvfKocEwsTmFgaTtEBAFWxSuoIRESSqda/se3s7GBqaoqMjIwyyzMyMtCoUSOJUhEREZGhqNZFx9zcHJ6enoiKitIsU6vViIqKgo+Pj06vrVQq4erqCi8vL11jEhERUTUl+aGr/Px8JCUlaR4nJycjMTER9erVQ7NmzRAWFoagoCB06dIFXbt2xaJFi1BQUIDg4GCd3jc0NBShoaGaWdtEREQkP5IXnbi4OPTu3VvzuHSicFBQECIjIzFy5EjcunULs2fPRnp6Ojw8PLB///5yE5SJiIiIHiV50fHz88PTLuUzadIkTJo06TklIiIiIrmo1nN0qhLn6BAREcmf0Rad0NBQnDt3DrGxsVJHISIioipitEWHiIiI5I9Fh4iIiGSLRYeIiIhky2iLDicjExERyZ/RFh1ORiYiIpI/oy06REREJH8sOkRERCRbLDpEREQkW0ZbdDgZmYiISP6MtuhwMjIREZH8GW3RISIiIvlj0SEiIiLZYtEhIiIi2WLRISIiItli0SEiIiLZMtqiw9PLiYiI5M9oiw5PLyciIpI/oy06REREJH8sOkRERCRbLDpEREQkWyw6REREJFssOkRERCRbLDpEREQkW0ZbdHgdHSIiIvkz2qLD6+gQERHJn9EWHSIiIpI/Fh0iIiKSLRYdIiIiki0WHSIiIpItFh0iIiKSLRYdIiIiki0WHSIiIpItFh0iIiKSLRYdIiIiki2jLTq8BQQREZH8GW3R4S0giIiI5M9oiw4RERHJH4sOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJltEWHaVSCVdXV3h5eUkdhYiIiKqI0Rad0NBQnDt3DrGxsVJHISIioipitEWHiIiI5I9Fh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGTL4IvO3bt30aVLF3h4eKB9+/ZYuXKl1JGIiIiomqghdQBd2djYIDo6GrVq1UJBQQHat2+PoUOHon79+lJHIyIiIokZ/B4dU1NT1KpVCwBQVFQEIQSEEBKnIiIioupA8qITHR2NQYMGwcHBAQqFAjt27Cj3HKVSCWdnZ1haWsLb2xsnT54ss/7u3btwd3dH06ZNMWPGDNjZ2T2n9ERERFSdSV50CgoK4O7uDqVSWeH6zZs3IywsDOHh4UhISIC7uzsCAgKQmZmpeU6dOnVw6tQpJCcn46effkJGRsbzik9ERETVmORFJzAwEHPnzsWQIUMqXP/tt98iJCQEwcHBcHV1xfLly1GrVi2sWbOm3HPt7e3h7u6OI0eOPPb9ioqKkJubW+aLiIiI5EnyovMkxcXFiI+Ph7+/v2aZiYkJ/P39ERMTAwDIyMhAXl4eACAnJwfR0dFo06bNY19z3rx5sLW11Xw5OjpW7SCIiIhIMtW66GRlZUGlUsHe3r7Mcnt7e6SnpwMAUlJS0LNnT7i7u6Nnz56YPHkyOnTo8NjXnDVrFnJycjRfqampVToGIiIiko7Bn17etWtXJCYmPvPzLSwsYGFhUXWBiIiIqNqo1nt07OzsYGpqWm5ycUZGBho1aiRRKiIiIjIU1bromJubw9PTE1FRUZplarUaUVFR8PHx0em1lUolXF1d4eXlpWtMIiIiqqYkP3SVn5+PpKQkzePk5GQkJiaiXr16aNasGcLCwhAUFIQuXbqga9euWLRoEQoKChAcHKzT+4aGhiI0NBS5ubmwtbXVdRhERERUDUledOLi4tC7d2/N47CwMABAUFAQIiMjMXLkSNy6dQuzZ89Geno6PDw8sH///nITlImIiIgeJXnR8fPze+otGyZNmoRJkyY9p0REREQkF9V6jk5V4hwdIiIi+TPaohMaGopz584hNjZW6ihERERURYy26BAREZH8segQERGRbLHoEBERkWwZbdHhZGQiIiL5M9qiw8nIRERE8me0RYeIiIjkj0WHiIiIZItFh4iIiGTLaIsOJyMTERHJn9EWHU5GJiIikj+jLTpEREQkfyw6REREJFssOkRERCRbLDpEREQkWyw6REREJFtGW3R4ejkREZH8GW3R4enlRERE8me0RYeIiIjkj0WHiIiIZItFh4iIiGSLRYeIiIhki0WHiIiIZItFh4iIiGTLaIsOr6NDREQkf0ZbdHgdHSIiIvnTuuisXbsWe/fu1Tx+//33UadOHXTv3h0pKSl6DUdERESkC62LzhdffIGaNWsCAGJiYqBUKvH111/Dzs4O06ZN03tAIiIiosqqoe1fSE1NRcuWLQEAO3bswLBhwzB+/Hj06NEDfn5++s5HREREVGla79GxtrbG7du3AQC//fYb+vbtCwCwtLTEvXv39JuOiIiISAda79Hp27cv3n77bXTq1AmXLl1C//79AQB///03nJ2d9Z2PiIiIqNK03qOjVCrh4+ODW7du4eeff0b9+vUBAPHx8Xj99df1HpCIiIiosrTeo1OnTh0sWbKk3PI5c+boJRARERGRvlTqOjpHjhzB6NGj0b17d6SlpQEAfvzxRxw9elSv4YiIiIh0oXXR+fnnnxEQEICaNWsiISEBRUVFAICcnBx88cUXeg9IREREVFlaF525c+di+fLlWLlyJczMzDTLe/TogYSEBL2Gq0q8BQQREZH8aV10Ll68iF69epVbbmtri7t37+oj03PBW0AQERHJn9ZFp1GjRkhKSiq3/OjRo3BxcdFLKCIiIiJ90LrohISEYOrUqfjf//4HhUKBGzduYMOGDZg+fTreeeedqshIREREVClan14+c+ZMqNVqvPjiiygsLESvXr1gYWGB6dOnY/LkyVWRkYiIiKhStC46CoUCH330EWbMmIGkpCTk5+fD1dUV1tbWVZGPiIiIqNK0LjqlzM3N4erqqs8sRERERHqlddHp3bs3FArFY9cfOnRIp0BERERE+qJ10fHw8Cjz+MGDB0hMTMTZs2cRFBSkr1xEREREOtO66CxcuLDC5Z988gny8/N1DkRERESkL5W611VFRo8ejTVr1ujr5YiIiIh0preiExMTA0tLS329HBEREZHOtD50NXTo0DKPhRC4efMm4uLi8PHHH+stGBEREZGutC46tra2ZR6bmJigTZs2+PTTT9GvXz+9BSMiIiLSldZFJyIioipyEBEREemd3uboGBqlUglXV1d4eXlJHYWIiIiqyDPt0albt+4TLxL4sOzsbJ0CPS+hoaEIDQ1Fbm5uucNxREREJA/PVHQWLVpUxTGIiIiI9O+Zig6veExERESGqNI39QSA+/fvo7i4uMyy2rVr6xSIiIiISF+0noxcUFCASZMmoWHDhrCyskLdunXLfBERERFVF1oXnffffx+HDh3CsmXLYGFhgVWrVmHOnDlwcHDAunXrqiIjERERUaVofehq9+7dWLduHfz8/BAcHIyePXuiZcuWcHJywoYNGzBq1KiqyElERESkNa336GRnZ8PFxQXAv/NxSk8nf+GFFxAdHa3fdEREREQ60LrouLi4IDk5GQDQtm1bbNmyBcC/e3rq1Kmj13BEREREutC66AQHB+PUqVMAgJkzZ0KpVMLS0hLTpk3DjBkz9B6QiIiIqLK0nqMzbdo0zZ/9/f1x4cIFxMfHo2XLlujYsaNewxERERHpQuuik5qaCkdHR81jJycnODk56TUUERERkT5ofejK2dkZvr6+WLlyJe7cuVMVmYiIiIj0QuuiExcXh65du+LTTz9F48aNMXjwYGzbtg1FRUVVkY+IiIio0rQuOp06dcL8+fPxzz//4Ndff0WDBg0wfvx42Nvb46233qqKjERERESVonXRKaVQKNC7d2+sXLkSBw8eRPPmzbF27Vp9ZiMiIiLSSaWLzvXr1/H111/Dw8MDXbt2hbW1NZRKpT6zEREREelE67OuVqxYgZ9++gnHjh1D27ZtMWrUKOzcuZNnXhEREVG1o/Uenblz58Lb2xvx8fE4e/YsZs2aJWnJSU1NhZ+fH1xdXdGxY0ds3bpVsixERERUvWi9R+eff/6BQqGoiiyVUqNGDSxatAgeHh5IT0+Hp6cn+vfvDysrK6mjERERkcS0LjrVqeQAQOPGjdG4cWMAQKNGjWBnZ4fs7GwWHSIiIqr8ZGR9iY6OxqBBg+Dg4ACFQoEdO3aUe45SqYSzszMsLS3h7e2NkydPVvha8fHxUKlUZa7cTERERMZL8qJTUFAAd3f3x56xtXnzZoSFhSE8PBwJCQlwd3dHQEAAMjMzyzwvOzsbb775Jn744YfnEZuIiIgMgNaHrvQtMDAQgYGBj13/7bffIiQkBMHBwQCA5cuXY+/evVizZg1mzpwJACgqKsLgwYMxc+ZMdO/e/YnvV1RUVOYqzrm5uXoYBREREVVHldqjU1JSgoMHD2LFihXIy8sDANy4cQP5+fl6DVdcXIz4+Hj4+/trlpmYmMDf3x8xMTEAACEExo4diz59+mDMmDFPfc158+bB1tZW88XDXERERPKlddFJSUlBhw4d8MorryA0NBS3bt0CAHz11VeYPn26XsNlZWVBpVLB3t6+zHJ7e3ukp6cDAI4dO4bNmzdjx44d8PDwgIeHB86cOfPY15w1axZycnI0X6mpqXrNTERERNWH1oeupk6dii5duuDUqVOoX7++ZvmQIUMQEhKi13DP4oUXXoBarX7m51tYWMDCwqIKExEREVF1oXXROXLkCI4fPw5zc/Myy52dnZGWlqa3YABgZ2cHU1NTZGRklFmekZGBRo0a6fW9iIiISH60PnSlVquhUqnKLb9+/TpsbGz0EqqUubk5PD09ERUVVeb9o6Ki4OPjo9NrK5VKuLq6wsvLS9eYREREVE1pXXT69euHRYsWaR4rFArk5+cjPDwc/fv31zpAfn4+EhMTkZiYCABITk5GYmIi/vnnHwBAWFgYVq5cibVr1+L8+fN45513UFBQoDkLq7JCQ0Nx7tw5xMbG6vQ6REREVH1pfehqwYIFCAgIgKurK+7fv4833ngDly9fhp2dHTZu3Kh1gLi4OPTu3VvzOCwsDAAQFBSEyMhIjBw5Erdu3cLs2bORnp4ODw8P7N+/v9wEZSIiIqJHKYQQQtu/VFJSgk2bNuH06dPIz89H586dMWrUKNSsWbMqMlap3Nxc2NraIicnB7Vr15Y6DlVTBRcKcDrgNExtTGFiKfl1Np+Z+r4aqjwVOh7oCKu2vC0KEcnHs/77XakLBtaoUQOjR4+udLjqQKlUQqlUVjjfiIiIiOThmYrOrl27nvkFX3755UqHeZ5CQ0MRGhqqaYREREQkP89UdAYPHvxML6ZQKLiHhIiIiKqNZyo62lyQj4iIiKi6MJxZlURERERaqlTRiYqKwsCBA9GiRQu0aNECAwcOxMGDB/WdrUrxgoFERETyp3XRWbp0KV566SXY2Nhg6tSpmDp1KmrXro3+/ftDqVRWRcYqwQsGEhERyZ/Wp5d/8cUXWLhwISZNmqRZNmXKFPTo0QNffPEFQkND9RqQiIiIqLK03qNz9+5dvPTSS+WW9+vXDzk5OXoJRURERKQPWhedl19+Gb/88ku55Tt37sTAgQP1EoqIiIhIH7Q+dOXq6orPP/8chw8f1txB/MSJEzh27Bjee+89fPfdd5rnTpkyRX9J9YxXRiYiIpI/re911bx582d7YYUCV69erVSo54n3uqJnwXtdERFVL1V2r6vk5GSdghERERE9L4bzX1MiIiIiLWm9R0cIgW3btuGPP/5AZmZmudtDbN++XW/hiIiIiHShddF59913sWLFCvTu3Rv29vZQKBRVkYuIiIhIZ1oXnR9//BHbt29H//79qyIPERERkd5oPUfH1tYWLi4uVZHlueK9roiIiORP66LzySefYM6cObh3715V5HlueK8rIiIi+dP60NWIESOwceNGNGzYEM7OzjAzMyuzPiEhQW/hiIiIiHShddEJCgpCfHw8Ro8ezcnIREREVK1pXXT27t2LAwcO4IUXXqiKPERERER6o/UcHUdHR94qgYiIiAyC1kVnwYIFeP/993Ht2rUqiENERESkP1ofuho9ejQKCwvRokUL1KpVq9xk5OzsbL2FIyIiItKF1kVn0aJFVRDj+VMqlVAqlVCpVFJHISIioiqiEEIIqUNI6Vlv807GreBCAU4HnIapjSlMLA3nXrjq+2qo8lToeKAjrNpaSR2HiEhvnvXfb6336Dzs/v37KC4uLrOMZYGIiIiqC63/a1pQUIBJkyahYcOGsLKyQt26dct8EREREVUXWhed999/H4cOHcKyZctgYWGBVatWYc6cOXBwcMC6deuqIiMRERFRpWh96Gr37t1Yt24d/Pz8EBwcjJ49e6Jly5ZwcnLChg0bMGrUqKrISURERKQ1rffoZGdna+5eXrt2bc3p5C+88AKio6P1m46IiIhIB1oXHRcXFyQnJwMA2rZtiy1btgD4d09PnTp19BqOiIiISBdaF53g4GCcOnUKADBz5kwolUpYWlpi2rRpmDFjht4DEhEREVWW1nN0pk2bpvmzv78/zp8/j4SEBLRs2RIdO3bUazgiIiIiXeh0HR0AcHZ2hrOzsx6iEBEREenXMx+6iomJwZ49e8osW7duHZo3b46GDRti/PjxKCoq0nvAqqJUKuHq6govLy+poxAREVEVeeai8+mnn+Lvv//WPD5z5gzGjRsHf39/zJw5E7t378a8efOqJGRVCA0Nxblz5xAbGyt1FCIiIqoiz1x0EhMT8eKLL2oeb9q0Cd7e3li5ciXCwsLw3Xffac7AIiIiIqoOnrno3LlzB/b29prHf/75JwIDAzWPvby8kJqaqt90RERERDp45qJjb2+vuX5OcXExEhIS0K1bN836vLw8mJmZ6T8hERERUSU9c9Hp378/Zs6ciSNHjmDWrFmoVasWevbsqVl/+vRptGjRokpCEhEREVXGM59e/tlnn2Ho0KHw9fWFtbU11q5dC3Nzc836NWvWoF+/flUSkoiIiKgynrno2NnZITo6Gjk5ObC2toapqWmZ9Vu3boW1tbXeAxIRERFVltYXDLS1ta1web169XQOQ0RERKRPWt/rioiIiMhQsOgQERGRbLHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWwZbdFRKpVwdXWFl5eX1FGIiIioihht0QkNDcW5c+cQGxsrdRQiIiKqIkZbdIiIiEj+WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItmRRdIYMGYK6deti+PDhUkchIiKiakQWRWfq1KlYt26d1DGIiIiompFF0fHz84ONjY3UMYiIiKiakbzoREdHY9CgQXBwcIBCocCOHTvKPUepVMLZ2RmWlpbw9vbGyZMnn39QIiIiMjiSF52CggK4u7tDqVRWuH7z5s0ICwtDeHg4EhIS4O7ujoCAAGRmZj7npERERGRoakgdIDAwEIGBgY9d/+233yIkJATBwcEAgOXLl2Pv3r1Ys2YNZs6cqfX7FRUVoaioSPM4NzdX+9BERERkECTfo/MkxcXFiI+Ph7+/v2aZiYkJ/P39ERMTU6nXnDdvHmxtbTVfjo6O+opLRERE1Uy1LjpZWVlQqVSwt7cvs9ze3h7p6emax/7+/nj11Vexb98+NG3a9IklaNasWcjJydF8paamVll+IiIikpbkh6704eDBg8/8XAsLC1hYWFRhGiIiIqouqvUeHTs7O5iamiIjI6PM8oyMDDRq1EiiVERERGQoqnXRMTc3h6enJ6KiojTL1Go1oqKi4OPjo9NrK5VKuLq6wsvLS9eYREREVE1JfugqPz8fSUlJmsfJyclITExEvXr10KxZM4SFhSEoKAhdunRB165dsWjRIhQUFGjOwqqs0NBQhIaGIjc3F7a2troOg4iIiKohyYtOXFwcevfurXkcFhYGAAgKCkJkZCRGjhyJW7duYfbs2UhPT4eHhwf2799fboIyERER0aMkLzp+fn4QQjzxOZMmTcKkSZOeUyIiIiKSi2o9R6cqcY4OERGR/Blt0QkNDcW5c+cQGxsrdRQiIiKqIkZbdIiIiEj+WHSIiIhItlh0iIiISLYkP+tKKkqlEkqlEiqVSuooRESkg0HfH5U6QqXsnvyC1BGMgtHu0eFkZCIiIvkz2qJDRERE8seiQ0RERLLFokNERESyxaJDREREssWzrqrwrCueCUBERCQto92jw7OuiIiI5M9oiw4RERHJH4sOERERyRaLDhEREckWiw4RERHJltEWHaVSCVdXV3h5eUkdhYiIiKqI0RYdnnVFREQkf0ZbdIiIiEj+WHSIiIhItlh0iIiISLZYdIiIiEi2WHSIiIhItlh0iIiISLaMtujwOjpERETyZ7RFh9fRISIikj+jLTpEREQkfyw6REREJFssOkRERCRbLDpEREQkWyw6REREJFssOkRERCRbLDpEREQkWyw6REREJFssOkRERCRbNaQOIBWlUgmlUgmVSiV1FCKiKjHo+6NSR6iU3ZNfkDoCyYjR7tHhLSCIiIjkz2iLDhEREckfiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOERERyRaLDhEREclWDakDSEWpVEKpVEKlUkkdhahaGfT9UakjVMruyS9IHYGIqiGj3aMTGhqKc+fOITY2VuooREREVEWMtugQERGR/LHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWyx6BAREZFssegQERGRbLHoEBERkWzJoujs2bMHbdq0QatWrbBq1Sqp4xAREVE1UUPqALoqKSlBWFgY/vjjD9ja2sLT0xNDhgxB/fr1pY5GREREEjP4PTonT56Em5sbmjRpAmtrawQGBuK3336TOhYRERFVA5IXnejoaAwaNAgODg5QKBTYsWNHuecolUo4OzvD0tIS3t7eOHnypGbdjRs30KRJE83jJk2aIC0t7XlEJyIiompO8qJTUFAAd3d3KJXKCtdv3rwZYWFhCA8PR0JCAtzd3REQEIDMzMznnJSIiIgMjeRFJzAwEHPnzsWQIUMqXP/tt98iJCQEwcHBcHV1xfLly1GrVi2sWbMGAODg4FBmD05aWhocHBwe+35FRUXIzc0t80VERETyVK0nIxcXFyM+Ph6zZs3SLDMxMYG/vz9iYmIAAF27dsXZs2eRlpYGW1tb/Prrr/j4448f+5rz5s3DnDlzqjy7MRn0/VGpI1TK7skvSB2BJGKon1mAn1tjZqifW6k/s5Lv0XmSrKwsqFQq2Nvbl1lub2+P9PR0AECNGjWwYMEC9O7dGx4eHnjvvfeeeMbVrFmzkJOTo/lKTU2t0jEQERGRdKr1Hp1n9fLLL+Pll19+pudaWFjAwsKiihMRERFRdVCt9+jY2dnB1NQUGRkZZZZnZGSgUaNGEqUiIiIiQ1Gti465uTk8PT0RFRWlWaZWqxEVFQUfHx+dXlupVMLV1RVeXl66xiQiIqJqSvJDV/n5+UhKStI8Tk5ORmJiIurVq4dmzZohLCwMQUFB6NKlC7p27YpFixahoKAAwcHBOr1vaGgoQkNDkZubC1tbW12HQURERNWQ5EUnLi4OvXv31jwOCwsDAAQFBSEyMhIjR47ErVu3MHv2bKSnp8PDwwP79+8vN0GZiIiI6FGSFx0/Pz8IIZ74nEmTJmHSpEnPKRERERHJRbWeo1OVOEeHiIhI/oy26ISGhuLcuXOIjY2VOgoRERFVEaMtOkRERCR/LDpEREQkWyw6REREJFtGW3Q4GZmIiEj+jLbocDIyERGR/Blt0SEiIiL5k/yCgVIrvVhhbm6u3l/7wb0Cvb/m86Dt98IYxlmQX4ACdQFMVCYwURnO/w/UKjXUajVy83OhylU9098xhu1pqGMEjGOc/B1UMWMZp7av+7SLDivE054hc9evX4ejo6PUMYiIiKgSUlNT0bRp08euN/qio1arcePGDdjY2EChUEgd55nk5ubC0dERqampqF27ttRxqgzHKS8cp3wYwxgBjrO6E0IgLy8PDg4OMDF5/J52oz90ZWJi8sQmWJ3Vrl3boD6UlcVxygvHKR/GMEaA46zObG1tn/ocw5lsQERERKQlFh0iIiKSLRYdA2RhYYHw8HBYWFhIHaVKcZzywnHKhzGMEeA45cLoJyMTERGRfHGPDhEREckWiw4RERHJFosOERERyRaLDhEREckWiw4RERHJFosOUTVmLCdFGsM4hRBQq9VSx3gujGF7kuHg6eUycOnSJaSmpiIlJQV2dnbo1asX6tSpI3UsvTOWcVZErVZDoVAYzP3YKkuu41SpVDA1NdU8lus4H2Us45QbIYSstpnR3+vK0C1cuBARERE4f/482rRpA0tLS6hUKvTs2RNjxoyBl5eXLD60xjJOlUqFqKgoZGdnIzMzE23atEHv3r1hbm4udTS9MpZxpqWl4eeff8alS5eQlpYGf39/jBw5EnZ2dlJH0ytj2Z7GovT3qBx+pwLco2PQ7t69C0dHR3z55Zd4++23cfnyZZw7dw7x8fGIj4+HSqXCF198AR8fH6mj6sRYxllQUICJEydi586dsLCwgIuLCwoLC1GzZk30798fo0ePhouLi8H/8jGWcebl5eGVV17BpUuX0KVLF5ibm+PYsWNIT09Hv3798OGHH6Jnz55Sx9SZMWxPtVqNwsJCWFtbSx2lShUXF+Po0aNo164dGjduXGadIW8/CDJYP/zwg/Dw8Ci3vLi4WBw7dkwEBgaKunXriitXrkiQTn+MZZxfffWVcHV1FSdOnBBCCBEXFyfWrVsnJk6cKHx8fMRrr70m7ty5I21IPTCWcX755ZeiU6dO4tatW0IIIe7cuSOSk5PFTz/9JAICAoS7u7uIioqSOKXujGF7Ll++XLi6uor58+eLixcvCpVKVe45+fn5IiYmRpSUlEiQUD+WLl0q6tatK958802xfPlycfLkyXLbLisrS6xevVoUFRVJE7ISWHQM2C+//CJat24tYmJiKlxfWFgovL29xQ8//PCck+mXsYzzhRdeEF9//XW55bm5uWLXrl2iadOmYsCAARX+kjUkxjLOV155RUyZMqXccpVKJZKTk8Urr7wi2rZtK7KysiRIpz/GsD09PT1Fq1athKOjozAzMxN9+vQRERERIi0tTfOcNWvWiBdffFHClLrr1auX6NWrl+jXr59o0qSJ8PDwEBMnThQ//fSTOHfunCgqKhKrV68WTk5OUkfVCs+6MmB9+vSBg4MDvv76a8TGxqKkpKTM+po1a8LMzAxpaWkSJdQPYxhnSUkJunTpgr179+LOnTtl1tnY2GDQoEGIjIzE9evXce7cOYlS6s5YxgkAAwcOxK5du5CSklJmuYmJCZydnbF48WJYWlrir7/+kiih7oxhe2ZmZsLExARz587FlStXsHfvXjRs2BBTpkxBu3bt8MYbb2Dv3r34/vvv0bp1a6njVtrt27cBAMHBwThw4ACOHTuGV199FSdPnsSMGTMwYcIEzJkzB5999hmGDBkicVotSd20SDdHjx4V7du3F1ZWVmLMmDFi9+7d4u+//xbx8fFi5cqVonbt2iIpKUnqmDqLjo4Wbm5uwsrKSowePVqW44yJiRFt2rQRM2fOrPAw3LVr14SVlZVISUmRIJ3+GMs4r1+/Lnr27Cm6d+8uIiIixOXLl0VxcbFm/ZUrV4SlpaW4evWqhCl1J/ftmZKSIubMmSP27dtXZnl2drZYt26d8PPzE2ZmZkKhUIhr165JlFJ3eXl5Ytu2beK3334rty42NlZMmTJFuLi4GOQ4ORlZJjZs2IBly5YhJiYGjRo1go2NDYqKijBt2jRMmTJF6nh6ExkZiRUrVuB///sfGjduLKtxqlQqrFy5Eh999BHMzMwwYsQIvPzyy2jQoAGuXbuGPXv2IC4uzqD3AADGM04AiIuLw7x583DmzBm4uLiga9eusLOzQ1FREaKiolBYWIjo6GipY+rEGLbn7du3YWZmhtq1a1c4KXfq1Kk4cuQIEhISJEqoH0VFRQAACwsLqFQqCCFQo8b/n5w9Z84c/PLLL0hMTJQoYeWw6BgwtVoNtVpd5oOYm5uLgwcPwtLSEp06dSo3c97QFBUVITk5GZmZmejVq5dmeU5ODg4cOAArKyt07tzZ4Mf5sOLiYixcuBDr1q3D+fPn0apVK+Tk5MDHxwczZ86Et7e31BH1wljGCQC//fYb1q9fjwsXLsDExATp6el47bXXMH78eLi4uEgdTy+MZXuWlJTAxMQEJib/zvy4f/8+XF1d8fbbb+PDDz+UOF3VUKvVuHfvHtq3b4/g4GDMnj1b6khaYdExQAUFBbCysiqzrPSKq6U/fHLw22+/Yd68eUhJSYEQArdv30bv3r0RGhqKfv36SR1PrzIzM3H27FlYWlqie/fumuXp6ek4fvw4WrRogXbt2hn8dUmMYZwPHjxAXFwc/vjjDzRr1gzu7u7o0KEDgH9/dq9du4Z27drJ4kJ6xrA979+/j5SUFKjVarRr106zXPx7Mg/y8/OxatUqTJw4EZaWlhImrTwhBC5fvqzZU961a1fUrl1bs16lUsHExATx8fHo2LGj4W1PaY6YkS6mTJki1qxZI06dOiXy8vLKrFOr1aKkpETk5uZKlE5/GjduLKZOnSo2bdok9u/fL5YvXy78/f1FzZo1hY+Pjzh+/LjUEfXiq6++EvXr1xdubm7Czs5O1KlTRwQFBYn4+Hipo+mVsYxz7Nixws7OTnTu3FnUq1dPmJiYiE6dOonvvvtO6mh6ZQzbMyIiQjg6OgpXV1fh5uYm2rdvLz744ANx6dKlMs97eO6VIZo+fbqoW7euaN++vbC2thaWlpZiwIAB5eYlGSoWHQOzZcsWoVAohJOTk+jWrZuYPn262L59u0hKStJc1+DevXuiT58+IiEhQeK0lbdlyxbh5ORU5heIWq0Wd+/eFQcOHBD9+/cXffv2Fbdv35Ywpe7Wr18vmjdvLpYuXSr+/PNPcejQITF//nzh5eUlatWqJcaMGSNu3rwpdUydGcs4165dK1q0aCGioqLE7du3RXFxsTh58qQIDg4WVlZWolWrViI6OlrqmDozhu25ceNG4eTkJMLDw8WWLVvE6tWrRVhYmOjUqZNo1aqVCA8PN6hryTzOunXrRIsWLcTWrVvFhQsXxJUrV8TWrVvFgAEDRI0aNYSvr684f/681DF1wqJjYEJCQkRwcLA4duyYmDlzpnBzcxOOjo6ib9++4vPPPxeHDh0SP/zwg7CwsJA6qk42bNggunbt+tgic+LECeHo6Cg2bdr0nJPpl7+/v3jvvffKLFOpVCI9PV1ERESIDh06iA8++ECidPpjLOMcNmyYmDhxoubxw9eOSUlJEQMHDhSDBg2SIppeGcP29PX1LTeGvLw88ddff4mPP/5YODk5iW+//VaidPrTv39/8e6771a47s8//xQ9evQQb7311nNOpV/ymdBhBFQqFZycnFC3bl10794d8+bNw9mzZxEREYEmTZpojhNPmzYNI0eOlDquTnx9fXH16lWEhITg77//LnfXZ29vb3h4eCA+Pl6ihLpTq9Vo3rw5cnNzyyw3MTGBvb09xo4di3HjxmHv3r1ISkqSKKXujGWcAODp6VnmjBQTExM8ePAADx48QLNmzTBlyhRcuHABUVFR0oXUkTFsz5KSEjRo0ABmZmZllltbW8PDwwOffvopRo8ejY0bNyIzM1OilLoTQsDNzQ3JyclllqvVaggh0KtXL4SGhuLYsWOIi4uTKKUeSN20SDvp6ema48OPHhe+d++eWLJkiVAoFLI4Tn7kyBHRrVs3MWDAALFgwQJx+PBhzbU4jhw5ImxtbcWxY8ckTqmbrVu3CoVCIT799FORnJxcbn1mZqaws7MTp0+ffv7h9MhYxvnXX38JGxsb0b9/f3H06NFy6+/fvy/q1asn4uLiJEinP8awPZcuXSrMzMzE2rVrK9yzfPXqVdGwYUNx7tw5CdLpz6FDh4RCoRAhISHi1KlT5dbn5uaK+vXri7/++uv5h9MTnnUlAyqVCgqFAiYmJlizZg2mTJmC/Px8qWPppPRjefjwYfzwww84fvw46tWrh7p16+LKlSswMTFB37598cMPP0icVHffffcdIiMj0aZNG/j6+qJ9+/Zo164dLCwssGDBAqxatQqpqalSx9TZd999h4iICLRp0wZ+fn6yHeeJEycwe/ZsZGVloVWrVvDx8UG/fv2gUqmwYMECHDt2DJcvX5Y6ps7kvj1LSkowc+ZM/Prrr/Dz88Mrr7yCli1bwt7eHmZmZli2bBm++uor3LhxQ+qoOtu+fTu++eYbWFtbo1OnTvD09ISPjw/Mzc0xb9487NixA9euXZM6ZqWx6BiQ0lP8HndKqhAC8+fPR0FBAebMmfOc0+mPWq0ud+rtjRs3sHfvXly7dg2Ojo5wdnaGv79/mWsIGar79+9j586dWL16NS5evIjGjRtDpVLh77//RufOnfHOO+9g1KhRUsesNLVaDRMTExQUFGDv3r1Ys2YNzp07h8aNG0OtVstmnMD/j/Xs2bPYs2cPTp48iZs3b+Ls2bMoKirCwIEDMWHCBAQEBEgdtdKMYXuWjjEnJweRkZFYsmQJkpOT0blzZzRt2hTHjx9H06ZNMWHCBISEhEgdVydCCKhUKkRHR2Pz5s04deoUFAoFrl+/jrS0NPj7++Odd94xvNs+PIRFxwCVbrInFR5Dvz4H8G+xU6lUqFGjhqyuD1SqokJ38eJFREVF4f79+3ByckKXLl3g5OQkYUrdFRYWoqCgAA0aNNAsS0pKwu+//47CwkI4OzvLYpwlJSXlindmZiaSk5NhZmYGMzMztGjRArVq1ZIooX4Yw/YUQiA3Nxe2traaZYmJidiyZQtycnLQrl079OrVC+3btzfo300qlQqmpqZllv3zzz9ISEiAWq2GnZ0d3NzcUL9+fYkS6geLjoFYt24d2rZti44dO5a5KFVF/1gasi+//BIdO3aEr69vmYsiPnjwAADKTQ6Ug5KSEqjVapiZmclmO5batm0bIiMj8ddff0GtVqN79+4YMmQIBg8eDGtra6njVYnS7WlwF1V7BsawPf/44w9ERETgwoULyM3Nhb+/P0aOHImePXtKHa3KyPl3EAAYbhU1IkePHkVwcDA++ugjTJ8+HWvXrsWFCxcAQHMoq6ioCJ9++qlBHy8+evQoPvzwQ3z22Wd47bXX8Nlnn+HEiRMAoPkf8b179zBlyhT8888/EqetvPPnz2PWrFk4evSoZo+Vubk5FAoFHjx4UO7u7IYqOjoa77//PmrWrImFCxfiv//9L7KzsxEUFIQOHTpg1apVUkfUixMnTqBnz5748ccfUVxcrNmewL+3RSi9Z1BWVhYM+f+VxrA9jx07htDQUKSkpGDYsGEYPHgwjh49Cj8/P3h4eOCXX34BAIPejsC/e6dGjRqFPXv2aPZElv4OKi4uls3vII3nPv2ZtDZ16lTh5eUlwsLChK+vr+jUqZMIDAwUH3zwgfjll1/E9evXRUxMjFAoFOWulGxIZsyYIXr27CkWLFgggoKCRM+ePYWPj48YMWKE+P7778XFixfFiRMnhEKhMOgrP7/55pvCyspK9OzZUwwbNkwsWLBAnDlzpsxzjh8/LsaPHy/UarVEKXU3fPhwERISUm75rVu3xPTp00WDBg3EwoULn38wPXvzzTeFmZmZcHJyEvXq1RPDhw8XBw4cKPOco0ePioCAAPHgwQOJUurOGLbn0KFDxbhx48osU6lUIjY2VowaNUq0aNFCbN++XaJ0+vPmm28KS0tL0bFjR+Hp6SmmTZtW7krzx44dEyNHjhQlJSUSpdQfw5/JaQRu376NHj16YMGCBSgqKsL+/fuxd+9eHD58GIcOHYKLiwvOnj2LPn36GPTu46ysLLRp0wZhYWEoKSnBiRMnEB0djYSEBPz000/Yvn07Ll++jICAANjY2Egdt9JOnz6Nd955Bw0aNEB8fDx+/vln/PLLL2jRogV69+6Nvn37Yt26dThy5IhB70YuKiqChYWF5nFxcTFMTExgZ2eH+fPnQwiB1atXY9iwYXB0dJQwqW6Sk5Px8ccfw9/fH7Gxsdi7dy9ef/111KpVC0OHDsX48eOxadMm3Lx506AnzxvD9szKyoK7u7vmcemk5C5dukCpVGL8+PGYN28eevXqZdDzVi5cuIDp06ejffv2iIuLQ1xcHPbt2wd7e3sEBgbi1VdfxYYNG3D27Nlyc3gMktRNi57uzJkzFd5z5MaNG2L16tVi2LBhQqFQiL1790qQTn9u3rwp/vjjj3LLb9++Lfbs2SOmTZtm8OO8dOmS8PX1FREREUIIobmlxUcffSQGDhwovL29ha+vr1AoFGLHjh3ShtXR+vXrRYMGDcTJkyfLLC+9WvCtW7dE8+bNxYkTJ6SIpxdpaWli3LhxYuXKlUKIf69tlZaWJqKiosTs2bNF165dRd26dYVCoRC7du2SOK1ujGF7Lly4UDRp0kQkJSWVWV66Z/Xq1auiVatWIjExUYp4enH16lUxYMAAsWLFCiGEEAUFBeKvv/4Sq1atEiEhIcLb21u0bdtWKBQKsXPnTonT6geLjoEovTigSqUSDx48KHNp+d27dwtbW1uJkunHo4doSsf58PJdu3YJKyur5x1N7y5fviwuXrxYbnlaWprYtm2b6Nevn6hTp44EyfRHrVaLvLw8MXLkSFGvXj3x5ptviu3bt4s7d+5onrNp0yZhbW0tXUg9ycrKEv/880+55ffv3xfXrl0T06dPl8XPZ15ennjttddkvT0zMzNF3759RZs2bUR4eLg4cuRImcPk27dvN/gxCvHvhWevXLlSbnl2drY4fvy4GDNmjMF/Zh/Gs64MmPi3qGL48OHIzc3FwYMHpY5UJUrHOXHiRGRnZ2PLli1SR9Kb0omqDx/SGDx4MKysrLBhwwYJk+lHfn4+IiMjsWvXLmRlZcHc3Bw2NjYQQiAtLQ0jRoww6Gs+VUQ8cnmHwYMHo0aNGti2bZuEqfQjLy8PERER+PXXX5GVlQVTU1PZbc9Lly5h2bJlOHr0KMzNzeHo6IhatWqhoKAA586dw0svvYT58+dLHVOvKvrMWlpaYtOmTRKm0h8WnWrswYMHOHfuHPbu3YvatWujU6dOcHZ2hr29PWrUqKG5BkJJSQlycnIM9phxYWEhTpw4ga1bt6JevXpo164dWrZsidatW6NevXqa56lUKuTl5aFOnTrSha1CQgjcvn0b3t7eWLduHXr06CF1JL25ceMGoqOjcf78eaSmpqKoqAihoaHw9PQsM+9DbvLz8zF58mRMnToVHh4eUsfRm4sXL+L48eO4du0arl+/jvv378tue545cwZ79uzBhQsXcOfOHRQWFuLdd99Fnz59DP5aSE9y9+5dDB06FF999RW8vLykjqMXLDrV2PTp07Fx40Y0bNgQ2dnZSE1NhYuLC15//XVMnToVdnZ2UkfUi7fffhu///47nJ2dkZWVhdTUVDRu3Bg9e/ZEaGgoOnXqJHVEvbhx4wY2btyI//3vf2jVqhXc3NzQrl07tGrVCtbW1mWuOPvwNYQMTVJSEpYsWYL4+HhNYe3WrRu8vb1l+Q/E065l9egkXkNVesrxw3sfSz+zhq70ZzMmJgYtW7aEh4cHvL290bx5c6hUKhQWFhr0CRCPevi2QRUpLCyU1c8qi041de7cOXTr1g2bNm1Cp06dYG9vj9TUVKxZswarV69Gfn4+lixZgtGjRxv0lZDPnTsHb29v7N+/H507d0bNmjVx584drF+/HitWrMD58+fx5Zdf4r333jPoCyNeu3YNI0aMQHZ2Njp37ozTp08jMzMTTZs2RWBgIGbMmFGmuBrqNr169SoGDhwIGxsbdOvWDX///TcyMjJgamoKd3d3TJ48GV26dJE6ps4yMzNx4sQJDBgwoMxZKY/+w3///n1YWloa7Pa8ffs2Ll68iO7du2uWqdVqlJSUwMTEBDVq1IAQAmq12mDPzqnoZzMjIwMODg4IDAzEhx9+WGbPsqGqaFuWbruHby1UeqTAUD+zFXp+04FIG3PnzhW9evXSPH74Wgb5+fli6tSpokOHDiIzM1OKeHqzYMEC8cILL2geFxUVlVk/f/580bx5c3H16tXnHU2vJkyYIAYMGCBSU1M1y5KTk0V4eLho0KCBaNy4cblrrxii//znP2LQoEFlJqmmpaWJpUuXik6dOgkbGxuxfv166QLqSWhoqFAoFMLOzk4EBQWJY8eOlVmvVqtFcnKymD9/vrh//75EKXU3ZcoUoVAoROvWrcWMGTPEhQsXyqwvKSkRKSkpYsuWLQZ7vZWn/Ww2atRI7N+/X8KE+vG0balSqTTb0pCv91QRw9/nKFPt2rVDZmam5grApXNxiouLYWVlhdDQUCgUCvz8888SJ9WNu7s7rl+/jqNHjwIAzM3NUVJSgnv37gEAxowZA3t7e4OfFHfmzBn06dMHTZs21Vz92NnZGZ988gnS09PRtWtXLFu2DIBhX3U1OTkZXbp0QZ06dTT3KnNwcMA777yDhIQEjBo1CitWrEBRUZFBjzMuLg7Tpk3D7NmzkZSUhF69esHJyQmzZs3C1atXoVAosHbtWixfvtygD1uVXpV96NChOHjwILy8vODp6YmFCxfizp07MDU1RWRkJD744AOD3aPztJ9Nb29vLF++HIBh/2w+bVuamJhotqUhX++pQhIXLXqMrKws0bZtW+Hq6iq2bdtW4f8KO3bsqLkWgqEqLCwU/v7+wsHBQfzwww+isLCw3HPc3d2FUqmUIJ3+zJ49W3Tp0qXMdiwuLtaMNyoqSrRs2bLcNUoMzbfffiuaN29e5tTVoqIizbgTExNF8+bNxZ9//ilVRJ1dv35dDB8+XHPtnNzcXBEXFyf++9//aq4/0rlzZ2FjY2PQVwq+du2aCAgIEOvWrRNFRUXi0qVLYuvWrSIkJES0aNFC1K5dWwQEBIj69euLb7/9Vuq4lWYMP5vGsi0fh0WnGktLSxMjR44UHTt2FP379xfh4eHi8OHDIjk5WYSFhYn69euL/Px8qWPqrLCwUEybNk00b95cuLm5ibfeekvs2LFDHD58WIwZM0Y4ODgY/DhjY2NFo0aNRJcuXcTu3bvLrb948aKwsLAQBQUFEqTTnytXrggPDw/h4uIiIiMjy60/e/asMDMzM+hxFhQUiF27dpW7MJ5KpRJZWVkiKipKDBw4UJiamlZY3A1FTk6OiIyMFIcPHy6z/O7duyIxMVGsXr1a9OzZ0+DHaQw/m8ayLR+Hk5GruczMTOzbtw8HDx5ESkoKLl68iKysLPj5+eHtt9/GG2+8IXVEnZROhMvLy8PRo0dx9OhRxMbGIi4uDiUlJfD398e4ceMwYMAAqaPqLCkpCR988AHi4uJQv3599OjRA/3798fFixexceNGODo6yuZaKzNnzsSmTZtQUlKCvn374qWXXsLZs2dx+PBhdOjQAT/++KPUMfVGVDBpc8yYMUhJSUF0dLREqfRLCKG5Ae3DRo4ciczMTPzxxx8SJdMPY/nZBOS/LSvColMNZWRkIDk5GRYWFqhZsyZcXFxgYmKCK1euoLCwEFZWVrCzszP4MwEqOjW1uLgYt27dQs2aNXH//n3Y2toa9KnWjyooKEBUVBQOHTqE2NhYnDlzBvXr18e4ceMwevRoODs7Sx1RJ6Xb9P79+zhz5gyio6Nx6NAhxMfHw8XFBaNGjcLQoUPRuHFjqaNW2tNOJ7937x5eeeUVvPPOOxgyZMhzTle1xENn6dy7dw++vr6YOXMmhg0bJnU0ncn9Z/NRct6Wj2LRqWZWrlyJiIgIJCQkoEaNGmjTpg3atWuHF198ES+//LLBXhTwSUr/h2Fqaiqf0xkfsm/fPty5cwcqlQqOjo7o2rUrrKysUFhYCFNTU+Tl5cnmmkiPUqlUmlNXc3JyYGtrK3WkKvfgwQPExcXBx8dH6ig6Kf3clpSUoEGDBvD29i7z+6eoqAgHDx406L2txvKzaQzb8klYdKqR27dvo1WrVggNDUVISAhyc3Oxb98+REVF4fLly2jfvj0WL16M5s2bG/Q1Du7cuYPOnTtj2LBhCA4Ohpubm2bdwxeyOn/+PJo2bWqwF+rKy8vDf/7zH/z+++8oKSlB48aNYWVlhfr166Nfv3549dVX0bRpUwCGfeG1kpISZGdno2HDhlJHqVLGMs5HP7cODg6wtrZG/fr14efnhxEjRsDJyUnqmDoxlp9NY9iWz0SCeUH0GIsXLxbe3t4Vrjt06JDw8vISrq6uZa5RYogWL14sFAqF6Nixo1AoFKJdu3bi66+/Funp6ZrnpKamCg8PjwpvPGco5s6dKzp06CCio6OFEP/ehX758uVi1KhRomPHjuLVV18Vd+/elTil7hYuXCjq1KkjJk2aJKKjoyuctJmTkyP27dunuTmtIXrWce7Zs6fc9aAMyZM+t+7u7mLEiBEG/7k1lp9NY9iWz4JFpxpZunSpcHNzE+fPnxdCCHHv3r0yvzDPnz8vWrduLbZs2SJVRL0YN26cGD9+vEhJSRHHjh0TkydPFo6OjsLExET4+vqKjRs3ikWLFolatWpJHVUnPXr0EIsWLSq3XKVSiQMHDohmzZqJwYMHS5BMv7p27Sq6d+8uvLy8hImJiWjbtq0IDw8XZ86c0VxEbunSpY8t8YbCWMZpDJ9bYxijEMYzzqcxzP1xMvXqq6/CxMQE33//vebS8ebm5lCr1QCAtm3bon79+khJSZE4aeUVFRXBzc0Nzs7OaNasGbp3746FCxfif//7H37++Wc0atQIkydPxrRp0/DBBx9IHbfSHjx4ADc3N/zyyy+4ffs2gH8PfZTOWenXrx+USiWSkpJw9uxZidNW3q1bt2Bubo533nkHJ0+exNmzZzFkyBBERkbCw8MDvr6+WL58OZYuXQpvb2+p41aasYzTGD63xjBGwHjG+Uykblr0L5VKJdRqtfj5559F06ZNRe3atUVISIhISEgQQghx48YN8dNPPwlra2uRnJwsbVgd3b9/X6SlpQkh/h33w4qLi8W+ffuEQqEQ169flyKe3sTExIiWLVuK//73vyIrK6vc+tTUVGFlZWXQ47xx44b49ttvy92+oqSkRERHR4uxY8cKW1tboVAoylxi39AYyziFMI7PrTGMUQjjGefTcDJyNVNUVIQrV67gzz//xM6dO3H06FEoFAo0adIEDx48wKhRo/Dpp59KHVNnN2/ehLm5eYVnkX322WeIiIjA1atXJUimH0IIlJSUICIiAh9++CFUKhVGjBiB119/HU5OTkhMTMSuXbtw5swZxMXFSR1XJ6W366hZs2aFk+SnT5+OQ4cOISEhQYp4emMM4zSGz60xjBEwnnE+CxadaiArKwubN2/G/PnzUb9+fdSrVw9169ZF165d0alTJxQWFuLq1asIDAxEq1atDPZsq9JxfvPNN2jQoAFq164NBwcHvPzyyxgwYABq1qwJtVqNVatWwcHBAQMHDpQ6sl7cvXsXkZGR+Omnn5CYmAhbW1tYWlqic+fOmDVrFrp16yZ1xCpz//59eHh4IDg42KAPRT6NHMdpDJ9bYxgjYDzjfBwWnWrgrbfewqlTpxAYGAhra2vcvn0bSUlJSEtLg5OTE+bMmQNXV1epY+rs4XHa2Njg9u3bOH/+PFJTU9GqVSuEhYUZ/LVHgH//51+zZs0yy4QQuHfvHvLz83HmzBlYW1sb9FwOoOJxVvScLVu24PXXX4e5uflzSqZfxjxOuX1ujWGMgPGM81mx6EhMCAFra2vs27cPvr6+mmVJSUk4cuQIVq1ahezsbGzbtg3t27eXOG3lPW6cV65cwZEjR7By5Urk5ORgy5YtZa6rY4jee+899OjRA56enmjUqFGFd6++c+cO6tata9DXQ3qWcd69exd16tR5/uH0iOP8f4b+uTWGMQLGM85n9rwmA1HFzp49K9q3by9iY2MrXF9YWCg6duwowsPDn28wPTOWcW7YsEEoFAphZmYmmjdvLqZNmyYOHTok0tPTNdeQycnJEa+88oo4ffq0xGkr73HjzMjIEA8ePBBCCJGfny8GDRokzpw5I3HayjP2ccrpc2sMYxTCeMapDRYdiRUWFoo+ffqIXr16iatXrwq1Wl3uOQsWLBCenp4SpNMfYxnnuHHjxDvvvCOuXLki5s6dK5ydnYVCoRCdO3cW8+bNEwkJCWLNmjWiRo0aUkfVCcfJcRoaYxijEMYzTm2w6FQDx48fFx4eHqJHjx5i/fr14saNG6KwsFAI8e+p2K+++qp44403JE6pO7mP88GDB+Lzzz8Xs2bNKrP81KlTYvz48cLW1lZYW1sLMzMzERwcLFFK3XGcHKehMYYxCmE849QWi041cfr0afHqq68KS0tLYWdnJwYPHiz+85//iObNmwsvLy9x6tQpqSPqhdzHeefOHXHhwgUhhBBFRUXl9lytX79eKBQKkZiYKEU8veE4/8VxGg5jGKMQxjNObdSQeo4Q/atDhw7YsmULMjMzsWfPHuzYsQPZ2dkIDg7G8OHD0a5dO6kj6oXcx1mnTh3NpNTSs2/UajWEEDA1NUVhYSEsLS3h7u4uYUrdcZwcp6ExhjECxjNObbDoVDMNGzbEW2+9hbfeesug75r7NMYyTgBlxpaXl4c5c+ZImKbqcJzyYgzjNIYxAsYzzsfh6eVEz9GDBw9gamoq62IHcJxyYwzjNIYxAsYzzoex6BAREZFsGU+lIyIiIqPDokNERESyxaJDREREssWiQ0RERLLFokNERESyxaJDREREssWiQ0RERLLFokNERESyxaJDREREsvV/DsB/6ZAUXcAAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plot_samples(data[2][1]) # plot the samples for (j,n) = 1,6, i.e. the configuration which contains both a physical and minimum-energy solution" ] }, { "attachments": {}, "cell_type": "markdown", "id": "cd31e406-cf82-4ef2-8df0-a17cf06a8d15", "metadata": {}, "source": [ "We can see that the best path is $b=(1,0,0,1,0,0)$, that is, the path $saf$." ] }, { "cell_type": "markdown", "id": "fbb8e055-474a-4456-823f-9348e06c6176", "metadata": {}, "source": [ "## Bibliography\n", "[1] Thomas Krauss and Joey McCollum. “Solving the Network Shortest Path Problem on a Quantum Annealer”. In: IEEE Transactions on Quantum Engineering 1 (2020), pp. 1–12. doi: 10.1109/TQE.2020.3021921.\n", "\n", "[2] Kamil Bradler and Hugo Wallner. Certain properties and applications of\n", "shallow bosonic circuits. 2021. doi: 10.48550/ARXIV.2112.09766. url:\n", "https://arxiv.org/abs/2112.09766." ] }, { "attachments": {}, "cell_type": "markdown", "id": "3065aa56-30e8-4aee-8b6c-f8aff1615c01", "metadata": {}, "source": [ "## Appendix\n", "If you want to try more complicated graphs, this code generates the objective function for any weighted, directed graph:" ] }, { "cell_type": "code", "execution_count": 18, "id": "dcfdc370", "metadata": {}, "outputs": [], "source": [ "def alpha_graph(graph): \n", " \"\"\"calculate alpha as sum of all the costs + 1\"\"\"\n", " alpha = 0 \n", " for i in range(0,len(graph)):\n", " for j in range(0,len(graph)):\n", " if graph[i][j] != 0 : alpha+=graph[i][j]\n", " return alpha+1\n", "\n", "def graph_to_s(graph, s, weight):\n", " \"\"\"calculate Hs\"\"\"\n", " ham = np.zeros((weight, weight))\n", " edges_coeff = [0]*weight\n", " constant_coeff = -1\n", " count = 0\n", " for i in range(0, len(graph)):\n", " for j in range(0, len(graph)):\n", " if graph[i][j] != 0 :\n", " if i == s:\n", " edges_coeff[count] = 1\n", " if j == s:\n", " edges_coeff[count] = -1\n", " count += 1\n", " \n", " for i in range(0, weight):\n", " for j in range(i+1, weight):\n", " ham[i][j] = edges_coeff[i]*edges_coeff[j]*2\n", " \n", " for i in range(0, weight):\n", " ham[i][i] += edges_coeff[i]*edges_coeff[i] + constant_coeff*edges_coeff[i]*2 \n", "\n", " return ham\n", "\n", "def graph_to_f(graph, t, weight):\n", " \"\"\"calculate Hf\"\"\"\n", " ham = np.zeros((weight, weight))\n", " edges_coeff = [0]*weight\n", " constant_coeff = +1\n", " count = 0\n", " for i in range(0,len(graph)):\n", " for j in range(0,len(graph)):\n", " if graph[i][j] != 0 :\n", " if i == t:\n", " edges_coeff[count] = 1\n", " if j == t:\n", " edges_coeff[count] = -1\n", " count += 1\n", " \n", " for i in range(0, weight):\n", " for j in range(i+1, weight):\n", " ham[i][j] = edges_coeff[i]*edges_coeff[j]*2\n", " \n", " for i in range(0,weight):\n", " ham[i][i] += edges_coeff[i]*edges_coeff[i] + constant_coeff*edges_coeff[i]*2 \n", "\n", " return ham\n", "\n", "def graph_to_i(graph, our_i, weight):\n", " \"\"\"calculate Hi\"\"\"\n", " ham = np.zeros((weight, weight))\n", " edges_coeff = [0]*weight\n", " count = 0\n", " for i in range(0, len(graph)):\n", " for j in range(0, len(graph)):\n", " if graph[i][j] != 0 :\n", " if i==our_i:\n", " edges_coeff[count] = 1\n", " if j==our_i:\n", " edges_coeff[count] = -1\n", " count += 1\n", " \n", " for i in range(0, weight):\n", " for j in range(i+1, weight):\n", " ham[i][j] = edges_coeff[i]*edges_coeff[j]*2 \n", " \n", " for i in range(0, weight):\n", " ham[i][i] += edges_coeff[i]*edges_coeff[i]\n", "\n", " return ham\n", "\n", "def graph_to_c(graph, weight):\n", " \"\"\"calculate Hc\"\"\"\n", " ham = np.zeros((weight, weight))\n", " count = 0\n", " for i in range(0,len(graph)):\n", " for j in range(0,len(graph)):\n", " if graph[i][j] != 0 :\n", " ham[count][count] = graph[i][j]\n", " count += 1\n", "\n", " return ham\n", "\n", "def graph_to_hamiltonian(graph,s,f): \n", " \"\"\"returns the graph adjacency matrix as a Hamiltonian\"\"\"\n", " weight = np.count_nonzero(graph)\n", " ham = np.zeros((weight, weight))\n", " alpha = alpha_graph(graph)\n", " ham += alpha*graph_to_s(graph, s, weight)\n", " ham += alpha*graph_to_f(graph, f, weight)\n", " for i in range(0,weight):\n", " if i!=s and i!=f:\n", " ham += alpha*graph_to_i(graph, i, weight)\n", " ham += graph_to_c(graph, weight)\n", "\n", " return ham" ] }, { "cell_type": "code", "execution_count": 19, "id": "7bbfb18f", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 2. 32. -32. -32. 32. 0.]\n", " [ 0. 1. 32. 0. -32. -32.]\n", " [ 0. 0. 35. 32. -64. -32.]\n", " [ 0. 0. 0. 2. -32. 32.]\n", " [ 0. 0. 0. 0. 35. 32.]\n", " [ 0. 0. 0. 0. 0. 4.]]\n" ] } ], "source": [ "\"\"\"example: our 5-edge graph\"\"\"\n", "# start,a,b,finish -> 0,1,2,3\n", "G1 = np.array([[0,2,1,0], [0,0,3,2], [0,3,0,4], [0,0,0,0]]) # weights associated to all possible paths\n", "G1_s = 0 # start index\n", "G1_f = 3 # finish index\n", "H_example = graph_to_hamiltonian(G1, G1_s, G1_f)\n", "print(H_example)" ] } ], "metadata": { "language_info": { "name": "python" } }, "nbformat": 4, "nbformat_minor": 5 }