11.3 消耗曲线仿真

# file: sim/ch11_consumption_curve.py
# tier: S
# verify: python3 sim/ch11_consumption_curve.py --config config/ --out out/
"""消耗曲线仿真：把 11.1 的养成消耗（建筑升级 / 科技）与活动型消耗（训练 / 行军 / 治疗）
汇总成 N 天曲线，并把 11.2 的市场回收单列。只做消耗侧，闭环（产出与消耗一起）留第 12 章。
产出 out/ch11_consumption.json / ch11_consumption.csv / ch11_consumption.svg"""
import argparse, csv, json, os, xml.etree.ElementTree as ET
from ch11_sink_slots import load, cost_at, tech_at, cval, out_at

COL = {"building": "#1f6feb", "tech": "#9333ea", "training": "#d97706", "healing": "#059669", "march": "#64748b"}
KEYS = ["building", "tech", "training", "healing", "march"]

def act_permille(c, day):
    return 1000 * c["activityGrowthNum"] ** (day - 1) // c["activityGrowthDen"] ** (day - 1)

def build_curve(c, sk, bl, jitter=None):
    val, tech, by_id = sk["resourceValue"], sk["tech"], {b["id"]: b for b in bl["buildings"]}
    levels, tl, order = dict(c["startLevels"]), c["techStartLevel"], c["upgradeOrder"]
    curve = []
    for day in range(1, c["days"] + 1):
        w = {k: 0 for k in KEYS}; by_res = {}
        eff = act_permille(c, day) * (jitter(day) if jitter else 1000) // 1000
        def add(cost, key=None):
            for r, v in cost.items():
                by_res[r] = by_res.get(r, 0) + v
            if key:
                w[key] += cost.get("wood", 0)
        if (day - 1) % c["upgradeEveryDays"] == 0:
            bid = order[((day - 1) // c["upgradeEveryDays"]) % len(order)]
            if levels[bid] < by_id[bid]["maxLevel"]:
                add(cost_at(by_id[bid], levels[bid]), "building"); levels[bid] += 1
        if (day - 1) % c["techEveryDays"] == 0 and tl < tech["maxLevel"]:
            add(tech_at(tech, tl), "tech"); tl += 1
        def scaled(base, count):
            return {r: v * count * eff // 1000 for r, v in base.items()}
        add(scaled(sk["unitCost"], c["trainUnitsPerDay"]), "training")
        add(scaled(sk["marchPerHour"], c["marchHoursPerDay"]), "march")
        add(scaled(sk["healPerWounded"], c["healWoundedPerDay"]), "healing")
        market = sk["referenceDay"]["tradesPerDay"] * eff // 1000 * sk["referenceDay"]["tradeAvgValue"] * c["taxPct"] // 100
        curve.append({"day": day, "wood": w, "totalWood": sum(w.values()),
                      "market": market, "totalValue": cval(by_res, val) + market})
    return curve

def draw(path, curve, res):
    W, H, pl, pr, pt, pb = 720, 360, 72, 24, 28, 52
    iw, ih = W - pl - pr, H - pt - pb
    hi = max(c["totalWood"] for c in curve) or 1; n = max(len(curve) - 1, 1)
    svg = ET.Element("svg", {"xmlns": "http://www.w3.org/2000/svg", "viewBox": f"0 0 {W} {H}"})
    ET.SubElement(svg, "rect", {"width": str(W), "height": str(H), "fill": "#fff"})
    ET.SubElement(svg, "text", {"x": str(pl), "y": "18", "font-size": "14"}).text = f"{res} 每日消耗曲线（{len(curve)} 天，仅消耗侧）"
    ET.SubElement(svg, "line", {"x1": str(pl), "y1": str(pt), "x2": str(pl), "y2": str(H - pb), "stroke": "#cbd5e1"})
    ET.SubElement(svg, "line", {"x1": str(pl), "y1": str(H - pb), "x2": str(W - pr), "y2": str(H - pb), "stroke": "#cbd5e1"})
    series = [("total", [c["totalWood"] for c in curve], "#111", 2.4)] + [(k, [c["wood"][k] for c in curve], COL[k], 1.4) for k in KEYS]
    for i, (name, ys, color, wd) in enumerate(series):
        pts = " ".join(f"{pl + iw * j / n:.1f},{pt + ih * (1 - v / hi):.1f}" for j, v in enumerate(ys))
        ET.SubElement(svg, "polyline", {"points": pts, "fill": "none", "stroke": color, "stroke-width": str(wd)})
        ET.SubElement(svg, "text", {"x": str(pl + 10 + i * 118), "y": str(H - pb + 22), "font-size": "11", "fill": color}).text = name
    os.makedirs(os.path.dirname(path), exist_ok=True)
    ET.ElementTree(svg).write(path, encoding="utf-8", xml_declaration=True)

def main():
    ap = argparse.ArgumentParser(); ap.add_argument("--config", default="config/"); ap.add_argument("--out", default="out/")
    a = ap.parse_args()
    c, sk, bl = load(a.config, "ch11_consumption"), load(a.config, "ch11_sinks"), load(a.config, "buildings")
    curve = build_curve(c, sk, bl)
    os.makedirs(a.out, exist_ok=True)
    with open(os.path.join(a.out, "ch11_consumption.csv"), "w", encoding="utf-8", newline="") as f:
        w = csv.writer(f); w.writerow(["day"] + KEYS + ["totalWood", "marketValue", "totalValue"])
        [w.writerow([x["day"]] + [x["wood"][k] for k in KEYS] + [x["totalWood"], x["market"], x["totalValue"]]) for x in curve]
    total = sum(x["totalWood"] for x in curve); first, last = curve[0]["totalWood"], curve[-1]["totalWood"]
    peak = max(curve, key=lambda x: x["totalWood"])
    summary = {"version": 1, "days": c["days"], "resource": c["resource"], "firstDayWood": first, "lastDayWood": last,
               "totalWood": total, "dayAvgWood": total // c["days"], "peakDay": peak["day"], "peakWood": peak["totalWood"],
               "growthPct": (last - first) * 100 // first if first else 0, "totalMarketValue": sum(x["market"] for x in curve)}
    with open(os.path.join(a.out, "ch11_consumption.json"), "w", encoding="utf-8") as f:
        json.dump(summary, f, ensure_ascii=False, indent=2)
    draw(os.path.join(a.out, "ch11_consumption.svg"), curve, c["resource"])
    print(json.dumps(summary, ensure_ascii=False))

if __name__ == "__main__":
    main()

# file: sim/ch11_consumption_bench.py
# tier: B
# verify: python3 sim/ch11_consumption_bench.py --config config/ --trials 2000
"""消耗曲线的不确定性基准（B 层，软告警不硬失败）。
确定性曲线假设每天活动强度按固定进度增长；真实玩家每天的训练量、行军时长、伤兵数会抖动。
对每日活动强度按 variancePct 做整数抖动，重复 trials 次，得到 N 天总消耗的分布，与点估计对照。
LCG 固定种子，结果可复现。"""
import argparse, json, os, sys
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
from ch11_sink_slots import load
from ch11_consumption_curve import build_curve

OUT = os.path.join(os.path.dirname(os.path.abspath(__file__)), "..", "out")

def lcg(seed):
    st = [seed & 0xFFFFFFFF]
    def nxt():
        st[0] = (st[0] * 1103515245 + 12345) % 2147483648
        return st[0]
    return nxt

def pct(xs, p):
    return xs[(p * (len(xs) - 1)) // 100] if xs else 0

def run(c, sk, bl, trials, seed, var):
    rng = lcg(seed); totals = []
    for _ in range(trials):
        jit = [(1000 - 10 * var) + rng() % (20 * var + 1) for _ in range(c["days"])]
        totals.append(sum(x["totalWood"] for x in build_curve(c, sk, bl, jitter=lambda d, j=jit: j[d - 1])))
    totals.sort(); det = sum(x["totalWood"] for x in build_curve(c, sk, bl)); p50 = pct(totals, 50)
    return {"trials": trials, "seed": seed, "variancePct": var, "days": c["days"], "deterministicTotal": det,
            "p10": pct(totals, 10), "p50": p50, "p90": pct(totals, 90), "p95": pct(totals, 95),
            "meanTotal": sum(totals) // len(totals), "minTotal": totals[0], "maxTotal": totals[-1],
            "spreadPctP10P90": (pct(totals, 90) - pct(totals, 10)) * 100 // (p50 or 1),
            "biasPct": (p50 - det) * 100 // (det or 1)}

def main():
    ap = argparse.ArgumentParser()
    ap.add_argument("--config", default=os.path.join(os.path.dirname(os.path.abspath(__file__)), "..", "config"))
    ap.add_argument("--trials", type=int, default=2000); ap.add_argument("--seed", type=int, default=20261012)
    ap.add_argument("--variance", type=int, default=25)
    a = ap.parse_args()
    c, sk, bl = load(a.config, "ch11_consumption"), load(a.config, "ch11_sinks"), load(a.config, "buildings")
    rep = run(c, sk, bl, a.trials, a.seed, a.variance)
    os.makedirs(OUT, exist_ok=True)
    with open(os.path.join(OUT, "ch11_consumption_bench.json"), "w", encoding="utf-8") as f:
        json.dump(rep, f, ensure_ascii=False, indent=2)
    print(f"trials={rep['trials']} seed={rep['seed']} variance=±{rep['variancePct']}% days={rep['days']}")
    print(f"deterministic={rep['deterministicTotal']}  p10={rep['p10']} p50={rep['p50']} p90={rep['p90']} p95={rep['p95']}")
    print(f"P10–P90 相对带宽={rep['spreadPctP10P90']}%  bias(P50-det)={rep['biasPct']}%")
    if rep["spreadPctP10P90"] > 60:
        print(f"WARN 消耗曲线不确定性偏大：P10–P90 带宽 {rep['spreadPctP10P90']}%，产出侧需留更多冗余")

if __name__ == "__main__":
    main()