viewof K      = Inputs.range([40, 200],   {value: 160, step: 10,   label: katexLabel("K", " (Home capital)")})
viewof L      = Inputs.range([40, 200],   {value: 80,  step: 10,   label: katexLabel("L", " (Home labor)")})
viewof Kf     = Inputs.range([40, 200],   {value: 100, step: 10,   label: katexLabel("K^{*}", " (Foreign capital)")})
viewof Lf     = Inputs.range([40, 200],   {value: 140, step: 10,   label: katexLabel("L^{*}", " (Foreign labor)")})
viewof alphaC = Inputs.range([0.1, 0.5],  {value: 0.15, step: 0.05, label: katexLabel("\\alpha_{C}", " (cloth K share)")})
viewof alphaW = Inputs.range([0.5, 0.9],  {value: 0.80, step: 0.05, label: katexLabel("\\alpha_{W}", " (wine K share)")})
viewof beta   = Inputs.range([0.1, 0.9],  {value: 0.5, step: 0.05, label: katexLabel("\\beta", " (cloth expenditure share)")})
viewof tariff = Inputs.range([0, 0.5],    {value: 0.3, step: 0.05, label: katexLabel("t", " (tariff rate)")})
viewof lambda = Inputs.range([1, 20],     {value: 2,  step: 0.5,   label: katexLabel("\\lambda", " (Foreign size)")})
viewof stage  = Inputs.range([1, 5],      {value: 5, step: 1,      label: "reveal stage"})

katexReady = {
  if (window.katex && typeof window.katex.renderToString === "function") {
    return window.katex;
  }
  return new Promise((resolve, reject) => {
    const check = () => {
      if (window.katex && typeof window.katex.renderToString === "function") {
        resolve(window.katex);
      } else {
        setTimeout(check, 30);
      }
    };
    check();
    setTimeout(() => reject(new Error("katex failed to load")), 8000);
  });
}

// KaTeX-rendered HTML label for sliders. Must be a block cell — Quarto OJS
// does not scan function bodies for reactive deps, so a bare `function`
// declaration referencing `katexReady` errors at call time.
katexLabel = {
  const ready = katexReady;
  return function (mathStr, trailing = "") {
    const span = document.createElement("span");
    const mathSpan = document.createElement("span");
    mathSpan.innerHTML = ready.renderToString(mathStr, { throwOnError: false });
    span.appendChild(mathSpan);
    if (trailing) span.appendChild(document.createTextNode(trailing));
    return span;
  };
}

// ── SVG text label with math-style sub/superscripts and italic math runs ──
// `_X` / `_{XX}` → subscript, `^X` / `^{XX}` → superscript, `$...$` → italic
// (math) run. `math: true` option makes the entire string italic by default.
function texLabel(parent, px, py, str, {
  color = "#000", fontSize = 12, anchor = "middle", bold = false, rotated = false, math = false
} = {}) {
  const textEl = parent.append("text")
    .attr("fill", color)
    .attr("font-size", fontSize)
    .attr("font-family", "serif")
    .attr("font-weight", bold ? "bold" : "normal")
    .attr("text-anchor",
      anchor === "start" ? "start" :
      anchor === "end" ? "end" : "middle")
    .attr("dominant-baseline", "middle")
    .attr("transform",
      rotated
        ? `translate(${px},${py}) rotate(-90)`
        : `translate(${px},${py})`);

  let i = 0;
  let inMath = math;
  while (i < str.length) {
    const c = str[i];
    if (c === "$") { inMath = !inMath; i++; continue; }
    if (c === "_" || c === "^") {
      i++;
      let content = "";
      if (str[i] === "{") {
        i++;
        while (i < str.length && str[i] !== "}") content += str[i++];
        i++;
      } else {
        content = str[i++] || "";
      }
      const tspan = textEl.append("tspan")
        .attr("baseline-shift", c === "_" ? "sub" : "super")
        .attr("font-size", `${fontSize * 0.72}px`)
        .text(content);
      if (inMath) tspan.attr("font-style", "italic");
    } else {
      let run = "";
      while (i < str.length && str[i] !== "_" && str[i] !== "^" && str[i] !== "$") {
        run += str[i++];
      }
      if (run) {
        const tspan = textEl.append("tspan").text(run);
        if (inMath) tspan.attr("font-style", "italic");
      }
    }
  }
  return textEl;
}

function svgAxisLabel(parent, px, py, text, {rotated = false, fontSize = 12} = {}) {
  return parent.append("text")
    .attr("text-anchor", "middle")
    .attr("dominant-baseline", "middle")
    .attr("font-size", fontSize)
    .attr("font-family", "serif")
    .attr("fill", "#000")
    .attr("transform",
      rotated
        ? `translate(${px},${py}) rotate(-90)`
        : `translate(${px},${py})`)
    .text(text);
}

// ── Collision-avoiding label placer (for future charts that need it) ──
function makeLabelPlacer({g, x, y, iw, ih}) {
  const placedBBs = [];
  function approxBB(str, px, py, anchor, fs = 11) {
    const hasSub = /_/.test(str);
    const w = str.length * (fs * 0.62);
    const h = hasSub ? fs * 2.3 : fs * 1.5;
    let x0 = px;
    if (anchor === "middle") x0 = px - w / 2;
    else if (anchor === "end") x0 = px - w;
    return {x0, y0: py - h / 2, x1: x0 + w, y1: py + h / 2};
  }
  function segHitsBB(ax, ay, bx, by, bb) {
    if ((ax < bb.x0 && bx < bb.x0) || (ax > bb.x1 && bx > bb.x1) ||
        (ay < bb.y0 && by < bb.y0) || (ay > bb.y1 && by > bb.y1)) return false;
    if (ax >= bb.x0 && ax <= bb.x1 && ay >= bb.y0 && ay <= bb.y1) return true;
    if (bx >= bb.x0 && bx <= bb.x1 && by >= bb.y0 && by <= bb.y1) return true;
    function si(p0x, p0y, p1x, p1y, p2x, p2y, p3x, p3y) {
      const s1x = p1x - p0x, s1y = p1y - p0y;
      const s2x = p3x - p2x, s2y = p3y - p2y;
      const denom = -s2x * s1y + s1x * s2y;
      if (denom === 0) return false;
      const s = (-s1y * (p0x - p2x) + s1x * (p0y - p2y)) / denom;
      const t = ( s2x * (p0y - p2y) - s2y * (p0x - p2x)) / denom;
      return s >= 0 && s <= 1 && t >= 0 && t <= 1;
    }
    return si(ax, ay, bx, by, bb.x0, bb.y0, bb.x1, bb.y0) ||
           si(ax, ay, bx, by, bb.x1, bb.y0, bb.x1, bb.y1) ||
           si(ax, ay, bx, by, bb.x1, bb.y1, bb.x0, bb.y1) ||
           si(ax, ay, bx, by, bb.x0, bb.y1, bb.x0, bb.y0);
  }
  function bbHitsCurve(bb, pts) {
    for (let k = 0; k < pts.length - 1; k++) {
      const [cx1, cy1] = pts[k];
      const [cx2, cy2] = pts[k + 1];
      if (segHitsBB(x(cx1), y(cy1), x(cx2), y(cy2), bb)) return true;
    }
    return false;
  }
  function bbHitsBB(a, b) {
    return !(a.x1 < b.x0 || a.x0 > b.x1 || a.y1 < b.y0 || a.y0 > b.y1);
  }
  function inPlot(bb) {
    return bb.x0 >= 2 && bb.y0 >= 2 && bb.x1 <= iw - 2 && bb.y1 <= ih - 2;
  }
  function placeLabel(str, candidates, avoid, color, bold, texStr) {
    const labelTeX = texStr || str;
    const pad = 3;
    for (const c of candidates) {
      const px = x(c.x), py = y(c.y);
      const raw = approxBB(str, px, py, c.anchor || "start");
      const bb = {x0: raw.x0 - pad, y0: raw.y0 - pad, x1: raw.x1 + pad, y1: raw.y1 + pad};
      if (!inPlot(bb)) continue;
      let hit = false;
      for (const curve of avoid) if (bbHitsCurve(bb, curve)) { hit = true; break; }
      if (hit) continue;
      for (const o of placedBBs) if (bbHitsBB(bb, o)) { hit = true; break; }
      if (hit) continue;
      placedBBs.push(bb);
      return texLabel(g, px, py, labelTeX, {color, fontSize: 11, anchor: c.anchor || "start", bold, math: true});
    }
    const c = candidates[0];
    return texLabel(g, x(c.x), y(c.y), labelTeX, {color, fontSize: 11, anchor: c.anchor || "start", bold, math: true});
  }
  return {placedBBs, approxBB, segHitsBB, bbHitsCurve, bbHitsBB, inPlot, placeLabel};
}

// ── Synced slider panel (subset selectable per slide) ──
function sliderPanel(keys = ["tariff", "stage"], opts = {}) {
  const all = {
    K:      Inputs.bind(Inputs.range([40, 200],   {value: 160, step: 10,   label: katexLabel("K", " (Home capital)")}), viewof K),
    L:      Inputs.bind(Inputs.range([40, 200],   {value: 80,  step: 10,   label: katexLabel("L", " (Home labor)")}), viewof L),
    Kf:     Inputs.bind(Inputs.range([40, 200],   {value: 100, step: 10,   label: katexLabel("K^{*}", " (Foreign capital)")}), viewof Kf),
    Lf:     Inputs.bind(Inputs.range([40, 200],   {value: 140, step: 10,   label: katexLabel("L^{*}", " (Foreign labor)")}), viewof Lf),
    alphaC: Inputs.bind(Inputs.range([0.1, 0.5],  {value: 0.15, step: 0.05, label: katexLabel("\\alpha_{C}", " (cloth K share)")}), viewof alphaC),
    alphaW: Inputs.bind(Inputs.range([0.5, 0.9],  {value: 0.80, step: 0.05, label: katexLabel("\\alpha_{W}", " (wine K share)")}), viewof alphaW),
    beta:   Inputs.bind(Inputs.range([0.1, 0.9],  {value: 0.5, step: 0.05, label: katexLabel("\\beta", " (cloth expenditure share)")}), viewof beta),
    tariff: Inputs.bind(Inputs.range([0, 0.5],    {value: 0.3, step: 0.05, label: katexLabel("t", " (tariff rate)")}), viewof tariff),
    lambda: Inputs.bind(Inputs.range([1, 20],     {value: 2,  step: 0.5,   label: katexLabel("\\lambda", " (Foreign size)")}), viewof lambda),
    stage:  Inputs.bind(Inputs.range([1, 5],      {value: 5, step: 1,      label: "reveal stage"}), viewof stage)
  };
  const div = document.createElement("div");
  div.style.fontSize = "13px";
  if (opts.columns && opts.columns > 1) {
    div.style.display = "grid";
    div.style.gridTemplateColumns = `repeat(${opts.columns}, 1fr)`;
    div.style.columnGap = "16px";
    div.style.rowGap = "2px";
  }
  for (const k of keys) all[k].style.maxWidth = "100%";
  for (const k of keys) div.appendChild(all[k]);
  return div;
}

ratio = ((1 - alphaC) / alphaC) * (alphaW / (1 - alphaW))

// Home autarky allocation (closed form: CD utility × CD production)
etaAut = (beta / (1 - beta)) * ((1 - alphaC) / (1 - alphaW))
lcAut = (etaAut * L) / (1 + etaAut)
kcAut = K * lcAut / (lcAut * (1 - ratio) + ratio * L)
qcAut = Math.pow(kcAut, alphaC) * Math.pow(lcAut, 1 - alphaC)
qwAut = Math.pow(K - kcAut, alphaW) * Math.pow(L - lcAut, 1 - alphaW)
priceAut = (beta * qwAut) / ((1 - beta) * qcAut)

// Foreign autarky — invariant to a common scaling of (Kf, Lf), so λ does not
// move it. Only the K/L ratio matters for Foreign's autarky price.
etaAutF   = (beta / (1 - beta)) * ((1 - alphaC) / (1 - alphaW))
lcAutF    = (etaAutF * Lf) / (1 + etaAutF)
kcAutF    = Kf * lcAutF / (lcAutF * (1 - ratio) + ratio * Lf)
qcAutF    = Math.pow(kcAutF, alphaC) * Math.pow(lcAutF, 1 - alphaC)
qwAutF    = Math.pow(Kf - kcAutF, alphaW) * Math.pow(Lf - lcAutF, 1 - alphaW)
priceAutF = (beta * qwAutF) / ((1 - beta) * qcAutF)

// Production at relative price p, handling diversification cone + corners.
// Convention: p = p_C / p_W with p_W ≡ 1.
function produce(p, Kend, Lend) {
  const logBC = -alphaC * Math.log(alphaC) - (1 - alphaC) * Math.log(1 - alphaC);
  const logBW = -alphaW * Math.log(alphaW) - (1 - alphaW) * Math.log(1 - alphaW);
  const det   = alphaW - alphaC;
  const rhsC  = Math.log(p) - logBC;
  const rhsW  = -logBW;
  const w = Math.exp((rhsC * alphaW - rhsW * alphaC) / det);
  const r = Math.exp(((1 - alphaC) * rhsW - (1 - alphaW) * rhsC) / det);
  const klC = (alphaC / (1 - alphaC)) * (w / r);
  const klW = (alphaW / (1 - alphaW)) * (w / r);
  const lc = (Kend - klW * Lend) / (klC - klW);
  if (lc < 0) {
    const qw = Math.pow(Kend, alphaW) * Math.pow(Lend, 1 - alphaW);
    return {qc: 0, qw, lc: 0, kc: 0, lw: Lend, kw: Kend, w, r, specialized: "wine"};
  }
  if (lc > Lend) {
    const qc = Math.pow(Kend, alphaC) * Math.pow(Lend, 1 - alphaC);
    return {qc, qw: 0, lc: Lend, kc: Kend, lw: 0, kw: 0, w, r, specialized: "cloth"};
  }
  const kc = klC * lc;
  const lw = Lend - lc;
  const kw = Kend - kc;
  const qc = Math.pow(kc, alphaC) * Math.pow(lc, 1 - alphaC);
  const qw = Math.pow(kw, alphaW) * Math.pow(lw, 1 - alphaW);
  return {qc, qw, lc, kc, lw, kw, w, r, specialized: null};
}

// ── SOE tariff setup ──
// World price pinned at Foreign's autarky price (λ-large limit).
pWorld = priceAutF

// Home domestic price under tariff t (SOE ⇒ full pass-through).
pHome = pWorld + tariff

// Reference income pins the demand curve in (Q, P) space so that welfare
// triangles have a stable frame. Use Home's GDP at free-trade world prices
// (standard Marshallian partial-eq convention).
prodHFree0 = produce(pWorld, K, L)
incomeRef = pWorld * prodHFree0.qc + prodHFree0.qw

// Home cloth supply and demand curves over a fine p grid.
//   S_H(p) = Q_C*(p) from the PPF (reduced form of CD production + PPF)
//   D_H(p) = β · I_ref / p  (CD preferences at fixed reference income)
cloth = {
  const pMin = 0.1;
  const pMax = Math.max(priceAut * 1.15, pWorld * 3, 2.5);
  const N = 200;
  const supply = [];
  const demand = [];
  for (let i = 0; i <= N; i++) {
    const p = pMin + (i / N) * (pMax - pMin);
    const pr = produce(p, K, L);
    supply.push({p, q: pr.qc});
    demand.push({p, q: beta * incomeRef / p});
  }
  return {supply, demand, pMin, pMax};
}

// Key equilibrium quantities at free trade (t=0) and under current tariff t.
qSFree = produce(pWorld, K, L).qc
qSTar  = produce(pHome,  K, L).qc
qDFree = beta * incomeRef / pWorld
qDTar  = beta * incomeRef / pHome
importsFree = qDFree - qSFree
importsTar  = qDTar  - qSTar
revenueTar  = tariff * importsTar

// Home production bundles at free trade and under tariff (needed for the PPF
// decomposition of the PS change).
prodHTar = produce(pHome, K, L)

// ── Wine-axis intercepts on the PPF (Graph 3 redesign) ──
//   I1      : line through A at slope -p^*       = q_W^A + p^* · q_C^A   (= incomeRef)
//   I2      : line through A at slope -(p^*+t)   = q_W^A + (p^*+t) · q_C^A
//   I3      : line through B at slope -(p^*+t)   = q_W^B + (p^*+t) · q_C^B   (BL3 — tangent at B)
//   I3tilde : line through B at slope -p^*       = q_W^B + p^* · q_C^B   (parallel to BL1, through B)
qcA = prodHFree0.qc
qwA = prodHFree0.qw
qcB = prodHTar.qc
qwB = prodHTar.qw

I1      = qwA + pWorld * qcA
I2      = qwA + pHome  * qcA
I3      = qwB + pHome  * qcB
I3tilde = qwB + pWorld * qcB

// ── Welfare decomposition (Marshallian — matches Graph 2's a/b/c/d) ──
// Demand is reduced-form Marshallian at fixed reference income I1, so the four
// pieces sum to Graph 2's CS-loss trapezoid by construction.
//   a = ΔPS gain                        (= I3 - I1; equals ∫_{p^*}^{p^*+t} S_H(p) dp)
//   b = production-distortion DWL       (= I1 - I3tilde; Harberger's left triangle)
//   c = tariff revenue                  (= t · imports under Marshallian demand)
//   d = consumption-distortion DWL      (= total CS loss − a − b − c)
//   csLossTotal = β · I1 · log(p_home / p_world)  for CD demand
aPSgain     = I3 - I1
bProdDWL    = I1 - I3tilde
cRevenue    = tariff * (qDTar - qcB)
csLossTotal = beta * I1 * Math.log(pHome / pWorld)
dConsDWL    = csLossTotal - aPSgain - bProdDWL - cRevenue

// ── Marshallian consumption point D (used by Graph 3's PPF panel) ──
// Lives off the PPF in the cloth direction (Home imports cloth: q_C^D > q_C^B).
// Wine residual derived from CD's fixed-income budget at consumer prices.
qcD = qDTar                       // β I1 / pHome (already named)
qwD = (1 - beta) * I1             // wine residual under fixed reference income I1

// Home PPF as a parametric sweep in L_C (same construction as week 4).
ppfPoints = {
  const pts = [];
  const N = 200;
  for (let i = 0; i <= N; i++) {
    const lc = (i / N) * L;
    const kc = K * lc / (lc * (1 - ratio) + ratio * L);
    const lw = L - lc;
    const kw = K - kc;
    const qc = Math.pow(kc, alphaC) * Math.pow(lc, 1 - alphaC);
    const qw = Math.pow(kw, alphaW) * Math.pow(lw, 1 - alphaW);
    pts.push([qc, qw]);
  }
  return pts;
}

KfL = Kf * lambda
LfL = Lf * lambda

// Free-trade world equilibrium price. Bisect on cloth-market clearing.
// CD utility ⇒ closed-form income; world demand = β·(GDP_W)/p.
ftWorldPrice = {
  const f = (p) => {
    const prH = produce(p, K, L);
    const prF = produce(p, KfL, LfL);
    const sCloth = prH.qc + prF.qc;
    const Ih = p * prH.qc + prH.qw;
    const If = p * prF.qc + prF.qw;
    const dCloth = beta * (Ih + If) / p;
    return sCloth - dCloth;
  };
  let lo = 0.05, hi = 5.0;
  for (let i = 0; i < 60; i++) {
    const mid = 0.5 * (lo + hi);
    if (f(mid) > 0) hi = mid; else lo = mid;
  }
  return 0.5 * (lo + hi);
}

// Reference incomes (pinned at FT levels — Marshallian convention).
prHFT = produce(ftWorldPrice, K, L)
prFFT = produce(ftWorldPrice, KfL, LfL)
IhFT  = ftWorldPrice * prHFT.qc + prHFT.qw
IfFT  = ftWorldPrice * prFFT.qc + prFFT.qw
IwFT  = IhFT + IfFT

// Equilibrium border price under tariff. Home producers face p+t, Foreign p.
// Demand uses fixed FT reference incomes (no rebate — keeps the chart Marshallian).
tarWorldPrice = {
  const t = tariff;
  const f = (p) => {
    const pH = p + t;
    const prH = produce(pH, K, L);
    const prF = produce(p,  KfL, LfL);
    const sCloth = prH.qc + prF.qc;
    const dCloth = beta * IfFT / p + beta * IhFT / pH;
    return sCloth - dCloth;
  };
  let lo = 0.05, hi = 5.0;
  for (let i = 0; i < 60; i++) {
    const mid = 0.5 * (lo + hi);
    if (f(mid) > 0) hi = mid; else lo = mid;
  }
  return 0.5 * (lo + hi);
}

// World RS curve (cloth/wine ratio) parameterised by border price p*.
// FT: both countries at p*. Tariff: Home at p*+t, Foreign at p*.
// `truncateAt`: stop the curve at this p* (used to cut tariff curves at the
// prohibitive boundary p*_prohib = p^H_aut − t, beyond which Home wouldn't
// import — tariff is irrelevant and the "Home at p*+t" assumption breaks).
function worldRSPoints(t, pMin, pMax, N = 100, truncateAt = null) {
  const pts = [];
  for (let i = 0; i <= N; i++) {
    const p = pMin + (i / N) * (pMax - pMin);
    if (truncateAt !== null && p > truncateAt) {
      pts.push({p, q: NaN});
      continue;
    }
    const prH = produce(p + t, K, L);
    const prF = produce(p,     KfL, LfL);
    const cloth = prH.qc + prF.qc;
    const wine  = prH.qw + prF.qw;
    pts.push({p, q: wine > 0 ? cloth / wine : NaN});
  }
  return pts;
}

// World RD curve. FT: β/((1-β)p) (homothetic, income drops out).
// Tariff: weighted by where each country's income is spent.
function worldRDPoints(t, pMin, pMax, N = 100, truncateAt = null) {
  const pts = [];
  for (let i = 0; i <= N; i++) {
    const p = pMin + (i / N) * (pMax - pMin);
    if (truncateAt !== null && p > truncateAt) {
      pts.push({p, q: NaN});
      continue;
    }
    const cloth = beta * IfFT / p + beta * IhFT / (p + t);
    const wine  = (1 - beta) * IwFT;
    pts.push({p, q: cloth / wine});
  }
  return pts;
}

// Foreign cloth exports as a function of border price p (= Foreign's price).
//   X_F(p) = Q_C^F(p)  −  β I^F(p) / p,  with I^F = p·Q_C^F + Q_W^F (CD)
// Positive when p > Foreign's autarky price (Foreign-side comparative advantage
// in cloth ⇒ Foreign is the cloth exporter).
function foreignExports(p) {
  const prF = produce(p, KfL, LfL);
  const If  = p * prF.qc + prF.qw;
  const dF  = beta * If / p;
  return prF.qc - dF;
}

// X_F curve sampled over a price band (only plotted where exports > 0).
xfCurve = {
  const pts = [];
  const N = 200;
  const pMin = 0.05, pMax = ftWorldPrice * 2.0;
  for (let i = 0; i <= N; i++) {
    const p = pMin + (i / N) * (pMax - pMin);
    const q = foreignExports(p);
    pts.push({p, q: q > 0 ? q : NaN});
  }
  return pts;
}

// ── Home welfare change as a function of tariff t (Marshallian, matches
// slide-1/2 convention with fixed reference incomes IhFT, IfFT) ──
// ΔW(t) = ΔPS + G − ΔCS  where
//   ΔCS(t) = β·IhFT·log(p_H_t / p^*_FT)              (CD demand at I_ref)
//   ΔPS(t) = ∫_{p^*_FT}^{p_H_t} S_H(p) dp           (trapezoidal, 30 pts)
//   G(t)   = t · M_t   with M_t = D_H(p_H_t) − S_H(p_H_t)
function homeWelfareChange(t) {
  // 1. Solve tariff border price (same bisection as tarWorldPrice).
  const f = (p) => {
    const pH = p + t;
    const prH = produce(pH, K, L);
    const prF = produce(p,  KfL, LfL);
    const sCloth = prH.qc + prF.qc;
    const dCloth = beta * IfFT / p + beta * IhFT / pH;
    return sCloth - dCloth;
  };
  let lo = 0.05, hi = 5.0;
  for (let i = 0; i < 60; i++) {
    const mid = 0.5 * (lo + hi);
    if (f(mid) > 0) hi = mid; else lo = mid;
  }
  const pStar = 0.5 * (lo + hi);
  const pH = pStar + t;
  const prH = produce(pH, K, L);

  // ΔCS, ΔPS, G
  const dCS = beta * IhFT * Math.log(pH / ftWorldPrice);   // loss (positive)
  let dPS = 0;
  const N = 30;
  for (let i = 0; i < N; i++) {
    const p1 = ftWorldPrice + (i     / N) * (pH - ftWorldPrice);
    const p2 = ftWorldPrice + ((i+1) / N) * (pH - ftWorldPrice);
    const s1 = produce(p1, K, L).qc;
    const s2 = produce(p2, K, L).qc;
    dPS += 0.5 * (s1 + s2) * (p2 - p1);
  }
  const dH_t = beta * IhFT / pH;
  const M_t  = Math.max(0, dH_t - prH.qc);
  const G    = t * M_t;
  return dPS + G - dCS;
}

// Welfare curve over a tariff grid.
welfareGrid = {
  const pts = [];
  const tMin = 0, tMax = 0.7;
  const N = 70;
  for (let i = 0; i <= N; i++) {
    const t = tMin + (i / N) * (tMax - tMin);
    pts.push({t, w: homeWelfareChange(t)});
  }
  // Normalize so ΔW(0) = 0 exactly (subtracts any bisection-precision noise).
  const w0 = pts[0].w;
  return pts.map(p => ({t: p.t, w: p.w - w0}));
}

// Optimal tariff: argmax of welfare grid.
optimalTariff = {
  let best = welfareGrid[0];
  for (const p of welfareGrid) {
    if (p.w > best.w) best = p;
  }
  return best;
}