guzhujushi
953e988983
用户本轮改了 markov_simulate / markov_steady_state 的入口校验，并写完 mm1_metrics，
三处算法本体都正确（mm1 的 5 个公式经探针与闭式解逐位对拍通过，含 W = Wq + 1/mu）。
Codex 只动工程侧：
Fixed
- mm1_metrics 的稳定性分界 "<" 应为 "<="（docstring 与既有测试的口径都是 lambda >= mu 不稳定）：
lambda == mu 原先掉进 try、rho == 1.0 触发 ZeroDivisionError，被 except 吞成 success=False。
两个分支的可观测输出完全相同（success=False / metrics={} / data=None），只有 message 不同，
故按决策 31 直接修并加回归测试，未改动任何口径。
- 三处报错文案的参数名抄错：markov_steady_state 的 max_iter <= 0 写成 "tol 必须是正数"，
mm1_metrics 的 service_rate 两条报错都写成 "tol"；半角逗号统一为全角。
Changed
- mm1_metrics 的 docstring 补 Raises 与 data 说明；call_params = {...} 补空格；
ruff 6 条（W293 x5 / F541 f-string 无占位符 / W292 文件末尾无换行）。
- markov_steady_state 的 try 前 3 个空行 -> 1；"异常: " -> "异常："。
- tests/test_solve.py：摘掉 2 条 mm1_metrics 的 skip，新增 4 条（lambda == mu 的回归、
非负校验回归、max_iter=3 未收敛判 success=False）。
- STATUS.md 决策 56 / 57 更新办结项 + 新增决策 58（mm1 的 4 条待拍板）；CHANGELOG；易错点 29。
校验：pytest 202 passed / 28 skipped（230 收集）；ruff 0；pyright 0
695 lines
26 KiB
Python
695 lines
26 KiB
Python
"""solve 包测试。"""
|
||
|
||
from __future__ import annotations
|
||
|
||
import numpy as np
|
||
import pytest
|
||
|
||
_NOT_IMPLEMENTED = "未实现：见 docs/03-roadmap.md L0/L1"
|
||
|
||
|
||
def _normal_batch(rng, size, **kwargs):
|
||
"""`monte_carlo` 的 sampler 契约：收下实现传来的 `size`，一次抽够。"""
|
||
return rng.normal(0, 1, size)
|
||
|
||
|
||
def test_monte_carlo_integral_converges_to_analytic_value():
|
||
"""∫₀¹ x² dx = 1/3。n=2e5 时误差应 < 0.01。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
r = monte_carlo_integral(lambda X: X[:, 0] ** 2, bounds=[(0, 1)], n=200_000, random_state=0)
|
||
assert r.success
|
||
assert abs(r.metrics["mean"] - 1 / 3) < 0.01 # 体积必须是正的
|
||
assert len(r.data) == 200_000 # data 是 n 个样本各自的贡献
|
||
|
||
|
||
def test_monte_carlo_integral_odd_dimension_sign():
|
||
"""∫∫∫ x dV over [0,1]³ = 1/2 —— 奇数维最容易暴露"体积符号写反"。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
r = monte_carlo_integral(lambda X: X[:, 0], bounds=[(0, 1)] * 3, n=200_000, random_state=0)
|
||
assert abs(r.metrics["mean"] - 0.5) < 0.01
|
||
|
||
|
||
def test_monte_carlo_integral_signed_bounds():
|
||
"""∫₀² x dx = 2 —— 边界带宽度时，体积算错和符号写反都会现形。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
r = monte_carlo_integral(lambda X: X[:, 0], bounds=[(0, 2)], n=200_000, random_state=0)
|
||
assert abs(r.metrics["mean"] - 2.0) < 0.02
|
||
|
||
|
||
def test_monte_carlo_integral_2d_and_params():
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
bounds = [(0, 1), (0, 1)]
|
||
r = monte_carlo_integral(lambda X: X[:, 0] * X[:, 1], bounds=bounds, n=200_000, random_state=0)
|
||
assert abs(r.metrics["mean"] - 0.25) < 0.01
|
||
assert r.name == "monte_carlo_integral"
|
||
assert r.params["dimension"] == 2
|
||
assert r.params["bounds"] == bounds
|
||
|
||
|
||
def test_monte_carlo_integral_rejects_bad_bounds():
|
||
"""四类非法 bounds 都要在入口抛 ValueError：区间反了 / 漏括号 / 空 / 参差不齐。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
for bad in ([(1, 0)], [0, 1], [], [(0, 1), (0, 2, 3)]):
|
||
with pytest.raises(ValueError):
|
||
monte_carlo_integral(lambda X: X[:, 0], bounds=bad, n=100, random_state=0)
|
||
|
||
|
||
def test_monte_carlo_integral_failure_message_stays_failure():
|
||
"""失败时 message 不能被包装层改写成"成功"（易错点 22）。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
r = monte_carlo_integral(lambda X: 1.0, bounds=[(0, 1)], n=100, random_state=0)
|
||
assert r.success is False
|
||
assert "成功" not in r.message
|
||
assert "失败" in r.message
|
||
assert r.params["dimension"] == 1 # 失败也要留下追溯信息
|
||
|
||
|
||
def test_monte_carlo_integral_accepts_numpy_integer_n():
|
||
"""§4.7：`np.int64` 不是 `int` 的子类，注解与运行期都要收（易错点 1）。"""
|
||
from mcm_toolkit.solve import monte_carlo_integral
|
||
|
||
r = monte_carlo_integral(lambda X: X[:, 0], bounds=[(0, 1)], n=np.int64(1000), random_state=0)
|
||
assert r.success
|
||
assert len(r.data) == 1000
|
||
|
||
|
||
def test_monte_carlo_probability_of_standard_normal():
|
||
"""P(Z > 0) = 0.5，n=5000 时误差应 < 0.05。"""
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
r = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=5000, random_state=0)
|
||
assert r.success
|
||
assert abs(r.metrics["probability"] - 0.5) < 0.05
|
||
|
||
|
||
def test_monte_carlo_is_reproducible():
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
a = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=1000, random_state=9)
|
||
b = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=1000, random_state=9)
|
||
assert a.metrics["probability"] == b.metrics["probability"]
|
||
|
||
|
||
def test_monte_carlo_returns_n_results_and_documented_metrics():
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
r = monte_carlo(_normal_batch, lambda x: x**2, n=2000, random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (2000,)
|
||
assert set(r.metrics) == {"mean", "std", "median", "p05", "p95", "ci_low", "ci_high"}
|
||
|
||
|
||
def test_monte_carlo_calls_simulator_only_once():
|
||
"""`simulator` 只被调用一次 —— 这是"必须向量化"的可执行版本。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
counter = {"calls": 0}
|
||
|
||
def simulator(params):
|
||
counter["calls"] += 1
|
||
return params
|
||
|
||
r = monte_carlo(_normal_batch, simulator, n=1000, random_state=0)
|
||
assert r.success
|
||
assert counter["calls"] == 1
|
||
|
||
|
||
def test_monte_carlo_data_is_reproducible_with_same_seed():
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
a = monte_carlo(_normal_batch, lambda x: x, n=2000, random_state=11)
|
||
b = monte_carlo(_normal_batch, lambda x: x, n=2000, random_state=11)
|
||
assert np.array_equal(a.data, b.data)
|
||
|
||
|
||
def test_monte_carlo_follows_set_seed_when_random_state_omitted():
|
||
from mcm_toolkit.core import set_seed
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
set_seed(5)
|
||
a = monte_carlo(_normal_batch, lambda x: x, n=500)
|
||
set_seed(5)
|
||
b = monte_carlo(_normal_batch, lambda x: x, n=500)
|
||
assert np.array_equal(a.data, b.data)
|
||
|
||
|
||
def test_monte_carlo_estimates_known_expectation():
|
||
"""E[X²] = 1（X ~ N(0, 1)）；n = 2e5 时误差应 < 0.02。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
r = monte_carlo(_normal_batch, lambda x: x**2, n=200_000, random_state=0)
|
||
assert abs(r.metrics["mean"] - 1.0) < 0.02
|
||
|
||
|
||
def test_monte_carlo_docstring_example_runs():
|
||
"""docstring 的 Examples 必须真的能跑通（此前那一版是坏的）。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
draw = lambda rng, size, **k: rng.normal(0, 1, size) # noqa: E731
|
||
r = monte_carlo(draw, lambda x: x**2, n=10000, random_state=0)
|
||
assert r.success
|
||
assert abs(r.metrics["mean"] - 1.0) < 0.1
|
||
|
||
|
||
def test_monte_carlo_accepts_numpy_integer_n():
|
||
"""§4.7：`np.int64` 不是 `int` 的子类，不能漏。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
r = monte_carlo(_normal_batch, lambda x: x, n=np.int64(100), random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (100,)
|
||
|
||
|
||
def test_monte_carlo_counts_boolean_outcomes():
|
||
"""布尔结果按 0/1 统计 —— 事件概率是 MC 最常见的用法。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
r = monte_carlo(
|
||
lambda rng, size, **k: rng.normal(0, 1, size) > 0,
|
||
lambda x: x,
|
||
n=5000,
|
||
random_state=0,
|
||
)
|
||
assert r.success
|
||
assert r.data.dtype == float
|
||
assert abs(r.metrics["mean"] - 0.5) < 0.05
|
||
|
||
|
||
def test_monte_carlo_rejects_nan_results():
|
||
"""结果里出现 NaN 就判失败，不静默给出一列 nan 指标。"""
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
r = monte_carlo(
|
||
lambda rng, size, **k: np.r_[np.zeros(size - 1), np.nan],
|
||
lambda x: x,
|
||
n=10,
|
||
random_state=0,
|
||
)
|
||
assert not r.success
|
||
assert "NaN" in r.message
|
||
|
||
|
||
def test_monte_carlo_rejects_bad_arguments():
|
||
from mcm_toolkit.solve import monte_carlo
|
||
|
||
with pytest.raises(TypeError):
|
||
monte_carlo(123, lambda x: x) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
monte_carlo(_normal_batch, 123) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
monte_carlo(_normal_batch, lambda x: x, n=True)
|
||
with pytest.raises(TypeError):
|
||
monte_carlo(_normal_batch, lambda x: x, n=10.0) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(ValueError):
|
||
monte_carlo(_normal_batch, lambda x: x, n=0)
|
||
|
||
|
||
def test_monte_carlo_probability_returns_documented_metrics():
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
r = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=20000, random_state=0)
|
||
assert r.success
|
||
assert set(r.metrics) == {"probability", "ci_low", "ci_high", "std_error"}
|
||
assert r.data.shape == (20000,)
|
||
assert 0.0 <= r.metrics["ci_low"] <= r.metrics["probability"]
|
||
assert r.metrics["probability"] <= r.metrics["ci_high"] <= 1.0
|
||
|
||
|
||
def test_monte_carlo_probability_std_error_matches_binomial_formula():
|
||
"""`std_error` 就是比例的标准误 √(p(1-p)/n)。"""
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
n = 20000
|
||
r = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=n, random_state=0)
|
||
p = r.metrics["probability"]
|
||
assert np.isclose(r.metrics["std_error"], np.sqrt(p * (1 - p) / n))
|
||
|
||
|
||
def test_monte_carlo_probability_ci_narrows_with_n():
|
||
"""区间宽度随 n 增大而收窄（O(1/√n)）。"""
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
def width(n):
|
||
r = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=n, random_state=0)
|
||
return r.metrics["ci_high"] - r.metrics["ci_low"]
|
||
|
||
assert width(100_000) < width(1_000)
|
||
|
||
|
||
def test_monte_carlo_probability_higher_confidence_is_wider():
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
def width(confidence):
|
||
r = monte_carlo_probability(
|
||
lambda z: z > 0, _normal_batch, n=5000, random_state=0, confidence=confidence
|
||
)
|
||
return r.metrics["ci_high"] - r.metrics["ci_low"]
|
||
|
||
assert width(0.99) > width(0.90)
|
||
|
||
|
||
def test_monte_carlo_probability_ci_matches_statsmodels_wilson():
|
||
"""区间必须对得上独立实现（statsmodels 的 wilson）—— `^` 写成异或那次就靠它拦住。"""
|
||
from statsmodels.stats.proportion import proportion_confint
|
||
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
n = 5000
|
||
r = monte_carlo_probability(lambda z: z > 0, _normal_batch, n=n, random_state=0)
|
||
count = int(r.data.sum())
|
||
lo, hi = proportion_confint(count, n, alpha=0.05, method="wilson")
|
||
assert np.isclose(r.metrics["ci_low"], lo)
|
||
assert np.isclose(r.metrics["ci_high"], hi)
|
||
|
||
|
||
def test_monte_carlo_probability_rejects_bad_confidence():
|
||
from mcm_toolkit.solve import monte_carlo_probability
|
||
|
||
for bad in (0.0, 1.0, -0.5, 1.5):
|
||
with pytest.raises(ValueError):
|
||
monte_carlo_probability(lambda z: z > 0, _normal_batch, n=100, confidence=bad)
|
||
|
||
|
||
def test_random_walk_shape_and_start():
|
||
"""`data` 形状 (n_paths, n_steps + 1)，第 0 列是起点。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(n_steps=100, n_paths=5, start=2.5, random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (5, 101)
|
||
assert np.all(r.data[:, 0] == 2.5)
|
||
|
||
|
||
def test_random_walk_default_steps_are_plus_minus_one():
|
||
"""默认是对称随机游走：每一步都恰好是 ±1。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(n_steps=200, n_paths=3, random_state=0)
|
||
assert set(np.unique(np.diff(r.data, axis=1))) == {-1.0, 1.0}
|
||
|
||
|
||
def test_random_walk_is_reproducible():
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
a = random_walk(n_steps=50, n_paths=4, random_state=7)
|
||
b = random_walk(n_steps=50, n_paths=4, random_state=7)
|
||
assert np.array_equal(a.data, b.data)
|
||
|
||
|
||
def test_random_walk_accepts_custom_step_sampler():
|
||
"""换成高斯增量就是带漂移的高斯游走 —— 建模差异全在 step_sampler 里。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(
|
||
n_steps=100,
|
||
n_paths=3,
|
||
step_sampler=lambda rng, size: rng.normal(0.5, 1.0, size),
|
||
random_state=0,
|
||
)
|
||
assert r.success
|
||
assert r.data.shape == (3, 101)
|
||
assert abs(r.data[:, -1].mean() - 50.0) < 20.0 # 100 步 × 漂移 0.5
|
||
|
||
|
||
def test_random_walk_rejects_bad_arguments():
|
||
"""§4.7 的分工：类型不对抛 TypeError，取值不对（<= 0）抛 ValueError。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
with pytest.raises(TypeError):
|
||
random_walk(n_steps=2.5) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
random_walk(n_steps=True)
|
||
with pytest.raises(TypeError):
|
||
random_walk(n_steps=10, n_paths=1.5) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(ValueError):
|
||
random_walk(n_steps=0)
|
||
with pytest.raises(ValueError):
|
||
random_walk(n_steps=10, n_paths=0)
|
||
with pytest.raises(TypeError):
|
||
random_walk(n_steps=10, start="a") # type: ignore[arg-type] # 就是要传错
|
||
|
||
|
||
def test_random_walk_rejects_2d_step_sampler():
|
||
"""`step_sampler` 必须返回**一维** —— 二维返回（含方阵）一律判失败（易错点 25 的根治）。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
def transposed(rng, size):
|
||
return rng.normal(0, 1, (4, 3)) # 无论要多少，一律给二维
|
||
|
||
r = random_walk(n_steps=4, n_paths=3, step_sampler=transposed, random_state=0)
|
||
assert r.success is False # 失败也要返回 Result，不能抛给调用方
|
||
assert "长度" in r.message
|
||
assert r.params["n_paths"] == 3 # 失败也要留下追溯信息
|
||
assert r.data is None
|
||
|
||
# 旧方案在 n_paths == n_steps 时对"转置"完全无感；现在同样拦住
|
||
assert random_walk(n_steps=3, n_paths=3, step_sampler=transposed, random_state=0).success is False
|
||
|
||
|
||
def test_random_walk_step_sampler_gets_flat_size_and_is_row_major():
|
||
"""契约：`size` 是整数 `n_paths * n_steps`，一维增量按**行主序**切成 n_paths 行。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
seen = {}
|
||
|
||
def sampler(rng, size):
|
||
seen["size"] = size
|
||
return np.arange(1, size + 1, dtype=float) # 1, 2, 3, 4, 5, 6
|
||
|
||
r = random_walk(n_steps=3, n_paths=2, step_sampler=sampler, random_state=0)
|
||
assert seen["size"] == 6 # 整数，不是 (n_paths, n_steps)
|
||
assert r.data.tolist() == [[0, 1, 3, 6], [0, 4, 9, 15]] # cumsum(1,2,3) / cumsum(4,5,6)
|
||
|
||
|
||
def test_random_walk_rejects_nan_steps():
|
||
"""增量含 NaN 判失败 —— 与 monte_carlo 的 NaN 口径（决策 48）一致。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(
|
||
n_steps=4,
|
||
n_paths=2,
|
||
step_sampler=lambda rng, size: np.full(size, np.nan),
|
||
random_state=0,
|
||
)
|
||
assert r.success is False
|
||
assert "NaN" in r.message
|
||
|
||
|
||
def test_random_walk_accepts_numpy_integer():
|
||
"""§4.7：`np.int64` 不是 `int` 的子类，注解与运行期都要收（易错点 1）。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(n_steps=np.int64(20), n_paths=np.int64(3), random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (3, 21)
|
||
|
||
|
||
def test_random_walk_metrics_are_terminal_positions():
|
||
"""终点位置的标准差按 √t 扩散（t = n_steps）—— 随机游走的核心性质。"""
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
n_steps = 100
|
||
r = random_walk(n_steps=n_steps, n_paths=2000, start=5.0, random_state=0)
|
||
assert abs(r.metrics["mean_terminal"] - 5.0) < 5.0 # 对称游走，期望还在起点附近
|
||
assert abs(r.metrics["std_terminal"] - np.sqrt(n_steps)) < 5.0 # √t 扩散
|
||
|
||
|
||
def test_random_walk_docstring_example_runs():
|
||
from mcm_toolkit.solve import random_walk
|
||
|
||
r = random_walk(n_steps=100, n_paths=5, random_state=0)
|
||
assert r.data.shape == (5, 101)
|
||
assert np.all(np.abs(np.diff(r.data, axis=1)) == 1)
|
||
|
||
def test_markov_simulate_shape_start_and_states():
|
||
"""`data` 形状 (n_paths, steps + 1)，第 0 列是初始状态，取值只能是合法状态编号。"""
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_simulate(P, initial=0, steps=100, n_paths=5, random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (5, 101)
|
||
assert np.all(r.data[:, 0] == 0)
|
||
assert set(np.unique(r.data)) <= {0, 1}
|
||
assert np.issubdtype(r.data.dtype, np.integer)
|
||
|
||
|
||
def test_markov_simulate_absorbing_states_never_move():
|
||
"""单位转移矩阵：状态自锁，整条路径恒等于初始状态（确定性自检）。"""
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
r = markov_simulate(np.eye(2), initial=1, steps=30, n_paths=4, random_state=0)
|
||
assert r.success
|
||
assert np.all(r.data == 1)
|
||
assert r.metrics["mean_terminal"] == 1.0
|
||
assert r.metrics["std_terminal"] == 0.0
|
||
|
||
|
||
def test_markov_simulate_is_reproducible():
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.7, 0.3], [0.4, 0.6]])
|
||
a = markov_simulate(P, initial=0, steps=40, n_paths=6, random_state=11)
|
||
b = markov_simulate(P, initial=0, steps=40, n_paths=6, random_state=11)
|
||
c = markov_simulate(P, initial=0, steps=40, n_paths=6, random_state=12)
|
||
assert np.array_equal(a.data, b.data)
|
||
assert not np.array_equal(a.data, c.data)
|
||
|
||
|
||
def test_markov_simulate_terminal_frequency_matches_stationary():
|
||
"""终点频率收敛到平稳分布 [2/3, 1/3] —— 马尔可夫链模拟的核心性质。"""
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_simulate(P, initial=0, steps=200, n_paths=20000, random_state=0)
|
||
freq = np.bincount(r.data[:, -1], minlength=2) / r.data.shape[0]
|
||
assert np.allclose(freq, [2 / 3, 1 / 3], atol=0.02)
|
||
|
||
|
||
def test_markov_simulate_metrics_are_terminal_states():
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_simulate(P, initial=0, steps=50, n_paths=500, random_state=3)
|
||
assert r.metrics["mean_terminal"] == float(np.mean(r.data[:, -1]))
|
||
assert r.metrics["std_terminal"] == float(np.std(r.data[:, -1]))
|
||
|
||
|
||
def test_markov_simulate_rejects_bad_transition():
|
||
"""非方阵、行和不为 1 都要在入口直接报错，不能悄悄模拟。"""
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
with pytest.raises(ValueError):
|
||
markov_simulate(np.array([[0.5, 0.5, 0.0], [0.2, 0.8, 0.0]]), initial=0)
|
||
with pytest.raises(ValueError):
|
||
markov_simulate(np.array([[0.5, 0.4], [0.2, 0.8]]), initial=0)
|
||
with pytest.raises(ValueError):
|
||
markov_simulate(np.array([[1.5, -0.5], [0.0, 1.0]]), initial=0) # 行和是 1，但含负元素
|
||
|
||
|
||
def test_markov_simulate_rejects_bad_arguments():
|
||
"""§4.7：步数 / 路径数的类型与正负都要卡住（bool 先判，易错点 1）。"""
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
with pytest.raises(ValueError):
|
||
markov_simulate(P, initial=0, steps=0)
|
||
with pytest.raises(ValueError):
|
||
markov_simulate(P, initial=0, n_paths=0)
|
||
with pytest.raises(TypeError):
|
||
markov_simulate(P, initial=0, steps=1.5) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
markov_simulate(P, initial=0, n_paths=1.5) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
markov_simulate(P, initial=0, steps=True) # type: ignore[arg-type] # 就是要传错
|
||
with pytest.raises(TypeError):
|
||
markov_simulate(P, initial=0, n_paths=True) # type: ignore[arg-type] # 就是要传错
|
||
|
||
|
||
def test_markov_simulate_accepts_numpy_integer():
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_simulate(P, initial=0, steps=np.int64(10), n_paths=np.int64(2), random_state=0)
|
||
assert r.success
|
||
assert r.data.shape == (2, 11)
|
||
|
||
|
||
def test_markov_simulate_docstring_example_runs():
|
||
from mcm_toolkit.solve import markov_simulate
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_simulate(P, initial=0, steps=50, n_paths=4, random_state=0)
|
||
assert r.data.shape == (4, 51)
|
||
assert np.all(r.data[:, 0] == 0)
|
||
|
||
|
||
def test_markov_steady_state_two_state_chain():
|
||
"""P=[[0.9,0.1],[0.2,0.8]] 的平稳分布是 [2/3, 1/3]。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
r = markov_steady_state(np.array([[0.9, 0.1], [0.2, 0.8]]))
|
||
assert r.success
|
||
assert np.allclose(r.data, [2 / 3, 1 / 3], atol=1e-6)
|
||
assert r.metrics["converged"] is True
|
||
assert r.metrics["iterations"] > 0
|
||
|
||
|
||
def test_markov_steady_state_matches_eigenvector():
|
||
"""与特征分解对拍：取特征值 1 对应的左特征向量（只是检验手段，实现仍走迭代）。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.5, 0.25, 0.25], [0.1, 0.8, 0.1], [0.2, 0.2, 0.6]])
|
||
eigenvalues, eigenvectors = np.linalg.eig(P.T)
|
||
ref = np.real(eigenvectors[:, int(np.argmin(np.abs(eigenvalues - 1.0)))])
|
||
ref = ref / ref.sum()
|
||
|
||
r = markov_steady_state(P)
|
||
assert np.allclose(r.data, ref, atol=1e-6)
|
||
|
||
|
||
def test_markov_steady_state_is_a_distribution():
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.7, 0.2, 0.1], [0.1, 0.6, 0.3], [0.05, 0.15, 0.8]])
|
||
r = markov_steady_state(P)
|
||
assert r.data.shape == (3,)
|
||
assert np.isclose(float(np.sum(r.data)), 1.0)
|
||
assert np.all(r.data >= 0.0)
|
||
assert r.data.dtype == np.float64
|
||
|
||
|
||
def test_markov_steady_state_slow_mixing_needs_more_iterations():
|
||
"""慢混合链：迭代次数显著变多，但结果仍要对得上特征向量（tol 是相邻差口径）。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.999, 0.001], [0.02, 0.98]])
|
||
r = markov_steady_state(P)
|
||
assert r.metrics["iterations"] > 100
|
||
assert np.allclose(r.data, [0.9523809524, 0.0476190476], atol=1e-6)
|
||
|
||
|
||
def test_markov_steady_state_accepts_numpy_integer():
|
||
"""§4.7：`np.int64` 不是 `int` 的子类，注解与运行期都要收（易错点 1）。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_steady_state(P, max_iter=np.int64(10000))
|
||
assert r.success
|
||
assert r.metrics["converged"] is True
|
||
|
||
|
||
def test_markov_steady_state_rejects_bad_transition():
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
with pytest.raises(ValueError):
|
||
markov_steady_state(np.ones((2, 3)) / 3)
|
||
with pytest.raises(ValueError):
|
||
markov_steady_state(np.array([[0.5, 0.4], [0.2, 0.8]]))
|
||
with pytest.raises(ValueError):
|
||
markov_steady_state(np.array([[1.5, -0.5], [0.0, 1.0]])) # 行和是 1，但含负元素
|
||
|
||
|
||
def test_markov_steady_state_non_convergence_is_success_false():
|
||
"""`max_iter` 不够时按"未收敛"处理：`success=False` 且 `converged=False`（口径 B）。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_steady_state(P, max_iter=3)
|
||
assert r.success is False
|
||
assert r.metrics["converged"] is False
|
||
assert r.metrics["iterations"] == 3
|
||
assert "未在最大迭代次数内收敛" in r.message
|
||
|
||
|
||
def test_markov_steady_state_rejects_bad_arguments():
|
||
"""取值错一律 `ValueError`（类型错要不要拆 `TypeError` 待定，见 STATUS 决策 57）。"""
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
with pytest.raises(ValueError):
|
||
markov_steady_state(P, max_iter=0)
|
||
with pytest.raises(ValueError):
|
||
markov_steady_state(P, tol=0)
|
||
|
||
|
||
def test_markov_steady_state_docstring_example_runs():
|
||
from mcm_toolkit.solve import markov_steady_state
|
||
|
||
P = np.array([[0.9, 0.1], [0.2, 0.8]])
|
||
r = markov_steady_state(P)
|
||
assert bool(np.allclose(r.data, [2 / 3, 1 / 3], atol=1e-6))
|
||
|
||
|
||
|
||
@pytest.mark.skip(reason=_NOT_IMPLEMENTED)
|
||
def test_knapsack_dp_known_optimum():
|
||
from mcm_toolkit.solve import knapsack_dp
|
||
|
||
r = knapsack_dp(weights=[1, 3, 4, 5], values=[1, 4, 5, 7], capacity=7)
|
||
assert np.isclose(r.metrics["total_value"], 9.0)
|
||
assert np.isclose(r.metrics["total_weight"], 7.0)
|
||
|
||
|
||
@pytest.mark.skip(reason=_NOT_IMPLEMENTED)
|
||
def test_linear_program_basic():
|
||
from mcm_toolkit.solve import linear_program
|
||
|
||
r = linear_program(c=[1, 2], A_ub=[[1, 1]], b_ub=[4], bounds=[(0, None)] * 2)
|
||
assert r.success
|
||
assert np.isclose(r.metrics["objective"], 4.0)
|
||
|
||
|
||
def test_mm1_metrics_closed_form():
|
||
"""与闭式公式逐项对拍（ρ = λ/μ，L = ρ/(1-ρ)，Wq = ρ/(μ-λ)，W = Wq + 1/μ）。"""
|
||
from mcm_toolkit.solve import mm1_metrics
|
||
|
||
arrival_rate, service_rate = 2, 3
|
||
rho = arrival_rate / service_rate
|
||
r = mm1_metrics(arrival_rate=arrival_rate, service_rate=service_rate)
|
||
|
||
assert r.success
|
||
assert np.isclose(r.metrics["rho"], rho, rtol=1e-12)
|
||
assert np.isclose(r.metrics["L"], rho / (1 - rho), rtol=1e-12)
|
||
assert np.isclose(r.metrics["Lq"], rho**2 / (1 - rho), rtol=1e-12)
|
||
assert np.isclose(r.metrics["W"], 1 / (service_rate - arrival_rate), rtol=1e-12)
|
||
assert np.isclose(r.metrics["Wq"], rho / (service_rate - arrival_rate), rtol=1e-12)
|
||
assert np.isclose(r.metrics["W"], r.metrics["Wq"] + 1 / service_rate, rtol=1e-12)
|
||
|
||
|
||
def test_mm1_rejects_unstable_system():
|
||
"""λ >= μ 时系统不稳定，应 success=False 而不是返回 inf（λ == μ 也要走这一支）。"""
|
||
from mcm_toolkit.solve import mm1_metrics
|
||
|
||
equal = mm1_metrics(arrival_rate=3, service_rate=3)
|
||
assert equal.success is False
|
||
assert "不稳定" in equal.message
|
||
assert equal.metrics == {}
|
||
|
||
assert mm1_metrics(arrival_rate=3, service_rate=2).success is False
|
||
|
||
|
||
def test_mm1_rejects_bad_arguments():
|
||
"""λ / μ 必须是正数：0、负数、bool、非数值都要在入口拦下。"""
|
||
from mcm_toolkit.solve import mm1_metrics
|
||
|
||
with pytest.raises(ValueError):
|
||
mm1_metrics(arrival_rate=0, service_rate=3)
|
||
with pytest.raises(ValueError):
|
||
mm1_metrics(arrival_rate=-1, service_rate=3)
|
||
with pytest.raises(TypeError):
|
||
mm1_metrics(arrival_rate=True, service_rate=3)
|
||
with pytest.raises(TypeError):
|
||
mm1_metrics(arrival_rate=2, service_rate="3") # type: ignore[arg-type] # 就是要传错
|
||
|
||
|
||
def test_mm1_docstring_example_runs():
|
||
from mcm_toolkit.solve import mm1_metrics
|
||
|
||
r = mm1_metrics(arrival_rate=2, service_rate=3)
|
||
assert round(r.metrics["rho"], 4) == 0.6667
|
||
|
||
|
||
def test_mm1_stability_boundary_is_success_not_exception():
|
||
"""λ == μ 不能靠 ZeroDivisionError 落到失败分支：必须是显式的"不稳定"口径。"""
|
||
from mcm_toolkit.solve import mm1_metrics
|
||
|
||
r = mm1_metrics(arrival_rate=5, service_rate=5)
|
||
assert r.success is False
|
||
assert "ZeroDivisionError" not in r.message
|
||
|
||
|
||
@pytest.mark.skip(reason=_NOT_IMPLEMENTED)
|
||
def test_simulated_annealing_finds_quadratic_minimum():
|
||
from mcm_toolkit.solve import simulated_annealing
|
||
|
||
f = lambda v: (v[0] - 3.0) ** 2 + (v[1] + 1.0) ** 2 # noqa: E731
|
||
r = simulated_annealing(f, x0=[0.0, 0.0], n_iter=5000, random_state=0)
|
||
assert np.allclose(r.params["x"], [3.0, -1.0], atol=0.1)
|