7.3 多重重要性采样

# snippet: sim/cg_trace_07_3_multiple_importance_sampling.py（节选：权重函数）
def weight(a, b, mode):
    if a <= 0 and b <= 0:
        return 0                      # 两条路径都采不到：权重取什么都不影响估计量
    if mode == "balance":
        return fdiv(a, a + b)
    a2, b2 = fmul(a, a), fmul(b, b)
    if a2 + b2 == 0:
        return 0                      # 幂启发式的定点下溢：p < 256 时 a² 被截断到 0
    return fdiv(a2, a2 + b2)

# snippet: 本机实跑（工作区 lantern/）
cd lantern
/opt/homebrew/bin/python3 sim/cg_trace_07_3_parity.py

#!/usr/bin/env python3
# file: lantern/sim/cg_trace_07_3_parity.py
# tier: S
# verify: /opt/homebrew/bin/python3 sim/cg_trace_07_3_parity.py
"""7.3 的紧凑驱动：调用 cg_trace_07_3_multiple_importance_sampling 的函数，打印正文数字。

三块输出：MIS 权重归一表（含幂启发式的定点下溢行）、三种光源尺寸下
{余弦, 光源, MIS-平衡, MIS-幂} 的方差对照、以及 P 层逐字节对照缓冲的 sha256。
"""

import hashlib
import os
import sys

HERE = os.path.dirname(os.path.abspath(__file__))
ROOT = os.path.dirname(HERE)
sys.path.insert(0, HERE)
import cg_trace_07_3_multiple_importance_sampling as m  # noqa: E402

ONE = 65536


def main():
    print("== 7.3 多重重要性采样 ==")
    print("权重归一（w_c + w_l 与 ONE 的残差，Q16.16）：")
    for r in m.weight_table():
        print(f"  p_c={r['pc']:>6} p_l={r['pl']:>6} : 平衡 {r['balanceSum']:<7} 幂 {r['powerSum']}")
    for label in ("small", "medium", "large"):
        st = m.study(os.path.join(ROOT, f"config/cg_trace_07_{label}.json"), label)
        print(f"[{st['label']:>6}] 半径 {st['radiusQ16'] / ONE:.4f}，每路径 {st['samplesPerTechnique']} 样本"
              f"（单路径 {st['singleRays']} 光线，MIS {st['misRays']} 光线）")
        for name in ("cosine", "light", "balance", "power"):
            d = st[name]
            print(f"    {name:>7}：相对标准差 {d['relStd']:.6f}｜相对偏差 {d['bias']:+.6f}")
        print(f"    MIS 平衡/较优单路径 = {st['balanceVsMin']:.4f}，"
              f"MIS 幂/较优单路径 = {st['powerVsMin']:.4f}，"
              f"残差 max（平衡/幂）= {st['weightResidualMaxQ16']['balance']}/"
              f"{st['weightResidualMaxQ16']['power']}")
    buf, _ = m.parity_values()
    print(f"P 层对照缓冲 sha256 : {hashlib.sha256(buf).hexdigest()}")
    return 0


if __name__ == "__main__":
    sys.exit(main())

== 7.3 多重重要性采样 ==
权重归一（w_c + w_l 与 ONE 的残差，Q16.16）：
  p_c= 65536 p_l=     0 : 平衡 65536   幂 65536
  p_c=     0 p_l= 65536 : 平衡 65536   幂 65536
  p_c= 65536 p_l= 65536 : 平衡 65536   幂 65536
  p_c= 32768 p_l= 65536 : 平衡 65535   幂 65535
  p_c= 32768 p_l=  9830 : 平衡 65535   幂 65535
  p_c=  9830 p_l= 32768 : 平衡 65535   幂 65535
  p_c=  2000 p_l= 40000 : 平衡 65535   幂 65535
  p_c=     1 p_l=     1 : 平衡 65536   幂 0
  p_c= 50000 p_l=     3 : 平衡 65535   幂 65536
[ small] 半径 0.1500，每路径 64 样本（单路径 64 光线，MIS 128 光线）
     cosine：相对标准差 1.074105｜相对偏差 -0.037993
      light：相对标准差 0.186531｜相对偏差 +0.001024
    balance：相对标准差 0.194863｜相对偏差 +0.000369
      power：相对标准差 0.198296｜相对偏差 -0.000502
    MIS 平衡/较优单路径 = 1.0447，MIS 幂/较优单路径 = 1.0631，残差 max（平衡/幂）= 1/1
[medium] 半径 0.2800，每路径 64 样本（单路径 64 光线，MIS 128 光线）
     cosine：相对标准差 0.575925｜相对偏差 -0.011372
      light：相对标准差 0.222106｜相对偏差 +0.001686
    balance：相对标准差 0.201072｜相对偏差 -0.001069
      power：相对标准差 0.206953｜相对偏差 -0.001244
    MIS 平衡/较优单路径 = 0.9053，MIS 幂/较优单路径 = 0.9318，残差 max（平衡/幂）= 1/1
[ large] 半径 0.5000，每路径 64 样本（单路径 64 光线，MIS 128 光线）
     cosine：相对标准差 0.308911｜相对偏差 +0.005613
      light：相对标准差 0.337317｜相对偏差 +0.007169
    balance：相对标准差 0.195892｜相对偏差 +0.001579
      power：相对标准差 0.201152｜相对偏差 +0.001230
    MIS 平衡/较优单路径 = 0.6341，MIS 幂/较优单路径 = 0.6512，残差 max（平衡/幂）= 1/1
P 层对照缓冲 sha256 : 70d0aef0c993bf23669edc33dc2232e111b5f0534efde455b3245d76396473d7

// file: lantern/core/test/cg_trace_07_3_parity.test.ts
// tier: P
// verify: node --test 'core/test/cg_trace_07_3_parity.test.ts'
// 7.3 的最小 P 复现：读 Python 侧 cg_trace_07_3_mis.json，把权重归一残差、
// 幂启发式的定点下溢、以及 MIS 相对两条单路径的比值钉住。

import assert from "node:assert/strict";
import { readFileSync } from "node:fs";
import { dirname, join } from "node:path";
import { test } from "node:test";
import { fileURLToPath } from "node:url";

import { ONE } from "../src/lantern/fixed.ts";
import { weight, weightTable } from "../src/lantern/cg_trace_07_3_multiple_importance_sampling.ts";

const root = join(dirname(fileURLToPath(import.meta.url)), "..", "..");
const mis = JSON.parse(
  readFileSync(join(root, "out", "cg_trace_07_3_mis.json"), "utf8"),
) as {
  small: { balanceVsMin: number; powerVsMin: number };
  medium: { balanceVsMin: number };
  large: { balanceVsMin: number };
  parity: {
    weightTable: { pc: number; pl: number; balanceSum: number; powerSum: number }[];
  };
};

test("【核心】权重表逐行与 Python 侧一致", () => {
  const rows = weightTable();
  assert.equal(rows.length, mis.parity.weightTable.length);
  rows.forEach((r, i) => {
    assert.equal(Number(r.pc), mis.parity.weightTable[i].pc);
    assert.equal(Number(r.pl), mis.parity.weightTable[i].pl);
    assert.equal(Number(r.balanceSum), mis.parity.weightTable[i].balanceSum);
    assert.equal(Number(r.powerSum), mis.parity.weightTable[i].powerSum);
  });
});

test("【核心】平衡启发式的权重和恒等于 ONE（残差 ≤ 1）", () => {
  for (const r of weightTable()) {
    if (r.pc === 0n && r.pl === 0n) continue;
    const d = r.balanceSum - ONE;
    assert.ok(d >= -1n && d <= 1n, `平衡权重和 ${r.balanceSum} 与 ONE 的残差 ${d} 应不超过 1`);
  }
});

test("【核心】MIS 在中等/大光源下严格优于两条单路径", () => {
  assert.ok(mis.medium.balanceVsMin < 0.95, `中等光源 ${mis.medium.balanceVsMin} 应 < 0.95`);
  assert.ok(mis.large.balanceVsMin < 0.7, `大光源 ${mis.large.balanceVsMin} 应 < 0.7`);
  assert.ok(mis.small.balanceVsMin > 0.95 && mis.small.balanceVsMin < 1.15);
});

test("【负例】幂启发式的定点下溢：两个 pdf 都小于 256 时权重退化成 0", () => {
  assert.equal(weight(1n, 1n, "power"), 0n);
  assert.equal(weight(255n, 255n, "power"), 0n);
  assert.ok(weight(256n, 256n, "power") > 0n);
  assert.equal(weight(1n, 1n, "balance"), 32768n);
});

✔ 【核心】权重表逐行与 Python 侧一致 (0.612208ms)
✔ 【核心】平衡启发式的权重和恒等于 ONE（残差 ≤ 1） (0.131916ms)
✔ 【核心】MIS 在中等/大光源下严格优于两条单路径 (0.072959ms)
✔ 【负例】幂启发式的定点下溢：两个 pdf 都小于 256 时权重退化成 0 (0.066875ms)
ℹ tests 4 / pass 4 / fail 0