{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "0",
   "metadata": {},
   "source": [
    "# Feynman rules of the Standard Model with `feynlag`\n",
    "\n",
    "This notebook walks through **feynlag**'s pipeline stage by stage: declare\n",
    "symmetries and fields, write the Lagrangian, check it, break the symmetry,\n",
    "diagonalize, and extract Feynman rules — using the electroweak sector of the\n",
    "Standard Model as the worked example. It ends with a guide to extending the\n",
    "model to a BSM scenario.\n",
    "\n",
    "Every cell below is executed live against the actual library — nothing here\n",
    "is hand-typed output. All conventions (metric signature, `P_L`, the sign of\n",
    "`D_μ`, the explicit `1/√2` in VEVs, vertex normalization) are fixed once in\n",
    "[`CONVENTIONS.md`](../CONVENTIONS.md); that file is the source of truth if a\n",
    "sign looks unfamiliar.\n",
    "\n",
    "**What we're building**: an `SU(2)_L × U(1)_Y` Higgs doublet with its\n",
    "potential and covariant kinetic term — enough to extract the Higgs\n",
    "self-couplings, gauge-boson masses, the Weinberg rotation, and every\n",
    "`hVV`/`hhVV`/Goldstone vertex."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1",
   "metadata": {},
   "outputs": [],
   "source": [
    "import sympy as sp\n",
    "from feynlag import (\n",
    "    Bilinear, DiracGamma, Dmu, ExternalParameter, InternalParameter,\n",
    "    Lagrangian, Model, Rotation, SU2, Scalar, U1, WeylFermion,\n",
    "    conjugate_pair, cubic_couplings, dag, diracPL, diracPR,\n",
    "    extract_fermion_vertices, fermion_feynman_rule, fermion_gauge_current,\n",
    "    fermion_mass_matrix, latex_feynman_table, quartic_couplings,\n",
    "    rotation_2x2,\n",
    ")\n",
    "\n",
    "sp.init_printing()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2",
   "metadata": {},
   "source": [
    "## 1. Declare the symmetries and parameters\n",
    "\n",
    "- `SU2`/`U1` are `GaugeGroup`s carrying explicit generator matrices (`σᵃ/2`\n",
    "  for SU(2)); they build their own gauge-boson fields on demand.\n",
    "- **`ExternalParameter`** is a number you feed in from experiment — it needs\n",
    "  a `value` for benchmarking and later UFO export.\n",
    "- **`InternalParameter`** is a symbol whose defining expression is *derived*\n",
    "  later (here, by the tadpole condition in §6). `unit_dim` is the mass\n",
    "  dimension, used by the dimension-≤-4 invariance check.\n",
    "- `positive=True` lets SymPy simplify `sqrt` expressions involving `v`,\n",
    "  `gw`, `g1` without spurious branch-cut assumptions (the \"positive dummy\n",
    "  symbols\" rule in `CONVENTIONS.md`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(ExternalParameter('gw'),\n",
       " ExternalParameter('g1'),\n",
       " ExternalParameter('v'),\n",
       " ExternalParameter('lam'),\n",
       " InternalParameter('mu2'))"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gw = ExternalParameter(\"gw\", 0.6535, positive=True)\n",
    "g1 = ExternalParameter(\"g1\", 0.3580, positive=True)\n",
    "SU2L, U1Y = SU2(\"SU2L\", coupling=gw), U1(\"U1Y\", coupling=g1)\n",
    "\n",
    "v   = ExternalParameter(\"v\", 246.0, positive=True, unit_dim=1)\n",
    "lam = ExternalParameter(\"lam\", 0.129)\n",
    "mu2 = InternalParameter(\"mu2\", unit_dim=2)   # will be *derived*, not set\n",
    "\n",
    "gw, g1, v, lam, mu2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4",
   "metadata": {},
   "source": [
    "## 2. Declare the fields\n",
    "\n",
    "`reps` assigns representations: `2` is the SU(2) doublet, `1/2` the\n",
    "hypercharge. This alone generates the two component symbols `Gp` (complex,\n",
    "the charged Goldstone-to-be) and `H0` (complex, the neutral component) plus a\n",
    "conjugate registry that `dag(H)` uses.\n",
    "\n",
    "Register the vacuum expectation value on the neutral component only — charge\n",
    "conservation forbids a VEV on `Gp`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAO8AAAAVCAYAAACuebnGAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAH8UlEQVR4Ae2b6ZHUOBSAm6kJgIUM2Aw4IljIgCOCGTJYav7xj4IMgAj2yACIYAoygAyAyYD9Po2lUtuW3e6W29A7r0qj48l6t/Qs91z78ePHKofnz5/fonyJY7T/pP206T+i/ynirupLDSytI+iv2ezKLr+eBsZs2Ic/ysVsnPB2Pkb7JuUpuN8pV4HbUo5d9PJK/dD8l3LDsX0BdN1c2zbbF/krOvU04AasLUvQwR/HmTz4kPZN6ldxbErNc9eZf0p50Dx30dTvwL2h3KLvJvCsGd97BW2Vc4+irMJ7ygXjj+xQR/7v01WegKc+ARfloftzADztZLOaUsDLweh2CVmg+Z5ym3JKedO2TR/+mmkzCB31A/WdnodeMmYA6si9AE7DnVFeUAzU5Oi0dTAD2oB4Sb/DGON7BXjwXeETdUdeGWH8O9U3ak/TjYH5o7raeLGRidAq2mzk0VnR8DWLbmdlurD4ErJA8yPs/EGdYihnL8fHtFmne51P2rTNYu+Y6zvxHdqmj2tE6ZtKCp68xQ0gzNjDH/iJKWYvL+Dl08Doxe+BxU1JbG2zTQlMnXdAul0tKItxqG1LkPAxeB/D7OQTkWcMXE9UAzddcvVQ/Ycx09OhOT2PzTL0pFlV3vtAeYQS/hK7/N+tbDYz24eiW9W0iCxNHGpbD5AO5PgjOqa1k4OK59wddHTfY9dO2w5FUlDG/u4ZX2IoBCc8l07W8P4LYyX8Ejyv0dzWZmuLzNP55XWbqWVJWYzHxxkv7WbAHzPq++gkR8V5TC19z/1Ce5MTW2IpLecZU1eV4+XRCcX17HuzbXvOCyJpD8l7F7xyjW1ITFsMijaDb3dsL97U5efcPrQ/M/aMOr7K0K0Kh6DbqJAlZdE/tXEptgLe4NVZU2DR3gTid9+NnsNZDIT8M9MTxnQi0+kPlBe0wy03tSe64zJfFVjbDULw2l0afaDz/yxZQh9/jg3Z7AzZ1K0Z1VtK7gA36JsOVg/eA9Ltam5ZWF8fu0FdynjdZKOv0uxAwBu8LmRaOwXiwkMnWO96MOyOdt4gPWU95XJnkjFP9Tkgbgi9PzaBD+nq9J33XXDxEuEreD+pPZuDwQ3X7LUZPOW6Vda2c5wwZrYzBxyKbtXN3LKcQUM/Kn2WNR6NjRIEvMHrbnxRmlUYDwvjLPlpmqYybhC0P7N4YaXDG6zxOZ3Nz0s56Jgr5lynTOUrX6evHTadjH57TjTa2qbEfNPQr9QxO7hP289ncX57nbn7JZvlun0ME23dtoO5Jp+Holt1Mqss+M3Yxq+dQhwUDBTwRwXk2LCRXwSdnBJTax3/Y2SYOgQkdVAQuLVAoW9AGOi1A5dlV24WbXqOR5Anv/+2aXvqpufA2zaAh3bHuObe6sg3tdmDxs9TZvmQ5zzLcawWHJJul5altDlHWwX8MT0DUUNPAU9O3xstQ7t5cG7mtJ1IWgZp38WQAdQ332e2BniIm0UnJXZR8BpMSEFqh3FlUD9tOQ1wn2mPMzQ7jNmspFsvsaoDOppNt6ztRqQ82k2n1RZmdf7gp7ruZ5ZFPekz96ATv2rQ7YAyDh2QAX/EJBUQgqyzRHkgHvsqthdgzjVldi0YsskdHM94ShsUcf0VYzJaA3QAYYgf8e3gLtFXuTrTZKgg05jN1L0bbAJoaqvOBWMFXqQxp27VsZeLbyUEv7663KUYBAkqyeF6s8jS8OfPH+X/4Qi/yjy0MQW8wauRJ11iQNiF3Tnc/Qy4NWBMxZpqOq8dDCvwBoRz0qbBmG2DNv00rJn3ndqfjO0KYaNhrTWnzhYNRgNfCu5samiqwFJgt+emfiWZxmwmvq3bTpZUiRdlm0W38KePaA9lMUuLGZmfElP6X1EOyMwjC+vqL/50WF31ZZwMJ/BnuyU/dVLAH9P4i1L6bOLEXlB5lN9AGsAG6AUl7hZ+X/RGV6V/o7TBU1fwk0a8WTale0A/rrGiHX+V5eXV5Assn2FNd2zr4MyMKes5dbx88jQSF3hiXFm8jAp42srVB66ZeO2b0DfGujvJ1Kw5aDNoqNfXFGUNmyftKE9ii7GteeHZ2XULjeDA1G6syhwgjmf9reVwDdbbhyzBV6DlXVAnA4qyNLW+GO+MWqjQDfj4jwl+nil9PvEEHfzHhL7Vh8YQQObvUruDjALz3K38r4tSII2uscsE6PqPCmYEaTek7Q/w/TfJFMC0N9YVc3eSieeLNpsq6668TKWXz4f2qG6bOfrnYFa0pBzK1PBZ9BPwbhLK66FnJpT8iX4AxjxIjLf215oO/ujykZDiDkV6M61a5c4xaIgWJV/wFwncho8X1OFktg8vts08UuA6PhF2lcmNopbNduVlouhr0wd1i451ZrOuTfxlSTkUalAW8PpN/JpRsp2vjtq2BAkfghfF+C4Rbo9LT9Qah5a7jwY532TNZr4ftBcDeDDl9IcZviKY5vu+f7ItQzVkYo0qNqvBy7Z68Dnoj+lWf0nvtyVaS8uxoSxuQL7v6kOdAGXcuDAO47u9yyZo448T5tIhfUeKt20Zqk4T4qaKZ81q/ozPn4j1MppRPGVO530tw++lCQ/pBrwCwVoyuYnsarNavGytliHdgjO1VM4xWFwOGRyRxexxSBZfJ0snssuv4cM7r6MChI18r7FTsNB2hzDohEf0O3n6Jer/+xeduJOqdG8U1dEmKV4VhUGrY7MqC18tslcNND5UfBXrw/8Hd3c97MWE15QAAAAASUVORK5CYII=",
      "text/latex": [
       "$\\displaystyle \\left( \\left[ Gp, \\  H_{0}\\right], \\  \\left( v, \\  H_{0 r}, \\  H_{0 i}\\right)\\right)$"
      ],
      "text/plain": [
       "([Gp, H₀], (v, H₀ ᵣ, H₀ ᵢ))"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "H = Scalar(\"H\", reps={SU2L: 2, U1Y: sp.Rational(1, 2)},\n",
    "           component_names=[\"Gp\", \"H0\"])\n",
    "H.expand_vev({H.components[1]: v})   # H0 -> (v + h + i*G0) / sqrt(2)\n",
    "\n",
    "H.components, H.vev_expansions[H.components[1]]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6",
   "metadata": {},
   "source": [
    "## 3. Write the Lagrangian — FeynRules style\n",
    "\n",
    "feynlag doesn't auto-generate the potential; you write it with library\n",
    "building blocks, and the library checks it's legal.\n",
    "\n",
    "- `dag(H)` returns the row-vector Hermitian conjugate; `H.mat` is the column\n",
    "  of components — their product is the standard invariant bilinear.\n",
    "- `Dmu(H)` reads `H.reps` and automatically pulls in **both** `SU2L` and\n",
    "  `U1Y` gauge bosons (created lazily the first time) with the correct\n",
    "  generator matrices — `D_μ = ∂_μ − i g Tᵃ Aᵃ_μ`.\n",
    "- Lagrangian terms are tagged by **sector** (`kinetic`, `potential`,\n",
    "  `yukawa`, `gauge`, `other`); hermiticity is checked per sector, and you can\n",
    "  extract vertices from one sector at a time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAAfCAYAAAAfi0EQAAAACXBIWXMAAA7EAAAOxAGVKw4bAAALvElEQVR4Ae2d63UVNxDHjY8LcEgHTgcEKgjpgEcFkA7I4RN84yQdECogoYMkFQTTQdJBwB2Q/0/WbLRa7d7nDnt9R+fo6jXSjEYzo+fatz5//nwSLjhwjBx4+fLlufr9PPf9IodPlH91jPyIPgcHNuHA2SbAARscuGEc+EkTxQ/WJ8VfK34p/43lRRgcCA60OXDazo7c4MBRcOCpJoz7RU9/UvxCeXeKvIgGB4IDDQ7E5NFgSmQdDQfYdbw/mt5GR4MDe+TArbjz2CM3o6mD5oB2HOw8HiiMY6uDHskg3oMDsfPw4HLgWDwH8lHVAxH67eKJDQKDAwvgwMFemEvZn4p//yj8w4uPwuX2NE24bnn1ywOPN++E77V8dxk+1UfB8dKKXce3ih/dSyv1OXRpSkAWVqbxcrVDwtfUpYM8tlJnWCF+r3At47CwsQ9yHDgg2WBCQOi/n0KX4X40WcrpE4X/TNW7KWXqZ+jSTRnMmfqRdWKgSwc3eagjvM3/U2EcL8wkLDelWckIK+pzhT+3+qT8NMGo7MeinAUJk8mN34Goj6FLxcBHdJwDkpWBLh3incdv6iLv8cMFByY5IIH/RQA/ZCPZguWbDp7qEprn+e6NnzgyM0KXMiMimOZAS5e6Ow8VppcmaoLVGGe/H6ab8y8VTdB2V+HkUYRRJjhWVsyYBm9G4XeV/ZLbw7iUK0+rvnbohWdtgnYAVF+eqfo9eY4zcNwpXSn/IQmFxk+MLvxN5Qo3+jLbC4/oYqHxRj7Rr7BzouGrLnFkEfU9dGnmMfeScS88YldPl3rHViICw8DZ1iIva0UXxGPIVhp7wWAEn8u/kmeisInjRPF0zqt8DCBfGbNC3cqprguerYjboZL6xaXcB4XN40Hlf1L5R4U7PWudG4/aZ4KD1q8U72RA6aN24kXokpMEzC3j1o258aj9ni7Vx1YYisXtOIw5Cpncfi/Szag6CQxn1+ygfpbvGQ2l3+WKrL62fq3lhSfTOhmIlvvy8Gdnp3bsC+smb1QO3xCkZvm6BHjgEQ7GHpl+tC5dRwIXujQy0JKZ0KUGb2pdOq1gWIm/rfIWkSwMzeQXwYJj4qAfTBxTL2Y472UXMwUz2ncvPKMEDAsw5vh9uMe5kbGJGv7ixsqvS1f/euFBZgbHVqvJu5kQkl1bHIQutYc4dKnNF3I7XeomDwkUq0n8TqtJWp/JYbD4roOVZNOpjHsb4LjHGIXLlT8q/DXHNwq88GxE1H6B0+Sgfo7JghnisfJ1qfHC87cIursuUUcAF7rkN8heMu6Fp9Ols4KHhrx5bCVDwsSCcWal/q/8O3n+lEP3DFJx2mBVc0/+ifxteY6PgOfSmotqjpGAsVUnce4dVhki2hydEFQf+rh/YIJZ5w6DfnDum1ymCfqNdtoj/bU88XQh7IVH+L6kY0ymxgNDPDmRr0m8Fx5kmie7+FEZWpPmWcGyfF0q7F3mK510S2GnbzsQErq0A/M2rOol4154Ol06KxiBcW8ajCy4HPN0R0HKS8dDCsvLaD7c4408ZW/k/yKt8EQhBv2TQpQXw2P5XF7Tdk9ZlK7duTLYLYw5JilcNyFcJ9u/mY5yonwMTfLQ8qf8K8WToipk0iQfHnnhESp/p75ipHD8dVn63HKMxVa7NmvMC0/GZ3LDIqAccyNnSSHyZfSWdCGDk8dMJfCKeOjSCgbto9hLxr3wZJ6YbF6UkwdG41XNNBGGoGFEMKzl/QDxbvWZO2Bn4CgpZeUq6bbycKz+ykmKfHCscsCV+Gt4M3pl2zVMMy16mLX/yoVGOzsrc2zV2NXgvPBcY/P/ZYLEPRRfBoZWefCBCd/GGtjkVIaBw7HT/FrptEBIOcMfLzxgvsroTQaH1CwnB/nqybD4iH4goz39zPm8KMQht7h1nkyHLl3zau5fLxn3wgO/Ol06IyUhRDAR0NJgUoTDIGDw66MghLyEfy+YK3naQZBtha5ocuDA1e2s+8KLdm3WSw1VP0l5hH9g8IBTPkavflYKvRg4JjqrN1BSlYObNqxvxA2eos4pfy941M5V12gRUT47K5vAipJ0RHii8prvwPDk9mEJPBFPbQu+2T/VM0GtDRyvd/5VPdut3VecY0qDr1F64QHvxxr5Lmn1CTlgd5rkYs22mpNxo25L/my8ezxX3W3/mRV0T/EkdEkM0jiHLg0FtJObs1yGcGJIWyv7RyrrCW1WHgSsW30qz4ydGYVeHcFiRDBiBqdkcrRfTyi5qBdQ73Yvp5+gUyhF0wmvGTWML4aOS/WE12hSOKak0A5/8C54mp1QpvC3BJp8dgMcNZW7vbFmpvIxXvXYlfDwqDWOLDK+M0DR8Yc8kwc0teTKCw8kmdx0gm90bhOqP8hi8/uXbdqzOmp3Sv66Xb7BK3yqOr/J23gxBuTdkf9QwNXR0CVxRDwKXdpMZ5GjTpdOs1RhHJMAiqGdAc5x0t0kkeEx+DDfhDZnp+CxfluCi2L04FWfPNpP9xQZn5JbuYRTbaRV00QLtqpqTVjwoaWk0PlrbtMLz0QX5inK40Hj9XgnhCrH4OPqcYSnjGM9SVwpz+ooeu288Bg+hdCGg54luyn5S7t88c76Qj8wfu9n6JCXjHvhmYFF0016ybgXnqK3Jn9XNnlgHM1g2BlqAT8wCgi5DTwvrkqD3ZokQAiM4bC2EX5WsRhsytOkZIVViGEywquilLTzdVbgTZdxDOgrgAdlqsMu5Ure2rdwbjwFWW5RxhXXmxyus9Iv/MHV4zg2Lh8FayuVVDH/eOExnEYD9CzZwd+kV0ZkltlSdzr9VFn5WIUq6FN5BGvN1GHoUs2R/ae9ZNwLj3Go06VTy1GI0GEEuOxMTmmMJoYE4U1OeRhNVpO24uGFVVpxKgSONt7Kl25spwKsGaruGKmsWMS50L5bpHvRTAPn+pwDY/B7TnnQzLYeWmvjd6JyaAGm7CtxJovvVA4vgKP+7HjA9QVcmhDVx54BK+hIgqpyG7OiqBlF0OBr7bzwGF7GMR07WsbSQvHU5A8ZLJ0tVmyBxeONgVN96sHXdY7TQpcGHNx7hpeMe+ExBnW6dJZzWLFgENNTW4PKIflvJJzPFPLNAxMDRgQjTV46clKIo+HWygfFaB0TYcy5SKzbUfbAYdDOBYtPhryGUP47eZ78QhsTBHBpYlP4t/LABY2tFaitqnlVBj04+ttNjilHPyp3wWP45gzVF8bmjTwhvKF/vK7jmXV5T0RZ4pHy4W36ZkchrjkeyqdNW1i44IGYhsOgrjvhNaq7ZN3NWLgrSvwnrTiLKiYMdLR7kECZOZUzNugST+nHxsLACUOXSm7sKS7eu8i4F54RtvyvS/wP80PxL168+CR/fw561e5r+cs52i7bnAuP2n0g/6zE5RkXbsbmTolT6c/yF2XernG1tzEe1bmUf7or7jnri76t5E/1LqhrtOX0Sp4LLnRpxPaJN6FL47zpdOl0ZHZZajaX1nbGt28aWVV7rE7nwsOKc51V5775Zu29UsR2bydaHRFnh2Y7P4PbNdwIT16lcaRjDx52xT9X/Y3lQn1jx8HOn7+EzQsr+skxV2tnreyeC13qsaOXCF3qseM6UevSWQNmyVkoCu/r7Rx4L7RmpqCInAXP5ubEo7Y9Jr5R3gg/f3aG40KOT7g345uaJ6MVtizYAs8joWIS+5IT62RvRdu5AJC/wV3cZMXrf2BF3csSTu1xxLXKhS6NcEj8C11q86anS73/59GGX1auBpbz4LcKyw8UtyZS7XDh9FyeVRtnwaziWvczKtreeeHZnsKbWVN8575gcG+1pN6KRnYdTByu/3NEeEOXliQIC6el1qVDnDxYacX/MF+4oC2BPAk7r+54YJEu/pdAU4sG0YdM8x8yXVe8GW/oUmtQIq/HAcnKQJcObvKgR+oIu4V7Cvd6fNXjViQOmgOSjXQfoHCuO7KD5o8RH7pknIhwjANjunQ6VmHJ+eoMR1Y8vWW7Hy440OIACwuemYeb4EDo0gRzosg40NSl/wALtecFF7kDTwAAAABJRU5ErkJggg==",
      "text/latex": [
       "$\\displaystyle lam \\left(Gp \\overline{Gp} + H_{0} \\overline{H_{0}}\\right)^{2} - \\mu_{2} \\left(Gp \\overline{Gp} + H_{0} \\overline{H_{0}}\\right)$"
      ],
      "text/plain": [
       "                   2                     \n",
       "    ⎛   __      __⎞       ⎛   __      __⎞\n",
       "lam⋅⎝Gp⋅Gp + H₀⋅H₀⎠  - μ₂⋅⎝Gp⋅Gp + H₀⋅H₀⎠"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "HdH = (dag(H) * H.mat)[0]          # H^dagger H, an SU(2)xU(1) invariant scalar\n",
    "V   = -mu2.s * HdH + lam.s * HdH**2\n",
    "\n",
    "DH = Dmu(H)                         # covariant derivative\n",
    "L = Lagrangian()\n",
    "L.add((dag(DH) * DH)[0], sector=\"kinetic\")\n",
    "L.add(-V, sector=\"potential\")       # L contains -V, per CONVENTIONS.md\n",
    "\n",
    "V"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8",
   "metadata": {},
   "source": [
    "## 4. Assemble the Model and check invariance\n",
    "\n",
    "`SU2L.bosons(\"W\")` / `U1Y.bosons(\"B\")` materialize the gauge-boson `Field`s\n",
    "(three real components `W_1, W_2, W_3` and one `B`) — the same objects\n",
    "`Dmu(H)` used internally, so they must be listed explicitly as model fields.\n",
    "\n",
    "`check_invariance()` is a real check: for every term it verifies (a) gauge\n",
    "invariance under each declared group by an infinitesimal variation\n",
    "`δφ = iαᵀφ`, requiring every generator's coefficient to vanish — including\n",
    "through the `Dmu`-built kinetic term, since `∂_μ` commutes with a global\n",
    "transformation; (b) hermiticity per sector; (c) mass dimension ≤ 4."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "InvarianceReport(8 checks, OK)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model = Model(\"SM-EW\", gauge_groups=[SU2L, U1Y],\n",
    "              fields=[H, SU2L.bosons(\"W\"), U1Y.bosons(\"B\")],\n",
    "              parameters=[gw, g1, v, lam, mu2], lagrangian=L)\n",
    "\n",
    "report = model.check_invariance()\n",
    "report"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "10",
   "metadata": {},
   "source": [
    "## 5. Break the symmetry: tadpoles and masses\n",
    "\n",
    "`solve_tadpoles` evaluates `V` on the vacuum, differentiates with respect to\n",
    "each registered VEV (`∂V/∂v = 0`), solves for `mu2`, and **registers the\n",
    "solution on the `InternalParameter`** — `mu2.expr` becomes `lam*v**2`. Every\n",
    "later stage (mass matrices, the physical Lagrangian, UFO export)\n",
    "automatically substitutes it.\n",
    "\n",
    "Scalar masses come from the Hessian of `V` at the vacuum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "11",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAHMAAAAaCAYAAACEuGN0AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAFlElEQVRoBe2a7VEbMRCGjYcCgHRAOgBSQaADCBUAHcDwj38Z6IBQQQY6gA4CdAAdwLgD8j6yVrPSnX1nG4ghtzOydKvVanff1ccdLLy8vPQ6+rcRODk5WZIFx9GK1VjviT+YxLLFOmEp2Rb/QOVc7as6mY73qhE4VZyJdyC1z9W4U/k65Ax/m3Dpe2HaGoASFO+o3QFJUN6e9hXrTTfNqdqr4q05Xi/isSPertoPKqzoRNnKVKcpWU8SXeM9IsDiuW0zkTBi62WhPahmaz6ycRmYYrJfP1pnV79PBATMr2ImwH0U/77g+0dwylZuCSbCz37ER2rLeXYWznuScr0hGHPpmmwGIHxo2h1Zodk2Wzkz59LDlkYpEGw5ANr7oECShNhPIgLWRFS3MidSMIfCZPS47WkOTQ7JB5BHAnELA1XzTN362PuMYG4qBlztPwxF4LAZMO0c5NxMl5s2znwqMGNQyOibNs7PkQyvg5x/1InkD4COoxXfWYKZHahe0Ldj0O5ULxd8VsWa+GeeP0tburiC36hucoxpmL8n2dptVnw7k9i6nlR4j94WP9mrdvBB/G8qeyoEjLmRZwu8Rl6FFbSrAtHmxT8lkdrEcl/liwrvhOnGqjY+sQrDe7zqLI7qm4r6NkoKmRxHLo03psa5ulsvh3f21WKMjsauaBMAbDQKDwUIdgqoHyNd+EbmE0QKALK1AYJP4q3YB+9CBbBNnhd25AFpJfLZCtFTxu046vnDGBVPJIglgudP0v4tYRYO8Qm0yK8YZBbGYHTKIPpGEIHJgqZxOI+en35M5B9Hnk3c6rujxg5Ulqm9zjFt7MrmR1bjsc388xcK2rzPBf2qGX+tAmErfWnV6hkQoCXxvf/wmSOQ+ogDIEIkmJ8THiuelT81aY4rFRYVO8WBys2ifsgysouvCrXbU82MFdAkQyAg7yTPrb47IlhHsqkVkJLDJgJa9wmSlQEAZaJis5e/ZT4V9AAmwfLEHFCpp7xB+xf+H5IvE6wENyid9Ad/VIg3gJ73YegBoy/VPmxSKJlRoIUMVH8Z/H03BvUEtvLdsWneFv3YBRB1gSKgWZJJzgCzldgTz2wPPuo5GyMd+Hjv5PQYKNNv/ap5+WeeEvxN9fkkGmqZ8Fc6WIjYz8o86zNeDYwmu+w8gD2KRoGWslz6cMCIRLm1hzessSsE388f29iTQIs2AEBP/SVgsDnP6nYpfMzkNR4e+jk30ed9HxUrLkUzkeYhUZiTMz7YFMBEqxgDVTA55McRxmeOaixbEsUCZmcketkK0G0EuH4bMv6sNXZV5ndKyxUbVhn9so9LDvYb1YEGSN5Hk8UfVis+0R+SJHbyXMbKQIgiU1ch4ZjXNPStEWsf9KIrZR3nhp0dJsOZC5lDXL0rpIkZhzPsAo0kec65F5Xs/WvMQOYn6LxGBNIzPpGkBDaQeNiALbZjkN0hKKqRQwe3RU+jVjKytlrZ7vyWCpB+XtocMQkAP8EU7Wc/ZtE/tGhvRBne+7gdhhub2jgBgGTpk9r+BhiGiIcjnJetvztqjJ2BFvSga8QPc7OrAIwll4nCvxD/UDVbHECxMjlW4IUtUjWEnXU7B6B5oJCF8InLY6mnJx5vB/yBn1iFXUPtSmxQ8hq04P9tJE7KasDRComP0xuqW60sUyB5AoRjBLwXn6lfK0Ntqv+mVuxIkAyrSVcmZ8lEt7AIHEkAmLY9A2q5esTqaJYItAZTQLDNsMLsktF2Xs47xmbnnvSFVdpWSSfXHIE6MFdGDLPzss35lVQItFf57pgUdg2LAAsko372NLzQsPrqCBC5XAzqOjveu0cAnLI7R3YBwhyBxXbIba7pfRPxjv5BBIQNlx/A/K52WlwVMLFNAryHcaZ1/zdLQOaEhAsXUC6OfHqtvCb9BbqMKAqV14ZmAAAAAElFTkSuQmCC",
      "text/latex": [
       "$\\displaystyle \\left\\{ \\mu_{2} : lam v^{2}\\right\\}$"
      ],
      "text/plain": [
       "⎧         2⎫\n",
       "⎨μ₂: lam⋅v ⎬\n",
       "⎩          ⎭"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "tadpole_solution = model.solve_tadpoles([mu2])\n",
    "tadpole_solution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "12",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAEAAAAATCAYAAAAgcwuHAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADjklEQVRYCeWWsVLUMRCHD7wHQHwDaK1QeoqjtAN8AsAn0LnuOgc7S6Czc6S0g4LGSuQN4A1krrQ7vy+TxP2HA688YGdy2fyy2exvs//kFiaTSe+xy2g0WoLjMPNcyf0u+Lj/2MlnfgeQ3S9c0Q/Rf9FWFwv4yPs9SA8CxwP0FbC1p5IAT/8iJKCqC0/hDqhss8LJWwFb9KvpDkDxYviQ51/T3zgGv8xY6spCBtq/auej7bzqxLxGbFvGb4zPaJI5ZGJnY2PjO+3o/Pz8JdgX+h+Mr9GToJ+B/WHwBvt3GX4wHTHL9TNtgD42cCvAcqg3pCCTnv4e6jfac7EgZq5TGWFublX4SF5emwaZx71FdG/HKwDfyihnDJaKYZjQ/msYz72aOfj0WelrNnQ/+RsrQKKCqSTQW6mJyY7MpGsekvjmy8O+Cnz273wFmLzC0rdyoaxA97MwixUrc/bgJsdPynvjN+2E5m37iT4JuhXkCazTdmnLND9B7S3PU+1p2rylKer+mamJR5eQ8bygWcFH9EnQjd1yd/97xU/glrDQDSVTXoZiY4A1iALas0ZiZtiNbZK27Ay8VhHjzTwndkwzQcV+O9tLbDnjxqAf76Mow+znJ6BJj2JSS/IifkufmgCs3OwkbxAXSfI0AurYScY1EqmvBmP1a7AxvXZxvQn2jqnVwdjAFfGYaPGaROY8IIkrHkrcU8zKajHxW9JvEZybbYP2NKrkTQ1iWll5AgZdyzAvlHC0v8BmTNOPCei8PowlprR+2pfH+MpLtIP9x7Tq389M5DXvVABOS+mZ1VYkY/DTnBtEPLEedoVkrRiwcXaqL206a4Dc9zLYaaZ0/Jd5ev/QuE+bsAFzMfH6mCo1AdnZKn09eXQvQU9KMbgUMJibJsm640o0Txl0j/mWpLDfZzlBx0VMTMee9WL6tzL1V/dmaEz1E3M+ixfjTJISgFNLb52+vfTM8E32ZCCF5DBjsWsrI52mBvj1oiuJFJpGVGLalD20U/xMrAqJOp8S6wTiuJNI9wJLydLgf9LPTr3AztDbhQOweEkZhIH6ZCVh7GfhqVWCjA3CpJbT9OZPJUmvnT7aP1N3VYy2xc8+6+MhSd5kJsm+rdqZyt9FfZoZNyi//1Zidj0JPw/JxCBcI34M/p7e8pOcFeATKBYT617xEmOYRKLtt+yEF+z2FD89MF8d/5d4gKly0OOBuf5e+Qv99WeimbUCpwAAAABJRU5ErkJggg==",
      "text/latex": [
       "$\\displaystyle 2 lam v^{2}$"
      ],
      "text/plain": [
       "       2\n",
       "2⋅lam⋅v "
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "h, G0 = sp.Symbol(\"H0_r\", real=True), sp.Symbol(\"H0_i\", real=True)\n",
    "m_h_squared = model.mass_matrix([h])[0, 0]\n",
    "m_h_squared"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "13",
   "metadata": {},
   "source": [
    "## 6. Gauge-boson masses and the Weinberg rotation\n",
    "\n",
    "Gauge-boson masses come from the *kinetic* sector evaluated on the vacuum\n",
    "(scalars are constant there, so all `∂_μ` terms vanish and only the\n",
    "`g²v²·(gauge field)²` piece survives). `W1, W2` decouple from `W3, B` at this\n",
    "order; the `(W3, B)` block is the familiar singular matrix whose non-zero\n",
    "eigenvalue is `m_Z²` and whose zero eigenvalue is the photon."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUsAAAB4CAYAAAB7LtApAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAXNUlEQVR4Ae2dXY7dthXHrw0/F84E6AKcHdjxCjzeQRqvIMkOEvgpeQviHdQ7COwduF2B0+zAzXuBjoNuYPr/cUSZo6F0qXspUaQOAV5+iDwk/4c89/BD1L3r6+tDSfPTTz9dqvyHso9kv5B9o7h/yE02OWgkF2YJDQFDYJcIPNhAqxGOn1EPuQjNj7L3CM8wOWjMKM6SGgKGwN4QuL+BBj8J6oB2+WcQTvXmoJFalqUzBAyBHSJQXFhKm/x3gPsv8v8tCCd5c9BIKsgSGQKGwG4RuPfjjz8+Vuv/NYLAWwmi2cJrhNZktMr5uxL8Xe7vkwknHuagMUHeHhkChkDDCEh+fFDzmN3eMXp2L9QsXynFdwOLAFvcqCJolD/I/V32+1MKzEHjlHItjyFgCDSDAHJoKAPf+taFGzxodeGU2Ke54yrdt4pkMwbDGuNDxSFsD3LZ3UZbfSq/00rlkvYXud91/j8UfiL/v2XRakn/vfxyDtTB0SKguBg9GsWzHzp3kgZpzBgChoAhMIWA5Mnr4XPFEfUVPw/4mWOUGc3vC7lI4IPcd3Lcpoz8CMXHcl/JXhOW5dnXsk69JSyLQCSMsAw3ZxT1yejZGL33SvXCp5yi4dOYawgYAobAOQjMEpYSSgg4tDp31Kcr+ELum8D/WumQxAhCJ0Tlfy6LUPXmV3lStFho36EnuqylPvXEWnXVRqdBq33/leUMKtp5Cm6tQnJWuwzPs+A7O3Pt+IdrlilgMCVmXdELQfIwhXaHyBnI3TO0vn6uL7/LJ7c3pO0DI54xeoqnzFD4jlCoN1ptZGnhvVzWclmWYMnhnfxOQ6+3ZWVqbniWwd2X2gL+c4Ul2s1vAQAILT+t9tG4vTATSEylmY6Hb+WEwjbMN+ZH2IbC8XJAbyxflfFqm1sTltv/4cgPZoRX2XSrEriRShueI8CsFN0K/nOFJWuFFwHGL+UPhaB/xMC+6gJf+kjcDrhYnjDZ0A8tR0/50axOPl40JLzRMBtjsTaCP38U/AGZSUfA8EzHaomUTeA/d82StcLnsv54D4KL9cehYcr4UukY3AhOdsFZf+McE1P1o1NwpQsNm0kvlM8JXrl3dq3CxA340aRjbfS48bzXOhto79JNMDyXRniafhP4P5hu492nElRuF5wn8kc1S8WjOc7VHu8W1sXkpjda0AYeqK0pWmOo3W+g1tutguFZljct4Z88DVejORKEZuiM/GiXaImx6WKXypwTEPCCEI18zKQI1LG8e4s3PMtyvBn852iWDF6EI5sPbPQwYJ/Jmlkfgc/XL7LpEg3PsuytAv85wtJrlaw/ePOthKf3b9ZVHede+VayLW4ja6QC/l+ac5dm0hAwPNNwWipVM/gnC8vKBM5SjF+crnDmKBblxKbaPs5v9Cxen9oLMDzLcrAl/JPXLMtCvrvS2RzjpMHQeM0y2+bZsIBGw4ZnWcY2gf9mhKX+gXhjxcwNArw+eut8agcM79EP36DqHpkzgYDhOQHOCo+awH8TwlKCkp31xyswrYoihAdnLK/k8o69M/IzBedCkm9uYuw3FQHDMxWpZdK1gn/ymuUyMLqzmpeibWtwdwFGi+TiDC4MYUMH95nCdlRLQJxgDM8TQMuYpXr8iwpLDXy0JTTK2NsqGflUHylh86dq3b8AUF8LtlVjw7MsP1rAv/Q0nKNH/UW/ZdlppRsChoAhMI5AMWEpIcl6nL3fPM4be2IIGAIbQqCIsJSgfCQMLuTaWuWGOoNVxRAwBMYRKLVmyVoln6bwN4ETPnRhLrw1jXOcZ/bEEDAECiBQRFhKGLKj2+/qKoymyTvnP8vPxoYZQ8AQMAQ2hUCRaXiIgIQjQtJrmByV6c8WhunMbwgYAoZASQSKaJZhgyUcOTZkR4dCUAK/8LlUkGUKtG9ue3qjuFmvO+agoXLNGAK7RqC4sNw1+mmNRzi6r2nKRWh+lJ17i1IOGmm1tVSGQKMIFJ+GN4przmbx5oM3aJenrOnmoOHrYK4hsEsETFhunO3SJsPjVazt8vGnWSYHjVkFWmJDoEEEbBpeCVMl8PgELt8Q708RzK16Dhpzy7T0hkArCJhmWQEnJeTQKJ2glN9/WXNWzXPQmFWgJTYEGkPAhOXGGSohxz2fCMiP8l/LnX25Rg4aG4fJqmcILI6ATcMXgFjCibOj7Fxj2JB5qDh3YYhcjgJx09JT+d36o1zScsaU76vj/0P2ifx8IC7cnFH0baPnMXru3Kqe8f32wzEatylayBAwBGIImGYZQ+WMOAkmtEAE3Sus/AhE7qI8KIwg5JPCxH/VhXnEpb6P8CgO4cqmjgsTN2a6/DF675XnaP4xuhZvCBgCdxEoqllqsPs3d7jclgPXaFfh7u/dGm84RnVHQNEmdy6yqyrfzeFafQz+10rHW0pojQhGzHPZd8538/OrnBQcovRE962sE9ABzd17hUlT/a02htaO//1SgAs41uK4NIONCzQtpozv5K9ZI2JKPPxGDlNu98aN2uYF5AvFhZeFuHyK6w1p+8CIZ4ye4ikzFL4jFPYTLUxa7G/VMLAF/BGWf+kQ9+7iDBBwbk1Pbi8w5EfLIswRmVoN2vFvvvJqE0KLT9sOBV8vzPSMqTlrmuErjF7j9KSOuQjbUDheDugdy9/0c2HRan+rgm+V4/9XD3IpzZJ1vNh5QdbaGOgIkBoN9b8IKv5S/lAI+kcIw6su8KWPxO06VixPmGzoh5ajp/xo5jFsh3n2FG61v9XCwybwZ83yfx3i3l2DAWhCscszvAbG817rXKNCOcqQoGKt8LksmzwYBBfrj0PDksNLpUO4IjjZBWc97YMsU3WPg4JJhuNEL5TPCV65MWyTCDWaqMn+VhGvasb/Px7nB96zlquBnKI1htrZWlXLUo7a15+DlD+qWSoezXGu9jhav9z0Rguq8IGwabq/bZ0lLeFfYhruBSEa1ZhJ6eBjeYvEq1NwhAfN0Bn50S7REm1K3GFSyGmyvxXC8pRim8F/dc0yEe3PE9NtKRnCH+HIZgIbPQj8Z7Jmto9Ajf1t+6im17AK/EsIS7cRMYKj/xfi3GVtxmuVrM94w6d+vX+zruo4937MzbYlUrFW+1ukqZuMagb/1YWlBiZHaeBqbKrt4+ZucBTvJY0LnOL4nlqBVvvbqXisna8l/EusWcIvNjfYKR4ar1lm2/wYFmDhXSJg/a0s25vAv5Sw5PW/W+cLO15yacTwDZiybN5I6fqH5g0UM6chYP3tNNxy5WoC/yLCUgOfc4BXcnlH2hn5mYJzocQ3NzH26xEQNuysP/Zhc+chYP1tHl65U7eC/+prlgEj0CK5OIMLH9jQwX2msB21ERDeCI9L+atbw/X135Br/a0sM6rHv5iwlBD4U7zrD3CX5eM2SxdGaNtolPZGzpkssv52JoBnZm8Bf6bhaC0Irakt/jOhsuwnIsDRI25kMmMIGAJlEGBzyil1CMtHstz043ei5TVTGgEJSdZzq3s/vjRuVr4hkBkBlsHcTWhFNngyN6Y5chKU/IFdyLW1yua4aw2qFYFia5a1ArZSvVmr/ELC0t/sTfjQhbkw2TTOlRhhxRgCHgETlh6JDbkShpwI6E8FKIymyTvnP8vPxpgZQ8AQWBkBm4avDPjc4iQcEZJew+SoVX82dS4tS28IGAKnI2Ca5enYrZJTwpFjQ3Z0KBPawpMFe5Y10Na5HeqN4ma9XpuDhso1sxICufhlwnIlhlkxm0EA4ei+vikXoflRdu6tSzlobAaQHVQkC79sGr6DnmJNvIUAb5J4g3Z5yhpwDhq+DuYuj0AWfpmwXJ5RVsKGEJA2GR7HYi2Yj2nNMjlozCrQEp+FQC5+2TT8LDZY5loR0ADioDHfrO9PHcxtSw4ac8u09KcjcC6/TLM8HXvLWSkCGjRolE5Qyu+/xDmrNTlozCrQEp+FQA5+oVn+pauFd8+qlGU2BHIjoI7O8Sk2YzCsMT5UnHtnXi6721w28lR+N6WWS1qOWfGJYfx/yD6Rn28kcS8o6b+XX467G6F//15xMXru6Jae8Qnjg9xJGqRp2aj9o/xIafdUfj2L4T/Fz1j6nPz6q2/TA3n+1wW865+ZawgUR0CDB82Pt5ncZQZy3ymMwDzIzyDiq5qvZK8Jy/Lsa1k2bw6EZVmnJIywDBf7SdIbPRujx/fdX/iEUzR8mlZdtX2UHyltnso/gX+Un0p/pTJj/M/Jr/674TYNT+GwpSmCgAYDAs5NmYMKXMiPwMTgf610HNRHEDohKv9zWZ9G3sOvsuHGDnExE6UnurxempI/RrOZuAR+HEgj+06WP55bhmeKyMnPVfmFZlnMCDynLqsCXP7LAWGmTrvvlJ4hhs+BKdbwMyOPFecOkfu+IhetL3xfnnzuphi5zvi0PhxzfZohPYUpMxS+sex7iDvGD3BCIGJj5lh+N/aH+IvQJD+H6Zfi1/1Yi9aIU4NY93kvl4V21oxYD+IfaQzoNaq1mTIMH8cK/kB/80wRJgxGP6320bi9MFMaNBqm4+FbOV7jDPNM+RmcoXC8HNCbytvys0l+CCP+2MI/rSEWk/mDxHP5uQq/ighLAeoWiENg5adDA/QtjSAAcDdew6dnNWtPF33ocHgpfygE/SP6DutXmC9vnJvfDstYnjDZ0A8tR0/5+fP+fZhgp+FUfozBk5p/Lj9X4deDsVYtHM+uZawDAia7lH6hfuFqbJa84SPWqB+8lX0uy6YCBsHF+uPQMCt5qXT0HwYau+As8XyQZS1z7tIOm0kvlM8JXrn2br4AEQ6p/FDqu2ZG/rn8XIVfpYQlanOsA/pOzfMpdf4uJ9qKMXw6fmqAuV1wgvJHNUvFoznO1R67Eu46uendLaHeGGFzlB9TrUvJPxf/uemn6jf1bPVpuBr2cKpC3bNw6pWQvJ0khs8NL4UDR0LQDJ2RH+0SLTE2I+lSmbMUAin86NLAJ8Yvmn5/nWBK/qXqnotuCc3SC0KmS2MmRaCO5a093vC54SD9A+HI+jYbA/SJZ7JmyiBwlB/iFX9k2P6Qf1DVo/mDtJv0lhCWKUB8npJox2max0cDjyUZzkua2QAC5/Lj3PwbgOBQQlheTTTca1Wcu9yrMXzEeQ2u61o7gOo+937MxZqquqCR/1MWN9X8TfnQEHtTMz/6Rox4Uvm1urBUxTgnR7VjzPNxfqNnpHntRhs+N7xN7cDt9oQ8LaM/idLoK56ppRg/DocHAoF/kLX/Cdm55BjI0HjNMtvO5rCASsKGTyWMsmq2jYDkY/9Zl/sKsOvIJQQx4bUUEm9E+Nbh4a4g/gGHr7ctVYct0zV8FuaO+jtvkJnZKAJb4Y/q8a2sWxJa/egQvFHhSOsrueHRAqbg3C7yjeyujeGzLPuFL8dbHi9bilE/FYGt8mf1NcsAQLRILs54KpcNHdxnCt9aWFbcXo3hswDn1b8uRXa3a+ILQJqV5Jb5U0xYChQWnvu3AbIi3gAxwyc/E4Upsxc0ytjbY/kLNIqzENg6f4pMw2chaIkNgXwIsP4UOzCdrwSjdA4Cm+aPCctzWGt5q0FAQpL18T3fN7BpXtXAHxOWm+5CVrkcCGggPhKdC7m2VpkD0Mw0auFPsTXLzHgbOUNgCgHWKvmOj7+Zn/ChC3MBtWmcU+gt/6wK/piwXL4jWAmFEZAw5IRFf8pCYTRNLuj4Wf4/C1dv98XXwh+bhu++q+4LAA1MhKTXMDm61p/13RcS22ztlvljmuU2+4zVaiEENBj719cWKsLICgHhfCmH6TVaPFfsvVHc0deYQ/54GnJ5iSCZhtIuYkxYLgKrETUEPiHgB71iZgmOTxSq9CEcP6PmchGaH2Xn3kGRgwZVyGJMWGaB0YgYApMIbGrQT9Y038PwpiP+JE5ZG85BI1uLbM0yG5RGyBAYRWBTg360lhkfSJsMj2mxRsxH+GaZHDRmFXgksWmWRwCyx4bAuQhsbdCf2545+dV2Pm39g9z+NMKc/KTNQWNumbH0plnGULE4Q2ABBIJBf3SjY4HiVyep9qJROkEpP5s0s00OGrMLHclgmuUIMBZtCOREIBj0fCnge9nq31FXGziGxeYNhjXJh75dcrkvlEtLaKscd9NT32bFXSqO50/ld1N0ue5Il1y+G36QO0mDNGsaE5Zrom1lZUFAg2h0kKYUMJVfz2KDGIHAmczvZPH/IftEfr4+GUu/6UGfgtGxNGq3O84j190cJved8iAwnVE4XKf10f4ZGHLp+CtZLh5HyJL3vewLl0g/ihul4dOs6do0fE20rayzEdAAYpAiqBhoaCpoJdyFmmSm8utZP4hF7KsuDN2vZdnRPSiOQc3mxaOJ9Ax6l14ueajvvcByZrBao3bQNjfFDhpxIT8CM8WQ9rXo8EIAfzhgepDLa6fhxhDRmzGmWW6GFVaRYwhoMPlB6s7vdekZeHyGw5kuDZsKfKHQDcLu0SEhf3QQKz+f5A0Fwa8KM6ij6VXOW9lkAS46tRm06eHnXx4rLmktVtg4gSgXLbJ/L19haIQ4K7gdU1RYChw3XREc3JTOvy1Tnc3+s6zNNsPnDuKTg7QbbAjUXqsbUDiWPzqIRYN8CODehP1U/qoGfd+I0z2M1d989g531mL7sSs/PIj+afl8csHVrVd2cZfK169rBuk24b1fqhYChcXb93LZLQMgFnXfyT/W0UtVtUi5hk8U9slBKszQdnpNJUJhMn+QvtdwRI+pOWtqodZ0S2PVcwZ9qBEx6MP0AekmvCwzoFV781Kevr1qO/hhj43lK6XBHpSHtCcfL4LG0qaIsBQwboFebt+x5acDEr71D740AFukb/iMcmVykI7m+vQgNT990Q1iuV9+yu4GNX23Fwzds6oGfdieU/zduOWDg+x0s4aMoANbZxR37E/LJ2Vz6IXSg+nm/2Ae+Fqv7KJ6x/5FABwG+N2xlau1meIMnwgr1C9YC3wuywDFMEhZP0wyM/Izy3mp9PRHBCe74CwZfZBlQ6KfbiqM8YPeCVY9b/4bP2pj//0s+W9pljeQHP9VPv50hn88xzMWSlFKWDJtiXUo3wl53mudhbApWazhM4L+uYM0Jf/cQTw3/UjTqohWW5le8647SxoHufxx8QcSU35I0oxZXVgK1IcJ6IXrIQnJ20li+MR5KVyODtIuDX809B+nGSrO/emm5I+XbLEDBNC0EY5MnRGYjOdnss2b1YWlEPWCENDHTIpAHctbe7zhE+fg0UGqAYx2g43tqB7NHy/WYkMEhDGzP45SjZrujyn6pzWaqYIHJYRlCiyfpyTacZrd4ZMySKf6w7n5p2jbs9sICOupP63biSsKlRCWVxP4eK2Kc5d7NU3ho4HDLOGfsnNmCxwoZ8D1RuHrPtCYR22beyluMQS2wocSmK0uLNVIDq/C7Njg8XF+o6dYpyhVcGv40B5hefY7vqJTjUAp1XfWKHfPfChyzlJM5bgAxz6GxmuW1RwnGDYgU9jwyQSkkTEEciFQSli+UQNuHfbtGoQGMnznNFdba6Jj+BTmljQo3jAzsxACNeJbRFgKKM5Y8gYAt444Iz9TcG53+eYmZr+/hk9Z3gt/zg4+LluLdkuvFd/V1yyDLoAWycUZ3M7Chg7uM4VvLewrbq/G8CnAefW/SxW72zXzpSGvGd9iwlKgsfDfvzK1NJNqo2/4rM8xYc7sBo0y9nbZ+hVqrMTa8S0yDW+sD1hz2kHgWw3o2IH2dlpYtiVV4xtqlh/UUYZQcnEBlzqYMQSaRkD9nPXzPd9HsCh/a8BXdeSilNgpHYcNwpL1mbHpsK3dOJjsp2UENEgYIBdyrb8vwOiK8PWXkUdRuHd93eyLEdEGW6QhMERAg5l1Sm4794a1Sy6KYErOBdWmcXpkTnBbwdeE5QnMtyxtI9BpQkzJPpOfjUgzGRGoFV/b4MnYCYxU/QhoIKNR+ukYR9v6s8D1t658C2rG9/8T0p95Z8cquQAAAABJRU5ErkJggg==",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}\\frac{gw^{2} v^{2}}{4} & 0 & 0 & 0\\\\0 & \\frac{gw^{2} v^{2}}{4} & 0 & 0\\\\0 & 0 & \\frac{gw^{2} v^{2}}{4} & - \\frac{g_{1} gw v^{2}}{4}\\\\0 & 0 & - \\frac{g_{1} gw v^{2}}{4} & \\frac{g_{1}^{2} v^{2}}{4}\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡  2  2                                ⎤\n",
       "⎢gw ⋅v                                 ⎥\n",
       "⎢──────    0         0           0     ⎥\n",
       "⎢  4                                   ⎥\n",
       "⎢                                      ⎥\n",
       "⎢          2  2                        ⎥\n",
       "⎢        gw ⋅v                         ⎥\n",
       "⎢  0     ──────      0           0     ⎥\n",
       "⎢          4                           ⎥\n",
       "⎢                                      ⎥\n",
       "⎢                    2  2            2 ⎥\n",
       "⎢                  gw ⋅v     -g₁⋅gw⋅v  ⎥\n",
       "⎢  0       0       ──────    ──────────⎥\n",
       "⎢                    4           4     ⎥\n",
       "⎢                                      ⎥\n",
       "⎢                        2       2  2  ⎥\n",
       "⎢                -g₁⋅gw⋅v      g₁ ⋅v   ⎥\n",
       "⎢  0       0     ──────────    ──────  ⎥\n",
       "⎣                    4           4     ⎦"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "W1, W2, W3 = SU2L.bosons().components\n",
    "B = U1Y.bosons().components[0]\n",
    "M2_gauge = model.gauge_mass_matrix([W1, W2, W3, B])\n",
    "M2_gauge"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "15",
   "metadata": {},
   "source": [
    "A `Rotation` is `new = R · old`; register it and every later\n",
    "Lagrangian extraction rewrites the weak-basis fields for you.\n",
    "\n",
    "- `rotation_2x2(θ)` is `[[cosθ, sinθ], [−sinθ, cosθ]]`; the minus sign on\n",
    "  `thetaW` gives the textbook convention `Z = cosθ_W W³ − sinθ_W B`.\n",
    "- `kind=\"unitary\"` tells `Rotation` to use `R†` (not `Rᵀ`) as the inverse —\n",
    "  needed for the complex combination `W± = (W1 ∓ iW2)/√2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "16",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "('rotations registered:', 2)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Z, A = sp.symbols(\"Z A\", real=True)\n",
    "thetaW = sp.atan(g1.s / gw.s)\n",
    "model.rotate(Rotation([W3, B], [Z, A], rotation_2x2(-thetaW)))\n",
    "\n",
    "Wp, Wm = sp.symbols(\"Wp Wm\")\n",
    "U_charged = sp.Matrix([[1, -sp.I], [1, sp.I]]) / sp.sqrt(2)\n",
    "model.rotate(Rotation([W1, W2], [Wp, Wm], U_charged, kind=\"unitary\"))\n",
    "\n",
    "\"rotations registered:\", len(model.rotations)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "17",
   "metadata": {},
   "source": [
    "## 7. Extract the Feynman rules\n",
    "\n",
    "What happens inside `feynman_rules`:\n",
    "\n",
    "1. `physical_lagrangian()` takes the full `L`, applies the vacuum shift\n",
    "   (`H0 -> (v+h+iG0)/sqrt(2)`), substitutes the tadpole solution for `mu2`,\n",
    "   then applies every registered `Rotation`'s substitution dict — all\n",
    "   lazily cached.\n",
    "2. Since the kinetic sector has `∂_μ` on every field, the result still\n",
    "   contains `PartialMu` heads. These get Leibniz-expanded and sent to\n",
    "   momentum space, `∂_μ φ → i·p(φ)·φ`, tagging each leg with its own\n",
    "   symbolic momentum `p(φ)`.\n",
    "3. `conjugate_map` rewrites `conjugate(Gp)` → the plain symbol `Gm` so the\n",
    "   commuting-symbol extractor can see it as an ordinary field leg (`Gp`\n",
    "   never got a VEV, so it stays a single complex physical field with its own\n",
    "   antiparticle).\n",
    "4. Every vertex gets the standard normalization\n",
    "   `i × coefficient × ∏(leg multiplicities)!`.\n",
    "\n",
    "**Dict keys are in SymPy's canonical sort order**, not the order you write\n",
    "them in — look them up with a small helper (this is exactly what the test\n",
    "suite does)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "18",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'38 vertices extracted'"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gp = H.components[0]\n",
    "Gm, cmap = conjugate_pair(Gp, \"Gm\")\n",
    "fields = [h, G0, Gp, Gm, Z, A, Wp, Wm]\n",
    "\n",
    "rules = model.feynman_rules(fields, conjugate_map=cmap, simplifier=sp.simplify)\n",
    "\n",
    "def rule(*flds):\n",
    "    key = tuple(sorted(flds, key=lambda s: s.sort_key()))\n",
    "    return rules[key]\n",
    "\n",
    "f\"{len(rules)} vertices extracted\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "hWW    = I*gw**2*v/2\n",
      "hhWW   = I*gw**2/2\n",
      "hZZ    = I*v*(g1**2 + gw**2)/2\n"
     ]
    }
   ],
   "source": [
    "print(\"hWW    =\", rule(h, Wp, Wm))       # i g^2 v/2 = i g m_W\n",
    "print(\"hhWW   =\", rule(h, h, Wp, Wm))    # i g^2/2\n",
    "print(\"hZZ    =\", rule(h, Z, Z))         # i (g^2+g'^2) v/2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "20",
   "metadata": {},
   "source": [
    "A momentum-dependent vertex — `A G⁺ G⁻` (the scalar-QED photon\n",
    "coupling) — carries the antisymmetric combination `p(Gp) − p(Gm)`, since it\n",
    "comes from the derivative (kinetic) term rather than the potential."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "21",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANQAAAA1CAYAAAA5xsU/AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAMyklEQVR4Ae2d4XHcthLHzxoV4CgVPLkDx67ASgeOXYGdDl7Gn+xvHrsDJRUkLx0kruDZ7sDpwBp14Px/EJYDgiAIHnkiT+DOUCCAxWKx2F0sQB61+/bt267G6/Xr1+dzjHsuOnPwEtJYK18hj0P3S4yhpM8czsmuQnjz5s1/NeyHU4c+F52pfPS0P/f89VSvu3hB2ZbIrRfnHl4iBD+QB0p/Dsvvyr3G9VRjeaz0lyljGkNHuPfV10tdP/o+r336l+p+1XWu/M9KJ/HkaTaJ6OE4rumjKTyCG/FbNEfCO4hcS+TWh5MyqC+S+Z9qMOvkrmEe/QR8UPrDFH7G0BEuSv1K11tdGI8Z0073KA5GdqHrnfKzK75ofhLtJ0qbfpVfLYhPjGRwjoR3ULmWyC2FcxpLVkgP4rI7lH+nsVzOMJ4iOpLlX+qL1ecH3f8T96syHBcGBc7fcf1MecYLv8cScQzK9pbkWiK3Dk5te6hnmow5VoFBOn7SWXmSxhQYy/90T1jWMbgAZ+9bP174xfMfA2Rle1tyLZFbCqcag9LgCa8mK20JHeHgZTEm9kVDodaV8P7QdUhg3M8O2cEctIdku4BcS+TWwnEhnxjFe7FpJtz7pHzSi/sBPxbOV13An7ou+vAdRuZPjp7qjKfvReJL2IfuOaFDWV0YoxRjQYlZDfoUmNAqGVZ5ehgAY3uhixCMPH1z/yKg20tHeDvhgU98/4/uk3IELwAmpAlD1YZ+GZ/xcqZ7xonM6ZuDjPe6wHmuC+CePVhyfKqjnLYl/AhtPvB8Tpat6EyVKzIq5cMEUCK3Fs6Jb/lKDL/XPeEHitkB1VMO3i8eF2Niw/tTB7mgoICe8fR/kYt5eqUyBGyAEZF/ZAWJlDoOXFLw3I+Jug+6OBZFaTmYQeGRi0GODji2V2mMxBqmUvVBuPc5qPuRvpXHofym6yl8+DJkjeHg/M58OTzSV8ijsi1g3KG8WpUHzswl26lyHcOHiaREbi2cE00KlovSAkwYCtQC4bAC4HWfWIXKwGPS2Xg7UNl9XS91YWi9oPosPdWHPOFZCYtCwNOE/eIlUMoO70EjeI3p7KK+UDpWFpyFAQKjP4MkHatUarjwNArUbzgueEGeGJfBmb+hPKRPOXz1AeO+dYMSj+E8TpXtFLmO5cPkWCK3Fs6pWqJAnz0F4uy3Ri1I8ZScSF1bmRcWWad8Pm+Dzk0ubYboXfXx5PuBfqhQ0PyoupxBoXQN/zTwEI4fwcfjd2MRbZSY9n10jJ5TXOGaTK3cpSrHMcUnqdcqZ6VhDNzTJ3TMK+vWAfwBcejGY4Bkfw77xtG4cfh8MvH9flDlIG5A4Ce16+t7TtlOketYPmx4tmhYPpW2cE6ZQLCUsmogyNZkqZxJpPx3XSFgPEy+U2KlCPWzpxPite5L6dHI04p5sn7jSWQlGQ3qw8ZvziA2VFZIxunwCjq4Eg48J0F03IqjlDCN0I29oJN50IfjRfkUL8g45gVH2Jq3qPMhJ+DQPd1Jz+jCfo1PpXPIdrJc9+CjRG4tnJNAAOxL3CqkjkOFYM8AxAqMosUT7hAH/oyhRx+xAnX6Fb8YfcxfzEZ2QoQMXTxZrKwowx8BsSE6jg/RcR41aBffmsdNGQKHDanxwEtL5uqHMubL7dmUD+dOxQ4og++lYA7ZziHXUj5MTiVya+E4g/KTgFK6SVGKcbVAOE045fGZyGYf00IuyBTSQ+mafj3ZVBkbzpaiJViAjlPiRB1FKWVlBbnWRThmMETHcFnxkyBe4aPTX4DcqfMyp10sc8JCnA7OgHpWqxjworEcY5xD5lPjGSvbQ8k1xYfJokRuLZwT35KJ2GlC/vaTEoZPzjv7cuuIPRAQeu6bkuG/Y+jhlfAADsQDE8MAHL8U+rJYyaiKAVocRXdANOjjoa6QLvdMYvzaTi8d4cIPihuexlHcgOrph1NL8Dp8q55+4ScOsZ2hqD52HOBaWRM+qiyEoT1WiDvrvfidRbaiM1WuY/gwGZTIrYVzSksxi4cj3MNayTdhiO7ZP7BUclzLaeD3uq50pcIjFedhDD3hclx8ySWq17p4FvMfXb8FZa1nVKrrAxS072gZQwXoj0MDgHFyhB179hwd11BtkOV3yiAzjAbejQ78spHHcJBjDJSHm2irRyGaebFCpRgn9ODbIoyg2t0yvviAI8Y5VH422WqMU+Q6hg+TRYnc2jhDv0lJ1ev3IJ90XfbUPVXdl1RdX1mOXl+bfcrhS9fDuK3KLuEhLu/L99Hpw1+yXLyew+9SPKjvVch2Dz4G5ZaS7amZYiqVRyA8YT/VvCngyyh/kWqTK5ubXq6vnjq8OZ469tZ4mfDZU0/zpriPToOwohvCVvhdCtYi27F8lMitg3MyIGWMCUYcyCAIO9g/ERqxl2hAed4uYOJo4+6Vt/DJ8IrpWYM5U/FDyARvhFUOdM+YyNvD7ZuKzN8UnQz6YlV+nIw3FSoenK+1yHYsHyVy68M5HZDqpa9/KQLsKVA+jMk2wU1zlbFHwGJzUEwvR2RiHQcG7KXYH3ESh5EDvOrE6zylytfQuWm+yr/IO16Nb4XRtch2Tz5K5JbE6fzA8FakvXAnEjIrEu/IvZ/Cylx0pvDQ11a8ER2wibfDkD7UVZYvJdsSueVwqjSoVWrQxtSdkMDQHupODHIbxCaB25LAZlC3JemtnyokcE9n6e3PHlUx7G2QmwQOI4FtD3UYuW5UK5XAFvJVOvHbsA8jgc2gDiPXjWqlEhh6sFupWPqHrWcQ256zXzxFNZLhvSLEI0TaDGrEpEkReIeRh6W8JbHBJoGOBDaD6ogkW8C7ikOvV7UIyPh4Xctebzr3lc3Lxi3kLXP0EtgMqnAK/eq0U9p6KbigOb+Jat6n0z3vgPFVqPhDLQWkjgNFY6zWiWyHEuU6yiozanXypHmx+CLohlWON8AJH+8q4ER4iZqL8JgfUmY/LXdXBLEZVMFMSikI1TCCsasT1FmdPnJTEdToRNz0biFfmZaP3jsZWRlh/HMQDCz1E3drchfSGp2Im7fNoAbUN1idOr8BG2jaqRYtwjx+g8WHPfYC0SB8XOxHgyVMV+pEnGi2kG9YQ9g3vR1Gy2N4w2Sly/1DgzyRm1o2/FxHAXM4kaMYqGdyW6EysyVlQHEfKSWE2Ru8MbFB5+tRO58nPcof/5UKwo9zDidS2uXieJtB5aeAk73s6iSlIYR7rItPnAF87OVC5W7v5JWKo3IMyk72MNB9TgzVbB7I8a06HAmflOOzB63PtPkxDP4rIeFxkFOdE9kMSrOeAq9U7l/JpOopEw7eF+NxeyKlKBHHw5zq2WEEeRS0dWws3EmrnujtDQV8278twlnwUR4bC33iZMKQ81p5xv1Il9tnejmszomIv4NDdQalyXbKoBRFyAEeGoNJgtqjbHy3gQ9aOlAZp3fQb74Iq3xT79EWTYb4Vj2rqH0BihCVZ0ghXCjTrNrC52vDPE4Iw9fVOZFwAIe8r8qgvLLw71r4HPTQCkFYk3ubAc/Ne32NYXr6zNeYb/yB3wHRwsOjvDGcUaD6FP98AXjoPcMhvq9Ew563PVNXofFgbDiM+MTzo9o0BqX7VTkR8XtrUI1BaZJZUQDCMVaflEJSvxMuKw8KnQTVm2L9HiFgAM2/+InqRmXVR5I/lTMOjs1Hf7FpDN++H4wnDPdsfGZwNqbwW/hWVmVak0G5VUOKworCvoen+aGyhAowtDqxXwBixSJEanlv9YFS4umhuffzJ7WfA4r5VmeMhRUPeRmkxodzieVg+NWlJ7WNWAqCwhOeJE/ZVM/qVRSyCTcMczAcPHi4f0LZoEcd1ypgiG/PJAcNzfgyZc+9TFcxtqWZqM6gvMA5bCBswgBiKHmQW/QveUQfD09oFitm3Odt5Yv49syw6jROwMvqTGUYmgNf1jgQK685rfYjLVIGfnnLCRVhjAPdsz95rDS5enk0w8UY2edwIsbzGrdCqW3nIMPT5Q3sTp3ajQJPa689FB2p/Ri+bR95raZfdREic6hxpYuy1jMq5auH04olgHKwj0I5bQV5pbInJTJRG0LHZr+kPEfFTb6ExhI4Y/gWbupgZOgUcYlhrabPmg2KsI/9DasRBwZ4bo5/8by9oHr2RRhe86tbX0b5i96G81XAX5bHVFcr4DvF1p0rq9agpGA8hGVFsSN0DKzE+2JMGJ8D0bivG8IgXrM5+GmX5/mm83F/F+V7HKvHi12tQfkpw4g4Quff22BgFvrlZvTSVxIu2t4JY1p7uHesfOfmYnV11R5K2EzIEHgoycnVg0KDsqbFqehy2DHLoURxpxviIhKofYVC6KxS/PO1ktVp1CSJJobKxp4QkcMP+vqqdPRbDmq3wRFI4F9elNYzRgqvbAAAAABJRU5ErkJggg==",
      "text/latex": [
       "$\\displaystyle \\frac{i g_{1} gw \\left(p{\\left(Gm \\right)} - p{\\left(Gp \\right)}\\right)}{\\sqrt{g_{1}^{2} + gw^{2}}}$"
      ],
      "text/plain": [
       "ⅈ⋅g₁⋅gw⋅(p(Gm) - p(Gp))\n",
       "───────────────────────\n",
       "       ___________     \n",
       "      ╱   2     2      \n",
       "    ╲╱  g₁  + gw       "
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rule(A, Gp, Gm)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22",
   "metadata": {},
   "source": [
    "If you want typed objects instead of a flat dict (useful for UFO\n",
    "export or for grouping by vertex topology), use `model.vertices(...)` — same\n",
    "arguments, returns `Vertex` objects classified into the closed catalog\n",
    "`SSS, SSSS, VSS, VVS, VVSS, VVV, VVVV, FFS, FFV` by inspecting the spin\n",
    "content of each leg."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "23",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['SSS', 'SSSS', 'VSS', 'VVS', 'VVSS']"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "verts = model.vertices(fields, conjugate_map=cmap, simplifier=sp.simplify)\n",
    "sorted({v.vertex_type for v in verts})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "24",
   "metadata": {},
   "source": [
    "## 8. Pure-gauge self-couplings\n",
    "\n",
    "Cubic/quartic pure-gauge self-couplings (`γW⁺W⁻`, `W⁺W⁻W⁺W⁻`, …) are *not*\n",
    "produced by `feynman_rules` — they come directly from the group's structure\n",
    "constants, since their Lorentz structure is universal. You supply the\n",
    "mixing matrix `U` from adjoint components `(W1, W2, W3)` to the physical\n",
    "bosons `(W+, W-, Z, A)`.\n",
    "\n",
    "`cubic_couplings` returns the vertex coefficient directly: the cubic carries\n",
    "one derivative, so the `i` from `∂_μ → i p_μ` cancels the Feynman-rule `i`.\n",
    "The `i` you see in `g(A,W⁺W⁻) = i·e` below therefore comes from the **complex\n",
    "W± rotation** `U`, not from a Feynman-rule factor — in a real basis the\n",
    "result is real, which is why the gluon comes out as `ggg = -g_s`. The\n",
    "quartic has no derivative, so there its `i` does survive and\n",
    "`assemble_vvvv(..., feynman_rule=True)` adds it.\n",
    "\n",
    "Both return **weak-basis** values that carry the group factor (`-g·f^{abc}`\n",
    "for the cubic). For a broken group in the physical basis that factor *is*\n",
    "the coupling, but an unbroken group exports as one particle plus a colour\n",
    "tensor, so there the coupling must be colour-stripped — see\n",
    "`export.ufo.vvvv.adjoint_vvv`/`adjoint_vvvv`.\n",
    "\n",
    "We print the magnitude to compare against the textbook value directly."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "g(A,W+,W-) = I*g1*gw/sqrt(g1**2 + gw**2)  -> |g(A,W+,W-)| - e  = 0\n",
      "g(Z,W+,W-) = I*gw**2/sqrt(g1**2 + gw**2)  -> |g(Z,W+,W-)| - g*cos(theta_W) = 0\n"
     ]
    }
   ],
   "source": [
    "cw, sw = gw.s / sp.sqrt(gw.s**2 + g1.s**2), g1.s / sp.sqrt(gw.s**2 + g1.s**2)\n",
    "U_mix = sp.Matrix([\n",
    "    [1 / sp.sqrt(2), 1 / sp.sqrt(2), 0, 0],\n",
    "    [sp.I / sp.sqrt(2), -sp.I / sp.sqrt(2), 0, 0],\n",
    "    [0, 0, cw, sw],\n",
    "])\n",
    "triple = cubic_couplings(SU2L, physical=[Wp, Wm, Z, A], U=U_mix)\n",
    "\n",
    "g_AWW = sp.simplify(triple[(A, Wp, Wm)])\n",
    "g_ZWW = sp.simplify(triple[(Z, Wp, Wm)])\n",
    "e = gw.s * sw   # = g sin(theta_W), the electric charge\n",
    "\n",
    "print(\"g(A,W+,W-) =\", g_AWW, \" -> |g(A,W+,W-)| - e  =\", sp.simplify(sp.Abs(g_AWW) - e))\n",
    "print(\"g(Z,W+,W-) =\", g_ZWW, \" -> |g(Z,W+,W-)| - g*cos(theta_W) =\",\n",
    "      sp.simplify(sp.Abs(g_ZWW) - gw.s * cw))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "26",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAACQAAAAuCAYAAABAm7v+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAADMElEQVRYCe2Y31HcQAzGOYYCSOjA6SAkFYR0QJIKgA7C8HT3liEdkFRAoIPQAYQOoIOQ6+Dy/eTV3t5i3x/b+F5OM7uStWvps6Rdez2YTCZbfdJoNNqVv7Pgswj8SPox8k5Q9MnO5fzEHUq+kPxH7Q26bbqe6VggDhKf55IL6d6iWwcgonOH8yoa9F1DOQhFhggdiq8tZRFTSNOhFPuuXEfKzLfAsMKIzr5kW2E2QMr6bsPhsFC7cL/huuC69xoKkWGpn1pEyo5CPyVSz/YhKcnpe7W/5dyta/ED6X9wLc7Gdqy2p/bgesmMsXRPxG2fEcdWnhb2HGzAI/k9M4Ck5GacW5GJk2duZJkaIPEz6XkanP1M9BJtB8aZE7WBjXdqNyh13yt4HW37QHDwVdcfEt2jZBz8Rqc5ROAWWfRR7cmkaceGZ3NRaT4g7tWwsxSlEeJpr2WEpzIKAJBJG/QoHQ6gz2rfTFIX5gLeIuF68TuNLQ3IIpQYu0wMIfLEYzcobmDFSRfOPY0S41wHjA56KNlyvaeMHEO5MdKSPzHz0N87QBSiZ3M1Topzmza5rnNANi4DMbSSicBMTSRGCslxbtBX6b7ITtUDJaZmRQf0C7VuxqgTNQXZWCnGnqcGsJHuA/hrtXh/0MUCL2cu7uPGGAywf7CK2GMsQtLbSy83JT2bGzRWY8+inniIJzV0M3uUrpeiCCifLYe2/4jbJpePv9T1jhxSeHxSxs/IoEN/9FKO6+xSQ4ChBowEhlQRenbjlVZIaaFdP9CbFjB57VwJzEqrox2M6d21NTSd0q/ky75fr3O8bQDNCY4NUdT9Hl0XINoU9YIAreXkOhfTZpXNDY8GO40QL2U1/yxZ5LtyvFNA8nClxodaY+oMkCLDEao1dQKIVAkJX4m0VtQJICHgYz49EjUG1RpQSFWrQk7RtwIkMIWMxYNkarip3AqQnPKno5NU+QM0BiQg/JLpLFWtAIVU7Yrnp1e325g3+vwQED8Y5I75CcHS54DAn5L0L1k+t/K6EaBKS1IKwD+xG/FPdXMW6RvXUI1hznS0xtRJhBQRipstwA+c/OC6lf77qsj+AwIocHsGUD5aAAAAAElFTkSuQmCC",
      "text/latex": [
       "$\\displaystyle \\frac{gw^{2}}{4}$"
      ],
      "text/plain": [
       "  2\n",
       "gw \n",
       "───\n",
       " 4 "
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "quartic = quartic_couplings(SU2L, physical=[Wp, Wm, Z, A], U=U_mix)\n",
    "sp.simplify(quartic[(Wp, Wm, Wp, Wm)])   # WWWW quartic"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "27",
   "metadata": {},
   "source": [
    "## 9. Adding fermions: the SM lepton sector\n",
    "\n",
    "Chiral fermions use a **separate track** from the scalar/gauge extractor\n",
    "above: a fermion sandwich `ψ̄ Γ ψ` is a c-number \"atom\", not a product of\n",
    "commuting symbols, so it can't go through the same `Poly`-based extractor.\n",
    "\n",
    "- `WeylFermion(name, reps=..., chirality=..., nflavors=...)` declares a\n",
    "  chiral multiplet; its gauge components are `IndexedBase` objects (one\n",
    "  `Indexed(base, flavor)` per flavor), not plain `Symbol`s.\n",
    "- `Bilinear(ψ̄[i], Γ, χ[j])` is the opaque atom representing a fermion\n",
    "  sandwich, where `Γ` is a Dirac structure built from `diracPL`/`diracPR`/\n",
    "  `DiracGamma(mu)*diracPL` (see [`CONVENTIONS.md`](../CONVENTIONS.md) for\n",
    "  the `P_L=(1-γ₅)/2` convention).\n",
    "- `fermion_gauge_current(field, flavor_index)` builds the interaction part\n",
    "  of `i ψ̄ γ^μ D_μ ψ` for you — the analogue of `Dmu` for fermions.\n",
    "\n",
    "**Two things this section does *not* verify** (documented gaps, not silent\n",
    "omissions): `check_invariance()`'s hermiticity check is *skipped* (not\n",
    "falsely failed) for any sector containing a `Bilinear`, since verifying\n",
    "`[ψ̄₁Γψ₂]† = ψ̄₂Γ̄ψ₁` needs Dirac-conjugation identities not yet\n",
    "implemented — write `+ h.c.` terms explicitly, as done below. And\n",
    "discrete-symmetry invariance checking (`S3`/`ZN`) doesn't yet cover fermion\n",
    "multiplets either. **Gauge invariance and mass-dimension checking of\n",
    "fermion terms, however, are fully implemented and exercised below.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "28",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJ8AAAAVCAYAAABYMDxeAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAFtElEQVRoBe2a7VEWMRDHT8YCGO0AOkDsADpArEDoQIdP8M3BDoQKUDrQDlQ6gA5EOsD/L2Zj7p7k7nL3vBwMO7NP3jab3eR/mxd4dn9/X8V0cnKyIb6xOuXfK3/oy29UvrK2p7T/DCx7Hpc9HjPRNqbaarhC/jk/Rr4zwAvgU/6l+FBt303uKS2fAc3fJ/X6pPRU6YtyDcU9lr5uHT4Cvj0v45wJ4KNBNS/jxmJ3F9BB9hB5X4uxDyLy8nF8U9sZFaukqds3Zm68b7vSseP1EIDuxOu+TPJZcpdROZmVzHfxlvhA7NbNgU8FlB0pfZXsucJK2UTEqJRyPriZmo1Tt2/M0nnfiNbMPeABiIFU3lDhWuml+E1oyGQkg65f4i/iuzUvx1bwOdNn5dUydMsb0fmFrcLYqds3Zk4i37419aiNHYhoyHYKEPsQOANvlYFvX51XvoW1WG5hf2YCWvoss2nq9o2ZC/Mtd+YvOr96nIG39TX9cJaKLxhjDF1U37colq25CVjUuH31Tt2+vn6k5My3mVcOrQc7EsyWXIIhZPc587GPFy2qBuJrYFAuAu/EoJ/nmN9i9HEZcGc15SvlOVMeiDfF7Pm1KKvyV9W/U8phNkWMNeN8SnAedd7eI+nCHwi7P0zFPmfRwJ+Bvs3gQ3rYZlm3Xue9hrno2wV82+LS896uBmcx2AbPxT8oK62UAqw/SjHKvgYuM8gTZZEP4FMdoOTMkDywqr4r7Kv7/MiPx3yEN03VAX7A6HyMR1u2ffHYpfkRvlXqy6sDxBMOwLsVhzlSvoSuJbwD+IhKKOpF3gE7e2EEN9AQ5VQmCkLordTGwv0gLyIqNseiri2q0Q7ZmP9KC/iVrfjDOLxrxjbxAeTGX5p9Y1we6Rvz4QKJUtb1XOzWfqBNYGAD8AGW3HaX0v1TBtx5IzDA/vphsoANsqgHOG0h91X/0bX+/2Fhm3X/W/0bk3TMhH4TwhZxiQ/WtZnaDoC++EvnGJEbH/urlnba5mVf096S8lDfWGtbS3xhnol4PL8kd4MeRqFvHfAVkR+cPrlJ54B6ZXJRijxfTbzlAl7qcgurJne2NPBSrpH0s23Tfx7gw0aOC3Ekr42XKPCxLcu+xPC9q4b6NvO8pflhzSDWbwi5gLemnoRAU1aiyIEs0YFznfvKIiMR40xHFIxBgmyluuTiqd4BXCIXyGWI82f4MjMyJdV2ROjssyL7Ou1qERjiW+q4YesCdoYQeLsFfCzcEARjQC1iaTEcmFT/RQwRlo0YowkyzktOh/oeiJsfgZ2nauOYQsmzNdYmJ6HDxPukzAUH6hmSXvMtbiu2L+6cy4/0Iae21De7ANpaxnp55YC4ODiSzXbcsqq2lMh3A/gAhClr6xDaNBBAAijNiOQApnbOCWyHcXstOnljtyVjgNykn8oxZSOjZHklh8PkqA6buGn/ipUU5NFnX7Xrhk4xkdzsjNUV2Rd3zOXn4ENOdW/fZAPryPk8RwF0CEieOaNPX+LPuFfP9XMh5r2mhBgovkhYX85zr2UMwLtVGi8YzxTnqmMSeD8DjBhxqjoiGHZUygMgu005h1QX20edfWWczwJgyYvRC2Dg0Ka6TpL8me+HjTbB6AmXKtrVNsi+TgMkIP2jfMiNIb2dvtFXcvgO2dmcZzL+fmtlZNDF2r1Viiw2h3YqOgiwHlb8P9/x8fG1eIt8k1V/Kt5p1k+5LHv3xOtTtLHvfI71oe84i5ijtrHVtiG+Ztw1j1DQHr7uDtQ+hGaib1HUm6BTj8GH1LSyA7ro6sDnQyb/7Oe2uVSPh1InH9gW7c9iD8Xsmp2PwYeaQ77g8QXO3BZtkY9mbjfuicTLPtSEW3PJO90U/XwMPqTmFXyFHfYZe6+RR2btX51VR4i0ZwZetuNLhHV9SjtmQPPGpYqJ55mBeUw+H3Wo6d28inVr89G3xX/vr/4CgUxCVVdkK9sAAAAASUVORK5CYII=",
      "text/latex": [
       "$\\displaystyle \\left( \\left[ nuL, \\  eL\\right], \\  \\left[ eR\\right]\\right)$"
      ],
      "text/plain": [
       "([nuL, eL], [eR])"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ll = WeylFermion(\"Ll\", reps={SU2L: 2, U1Y: -sp.Rational(1, 2)},\n",
    "                 chirality=\"L\", nflavors=3, component_names=[\"nuL\", \"eL\"])\n",
    "eR = WeylFermion(\"eR\", reps={U1Y: -1}, chirality=\"R\", nflavors=3,\n",
    "                 component_names=[\"eR\"])\n",
    "\n",
    "Ll.components, eR.components"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "29",
   "metadata": {},
   "source": [
    "The Yukawa Lagrangian, written by hand exactly like the potential\n",
    "was — a scalar factor times a `Bilinear`, plus its hermitian conjugate\n",
    "(built explicitly, since it can't be auto-derived yet):\n",
    "`L_Yuk = -Ye[i,j] L̄_i H e_R[j] + h.c.`, with the `P_R` sandwich carrying the\n",
    "right-handed projection and `Ye` a general flavor matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "30",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'Yukawa and gauge-current terms built'"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "i, j = sp.symbols(\"i j\", integer=True)\n",
    "Ye = sp.IndexedBase(\"Ye\")\n",
    "nuLbar, eLbar = Ll.bar_components\n",
    "eRc, eRbar = eR.components[0], eR.bar_components[0]\n",
    "Gp, H0 = H.components   # (re)bind: Gp was set in §7, H0 not yet named\n",
    "\n",
    "LYuk = -(Ye[i, j] * Gp * Bilinear(nuLbar[i], diracPR, eRc[j])\n",
    "         + Ye[i, j] * H0 * Bilinear(eLbar[i], diracPR, eRc[j]))\n",
    "LYuk += -(sp.conjugate(Ye[i, j]) * sp.conjugate(Gp)\n",
    "          * Bilinear(eRbar[j], diracPL, Ll.components[0][i])\n",
    "          + sp.conjugate(Ye[i, j]) * sp.conjugate(H0)\n",
    "          * Bilinear(eRbar[j], diracPL, Ll.components[1][i]))\n",
    "\n",
    "current = fermion_gauge_current(Ll, i) + fermion_gauge_current(eR, i)\n",
    "\"Yukawa and gauge-current terms built\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "31",
   "metadata": {},
   "source": [
    "To add fermions to the model, build the **complete** Lagrangian\n",
    "(reusing the exact same `DH`, `V` expressions from §3) and construct one\n",
    "fresh `Model` with every field — this mirrors\n",
    "[`examples/sm_scalar_gauge.py`](sm_scalar_gauge.py), which does the same\n",
    "thing from the start rather than mutating the bosonic-only `model` from\n",
    "§4 in place. The already-found tadpole solution and `Z`/`A`/`W±` rotations\n",
    "are plain Python objects — reuse them directly."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "32",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "InvarianceReport(14 checks, OK)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "L_full = Lagrangian()\n",
    "L_full.add((dag(DH) * DH)[0], sector=\"kinetic\")\n",
    "L_full.add(-V, sector=\"potential\")\n",
    "L_full.add(LYuk, sector=\"yukawa\")\n",
    "L_full.add(current, sector=\"yukawa\")\n",
    "\n",
    "model_leptons = Model(\"SM-EW-leptons\", gauge_groups=[SU2L, U1Y],\n",
    "                      fields=[H, Ll, eR, SU2L.bosons(), U1Y.bosons()],\n",
    "                      parameters=[gw, g1, v, lam, mu2], lagrangian=L_full)\n",
    "\n",
    "report = model_leptons.check_invariance()\n",
    "report"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "33",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'tadpoles solved, rotations applied'"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model_leptons.solve_tadpoles([mu2])\n",
    "for rot in model.rotations:      # reuse the Z/A and W+/W- rotations from §6\n",
    "    model_leptons.rotate(rot)\n",
    "\"tadpoles solved, rotations applied\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "34",
   "metadata": {},
   "source": [
    "### Lepton mass matrix\n",
    "\n",
    "`fermion_mass_matrix` collects the coefficient of the\n",
    "`Bilinear(ψ̄_a, Γ, χ_b)` sandwich at the vacuum — the Lagrangian mass term\n",
    "is `−ψ̄ M χ`, so the returned matrix already has the right sign:\n",
    "`M_e = Y_e v/√2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "35",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAASAAAABeCAYAAACD3MmFAAAACXBIWXMAAA7EAAAOxAGVKw4bAAATzklEQVR4Ae2dW47cNhaG20YW0NNLcIBZgBOvYNo7cGBg5tmdHTjwwyB5azg7sP0+gO3swJkVuOMFDOAsYewVTM//0aTMUkkssnRhSToEVKTIw9v/65w61PXO7e3t2VThl19+eau2Hx1qX3J3DslY+fEIGA/HYzdmTfFwUNm2pgt3pjJAAvIK8hS/HJNEa6sMAeOhDK+ppI2HbmS/6c4eJfcHgf5wlJaskSEIGA9D0BuvrvHQgeU3MhL3lf9HRxlZv6n8h56y3mzVearCF0FA++dKP/P793z8RPmfg8yhWLLPJeO8KsXX2v81rqN9+qSPN9qea//PuHyLaY+J8VCZ/JgHpTelC5rvR8EfdH6HCZXdiT0gFBrhOByrxA/VeOz9YBB+DA0rjVJg9L4NeYdi1flJG2IYmr5l3d8k8wEhCw4B4+E0DoSYh63pAo5DO2Ab3Lnhu1HJCynvy9b2e1SelVT9He/HV7pS/mXUAIO6pzy8r5Jw7YWDN+V21Q6eEd6aGR8PkLAwHjwWNaMOHjalC5p/26bgPLwLnMQeUMg7OlZnuJextQ9t4f3chJ12rHoYJ4zRA6Xdkk+xs5yKfwrySn/WxgRQLpevfYzPjeJOb035weXFuyMd2lByncHP2XioTG8PD9V0AThOTR9iDyiLLk0AzyWci2nXwTPZc7kkjxX8HAlDwp/K+6ANo3BfMUvAR34f0ffautaOwTA9lSzjcO1QoR18Wyz1XivNGOgjXhq2qyxmX3MxHk6ArVIe/HE4uy4AlfpG105KH4oMkCaAp4LxYLnGZJrg9zEkyWWbyvF0WP995ytfKMY4kIcxceQo/o19L9NEyiePPjB2yKf6eyUZjBwbCovxatw/pRcZNA/j4QSYG8qD6s+pCyB2cvpQtAQTYCj77x54DECzPFIa5Y73tbsbVA+PBrnvlA6GxhkZ7T9WPkbHBe1DTp+xCHVSxod2UNQ3agtP6b/auHrm+lX6TGnG4wyfYtKcR9ozesp3oVQ+1Bs71jg2zQN4ei5eKMnl7YbTsbFOtTeEBz/+OXWBqfTqgx9PiS7Q1rk29IaLSW89HkrmhyIDFDV7rTQdOoVW7AaiuPcksMoYKFey3BLI758pDgrPhOJL/pcqY8nUFZDl4DsUGFfqkjx3anPlzB3AinFPg2fW1XapfFcbY+ZtkgfxxJ8TxxPbKYQiHjT+GroATil9KD220f+/0Khi2v2krfiJhrs0UBrUYVge4QURkt6P5AEcg8HSjWUaBxDeEl5JCKTdvpfvNGYqY7K0d8j7oV1k6KsJqs9Yz0I7ip3x8QIs02h/L5TK7zUwQYbGtDkegFHzZkndeMsTQFvUZAkPkq2lC8ypUx80JqdTirN0wYMT/1Ezp7iuFzkcHesB0TIG5JUG/VrxheJOg4GgAp4FkyRugupwPikE0o+V9z0Zivfu9VEeS6ngJSHLFa3gQVGtHZDFAwoAMQb+rQgX2pzBc3tffthHtmsupfJRs5Mmt8bDpGAOaDyXh1q6wNT69KH42G7pHX/qQS+LIDzaAGkAnC+h439ri63h3gAk51y1vYIoQzLuvEaUtZeUDEZpzzDtCfoMyWOVYyMXi2JsAD4O7PdZ8lL5uN3J0prj1niYDMshDefyILkqusDc1HenPiif4hJdQN4F1WVlw6mVrj9tL9Uf3e0vyirB6nNSOuWFZDU0t5DGDBntcXNVrclT+jyMK0c+yFaIN8NDBWxLulwkDznHdqwLARDl4YA446M09+YVh6M9IHpSp6zDT2YtXjx7uY0eOIwO/wCNG6l8lmK4y/E/Vq+85KqFrfGg+d4X2Jfa4OyZ9t97DKpxQMd+DEvVh95jW/Pa0wXloRvwwP14ityfed9FI8q7w88//3xf2622e7yaw7avGAiTR4bHVzxqYWE8rIsD8XmlTYfT7dnQJVi3VVtPLv+wFuojYDyslIPsJZjcrINvc6uP0fgj0Lxz7jc6k1zxPRDHjNZ4SKNmPKTxGVIqbLN0gT5yeSgxQLMo2BCAtlA3l9gtYFFzjsbDOOjbEmwcHK0VQ8AQOAKBbA/oiLZHr6J/HS6LP/MNc2ae8ET5XFK3MBMCxsNMQCe6WQsHizJA4mPw2+QSnFpRPgLGQz5WU0mugoOlLcGuZPkvI0a5EYrnt7gfwcJ8CBgP82Hd19MqOFiaAeKxips+Rix/NgSMh9mg7u1oFRwsagkmT6f9HBgk9L4RsZc6KxiEgPEwCL5RKq+Fg6V5QA15ftnFC5SSD8I2FSwxCQLGwySwFjW6ZA4WaYAE+D0xxPmf5s2KRYyZ8CgIGA+jwDiokaVzsKglGEx5wHkCt+/NioMItcp5CBgPeThNKbUGDhZlgDzg4f0j4coX54F4DYKFmRAwHmYCOtHNWjhYlAESHzlvk0vQZkUjIWA8jATkgGZWwcGiDJCsfvxungHcWdUhCBgPQ9Abp+5aOJjUAAkk3rTPlapkkJw96JpEaFih8TAMv7Fqi4eDb5TYmi5MZoAEJC+Qf6e4ecvgWERaO/kIGA/5WE0p6Xn4UXH7XrYpuz35ticzQJo5r3hcxWeQT57F9ACNhzQ+c5UaDx1IT2KAZHh4QXXz8iLtD36KXW1w3w9eFYEPIu68f9b3yZPyb7TxoF7zcnkqbDF4TIyHyuTHPChtuhDxMYkBUvsPBXTs/Qx+clftce8PQ8e49bmxfOX0qM+D0PAKg/FwGqTGPJguRJyMfie0DMCO9+P7GuvJ3WvfHp5OE9QnnhHfxzLj41ExHprDo2qigwfThYiRUT0ggY17GVv70FXyyV3Vu5QgNxY+UNqdtFbMkutMcXOTodJ8CRXvByPn8rWP8blR3LnkUn5weT9KjnRoQ8l1Bj9n46EyvT08VNMF4Dg1fSj2gDQB3r8TzsW0KcYzcYYjLpD8S22fozxIcE+xKx+jwPfiOafzyO8j+l4bz3y1QzBMfI+IcfQ+De/b4oat10ozBvqIl4btthezr7kYDyfAVikP/jicXReASn2jayelD3hAfHL4X9r+py0ZNAE8FRQY5X+jrQHSTw5D0ngsXY2pHE8nfor9QvsYB/IwJq5NxSypHrTbUB4yfMYZY8eVBdJ94ZUKPkiGDWOG4XsXC/t8TtTSVjOfWCakvawbp/JojzF2el6hzhSx+tw0D2DqucjibQoO/BgG8aA5zKkLDDmpDyWYSpa5Y9DQg2+1vVVeShcl0oT/KIXNOftGGwbg79r+qS0ZfAd8ipnOMQCxscEziff32lI9Botc8xS78pwCK36s/OarktqHnB1jof0QQp1DE2acGEo8JQwtV88aI+P7YExsOYEbKznRHYwk/yazvw5E/TPvzfJwBG853BbLDOFBdefWBebXqw9HYIrBcU8mKMYQfdKWe0PxXyWLzfkHBuiYcK1KDMAptB8AS4Lek8AqA/DUU+yAE9+0eKk6O5fao4Eiy7/foQAwvZfk/XjxjvaWje2GJUNbzNEZH1/O/nkrr111yv3N8QCYwpvjLIu3KcGP2i7iQeOvoQsMt1cfjsA0/uNlPrFeRNCkk3fTxd2lGiyeCl4IXhAh6f14wDEYL5RmmYZ3g7eEVxICabfv5TnI9oLKnCFQwSHvh7rI0FcTVP+gsWmEdxN4ivF4KWUf8KsEzWWLPFTBOtVpCQ+S5XipoQtMYTR90DziUw/oVOw80FdWONYDonEMyCsN5LXiC8WdBgNBhZwndzk/81jtfE8FxVzt2gnKYykVJoosV7RiIHbkvSweULDQGC/+rY4JGBuMUBzYP8ryx40MTG+Nh4FwTVY9l4dausDE0Z2x9MEBKd3CmLKySem/k+380Ufi7/OheG33+Fh8yaY6H7V9OqZuST9Ty/p5nLf7Uf5Onvb/iGXa+3HZnGk//s3wELD1897hKJTViNfAQx+myt/DWXnPQ77ip7mYS/ZKm8Rv3UnoTsOUmYnVxxNJeSGZTc0vpnGzPON8Ep7MM+2/1+ZOhCvGa+LfKn4FCFfKuAeJ+VIneGNKVg2b4kEc9PJWlYUvq4JF6kMK0y5dUB66AQ9cEQd2dKLvnC3l3UGW6GgPKNfiLVVO2Dxa6tjXNG7joWxlMgX3Y3KgthoP6G63WbJcjwBejoX6CBgPK+XADFCCWLmWeyfCE+JWNBECxsNEwBY0OxUH2VfBNICDb3MrmM/qRIVP7k1Yg+ZuPKThMx7S+MxVmstDiQGaRcHmAmip/eQSu9T5LWXcxsM4TNkSbBwcrRVDwBA4AoFsD+iItkevon8dbiQMd19zmZzwRPm1bwb8MpKN/BoP9YleCweLMkCinbs4uWPaBaW5C5P7EXga18J8CBgP82Hd19MqOFjaEuxKRucyYoRnUHgglBuiLMyHgPEwH9Z9Pa2Cg6UZILyfmz5GLH82BIyH2aDu7WgVHCxqCSZPp31fDiT0vhGxlzorGISA8TAIvlEqr4WDpXlADXl+2RW/WbEps8R8CBgP82Hd19OSOVikARLg90QG53+aNyv2kWP50yFgPEyHbW7LS+dgUUswSPGA8/4R93J5v0/+Ip/Izz3QTk3O4248VCRmDRwsygB5wMMLkMKVL84D8ToKCzMhYDzMBHSim7VwsCgDJD5y3iaXoM2KRkLAeBgJyAHNrIKDSQ2QrDRfkeBEcTJILus5M8nFLwdLtmmFXxEwHr5iUTMlHg4+0L01XZjMAAlI3t/8TvGpvDWw5rFXrW/joRr0Ox17Hn5U3L6VZEduazuTGSAByetLV/EV0oUfFMbDaRBoPHTwMIkBkuHhvcmcLHZB+4MfIlUbXHbHqyLwPbKd98/6PnlQ9Y02npPZ/FUxj4nxwBFTMcQ8KG26EHExiQFS+w8FdOz9DH5wTu1xyZehY9z63Fi+Wnrc50FoeX3BeDgNTmMeTBciTka/EVEGYMf78X2N9eDctW8PT6cJ6hPPiO+0m/HxqBgPzeFRNdHBg+lCxMioHpDAxr2MrX3oKvngnOpdSpD7eh4o7U5aK2bJdaa4ucdHaT5EiPeDkXP52sf43CjuXHIpP7i8HyVHOrSh5DqDn7PxUJneHh6q6QJwnJo+FHtAmgCvvwjnYtoU45k4wxEXSP6lts9RHiS4h0iVj1Hgc82c03nk9xF9r41HLtohGCa+R8Q4eh9G9W1xv8RrpRkDfcRLw3bbi9nXXIyHE2CrlAd/HM6uC0ClvtG1k9KHIg9IE8BTQYFR/jfaGiD95DAkjcciub2gcjyd+CHSC+1jHMjDmLg2FbOketBuQHnI8I1rjB1XFlLfiH8lmQ+SYcOYYfjeaXNBecwHUijjpWZvU+35Ntw4fR3G2Ol5qXyy4Me9WR4A1nPBCXaOgeY4nAz0joaH8qD6c+oCM+jVBz+XEl0o0p0O+FxWqQFC2X/3g8UAxMYGzyTedx3EP6qHoiPXPESqPKfAih8r332VlDrah5zGWJAXhVAnZXwQByQMJZ4S33bn6ll8sGJw3M2NigH/k7bUTZHcWMmJbteGYv5NvtM2a1C/m+ZB8+fY4FhiqxaG8KC6c+sCOKX0oVgXNIcS3enkqcgARS1cK82AnUIrdpZTce9JYJUBeOrhRcCJb1q8VJ2dS+1R/8jy73coMK7UJfnYeDA+Z1i6Go3mGMuwDDrXFud1VZ8qb3M8AKTw5jjDq91b7k8F9IF2i3jQuGvoAlNI6UO2LngsSuU7IbzbmXsgUwDiqeCF4AURkt6PBxyD8UJplmn8g+Et4ZWEQNrte3kOsr2gMkCEwEPeD3WRoa8mqH5z0CodL5/Ijw1gU8cnLhTH4yWbfcZSJWj8W+ShCtapTkt4kCzHSw1dYAq9+lCoC2el8n34HesB0R4G5JUG8lrxheJOg4GgQs6Dc5yfeax2vqeC4r17fZTHUioYCWS5ohUbEarGAVk8IEjHS8F48W+1E1TOAYF3lpoDxgYjFAf2abdm2BoPNbFO9Z3LQy1dYOwH9SFTFxocSuWbiiGhD8Xf13ar7V7pR+1V56O2T8fULe1rKnmN/bm2c9pX/DTuJ+SHPO3/EdJefmc/LpszrXFtioeArZ+34y7k1YyXzoPGn60L4JyST/GgelfaJHJ7NsQDwoZh9fFEUl5IsHUnF2vc/BuxROOqHuNjHu68k/bxmiiPn8Dnigv3ICGH9xO8MSWrhk3xIA7gjPOAcPBM+++1NRcwKjKxWB6EX5EupOSL8JclOtoDSlm5NZQJm0drmMfS52A83J7V5nBMDtRW4wHdLbJW2xPmH9ZCfQSMh5VyYAYoQazczL0T4QlxK5oIAeNhImALmp2Kg+xzQBrAwbe5FcxndaLCJ3UD42jzNR7SUBoPaXzmKs3locQAzaJgcwG01H5yiV3q/JYybuNhHKZsCTYOjtaKIWAIHIFAtgd0RNujV9G/DjcShruvuUxOeKL82jcDfhnJRn6Nh/pEr4WDRRkg0T74bXL1D51VjMB4qE/jKjhY2hLsSpb/MuKe57d4IJQb0yzMh4DxMB/WfT2tgoOlGSCeF7vpY8TyZ0PAeJgN6t6OVsHBopZg8nTa9+VAQu8bEXups4JBCBgPg+AbpfJaOIgN0EdNqg0Ob/w7leeddsamcbHsit+suFNuO/MgYDzMg3Oql1PmQGPjXezhgtHeNDBAPFiJJ9EVTvIhU02KCXH+p3mzYtfgLW9aBIyHafHNaX0BHKCnveEOD7ktKXjAeXePM5p+/0zxSRrLJWFbMlaPu/FQAtrIsmvgYFEGyAPuXh4WcYkhQhHsXqAIlCmTxsOU6Oa1vRYO4nNAeTOvK5XzNrm6I9xG78ZDfZ5XwcH/AboB2lPlgwYtAAAAAElFTkSuQmCC",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}\\frac{\\sqrt{2} v {Ye}_{0,0}}{2} & \\frac{\\sqrt{2} v {Ye}_{0,1}}{2} & \\frac{\\sqrt{2} v {Ye}_{0,2}}{2}\\\\\\frac{\\sqrt{2} v {Ye}_{1,0}}{2} & \\frac{\\sqrt{2} v {Ye}_{1,1}}{2} & \\frac{\\sqrt{2} v {Ye}_{1,2}}{2}\\\\\\frac{\\sqrt{2} v {Ye}_{2,0}}{2} & \\frac{\\sqrt{2} v {Ye}_{2,1}}{2} & \\frac{\\sqrt{2} v {Ye}_{2,2}}{2}\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡√2⋅v⋅Ye[0, 0]  √2⋅v⋅Ye[0, 1]  √2⋅v⋅Ye[0, 2]⎤\n",
       "⎢─────────────  ─────────────  ─────────────⎥\n",
       "⎢      2              2              2      ⎥\n",
       "⎢                                           ⎥\n",
       "⎢√2⋅v⋅Ye[1, 0]  √2⋅v⋅Ye[1, 1]  √2⋅v⋅Ye[1, 2]⎥\n",
       "⎢─────────────  ─────────────  ─────────────⎥\n",
       "⎢      2              2              2      ⎥\n",
       "⎢                                           ⎥\n",
       "⎢√2⋅v⋅Ye[2, 0]  √2⋅v⋅Ye[2, 1]  √2⋅v⋅Ye[2, 2]⎥\n",
       "⎢─────────────  ─────────────  ─────────────⎥\n",
       "⎣      2              2              2      ⎦"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "M_e = fermion_mass_matrix(LYuk, eLbar, eRc, model_leptons.vacuum, 3, (i, j),\n",
    "                          gamma=diracPR)\n",
    "M_e"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "36",
   "metadata": {},
   "source": [
    "### Fermion Feynman rules\n",
    "\n",
    "`model.physical_lagrangian(sector=...)` applies the vacuum shift, tadpole\n",
    "solution, and every registered rotation to *any* sector — including\n",
    "`\"yukawa\"` — so the h/G0-lepton couplings and the already-rotated gauge\n",
    "currents come out of a single `extract_fermion_vertices` pass. Note the\n",
    "Yukawa term carries **two** flavor indices (`Ye[i,j]` connects the\n",
    "doublet's flavor `i` to `eR`'s flavor `j`), while the gauge currents are\n",
    "flavor-diagonal (same index on both legs) — the lookup key must match."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "37",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "h  e_L e_R -> -sqrt(2)*I*Ye[i, j]*PR/2\n",
      "W+ nu  e   -> sqrt(2)*I*gw*gamma(mu)*PL/2\n",
      "Z  e_L e_L -> I*(g1**2 - gw**2)*gamma(mu)*PL/(2*sqrt(g1**2 + gw**2))\n",
      "Z  e_R e_R -> I*g1**2*gamma(mu)*PR/sqrt(g1**2 + gw**2)\n",
      "A  e_L e_L -> -I*g1*gw*gamma(mu)*PL/sqrt(g1**2 + gw**2)\n",
      "A  e_R e_R -> -I*g1*gw*gamma(mu)*PR/sqrt(g1**2 + gw**2)\n"
     ]
    }
   ],
   "source": [
    "Gm, cmap = conjugate_pair(Gp, \"Gm\")\n",
    "boson_fields = [h, G0, Gp, Gm, Z, A, Wp, Wm]\n",
    "\n",
    "LYuk_phys = model_leptons.physical_lagrangian(sector=\"yukawa\").xreplace(cmap)\n",
    "table = extract_fermion_vertices(LYuk_phys, boson_fields)\n",
    "\n",
    "mu = sp.Symbol(\"mu\", integer=True)\n",
    "gammaL = DiracGamma(mu) * diracPL\n",
    "gammaR = DiracGamma(mu) * diracPR\n",
    "nuL, eL = Ll.components\n",
    "\n",
    "def fermion_rule(bar, gamma, field, boson):\n",
    "    coeff = table[(bar, gamma, field)][1][(boson,)]\n",
    "    return sp.simplify(fermion_feynman_rule(coeff, gamma, (boson,)))\n",
    "\n",
    "print(\"h  e_L e_R ->\", fermion_rule(eLbar[i], diracPR, eRc[j], h))\n",
    "print(\"W+ nu  e   ->\", fermion_rule(nuLbar[i], gammaL, eL[i], Wp))\n",
    "print(\"Z  e_L e_L ->\", fermion_rule(eLbar[i], gammaL, eL[i], Z))\n",
    "print(\"Z  e_R e_R ->\", fermion_rule(eRbar[i], gammaR, eRc[i], Z))\n",
    "print(\"A  e_L e_L ->\", fermion_rule(eLbar[i], gammaL, eL[i], A))\n",
    "print(\"A  e_R e_R ->\", fermion_rule(eRbar[i], gammaR, eRc[i], A))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "38",
   "metadata": {},
   "source": [
    "The `h e_L e_R` rule, specialized to a single flavor with\n",
    "`Y_e = √2 m_ℓ/v` (the diagonal-mass convention), reduces to the pinned\n",
    "value `−i m_ℓ/v · P_R` (`+ h.c.`) — and the `W⁺νe`/`Z`/`γ` couplings above\n",
    "match `i g/√2 · γ^μP_L`, the `g_Z(T³−Qs²_W)` pattern, and `∓ie` respectively,\n",
    "all cross-checked in\n",
    "[`tests/test_fermion_sector.py`](../tests/test_fermion_sector.py) (which\n",
    "also covers `diagonalize_svd`/`diagonalize_takagi` for 3-generation mass\n",
    "matrices and a toy type-I seesaw — not repeated here)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "39",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAE8AAAAqCAYAAADyDQZvAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAEqUlEQVRoBe2a7VHcMBCGDUMBhFQQ6CAJFQAd8FFBQgdh+AX/GOiAUAEDHUAqCKED6CAJHZD30e1qZCNzd9i+4MSa0Wm1Wu2uX+1KPttzj4+PRRfl8PDwi/SuqN3tQv8kOs2HDcmum/y12gejF9Uuqd6o7knW+TY8vpnrELw7mb+UU3vj3ehWQj4QIbdqP1QtiXcl3kfVd6KnAnChqqytvhxZaUtXEz3y473NP6/Rcyo+kflZ9aRGJsuez3L/LWaasrkrWzbm29zgc7z/ATz2vEIReFsDRBjXWF1k1kwrilbTVg6yCRP+pOwP9b+6ZdFEACm0qvpJlc2aw+SnKhdwJZkTVWR2VCnQx+Kx0b+0YDc732wxzoFRB26t3bYjb19OsG9cqB5XrG7YGACfqW6qj9PIb6kCEsAvGZ+Dhv0IXS8q0gP4FA6FWMRfVsU/dG+Jnmqvc0WtRZ45+t0UA8a9G9EYq+sXwB5zX3GYKKQsip9GCXzAzhbJ+gKdio72EmHsUlYTWfrcCRDppTsB9UMUaoyWrOH0xT5+4HNJfkHMtgrKPfS3pfQoUXyjsQdVHAG86r2fR0hMc5vLrYXrNNaokS6AY7EYDwsyGin9+n7HYo4t0nmtii9EZslH9dmGWKTIby1tpZRVKtRuqgGkCISPiRciQf00usQOex73YUEHDCssQlW2kBxgkfaXatfVPpGx+djLgm/juYY9Oafvl/jcD8bSGnhRY1Hsi+bm2CMtGQoHQe5iuMiSw5oPj0Vg3yvUh/bC6gMcJeWPOPqVvEdzSW8UqCew61tMkJIueADHQRfLQqRaIOwCcdr3BoB0Ggs4ESMShs0hkkoOqw9ARCPbAePpXA6WNfH5C1jSp74X5ClVvSNu5tfssBikLTag/c7hTXVKq+BhFAMyzN4BzcYcivVx5txY3pCaYY4zrEXWo2ZX88MiJHqIBOw8VOZ5dwdC467D+c+1AB4PM83FB07krI355zRNOyZjpCQpy6qtq02jAjDTQ8XV42Aq53wOBE5foiukrg0Q2VxMTlcheTZ1oi2kLbQqOiYpHDARbM3DDrbz83kw0Kd6cHBwqnrRhc/S+1t1M9Wt/jH8lOd0q5E3ydI2kVEkkFZpFBbGa6I2zJUeMoMsiJFnSuFz0gZbJheGOnskFbR39KMLSNOIfS93gk9s3fSRsiwO/zaOxCNlC7Wk/5lq2KvVj/9GegkeF/UaSq/S9jUAlvowgJeiMSU9gDclYKn4AF6KxpT0nO5hOJ6/qdJOWngGVjrh1O/mNdykHv0FueG0bQD6kLYDeA0QaDB1iLwBvAYINJg6RN4AXgMEGkwdIm8ArwECDaYOkdcAvLZfADVwZfxU/QXkLyTvR/ii6U79+O5DNC+b+HzDX0mq223pW+T5tzB8KcALorQsqRPemKXMLunegKeI4nG4fwvDI/Pqtym8kK7yusSu6FPapq8at4XKUQWZmQKH7T5Fnr+QefItjIHIe+KZ7Xe9As8AogkpK6ACmAl/6s9ik7kvInsTecnVLYuuPoglGkvvcxP5zsg+ggdwABiKIhCaD3Nmvuf18kmygCLKuDUJX0CpH+/3AqIz+vkDefhdJEgoB7wAAAAASUVORK5CYII=",
      "text/latex": [
       "$\\displaystyle - \\frac{i m_{\\ell} P_R}{v}$"
      ],
      "text/plain": [
       "-ⅈ⋅mₑₗₗ⋅PR \n",
       "───────────\n",
       "     v     "
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "m_ell = sp.Symbol(\"m_ell\", positive=True)\n",
    "hll_coeff = table[(eLbar[i], diracPR, eRc[j])][1][(h,)]\n",
    "hll_rule = fermion_feynman_rule(hll_coeff, diracPR, (h,))\n",
    "sp.simplify(hll_rule.subs(Ye[i, j], sp.sqrt(2) * m_ell / v.s))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "40",
   "metadata": {},
   "source": [
    "## 10. Output: LaTeX tables and UFO export\n",
    "\n",
    "`latex_feynman_table` produces a ready-to-paste `\\begin{array}` table."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "41",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\\begin{array}{|c|c|}\n",
      "\\hline\n",
      "\\textbf{Interaction} & \\textbf{Coefficient} \\\\\n",
      "\\hline\n",
      "$A A Gm Gp$ & $\\frac{2 i g_{1}^{2} gw^{2}}{g_{1}^{2} + gw^{2}}$ \\\\\n",
      "\\hline\n",
      "$A Gm Gp Z$ & $\\frac{i g_{1} gw \\left(- g_{1}^{2} + gw^{2}\\right)}{g_{1}^{2} + gw^{2}}$ \\\\\n",
      "\\hline\n",
      "$A Gm H_{0 i} Wp$ & $- \\frac{g_{1} gw^{2}}{2 \\sqrt{g_{1}^{2} + gw^{2}}}$ \\\\\n",
      "\\hline\n",
      "$A Gm H_{0 r} Wp$ & $\\frac{i g_{1} gw^{2}}{2 \\sqrt{g_{1}^{2} + gw^{2}}}$ \\\\\n",
      "\\hline\n",
      "\\end{array}\n"
     ]
    }
   ],
   "source": [
    "preview = {k: rules[k] for k in list(rules)[:4]}\n",
    "print(latex_feynman_table(preview))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "42",
   "metadata": {},
   "source": [
    "For a MadGraph-importable model directory, use\n",
    "`feynlag.export.ufo.write_ufo(...)` — it needs a `UFOParticle` spec per\n",
    "physical field (PDG code, spin, mass parameter name, color, charge). See\n",
    "[`tests/test_ufo_export.py`](../tests/test_ufo_export.py) for a complete\n",
    "worked call, including a numeric round-trip check that the exported `hWW`\n",
    "coupling evaluates to `i g² v/2`."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43",
   "metadata": {},
   "source": [
    "## 11. Extending the SM: writing a BSM model\n",
    "\n",
    "Every BSM model follows exactly the pipeline above; the only things that\n",
    "change are *what you declare in §1–§3*. Worked templates already live in\n",
    "[`examples/thdm.py`](thdm.py) (a second Higgs doublet, softly-broken `Z₂`)\n",
    "and [`examples/thdm_s3.py`](thdm_s3.py) (a discrete `S₃` flavor symmetry\n",
    "across three doublets) — read them alongside this checklist:\n",
    "\n",
    "1. **New fields**: declare extra `Scalar`/`Fermion`/`GaugeBoson`s with their\n",
    "   representations under your (possibly extended) symmetry groups. Extra\n",
    "   gauge groups are just more `SU2`/`U1`/`SU3` objects passed to\n",
    "   `gauge_groups=[...]`; a gauged extra `U(1)'` is created exactly like\n",
    "   `U1Y` above.\n",
    "2. **New symmetries**: for a discrete symmetry (softly-broken `Z₂` picking\n",
    "   out Yukawa alignment in a 2HDM, or an `S₃` flavor symmetry across three\n",
    "   doublets), instantiate `ZN(name, N)` or `S3()`, then\n",
    "   `group.assign(irrep, *fields)`. `S3.doublet_product(...)` gives you the\n",
    "   `2⊗2 = 1⊕1'⊕2` Clebsch–Gordan contractions to build invariant terms\n",
    "   without doing the group theory by hand — `thdm_s3.py` builds its entire\n",
    "   potential this way.\n",
    "3. **Extend the potential/Yukawas**: write the new invariant contractions\n",
    "   (`dag(H1)*H2.mat`, CG products, …) and add them to a `Lagrangian` sector.\n",
    "   Call `model.check_invariance()` immediately — a forbidden term (wrong\n",
    "   representation, missing `+ h.c.`, wrong discrete charge) fails loudly\n",
    "   instead of silently producing wrong Feynman rules downstream.\n",
    "4. **VEVs on every scalar that gets one**, then `solve_tadpoles([...])` for\n",
    "   as many mass parameters as you have independent tadpole conditions. If\n",
    "   the system is *over*-constrained (more tadpole equations than free mass\n",
    "   parameters — the `thdm_s3.py` case), that's physical: it's forcing a\n",
    "   vacuum alignment, visible as a residual condition on the VEV ratios after\n",
    "   solving the others (exactly how `thdm_s3.py` recovers the `√3`\n",
    "   alignment).\n",
    "5. **Mass matrices grow with the field content**: `model.mass_matrix([...])`\n",
    "   accepts any list of same-block fluctuation symbols — a 3×3 CP-even block\n",
    "   for a 3HDM works identically to the SM's 1×1 above. For blocks too large\n",
    "   to diagonalize in closed form (≥3×3 symbolic), don't call\n",
    "   `sp.Matrix.eigenvects` directly; either parametrize a rotation and solve\n",
    "   the off-diagonal conditions numerically at benchmark points, or export to\n",
    "   UFO and diagonalize numerically downstream.\n",
    "6. **Register every new rotation** the same way as `alpha`/`beta` in\n",
    "   `thdm.py`'s CP-even block.\n",
    "7. **Extract vertices exactly as in §7**, just with a longer `fields` list.\n",
    "   Chiral fermions use a parallel \"bilinear\" track\n",
    "   (`feynlag.Bilinear`, `fermion_gauge_current`,\n",
    "   `extract_fermion_vertices`) since a fermion sandwich is a c-number atom,\n",
    "   not a product of commuting symbols — §9 above works through the full\n",
    "   SM-lepton case (Yukawa, gauge currents, mass matrix, Feynman rules,\n",
    "   including how `check_invariance()` now covers fermion gauge invariance\n",
    "   and mass dimension); see\n",
    "   [`tests/test_fermion_sector.py`](../tests/test_fermion_sector.py) for the\n",
    "   3-generation `diagonalize_svd`/`diagonalize_takagi` and seesaw pieces not\n",
    "   repeated there.\n",
    "\n",
    "The invariant is: nothing about `Model`, `check_invariance`,\n",
    "`solve_tadpoles`, `mass_matrix`, `rotate`, or `feynman_rules` is\n",
    "SM-specific — they operate purely on whatever fields, symmetries, and\n",
    "Lagrangian terms you declared. Growing the model is additive, not a\n",
    "rewrite. As a taste of step 1–2, here is the *declaration* half of a second\n",
    "Higgs doublet (the full potential and pipeline is in `thdm.py`):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "44",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJkAAAAVCAYAAABVLkwZAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAF0klEQVRoBe2a7XEUORCG15QD8EEGJgNzRABkAEUEQAZc8cv+54IMDBkAGXAXAQcZQAYYMvA9jzw9pdHOrGe1WrNV566S9dXqD+lVSxrv3sXFxWJX6OTk5NqMQdfervj9O+y4zrne2yWQ/Y7JvtG5/Rm4tX0VNxr+7zOwX04AYfQlbfdJj7u+v8l/0f7EOvlzskekh6QDUuonf0bfL/JmhLydsaWFUy38QcYhtvzV2fMn+bl12r92bRtlyHFNY42VFWv6ib63nf4X5GHDlfomj0uEeD/6Sn5vTArtP2k/J7871t+ybZdsaeFXrT+ME2Bn5G7yRJRfU3AzPqLshq8mxivnFemUJKACYAvKBp0ILq+pv6U+i0aPSwQcdaNHjaZfZ0X8aP8szTOZdsEWbHhIcndvTBv6I6Be5EYgz4giGD7k7euWkfOJMcq+R/kNqQeYsqh/NIdc+7XWfRRkCHmqNEjFY+RRKU31X/a2+bsLtrihTC1oE3+c928seGmLi35AuwBYmxjnOipbgH1fIUAge3VaxbM0fApkCUQIm0Jsup8hbap/SdEGDbtkywZu9EM38cf5/s66DKJML7liIyDL6KhN3rOm5IYK73/vozI3359g9LhcBSAvnKucnRBb1bxLtlQ5UAyq9gcQxOYuRC7S9Yb+dPknty5wfMA9IxnhrN/pyumRBp/t3sNcyzl3LCPYGSnRXD1LIGOgxkiHlKfOecP12ohOUtf4s0u2rGH2JOs2/EGmgBIs3s2CntLui9P1+4d0SvmNneRGLtu9xMf9rgcObZPEWCNd/oqdpWcJZAiJl8sThOYCk3LaRL4vjcF9jHaB58tE0mlp8FkDHh2UfpDuUM8nJnUUf6psUUaFrkL1VqrV/qywRsB8xN8AkaD73PG7DkapuLTb/I3kGkoRUFadWpecxV9kztYzBrKkGCFLAOv0xESVhvmsjZ2xoOzu+EJKnzio+zr7QR6T4YvNby8hrxM/yKpsqdS1YJw2x8Tnhty2Qn/vX9bpZ56pYyxjS8Uqf0ohUe/sFUS5fuuxdgLhNPi73GCwgMc8BYOMv2O5zGgXjOUnKi/+BofZevYHUi8rGlYCKGdzopxYQ2dOz2n7QIqxRi3bjkg6bf1BDJCPJMg8lqdeK7W21OhaYMcYiGw3cmtn2iDhQ0Ve68+SKmxx094mH2xS6mldyGOzxHqEDPkFiumccgJddOY5/REQ3Hzq83GQ7m6Ol5f8Sj23ZAzKBgyOwqzfSZJKw21zgf61UBJy3TE6U4JJQ0PmYFitLTW6Boq3VKn1Z8wcZAn6u+R9BNNvU8YvmIw2CQxZu6B439VTxCvGZax9MSLe2OPgSj0DkCEydsUYiNQYqF0CIYYOvhDDK+gipAqwMXIn3R7roK3WlhpdEyY0ba71Z2AE8+ymvE9e3mcFnvMZ5FoN1pExRiNBF2Mjd+woMUaALcnKmJf6Sj37GbPFpAymhPCiz2qaKPoHxpd89DsRyhr9l1TGL8CmQNHElpm6MratFTf2h3l1wb3oe9XwCMvJO24cb86pa9CDrhsrqB5QFmgLcoOA0dBrjsfnIFJRV8YrkifQWGCZpWcfQTK+I5lHWNSRz/SF0ToUiF7QrkLvU6mfck+0yeedyK/HyRnKkfd8XUGd/REKv/VNbZmlqzRkG/VG/uSmOe/OrxGppDwwGF0kP2PES9JvZP5/s59vGaj7Mv2Dog835Tt/weN/F/zKoM4esJSD5unx92St0vHx8SHpLOR19UPrlH+SjqKva7uQJ29rUW6tC3mPSS9b2HYdMrD1jPRl27rm6imPy0Do2nmHdiOeu8cwK3kvi3P/lLLITzsOHsvuotg1VJtRa13u7qkI2czohoLS3DaUNyVqlp7Jn/pMSZ1qByz+9Oeg7Ke9/5kzZY9RyY+xfn8RkFtZvOvUpUO7QvjtGrgWHnP5R9imJq6jpxnImnpwI6xqBlh4Hxde1D1JPDH87dngMk/bxrSunv8Al0Ki9hzV5nAAAAAASUVORK5CYII=",
      "text/latex": [
       "$\\displaystyle H_{20} \\overline{H_{0}} + H2p \\overline{Gp}$"
      ],
      "text/plain": [
       "    __       __\n",
       "H₂₀⋅H₀ + H2p⋅Gp"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "H2 = Scalar(\"H2\", reps={SU2L: 2, U1Y: sp.Rational(1, 2)},\n",
    "            component_names=[\"H2p\", \"H20\"])\n",
    "v2 = ExternalParameter(\"v2\", 246.0 * 0.31, positive=True, unit_dim=1)\n",
    "H2.expand_vev({H2.components[1]: v2})\n",
    "\n",
    "# a new SU(2)xU(1)-invariant mixed bilinear, available for the potential:\n",
    "H1dH2 = (dag(H) * H2.mat)[0]\n",
    "H1dH2"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.13.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}