examples/braket_features/Noise_models/Noise_models_on_Rigetti.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "74154c03-d48f-4f41-867b-5e30254ce31a",
   "metadata": {},
   "source": [
    "# Noise models on Rigetti\n",
    "\n",
    "This notebook shows how to construct a noise model from device calibration data for Rigetti Aspen-M-3. We compare the measurement outcomes of circuits run on a noisy simulator with the same circuits run on quantum processing units (QPUs), to show that simulating circuits with noise models more closely mimics QPUs.\n",
    "\n",
    "**Before you begin**: We recommend being familiar with [Noise models on Amazon Braket.](https://github.com/amazon-braket/amazon-braket-examples/blob/main/examples/braket_features/Noise_models/Noise_models_on_Amazon_Braket.ipynb)\n",
    "Additionally, users should be familiar with [Running quantum circuits on QPU devices](https://github.com/amazon-braket/amazon-braket-examples/blob/main/examples/getting_started/2_Running_quantum_circuits_on_QPU_devices/2_Running_quantum_circuits_on_QPU_devices.ipynb). \n",
    "\n",
    "### Table of Contents\n",
    "\n",
    "- Noise model for Rigetti\n",
    "  - Loading device calibration data\n",
    "  - Comparing noisy simulator results to QPU results\n",
    "  - Smaller noise models compared to QPU results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "e5e15fcd-55a2-4168-acc6-a38100f94511",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "from braket.aws import AwsDevice\n",
    "from braket.circuits import Circuit, Gate\n",
    "from braket.circuits.noise_model import GateCriteria, NoiseModel, ObservableCriteria\n",
    "from braket.circuits.noises import AmplitudeDamping, BitFlip, Depolarizing, PhaseDamping, TwoQubitDepolarizing\n",
    "from braket.devices import Devices, LocalSimulator"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4da3e37-e93c-4273-98bc-d67b26b8c822",
   "metadata": {},
   "source": [
    "Braket provides access to hardware providers' reported calibration data. \n",
    "This can be used to construct noise models to approximate the behavior of the QPU when running circuits on a noisy simulator.\n",
    "In this tutorial, we focus on local noise models with no crosstalk interactions. Real devices can have crosstalk and unexpected effects that can further degrade the results.\n",
    "\n",
    "The Aspen-M-3 calibration data is available on the Braket devices page. Under qubit specs, the calibration data include the qubit index, with corresponding values for the $T_1$, $T_2$, fidelity from randomized benchmarking (fRB),  fidelity from simultaneous randomized benchmarking (fsRB),  readout fidelity (fRO), and active reset fidelity.\n",
    "Under \"edge specs\", the data includes the RB fidelity for C-Phase, XY, and CZ gates for each connected edge in the device topology."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "38001b96-faf7-49f0-b5a5-2d52820fcfc9",
   "metadata": {},
   "source": [
    "**One-qubit calibration data (Qubit specs)**\n",
    "<div align=center><img src=\"aspenm_qubit_specs.png\"/></div>\n",
    "\n",
    "\n",
    "**Two-qubit calibration data (Edge specs)**\n",
    "<div align=center><img src=\"aspenm_edge_specs.png\"/></div>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "276d680f-887d-4929-bcaa-e276a6303496",
   "metadata": {},
   "source": [
    "We can programmatically access all the calibration data with the Braket SDK. First we load the AwsDevice using the ARN for Rigetti Aspen-M-3."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "a55ea9d5-05f3-4da5-a305-a391e914d012",
   "metadata": {},
   "outputs": [],
   "source": [
    "aspen_m = AwsDevice(Devices.Rigetti.AspenM3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a8a42ae1-51ac-4bf1-bd6a-81ffd982b5f8",
   "metadata": {},
   "source": [
    "The properties dictionary contains keys \"provider\" and \"specs\", with \"1Q\" and \"2Q\" keys for the one- and two-qubit calibration data. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "153b5b63-8205-4a7e-bb18-0305ef428f36",
   "metadata": {},
   "outputs": [],
   "source": [
    "aspen_m_specs = aspen_m.properties.dict()[\"provider\"][\"specs\"]\n",
    "\n",
    "one_qubit_data = aspen_m_specs[\"1Q\"]\n",
    "two_qubit_data = aspen_m_specs[\"2Q\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "661cca17-01fa-40be-bd3e-741bacc74ce4",
   "metadata": {},
   "source": [
    "The keys of the \"1Q\" dictionary is the qubit number. \n",
    "For Aspen-M-3, there are 79 qubits indexed by the first digit representing the octagon (0 to 4, and 10 to 14)\n",
    "and the second representing the qubit in that octagon (0 to 7). We can get all qubit indices with `one_qubit_data.keys()` or with `aspen_m.topology_graph.nodes`.\n",
    "\n",
    "The keys of the \"2Q\" dictionary are the connected qubit pairs separated by a hyphen. \n",
    "For example, if qubit 0 and 1 are connected the key is \"0-1\".\n",
    "\n",
    "\n",
    "#### One-qubit noise \n",
    "\n",
    "Let's look at the one qubit calibration data for qubit 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "8d4ecc35-f4dc-49d8-b6ff-adaee560ac4b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'T1': 2.4443429874176914e-05,\n",
       " 'T2': 1.627681336532139e-05,\n",
       " 'f1QRB': 0.9985669301541129,\n",
       " 'f1QRB_std_err': 0.0004092835114091019,\n",
       " 'f1Q_simultaneous_RB': 0.9977676229871197,\n",
       " 'f1Q_simultaneous_RB_std_err': 0.00021522363993495725,\n",
       " 'fActiveReset': 0.9905000000000002,\n",
       " 'fRO': 0.9795}"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "one_qubit_data[\"0\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4b676762-e642-4a98-89e3-86d79367071a",
   "metadata": {},
   "source": [
    "For each qubit, there are various metrics of the quality:\n",
    "\n",
    "- **T1**: Thermal relaxation time is related to the time it takes for the excited state, |1⟩, to decay into the ground state, |0⟩. The probability of remaining in the excited state is $p(|1⟩)\\sim e^{-t/T_1}$\n",
    "\n",
    "- **T2**: The dephasing time, is the decay constant for the scale for a |+⟩ state to decohere into the completely mixed state. $p(|+⟩)\\sim e^{-t/T_2}$ \n",
    "\n",
    "- **Fidelity (RB)**: Single-qubit randomized benchmarking fidelities. RB fidelity quantifies the average gate fidelity where the average is over all Clifford gates. RB describes an *effective* noise model with gate-independent depolarizing noise on each Clifford gate.\n",
    "\n",
    "- **Fidelity (sRB)**: Single-qubit simultaneous randomized benchmarking fidelities. These are extracted by running single-qubit RB on all qubits simultaneously. Note that we expect the sRB fidelity to be lower than standard RB fidelity due to non-local crosstalk type noise on the device. \n",
    "\n",
    "- **Active reset fidelity**: Single-qubit active reset fidelities represents the accuracy to which qubits are reinitalized into the ground state |0⟩.\n",
    "\n",
    "- **Readout fidelity**: Single-qubit readout fidelities describes the probability of a bit flip error before readout in the computational basis. The readout fidelity is related to the probability of correctly measuring the ground state and excited states respectively, e.g. $f_{RO} =\\frac{p(0|0)+p(1|1)}{2}$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "296708a3-771e-42b4-9ba8-f3db1a02970b",
   "metadata": {},
   "source": [
    "Now that we know how to extract and use the calibration data, we can build a simple noise model. For every qubit we will add:\n",
    "- amplitude dampening noise with probability $p= 1-e^{-t/T_1}$ for every gate\n",
    "- phase dampening noise with probability $p= 0.5(1-e^{-t/T_2})$ for every gate\n",
    "- depolarizing noise with probability $p=1-f_{sRB}$ (from simultaneous RB fidelity) for  every gate\n",
    "- readout bit flip noise with probability $p=1-f_{RO}$ to measurements \n",
    "\n",
    "Technically, the sRB fidelity already includes effects from $T_1$/$T_2$, however to be explicit we add these as separate terms. In a sense, this model might overestimate the noise on the QPU.  \n",
    "\n",
    "For the dampening noises, we first need the gate times to model $T_1$ and $T_2$ errors.\n",
    "From the Braket Aspen-M-3 device page: \"These gates offer fast (40ns and 180ns) 1Q and 2Q gate times.\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "cd27adfb-f7d6-4583-b854-794e38152c90",
   "metadata": {},
   "outputs": [],
   "source": [
    "# median 1q and 2q gate times\n",
    "gate_time_1_qubit = 40e-9\n",
    "gate_time_2_qubit = 240e-09"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5769f2f4-de1c-4b47-9e4e-fd8d39b9e975",
   "metadata": {},
   "source": [
    "To create the noise model, we iterate over all qubits keys in `one_qubit_data`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "508070bc-b1d2-4737-a2fc-d7711e419996",
   "metadata": {},
   "outputs": [],
   "source": [
    "noise_model = NoiseModel()\n",
    "\n",
    "# Single-qubit noise\n",
    "for q, data in one_qubit_data.items():  # iterate over qubits\n",
    "\n",
    "    # T1 dampening\n",
    "    t1 = data[\"T1\"]\n",
    "    damping_prob = 1 - np.exp(-(gate_time_1_qubit / t1))\n",
    "    noise_model.add_noise(AmplitudeDamping(damping_prob), GateCriteria(qubits=int(q)))\n",
    "\n",
    "    # T2 phase flip\n",
    "    t2 = data[\"T2\"]\n",
    "    dephasing_prob = 0.5 * (1 - np.exp(-(gate_time_1_qubit / t2)))\n",
    "    noise_model.add_noise(PhaseDamping(dephasing_prob), GateCriteria(qubits=int(q)))\n",
    "\n",
    "    # 1q RB depolarizing rate from simultaneous RB\n",
    "    depolar_rate = 1 - data[\"f1Q_simultaneous_RB\"]\n",
    "    noise_model.add_noise(Depolarizing(depolar_rate), GateCriteria(qubits=int(q)))\n",
    "\n",
    "    # 1q RB readout\n",
    "    readout_value = 1 - data[\"fRO\"]\n",
    "    noise_model.add_noise(BitFlip(readout_value), ObservableCriteria(qubits=int(q)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "9d20d24c-17e2-4bba-b317-fb2d156bedb0",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of terms in noise model is: 320\n",
      "Number of parameters in noise model is: 320\n"
     ]
    }
   ],
   "source": [
    "num_params = sum(len(item.noise.parameters) for item in noise_model.instructions)\n",
    "print(f\"Number of terms in noise model is: {len(noise_model.instructions)}\")\n",
    "print(f\"Number of parameters in noise model is: {num_params}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "790c62f7-776a-4e4c-9405-7b7a14eee5b0",
   "metadata": {},
   "source": [
    "#### Two-qubit noise \n",
    "Next we consider adding two-qubit noise to the model. \n",
    "\n",
    "Let's first look at the data provided in the Aspen-M-3 device calibration data. On the first connect, \"0-1\", the properties are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "626cbe56-3649-48eb-a2a1-f669eb790b47",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'fCPHASE': 0.8991973184700036,\n",
       " 'fCPHASE_std_err': 0.019391126441846068,\n",
       " 'fCZ': 0.9406104552003762,\n",
       " 'fCZ_std_err': 0.011069113247524152,\n",
       " 'fXY': 0.8365560014004874,\n",
       " 'fXY_std_err': 0.014585224789608713}"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "two_qubit_data[\"0-1\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4bcd3be4-a6e4-4cf3-841e-1a451e01803a",
   "metadata": {},
   "source": [
    "Here, we see the fidelity per gate (CPhase, CZ, or XY) and the associated standard error. \n",
    "\n",
    "Next we loop over the entries in the `two_qubit_data` dictionary and add two-qubit depolarizing noise to the model. Notice that Aspen-M-3 has symmetric connections (\"0-1\" and \"1-0\") so we need to add noise in both directions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "80f57550-6549-46b6-9a6d-ff7f0f03e90d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Two-qubit noise\n",
    "for pair, data in two_qubit_data.items():  # iterate over qubit connections\n",
    "\n",
    "    # parse strings \"0-1\" to integers [0, 1]\n",
    "    q0, q1 = (int(s) for s in pair.split(\"-\"))\n",
    "\n",
    "    if \"fCPHASE\" in data:\n",
    "        phase_rate = 1 - data[\"fCPHASE\"]\n",
    "        noise_model.add_noise(\n",
    "            TwoQubitDepolarizing(phase_rate),\n",
    "            GateCriteria(\n",
    "                Gate.CPhaseShift, [(q0, q1), (q1, q0)]\n",
    "            ),  # symmetric connections\n",
    "        )\n",
    "\n",
    "    if \"fXY\" in data:\n",
    "        xy_rate = 1 - data[\"fXY\"]\n",
    "        noise_model.add_noise(\n",
    "            TwoQubitDepolarizing(xy_rate), GateCriteria(Gate.XY, [(q0, q1), (q1, q0)])\n",
    "        )\n",
    "\n",
    "    if \"fCZ\" in data:\n",
    "        cz_rate = 1 - data[\"fCZ\"]\n",
    "        noise_model.add_noise(\n",
    "            TwoQubitDepolarizing(cz_rate), GateCriteria(Gate.CZ, [(q0, q1), (q1, q0)])\n",
    "        )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "82376d08-84bc-4249-9c8d-d7aacdb73113",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of terms in noise model is: 612\n",
      "Number of parameters in noise model is: 612\n"
     ]
    }
   ],
   "source": [
    "num_params = sum(len(item.noise.parameters) for item in noise_model.instructions)\n",
    "print(f\"Number of terms in noise model is: {len(noise_model.instructions)}\")\n",
    "print(f\"Number of parameters in noise model is: {num_params}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5913a8a2-a9de-4e32-a35d-c48191891b1f",
   "metadata": {},
   "source": [
    "### Compare circuits run on device vs simulator with a noise model\n",
    "\n",
    "Let's just look at the first 5 qubits. Note that to ensure the noise model applied T1 and T2 noise during the time between gate, we manually add identity gates to each moment."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "6fabbd05-b480-43c6-b43a-d396d81bd687",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "T  : |   0    |   1    |2|\n",
      "                          \n",
      "q0 : -Rx(0.50)-Rx(3.14)-C-\n",
      "                        | \n",
      "q1 : -Rz(0.50)-Rx(3.14)-Z-\n",
      "                          \n",
      "q2 : -Rz(0.50)-Rx(3.14)---\n",
      "\n",
      "T  : |   0    |   1    |2|\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(42)\n",
    "\n",
    "circ = Circuit().rx(0, 0.5).rz(1, 0.5).rz(2, 0.5).rx(0, np.pi).rx(1, np.pi).rx(2, np.pi).cz(0, 1)\n",
    "print(circ)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "e2af07ed-33de-4709-9e7c-5e3c3caaa751",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "T  : |                     0                      |                     1                      |      2      |\n",
      "                                                                                                              \n",
      "q0 : -Rx(0.50)-AD(0.0016)--PD(0.0012)-DEPO(0.0022)-Rx(3.14)-AD(0.0016)--PD(0.0012)-DEPO(0.0022)-C-DEPO(0.059)-\n",
      "                                                                                                | |           \n",
      "q1 : -Rz(0.50)-AD(0.00086)-PD(0.0015)-DEPO(0.025)--Rx(3.14)-AD(0.00086)-PD(0.0015)-DEPO(0.025)--Z-DEPO(0.059)-\n",
      "                                                                                                              \n",
      "q2 : -Rz(0.50)-AD(0.0018)--PD(0.0011)-DEPO(0.0041)-Rx(3.14)-AD(0.0018)--PD(0.0011)-DEPO(0.0041)---------------\n",
      "\n",
      "T  : |                     0                      |                     1                      |      2      |\n"
     ]
    }
   ],
   "source": [
    "noisy_circ = noise_model.apply(circ)\n",
    "\n",
    "print(noisy_circ)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "9304dcf1-a1e8-4e41-a6b5-769841f4eefe",
   "metadata": {},
   "outputs": [],
   "source": [
    "simulator = LocalSimulator()  # noise free simulator\n",
    "task = simulator.run(circ, shots=100_000)\n",
    "free_probs = task.result().measurement_probabilities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "0e1d9cfe-03d6-406b-b4ab-89166a907681",
   "metadata": {},
   "outputs": [],
   "source": [
    "noisy_simulator = LocalSimulator(\"braket_dm\")\n",
    "noisy_task = noisy_simulator.run(noisy_circ, shots=100_000)\n",
    "noisy_probs = noisy_task.result().measurement_probabilities"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b3f370bd-a8b3-4d5b-a724-da549fc4c81f",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Note:</b> The below section runs tasks on the Rigetti Aspen-M-3 device. When you run this notebook, make sure the device is currently available. You can find QPU availability windows on the <a href=\"https://us-west-1.console.aws.amazon.com/braket/home?region=us-west-1#/devices\">Devices page</a> in the Amazon Braket Console.\n",
    "</div>\n",
    "\n",
    "<div class=\"alert alert-block alert-info\">\n",
    "<b>Note:</b> Running the circuit below will result in charges on your AWS account.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "c8f4664a-0527-46fb-87f2-2a0c85f0bd79",
   "metadata": {},
   "outputs": [],
   "source": [
    "aspen_m_task = aspen_m.run(circ, shots=100_000, disable_qubit_rewiring=True)\n",
    "aspen_m_result = aspen_m_task.result()\n",
    "aspen_m_probs = aspen_m_result.measurement_probabilities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "1676a45c-5ef8-49cc-b227-69f69fa2004a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>aspen-m</th>\n",
       "      <th>noisy_sim</th>\n",
       "      <th>free_sim</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>000</th>\n",
       "      <td>0.00136</td>\n",
       "      <td>0.00016</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>111</th>\n",
       "      <td>0.68340</td>\n",
       "      <td>0.85721</td>\n",
       "      <td>0.93927</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>110</th>\n",
       "      <td>0.07619</td>\n",
       "      <td>0.00635</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>011</th>\n",
       "      <td>0.12967</td>\n",
       "      <td>0.07411</td>\n",
       "      <td>0.06073</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>101</th>\n",
       "      <td>0.07366</td>\n",
       "      <td>0.04401</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>100</th>\n",
       "      <td>0.00777</td>\n",
       "      <td>0.00023</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>010</th>\n",
       "      <td>0.01303</td>\n",
       "      <td>0.00045</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>001</th>\n",
       "      <td>0.01492</td>\n",
       "      <td>0.01748</td>\n",
       "      <td>NaN</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     aspen-m    noisy_sim  free_sim\n",
       "000    0.00136    0.00016       NaN\n",
       "111    0.68340    0.85721   0.93927\n",
       "110    0.07619    0.00635       NaN\n",
       "011    0.12967    0.07411   0.06073\n",
       "101    0.07366    0.04401       NaN\n",
       "100    0.00777    0.00023       NaN\n",
       "010    0.01303    0.00045       NaN\n",
       "001    0.01492    0.01748       NaN"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "free_sim = pd.DataFrame.from_dict(free_probs, orient=\"index\").rename(\n",
    "    columns={0: \"free_sim\"}\n",
    ")\n",
    "noisy_sim = pd.DataFrame.from_dict(noisy_probs, orient=\"index\").rename(\n",
    "    columns={0: \"noisy_sim\"}\n",
    ")\n",
    "aspenm = pd.DataFrame.from_dict(aspen_m_probs, orient=\"index\").rename(\n",
    "    columns={0: \"aspen-m\"}\n",
    ")\n",
    "df = aspenm.join(noisy_sim).join(free_sim)\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "472e6850-2603-4174-bf92-b04eb0c38c2f",
   "metadata": {},
   "source": [
    "We can compute the fidelity between the free simulation and Rigetti, as well as the noisy simulation and Rigetti. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "5b7d1740-cb46-433e-8c2b-7e8b7c10fb4e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Total fidelity between Aspen-M and noise-free simulator is 0.8899252242353248\n",
      "\n",
      "Total fidelity between Aspen-M and noisy simulator is 0.9627231692802513\n"
     ]
    }
   ],
   "source": [
    "def fidelity(p, q):\n",
    "    return np.sum(np.sqrt(p * q))\n",
    "\n",
    "\n",
    "f_free_aspen = fidelity(df[\"free_sim\"], df[\"aspen-m\"])\n",
    "f_noisy_aspen = fidelity(df[\"noisy_sim\"], df[\"aspen-m\"])\n",
    "\n",
    "print(f\"\\nTotal fidelity between Aspen-M and noise-free simulator is {f_free_aspen}\")\n",
    "print(f\"\\nTotal fidelity between Aspen-M and noisy simulator is {f_noisy_aspen}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fd4d88b7-015e-445b-ba26-1bdbc6a2bd16",
   "metadata": {},
   "source": [
    "To better visualize, we can also plot the output probability distributions from each circuit:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "9329ed72-2f30-4383-b48f-0d50e1a5dd11",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsIAAAF+CAYAAACI8nxKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAABHBklEQVR4nO3dd5hU5fnG8ftZwF167x1BOgKhCQqrIkUEDGoUCQiKirGAMSooEsSoYEURgyY/qgVULIiIRCJWECkrCBsiEKqAiLI0qfv+/jhnx9lemN1lOd/Pde117Zz6zDtnZu555z1nzDknAAAAIGii8rsAAAAAID8QhAEAABBIBGEAAAAEEkEYAAAAgUQQBgAAQCARhAEAABBIBGEgD5jZADNblN91pGRmh8ysXoS3aWY2zcx+MbPlkdx2Nut4wMz+mUvbXmJmQ3Nj23kpL49LM7vIzDbk0rbHmtkrubHt/GZm083sb/ldB3C2IgijQDGz681shR/gdpnZh2Z2YX7XlRnn3KvOuW75XUdKzrkSzrnNEd7shZIuk1TDOdcuwtvOMufcY865fA+rZrbFzLrmdx1pycvj0jn3uXOuYV7sKyNnY7A0s1gzc2Z2f37Xkhb/g6Mzs/NTTH/Hnx6bwbovm9kGM0s0s8G5XCoCiCCMAsPM/ixpoqTHJFWWVEvSi5L65mNZmTKzwvldQx6rLWmLc+5wWjMD2B45RludecysUH7XkIYbJP0saVB+F5KB/yqsPjMrL+kCSXszWe9bSX+StCr3SkOQEYRRIJhZaUnjJN3unHvbOXfYOXfCOfe+c+5ef5loM5toZj/4fxPNLNqfF2tmO8zsPjP70e9NvtLMLjez/5rZz2b2QNj+xprZW2Y2x8wOmtmq8N4MMxtpZpv8eevN7Pdh8wab2Zdm9qyZ7ZM01p/2RdgyzsyGmdn3ZrbfzCabmfnzCpnZ02b2k5n9z8zu8JdPMxT5PY5/MbM1Zpbg1xwTNv9mM9vo38d5ZlYtRR31/f8v9+/LQTPbaWZ/CVvuCjOL82v9ysxapFPLTZL+KekCv9f+4bC2v9/MdkuaZmZRYW24z8zeMLNyYdvp4O9nv5l9m0mP0f1+vQf9nqNLwx7DV/z/6/j3dYiZbTdv2MYwM2vrt9t+M3shxeP/StjtpPVTPQZmdq6Z/du/Hz+Z2atmVsafN0veB7b3/fa4z5/ex8zW+ftdYmaNUzye95vZGkmH03vcw5ZPat97wo7tIWHzS5vZTDPba2ZbzWy0mUX580LHpXme9bdxwMzWmlkzv432WFgANLN+ZvZtOvWkeRwl1Znift7rt/9hM/s/M6ts3rc8B83sYzMrm9a6Yeun2dNuZm+a2W7zng+fmVlTf/otkgZIus9/PN73pzf2H4f9/uPSJ2xb083s72a2wMwOS7o4o8fDX+c5/zg7YGYrzeyisHlj/eN9pn8/15lZm7D5rcx7vTloZnMkxaS5k9+WLy7pakm3S2qQYlsxZvaKf2zuN7NvzKyyP2+JmT1uZsv9Ot+zLD4H/XUfMe917qCZLTKzCpk0y6uSrg07jvpLekfS8YxWcs5Nds4tlnQ0k+0DOeOc44+/M/5PUg9JJyUVzmCZcZKWSaokqaKkryQ94s+L9dcfI6mIpJvl9US8JqmkpKaSfpVU119+rKQT8t5gikj6i6T/SSriz79GUjV5HyavlXRYUlV/3mB/X3dKKiypqD/ti7BanaT5ksrIC0p7JfXw5w2TtF5SDUllJX3sL5/mfZe0RdJyv55ykuIlDfPnXSLpJ0mtJUVLmiTpsxR11Pf/3yXpIv//spJa+/+3kvSjpPaSCsnrfdoiKTqdelLe16S2n+DXUFTScP+xquFPe0nS6/7y1SXtk3S5376X+bcrprGvhpK2S6rm364j6dywx/CVsOlO0hR5waKbvDfWd+UdL9X9+9gl5bop1i/s314iaaj/f32/xmh5x91nkiameHy6ht0+T97xcpm8Y+s+SRslnRO2fJykmpKKZuG5kdS+4/ztXS7piKSy/vyZkt6Td5zXkdczd1PKx0pSd0kr5R2TJqmxfjum10vqGbbPdyTdk0496R1HsZJ2pGiXZfK+3Ulq/1XyjrcYSf+W9Ne01k3Zrmk8Xjf69zda3rdIcWHzpkv6W9jtIn77PyDpHHnPmYOSGoYtnyCpk7zjMSYLj8kfJZWX9/y/R9LupPX8Wo/6j1MhSY9LWubPO0fSVkl3+3VdLe916G8Z7Gug3+aFJL0vaVLYvFv9acX8+b+TVCrsGN4pqZmk4pLm6rfnS4bPQX/dTfKO5aL+7fEZ1LhE0lBJi+QfR/Jesy6QtENSbBba9AtJgzNbjj/+svtHjzAKivKSfnLOncxgmQGSxjnnfnTO7ZX0sLw3iSQnJD3qnDshabakCpKec84ddM6tk/dmHz6GbaVz7i1/+WfkvTl3kCTn3JvOuR+cc4nOuTmSvpcUPh72B+fcJOfcSefcr+nUO945t985t03SJ5Ja+tP/4Ne1wzn3i6TxGTeNJOl5v56f5b3xJW1rgKSpzrlVzrljkkbJ662tk8Y2TkhqYmalnHO/OOeSvoq8RdJLzrmvnXOnnHMzJB1LaossSpQXao757TFM0oP+fTwmLxxc7fd+/lHSAufcAr99/yVphbw35ZROyQs7TcysiHNui3NuUwZ1POKcO+qcWyQvjL7uHy87JX0uL4Rli3Nuo3PuX/592yvvWOmSwSrXSvrAX+eEpKfkhYmOYcs875zbnsGxk9IJecf+CefcAkmHJDX0e9+ukzTKP863SHpayZ8X4dsoKamRJHPOxTvndvnzZsh7XOT3GnaX9yEyvVrSOo7SMsk5tyes/b92zq12zh2VF7az/XhIknNuqn9/k46t8837ViktHSSVkPd8PO6c+7e8D6n9w5Z5zzn3pX88Ztoz6Zx7xTm3z3/+Py3vGA0fH/2Ff3yfkjRLv73udJAXgCf6j+Vbkr7JZHc3SJrjb+s1SdeZWRF/3gl5r531/efuSufcgbB1ZznnvnPeMKaHJP3BP2ay8hyc5pz7r3+MvqHfXnMyMlPSIDNrJKmMc25pFtYBchVBGAXFPkkVMvmauJq83pQkW/1poW34bxaS1/srSXvC5v8q7w0xyfakf5xzifJ6LqpJkpkNst+GCuyX16tSIa11M7A77P8jYfuulmL9091WqE2cc4fktWX1NLZxlbw3uq1m9qmZXeBPry3pnqT76t/fmpKqmXfVgUP+34cZ1Lc3RYCoLemdsO3Fywu1lf1516TY34WSqqbcqHNuo6QR8sLOj2Y228KGfqQh5eOd0eOfJf7X+bP9YQAHJL2i5MdCSikfk0R5j3H4Y5KVxzzcvhQfEpOOgQryglXK50Wqx98PgC9ImiyvLV82s1L+7Fck9fa/hv+DpM/DQnJK6R1HacmNx6OQmY03b9jNAXk9x1L6j0k1Sdv9xyFJyjbK1uNh3lCleH9oxn5JpVPsP+XzNcZ/basmaadzzqWoJb391JQ3VONVf9J78j6w9/Jvz5L0kaTZ5g0XeyIsJKe8X1vlHSsVlLXnYJqvOWY2Jew14QEl97a8Hvc7/NpS3p9DYX+10rvfQCQRhFFQLJXXC3llBsv8IO8FPEktf1pO1Uz6x7wxlTUk/WBmtSX9Q96LeXnnXBlJ38n7OjlJ+BtZdu3y95WqjhxI1iZ+kCkv7yvRZJxz3zjn+sobKvCuvF4eyXuzfNQ5Vybsr5hz7nXnXXWghP/XM4M6UrbHdnlfkYZvM8bvGdwur6cqfF5x51yaPePOudeccxf699PJG4Jxug7L+zo5SZUMln3M329z51wpeb1pGR0LKR8Tk/cY78xgnZz6SV6vYMrnRarHX5Kcc887534nqYm8r73v9afvlPcc7CevNzlViAnbRnrH0elI9nj4vZYV01n2enkn0HaVF0DrJK2WVGKK5X+QVNN/jidJ2UZZfjz88cD3yfvAUNZ/fUhQ8mMiPbskVfePifBa0jNQ3vv4++aNv98sLwjfIEl+r/LDzrkm8r5xuELJT6gLf22pJe9Y+UnZfA6Gc84NC3tNeCzFvCOSPpR0m9I4hsLWK+G8b8qAXEcQRoHgnEuQN753snknuRUzsyJm1tPMnvAXe13SaDOr6J+4MUZeT1ZO/c68k4IKy+t1PCZvTGNxeW+MeyXJvBOTmp3GflJ6Q9JwM6tu3klXp3NJpNclDTGzluadOPiYvK+ft4QvZGbn+L27pf2v6w/IG84geaF/mJm1N09xM+tlZiVPo64pkh71P1TIf8ySrv6R1PvY3e/dizHvZKkaKTdiZg3N7BL/vh2V14uYmHK5HIiT1NnMavlfqY/KYNmS8oYiJJhZdfnhMcweSeHXan5DUi8zu9TvnbtH3rH1VXo7MO+ErenZvRP+NyBvyGvrkn57/1lpPC/MOymuvV/TYXntGd6WM+UFvObyevbSqjOj4+h0/Fder2kvv77R8oYbpKWkvPbcJy88P5ZifsrH42t5PZr3+a8psZJ6yxs+lSbzTjLcksH+T8p7fShsZmMklUpn2ZSW+uve5dfST8mHXKV0g7whYC3D/q6SdLmZlTezi82suf/B4YC8oBv+ePzRzJqYWTF5Y8zf8o+ZLD8Hc+ABeWPxt2RlYf+YipH3QaKIXwvZBRHDwYQCwx9r92d5b4J75fVa3CGv10mS/iZvHNsaSWvlnXhzOtcLfU/eeM5f5PW89PN7WNbLG2e5VN6banNJX57GflL6h7yTStZIWi1pgbw3x1MZrZQW59zH8sb+zZXX23SuvDGjaRkoaYv/dfIweeOL5ZxbIe/kwhfktcVGeSdZnY7nJM2TtMjMDsr7gNHe3992eT16D+i3x/lepf16FS1vDPVP8r6qraSMQ2uW+GMi58h7DFbKGzOanoflnYyYIOkDpQ6Jj8v7gLbfzP7inNsgr9d4kl93b0m9nXMZnT1fUzk/xu6UF2w3yzvh6DVJU9NYrpS8Y+8XeV+T75P0ZNj8d+QPafF79tKT5nF0OvwPwn+Sd0WSnfLuz450Fp8pr/6d8sb9L0sx///kjWHeb2bv+u3eW1JPeY/Hi5IGOef+k0FJGT0eH0laKC+8b5X3gSJLQyv8WvrJe379LO/1J70PHR3kPR6TnXO7w/7myXuO9pf3TcZb8kJwvKRPlbwndpa8kwF3y+tJvsuvIzvPwWxx3rkMX2S+ZMgieR9wO0p62f+/8+nWASSx5EORAEjeJY7knWDyxzOglp6Spjjname6MM46ZnaOvGuptvB7WfOzlk2SbvU/YAWWeb/GN9w5F5/fteSUmS2Rd5WIXPn1RaCg4GLtwBnGzIrKOwFmkbyTx/4qrzcOAeT3EjbOdMFcZmZXyRsS9O/8riW/uTPwVyIB5AxBGDjzmLyv2+fI+xrwA3njnYF84fceNpE0MMXVFQCgQGNoBAAAAAKJk+UAAAAQSARhAAAABFK+jRGuUKGCq1OnTn7tHgAAAAGxcuXKn5xzqX6IJ9+CcJ06dbRixYr82j0AAAACwszS/LlyhkYAAAAgkAjCAAAACCSCMAAAAAKJH9QAAADIAydOnNCOHTt09OjR/C7lrBUTE6MaNWqoSJEiWVqeIAwAAJAHduzYoZIlS6pOnToys/wu56zjnNO+ffu0Y8cO1a1bN0vrMDQCAAAgDxw9elTly5cnBOcSM1P58uWz1eNOEAYAAMgjhODcld32JQgDAADgjHfjjTeqUqVKatasWcS2yRhhAACAfFBn5AcR3d6W8b0iur0zzeDBg3XHHXdo0KBBEdsmPcIAAAABcuWVV+p3v/udmjZtqpdfflmnTp3S4MGD1axZMzVv3lzPPvusJCk2NlbDhw9Xy5Yt1axZMy1fvlySdPjwYd14441q166dWrVqpffee0+SNH36dPXr1089evRQgwYNdN9996W5/+nTp+vKK6/UZZddpjp16uiFF17QM888o1atWqlDhw76+eef01yvc+fOKleuXETbgh5hAACAAJk6darKlSunX3/9VW3bttXvfvc77dy5U999950kaf/+/aFljxw5ori4OH322We68cYb9d133+nRRx/VJZdcoqlTp2r//v1q166dunbtKkmKi4vT6tWrFR0drYYNG+rOO+9UzZo1U9Xw3XffafXq1Tp69Kjq16+vCRMmaPXq1br77rs1c+ZMjRgxIi+agh5hAACAIHn++ed1/vnnq0OHDtq+fbuOHz+uzZs3684779TChQtVqlSp0LL9+/eX5PXGHjhwQPv379eiRYs0fvx4tWzZUrGxsTp69Ki2bdsmSbr00ktVunRpxcTEqEmTJtq6dWuaNVx88cUqWbKkKlasqNKlS6t3796SpObNm2vLli252wBhCMIAAAABsWTJEn388cdaunSpvv32W7Vq1UrHjh3Tt99+q9jYWE2ZMkVDhw4NLZ/yKgxmJuec5s6dq7i4OMXFxWnbtm1q3LixJCk6Ojq0bKFChXTy5Em98847atmypVq2bKkVK1akWi4qKip0OyoqSidPntT27dtD60yZMiXX2oOhEThjNJ/RPEfrrb1hbYQrAQDg7JSQkKCyZcuqWLFi+s9//qNly5bpp59+UmJioq666io1bNhQf/zjH0PLz5kzRxdffLG++OILlS5dWqVLl1b37t01adIkTZo0SWam1atXq1WrVunu8/e//71+//vfh24nDcHISM2aNRUXF3da9zUrCMIAAAAB0aNHD02ZMkWNGzdWw4YN1aFDB+3cuVOxsbFKTEyUJD3++OOh5WNiYtSqVSudOHFCU6dOlSQ99NBDGjFihFq0aKHExETVrVtX8+fPz/Xa+/fvryVLluinn35SjRo19PDDD+umm246rW2acy5C5WVPmzZtXFL3OCDRIwwAOLvFx8eHhhAUBLGxsXrqqafUpk2b/C4lW9JqZzNb6ZxLdUcYIwwAAIBAYmgEAAAAUlmyZEl+l5Dr6BEGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIHHVCAAAgPwwtnSEt5cQ2e1JGjNmjDp37qyuXbtGfNuSNG/ePK1fv14jR47Mle1nhiAMAACANI0bNy5Xt9+nTx/16dMnV/eREYZGAAAABMSWLVvUuHFj3XzzzWratKm6deumX3/9VXFxcerQoYNatGih3//+9/rll18kSYMHD9Zbb70lSRo5cqSaNGmiFi1a6C9/+YsOHjyounXr6sSJE5KkAwcOJLud0vPPPx9a/7rrrpMkTZ8+XXfccUdoX7fddps6dOigevXqacmSJbrxxhvVuHFjDR48OFfagyAMAAAQIN9//71uv/12rVu3TmXKlNHcuXM1aNAgTZgwQWvWrFHz5s318MMPJ1tn3759euedd7Ru3TqtWbNGo0ePVsmSJRUbG6sPPvhAkjR79mz169dPRYoUSXO/48eP1+rVq7VmzRpNmTIlzWV++eUXLV26VM8++6z69Omju+++W+vWrdPatWsVFxcX0XaQCMIAAACBUrduXbVs2VKS9Lvf/U6bNm3S/v371aVLF0nSDTfcoM8++yzZOqVLl1ZMTIxuuukmvf322ypWrJgkaejQoZo2bZokadq0aRoyZEi6+23RooUGDBigV155RYULpz06t3fv3jIzNW/eXJUrV1bz5s0VFRWlpk2basuWLad5z1MjCAMAAARIdHR06P9ChQpp//79ma5TuHBhLV++XFdffbXmz5+vHj16SJI6deqkLVu2aMmSJTp16pSaNWuW7jY++OAD3X777Vq1apXatm2rkydPpltbVFRUsjqjoqLSXP50EYQBAAACrHTp0ipbtqw+//xzSdKsWbNCvcNJDh06pISEBF1++eV69tln9e2334bmDRo0SNdff32GvcGJiYnavn27Lr74Yk2YMEEJCQk6dOhQ7tyhbOCqEQAAAPkhFy53llMzZszQsGHDdOTIEdWrVy803CHJwYMH1bdvXx09elTOOT3zzDOheQMGDNDo0aPVv3//dLd/6tQp/fGPf1RCQoKcc7rrrrtUpkyZ3Lo7WWbOuXzZcZs2bdyKFSvyZd84MzWf0TxH6629YW2EKwEAIPLi4+PVuHHj/C4j4t566y299957mjVrVn6XIintdjazlc65NimXpUcYAAAAOXLnnXfqww8/1IIFC/K7lBwhCAMAACBHJk2alGra7bffri+//DLZtOHDh2c4hji/EIQBAAAQMZMnT87vErKMq0YAAAAgkAjCAAAACCSCMAAAAAKJMcIAAAD5IKeXDU0PlxPNPnqEAQAAAuT5559X48aNNWDAgDzbZ8eOHfNsX9lBjzAib2zpnK1Xt1Zk6wAAAKm8+OKL+vjjj1WjRo3QtJMnT6pw4dyLhV999VWubft00CMMAAAQEMOGDdPmzZvVs2dPlS5dWgMHDlSnTp00cOBA7d27V1dddZXatm2rtm3bhq4FfPjwYd14441q166dWrVqpffeey/d7a9bt07t2rVTy5Yt1aJFC33//feSpBIlSkiSlixZoi5duqhv376qV6+eRo4cqVdffVXt2rVT8+bNtWnTptxvhDD0CAMAAATElClTtHDhQn3yySd64YUX9P777+uLL75Q0aJFdf311+vuu+/WhRdeqG3btql79+6Kj4/Xo48+qksuuURTp07V/v371a5dO3Xt2lXFixdPc/vDhw/XgAEDdPz4cZ06dSrVMt9++63i4+NVrlw51atXT0OHDtXy5cv13HPPadKkSZo4cWIetISHIAwAABBQffr0UdGiRSVJH3/8sdavXx+ad+DAAR06dEiLFi3SvHnz9NRTT0mSjh49qm3btqlx48aptnfBBRfo0Ucf1Y4dO9SvXz81aNAg1TJt27ZV1apVJUnnnnuuunXrJklq3ry5Pvnkk4jfx4wQhAEAAAIqvFc3MTFRy5YtU0xMTLJlnHOaO3euGjZsmOn2rr/+erVv314ffPCBLr/8cr300ku65JJLki0THR0d+j8qKip0OyoqSidPnjydu5NtBGEAAIB8cKZd7qxbt26aNGmS7r33XklSXFycWrZsqe7du2vSpEmaNGmSzEyrV69Wq1at0tzG5s2bVa9ePd11113atm2b1qxZkyoIn0k4WQ4AAAB6/vnntWLFCrVo0UJNmjTRlClTJEkPPfSQTpw4oRYtWqhp06Z66KGH0t3GG2+8oWbNmqlly5b67rvvNGjQoLwqP0fMOZcvO27Tpo1bsWJFvuwbuSyHl09rnsPLp51pn6gBAEhLfHx8muNqEVlptbOZrXTOtUm5LD3CAAAACCTGCAMAACBbPvroI91///3JptWtW1fvvPNOPlWUMwRhAACAPOKck5nldxmnrXv37urevXt+l5FKdof8MjQCAAAgD8TExGjfvn3ZDmvIGuec9u3bl+rybxmhRxgAACAP1KhRQzt27NDevXvzu5SzVkxMjGrUqJHl5QnCAAAAeaBIkSKqW7dufpeBMAyNAAAAQCARhAEAABBIWQrCZtbDzDaY2UYzG5nG/Fpm9omZrTazNWZ2eeRLBQAAACIn0yBsZoUkTZbUU1ITSf3NrEmKxUZLesM510rSdZJejHShAAAAQCRlpUe4naSNzrnNzrnjkmZL6ptiGSeplP9/aUk/RK5EAAAAIPKyEoSrS9oednuHPy3cWEl/NLMdkhZIujOtDZnZLWa2wsxWcOkQAAAA5KdInSzXX9J051wNSZdLmmVmqbbtnHvZOdfGOdemYsWKEdo1AAAAkH1ZCcI7JdUMu13DnxbuJklvSJJzbqmkGEkVIlEgAAAAkBuyEoS/kdTAzOqa2TnyToabl2KZbZIulSQzaywvCDP2AQAAAGesTIOwc+6kpDskfSQpXt7VIdaZ2Tgz6+Mvdo+km83sW0mvSxrs+CFtAAAAnMGy9BPLzrkF8k6CC582Juz/9ZI6RbY0AAAAIPfwy3IAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIpML5XQDOXHVGfpCj9bbERLgQAACAXECPMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCSCMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCSCMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCSCMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCSCMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCSCMAAAAAKJIAwAAIBAIggDAAAgkLIUhM2sh5ltMLONZjYynWX+YGbrzWydmb0W2TIBAACAyCqc2QJmVkjSZEmXSdoh6Rszm+ecWx+2TANJoyR1cs79YmaVcqtgAAAAIBKy0iPcTtJG59xm59xxSbMl9U2xzM2SJjvnfpEk59yPkS0TAAAAiKysBOHqkraH3d7hTwt3nqTzzOxLM1tmZj3S2pCZ3WJmK8xsxd69e3NWMQAAABABkTpZrrCkBpJiJfWX9A8zK5NyIefcy865Ns65NhUrVozQrgEAAIDsy0oQ3impZtjtGv60cDskzXPOnXDO/U/Sf+UFYwAAAOCMlJUg/I2kBmZW18zOkXSdpHkplnlXXm+wzKyCvKESmyNXJgAAABBZmQZh59xJSXdI+khSvKQ3nHPrzGycmfXxF/tI0j4zWy/pE0n3Ouf25VbRAAAAwOnK9PJpkuScWyBpQYppY8L+d5L+7P8BAAAAZzx+WQ4AAACBlKUeYSCSEp3TbfOP6q34k/r5V6dPbiim2DrpH4ou0emHmT/owDcHdOrwKdW5v45KNC6RhxUDAICzET3CyHMLvj+paXEn9H7/otp1Twl1rFkow+UPrjmo/Z/vV+0RtdVwYkMVa1AsjypFOOecxo4dq2rVqqlo0aKKjY3VunXrMl3vueeeU6NGjVS0aFHVqFFDt99+uw4dOhSaf/DgQY0YMUK1a9dW0aJF1bFjR33zzTcR2TcAABkhCCPPbfw5UVVLmjrWLKwqJaJ0TiHLcPnje46rcJnCKtagmIqUKaKowqkP2+PHj+dWufA98cQTevrppzVp0iR98803qlSpki677DIdPHgw3XVee+013XfffXrwwQcVHx+vmTNnasGCBRo+fHhomaFDh+qjjz7SjBkztHbtWnXr1k1du3bVzp2/XaUxJ/sGACAzBGHkqcHv/qq7PzqmbQlO9vAB1ZmYcZDZ8Y8d2v36bp3Yd0LfDf5OG+7ZIEna/Phm/TDjB+2avUsVK1ZUp06dJEnr169Xr169VLJkSVWqVEn9+/fX7t27k21z2rRpatKkiWJiYnTeeefp2WefVWJiYpbqf+aZZ9SiRQsVL15c1atX19ChQ7V///7Q/ISEBA0cOFCVKlVSTEyM6tWrp4kTJ4bmm5leeOEF9erVS8WKFVPt2rX1yiuvJNvHzp07dd1116ls2bIqW7asevXqpe+//z40f+zYsWrWrJlmz56tc889VyVLltSVV16pn376KUv3ISecc5o4caJGjhypq666Ss2aNdOMGTN08OBBvfbaa+mu99VXX6lDhw4aOHCg6tSpo0suuUSDBg3S119/LUn69ddfNXfuXI0fP16xsbGqX7++xo4dq/r16+vvf//7ae0bAIDMEISRp57rEaMxnc9RjVKmXfeU0Dc3F89w+aoDqqpi34oqXK6wGk5sqHP/em5o3v6l+yUnff7555o5c6Z27dqlzp07q1mzZlq+fLk+/vhjHTp0SH379g0F3X/84x964IEHNG7cOMXHx+vpp5/WhAkT9OKLL2ap/qioKE2cOFHr1q3Ta6+9puXLl+vOO+8MzR89erTWrl2r+fPna8OGDZo6daqqV0/+i+R//etf1adPH8XFxemWW27RoEGDtGLFCknSkSNHdPHFFysmJkaffvqpli5dqqpVq6pr1646cuRIaBtbtmzRnDlz9M4772jRokVavXq1HnzwwQxrb9q0qUqUKJHuX9OmTdNd93//+592796tbt26haYVLVpUnTt31ldffZXuehdeeKHi4uK0bNkySdK2bds0b948XX755ZKkkydP6tSpU4qJiUm2XtGiRfXFF1+c1r4BAMgMJ8shT5WOMZWMNhUyqUqJzD+HFSpWSIViCsnMVKRMkWTzilQooqr9q6pRo0aSpDFjxuj888/XhAkTQsvMnDlT5cqV04oVK9SuXTs98sgjeuKJJ3T11VdLkurWrauRI0fqxRdf1B133JFpPSNGjAj9X6dOHT3xxBPq27evZsyYoaioKG3dulWtW7dWu3btJEm1a9dOtY1+/frp1ltvlSQ9+OCD+uSTTzRx4kS98sormj17tpxzmjZtmsy8ISMvvfSSKlWqpPnz5+sPf/iDJC9ATp8+XaVLl5Yk3XLLLZo2bVqGtS9YsEAnTpxId36RIkXSnZfUq165cuVk0ytXrpxsCENK1113nfbt26fOnTvLOaeTJ09q4MCBoceoZMmSuuCCC/S3v/1NzZo1U5UqVfT6669r6dKlql+//mntGwCAzBCEUWAVrVM02e2VK1fqs88+U4kSqa8osWnTJtWtW1fbt2/Xrbfeqttuuy007+TJk/IuhZ25f//733r88ccVHx+vhIQEnTp1SsePH9fu3btVrVo13Xbbbbr66qu1cuVKXXbZZerdu7e6dOmSbBsXXHBBqtsffPBB6D7873//U8mSJZMtc+TIEW3atCl0u3bt2qEQLEnVqlXTjz/+mGHtaYXy3Pbpp5/qkUce0Ysvvqj27dtr48aNGj58uP76179q3LhxkqRZs2bpxhtvVI0aNVSoUCG1bt1a/fv318qVK/O8XgBAsBCEUWBFRSfvUU5MTFSvXr301FNPpVq2cuXKoaEFU6ZMUceOHbO9v61bt6pXr166+eabNW7cOJUvX16rVq1S//79Qyfr9ezZU1u3btWHH36oxYsXq1evXrrmmmsy7a0Nvw8tW7bU7NmzU80rV65c6P+Uvbdmluk456ZNm2rr1q3pzq9du3a6V2KoUqWKJGnPnj2qVatWaPqePXtC89IyevRo9e/fX0OHDpUkNW/eXIcPH9bQoUM1ZswYFS5cWOeee64+/fRTHT58WAcOHFDVqlV17bXXql69eqe1bwAAMkMQxlmjdevWeuONN1S7du00v+YvWbKkqlWrpk2bNmnQoEHZ3v6KFSt0/PhxPfvssypUyLvk2/z581MtV6FCBQ0cOFADBw5Uz5491b9/f02ZMkXR0dGSpGXLlunGG28MLb9s2TI1btw4dB9ef/11VahQQWXKlMl2jRk5naERdevWVZUqVfSvf/1Lbdu2lSQdPXpUn3/+uZ588sl01zty5EiorZIUKlQozR744sWLq3jx4vrll1/00Ucf6YknnjitfQMAkBmCMM4at99+u/7xj3/o2muv1f3336+KFStq8+bNeuONN/T000+rZMmSevjhh3XnnXeqTJkyuvzyy3XixAmtWrVKO3fu1KhRozLcfoMGDZSYmKiJEyeqX79+WrZsWbIrQkjeOOXWrVuradOmOnnypN5++23Vq1cvFIIl6e2331bbtm0VGxurt956S4sXLw5dRWHAgAF66qmn1LdvX40bN061atXS9u3b9d5772nYsGFq0KBBjtvndIZGmJlGjBihxx57TI0aNdJ5552nv/3tbypRooSuv/760HKXXnqp2rVrp8cff1yS1Lt3bz3zzDNq06ZNaGjEQw89pCuuuEKFC3svPx999JESExPVqFEjbdy4Uffee68aNWqkIUOGZGvfAABkF0EYZ41q1arpyy+/1KhRo9SjRw8dPXpUtWrVUrdu3UJBdOjQoSpevLiefPJJjRo1SkWLFlXTpk2zdKJcixYt9Nxzz2nChAkaPXq0OnbsqKeeekrXXnttaJno6Gg9+OCD+t///qeYmBh16NBB77//frLtjB07VnPnztVdd92lihUratq0aaGezmLFiumzzz7TyJEjdc011yghIUHVqlXTxRdfrLJly0awtbLvvvvu06+//qrbb79dv/zyi9q3b69FixYlG8+8adMm1axZM3R79OjRMjM99NBD2rFjhypUqKDevXvr0UcfDS2TkJCgUaNGaceOHSpXrpyuuuoqPfroo8l6qLOybwAAssuyepJQpLVp08YlXTIKZ6Y6Iz/I0XpbYnLWS9e8bq3MF0rD2hvW5mi9/GBmevPNN0NXrQAAALnPzFY659qknM51hAEAABBIDI1AvtqWkKgmkw9Jkn6NWp9qfv3H6uuc8ufkSS2vvvpq6Pq+KWV0RQUAAFAwEYSRr6qVNMUN8677e3mNqqnmp/wRjdzUp08ftW/fPs15GV1RITvyaygSAABIjSCMfFU4ylS/nPcLatGVozNZOneVLFmSk68AAAgQxggjzyU6p1vf/1Xlnzgoe/iAlmw5mev7XLJkicxMP/30U67vKy+tWLFCZqYtW7ZkeZ3BgwfriiuuyL2iAAAoIAjCyHMLvj+paXEn9H7/otp1Twl1rFko85VOU8eOHbVr1y6VL18+1/eFs8fcuXPVpEkTRUdHq0mTJnrnnXcyXeeNN95Qy5YtVaxYMdWuXTvNH/04fvy4xowZo7p16yo6Olq1atXS888/H5r/5ptvqk2bNipTpoyKFy+uli1basaMGRG9bwAAhkYgH2z8OVFVS5o61sy7w++cc87h53iRLUuXLtW1116rhx9+WP369dPbb7+ta665Rl9++WW6Y8k//PBDXX/99Xr++efVo0cPxcfH6+abb1bRokWTXav6uuuu044dO/Tyyy+rQYMG2rNnj3799dfQ/PLly2v06NFq1KiRihQpovnz5+umm25SxYoVdfnll+f6fQeAoKBHGHlq8Lu/6u6PjmlbgpM9fEB1Jh7MdJ3Nj2/WDzN/0O63div+jnjF3xmvXbN3ySX+duLZL7/8ohtuuEFly5ZV0aJF1bVr12RXeUg5NCIhIUEDBw5UpUqVFBMTo3r16iX7lbiEhATdcsstqlSpkkqWLKkuXbooq9e9nj59ukqUKKEPP/xQjRo1UrFixdSnTx8lJCTorbfeUoMGDVS6dGkNHDgwWfg5duyYRowYocqVK4d+jOOLL75Itu2FCxeqUaNGiomJ0UUXXaT//ve/qfb/1VdfqUuXLipWrJiqV6+u2267TQcOHMhS7SktXLhQF110kcqWLaty5cqpe/fuio+PT7bMuHHjVLt2bUVHR6tKlSrJfr46NjZWw4YN0/Dhw1W2bFmVLVtW9957rxITE0PLHD9+XPfff79q1KihYsWKqW3btvroo49C85Meu8WLF6t9+/YqVqyY2rRpo1WrVuXoPmXVxIkTdfHFF+vBBx9U48aN9eCDDyo2NjbVrwmGmzVrlnr37q0//elPqlevnnr16qVRo0ZpwoQJoRMlFy1apMWLF2vBggW67LLLVKdOHbVv316xsbGh7VxyySW68sor1ahRI5177rkaPny4WrRooc8//zxX7zMABA1BGHnquR4xGtP5HNUoZdp1Twl9c3PxLK23f+l+WZSp3uh6qjawmvYt2qeE5Qmh+YMHD9bXX3+t9957T8uXL1exYsXUo0ePZEEz3OjRo7V27VrNnz9fGzZs0NSpU1W9enVJ3pUdevXqpZ07d2r+/PlavXq1OnfurEsuuUS7du3KUr3Hjh3T008/rVdffVWLFy/WihUrdNVVV2nGjBmaO3eu3n33Xc2fP18vvvhiaJ377rtPc+bM0dSpU7V69Wo1b95cPXr0CO1z+/btuvLKK3XZZZcpLi5Od955p+67775k+127dq26deumPn366Ntvv9Xbb7+tuLg43XjjjVmqO6XDhw9rxIgRWr58uZYsWaLSpUurd+/eOn78uCRv6MBTTz2lF198Ud9//73mz5+vdu3aJdvGq6++qsTERC1dulQvvfSSXn755WRhcsiQIfr000/12muv6bvvvtMNN9yg3r1769tvv022nVGjRmn8+PFatWqVypcvrwEDBmR4FY5hw4apRIkSGf5t27Yt3fWXLl2qbt26JZvWvXt3ffXVV+muc+zYMcXExCSbVrRoUe3YsUNbt26VJL377rtq27atnnnmGdWoUUMNGjTQXXfdpUOHDqW5TeecFi9erA0bNqhz587p7hsAkH0MjUCeKh1jKhltKmRSlRJZ/xwWXS1alftV9v6vEq2fP/1Zh9cfVpkOZfT9999r3rx5+vTTT0NBYdasWapVq5ZeffVVDR06NNX2tm7dqtatW4dCW+3atUPzPvnkE8XFxWnv3r0qWrSoJOmRRx7R+++/r1mzZqUKn2k5efKkJk+erIYNG0qSrr/+ej377LPas2ePKlSoIEnq27evPvnkE91zzz06fPiw/v73v+uf//ynevXqJUmaMmWK/v3vf2vy5Mn629/+pr///e+hsaRmpkaNGum///2vHnroodB+n3zySV177bW65557QtP+/ve/q1WrVvrxxx9VqVKlLLe5JF111VXJbk+bNk2lSpXS8uXLdeGFF2rr1q2qWrWqunXrpiJFiqhWrVpq0yb5D/dUrVo1Vc3PPPOM/vznP2vTpk16/fXXtWXLFtWq5f2y4B133KGPP/5YL730UrIPCo888oguvvhiSdKYMWN04YUXaufOnapRo0aatY8bN05/+ctfMrx/1apVS3fe7t27Vbly5WTTKleurN27d6e7Tvfu3TVixAgtWrRIXbt21caNG/X0009Lknbt2qU6depo8+bN+uKLLxQdHa25c+dq//79uvPOO/XDDz/orbfeCm0rISFB1atX17Fjx1SoUCFNnjxZPXv2zPD+AACyhyCMAiGmZvJetiJliujkAe9qE/Hx8YqKitIFF1wQml+6dGk1b95c69en/pEOSbrtttt09dVXa+XKlbrsssvUu3dvdenSRZK0cuVKHTlyRBUrVky2ztGjR7Vp06Ys1RsdHR0KwZIXoKpUqRIKwUnTkurbtGmTTpw4oU6dOoXmFypUSBdccEFomfj4eHXo0EFmFlom/D4n1b5x40bNmTMnNC2p13TTpk3ZDsKbNm3SQw89pK+//lp79+5VYmKiEhMTQz2p11xzjZ577jnVrVtX3bt3V48ePdSnTx9FR/92Kby0an7ooYd04MABrVq1Ss45NWnSJNl+jx07pksuuSTZtBYtWoT+TwqwP/74Y7pBuFKlStm+v6fr5ptv1qZNm9S3b1+dOHFCpUqV0vDhwzV27FhFRXkf/BITE2Vmeu2111S6dGlJ0gsvvKDu3btrz549ofBdsmRJxcXF6dChQ1q8eLH+/Oc/q06dOrr00kvz9D4BwNmMIIwCwQpZ6olZ+G2K8AAWrmfPntq6das+/PBDLV68WL169dI111yjadOmKTExUZUrV05zPGapUqWyVG/hwsmfWmaW6kc5zCzZWNns3oe0JCYmaujQobr77rtTzUsa+pEdV1xxhWrUqKGXXnpJ1atXV+HChdWkSZPQ0IiaNWtqw4YNWrx4sT7++GPdc889evjhh/X111+rePHMh70khcJvvvkmVfsk9cYnCZ+f1CYZtd+wYcP0yiuvZLj/9evXh3qiU6pSpYr27NmTbNqePXsyPOnSzDRhwgQ99thj2r17typWrKjFixdLkurVqyfJ6yGvXr16KARLUuPGjSVJ27ZtCwXhqKgo1a9fX5LUsmVLxcfH67HHHiMIA0AEEYRR4DVu3Dg0BjVpaMSBAwe0du1aDRkyJN31KlSooIEDB2rgwIHq2bOn+vfvrylTpqh169bas2ePoqKiQuElt5177rk655xz9OWXX+rcc8+VJJ06dUpLly7V9ddfL8m7n3PnzpVzLhQEly1blmw7rVu31rp160IB6nTs27dP//nPf/Tiiy+GhiSsWrVKJ08mv+5zTEyMevXqpV69emnkyJGqUqWKvvzyy9D42q+//jpVzdWqVVOpUqXUqlUrOee0e/fu0D4i5XSHRlxwwQX617/+pXvvvTc07V//+pc6duyY6b4LFSoU+uDx+uuv64ILLgh9w9CpUye9+eabOnTokEqU8H5VMemkx/AhOiklJibq2LFjme4bAJB1BGEUeA0aNFDfvn1166236uWXX1aZMmX04IMPqlSpUqEQmdKYMWPUunVrNW3aVCdPntTbb7+tevXqKTo6Wl27dlWnTp3Ut29fPfHEE2rUqJF2796thQsXqmvXrrrooosifh+KFy+u2267Tffff78qVKigunXrhsYU/+lPf5Lk9XA+/fTTGjFihP70pz9p7dq1mjJlSrLt3H///erQoYOGDRumW2+9VSVLltR//vMfvf/++3rppZeyVVPZsmVVoUIF/eMf/1DNmjW1c+dO3Xvvvcl6u6dPn66TJ0+qffv2KlGihObMmaMiRYqoQYMGoWV++OGHZDU/+eSTGj16tCTpvPPO04ABAzR48GA9/fTTat26tX7++WctWbJE9erVU79+/XLapKc9NGL48OHq3Lmzxo8fryuvvFLvvPOOPvnkk2RX8hg1apSWL18e6vX96aef9Oabbyo2NlbHjh3TtGnT9Oabb+rTTz8NrXP99dfrkUce0ZAhQzR27Fjt379fw4cP19VXXx2q99FHH1X79u1Vr149HTt2TAsWLNCsWbM0adKkHN8fAEBqXDUCZ4Vp06apXbt26tOnj9q1a6cjR45o4cKFqb5eTxIdHa0HH3xQ559/vjp16qSDBw/q/fffl+R9vb1gwQJdcskluvnmm9WwYUP94Q9/0IYNGzLsQTxdEyZM0LXXXqshQ4aoZcuWWrNmjRYuXKiqVatKkmrVqqW3335bCxcu1Pnnn69nn31W48ePT7aNFi1a6LPPPtOWLVvUpUsXnX/++Ro1alSqk76yIioqSnPmzNGaNWvUrFkz3X777XrkkUeSjf8tU6aM/u///k8XXXSRmjVrprlz5+rtt99W3bp1Q8sMGDBAp06dUvv27XXzzTfrpptuSjZ0Y9q0aRoyZIjuu+8+NWrUSFdccYU+++yzDHtH80LHjh01e/ZsTZ8+XS1atNDMmTM1Z86cZNcQ3rVrV6px4zNnzlTbtm3VqVMnrVu3TkuWLEl2JY0SJUro448/VkJCgtq2bas//OEP6tKli6ZOnRpa5tChQ7rtttvUtGlTderUSXPnztXMmTM1bNiw3L/jABAgltHlh3JTmzZtXFavy4r8UWfkBzlab0tM2r2wmWleN+2xmplZe8PaHK2H3BcbG6tmzZrphRdeyO9SAAABZmYrnXNtUk6nRxgAAACBxBhh5KttCYlqMtn7IYFfo1Jf6qz+Y/V1Tvlz8rqsDPXs2TPdX/h64IEH9MADD+RxRTmzbdu2VJctC5fRFRUAADgbEISRr6qVNMUN886cv7xG1VTzi5QpkmpafvvnP/+Z7i/WlStXLo+ryblq1aopLi4uw/mna8mSJae9DQAAcgtBGPmqcJSpfjnvslrRlaMzWfrMkJPr8Z6JChcuHJHLrAEAUFAxRhgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAARSloKwmfUwsw1mttHMRmaw3FVm5sysTeRKBAAAACIv0yBsZoUkTZbUU1ITSf3NrEkay5WUNFzS15EuEgAAAIi0rPQIt5O00Tm32Tl3XNJsSX3TWO4RSRMkHY1gfQAAAECuyEoQri5pe9jtHf60EDNrLammc+6DjDZkZreY2QozW7F3795sFwsAAABEymmfLGdmUZKekXRPZss65152zrVxzrWpWLHi6e4aAAAAyLGsBOGdkmqG3a7hT0tSUlIzSUvMbIukDpLmccIcAAAAzmRZCcLfSGpgZnXN7BxJ10malzTTOZfgnKvgnKvjnKsjaZmkPs65FblSMQAAABABmQZh59xJSXdI+khSvKQ3nHPrzGycmfXJ7QIBAACA3FA4Kws55xZIWpBi2ph0lo09/bIAAACA3MUvywEAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAokgDAAAgEAiCAMAACCQshSEzayHmW0ws41mNjKN+X82s/VmtsbMFptZ7ciXCgAAAEROpkHYzApJmiypp6QmkvqbWZMUi62W1MY510LSW5KeiHShAAAAQCRlpUe4naSNzrnNzrnjkmZL6hu+gHPuE+fcEf/mMkk1IlsmAAAAEFlZCcLVJW0Pu73Dn5aemyR9eDpFAQAAALmtcCQ3ZmZ/lNRGUpd05t8i6RZJqlWrViR3DQAAAGRLVnqEd0qqGXa7hj8tGTPrKulBSX2cc8fS2pBz7mXnXBvnXJuKFSvmpF4AAAAgIrIShL+R1MDM6prZOZKukzQvfAEzayXpJXkh+MfIlwkAAABEVqZB2Dl3UtIdkj6SFC/pDefcOjMbZ2Z9/MWelFRC0ptmFmdm89LZHAAAAHBGyNIYYefcAkkLUkwbE/Z/1wjXBQAAAOQqflkOAAAAgUQQBgAAQCARhAEAABBIBGEAAAAEEkEYAAAAgUQQBgAAQCARhAEAABBIBGEAAAAEEkEYAAAAgUQQBgAAQCARhAEAABBIBGEAAAAEEkEYAAAAgUQQBgAAQCARhAEAABBIBGEAAAAEEkEYAAAAgUQQBgAAQCAVzu8CACC76oz8IEfrbRnfK8KVAAAKMnqEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBGEAQAAEEgEYQAAAAQSQRgAAACBRBAGAABAIBXO7wKAIKgz8oMcrbdlfK8IVwIAAJIQhAEEx9jSOVwvIbJ1nGH4oAYgqBgaAQAAgEAiCAMAACCQCMIAAAAIJIIwAAAAAomT5QAgE81nNM/RemtvWBvhSgAAkUSPMAAAAAKJIAwAAIBAIggDAAAgkAjCAAAACCROlgMA5Ay/1AeggKNHGAAAAIFEEAYAAEAgEYQBAAAQSARhAAAABBJBGAAAAIFEEAYAAEAgEYQBAAAQSARhAAAABBI/qAEAAID0ncU/nkMQBgAAZ6Q6Iz/I0XpbxveKcCVnhxy3Z0yECzmDZGlohJn1MLMNZrbRzEamMT/azOb48782szoRrxQAAACIoEx7hM2skKTJki6TtEPSN2Y2zzm3PmyxmyT94pyrb2bXSZog6drcKBgAgDNZTnrd6MEE8kdWhka0k7TRObdZksxstqS+ksKDcF9JY/3/35L0gpmZc85FsFYgeM7icVkAAOS3rATh6pK2h93eIal9ess4506aWYKk8pJ+ikSRWcVYIgAAAGSVZdZpa2ZXS+rhnBvq3x4oqb1z7o6wZb7zl9nh397kL/NTim3dIukW/2ZDSRsidUdyUQXlcaA/y9GekUNbRhbtGVm0Z+TQlpFFe0ZWQWnP2s65iiknZqVHeKekmmG3a/jT0lpmh5kVllRa0r6UG3LOvSzp5axWfCYwsxXOuTb5XcfZgvaMHNoysmjPyKI9I4e2jCzaM7IKentm5aoR30hqYGZ1zewcSddJmpdimXmSbvD/v1rSvxkfDAAAgDNZpj3C/pjfOyR9JKmQpKnOuXVmNk7SCufcPEn/J2mWmW2U9LO8sAwAAACcsbL0gxrOuQWSFqSYNibs/6OSrolsaWeMAjWUowCgPSOHtows2jOyaM/IoS0ji/aMrALdnpmeLAcAAACcjbL0y3IAAADA2YYgDAAAgEAiCAMAACCQsnSyXBCZWTlJcs79nN+1AMg7ZlbCOXcov+sAEDlmVlner+BK0k7n3J78rAdnDk6WC2NmtSQ9IelSSfslmaRSkv4taaRzbku+FXeWMbNGzrn/5HcdBYmZmaR2Cnsxl7Sca3ZHlpltc87Vyu86zhZ8sMg5wtvpM7OWkqbI+6GvpB8DqyHvPf5PzrlV+VPZ2cXM1jrnmud3HTlBj3BycyRNlDTAOXdKksyskLxLw82W1CH/SjvrLJJE2MgiM+sm6UVJ3yv5i3l9M/uTc25RvhVXAJnZn9ObJalEXtYSAOvFcz1b0gtvZrZfhLfsmi7pVufc1+ETzayDpGmSzs+PogoiM+uX3ixJVfKylkgiCCdXwTk3J3yCH4hnm9kj+VRTgWVmz6c3S1KZPCzlbPCcpK4pv5Uws7ryrvHdOD+KKsAek/SkpJNpzOPciWzig0XETRfhLVKKp2xHSXLOLTOz4vlRUAE2R9KrktL6FjImj2uJGIJwcivN7EVJMyRt96fVlPfz0avzraqCa4ikeyQdS2Ne/zyupaArLGlHGtN3SiqSx7WcDVZJetc5tzLlDDMbmg/1FHR8sIgswlvkfGhmH0iaqeTv64MkLcy3qgqmNZKecs59l3KGmXXNh3oigiCc3CBJN0l6WMnHYSb9jDSy5xtJ3znnvko5w8zG5n05BdpUSd+Y2WwlfzG/ThybOTFE0r505rXJy0LOEnywiCzCW4Q45+4ys56S+ir5+/pk/1dzkXUjJB1IZ97v87COiOJkOeQa/8obR51zR/K7lrOBmTVW6hfzec659flXFSCZWUNJPzvn9qYxrzIneWVfOuFtHuENiCyCcBgzKyyvR/hKJX/xeU/S/znnTuRTaQDyiJm97Jy7Jb/rAJC7eK5nT1hG+r2kav7kAp+RCMJhzOx1eZdUmaHfxmPWkDdGuJxz7tp8Ku2sY2YfOud65ncdZwPaMvuSrhOe1ixJ3zrnauRlPQWdmZWWNEpeJ0IleSfT/CjvDXK8c25/vhVXAIW1Z19JlUV75hjP9cg5WzMSY4ST+51z7rwU03ZIWmZm/82PggoyM2ud3ixJLfOwlAKPtoy4vZK2ymu/JM6/XSlfKirY3pB3vfVY59xuSTKzKvLeIN+Q1C0fayuIktrz4hTtOVi0Z3bxXI+cszIjEYST+9nMrpE01zmXKElmFiXvOsK/5GtlBdM3kj5V8hegJGXytpQCj7aMrM2SLnXObUs5w8y2p7E8MlbHOTchfIIf4CaY2Y35VFNBll57jjezIflUU0HFcz1yzsqMRBBO7jpJEyRN9i9cLnkh4xN/HrInXt61ML9POYMXoGyjLSNroqSyklK9Ocr7dUlkz1Yzu0/SjKQT4/xfRRus3656gKyjPSNnoniuR8pZmZEYI5xCOmfmv+eci8+/qgomM7ta0lrn3IY05l3pnHs376sqmGjLyDOzRkr7rHye69lkZmUljZTXnklfN++Rd+nJ8c65AttblB9oz8jiuR45Z2NG4kLnYczsfkmvyRs/9LX/J0mvm9nIfCusgHLOvZVWcPOVzdNiCjjaMrL83rbZ8oaaLPf/TDzXc8Q594tz7n7nXCPnXDn/r7Fz7n55J9AhG2jPyOG5Hjlna0aiRziMP9i7acpLgJjZOZLWOeca5E9lZx8z2+acq5XfdZwNaMvs47medzg+I4v2zB6e65FztrYlY4STS5R3bbytKaZX9echG8xsTXqz5F0SCFlEW0Ycz/UI4viMLNozoniuR85Z2ZYE4eRGSFpsZt/rtxMSakmqL+mO/CqqAKssqbtSn01qklL97DIyRFtG1gjxXI8kjs/Ioj0jZ4R4rkfKCJ2FbUkQDuOcW2hm50lqp+QDwb9xzp3Kv8oKrPmSSjjn4lLOMLMleV5NwUZbRhDP9Yjj+Iws2jNCeK5HztnalowRBgAAQCBx1QgAAAAEEkEYAAAAgUQQBgAAQCARhAEAABBIBGEAAAAE0v8DCU8ZI5UhbHIAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 864x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "df.plot.bar(\n",
    "    title=\"Comparing noise-free simulator, noisy simulator, and Aspen-M\",\n",
    "    figsize=(12, 6),\n",
    ")\n",
    "\n",
    "text = f\"f_free_aspen = {f_free_aspen:.3f} \\nf_noise_model_aspen = {f_noisy_aspen:.3f}\"\n",
    "plt.text(1, 0.5, text, fontsize=14)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aa8ffab8-e062-4c33-9c9d-4a312cd76b63",
   "metadata": {},
   "source": [
    "We confirm that the simulator with a noise model is closer to the distribution produced by Aspen-M-3."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c4571a5e-2f8d-4767-b28e-a527d2490830",
   "metadata": {},
   "source": [
    "### Smaller, reduced noise models\n",
    "\n",
    "The full Rigetti Aspen-M-3 noise model contains due to non-uniform qubit noise. We can obtain simpler, smaller noise models by coarse graining the model above.\n",
    "\n",
    "Here, we consider taking the average over all qubits for the $T_1$, $T_2$, depolarizing, and readout depolarizing rates. This is a substantially smaller noise model, but may be less accurate. \n",
    "\n",
    "We use the `from_filter` function to extract all instructions with amplitude dampening noise in the model. We then compute the mean of the error probabilities. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "eac8479f-4b29-42d8-a15a-ba10aae82d3f",
   "metadata": {},
   "outputs": [],
   "source": [
    "avg_t1 = np.mean(\n",
    "    [\n",
    "        n.noise.parameters\n",
    "        for n in noise_model.from_filter(noise=AmplitudeDamping).instructions\n",
    "    ]\n",
    ")\n",
    "avg_t2 = np.mean(\n",
    "    [\n",
    "        n.noise.parameters\n",
    "        for n in noise_model.from_filter(noise=PhaseDamping).instructions\n",
    "    ]\n",
    ")\n",
    "avg_depo = np.mean(\n",
    "    [\n",
    "        n.noise.parameters\n",
    "        for n in noise_model.from_filter(noise=Depolarizing).instructions\n",
    "    ]\n",
    ")\n",
    "avg_readout = np.mean(\n",
    "    [n.noise.parameters for n in noise_model.from_filter(noise=BitFlip).instructions]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d1a64cae-226c-45d2-9af3-4b5b7a36c7a4",
   "metadata": {},
   "source": [
    "Now we construct a new noise model with the mean values above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "e49eb24c-5f24-4084-9d4c-c2e7e82fc259",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Gate Noise:\n",
      "  AmplitudeDamping(0.0019483910859341453), GateCriteria(None, None)\n",
      "  PhaseDamping(0.0018677309884621795), GateCriteria(None, None)\n",
      "  Depolarizing(0.007153677053977023), GateCriteria(None, None)\n",
      "Readout Noise:\n",
      "  BitFlip(0.053481250000000015), ObservableCriteria(None, None)\n"
     ]
    }
   ],
   "source": [
    "simple_noise_model = NoiseModel()\n",
    "simple_noise_model.add_noise(AmplitudeDamping(avg_t1), GateCriteria())\n",
    "simple_noise_model.add_noise(PhaseDamping(avg_t2), GateCriteria())\n",
    "simple_noise_model.add_noise(Depolarizing(avg_depo), GateCriteria())\n",
    "simple_noise_model.add_noise(BitFlip(avg_readout), ObservableCriteria())\n",
    "\n",
    "print(simple_noise_model)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b44bbbd1-fb65-4b8f-8dfa-0ea5aa7979b2",
   "metadata": {},
   "source": [
    "We can see the resultant circuits contain qubit-independent noise:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "c83e96f6-54b0-445b-bcdd-8a7a16691b62",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "T  : |                     0                     |                 1                  |                     2                     |\n",
      "                                                                                                                                   \n",
      "q0 : -Rx(0.75)-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-C-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-I--------AD(0.0019)-PD(0.0019)-DEPO(0.0072)-\n",
      "                                                  |                                                                                \n",
      "q1 : -I--------AD(0.0019)-PD(0.0019)-DEPO(0.0072)-Z-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-Rx(0.50)-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-\n",
      "                                                                                                                                   \n",
      "q2 : -Rz(0.50)-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-I-AD(0.0019)-PD(0.0019)-DEPO(0.0072)-I--------AD(0.0019)-PD(0.0019)-DEPO(0.0072)-\n",
      "\n",
      "T  : |                     0                     |                 1                  |                     2                     |\n"
     ]
    }
   ],
   "source": [
    "simple_noisy_circ = simple_noise_model.apply(circ)\n",
    "print(simple_noisy_circ)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e2c20b79-cbf2-4897-80b3-452c78d56a03",
   "metadata": {},
   "source": [
    "We run the circuit on a noisy simulator below"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "a1bb19e8-6824-4a50-b9fb-ae06c14359b6",
   "metadata": {},
   "outputs": [],
   "source": [
    "simple_noisy_task = noisy_simulator.run(simple_noisy_circ, shots=100_000)\n",
    "simple_noisy_probs = simple_noisy_task.result().measurement_probabilities"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "76418e5c-be5b-4b3a-8916-c0f56f003cd7",
   "metadata": {},
   "source": [
    "and add it to the previous dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "b4d27a95-612d-4dea-9f6a-98d4525b0a70",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>aspen-m</th>\n",
       "      <th>noisy_sim</th>\n",
       "      <th>free_sim</th>\n",
       "      <th>simple_noisy_sim</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>000</th>\n",
       "      <td>0.53302</td>\n",
       "      <td>0.74544</td>\n",
       "      <td>0.81380</td>\n",
       "      <td>0.78149</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>100</th>\n",
       "      <td>0.11748</td>\n",
       "      <td>0.13151</td>\n",
       "      <td>0.12486</td>\n",
       "      <td>0.13176</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>010</th>\n",
       "      <td>0.24681</td>\n",
       "      <td>0.08706</td>\n",
       "      <td>0.05338</td>\n",
       "      <td>0.06130</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>110</th>\n",
       "      <td>0.05687</td>\n",
       "      <td>0.02742</td>\n",
       "      <td>0.00796</td>\n",
       "      <td>0.01061</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>011</th>\n",
       "      <td>0.01819</td>\n",
       "      <td>0.00055</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.00096</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>111</th>\n",
       "      <td>0.00477</td>\n",
       "      <td>0.00021</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.00015</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>001</th>\n",
       "      <td>0.01892</td>\n",
       "      <td>0.00650</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.01179</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>101</th>\n",
       "      <td>0.00394</td>\n",
       "      <td>0.00131</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.00194</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       aspen-m  noisy_sim  free_sim  simple_noisy_sim\n",
       "000    0.53302    0.74544   0.81380           0.78149\n",
       "100    0.11748    0.13151   0.12486           0.13176\n",
       "010    0.24681    0.08706   0.05338           0.06130\n",
       "110    0.05687    0.02742   0.00796           0.01061\n",
       "011    0.01819    0.00055       NaN           0.00096\n",
       "111    0.00477    0.00021       NaN           0.00015\n",
       "001    0.01892    0.00650       NaN           0.01179\n",
       "101    0.00394    0.00131       NaN           0.00194"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "simple_noisy_sim = pd.DataFrame.from_dict(simple_noisy_probs, orient=\"index\").rename(\n",
    "    columns={0: \"simple_noisy_sim\"}\n",
    ")\n",
    "df = df.join(simple_noisy_sim)\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b1699bc2-712f-4a99-b9d6-0a0920114a13",
   "metadata": {},
   "source": [
    "We compute the fidelity between the simple noise model and the QPU:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "d7710e0b-29ac-4968-aa00-e23cb8b32269",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Total fidelity between Aspen-M and full noise model is: 0.9582415634979304\n",
      "\n",
      "Total fidelity between Aspen-M and simple noise model is: 0.9401126418474742\n",
      "\n",
      "Total fidelity between Aspen-M and noise-free is: 0.915784824249657\n"
     ]
    }
   ],
   "source": [
    "f_simple = fidelity(df[\"simple_noisy_sim\"], df[\"aspen-m\"])\n",
    "\n",
    "print(f\"\\nTotal fidelity between Aspen-M and full noise model is: {f_noisy_aspen}\")\n",
    "print(f\"\\nTotal fidelity between Aspen-M and simple noise model is: {f_simple}\")\n",
    "print(f\"\\nTotal fidelity between Aspen-M and noise-free is: {f_free_aspen}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "e36e19d0-664a-4cb4-a2e3-0061cd60a6c1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsIAAAF+CAYAAACI8nxKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAABPN0lEQVR4nO3deZxOdf/H8ddnDMbW2PdtZBuMdchSUcmSEG2WMJWKUNwlKql0Ky1KpNTdbUnu4qZCFG0iKUbGOnVbmixRKEtkne/vj+ty/WbGzLhwzYyZ6/18PObxmOuc7znnc84crvf1vb7nHHPOISIiIiISbEKyugARERERkaygICwiIiIiQUlBWERERESCkoKwiIiIiAQlBWERERERCUoKwiIiIiISlBSERSRbMLOeZrY4k7f5mJm9nUHrXmJmfTNi3ZcKM2tlZjuzug4RkbQoCIsEGTPrYWaxZvaXme02s0/M7MqsrutcnHMznHNtMnmbzzrnsjysmlmCmbXO6joygnlsM7NNWV1Laswsxsycmb2SYnpn7/Sp6Sx7m5l9a2ZHzWxJRtcqIudPQVgkiJjZP4BxwLNAKaAi8DrQOQvLOiczC83qGrKrbHDsrgZKAlXMrHFWF5OGrcBtKY5lH+B/51juDzz/3sZkUF0icpEUhEWChJmFA6OAAc65D5xzR5xzJ51z851zQ71t8prZODP71fszzszyeue1MrOdZvaImf3u7U2+ycxuMLP/mdkfZvZYku09ZWazzWymmR02sx/MrF6S+cPNbKt33iYz65JkXoyZLTezV8xsP/CUd9o3Sdo4M+tnZpvN7ICZTTQz887LZWZjzWyfmf1sZgO97VMNhWY2zMx2eWv5ycyuS7IP73p/r+xdx51mtsPM/vRuv7GZrfPW8FqK/X83yevKadVgZpeb2Zdmtt9b8wwzK+ydNx3PB5b53l78R7zTO5nZRu92l5hZZJL1JXj3aR1wxMxC09rHVGrpYGZrzOyQdz+fSmUf+pjZdm+tjyeZn8/MpnqPzSbAn2DbB5gLLPT+nrSWGG9v8WHv37FnkunLzew1MztoZj8m3R8zCzezf3vP0V1m9k8zy5Vk2W/M7CVvnT+bWftz1LgHWA+09a6jKNAcmJfeQs65z51zs4Bf/TgOIpIFFIRFgkczIAz4MJ02jwNNgfpAPaAJMCLJ/NLedZQDRgL/Au4AGgFXAU+YWUSS9p2B/wJFgf8AH5lZbu+8rd5lwoGngXfNrEySZa8AtuHpuR6dRr034glbdYHb8AYV4B6gvXc/GgI3pbXDZlYDGAg0ds4V8q4jIa323rqqAbfj6e17HGgN1MbTa9gynWXTLAN4DigLRAIVgKcAnHO9gO1AR+dcQefcC2ZWHXgPGAyUwBMi55tZniTr7A50AAoDl5/HPh4BenuX6wD0N7ObUrS5EqgBXAeMTBLCn/Ru63LvNvqQDjPLD9wCzPD+dDuzD2ZWABgPtPfW3ByIS7L4FXjOoeLe7X7gDagAU4FTQFWgAdAG6Jti2Z+8y74A/PvMh6h0vIPnuAB0wxPej59jGRG5xCkIiwSPYsA+59ypdNr0BEY55353zu3FE1B7JZl/EhjtnDsJvI8nSLzqnDvsnNsIbMIToM9Y7Zyb7W3/Mp4Q3RTAOfdf59yvzrlE59xMYDOe4H3Gr865Cc65U865v9Ood4xz7oBzbjvwFZ7gC55Q/Kpzbqdz7k/S/2r6NJAXqGVmuZ1zCc65rem0f8Y5d8w5txhPaHzPe7x2AcvwBK/z4pzb4pz7zDl33HvcXwbSC9S3Awu8y5wEXgLy4QmLZ4x3zu3wHju/99E5t8Q5t977d1mHJ3CnrOVp59zfzrm1wFr+/29+G57z4w/n3A48QTY9XfGEycXAAiA3nvB9RiJQx8zyOed2e8+xM34Hxnm/1ZiJJ9h2MLNSwA3AYO+3Hr8Dr+AJr2f84pz7l3PuNDANKIPnA1d6PgRameebld54grGIZHMKwiLBYz9QPK3hAV5lgV+SvP7FO823Dm94ADgTTn9LMv9voGCS1zvO/OKcSwR2nlmfmfU2szjvV/sHgDp4gvVZy6ZjT5LfjybZdtkUy6e5LufcFjw9q08Bv5vZ+2ZWNq32nL2/6e2/X8yslHe7u8zsEPAuyY9FSsn+Tt5juwNPT/0ZSY+93/toZleY2VdmttfMDgL9UqnF3+Oe9FxKTR9glvfDzjFgjncazrkjeAJ/P2C3mS0ws5pJlt3lnHMptlUWqIQnUO9Ocm69iWcc8ln1O+eOen8taGZXeYef/GVmSUM33g8UC/B8Q1LMObc86Xwzm5Rk2ccQkWxBQVgkeKzA0/t2UzptfsUTJM6oyMWNb6xw5hczCwHKA7+aWSU8wyoG4gkVhYENeIYInJE05Jyv3d5tnVVHapxz/3HOXYln3x3w/EVs+4wjQP4kr0un0/ZZ73ajnHOX4Rlukt6xSPZ38n6tXwHYldYy57GP/8Ez9rWCcy4cmJSilvTsJvmxrphWQzMrD1wL3GFme8xsD55hEjeYWXFvzYucc9fj6bH9Ec85c0a5FMMZzpyrO/Cc58Wdc4W9P5c552qfq3jn3DLv8JOCabR/B3gIzweVlMv2S7Lss+falohcGhSERYKEc+4gnnG9E81zkVt+M8ttZu3N7AVvs/eAEWZWwhtGRpLKm/55aGRmXb290IPxBJTvgAJ4wtheADO7E0+PcKDMAh40s3Lei86GpdXQzGqY2bXmuSjwGJ5e3cQA1BAHXG1mFb1fpz+aTttCwF/AQTMrBwxNMf83oEqS17PwDAO4zjvm+iE8x/bb1FZ+nvtYCPjDOXfMzJoAPdLbyRRmAY+aWRFv0B2UTtteeO66UAPPkJb6QHU83xp09/aSd/aOFT6O5/gkrbkk8ID3HL4Vz9jqhc653XiGWow1s8vMLMQ8FyNeyNjtlL4Grgcm+NPYPBdthgGhQIiZhSUZIy8ilwAFYZEg4pwbC/wDz9e7e/H0ng0EPvI2+ScQC6zDc5X8D95pF2ounq+3/8QTfLp6x3RuAsbi6aX+DYgClqe5lvP3LzxhaB2wBs/FZKfwjJVNKS+eMcT78HxlXpL0Q6tfnHOfATO9NawGPk6n+dN4Luo7iOfr9w9SzH8OzweUA2b2sHPuJzy9xhO8dXfEczHdiTTWfz77eD8wyswO4/kgNCu9/UxlP34BfsZz/Ken07YP8Lpzbk/SHzw90H3wvD/9A08v7x94xin3T7L893guWtyH52LKW5xz+73zegN58IxZ/xOYjadX+aI4jy+cc3/4uUgvPB863sBzYejfJO/VFpEsZsmHWImIBIZ5brtV1Tl3xyVQS3tgknOu0jkbyyXPzGKAvt6hHiIiF0w9wiKS45jnfrY3mOf+ueXw3F4rvdvGiYhIEFIQFpGcyPB8Tf8nnqER8Xi+5hcREfHR0AgRERERCUrqERYRERGRoKQgLCIiIiJBKb0nTGWo4sWLu8qVK2fV5kVEREQkSKxevXqfc65EyulZFoQrV65MbGxsVm1eRERERIKEmaX6yHcNjRARERGRoKQgLCIiIiJBSUFYRERERIJSlo0RFhEREcksJ0+eZOfOnRw7diyrS5EMFBYWRvny5cmdO7df7RWERUREJMfbuXMnhQoVonLlyphZVpcjGcA5x/79+9m5cycRERF+LaOhESIiIpLjHTt2jGLFiikE52BmRrFixc6r119BWERERIKCQnDOd75/YwVhEREREUnTXXfdRcmSJalTp05WlxJwGiMsIiIiQafy8AUBXV/CmA4BXd+lJCYmhoEDB9K7d++sLiXg1CMsIiIikgluuukmGjVqRO3atXnrrbc4ffo0MTEx1KlTh6ioKF555RUAWrVqxYMPPkj9+vWpU6cOK1euBODIkSPcddddNGnShAYNGjB37lwApk6dSteuXWnXrh3VqlXjkUceSXX7U6dO5aabbuL666+ncuXKvPbaa7z88ss0aNCApk2b8scff6S63NVXX03RokUz4IhkPfUIi4iIiGSCyZMnU7RoUf7++28aN25Mo0aN2LVrFxs2bADgwIEDvrZHjx4lLi6OpUuXctddd7FhwwZGjx7Ntddey+TJkzlw4ABNmjShdevWAMTFxbFmzRry5s1LjRo1GDRoEBUqVDirhg0bNrBmzRqOHTtG1apVef7551mzZg1DhgzhnXfeYfDgwZlxKC4Z6hEWERERyQTjx4+nXr16NG3alB07dnDixAm2bdvGoEGD+PTTT7nssst8bbt37w54emMPHTrEgQMHWLx4MWPGjKF+/fq0atWKY8eOsX37dgCuu+46wsPDCQsLo1atWvzyyy+p1nDNNddQqFAhSpQoQXh4OB07dgQgKiqKhISEjD0AlyAFYREREZEMtmTJEj7//HNWrFjB2rVradCgAcePH2ft2rW0atWKSZMm0bdvX1/7lHc/MDOcc8yZM4e4uDji4uLYvn07kZGRAOTNm9fXNleuXJw6dYoPP/yQ+vXrU79+fWJjY89qFxIS4nsdEhLCqVOn2LFjh2+ZSZMmZdjxuFRoaMQ5RE2L8qvd+j7rM7gSERERya4OHjxIkSJFyJ8/Pz/++CPfffcd+/btIzExkZtvvpkaNWpwxx13+NrPnDmTa665hm+++Ybw8HDCw8Np27YtEyZMYMKECZgZa9asoUGDBmlus0uXLnTp0sX3+swQjPRUqFCBuLi4i9rX7ERBWERERCSDtWvXjkmTJhEZGUmNGjVo2rQpu3btolWrViQmJgLw3HPP+dqHhYXRoEEDTp48yeTJkwF44oknGDx4MHXr1iUxMZGIiAg+/vjjDK+9e/fuLFmyhH379lG+fHmefvpp7r777gzfbmYw59y5G5m1A14FcgFvO+fGpJhfEZgGFPa2Ge6cW5jeOqOjo92ZbvpLmXqERUREsr/4+HjfMIJLXatWrXjppZeIjo7O6lKypdT+1ma22jl31gE95xhhM8sFTATaA7WA7mZWK0WzEcAs51wDoBvw+gXWLiIiIiKSKfwZGtEE2OKc2wZgZu8DnYFNSdo44MyljuHAr4EsUkRERCRYLFmyJKtLCBr+3DWiHLAjyeud3mlJPQXcYWY7gYXAoNRWZGb3mlmsmcXu3bv3AsoVEREREQmMQN0+rTsw1TlXHrgBmG5mZ63bOfeWcy7aORddokSJAG1aREREROT8+ROEdwFJH01S3jstqbuBWQDOuRVAGFA8EAWKiIiIiGQEf4LwKqCamUWYWR48F8PNS9FmO3AdgJlF4gnCGvsgIiIiIpescwZh59wpYCCwCIjHc3eIjWY2ysw6eZs9BNxjZmuB94AY58992UREREREsohfD9Tw3hN4YYppI5P8vgloEdjSRERERDLIU+EBXt/BgK5u5MiRXH311bRu3Tqg6z1j3rx5bNq0ieHDh2fI+rMLPVlORERE5BIzatSoDF1/p06d6NSp07kb5nCBumtE0IuvGen3j4iIiASXhIQEIiMjueeee6hduzZt2rTh77//Ji4ujqZNm1K3bl26dOnCn3/+CUBMTAyzZ88GYPjw4dSqVYu6devy8MMPc/jwYSIiIjh58iQAhw4dSvY6pfHjx/uW79atGwBTp05l4MCBvm3179+fpk2bUqVKFZYsWcJdd91FZGQkMTExGXxkspaCsIiIiEgm2Lx5MwMGDGDjxo0ULlyYOXPm0Lt3b55//nnWrVtHVFQUTz/9dLJl9u/fz4cffsjGjRtZt24dI0aMoFChQrRq1YoFCxYA8P7779O1a1dy586d6nbHjBnDmjVrWLduHZMmTUq1zZ9//smKFSt45ZVX6NSpE0OGDGHjxo2sX7+euLi4gB6HS4mCsIiIiEgmiIiIoH79+gA0atSIrVu3cuDAAVq2bAlAnz59WLp0abJlwsPDCQsL4+677+aDDz4gf/78APTt25cpU6YAMGXKFO688840t1u3bl169uzJu+++S2ho6qNiO3bsiJkRFRVFqVKliIqKIiQkhNq1a5OQkHCRe37pUhAWERERyQR58+b1/Z4rVy4OHDhwzmVCQ0NZuXIlt9xyCx9//DHt2rUDoEWLFiQkJLBkyRJOnz5NnTp10lzHggULGDBgAD/88AONGzfm1KlTadYWEhKSrM6QkJBU2+cUCsIiIiIiWSA8PJwiRYqwbNkyAKZPn+7rHT7jr7/+4uDBg9xwww288sorrF271jevd+/e9OjRI93e4MTERHbs2ME111zD888/z8GDB/nrr78yZoeyId01QkRERIJPgG93dqGmTZtGv379OHr0KFWqVPENdzjj8OHDdO7cmWPHjuGc4+WXX/bN69mzJyNGjKB79+5prv/06dPccccdHDx4EOccDzzwAIULF86o3cl2LKueexEdHe1iY2OzZNvnI2palF/tZj3n/9cGkT/GX2g5IiIicgHi4+OJjMxZd26aPXs2c+fOZfr06VldyiUltb+1ma12zkWnbKseYREREZFsZtCgQXzyyScsXLjw3I0lTQrCIiIiItnMhAkTzpo2YMAAli9fnmzagw8+mO4Y4mCnICwiIiKSA0ycODGrS8h2dNcIEREREQlKCsIiIiIiEpQUhEVEREQkKGmMsIiIiAQdf2+P6q/1fdYHdH2SOYIzCD8V7n/biIoZV4eIiIgEjfHjx/PGG2/QsGFDZsyYkSnbbN68Od9++22mbCs7Cs4gLCIiIpLJXn/9dT7//HPKly/vm3bq1ClCQzMujikEp09jhEVEREQyWL9+/di2bRvt27cnPDycXr160aJFC3r16sXevXu5+eabady4MY0bN/bdC/jIkSPcddddNGnShAYNGjB37tw0179x40aaNGlC/fr1qVu3Lps3bwagYMGCACxZsoSWLVvSuXNnqlSpwvDhw5kxYwZNmjQhKiqKrVu3ZvxBuASpR1hEREQkg02aNIlPP/2Ur776itdee4358+fzzTffkC9fPnr06MGQIUO48sor2b59O23btiU+Pp7Ro0dz7bXXMnnyZA4cOECTJk1o3bo1BQoUSHX9Dz74ID179uTEiROcPn36rDZr164lPj6eokWLUqVKFfr27cvKlSt59dVXmTBhAuPGjcuEI3FpURAWERERyWSdOnUiX758AHz++eds2rTJN+/QoUP89ddfLF68mHnz5vHSSy8BcOzYMbZv305kZORZ62vWrBmjR49m586ddO3alWrVqp3VpnHjxpQpUwaAyy+/nDZt2gAQFRXFV199FfB9zA4UhEVEREQyWdJe3cTERL777jvCwsKStXHOMWfOHGrUqHHO9fXo0YMrrriCBQsWcMMNN/Dmm29y7bXXJmuTN29e3+8hISG+1yEhIZw6depidifbUhAWERGRoHMp3e6sTZs2TJgwgaFDhwIQFxdH/fr1adu2LRMmTGDChAmYGWvWrKFBgwaprmPbtm1UqVKFBx54gO3bt7Nu3bqzgrCcTRfLiYiIiGSh8ePHExsbS926dalVqxaTJk0C4IknnuDkyZPUrVuX2rVr88QTT6S5jlmzZlGnTh3q16/Phg0b6N27d2aVn62Zcy5LNhwdHe1iY2OzZNvncx/hKD/vIzzrOf+/Uoj8Md7vtiIiInLx4uPjUx1bKzlPan9rM1vtnItO2VY9wiIiIiISlDRGWERERCSbWLRoEcOGDUs2LSIigg8//DCLKsreFIRFREREsom2bdvStm3brC4jx9DQCBEREREJSgrCIiIiIhKU/ArCZtbOzH4ysy1mNjyV+a+YWZz3539mdiDglYqIiIiIBNA5xwibWS5gInA9sBNYZWbznHO+ZwE654YkaT8ISP1uzyIiIiKXgPiagb2Vmm6Nmj350yPcBNjinNvmnDsBvA90Tqd9d+C9QBQnIiIiklP17duXTZs2nbuhHwoWLBiQ9ZzLyJEj+fzzzzNs/fPmzWPMmDEZtv6U/LlrRDlgR5LXO4ErUmtoZpWACODLNObfC9wLULGifw+qEBEREcmJ3n777awu4byNGjUqQ9ffqVMnOnXqlKHbSCrQF8t1A2Y7506nNtM595ZzLto5F12iRIkAb1pERETk0nTkyBE6dOhAvXr1qFOnDjNnzqRVq1acecpuwYIFGTp0KLVr16Z169asXLmSVq1aUaVKFebNmwfA1KlT6dy5M61ataJatWo8/fTTqW7rxRdfpHHjxtStW5cnn3wyzZoSEhKIjIzknnvuoXbt2rRp04a///4bgLi4OJo2bUrdunXp0qULf/75JwAxMTHMnj0bgOHDh1OrVi3q1q3Lww8/zOHDh4mIiODkyZMAHDp0KNnrlMaPH+9bvlu3br59HDhwoG9b/fv3p2nTplSpUoUlS5Zw1113ERkZSUxMzPkc/jT5E4R3ARWSvC7vnZaabmhYhIiIiEgyn376KWXLlmXt2rVs2LCBdu3aJZt/5MgRrr32WjZu3EihQoUYMWIEn332GR9++CEjR470tVu5ciVz5sxh3bp1/Pe///UF6TMWL17M5s2bWblyJXFxcaxevZqlS5emWdfmzZsZMGAAGzdupHDhwsyZMweA3r178/zzz7Nu3TqioqLOCt379+/nww8/ZOPGjaxbt44RI0ZQqFAhWrVqxYIFCwB4//336dq1K7lz505122PGjGHNmjWsW7eOSZMmpdrmzz//ZMWKFbzyyit06tSJIUOGsHHjRtavX09cXFya++Uvf4LwKqCamUWYWR48YXdeykZmVhMoAqy46KpEREREcpCoqCg+++wzhg0bxrJlywgPD082P0+ePL5wHBUVRcuWLcmdOzdRUVEkJCT42l1//fUUK1aMfPny0bVrV7755ptk61m8eDGLFy+mQYMGNGzYkB9//JHNmzenWVdERAT169cHoFGjRiQkJHDw4EEOHDhAy5YtAejTp89ZYTo8PJywsDDuvvtuPvjgA/Lnzw94xj1PmTIFgClTpnDnnXemue26devSs2dP3n33XUJDUx+t27FjR8yMqKgoSpUqRVRUFCEhIdSuXTvZcblQ5wzCzrlTwEBgERAPzHLObTSzUWaWdBBHN+B955y76KpEREREcpDq1avzww8/EBUVxYgRI84aa5s7d27MDICQkBDy5s3r+/3UqVO+dmfapPXaOcejjz5KXFwccXFxbNmyhbvvvjvNus5sByBXrlzJtpWe0NBQVq5cyS233MLHH3/sC/EtWrQgISGBJUuWcPr0aerUqZPmOhYsWMCAAQP44YcfaNy4carbTnocktaa8rhcKL8eseycWwgsTDFtZIrXT110NSIiIiKZILNvd/brr79StGhR7rjjDgoXLnzBF8p99tln/PHHH+TLl4+PPvqIyZMnJ5vftm1bnnjiCXr27EnBggXZtWsXuXPnpmTJkn5vIzw8nCJFirBs2TKuuuoqpk+f7usdPuOvv/7i6NGj3HDDDbRo0YIqVar45vXu3ZsePXrwxBNPpLmNxMREduzYwTXXXMOVV17J+++/z19//eV3jYHiVxAWERERkQu3fv16hg4dSkhICLlz5+aNN97g4YcfPu/1NGnShJtvvpmdO3dyxx13EB0dnWx+mzZtiI+Pp1mzZoDnIrx33333vIIwwLRp0+jXrx9Hjx6lSpUqvuEOZxw+fJjOnTtz7NgxnHO8/PLLvnk9e/ZkxIgRdO/ePc31nz59mjvuuIODBw/inOOBBx6gcOHC51VjIFhWjWSIjo52KQd4Z5qnws/dxisqwr/bvM16zv/ued10W0REJHPFx8cTGRnYh2hktqlTpxIbG8trr72W1aWka/bs2cydO5fp06dnyfZT+1ub2WrnXHTKtuoRFhEREZGAGDRoEJ988gkLFy48d+NLgIKwiIiISDYQExNzQffP3b9/P9ddd91Z07/44guKFSsWgMr+34QJE86aNmDAAJYvX55s2oMPPpjuHSUyi4KwiIiIBAXn3Fl3WQgGxYoVC8g9dy/UxIkTM21b5zvkN9BPlhMRERG55ISFhbF///7zDkqSfTjn2L9/P2FhYX4vox5hERERyfHKly/Pzp072bt3b1aXIhkoLCyM8uXL+91eQVhERERyvNy5cxMREZHVZcglRkMjRERERCQoKQiLiIiISFBSEBYRERGRoKQgLCIiIiJBSUFYRERERIKSgrCIiIiIBCUFYREREREJSgrCIiIiIhKUFIRFREREJCgpCIuIiIhIUFIQFhEREZGgpCAsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQUlBWDJdYmIi9913H8WKFcPMWLJkSUDbi4iIiPjDryBsZu3M7Ccz22Jmw9Noc5uZbTKzjWb2n8CWKTnJwoULmTJlCvPnz2f37t00b948oO0vRdu3b6djx44UKFCA4sWL88ADD3DixIl0l9m6dStdunShRIkSXHbZZdx222389ttvydqMHj2aFi1aUKBAAcwszXW9++671K9fn7CwMIoXL07v3r0Dsl8iIiLZ2TmDsJnlAiYC7YFaQHczq5WiTTXgUaCFc642MDjwpUpOsWXLFsqUKUPz5s0pXbo0efLkuej25wqVWen06dN06NCBw4cPs2zZMt577z1mz57NQw89lOYyR44coU2bNjjn+PLLL1m+fDknTpygY8eOJCYm+todP36crl27Mnjw4DTXNX78eIYOHcrDDz/Mhg0b+Oqrr+jcuXMgd1FERCRbCvWjTRNgi3NuG4CZvQ90BjYlaXMPMNE59yeAc+73QBcqOUNMTAzTpk0DwMyoVKkSCQkJ592+VatWREZGUqBAAaZNm0blypVZtWoVmzZtYujQoSxdupR8+fJx3XXX8corr1C6dGnfOqdMmcKLL77Itm3bqFixIv379+fBBx8kJCRjRgotXryYjRs38ssvv1ChQgUAXnjhBfr27cvo0aO57LLLzlpm+fLl/Pzzz8TGxlKkSBEApk2bRpEiRfjyyy9p3bo1AKNGjQJg9uzZqW77wIEDPProo3z00Udcf/31vulRUVEB3UcREZHsyJ93/nLAjiSvd3qnJVUdqG5my83sOzNrF6gCJWd59dVXGTlyJOXLl2f37t2sWrXqgtu/++67OOdYtmwZ77zzDrt37+bqq6+mTp06rFy5ks8//5y//vqLzp07+3pR//Wvf/HYY48xatQo4uPjGTt2LM8//zyvv/56mjUsW7aMggULpvvz7LPPprn8ihUriIyM9IVggLZt23L8+HFWr16d6jLHjx/HzAgLC/NNCwsLIyQkhG+++SbdY5bU4sWLOX36NL/99hu1atWiXLlydOnShW3btvm9DhERkZzKnx5hf9dTDWgFlAeWmlmUc+5A0kZmdi9wL0DFihUDtGnJTsLDwylUqBC5cuVK1kt7Ie0jIiIYO3as7/XIkSOpV68ezz//vG/aO++8Q9GiRYmNjaVJkyY888wzvPDCC9xyyy2+dQwfPpzXX3+dgQMHplpDdHQ0cXFx6dZZtGjRNOft2bOHUqVKJZtWvHhxcuXKxZ49e1JdpmnTphQsWJChQ4f69mf48OGcPn2a3bt3p1tLUtu2bSMxMZF//vOfjBs3jqJFizJq1CiuueYa4uPjyZ8/v9/rEhERyWn8CcK7gApJXpf3TktqJ/C9c+4k8LOZ/Q9PME7W3eecewt4CyA6OtpdaNEiAI0aNUr2evXq1SxdupSCBQue1Xbr1q1ERESwY8cO7rvvPvr37++bd+rUKZxL+3TMly8fVatWDVzhfihRogT//e9/6d+/P6+//johISF0796dhg0bntcQjsTERE6ePMn48eNp06YNADNmzKB06dLMnz+f22+/PaN2QURE5JLnTxBeBVQzswg8Abgb0CNFm4+A7sAUMyuOZ6iEvnuVDFWgQIFkrxMTE+nQoQMvvfTSWW1LlSrF0aNHAZg0adJ53Xli2bJltG/fPt02jz32GI899liq80qXLs3y5cuTTdu3bx+nT59Ot1e8TZs2bN26lX379hEaGkrhwoUpXbo0VapU8bv2MmXKAFCr1v9f3xoeHk7ZsmXZvn273+sRERHJic4ZhJ1zp8xsILAIyAVMds5tNLNRQKxzbp53Xhsz2wScBoY65/ZnZOEiKTVs2JBZs2ZRqVIlcufOfdb8QoUKUbZsWbZu3Xpetw+72KERzZo145///Cc7d+6kfPnyAHz22WfkzZv3rF7t1BQvXhyAL7/8kt9//51OnTr5XXuLFi0A+Omnn3zb/uuvv9i9ezeVKlXyez0iIiI5kV9jhJ1zC4GFKaaNTPK7A/7h/RHJEgMGDOBf//oXt99+O8OGDaNEiRJs27aNWbNmMXbsWAoVKsTTTz/NoEGDKFy4MDfccAMnT57khx9+YNeuXTz66KOprvdih0a0adOG2rVr07t3b8aOHcv+/fsZOnQo99xzj++OEStXrqR379688847NGnSBPDc3aJmzZqULFmSFStW8OCDDzJkyBBq1KjhW/f27dv5448/fHfeOBPYq1atSsGCBalevTqdO3fmwQcf5M0336RIkSI8+eSTlCxZkhtvvPGC90lERCQn0JPlJMcoW7Ysy5cvJyQkhHbt2lG7dm0GDBhA3rx5yZs3LwB9+/Zl8uTJTJ8+nXr16nHVVVfx1ltvERERkWF15cqViwULFpA/f35atGjB7bffzs0335xsCMfRo0f56aeffMM3wNOL26VLFyIjIxk1ahSPP/74WcM+Ro4cSYMGDRg6dCgADRo0oEGDBsTGxvraTJ8+nWbNmtGxY0datGjBsWPH+OKLL3ShnIiIBD1L7yKhjBQdHe2SvllnqqfC/W4aFeHf3S1mPXfK73VG/hjvd1sRERERuThmtto5F51yunqERURERCQoKQhLltq+fXu6D6rQnQ1EREQkowTqgRoiF6Rs2bLp3pGhbNmymVeMiIiIBBUFYclSoaGhmf6wChERERHQ0AjJAomJidx3330UK1YMM2PJkiUXvK4lS5ZgZuzbty9wBaahVatWaT6GOatcijUFwo033khMTIzf7RMSEjAzsuwCXBERyZYUhCXTLVy4kClTpjB//nx27959Xk95S6l58+bs3r2bYsWKBbDC7OODDz7gueeey+oygsbXX39No0aNCAsLo0qVKkyaNOmcy3zxxRc0b96cQoUKUbp0aYYNG8apU6nfZWbz5s0UKlQo1ceEX8i2RUQkfQrCkum2bNlCmTJlaN68OaVLlyZPnjwXvK48efJQunRpzCyAFWYfRYsWpVChQlldRlD4+eefueGGG2jevDlr1qzh0UcfZdCgQcyZMyfNZdauXcsNN9zA9ddfz5o1a5g5cybz5s1j+PDhZ7U9ceIE3bp14+qrrw7ItkVE5NwUhCVTxcTEMGTIELZv346ZUbly5XMus3TpUpo2bUrBggUJDw+nSZMmbNiwATh7aMTUqVMpWLAgn3zyCTVr1iR//vx06tSJgwcPMnv2bKpVq0Z4eDi9evXi77//9m2jVatW9OvXjwcffJAiRYpQpEgRhg4dSmJiYpp1nThxgmHDhlG+fHny589P48aNWbRokV/H4UzdX3zxBVdccQX58+cnOjqaH374IVm7Dz74gKioKPLmzUuFChUYPXo0Se/9nXJoxAcffEDdunXJly8fRYsWpWXLlvz222+++fPnz/f1KkZERPD4449z4sQJv2quXLkyo0aNIiYmhkKFClGhQgVmzpzJgQMH6NatGwULFqRatWosXrw42XJLly7liiuuICwsjFKlSjFkyJBk2zx69CgxMTEULFiQUqVK8eyzz5617Ys51oEyadIkypYty4QJE4iMjOSee+6hT58+Zz3kJKmZM2dSq1Ytnn76aapWrUrLli154YUXmDhxIocPH07WdtiwYdStW5dbb701INsWEZFzUxCWTPXqq68ycuRIypcvz+7du1m1alW67U+dOkXnzp258sorWbt2Ld9//z2DBw8mV65caS5z/Phxxo4dy4wZM/jiiy+IjY3l5ptvZtq0acyZM4ePPvqIjz/+mNdffz3ZcjNmzCAxMZEVK1bw5ptv8tZbbzFu3Lg0t3PnnXfy9ddf85///IcNGzbQp08fOnbsyNq1a/0+Ho8++ihjxozhhx9+oFixYvTs2dMXdFevXs2tt95K165dWb9+PWPGjOG5557jtddeS3Vde/bsoVu3bvTp04f4+HiWLl1Kr169fPMXLVpEz549GThwIBs3bmTy5MnMnj2bxx57zO96x40bR5MmTfjhhx+47bbb6NOnDz169OCGG24gLi6Oq6++mjvuuINjx44BsGvXLtq3b0+DBg1Ys2YN//73v3nvvfeSPc764Ycf5rPPPmPOnDl88cUXrFmzhqVLlwb8WD/77LPp3qqvYMGCLFu2LM3lV6xYQZs2bZJNa9u2LbGxsZw8eTLVZY4fP05YWFiyafny5ePYsWOsXr3aN23BggV8/PHHTJgwIWDbFhERPzjnsuSnUaNGLss8eZnfP3Wm1vHrZ1ONmn7/BLsXX3zRVapUya+2+/fvd4BbsmRJqvO/+uorB7i9e/c655ybMmWKA9yPP/7oa/PQQw+5kJAQXxvnnOvTp4/r0KGD73XLli1dtWrVXGJiom/aM88848qVK5eszYABA5xzzm3ZssWZmfvll1+S1dO5c2fXv3//c+7Xmbo//fRT37RvvvnGAW7Hjh3OOed69OjhrrnmmmTLPfnkk2nWtHr1age4hISEVLd51VVXuVGjRiWb9uGHH7oCBQok2++0VKpUyXXr1s33+vDhww5wgwYN8k37+eefHeBWrVrlnHPusccec1WrVnWnT5/2tZkyZYrLkyePO3LkiDt8+LDLkyePe/fdd5OtNzw83PXp08c559+xTrnd1Ozfv99t3rw53Z+jR4+muXy1atXc008/nWza119/7QD366+/prrMokWLnJm56dOnu5MnT7qdO3e6q666ygHuP//5j3POuV27drkyZcq47777znd8ChQocNHbFhGR/wfEulTyqHqE5ZJWtGhRYmJiaNu2LR06dODll18+50M28ubNS40aNXyvS5UqRenSpSlevHiyab///nuy5Zo2bZpsrHGzZs3YtWsXhw4dOmsbP/zwA845atWqlaxHccGCBWzdutXv/atbt67v9zP3TD5TV3x8PC1atEjW/sorr0yzpnr16tG6dWvq1KnDzTffzBtvvMHevXt981evXs3o0aOT1dujRw+OHDnCnj17zrveggULkj9/fqKionzTSpUqddY+NG3alJCQ//+v5sorr+TEiRNs2bKFrVu3cuLECZo1a5ZsvUnXGahjXbRoUapWrZruT758+fxenz/atGnDSy+9xIABAwgLC6N69erccMMNAL5j0qtXL/r3788VV1wR0G2LiMi5KQjLJW/KlCl8//33XH311cybN48aNWqkOz40NDT57bHNjNy5c581Lb3xv+eSmJiImbFq1Sri4uJ8P/Hx8UyePNnv9SSt60wI96eu1C4OzJUrF4sXL2bx4sXUrVuXf//731SrVs03fCAxMZEnn3wyWb3r1q1j8+bNlChR4rzrPVNHIPchNYE61hc7NKJ06dLJxlsD/Pbbb4SGhib7kJXSP/7xDw4cOMD27dvZt28fnTt3BqBKlSoAfPnllzz99NOEhoYSGhrK3XffzZEjRwgNDeWtt966qG2LiEj69EANyRbq1atHvXr1GDZsGO3bt2fatGm0bds2oNv4/vvvcc75Atp3331H2bJlueyyy85q26BBA5xz7Nmzh2uuuSagdZwRGRnJ8uXLk0375ptvKF++fJp3ijAzmjVrRrNmzRg5ciS1a9dm5syZ1KtXj4YNG/Ljjz9m6gNMIiMjmTVrFomJib4e0G+++YY8efJw+eWXk5iYSO7cufnuu+98wfDIkSNs2LCByy+/HAjcse7Xrx+33XZbum3KlSuX5rxmzZrx4YcfJpv22WefER0dfdYHhJTMzNfj/95771GhQgUaNmwIwPr165O1nTt3LqNHj2blypW+ei5m2yIikjYFYbmk/fzzz7z55pt06tSJcuXKsW3bNtatW0f//v0Dvq1ff/2VwYMHc//997N+/XpefPFFRowYkWrb6tWr07NnT2JiYhg7diwNGzbkjz/+YMmSJVSpUoWuXbtedD0PPfQQjRs35qmnnqJHjx6sWrWKsWPHpnpXBfAE988//5y2bdtSqlQp1qxZw44dO6hVqxYAI0eO5MYbb6RSpUrcdttthIaGsmHDBlauXMkLL7xw0fWm5v7772fcuHHcf//9PPjgg2zbto3hw4czcOBA8ufPD8Ddd9/NsGHDKFGiBGXLlmXUqFGcPn3at45AHeuiRYtStGjRC96Xfv368dprrzF48GDuu+8+li9fztSpU3nvvfd8bV577TVee+01fvzxR9+0F198kXbt2hESEsIHH3zAmDFjmDVrlu+Czzp16iTbTmxsLCEhIcmm+7NtERE5fwrCcknLnz8///vf/7j11lvZt28fpUqVomfPngwbNizg2+rZsyenT5/miiuuwMy4++67GTJkSJrtp0yZwujRo3nkkUfYuXMnRYsWpUmTJgHrIW7YsCH//e9/efLJJ3n22WcpVaqUL0SmJjw8nOXLlzNhwgQOHDhAhQoVeOKJJ7jjjjsAz10GFixYwDPPPMNLL71EaGgo1atXP68nuJ2vcuXK8cknnzB06FDq169P4cKF6dGjR7Iw/9JLL3HkyBG6dOlC/vz5GTRoEEeOHEm2now+1v6IiIhg4cKFDBkyhDfeeIOyZcsyfvx4br75Zl+bffv28dNPPyVb7pNPPmH06NEcP36cevXqMXfuXNq3bx/wbYuIyPkzl+SepJkpOjraZdnjUJ8K97tpVERFv9rNei71J0WlJvLHeL/bSuZo1aoVderUSfPWZCIiIpJ9mdlq51x0yum6WE5EREREgpKCsGSp7du3p3sV/7lulXap6tevX5r71K9fv6wu7yzLli1L9+8gIiKSE+WoMcKVhy/wq11C2LnbSOYoW7YscXFx6c7PDEuWLAno+kaNGsXDDz+c6rzU7kKR1aKjo9P9O4iIiOREOSoIS/YTGhqaqbfzyiwlS5akZMmSWV2G3/Lly5cj/w4iIiLp0dAIEREREQlKCsKS6RITE7nvvvsoVqwYZhbwYQmpWbJkCWbGvn37MnxbmSk2NhYzIyEhwe9lYmJiuPHGGzOuKBERkWxCQVgy3cKFC5kyZQrz589n9+7dNG/ePMO32bx5c3bv3k2xYsUyfFsCx48fZ9CgQRQvXpwCBQrQqVMndu7cme4yhw8fZvDgwVSqVIl8+fLRvHlzVq1alaxNTEwMZpbsp2nTpsna7Nmzh169elG6dGny589PvXr1mDFjRsD3UUREsj8FYcl0W7ZsoUyZMjRv3pzSpUuTJ0+eDN9mnjx5KF26tO/xyZKxBg8ezJw5c3jvvfdYtmwZhw4d4sYbb0z2xLiU+vbty6JFi5g2bRrr16+nTZs2tG7dml27diVr17p1a3bv3u37WbhwYbL5vXv3Jj4+nrlz57JhwwZ69+5Nr169WLp0aYbsq4iIZF8KwpKpYmJiGDJkCNu3b8fMqFy58jmXadWqFffffz+PPfYYxYsXp2TJkjz88MMkJib62vz555/06dOHIkWKkC9fPlq3bs3GjRt981MOjTh48CC9evWiZMmShIWFUaVKFcaNG+drf/DgQe69915KlixJoUKFaNmyJf4+AGbq1KkULFiQTz75hJo1a5I/f346derEwYMHmT17NtWqVSM8PJxevXrx999/+5Y7fvw4gwcPplSpUoSFhdG0aVO++eabZOv+9NNPqVmzJmFhYVx11VX873//O2v73377LS1btiR//vyUK1eO/v37c+jQIb9qD4SDBw/y73//mxdffJHrr7+ehg0bMn36dNatW8fnn3+e6jJ///03c+bMYcyYMbRq1YqqVavy1FNPUbVqVd54441kbfPmzUvp0qV9Pykfm/ztt98yYMAArrjiCqpUqcJDDz1EhQoVWLlyZYbts4iIZE8KwpKpXn31VUaOHEn58uXZvXv3WV99p2XGjBmEhoby7bff8tprrzFu3Dhmzpzpmx8TE8P333/P3LlzWblyJfnz56ddu3bJgmZSI0aMYP369Xz88cf89NNPTJ48mXLlygHgnKNDhw7s2rWLjz/+mDVr1nD11Vdz7bXXsnv3br/qPX78OGPHjmXGjBl88cUXxMbGcvPNNzNt2jTmzJnDRx99xMcff8zrr7/uW+aRRx5h5syZTJ48mTVr1hAVFUW7du1829yxYwc33XQT119/PXFxcQwaNIhHHnkk2XbP9KR26tSJtWvX8sEHHxAXF8ddd93lV91npHdP4YIFC6b7iODVq1dz8uRJ2rRp45tWoUIFIiMj+fbbb1Nd5tSpU5w+fZqwsOT3NsyXL99ZHwa++eYbSpYsSfXq1bnnnnv4/fffk82/8sormTVrFvv37ycxMZG5c+eyd+9eWrdufV7HQEREcj7dPk0yVXh4OIUKFSJXrlyULl3a7+Vq1arFqFGjAKhevTr/+te/+OKLL+jevTubN29m3rx5fP3111x99dUATJ8+nYoVKzJjxgz69u171vp++eUXGjZsSJMmTQCoVKmSb95XX31FXFwce/fuJV++fAA888wzzJ8/n+nTp58VPlNz6tQpJk6cSI0aNQDo0aMHr7zyCr/99hvFixcHoHPnznz11Vc89NBDHDlyhDfeeIO3336bDh06ADBp0iS+/PJLJk6cyD//+U/eeOMNKlasyPjx4zEzatasyf/+9z+eeOIJ33ZffPFFbr/9dh566CHftDfeeIMGDRrw+++/+31Lt3PdU/jMcUnNnj17yJUrl28/zyhVqhR79uxJdZlChQrRrFkz/vnPf1KnTh1Kly7Ne++9x4oVK5Ld1q1du3Z07dqViIgIEhISGDFiBNdeey2rV68mb968AMyaNYtu3bpRvHhxQkNDyZs3L++99x7169f3a99FRCR4KAhLtlC3bt1kr8uWLevrCYyPjyckJIRmzZr55oeHhxMVFcWmTZtSXV///v255ZZbWL16Nddffz0dO3akZcuWgKdH8+jRo5QoUSLZMseOHWPr1q1+1Zs3b15fCAZPCCxdunSycFiqVClffVu3buXkyZO0aNHCNz9Xrlw0a9bM1yY+Pp6mTZsmG+ecdJ/P1L5ly5ZkveXOOd82/A3CWXFP4enTp3PXXXdRvnx5cuXKRcOGDenevTurV6/2tenWrZvv96ioKBo1akSlSpVYsGABXbt2BTy9/fv27ePzzz+nePHifPTRR/Tu3ZulS5dSr169TN8vERG5dPkVhM2sHfAqkAt42zk3JsX8GOBF4MxVLa85594OYJ0S5HLnzp3stZklGyOclrQujmvfvj2//PILn3zyCV988QUdOnTg1ltvZcqUKSQmJlKqVCmWLVt21nL+PhUuNDT5Py0zC/g+pCYxMZG+ffsyZMiQs+adGfrhj3M9Vvmqq67ik08+SXVe6dKlOX36NPv27Uv2YeK3337jqquuSnOdl19+OV9//TVHjhzh0KFDlClThttvv50qVaqkuUzZsmUpX748mzdvBjxhf8KECcTFxflCb7169Vi2bBkTJkzg7bf135KIiPy/cwZhM8sFTASuB3YCq8xsnnMuZVfbTOfcwAyoUSRdkZGRJCYmsmLFCt/QiEOHDrF+/XruvPPONJcrXrw4vXr1olevXrRv357u3bszadIkGjZsyG+//UZISEi6ISyQLr/8cvLkycPy5cu5/PLLATh9+jQrVqygR48egGc/58yZg3POF46/++67ZOtp2LAhGzduvOge3YsZGtGoUSNy587NZ5995qt9586dxMfH+3WrvAIFClCgQAH+/PNPFi1axAsvvJBm23379rFr1y7KlCkDwNGjRwFPb3pSuXLl8utDh4iIBBd/eoSbAFucc9sAzOx9oDOQ+nfOIpmsWrVqdO7cmfvuu4+33nqLwoUL8/jjj3PZZZf5glhKI0eOpGHDhtSuXZtTp07xwQcfUKVKFfLmzUvr1q1p0aIFnTt35oUXXqBmzZrs2bOHTz/9lNatW6fbq3mhChQoQP/+/Rk2bBjFixcnIiLCN6b4/vvvB6Bfv36MHTuWwYMHc//997N+/XomTZqUbD3Dhg2jadOm9OvXj/vuu49ChQrx448/Mn/+fN58802/67mYIB0eHs7dd9/NI488QsmSJSlWrBj/+Mc/qFu3brIL1mrWrMnAgQMZONDz+XnRokUkJiZSs2ZNtmzZwtChQ6lZs6bvw8xff/3FU089xc0330yZMmVISEjg0UcfpWTJknTp0sW3zqpVq3L//ffz0ksvUaxYMT766CM+++wz5s6de8H7JCIiOZM/d40oB+xI8nqnd1pKN5vZOjObbWYVUluRmd1rZrFmFrt3794LKFckdVOmTKFJkyZ06tSJJk2acPToUT799NM0ey7z5s3L448/Tr169WjRogWHDx9m/vz5gGcowsKFC7n22mu55557qFGjBrfddhs//fQTZcuWzbB9eP7557n99tu58847qV+/PuvWrePTTz/19XZWrFiRDz74gE8//ZR69erxyiuvMGZMslFK1K1bl6VLl5KQkEDLli2pV68ejz76KKVKlcqwulMzbtw4unTpwu23306LFi0oWLAg8+fPT9ZT+9NPPyV70t/BgwcZOHAgNWvWpHfv3lx55ZUsWrTIN6QkV65crF+/ns6dO1O9enX69OlDjRo1WLFiBYUKFQI8Q2gWLlxIiRIl6NixI3Xr1uWdd95hypQpdOzYMVOPgYiIXPrszIU0aTYwuwVo55zr633dC7gi6TAIMysG/OWcO25m9wG3O+euTW+90dHRzt/7svqr8vAFfrVLCEu9lzA1UREV/Wo367lTfq8z8sd4v9uKiIiIyMUxs9XOueiU0/3pEd4FJO3hLc//XxQHgHNuv3PuuPfl20CjCy1URERERCQz+BOEVwHVzCzCzPIA3YB5SRuYWZkkLzsB6vIUv2zfvj3dBzds3749q0s8S/v27dOs99lnn83q8kRERMRP57xYzjl3yswGAovw3D5tsnNuo5mNAmKdc/OAB8ysE3AK+AOIycCaJQcpW7ZsuncoyMgxuRfq7bffTvOJdSkf9ysiIiKXLr/uI+ycWwgsTDFtZJLfHwUeDWxpEgxCQ0Oz5OENF+N87scrIiIily5/hkaIiIiIiOQ4CsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQUlBWERERESCkoKwiIiIiAQlBWERERERCUoKwiIiIiISlBSERURERCQoKQiLiIiISFBSEBYRERGRoKQgLCIiIiJBSUFYRERERIKSgrCIiIiIBCUFYREREREJSgrCIiIiIhKUFIRFREREJCgpCIuIiIhIUFIQFhEREZGgpCAsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaDkVxA2s3Zm9pOZbTGz4em0u9nMnJlFB65EEREREZHAO2cQNrNcwESgPVAL6G5mtVJpVwh4EPg+0EWKiIiIiASaPz3CTYAtzrltzrkTwPtA51TaPQM8DxwLYH0iIiIiIhnCnyBcDtiR5PVO7zQfM2sIVHDOLUhvRWZ2r5nFmlns3r17z7tYEREREZFAueiL5cwsBHgZeOhcbZ1zbznnop1z0SVKlLjYTYuIiIiIXDB/gvAuoEKS1+W9084oBNQBlphZAtAUmKcL5kRERETkUuZPEF4FVDOzCDPLA3QD5p2Z6Zw76Jwr7pyr7JyrDHwHdHLOxWZIxSIiIiIiAXDOIOycOwUMBBYB8cAs59xGMxtlZp0yukARERERkYwQ6k8j59xCYGGKaSPTaNvq4ssSEREREclYerKciIiIiAQlBWERERERCUoKwiIiIiISlBSERURERCQoKQiLiIiISFBSEBYRERGRoKQgLCIiIiJBSUFYRERERIKSgrCIiIiIBCUFYREREREJSgrCIiIiIhKUFIRFREREJCgpCIuIiIhIUFIQFhEREZGgpCAsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQUlBWERERESCkoKwiIiIiAQlBWERERERCUoKwiIiIiISlBSERURERCQoKQiLiIiISFDyKwibWTsz+8nMtpjZ8FTm9zOz9WYWZ2bfmFmtwJcqIiIiIhI45wzCZpYLmAi0B2oB3VMJuv9xzkU55+oDLwAvB7pQEREREZFA8qdHuAmwxTm3zTl3Angf6Jy0gXPuUJKXBQAXuBJFRERERAIv1I825YAdSV7vBK5I2cjMBgD/APIA16a2IjO7F7gXoGLFiudbq4iIiIhIwATsYjnn3ETn3OXAMGBEGm3ecs5FO+eiS5QoEahNi4iIiIicN3+C8C6gQpLX5b3T0vI+cNNF1CQiIiIikuH8CcKrgGpmFmFmeYBuwLykDcysWpKXHYDNgStRRERERCTwzjlG2Dl3yswGAouAXMBk59xGMxsFxDrn5gEDzaw1cBL4E+iTkUWLiIiIiFwsfy6Wwzm3EFiYYtrIJL8/GOC6REREREQylJ4sJyIiIiJBSUFYRERERIKSgrCIiIiIBCUFYREREREJSgrCIiIiIhKUFIRFREREJCgpCIuIiIhIUFIQFhEREZGgpCAsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKIVmdQEiwaDy8AV+t00Y0yEDKxEREZEz1CMsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQUlBWERERESCkoKwiIiIiAQlBWERERERCUoKwiIiIiISlBSERURERCQo+RWEzaydmf1kZlvMbHgq8/9hZpvMbJ2ZfWFmlQJfqoiIiIhI4JwzCJtZLmAi0B6oBXQ3s1opmq0Bop1zdYHZwAuBLlREREREJJD86RFuAmxxzm1zzp0A3gc6J23gnPvKOXfU+/I7oHxgyxQRERERCSx/gnA5YEeS1zu909JyN/BJajPM7F4zizWz2L179/pfpYiIiIhIgAX0YjkzuwOIBl5Mbb5z7i3nXLRzLrpEiRKB3LSIiIiIyHkJ9aPNLqBCktflvdOSMbPWwONAS+fc8cCUJyIiIiKSMfzpEV4FVDOzCDPLA3QD5iVtYGYNgDeBTs653wNfpoiIiIhIYJ0zCDvnTgEDgUVAPDDLObfRzEaZWSdvsxeBgsB/zSzOzOalsToRERERkUuCP0MjcM4tBBammDYyye+tA1yXiIiIiEiG0pPlRERERCQoKQiLiIiISFBSEBYRERGRoKQgLCIiIiJBSUFYRERERIKSgrCIiIiIBCUFYREREREJSgrCIiIiIhKUFIRFREREJCgpCIuIiIhIUFIQFhEREZGgpCAsIiIiIkFJQVhEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQSk0qwuQHOCp8PNoe9CvZvE1I/1eZeSP8f5vX0RERMRLPcIiIiIiEpQUhEVEREQkKGlohGSqqGlRfrWblcF1iIiIiKhHWERERESCkoKwiIiIiAQlDY2QNFUevsCvdglhGVyIiIiISAZQj7CIiIiIBCUFYREREREJSn4FYTNrZ2Y/mdkWMxueyvyrzewHMztlZrcEvkwRERERkcA6ZxA2s1zARKA9UAvobma1UjTbDsQA/wl0gSIiIiIiGcGfi+WaAFucc9sAzOx9oDOw6UwD51yCd15iBtQoIiIiIhJw/gyNKAfsSPJ6p3eaiIiIiEi2lakXy5nZvWYWa2axe/fuzcxNi4iIiIgk408Q3gVUSPK6vHfaeXPOveWci3bORZcoUeJCViEiIiIiEhD+BOFVQDUzizCzPEA3YF7GliUiIiIikrHOGYSdc6eAgcAiIB6Y5ZzbaGajzKwTgJk1NrOdwK3Am2a2MSOLFhERERG5WH49Ytk5txBYmGLayCS/r8IzZEJELtZT4X62O5ixdYiIiORwerKciIiIiAQlBWERERERCUoKwiIiIiISlBSERURERCQoKQiLiIiISFBSEBYRERGRoKQgLCIiIiJBSUFYRERERIKSgrCIiIiIBCW/niwnItlbfM1Iv9pF/hifwZWIiIhcOtQjLCIiIiJBST3CIpLtVB6+wO+2CWM6ZGAlIiKSnSkIi2RTUdOi/G47KwPrEBERya40NEJEREREgpKCsIiIiIgEJQVhEREREQlKCsIiIiIiEpQUhEVEREQkKCkIi4iIiEhQUhAWERERkaCkICwiIiIiQUlBWERERESCkp4sJyI521PhfrY7mLF1iIjIJUc9wiIiIiISlBSERURERCQoaWiEiIiIXJIqD1/gd9uEMR0ysBLJqRSERUSCnL9hQ0FDRHIaDY0QERERkaCkHmERERGRVMTXjPSrXeSP8RlciWQUv4KwmbUDXgVyAW8758akmJ8XeAdoBOwHbnfOJQS2VBGRrOfvGyPozVFE5FJ3ziBsZrmAicD1wE5glZnNc85tStLsbuBP51xVM+sGPA/cnhEFi4hkhKhpUX61m5XBdYjIBdI9w+UC+NMj3ATY4pzbBmBm7wOdgaRBuDPwlPf32cBrZmbOORfAWkVERC55uvhQLlU6N8/mTxAuB+xI8noncEVabZxzp8zsIFAM2BeIIkVE5BLgb48b+N3rpqEmIpegIPq3bufqtDWzW4B2zrm+3te9gCuccwOTtNngbbPT+3qrt82+FOu6F7jX+7IG8FOgdiQDFUeBPpB0PANHxzKwdDwDS8czcHQsA0vHM7Cyy/Gs5JwrkXKiPz3Cu4AKSV6X905Lrc1OMwsFwvFcNJeMc+4t4C1/K74UmFmscy46q+vIKXQ8A0fHMrB0PANLxzNwdCwDS8czsLL78fTnPsKrgGpmFmFmeYBuwLwUbeYBfby/3wJ8qfHBIiIiInIpO2ePsHfM70BgEZ7bp012zm00s1FArHNuHvBvYLqZbQH+wBOWRUREREQuWX7dR9g5txBYmGLayCS/HwNuDWxpl4xsNZQjG9DxDBwdy8DS8QwsHc/A0bEMLB3PwMrWx/OcF8uJiIiIiORE/owRFhERERHJcRSERURERCQoKQiLiIiISFDy62K5YGRmRQGcc39kdS0iZ5hZKTxPcgTY5Zz7LSvrEUmPmRV0zv2V1XWIiKRFF8slYWYVgReA64ADgAGXAV8Cw51zCVlWXDam8HbxzKw+MAnPw2rOPNCmPJ7z9H7n3A9ZU1n2ZWYGNCHJuQms1D3QA8fMtjvnKmZ1HTmJmdV0zv2Y1XXkBGa23jkXldV15BTZ9YOveoSTmwmMA3o6504DmFkuPLeGex9omnWlZT9phTczO4DC2/maCtznnPs+6UQzawpMAeplRVHZlZm1AV4HNpP8g0VVM7vfObc4y4rLZszsH2nNAgpmZi1BYjGgDxd+MrOuac0CSmdmLUFgE9nw3FQQTq64c25m0gneQPy+mT2TRTVlZ1NReAuUAimPI4Bz7jszK5AVBWVzrwKtU37LY2YReO6ZHpkVRWVTzwIvAqdSmafrUC6AmY1PaxZQOBNLyQlmAjOA1L7pCcvkWrK9nPjBV0E4udVm9jowDdjhnVYBz+Oj12RZVdmXwlvgfGJmC4B3SH5u9gY+zbKqsq9QYGcq03cBuTO5luzuB+Aj59zqlDPMrG8W1JMT3Ak8BBxPZV73TK4lu1sHvOSc25Byhpm1zoJ6srsc98FXQTi53sDdwNMkHzd45jHScn4U3gLEOfeAmbUHOpP83JzoffKjnJ/JwCoze5/k52Y39G/9fN0J7E9jXnRmFpKDrAI2OOe+TTnDzJ7K/HKytcHAoTTmdcnEOnKKHPfBVxfLSYZKI7zNU3iTrGZmkaR+bm7KuqpEfHctOuacO5rVtYgkZWY1gD+cc3tTmVcqO14MryCchJmF4ukRvonkb45zgX87505mUWkiaTKzt5xz92Z1HSIp6dyUrJbkfb0LUNY7We/r4qMgnISZvYfndlTT+P/xg+XxjBEu6py7PYtKy5bMLBx4FE+vWyk8Fyv8juc/oDHOuQNZV132cua+1qnNAtY658pnZj05mZl94pxrn9V1ZBc6NzOXzs/zo/f1wEryvn4TUJIc8L6uMcLJNXLOVU8xbSfwnZn9LysKyuZm4bkH8zXOuT0AZlYaiPHOa5N1pWU7e4Ff8ISLM5z3dcksqSgbM7OGac0C6mdiKTmBzs0A0/kZUHpfD6wz7+utUryv9yGbvq8rCCf3h5ndCsxxziUCmFkInvsI/5mllWVPlZ1zzyed4P2HM8bM7syimrKrbcB1zrntKWeY2Y5U2kv6VgFfkzy8nVE4c0vJ9nRuBp7Oz8DR+3pgpfW+/ryZ3ZVFNV0UBeHkugHPAxO9D30Az386X3nnyfn5xcweAaadGUDvfcpcDP9/pb74ZxxQBDgrbOB5GqKcn3g897jenHKGwtt5G4fOzUDT+Rk4el8PrBz3vq4xwimkcSX5XOdcfNZVlT2ZWRFgOJ7jeeYr0t/w3I5ujHNOn8bPg5nVJPW7HOjcPE9mdguw3jn3UyrzbnLOfZT5VWVfOjcDS+dnYOl9PXBy4vt6trz5cUYxs2HAf/CMb/ve+wPwnpkNz7LCsinn3J/OuWHOuZrOuaLen0jn3DA8A+3FT95P4O/j+ap0pffH0Ll5QZxzs1MLGV5FMrWYbE7nZuDp/Awcva8HVk58X1ePcBLegfO1U95OxczyABudc9WyprKcx8y2O+ey3TPJs4rOzcyjc/P86NzMXDo/z4/Oz8yTXc9NjRFOLhHPfQZ/STG9jHeenAczW5fWLDy3UxP/6dwMIJ2bAaVzM8B0fgaUzs8AyonnpoJwcoOBL8xsM/8/6LsiUBUYmFVFZWOlgLacfWWuAWc9OlTSNRidm4GkczNwBqNzM9B0fgbOYHR+BlKOOzcVhJNwzn1qZtWBJiQfVL/KOXc66yrLtj4GCjrn4lLOMLMlmV5NNqZzM+B0bgaIzs0MofMzQHR+BlyOOzc1RlhEREREgpLuGiEiIiIiQUlBWERERESCkoKwiIiIiAQlBWERERERCUoKwiIiIiISlP4P7V4/b/1UwdwAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 864x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "df.plot.bar(title=\"Comparing simulators and Aspen-M\", figsize=(12, 6))\n",
    "text = f\"f_free = {f_free_aspen:.3f} \\nf_simple_noise_model = {f_simple:.3f} \\nf_noise_model = {f_noisy_aspen:.3f}\"\n",
    "plt.text(1, 0.5, text, fontsize=14)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6da72e22-a594-41e7-a4e2-cb82cfe0726a",
   "metadata": {},
   "source": [
    "We see that compared to the full noise model, the simple model is less accurate, however, it is still a significant improvement over the noise-free case and has far fewer parameters in the model. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b2c6b2ba-cb54-43f6-9a89-e29066b5192a",
   "metadata": {},
   "source": [
    "## Summary\n",
    "\n",
    "In this notebook, we showed how to construct a noise model for Rigetti Aspen-M-3 based only on the available calibration data.\n",
    "We used a coarse assumption of gate-independent single-qubit depolarizing noise and gate-dependant two-qubit noise. Our qubit-dependent model could be improved in many ways. We could add gate-dependence noise, or change the depolarizing channel to Pauli channels. "
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "a86a9ff96c3b20b9c94872c880ce00e370ced2a9959f9303c435918fb220e358"
  },
  "kernelspec": {
   "display_name": "conda_braket",
   "language": "python",
   "name": "conda_braket"
  },
  "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.5"
  },
  "toc-autonumbering": false,
  "toc-showcode": false,
  "toc-showmarkdowntxt": false,
  "toc-showtags": true
 },
 "nbformat": 4,
 "nbformat_minor": 5
}