OpenJij/OpenJij

View on GitHub
docs/tutorial/ja/physics/ClassicalSystem.ipynb

Summary

Maintainability
Test Coverage
{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# core interfaceを用いた古典イジング模型の数値シミュレーション"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "このチュートリアルではOpenJijのcore interface (core python interface)を用いて、ランダム相互作用およびランダムな縦磁場をもつイジング模型の数値シミュレーションを行います。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "まずはGraphを定義し、今回数値シミュレーションを行いたい系の$J_{ij}, h_i$を定義します。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openjij.cxxjij.graph as G\n",
    "#問題サイズを100とします。\n",
    "N = 100\n",
    "\n",
    "graph = G.Dense(N)\n",
    "#sparseの場合は以下を使用します。\n",
    "#graph = G.Sparse(N)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "次に、$J_{ij}, h_i$を設定します。\n",
    "今回は平均0、標準偏差1のGauss分布から生成される値を設定します。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pip install numpy #乱数生成にnumpyを使います。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "mu, sigma = 0, 1\n",
    "\n",
    "for i in range(N):\n",
    "    for j in range(N):\n",
    "        #Jijの値が大きくなりすぎてしまうので、1/Nで規格化を行なっています。\n",
    "        graph[i,j] = 0 if i == j else np.random.normal()/N\n",
    "\n",
    "for i in range(N):\n",
    "    graph[i] = np.random.normal()/N"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "縦磁場に関しては、`graph[i]`でも、`graph[i,i]`でもどちらでもアクセスできます。また、イジングモデルの定義上、$J_{ij}$と$J_{ji}$は自動で同じ値となります。試しに以下のように出力を行なってみましょう。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5\n",
      "0.5\n",
      "-0.6\n",
      "-0.6\n"
     ]
    }
   ],
   "source": [
    "graph[20] = 0.5\n",
    "print(graph[20,20])\n",
    "print(graph[20])\n",
    "graph[12,34] = -0.6\n",
    "print(graph[12,34])\n",
    "print(graph[34,12])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## システムの設定 - System -\n",
    "\n",
    "続いて計算を行うためのシステムを定義します。  \n",
    "\n",
    "ここでは古典イジング模型の数値シミュレーションを行いたいので古典イジング模型のシステムを作成してみます。`system.make_classical_ising`で作成できます。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openjij.cxxjij.system as S\n",
    "\n",
    "mysystem = S.make_classical_ising(graph.gen_spin(), graph)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ここで、1つ目の引数にはランダムに生成したスピン、2つめにはGraphそのものを代入します。これにより初期のスピン配位が`graph.gen_spin()`となる古典イジング模型のシステムの作成ができます。\n",
    "\n",
    "システムに直接アクセスして、値を読むことも可能です。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 1. -1.  1.  1. -1.  1.  1.  1. -1. -1. -1. -1. -1.  1. -1. -1.  1. -1.\n",
      "  1.  1. -1. -1.  1. -1.  1. -1. -1. -1. -1. -1.  1. -1. -1.  1. -1. -1.\n",
      "  1. -1.  1. -1.  1.  1.  1.  1.  1. -1. -1.  1. -1. -1.  1. -1.  1. -1.\n",
      "  1. -1.  1.  1.  1.  1.  1. -1. -1.  1.  1.  1. -1.  1.  1. -1.  1.  1.\n",
      " -1. -1. -1. -1.  1. -1.  1. -1. -1.  1. -1.  1. -1. -1.  1.  1. -1.  1.\n",
      " -1. -1.  1. -1. -1.  1.  1.  1.  1.  1.  1.]\n"
     ]
    }
   ],
   "source": [
    "print(mysystem.spin)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## アルゴリズムの実行 -Updater, Algorithm-\n",
    "\n",
    "Systemを定義した後はUpdaterを選択し、Algorithmを実行していきます。\n",
    "\n",
    "### Updater\n",
    "\n",
    "Systemに対して使用できるUpdaterは決められており、古典イジング模型に対するUpdaterは主に\n",
    "\n",
    "- [SingleSpinFlip](https://github.com/OpenJij/OpenJij/blob/ec41aecfbac7e4c895e1e7a1718f06eb7ffae0ba/src/updater/single_spin_flip.hpp#L40) (メトロポリス・ヘイスティング法によるスピン1つずつのアップデート)\n",
    "- [SwendsenWang](https://github.com/OpenJij/OpenJij/blob/ec41aecfbac7e4c895e1e7a1718f06eb7ffae0ba/src/updater/swendsen_wang.hpp#L45) (SwendsenWang法によるクラスターアップデート)\n",
    "\n",
    "の2つが用意されています。\n",
    "Algorithmを実行するには、**スケジュールリスト**が必要になります。まずスケジュールリストを作成するところから始めましょう。\n",
    "\n",
    "### Algorithm\n",
    "\n",
    "#### スケジュールリスト\n",
    "\n",
    "スケジュールリストは`(パラメータ, モンテカルロステップ数)`のリストで与えられるものです。パラメータ部分に入力する値はSystemによって異なります。例えば、古典イジング模型ならばパラメータとして温度の逆数である逆温度$\\beta$を設定します。\n",
    "なので、古典イジング模型の場合には、温度スケジュールのリストを指しています。\n",
    "ここでは、例として以下のように設定してみましょう。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "schedule_list = [(0.01, 10),(10, 80),(0.1, 30)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "この場合、逆温度$\\beta=0.01$で10モンテカルロステップ、$\\beta=10$で80ステップ、$\\beta=0.1$で30ステップの計120モンテカルロステップを実行することを意味します。  \n",
    "アニーリングを実行するにあたっては、逆温度は等比級数で増加させていくことが多いため、以下のように`utility`にある`make_classical_schedule_list`を使うと逆温度は等比級数で増加させていくスケジュールを簡単に作成することができます。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[((beta: 0.100000) mcs: 20), ((beta: 0.199474) mcs: 20), ((beta: 0.397897) mcs: 20), ((beta: 0.793701) mcs: 20), ((beta: 1.583223) mcs: 20), ((beta: 3.158114) mcs: 20), ((beta: 6.299605) mcs: 20), ((beta: 12.566053) mcs: 20), ((beta: 25.065966) mcs: 20), ((beta: 50.000000) mcs: 20)]\n"
     ]
    }
   ],
   "source": [
    "import openjij.cxxjij.utility as U\n",
    "schedule_list = U.make_classical_schedule_list(0.1, 50, 20, 10)\n",
    "print(schedule_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "上の例では$\\beta=0.1$から$\\beta=50$まで、各温度で20モンテカルロステップ計算しながら10段階で温度を変えていく設定例です。計200モンテカルロステップの計算を行います。\n",
    "\n",
    "#### Algorithmの実行\n",
    "\n",
    "続いて、Algorithmを実行します。`Algorithm_[Updater]_run`のように書くことで、指定したUpdaterで計算を行います。次の例ではSingleSpinFlipを実行しています。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openjij.cxxjij.algorithm as A\n",
    "A.Algorithm_SingleSpinFlip_run(mysystem, schedule_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "一瞬で処理が終わりましたが、この間に計200モンテカルロステップの計算が行われています。\n",
    "> `A.Algorithm_SingleSpinFlip_run(mysystem, seed, schedule_list)`とすることで、seedを固定したまま計算を行うことができます。結果に再現性をもたせたい際に使うことができます。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "callbackを使用することで、Algorithmの実行中に1モンテカルロステップごとのシステムを取得することができます。古典イジング模型の場合は、システムとパラメータ (逆温度)を引数を持つ関数を作成すれば良いです。  \n",
    "例として、以下ではシステムのエネルギーの値を記録するcallbackを作っています。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "energies = []\n",
    "\n",
    "def callback_log_energy(system, beta):\n",
    "    #graphは以前にGraphモジュールにて定義したオブジェクトです\n",
    "    energies.append(graph.calc_energy(system.spin))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "このcallbackを用いて同じAlgorithmを実行します。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "#スケジュールをもっと長く取ります (計20000モンテカルロステップ)\n",
    "schedule_list = U.make_classical_schedule_list(0.1, 50, 200, 100)\n",
    "A.Algorithm_SingleSpinFlip_run(mysystem, schedule_list, callback_log_energy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "記録したシステムのエネルギーを、横軸をモンテカルロステップ、縦軸をエネルギーでプロットすると次のようになります。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Requirement already satisfied: matplotlib in /opt/conda/lib/python3.9/site-packages (3.5.2)\n",
      "Requirement already satisfied: kiwisolver>=1.0.1 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (1.4.4)\n",
      "Requirement already satisfied: python-dateutil>=2.7 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (2.8.2)\n",
      "Requirement already satisfied: pyparsing>=2.2.1 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (2.4.7)\n",
      "Requirement already satisfied: numpy>=1.17 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (1.23.2)\n",
      "Requirement already satisfied: pillow>=6.2.0 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (9.2.0)\n",
      "Requirement already satisfied: cycler>=0.10 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (0.11.0)\n",
      "Requirement already satisfied: fonttools>=4.22.0 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (4.37.1)\n",
      "Requirement already satisfied: packaging>=20.0 in /opt/conda/lib/python3.9/site-packages (from matplotlib) (21.3)\n",
      "Requirement already satisfied: six>=1.5 in /opt/conda/lib/python3.9/site-packages (from python-dateutil>=2.7->matplotlib) (1.16.0)\n"
     ]
    }
   ],
   "source": [
    "!pip install matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAGwCAYAAABRgJRuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAABiU0lEQVR4nO3dd3gUVdsG8HvTC6RAEkJIQgi9BghFeosUUcFeEAERBUFRQAWVYgNExfaioChgA+UTROkd6TV0CCQQSkiBQBqk7/n+CFmyyfad3dly/64rF8nMmZlnsmTn2VMVQggBIiIiIjvnIncARERERFJgUkNEREQOgUkNEREROQQmNUREROQQmNQQERGRQ2BSQ0RERA6BSQ0RERE5BDe5A7AmpVKJa9euoXr16lAoFHKHQ0RERAYQQiA3NxdhYWFwcdFeH+NUSc21a9cQEREhdxhERERkgitXriA8PFzrfqdKaqpXrw6g7Jfi5+cnczRERERkiJycHERERKie49o4VVJT3uTk5+fHpIaIiMjO6Os6wo7CRERE5BCY1BAREZFDYFJDREREDoFJDRERETkEJjVERETkEJjUEBERkUNgUkNEREQOgUkNEREROQQmNUREROQQmNQQERGRQ2BSQ0RERA6BSQ0RERE5BCY1DkAIgYLiUrnDICIikhWTGgcw7vd4NJm6Hpcz78gdChERkWyY1DiANSdSAQC/7b8kcyRERETyYVLjQEqVQu4QiIiIZMOkxoGUCudJakpKlXKHQERENoZJjQNxlpzmmy3n0WzaBpxMyZY7FCIisiFMahyI0kmyms83nUNRqRIfrD4tdyhERGRDmNRIrFQpEDV5DXp/tt3q11ZY/YrSKilVQjhJYkZERNJjUiOxRbsvAgAu3Lht9WvbQjqgVAocu5KFohLj+rzkF5Wi48wteO7H/RaKjIiIHB2TGoml5xRIcp51J1Ltss/IdzuSMGjebry69IhRx+29cAOZt4uwOzHTQpEREZGjc5M7AEejUJjfCHT8ahbG/FaWFCTPHmjwcSm38s2+trl+2HkBALDhVLrMkRARkbNhTY3EpOjXci49z6TjtpzNkODqRERE9olJjdQkyGr+jk8x/yTOwBY6ERERkc1gUiMxRYWsxtSRPLuTbmjd9/v+yxi+6ADuFJUYdc5f9ibjz4NXTIrHGPY+AouIiOwX+9TYmXdWngAALNlzCWN61jfomOu5hZi66hQA4JG2deDu6iC5LDMoIiKqwEGebrZDgn7CBskrLDa4bMVaHYeaoM+BboWIiMzHpEZiFXMaU/MHe66AkGL0V0WFJaWY+OcxrD5+TdLzEhGR42FSY4OkTgysSVM/opJSJUYsOoC5GxOQX1SKpxbsxYIdSQad79d9l/HXkasY93u81KHqte1sBp5asBdXbt6x+rWJiMh47FMjMZcKCYklW0cUdlSfsz3hOrbd/Qr09cD+izex/+JNvNxDf5+gG3mFZl27pFQJF4UCLi7G/75GLD4IAJi0/Bj+eLmTWXEQEZHlsaZGYg7VZ8UEmmqZCissmZBfXGq1WApLStHg3XWIfmetWee5ebtIooiIiMiSmNRI7K8jV80+hyMv6mjNGqaTKTmq7y9nsgmJiMjRMamRWHrOveaS09dydJTUTmnHOY2+lEWu7kKZt01vxrJ0zJl5hUYlsrduF2Fl/FXkF1mv1ouIyB4wqbGgf21wxE5egXGT9hmr4qNZqRRIzVZfj6pi15Z1J1ItGkvFaEpsNFP84+BlxH60GZ9tTDD4mGGLDuCNP47h/X9PWTAyIiL7w6TGgkpKBTJypVm1WyqxH2226Pkr9j+ZtPwYOs3ain+P3UvuKjY/lS/aKbWiEiX2JN1Q68vzxPy9SMu2rdcCAN7+q2wyxXnbDBsNBgDHr5at3l7x92orCktKsTvxBgqs2HeKiKgckxoL+mn3RXT4eIvJzVC2oKRUidwCwyf6q2jF3TWs1p9KU22zRvPTB6tP4dkf9uPtv46rbf9t/yXLX9wAV27eMfl3auum/X0KQxbuxzsrTsgdChE5ISY1VvD3UftdoPKh/+1GyxkbkZFje7Uc2vy67zIA4MrNfD0lrS/5xm10m7MN7T/WXWOWV1hilx3G/zhUtr7YCi7KSkQyYFJjY3IM/ARvrQ63Z1LLapm2nM2Q5Hy2MrHgxlNpeP/fUygpVeote/WWdMlR+WKlBcXar3suPRctpm/A6F8PS3ZdIiJnwKTGhsRfvoVWMzbKHYZKeoXaGTusNNDppV8OY9HuZINqFO4UlWKtxTs137NkTzIAYMOpdKtdk4jIETCpsSGfrD8rdwhq3vq/4/oLWVnFJQuEBHM2G9qs9uOui2ZfC3C85JCIyJYwqZGYOa0rpTY27PhyhQTCRlqN0OPTbSYfayO3oJdt/S8gIrIfTGokZs6DU9tcKkIIJKTlorhC/49cC883A9hGElA5hoq/IjnXv1q0+yL6frFD1iH7ttI/iYjIVjCpkZi2T9nFpUr8e+yaSQ/BPw5eQb8v/8NLPx9SbVt8t99FZd9sOY9jV7KMvoZGGp6Z2XeKsel0OopK9Hew1UTbiB5tzTK6ai2kaH4yNTF4/9/TOJeehy82nTPyeiZdziwHLt7E5xsT1JJiIiJHxKRGYtqeWT/svIBXl8bjga92aj1W04N9y5l0LPjvAgBgW8J1vdf/fNM5DJq325BQ9dJ0L88u3IdRPx/C11vOm3ROS8/s+42JcZmqUEdydz23EOtOpKKkVImNp9KwMv6qpH1qDB3y/eSCvfhmayJ+2Wsb8/QQEVmK3SQ1s2bNQvv27VG9enWEhIRg8ODBSEgwfGp5a9H2yX/z6bKRLDfyjFvxeeSSQ7h447bZcZlC072cujuR4EoT5yHRlvRpq8HQV7EhhMA7K0/g2+2JKCguxee6ak7MqCYxZc6YB77eiTG/HcGPuy7ipV8O440/jhnUMdmYS+06fwOLdxvWiVmu/0f6FJUocafI8s2pROT43OQOwFA7duzA2LFj0b59e5SUlOCdd95B3759cfr0afj6+sodnmwka2qq5ExqDhIz8iQ/r9TNL8evZuP3/WWT7c1ZL12SW3G5B510JCDXc8sW0Zy17t6otux8aWcSfu7H/QCAprX90DG6pqTntpZOs7Yg83YRTn/QDz4edvOWREQ2yG7eQdavX6/28+LFixESEoLDhw+je/fuGo8pLCxEYeG91Zlzciy/XEGAtzsyKz0QKyceS/Yk48l2EfD2cDX7elI1NVU2YtFBi5xXavlGrjH09ZbziKzhU2V7TkEx5m48h37NQ/Hq0njcyDN9VW8pZeQUIMTPS2+5lCzbmz3ZUOV/LwlpuWgTGShzNERkz+ym+amy7OyyRf1q1KihtcysWbPg7++v+oqIiLB4XPVDqlXZVrnaf/o/pzBng3xz0twpKsGKI1dxS0dtROZt/Q/14lKlqjZCDpcy7+gvVMHxq1mYu+kcXv/jaJV9n6w7i8V7kvHMD/s0JjRyjTTalqB9JufbRVw0koioIrtMapRKJV5//XV06dIFLVq00FpuypQpyM7OVn1duXLF8sFpaI5QKNSHIgPA7sSy6fLPp+diwFc7sf5kWtUDLWT6qlOY8OcxDF+svTbGkH4dD32zC+0/3ozz6bkSRgfsSbyBBTuS9PZjycgtNKr/ya07VZt+yq9x3gJNbaazTGdqKUaLERHZMrtpfqpo7NixOHnyJHbt2qWznKenJzw9Pa0UlXYKKKo8oMt/HL/sKM6k5mD0r4fROiLAKvH8c+waAPP745xNK0tm/j12DRP6NjY3LJVnF5b1E2lQqdbL7M6kGjKgzzaeQ3pOIQ5cvGnSKTefMW4pA6YVRESWY3c1NePGjcPq1auxbds2hIeHyx2OycofbnmFJVW22ZvyuA0ZITRzreHNbhWXRACkXViyol/26R/qrO3ecgpKZG2CIyKie+wmqRFCYNy4cVi5ciW2bt2KevXqyR2SwVwUVRMW5d2HpKuL6X01MnV0Zv3v3HWcvpZj0lBkQH+CVfG85d+aOsxbm8p9RjT3dZH0kloJAPsvZGocFXXrjuHD9CuG+9+5qvMOCSGw9MCVCuU5azARkaHspvlp7Nix+P3337Fq1SpUr14daWllfVD8/f3h7e0tc3S6KRSKqi0fd38255E149/TWvc9/9MBAECn6JqY/VhLM66iX3lfDXP6BWnKvXLyi9WSvol/HjPoOEuIv5yFp77fh+pebjgxo5/e8tpmXK4YbvlrVO5aVj683M0fEUdE5Kzspqbmu+++Q3Z2Nnr27InatWurvv744w+5Q9NLoajaSVNU2Geqf+/2jdFl74XMKsOeTbnm/B1JeHN51aQCsGxikZx5b+RYanbVieu2GNOnRYJqHUPW3MotKEbM+xuNPnfn2Vux9az20U66WPI1+O/c9SpNgUREtshuampMbUaxBQoFcOu2+sib8uYnFyu0n1RuwtD3qyxViiorhl/KvKM2hLri3m+3J+H1uEYmxZZbUAxvd1et+cbaE7prfxbuMmw2XUuqHPqOc9eNnj+n3OI95t/Pkcu34OflXqWjtSl/QvsuZKpqlJJnDzQ7NiIiS7Kbmhq7oeHh7KJQVJkczZo9JYzNm+q/s9boazR6b50qUTNGyxkbMfDrXRofuPaaxkqaqFY41aHkm3h1abzmYnfLpWbn49Fv9yBu7o4qZQSAl385hNeXaT6HJocv3dK5/3puoV1/4CAix8KkRmqa5qnRUEzTJ3mlhRd7NMSJq9kmH6upecgQCVrmubHXh6UZfb91dgx+fP5evU2OuiYkTM3Kx4ZT6fj76DUUmFiTVNEv+y6h/ceb8dlG21uDjYicE5Maibm7GfZES8+x3jDgyhHpWlk6K9+wkTyaEpjyxS5NoalyQ1tnW5OvIenZdFxHhtmHDcn/pM6Zp686CQCYty1J2hMTEZmISY3EHo4Jq7LN0IecpZ6FD//PMutDWdqSvfrnj7FFupqfjK18yi0owfhl8TqXSzCFFJVgNlCxSESkxm46CtsLd9eqeaKhuYqlWlsM7bR6KfM2hv54QH9BC7CnlqbKI64q5zAa/guYbO7GBNwuKsWqo7qbnQxJiGVavoqIyGpYU2MFuh4mtvSgmaRlyDapG7nkkNZ92XeKMf2fU1r3G/t6G7popSFJoT0ljkREpmBSIzFNDw5dzU/XskzrXGsJmmbLtZYSB2nLmLn2DK7cNH05hxMppnfUNoahi1vaUtJNRKQPkxqJaXpU6HouVFz7yZm9/MthuUMw2aLdyarvz2foXrHc0rUlus5/0koJk6nk6GBNRI6FSY3ENA1D5nu1bThq5qrk2vy2/7Lqe1uub8rUUROXW1CsdR8Rkb1gUiMxDzdNHYU1ZzVn00wfAk32yVYS3Iq594+7LqLljI1YeuCy9gOIiOwAkxqJaVp1W9vkcv2/3Kn2s6H9HMj2ZN8pq+mIv5yls9zPNjhM/cPVZQujTllxQuZIiIjMw6SGSAJTVh6XOwSN0rTM8qwtfS4sKa2y7hcRkb1gUiMxczqCnkyRrznKXpcksBX7LtyUOwSN3vvbuNqXljM24oGvduot97kBSyPYwrIfRORcmNRI7E6R/Y1mOn41C60/2ISk67flDsVu2WpSmJajuaZGW9eeohKl1ubScqnZ+fhma6Lea3eYuRnvrDQ8qbLV3yER2Q8mNRJ7+y/76peQdD0Pry87iux8jn4xx607tvf7KylVaq39Myd9KCg2bE2uG3lF+H0/Ox8TkfUwqXFyfT7fIXcIJJHEjFx8tOa06uf/bdNfm0JE5Ei49hPhwg02O0nh2+3yJRE/703GseXqk+t9ufm8TNEQEcmDNTVEEpmzXn/nWUs5dtW42YLN6b/Cvi9EZKuY1BCRVQkhsOzAZYvN8ExEzovNTxLJulOEq7dMX8iQyJoMrWupOBt2qVJonFzSWNvPXcfkuxP9Jc8eaPb5iIjKsaZGIl0/2YYHv9kldxhEkjlxNRufrD+r+vmJ+XskOW9iep4k5yEiqoxJjUS42jbZkzXHU/WWeeh/6kn6kbtLQHA1bSKyVUxqiJzQlBUnsOzAZZNm/dXUUfiRb3ejqET7/DVfbDqHx7/bg4LiUqOvR0RkKPapIXJSk1ecQLFESxnEX85Co/fWad3/1Zay4eUr41O0lmENEBGZizU1RE5s6t8nNW6fssIyC3QWlxo2GzERkSmY1BBRFUsPXJE7BCIiozGpISKjmNtgVbGV6c9DTJ6ISDpMaohINm/9n2WauYjIOTGpISKjmNudV+pVFg4m38S1LE58SUQc/UREduzYlSw8MX8vAM5OTESsqSEiO3b40i25QyAiG8KkhoiMwjW6ichWMakhIqviHHtEZClMaojIaq7eysfSA5flDgNFJUocu5Jl0jIRRGS7mNQQkdV8/98FJF2/bfLxm0+nY9VR7UstGOq1pfEYNG83vt2eaPa5iMh2MKkhIqMY23qka6HLim7dLsKF63m4dbtI4/49STfw4s+HMH7ZUaRlFxgZhbr1p9IAAD/uumjWeYjItnBINxEZLDEjF19tMa524/LNOwaVG7H4oOr7r55ujUGt66jtf/aH/arvs/KLEOrvZVQcROT4mNQQkcHi5v5nwjE7jD5m/LKjaBdVA3UCvI0+loicF5ufiMgmZeYVat2XkJZrxUiIyF4wqSEim6RrYNL4ZUcluYaC48uJHAqTGiKyScKMRaK2JWSg7xc7cPxqlnQBEZHNY1IjgcuZhnWEJCLDmTODzIhFB3EuPQ8jFh3UX5iIHAaTGglcumn6vBtEpJkUq3nnFZbo3M/GJyLHwqSGiGwUZ/slIuMwqSEim7TlTAbeXXkChSWlcodCRHaCSQ0R2aRvtyfht/2XsXh3ssHHFBSXYu7GBIPLZ+UX48Ulh7D+ZJrGc93W03xFRLaFSY0EFGyZJ7KYrWczcDYtx6Cy321Pwtdb7814XD5iWwiB+Mu3kJ1frFa+VCmw+Uw6Rv96uMq52n20Gc2nb8CdIiY2RPbC7pKaefPmISoqCl5eXujYsSMOHDggd0hEZEH7L95E/y93VtleUFy1WarypHzlnY03nk7HI9/uQb8vDJ8RubyTcWJGnhHREpGc7Cqp+eOPPzBhwgRMnz4dR44cQUxMDPr164eMjAy5QyMiK2s1YyMKDVwss7x5KS3HvIUwici22VVSM3fuXIwaNQojRoxAs2bNMH/+fPj4+OCnn36SOzQisrKiUiX2XsiUOwwisiF2k9QUFRXh8OHDiIuLU21zcXFBXFwc9u7dq/GYwsJC5OTkqH1ZAmdaJ5LHf+euq76/ebuoyv6KfWr0ySssQWIG15Qismd2k9TcuHEDpaWlqFWrltr2WrVqIS2t6sgFAJg1axb8/f1VXxEREdYIlYhkMPTH/WYd3+aDjYib+x+OXL6ltl2KSQCJyDrsJqkxxZQpU5Cdna36unLlitwhEZGFnLqWY1ataXFpWfay6XS6RBERkbW5yR2AoYKCguDq6or0dPU3nPT0dISGhmo8xtPTE56entYIj4hsQOValYJiwzoSE5FjsJuaGg8PD8TGxmLLli2qbUqlElu2bEGnTp1kjIzrxxDZssOXbukvREQOwW5qagBgwoQJGDZsGNq1a4cOHTrgyy+/xO3btzFixAi5QyMiG6Cp+el6bqFR56hc28MuNUT2w66SmqeeegrXr1/HtGnTkJaWhtatW2P9+vVVOg8TkXNip14i52ZXSQ0AjBs3DuPGjZM7DCJyUJVrewwZDk5EtsFu+tQQEemjqfnp8s3b+PvoNYPPIYTmJRg0qbyWFBHJi0kNETm0mWvPGlV+/o4kNJm6XmcZIQTmbkxAzPsb8c8xwxMmIrIsJjVS4PAnIrtRWGJYLUy5yo1PF67nod1Hm1WrgU9bdVKiyDS7XViCxbsv4lpWvkWvQ+QImNQQkcMoUerv/3I+3bxVtz9YfRqZGpZksJQPV5/GjH9PY/C83Va7JpG9YlJDRA6jpNTyk+1Zu9/wjrvrW2UYOTSdyBkxqZGAkpOWEjktDo4ish1MaiRwI4+foIhswdk0/atsbz5j3NpOTFqI7AeTGgl4ubvKHQIRAUjNLtBb5svN5yW9pjmLaBKRtJjUSKC2v5fcIRCRjThy+RYW7rwApQGdlolIWnY3o7At4ic1IudR3nG3XOXmqUe/3QMACKrmicFt6lgrLCICkxoiIj0Mq3H579x1bKnQXyfpunlDx4nIeExqiIh0MLSj8PM/HTDinAI38ooQXN1TsusTEfvUEBFZhK5kZOKfx9D+481qNTtEZD4mNUREOhSWSD8R1Yr4FADA3E3nJD83kTNjUmNlri7sVUxkT4Ys3I+jV7IghMCqoykay8zblmjSuU9dy0FuAVf6JpIKkxor+/uVLnKHQERGGjxvN55asA/jlx3VuP/TDQkmn3tvUqbJxxKROiY1ElAYsUw3h38T2acDyTeNKi8MHDVFRNJhUkNEZAGJGXmYte4MblpxRW8iZ8ch3UREFrDhVNnIpuQbt7FgaDuTz8PaXSLDsaaGiMgMQs9EMsevZpt5frMOJ3IqTGqIiKzkTlGJ3CEQOTQmNVbGqmQix6Iw8I96wY4kNJu2AauPX7NwRETOi0mNlUXV9JU7BCKyopK7q3XPWncWAPDm8uNyhkPk0JjUWJm3u6vcIRCRhPT1qbmeW2ilSIiISQ0RkRlyCthPhshWMKmRgDH9ZNinhsi55ReXqv1saJ8cItKPSQ0RkYz0NV9xZmIiwzGpsTJ+KiMiIrIMJjVERBb27zHtw7iN/aAz+a/jeG7hfiiVrMEhqozLJBARWdirS+ONPia3oBi/7LuE9Bz10VPLDl4BABy7moU2kYGSxEfkKJjUEBHZoI9Wn8Efh65o3c+KGqKq2PxERCSjv+NTVN+vPn4Ng+btxpWbd7D/YqaMURHZJyY1REQyWnMiFbsTbwAAxv0ej2NXsjB11UmZoyKyT0xqiIhkduTSLbWfczmhH5FJmNQQEcksp6C4yjZO/0BkPCY1RERE5BCY1BARyYy1MkTSYFIjAT2znBMRGUVx94uIjMOkhohIZov3JKv9zM9JRKZhUkNEJLOiEiX2JN1Q36ihquZceq51AiKyU0xqJODIq+h2qFdD7hCInMLZVP0Jy3t/c/4aIl2Y1JBOXz3dWu4QiJyOIf1p9PUtLiguxb/HriHrTpEkMRHZAyY1DqJvs1ro26wWvNylfUl9PNwwqW8jSc9JRPplVFrI0lgz157Bq0vjMfTHAxJFRGT7mNTIYEjHSMnP2al+TXz/fDvU9veW/NzjejeU/JxEpK5yI3ZeoXmzCpevKXUiJdus8xDZEyY1EjB2SPfHj7SUPIbymmgh8fhyTp9BZB3HrmQZVV7fn7rj9vQj0o5JjUxc7CRZsJMwiezeP8eu6S2j7e9RCIErN+9I/qGGyN4wqXEQ5TOS8i2NyP5l5GruT1Px7zu/qFT1/VdbzqPbnG2Yu+mchSMjsm0mJTU9evTAzz//jPz8fKnj0Sg5ORkjR45EvXr14O3tjfr162P69OkoKmKv/sqk/qDG6duJrO/yzTt6y8xce0b1/ZebzwMAvtmaaLGYiOyBSUlNmzZtMGnSJISGhmLUqFHYt2+f1HGpOXv2LJRKJRYsWIBTp07hiy++wPz58/HOO+9Y9LqGsoXaEUvlHkxpiGzT6dQcuUMgsjkmJTVffvklrl27hkWLFiEjIwPdu3dHs2bN8NlnnyE9PV3qGNG/f38sWrQIffv2RXR0NB5++GFMmjQJK1askPxaRES2qvKHjI2n0vDE/D2yxEJki0zuU+Pm5oZHH30Uq1atwtWrV/Hss89i6tSpiIiIwODBg7F161Yp46wiOzsbNWronu22sLAQOTk5al/25L2BTQ0ua6kaFbY+Edmul345jIPJt9S2jf7lMEpKlTJFRCQvszsKHzhwANOnT8fnn3+OkJAQTJkyBUFBQXjwwQcxadIkKWKsIjExEd988w1efvllneVmzZoFf39/1VdERIRF4jGFIX1V/LzcjT6vIy/ZQET6rT+VhnUn02yjXZzIykxKajIyMvD555+jRYsW6NatG65fv46lS5ciOTkZ77//PhYuXIiNGzdi/vz5Os8zefJkKBQKnV9nz55VOyYlJQX9+/fHE088gVGjRuk8/5QpU5Cdna36unLliim3q5fFhlEaU0tSPvpJ6o7C7FVDZDMKSgyrgTF34j4ie+VmykHh4eGoX78+XnjhBQwfPhzBwcFVyrRq1Qrt27fXeZ6JEydi+PDhOstER0ervr927Rp69eqFzp074/vvv9cbp6enJzw9PfWWk4OLAijVU0YBIKiaJ27kGTBd+t1sZlDrMMzblmR2fLZs6H118cu+S3KHQWR1hn6AUnK+GnJSJiU1W7ZsQbdu3XSW8fPzw7Zt23SWCQ4O1pgQaZKSkoJevXohNjYWixYtgouLfU+x46JQQF/9sEKhwAeDmuOV347oPV95c9b4Po3QKjwA7aNqYG9SJtpEBqDzbMv2b7KmDwY1Z1JDTsvQXEWpZFJDzsmkzEBfQiO1lJQU9OzZE5GRkfjss89w/fp1pKWlIS0tzapxSMnVgCmFFQAeaFnboPO53E1qPNxc0K95KGr4emBgq9oIC1BfC2rV2C54/+HmBscpRUdhTzfzEtDfX+yo+j7M35tz55DTMnQdp1Ile9eRczKppqZNmzYaHywKhQJeXl5o0KABhg8fjl69epkdIABs2rQJiYmJSExMRHh4uNo+W5gW3JQIXAx4MBvz7Da0bExEAGIiAjD9n1OGn9wM7aMCcfyqaQvqvd2/CZJv3Ean+jUljorIsUldUXP40i1cvnkbj7QJ11+YSEYmfYTu378/Lly4AF9fX/Tq1Qu9evVCtWrVkJSUhPbt2yM1NRVxcXFYtWqVJEEOHz4cQgiNX/bKkLWfjElq7GUtKWOM6VkfnzzeSi2BrlzzRERVKYWQtIv/Y9/twRt/HMNRIxfdJLI2k2pqbty4gYkTJ2Lq1Klq2z/66CNcunQJGzduxPTp0/Hhhx9i0KBBkgTqaOYNaYuhPx7QWcaYkUeO3iSzdNR9uJaVj2ZhfnKHQmQXLPGR71LmbbSOCLDAmYmkYVJNzZ9//olnnnmmyvann34af/75JwDgmWeeQUJCgnnRObBuDfV3kDaupkZ/4Zq+HoafUEKje9Q3+xyd6tfEY7Gs+iYyhLZK7ILiUjz0zS58uPq0dQMishKTkhovLy/s2VN1au49e/bAy8sLAKBUKlXfOzpbaAWTujnLnGMqG9+nIf4Z1wVNQqsbfMzjTGCITKatm/DaE6k4kZKNH3ddtHJERNZhUvPTq6++itGjR+Pw4cOquWgOHjyIhQsXqhaZ3LBhA1q3bi1ZoM7ImCYlQ2pqTIrhbhPYwFa1seZ4qknHu7go0Co8AG6uttlENvfJGBy4eBPLDlpmckYiW1HCod7k4ExKat577z3Uq1cP//vf//DLL78AABo3bowffvgBzz77LABg9OjRGDNmjHSROqHyFOCNuEb49/g1JGbkaS9r4Xzh8ydiTEpqTGXN9OfRtuFoUcefSQ05jOTMO5pnFWZOQw7O6KSmpKQEM2fOxAsvvIAhQ4ZoLeft7TyjVE5czTLr+Cah1eHt4Yr4y+rnKU9Uxsc1xPi4hoiavEbrOSxWU3P3tF7urhjYsjbWnDA9sTGmmY7vvUSm+33/ZblDIJKF0X1q3NzcMGfOHJSUcG2Rcr+Z+QaiUCjQIsy/6nYj6iuMTWo+GGTYBHxuFTrrfPV0a4zsWs+o69iLir+96CBfq133RQf9fZJt4pR85OhM6ijcp08f7NixQ+pY7Jabq/lLNrzYrezh9kSFDrKWmHyv3POdogw8770Tu7m6oH1UDaOu83g70zr8ytn7ZsMb3fFAy1CLnb9Dhd/huN4NLHYdIiJnY1KfmgEDBmDy5Mk4ceIEYmNj4eur/sn24YcfliQ4ZyGEQN2avkj4qD883Vyx/PBVAIaNaOrRKBjHr2ahRyPD1tAyV1zTEPRqHIymtf3w7XbdC2c2CKmmlqTZYvPTouHtq1zPRaHAo23CsfaEZZbh8HS/lwS7S5AQExnKFkZqElmSSUnNK6+8AgCYO3dulX0KhQKlpfrWn3YsUs1s7OnmqvazIQ+8xSPao1QpJKktMoSbqwsWjegAAHqTmlbh/jY/KWBzDZP5KQB0iFavkZp4fyN8vumc2rZD78Wh3UebtZ67Z+NgbE+4XvX8Nv47ISKyVyY9CZVKpdYvZ0toLOHl7tHoXL+mQbUvCoVC0oTmidhwvBHXSLLzVaTrWT7/ubbqZfWca0ALaZqHtCUYlbf6ebtXKRNUzVPnwqTa9tjoyHZyAuZ+/GJCTrbO7KdhQUGBFHFQBVMeaIrfR91ngdoX/W9Ift7uaBVRtdOySSq9g1as0GoYUk1tX/8WtbH61a4Gn/qbZ9pgw+vdzYkOgOFrZlWujXvCjMkBfTzuVZCyNYCISDomPTVLS0vx4Ycfok6dOqhWrRouXLgAAJg6dSp+/PFHSQMk26WrlsIULercS6b0fSB0c3VB40ozFPtrqE0pN75PQ43bvT3KmvwqX67iJ9IAH3fEVFrvZtajLTUep+0cFfVpGoLeTULwmpaYKmofFYiPBrfQW45ICsWlSrlDIDKLSUnNxx9/jMWLF2POnDnw8Li3nlCLFi2wcOFCyYJzdOVrMfVsHCJzJPfYSkfCdkaOsgJ017r0bV5L4/byWhNdt73mtW5wc7n3pxJV00dVi6Yr+Qqu5qlxu7urC34a3h4T7tffzBfi54Xn7qurtxyRLt//l4SoyWuQmpWv2jbhj6NqZc6k5qDxe+vw6YazVo6OSDomJTU///wzvv/+ewwZMgSurvc6t8bExODsWef7gzA1EVjzWjfMeawVXo/T/4ldE18PV/2FjGTJeSyMOXOdAPkmb9SUqFT8vXSody/h0jSX0PdDYxHXNARvD2ii8fzmrnKc8FF/s44n5zNzbdn78tdbE1XbVsSn4KvN55F1pwgAMHvdWSgFMG+b7gEARLbMpKQmJSUFDRpUnV9DqVSiuLjY7KDsXaCP9maQikL9vfBk+wh4uZuWnJhSm1HRy92jy/7tEa3aZq2amoqXkWrxSik7Mfp6uKJjvRpoHRGA2n7qC7NOfbBZhYtWPbZv81AsHNYeNTSsir7pje6IqOFjcBzdGwZV2VZ5lByRqb7YfA5v/d9xucMgkoxJSU2zZs2wc+fOKtv/7//+D23atDE7KHvn6mIfc4+83b8JNr3RHZP7q9coxIQHAADczRymoys/qu51r7NsPQvM4NsmMsCs4xUKBZa9dB9WvtIZLi4KtWSvute9pPWlbtEajtauYS31fkD6pgN4IjbCqPMTGetA8k25QyCSjEnz1EybNg3Dhg1DSkoKlEolVqxYgYSEBPz8889YvXq11DHaHXsZ9ejioqjykAWAGr4eOPReHHz0NG+Zc5sd6tWostaVlFa+0kXnWlmA/iauijU/YVrK1gm0bDOZi8SdsS0lsoYPLt+8I3cYZAKjJsW0lU53RFqYlNQMGjQI//77Lz744AP4+vpi2rRpaNu2Lf7991/cf//9Usdo8yr3Q7HWY8gSyVP5m1aQlk6uUjFmXStDVX7DrebpplqpuPL1BrasjQl9DZ+PJ7i6J/4e28WsfkyebpavwVvyQge4KIChPx4w6figap64kVdo9HFcU8j+Xc81/nUnsjUmJTUA0K1bN2zatEnKWMgGaKuRsHV+Xsb9V543pK3Wfdr65pjTwXfFK53RWEOtmNTMWS7jkTZ1MPfJGNSbslZte3SwL9xdXJCQnmtueGSjlEqB06k5codBZDazPjoWFRXh6tWruHz5stoXWcbiEe3Rok7Vaf0NpatmZ/GI9ni+U10M7xJl8vkrK1Wqf3qvWJOiLZbqnmXJScU5awzxmIbOxjMeLluJfEzP+kadS0r73+mDjW90R9vIQPh6Vk28bKl+QwihMaHbMqEH2tYN0HOshYIiixNCoJQvIDkIk2pqzp8/jxdeeAF79uxR217+puhsSyUY+n7Q0sgHdWU9G4egZ+MQvX1FtNHVPaP83FL659g1fP2M5o7j/ZuH4jsNa0cdfC8OhcVKnRPpVbR8dCesPnYNE/s2xt/xKWr7Ho8NR+8mIajh64HT1yzzKbSahkSlolp+XqhVafSUvTFkVBmfifbtjIG1NBtPp6OwRIkn27EDO9kmk5Ka4cOHw83NDatXr0bt2rW5HogBjk67X+8D0FK+eCoGH685g/nPxcpyfU0qz9Bbzsvd1agh7u2jaqC9jqHtmoZVS2lAi1AMbFUbSqXAupPGr+pd3cD/Ex8Oao6pq04ZfX5jaPo7/ml4O4tek+QnADz8v90GlV1zPBVrjqeiS4MgWeeSItLGpKfs0aNHcfjwYTRponlyMWdT+UNqo1rVkVGp012Aj/QPV0NTyUfahGNw6zqyJ5+v9GqA15bGY2DL2hY5vzmVBdFBvois4YMAA+cYKufm6oJ5z7bFvguZJiU1cr8mFZU/pN5/uDmm/3MKc5+MQe8mmmdiroyjYuxXbkGJ0cdczrzDpIZsksnz1Ny4cUPqWBzCmQ/66x0KLQdLPDwn9WsMABjSMdKg8g/HhGHX273wjZYmKUuqePua3ozdXF2wbVJP/P1KFytGZboXu9aT/Jyv9CrrezSscxROf9APj7bVPilir8bqHZKZ0jiXmWvPyB0CkUYmJTWffPIJ3nrrLWzfvh2ZmZnIyclR+3Jm3h6ueHdgUwRX90S4hecwkdvL3aOxbVJPfDhI84KLwdWrDgsPD/Sx2Nwrus5asSLhl5EdNJZxdVHY7rwwlZLS9x5shjWvGb6quSEqrh5e8XtDaJrvyBDmjNYi+SRdz5M7BCKNTEpq4uLisG/fPvTu3RshISEIDAxEYGAgAgICEBgYKHWMNq9y1Xvdmr448E4f1TIEjkqhUKBekG+VRGDK3TWPlr10nxxh6RVVU/oZjKU0oEUoAOCZDrbbGbNi0jP9oWZ4OCZM8mv0bmI7C72SOhtN/YlM61Ozbds2qeNwONboK+Fqo7UKL/eoj5e6R9tUfxFLhyLl6b94qjWe7nAT90XrXtvL3G4sz3SIwNIDV0w69r0Hm+LCjdsY1qkunu4QiRVHrmos16FeDRy4aNo0/L0aB2Pr2QyTjiXLsqW/baKKTKqp6dGjB1xcXPDDDz9g8uTJaNCgAXr06IHLly+rrdpNlmXLbyy2HJslSNmnxMvdFT0aBVt84cpO9asulmmo2v7eWDe+G57uoL0/1fznYvHBoOYmX+MZHecmeTnXXzfZE5OSmr/++gv9+vWDt7c34uPjUVhYNtInOzsbM2fOlDRAe5CaXSDLdW20okY2hiYW9phv2VLIwztHGVSuf4tQNAn1w5kP+iPUhLl63FztY2FYp2RL/yGJKjDpXeOjjz7C/Pnz8cMPP8Dd/d4Q2C5duuDIkSOSBWcv7hTJM9mgiz0+nUkyco2iLp+p2VDeHq7YNqmnxn38L2yfTBkGTmQNJiU1CQkJ6N69e5Xt/v7+yMrKMjcmx2Hhd2wmNc5N1yKS0cH6O0Mb97/HvP9r3h6uiKzhY9Y5iIj0MSmpCQ0NRWJiYpXtu3btQnS0Y4/4MUan6JoALJfbMKcxja319zF13h5dNTWT+jZWfd8mMgBAWcfyga3uTXxo7YoeU1byDtEwLQARkTYmJTWjRo3C+PHjsX//figUCly7dg2//fYbJk2ahDFjxkgdo91qEFINWyb2wNGpfS1yfltcU8iUvhPWYGN5jJqHDBgObWytXMXS85+Lxciu9bBlQg/Me1Z9dXJjVzfXpn+LUJNmmNU3u/Ty0Z0M7sMDAP2aGzYDMhE5JpOSmsmTJ+PZZ59Fnz59kJeXh+7du+PFF1/Eyy+/jFdffVXqGO1a/eBq8Ddy6n195j8XiwEtQjE+rqGk5zXHX2M6o2uDICx5QfPEdrpINcX+g3drIZrVNn0lc1s1qHUYomr64Ln77o0IijCwOaeWnxemPtgMUUGWm5/Hx8MNO9/qpXNuGYWGJixtC5d+NLhsQse6NX2N6sPj5yXt35qzKilVAgDOpeci/vItmaMhMpxJH9MUCgXeffddvPnmm0hMTEReXh6aNWuGatWqSR0fadC/RSj6352gzVbE1g3Ery92lDWG9wY2Q/uoGuje0PFmqfX1dMO2ST3Vms5q+Hpg/evdMGn5MZxMUZ/Ju5a//hozqSuvTJmNWVtToBuH9slq8Z5kxEQE4In5e+UOhcgoZtU9e3h4oFmzZlLFQmQWL3dXDGpdR+4wLEZTAtAk1A8RgT5Vkpq2kYGY/lAz1K1pO51z20QG4PLNO2rbtKUunc2YQ4fM99Earu1E9kmaBnUiG9cwpDqah/mhhq/0q6XLTVvr3Ygu+he9tGZn4Q8eboFVR68ZVDbShpIxIrIfTGrIKbi6KLD6VWkXgLQVpowqAozrPN0mMgBLD1w26Trl/H3cEejjjlt3is06DxGRNpyyk5yGQqGw2HBuTSuS2zohDO9X83jbcMx5vBU2T+hh3jUr/axQAH3udi5+tI3jNh06os82JKg6FBPZCtbUEEmgfnA1fPJYS7tLbibc3wgz/j2NJ9uF6yzn4qLAk+2kXzXc3dUFC4bG4vLNO0jJyseK+BSzzifTJMtO6X/bElHD1wMvdNXfzElkLUxqiCTyVHt5FmA0dUS8QgEM6xyFno1DZJvtNzrYF26uLogOroaUrHxJz90gpBoSM/IkPSep25OUyaSGbAqbn0h2lpw/xRmYUzuhUCgQFeRr0nBsjbEYmWGFB0qfTC0a0R4T72+ETW9UXcqFiBwbkxqSzfLRnfDuA031zipLjsOjwsrbMx9pafBx859rixZ1/PB0e/1NYL0ah+DVPg019p/q39y25neyd5vPpGPLmXS5wyBSYVJDsmkfVQOjukfb3FpM9sbk5ifJp9/Tv67WD8+3Q50Ab8x7ti2e7Wh4c13/FrWx+tVuqKehVu+vMZ21HvfLyHszXH/5VGu8cX8jg69Jhhm55JDcIRCpsE8NkZ2LsqM5XWIiArB7cm+Tj28Z7l9lW2zdQK3lu1WYXbpBSDWbXgPM0RUUl2LV0RT0ahyCEBtdI47sH2tqiOzc6/c3wnP3RWLZS/cZdZwlHvBSreOlTef6QVgwNBb3N5N24crO9WtKej6qava6s3j7rxN45Ns9codCDszukprCwkK0bt0aCoUCR48elTscItlV83TDR4Nb4r5owx7MT7ePQONa1dGnqfbFJ21Zv+ahiA7W3Lncx8NV63HhgdpXEe/Z2PHWC7M1v+8vm7xR6lFuRBXZXVLz1ltvISwsTO4wiOzW7MdaYf3r3eDppj0BkENTI1ZXr9wfaOYjLdEmMgDj+1RduX7X272weUJ3BPh4aO1/NKxzlDGhUiVfbT6vt0wRJ+ojK7CrpGbdunXYuHEjPvvsM7lDIbJrluqcbc55g6p5Yvfk3jg67X6jj322YyRWvtIFNatVnfwwPNAHDUKq6zze080Vr/ZuYPR1qcwXm89h2d1lNOZtS8STC/aioLhU5qjIGdlNUpOeno5Ro0bhl19+gY+PYR0jCwsLkZOTo/ZFRJbz3sCmqOnrgSkDmph0fJ0AbwT4yLPo6Itdo2W5rqOYvOIEDl+6iU83JODAxZtoMnU99iTeAADkFzHBIeuwi6RGCIHhw4dj9OjRaNeuncHHzZo1C/7+/qqviAjpp3knonuig6vh0HtxeLlHfblD0cjUxT/JMFduqveXeXbhfgD8vZP1yJrUTJ48WbXIoLavs2fP4ptvvkFubi6mTJli1PmnTJmC7Oxs1deVK1csdCdEVM4a8w7xIWmb9l/MrLJt5OKDmL3urAzRkDOSdZ6aiRMnYvjw4TrLREdHY+vWrdi7dy88PdXby9u1a4chQ4ZgyZIlGo/19PSscgyRLSpfp6hBSDW5Q3FoOicc5Bw2Zlt6oOoHxy1nM2SIhJyVrElNcHAwgoP1D6X8+uuv8dFHH6l+vnbtGvr164c//vgDHTt2tGSIRFax5IUOWLInmaNwLKwhk0Yih2YXMwpHRqpPp16tWtkbU/369REeHi5HSESSqhPgjXceaCp3GA5PioU7x/Ssj++2J0kQDRFJzS46ChMRWZqhXYFcudYCkc2yi5qayqKioiw+HTsRERHZF9bUEJHdeaZ9WZN0LzOWN2gSWjYh34AWoQCAah5uiKrpgzoB2pdTICLbZpc1NUTk3KKCfHHy/X7w1bHWkz7Navth+ehOqOZZ9jbo4qLAlok9IYRAg3fXqZV9sl04/jx01ayYicjyWFNDRHapmqebeXPiKIDqXu5q53B1UcDNVf1tcdHw9pjzeIzaNi93w94632XnbyKrYlJDRFRJ+dDvTtE10auJ+mrmCgVwdFpf1c99moQgwMdd43lGdq2HMH8vywVKRGqY1BARVfLrix3xVv/GmDekrcb9Xu73mr0ah1bHMx0iNZZzcVHgv7d6WSRGIqqKSQ0RUSW1/LzwSs8GqOFr2OKaE+5vhG4NgzTuq9ycRdqVlCrlDoHsHP/aiMgpxdYNlOxc7q4uWpMaqirpeh6u3LyDV5fGY+zvR6BUChy5fAsN31vHiQ3JLBz9REROZduknjiUfBOPtjVtNnJTuiYH+rjj1p1inWWG3lcXv+y7ZFJM9qbP5zvUfn4iNhyz152FEMAn689iTE/bXOWdbB9raojIqdQL8sUT7SLgKsGSCRXpWixT3yitNpEB6FS/pqTx2JOCYjY7kTSY1BARGaHyXOaGjCp/sl2Ezv1ebqbPt0NE9zCpISKysIl9G+HRtnXkDsOmVazN4jI4ZComNUREEnigVW2t+9xdXdCzcYjW/QICzv4cr5jIjF92VL5AyK4xqSEiMkLl1qZGtcrWkKoT4I3jM/qifrCv9YNyAGfTclXf/3PsmoyRkD1jUiOBJ9uVjaKIqMGF8Iicxb/juuKDQc3xUKsw1TY/L3fJOyATkeE4pFsC4YE+AICuDUxfMZiI7EvLcH+0DPevst3Zm5GI5MSaGgnwTYyI9NHV+dXPS/PaUaQuO78YG0+loaiEQ8BJMyY1EjJnwWAishN6/tBNeR94f1BziCqDxZ3H/x2+YlC55xbux0u/HMbcTecsHBHZKyY1REQSMqXmtra/c/fH23wmw6ByJ1KyAQCrjqYgIS0Xyw9d4fBvUsM+NRJw5k9YRCQNTTMSu7koUKLk+4sm/b78DwDg4+GGgTqG05NzYU2NhNj6RETaUhBPN81vt3UCtNfSrBrXRYKI7FNiRtmil5pUrJw5dS3bShGRPWBSIwHWfhKRPn2a1kLXBuoreY/sWg9rx3fTekzzsKqjq5xF3Nwd6DZnG45eyZI7FLIjTGokxI7CRKStj4e7qwt+fbEjPnsiRrWtXd1A+HuXjXxiM7Zm01ad1Lmf77tUEZMaIiIrejw2XPW9h5YmKTKcrtXRyfnwL0oC/HxFROUah1bXW+aNuEaIa1pL53pQ2iwa0d6UsOzW8avZWLjzgtb9rKmhijj6SUL8xEBEHw5qgRq+HniqXaTWMuPjGpp8/l6NQzDv2bYY+/sRk89hbz5acwYt6zhv/yIyHGtqiIiMoO+jS81qnvhocEuNSyjoYsiAg9r+XgDglEOY03IKNG7nR0mqiEmNFDj8iYisYOqDzQwqN6pbPQtHIi+1TtVsf6IKmNRIiH9bRGRJmhYA9/VwrbKtT9NaVojGurR9duTbLlXEpEYCrKchIrloev+pF+Rr9Tisif0XSRsmNRLinxkRyWlw6zCsfrUravl54cunWssdjqQqNjlV/J415FQRkxoiIiNY6iFas5qHUeXff7g5AOCrp9uotkXW9EULJxglVHE5rC83n8eepBvyBUM2hUmNBNhPmIjM1Sm6Jl7t3QBPVJicDwB+e7GjxvLDOkfh3EcDcH8zzf1nHG2G4sl/nVB9fz23UG3fcwv3WzscslFMaiSkYD0oEZlIoVBgYt/GeCgmTG17F7X1otTfY3TNSOxoH7YKS5Ra93EhcyrHpEYCjvaJiIi0q+XnJXcIBnG0pIbIEJxRmIjIAAuGxmJvUmaV5iFr8tYwfFtNhUxGyayGnBCTGiIiA/RrHop+zUNlufbkAU1wJjUH3dSaonRjSkPOiEmNBPiBiIikoqlr3uge9Y0+j+AbEzkh9qmREPsJE5GtuL9ZWa1S64gAeQMhsiLW1EjgUuYdAEC6lgXXiIgMZdZsuRU+WdXw9cDZD/vD1UWBhu+ukyAyItvHmhoJrDmRCgBYeyJN5kiIiO7xcneFC6uQyYkwqSEisiHMQYhMx6SGiIjs3sKdF/D7/svILypVbdt5/jqm/n0SBcWlOo4kR8I+NUREZPc+WnMGAPDOyhP4781eiKzpg6E/HgAAhFT3xKt9GsoZHlkJa2qIiGyIWa1PGoZxO2Nr1m8HLqn9nJKVX7Z9/yW8/X/HoeS6Cg6LNTVERA7MKR/fWm763ZUnAQBxzWppXQiU7BtraoiIyKHoS+RyC4qtEgdZH5MaIiJHwaFT5OTsKqlZs2YNOnbsCG9vbwQGBmLw4MFyh0REJC2J8xIXJ81zuEyEc7KbPjV//fUXRo0ahZkzZ6J3794oKSnByZMn5Q6LiEhSZs0orOl8CgUOvReHklKB+2ZtkfTctur7/y7A3fXe73HZwSuY/VgrGSMia7GLpKakpATjx4/Hp59+ipEjR6q2N2vWTMaoiIjsQ1A1T7Wf20cF4mDyLZmisY5525K07ltzPBV9mtSCv4+7FSMia7CL5qcjR44gJSUFLi4uaNOmDWrXro0BAwborakpLCxETk6O2hcRkbOTujbI3mw5m4ERiw/IHQZZgF0kNRcuXAAAzJgxA++99x5Wr16NwMBA9OzZEzdv3tR63KxZs+Dv76/6ioiIsFbIREQmYV9f6zhyOUvuEMgCZE1qJk+eDIVCofPr7NmzUCqVAIB3330Xjz32GGJjY7Fo0SIoFAosX75c6/mnTJmC7Oxs1deVK1esdWtERFZXzdNV7hCIZCVrn5qJEydi+PDhOstER0cjNbVsFeyKfWg8PT0RHR2Ny5cvaz3W09MTnp6eWvcTETmCDwe3wKbT6Rh6X5TBx8x/Lhajfz0MAKgX5IuLN25bKDoi65E1qQkODkZwcLDecrGxsfD09ERCQgK6du0KACguLkZycjLq1q1r6TCJiKzGlNanoffVxdD7jHsv7N8i1KxrEtkiu+hT4+fnh9GjR2P69OnYuHEjEhISMGbMGADAE088IXN0RER2plIWM6h1HQBAdS+7GBBrEs5b4xzs5n/wp59+Cjc3NwwdOhT5+fno2LEjtm7disDAQLlDIyKya53q18SDMT1QJ8AbJ1Ky8cm6szh0ybGGfC/Zkyx3CGQFdlFTAwDu7u747LPPkJ6ejpycHGzatAnNmzeXOywiIkkFVbdcP8CQu+fuq2Exx/rB1eDl7or2UTUw42HHe2+d8e9puUMgK7CbmhoiImdQP7gaZj/aEjWrSZ/crB3fDUcu3ULvJiFq2ys3zXBYOdkrJjVERDbm6Q6RFjlvUDVP9G0eqrecs0/OR/bLbpqfiIiIiHRhUkNERGrY/ET2ikkNEZGTU1TKYjQlNXUCvK0UDZHpmNQQEZGaujV8q2zrGF1DY9l6QVXLEsmFSQ0REanx9nDFsel9cfqDfjrLLXvpPitFZBkJablyh0ASY1JDRERV+Hu7w8dD9wDZ+6JrWikay+j35X+4cD1P7jBIQkxqiIjIIB8ObiF3CJKbufYM3vjjKAqKS3WWKylVWikiMgeTGiIiJxcV5GNQOWMXzbQHm89kYGV8Cn7cdVFrmb/jU9Dg3XVYdyLVipGRKZjUEBE5qd2Te2PD690RUt1Lf2Et60E6yujv67mFWve9/sdRAMCY345YKRoyFZMaIiInVSfAG41Dqxt1zCs966v9PHlAEylDko0QAkqlwN6kTGTfKZY7HDIRkxoiIjJYo1rqSVDf5qEY2Kq2TNFIRwD449AVPPPDPsR8sBF7Em/gys07codFRuLaT0REZDBNE/P5edn/o0QpBKasOKH6+dmF+wEAybMHyhUSmYA1NUREZLDKsw8DgNDS38aeFBZzdJMjYFIjoZiIALlDICKyKE0dg0uV9p/VLD98Ve4QSAJMaoiIyGDe7q5VtjlATkMOwv4bQomIyGp6Ng5G7yYhaFHHX7VNOEL7k5Gu3rqDf4+l4tmOkfD3dpc7HLqLNTUScpT5GoiIKitPW9xcXfDT8PaYcH8j1b6hnRxvUj5tXl8Wj9uFJXjk2z34ZP1ZvLvyhP6DyGpYUyMh5/usQkQEtIkMlDsEq/n76DVE1vBRTda3JylT5oioItbUEBERaaGpaS2jwuzDztj0ZsuY1EiIzU9ERMCQjpFwd1XgwVa1seH17tj4Rnc8b6dNVGk5BTr3O8LIL0fC5iciItJLX2fYmY+0xOnUbPRoFIKejYMx4+HmcHe997n5yXYR+HnvJUuHKbmSUt1JS05BCa5l5SMswNtKEZEuTGqIiEirr55ujRVHUvB6XEOd5Z7tGKlzv5urfdZld5uzTW+Zsb8fwc8vdMDB5Jvo2iAYHm5sBJELkxoiItJqUOs6GNS6jtnnaVyrOuKahqCgWIldiTckiEw+17LVm6QS0nLR9ZNtyM4vxvDOUZjxcHOZIiOmk0REZHEKhQILh7XHry92xF9jOssdjln+O3e9yrbs/LKVvRfvSbZyNFQRkxoiIrKq2LqONQT8TlGp3CHQXUxqJKRp9VoiInJOQghkaBg9dT23EE/M34MVR7jelNTYp4aIiMgCPlx9Bj/tvggASJ49ULX90w1ncTD5Fg4m38KjbcPlCs8hMakhIiKS2Fv/dwx/HtJcE5OTX2LlaJwHm5+IiIgkVjmh2V1hxBe7KlgOkxoiIpJVDV8PuUOwuCEL92N34g0UlpSqJTVFJUr5gnJATGqIiEhWfl73ekI0Ca0uYyTSuHm7SOP2IQv3Y/C8PVBUWFRn3rZEa4XlFJjUSIg1ikRE5nkoJkzuEMzW7ZOtWvedSc1BQfG9IeCbz6RbIySnwaRGQlzWjIjIOE/EhkPhYJ1MbuuZt+Z2kfaOwsevZuFaVr7UITkNJjVERGR1/4zrgpe6R2O6Ey4pcOxKtup7UeHTcGJGHh7+3250nq29pod045BuIiKyulbhAWgVHiB3GLLIL9Zck3P8apZ1A3FArKmRkGNVoBIRWZ8Q+hvyfxzWzgqRWIe2uy3QkviQbkxqiIhIVsZ8IIwJ90efprXQvVGwxeKRw/JDVzDhz2Oqn1fGp8gYjf1iUkNERLIKNGKemrf7NwEALBnRHv97to2lQrKqH3ddxJv/d1xtWz4XyTQJkxoiIpLV3Cdj0D4q0KBmJV/Psq6gCoUCD7bSP/z7w0G23RH51u0ifLj6tMZ9q49fQ6dZW/DJ+rPILSi2cmT2iUkNERHJqm5NXywf3Rl9mtaCpi41q1/tqvq+8ujvta9103lud1fbfsylaVjFGwBOpmRj3O/xSM0uwHfbk1RNU4UlpXhywV60/XATrt66Y81Q7YJtv9pEROT0WtTxV32vqNQDp1mYn7XDsYoVlfrUbDqdjtTsfMzbmogDF2/i5u0izPjnlEzR2S4O6ZaQo00gRUQkp+ZhfnigZW21bc78Nvvm8uPYVWFhzNuFZf1u0rILUMPXAx5urKfgb4CIiGxGxdanNa91w9heDfQe07FeDQBAw5BqOsv1a17LnNBkVzGhKXcmNQf3zdqCh77ZBaBsSHxmXqHq+1d+O4yxvx8BAJSUKjH2tyP4cddF6wVtZUxqJNC/eSgAYFS3aJkjISJybJpqauY/F4sPBzXH8tGd4KGjD03riEALRmZ9ey9kYsBXOwEACem5AIB6U9Yi9qPN+OG/C7ieV4i1J9Kw5ngqMnIL8N7fJ7HmRKrWjsmOgEmNBL4d0hZ7JvdG/xahcodCROSQujSoiUa1qqFxraqreAf6emBopygE+Hjgl5EdVNtbRwSolRtYqSnL0SSk5aq+/3jtGXT4eIvq59nrzmLZwStyhGVV7FMjARcXBcICvOUOg4jI7kXW8NG4/deRHSFE2futLhWbr/4a0xl/Hb6q+tnN1bE75PT78j+t+1YcUe94nHWnCAE+hs8PVFF+USk83Fzgque1kIPd1NScO3cOgwYNQlBQEPz8/NC1a1ds27ZN7rCIiEhCD8eE4Y24Rvj9xY5q2xUKhd6EpjJXF4XadMUuztzLuJLWH2wy6bisO0VoOm09Hrzbh8fW2E1S8+CDD6KkpARbt27F4cOHERMTgwcffBBpaWlyh0ZERBJxcVFgfFxDdG4QZNLxoX5e2s/NnEbNK78dRmGJ+szFBcWlGPbTASzZk6zlmLJOx2dScywdnknsovnpxo0bOH/+PH788Ue0atUKADB79mx8++23OHnyJEJDNfdlKSwsRGFhoernnBzbfBGIiEgaUUG++PqZNgi6u/RCxTyG026oW3siDXUCErAyPgUPtgrDjIebY/mhK9hx7jp2nLuOPw5ewQMtQ/HcfXXx066LKCxVYk9Spup4IYTN/U7toqamZs2aaNy4MX7++Wfcvn0bJSUlWLBgAUJCQhAbG6v1uFmzZsHf31/1FRERYcWoiYhIDg/HhKlqerzcXVXbvdzt4pFnVT/svIgbeUVYvCcZ2XeKcSIlW7XvdGoOPtt4Dq0/2ISvtyZiwY4Lase2/3gz4i/fAgCcTctBj0+3YdVReRfitItXWKFQYPPmzYiPj0f16tXh5eWFuXPnYv369QgM1D5Eb8qUKcjOzlZ9Xbni+D2/iYjonupebhW+d0dtf+3NU85u1C+H8Oehq/oL3nUjrwiPfLsHJ1Oy8cpvR3Ap8w7GLzuKklKlBaPUTdakZvLkyVAoFDq/zp49CyEExo4di5CQEOzcuRMHDhzA4MGD8dBDDyE1NVXr+T09PeHn56f2RUREzqOmr6faz0tH3QcACPRxx8n3+8kRks06cPGmScc9+M0uXLh+W/Vz5u0iqUIymkIITcuHWcf169eRmZmps0x0dDR27tyJvn374tatW2qJScOGDTFy5EhMnjzZoOvl5OTA398f2dnZTHCIiJzEF5vOITzQG0+0K+uCcC0rHwE+7vDxcEPU5DUyR+d4Nk/ojgYhVecTMoehz29ZOwoHBwcjODhYb7k7d8pWInVxUa9YcnFxgVIpXzUXERHZvjfub6T2c8V5xd7s1xifbkjQe47XejfA11sTVT+vea0rBn5tm8Oa5fbFpvOYN6StLNe2iz41nTp1QmBgIIYNG4Zjx47h3LlzePPNN3Hx4kUMHDhQ7vCIiMhORdX0VX2vqznq5R71Vd+/3b8Jmof5ay3r7OoF+eovZCF2kdQEBQVh/fr1yMvLQ+/evdGuXTvs2rULq1atQkxMjNzhERGRnRrQIhSv9WmIn4a3QzVP7Y0XvhX2Ke/22ogOlu/hbcs2npZv/ji7mKcGANq1a4cNGzbIHQYRETkQFxcFJlRqntKmupcbcgtK0LNxWbeJFWM6q83MGx7ojau38i0Spz1RQL65a+yipoaIiMgaPhjUHADwYKuqi1/umdwbWyb2UDU9VV476fcX76tyTDVPNzQJ1d1pduh9dU0N1yZNeaCJbNe2m5oaIiIiS3u+UxSG3lcXCoUCq4+XjYyq5Vc2LLy6lzuqe7lrPTaypg9WvNIZV27eQU1fT+QVlqB9VCCKSpXoNGsrAGDrxB5YuOsiHmtbByuOpOC5++qiaW0/jO5ZH11mb1U7n6uLAqVKgVmPtkSgjztG/3pEb/xTH2yGD1efNvX2JdGzcYhs15Z1SLe1cUg3EREZqny4d92aPtjxZi+NZd744yhWxqdgSMdIfPxIS2uG51QMfX4zqSEiItJgzfFUzNlwFvOebYsWdbSPdsrMK0TNap5a95P57GKeGiIiIls1sFVtDNTQt6YyJjS2gx2FiYiIyCEwqSEiIiKHwKSGiIiIHAKTGiIiInIITGqIiIjIITCpISIiIofApIaIiIgcApMaIiIicghMaoiIiMghMKkhIiIih8CkhoiIiBwCkxoiIiJyCExqiIiIyCEwqSEiIiKH4CZ3ANYkhAAA5OTkyBwJERERGar8uV3+HNfGqZKa3NxcAEBERITMkRAREZGxcnNz4e/vr3W/QuhLexyIUqnEtWvXUL16dSgUCsnOm5OTg4iICFy5cgV+fn6SndeWOPo98v7sn6Pfo6PfH+D498j7M50QArm5uQgLC4OLi/aeM05VU+Pi4oLw8HCLnd/Pz88h/6NW5Oj3yPuzf45+j45+f4Dj3yPvzzS6amjKsaMwEREROQQmNUREROQQmNRIwNPTE9OnT4enp6fcoViMo98j78/+Ofo9Ovr9AY5/j7w/y3OqjsJERETkuFhTQ0RERA6BSQ0RERE5BCY1RERE5BCY1BAREZFDYFIjgXnz5iEqKgpeXl7o2LEjDhw4IHdIVcyaNQvt27dH9erVERISgsGDByMhIUGtTM+ePaFQKNS+Ro8erVbm8uXLGDhwIHx8fBASEoI333wTJSUlamW2b9+Otm3bwtPTEw0aNMDixYstfXsAgBkzZlSJv0mTJqr9BQUFGDt2LGrWrIlq1arhscceQ3p6uto5bPn+oqKiqtyfQqHA2LFjAdjf6/fff//hoYceQlhYGBQKBf7++2+1/UIITJs2DbVr14a3tzfi4uJw/vx5tTI3b97EkCFD4Ofnh4CAAIwcORJ5eXlqZY4fP45u3brBy8sLERERmDNnTpVYli9fjiZNmsDLywstW7bE2rVrLX6PxcXFePvtt9GyZUv4+voiLCwMzz//PK5du6Z2Dk2v++zZs23iHvW9hsOHD68Se//+/dXK2PNrCEDj36RCocCnn36qKmOrr6EhzwVrvm9K8iwVZJZly5YJDw8P8dNPP4lTp06JUaNGiYCAAJGeni53aGr69esnFi1aJE6ePCmOHj0qHnjgAREZGSny8vJUZXr06CFGjRolUlNTVV/Z2dmq/SUlJaJFixYiLi5OxMfHi7Vr14qgoCAxZcoUVZkLFy4IHx8fMWHCBHH69GnxzTffCFdXV7F+/XqL3+P06dNF8+bN1eK/fv26av/o0aNFRESE2LJlizh06JC47777ROfOne3m/jIyMtTubdOmTQKA2LZtmxDC/l6/tWvXinfffVesWLFCABArV65U2z979mzh7+8v/v77b3Hs2DHx8MMPi3r16on8/HxVmf79+4uYmBixb98+sXPnTtGgQQPxzDPPqPZnZ2eLWrVqiSFDhoiTJ0+KpUuXCm9vb7FgwQJVmd27dwtXV1cxZ84ccfr0afHee+8Jd3d3ceLECYveY1ZWloiLixN//PGHOHv2rNi7d6/o0KGDiI2NVTtH3bp1xQcffKD2ulb8u5XzHvW9hsOGDRP9+/dXi/3mzZtqZez5NRRCqN1bamqq+Omnn4RCoRBJSUmqMrb6GhryXLDW+6ZUz1ImNWbq0KGDGDt2rOrn0tJSERYWJmbNmiVjVPplZGQIAGLHjh2qbT169BDjx4/XeszatWuFi4uLSEtLU2377rvvhJ+fnygsLBRCCPHWW2+J5s2bqx331FNPiX79+kl7AxpMnz5dxMTEaNyXlZUl3N3dxfLly1Xbzpw5IwCIvXv3CiFs//4qGz9+vKhfv75QKpVCCPt+/So/LJRKpQgNDRWffvqpaltWVpbw9PQUS5cuFUIIcfr0aQFAHDx4UFVm3bp1QqFQiJSUFCGEEN9++60IDAxU3Z8QQrz99tuicePGqp+ffPJJMXDgQLV4OnbsKF5++WWL3qMmBw4cEADEpUuXVNvq1q0rvvjiC63H2Mo9aktqBg0apPUYR3wNBw0aJHr37q22zV5ew8rPBWu+b0r1LGXzkxmKiopw+PBhxMXFqba5uLggLi4Oe/fulTEy/bKzswEANWrUUNv+22+/ISgoCC1atMCUKVNw584d1b69e/eiZcuWqFWrlmpbv379kJOTg1OnTqnKVPx9lJex1u/j/PnzCAsLQ3R0NIYMGYLLly8DAA4fPozi4mK12Jo0aYLIyEhVbPZwf+WKiorw66+/4oUXXlBbnNXeX79yFy9eRFpamlos/v7+6Nixo9rrFRAQgHbt2qnKxMXFwcXFBfv371eV6d69Ozw8PFRl+vXrh4SEBNy6dUtVxhbuGSj7u1QoFAgICFDbPnv2bNSsWRNt2rTBp59+qla1b+v3uH37doSEhKBx48YYM2YMMjMz1WJ3pNcwPT0da9aswciRI6vss4fXsPJzwVrvm1I+S51qQUup3bhxA6WlpWovJgDUqlULZ8+elSkq/ZRKJV5//XV06dIFLVq0UG1/9tlnUbduXYSFheH48eN4++23kZCQgBUrVgAA0tLSNN5r+T5dZXJycpCfnw9vb2+L3VfHjh2xePFiNG7cGKmpqXj//ffRrVs3nDx5EmlpafDw8KjysKhVq5be2Mv36Spjjfur6O+//0ZWVhaGDx+u2mbvr19F5fFoiqVirCEhIWr73dzcUKNGDbUy9erVq3KO8n2BgYFa77n8HNZSUFCAt99+G88884zaYoCvvfYa2rZtixo1amDPnj2YMmUKUlNTMXfuXNV92Oo99u/fH48++ijq1auHpKQkvPPOOxgwYAD27t0LV1dXh3sNlyxZgurVq+PRRx9V224Pr6Gm54K13jdv3bol2bOUSY0TGjt2LE6ePIldu3apbX/ppZdU37ds2RK1a9dGnz59kJSUhPr161s7TKMNGDBA9X2rVq3QsWNH1K1bF3/++afVHsbW8uOPP2LAgAEICwtTbbP318+ZFRcX48knn4QQAt99953avgkTJqi+b9WqFTw8PPDyyy9j1qxZNj/d/tNPP636vmXLlmjVqhXq16+P7du3o0+fPjJGZhk//fQThgwZAi8vL7Xt9vAaansu2Bs2P5khKCgIrq6uVXqCp6enIzQ0VKaodBs3bhxWr16Nbdu2ITw8XGfZjh07AgASExMBAKGhoRrvtXyfrjJ+fn5WTywCAgLQqFEjJCYmIjQ0FEVFRcjKyqoSm77Yy/fpKmPN+7t06RI2b96MF198UWc5e379yuPR9bcVGhqKjIwMtf0lJSW4efOmJK+ptf6GyxOaS5cuYdOmTWq1NJp07NgRJSUlSE5OBmAf91guOjoaQUFBav8nHeE1BICdO3ciISFB798lYHuvobbngrXeN6V8ljKpMYOHhwdiY2OxZcsW1TalUoktW7agU6dOMkZWlRAC48aNw8qVK7F169YqVZ2aHD16FABQu3ZtAECnTp1w4sQJtTeh8jfhZs2aqcpU/H2Ul5Hj95GXl4ekpCTUrl0bsbGxcHd3V4stISEBly9fVsVmL/e3aNEihISEYODAgTrL2fPrV69ePYSGhqrFkpOTg/3796u9XllZWTh8+LCqzNatW6FUKlUJXadOnfDff/+huLhYVWbTpk1o3LgxAgMDVWXkuufyhOb8+fPYvHkzatasqfeYo0ePwsXFRdVsY+v3WNHVq1eRmZmp9n/S3l/Dcj/++CNiY2MRExOjt6ytvIb6ngvWet+U9FlqVLdiqmLZsmXC09NTLF68WJw+fVq89NJLIiAgQK0nuC0YM2aM8Pf3F9u3b1cbVnjnzh0hhBCJiYnigw8+EIcOHRIXL14Uq1atEtHR0aJ79+6qc5QP3evbt684evSoWL9+vQgODtY4dO/NN98UZ86cEfPmzbPakOeJEyeK7du3i4sXL4rdu3eLuLg4ERQUJDIyMoQQZUMTIyMjxdatW8WhQ4dEp06dRKdOnezm/oQoGxEQGRkp3n77bbXt9vj65ebmivj4eBEfHy8AiLlz54r4+HjVyJ/Zs2eLgIAAsWrVKnH8+HExaNAgjUO627RpI/bv3y927dolGjZsqDYcOCsrS9SqVUsMHTpUnDx5Uixbtkz4+PhUGSrr5uYmPvvsM3HmzBkxffp0yYYD67rHoqIi8fDDD4vw8HBx9OhRtb/L8lEje/bsEV988YU4evSoSEpKEr/++qsIDg4Wzz//vE3co677y83NFZMmTRJ79+4VFy9eFJs3bxZt27YVDRs2FAUFBapz2PNrWC47O1v4+PiI7777rsrxtvwa6nsuCGG9902pnqVMaiTwzTffiMjISOHh4SE6dOgg9u3bJ3dIVQDQ+LVo0SIhhBCXL18W3bt3FzVq1BCenp6iQYMG4s0331Sb50QIIZKTk8WAAQOEt7e3CAoKEhMnThTFxcVqZbZt2yZat24tPDw8RHR0tOoalvbUU0+J2rVrCw8PD1GnTh3x1FNPicTERNX+/Px88corr4jAwEDh4+MjHnnkEZGamqp2Dlu+PyGE2LBhgwAgEhIS1Lbb4+u3bds2jf8nhw0bJoQoG9Y9depUUatWLeHp6Sn69OlT5b4zMzPFM888I6pVqyb8/PzEiBEjRG5urlqZY8eOia5duwpPT09Rp04dMXv27Cqx/Pnnn6JRo0bCw8NDNG/eXKxZs8bi93jx4kWtf5flcw8dPnxYdOzYUfj7+wsvLy/RtGlTMXPmTLWkQM571HV/d+7cEX379hXBwcHC3d1d1K1bV4waNarKQ8qeX8NyCxYsEN7e3iIrK6vK8bb8Gup7Lghh3fdNKZ6lirs3RkRERGTX2KeGiIiIHAKTGiIiInIITGqIiIjIITCpISIiIofApIaIiIgcApMaIiIicghMaoiIiMghMKkhIiIih8CkhoicnkKhwN9//y13GERkJiY1RE5g+PDhUCgUGD16dJV9Y8eOhUKhwPDhwyW95owZM9C6dWvJzrdt2zY88MADqFmzJnx8fNCsWTNMnDgRKSkpkl3DEpKTk6FQKFQLjBKR5TCpIXISERERWLZsGfLz81XbCgoK8PvvvyMyMlLGyPRbsGAB4uLiEBoair/++gunT5/G/PnzkZ2djc8//9zk8xYVFUkYJRHJjUkNkZNo27YtIiIisGLFCtW2FStWIDIyEm3atFErW1hYiNdeew0hISHw8vJC165dcfDgQdX+7du3Q6FQYMuWLWjXrh18fHzQuXNnJCQkAAAWL16M999/H8eOHYNCoYBCocDixYsBAFlZWXjxxRcRHBwMPz8/9O7dG8eOHdMa99WrV/Haa6/htddew08//YSePXsiKioK3bt3x8KFCzFt2jQAQGZmJp555hnUqVMHPj4+aNmyJZYuXap2rp49e2LcuHF4/fXXERQUhH79+mm85okTJ9C7d294e3ujZs2aeOmll5CXl6c1xlu3bmHIkCEIDg6Gt7c3GjZsiEWLFgEA6tWrBwBo06YNFAoFevbsqTpu4cKFaNq0Kby8vNCkSRN8++23qn3lNTzLli1D586d4eXlhRYtWmDHjh1a4yBydkxqiJzICy+8oHrYAsBPP/2EESNGVCn31ltv4a+//sKSJUtw5MgRNGjQAP369cPNmzfVyr377rv4/PPPcejQIbi5ueGFF14AADz11FOYOHEimjdvjtTUVKSmpuKpp54CADzxxBPIyMjAunXrcPjwYbRt2xZ9+vSpcu5yy5cvR1FREd566y2N+wMCAgCU1TrFxsZizZo1OHnyJF566SUMHToUBw4cUCu/ZMkSeHh4YPfu3Zg/f36V892+fRv9+vVDYGAgDh48iOXLl2Pz5s0YN26clt8qMHXqVJw+fRrr1q3DmTNn8N133yEoKAgAVNffvHkzUlNTVUnlb7/9hmnTpuHjjz/GmTNnMHPmTEydOhVLlixRO/ebb76JiRMnIj4+Hp06dcJDDz2EzMxMrbEQOTWj1/UmIrszbNgwMWjQIJGRkSE8PT1FcnKySE5OFl5eXuL69eti0KBBYtiwYUIIIfLy8oS7u7v47bffVMcXFRWJsLAwMWfOHCGEENu2bRMAxObNm1Vl1qxZIwCI/Px8IYQQ06dPFzExMWpx7Ny5U/j5+YmCggK17fXr1xcLFizQGPuYMWOEn5+fSfc9cOBAMXHiRNXPPXr0EG3atKlSDoBYuXKlEEKI77//XgQGBoq8vDzV/jVr1ggXFxeRlpam8ToPPfSQGDFihMZ9Fy9eFABEfHy82vb69euL33//XW3bhx9+KDp16qR23OzZs1X7i4uLRXh4uPjkk0+03zSRE3OTNaMiIqsKDg7GwIEDsXjxYgghMHDgQFWNQrmkpCQUFxejS5cuqm3u7u7o0KEDzpw5o1a2VatWqu9r164NAMjIyNDaR+fYsWPIy8tDzZo11bbn5+cjKSlJ4zFCCCgUCr33VlpaipkzZ+LPP/9ESkoKioqKUFhYCB8fH7VysbGxOs9z5swZxMTEwNfXV7WtS5cuUCqVSEhIQK1ataocM2bMGDz22GM4cuQI+vbti8GDB6Nz585ar3H79m0kJSVh5MiRGDVqlGp7SUkJ/P391cp26tRJ9b2bmxvatWtX5XUgojJMaoiczAsvvKBqSpk3b55Z53J3d1d9X554KJVKreXz8vJQu3ZtbN++vcq+8makyho1aoTs7GykpqaqEidNPv30U3z11Vf48ssv0bJlS/j6+uL111+v0hm4YrIilQEDBuDSpUtYu3YtNm3ahD59+mDs2LH47LPPNJYv75/zww8/oGPHjmr7XF1dJY+PyFmwTw2Rk+nfvz+KiopQXFyssaNs/fr1VX1OyhUXF+PgwYNo1qyZwdfx8PBAaWmp2ra2bdsiLS0Nbm5uaNCggdpX5Rqjco8//jg8PDwwZ84cjfuzsrIAALt378agQYPw3HPPISYmBtHR0Th37pzB8ZZr2rQpjh07htu3b6u27d69Gy4uLmjcuLHW44KDgzFs2DD8+uuv+PLLL/H9998DKPs9AFD7XdSqVQthYWG4cOFCld9Decficvv27VN9X1JSgsOHD6Np06ZG3xeRM2BNDZGTcXV1VTVfaKoV8PX1xZgxY/Dmm2+iRo0aiIyMxJw5c3Dnzh2MHDnS4OtERUXh4sWLOHr0KMLDw1G9enXExcWhU6dOGDx4MObMmYNGjRrh2rVrWLNmDR555BG0a9euynkiIiLwxRdfYNy4ccjJycHzzz+PqKgoXL16FT///DOqVauGzz//HA0bNsT//d//Yc+ePQgMDMTcuXORnp5uVCIGAEOGDMH06dMxbNgwzJgxA9evX8err76KoUOHamx6AoBp06YhNjYWzZs3R2FhIVavXq1KPEJCQuDt7Y3169cjPDwcXl5e8Pf3x/vvv4/XXnsN/v7+6N+/PwoLC3Ho0CHcunULEyZMUJ173rx5aNiwIZo2bYovvvgCt27dUnXIJiJ1rKkhckJ+fn7w8/PTun/27Nl47LHHMHToULRt2xaJiYnYsGEDAgMDDb7GY489hv79+6NXr14IDg7G0qVLoVAosHbtWnTv3h0jRoxAo0aN8PTTT+PSpUtaEwYAeOWVV7Bx40akpKTgkUceQZMmTfDiiy/Cz88PkyZNAgC89957aNu2Lfr164eePXsiNDQUgwcPNjjecj4+PtiwYQNu3ryJ9u3b4/HHH0efPn3wv//9T+sxHh4emDJlClq1aoXu3bvD1dUVy5YtA1DWD+brr7/GggULEBYWhkGDBgEAXnzxRSxcuBCLFi1Cy5Yt0aNHDyxevLhKTc3s2bMxe/ZsxMTEYNeuXfjnn3+01moROTuFEELIHQQREalLTk5GvXr1EB8fL+nMzESOjDU1RERE5BCY1BAREZFDYPMTEREROQTW1BAREZFDYFJDREREDoFJDRERETkEJjVERETkEJjUEBERkUNgUkNEREQOgUkNEREROQQmNUREROQQ/h97ptI5CKBdOwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "plt.plot(range(len(energies)), energies)\n",
    "plt.xlabel('Monte Carlo step')\n",
    "plt.ylabel('energy')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "アニーリングが進むに連れ徐々にエネルギーが低くなっているのが分かります。\n",
    "\n",
    "このようにcallbackは、Algorithmの動作中にシステムの様子を知りたい時に有用です。\n",
    "\n",
    "## 結果の取得 -Result-\n",
    "\n",
    "`result.get_solution`で計算結果であるスピン列を取得できます。古典イジング模型の場合は直接`mysystem.spin`を参照することで、スピン列を取得も可能です。しかし、`result.get_solution`はそれ以外のシステムについてもスピン列を取得できる便利なメソッドです。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1, 1, -1, 1, -1, -1, 1, -1, 1, 1, 1, -1, 1, -1, 1, 1, -1, 1, -1, -1, -1, -1, -1, 1, 1, -1, 1, -1, 1, -1, -1, -1, 1, -1, 1, -1, 1, 1, -1, 1, -1, -1, -1, -1, 1, -1, 1, -1, 1, -1, 1, 1, -1, 1, 1, -1, -1, 1, 1, 1, 1, -1, 1, -1, -1, -1, 1, 1, -1, 1, -1, -1, 1, 1, -1, 1, -1, -1, -1, -1, -1, -1, -1, 1, -1, -1, 1, -1, -1, -1, -1, -1, -1, -1, -1, 1, -1, -1, 1, -1]\n"
     ]
    }
   ],
   "source": [
    "import openjij.cxxjij.result as R\n",
    "print(R.get_solution(mysystem))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "このスピン列がアニーリングによって得られた答えです。ハミルトニアンの基底状態 (に近い状態)であることが期待されます。"
   ]
  }
 ],
 "metadata": {
  "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}