# smooth-2d: cubic tensor-product B-spline fit of Taylor's plain tables (pol_mid.json)

**Result:** 2,665 bytes of table data. The loss on the report set is **+0.0097 (se 0.0090, n=2048)**, which is not significantly different from zero. Targets were at most 12 KB and a loss under 0.2.

I did all the work on Taylor's pure-Python replica (gtmcknight/amm-sim @ 4812a638, `research/solved`), run with plain `python3 -I`. I did not run anything on the Rust engine or in Solidity, so gate F still needs a contract.

## Gate P (before any change)
- `python3 -I gates.py` reproduced the README to every printed digit: gaussian+plain gave 443.455598 / 571.463409 / 525.168028 on seeds 10000/10001/10002. The other three columns matched too.
- `python3 -I real.py pol_mid.json plain:soft=1` printed `edge 532.386 ... vs HA +2.828 se 0.041`, the same as the README's +2.828.

## Step 1: the plain tables in pol_mid.json
| table | shape | rows (axis 1) | columns (axis 2) | value range (bps) |
|---|---|---|---|---|
| A0, B0 | 321 x 13 | w = (row-160) x 5 bps, -800..+800, step 5 | phase 0: R = (col-6) x 5 bps, -30..+30, step 5 (real.py clips R to +-30) | -24.262 .. 800.0 |
| A1, B1 | 321 x 49 | same | R = W2L[col] x 5 bps. W2L is in 5-bps lattice units: -160..160, non-uniform (steps of 1 near 0, 30 at the ends), so R spans -800..+800 bps | -31.323 .. 800.0 |
| A2, B2 | 321 x 49 | same | same as A1 | -31.333 .. 800.0 |

The values are ask (A) and bid (B) offsets from the price in bps. Other keys: g, lam 0.8, sigma 9.45, size 20, H_bps 5, NW 321, W0 160, N2 49, C2 24, FB 6, W2L, and VD0/VD1 (value differences, which real.py does not use). Together the six tables hold 321x(13+49+49)x2 = 71,262 values.

Two facts make the tables compressible:
1. **The tables are exactly mirror-symmetric:** B_p[i][j] == A_p[320-i][n2-1-j] for every cell in all three phases (max |diff| = 0.0). Only the three ask tables need to be stored.
2. **The ask fee A - w is never negative, and in 47% of cells it is exactly 0.** Those are cells where the DP's constraint a >= w binds. I fit the fee surface f = A - w, not A.

## Method
- For each phase p (0, 1, 2) there is one surface f_p(i, j) = sum_ab C_p[a][b] Bu_a(i) Bv_b(j). It uses cubic (degree 3) clamped B-splines in the row index i (0..320, which is w) and the column index j (0..12 or 0..48, which is R). The column index is a monotone reparametrization of R: it is the grid position real.py already interpolates in.
- **Decoding:** fee = max(0, f_p) at every grid node, A_p = w + fee, and B_p is A_p mirrored. The decoded tables then go through real.py's own bilinear lookup unchanged.
- **Weights:** I instrumented real.py's lookup (plain:soft=1) on choose seeds 60000..60255 and accumulated bilinear cell weights for each phase. The ask table gets h + mirror(h), because the bid lookups land on the mirrored ask cells. The least-squares weight is h/mean(h) + 1e-3.
- **Knots:** interior knots sit at quantiles of 0.8 x (visit marginal) + 0.2 x uniform, rounded to integer grid indices, with duplicates dropped. Rows get dense knots for |w| <= ~100 bps and sparse knots out to +-800.
- **Censoring:** cells whose fee is exactly 0 are one-sided. The fit is free to go below 0 there and is penalised only above 0. This is iteratively reweighted least squares: the target is min(fit, 0), with up to 30 iterations and ridge 1e-6.
- **Quantisation:** coefficients are stored as int16 in units of 1/8 bps.
- **Scoring:** `real.py` itself only prints a mean vs HA over seeds 10000..10000+NSEED. To get per-seed edges on any seed range, I wrote a small wrapper that calls `real.ev((path, {'soft': 1.0}, seeds))`, the same code path as `real.py <tables> plain:soft=1`, and pairs per seed. The README uses `real.py pol_mid.json plain:soft=1` too. As a cross-check, I also ran `python3 -I real.py <decoded.json> plain:soft=1` directly. It printed `edge 532.376 ... vs HA +2.819 se 0.041` against +2.828 for the full tables, a difference of 0.009 that matches the paired loss.
- loss = mean over seeds of (edge with full tables - edge with decoded tables), in real.py's edge units.

## Encoding (exact byte layout, big-endian)
```
byte 0      u8  k        spline degree (3)
byte 1      u8  s        coefficient shift: 1 unit = 2^-s bps (s = 3, so 1/8 bps)
then for phase p = 0, 1, 2 in order:
  u8  nu                 number of interior row knots
  u8  nv                 number of interior column knots
  nu x u16  KU[]         interior row knots, row-index units (0 < KU < 320), increasing
  nv x u8   KV[]         interior column knots, column-index units (0 < KV < n2-1; n2 = 13 for p=0, 49 otherwise)
  (nu+k+1)*(nv+k+1) x i16  C[a][b]   row-basis index a major, column-basis index b minor
```
- The full knot vectors are clamped: row [0,0,0,0, KU..., 320,320,320,320] and column [0 x4, KV..., (n2-1) x4].
- Chosen blob: phase 0 has nu=22, nv=5, 26x9 coefficients (519 B). Phases 1 and 2 each have nu=24, nv=14, 28x18 coefficients (1072 B). Total = 2 + 519 + 1072 + 1072 = **2,665 bytes**.
- **Value at grid node (i, j) of phase p:** ask fee F = max(0, sum_ab C[a][b] Bu_a(i) Bv_b(j)) / 2^s. Then A = w_i + F and the bid fee at (i, j) is the ask fee at (320-i, n2-1-j), with w_i = (i-160) x 5.
- **Solidity lookup equivalent to the reference:** compute pw, p2, i, j, f, g exactly as `DPPolicy._fees` does. Get the ask fee at the 4 nodes (i|i+1, j|j+1) and the bid fee at their 4 mirrored nodes. Interpolate bilinearly in fee space and apply the a > w / b > -w rule: because A - w is linear in w along a row, bilinear on A and then subtracting w gives the same result as bilinear on the fee.

## Tradeoff (choose set: seeds 60000..61023, n=1024, paired vs the full tables)
"Knots requested" is interior knots for phase 0 (rows x cols) / phases 1 and 2 (rows x cols). Rounding to integer indices can drop a few.

| setting | degree | knots requested | bytes | choose-set loss | se |
|---|---|---|---|---|---|
| s6x4 | 3 | 6x4 / 6x4 | 515 | +0.2087 | 0.0172 |
| q8x6 | 2 | 8x6 / 8x6 | 599 | +0.0968 | 0.0177 |
| s8x6 | 3 | 8x6 / 8x6 | 719 | +0.0867 | 0.0162 |
| s12x8 | 3 | 12x8 / 12x8 | 1,091 | +0.0519 | 0.0155 |
| s16x10 | 3 | 16x10 / 16x10 | 1,609 | +0.0366 | 0.0133 |
| **s24x14** | 3 | 24x10 / 24x14 | **2,665** | **+0.0098** | 0.0127 |
| s32x20 | 3 | 32x10 / 32x20 | 4,163 | +0.0102 | 0.0110 |

I picked s24x14 because it is the smallest setting whose choose-set loss is at the noise floor; s32x20 costs 56% more bytes and gains nothing measurable. Worst 5 seeds on the choose set for s24x14: 60822 +2.062, 60205 +1.266, 60916 +1.175, 60365 +1.140, 60455 +1.074.

## Chosen
- **Setting:** s24x14, cubic, s = 3 (1/8 bps). Fit command (scratch script, not shipped): `python fit.py s24x14 3 3 24 10 24 14`.
- **Bytes:** 2,665. sha256 6749a7d91c465a96ad29fa21a4f1e116ac00aaf53233688d27e78f4b42925a18.
- **Report set** (seeds 10000..12047, scored once): **loss +0.0097, se 0.0090, n = 2048**.
- **Worst 5 seeds** (full edge, decoded edge, full - decoded):

| seed | full | decoded | diff |
|---|---|---|---|
| 10359 | 666.107323 | 664.226064 | +1.881 |
| 11826 | 631.554259 | 629.937974 | +1.616 |
| 10652 | 709.461070 | 707.957388 | +1.504 |
| 11741 | 698.990976 | 697.538474 | +1.453 |
| 11696 | 717.264942 | 715.835881 | +1.429 |

- **Fit quality on the grid, visit-weighted:** mean |fee error| 0.17 / 0.26 / 0.26 bps and rms 0.50 / 0.85 / 0.71 bps for phases 0 / 1 / 2.
- **Worst single cells:** 23 bps (phase 0) and 59 bps (phases 1 and 2), all at the fee jumps for |w| <= 100 bps. The policy is insensitive to them on average.

## Lookup cost (one afterSwap lookup returns both fees; not measured in gas)
Per table (ask, then bid at the mirrored nodes):
- **Knot spans:** binary search over at most 24 u16 row knots and at most 14 u8 column knots. That is about 5 + 4 comparisons and reads for each of the 2 rows and 2 columns, or about 18 small reads.
- **Basis values:** degree-3 Cox-de Boor at 2 row nodes and 2 column nodes. Each point takes about 12 divisions and 24 multiply/add operations, so about 48 divisions and 96 multiply/add operations in total.
- **Coefficients:** at most 5 x 5 = 25 i16 reads (16 when both rows and both columns share a span).
- **Contraction:** 4 nodes x (16 multiply-adds + 4) = about 80 multiply-adds, then 4 max(0, ·).
- **Bilinear:** 8 multiplies.

A lookup for both tables is therefore about 2 x (about 18 knot reads + at most 25 coefficient reads) = **at most about 86 two-byte reads** and about **2 x 250 = about 500 arithmetic operations, about 100 of them divisions**. That is around 10-20k gas in a hand-written fixed-point implementation; this is an estimate, not a measurement.

Fixed point at 1e18 reproduces the float decoder to far below 0.01 bps. Because fee = max(0, ·) and bilinear interpolation are continuous, rounding cannot move a fee by more than the rounding itself.

A cheaper option is to evaluate the spline once at (pw, p2) directly instead of at 4 nodes. That would cost about a quarter of the arithmetic, but it would no longer equal the decoded-table reference scored here, so it would need its own gate Q run.

## Decoder (standalone; standard library + numpy; 66 lines)
```python
import struct
import numpy as np

W2L = [-160.0, -130.0, -100.0, -80.0, -60.0, -50.0, -40.0, -30.0, -25.0, -20.0, -17.0, -14.0, -12.0, -11.0, -10.0,
       -9.0, -8.0, -7.0, -6.0, -5.0, -4.0, -3.0, -2.0, -1.0, 0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0,
       11.0, 12.0, 14.0, 17.0, 20.0, 25.0, 30.0, 40.0, 50.0, 60.0, 80.0, 100.0, 130.0, 160.0]


def bspline_basis(x, inner, lo, hi, k):
    """Clamped B-spline basis of degree k on [lo, hi] with interior knots `inner`, at points x (any floats).
    Returns a (len(x), len(inner) + k + 1) matrix. Cox-de Boor; the right end x == hi uses the last interval."""
    t = np.array([lo] * (k + 1) + list(inner) + [hi] * (k + 1), dtype=float)
    x = np.asarray(x, dtype=float)
    nb = len(t) - 1
    B = np.zeros((len(x), nb))
    for m in range(nb):
        if t[m] < t[m + 1]:
            B[:, m] = (x >= t[m]) & (x < t[m + 1])
    last = max(m for m in range(nb) if t[m] < t[m + 1])
    B[x >= hi, :] = 0.0
    B[x >= hi, last] = 1.0
    for d in range(1, k + 1):
        Bn = np.zeros((len(x), nb - d))
        for m in range(nb - d):
            a = t[m + d] - t[m]
            b = t[m + d + 1] - t[m + 1]
            if a > 0:
                Bn[:, m] += (x - t[m]) / a * B[:, m]
            if b > 0:
                Bn[:, m] += (t[m + d + 1] - x) / b * B[:, m + 1]
        B = Bn
    return B


def decode(blob):
    """blob (bytes) -> dict in pol_mid.json's table format: A0 B0 A1 B1 A2 B2 plus the grid constants.
    Layout (big-endian): u8 degree k, u8 shift s; then for phase 0, 1, 2: u8 nu, u8 nv,
    nu x u16 interior row knots (row index units), nv x u8 interior column knots (column index units),
    (nu+k+1)*(nv+k+1) x i16 coefficients, row-basis major; one unit = 2**-s bps of ask fee."""
    k, s = blob[0], blob[1]
    pos = 2
    NW, W0, H = 321, 160, 5.0
    out = {'H_bps': H, 'NW': NW, 'W0': W0, 'N2': 49, 'C2': 24, 'FB': 6, 'W2L': list(W2L),
           'g': 0.05368206098114493, 'lam': 0.8, 'sigma': 9.45, 'size': 20.0}
    w = (np.arange(NW) - W0) * H
    for p in range(3):
        nu, nv = blob[pos], blob[pos + 1]
        pos += 2
        ku = struct.unpack('>%dH' % nu, blob[pos:pos + 2 * nu])
        pos += 2 * nu
        kv = struct.unpack('>%dB' % nv, blob[pos:pos + nv])
        pos += nv
        mu, mv = nu + k + 1, nv + k + 1
        C = np.array(struct.unpack('>%dh' % (mu * mv), blob[pos:pos + 2 * mu * mv]), dtype=float)
        pos += 2 * mu * mv
        C = C.reshape(mu, mv) / float(1 << s)
        n2 = 13 if p == 0 else 49
        Bu = bspline_basis(np.arange(NW), ku, 0, NW - 1, k)
        Bv = bspline_basis(np.arange(n2), kv, 0, n2 - 1, k)
        fee = np.maximum(Bu @ C @ Bv.T, 0.0)          # ask fee in bps at every grid node, never negative
        A = w[:, None] + fee                            # ask offset from the price, bps
        Bt = -w[:, None] + fee[::-1, ::-1]              # bid table is the ask table mirrored in (w, R)
        out['A%d' % p] = [[round(float(v), 6) for v in r] for r in A]
        out['B%d' % p] = [[round(float(v), 6) for v in r] for r in Bt]
    assert pos == len(blob)
    return out
```
Use: `import base64, json; json.dump(decode(base64.b64decode(DATA)), open('decoded.json', 'w'))`. The output is loadable by `python3 -I real.py decoded.json plain:soft=1` (the path is relative to the cwd). It contains A0 B0 A1 B1 A2 B2 and the grid constants, but omits VD0/VD1, which real.py does not read. I checked that decoding this base64 round-trips bit-identically to the tables I scored.

## Data (2,665 bytes, base64)
```
AwMWBQA7AGwAggCKAI8AkgCTAJUAlgCYAJsApACnAKkAqgCsAK0AsAC1AL0A0wEEAQMICgsYPhhCGEEYVRhBGHUYdxikGIEVLBUiFUEVFxW5FMkVCRQsFNgPlA+pD4IP3w70EIMQDxF4EHcI2wjKCPcItwmZCFgInweLCE0E4QT7BNAFPARnBdEGFQccBqUDLgMJA2kC/ASlBQ0EQgVMBHYCcQKUAmcC+wQEAp0FQADkBCsCDwIrAhEC1ALZBCEC+wR0A1MB2gHZAfUCbwMKAksDiwKDA2MBxAHeAccCVgJ4A1cCwQOCAvkBoAGPAcgB8QKKAoEC0wJvAtoBrwGkAXwB3gIEAjgCxQI1AocBAQH+ARIA5QG0AakBzAG1Ad4BOP6TAWv/6QDGAOEBKgECASH/LgPp/tAA8P+lAIEAIQEBAIIByvwRAcb+cwBz/5YA7/6rAD78WwLq/fUAzP+zADL+9wIrAA4Alf8ZAGb/rgAc/+wAP/92/+3/pQAZ/+QACv/8AAP/5wAyAAQAA//4AAT//AAC//4ACP/v//v//gAC//4AAf//AAH/+gAMAAEAAP//AAD//wAAAAAAA//6//4AAAAB//8AAAAAAAD/+wAJAAEAAP//AAH//wAA//8ABf/2//0AAAABAAAAAAAAAAD/+wAJAAEAAAAAAAAAAAAAAAAAAAAA//8YDgA+AHUAiwCRAJQAlgCYAJkAmgCcAJ0AnwCgAKIAowClAKYApwCpAKsArgC0AMoBAQoPEhQVFhcYGRobHSAlGAUX/BgcGBwYKBgpGDAYMBgXGCcYGRg1GEcYpRh0GFMX5RgyFNIU2RSjFOsU3BT5FNQU8xV8FZwVShUbFEYTxBQTFQ4WDRVrDqsOtQ7xDtgO7w7WDwsO5A4/DfYOXw7FD84QkRC+Dv0NSg5cB64HiQeMB44Hpge8B50Huwg/CHwIOgfuBywGkAYhB0YJ3gi0AzcDhwNmA5IDlAOTA7QDpwNEAwsDMwNgBBEEowVtBgMFtQXQAiYBnwISAfkCCwIxAhkCJAK2Az0DlgPNA9UD8gO9A6MD7APQAcQBEwGhAbIBuwHFAgICpAKwAp0CmAKgAssC+ANPA3oD1wOxAVkAsQE8AXsBjwG2AecB8wIPAlICfwKlAsIC6wL1AvkDSAMcALwAsAEAAUQBTgGaAbwB2gICAhACKAI7Al0CiAK7AuADNAMjAXH/awEZAQwBTAFWAZoBpwHGAeUCAQIbAj8CZwKKApsDAQLDAJr/7QDKANgBXgFYAV4BhgGVAbYBzgHrAf8CNAJaAnYC6wL1Adn+9QDKAKsBDQEnAS0BVAFiAYMBnQG1Ad4B/wIrAjcCvgI4/YAAE//tAMIBOAExAPcBFgEuAUEBYQGAAZIBxAHyAgIC5QLHAB4BX/7vAQQA1gD/ANwA4AD9AQ8BLQFIAWcBmwGrAhoCQQdlACr71f+kAQ0AvQDTAJQApwDCANAA8AEIATQBSwGiAR4GpARWAET5JgBqAKYAowChAJ4AfwCOAKgAuQDcAPgBJwFRAVQGzwRg/+7+yADXAFoAmgBqAIwARABTAGgAfwCYAMIA3gFBAQ8GTQQ8ACf/VAAqAGQAWQBOAEUAQAAnAEgASABqAJQAtwEVAQ0F5AQv/0gAAgAoAEYAMwA6ADYAFgAg//UANgAyAGIAfgD9AJwGEAPdABr/igAbACYAMgAkAAQALQAJ/2T//wAUADgAbgCdAPEE8gNt////cQAa/9AADP+r////sQAaAAgAM/9z/9IALACz//QGSgRk//4AFP/9AAj//AAP//8ADv/4//H/6gAXABf/sAAlAIUEHQIJAAD/9AAC//gABP/xAAH/8gAF//8AD//Q//wADP/r/wkF7wQ5AAAAAgAAAAP//wAFAAAABAAC/+wAFf+sAF/++wJl/EgDK/x/AAD//QAB//sAA//1AAD/9gABAA3/8QBA/2QBC/06A9z8igJPAAAAAv//AAT//QAKAAAACf/9//L///+7AHn+2gJE/GUD8v0qAAD//gAB//wAA//0AAD/9QADAAQAAgAd/68Aiv7XAF/+RwC1AAAAAAAAAAAAAAAAAAAAAP//////+//4//7/4gAk/9gAXv/SGA4APgB1AIoAkACTAJUAlwCZAJoAnACdAJ8AoACiAKMApQCmAKgAqgCsAK8AtQDKAQEJDRETFBYXGBkbHB4iJhgFF/8YFBgYGB4YGhghGBUYMRgiF+8YExf3GCsYMxh/GDwYVxTMFNoUrBTOFMsU6xTZFRkUkBUOFVYVlRWyFNEU4BQSFN0Ulg66DqwO5A7jDuoOyg7kDpcPSA6nDkAN2w3dDvIPARBPD1wPrQepB6EHmQeLB6EHugevB+gHbwfiCEQIqgh/B94HqwaJBxcG6QNuA4wDjAOVA6MDnAOvA4kD7wOXA0IC1gL/A18DrgTHBb4FwQIzAdwCHQIfAhoCPAI2AnEB4AJIAr0DwgQaBHwESwSGA+EECgHDAREBtAG6AcMBwQHaAVQCpQNEAz8DEgMAAwUDOwMrA34DoQFHASYBPAGdAZEBrgGgAkECcAIpAm8CsgL6Ay4DJwNZA0YDMwFwAEQBDgFEAUcBcQHGAckCBQJAAm4CcgKIAqUCwgK9AxQDAADVADwA3wESATYBOAGGAbEBzwHvAgMCNgJfAowCkQK4As4CqwEO/rQA/ADLAPsBIgFNAXYBlQG2AdwB+QIXAkACTgJjArICqACcABAAPwDVAMMBMwD/AVgBXQGEAaIByQHvAhcCHgJNAkoCIP/c/f0AfQB3ANEBKAEPAPUBNQE9AXEBiwGpAdkB6gHrAqICyP/SAD7/tQCFAPIA4QDwANEA/gEOATkBVQGBAaoBpAIkAhoAfwAX/+X/QwCVAO0AtQC/AIIAyQDNAPkBJQE+AXEBcQIqAnsHYwAW/JT/lADFAJsAlwCfAH4AjQCmAMQA6QERATcBcgEoCGsDUgAC/1gA4gBCAKcAaACjAEwATQBrAIEAtgDQAQcBPQFsBfYDCf5tAGj/1ACJADwAUQBBADMALwBCAFQAdQCzAMoBMAC5CMYCrP/bAND/2wBQACwAVQAMAGH/9v/zACoARQBsALYA4gD5BjsDiwAc/Z8AOgAL//X//AAZ/8kAGv+A//3/4gBKAGoAqgDNBxkB7v/1ADUALf/LAAT/8//2AAAABv/6AC//lP/XADwAkABIB1cBOQAD/67/8wAL//8AAwAAAAH/+QAE/80APf/r/83/+ADLBLsBn/////8ACv/2AAH//QAB//wABv/tABn/tAAQABgABv9gBxkCKwAA//z/+wADAAAAAAAC//gADf/OACj/nwBH/sYA3f23Ahf9TwAAAAAABf/8//8AAP/9AAb/9gAi/+EAR/+rAQn+0AMZ/SACGQAAAAD/+gAEAAEAAQAC//cADf/HAB7/mgBO/nAAj/v/A4/9ZwAA//8AA//9AAD/////AAEAAAAIAAYAC//2AF4AXf/4/usAmwAAAAD//wABAAAAAAAA//8AAf/4AAH/6QAI/9D/mwAMAED/2g==
```

## Problems
- **The real.py CLI does not fit the scoring protocol.** It only runs seeds 10000..10000+NSEED and prints a mean vs HA, with no per-seed values and no way to choose seeds 60000..61023. I scored through `real.ev`, the same function the CLI calls, with options `{'soft': 1.0}`. I also ran the CLI directly on the full and decoded tables over 10000..12047: +2.828 and +2.819 vs HA.
- **The choose set is also the fitting data.** The visit weights (seeds 60000..60255) and the knot quantiles come from it. The report set is untouched and agrees (+0.0097 vs +0.0098).
- **The column axis of phases 1 and 2 is the W2L position**, i.e. a non-uniform coordinate in R, not R in bps. A Solidity port needs W2L (49 values, a constant) for `_pos2` anyway.
- **Large errors at a few cells:** up to 59 bps at the fee discontinuities. Paired loss is the measure that matters and it is at the noise floor, but a few seeds lose up to about 1.9.
- **Not done here:** no Solidity, no engine run, no gate F, and no gas measurement. The lookup cost above is an operation count.