Solution 6.12.ipynb

solution-6-12.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": [],
      "authorship_tag": "ABX9TyND83Ii6Uo/MiDum0F6EbXv",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/gist/enakai00/f091e7713c8966c35bd88f56565fd0b5/solution-6-12.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "id": "UVT1PUKoBYCV"
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "\n",
        "plt.rcParams.update({'font.size': 12})"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "- $\\Phi_0 = 2.07\\times 10^{-15}$\n",
        "\n",
        "- $\\displaystyle\\varphi_0 = \\frac{\\Phi_0}{2\\pi}$\n",
        "\n",
        "- $\\alpha = 0.82$\n",
        "\n",
        "- $I_C = 570 \\times 10^{-9}$\n",
        "\n",
        "- $S = 10 \\times 10^{-6}$\n",
        "\n",
        "- $E_J = I_C \\varphi_0$ : p.55\n",
        "\n",
        "- $\\displaystyle\\delta\\varphi = \\frac{2\\pi\\phi}{\\Phi_0}$ : Relationship between the phase jump and branch flux.\n",
        "\n",
        "- $I=I_C\\sin(\\delta\\varphi)$ : Relationship between the circuit current and  critial current\n",
        "\n",
        "At the operation point $\\Phi = \\Phi_0 / 2$ and $\\phi_- = 0$\n",
        "\n",
        "- $\\displaystyle V = -\\alpha E_J \\cos\\left(\\pi + \\frac{2\\pi\\phi_+}{\\Phi_0}\\right) - 2E_J\\cos\\left(\\frac{2\\pi\\phi_+}{2\\Phi_0}\\right)$\n",
        "\n",
        "In the following code, `plus` stands for $\\phi_+/\\Phi_0$."
      ],
      "metadata": {
        "id": "UW-L39tEkT82"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "Phi0 = 2.07 * 1e-15\n",
        "Alpha = 0.82\n",
        "IC = 570 * 1e-9\n",
        "S = 10 * 1e-6\n",
        "EJ = IC * Phi0 / (2 * np.pi)\n",
        "\n",
        "def v(plus):\n",
        "    V = -Alpha * EJ * np.cos(np.pi+plus*2*np.pi)\\\n",
        "            - 2 * EJ * np.cos(plus*2*np.pi / 2)\n",
        "    return V"
      ],
      "metadata": {
        "id": "0y07g0oOAbng"
      },
      "execution_count": 2,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "fig = plt.figure(figsize=(5, 8))\n",
        "subplot = fig.add_subplot(1, 1, 1)\n",
        "\n",
        "xs = np.linspace(-1.2, 1.2, 100)\n",
        "ys = v(xs)\n",
        "bottom = np.min(ys)\n",
        "top = np.max(ys)\n",
        "barrier = v(0)\n",
        "min_phi = np.abs(xs[np.argmin(ys)])\n",
        "max_phi = np.abs(xs[np.argmax(ys)])\n",
        "subplot.plot(xs, ys)\n",
        "subplot.plot([-1.2, 1.2], [bottom, bottom], linestyle='--', color='gray')\n",
        "subplot.plot([-.8, .8], [barrier, barrier], linestyle='--', color='gray')\n",
        "subplot.plot([-1.2, 1.2], [top, top], linestyle='--', color='gray')\n",
        "subplot.plot([min_phi, min_phi], [bottom-(barrier-bottom)/10, barrier+(barrier-bottom)/10], label='Energy barrier')\n",
        "subplot.plot([max_phi, max_phi], [bottom-(top-bottom)/50, top+(top-bottom)/50], label='Total depth')\n",
        "subplot.set_xlabel('φ/Φ0')\n",
        "subplot.set_ylabel('V')\n",
        "subplot.legend()\n",
        "\n",
        "print('(1)')\n",
        "print(' Energy barrier : {}'.format(barrier - bottom))\n",
        "print(' Energy barrier / Total depth : {}'.format((barrier-bottom) / (top - bottom)))\n",
        "print()\n",
        "print('(2)')\n",
        "print(' Phase gap : {}'.format(min_phi*2*np.pi))\n",
        "print(' I = IC * {}'.format(np.sin(min_phi*2*np.pi)))"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 645
        },
        "id": "mtm8D4gyA1ah",
        "outputId": "219e2d16-0e6d-4129-858b-cbe7a040a7da"
      },
      "execution_count": 3,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "(1)\n",
            " Energy barrier : 4.6629284768074576e-23\n",
            " Energy barrier / Total depth : 0.05846219366682481\n",
            "\n",
            "(2)\n",
            " Phase gap : 1.9039955476301773\n",
            " I = IC * 0.9450008187146686\n"
          ]
        },
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 360x576 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAU8AAAHzCAYAAACgxZ/8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdeZRcZZ3/8fdTvVfve6eXdHcSspBAAjQIAQMKKAIRBAdFQFAWdYZRFP3pjPwkDoxnlJmjzjjjiKIYwAV/AgoCOuKwyZYACSHQhGydXpPe97Xq+f3RXZ1O0t3prlTVvVX1eZ3TR+xby7dSXZ96tvtcY61FRETmx+N0ASIi0UjhKSISBIWniEgQFJ4iIkFQeIqIBEHhKSISBIWniEgQYiY8jTE3G2M2G2OGjTH3hvixU4wx9xhj6owxvcaYLcaYD005frox5n+MMR3GmFZjzG+MMQtCWYOIuEvMhCfQBNwJ/DQMj50I1ANnA9nAbcCDxpiqieO5wN1AFVAJ9AI/C0MdIuISJtbOMDLG3AmUW2uvm/K7ixkP1irgLeCz1to3jvF53gC+aa397TTHTgaesdZmHstziIh7xVLLc1rGmJMYb41+BsgHfgT83hiTcgyPWQwsBbbPcJN1sxwTkRgQ8+EJ3AT8yFr7srXWZ639OTAMnB7MgxljkoAHgJ9ba2unOX4i8A3gK8dQs4i4XDyEZyVwqzGmK/ADVAClAMaYp40xdoaf56c+kDHGA9wHjAA3H/5ExpglwBPAF6y1z4X7hYmIcxKdLiAC6oF/ttb+83QHrbXnzOVBjDEGuAcoBi601o4edrwS+DNwh7X2vmOqWERcL2ZansaYRGNMKpAAJBhjUo0xicCPgc8aY95jxqUbYy4yxsx3MueHwApgvbV28LDnLgP+AvzAWvvfIXg5IuJyMTPbbozZANx+2K+/aa3dYIy5ALgDOA4YBJ4HPm2t7Z3jY1cCexkfKx2bcugz1toHjDG3AxuA/qn3s9ZmzP+ViEg0iJnwFBGJpJjptouIRJLCU0QkCDEx215QUGCrqqqcLkNEYsyrr77aZq0tnO5YTIRnVVUVmzdvdroMEYkxxpi6mY6p2y4iEgSFp4hIEBSeIiJBUHiKiARB4SkiEgSFp4hIEBSeIiJBUHiKiARB4SkiEgSFp4hIEBSeIiJBUHiKiARB4SkiEgSFp4hIEBSeIiJBUHiKiARB4SkiEgSFp4i4zrdf+TbffuXbTpcxq5i4DIeIxJbajlqnSzgqtTxFRIKg8BQRCYLCU0QkCApPEZEgKDxFRIIQt7Pt99577xG/W7lyJaeeeiqjo6M88MADRxxfs2YNa9asYWBggAcffHDy933DYwyN+lh7+mmsrTmJ7u5uHn744SPuf8YZZ7Bs2TLa2tp47LHHjji+bt06Fi1aREtLC08++eQRx88991wqKiqor6/nqaeeOuL4BRdcQElJCbt37+bZZ5894vjFF19MQUEB77zzDi+++OIRxz/ykY+QnZ3Nm2++yebNm484fsUVV+D1etmyZQtbtmw54vhVV11FUlISmzZtYvv27Uccv+666wB44YUX2LFjxyHHkpKSuOqqqwB45pln2LNnzyHHvV4vV1xxBQB//vOfaWhoOOR4VlYWl112GQBPPvkkLS0thxzPz89n/fr1ADz66KO0t7cfcrykpIQLLrgAgIceeoienp5DjpeXl3PeeecB8OCDDzIwMHDI8erqas4++2wAHnjgAUZHRw85vnTpUtauXQuE9m8voKamhlWrVkX93962bW/y3IsvU5/QRGKCZ/LfKlR/e6EUt+EZKo1dgzR0DmAt/OHhN0n7ax93XrTI6bJEos7PX9jLb/+8jdKxLgZLfXjG/PQxRkaKO2PKWGudruGY1dTU2OlaSuE0NOrjSw9u4fFtLaxfXcqnzqzi5d0d3PP8HrJSE/n935/l2jddxG2ef7eNa376Mu+pzuOjp1Twy31fY1dbPwN7b+LOj6ziipoKR+oyxrxqra2Z7pg+3UH64dO7eHxbC/944XJufO8ijDGcvDCXkxbm8Ikfv8Q/PrSN7398DcYYp0sVcbUDPUPc8uvXWVyYwU+vOxVvciJ/aE3ihNJs/CaXr/32DU6qyOG44kynSz2EJoyC0No7zI+f281FJyzgpnWLDwnI0xflc+sHlvH7rU384pV9DlYp4n4+v+Xzv3qdvuEx/uuqk/EmH2zPJSYY/uPK8d/d9cd3HKxyegrPIPzHX95leMzPlz+4bNrjnzt7MWctKeBfHq+lf3gswtWJRI/HtzXz0u4OvvnhlSydpmWZl57MZ89exJ/e2s+rdR0OVDgzhec87W3r5xcv7+PK0yqoLkif9jYej+FLH1hK7/AYD7/eGOEKRaLHxhf3sjDPy9+cMvOY5qfPqqYwM4V/eaIWN83RKDzn6d/+ZwdJCR4+f+5xs97upIocTijLZuOLe131hou4xfambjbt7eSTZ1Ti8cw8N+BNTuSW845j095Onnr7QAQrnJ3Ccx46+kd4YlszV5++kKLM1Flva4zhk2dUsmN/Hy/tdld3Q8QN7nuxjtQkz6ytzoAraioozU7lgZfrIlDZ3Cg85+Hxbc2M+S0fOal8Trdfv7qUXG8SG1/cG9a6RKJN18AIj2xp5CMnlZHtTTrq7ZMSPKxfXcpz77bR2T8SgQqPTuE5D7/f2sSSogxWLJjbkonUpASuOLWCP721n6auwTBXJxI9frO5gaFRP9ecXjXn+6xfXcqY3/L4m83hK2weFJ5z1Nw9yKa9HXx4dem81m5e/Z5KfH7LQ681HP3GInHiwc311FTmcnxp1pzvs7I0i0WF6Ty6tSmMlc2dwnOOHtvajLXw4dWl87pfRZ6X1eXZPFXrnoFuESfVtffz7oE+LjxhwbzuZ4zhw6tLeXlPBy3dQ2Gqbu4UnnP0+61NnFieTdUMy5Nm8/7lxWyp76KtbzgMlYlEl79MNCTOXVE07/uuX12KtfDYG863PhWec7CnrZ9tjd2sP3F+rc6Ac1cUYS08/U5riCsTiT5/qT3A4sJ0KvPn3xBZXJjBytIsV3TdFZ5z8OSb49ubXbx6ft2MgJWlWRRnpfC/6rpLnOsbHuPl3R2cu6I46MdYv7qUrQ3djk/CKjzn4IVdbSwrzmRBdlpQ9zfG8L5lRTy7o5WRMX+IqxOJHs+/28aIz8/7ls2/yx5w9tJCAF7Y1X6UW4aXwvMohsd8bNrbwdol+cf0OO9fXkTv8Bib92rBvMSvv9TuJzM1kZqq3KAfY1lxJvnpybywsy2Elc2fwvMoXt/XxdCon7WLC47pcc5cUkByokez7hK3/H7LX2pbOXtpIUkJwUePx2M4fXE+L+xqd/TUZ4XnUbywsw2Pgfcsyjumx0lPSeSMRfk89fb+EFUmEl22NXbT1jcc1Cz74c5cXEBLzxC72/pDUFlwFJ5H8cKudk4ozyEr9einkB3N2UsL2ds+4PhAt4gT/rprvJv93uMKj/mx1i4eH0ZzctxT4TmL/uExttR3Tb5RxyowzrO5rjMkjycSTV7d28mignQKMlKO+bEq872UZqc6Ou6p8JzFK3s7GPNbzjzG8c6A4xdk4U1O4FVNGkmc8fstr+7rPKaJoqmMMaxdUsCLu9vx+50Z91R4zuLFXe0kJ3g4pTI0b3higoc1FTlqeUrc2d3WR9fAKDWVxzZ3MNXaxfl0DYzydkvP0W8cBgrPWfx1ZxsnV+aQlpwQssesqczl7eYe+nR5Dokjm/eONxhOCVHLE5hcAfPCTmfGPRWeM+gdGuWt5h5OXxSa8c6AU6ry8Ft4fZ9anxI/Ntd1kpeezKIg9oaYSUl2KtUF6by8R+HpKtsaurEWTloYum9KgJMX5uAxB7+JReLB5r0dnFKZG/JLcZ9UkcOW+m5H1nsqPGewtaEbgNXl2SF93MzUJJaVZPGqxj0lTrT2DrO3fYCaEM0dTLW6Ioe2vmGaHdiiTuE5g631XVTme8nxJof8sWsqc3ltXydjPp3nLrEv0FAI1Uz7VKsrcgB4o6Er5I99NArPGWxt6GJ1eU5YHrumKpeBER+1Lb1heXwRN3m1roPkRA+rykLbiwNYsSCTpATDlvrukD/20Sg8p3GgZ4jm7qHJb7VQq6kaX66hTUIkHmyu62R1eTYpiaFbtRKQkpjAigVZbK1Xy9MVwjXeGVCWk0ZxVgpbHHjDRSJp1Odne2NPyCdep1pdnsO2xu6IL5ZXeE5ja30XCR7DytLwhCfAqtJstjc5s7hXJFLe3d/HiM8fli57wOqKHPqGx9jd1he255iOwnMaWxu6WFacGdLF8YdbWZrFrtY+Bkd8YXsOEadtbxrvxa2cx1Uy5yvQQ4z0uKfC8zDWWrbWd7G6InzflADHl2bjt1Dr0KllIpGwvakHb3IC1UFcr2iuFhVmkJGSGPFxT4XnYfa2D9AzNBa2mfaAVWXj38RvqusuMWx7UzfHL8jC4wnt4vipEjyGE8qyI75cSeF5mMAbEK6Z9oCynDSy05J4qynySyxEIsHvt7zV1BPWLnvA6ooc3mruYXgscsNgCs/DbK3vJjXJw3FFGWF9HmMMK0uzNGkkMauuY4D+EV9YJ14DVpdnM+qz1DZHbu20wvMwtS09LCvJIvEYrrEyVytLs6ht6WVUZxpJDApMFh0fgZZn4DkiOYeg8JzCWkttSy/LizMj8nyryrIZGfOz80Bkl1iIRMKbjT0kJRiWRuDzVJHrxZucwNtqeTqjtW+Yjv4Rli+ITHgGxoLUdZdYtL2pm6XFmSQnhj9mPJ7xkH4ngqc8KzynCIyXLC8JfzcDoLogg7SkhMnujUissDZyk0UBKxZkUtvSE7Ht6RSeUwTGS5aXRKblmeAxrFiQqZanxJz9PcO094+E9cyiwy0rzqRzYJTW3uGIPJ/Cc4rall6Ks1LITQ/9NnQzWVmazVtNPY5dxEokHN5sDP+ZRYdbvmD8ud6OUNdd4TlFbXNvxLrsAStLs+gbHqO+cyCizysSTm8392BM5IbA4GCPsbY5Mj05heeEUd/4rHekuuwBSyeeb8d+zbhL7Hhnfy8VuV7SUxIj9pw53mRKslIjNmmk8Jywp62fEZ8/YjPtAUsmFuO/e0AbI0vs2HmgL+wnmkxn+YJMddsjLbCre6S77VmpSSzITuVdtTwlRoz5/Oxu7ee4CK2Xnmp5SRY7D0TmxBOF54Ta5h4SPYbFhZH/tjyuOJMd+9XylNiwt32AEZ+fpcUOtDxLMhn1Wfa09Yf9uRSeE2pbellcmBGRBb2HW1qUwc4Dffg04y4x4N2JhkAkziw6XGDY7e0ITBopPCe809Ib8fHOgOOKMxge89OgGXeJATv292EMjvTiFhVkkJRgInJxRYUn0D04SmPXIMsiPNMeEBgb0oy7xIJ3D4zPtIfzSgwzSU70sLgwIyLLlRSeMLkxx9Iih8JzYlZS454SC97d3+fIeGfA0uJMdraGvyGi8AR2TfxDL3FgaQVAZmoSpdmp2l1Jot6oz8/utj5HZtoDFhdm0NA5yNBoeDdGVngCuw70kZzgoTw3zbEalmjGXWJAXXs/oz7ryBrPgMVF6VgLu1vDO+Ou8GS85VldkB6RDZBnohl3iQWBcXsnZtoDAj3IXWHuuis8gV2t/SwuCt/V/eZiaXEmw2N+6js04y7R610HZ9oDqvLTMUbhGXbDYz72dQw4+mYDLCkOnKapcU+JXjscnGkPSE1KoCLXG/Y5hLgPz7r2AXx+69hkUYBm3CUWvLu/19GZ9oDFhensitcxT2PMccaYIWPM/eF8nl0T305OtzwDM+7vKjwlSo36/Oxpc+ac9sMtKcpgd2t45xBcG57AfwKbwv0kgXGRRYXOjnkCLC7KYHcEzskVCYf6jgFGfdbxhgiMN4aGx/w0dQ2G7TlcGZ7GmI8DXcBT4X6unQf6KMtJw5scuX0HZ7KoIJ09rf0RuwaLSCgFNuNwS0MECOtiedeFpzEmC/gn4EuReL5drf2ueLMBqgvS6R0eo61vxOlSROZtMjwLnP88LZlo/e4K46SR68ITuAO4x1rbMNuNjDE3GWM2G2M2t7a2BvVE1lp2tfY5PlkUUD3xhu+OwKllIqG2q7WfXG8SOd7IXQNsJrnpyeSlJ4d1uZKrwtMYswY4D/ju0W5rrb3bWltjra0pLCwM6vlaeoYYGPG5YowGDn5jR2IvQpFQ29PWxyKXfJZgYsb9QPg+S84P9B3qHKAK2GeMAcgAEowxx1trTw71k+10yUx7QGlOGsmJHoWnRKU9bf2897jgGjLhsLgwgz+9tT9sj++qlidwN7AYWDPx89/AH4APhuPJAuMhbum2J3gMVfnesK9PEwm1/uEx9vcMU+2C8c6AJUUZdPSP0NEfnjkEV7U8rbUDwOT5icaYPmDIWhvcoOZR7GrtJys1kYIM58doAhYVZOhicBJ13DRZFLB4yhxCXnpeyB/fbS3PQ1hrN1hrrw7X4++eGKOZGCJwherCdPZ1DDAWgQtYiYTK7sllSu7oxcHBJVPhmjRydXiG2962AVd1M2B8udKoz9LQGb7FvSKhtqe1H2OgMt/rdCmTynLSSPQY9raHZ7OduA3PoVEfTd2DVOW7KzwXF2rGXaLPnrY+SrPTSE1ybkOQwyUmeFiY56WuPTyfpbgNz/qOAayFqgL3fFMCVBdMjNMoPCWK7G5zz8kmU1Xme9nTppZnSAVadm5reeZ6k8hOS9JCeYka1lr2tPa7arIooKognbr28JzyHLfhubfdneFpjGFRYbq67RI12vpG6B0ec938AYzPIQyM+GjtHQ75Y8dxeA6Mt/K8SU6XcoTqgvSwX39FJFQCvaRqF820B1Tmh28OIX7Ds61/8h/WbRYVpNPSM0T/8JjTpYgclRvXeAZUT3zG68Iw4x7X4enGbgYcXCunrrtEgz1t/SQneijNce7qszMpzUkl0WPYE4YZ97gMz/FlSkOuG+8MCKyVC8e3pUio7WnrZ2GelwSPe042CQgsV9qrbnto7Ju4QqXblikFBIYT6jrU8hT329cxQJWLFscfrqogPSwL5eMyPN26TCkgI2X8fPt9anmKy1lr2dcxwMI8d36WYLwnF47lSnEZnnUuXaY01fiZEQpPcbfWvmEGRnyuOi3zcOFarhSX4bmnzb3LlAIq89PDdlqZSKgEvuAXujg8q8K0XCkuw7Ou3b3LlAIW5nlp7hlieMzndCkiMwqEp5t7cYHa9oa4MRKX4enmZUoBlflerIX6Du2uJO61r70fjxnfwcitSnNSSUoI/e5KcReebl+mFBBoGe/TjLu4WF3HwOTlY9wqMcFDRW7olyu59xWHiduXKQVoradEg7r2AVdPFgWEY7lS3IVn79AoVfleFhW47zzcqfLTk0lPTlB4iqtFw/wBwPELskhPDu1eo666hlEknFKZx9NfeZ/TZRyVMYaFmnEXF+sZGqVzYJTKPPe3PL/8wWUhf8y4a3lGk8o8L3UdanmKOwVO4oiGbns4KDxdrDLfS0PHID5/6DdyFTlWk2s8XXx2UTgpPF2sMj+dEZ+flp4hp0sROUJg7wU3L5APJ4Wni03OuGtrOnGhurYBCjJSyEiJu6kTQOHpagsnBuI17iluVNfRH7fjnaDwdLXSnDSSEoyWK4kr7WsfiIqZ9nBReLpYgsdQnuvVWUbiOsNjPpp7huJ2vBMUnq6nrenEjeo7BrE2fpcpgcLT9RbmeanXmKe4TH1nYJmSwlNcqiIvjZ6hMboHR50uRWRSw8QXekWuwlNcKvDHqdanuEl95yApiR4KM1OcLsUxCk+Xq8hTeIr77GsfoDw3DWPcd8XMSFF4utxky7NT4SnuUd85MPnFHq8Uni6X7U0iMzVRO8qLq9R3DMT1eCcoPKNCRa5XLU9xje7BUXqGxqjIc++lNyJB4RkFKvLSNOYprlGvmXZA4RkVKnK9NHQOYq22phPnNUz0gjTmKa5XkedleMxPa++w06WITI6/q+UprhcYW9K4p7hBfecAmamJZHuTnC7FUQrPKLBwcq2nZtzFeZppH6fwjALlOstIXKS+czDuZ9pB4RkVUpMSKMxMUbddHGetpaFTLU9QeEaNitw0ddvFca19wwyN+uN+ph0UnlGjIk8L5cV5kzPt6rYrPKNFRa6X5u4hxnx+p0uRODa5xlPddoVntKjIS8PntzR36zLE4pzApGW5wlPhGS20r6e4QX3HIAUZKaQlJzhdiuMUnlFicl9PjXuKgxq6BjTeOUHhGSVKslPxGGjs1Iy7OKehc1Bd9gkKzyiRlOBhQXYaDQpPcYjfb2nqGqQ8Vy1PUHhGlbKcNBq6FJ7ijAO9w4z6LGU5Ck9QeEaV8tw0ddvFMYFlSmp5jlN4RpGy3DSauwcZ1VpPcUBgyEjhOU7hGUXKc9PwW2jRWk9xQOPEkFFZjiaMQOEZVQKznJo0Eic0dA5QkJGsNZ4TFJ5RJDBQ36hJI3FAQ+egJoumUHhGkQU5qRhzcOBeJJIatcbzEArPKJKSmEBRZoq67RJxfr+loWuQMk0WTVJ4RpnyXK+WK0nEtfUPMzLm10z7FArPKFOem0ZDl7rtEllapnQkhWeUKctJo7lrCJ9f13CXyAn0drRM6SCFZ5Qpz/Uy5rfs79FaT4mcQMtTY54HKTyjTOCPV5NGEkkNnQPkeJPISEl0uhTXUHhGmcCYU6PGPSWCGrWb0hEUnlEmsEi5QVfSlAhq6BykXOOdh1B4RpnUpAQKMlJ0lpFEjLWWxk6t8TycwjMKledqU2SJnI7+EQZHfeq2H0bhGYXKctN0iqZEzORMu85rP4TCMwqV56bR1DWEX2s9JQImt6JTy/MQCs8oVJaTxojPT1vfsNOlSBwILJDXhNGhFJ5RSFvTSSQ1dg2SkZJIVprWeE6l8IxCZbkKT4mcxq7xfTyNMU6X4ioKzyhUGmh5asZdIkDLlKan8IxCWalJZKYmquUpEdHYNUhpTqrTZbiOwjNKleWk0aTwlDDrGx6je3BUuylNQ+EZpcpytFBewq9Jy5RmpPCMUmW5aeq2S9gd3MdT3fbDKTyjVFlOGr1DY/QMjTpdisQwXat9ZgrPKBWYcde4p4RTY9cgSQmGoswUp0txHYVnlJpc66lxTwmjxs5BFmSn4fFojefhFJ5RqlxnGUkEaJnSzBSeUaogI4XkBI/CU8KqqWtQ450zUHhGKY/HsCAnVd12CZtRn5/9PUNapjQDhWcUK8vRciUJn5buIfxWy5RmovCMYqU6y0jCSMuUZqfwjGJlOWkc6B1mZMzvdCkSgxp1rfZZKTyjWFluGtZCc7danxJ6gZbngmx126fjqvA0xqQYY+4xxtQZY3qNMVuMMR9yui63KtPWdBJGjZ2DFGSkkJqU4HQpruSq8AQSgXrgbCAbuA140BhT5WBNrqUd5SWcmrq1j+dsXBWe1tp+a+0Ga+1ea63fWvsYsAc4xena3GjBxCyowlPCobFzUDPts3BVeB7OGFMMLAW2T3PsJmPMZmPM5tbW1sgX5wIpiQkUZqZoxl1Czlo7efkNmZ5rw9MYkwQ8APzcWlt7+HFr7d3W2hprbU1hYWHkC3SJ8U2Rh5wuQ2JMR/8Iw2N+hecsXBmexhgPcB8wAtzscDmupoXyEg6Bv6lSheeMXBeeZvwSffcAxcDl1lptWDmL0pxUGrsGsdY6XYrEkCaF51G5LjyBHwIrgPXWWjWpjqIsJ42RMT/t/SNOlyIxpHFiKKhcs+0zclV4GmMqgc8Aa4AWY0zfxM9VDpfmWroMsYRDY+cg3uQEstOSnC7FtRKdLmAqa20doF1X52HqjvKrK3IcrkZiRVPXIKU5aYyPosl0XNXylPkLdKs0aSSh1NStZUpHo/CMctlpSXiTExSeElKNnYOaLDoKhWeUM8ZoazoJqaFRH+39Izq76CgUnjFAC+UllAJfxDqvfXYKzxhQqoXyEkKTC+SzFZ6zUXjGgLKcVDr6Rxgc8TldisQALZCfG1ctVZLgBLpXTd2DLC7MmPF2fr+ftrY2urq68PkUtPEoNTWV8vJykpJmXr/Z2DWEx0CJNkGelcIzBgS6V42ds4dnQ0MDxhiqqqpISkrSGr44Y62lvb2dhoYGqqurZ7xdY+cgxVmpJCWoYzob/evEgKkL5WfT399PWVkZycnJCs44ZIwhPz+foaHZJxcDC+RldgrPGFCSnYrHHD08ATweveXxbC5fmlogPzf6JMWApAQPxVmpNGjGXY6R329p7hpSy3MOFJ4xQgvlnWOMYefOnWF7/H379pGRkRGRSb62vmFGfH4tkJ8DhWeMiIWF8lVVVaSlpZGRkTH5c/PN2gt74cKF9PX1kZAQ/qtYNmqB/Jxptj1GlOak8cSbzfj9Fo8neieDHn30Uc4777ywPsfY2BiJie7705+urmOtdb731w7yc6eWZ4woy0ll1Gc50DvsdClhce+993LWWWfx5S9/mdzcXKqrq3niiScmj3d3d3P99dezYMECysrKuO222ya7uffeey9nnnkmX/ziF8nPz2fDhg20t7ezfv16srKyOPXUU7nttts466yzAPi7v/s7br311kOe/8Mf/jDf/e53Z6zv8ccfZ9GiRRQUFPCVr3wFv98PwK5du3j/+99Pfn4+BQUFXHXVVXR1dU3er6qqim9/+9uceOKJpKens3PnTowx3HPPPSxcuJD3v//97N27F2MMY2NjQb3W+dAC+blz39evBGXqQvk5L25+4mvQsi2MVQElJ8CH/iUkD/Xyyy9z7bXX0tbWxt133831119PY2Mjxhiuu+46ioqK2LlzJ/39/Vx88cVUVFTwmc98ZvK+H//4x9m/fz+jo6N86lOfIj09nZaWFvbu3csHP/hBKisrAbj22mu59NJLueuuu/B4PLS1tfHnP/+ZH//4xzPW9vDDD7N582b6+vo477zzWLZsGTfccAPWWv7hH/6BdevW0dPTw+WXX86GDRv43ve+N3nfX/7yl/zhD3+goKCA/fv3A/DMM8/w9ttv4/F4Jn8XMN/XOh9NXUNkpiaSlapNkI9GLc8YMde1nm536aWXkpOTM/kzNbAqKyu58cYbSUhI4Nprr6W5uZn9+/ezf/9+Hn/8cb73ve+Rnp5OUVERX/ziF/nVr341ed/S0qE2rmcAACAASURBVFL+/u//nsTERJKTk/ntb3/LN7/5TbxeL8cffzzXXnvt5G1PO+00srOzeeqppwD41a9+xTnnnENxcfGMdX/1q18lLy+PhQsXcsstt/DLX/4SgCVLlnD++eeTkpJCYWEhX/rSl3jmmWcOue/nP/95KioqSEs72NrbsGED6enph/wOmPdrPfz+R9PQqWVKc6WWZ4wI6nIcIWoRhtIjjzwy45hnSUnJ5H97vV4A+vr66OjoYHR0lAULFkwe9/v9VFRUTP7/qf/d2trK2NjYjMdhvPV5//33c/7553P//ffzhS98Yda6p96/srKSpqYmYDzsvvCFL/Dcc8/R29uL3+8nNzd3xvvO9juAurq6eb3W+dIC+blTeMaIrNQkMlMTo77lGYyKigpSUlJoa2ubcXJk6uLwwsJCEhMTaWhoYOnSpQDU19cfcvurr76aVatWsXXrVt5++20uvfTSWWuor69n5cqVwPjSotLSUgD+8R//EWMM27ZtIy8vj0ceeeSIFQTTLVyfaTH7fF/rfDV1D3JKZe7RbyjqtseS8Wu4R/dypWAsWLCAD3zgA9x666309PTg9/vZtWvXEd3jgISEBC677DI2bNjAwMAAtbW1bNy48ZDblJeXc+qpp3LNNddw+eWXH7X7e9ddd9HZ2Ul9fT3f//73+djHPgZAb28vGRkZZGdn09jYyF133RXR1zof/cNjdA2MquU5RwrPGBIL+3quX7/+kHWeH/nIR+Z0v40bNzIyMsLxxx9Pbm4uH/3oR2lubp7x9j/4wQ/o7u6mpKSEa665hiuvvJKUlJRDbnPttdeybds2rrnmmqM+/yWXXMIpp5zCmjVruOiii7j++usBuP3223nttdfIzs7moosu4rLLLpvT65nNfF/rXB2cadcC+bkw1lqnazhmNTU1dvPmzU6X4bj/+8ib/H5rE1tv/8C0x99++21WrFgR4aqiw1e/+lVaWlr4+c9/Pvm7Z599lquvvpq6urqY2khlpr+Dp985wHU/28RvP3cGp1TmOVDZQZ968lMA/OyCnzlahzHmVWttzXTH1PKMIaU5aXQPjtI3POZ0Ka5XW1vLG2+8gbWWV155hXvuueeQVu7o6Cjf//73ueGGG2IqOGejBfLzM2t4GmMUrlEk0N2Kx0mj+ert7eWyyy4jPT2dj33sY9x6661ccsklwHjLLCcnh+bmZm655RaHK42cpq5BEj2Gokx12+fiaLPtjcaY+4CN1to3I1GQBG/qNdyXFmc6XI27nXrqqTNu5rFixQr6+/sjXJHzGjvHT7BIiOLTeyPpaC3LzwLVwCZjzGvGmC8YYwojUJcEIVYWyoszmrQV3bzMGp7W2t9Za/8GWAD8CPgboMEY83tjzOXGGJ3D5SJFmakkesz8FsqLTGjsGqRc4TlncxrTtNZ2WWt/ZK09C1gBbAa+Cxz7+ggJmQSPoSQ7VS1Pmbcxn5+WHrU852NeE0LGmGSgBngPUAyEeVcJma/SGNjXUyLvQO8wPr9VeM7DnMLTGHOWMeZuYD9wJ/ASsNRa+75wFifzVx4DC+Ul8pq0CfK8HW2p0gZjzE7g0YlfXWytXWqtvcNaWxf+8mS+SnPSaOkZYsznd7oUV5rPJTM2bNjA1VdfHZY6AvuTusXkDvI6u2jOjtbyfA9wG7DAWnuTtfavEahJjkFpTho+f/Rtijz1lEyPx3PI5TgeeOCBae/z9NNPU15eHuFK5+/wzYzdSAvk52/WdZ7W2g9FqhAJjbIpaz2j6YPQ19c3+d9VVVX85Cc/CfvlOOSgpq5Bcr1JeJO10dpc6QyiGFMWY2cZDQ8Pc8stt1BaWkppaSm33HILw8PD9Pf386EPfYimpqbJFmpTUxOvvPIKZ5xxBjk5OSxYsICbb76ZkZGROT3Xnj17OPvss8nMzOT888+nra3tkOMvvfQSa9euJScnh9WrV/P0009PHjvnnHP4h3/4B0477TSysrK45JJL6OjoAGDdunUA5OTkkJGRwYsvvjh5v5kuKxJpjZ3R9WXrBvqaiTGTmyLPITy//cq3qe2oDWs9y/OW89XTvhr0/f/5n/+Zl156iS1btmCM4ZJLLuHOO+/kjjvu4IknnuDqq6+moaFh8vbNzc1897vfpaamhoaGBj70oQ/xX//1X3M6zfITn/gEZ5xxBn/60594+eWXueiiiyZP2WxsbOSiiy7ivvvu44ILLuCpp57i8ssvp7a2lsLC8fNGNm7cyB//+Eeqq6v55Cc/yec//3nuv/9+nn32Waqrq+nq6prcg/Odd96Z9bIikdbUNcTCfG/EnzeaqeUZY7zJieR6k2JmofwDDzzAN77xDYqKiigsLOT222/nvvvum/H2p5xyCqeffjqJiYlUVVXxmc98Zk57Xe7bt49NmzZxxx13kJKSwrp161i/fv3k8fvvv58LL7yQCy+8EI/Hw/nnn09NTQ2PP/745G2uueYaVq1aRXp6OnfccQcPPvjgrNdan+myIpFmraWxS5ffmC+1PGPQ+FrPo4fnsbQII6WpqWnywmxw6CUuprNjxw6+9KUvsXnzZgYGBhgbG+OUU06Z0/Pk5uaSnp5+yHMFdpivq6vjN7/5DY8++ujk8dHRUd73voOr9Q6/FMfo6OgRXf+pZrqsSKT1DI3RNzym8JwntTxjUFkMrfUsLS2lru7gqripl7iYrnv7uc99juXLl/Puu+/S09PDt771LeayZ+2CBQvo7Ow8ZEOQffv2Tf53RUUF11xzDV1dXZM//f39fO1rX5u8zdRLeezbt4+kpCQKCgpcv6VdoJeiNZ7zo/CMQWW5aTR2Ds4pNNzuyiuv5M4776S1tZW2tjb+6Z/+aXLtZXFxMe3t7XR3d0/evre3l6ysLDIyMqitreWHP/zhnJ6nsrKSmpoabr/9dkZGRnj++ecPaWVeffXVPProo/zxj3/E5/MxNDTE008/fch46/33389bb73FwMAA3/jGN/joRz9KQkIChYWFeDwedu/eHaJ/ldA6uMZT4TkfCs8YVJaTRv+Ij+7B+V2z241uu+02ampqOPHEEznhhBM4+eSTue222wBYvnw5V155JYsWLSInJ4empib+9V//lV/84hdkZmZy4403Tl5LaC5+8Ytf8PLLL5OXl8c3v/lNPvnJT04eq6io4He/+x3f+ta3KCwspKKigrvuugu//+DJCNdccw3XXXcdJSUlDA0N8e///u/AeJf861//OmeeeSY5OTm89NJLIfrXCY3GzgFALc/50mU4YtCTbzbz2ftf47G/P4tVZdmTv9dlOMLnnHPO4eqrr+aGG25wupSjOvzv4J//8BYbX6yj9o4LXDPEoMtwiCPms1xJJDDT7pbgjBYKzxgUGLuKleVKEl6NnYPqsgdBS5ViUF56MqlJnpg5yygaTD3bKNo0dg2xYkGW02VEHbU8Y5AxJqaWK0n4DI36aOsb1kx7EBSeMaos1ztteE6dHZb4c/gEsfbxDJ7CM0aV5aQdMeaZnp5OY2MjIyMjMbEGVObHWkt7ezupqQf37NQaz+BpzDNGleWk0t4/wuCIj7TkBADKy8tpa2ujrq7O1XtLSvikpqYesgdq4AtWOyrNn8IzRk3d13NJUQYAHo+HoqIiioqKnCxNXKSxaxCPgZJs7SA/X+q2x6iynPGNJjTjLrNp7BqkJCuVpARFwXzpXyxGTW15isxEazyDp/CMUcWZKSR4jBbKy6y0j2fwFJ4xKjHBQ0lWqlqeMiOf39LSPaSWZ5AUnjFsuuVKIgH7e4YY81vNtAdJ4RnDynJ1lpHMTGs8j43CM4aV5aTR0jPEmE9nFcmRAr2ScnXbg6LwjGFluWnj41o9Q06XIi4UaHmq2x4chWcMC3THmroUnnKkxq5B8tKT8SbrXJlgKDxjWKA71jBxmQWRqRo6B9VlPwYKzxgW6I41aMZdptHQOaDwPAYKzxiWmpRAUWaKWp5yBGstjZ2DlOd6nS4laik8Y5yWK8l02vpGGB7za5nSMVB4xrjyXK+67XKEQG9E3fbgKTxjXHluGk1dg/j82vxYDmqYXOOpbnuwFJ4xrjw3jVGf5UCvlivJQYHw1HntwVN4xjhdhlim09g1QI43iYwUrfEMlsIzxgW6ZRr3lKm0xvPYKTxjnBbKy3QaOgcpz9F457FQeMa41KQECjJS1PKUSdZaLZAPAYVnHNBaT5mqo3+EoVG/JouOkcIzDpTnpqnlKZO0TCk0FJ5xoDx3fEd5v9Z6ClPDUy3PY6HwjAPluV5GfH5a+4adLkVcIDB5qG77sVF4xoFy7a4kUzR2DZKVmkhWapLTpUQ1hWcc0HIlmapBuymFhMIzDpTlquUpB2mZUmgoPOOANzmR/PRkhadMrPFUyzMUFJ5xoiw3Td12oXNglIERnyaLQkDhGScqtK+nAPUd41+gFQrPY6bwjBPleVrrKVA/0fuoyFO3/VgpPONExcRaz/3a1zOu1XeM9z4UnsdO4RknAh+WwIdH4lN95wC52sczJBSecSIwxhUY85L4VN8xoFZniCg840RZbhrGwD6FZ1xr6BykQsuUQkLhGSdSEhMoyUqdnDCQ+OPzj+/jqZZnaLguPI0xecaYh40x/caYOmPMJ5yuKVZU5Hpp0Jhn3NrfM8Soz1KRp2VKoeC68AT+ExgBioGrgB8aY1Y6W1JsKM9LU8szjh1c46mWZyi4KjyNMenA5cD/tdb2WWufB34PXONsZbGhItdLS88Qw2M+p0sRB9R3aplSKLkqPIGlwJi1dseU320F1PIMgYo8L9bqMsTxqr5jAGOgNCfV6VJigtvCMwPoOex33UDm4Tc0xtxkjNlsjNnc2toakeKi3eRyJYVnXKrvHKAkK5WUxASnS4kJbgvPPiDrsN9lAb2H39Bae7e1tsZaW1NYWBiR4qLdwvzAQnmNe8ajho5BddlDyG3huQNINMYcN+V3q4HtDtUTU4ozU0lO8GjSKE7Vdw5osiiEXBWe1tp+4CHgn4wx6caYM4FLgPucrSw2eDxmfGs6LVeKO8NjPlp6hrRMKYRcFZ4T/hZIAw4AvwQ+Z61VyzNEynO1XCkeNXYOYq2WKYWS63YHsNZ2AJc6XUesqsjzsm1bs9NlSIRpmVLoubHlKWFUkeula2CU3qFRp0uRCJpcIK9ue8goPONM4MOjreniS33nAMkJHooztcYzVBSecWbhRLdNuyvFl4aOQcpy0/B4jNOlxAyFZ5xZmKe1nvGorqN/8r2X0FB4xpkcbzLZaUnUdfQ7XYpEiLWWurYBqvIVnqGk8IxDlfle6trV8owXnQOj9A6PsTA/3elSYorCMw4tzFN4xpO69vFeRqW67SGl8IxDlfleGrsGGfX5nS5FIiAwOVipbntIKTzjUGVeOj6/palLy5XiQaCXoQXyoaXwjEOB3ZXUdY8Pde3jW9GlJmkrulBSeMahQPetTsuV4sK+jv7JL0wJHYVnHCrOTCU50cO+di1Xigd17QOaLAoDhWcc8ngMlZpxjwsDI2Mc6B2mqkDLlEJN4RmnKvO9OkUzDgTeY51dFHoKzzi1MC+dfR0DWGudLkXCKNC70DKl0FN4xqnKfC8DIz5a+4adLkXCaF8gPPPUbQ81hWecCsy+7tO4Z0yr6+gnOy2JbG+S06XEHIVnnArMvu5VeMa0uvYBddnDROEZp8pzvXgMWq4U4/Z1DGiyKEwUnnEqOdHDguw0LZSPYaM+Pw2dg2p5honCM45pa7rY1tQ1iM9vNVkUJgrPOFZVkM5eddtj1l4tUworhWccq85Pp2tglK6BEadLkTDY09oHQHWhWp7hoPCMY9UTp+ztaVPrMxbtbR8gIyWRwowUp0uJSQrPOBZokSg8Y9Putn6qC9IxRlfMDAeFZxyryPWS4DEKzxi1p61vsnchoafwjGPJiR7Kc9PYrfCMOcNjPho6B7WbUhgpPONcdUE6e1oVnrFmX/sA1sIihWfYKDzjXPXEciXtrhRbAkMx6raHj8Izzi0qSGdgxMeBXu2uFEsC4alue/goPONc4MO1W133mLKnrZ/89GSy07SbUrgoPOOc1nrGpsAyJQkfhWecK81OIznRw562PqdLkRDao/AMO4VnnPN4DNX56exp0wYhsaJveIzW3mGdlhlmCk8ZX66klmfM2DsxBKNlSuGl8BSqCsYvBjfm8ztdioTAbs20R4TCU1hUkM6oz9LYNeh0KRICgZMeqvIVnuGk8JTJsTGdphkb9rT1UZaTRmpSgtOlxDSFpxxcrqS1njFhT1s/VQXaADncFJ4yuZh6Z6smjaKdtZZdrf0sKcxwupSYp/AUjDEcV5TBzv0Kz2jX0jNE3/AYS4oUnuGm8BQAlhRlqOUZA96d+AJcUpTpcCWxT+EpwHh4dvSP0NGv6xlFs50HAuGplme4KTwFOPhhC3z4JDrtbO0jOy2Jgoxkp0uJeQpPAQ6G57sHeh2uRI7Fzv19HFeUoesWRYDCU4DxDULSkhLU8oxyO1v71GWPEIWnAOMbhCwpylB4RrH2vmE6+kcUnhGi8JRJCs/oFnjvFis8I0LhKZOWFGXQ3D2+TlCiT2Cp2XEKz4hQeMqkQHdvl1qfUWnngT7SkhIozU5zupS4oPCUSVquFN12HuhjcVE6Ho9m2iNB4SmTKvO8JCUY3lV4RqWdB/o4TmcWRYzCUyYlJnioLkhXyzMK9Q2P0dw9pJn2CFJ4yiGWFGWwS+e4R53AOPVi7aYUMQpPOcSSokzq2vsZGvU5XYrMw47942eGHVes8IwUhaccYnlJJn6rSaNo805LLymJHl16I4IUnnKI5SXjEw5vN/c4XInMR21LL8tKMknQTHvEKDzlEJX56aQmeaht0QYh0aS2pWfyi08iQ+Eph0jwGJYVZ1LbopZntGjtHaatb4TlJVlOlxJXFJ5yhOUlWbzd3Iu11ulSZA4CX3TLF6jlGUkKTznC8gWZdPSP0No37HQpMge1zeNDLGp5RpbCU44Q+BAGPpTibm+39FCclUJeunaPjySFpxwhMPGgcc/oUNvcq1anAxSecoTc9GRKslLV8owCoz4/Ow/0abzTAQpPmdbyBZm8reVKrrenrZ8Rn58VanlGnMJTprW8JIudB3oZ9fmdLkVmETiZQS3PyFN4yrRWLMhk1GfZ3drvdCkyi9qWXpISDIsKdE57pCk8ZVqTM+6aNHK12uYeFhdmkJyoj3Kk6V9cprWoMJ3kBA9v6Rx3V3u7uVenZTpE4SnTSkrwsKwkkzcbu50uRWZwoHeIlp4hVpVlO11KXFJ4yoxOKM/mjYZunabpUoEvthPLcxyuJD4pPGVGq8uz6R0aY2/7gNOlyDS21nfjMbCyVMuUnKDwlBmdUDbeonmjocvhSmQ62xq7WVKUQXpKotOlxCWFp8zouOIMUhI9bGvQuKfbWGt5o6F78gtOIk/hKTNKSvBwfGkWb2jSyHVaeoZo6xvmxHJNFjlF4SmzWl2ew5uN3fj8mjRyk631gckihadTFJ4yqxPKshkY8bFblyN2lW2NXSR6DCsWaLLIKQpPmVWgZfOGxj1d5Y2GbpYWZ5KalOB0KXFL4SmzWlSYgTc5gW0a93QNay3bGrvVZXeYwlNmleAxrCrLZquWK7lGfccgXQOjWhzvMIWnHNWJZdm81dSj7elc4o3G8S8ytTydpfCUo1pdkcPwmH9y70hx1uv7ukhJ9LC0WBuCOEnhKUdVU5ULwKa9nQ5XIgCb9nawpiJH29A5TP/6clQLstMoz01j054Op0uJe/3DY2xv6uHUqjynS4l7Ck+Zk9Oq8ti0t0M7LDnstX2d+PyWU6sVnk5zRXgaY1KMMfcYY+qMMb3GmC3GmA85XZccdGp1Hu39I+xp02U5nLRpbyceAycv1Ey701wRnkAiUA+cDWQDtwEPGmOqHKxJpgh0EzftVdfdSZv2dLBiQRaZqUlOlxL3XBGe1tp+a+0Ga+1ea63fWvsYsAc4xenaZNziwnTy05N5ZY8mjZwyMubn9fpOjXe6hCvC83DGmGJgKbDd6VpknDGGmqpctTwd9GZTN0Ojfk7TeKcruC48jTFJwAPAz621tbPc7iZjzGZjzObW1tbIFRjHTq3KY1/HAPt7hpwuJS5tnvjiCiwdE2dFJDyNMU8bY+wMP89PuZ0HuA8YAW6e7TGttXdba2ustTWFhYVhfgUCGvd02it7OqnK91KUmep0KUKEwtNae4611szwcxaAMcYA9wDFwOXW2tFI1CZzt7I0C29ygtZ7OsDvt2yu69B4p4u46eInPwRWAOdZawedLkaOlJjg4dSqPJ7b2eZ0KXFne1MPXQOjnL4o3+lSZIIrxjyNMZXAZ4A1QIsxpm/i5yqHS5PDnL20kN2t/dR36IqakfT0OwcAWLdUQ1Ru4YrwtNbWTXThU621GVN+HnC6NjnUOcvGP7yBD7NExjM7WjmhLJvCzBSnS5EJrghPiR7VBekszPPy9Dta4RAp3QOjvLavc/KLS9xB4SnzYozh7KWFvLCrneExn9PlRJcnvjb+M0/P7WzFb1F4uozCU+btnGWFDI762KSzjeanZdv4zzw9/U4rWamJrNbO8a7iptl2mca99957xO9WrlzJqaeeyujoKA88cOSw8Jo1a1izZg0DAwM8+OCDRxyvqalh1apVdHd38/DDDx9x/IwzzmDZsmW0tbXx2GOPHXH8tDPOJDnBw/++voOdfz3y+LnnnktFRQX19fU89dRTRxy/4IILKCkpYffu3Tz77LNHHL/44ospKCjgnXfe4cUXXzzi+Ec+8hGys7N588032bx58xHHr7jiCrxeL1u2bGHLli2M+vx0D44yMja+E37+6nM5oSIP2naxo/btI+5/3XXXAfDCCy+wY8eOQ44lJSVx1VXj85jPPPMMe/bsOeS41+vliiuuAODPf/4zDQ0NB193SwsJiQkE2o9PPvkkLS0th9w/Pz+f9evXA/Doo4/S3t7OSF0n671J3H9fMyUlJVxwwQUAPPTQQ/T0HLpBdXl5Oeeddx4ADz74IAMDh07sVVdXc/bZZwPwwAMPMDp66IrApUuXsnbtWmDuf3uBf694o/CUeUtNSuC06jw217VS4dKe5MiYnxd3tfFuQxcDI4cOL9zzx1p8JLAyuZUTvb2U5qSRmerOj0L/yBgjPj85adoIxG1MLOzPWFNTY6drgUj4/OS53dz5h7d57v+8j4o8r9PlTLLW8tvXGvnu/+ygsWuQ1RU5XLiqhDMW57OkKAODYXjMx2v7OnlhZzuPbGmkrW+E9x5XwNcvWsHykjBeB/1nF43/76f+MOe7/Of/7uSuP77DK18/N67OLPrUk58C4GcX/MzROowxr1pra6Y7pjFPCcp5K4oBeOLNZocrOahrYITP3f8aX/7NVgozU9j46dN45G/X8pmzF3NieQ7e5ETSkhPI8Sbz/uXF3Hbx8Tz7f97H1y9cwfamHj78g7/y8xf2umrD58e3NbO6PDuugjNaKDwlKFUF6awuz+Z3W5qcLgWA1/d1csH3nuOp2v18/cIVPPS5taxbWsj4Wb8z8yYncuO6Rfzpi+s4a0kBt/9+Ozf8fDO9Q86fHbzzQC/bm3q4ZE2Z06XINBSeErQPrylje1MPOw/0OlrHX2r384kfv0xyooeH//ZMbly3CI9n9tA8XEFGCvdcW8OG9cfz9I5WPn73S7T2Doep4rn53ZYmPAYuPnGBo3XI9BSeErT1Jy7AY3C09fn/Xm3gxo2vsqQog99+bi2ryoK/lrkxhuvOrOYn19awu7Wfy3/4AnsduuyItZbfbWli7eICirLUZXcjhacErSgrlbWLC/jdliZHxgkf3FTPl3+zlbWL8/nlTaeH7NTF9y0r4pc3nU7v0Cgfv/sl6tojH6Cv13exr2OAS9aURvy5ZW4UnnJMPrymlH0dA2yp74ro8z7yeiNffegN1i0t5CfX1pCREtqlRmsqcvjFjaczNObjEz9+mYbOyG6E8vstTSQnevjgqpKIPq/MncJTjskFq0pITvREtOv+xLZmvvTgFk6vzufua04hJTEhLM+zYkEW91//HnqGRvnEj1/mQIR20B/z+XnsjSbOW1FEli705loKTzkmWalJnH98MQ+91kD/8FjYn++FXW184VdbWFORwz3X1ZCaFJ7gDFhVls3GT59GW98w1/5sEz0RmIX/01v7aesb4VLNsruawlOO2Q1nVdMzNMavNtWH9XnebOzmpo2vUpnv5afXnYo3OTJnBZ20MJcfXn0K7+7v5aaNmxkaDd+GKNZafvTMLqoL0jl3Yi2tuJPCU47ZSQtzOa06j3ue282ozx+W56jvGOC6n20iMzWRn3/6NHK8yWF5npmcvbSQf/2b1by0u4Mv/noLPn94Jshe2t3B1oZubnzvIhLmudxKIkvhKSHx2bMX0dQ9xGNvhH7ss71vmE/+9BVGfX42fvo0SnPSQv4cc3HpSWXcdtEKnnizhW8+uj0sKwx+9OwuCjKSuexkddndTuEpIXHO0iKWFmfwo2d2hzRU+ofH+PS9m2jqGuSea2s4rjgzZI8djBveu4ib1i1i44t1/NfTu0L62G839/D0O61ct7Yq7GO5cuwUnhISHo/hpnWLqW3p5X/e2h+Sxxwe8/HZ+19lW2M3P/jEydS45MqRX7tgOZeuKeWuP77DL17eF7LH/cH/7sSbnMDVp1eG7DElfBSeEjIfXl3KcUUZbPj9dvqOceZ91Ofn5l+8znPvtvEvl53I+ce7Z/LE4zF856OrOWdZIV9/ZBuPvN54zI/5v7UH+MMbzdz43kURH8+V4Cg8JWSSEz18+6Mn0twzxHeerA36cXx+y5d/s5X/eWs/3/zwSq44tSKEVYZGcqKH/776FE6vzufW32zlyWPYXapveIyvP7yNpcUZ/O37FoewSgknhaeE1MkLc7lubRUbX6xj096Oed9/aNTH3z7wKr/b0sRXL1jOtWurQl9kiKQmJfCTa2tYXZ7N3/3idX6zzTNTfwAACNlJREFUObilWt95spbmniH+5fITw7bgX0JP4Skh9+UPLKM8N41bH9xKS/fcz8rpHhzlk/e8wp/e2s/t64/nc+e4vxWWnpLIxuvfw9rF+Xzl/73Bfz29c14TZk9sa2bji3V8am01Jy/MDWOlEmoKTwm59JRE/v3Kk+joH+Fjd784p/PCX9vXyaX/+Vder+/k3z9+Ep86szoClYZGRkoi91x7Kh9eXcp3nnyHL/xqC10DI0e93++2NHLzL1/n5IU5fPmDSyNQqYSSwlPC4uSFudx3/WnjAfqjl2bcOKR3aJTvPFnLR3/4AiNjfu6//j2sXx19OwklJ3r43sfWcOv5S3l8WzMf/N6z/Gl7C/5pFtOPjPn52V/3cMuvt1BTmcvG698TsbOlJHTi9hpGbrwq5bp161i0aBEtLS08+eSTRxx321UpD3fVVVeRlJTEpk2b2L59OzC+TrO2pZdRn5+OsrVcUVNB99436WreR9fgCB39I/j8lsIsL1+5+QayUpPmfVVKgKysLC677DJg7lelnCqUV6X873vuZUfz+IXnUhITKMxMYcHCai5t/leGRv18u3c9I2N+crzJLC3OwGOM/vYO+9v7Nb8G4GN8DAjub2+qYK/wOds1jPR1J2GVnpLImoocDvQOs61zkC//ZiurElso9/SS6DEUZKRQmJlCbkZazOwglJmaxKqybDr6R2jtHaahc4CX2htZ6ekDIC0tgUUFGeR4Y+P1xqu4bXlK5I35/Oxq7cdiMRgW5nlJS4792eUDPUN0DoxS8buP4vEYUm88smUnh4qGq2eq5SkRk5jgYVmJs6dXOqEoK3X8Uhoa14wpejdFIqXkBKcrkBBSeIpEyof+xekKJIS0VElEJAgKTxGRICg8RUSCoDFPEXGd5XnLnS7hqBSeIuI6Xz3tq06XcFTqtouIBEHhKSISBIWniEgQFJ4iIkFQeIqIBEHhKSISBIWniEgQFJ4iIkFQeIqIBEHhKSISBIWniEgQFJ4iIkFQeIqIBEHhKSISBIWniEgQFJ4iIkFQeIqIBEHhKSISBGOtdbqGY2aMaQXqnK4DKADanC4ihGLt9YBeU7Rwy2uqtNYWTncgJsLTLYwxm621NU7XESqx9npArylaRMNrUrddRCQICk8RkSAoPEPrbqcLCLFYez2g1xQtXP+aNOYpIhIEtTxFRIKg8BQRCYLCM0jGmJuNMZuNMcPGmHvncPsvGmNajDE9xpifGmNSIlDmvBhj8owxDxtj+o0xdcaYT8xy2w3GmFFjTN+Un0WRrHcmc30dZty3jTHtEz/fNsaYSNc7F/N4Ta59X6aaz+fHrZ8dhWfwmoA7gZ8e7YbGmA8CXwPOBSqBRcA3w1pdcP4TGAGKgauAHxpjVs5y+19bazOm/OyOSJVHN9fXcRNwKbAaOBFYD3wmUkXO03zeG7e+L1PN6fPj5s+OwjNI1tqHrLWPAO1zuPm1wD3W2u3W2k7gDuC6cNY3X8aYdOBy4P9aa/ustc8Dvweucbay+Znn67gW+DdrbYO1thH4N1z2vkDsvDdTzePz49rPjsIzMlYCW6f8/61AsTEm36F6prMUGLPW7pjyu62M1z6T9caYDmPMdmPM58Jb3pzN53VM977M9nqdMt/3xo3vS7Bc+9lReEZGBtA95f8H/jvTgVpmkgH0HPa7bmau8UFgBVAI3Ah8wxhzZfjKm7P5vI7p3pcMF457zuc1ufV9CZZrPzsKz2kYY542xtgZfp4P4iH7gKwp/z/w373HXu3czOE1HV5joM5pa7TWvmWtbbLW+qy1LwDfBz4a3lcxJ/N5HdO9L33WfYuf5/yaXPy+BMvxz85MFJ7TsNaeY601M/ycFcRDbmd8UiJgNbDfWjuX8dKQmMNr2gEkGmOOO6zO7XN9CsANLbb5vI7p3pe5vt5IOpb3xi3vS7Ac/+zMROEZJGNMojEmFUgAEowxqcaYxBluvhG43hhzvDEmB7gNuDdCpc6JtbYfeAj4J2NMujHmTOAS4L7pbm+MucQYkzux3Oc04PPA7yJX8fTm+To2Al8yxpQZY0qBW3HZ+wLze01ufV8ON4/Pj3s/O9Za/QTxA2xg/Ft96s+GiWMLGe9uLJxy+y8B+xkfu/oZkOL0a5jmNeUBjwD9wD7gE1OOvZfxLm3g//+S8ZnSPqAW+LzT9R/tdUzzGgzwHaBj4uc7TJyy7Lafebwm174vh72eaT8/0fTZ0bntIiJBULddRCQICk8RkSAoPEVEgqDwFBEJgsJTRCQICk8RkSAoPEVEgqDwlJhmjCk1xjRM8/sGY8z7ZriPKzffFXdReEqsuxB4EsAYk2WM+Z4xpgUoA/5gjPmjMebkwI3dvPmuuIvCU6KaMWatMWarMWadMabVGPOqMeb4KTe5EHh8Ypu5x4DTgHOARuBK4CXgeWPMCRO3d+3mu+IuOj1TopYxxgu8y3jL8BRgkPFQ/ARwMpAINAPVjF9m43mgwlrbMNGVv9pa+7Qx5s/AbmvtTcaYrcC3rLW/nniOAqAVKLAu2MlH3EMtT4lmpzG+ucePgXSgC/gP4ATGN5hYB2y11vYy3gVvs9YeMf4JvP7/27tflgiiKAzjz0UQFpMYVDAYTX4B/QaC38BiEjRZxGY2itlkEMS6Vcsu2AWLQWE3WzXIMdwbVtw1DIp3lueX5g8zzJSXOXPgXGC1bFc7fFd1MTzVZovAMEbKp4h4A16BJUrJXk69AAsppeUx91kHnst2tcN3VRfDU202BFZGl81IKXWAeWDA1/DsAz3gamQp3rmU0jGwCZyVY9UO31VdDE+12T15vuU+uXyfIf//7AGz5LmPjwDl63SLXKL3yd32a3LzaCMiHso96x2+q6rYMFKrlbXLL8gNog9yU2gX2AbWIuJgwnUDYCcibsecOwSOgA5wA+xFxPvfvIHayvDUVEgpXQJPEXFS9rvAeUR0f7xQasiyXdPqDvj2VSn9lkkLlkmtFhGn//0Mmm6W7ZLUgGW7JDVgeEpSA4anJDVgeEpSA4anJDVgeEpSA59P2QnZM0x5BQAAAABJRU5ErkJggg==\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "$\\displaystyle V = -\\alpha E_J \\cos\\left(\\pi + \\frac{2\\pi\\Phi_B}{\\Phi_0} + \\frac{2\\pi\\phi_+}{\\Phi_0}\\right) - 2E_J\\cos\\left(\\frac{2\\pi\\phi_+}{2\\Phi_0}\\right)$\n",
        "\n",
        "$\\Phi_B = BS$\n",
        "\n",
        "$\\mu = \\displaystyle\\left.\\frac{\\partial V}{\\partial B}\\right|_{B=0} = \\left.\\frac{\\partial V}{\\partial \\Phi_B}\\right|_{\\Phi_B=0}\\times S = \\alpha E_J \\sin\\left(\\pi + \\frac{2\\pi\\phi_+}{\\Phi_0}\\right)\\frac{2\\pi S}{\\Phi_0}$\n"
      ],
      "metadata": {
        "id": "N2UlIzZSq4wF"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "mu = Alpha * EJ * np.sin(np.pi+min_phi*2*np.pi) * 2 * np.pi * S / Phi0 \n",
        "print('Mangnetic moment : {}'.format(mu))"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "RlBR_ZNEqomB",
        "outputId": "78a2ae44-2969-42aa-bbc7-52c65ff2fec8"
      },
      "execution_count": 4,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Mangnetic moment : -4.416933826672361e-12\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "The order of magnetic moment for atoms are given by\n",
        "\n",
        "$\\displaystyle m_B = \\frac{e\\hbar}{2m_e} \\sim 10^{-23}$\n",
        "\n",
        "The magnetic moment of qubit is $10^{11}$ times larger than atoms. This comes from the large current $I$ on the circuit (compared to the spin of atoms.)"
      ],
      "metadata": {
        "id": "3Dunb4JktFTu"
      }
    }
  ]
}