{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "226dfb5f",
   "metadata": {},
   "source": [
    "# Exploring Hypothetical Hardware\n",
    "\n",
    "The [previous notebooks](0_getting_started.ipynb) showed how to use the quantum\n",
    "resource estimator with built-in models like `GateBased`, `SurfaceCode`, and\n",
    "`Litinski19Factory`.  These models capture well-studied hardware assumptions,\n",
    "but what if your hardware doesn't match them?\n",
    "\n",
    "In this notebook we build **custom models** to answer three \"what if\" questions:\n",
    "\n",
    "- **Architecture model:** What if my qubits had different gate speeds and error rates?\n",
    "- **QEC ISA transform:** What if I had a QEC code with a higher threshold or lower overhead?\n",
    "- **Factory transform:** What if my magic state distillation protocol were different?\n",
    "\n",
    "Each model plugs into the estimator's layered architecture.  You implement a\n",
    "small Python class and the estimator takes care of the rest: enumerating\n",
    "configurations, pruning infeasible ones, and finding the Pareto-optimal results.\n",
    "\n",
    "We finish by composing all three custom layers into a full estimation stack and\n",
    "comparing it against the built-in models."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5dafb171",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/javascript": "// Copyright (c) Microsoft Corporation.\n// Licensed under the MIT License.\n\n// This file provides CodeMirror syntax highlighting for Q# magic cells\n// in classic Jupyter Notebooks. It does nothing in other (Jupyter Notebook 7,\n// VS Code, Azure Notebooks, etc.) environments.\n\n// Detect the prerequisites and do nothing if they don't exist.\nif (window.require && window.CodeMirror && window.Jupyter) {\n  // The simple mode plugin for CodeMirror is not loaded by default, so require it.\n  window.require([\"codemirror/addon/mode/simple\"], function defineMode() {\n    let rules = [\n      {\n        token: \"comment\",\n        regex: /(\\/\\/).*/,\n        beginWord: false,\n      },\n      {\n        token: \"string\",\n        regex: String.raw`^\\\"(?:[^\\\"\\\\]|\\\\[\\s\\S])*(?:\\\"|$)`,\n        beginWord: false,\n      },\n      {\n        token: \"keyword\",\n        regex: String.raw`(namespace|open|as|operation|function|body|adjoint|newtype|controlled|internal)\\b`,\n        beginWord: true,\n      },\n      {\n        token: \"keyword\",\n        regex: String.raw`(if|elif|else|repeat|until|fixup|for|in|return|fail|within|apply)\\b`,\n        beginWord: true,\n      },\n      {\n        token: \"keyword\",\n        regex: String.raw`(Adjoint|Controlled|Adj|Ctl|is|self|auto|distribute|invert|intrinsic)\\b`,\n        beginWord: true,\n      },\n      {\n        token: \"keyword\",\n        regex: String.raw`(let|set|use|borrow|mutable)\\b`,\n        beginWord: true,\n      },\n      {\n        token: \"operatorKeyword\",\n        regex: String.raw`(not|and|or)\\b|(w/)`,\n        beginWord: true,\n      },\n      {\n        token: \"operatorKeyword\",\n        regex: String.raw`(=)|(!)|(<)|(>)|(\\+)|(-)|(\\*)|(/)|(\\^)|(%)|(\\|)|(&&&)|(~~~)|(\\.\\.\\.)|(\\.\\.)|(\\?)`,\n        beginWord: false,\n      },\n      {\n        token: \"meta\",\n        regex: String.raw`(Int|BigInt|Double|Bool|Qubit|Pauli|Result|Range|String|Unit)\\b`,\n        beginWord: true,\n      },\n      {\n        token: \"atom\",\n        regex: String.raw`(true|false|Pauli(I|X|Y|Z)|One|Zero)\\b`,\n        beginWord: true,\n      },\n    ];\n    let simpleRules = [];\n    for (let rule of rules) {\n      simpleRules.push({\n        token: rule.token,\n        regex: new RegExp(rule.regex, \"g\"),\n        sol: rule.beginWord,\n      });\n      if (rule.beginWord) {\n        // Need an additional rule due to the fact that CodeMirror simple mode doesn't work with ^ token\n        simpleRules.push({\n          token: rule.token,\n          regex: new RegExp(String.raw`\\W` + rule.regex, \"g\"),\n          sol: false,\n        });\n      }\n    }\n\n    // Register the mode defined above with CodeMirror\n    window.CodeMirror.defineSimpleMode(\"qsharp\", { start: simpleRules });\n    window.CodeMirror.defineMIME(\"text/x-qsharp\", \"qsharp\");\n\n    // Tell Jupyter to associate %%qsharp magic cells with the qsharp mode\n    window.Jupyter.CodeCell.options_default.highlight_modes[\"qsharp\"] = {\n      reg: [/^%%qsharp/],\n    };\n\n    // Force re-highlighting of all cells the first time this code runs\n    for (const cell of window.Jupyter.notebook.get_cells()) {\n      cell.auto_highlight();\n    }\n  });\n}\n",
      "text/plain": []
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from dataclasses import KW_ONLY, dataclass, field\n",
    "from math import ceil\n",
    "from typing import Generator\n",
    "\n",
    "import qdk\n",
    "from qdk.qre import (\n",
    "    Architecture,\n",
    "    ISAContext,\n",
    "    ConstraintBound,\n",
    "    ISA,\n",
    "    ISARequirements,\n",
    "    ISATransform,\n",
    "    LOGICAL,\n",
    "    PHYSICAL,\n",
    "    constraint,\n",
    "    estimate,\n",
    "    linear_function,\n",
    "    plot_estimates,\n",
    ")\n",
    "from qdk.qre.models import GateBased, SurfaceCode, Litinski19Factory\n",
    "from qdk.qre.application import QSharpApplication\n",
    "from qdk.qre.instruction_ids import (\n",
    "    H, S, S_DAG, T, T_DAG, CNOT, CZ, MEAS_X, MEAS_Y, MEAS_Z, LATTICE_SURGERY,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4ed9261a",
   "metadata": {},
   "source": [
    "## Setup and baseline\n",
    "\n",
    "We start with a simple Q# program and a baseline estimation using the built-in\n",
    "models.  This gives us a reference point to compare against later."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ccd60531",
   "metadata": {
    "vscode": {
     "languageId": "qsharp"
    }
   },
   "outputs": [],
   "source": [
    "%%qsharp\n",
    "\n",
    "import Std.Arithmetic.*;\n",
    "import Std.Math.*;\n",
    "import Std.Convert.*;\n",
    "\n",
    "/// An 8-bit adder with single-qubit rotations.\n",
    "/// The rotations introduce T gates via synthesis, which exercises\n",
    "/// the magic state factory during resource estimation.\n",
    "operation AdderWithRotations() : Unit {\n",
    "    use a = Qubit[8];\n",
    "    use b = Qubit[8];\n",
    "    for i in 0..7 {\n",
    "        Ry(PI() / IntAsDouble(i + 2), a[i]);\n",
    "    }\n",
    "    RippleCarryCGIncByLE(a, b);\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "41bfdcd6",
   "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>name</th>\n",
       "      <th>qubits</th>\n",
       "      <th>runtime</th>\n",
       "      <th>error</th>\n",
       "      <th>physical_compute_qubits</th>\n",
       "      <th>physical_factory_qubits</th>\n",
       "      <th>physical_memory_qubits</th>\n",
       "      <th>factories</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Built-in models</td>\n",
       "      <td>9567</td>\n",
       "      <td>0 days 00:00:00.000225</td>\n",
       "      <td>0.007312</td>\n",
       "      <td>3087</td>\n",
       "      <td>6480</td>\n",
       "      <td>0</td>\n",
       "      <td>8×T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Built-in models</td>\n",
       "      <td>13221</td>\n",
       "      <td>0 days 00:00:00.000141</td>\n",
       "      <td>0.009416</td>\n",
       "      <td>1071</td>\n",
       "      <td>12150</td>\n",
       "      <td>0</td>\n",
       "      <td>15×T</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "              name  qubits                runtime     error  \\\n",
       "0  Built-in models    9567 0 days 00:00:00.000225  0.007312   \n",
       "1  Built-in models   13221 0 days 00:00:00.000141  0.009416   \n",
       "\n",
       "   physical_compute_qubits  physical_factory_qubits  physical_memory_qubits  \\\n",
       "0                     3087                     6480                       0   \n",
       "1                     1071                    12150                       0   \n",
       "\n",
       "  factories  \n",
       "0       8×T  \n",
       "1      15×T  "
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "app = QSharpApplication(qdk.code.AdderWithRotations)\n",
    "\n",
    "# Baseline: built-in architecture, QEC, and factory\n",
    "baseline_arch = GateBased(error_rate=1e-4, gate_time=100, measurement_time=500)\n",
    "baseline = estimate(\n",
    "    app, baseline_arch,\n",
    "    isa_query=SurfaceCode.q() * Litinski19Factory.q(),\n",
    "    max_error=0.01,\n",
    "    name=\"Built-in models\",\n",
    ")\n",
    "baseline.add_qubit_partition_column()\n",
    "baseline.add_factory_summary_column()\n",
    "baseline.as_frame()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef447c3b",
   "metadata": {},
   "source": [
    "## Building a custom architecture model\n",
    "\n",
    "An `Architecture` defines the physical instruction set of a quantum computer:\n",
    "the native operations the hardware can perform, each with a gate time and error\n",
    "rate.\n",
    "\n",
    "The built-in `GateBased` architecture uses a single error rate for all gates.\n",
    "In practice, two-qubit gates are typically noisier and slower than single-qubit\n",
    "gates, and measurement fidelities may differ from gate fidelities.  Let's build\n",
    "a more flexible architecture that models these differences.\n",
    "\n",
    "To implement a custom architecture, subclass `Architecture` and implement one\n",
    "method: `provided_isa`.  This method receives an `ISAContext` and uses it to\n",
    "register each physical instruction."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "5e98f44a",
   "metadata": {},
   "outputs": [],
   "source": [
    "@dataclass\n",
    "class ParameterizedArchitecture(Architecture):\n",
    "    \"\"\"\n",
    "    A configurable architecture with separate timing and error parameters for\n",
    "    single-qubit gates, two-qubit gates, and measurements.\n",
    "\n",
    "    Unlike the built-in ``GateBased`` model (which uses one error rate for all\n",
    "    operations), this model lets you specify different error rates for each gate\n",
    "    class, reflecting hardware where two-qubit gates are noisier than\n",
    "    single-qubit gates.\n",
    "    \"\"\"\n",
    "\n",
    "    _: KW_ONLY\n",
    "    single_qubit_gate_time: int\n",
    "    two_qubit_gate_time: int\n",
    "    measurement_time: int\n",
    "    single_qubit_error_rate: float\n",
    "    two_qubit_error_rate: float\n",
    "    measurement_error_rate: float\n",
    "\n",
    "    def provided_isa(self, ctx: ISAContext) -> ISA:\n",
    "        instructions = []\n",
    "\n",
    "        # Single-qubit gates: Clifford gates and T gates\n",
    "        for gate_id in [H, S, S_DAG, T, T_DAG]:\n",
    "            instructions.append(\n",
    "                ctx.add_instruction(\n",
    "                    gate_id,\n",
    "                    encoding=PHYSICAL,\n",
    "                    arity=1,\n",
    "                    time=self.single_qubit_gate_time,\n",
    "                    error_rate=self.single_qubit_error_rate,\n",
    "                )\n",
    "            )\n",
    "\n",
    "        # Measurements\n",
    "        for meas_id in [MEAS_X, MEAS_Y, MEAS_Z]:\n",
    "            instructions.append(\n",
    "                ctx.add_instruction(\n",
    "                    meas_id,\n",
    "                    encoding=PHYSICAL,\n",
    "                    arity=1,\n",
    "                    time=self.measurement_time,\n",
    "                    error_rate=self.measurement_error_rate,\n",
    "                )\n",
    "            )\n",
    "\n",
    "        # Two-qubit entangling gates\n",
    "        for gate_id in [CNOT, CZ]:\n",
    "            instructions.append(\n",
    "                ctx.add_instruction(\n",
    "                    gate_id,\n",
    "                    encoding=PHYSICAL,\n",
    "                    arity=2,\n",
    "                    time=self.two_qubit_gate_time,\n",
    "                    error_rate=self.two_qubit_error_rate,\n",
    "                )\n",
    "            )\n",
    "\n",
    "        return ctx.make_isa(*instructions)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee910f8e",
   "metadata": {},
   "source": [
    "### What if my qubits had different characteristics?\n",
    "\n",
    "Let's explore two hypothetical hardware scenarios:\n",
    "\n",
    "- **Fast but noisy:** short gate times, higher error rates (reminiscent of\n",
    "  some superconducting platforms)\n",
    "- **Slow but clean:** longer gate times, lower error rates (reminiscent of\n",
    "  some trapped-ion platforms)\n",
    "\n",
    "We use the built-in `SurfaceCode` and `Litinski19Factory` transforms so that we\n",
    "isolate the effect of the architecture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b7f6079c",
   "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>encoding</th>\n",
       "      <th>arity</th>\n",
       "      <th>space</th>\n",
       "      <th>time</th>\n",
       "      <th>error</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>id</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>H</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>50</td>\n",
       "      <td>0.0001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>MEAS_X</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>200</td>\n",
       "      <td>0.0005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>S</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>50</td>\n",
       "      <td>0.0001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>T</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>50</td>\n",
       "      <td>0.0001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>MEAS_Y</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>200</td>\n",
       "      <td>0.0005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>CZ</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>100</td>\n",
       "      <td>0.0005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>S_DAG</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>50</td>\n",
       "      <td>0.0001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>T_DAG</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>50</td>\n",
       "      <td>0.0001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>MEAS_Z</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>200</td>\n",
       "      <td>0.0005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>CNOT</th>\n",
       "      <td>PHYSICAL</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>100</td>\n",
       "      <td>0.0005</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        encoding  arity  space  time   error\n",
       "id                                          \n",
       "H       PHYSICAL      1      1    50  0.0001\n",
       "MEAS_X  PHYSICAL      1      1   200  0.0005\n",
       "S       PHYSICAL      1      1    50  0.0001\n",
       "T       PHYSICAL      1      1    50  0.0001\n",
       "MEAS_Y  PHYSICAL      1      1   200  0.0005\n",
       "CZ      PHYSICAL      2      2   100  0.0005\n",
       "S_DAG   PHYSICAL      1      1    50  0.0001\n",
       "T_DAG   PHYSICAL      1      1    50  0.0001\n",
       "MEAS_Z  PHYSICAL      1      1   200  0.0005\n",
       "CNOT    PHYSICAL      2      2   100  0.0005"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fast_noisy = ParameterizedArchitecture(\n",
    "    single_qubit_gate_time=50,\n",
    "    two_qubit_gate_time=100,\n",
    "    measurement_time=200,\n",
    "    single_qubit_error_rate=1e-4,\n",
    "    two_qubit_error_rate=5e-4,\n",
    "    measurement_error_rate=5e-4,\n",
    ")\n",
    "\n",
    "slow_clean = ParameterizedArchitecture(\n",
    "    single_qubit_gate_time=500,\n",
    "    two_qubit_gate_time=1000,\n",
    "    measurement_time=1000,\n",
    "    single_qubit_error_rate=1e-6,\n",
    "    two_qubit_error_rate=5e-6,\n",
    "    measurement_error_rate=1e-5,\n",
    ")\n",
    "\n",
    "# Inspect the physical ISA provided by the fast-but-noisy architecture\n",
    "fast_noisy.provided_isa(fast_noisy.context()).as_frame()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "de4ec193",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 23608 characters]",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isa_query = SurfaceCode.q() * Litinski19Factory.q()\n",
    "\n",
    "results = [\n",
    "    estimate(app, fast_noisy, isa_query, max_error=0.01, name=\"Fast but noisy\"),\n",
    "    estimate(app, slow_clean, isa_query, max_error=0.01, name=\"Slow but clean\"),\n",
    "    baseline,\n",
    "]\n",
    "\n",
    "plot_estimates(results, figsize=(8, 5), runtime_unit=\"ms\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f4cd8012",
   "metadata": {},
   "source": [
    "## Building a custom QEC ISA transform\n",
    "\n",
    "An `ISATransform` maps one instruction set to another.  A QEC transform\n",
    "consumes physical instructions and produces logical instructions by encoding\n",
    "qubits in an error-correcting code.\n",
    "\n",
    "The key relationships in any QEC code are:\n",
    "\n",
    "- **Space**: How many physical qubits per logical qubit?  (Depends on code\n",
    "  distance.)\n",
    "- **Time**: How long does one logical operation take?  (Depends on syndrome\n",
    "  extraction cycles.)\n",
    "- **Error**: What is the logical error rate?  (Depends on physical error rate,\n",
    "  code distance, and the code's threshold.)\n",
    "\n",
    "To implement a QEC transform, subclass `ISATransform` and implement two\n",
    "methods:\n",
    "\n",
    "1. `required_isa`: a **static method** that declares what physical\n",
    "   instructions the transform needs.  Because it is static, it cannot reference\n",
    "   instance fields (such as the error-correction threshold); instead, you\n",
    "   specify fixed bounds using `ConstraintBound`.  `ConstraintBound` supports\n",
    "   comparisons like `.lt(value)` (less than), `.le(value)`, `.eq(value)`,\n",
    "   `.gt(value)`, and `.ge(value)`.\n",
    "2. `provided_isa`: computes the logical instruction properties from the\n",
    "   physical ones and yields one or more ISAs.\n",
    "\n",
    "The `distance` field uses `metadata={\"domain\": range(3, 22, 2)}` to tell the\n",
    "estimator to automatically search over odd distances from 3 to 21.  The\n",
    "estimator evaluates all distances and keeps only the Pareto-optimal ones."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "3f829830",
   "metadata": {},
   "outputs": [],
   "source": [
    "@dataclass\n",
    "class GenericQEC(ISATransform):\n",
    "    \"\"\"\n",
    "    A configurable QEC code model with tunable threshold, overhead,\n",
    "    and error suppression.\n",
    "\n",
    "    This transform consumes physical H, CNOT, and MEAS_Z instructions and\n",
    "    produces a logical LATTICE_SURGERY instruction.  The resource formulas\n",
    "    are parameterized so you can model different QEC code families.\n",
    "\n",
    "    The logical error rate follows:\n",
    "\n",
    "        error = crossing_prefactor × (p_physical / threshold) ^ ⌊(d+1)/2⌋\n",
    "\n",
    "    Space per data qubit is:\n",
    "\n",
    "        qubits_per_data_qubit × d²\n",
    "    \"\"\"\n",
    "\n",
    "    crossing_prefactor: float = 0.03\n",
    "    error_correction_threshold: float = 0.01\n",
    "    qubits_per_data_qubit: int = 2\n",
    "    syndrome_extraction_depth: int = 4\n",
    "    _: KW_ONLY\n",
    "    distance: int = field(default=3, metadata={\"domain\": range(3, 22, 2)})\n",
    "\n",
    "    @staticmethod\n",
    "    def required_isa() -> ISARequirements:\n",
    "        # required_isa is static, so we cannot reference instance fields here.\n",
    "        # Instead, we use ConstraintBound to express fixed constraints on the\n",
    "        # implementation ISA.  ConstraintBound supports .lt(), .le(), .eq(),\n",
    "        # .gt(), and .ge() comparisons.\n",
    "        return ISARequirements(\n",
    "            constraint(H, error_rate=ConstraintBound.lt(0.01)),\n",
    "            constraint(CNOT, arity=2, error_rate=ConstraintBound.lt(0.01)),\n",
    "            constraint(MEAS_Z, error_rate=ConstraintBound.lt(0.01)),\n",
    "        )\n",
    "\n",
    "    def provided_isa(\n",
    "        self, impl_isa: ISA, ctx: ISAContext\n",
    "    ) -> Generator[ISA, None, None]:\n",
    "        cnot = impl_isa[CNOT]\n",
    "        h = impl_isa[H]\n",
    "        meas_z = impl_isa[MEAS_Z]\n",
    "\n",
    "        # Physical error rate is the worst case across all required gates\n",
    "        physical_error_rate = max(\n",
    "            cnot.expect_error_rate(),\n",
    "            h.expect_error_rate(),\n",
    "            meas_z.expect_error_rate(),\n",
    "        )\n",
    "\n",
    "        d = self.distance\n",
    "\n",
    "        # Space: physical qubits per data qubit, scaled by distance².\n",
    "        #\n",
    "        # Because LATTICE_SURGERY has variable arity (it can operate on any\n",
    "        # number of logical qubits), space and error_rate must be provided as\n",
    "        # functions of arity rather than as fixed numbers.  The estimator\n",
    "        # provides several helpers for this:\n",
    "        #\n",
    "        #   linear_function(slope)          → f(n) = slope × n\n",
    "        #   constant_function(value)        → f(n) = value\n",
    "        #   block_linear_function(k, s, o)  → f(n) = s × ⌈n/k⌉ + o\n",
    "        #   generic_function(callable)      → f(n) = callable(n)\n",
    "        #\n",
    "        # For QEC codes, resources scale linearly with the number of logical\n",
    "        # qubits, so linear_function is the natural choice.\n",
    "        space = linear_function(self.qubits_per_data_qubit * d**2)\n",
    "\n",
    "        # Time: syndrome extraction cycle × code distance\n",
    "        code_cycle_time = (\n",
    "            h.expect_time()\n",
    "            + self.syndrome_extraction_depth * cnot.expect_time()\n",
    "            + meas_z.expect_time()\n",
    "        )\n",
    "        time = code_cycle_time * d\n",
    "\n",
    "        # Error: exponential suppression below threshold\n",
    "        logical_error = self.crossing_prefactor * (\n",
    "            (physical_error_rate / self.error_correction_threshold)\n",
    "            ** ((d + 1) // 2)\n",
    "        )\n",
    "        error = linear_function(logical_error)\n",
    "\n",
    "        yield ctx.make_isa(\n",
    "            ctx.add_instruction(\n",
    "                LATTICE_SURGERY,\n",
    "                encoding=LOGICAL,\n",
    "                arity=None,\n",
    "                space=space,\n",
    "                time=time,\n",
    "                error_rate=error,\n",
    "                # transform=self records that this instruction was produced by\n",
    "                # this QEC transform instance, and source=[...] records which\n",
    "                # physical instructions it was derived from.  Together, they\n",
    "                # form the provenance chain that you can inspect in the\n",
    "                # estimation results (see notebook 2).\n",
    "                transform=self,\n",
    "                source=[cnot, h, meas_z],\n",
    "                # Extra keyword arguments (like distance) are stored as\n",
    "                # properties on the instruction.  This lets you retrieve the\n",
    "                # code distance used for each Pareto-optimal result when\n",
    "                # analysing the estimation output.  Further, the property can be\n",
    "                # required and read by a parent ISATransform to create new\n",
    "                # logical instructions based on this one.\n",
    "                distance=d,\n",
    "            ),\n",
    "        )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2981ec35",
   "metadata": {},
   "source": [
    "### What if we had a better QEC code?\n",
    "\n",
    "Let's compare two configurations:\n",
    "\n",
    "- **Surface-code-like:** threshold 1%, prefactor 0.03, ~2d² physical qubits\n",
    "  per logical qubit (the `GenericQEC` defaults).\n",
    "- **Optimistic code:** threshold 3%, prefactor 0.01, ~d² qubits\n",
    "  per logical qubit (a hypothetical code with better parameters).\n",
    "\n",
    "We use the built-in `GateBased` architecture and `Litinski19Factory` to isolate\n",
    "the effect of the QEC code.  Because `GenericQEC` is a dataclass, we can pass\n",
    "the non-default parameters directly to `.q()`. The estimator fixes those values\n",
    "while still searching over the `distance` hyperparameter."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "175568cb",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 26556 characters]",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "arch = GateBased(error_rate=1e-4, gate_time=100, measurement_time=500)\n",
    "\n",
    "# Surface-code-like parameters (the defaults)\n",
    "surface_like = GenericQEC.q() * Litinski19Factory.q()\n",
    "\n",
    "# Optimistic hypothetical code: higher threshold, lower overhead\n",
    "optimistic = GenericQEC.q(\n",
    "    crossing_prefactor=0.01,\n",
    "    error_correction_threshold=0.03,\n",
    "    qubits_per_data_qubit=1,\n",
    ") * Litinski19Factory.q()\n",
    "\n",
    "results = [\n",
    "    estimate(app, arch, surface_like, max_error=0.01, name=\"Surface-code-like QEC\"),\n",
    "    estimate(app, arch, optimistic, max_error=0.01, name=\"Optimistic QEC\"),\n",
    "    baseline,\n",
    "]\n",
    "\n",
    "plot_estimates(results, figsize=(8, 5), runtime_unit=\"ms\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6229d30a",
   "metadata": {},
   "source": [
    "## Building a custom factory transform\n",
    "\n",
    "A factory transform produces high-fidelity magic states (typically T states) by\n",
    "consuming lower-level instructions.  It models the distillation process: take\n",
    "many noisy copies of a state, run a distillation protocol, and output fewer\n",
    "copies with lower error.\n",
    "\n",
    "The built-in `Litinski19Factory` uses pre-computed protocol data from\n",
    "[arXiv:1905.06903](https://arxiv.org/abs/1905.06903).  Each entry in the table\n",
    "specifies an output error rate, a space cost in physical qubits (which already\n",
    "accounts for the QEC overhead inside the factory), and a time cost in syndrome\n",
    "extraction cycles.\n",
    "\n",
    "We can build a custom factory that starts from this data and explores \"what if\"\n",
    "scenarios by scaling the space and cycle costs:\n",
    "\n",
    "- `space_factor < 1`: What if the factory were more compact?\n",
    "- `cycle_factor < 1`: What if distillation were faster?\n",
    "\n",
    "This lets you study how improvements in factory design propagate through the\n",
    "full resource estimate, without having to derive the protocols from scratch."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "3f4951e1",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Reference data from Litinski (arXiv:1905.06903), Table 1.\n",
    "# Assumes Clifford and T error rates at most 1e-4.\n",
    "# Each tuple: (output_error_rate, space_in_physical_qubits, cycles)\n",
    "_LITINSKI_TABLE = [\n",
    "    (4.4e-8,    810, 18.1),\n",
    "    (1.5e-9,    762, 36.2),\n",
    "    (9.3e-10,  1150, 18.1),\n",
    "    (1.9e-11,  2070, 30.0),\n",
    "    (2.4e-15, 16400, 90.3),\n",
    "    (6.3e-25, 18600, 67.8),\n",
    "]\n",
    "\n",
    "\n",
    "@dataclass\n",
    "class ScaledLitinskiFactory(ISATransform):\n",
    "    \"\"\"\n",
    "    A factory transform based on the Litinski (2019) distillation protocols,\n",
    "    with configurable scaling factors for space and cycle costs.\n",
    "\n",
    "    The base data is taken directly from Table 1 in arXiv:1905.06903.\n",
    "    Each entry specifies a distillation protocol's output error rate, space\n",
    "    footprint (in physical qubits), and time (in syndrome extraction cycles).\n",
    "    The space already includes the QEC overhead within the factory itself.\n",
    "\n",
    "    Scaling factors let you explore hypothetical improvements:\n",
    "\n",
    "    - ``space_factor=0.5`` models a factory that uses half the qubits.\n",
    "    - ``cycle_factor=0.5`` models a factory that runs twice as fast.\n",
    "\n",
    "    Note:\n",
    "        The error rates are kept fixed when scaling.  In practice, changing\n",
    "        space or cycles would affect error rates too, but modelling that\n",
    "        relationship requires protocol-specific analysis.\n",
    "    \"\"\"\n",
    "\n",
    "    space_factor: float = 1.0\n",
    "    cycle_factor: float = 1.0\n",
    "\n",
    "    @staticmethod\n",
    "    def required_isa() -> ISARequirements:\n",
    "        return ISARequirements(\n",
    "            constraint(T, error_rate=ConstraintBound.le(1e-4)),\n",
    "            constraint(H, error_rate=ConstraintBound.le(1e-4)),\n",
    "            constraint(CNOT, arity=2, error_rate=ConstraintBound.le(1e-4)),\n",
    "            constraint(MEAS_Z, error_rate=ConstraintBound.le(1e-4)),\n",
    "        )\n",
    "\n",
    "    def provided_isa(\n",
    "        self, impl_isa: ISA, ctx: ISAContext\n",
    "    ) -> Generator[ISA, None, None]:\n",
    "        cnot = impl_isa[CNOT]\n",
    "        h = impl_isa[H]\n",
    "        meas_z = impl_isa[MEAS_Z]\n",
    "        t = impl_isa[T]\n",
    "\n",
    "        # Syndrome extraction time from the physical ISA\n",
    "        syndrome_extraction_time = (\n",
    "            4 * cnot.expect_time() + h.expect_time() + meas_z.expect_time()\n",
    "        )\n",
    "\n",
    "        for error_rate, base_space, base_cycles in _LITINSKI_TABLE:\n",
    "            # Apply scaling factors to space and time\n",
    "            space = ceil(base_space * self.space_factor)\n",
    "            time = ceil(syndrome_extraction_time * base_cycles * self.cycle_factor)\n",
    "\n",
    "            yield ctx.make_isa(\n",
    "                ctx.add_instruction(\n",
    "                    T,\n",
    "                    encoding=LOGICAL,\n",
    "                    arity=1,\n",
    "                    space=space,\n",
    "                    time=time,\n",
    "                    error_rate=error_rate,\n",
    "                    transform=self,\n",
    "                    source=[cnot, h, meas_z, t],\n",
    "                ),\n",
    "            )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "739582e7",
   "metadata": {},
   "source": [
    "### What if the factory were more compact or faster?\n",
    "\n",
    "Let's compare three scenarios, all using the built-in `SurfaceCode` for QEC:\n",
    "\n",
    "- **Litinski baseline:** the built-in `Litinski19Factory` (our reference).\n",
    "- **Compact factory:** `space_factor=0.5`: the factory uses half the physical\n",
    "  qubits (perhaps due to better routing or a denser layout).\n",
    "- **Fast factory:** `cycle_factor=0.5`: the factory completes in half the\n",
    "  cycles (perhaps due to a faster syndrome extraction scheme)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "1de231fb",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 34124 characters]",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "arch = GateBased(error_rate=1e-4, gate_time=100, measurement_time=500)\n",
    "\n",
    "# Litinski baseline (unscaled, should closely match built-in Litinski19Factory)\n",
    "litinski_baseline = SurfaceCode.q() * ScaledLitinskiFactory.q()\n",
    "\n",
    "# Compact factory: half the physical qubits\n",
    "compact = SurfaceCode.q() * ScaledLitinskiFactory.q(space_factor=0.5)\n",
    "\n",
    "# Fast factory: half the cycles\n",
    "fast = SurfaceCode.q() * ScaledLitinskiFactory.q(cycle_factor=0.5)\n",
    "\n",
    "results = [\n",
    "    estimate(app, arch, litinski_baseline, max_error=0.01, name=\"Litinski baseline\"),\n",
    "    estimate(app, arch, compact, max_error=0.01, name=\"Compact factory (0.5× space)\"),\n",
    "    estimate(app, arch, fast, max_error=0.01, name=\"Fast factory (0.5× cycles)\"),\n",
    "    baseline,\n",
    "]\n",
    "\n",
    "plot_estimates(results, figsize=(8, 5), runtime_unit=\"ms\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c5da8f59",
   "metadata": {},
   "source": [
    "## Composing a full custom stack\n",
    "\n",
    "Now that we have all three custom layers, let's compose them into a complete\n",
    "estimation stack.  The `*` operator creates a Cartesian product over all\n",
    "hyperparameter configurations:\n",
    "\n",
    "- `GenericQEC` searches over code distances 3–21\n",
    "- `ScaledLitinskiFactory` produces multiple factory configurations from\n",
    "  the Litinski protocol table\n",
    "\n",
    "The estimator evaluates all combinations and returns only the Pareto-optimal\n",
    "results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "bba890f4",
   "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>name</th>\n",
       "      <th>qubits</th>\n",
       "      <th>runtime</th>\n",
       "      <th>error</th>\n",
       "      <th>physical_compute_qubits</th>\n",
       "      <th>physical_factory_qubits</th>\n",
       "      <th>physical_memory_qubits</th>\n",
       "      <th>factories</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Full custom stack</td>\n",
       "      <td>3151</td>\n",
       "      <td>0 days 00:00:00.000299</td>\n",
       "      <td>0.006790</td>\n",
       "      <td>3150</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1×T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Full custom stack</td>\n",
       "      <td>9630</td>\n",
       "      <td>0 days 00:00:00.000292500</td>\n",
       "      <td>0.007312</td>\n",
       "      <td>3150</td>\n",
       "      <td>6480</td>\n",
       "      <td>0</td>\n",
       "      <td>8×T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Full custom stack</td>\n",
       "      <td>13284</td>\n",
       "      <td>0 days 00:00:00.000183300</td>\n",
       "      <td>0.009416</td>\n",
       "      <td>1134</td>\n",
       "      <td>12150</td>\n",
       "      <td>0</td>\n",
       "      <td>15×T</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                name  qubits                   runtime     error  \\\n",
       "0  Full custom stack    3151    0 days 00:00:00.000299  0.006790   \n",
       "1  Full custom stack    9630 0 days 00:00:00.000292500  0.007312   \n",
       "2  Full custom stack   13284 0 days 00:00:00.000183300  0.009416   \n",
       "\n",
       "   physical_compute_qubits  physical_factory_qubits  physical_memory_qubits  \\\n",
       "0                     3150                        1                       0   \n",
       "1                     3150                     6480                       0   \n",
       "2                     1134                    12150                       0   \n",
       "\n",
       "  factories  \n",
       "0       1×T  \n",
       "1       8×T  \n",
       "2      15×T  "
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Full custom stack: custom architecture + custom QEC + custom factory\n",
    "custom_arch = ParameterizedArchitecture(\n",
    "    single_qubit_gate_time=100,\n",
    "    two_qubit_gate_time=200,\n",
    "    measurement_time=400,\n",
    "    single_qubit_error_rate=5e-5,\n",
    "    two_qubit_error_rate=1e-4,\n",
    "    measurement_error_rate=1e-4,\n",
    ")\n",
    "\n",
    "custom_results = estimate(\n",
    "    app, custom_arch,\n",
    "    isa_query=GenericQEC.q() * ScaledLitinskiFactory.q(),\n",
    "    max_error=0.01,\n",
    "    name=\"Full custom stack\",\n",
    ")\n",
    "\n",
    "custom_results.add_qubit_partition_column()\n",
    "custom_results.add_factory_summary_column()\n",
    "custom_results.as_frame()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ce953fa",
   "metadata": {},
   "source": [
    "> **Why does the first result show only 1 factory qubit?**\n",
    ">\n",
    "> When the physical T-gate error rate is low enough, the estimator may find that\n",
    "> it can meet the overall error budget *without* distillation, using the raw\n",
    "> physical T gate directly.  In that case the \"factory\" is just a single physical\n",
    "> qubit producing T states one at a time, which is why the qubit count is so\n",
    "> small.  The trade-off is a much longer runtime (thousands of serial T-gate\n",
    "> executions instead of a few batched factory runs).  Subsequent rows in the\n",
    "> table use the `ScaledLitinskiFactory` and show the more typical pattern of\n",
    "> dedicating thousands of physical qubits to a distillation factory in exchange\n",
    "> for faster execution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "c65d0684",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 22988 characters]",
      "text/plain": [
       "<Figure size 800x500 with 1 Axes>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "plot_estimates([custom_results, baseline], figsize=(8, 5), runtime_unit=\"ms\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ae29b99a",
   "metadata": {},
   "source": [
    "## Key takeaways\n",
    "\n",
    "1. **Three extension points.** The estimator's model stack has three layers:\n",
    "   `Architecture`, QEC `ISATransform`, and factory `ISATransform`, each of\n",
    "   which can be replaced independently.\n",
    "\n",
    "2. **Implement two methods.** Every custom QEC or factory model is a dataclass\n",
    "   that implements `required_isa` (what physical instructions it needs) and\n",
    "   `provided_isa` (what logical instructions it provides, with resource\n",
    "   formulas).  A custom architecture implements only `provided_isa`.\n",
    "\n",
    "3. **Automatic search.** Mark hyperparameters as `KW_ONLY` fields with a\n",
    "   `metadata={\"domain\": ...}` and the estimator searches over them\n",
    "   automatically, keeping only the Pareto-optimal configurations.\n",
    "\n",
    "4. **Mix and match.** Custom models compose freely with built-in ones using the\n",
    "   `*` operator.  You can replace one layer while keeping the others at their\n",
    "   defaults, which is useful for isolating the effect of a single modeling choice.\n",
    "\n",
    "5. **Provenance tracking.** Pass `transform=self` and `source=[...]` when\n",
    "   adding instructions to preserve the link between physical and logical\n",
    "   instructions.  This lets you inspect the full derivation chain in the\n",
    "   results (see [Analysing Results](2_analysing_results.ipynb))."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.11.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
