tutorials/combinatorial_optimization/ArbitrageOpportunityOptimation_CN.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 量子金融应用：最佳套利机会\n",
    "\n",
    "<em> Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved. </em>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 概览\n",
    "\n",
    "当前量子计算应用到金融问题上的解决方案通常可分为三类量子算法，即量子模拟，量子优化以及量子机器学习 [1,2]。许多的金融问题本质上是一个组合优化问题，解决这些问题的算法通常具有较高的时间复杂度，实现难度较大。得益于量子计算强大的计算性能，未来有望通过量子算法解决这些复杂问题。\n",
    "\n",
    "量桨的 Quantum Finance 模块主要讨论的是量子优化部分的内容，即如何通过一些量子算法解决实际金融应用中的优化问题。本文主要介绍如何使用量子算法求解最佳套利机会问题 [3]。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 最佳套利机会问题\n",
    "\n",
    "套利（Arbitrage）一般指相同的资产在不同的市场中有不同价格，可以通过在多个市场间的转换，获得正向收益。也就是说，如果给定一组资产和交易成本，那么就可能在不同市场间存在一个交易循环产生正向收益。\n",
    "\n",
    "这个问题可以用图论的语言来描述：给定一个带权重的有向图 $G$ ，顶点 $i$ 表示一种市场货币，从顶点 $i$ 到顶点 $j$ 的边的权重表示将货币 $i$ 转换为货币 $j$的转换率 $c_{ij} $。也就是说如果拥有数量为 $x_i$ 的货币 $i$，那么通过交易转换我们将得到 $c_{ij}x_i$ 的货币 $j$。通常来讲转换率不是对称的，即 $c_{ij} \\neq c_{ji}$，并且我们假设交易成本（服务供应商收取的服务费，买入和卖出时的手续费等）是已经包含在该转换率中的。\n",
    "\n",
    "我们固定可交易的货币种类数量为 $n$，则该优化问题就是在给定的加权有向图上寻找一个获利最大回路，该回路中包含的顶点数量为 $K (K \\leq n)$。找到的最佳套利回路就是从某一市场出发，使得资产转换一圈后收益最大的交易路径。在该问题中，用户可以自行定义套利回路中包含的顶点数 $K$。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 编码最佳套利机会问题\n",
    "\n",
    "为了将最佳套利机会问题转化成一个参数化量子电路（parameterized quantum circuits, PQC）可解的问题，我们首先需要编码该问题的哈密顿量。第一步是将该问题转化为一个整数规划问题：给定一个整数 $n$, 即最佳套利机会问题涉及的货币种类数量, 也就是图形 $G$ 的顶点数量为 $|V| = n$ 个，那么对于每个顶点 $i \\in V$，我们定义 $n$ 个二进制变量 $x_{i,k}$，$k \\in [0,K-1]$：\n",
    "\n",
    "$$\n",
    "x_{i, k}=\\begin{cases}\n",
    "1, & \\text{资产转换到市场 $i$ 的次序为 $k$}\\\\\n",
    "0, & \\text{其他情况}\n",
    "\\end{cases}.\n",
    "\\tag{1}\n",
    "$$\n",
    "\n",
    "参考旅行商问题（[TSP 教程](./TSP_CN.ipynb)），因为 $G$ 有 $n$ 个顶点，所以我们共有 $n^2$ 个变量 $x_{i,k}$，所有这些变量的取值我们用 $x = x_{0,0} x_{0,1}\\dots x_{n-1,K-1}$ 来表示。在这里我们假设 $x$ 对应一个套利回路，那么对于图中的每一条边 $(i,j,w_{i,j})\\in E$，条件 $x_{i,k}·x_{j,k+1} = 1$ 成立当且仅当套利回路中第 $k$ 个顶点为顶点 $i$ 且 第 $k+1$ 个顶点是顶点 $j$。可以通过下面的式子计算资产在套利回路转换一圈后所获利的对数：\n",
    "\n",
    "$$\n",
    "P(x) = - \\sum_{i,j\\in V} \\log(c_{ij}) \\sum_{k=0}^{K-1}  x_{i,k}x_{j,k+1}, \\tag{2}\n",
    "$$\n",
    "\n",
    "根据套利回路的定义，$x$ 如果对应一条套利回路，那么就需要满足如下的限制：\n",
    "\n",
    "$$\n",
    "\\sum_{i=0}^{n-1} x_{i,k} = 1 \\quad \\text{和} \\quad \\sum_{k=0}^{K-1}\\sum_{(i,j)\\notin E}x_{i,k}x_{j, k+1}, \\tag{3}\n",
    "$$\n",
    "\n",
    "其中第一个式子用来保证找到的 $x$ 所代表的回路中同一时刻只可能出现一个顶点，第二个保证不存在的边不会出现在最后的套利回路中。这两个式子共同保证了参数化量子电路找到的 $x$ 是个简单回路。因为后续要用梯度下降方法对损失函数进行优化求最小值，对上述式子稍作调整给出损失函数：\n",
    "\n",
    "$$\n",
    "C_x = - P(x) + A\\sum_{k=0}^{K-1} \\left(1 - \\sum_{i=0}^{n-1} x_{i,k}\\right)^2 + A\\sum_{k=0}^{K-1}\\sum_{(i,j)\\notin E}\n",
    "x_{i,k}x_{j,k+1}.\n",
    "\\tag{4}\n",
    "$$\n",
    "\n",
    "在这个式子中 $V$ 是图中顶点个数，$E$ 表示图中边的集合，$K$ 是套利圈中顶点的个数，$A$ 是惩罚参数，它保证了上述的限制被遵守。因为我们想要在找套利回路 $P(x)$ 最大值的同时保证 $x$ 确实表示一个回路，所以我们需要设置一个大一点的 $A$，最起码大过图 $G$ 中边的最大的权重，从而保证不遵守限制的路线不会成为最终的路线。\n",
    "\n",
    "我们现在需要将损失函数 $C_x$ 转化为一个哈密顿量从而完成最佳套利问题的编码。每一个二进制变量可以取 $0$ 和 $1$ 两个值，分别对应量子态 $|0\\rangle$ 和 $|1\\rangle$。每个二进制变量都对应一个量子比特，所以我们需要 $n^2$ 个量子比特来解决最佳套利机会问题。因为我们的变量 $x_{i,k}$ 的值为 $0$ 和 $1$，所以我们要构造一个本征值和它对应的哈密顿量。泡利 $Z$ 的本征值为 $\\pm 1$，于是我们构造的哈密顿量为 $\\frac{I-Z}{2}$， 对应的本征值即为 $0$ 和 $1$。我们现在将二进制变量映射到这个哈密顿矩阵上，从而使 $C_x$ 转化成哈密顿矩阵：\n",
    "\n",
    "$$\n",
    "x_{i,k} \\mapsto \\frac{I-Z_{i,k}}{2}, \\tag{5}\n",
    "$$\n",
    "\n",
    "这里 $Z_{i,k} = I \\otimes I \\otimes \\ldots \\otimes Z \\otimes \\ldots \\otimes I$，也就是说 $Z$ 作用在位置为 $(i,k)$ 的量子比特上。通过这个映射，如果一个编号为 $(i,k)$ 的量子比特的量子态为 $|1\\rangle$，那么对应的二进制变量的取值为 $x_{i,k} |1\\rangle = \\frac{I-Z_{i,k}}{2} |1\\rangle = 1|1\\rangle $，也就是说顶点 $i$ 在最佳套利回路中的交易顺序是 $k$。同样地，对于量子态为 $|0\\rangle$的量子比特 $(i,k)$，它所对应的二进制变量的取值为 $x_{i,k}|0\\rangle  = \\frac{I-Z_{i,k}}{2} |0\\rangle = 0 |0\\rangle $。\n",
    "\n",
    "我们用上述映射将 $C_x$ 转化成量子比特数为 $n^2$ 的系统的哈密顿矩阵 $H_C$，从而实现了最佳套利机会问题的量子化。这个哈密顿矩阵 $H_C$ 的基态即为最佳套利机会问题的最优解。在接下来的部分，我们将展示如何用参数化量子电路找到这个矩阵的基态，也就是对应最小本征值的本征态。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Paddle Quantum 实现\n",
    "\n",
    "要在量桨上实现用参数化量子电路解决量子金融中的最佳套利机会问题，首先要做的便是加载需要用到的包。其中 `networkx` 包可以帮助我们方便地查看和处理图。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:00:15.901429Z",
     "start_time": "2021-05-17T08:00:12.708945Z"
    }
   },
   "outputs": [],
   "source": [
    "#加载需要的包\n",
    "import numpy as np\n",
    "import networkx as nx\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "#加载飞桨，量桨相关的模块\n",
    "import paddle\n",
    "import paddle_quantum\n",
    "from paddle_quantum.ansatz import Circuit\n",
    "from paddle_quantum.finance  import arbitrage_opportunities_hamiltonian"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来，我们先生成用于最佳套利机会问题的图 $G$。为了运算方便，图中的顶点从0开始计数。\n",
    "\n",
    "同时，为了方便运算验证，这里构造的加权有向图中对权值做了特殊设定。在实际使用过程中用户可以自行构造需要计算的货币类型，设置真实的货币汇率。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:00:16.212260Z",
     "start_time": "2021-05-17T08:00:15.918792Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1MAAADnCAYAAAD7CwxiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAA2QUlEQVR4nO3deXyU9b3+/2vuzGSyESBsCoRAZAcVMSiyqECoJMX14I76KIpo9RztofagPWorHvVnqdoecYFqF2xVpP1p1SQKhCMISlhUdkQSIAFkixCyz2Tu7x80I5FlJmEm98w9r+c/JZnJ5DP46Ly57utz37fDNE1TAAAAAIBmMaxeAAAAAABEI8IUAAAAALQAYQoAAAAAWoAwBQAAAAAtQJgCAAAAgBYgTAEAAABACxCmAAAAAKAFCFMAAAAA0AKEKQAAAABoAcIUAAAAALQAYQoAAAAAWoAwBQAAAAAtQJgCAAAAgBYgTAEAAABACxCmAAAAAKAFCFMAAAAA0AKEKQAAAABoAcIUAAAAALQAYQoAAAAAWoAwBQAAAAAt4LR6AZHK5zO1s7xaJQcrVevxydPgkyvOUILLUK+OKcpIS5JhOKxeJgAArYr5CADfI0z9i89navn2gyrcsl+rdpRr2/5KGQ6HnIZDpkyZpuRwSA455PWZ8pmm+nRO0bCeaRrbv7NGntOR4QEAsB3mIwCcmsM0TdPqRVjpSI1H81eXau6yYlXVeVVd36Dm/IU4JCXFxynZ7dTU0Zm6IStdbRNd4VouAACtgvkIAIHFbJiqqW/Q0/mb9fbqUjkcUq3Hd8avmegy5DOlG7PS9XDOACXGx4VgpQAAtB7mIwAELybDVFFJue5/c60qajyq9Z75kPihBKeh1ESXZt8yVMN6poX89QEACAfmIwA0T0yFqTpvg2Z+sEkL1paF5EhbIAkuQ5OGdtejEwfK7eQoHAAgMjEfAaBlYiZMVdV5Nfm1ldq8tyIsR9tOJcFpaGDXVM2bcrGS3VzvAwAQWZiPANByMRGmquq8mvTqChUfqFJdKw6KRm6nocxOyVowbQQDAwAQMZiPAHBmbH/T3jpvgya/ttKyQXFsDT4VH6jSba+vVJ23wZI1AABwPOYjAJw524epmR9s0ua9FZYNikZ1Xp827anQzA82W7oOAAAk5iMAhIKtw1RRSfmxk2ktHhSNar0+LVhbqlU7yq1eCgAghjEfASA0bBumauobdP+ba1vlqkTNUevx6b6/rVVNPdsZAACtj/kIAKFj2zD1VP5mVdR4rF7GSVXUePR0AdsZAACtj/kIAKFjyzB1pMaj+atLI2b7wg/Ven16e1WpjkToMAMA2BPzEQBCy5Zhav7qUjkcVq/i9AyH9M7qUquXAQCIIcxHAAgt24Upn8/U3GXFEbcX/IdqPD7NWVYsn8/2t/kCAEQA5iMAhJ7twtTy7QdVVee1ehlBqazzakXxIauXAQCIAcxHAAg924Wpwi37VR0lVwKq8TSocMt+q5cBAIgBzEcACD3bhalVO8oVLRsDTFNatYMjbwCA8GM+AkDoOa1eQCj5fKa27a8M+vkXZrTXvZedowsz2ivJFadvK2q1ZOsB/U/eJm17MleStHHPEf34fz+VJP1oYBfNuS1LZd9Va9SzS/T21OG6OLODFm7ap6nzVkuS2ie5VDj9crVPitfMDzfptU9LTruGr/dVyjRNOSL9jGAAQNRiPgJAeNiqmdpZXi0jyA/dK887W29PHa7sAV2093CN/vHFbu0qr9atF/dQoivO/7xBXdvqikFnnfQ1Hv7/16vO06DxA7soe0BnSdIjOQPUPileX5Ye1h+Xn35QSJLhcGjnoeqg1gwAQEswHwEgPGzVTJUcrJTTCDwsElyGnrh6sJxxhv7xRZmmv/OVzH/tfeiRlqQaz/d7yn0+Uw+O66OPNn57wusUH6zS7wu/0UNX9NPjEwep1uPTvw3tLk+DTzP+sU7BXIjIaThUcrBKPTsmB/0+AQBoDuYjAISHrZqpWo9PZhA7wrMy0tQ+KV6S9GLhN/5BIUm7yqvlafj+G3kb9mrA2anKHXzyo2+vLt2uzXsrlJ6WpNfuyJJhODRnabG2fHs0qDWbkmq90XFCMAAgOjEfASA8bBWmPA2+Jh/8p9IhJd7/592Ha0773D+t2KFDlXV6YFzfk+7b9vpMzfjHOjX4TLmdcSo+WKnfF24Les2maao+Qu9EDwCwB+YjAISHrcKUK84I6s7uhyrr/X/u1i7xtM+trm/QnGXF6ndWG0089+yTPuersiNas/M7SdLfVu5SXTM+/B0Oh+KdtvrPAACIMMxHAAgPW31KJbgMORR4WqzZ+Z0OVx8bGPeP7d1kwHRrl3jCvvI/f7ZDB47W6cenGBaS5PvXIT9fMIf+juOQlOCMC/g8AABaivkIAOFhqzDVq2OKvEGc1VrjadDj729Ug8/UdRd014f3j9JT156r1+7IUuH0y5QU3/TDu9bj06tLt8sI4uTd5vL6TPXi5FoAQBgxHwEgPGwVpjLSkoI+8vXel3t089zPVbhlv7q2S9S/De2m3p1S9Naq0iZXK2o07/Od2l9RG+oly2eayuiQFPLXBQCgEfMRAMLDYZrN7N0j3MT/XaYNeyqsXkbQzu2WqvfvHy2fz6ft27fryy+/1Lp165SVlaWrr77a6uUBAGwi2uZjz1RDH08fp/j4+MBPBgCL2Oo+U5I0rGeaNu6pCOICsNZzOKRu8bXq3LmzDh8+LI/H438sPT2dMAUACJlomo+mz6d1i95Vh6du0Lhx45STk6OcnBz16NHD6qUBQBO22uYnSWP7dz5hT3ekSnTFqa5kjQ4cONAkSElSx44dtXv3botWBgCwm2iaj/FxUhffQVVWVuq9997TPffco4yMDA0aNEgPPfSQCgsLVV9fH/iFACDMbLfNz+czNfyZxdp/tM7qpQTUpY1by38xRrfeeovmz59/0uecd955ys3NVU5Oji655BK5XK5WXiUAwA6ibT5+NmOcdu8uU0FBgfLz87Vo0SIdPfr9DX9TUlJorQBYznZhSpLmLivWbxduVa0ncm/2l+gyNH18P901OlOmaerZZ5/VjBkz/I9nZ2drxYoVqq6u9n8vNTVV48ePV05OjiZMmKBu3bpZsXQAQJSKtvl4vPr6ei1fvlz5+fnKz8/Xhg0bmjw+cOBA/8HHUaNGca4VgFZhyzB1pMaji55a1KybA7Y2t9NQ0SPZapv4fdP02muv6e6779ZFF12kzz77THV1dVq6dKl/cGzZsqXJa9BaAQCaI1rn48mUlpbSWgGwnC3DlCQ99t4GzV9dqtoIHBgJTkM3DEvXE1cNPuGx4uJipaWlqV27dic8tmPHDn+wWrx4Ma0VAKDZonU+ng6tFQCr2DZM1dQ36LJZSyJyb3iXVLf+b/oYJZ7BicC0VgCAlrD7fJSatlYLFy5UZWWl/zFaKwChZNswJUmrdpTrttdXRtTe8ASXoTemXKysnmkhfV1aKwBAsGJpPtJaAQgnW4cpSfrvd9drwZqyiNjOkOA0NOnCdD15TfO2LzRXbW2tli1bRmsFADilWJyP0vetVV5enhYtWkRrBeCM2D5M1XkbdPPcz7VxT4WlJ9y6nYYGdU3Vm1OHy+1s3ft8lJSU+AdHYWEhrRUAgPkoWisAZ872YUqSquq8mvTqChUfqLJkYLidhjI7JWvBtBFKdjtb/fcfj9YKANCI+dgUrRWA5oqJMCUdGxi3vb5Sm/ZUtOqWhgSnoYFdUzVvysURMSh+iNYKAGIb8/HkArVWgwYN8gcrWisgdsVMmJKObWmY+cFmLVhb2ion3Sa4DE0amq5HJw5o9a0LLUFrBQCxifkYWGlpqX8+0loBaBRTYarRqh3luu9va1VR4wnLUbgEp6HURJdm3zJUw0J8VaLWRGsFALGF+RgcWisAjWIyTEnH7rPxdMFmvb2qVIZDqgnBkbhElyGfKd04LF0PTxhwxvfJiCS0VgAQG5iPzUdrBcSumA1TjY7UePTO6lLNWVasyjqvajwNas7fiMMhJbrilOJ26u7Rmbo+K11tE+0fJEpKSvyDg9YKAOyH+dgyja1VXl6e8vPztXHjxiaP01oB9hLzYaqRz2dq+faDWrL1gIpKDmnb/koZDofiDOno0aMyjDglJyfLIcnrM+UzTfXtkqJhPTtobP/OGpHZQYbhsPptWILWCgDs61Tz0Wk4ZEoyTVMOh4P5eAq7du1SQUEBrRVgU4SpU/D5TO0qr9aabaW6/Sd3ql1aR/1hzitKcMapV8dkZXRIksMRu8PhdGitAMC+GudjycEq1XobVO/1Kd5pMB+DQGsF2A9hKoC9e/eqa9euOuuss7R3716rlxN1aK0AwL7+53/+RzNnztScOXN0++23W72cqBNMa9U4I9PT0y1cKYBTIUwFQJgKrWBbq5ycHHXt2tXClQIATsfj8ahjx46qqKhQRkaGSkpKaKTOAK0VEJ0IUwEQpsKnsbVqHBxbt25t8vj555/vHxy0VgAQWV5//XXdeeed/q8LCgp0xRVXWLgiewnUWmVnZ/tnJK0VYB3CVACEqdZDawUA0cHj8ahfv34qKSnxf2/48OFasWIF7VQY1NfX69NPP/XPSForIHIQpgIgTFmD1goAItcbb7yh22677YTvFxYWasyYMRasKLbQWgGRgzAVAGEqMoS7taqtrVVCQkIolwwAtvXXv/5Vt912m47/J4TT6VRBQYHGjRtn4cpiD60VYC3CVACEqcgT6tbqD3/4g/76179q//79euaZZzRx4kS2qQBAABUVFZoxY4Zefvll/frXv9Z//ud/KiUlxeplxbxdu3b5g9XixYtprYAwI0wFQJiKfMG0Vj/5yU+Uk5MjwzCa/GxRUZFuuukmff7559q4caN+85vf6O9//7sSExP9z2loaFBcXFyrvR8AiBbTp0/Xc889p1mzZmn69OlWLwc/EI7WasOGDerWrZvat28frmUDUcUI/BQgsvXq1Us//elP9f777+vQoUP66KOP9OCDD6pfv36qqKjQ3//+d23dulU+n6/Jz1VXV+vDDz/UT37yE3Xu3FnnnHOOJKm8vNz/nIaGBs2ZM0cjRozQzJkztXbt2lZ9bwAAtFR8fLzGjh2r3/zmN9qwYYN27typV155RVdffbVSUlK0ceNGzZo1S+PGjVOHDh107bXXas6cOSorKzvp6x0+fFgvvfSSRo0apeHDh+ujjz5q5XcERB7CFGwlISFBP/rRj/T8889ry5YtKi4u1uzZs3XrrbfK6XQ2eW5ZWZm+/fZb//7+8vJyDRgwQF9//bX/Ofv379fevXs1cOBAeb1ePfbYY/rggw9a9T0BABAKPXr00LRp0/Tuu+/q0KFDWrx4sX7+859r0KBBqqys1Lvvvqtp06bpZz/7WZPtgY3cbreefvppbdy4UU888YReeumlJrtBgFjkDPwUIHo1tlYns3fvXlVVVen888+XJO3bt08VFRXq2bOn/zllZWU6dOiQHn30UWVkZEiSvF5v2NcNAEA4NbZWjc3V8eda3XTTTSc9/y0xMVFut1uS1L17dx06dEi1tbVKSkpq7eUDEYMwhZhVUVGhPXv2KDk5WRUVFdq2bZvS0tLUq1cv/3PKy8v14YcfqqSkRPfee69ycnJOaLgAAIh2ja3VtGnTdKrT6evr6/3nVc2dO1dDhw7lSriIeWzzQ8xq06aNOnbsKEn6+OOP9cknn+iaa66RJP/5VZdffrlWr16tm266SX/5y1+0a9cuq5YLAECrONUVbePj47Vnzx5NnDhRhw4d0v333+9vpfbt26fS0tLWXCYQEQhTiFmjRo1Senq6zjrrLP3xj3/Ufffdp8GDB0uS/6p/brdbHTt21O23367evXvrueeea/IabPkDAMSCqqoqvfXWW5o8ebIuvfRSvfTSS+rbt6+kYxd0evrpp9WjRw8NHjxYDz30kJYsWaL6+nqLVw2EH2EKMcvpdOq3v/2tvv32W/3pT3/S5Zdfrj//+c+aN2+eJGn58uXas2eP//lfffWV/3wq0zTl9Xr14osvasiQIXr44Ye1dOlSeTweK94KAABhU19frzvuuEMPPvigHnjgAf3iF79ock6Vy+VSfX19kysEjh07tskVAmmtYFfcZyoA7jMVezwej1wul1588UXNmTNHnTp10tlnn636+nr97ne/09lnny3p2FG6SZMmqaCgwP+zbdu21fjx45WTk6MJEyaoa9euVr0NAAg77jMVGxoaGjR//nzl5+dr+/btMk1TM2bM0JVXXtlkS+Dx97XKy8vTpk2bmrzOoEGDlJubq5ycHI0cOTKo+1oBkY4wFQBhCitXrtT27dt1/fXXy+VyNXmstrZWS5cu9V8BaevWrU0eP//88/03RLzkkktO+HkAiGaEqdi0cuVKOZ1OXXjhhad93s6dO1VQUKD8/HwtWrRIVVVV/sdSUlKUnZ3tn5Hp6enhXjYQFoSpAAhTaI7i4mJ/sCosLFRNTY3/MVorAHZDmEKw6urq/K1Vfn7+Ca3V4MGD/cGK1grRhDAVAGEKLRVsa5Wbm6tLLrmES64DiDqEKbQUrRXsgjAVAGEKoUJrBcBuCFMIBVorRDPCVACEKYRDMK1V40m6tFYAIhVhCuFAa4VoQpgKgDCF1kBrBSAaEaYQbrRWiHSEqQAIU2htgVqrIUOGNLlCIK0VAKsQptDaGlurvLw8LV68uElr1aZNG40bN86/s6N79+4WrhSxgjAVAGEKVqO1AhCpCFOwEq0VIgFhKgDCFCLJ8a1VXl6evv766yaP01oF5vOZ2llerZKDlar1+ORp8MkVZyjBZahXxxRlpCXJMByBXwgAYQoRhdbqzDAfW4YwFQBhCpEsmNYqNzdXEyZM0Nlnn23hSq3j85lavv2gCrfs16od5dq2v1KGwyGn4ZApU6YpORySQw55faZ8pqk+nVM0rGeaxvbvrJHndGR4AKdAmEKkorUKjPkYGoSpAAhTiBa0Vk0dqfFo/upSzV1WrKo6r6rrG9ScDzuHpKT4OCW7nZo6OlM3ZKWrbaIrXMsFohJhCtEiUGt1/BUC7d5aMR9DizAVAGEK0SpWW6ua+gY9nb9Zb68ulcMh1Xp8Z/yaiS5DPlO6MStdD+cMUGJ8XAhWCkQ/whSiUay2VszH8CBMBUCYgh3ESmtVVFKu+99cq4oaj2q9Zz4kfijBaSg10aXZtwzVsJ5pIX99INoQpmAHO3fu9Acru7ZWzMfwIUwFQJiCHZ2utWrXrl2TKwRGQ2tV523QzA82acHaspAcaQskwWVo0tDuenTiQLmdsXcUDmhEmILdBNta5ebmauTIkXK5Int7G/Mx/AhTARCmYHfR3lpV1Xk1+bWV2ry3IixH204lwWloYNdUzZtysZLdkfV3ArQWwhTsLppbK+Zj6yBMBUCYQqyJptaqqs6rSa+uUPGBKtW14qBo5HYayuyUrAXTRsTEwAB+iDCFWBJMa9V46XWrWyvmY+shTAVw5MgR/fKXv1Rqaqqeeuopq5cDtKpgW6vc3FwNHz68VVurOm+DbprzuTbtrbBkUDRyOw0N6pqqN6cOj5ktDUAjwhRiWaS2VszH1kWYCsA0TdXW1srhcCghIcHq5QCWiqTW6r/fXa8Fa8padevCqSQ4DU26MF1PXjPY6qUArYowBRzT2Frl5eUpPz9fmzdvbvJ4a7ZWzMfWRZgC0CK1tbX65JNP/OHqZK1V4+AIdWtVVFKu2/+4slVOpg1WgsvQvCkXx9xVjBDbCFPAyVnVWjEfWx9hCkBItFZrVVPfoMtmLdH+o3WhWHZIdW7j1ic/HxOT99lAbCJMAYEFaq3OPffcJve1amlrxXy0hmH1AgDYQ2Zmpu677z598MEHKi8vV0FBgR544AH17dtXhw8f1jvvvKMpU6aoa9euuuCCC/TLX/5Sn376qbxeb7N+z1P5m1VR4wnTuzgzFTUePV2wOfATAQAxw+12a9y4cfrtb3+rTZs2aceOHXr55Zd11VVXKTk5WevXr9ezzz6rMWPGqEOHDrruuus0d+5clZWVNXmdsrIyTZkyRStXrjzp72E+WoNmCkDYhaq1OlLj0UVPLbL0hNpA3E5DRY9kq21iZN97BAgFmingzDSntXrvvff0wgsvKDExUe+8845+/OMf+5/HfLQOYSpIDQ0NWrVqlXbv3i1J6tatmy666CIZBuUe0ByBzrW64IIL/IPjh+dazV1WrN8u3BpRe8F/KNFlaPr4frprdKbVSwHCjjAFhNbpzrUyDEM+37H5FxcXpz/96U+aPHmyJOajlQhTQSgoKND999+v3r17q0ePHpKk0tJSbdu2TS+++KImTJhg8QqB6BVsa/WjK67QtX/aFJF7wX+ocxu3Pp8xTobhsHopQFgRpoDwqaur07Jly5Sfn6/33ntP27dvP+E59957r158cbaGP7OY+WgRwlQQ+vXrp48//lgZGRlNvr9z505dccUV2rJli0UrA+zldK1VSu8sdbvxV6ptsHCBQUqKj9Oc27I0qndHq5cChNXSpUv1xRdfaPTo0Ro6dKjVywFs65VXXtG99957wvfbt2+v91Zu1T1vrFFVfeQPSDvOR3vfkjhEfD7fSc/j6Nq1q79uBXDmEhISdMUVV+iKK67QCy+80KS12pgwUHWRPyckSTWeBhVu2W+rYQGczKWXXqpLL73U6mUAthcfHy9JSkpK0nnnnadzzz1Xbdu21R133KG/b9mv6igIUpI95yPNVBCeeeYZvfXWW7rpppv89wIoKyvTW2+9pZtvvln/9V//ZfEKAfub+L/LtGFPRdDP//QXY9S9fZLunrdaA89O1YPZff2P1XoatKu8Wn9duUt//myHfpbdVw+M66O9R2qU/dwn/qN7b08droszO2jhpn2aOm91s9Z7brdUvX//6Gb9DAAAJ2Oapg4cOKAOHTooLq7p5cWbOx8la2ek3eYjV08IwowZM/Tmm2/K4XCoqKhIRUVFcjgceuuttwhSQCvw+Uxt2195xq+z81CVXl9eosIt+9W7U4p+fdUg3TQsXbOXfKNv9lfq7LaJ+tn4YwPl34Z208WZHVRR69Gj721o9u/6el+lOFYFAAgFh8Ohzp07nxCkQjUfpdabkXabj2zzC8Ljjz+uX//61xowYIDVSwFi0s7yahmOMz9Zdeu3R/XEB5skSX+4PUvZA7poTL/OemtVqWb8Y53m332J7rikpz7euE8zJhz7//uzH23VtxW1zf5dhsOhnYeq1bNj8hmvGwCAkwnVfJRab0babT7STAVh2LBhVi8BiGklByvlDOGVf7qkupXZ6diHeHlVvSRp9c7v9LeiXXLFGZp350Xq1MatVTvK9cbnO1v0O5yGQyUHqwI/EQCAFgr1fJTCPyPtNh9ppoIwceJEq5cAxLRaj0+mznxLwI8GnaUdT39/k8Pvqus1d1mx/+v/r2CLsgd00VltE1TnbdCMf6xv8e8yJdV6o+OEYKCl1q9fr7PPPlsdO554MnldXZ3cbrcFqwJiR6jmo9R6M9Ju85FmKggPPvigvvnmmxO+X1lZqVdffdWCFQGxxdPgUyi2VzfuB5+95Bv919/X6dJnl6j4uKNjR+u8mr+mVJK0bNtBbT/Q8n3opmmqPoLvRA+Ewh133CGXy+X/uvEKt9XV1VzlD2gFoZqPUuvNSLvNR8JUEBYtWqTevXv7v37ppZckSSkpKZo9e7ZVywJihivOUCi2hDfuB//Nx1v19upSHa3znvAcn89s8r8t5XA4FO/kIxb25vP51LZtW//XF154oaRjl2+ur6+3allAzAjVfJRab0babT7a552EUeO1/Rs999xz/j//8KoqAEIvwWXIoei6W7pDUoKTzwfYm8/nU01NjSSpoqJC27dvV3V1tRoaGuT1nvgPMQChxXy0HudMBaF9+/Z6//33NXHiRL399tvyeDx68803lZKSojZt2li9PMD2enVMkbeZR8Ear27kabDm8qten6leNrlSEXAqN954o6644gqNHz9eS5Ys0X333adLL71UTqdTkyZNsnp5gO21ZD5K1s5Iu81HbtobhM2bN+uOO+7Q119/rWuuuUaPPvqo/uM//kM+n0/PP/+8+vfvb/USAVvz+UwNeLxAdUHuse7ePlH/N/1yOeMMjX/+k5Ddg6M53E5DW56YIEeo9l8AESovL09btmzRhAkTNHDgQK1atUqmaeqiiy6yemmA7TV3PkrWz0i7zUeaqSAMGDBARUVFTb734YcfWrQawJ4aGho0fPhw7dmzR0OGDNGQIUM0aNAgnXfeeRo8eLD6dE4J6g7v/zG2t+4anSlnnKH1u4/omzO4iMSZ6NslxTaDAjid3Nxc5ebm+r/mdiJA6ygrK1N+fr7iquIl94lX1DyZSJiRdpuPhKkgvPHGGxo1apR69uzZ5Pu1tbVav349gwMIAZ/Ppw0bNqi2tlZ79uxRXl6e/7G8vDwN69lDG/dUBLwAbHr7JNXWN2jp1wf0VN7mkF3lqDkcDmlYzw6t/4uBVlZcXKzdu3erV69e6t69u9XLAWzN4/Fo+fLlys/PV35+vtavP3Zp8vbjpio160rJEfhSCFbPSDvOR8JUEJ555hl99dVX/q+Liop00UUXKT4+Xnfffbe++OILC1cHRLe6ujp9+umnysvLO+kFXVwul3r37q0UtdP81aWqqj/9vSke+vu6cC01aImuOI3t39nqZQBhN3v2bD333HOaNWuWpk+fbvVyANtpbJ/y8/O1aNEiHT161P9YcnKyxo0bpwFjRuif5U5VewJv9bN6RtpxPhKmghAfH9/kH3lTpkzRhg0bZBhcDBFoiZ07d/qHw+LFi1VVdfI7oWdmZqqwsFAZGRk6x2cq2e0MGKYiQRu3UyMy7XXkDQAQfqdqnxoNHDhQOTk5ysnJ0ahRo+R2u+XzmSp8ZrGqPXUWrTp4dpyPhKkgGIahnTt3KiMjQ1u3btWuXbtUUlKixMREq5cGRIXj26f8/Hxt3ry5yeODBw9Wbm6uzjvvPE2ePFmSlJWVpby8PHXq1EmSZBgOTR2dqd8u3KraII6+WSXRZWjq6EwZhn32gwMAwieY9ik3N1cTJkxQRkbGCT/PfLQWYSoIjz/+uEaOHKl+/fqpsrJS8+bN04gRI+Tz+fTqq69avTwgIp2ufWrTpo2ys7OVk5OjCRMmKD093f/YsmXLVFFRoVdfffWEWw/ckJWuWR9vbbX30BI+U7o+Kz3wEwEAMakl7VMgzEfrEKaCcOWVV2rkyJHasWOHBg8erPj4eF155ZWSxFY/4F8a26f8/Hzl5eWdsn3KycnRiBEjTrgZdqNXXnnllL+jbaJLN2ala/7qUtU24zKwrSXBaeiGYelqm+iyeikAgAhSVlamgoIC5efna+HChSdtnxoD1Mnap0CYj9YhTAUpLS1NaWlp/q8JUUDL26cz8XDOABVs/Fa1RyNvb3jbJJcenjDA6mUAACwWjvYpEOajNQhTAIJ2fPuUn5+vTZs2NXk82PbpTCTGx2n2LUN12+srI2pveILL0Oybhyox/sQrEgIA7C/c7VMgzEdrEKYCqKqq0ptvvqmkpCTdcsstVi8HaHVWtE+BDOuZpklDu2vBmrKI2M6Q4DQ0aWi6snqmBX4yAMAWArVPAwYMUE5OjnJzc0PWPgXCfGx9hKkAKioqNHXqVJ111lmEKcSESGifgvHoxIHatLdCG/dUqM7CgeF2GhrYNVWPTrTn9gUAwPesbp+CwXxsXYQpABHZPgXidsZp3pSLNenVFSo+UGXJwHA7DWV2Sta8KRfL7bTn9gUAiGWR2D4FwnxsXYQpIAZFS/sUSLLbqQXTRui211dq056KVt3SkPCvI27zplysZDcfpQBgF9HQPgXCfGw99n+HACRJu3bt8l+2PFrap2Aku516c+pwzfxgsxasLW2Vk24TXMf2gD86cYDtj7gBgN1FY/sUDOZj6yBMATYVTPvUeGRt5MiREds+BcPtjNOT1wzW1UO66r6/rVVFjScsR+ESnIZSE12afctQDbPxybQAYHe7d+/2z8fTtU8TJkxQz549rVvoGWI+hh9hCrARu7ZPwRrWM02f/HyMni7YrLdXlcpwSDUhOBKX6DLkM6UbhqXr4QkDbHt5VwCwK4/HoxUrVigvL++07VNOTo5Gjx4dNe1TsJiP4UOYAqJYLLVPwUqMj9MTVw3W9PH99M7qUs1ZVqzKOq9qPA0yzeBfx+GQEl1xSnE7dffoTF2fZc87twOAXR3fPi1atEgVFRX+x+zUPgWL+RgehCkgysR6+xSstoku3TU6U1NG9tLy7Qe1ZOsBFZUc0rb9lTIcDjkNh0xJpmnK4XDIIcnrM+UzTfXtkqJhPTtobP/OGpHZQYbhsPrtAAACiPX2KVjMx9AiTAERrr6+XsuWLaN9aiHDcGh0n04a3aeTJMnnM7WrvFolB6tU621QvdeneKehBGecenVMVkaHJDkcDAcAiAana5+SkpKaXHkvFtqn5mA+hgZhCohAje1T432fKisr/Y/RPp0Zw3CoZ8dk9eyYbPVSAADNRPsUPszHliFMARGA9gkAgJOjfUIkI0wBFqF9AgDgRLRPiCaEKaCV0D4BAHBytE+IVoQpIIxonwAAsaKmpkaJiYn+q8CdTqD2qX///srNzaV9QsQjTAEhVF9f77/vU15eHu0TAMDWioqKtHDhQi1ZskRnnXWW3njjjVMGKdon2BFhCjhDtE8AgFhRWVmplJQU/9fDhw/XXXfdpWeeeUZZWVkn/Zn9+/dr/PjxWrduXZPv0z7BDghTwHGWL1+uX/3qV/rpT3+qa6+9Vj6fT4ZhNHkO7RMAIJYcPHhQc+bM0ccff6yjR49q6tSp+vGPf6z09HSNHTtWV111lbKysk65vS81NVV79uyhfYItEaYAyR+avv76a+3evVuffvqprr32Wv9QME1Ta9as0ZNPPkn7BACICY3h6Msvv9R3332n559/XgkJCZo1a5Z2796tmTNn6rrrrtPkyZN144036uDBg5o0aZJuvvnmE16rsLBQffv2pX2C7RCmAEmGYcjj8eijjz7Sr371K73xxhtNtjI4HA4lJSXpvffek0T7BACwl+LiYv3ud79TXFycbrnlliZN0/nnn68LLrhAHTp0kCSNHj1an3/+uSTp5ptv1vvvv6/s7Gz17dtXN9xwgxITE3XllVcqLi5OkpSQkKBzzz3XsvcGhBNhCtCxo2+HDx9W27Zt1atXLyUnJ2v9+vW65JJL/MOkT58++vOf/6wxY8bQPgEAbKOmpkbPP/+8DMNQr1699LOf/Uyvvfaa+vbtK0nq1KmTpO93cRQUFGjy5MmSpPbt2ys/P9//WlOnTtXKlSuVk5PjD1OAnRGmELOO39vtcDj01ltvacSIERo2bJj69u2r+fPnq02bNho8eLAkyeVy6fbbb7dyyQAAhFxlZaXef/997dixQ9KxcPXaa6/pscceU3JysqRjM9MwDG3btk0HDx5U//79T/paXq9X1dXVcrvdQV0iHYh2hCnElMYr73311VeaNWuWkpKSJB0bEj6fTwsXLtTKlSv1z3/+U9XV1erTp48/TAEAYEcNDQ0aMmSItm3bpj59+mjkyJGaN2+e1q1b12SHhiTNnTtXw4YNU+/evf0/u337dv3hD3/Q2rVrVVVVpTlz5kgSQQoxwQj8FCB61dfXq7CwUA899JAGDx6sjIwM3XPPPXr55Zeb3N/C4XDo8OHDOnr0qHr16qWXXnpJo0aN0k033WTh6gEACD+n06mMjAz/jXO7d++ujh07auPGjZKOHXCUpJ07d6qoqEj333+/3njjDT3yyCOKi4uTx+NRSkqKnnzySX322WecH4WYQjMF2zndfZ9SUlKUnZ2t3NzcJvfJkKTHHnvMfxTt6NGjWrRokdauXavs7OxWXT8AAOHg8Xi0YsUK7dmzR1dddZV/C1+7du3Uo0cPFRUV6brrrlP79u1VW1vrv+BE4y1C/vKXv2jp0qWaOHGiBg0apNzcXDU0NGjQoEEaNGiQZe8LsBJhClHv+Ps+5efn+4+kNQr2ynsOh0MNDQ3+K/f9/ve/b43lAwAQNrt371ZBQYHy8vK0aNEiVVRUKC0tTdddd53/OU6nU5dddpnuueceHTlyRG3bttXnn3+um2++WXl5edq+fbv+/d//XcOHD1dRUdEpb84LxCLCFKJSaWmpPzwtWrTolO1Tc+/7xJWHAADRrLF9apyR69ata/J4//79lZOTI4/H0+SeT1lZWRo+fLimTZumTZs2qX///urfv7+6deum0aNHS5LGjx/fqu8FiAaEqVPwer3Kzc3VmjVrJEn79u1Tp06ddOutt+qFF16wdnExKFTtEwAAdnOy9qlRUlKSxo4d6z/A2KtXr1O+zvPPP6+VK1fK6/Xq8ssvlyT/VkAAJ0eYOo2VK1f6P5BM09TBgwf94QrhF0z71BiguO8TACBWBNs+5eTkaPTo0UpISAjqdV0ul0aNGhWOJQO2RZg6BafTqQceeEAzZ85s8v3p06dbtCL7o30CAODkGtun/Px8LVy4sMXtE4DQIkydxoMPPqgXXnhBR48elSQNGTJEV199tcWrshfaJwAAThSu9glAaBGmTiMtLU133nmn/xypxx9/nBvQnaH6+notX75ceXl5J22fGi+1SvsEAIg1tE9A9CFMBTBjxgy9+OKLateuHa1UCwXbPk2YMEE9evSwcKUAALQej8ejzz77zH+AkfYJiD6EqQC6dOmiqqoqxcXF0UoFqbF9ys/PV15e3knbp5ycHOXm5tI+AQBiCu0TYC+EqVPw+UztLK9WycFK1Xp88jT45IozlOAy1KtjijLSkmQYhKtGtE8AEBuOn4+VaX2UNPAyHUxMV+GWfczHk2hsnxoPMNI+AfbiME3TtHoRkcDnM7V8+0EVbtmvVTvKtW1/pQyHQ07DIVOmTFNyOCSHHPL6TPlMU306p2hYzzSN7d9ZI8/pGFPD4/j2KT8/Xxs2bGjyOO0TANhDoPno9TbI6YxjPh4nmPapMUDRPgHRLebD1JEaj+avLtXcZcWqqvOqur5BzfkLcUhKio9TstupqaMzdUNWutomusK1XEvRPgFA7GA+Bo/2CYhdMRumauob9HT+Zr29ulQOh1Tr8Z3xaya6DPlM6casdD2cM0CJ8XEhWKl1gm2fcnJyNGrUKNonALAB5mNw9uzZ45+PtE9A7IrJMFVUUq7731yrihqPar1nPiR+KMFpKDXRpdm3DNWwnmkhf/1won0CgNjFfDy149un/Px8ffXVV00ep30CYlNMhak6b4NmfrBJC9aWheRIWyAJLkOThnbXoxMHyu2MzKNwtE8AAObjydE+AQgkZsJUVZ1Xk19bqc17K8JytO1UEpyGBnZN1bwpFyvZHRkXTywrK/Pv6168eLGOHj3qfywlJUXjxo3zX5aV9gkA7I35+D3aJwDNFRNhqqrOq0mvrlDxgSrVteKgaOR2GsrslKwF00ZYMjBonwAAJxPr81E61j4VFBQoLy+P9glAs9k+TNV5G3TTnM+1aW+FJYOikdtpaFDXVL05dXirbGkIpn1qHA60TwAQe2J1PtI+AQgl24ep/353vRasKWvVrQunkuA0NOnCdD15zeCQvzbtEwCgOWJlPkrft0+N5z4dOXLE/xjtE4AzERmblMOkqKT82Mm0ETAoJKnW69OCtaW6ekjXkFzFqLF9arzyHu0TACAYdp+PtE8AWottm6ma+gZdNmuJ9h+ts3opJ+jcxq1Pfj6m2ffZCNQ+DRw4ULm5ubRPAIBTiub5WFtbe8rgQ/sEwAq2baaeyt+sihqP1cs4qYoaj54u2Kwnrgq8nYH2CQAQStE4Hz0ej+6++2795S9/UWFhoS677LKA7VO/fv2Uk5Oj3Nxc2icAYWPLZupIjUcXPbXI0hNqA3E7DRU9kq22ia4m3/d4PFq+fLny8vJonwAAIRWN87G6ulrXX3+98vLyJEmjR49Wly5daJ8ARARbNlPzV5fK4bB6FadnOKR3VpfqrtGZtE8AgFYRbfNx//79GjduXJMDi8uWLfP/mfYJgNVs10z5fKaGP7M4IveC/1CCWSfnh7/Sxg3rm3yf9gkAEGrRNB87t3Hr8xnj1LXr2dq3b98Jjz/yyCO66667aJ8AWM52zdTy7QdVVee1ehlBqfb4dOCooeTkZGVnZysnJ0cTJkxQRkaG1UsDANhMNM3HyjqvVhQfUrdu3U4aprp27UqQAhARbBemCrfsV3V9g9XLCIoRn6Bbf/GMZt85Rm632+rlAABsLJrmY42nQYVb9mvNmjWqra3Vxo0b9eWXX+rLL79UWVmZsrOzrV4iAEiyYZhataNcgfYtfvqLMerePkl3z1utjzcdO+L12MSBmjKylxasKdXPF6zTiHM66D+z+6rfWW1kOBw6UFmnDbuP6P43v5AkTRraXbOuP1+SVO/1qbreqx2HqvXxpm/12qclQZ7c69A+XzJBCgAQdtE0H01TWrXjkCQpISFBF154oS688MKWv3kACBNbhSmfz9S2/ZVn/DpdUt36w+1ZcsUZKtjwrY7WeZXZMVnZA7qc8NxDlXX657o9Sm+fpEv7dNKQ9HbKHtBFN8/9PKiB8fW+SpmmKUeknxEMAIhazEcACA9bhamd5dUyQvChOyS9vZLinVq46Vv9+1tf+L//w8uYS9K+ilr9+v1NkqRBXVP1j3tHaGiP9poyspde/mR7wN9lOBzaeahaPTsmn/G6AQA4GeYjAISHYfUCQqnkYKWcxpkPiwNHayVJY/t30bwpF+ln2X00rGd7VdSe/iaHG/dU6OONx7ZFjOvfOajf5TQcKjlYdWYLBgDgNJiPABAetgpTtR6fzIA7wgNbu+uwXl26XaZpanSfTnpgXF+9M22E/nnfKKUmnL7M2324RpLUISW486BMSbXe6DghGAAQnZiPABAetgpTngafgrlrVq3n2F5tV9z3bz/+X3+u+ddjT+dvUdZTi3Tf39bqbyt3qt7r07nd2ur6rPTTvna3domSju0VD4ZpmqqP4DvRAwCiH/MRAMLDVmHKFWcEdWf3neXHtg0M7dFO0rG7rQ9JP/bnHYeq1LVtgnqkJelwtUcfrt+rR97doP/bul+SlOI+9ZG3QV1T9aNBx07CXbxlf1Brdjgcinfa6j8DACDCMB8BIDxsdQGKBJchhwJPi7eKSjWufxfdOSpTF2akqW2iS706Jqui1qN/frVH53Zrq7m3ZemLXd9p+4EqJcYbGtO/s3w+U8u/OdjktbqkJujxKweqe7skXda3k+Kdhtbu+k6vLy8Jas0OSQnOuJa8XQAAgsJ8BIDwsFWY6tUxRV5f4H0MCzfv0wNvf6G7RmXqnE7JqvP6tHTbAf3mo606cLROX+87qn+sLVNWzzQNODtVpqSt3x7VnGXFWr3zuyav1SHFrckXZ6iqzqtNe4/oo4379NqnJapvCG5rgtdnqhdXKgIAhBHzEQDCw2Gaweyijg4+n6kBjxcEecPcyOB2GtryxATuowEACBvmIwCEh602IxuGQ306p1i9jGbp2yWFQQEACCvmIwCEh63ClCQN65kWxK7wyOBwSMN6drB6GQCAGMB8BIDQs12YGtu/s5Lio+OE1URXnMYGefNCAADOBPMRAELPdmFq5DkdlXyay7NGkjZup0ZkcuQNABB+zEcACD3bhSnDcGjq6EwluCL7rSW6DE0dnSnDiJZNFwCAaMZ8BIDQi+xP1Ba6ISs9qDu9W8lnKuDd4gEACCXmIwCEli3DVNtEl27MSldChN45PcFp6MZh6Wqb6LJ6KQCAGMJ8BIDQisxP0xB4OGeAUiP0w7htkksPTxhg9TIAADGI+QgAoWPbMJUYH6fZtwyNuL3hCS5Ds28eqsQouaISAMBemI8AEDqR9UkaYsN6pmnS0O4Rs50hwWlo0tB0ZfVMs3opAIAYxnwEgNCIjE/RMHp04kAN7Joqt8UDw+00NLBrqh6dyPYFAID1mI8AcOZsH6bczjjNm3KxMjslWzYw3E5DmZ2SNW/KxXI72b4AALAe8xEAzpztw5QkJbudWjBthAZ1TW31LQ0JTkODuqZqwbQRUXOzRABAbGA+AsCZcZhmpN9xInTqvA2a+cFmLVhbqlqPL+y/L8F1bA/4oxMHcMQNABCxmI8A0DIxFaYardpRrvv+tlYVNR7VekM/NBKchlITXZp9y1AN42RaAECUYD4CQPPEZJiSpJr6Bj1dsFlvryqV4ZBqQnAkLtFlyGdKNw5L18MTBnB5VwBA1GE+AkDwYjZMNTpS49E7q0s1Z1mxKuu8qvE0qDl/Iw6HlOiKU4rbqbtHZ+r6LO7cDgCIfsxHAAgs5sNUI5/P1PLtB7Vk6wEVlRzStv2VMhwOOQ2HTEmmacrhcMghyesz5TNN9e2SomE9O2hs/84akdlBhuGw+m0AABBSzEcAODXC1Cn4fKZ2lVer5GCVar0Nqvf6FO80lOCMU6+OycrokCSHg+EAAIgtzEcA+B5hCgAAAABaICbuMwUAAAAAoUaYAgAAAIAWIEwBAAAAQAsQpgAAAACgBQhTAAAAANAChCkAAAAAaAHCFAAAAAC0AGEKAAAAAFqAMAUAAAAALUCYAgAAAIAWIEwBAAAAQAsQpgAAAACgBQhTAAAAANAChCkAAAAAaAHCFAAAAAC0AGEKAAAAAFqAMAUAAAAALUCYAgAAAIAWIEwBAAAAQAsQpgAAAACgBf4fIe67fW4rSBoAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 1080x288 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# n 代表图形 G 的顶点数量————指不同国家的货币\n",
    "n = 3\n",
    "nodes = [ \"JPY\", \"CNY\", \"USD\"]\n",
    "G = nx.DiGraph()\n",
    "G.add_nodes_from(nodes)\n",
    "edges = [(\"JPY\",\"CNY\", 0.5), (\"CNY\",\"JPY\",2), (\"CNY\",\"USD\", 0.33), (\"USD\",\"CNY\",3),(\"JPY\",\"USD\", 0.25), (\"USD\",\"JPY\",4)]\n",
    "G.add_weighted_edges_from(edges)\n",
    "\n",
    "# 两幅图分别表示不同方向的转换率\n",
    "G1 = nx.DiGraph()\n",
    "G1.add_nodes_from(nodes)\n",
    "edges1 = [(\"JPY\",\"CNY\", 0.5),  (\"CNY\",\"USD\", 0.33), (\"USD\",\"JPY\",4)]\n",
    "G1.add_weighted_edges_from(edges1)\n",
    "\n",
    "G2 = nx.DiGraph()\n",
    "G2.add_nodes_from(nodes)\n",
    "edges2 = [(\"CNY\",\"JPY\",2), (\"USD\",\"CNY\",3),(\"JPY\",\"USD\", 0.25)]\n",
    "G2.add_weighted_edges_from(edges2)\n",
    "\n",
    "options = {\n",
    "    \"with_labels\": True,\n",
    "    \"font_weight\": \"bold\",\n",
    "    \"font_color\": \"white\",\n",
    "    \"node_size\": 2000,\n",
    "    \"width\": 2\n",
    "}\n",
    "fig, ax = plt.subplots(1, 2, figsize=(15, 4))\n",
    "for i, a in enumerate(ax):\n",
    "    a.axis('off')\n",
    "    a.margins(0.20)\n",
    "nx.draw(G1, pos=nx.circular_layout(G1), ax=ax[0], **options)\n",
    "nx.drawing.nx_pylab.draw_networkx_edge_labels(G1, pos=nx.circular_layout(G1), ax=ax[0], edge_labels=nx.get_edge_attributes(G1, 'weight'))\n",
    "nx.draw(G2, pos=nx.circular_layout(G2), ax=ax[1], **options)\n",
    "nx.drawing.nx_pylab.draw_networkx_edge_labels(G2, pos=nx.circular_layout(G2), ax=ax[1], edge_labels=nx.get_edge_attributes(G2, 'weight'))\n",
    "plt.axis(\"off\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 编码哈密顿量\n",
    "\n",
    "这里我们将式（4）中的二进制变量用式（5）替换，从而构建哈密顿量 $H_C$。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:00:16.237497Z",
     "start_time": "2021-05-17T08:00:16.219567Z"
    }
   },
   "outputs": [],
   "source": [
    "# 惩罚参数 等于构造的有向图中权值最大的边\n",
    "penalty = 4 \n",
    "# 在本样例中，设定交易回路内货币数量为全部的货币种类\n",
    "K = n \n",
    "# 构建最佳套利问题的哈密顿量\n",
    "hamiltonian = arbitrage_opportunities_hamiltonian(G, penalty, n, K)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 计算损失函数\n",
    "\n",
    "我们使用 $U_3(\\vec{\\theta})$ 和 $\\text{CNOT}$ 门构造的参数化量子电路，通过调用量桨内部的 [`complex_entangled_layer()`](https://qml.baidu.com/api/paddle_quantum.ansatz.circuit.html#Circuit.complex_entangled_layer) 来实现。\n",
    "\n",
    "上述电路会给出一个输出态 $|\\vec{\\theta}\\rangle$，由此输出态，我们可以定义最佳套利机会问题在经典-量子混合模型下的损失函数：\n",
    "\n",
    "$$\n",
    "L(\\vec{\\theta}) =  \\langle\\vec{\\theta}|H_C|\\vec{\\theta}\\rangle.\n",
    "\\tag{6}\n",
    "$$\n",
    "\n",
    "之后我们利用经典的优化算法寻找最优参数 $\\vec{\\theta}^*$。下面的代码给出了通过量桨和飞桨搭建的完整网络："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:00:16.258893Z",
     "start_time": "2021-05-17T08:00:16.241066Z"
    }
   },
   "outputs": [],
   "source": [
    "class AONet(paddle.nn.Layer):\n",
    "\n",
    "    def __init__(self, num_qubits, p, dtype=\"float64\"):\n",
    "        super(AONet, self).__init__()\n",
    "\n",
    "        self.depth = p\n",
    "        self.num_qubits = num_qubits\n",
    "        self.cir = Circuit(self.num_qubits)\n",
    "        self.cir.complex_entangled_layer(depth=self.depth)\n",
    "        \n",
    "\n",
    "    def forward(self):\n",
    "        \"\"\"\n",
    "        前向传播\n",
    "        \"\"\"\n",
    "        # 运行电路\n",
    "        state = self.cir(init_state)\n",
    "        # 计算损失函数\n",
    "        loss = loss_func(state)\n",
    "        \n",
    "        return loss, self.cir"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 训练量子神经网络\n",
    "\n",
    "定义好了量子神经网络后，我们使用梯度下降的方法来优化其中的参数，使得式（6）的期望值最小。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "SEED = 100   # 随机数种子\n",
    "p = 1        # 量子电路的层数\n",
    "ITR = 120    # 迭代次数\n",
    "LR = 0.4     # 梯度下降优化速率 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "使用飞桨，优化上面定义的网络。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "循环数： 10     损失: 0.4704\n",
      "循环数： 20     损失: 0.1302\n",
      "循环数： 30     损失: -0.2744\n",
      "循环数： 40     损失: -0.4700\n",
      "循环数： 50     损失: -0.5512\n",
      "循环数： 60     损失: -0.5684\n",
      "循环数： 70     损失: -0.5821\n",
      "循环数： 80     损失: -0.5833\n",
      "循环数： 90     损失: -0.5843\n",
      "循环数： 100     损失: -0.5847\n",
      "循环数： 110     损失: -0.5849\n",
      "循环数： 120     损失: -0.5849\n"
     ]
    }
   ],
   "source": [
    "# 固定 paddle 随机种子\n",
    "paddle.seed(SEED)\n",
    "# 定义电路所需比特数\n",
    "num_qubits = n * K\n",
    "# 定义量子神经网络\n",
    "net = AONet(num_qubits, p)\n",
    "# 定义初始态\n",
    "init_state = paddle_quantum.state.zero_state(num_qubits)\n",
    "# 使用 Adam 优化器\n",
    "opt = paddle.optimizer.Adam(learning_rate=LR, parameters=net.parameters())\n",
    "# 定义损失函数\n",
    "loss_func = paddle_quantum.loss.ExpecVal(hamiltonian)\n",
    "# 梯度下降优化循环\n",
    "for itr in range(1, ITR + 1):\n",
    "    # 运行上面定义的网络\n",
    "    loss, cir = net()\n",
    "    #计算梯度并优化\n",
    "    loss.backward()\n",
    "    opt.minimize(loss)\n",
    "    opt.clear_grad()\n",
    "    if itr % 10 == 0:\n",
    "        print(\"循环数：\", itr, \"    损失:\", \"%.4f\"% loss.numpy())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 解码量子答案\n",
    "\n",
    "当求得损失函数的最小值以及对应的一组参数 $\\vec{\\mathbf{\\theta}}^*$后，我们的任务还没有完成。为了进一步求得最佳套利机会问题的近似解，需要从电路输出的量子态 $|\\vec{\\mathbf{\\theta}}^*\\rangle$ 中解码出最佳套利问题的答案。物理上，解码量子态需要对量子态进行测量，然后统计测量结果的概率分布（我们的测量结果是表示最佳套利机会问题答案的比特串）：\n",
    "\n",
    "$$\n",
    "p(z) = |\\langle z|\\vec{\\mathbf{\\theta}}^*\\rangle|^2.\n",
    "\\tag{6}\n",
    "$$\n",
    "\n",
    "在量子参数化电路表达能力足够的情况下，某个比特串出现的概率越大，意味着其是最佳套利机会问题最优解的可能性越大。\n",
    "\n",
    "量桨提供了查看参数化量子电路输出状态的测量结果概率分布的函数："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "利用哈密顿量找到的解的比特串形式： 100001010\n"
     ]
    }
   ],
   "source": [
    "# 模拟重复测量电路输出态 1024 次\n",
    "final_state = cir(init_state)\n",
    "prob_measure = final_state.measure(shots=1024)\n",
    "arbitrage_opportunity_route = max(prob_measure, key=prob_measure.get)\n",
    "print(\"利用哈密顿量找到的解的比特串形式：\", arbitrage_opportunity_route)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "概率最大的测量结果是表示最佳套利机会问题答案的比特串。返回的二进制字符串每 $n$ 位一组，每组中的第 $k$ 位为 $1$ 表示该笔资产交易到该市场时的顺序为 $k$ 。如果最后的结果不是类似的有效解，读者可以通过调整参数化量子电路的参数值，即随机种子 `SEED`、层数 `p`、迭代次数 `ITR` 和梯度下降优化速率 `LR`，来获得更好的训练效果。\n",
    "\n",
    "我们将找到的比特串映射回经典解，即转化成了 ``dictionary`` 的形式。其中 ``key`` 代表货币名称，``value`` 代表顶点在套利回路中的顺序，即一笔资产交易到该市场的顺序。在下面的代码中，我们还计算了该笔资产经找到套利回路交易一圈后的正向回报率。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:02:14.554317Z",
     "start_time": "2021-05-17T08:02:14.500593Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'JPY': 0, 'CNY': 2, 'USD': 1}\n",
      "正向回报率:  1.5\n"
     ]
    }
   ],
   "source": [
    "solution = {nodes[i]:t for i in range(n) for t in range(n) if arbitrage_opportunity_route[i * n + t] == '1'}\n",
    "print(solution)\n",
    "rate = sum([np.log2(G[u][v][\"weight\"]) if solution[v] == (solution[u] + 1) % n else 0 for (u, v) in G.edges])\n",
    "print(\"正向回报率: \", 2**rate)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "为了清晰地表达资产的交易路径，依然选择以图的形式呈现：\n",
    "* 顶点中的数字代表交易到该货币市场的顺序，字母代表不同的货币类型\n",
    "* 红色的有向边表示找到的最佳路线"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-17T08:02:14.571954Z",
     "start_time": "2021-05-17T08:02:14.559634Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1MAAADnCAYAAAD7CwxiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAA4PUlEQVR4nO3deXiU9b3//+d9ZyaTjQAh7IRA2BeVYgLI4sJSScT14FrUywWX6jnaY63oOWgrHvXXWrUesRWqbYW6Yn+1aoKCcARRCIsLu0hCSNhDlJB9lvv7x5BADGGGMMlsr8d1cZXMfefOZ+LVefO635/78zEsy7IQERERERGR02IGewAiIiIiIiLhSGFKRERERESkBRSmREREREREWkBhSkREREREpAUUpkRERERERFpAYUpERERERKQFFKZERERERERaQGFKRERERESkBRSmREREREREWkBhSkREREREpAUUpkRERERERFpAYUpERERERKQFFKZERERERERaQGFKRERERESkBRSmREREREREWkBhSkREREREpAUUpkRERERERFpAYUpERERERKQFFKZERERERERaQGFKRERERESkBWzBHkCo8ngsisqqKCytoMbpwen2YI8xibOb9E1NIj0lAdM0gj1MERGRNqX6KCJynMLUMR6PxaqdpSzbdpC1u8rYcbAC0zCwmQYWFpYFhgEGBi6PhceyGNAliaw+KUwc3IVx/VJVPEREJOKoPoqINM+wLMsK9iCC6Ui1k7fXFTN/ZQGVtS6q6tyczi/EABJiY0h02Jg5IYNrMtNoH29vreGKiIi0CdVHERHfojZMVde5eSpvK2+tK8YwoMbpOeNrxttNPBZcm5nGw9lDiI+NCcBIRURE2o7qo4iI/6IyTOUXlnHvGxsor3ZS4zrzIvFjcTaT5Hg7c28YSVaflIBfX0REpDWoPoqInJ6oClO1LjdzPtjCog0lAbnT5kuc3WT6yF7MnjYUh0134UREJDSpPoqItEzUhKnKWhczXlnD1n3lrXK3rTlxNpOhPZJZcOtoEh1a70NEREKL6qOISMtFRZiqrHUx/eXPKThUSW0bFop6DptJRudEFt05VgVDRERChuqjiMiZifhNe2tdbma8siZohcI7Bg8Fhyq58dU11LrcQRmDiIjIiVQfRUTOXMSHqTkfbGHrvvKgFYp6tS4PW/aWM+eDrUEdh4iICKg+iogEQkSHqfzCMu/DtEEuFPVqXB4WbShm7a6yYA9FRESimOqjiEhgRGyYqq5zc+8bG9pkVaLTUeP0cM/rG6iu03QGERFpe6qPIiKBE7Fh6sm8rZRXO4M9jJMqr3by1GJNZxARkban+igiEjgRGaaOVDt5e11xyExf+LEal4e31hZzJESLmYiIRCbVRxGRwIrIMPX2umIMI9ijODXTgHfWFQd7GCIiEkVUH0VEAiviwpTHYzF/ZUHIzQX/sWqnh3krC/B4In6bLxERCQGqjyIigRdxYWrVzlIqa13BHoZfKmpdfF5wONjDEBGRKKD6KCISeBEXppZtO0hVmKwEVO10s2zbwWAPQ0REooDqo4hI4EVcmFq7q4xwmRhgWbB2l+68iYhI61N9FBEJPFuwBxBIHo/FjoMVDV+fm96Ruy/ox7npHUmwx7C/vIbl2w/xP7lbcLotdj11CQCb9x7hkv/9DICfDu3KvBszKfm+ivG/Xc5bM8cwOqMTS7YcYOaCdQB0TLCz7IEL6ZgQy5wPt/DKZ4XNjumyc3owY3Q6/TonkhBrY3dZFa+sKuTtYw/XfnugAsuyMEL9iWAREQlbqo8iIq0jojpTRWVVmMc+dC89uztvzRzD5CFd2fdDNf/4cg+7y6r42ejexNtjGn3fsB7tuXhYt5Ne8+H/fyO1TjdThnZl8pAuADySPYSOCbF8VfwDf1nVfKEAOH9AKr06xrNiRylri8oY1K0dv/23s5k02Hst0zAoOlx1pm9dRESkWaqPIiKtI6I6U4WlFdhMgzi7yeOXD8cWY/KPL0t44J2vsY7NbeidkkC1s/GccY/H4v5JA/ho8/4m1yworeSFZd/x4MWDeGzaMGqcHv5tZC+cbg+z/vENvhYbenXVLn717vHz3pw5hjEZnRg/IJVPth3EZhoUllbSJzUxEL8CERGRJlQfRURaR0R1pmqcHiwsMtNT6JgQC8CLy75rKBQAu8uqcLobf8LnbtrHkO7J5Aw/+d23l1fsZOu+ctJSEnjl5kxM02DeigK27T/qc0xb9pU3Kij2GO+vfP+RGgAsoMYVHg8Ei4hIeFJ9FBFpHREVppxuD5YFnZJiG17b80O1z+/76+e7OFxRy32TBp50brbLYzHrH9/g9lg4bDEUlFbwwrIdpz2+28b35dz0jhSWVvL3NbsBsCyLuhDdiV5ERCKD6qOISOuIqDBljzExDDhcUdfwWs8O8T6/r6rOzbyVBQzq1o5pZ3U/6TlflxxhfdH3ALy+Zje1p/kBf/+kAcy+ZChFhyv52Z9XU3Fsrw/DMIi1RdR/BhERCTGqjyIirSOiPqXi7CYGBuuLvueHKm/BuHdif068mdazQzw2s+ndtb99sYtDR2u5pJliAeA5Nh/CY/m/uKxhwJzLh3P/5IFs2nOE6X/6gr3HpjAAGECcLab5C4iIiJwh1UcRkdYRUQtQ9E1NwuWxqHa6eez9zTx79Qiu+kkvBndtx1clR+ia7GB8/1Sy/mcp5TWNd4GvcXp4ecVO/vuSoX7/vF4d4vnsoYkAnP2bj5pcE+CBKYO4cUw6bo/F5r3l3H1hPwB2Ha7ktS+KcHks+urhWhERaUWqjyIirSOiwlR6SkLDXbH3vtrLvh9quOuCfozs3YF+XZLYf6SGN9cWN1mtqN6C1UXcMSGDLslx/v3AYzfwPB4LdzPLFnVLdgAQYxpcm5XW8PrqgsO89kURHssivVOCn+9QRETk9Kk+ioi0DsOyTqMnHwam/e9KNu0tb5OfNWVIV+bflMlrX+zi0X9tbtE1zuqZzPv3TsDj8bBz506++uorvvnmGzIzM7n88ssDPGIREYlWYVcfk03ef2ASxMb6PllEJEgiqjMFkNUnhc17y2mLhDi2fyd2Ha7kqbxtLfp+w4CesTV06dKFH374AafT2XAsLS1NYUpERAImrOqjx0PW0n/Ck9fApEmQne3907t3YAcqInKGImoBCoCJg7uQENs2D6z+5v0tXPjM/zU7LcKXeHsMtYXrOXToUKMgBZCamsqePXsCMUwREZHwqo8xMNFTChUV8N57cNddkJ4Ow4bBgw/CsmVQV+f7QiIirSziwtS4fqkkOsKj4dbOYWP+nAe55pprmhz78ssv6dWrF+eccw4PP/wwK1asaBK4RERE/BVW9TEpnrGf/AN274Z58+DKK6FdO9iyBZ55xtut6tQJrrgCXn7Ze56ISBBE3DNTAPNXFvD7JdupcYbuZn/xdpMHpgzi9gkZWJbFb3/7W2bNmtVwfPLkyXz++edUVVU1vJacnMyUKVPIzs5m6tSp9OzZMxhDFxGRMBVu9bGRujpYtQry8rx/Nm1qfHzoUMjJ8U4HHD9ez1qJSJuIyDB1pNrJqCeXnvbGgW3JYTPJf2Qy7ePtDa+98sor3HHHHYwaNYovvviC2tpaVqxYQV5eHnl5eWzb1nju+dlnn01OTg7Z2dmcd9552O32H/8YERGRBuFaH0+quBgWL/YGq6VL4ejR48eSkvSslYi0iYgMUwCPvreJt9cVUxOCBSPOZnJNVhqPXza8ybGCggJSUlLo0KFDk2O7du1qCFaffPKJulYiInLawrU+npK6ViISJBEbpqrr3FzwzHIOHq0N9lCa6Jrs4P8euIj4M3gQWF0rERFpiUivj0DjrtWSJd6FLOqpayUiARSxYQpg7a4ybnx1TUjNDY+zmyy8dTSZfVICel11rURExF/RVB/VtRKR1hTRYQrgv/+5kUXrS0JiOkOczWT6uWk8ccVpTl84TTU1NaxcuVJdKxERaVY01kfgeNcqN9f7rJW6ViJyBiI+TNW63Fw/fzWb95YH9YFbh81kWI9k3pg5Boetbfb5qFdYWMjixYvJzc1l2bJl6lqJiIjqI6hrJSJnLOLDFEBlrYvpL39OwaHKoBQMh80ko3Mii+4cG/Q9PtS1EhGReqqPP6KulYicpqgIU+AtGDe+uoYte8vbdEpDnM1kaI9kFtw6OjQKxY+oayUiEt1UH5vhq2s1bNjxYKWulUjUipowBd4pDXM+2MqiDcVt8tBtnN1k+sg0Zk8b0vZTF1pAXSsRkeik+uiH4uLjwUpdKxE5JqrCVL21u8q45/UNlFc7W+UuXJzNJDneztwbRpIV6FWJ2pC6ViIi0UX10U/qWonIMVEZpsC7z8ZTi7fy1tpiTAOqA3AnLt5u4rHg2qw0Hp465Mz3yQgh6lqJiEQH1ccWUNdKJGpFbZiqd6TayTvripm3soCKWhfVTjen8xsxDIi3x5DksHHHhAyuzkyjfXzkB4nCwsKGYKWulYhI5FF9bKH6rlVurjdcbd7c+Li6ViIRJerDVD2Px2LVzlKWbz9EfuFhdhyswDQMYkw4evQophlDYmIiBuDyWHgsi4Fdk8jq04mJg7swNqMTpmkE+20EhbpWIiKRq7n6aDMNLMCyLAzDUH1szu7d3hUC1bUSiUgKU83weCx2l1WxfkcxN91yGx1SUvnzvD8RZ4uhb2oi6Z0SMIwoLg6noK6ViEjkqq+PhaWV1Ljc1Lk8xNpM1Ud/qGslEnEUpnzYt28fPXr0oFu3buzbty/Ywwk76lqJiESw//kfmDMH5s2Dm24K9mjCjz9dq/pNg9PSgjdOEWmWwpQPClOB5W/XKjs7mx49egRxpCIickpOJ6SmQnk5pKdDYaH3QSlpGXWtRMKSwpQPClOtp75rlZubS15eHtu3b290/JxzzmkIVupaiYiEmFdfhdtuO/714sVw8cXBG0+k8dW1mjz5eLhS10okaBSmfFCYajvqWomIhAmnEwYN8naj6o0ZA59/ru5Ua6irg88+O778urpWIiFDYcoHhangUNdKRCSELVwIN97Y9PVly+Cii9p+PNFGXSuRkKEw5YPCVGho7a5VTU0NcXFxgRyyiEjk+vvfvWHqxH9C2Gzef+BPmhS8cUUjda1EgkphygeFqdAT6K7Vn//8Z/7+979z8OBBnn76aaZNm6ZlfUVEfCkvh1mz4I9/hN/8Bv7zP71dEQmu3buPB6tPPlHXSqSVKUz5oDAV+vzpWt1yyy1kZ2djmmaj783Pz+e6665j9erVbN68md/97ne8++67xMfHN5zjdruJiYlps/cjIhI2HngAnn0WnnnG+3cJLa3Qtdq0aRM9e/akY8eOrTRokfBi+j5FJLT17duXn//857z//vscPnyYjz76iPvvv59BgwZRXl7Ou+++y/bt2/F4PI2+r6qqig8//JBbbrmFLl260K9fPwDKysoaznG73cybN4+xY8cyZ84cNmzY0KbvTUREpMViY2HiRPjd72DTJigqgj/9CS6/3Nul2rzZG4QnTYJOneDKK717hpWUnPRyP/zwAy+99BLjx49nzJgxfPTRR238hkRCj8KURJS4uDh++tOf8txzz7Ft2zYKCgqYO3cuP/vZz7DZbI3OLSkpYf/+/Uw6Nr+/rKyMIUOG8O233zacc/DgQfbt28fQoUNxuVw8+uijfPDBB236nkRERAKid2+480745z/h8GHvNMBf/tLboaqo8L5+553wi180nh54jMPh4KmnnmLz5s08/vjjvPTSS41mg4hEI5vvU0TCV33X6mT27dtHZWUl55xzDgAHDhygvLycPn36NJxTUlLC4cOHmT17Nunp6QC4XK5WH7eIiEirqu9a1XeuTnzW6rrrTvr8W3x8PA6HA4BevXpx+PBhampqSEhIaOvRi4QMhSmJWuXl5ezdu5fExETKy8vZsWMHKSkp9O3bt+GcsrIyPvzwQwoLC7n77rvJzs5u0uESEREJe/VdqzvvbLxK4wnq6uqIPfZc1fz58xk5cqRWwpWop2l+ErXatWtHamoqAB9//DGffvopV1xxBUDD81UXXngh69at47rrruO1115j9+7dwRquiIhI22hmRdvY2Fj27t3LtGnTOHz4MPfee+/xrtSBA1Bc3IaDFAkNClMStcaPH09aWhrdunXjL3/5C/fccw/Dhw8HaFj1z+FwkJqayk033UT//v159tlnG11DU/5ERCQaVFZW8uabbzJjxgzOP/98XnrpJQYOHOg9WFUFTz3l7W4NHw4PPgjLl3tXExSJcApTErVsNhu///3v2b9/P3/961+58MIL+dvf/saCBQsAWLVqFXv37m04/+uvv254nsqyLFwuFy+++CIjRozg4YcfZsWKFTidzmC8FRERkVZTV1fHzTffzP333899993Hr371K5JOfKbKbvcGpxNXCJw4sfEKgepaSYTSPlM+aJ+p6ON0OrHb7bz44ovMmzePzp070717d+rq6vjDH/5A9+7dAe9duunTp7N48eKG723fvj1TpkwhOzubqVOn0qNHj2C9DRGR1qd9pqKC2+3m7bffJi8vj507d2JZFrNmzeLSSy9tvMn9ifta5ebCli2NLzRsGOTkePe1GjfOr32tREKdwpQPClOyZs0adu7cydVXX43dbm90rKamhhUrVjRsGrx9+/ZGx8855xyys7PJzs7mvPPOa/L9IiJhTWEqKq1Zswabzca555576hOLimDxYm+4WroUKiuPH0tKgsmTj28anJbWuoMWaSUKUz4oTMnpKCgoaAhWy5Yto7q6uuGYulYiEnEUpsRftbXHu1Z5eU27VsOHHw9W6lpJGFGY8kFhSlrK365VTk4O5513npZcF5HwozAlLaWulUQIhSkfFKYkUNS1EpGIozAlgaCulYQxhSkfFKakNfjTtcrJyWl41kpdKxEJSQpT0hrUtZIwojDlg8KUtAV1rUQkLClMSWtT10pCnMKUDwpT0tZ8da1GjBjRaIVAda1EJGgUpqSt1XetcnPhk08ad63atYNJk44vv96rV/DGKVFDYcoHhSkJNnWtRCRkKUxJMKlrJSFAYcoHhSkJJSd2rXJzc/n2228bHVfXyjePx6KorIrC0gpqnB6cbg/2GJM4u0nf1CTSUxIwTcP3hUREYUpCi7pWZ0T1sWUUpnxQmJJQ5k/XKicnh6lTp9K9e/cgjjR4PB6LVTtLWbbtIGt3lbHjYAWmYWAzDSwsLAsMAwwMXB4Lj2UxoEsSWX1SmDi4C+P6pap4iDRHYUpClbpWPqk+BobClA8KUxIu1LVq7Ei1k7fXFTN/ZQGVtS6q6tyczoedASTExpDosDFzQgbXZKbRPt7eWsMVCU8KUxIufHWtTlwhMMK7VqqPgaUw5YPClISraO1aVde5eSpvK2+tK8YwoMbpOeNrxttNPBZcm5nGw9lDiI+NCcBIRSKAwpSEoyjtWqk+tg6FKR8UpiQSREvXKr+wjHvf2EB5tZMa15kXiR+Ls5kkx9uZe8NIsvqkBPz6ImFHYUoiQVHR8WAVoV0r1cfWozDlg8KURKJTda06dOjQaIXAcOha1brczPlgC4s2lATkTpsvcXaT6SN7MXvaUBy26LsLJ9JAYUoijb9dq5wcb9fKHtrT21QfW5/ClA8KUxLpwr1rVVnrYsYra9i6r7xV7rY1J85mMrRHMgtuHU2iI7R+JyJtRmFKIl0Yd61UH9uGwpQPClMSbcKpa1VZ62L6y59TcKiS2jYsFPUcNpOMzoksunNsVBQMkSYUpiSa+NO1ql96PchdK9XHtqMw5cORI0f4r//6L5KTk3nyySeDPRyRNuVv1yonJ4cxY8a0adeq1uXmunmr2bKvPCiFop7DZjKsRzJvzBwTNVMaRBooTEk0C9Gulepj21KY8sGyLGpqajAMg7i4uGAPRySoQqlr9d//3Mii9SVtOnWhOXE2k+nnpvHEFcODPRSRtqUwJeJV37XKzfWGq61bGx9vw66V6mPbUpgSkRapqanh008/bQhXJ+ta5eTkkJ2dHfCuVX5hGTf9ZU2bPEzrrzi7yYJbR0fdKkYS5RSmRE4uSF0r1ce2pzAlIgHRVl2r6jo3FzyznINHawMx7IDq0s7Bp7+8KCr32ZAopTAl4puvrtVZZzXe16qFXSvVx+BQmBKRgGvNrtXs9zbxzrrikJi+8GNxNpNrstJ4/LLInc4g0ojClMjpa2nXqqQEHn0U7rwTRo9uclnVx+BQmBKRVheortWRaiejnlwa1AdqfXHYTPIfmUz7+NDee0QkIBSmRM7M6XSt3nsPnn8e4uPhnXfgkksaTlN9DB6FKT+53W7Wrl3Lnj17AOjZsyejRo3CNM0gj0wkvPjqWv3kJz9p2Nfqx12r+SsL+P2S7SE1F/zH4u0mD0wZxO0TMoI9FJHWpzAlElin6lqZJniO1b+YGPjrX2HGDED1MZgUpvywePFi7r33Xvr370/v3r0BKC4uZseOHbz44otMnTo1yCMUCV/+dq1+evHFXPnXLSE5F/zHurRzsHrWJEzTCPZQRFqXwpRI66mthZUrvcHqvfdg586m59x9N54X5zLm6U9UH4NEYcoPgwYN4uOPPyY9Pb3R60VFRVx88cVs27YtSCMTiSyn6lol9c+k57W/psYdxAH6KSE2hnk3ZjK+f2qwhyLSulasgC+/hAkTYOTIYI9GJHL96U9w991NX+/YkZVrtnPXwvVU1oV+gYzE+hjZWxIHiMfjOelzHD169MDjCd12qki4iYuL4+KLL+biiy/m+eefb9S12hw3lNrQrxMAVDvdLNt2MKKKhchJnX++94+ItK7YWO//JiTA2Wd7n6Vq3x5uvpll2w5SFQZBCiKzPipM+eG2225j1KhRXHfddfQ6tqpKSUkJb775JrfddluQRycSuTIyMrjnnnu45557mPa/K9m0t7zR8c9+dRG9OiZwx4J1DO2ezP2TBzYcq3G62V1Wxd/X7OZvX+ziF5MHct+kAew7Us3kZz9tuIP31swxjM7oxJItB5i5YN0pxzOkWzsenTaMn/TuQLXTzeJN+3niwy1N7gZaFqzddThAvwUREYl6t9wC06ZBp07e56VOsHb5Sn48zUz1se1o9QQ/zJo1izfeeAPDMMjPzyc/Px/DMHjzzTd56KGHgj08kYjn8VjsOFjh17lFhyt5dVUhy7YdpH/nJH5z2TCuy0pj7vLv+O5gBd3bx/OLKd6i8m8jezI6oxPlNU5mv7fplNdNjI1hwW2jOa9fJ5ZtO0jJ99VcP6o3T1111knP//ZABZpFLSIiAWEY0KVLkyCl+hh8ClN+eOyxxxgyZAgPPfQQL7zwAi+88AIPPfQQgwcPDvbQRKJCUVkVpuHfw6rb9x/l8Q+28PPXN7Bs+0EALhrUhTq3h1n/+AaPx+Lm8/owqk8Ks6YOAeC3H21nf3nNKa97bVYaqUkOlm49wM9f38C1876gxunmkrN6kNYxvsn5pmFQdLjqNN+piIiI/1Qfg09hyg9ZWVnBHoJIVCssrcB2miv/dE12kNE5EYCyyjoA1hV9z+v5u7HHmCy4bRSd2zlYu6uMhauLfF5vWI/2AGzccwSAqjo3Ow9VEGMaDOme3OR8m2lQWFrZ5HUREZFAUX0MPj0z5Ydp06YFewgiUa3G6cFqMiP85H46rBu7njq+keH3VXXMX1nQ8PX/t3gbk4d0pVv7OGpdbmb9Y6Nf101N8j78W1nranit/oHfzkmOJudbQI0rPB4IFmmpjRs30r17d1JTmz5MXltbi8PR9P8bIhI4qo/Bp86UH+6//36+++67Jq9XVFTw8ssvB2FEItHF6fbg7/Tq+jnhc5d/x0PvfsP5v11OwQl3wI7Wunh7fTEAK3eUsvOQf3PNSyu8d+8SHcfvQSXGeueuH6poureHZVnUhfBO9CKBcPPNN2O32xu+rl/htqqqivO1yp9Iq1N9DD6FKT8sXbqU/v37N3z90ksvAZCUlMTcuXODNSyRqGGPMfFzSnjDnPDffbydt9YVc/SEO2X1PB6r0f/6Y8uxlQTP6dUB8BaKfp2T8Hgstu0/2uR8wzCItekjViKbx+Ohffv2DV+fe+65ACQkJFBXVxesYYlEDdXH4Iucd9KKYuvX9j/m2Wefbfh7zI9WVRGRwIuzmxi07m7pn/3qInY9dQk/Hdr1pMffXLubsso6Jg7uwks3jOStO87DYY8hd9M+dpc1fZDWAOJs+nyQyObxeKiurgagvLycnTt3UlVVhdvtxuVq+g81EQks1cfg0zNTfujYsSPvv/8+06ZN46233sLpdPLGG2+QlJREu3btgj08kYjXNzUJ10nuktWvYOR0B26J1ZP9HIDKOjc3vrKG2dOGMnFwF2qcbt5aW8ycD7c0e52+qYkBG5dIKLr22mu5+OKLmTJlCsuXL+eee+7h/PPPx2azMX369GAPTyTiqT4Gn2FF0kLvrWTr1q3cfPPNfPvtt1xxxRXMnj2b//iP/8Dj8fDcc89piXSRVubxWAx5bDG1J8yx7tUxnv974EJsMSZTnvvU7302TiY5zsaXs3/Ktv3lXD53VbMF43Q4bCbbHp+K4e/8C5EwlZuby7Zt25g6dSpDhw5l7dq1WJbFqFGjgj00kYin+hh86kz5YciQIeTn5zd67cMPPwzSaEQik9vtZsyYMezdu5cRI0YwYsQIhg0bxtlnn83w4cMZ0CWJTcfmZf/HxP7cPiEDW4zJxj1H+M7Ph2Sbc16/VFweD//59tcBKRQAA7smRUyhEDmVnJwccnJyGr7WdiIibaSkBDMvjwGVsWxyeFfUVH1sewpTfli4cCHjx4+nT58+jV6vqalh48aNKhwiAeDxeNi0aRM1NTXs3buX3NzchmO5ublk9enN5r3lWEBaxwRq6tys+PYQT+Zu9Xslo+Z8tHk/g2YvPrOLnMAwIKtPp4BdTyRkFRTAnj3Qty/06hXs0YhENqcTVq2CvDzvn43epcuzJs1kc+alWIap+hgEClN+ePrpp/n6668bvs7Pz2fUqFHExsZyxx138OWXXwZxdCLhrba2ls8++4zc3NyTLuhit9vp378/SXTg7XXFVNa5efDdb4IwUv/F22OYOLhLsIch0vrmzoVnn4VnnoEHHgj2aEQiT0nJ8fC0dCkcPWF1vMREmDSJiReN5e0yG5VOj+pjEChM+SE2NrbRP/JuvfVWNm3ahGlqMUSRligqKiIvL4+8vDw++eQTKitPvhN6RkYGy5YtIz09nX4ei0SHjcq60N/or53DxtiMyLrzJiIibaCZ7lODoUMhO9v7Z/x4cDgY57FIfPoTKp1N93QKNZFYHxWm/GCaJkVFRaSnp7N9+3Z2795NYWEh8fHxwR6aSFg4sfuUl5fH1q1bGx0fPnw4OTk5nH322cyYMQOAzMxMcnNz6dy5MwCmaTBzQga/X7KdGmfobvYXbzeZOSED04yc+eAiItKK/Og+kZMDU6dCenqTb1d9DC6FKT889thjjBs3jkGDBlFRUcGCBQsYO3YsHo+Hl19+OdjDEwlJp+o+tWvXjsmTJ5Odnc3UqVNJS0trOLZy5UrKy8t5+eWXm2w9cE1mGs98vL3N3kNLeCy4OjPN94kiIhKdWtB98kX1MXgUpvxw6aWXMm7cOHbt2sXw4cOJjY3l0ksvBdBUP5Fj6rtPeXl55ObmNtt9ys7OZuzYsU02w673pz/9qdmf0T7ezrWZaby9rpgaV+jdfYuzmVyTlUb7eHuwhyIiIqGkpAQWL/aGpyVLTt59qg9QJ+k++aL6GDwKU35KSUkhJSWl4WuFKJGWd5/OxMPZQ1i8eT81R0Nvbnj7BDsPTx0S7GGIiEiwtUL3yRfVx+BQmBIRv53YfcrLy2PLlsa7m/vbfToT8bExzL1hJDe+uiak5obH2U3mXj+S+NimKxKKiEgUaOXuky+qj8GhMOVDZWUlb7zxBgkJCdxwww3BHo5ImwtG98mXrD4pTB/Zi0XrS0JiOkOczWT6yDQy+6T4PllERCKDr+7TkCHe4JSTE7Duky+qj21PYcqH8vJyZs6cSbdu3RSmJCqEQvfJH7OnDWXLvnI27y2nNogFw2EzGdojmdnTInP6goiInCDI3Sd/qD62LYUpEQnJ7pMvDlsMC24dzfSXP6fgUGVQCobDZpLROZEFt47GYYvM6QsiIlEtBLtPvqg+ti2FKZEoFC7dJ18SHTYW3TmWG19dw5a95W06pSHu2B23BbeOJtGhj1IRkYgRBt0nX1Qf207kv0MRAWD37t0Ny5aHS/fJH4kOG2/MHMOcD7ayaENxmzx0G2f3zgGfPW1IxN9xExGJeGHYffKH6mPbUJgSiVD+dJ+ys7PJzs5m3LhxIdt98ofDFsMTVwzn8hE9uOf1DZRXO1vlLlyczSQ53s7cG0aSFcEP04qIRLw9e46Hp1N1n6ZOhT59gjbMM6X62PoUpkQiSKR2n/yV1SeFT395EU8t3spba4sxDagOwJ24eLuJx4JrstJ4eOqQiF3eVUQkYjmd8PnnkJt76u5TdjZMmBA23Sd/qT62HoUpkTAWTd0nf8XHxvD4ZcN5YMog3llXzLyVBVTUuqh2urEs/69jGBBvjyHJYeOOCRlcnRmZO7eLiESsE7tPS5dCefnxYxHUffKX6mPrUJgSCTPR3n3yV/t4O7dPyODWcX1ZtbOU5dsPkV94mB0HKzANA5tpYAGWZWEYBgbg8lh4LIuBXZPI6tOJiYO7MDajE6ZpBPvtiIiIL1HeffKX6mNgKUyJhLi6ujpWrlyp7lMLmabBhAGdmTCgMwAej8XusioKSyupcbmpc3mItZnE2WLom5pIeqcEDEPFQUQkLJyq+5SQ0HjlvSjoPp0O1cfAUJgSCUH13af6fZ8qKioajqn7dGZM06BPaiJ9UhODPRQRETld6j61GtXHllGYEgkB6j6JiIg0Q90nCWEKUyJBou6TiIjISaj7JGFEYUqkjaj7JCIi0gx1nyRMKUyJtCJ1n0REJFpUV1cTHx/fsArcKfnqPg0eDDk56j5JyFOYEgmgurq6hn2fcnNz1X0SEZGIlp+fz5IlS1i+fDndunVj4cKFzQcpdZ8kAilMiZwhdZ9ERCRaVFRUkJSU1PD1mDFjuP3223n66afJzMw8+TcdPAhTpsA33zR+Xd0niQAKUyInWLVqFb/+9a/5+c9/zpVXXonH48E0zUbnqPskIiLRpLS0lHnz5vHxxx9z9OhRZs6cySWXXEJaWhoTJ07ksssuIzMzs/npfcnJsHevuk8SkRSmRKAhNH377bfs2bOHzz77jCuvvLKhKFiWxfr163niiSfUfRIRkahQH46++uorvv/+e5577jni4uJ45pln2LNnD3PmzOGqq65ixowZXHvttZSWljJ9+nSuv/76phdbtgwGDlT3SSKOwpQIYJomTqeTjz76iF//+tcsXLiw0VQGwzBISEjgvffeA9R9EhGRyFJQUMAf/vAHYmJiuOGGGxp1ms455xx+8pOf0KlTJwAmTJjA6tWrAbj++ut5//33mTx5MgMHDuSaa64hPj6eSy+9lJiYGO/F4+LgrLOC9dZEWpXClAjeu28//PAD7du3p2/fviQmJrJx40bOO++8hmIyYMAA/va3v3HRRRep+yQiIhGjurqa5557DtM06du3L7/4xS945ZVXGDhwIACdO3cGjs/iWLx4MTNmzACgY8eO5OXlNVxr5syZrFmzhuzs7ONhSiSCKUxJ1DpxbrdhGLz55puMHTuWrKwsBg4cyNtvv027du0YPnw4AHa7nZtuuimYQxYREQm4iooK3n//fXbt2gV4w9Urr7zCo48+SmJiIuCtmaZpsmPHDkpLSxk8ePBJr+VyuaiqqsLhcPi3RLpImFOYkqhSv/Le119/zTPPPENCQgLgLRIej4clS5awZs0a/vWvf1FVVcWAAQMawpSIiEgkcrvdjBgxgh07djBgwADGjRvHggUL+OabbxrN0ACYP38+WVlZ9O/fv+F7d+7cyZ///Gc2bNhAZWUl8+bNA1CQkqhg+j5FJHzV1dWxbNkyHnzwQYYPH056ejp33XUXf/zjHyk/YX8LwzD44YcfOHr0KH379uWll15i/PjxXHfddUEcvYiISOuz2Wykp6ez8djGub169SI1NZXNmzcD3huOAEVFReTn53PvvfeycOFCHnnkEWJiYnA6nSQlJfHEE0/wxRdfcJaej5Ioos6URJxT7fuUlJTE5MmTycnJabRPBsCjjz7acBft6NGjLF26lA0bNjB58uQ2Hb+IiEircDrh88+9y5Rfdhkcm8LXoUMHevfuTX5+PldddRUdO3akpqamYcGJ+i1CXnvtNVasWMG0adMYNmwYOTk5uN1uhg0bxrBhw4L2tkSCSWFKwt6J+z7l5eU13Emr5+/Ke4Zh4Ha7G1bue+GFF9pi+CIiIq1nzx5YvBhyc2HpUigvh5QUuOqqhlNsNhsXXHABd911F0eOHKF9+/asXr2a66+/ntzcXHbu3Mm///u/M2bMGPLz85vfnFckCilMSVgqLi5uCE9Lly5ttvt0uvs+aeUhEREJa/Xdp7w8759vvml8fPBg74a5TmejPZ8yMzMZM2YMd955J1u2bGHw4MEMHjyYnj17MmHCBACmTJnSlu9EJCwoTDXD5XKRk5PD+vXrAThw4ACdO3fmZz/7Gc8//3xwBxeFAtV9EhERiTgn6z7VS0iAiRMhJwemToW+fZu9zHPPPceaNWtwuVxceOGFAA2r+YnIySlMncKaNWsaFimwLIvS0tKGcCWtz5/uU32A0r5PIiISNfztPmVnw4QJ3k1z/WC32xk/fnwrDFgkcilMNcNms3HfffcxZ86cRq8/8MADQRpR5FP3SUREpBn13ae8PFiypMXdJxEJLIWpU7j//vt5/vnnOXr0KAAjRozg8ssvD/KoIou6TyIiIifRSt0nEQkshalTSElJ4bbbbmt4Ruqxxx7TBnRnqK6ujlWrVpGbm3vS7lP9UqvqPomISNRR90kk7ChM+TBr1ixefPFFOnTooK5UC/nbfZo6dSq9e/cO4khFRETakNMJX3zhXThC3SeRsKQw5UPXrl2prKwkJiZGXSk/1Xef8vLyyM3NPWn3KTs7m5ycHHWfREQkuqj7JBJRFKaa4fFYFJVVUVhaQY3Tg9PtwR5jEmc36ZuaRHpKAqapcFVP3ScRkejQqD6mDMA59ALs8WnEbTug+ngy9d2nvDxvB0rdJ5GIYliWZQV7EKHA47FYtbOUZdsOsnZXGTsOVmAaBjbTwMLCssAwwMDA5bHwWBYDuiSR1SeFiYO7MK5falQVjxO7T3l5eWzatKnRcXWfREQig8/66HJj2GJUH0/kT/epPkCp+yQS1qI+TB2pdvL2umLmryygstZFVZ2b0/mFGEBCbAyJDhszJ2RwTWYa7ePtrTXcoFL3SUQkeqg+ngZ1n0SiVtSGqeo6N0/lbeWtdcUYBtQ4PWd8zXi7iceCazPTeDh7CPGxMQEYafD4233Kzs5m/Pjx6j6JiEQA1Uc/7d17fNlydZ9EolZUhqn8wjLufWMD5dVOalxnXiR+LM5mkhxvZ+4NI8nqkxLw67cmdZ9ERKKX6uMpnNh9ysuDr79ufFzdJ5GoFFVhqtblZs4HW1i0oSQgd9p8ibObTB/Zi9nThuKwheZdOHWfRERE9bEZ6j6JiA9RE6Yqa13MeGUNW/eVt8rdtubE2UyG9khmwa2jSXSExuKJJSUlDcuWf/LJJxw9erThWFJSEpMmTSInJ0fdJxGRKKD6eAJ1n0TkNEVFmKqsdTH95c8pOFRJbRsWinoOm0lG50QW3Tk2KAVD3ScRETmZaK+PgLf7tHixd+EIdZ9E5DRFfJiqdbm5bt5qtuwrD0qhqOewmQzrkcwbM8e0yZQGf7pP9QFK3ScRkegTrfVR3ScRCaSID1P//c+NLFpf0qZTF5oTZzOZfm4aT1wxPODXVvdJREROR7TUR+B496n+2acjR44fU/dJRM5AiExSbh35hWXeh2lDoFAA1Lg8LNpQzOUjegRkFaP67lP9ynvqPomIiD8ivT6q+yQibSViO1PVdW4ueGY5B4/WBnsoTXRp5+DTX1502vts+Oo+DR06lJycHHWfRESkWWFdH2tqmg8+6j6JSBBEbGfqybytlFc7gz2MkyqvdvLU4q08fpnv6QzqPomISCCFZX10OuGOO+C112DZMrjgAt/dp0GDvMEpJ0fdJxFpNRHZmTpS7WTUk0uD+kCtLw6bSf4jk2kfb2/0utPpZNWqVeTm5qr7JCIiARWW9bGqCq6+2rvaHniDUdeu6j6JSEiIyM7U2+uKMYxgj+LUTAPeWVfM7RMy1H0SEZE2EW71kYMHYdIkOPHG4sqVx/+u7pOIBFnEdaY8HosxT38SknPBfyzOqsX24a/ZvGljo9fVfRIRkUALp/rYpZ2D1bMmYfboDgcOND3hkUfg9tvVfRKRoIu4ztSqnaVU1rqCPQy/VDk9HDpqkpiYyOTJk8nOzmbq1Kmkp6cHe2giIhJhwqk+VtS6+LzgMON79jx5mOrRQ0FKREJCxIWpZdsOUlXnDvYw/GLGxvGzXz3N3NsuwuFwBHs4IiISwcKpPlY73SzbdpDx69d7V/DbvBm++sr7p6QEJk8O9hBFRIAIDFN3XdCPxy4dxh0L1vHxFu/drEenDeXWcX1ZtL6YXy76BoCx/Trxn5MHMqhbO0zD4FBFLZv2HOHeN74EYPrIXjxz9TkA1Lk8VNW52HW4io+37OeVzwpP6+Hdcf06cf/kgZzVsz1x9hhWFxzmuvmrAYMDnkQFKRERaXVrd5Wx8lcX0atjQljUyLW7DntPiIuDc8/1/hERCTERFaY8Hgun2/cHeNdkB3++KRN7jMniTfs5WusiIzWRyUO6Njn3cEUt//pmL2kdEzh/QGdGpHVg8pCuXD9/td/Fom/nJOLtMWw/cJRzenVodOzbAxVYloUR6k8Ei4hI2PJ4LHYcrPDr3FCpkaqPIhIOIipMFZVVERtj+jxvRFpHEmJtLNmyn39/88uG13+8TDnAgfIafvP+FgCG9UjmH3ePZWTvjtw6ri9//HSnX+NauLqIhauLuHVcnyZhyjQMig5X0Sc10a9riYiInK6isipMP0NJqNRI1UcRCQe+k0cYKSytAD9qxaGjNQBMHNyVBbeO4heTB5DVpyPlNafexHDz3nI+3uydFjFpcJczHi+AzTQoLK0MyLVEREROprC0ApvpX5gKlRqp+igi4SCiwlSN078pBRt2/8DLK3ZiWRYTBnTmvkkDeefOsfzrnvEkx526Wbfnh2oAOiUF5jknC6hxhccDwSIiEp5qnB4s/NsJJVRqpOqjiISDiApTTreHGqf3g9d+wnS/+ql/1SeErafytpH55FLueX0Dr68pos7l4aye7bk6M+2UP6Nnh3jAO088ECzLoi6Ed6IXEZHw53R7sKzjNx3DoUaqPopIOIioMGWPMSn5vgqAkb07AN6d1Eekef++67B3ukCP9nH0TknghyonH27cxyP/3MT/bT8IQJKj+btuw3ok89Nh3gdwP9l2MCBjNgyDWFtE/WcQEZEQY48xMQwoKvPWwXCokaqPIhIOImoBiji7ybvrS7hgYBduG5/BuekptI+30zc1kfIaJ//6ei8Ag7snM//GTL7c/T07D1USH2ty0eAueDwWq74rbXTNrslxPHbpUHp1SOCCgZ2JtZls2P09r64qBGBM3xTevOM8yqudnP34xycdV2Z6R67LSqN/l3YA9OucyDPTz2bnoUoWri4izhbTir8VERGJdnF2EwODN/OLmTS4a1jUyJLvq1UfRSTkRVSY6puaxJKtB7nvrS+5fXwG/TonUuvysGLHIX730XYOHfVOO/j2wFH+saGEzD4pDOmejAVs33+UeSsLWFf0faNrdkpyMGN0OpW1LrbsO8JHmw/wymeF1B1bgr1+yVaXp/m56H06JTL93ONTIzq3i2P6uWmsLjjMXz/fRV+tVCQiIq2ob2oSLo/Fkq0HwqZG5heW0aWd9mEUkdBmWJbl3xOpYcDjsRjy2OLT2izwTN0ytg+PXTqMJ3O3Mm9lwWl/v8Nmsu3xqdpHQ0REWk0w6iOcWY1UfRSRcBBRk5FN02BAl6Q2/Zlj+6WydlcZ8z87/SAFMLBrkgqFiIi0qmDURzizGqn6KCLhIKKm+QFk9Ulh895yPxeAPXMzF6xr8fcaBmT16RTA0YiIiJxcW9dHaHmNVH0UkXARUZ0pgImDu5AQGx4PrMbbY5gYoM1/RURETkX1UUQk8CIuTI3rl0riKZZuDSXtHDbGZujOm4iItD7VRxGRwIu4MGWaBjMnZBBnD+23Fm83mTkhA9PUfHAREWl9qo8iIoEX2p+oLXRNZhqhvkahx8LnTvIiIiKBpPooIhJYERmm2sfbuTYzjbgQ3Tk9zmZybVYa7ePtwR6KiIhEEdVHEZHACs1P0wB4OHsIySH6Ydw+wc7DU4cEexgiIhKFVB9FRAInYsNUfGwMc28YGXJzw+PsJnOvH0l8mKyoJCIikUX1UUQkcELrkzTAsvqkMH1kr5CZzhBnM5k+Mo3MPinBHoqIiEQx1UcRkcAIjU/RVjR72lCG9kjGEeSC4bCZDO2RzOxpmr4gIiLBp/ooInLmIj5MOWwxLLh1NBmdE4NWMBw2k4zOiSy4dTQOm6YviIhI8Kk+ioicuYgPUwCJDhuL7hzLsB7JbT6lIc5mMqxHMovuHBs2myWKiEh0UH0UETkzhmWF+o4TgVPrcjPng60s2lBMjdPT6j8vzu6dAz572hDdcRMRkZCl+igi0jJRFabqrd1Vxj2vb6C82kmNK/BFI85mkhxvZ+4NI8nSw7QiIhImVB9FRE5PVIYpgOo6N08t3spba4sxDagOwJ24eLuJx4Jrs9J4eOoQLe8qIiJhR/VRRMR/URum6h2pdvLOumLmrSygotZFtdPN6fxGDAPi7TEkOWzcMSGDqzO1c7uIiIQ/1UcREd+iPkzV83gsVu0sZfn2Q+QXHmbHwQpMw8BmGliAZVkYhoEBuDwWHstiYNcksvp0YuLgLozN6IRpGsF+GyIiIgGl+igi0jyFqWZ4PBa7y6ooLK2kxuWmzuUh1mYSZ4uhb2oi6Z0SMAwVBxERiS6qjyIixylMiYiIiIiItEBU7DMlIiIiIiISaApTIiIiIiIiLaAwJSIiIiIi0gIKUyIiIiIiIi2gMCUiIiIiItICClMiIiIiIiItoDAlIiIiIiLSAgpTIiIiIiIiLaAwJSIiIiIi0gIKUyIiIiIiIi2gMCUiIiIiItICClMiIiIiIiItoDAlIiIiIiLSAgpTIiIiIiIiLaAwJSIiIiIi0gIKUyIiIiIiIi2gMCUiIiIiItICClMiIiIiIiItoDAlIiIiIiLSAgpTIiIiIiIiLfD/AC1EX2oXctz7AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1080x288 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "label_dict = {i: str(i) + \", \" + str(t) for i, t in solution.items()}\n",
    "edge_color1 = [\"red\" if solution[v] == (solution[u] + 1) % n else \"black\"\n",
    "              for (u, v) in G1.edges]\n",
    "edge_color2 = [\"red\" if solution[v] == (solution[u] + 1) % n else \"black\"\n",
    "              for (u, v) in G2.edges]\n",
    "\n",
    "# 在图上画出找到的最优路线\n",
    "fig, ax = plt.subplots(1, 2, figsize=(15, 4))\n",
    "for i, a in enumerate(ax):\n",
    "    a.axis('off')\n",
    "    a.margins(0.20)\n",
    "nx.draw(G1, pos=nx.circular_layout(G1), labels=label_dict, edge_color=edge_color1, ax=ax[0], **options)\n",
    "nx.drawing.nx_pylab.draw_networkx_edge_labels(G1, pos=nx.circular_layout(G1), ax=ax[0], edge_labels=nx.get_edge_attributes(G1, 'weight'))\n",
    "nx.draw(G2, pos=nx.circular_layout(G2), labels=label_dict, edge_color=edge_color2, ax=ax[1], **options)\n",
    "nx.drawing.nx_pylab.draw_networkx_edge_labels(G2, pos=nx.circular_layout(G2), ax=ax[1], edge_labels=nx.get_edge_attributes(G2, 'weight'))\n",
    "plt.axis(\"off\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "上面给出的左图展示了未被选中的交易路线，右图展示了该算法找到的最优套利回路。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 结语\n",
    "\n",
    "在本教程中，针对套利问题，给出了一种寻找最优套利回路的经典量子混合算法。对于 $n$ 个货币市场，在该解决方案中，我们需要 $n^2$ 个量子比特。在样例中，我们设定的不同市场间的汇率比较特殊且把交易成本包含在内，所以会产生高达 $150\\%$ 的正向收益率。用户在使用时可以根据需求，自行设定真实的货币汇率。\n",
    "\n",
    "在实际的金融市场中，一般不会有如此高的投资回报率，而且考虑的货币市场数量比较大，需要考虑的因素也更多。因此，随着货币数量的增加，所需的量子比特数也急速增加，实现难度也就越大。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "_______\n",
    "\n",
    "## 参考文献\n",
    "\n",
    "[1] Orus, Roman, Samuel Mugel, and Enrique Lizaso. \"Quantum computing for finance: Overview and prospects.\" [Reviews in Physics 4 (2019): 100028.](https://arxiv.org/abs/1807.03890)\n",
    "\n",
    "[2] Egger, Daniel J., et al. \"Quantum computing for Finance: state of the art and future prospects.\" [IEEE Transactions on Quantum Engineering (2020).](https://arxiv.org/abs/2006.14510)\n",
    "\n",
    "[3] Rosenberg, G. \"Finding optimal arbitrage opportunities using a quantum annealer.\" [1QB Information Technologies Write Paper (2016): 1-7.](https://1qbit.com/whitepaper/arbitrage/)"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "3b61f83e8397e1c9fcea57a3d9915794102e67724879b24295f8014f41a14d85"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.7"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}