Evolution_in_Closed_and_Open_Systems.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "colab_type": "text",
        "id": "view-in-github"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/MonitSharma/Learn-Quantum-Machine-Learning/blob/main/Evolution_in_Closed_and_Open_Systems.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4Bmw8IqFa1CV"
      },
      "source": [
        "Classical probability distributions can be written as a stochastic vector, which can be transformed to another stochastic vector by applying a stochastic matrix. In other words, the evolution of stochastic vectors can be described by a stochastic matrix.\n",
        "\n",
        "Quantum states also evolve and their evolution is described by unitary matrices. This leads to some interesting properties in quantum computing. Unitary evolution is true for a closed system, that is, a quantum system perfectly isolated from the environment. This is not the case in the quantum computers we have today: these are open quantum systems that evolve differently due to to uncontrolled interactions with the environment. In this notebook, we take a glimpse at both types of evolution.\n",
        "\n",
        "\n",
        "# Unitary evolution\n",
        "\n",
        "A unitary matrix has the property that its conjugate transpose is its inverse. Formally, it means that a matrix $U$ is unitary if $UU^\\dagger=U^\\dagger U=\\mathbb{1}$, where $^\\dagger$ stands for conjugate transpose, and $\\mathbb{1}$ is the identity matrix. A quantum computer is a machine that implements unitary operations.\n",
        "\n",
        "As an example, we have seen the NOT operation before, which is performed by the X gate in a quantum computer. While the generic discussion on gates will only occur in a subsequent notebook, we can study the properties of the X gate. Its matrix representation is $X = \\begin{bmatrix} 0 & 1\\\\ 1 & 0\\end{bmatrix}$. Let's check if it is indeed unitary:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "ExecuteTime": {
          "end_time": "2018-11-19T19:49:07.614968Z",
          "start_time": "2018-11-19T19:49:07.530927Z"
        },
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "8hbMef9ra1Ck",
        "outputId": "ec9b6328-e1c1-4652-ae0d-36ffa52ac5c9"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "XX^dagger\n",
            "[[1 0]\n",
            " [0 1]]\n",
            "X^daggerX\n",
            "[[1 0]\n",
            " [0 1]]\n"
          ]
        }
      ],
      "source": [
        "import numpy as np\n",
        "X = np.array([[0, 1], [1, 0]])\n",
        "print(\"XX^dagger\")\n",
        "print(X @ X.T.conj())\n",
        "print(\"X^daggerX\")\n",
        "print(X.T.conj() @ X)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "z0s5kKpca1Cp"
      },
      "source": [
        "It looks like a legitimate unitary operation. The unitary nature ensures that the $l_2$ norm is preserved, that is, quantum states are mapped to quantum states."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "ExecuteTime": {
          "end_time": "2018-11-19T19:49:07.626531Z",
          "start_time": "2018-11-19T19:49:07.618259Z"
        },
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "aaViof-Qa1Cq",
        "outputId": "3039c06f-1117-4ef6-aa01-472dd84ff434"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "The norm of the state |0> before applying X\n",
            "1.0\n",
            "The norm of the state after applying X\n",
            "1.0\n"
          ]
        }
      ],
      "source": [
        "print(\"The norm of the state |0> before applying X\")\n",
        "zero_ket = np.array([[1], [0]])\n",
        "print(np.linalg.norm(zero_ket))\n",
        "print(\"The norm of the state after applying X\")\n",
        "print(np.linalg.norm(X @ zero_ket))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "8GgFIPwDa1Cr"
      },
      "source": [
        "Furthermore, since the unitary operation is a matrix, it is linear. Measurements are also represented by matrices. These two observations imply that everything a quantum computer implements is actually linear. If we want to see some form of nonlinearity, that must involve some classical intervention.\n",
        "\n",
        "Another consequence of the unitary operations is reversibility. Any unitary operation can be reversed. Quantum computing libraries often provide a function to reverse entire circuits. Reversing the X gate is simple: we just apply it again (its conjugate transpose is itself, therefore $X^2=\\mathbb{1}$)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "ExecuteTime": {
          "end_time": "2018-11-19T19:49:08.710098Z",
          "start_time": "2018-11-19T19:49:07.629733Z"
        },
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "drJJ58sta1Cs",
        "outputId": "96192518-8d9a-4120-ab0d-4c6b4dd23b6d"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "[1.+0.j 0.+0.j]\n"
          ]
        }
      ],
      "source": [
        "import numpy as np\n",
        "from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister\n",
        "from qiskit_aer import AerSimulator\n",
        "from qiskit.quantum_info import Statevector\n",
        "np.set_printoptions(precision=3, suppress=True)\n",
        "\n",
        "backend_statevector = AerSimulator(method='statevector')\n",
        "q = QuantumRegister(1)\n",
        "c = ClassicalRegister(1)\n",
        "circuit = QuantumCircuit(q, c)\n",
        "circuit.x(q[0])\n",
        "circuit.x(q[0])\n",
        "\n",
        "# get the statevector\n",
        "state = Statevector.from_instruction(circuit)\n",
        "print(state.data)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Xed2txNEa1Cu"
      },
      "source": [
        "which is exactly $|0\\rangle$ as we would expect.\n",
        "\n",
        "In the next notebook, you will learn about classical and quantum many-body systems and the Hamiltonian. In the notebook on adiabatic quantum computing, you will learn that a unitary operation is in fact the Schrödinger equation solved for a Hamiltonian for some duration of time. This connects the computer science way of thinking about gates and unitary operations to actual physics, but there is some learning to be done before we can make that connection. Before that, let us take another look at the interaction with the environment."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "D0J9W8Xxa1C-"
      },
      "source": [
        "# Interaction with the environment: open systems\n",
        "\n",
        "Actual quantum systems are seldom closed: they constantly interact with their environment in a largely uncontrolled fashion, which causes them to lose coherence. This is true for current and near-term quantum computers too.\n",
        "\n",
        "![](https://raw.githubusercontent.com/qosf/qml-mooc/9143b9aba0e58a049e2e4b98efd4b7caf471bd6e/figures/open_system.svg)\n",
        "\n",
        "This also means that their actual time evolution is not described by a unitary matrix as we would want it, but some other operator (the technical name for it is a completely positive trace-preserving map).\n",
        "\n",
        "Quantum computing libraries often offer a variety of noise models that mimic different types of interaction, and increasing the strength of the interaction with the environment leads to faster decoherence. The timescale for decoherence is often called $T_2$ time. Among a couple of other parameters, $T_2$ time is critically important for the number of gates or the duration of the quantum computation we can perform.\n",
        "\n",
        "A very cheap way of studying the effects of decoherence is mixing a pure state with the maximally mixed state $\\mathbb{1}/2^d$, where $d$ is the number of qubits, with some visibility parameter in $[0,1]$. This way we do not have to specify noise models or any other map modelling decoherence. For instance, we can mix the $|\\phi^+\\rangle$ state with the maximally mixed state:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "ExecuteTime": {
          "end_time": "2018-11-19T19:49:08.730307Z",
          "start_time": "2018-11-19T19:49:08.714253Z"
        },
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "iK1FfvoFa1C_",
        "outputId": "145cf394-3dea-47b6-bd83-3fbf98f7940c"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Maximum visibility is a pure state:\n",
            "[[0.5 0.  0.  0.5]\n",
            " [0.  0.  0.  0. ]\n",
            " [0.  0.  0.  0. ]\n",
            " [0.5 0.  0.  0.5]]\n",
            "The state is still entangled with visibility 0.8:\n",
            "[[0.45 0.   0.   0.4 ]\n",
            " [0.   0.05 0.   0.  ]\n",
            " [0.   0.   0.05 0.  ]\n",
            " [0.4  0.   0.   0.45]]\n",
            "Entanglement is lost by 0.6:\n",
            "[[0.4 0.  0.  0.3]\n",
            " [0.  0.1 0.  0. ]\n",
            " [0.  0.  0.1 0. ]\n",
            " [0.3 0.  0.  0.4]]\n",
            "Barely any coherence remains by 0.2:\n",
            "[[0.3 0.  0.  0.1]\n",
            " [0.  0.2 0.  0. ]\n",
            " [0.  0.  0.2 0. ]\n",
            " [0.1 0.  0.  0.3]]\n"
          ]
        }
      ],
      "source": [
        "def mixed_state(pure_state, visibility):\n",
        "    density_matrix = pure_state @ pure_state.T.conj()\n",
        "    maximally_mixed_state = np.eye(4)/2**2\n",
        "    return visibility*density_matrix + (1-visibility)*maximally_mixed_state\n",
        "\n",
        "ϕ = np.array([[1],[0],[0],[1]])/np.sqrt(2)\n",
        "print(\"Maximum visibility is a pure state:\")\n",
        "print(mixed_state(ϕ, 1.0))\n",
        "print(\"The state is still entangled with visibility 0.8:\")\n",
        "print(mixed_state(ϕ, 0.8))\n",
        "print(\"Entanglement is lost by 0.6:\")\n",
        "print(mixed_state(ϕ, 0.6))\n",
        "print(\"Barely any coherence remains by 0.2:\")\n",
        "print(mixed_state(ϕ, 0.2))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ifZmi3oaa1DA"
      },
      "source": [
        "Another way to look at what happens to a quantum state in an open system is through equilibrium processes. Think of a cup of coffee: left alone, it will equilibrate with the environment, eventually reaching the temperature of the environment. This includes energy exchange. A quantum state does the same thing and the environment has a defined temperature, just like the environment of a cup of coffee.\n",
        "\n",
        "The equilibrium state is called the thermal state. It has a very specific structure and we will revisit it, but for now, suffice to say that the energy of the samples pulled out of a thermal state follows a Boltzmann distribution. The Boltzmann -- also called Gibbs -- distribution is described as $P(E_i) = \\frac {e^{-E_{i}/T}}{\\sum _{j=1}^{M}{e^{-E_{j}/T}}}$, where $E_i$ is an energy, and $M$ is the total number of possible energy levels. Temperature enters the definition: the higher the temperature, the closer we are to the uniform distribution. In the infinite temperature limit, it recovers the uniform distribution. At high temperatures, all energy levels have an equal probability. In contrast, at zero temperature, the entire probability mass is concentrated on the lowest energy level, the ground state energy. To get a sense of this, let's plot the Boltzmann distribution with vastly different temperatures:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "ExecuteTime": {
          "end_time": "2018-11-19T19:49:09.226294Z",
          "start_time": "2018-11-19T19:49:08.733112Z"
        },
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 283
        },
        "id": "RIjKKK3la1DB",
        "outputId": "2284da30-c234-4606-c1c2-0f2473fc08f0"
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<matplotlib.legend.Legend at 0x172c55d90>"
            ]
          },
          "execution_count": 5,
          "metadata": {},
          "output_type": "execute_result"
        },
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhcAAAGZCAYAAAA6ixN9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAABGQElEQVR4nO3dd3hUZd7G8e/MpHdIaJHQEUU6iooUQRaRXUURsCu2VQGxL7rq6lpWF913XSy4qy5rb6hYVkVFFAUUVqo0BVEinVDSIGVm3j8Gkpwzk2SSzMyZcn+uay5ynnNm5pcAyZ3nPMXmdrvdiIiIiASI3eoCREREJLooXIiIiEhAKVyIiIhIQClciIiISEApXIiIiEhAKVyIiIhIQClciIiISEApXIiIiEhAxVnxpi6Xi23btpGeno7NZrOiBBEREWkgt9tNUVERubm52O21909YEi62bdtGXl6eFW8tIiIiTZSfn0/btm1rPW9JuEhPTwc8xWVkZFhRQq1+3FnEOU8tMrStvnekeljC2eKn4Iu/VB8nNYOpy8ARb11NIiJRqLCwkLy8vKqf47WxJFwc+UGdkZERduEip8KBPTHF0JaYkkZSvMOiiqReA86HxQ9VH7v3Q8Fy6DLCspJERKJZfb9wa0CnSbKPEHGowmlBJeK3Zh3gqP7Gtu/ftqQUERFRuPCS6CNcHFS4CH89zjUer3sfKg5ZU4uISIxTuDDx3XPhsqASaZDjzgFqdNOVFcLGzywrR0QklilcmMQ7bDjsxntJB8vVcxH2MnKh/UBj2/dvWVOLiEiMU7gwsdlsJMUZvyy6LRIheow1Hv/wMZSXWFOLiEgMU7jwITnBeGukTOEiMnQ/G2w1/u4qSmHDR5aVIyISqxQufEiMM4YL9VxEiNQc6HSqsU2zRkREQk7hwgdzz4UGdEYQ86yRjZ/Cwf2WlCIiEqsULnwwzxhRz0UEOea34EioPnaWw/oPrKtHRCQGKVz4kBSvAZ0RKzkLuvzG2KZZIyIiIaVw4YN5qW8N6Iww5lkjP30JxbutqUVEJAYpXPhgDhda5yLCdDsD4mvsD+N2wrp3ratHRCTGKFz4YB5zcahS4SKiJKR6AkZNq3VrREQkVBQufPAa0Fmu2SIRxzxrZMsi2J9vTS0iIjFG4cIHDeiMAl1GQFKmsW31m9bUIiISYxQufEjSCp2RLy7Rs2JnTaveALfbknJERGKJwoUPSVqhMzr0Os94vHsd7PzemlpERGKIwoUP3it0KlxEpHYnQ0ZbY9uq162pRUQkhihc+KAVOqOE3Q69xhvbVs8Gl/4+RUSCSeHCB+8BnZotErHMt0aKtsPPX1lTi4hIjFC48EErdEaRlsdC657GtlWaNSIiEkwKFz54rdCpcBHZek4wHq99FyoOWlOLiEgMULjwwWuFToWLyNZzHGCrPi4vgg0fWVaOiEi0U7jwwTxbRHuLRLiMXOg4xNi26g1rahERiQEKFz6Y17k4pAGdkc88sHPjp1BSYE0tIiJRTuHCh+QE45el3OnC6dLKjhHt2DMhLqn62FUJa9+xrh4RkSimcOFDoqnnAjTuIuIlZXjvlLpSC2qJiASDwoUP5jEXoHARFcy3Rn5dAgWbrKlFRCSKKVz4YJ4tApqOGhU6nwYp2ca2la9aU4uISBRTuPDBvM4FqOciKsQlQE/TcuArXwOXBuyKiASSwoUPDruNBIfxS6MZI1Giz4XG4wP5Wg5cRCTAFC5qkei1v4h6LqJC617Q8jhj24pXrKlFRCRKKVzUQqt0Rimbzbv3Yt17UFZkTT0iIlFI4aIWWqUzivWaALYaf78VpZ79RkREJCAULmphXqVTt0WiSFpL6PobY9sKzRoREQkUhYtaJCWYt13XgM6oYr418svXsHezNbWIiEQZhYtaJMVpQGdUO3oUJGUZ21a+ZkkpIiLRRuGiFuYxFxrQGWXiEn2sefGq1rwQEQkAhYtamGeLqOciCvW5wHi8/xfYssiaWkREoojCRS3Mq3QqXESh3H7Q4hhjmwZ2iog0mcJFLczhQgM6o5DNBr1NvRdr3tGaFyIiTaRwUYsk8wqdWuciOvU6D2w1/q4rSjwBQ0REGk3hohZeK3RWKlxEpYw20HWksW3Zi9bUIiISJRQuauE1oFM9F9Gr36XG41+XwK711tQiIhIFFC5qoQGdMaTrSEhtaWxbrt4LEZHGUriohVbojCGOeO9pqStfhcpya+oREYlwChe10AqdMaav6dZIaQFs+NCaWkREIpzCRS20QmeMyekC7QYa23RrRESkURQuaqEVOmOQeWDnxnmwP9+aWkREIpjCRS3MAzrVcxEDuo+BxIwaDW5Y8Ypl5YiIRCqFi1p4hwsN6Ix6CSnQc5yxbflL2sxMRKSBFC5q4bVCp3ouYkPfS4zHB7bA5i8sKUVEJFIpXNTCPObC6XJT4dRvsFEvty+06mFsW/aCNbWIiEQohYtamGeLgHovYoLN5j2wc90HULzbmnpERCKQwkUtkuK8w8UhLQEeG3qOB0di9bGrAla8bF09IiIRRuGiFr56LjSoM0akNIfjzjG2fTdLAztFRPwUZ+Wb37f4PhJTE7Hb7NiwYT+89bXdZjd8fOSczWar/hgbNpvx47quO3Jsft2ar2E+F5e5GrB5Hm4bc385SO6B5Hrf39/38Ld2O4fba74HdrCZvj413qvRX7vDxzHv+Ctg1WvVx/t+9gzs7DzcqopERCKGpeHiw80f4kj27iEIF8m5xuMZq62pwwrm4FFncKsZgPwIMAEPX3WEpIa+huH923aCg/uwuz0R0774fmx7/9fk96jtWnOQ9Pqamq6z2w4HzMPXVQXOw+cC9bo1P4+a5+r9N1DHvxXDvwFT7TXrEZHIZGm4kPDlch++BeC2tg5LxQPx6dXHzj2wVkuCh5K/4TRowdVHr6GvUFln7+KR8Ge6Hhve4auBn4ev2hv8Ofrx/g16j3rCM4d7XM1Bsravp9d1Ps4F5HV9/J0p5DaewoWIhC2FXAkH/gSyUAY3f96fAPZq1ny/8lL/dou2NFxcfMzFJKQm4MZd9U3E5XZVHbvdbq9zR87XPGf4+PA5X9e5cIEb43W4qq6veo/D163dvp+yShfgBpuLts2SSUtyeF6zxnN9vWfVOV/X4fk8cGM4V/NzMX9OIiJiDTdunG7nkYOY5jzo36xJS8PFlH5TyMjIqP9Ci5zxj69Yt72w6vi6Ab05p29bS2rxN0wFLIgdfp2az6u6rsbzaoafqmN31VV1P/dwiKv5njXrqHmd2+2u9/PyFfbqfY9avgZu9+FvJiW7ca99F0/EtOGygbvTMNwZuYbr3J43MYTVmmESN1XXVb2P6etS29eqqpaafx+4cblM1+HG6XIa/j7M1x35e3W6nb6DbB3/PkRE/KXbInVINi0BbuVUVHO3mYTQxhWwfUX1cfMSGP2AVdVYxtz7VleYamzvXs3eRV89fr5e19zjV1/PYM3aq4JTbe9Zsyf1cBsYQ67h2N9g3ZTw7mew9uc6N4eD6+EAXtvX0xCOzQH88L8NX/8W/ArgR0Jvjfc48vn6vE4igsJFHcxrXRzUIlqx6fjL4f0bqo9//MSzFXtWnnU1WeBIwMUGDsJ3lpdEN5+9e7WEoiNh1BwQ/emp89WTG9TeYtNte18B1OvzqBmmwed15nO1BXNfr+vruoPFB1nHunr/nhQu6mBepVPLf8eoHuNg7l1QXuQ5drs8+40Mv9PaukRikM1mw2FTuLVKYWEhf+Ev9V6nPvY6JJl6LsoULmJTYhr0mmBsW/YCOCusqUdEJMwpXNRBPRdS5fgrjMfFO2Dd+9bUIiIS5hQu6pCcED4DOsVirXtA3knGtiXPWFOLiEiYU7ioQ3K8ei6khgFXG4+3LIId31tTi4hIGFO4qEOSwoXUdOxZkNbK2LZUvRciImYKF3UwhwsN6IxxcQnQ/3Jj26o34OA+a+oREQlTChd1UM+FeOk/Eew1ZnBXlMKKVywrR0QkHClc1ME85kIDOoWMNnDsmca2Jc+AS/82RESOULiog3m2iFboFAAG/N54vG8zbJpnTS0iImFI4aIO5nUuDum2iAC0Oxla9TC2LfmXNbWIiIQhhYs6mFfoVLgQAGw272mpP34Ke3+yph4RkTCjcFEHrdApteo5HhIzazS4YelzlpUjIhJOFC7qYN4VVQM6pUpCKvS92Ni2/EUoK7amHhGRMKJwUQdfK3S63W6LqpGwc8KVxuNDB2Dlq9bUIiISRhQu6pAU7/3lKatU74Uclt0Zjh5lbPv2aU1LFZGYp3BRB3PPBWhQp5icdJ3xuGAjbPzUmlpERMKEwkUdEn2ECw3qFIOOQ6Hlcca2b56yphYRkTChcFEH3z0X6vKWGmw2796Ln76AnWstKUdEJBwoXNQh3mHDYbcZ2rRKp3jpOR5Sso1t6r0QkRimcFEHm81GUpxpCXDdFhGz+CQ43jRzZNUbULLHmnpERCymcFEP81oX2nZdfDrhKrDHVx87y+B/s6yrR0TEQgoX9UjUKp3ij/RW0HOcsW3pM1BZZk09IiIWUrioh1bpFL+ZB3YW74Q171hTi4iIhRQu6uFrlU4Rn9r0hvaDjG2LnwSt6ioiMUbhoh7mVToVLqRO5t6LHatg8wJrahERsYjCRT2S4jWgUxqg2xnQrKOxbdEMa2oREbGIwkU9zOFC61xInewOOHmysW3jZ7BzjTX1iIhYQOGiHuYxF4cqFS6kHn0u8l5Ua9Hj1tQiImIBhYt6eA3oLNdsEalHQgqccLWxbfWbcGCrNfWIiISYwkU9NKBTGmXA1RCXVH3sqoRvZ1pXj4hICClc1CNJK3RKY6TmeG6P1PS//8ChA5aUIyISSgoX9UjSCp3SWCdPBmpsfFdeBN89b1k5IiKhonBRD/MKnQoX4rfsznDs74xt38yEynJr6hERCRGFi3qkmMJFqaaiSkMMvMF4XLQNvn/LmlpEREJE4aIe6UlxhuOiQ5UWVSIRKe8EaHeysW3RDC0JLiJRTeGiHmmJ8Ybj4rIKiyqRiDVwqvF411r4Ya41tYiIhIDCRT3SEtVzIU109CjI6WZs++pR9V6ISNRSuKiH+bZI8aFK3PqhIA1ht8Ogm4xtvy6Fn7+2ph4RkSBTuKiHOVxUutyUVWqVTmmgnuMgs52x7au/WVOLiEiQKVzUIz0p3qut8JDGXUgDOeLhFNPYi5/mw9Zl1tQjIhJEChf1SE10eLUVa9yFNEbfiyG1hbHt6/+zphYRkSBSuKhHYpyDhDjjl0mDOqVR4pO9t2Nf9z7s3mBNPSIiQaJw4YcM86DOMoULaaTjr4TETGPb149ZUoqISLAoXPhB01ElYJIy4MTfG9tWvQ77frGmHhGRIFC48EOa1yqdGtApTXDitRCXXH3sdsKix62rR0QkwBQu/JDutUqnei6kCVJzoP9EY9uyF6BohyXliIgEmsKFH8w9F5otIk02cArYa4RWZxks/Id19YiIBJDChR+8Ni9Tz4U0VWZb6HuRse1//4aindbUIyISQAoXfkjXgE4JhkE3g73Gv63KQ54dU0VEIpzChR80oFOColl76H2BsW3pc1C8y5p6REQCROHCD+YlwDWgUwJm8C1gq7EKbOVBzRwRkYincOEH8zoXGtApAdO8o4/ei2ehZI819YiIBIDChR+8BnQqXEggDTH1XlSUqvdCRCKawoUfzOFCt0UkoJp3gl7nGduWPAMlBdbUIyLSRAoXfkgzLaKlAZ0ScENuBVuN/44VJbD4CevqERFpAoULP/jquXC73RZVI1EpuzP0HG9sW/Iv9V6ISERSuPCDeUCnyw2l5U6LqpGoNeQ2Y+9FeTEsfMyyckREGkvhwg/mngvQoE4Jgpyu0HOCsW3JM9pzREQijsKFH8w9FwDFZRp3IUFw6jTvdS+++j/r6hERaQSFCz/EOewkxzsMbeq5kKBo3gn6Xmxs+24W7M+3ph4RkUZQuPCT9xLgChcSJENuA0dC9bGzHBY8Yl09IiINpHDhJ611ISGTlQf9Lze2LX8JCjZZU4+ISAMpXPjJvDOqlgCXoBp8M8QlVR+7nfDldOvqERFpAIULP5k3LyvUQloSTOmtYcDVxrbVb8DuDdbUIyLSAAoXfvLavEy3RSTYTrkREtKqj90umP8Xy8oREfGXwoWfNKBTQi41B066zti2dg5sW25JOSIi/lK48JPXgE6FCwmFk6dAUqax7bN7LSlFRMRfChd+8hrQqdsiEgrJWTDoJmPbT1/Aps+tqEZExC8KF34y3xbRgE4JmQHXQHobY9tn94LLZUk5IiL1Ubjwk3m2iHouJGQSUuDUO4xt21fCmretqUdEpB4KF34yzxbRgE4JqT4XQc7RxrbPH4DKcmvqERGpg8KFnzSgUyzliIPT/mRs27cZlj1vTT0iInVQuPCTlv8Wyx3zO2h7grHty79CWbE19YiI1ELhwk9pid5jLpwut0XVSEyy2WDEvca2kt2w+ElLyhERqY3ChZ/MPRcAJeXqvZAQ6zAIuo40ti2aAUU7ralHRMQHhQs/maeigsZdiEVOuwewVR+XF8MXWhZcRMKHwoWfUhO8w4VmjIglWveA3hcY25a9ADvXWlOPiIiJwoWfHHabj83LtJCWWOS0uyEuufrY7YJP77auHhGRGhQuGsAcLgrVcyFWyciFgdcb2zZ+BhvnWVOPiEgN3n39Uqv0pDh2FFYfa8yFWOqUGzzrXBTXGMz5yd3Q6VSwOywrSyTSud1uKisrcTqdVpcScg6Hg7i4OGw2W/0X10HhogHMgzq11oVYKjENht8F79Xowdi1Bpa/BP0vs64ukQhWXl7O9u3bKS0ttboUy6SkpNCmTRsSEhIa/RoKFw3gvQS4xlyIxfpcBN/+E3Z+X902/0Hoca4nfIiI31wuF5s3b8bhcJCbm0tCQkKTf4OPJG63m/Lycnbv3s3mzZvp2rUrdnvjRk8oXDRAhnnzMt0WEavZHTDyfnjxnOq24p2w8B8w/E7r6hKJQOXl5bhcLvLy8khJSbG6HEskJycTHx/PL7/8Qnl5OUlJSY16nUZFkvnz5zfqzSKdV8+FbotIOOg8HLr8xti2aAbsz7emHpEI19jf1qNFID7/Rr3CqFGj6Ny5Mw888AD5+bHzDcw85kLrXEjYGHk/2Gr8d648pKmpImKZRoWLrVu3MmXKFGbPnk2nTp04/fTTeeONNygvj+7tn7UzqoStlsfC8VcY29a8Az9/bU09IhLTGhUucnJyuOmmm1ixYgXffvstRx99NJMmTSI3N5epU6eycuXKQNcZFrxvi2hAp4SRYXdCUpax7aPbwRV70+lExFpNvrHSr18/7rjjDqZMmUJxcTH//ve/6d+/P4MHD2bNmjWBqDFsaECnhLWU5p6pqTXtXA3f/ceSckQkNObOnYvNZqvz8cknn4S0pkaHi4qKCmbPns3o0aNp3749c+fO5YknnmDnzp1s3LiR9u3bM378+EDWajmvMRca0Cnhpv/l0LK7se3zB+DgPmvqEZGgGzJkCNu3b696ZGdnc/fddxvaTjvttJDW1KipqNdffz2vvvoqbrebSy65hOnTp9OjR4+q86mpqTz66KPk5uYGrNBw4L3OhcKFhBlHHIx6GF44q7rt4F6Y/xCMnm5dXSIRxuVys6/U2nGEzVISsNvrX2cjOTmZ5GTPXkNbt26loKCAwYMH07p162CXWKtGhYu1a9fy+OOPM3bsWBITE31ek5OTE3VTVjWgUyJCp6Fw7Jmw7v3qtqXPQv+J0Kp7rU8TkWr7Ssvp/8Bnltbw3V0jyE7z/TO2NsuXLwc8Qxas1KjbIvfccw/jx4/3ChaVlZUsWLAAgLi4OIYOHdr0CsOIOVwcrHBS6XRZVI1IHUY+AI4a/z/dTvj4dnC7ratJRIJu2bJl5OXlkZ2dbWg/55xzaNasGePGjQtJHY0KF8OGDWPv3r1e7QcOHGDYsGFNLipcpSXGe7VpfxEJS806wClTjW2bv/RMTxWRqLVs2TKfvRY33HADL7zwQsjqaFS4cLvdPtdbLygoIDU1tclFhStzzwVo3IWEsUE3QcZRxra5f4SyImvqEZGgqy1cnHrqqaSnp4esjgaNuRg7diwANpuNiRMnGm6LOJ1OVq1axcCBAwNbYRhJSXBgt4GrRs+ywoWErYRUOP0v8GaNHVKLtsMXD8PpD1pXl0gEaJaSwHd3jbC8hobYs2cP+fn5lo+3gAaGi8zMTMDTc5Genl41OhUgISGBk046iauvvjqwFYYRm81GWmIchTUChW6LSFjrPsaz98imz6vbvpkJfS6EVsdZV5dImLPbbQ0eTGm1ZcuWAdYP5oQGhotZs2YB0KFDB2699daovgVSm/SkeFO40CqdEsZsNhj9KDx1EjgPT6tzO+G/t8LlH3rOi0hUWL58Oa1atQqLZSAaPVskFoMFaK0LiUDZneGUG41tWxbBytcsKUdEgmPatGns2LHD6jKABvRc9OvXj3nz5tGsWTP69u3rc0DnEUe6ZqKReVCnwoVEhME3w6rXYf8v1W2f3g3dzoDkLMvKEpHQGDFiBCtXrqSkpIS2bdvy5ptvcvLJJwft/fwOF2PGjKkawHn22WcHq56wZ14CXGMuJCLEJ8MZ0+HV86rbSnZ7lgb/7aPW1SUiIfHZZ6FdEMzvcHHPPff4/DjWeN8W0ZgLiRDdRkG30bDhw+q2pc9C7wugbX/r6hKRqNPkXVFjTbp2RpVINuphiEuu0eCG96eCUyFZRALH756LZs2a1TnOoiZfq3dGC425kIjWrD2cOg0+u7e6bef3sPhJGHSjVVWJSJTxO1w89thjQSwjcqSbb4tozIVEmpOnwOrZnlBxxBcPQ/ezoHkn6+oSkajhd7i47LLL6r8oBngN6FTPhUQaRzyc+Q94dgRweLnZyoPwwc1wyTta+0JEmszvMReFhYWGj+t6RDOvAZ1aREsiUdvjYYBpNd2f5sPqN62pR0SiSoPGXGzfvp2WLVuSlZXlc/zFkQ3NnE5nQIsMJxrQKVFj+N2w7gMo2lbd9vHt0GUEpDS3ri4RiXh+h4vPP/+c5s0933Dmz58ftILCnQZ0StRIyoDRj8DrF1W3lRbA3DvhnJnW1SUiEc/vcDF06FCfH8ca79siChcSwY79HRzzO1j/QXXbyleg5zjocpp1dYlIRGvQxmU17du3j+eee45169YB0L17dy6//PKq3o1oZe65KK90UVbpJDHOYVFFIk00+hH46UsoL6pue/8GmLQYEtOtq0tEIlajFtFasGABHTp0YMaMGezbt499+/YxY8YMOnbsyIIFCwJdY1gxzxYBjbuQCJeRCyPvM7YdyDeuhSEiYWvu3LnYbLY6H5988klIa2pUz8XkyZM577zzmDlzJg6H5zd2p9PJpEmTmDx5MqtXrw5okeEkwzSgEzz7i2SnJVpQjUiA9JsI378NP39V3bb0Weh+NnQcbFVVIuKHIUOGsH379qrjHj16MGnSJCZNmlTV1qJFi5DW1Kiei40bN3LLLbdUBQsAh8PBzTffzMaNGwNWXDhKjLMTZzfOlNGgTol4djucNQPiU4zt702B8hJrahIRvyQnJ9O6dWtat26N0+mkoKCAwYMHV7W1bt3a8PM6FBrVc9GvXz/WrVtHt27dDO3r1q2jd+/eASksXNlsNtKS4thfWr2+hcKFRIXmneC0P3mmox6x72f4/EEY9RfLyhKxhMsFBy3eyiK5uSf4N8Dy5csBz89pK/kdLlatWlX18dSpU7nhhhvYuHEjJ510EgDffPMNTz75JA8//HDgqwwz6aZwoW3XJWoM+D2seQfyv61u++YpOO5syBtgWVkiIXdwLzzS2doabtsEqTkNesqyZcvIy8sjOzu7qi0/P59LLrmEXbt2ERcXx91338348eMDXa2B3+GiT58+2Gw23G53Vdsf/vAHr+suvPBCzjvvvMBUF6bSEuOBg1XH2nZdoobdAWOehJmngLPscKMb5kyCa7+C+OQ6ny4i1lq2bJlXr0VcXByPPfYYffr0YceOHfTv35/Ro0eTmpoatDr8DhebN28OWhGRxrx5mXouJKrkdIVhf4TP7qluK/gR5t0Hox6yri4RqdeyZcu46qqrDG1t2rShTZs2ALRu3ZqcnBz27t0bHuGiffv2QSsi0miVTol6J0+Bte/CtmXVbd88Bd1Ga/aISJjas2cP+fn5dY63+O6773A6neTl5QW1lkYvogWwdu1atmzZQnl5uaH9rLPOalJR4c681oXChUQdRxyc8zQ8PbjG7RE8t0euW+hZOlwkmiU394x5sLqGBli2zPPLQG3hYu/evVx66aU888wzTS6tPo0KFz/99BPnnHMOq1evNozDOLKZWTRvXAbePRfF2hlVolGLbjDiHpj7x+q2A1s8x2OesK4ukVCw2xs8mNJqy5cvp1WrVuTm5nqdKysr4+yzz+b2229n4MCBQa+lUetc3HDDDXTs2JFdu3aRkpLCmjVrWLBgAccffzxffPFFgEsMP54BndW0QqdErROvg/aDjG3LX4Qf5lpTj4jUatq0aezYscOr3e12M3HiRIYPH84ll1wSkloaFS4WL17MfffdR05ODna7HbvdzqBBg3jooYeYOnVqoGsMOxpzITHDboezn4SENGP7e9dDqcVrAIiIXxYuXMjrr7/OnDlz6NOnD3369An6StqNui3idDpJT/dsaJSTk8O2bdvo1q0b7du3Z8OGDQEtMBxlmMLF/oO6LSJRrFkHOP0v8H6NXxyKd8IHN8H4/4DNVtszRSQMDBo0CJfLFdL3bFTPRY8ePVi5ciUAJ554ItOnT2fhwoXcd999dOrUKaAFhqPmqcZ9RPaWlNdypUiU6HcpdD3d2LZ2Dqx8zZJyRCS8NSpc3HXXXVUp6L777mPz5s0MHjyYDz/8kBkzZgS0wHDUPDXBcFxQXFbLlSJRwmbz7D2S3MzY/uGtsPcna2oSkbDVqNsip59e/RtMly5dWL9+PXv37qVZs2ZVM0aiWU6aMVwUHqqkvNJFQlyjsppIZEhvDWf+A964tLqtvBjeuhqu+Bgc3jsGi0hsavJPw/z8fPLz82nevHlMBAvw7rkA2FeqWyMSA7qPgb6m0eZb/wcLHrGmHhEJS40KF5WVldx9991kZmbSoUMHOnToQGZmJnfddRcVFdE/uDErJQHTruvs0a0RiRWjHobmpg2dFjwCvyy2ph4RCTuNChfXX389//rXv5g+fTrLly9n+fLlTJ8+neeeey4mpqI67DaapRh7LzSoU2JGYhqc+wzYa9xVdbvg7d/DoQPW1SUiYaNRYy5eeeUVXnvtNc4444yqtl69epGXl8cFF1zAzJkzA1ZguMpOS6CgRqAoKFa4kBhyVH/P5mbz7qtuO7AF/nsLjH1G01NFYlyjei4SExPp0KGDV3vHjh1JSPAejxCNsk3TUQvUcyGx5pQbvVfvXP0mrHjZknJEJHw0KlxMmTKF+++/n7Ky6nEGZWVlPPjgg0yZMiVgxYWz5mmajioxzu6Asf+EpExj+39vhV3rralJRMKC37dFxo4dazj+7LPPaNu2Lb179wZg5cqVlJeXc9pppwW2wjCVk6oxFyJktoWzHjdOT608CG9OhKs/h4QUy0oTEev4HS4yM42/nZx77rmG42DvDR9uzKt07tGYC4lV3cfACVfB0mer23avg4+neYKHiATV3LlzGTVqVL3XjBw5MkQVNSBczJo1K5h1RJzsNHPPhW6LSAwb+SBs+RZ21tgMadkL0GEI9BpvXV0iMWDIkCFs37696rhHjx5MmjSJSZMmVbW1aNEipDU1arbIEbt3767aqKxbt24hL95K2eYlwHVbRGJZfJJnE7N/DoGKkur2D26Eo/pBdufanikiTZScnExycjIAW7dupaCggMGDB9O6dWvLampUuCgpKeH666/nhRdeqNpjxOFwcOmll/L444+TkhL991mz00ybl+m2iMS6nC7wu7/DO7+vbisvhjcvgys/hfhk62oTaSCX28X+sv2W1pCVmIXd1rB5F8uXLwegX79+wSjJb40KFzfffDNffvkl77//PqeccgoAX3/9NVOnTuWWW26JiXUuzEuAF5VVcqjCSVK8w6KKRMJA7/Ng8wJY8VJ1247V8OFtMOYJ6+oSaaD9ZfsZ+vpQS2v48rwvaZ7UvEHPWbZsGXl5eWRnZ1e17d+/nxEjRlBZWUllZSU33HADV199daDLNWhUuHjrrbeYPXs2p556alXb6NGjSU5OZsKECTERLsybl4Fnxkhuln47kxg3ejr8uhT2bKhuW/4i5J0I/S6p/Xki0mTLli3z6rVIT09nwYIFpKSkUFJSQo8ePRg7dqwhgARao9a5KC0tpVWrVl7tLVu2pLS0tMlFRYKMpHgcpg1GNB1VBEhIhQkvQHyqsf3DW2H7KmtqEokRvsKFw+GoGq5QVlaG2+3G7XYHtY5GhYuTTz6Ze+65h0OHDlW1HTx4kD//+c+cfPLJASsunNntNq9bI9q8TOSwlsfAWTOMbZWH4I1L4OB+S0oSiXZ79uwhPz/f53iL/fv307t3b9q2bcttt91GTk5OUGtp1G2Rxx57jFGjRnktopWUlMTcuXMDWmA4y05NYHdRdaBQz4VIDT3HQf4SWPLP6rZ9P8M718L5r4C9Ub/biIREVmIWX573peU1NMSyZcsA34M5s7KyWLlyJTt37mTs2LGMGzfO5x2IQGlUuOjZsyc//vgjL7/8MuvXe5b5veCCC7jooouqpsPEAvNaF9q8TMRk5AOwbZlnDMYRP3wEC/8Og2+xri6Retht9gYPprTa8uXLadWqFbm5ubVe06pVK3r37s1XX33FuHHjglZLg8NFRUUFxxxzDB988EHQR5uGO/MqnVrrQsQkLgHGPw//HAylBdXtnz8AbXpDlxHW1SYSZaZNm8a0adO82nfu3ElKSgrp6ekcOHCABQsWcN111wW1lgb3S8bHxxvGWsQyr4W0NOZCxFvmUXDuc0CNAdBuF8y+Ago2WVaWSKz45ZdfGDx4ML1792bw4MFcf/319OzZM6jv2aibnpMnT+avf/0rlZWVga4nopino2rMhUgtOg+D4XcZ2w4dgNcugrJia2oSiREDBgxgxYoVrFy5klWrVnHNNdcE/T0bNeZi6dKlzJs3j08++YSePXuSmmqccvb2228HpLhw57V5mcKFSO0G3wI7VsHad6vbdq+DOdd5pq7abLU/V0QiSqPCRVZWlteuqLFIm5eJNIDNBmOegj0/wq611e3r3oOv/gZDbrWuNhEJqAaFC5fLxSOPPMIPP/xAeXk5w4cP5957742pGSI1eY+5UM+FSJ0S0+D8l+Ffp3puixzx+QPQuiccfbplpYlI4DRozMWDDz7IH//4R9LS0jjqqKOYMWMGkydPDlZtYc+8eVlpuZOD5U6LqhGJEM07wbh/g2FDJje8dRXsWm9ZWSISOA0KFy+88AJPPfUUc+fOZc6cObz//vu8/PLLVTujxhrzCp0ABbo1IlK/LiNgxL3GtrJCePU8KCnw+RQRiRwNChdbtmxh9OjRVccjRozAZrOxbdu2gBcWCTKS4oh3aH8RkUYZOBV6mBbx2fczvH4xVOr/kVgn2PtuhLtAfP4NCheVlZUkJSUZ2uLj46moqGhyIZHIZvPeX0TjLkT8ZLN5tmHPNS1VvGURfHATxPg3eAm9+Ph4gJjZgLM2Rz7/I1+PxmjQgE63283EiRNJTKwea3Do0CGuvfZaw3TUWJmKCpCdmsjOwupbIdq8TKQB4pPhglfhmeFQuLW6fcVL0KIbnDLVutok5jgcDrKysti1axcAKSkp2GJoirTb7aa0tJRdu3aRlZWFw+Fo9Gs1KFxcdtllXm0XX3xxo988GnhPR1XPhUiDpLf2BIx/j4KKGr8xfvonyO4Cx4yu/bkiAda6dWuAqoARi7Kysqq+Do3VoHAxa9asJr1ZNPKajqpwIdJwbXrD2H95xltUOTyD5PL/Qm5fy0qT2GKz2WjTpg0tW7aMyVv+8fHxTeqxOKJRi2hJNa/NyzTmQqRxjj0TTrsH5v25uq2iBF45D676DLLaWVebxByHwxGQH7KxqlF7i0g1r23XNRVVpPEG3QS9LzC2Fe+El8fDwf2WlCQiDadw0UTm2yIacyHSBDYbnDkDOgw2tu9erymqIhFE4aKJzKt06raISBPFJcB5L0JON2P7z1/Be9driqpIBFC4aCJft0VifQEWkSZLbgYXvQmpLY3tq16D+X+xpiYR8ZvCRROZb4scqnBRqv1FRJquWXu48HWITzG2L5gOS5+zpiYR8YvCRROZb4uAxl2IBMxR/WDcLNMmZ8B/b4G171pTk4jUS+GiiVITHCTEGb+MWqVTJIC6jYLRj5oaD6+BsfkrS0oSkbopXDSRzWYjRzNGRILrhCth6DRjm7McXrsQtq+ypiYRqZXCRQA0Nw/q1IwRkcA79Q7of7mxrawQXjoX9m62piYR8UnhIgCyTat07tFCWiKBZ7PBb//mWcmzppJd8OLZULjdkrJExJvCRQB4LaSlnguR4LA7YOyz0H6QsX3fz56AUVJgRVUiYqJwEQDea10oXIgETXwSXPAKtOppbN+9Hl4aC4cKralLRKooXASA1+ZlChciwZWUCRe/Bc07G9u3r/BsdFZe6vNpIhIaChcB4NVzoamoIsGX3goufRcy2hrbtyw6vA+J/h+KWEXhIgC0eZmIRbLy4LL3vJcJ3zQP3roSnBXW1CUS4xQuAsDX5mXaX0QkRLI7wyXvQFKWsX3d+/D278FZaUlZIrFM4SIAzD0X5U4XxWX6hiYSMq17eMZgJKQZ29e8DXOuBZf2+xEJJYWLADCPuQAtpCUScm2P92x0FpdsbF/9JsyZpIAhEkIKFwGQkhBHcrzD0KYZIyIW6DAILnwN4pKM7ateg/emgstlTV0iMUbhIkCap2rGiEhY6HQqnP8KOEw7Fq94Cd5XwBAJBYWLAMlJ04wRkbDR5TQ4/2VwmG5ZLn8R3p2sWyQiQaZwESBePRcKFyLW6vobmPAC2OON7StfgXeu0SwSkSBSuAgQX9NRRcRi3c6ACc97B4zVb2odDJEgUrgIEPN01ALtjCoSHo75re9bJGvnwJsToVK/CIgEmsJFgJino+7RgE6R8HH06XCBj1kk6z+A1y+CioPW1CUSpRQuAqRVhvGb1rb9hyyqRER86nIaXPgGxKcY23/8BF46V7upigSQwkWAtG1m/Ia1dd9BXC4tAS4SVjoN9b2S5y8L4fkzoaTAmrpEoozCRYDkNTOuCljudLFbt0ZEwk/7gYf3Isk0tm9fAbPOgMJtlpQlEk0ULgIkJy2RhDjjl/PXfaUWVSMidcobABP/C6ktjO17NsC/T4eCTdbUJRIlFC4CxG630TbL2HuRv1eDxETCVuuecMVcyMwztu/fAv8eBdtXWlOXSBRQuAigts2N4y7UcyES5rI7wxUfQ3ZXY3vJLpj1W/jpC0vKEol0ChcB1NY07uLXfeq5EAl7mW3h8o88PRk1lRfBS+Ng9Wxr6hKJYAoXAaRwIRKh0lp4xmB0GGxsd1V4VvL8ZqY1dYlEKIWLADJPR83XbRGRyJGUCRfNhu5jvM99fDt8+iftqCriJ4WLADJPR922/yBOrXUhEjnik2DcLBjwe+9zC/8Bb12h1TxF/KBwEUDmnosKp5tdRVqpUySi2B1wxnQ47U/e59a8A8+fBSV7Ql+XSARRuAignLQEEr3WutBvOSIRx2aDwbfAmKfAHmc89+sSeHYE7NloTW0iEUDhIoBsNpvXoM78vRp3IRKx+l4EF70JiRnG9n2b4bkR8Msia+oSCXMKFwGW57XWhXouRCJa5+GexbYy2hrbD+7z3CJZ/pI1dYmEMYWLAPOejqqeC5GI16o7XD0P2vQ2trsq4N3JMPdOcDmtqU0kDClcBJh5UKd6LkSiRHprz2JbR5/hfW7xE/Dq+dq2XeQwhYsA8xpzoZ4LkeiRkArnvwwnT/E+9+Mn8NxvYO9Poa9LJMwoXARYnqnnYvv+Q1Q6tfCOSNSwO+D0B2HMk2CPN57bvR6eGQ4b51lTm0iYULgIMHPPRaXLzc6iMouqEZGg6XsxXPY+pGQb2w/ug5fHwdd/B7cW0ZPYpHARYM1TE0iOdxjaftV0VJHo1P5kuHo+tDzO2O52wWf3wpsToazYispELKVwEWA+17rQoE6R6NWsPVw5F4490/vc2jmecRgFm0JeloiVFC6CwHutC/VciES1xHSY8OLhJcNtxnO71sK/hsG6DywpTcQKChdBoK3XRWLQkSXDL5rt2WG1prID8PpF8Mld4Kywpj6REFK4CAItpCUSw7qOgN9/AS27e59b9LhnVc/C7SEvSySUFC6CwLyQVv5e9VyIxJTmneDKT6HneO9zWxbBPwfDT1+Gvi6REFG4CALzWhc7CrXWhUjMSUyDsc/Ab/8GjgTjuZLd8MIY+PxBcFZaU59IEClcBIH5tojT5Wb7gUMWVSMilrHZ4ISr4IqPIbOd6aQbFkyH58+EA79aUp5IsChcBEFWSjypCaa1LjSoUyR2HdUfrvkSup7ufW7LInh6EKz/MPR1iQSJwkUQeNa6MI270KBOkdiW0hwueA1+cx/Y44znDu6D1y6AD2+DCv0iIpFP4SJI8pprOqqImNjtcMoNcPnHkGW+TQIs+ZdnTYwdq0Nfm0gAKVwEiffW6+q5EJHD8k6Aa76C7mO8z+1e59n8bOEMcGkguEQmhYsg0UJaIlKn5CwY/zz87u8Ql2Q85yyHT++GF8fAga2WlCfSFAoXQeIVLrR5mYiY2Wxw/BVwzQJo3cv7/OYFMPNkWPm6dliViKJwESTm2yI7Cg9RXqkuThHxoUU3uGqeZzyGeW+SQwfgnd/D6xdD8S5LyhNpKIWLIDEvpOVyww6tdSEitYlL8Mwkuex9yDjK+/z6D+Cpk2DNnJCXJtJQChdBkpEcR3qicbqZBnWKSL06DobrFkLPCd7nSgvgzctg9hVQsif0tYn4SeEiSGw2G22ba60LEWmE5GZw7jOebdxTcrzPf/8WPDkAVs/WWAwJSwoXQaQZIyLSJN3PgknfwLFnep8rLYC3roRXL4DCbaGvTaQOChdBpHAhIk2W1sLTgzH2GUjK8j7/w0fw5Inwv1laF0PChsJFEJlnjGzeU2JRJSIS0Ww26DUBJi/x3YtRVggf3AizzoBd60JenoiZwkUQdW2ZZjjesKMIp0v3R0WkkdJbwXkvwYQXILWl9/n8bzyboM27T3uUiKUULoKoe26G4fhghVO9FyLSdN3HwORvofcF3udclfDV3zzTVjfOC31tIihcBFVOWiKtMhINbWu3F1pUjYhElZTmcM7TcMk70Kyj9/l9P8NLY+H1S2B/fsjLk9imcBFk3dsYey/WbDtgUSUiEpU6D4dJi2HwLd5buQOse88zbfWr/4PK8tDXJzFJ4SLIzLdG1m5Tz4WIBFh8Mpz2J7j2a8g70ft8RSnM+7NnnxLdKpEQULgIsuNyMw3Ha7cV4taiNyISDC2Phcs/hjNnQHJz7/MFGz23Sl45Hwo2hb4+iRkKF0Fmvi1SUFLOrqIyi6oRkahnt0P/y+D67zw7rpo3QoPqtTE+uRsOqTdVAk/hIsjaNU8hzbTHiG6NiEjQpTSH3/0drv4cjurvfd5VAYtmwOP94LvnweUMfY0StRQugsxut3Fsm3RDmwZ1ikjIHNUPrvwMznocUlt4ny/ZDe9P9ayP8eNn2qtEAkLhIgTMt0Y0HVVEQspuh36XwvXLYOBUsMd7X7NrLbx8Lrx4DuxYHfoaJaooXISAZoyISFhIyoCR93sW4Dr6DN/X/DQfnh4M71wL+7eEtj6JGgoXIWCeMfJzQSlFhyosqkZEYl52Z7jwNbhkDrTq4eMCN6x8FR7vDx/fASV7Ql2hRDiFixDo0jKNOLtxxPb6HUUWVSMicljnYXDNAhjzFKS38T7vLIdvnoJ/9IYvHoYyfd8S/yhchEBSvIMupk3MdGtERMKC3QF9L/KMxxh2FySkeV9TXgxfPOQJGQv/AeWloa9TIorCRYiYx11oxoiIhJWEFBh6G0xdASde63vQZ2kBfPonT8j4ZiZUHAp5mRIZFC5CRDNGRCQipLWAM/7qWYSr1/n4XISrZBd8fDvM6ANLnlHIEC8KFyFi7rn4YUcxFU6XRdWIiNSjWXsY+0+4bmHtM0uKtsOHt3pCxjczoeJgSEuU8KVwESLHtTHOGCl3uti4q9iiakRE/NTqOM/Mkqs/hy4jfF9TtN3Tk/FYL1j0OJSXhLZGCTsKFyGSmRLPUVnJhjYN6hSRiHFUf7j4LbjiE+g41Pc1Jbvgk7vg78d5ZpeU7g1tjRI2FC5CyGsxLY27EJFI0+5EuOw9mPhh7SHj4D7P7JK/94C5d8KBraGtUSyncBFCx2nGiIhEiw6neELGFZ/UfrukogQWP+GZXTJnMuxaF9oaxTIKFyHkNWNkWyFubRIkIpGs3Yme2yVXfQ5Hj/J9jasCVrwET50EL4+HzQu0QVqUU7gIIfNtkcJDlWzdr9HVIhIF2vaHC1+H6xZBzwlgc/i+7sdP4Pkz4V9DYdUbUFke2jolJBQuQuiorGQyk40L06zRoE4RiSatjoNzn/Gsk3HCVeBI9H3d9pXw9tXwWE/4cjoU7w5tnRJUChchZLPZfN4aERGJOs07wm//Bjd9D0P+AMnNfF9XvAPmP+iZYTJnMmxfFdo6JSgULkLMfGvk+60a1CkiUSytJQy/E25aA6MfhWYdfV/nLPOMy/jnYHj2N7DqTagsC22tEjAKFyHW8yjjYlpLNu/VSp0iEv0SUmHA1Z7bJRNehHYDa7/21yXw9lWe3ox598H+/NDVKQGhcBFiAztnG46LyipZ9ss+i6oREQkxuwO6nwVXfOTZ7r33heBI8H1tyW746m/wj17w8gTY8BE4K0NbrzSKwkWItcxI8hp38eUPGsgkIjGoTW84Z6bnlsmpd0Baa9/XuV3w41x49XxP0Jj/EBz4NbS1SoMoXFhgaLcWhmOFCxGJaWkt4dTbPYM/x/8H2g+q/drCrfDlw57VP18cC2ve0diMMKRwYYFTjzaGizXbCtlVpC2LRSTGOeLhuHPg8v/CpG88U1kT0mu52A2b5sGbE+Fvx8BHt8OO70NZrdRB4cIC/do3Iz0xztC24Ic9FlUjIhKGWh7rmcp66wY463HI7Vf7tQf3wrcz4elTYOYgWPwkFO8KXa3iReHCAvEOO6d0yTG0fbFB/xFERLwkpEK/S+H38z0DQPtfDokZtV+/czXM/aOnN+PlCfD9W1BeGrp6BVC4sIx53MVXP+7B6dJa+yIitWrTG858DG7ZAGc/XffYDLfTMwh09hXwaFd451rYOE+zTUIkrv5LJBiGmsZdHDhYwcpf99OvXS2r2ImIiEdCCvS5wPMo2ATLX4JVr3sGe/pSXgwrX/U8Ult6xnX0OBfangB2/Y4dDPqqWiQ3K5mjW6UZ2r7YoFkjIiINkt0ZRtwDN66GS+ZAr/MhPqX260t2wZJ/wr9HevY1mXsnbF2mXVoDTOHCQubeC01JFRFpJLsDOg+Dsf+EW3+AMU9Bp1PBVsePucJfYfET8Mww+Edv+ORu+PU7BY0AsLndof8qFhYWkpmZyYEDB8jIqGNgTpT7+sc9XPzct1XHNht8d9dvaJ5ay2p1IiLSMEU7PIM6V70B21f495zMPDj2LDj2TMgb4AkuAvj/81vhwkJllU76/PlTDlY4q9r+cX4fxvQ5ysKqRESi1O4fYM3bnrCx5wf/npPaArqdAcecCR2HQHxScGsMcwoXEeLK/yxl3vrqaahj+x7F/53Xx7qCRESindsNO9dUB419P/v3vIQ06HIadBsNXUdCSvOglhmO/P35rdkiFhvarYUhXCz4cTculxu73WZhVSIiUcxmg9Y9PI/hd8OOVbD2XVgzB/Zuqv155cWe69a+6xnLkXcSdBsFR4+CnKM9ryuAei4st6WglCGPzDe0vT9lED3bZtbyDBERCQq3G3at9YSHdR/ArjX+Pzervac3o+tI6DDIM102CqnnIkK0y06hY04qm/eUVLV9+cMuhQsRkVCz2aDVcZ7HsD961tBY/1/PI/9boI7fxff/Akuf8TzikjwBo8sI6Hwa5HSNuV4N9VyEgXvfW8N/Fv1cdXxM63Q+umEwthj7xygiEraKdnpW/NzwEWyaD5UH/X9uZjvoMhw6D/cMCk2O3MUSNaAzgnz1424ueW6Joe3Na0/mhA6xN1hIRCTsVRyEzQtgw4fw46e1rwzqi80Obfp41uTodCrknQhxicGqNOAULiKI0+Vm2KNfsGVv9eY6v+3VhicvrGMXQBERsd6RmSc/fuJ55H8Lbpf/z49LgnYnQYfB0HEo5PYFR/iOWFC4iDDPfvUTD/x3XdVxnN3G19OG0zoztudUi4hElIP74KcvYONnno3SirY37PkJ6YfDximejdly+4AjPhiVNorCRYQ5UFrBSQ/NMyyoNfW0rtz8m6MtrEpERBrN7YZd6zxB46f58MsiqDzUsNeIT4V2J0L7gdBuIBzVD+KTg1OvHxQuItAdb6/m1SVbqo5z0hJZdPtwEuK0BYyISMSrOOS5bfLTF56wsW0Fdc5A8cWR4Ll10u5kTw9H3okhXcxL4SICrd9RyKjHvjK0aTlwEZEodXAf/LzQMzh08wLYva7+5/iS3dXTu5F3+JHdNWhbyStcRKjz/rmYbzfvrTru1y6LtyedYmFFIiISEkU74ZevPYHjl4Wwe33jXicx03P7pO0JnsdR/SE1OyAlahGtCHXZwA6GcLFsy35W/3pAi2qJiES79FbQ41zPA6B4tydk/LIItiyGnd/7NxOl7IDntstPNVZ/btbBEzKO6g+5/aBNL0hIDcqnAeq5CDsVTheD/zqfHYXVg37G9W/Lo+N7W1iViIhY7tAB+HUp/LIYtnwDW79r2GJeNdns0OIYz5obuX09s1Ja9ah32XLdFolgT3z+I49+Ur0dcEKcnW/uOI3mqQkWViUiImHFWeHZdC1/iSds5C+Bom2Nfz2bHXK6eXo12vT2PFr1gOSsqksULiLYnuIyBj70OeXO6u6v207vxuRhXSysSkREwt6BrbD1f/Dr4ce25Y3v3Tgiqx206gmte1KY1pnMAecpXESqm19fwdvLq5eUTU1w8OnNQ8nNsm5+s4iIRBhnpWcWytbvDj+WeXZ+bcgqojUUlrnJfLhI4SJSrf71AGc+8bWhbcSxLXnm0uO1oZmIiDReeSnsWA3bV3h6NratgD0b/Aoc/oYLzRYJUz3bZnLe8Xm8/r/8qrbP1u3io+93MLpnGwsrExGRiJaQ4lkXo92J1W3lJZ49UravrH7sWgeuika9hcJFGPvj6GOZt34Xe4rLqtrueW8Np3TJITM5fNaaFxGRCJeQCnkDPI8jKss9PRo7vvf0dOxcDT+vBIrqfTndFglz76/cxvWvLje0XXhiO/5yTk+LKhIRkVhVeOAAmVlZ9f781qYVYe53vdowrFsLQ9sr325hSY2FtkRERELCzzF/Chdhzmazcf/ZPUhJcBja73h7FWWVzlqeJSIiYh2FiwjQtlkKt4zsZmjbtLuEB/+7DgvuaomIiNRJ4SJCTBzYgV6m/UVeWPwLf6uxkqeIiEg4ULiIEA67jYfG9iTeYbzf9cT8jcz8YpNFVYmIiHhTuIggx+Vm8uj43l7jaf768XpeWPyzJTWJiIiYKVxEmDF9juIhH9NQ//TuGt6sseCWiIiIVRQuItD5A9px9++6e7VPe2sV/1qwCZdLgzxFRMQ6ChcR6spBHbnlN0cb2lxu+MuH67no2W/Ztr+Ju+CJiIg0ksJFBJsyvAvXDO3k1b74pwJGPbaA91Zus6AqERGJdQoXEcxms3H7qGO4fngXr0GehYcqmfrqcqa+upz8vaXWFCgiIjFJe4tEicWbCrjljRVsO3DI65zDbmNM71yuGdqZbq3TLahORESigb8/vxUuosiBgxX86d3veXdF7bdDRhzbkqsGd2JAh+bY7f6tES8iIgIKFzHt3RVbuWvO9xQdqqz1mpy0RH7TvSUjj2vNwM7ZJMY5ar1WREQEFC5i3s7CQzz95SZeW5LPwYq6NzhLTXAwsEsOffKy6N02i55tM8lMjg9RpSIiEikULgSAvSXlPL/oZ55f/DP7Syv8fl6nnFSObZNB++yUw49UOmSn0jI9UbdTRERilMKFGJSUVfLa0nxeX7qFH3YWN/p1HHYbOWkJtEhPpEVaIi3SE2mWkkBGcjwZSXGH/4wnJcFBSkIcyQkOkhMcpMQ7SIy3k+CwE+fQJCURkUjk78/vuBDWJBZKTYzjykEduXJQRzbvKeHTtTv4ZM1Ovtuyj4bES6fLzc7CMnYWljW6FrsNEuI8QSMhzk6c3U6cw0a8w47DbiPObsNusxHn8PzpsNtw2GzYbGC32bDbPX+CZzqu3QY2PG2eZs+fNjj855FrqfrzSJunwfBH1eua26rP+fd5qn9HRKJNWal/v5xaEi6OdJYUFhZa8fYxLzsBzu/TgvP7tGB30SEWbyrg+60HWLW1kA07iqhwuoL6/i6gEtDqGyIikcVV5vnOXd9ND0vCRVFREQB5eXlWvL2IiIg0QVFREZmZmbWet2TMhcvlYtu2baSnp1d1P4uIiEh4c7vdFBUVkZubi91e+/g5S8KFiIiIRC8N2xcREZGAUrgQERGRgFK4EBERkYBSuBAREZGAUrgQiVETJ07EZrN5PUaNGmV1aSIS4bRCp0gMGzVqFLNmzTK0JSYmBu39ysvLSUhICNrri0h4UM+FSAxLTEykdevWhkezZs0AzxLozz77LOeccw4pKSl07dqV9957z/D877//njPOOIO0tDRatWrFJZdcwp49e6rOn3rqqUyZMoUbb7yRnJwcTj/9dADee+89unbtSlJSEsOGDeP555/HZrOxf/9+SkpKyMjIYPbs2Yb3mjNnDqmpqVWL8IlI+FK4EJFa/fnPf2bChAmsWrWK0aNHc9FFF7F3714A9u/fz/Dhw+nbty//+9//+Pjjj9m5cycTJkwwvMbzzz9PQkICCxcu5Omnn2bz5s2MGzeOs88+m5UrV3LNNddw5513Vl2fmprK+eef79WjMmvWLMaNG0d6enrwP3ERaRq3iMSkyy67zO1wONypqamGx4MPPuh2u91uwH3XXXdVXV9cXOwG3B999JHb7Xa777//fvfIkSMNr5mfn+8G3Bs2bHC73W730KFD3X379jVcM23aNHePHj0MbXfeeacbcO/bt8/tdrvd3377rdvhcLi3bdvmdrvd7p07d7rj4uLcX3zxReC+ACISNBpzIRLDhg0bxsyZMw1tzZs3r/q4V69eVR+npqaSkZHBrl27AFi5ciXz588nLS3N63U3bdrE0UcfDUD//v0N5zZs2MAJJ5xgaBswYIDX8XHHHcfzzz/P7bffzksvvUT79u0ZMmRIIz5LEQk1hQuRGJaamkqXLl1qPR8fH284ttlsuFyeXXOLi4s588wz+etf/+r1vDZt2hjeozGuuuoqnnzySW6//XZmzZrF5Zdfrr2IRCKEwoWINEq/fv1466236NChA3Fx/n8r6datGx9++KGhbenSpV7XXXzxxfzhD39gxowZrF27lssuu6zJNYtIaGhAp0gMKysrY8eOHYZHzdkedZk8eTJ79+7lggsuYOnSpWzatIm5c+dy+eWX43Q6a33eNddcw/r165k2bRo//PADb7zxBv/5z38ADD0TzZo1Y+zYsdx2222MHDmStm3bNulzFZHQUbgQiWEff/wxbdq0MTwGDRrk13Nzc3NZuHAhTqeTkSNH0rNnT2688UaysrLq3Iq5Y8eOzJ49m7fffptevXoxc+bMqtki5jU2rrzySsrLy7niiisa/0mKSMhpy3URsdyDDz7I008/TX5+vqH9xRdf5KabbmLbtm1afEskgmjMhYiE3FNPPcUJJ5xAdnY2Cxcu5JFHHmHKlClV50tLS9m+fTsPP/ww11xzjYKFSITRbRERCbkff/yRMWPG0L17d+6//35uueUW7r333qrz06dP55hjjqF169bccccd1hUqIo2i2yIiIiISUOq5EBERkYBSuBAREZGAUrgQERGRgFK4EBERkYBSuBAREZGAUrgQERGRgFK4EBERkYBSuBAREZGAUrgQERGRgPp/UBEHr6SFUM8AAAAASUVORK5CYII=",
            "text/plain": [
              "<Figure size 640x480 with 1 Axes>"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "\n",
        "temperatures = [.5, 5, 2000]\n",
        "energies = np.linspace(0, 20, 100)\n",
        "fig, ax = plt.subplots()\n",
        "for i, T in enumerate(temperatures):\n",
        "    probabilities = np.exp(-energies/T)\n",
        "    Z = probabilities.sum()\n",
        "    probabilities /= Z\n",
        "    ax.plot(energies, probabilities, linewidth=3, label = \"$T_\" + str(i+1)+\"$\")\n",
        "ax.set_xlim(0, 20)\n",
        "ax.set_ylim(0, 1.2*probabilities.max())\n",
        "ax.set_xticks([])\n",
        "ax.set_yticks([])\n",
        "ax.set_xlabel('Energy')\n",
        "ax.set_ylabel('Probability')\n",
        "ax.legend()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "PCpbmGLha1DC"
      },
      "source": [
        "Here $T_1<T_2<T_3$. Notice that $T_1$ is a low temperature, and therefore it is highly peaked at low energy levels. In contrast, $T_3$ is a very high temperature and the probability distribution is almost completely flat."
      ]
    }
  ],
  "metadata": {
    "colab": {
      "include_colab_link": true,
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.12.2"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}