diff --git a/E03_StabSection/E03_StabSection_Engine_Prob.py b/E03_StabSection/E03_StabSection_Engine_Prob.py
new file mode 100644
index 000000000..fc6e913c0
--- /dev/null
+++ b/E03_StabSection/E03_StabSection_Engine_Prob.py
@@ -0,0 +1,441 @@
+"""E03 2D 확률 해석 · 민감도 (요구 351 · PLAN 537-1) — Slide2 · SLOPE/W Probabilistic 과 같은 생각.
+
+확률 변수 = 물성(c · φ · γ) · 지하수위(수위선 위아래 이동 m) · r_u — 정규 · 대수정규 · 절단(평균 아래 · 위
+거리 = Slide2 Relative Min · Max · 물리 한계 c ≥ 0 · 0 ≤ φ < 89 · γ > 0 · 0 ≤ r_u < 1 도 절단) · 변수 사이 상관
+(Slide2 c–φ Correlation · 표준정규에서 촐레스키) → 몬테카를로 · 라틴 하이퍼큐브 N 회
+→ 안전율 표본 → 평균 · 표준편차 · Pf = P(F < 1) · β 정규 (μ − 1)/σ · β 대수정규
+ln(μ/√(1+V²))/√ln(1+V²) · 히스토그램 · 누적분포 · P(F < 기준안전율).
+갈래: fixed = 결정론 임계면 고정(절편 기하 한 번 · 물성만 바꿔 묶음 평가 · 빠름 · Slide2 Global Minimum) ·
+search = 표본마다 임계면 다시 찾음(`Engine_Search.analyze` · 느림 · 프로세스 병렬 · Slide2 Overall Slope).
+민감도 = 변수 하나씩 범위(절단 끝 · 없으면 ±3σ) 50 등분 훑기(나머지 평균 · 임계면 고정 · Slide2 Sensitivity 와 같음)
+· 토네이도 = 끝값 F 폭 순.
+γ 변수 = γt 표본 · γsat 는 같은 차만큼 옮김(γsat − γt 그대로).
+search 일꾼 = `python -m E03_StabSection.E03_StabSection_Engine_Prob <job.json>`(서버 프로세스에서 풀을 열지 않음 ·
+`Worker` 머리글) · 한 줄마다 {"done": n} · 끝에 job["out"] 파일.
+NumPy · SciPy(stats) 만 · 계산 자리 = ③ 서버 단독.
+"""
+
+from __future__ import annotations
+
+import copy
+import json
+import math
+import sys
+import time
+from pathlib import Path
+
+import numpy as np
+from scipy.stats import norm
+
+from E03_StabSection import E03_StabSection_Engine_LEM as lem
+from E03_StabSection.E03_StabSection_Engine_Search import analyze
+from E03_StabSection.E03_StabSection_Engine_Station import MIN_DEPTH_M
+from E03_StabSection.E03_StabSection_Engine_Surface import (
+    circle_slices,
+    evaluate,
+    poly_slices,
+    three_point_center,
+)
+
+PROPS = ("c", "phi", "gamma", "water", "ru")
+#: 물리 한계(절단) — 지하수 이동은 한계 없음
+BOUNDS = {"c": (0.0, math.inf), "phi": (0.0, 89.0), "gamma": (1e-3, math.inf), "ru": (0.0, 0.99)}
+#: 절단 없을 때 민감도 훑기 범위(±σ 배)
+SWEEP_SIGMA = 3.0
+SWEEP_POINTS = 51
+HIST_BINS = 30
+CDF_POINTS = 101
+N_SLICES = 50
+
+
+# ---- 분포 · 표본 --------------------------------------------------------------
+def _bounds(v: dict) -> tuple[float, float]:
+    """절단 구간 [a, b] — 평균 ∓ 아래 · 위 거리 ∩ 물리 한계."""
+    lo_phys, hi_phys = BOUNDS.get(v["prop"], (-math.inf, math.inf))
+    mean = float(v["mean"])
+    a = mean - float(v["lo"]) if v.get("lo") is not None else -math.inf
+    b = mean + float(v["hi"]) if v.get("hi") is not None else math.inf
+    return max(a, lo_phys), min(b, hi_phys)
+
+
+def _lognormal(v: dict) -> tuple[float, float]:
+    """대수정규 (μ_ln, σ_ln) — 평균 · 표준편차 그대로 맞춤."""
+    mean, std = float(v["mean"]), float(v["std"])
+    s2 = math.log1p((std / mean) ** 2)
+    return math.log(mean) - 0.5 * s2, math.sqrt(s2)
+
+
+def ppf(v: dict, u: np.ndarray) -> np.ndarray:
+    """균등 u(0~1) → 변수 값(절단 분포 역함수)."""
+    u = np.asarray(u, float)
+    mean, std = float(v["mean"]), float(v.get("std") or 0.0)
+    if std <= 0:
+        return np.full(u.shape, mean)
+    a, b = _bounds(v)
+    if v.get("dist") == "lognormal" and v["prop"] != "water":
+        if mean <= 0:
+            raise ValueError(f"대수정규는 평균 > 0: {v}")
+        mu, s = _lognormal(v)
+        za = (math.log(a) - mu) / s if a > 0 else -math.inf
+        zb = (math.log(b) - mu) / s if b < math.inf else math.inf
+        pa, pb = norm.cdf(za), norm.cdf(zb)
+        return np.exp(mu + s * norm.ppf(pa + u * (pb - pa)))
+    pa, pb = norm.cdf((a - mean) / std), norm.cdf((b - mean) / std)
+    return mean + std * norm.ppf(pa + u * (pb - pa))
+
+
+def sample(
+    variables: list[dict],
+    n: int,
+    *,
+    sampling: str = "lhs",
+    seed: int = 0,
+    correlations: list | None = None,
+) -> np.ndarray:
+    """표본 (n, 변수 수) — lhs = 라틴 하이퍼큐브(변수마다 n 칸 하나씩 · 칸 차례 섞음) · mc = 단순 무작위.
+
+    correlations = [(i, j, ρ), …] — 표준정규 칸에서 촐레스키로 섞은 뒤 분포로(LHS 층은 섞인 변수만 흐려짐)."""
+    rng = np.random.default_rng(int(seed))
+    k = len(variables)
+    if sampling == "lhs":
+        order = np.argsort(rng.random((n, k)), axis=0)
+        u = (order + rng.random((n, k))) / n
+    else:
+        u = rng.random((n, k))
+    u = np.clip(u, 1e-12, 1 - 1e-12)
+    if correlations:
+        corr = np.eye(k)
+        for i, j, rho in correlations:
+            corr[i, j] = corr[j, i] = float(rho)
+        u = norm.cdf(norm.ppf(u) @ np.linalg.cholesky(corr).T)
+        u = np.clip(u, 1e-12, 1 - 1e-12)
+    return np.column_stack([ppf(v, u[:, j]) for j, v in enumerate(variables)]) if k else u
+
+
+def sweep_range(v: dict) -> tuple[float, float]:
+    """민감도 훑기 범위 — 절단 끝 · 없으면 평균 ±3σ(물리 한계 안)."""
+    mean, std = float(v["mean"]), float(v.get("std") or 0.0)
+    a, b = _bounds(v)
+    lo = a if math.isfinite(a) and v.get("lo") is not None else mean - SWEEP_SIGMA * std
+    hi = b if math.isfinite(b) and v.get("hi") is not None else mean + SWEEP_SIGMA * std
+    lo_phys, hi_phys = BOUNDS.get(v["prop"], (-math.inf, math.inf))
+    return max(lo, lo_phys), min(hi, hi_phys)
+
+
+# ---- 통계 -------------------------------------------------------------------
+def statistics(F: np.ndarray, required: float | None = None, bins: int = HIST_BINS) -> dict:
+    """안전율 표본 → 평균 · 표준편차 · Pf · β 정규 · 대수정규 · 히스토그램 · 누적분포(해 없는 표본 = 셈 밖)."""
+    F = np.asarray(F, float)
+    ok = F[np.isfinite(F)]
+    out = {"n": int(len(F)), "n_valid": int(len(ok))}
+    if len(ok) < 2:
+        return {**out, "mean": None, "std": None, "pf": None, "beta_normal": None}
+    mean, std = float(ok.mean()), float(ok.std(ddof=1))
+    V = std / mean if mean > 0 else math.nan
+    beta_ln = (
+        math.log(mean / math.sqrt(1 + V * V)) / math.sqrt(math.log1p(V * V))
+        if mean > 0 and V > 0
+        else None
+    )
+    counts, edges = np.histogram(ok, bins=bins)
+    p = np.linspace(0.0, 1.0, CDF_POINTS)
+    out.update(
+        mean=mean,
+        std=std,
+        cov=V if math.isfinite(V) else None,
+        min=float(ok.min()),
+        max=float(ok.max()),
+        pf=float(np.mean(ok < 1.0)),
+        n_fail=int(np.sum(ok < 1.0)),
+        beta_normal=(mean - 1.0) / std if std > 0 else None,
+        beta_lognormal=beta_ln,
+        required=required,
+        p_below_required=float(np.mean(ok < required)) if required else None,
+        histogram={"edges": edges.tolist(), "counts": counts.tolist()},
+        cdf={"p": p.tolist(), "F": np.quantile(ok, p).tolist()},
+    )
+    return out
+
+
+# ---- 변수 → 모형 ---------------------------------------------------------------
+def perturb(model: dict, variables: list[dict], row) -> dict:
+    """모형 사본 — 표본 한 줄을 물성 · 지하수위에 넣음(search 갈래)."""
+    out = dict(model)
+    mats = {m["id"]: dict(m) for m in model["materials"]}
+    dz = 0.0
+    for v, x in zip(variables, row):
+        if v["prop"] == "water":
+            dz += float(x)
+            continue
+        if v["prop"] == "ru":
+            gw = model.get("groundwater") or {}
+            if gw.get("mode") == "ru":
+                out["groundwater"] = {**gw, "ru": float(x)}
+            continue
+        m = mats[v["material"]]
+        if v["prop"] == "gamma":
+            gt = float(m["gamma_t"])
+            m["gamma_sat"] = float(m.get("gamma_sat", gt)) + float(x) - gt
+            m["gamma_t"] = float(x)
+        else:
+            m[v["prop"]] = float(x)
+    out["materials"] = list(mats.values())
+    gw = out.get("groundwater") or {}
+    if dz and gw.get("mode") == "line":
+        out["groundwater"] = {**gw, "line": [[x, z + dz] for x, z in gw["line"]]}
+    return out
+
+
+class Fixed:
+    """임계면 고정 — 모형 · 면 하나를 한 번 준비해 표본마다 물성 칸만 갈아 끼우고 절편 묶음으로 평가."""
+
+    def __init__(self, model: dict, surface: dict, direction: int, n_slices: int = N_SLICES):
+        self.prep = lem.Prepared(model, direction)
+        self.n = n_slices
+        s = self.prep.direction
+        self.index = {m["id"]: i for i, m in enumerate(self.prep.materials)}
+        if surface["kind"] == "circle":
+            xs = sorted((s * surface["entry"][0], s * surface["exit"][0]))
+            self.geom = (
+                "circle",
+                s * surface["center"][0],
+                surface["center"][1],
+                surface["radius"],
+            )
+            self.geom += tuple(xs)
+        else:
+            pts = sorted((s * x, z) for x, z in surface["points"])
+            self.geom = ("polyline", np.array([p[0] for p in pts]), np.array([p[1] for p in pts]))
+
+    def _prep(self, variables: list[dict], row) -> lem.Prepared:
+        base = self.prep
+        p = copy.copy(base)
+        gt, gs = base.gamma_t.copy(), base.gamma_sat.copy()
+        c, tp = base.c.copy(), base.tanphi.copy()
+        dz = 0.0
+        for v, x in zip(variables, row):
+            if v["prop"] == "water":
+                dz += float(x)
+                continue
+            if v["prop"] == "ru":
+                if base.gw_mode == "ru":
+                    p.ru = float(x)
+                continue
+            i = self.index[v["material"]]
+            if v["prop"] == "c":
+                c[i] = x
+            elif v["prop"] == "phi":
+                tp[i] = math.tan(math.radians(x))
+            else:
+                gs[i] += x - gt[i]
+                gt[i] = x
+        p.gamma_t, p.gamma_sat, p.c, p.tanphi = gt, gs, c, tp
+        p.table = base.table.copy()
+        p.table[:, 4, :] = gt[base.mat]
+        p.table[:, 5, :] = gs[base.mat] - gt[base.mat]
+        if dz and base.gw_mode == "line":
+            p.water_z = base.water_z + dz
+        return p
+
+    def _slices(self, prep: lem.Prepared) -> tuple[dict, np.ndarray]:
+        g = self.geom
+        if g[0] == "circle":
+            _, xc, zc, R, xE, xX = g
+            sl = circle_slices(prep, [xc], [zc], [R], [xE], [xX], self.n)
+            return sl, np.array([[xc, zc]])
+        vx, vz = g[1][None, :], g[2][None, :]
+        return poly_slices(prep, vx, vz, self.n), three_point_center(vx, vz)
+
+    def F(self, variables: list[dict], X: np.ndarray, method: str) -> np.ndarray:
+        """표본 줄마다 안전율(해 없음 = nan)."""
+        X = np.atleast_2d(np.asarray(X, float))
+        parts = [self._slices(self._prep(variables, row)) for row in X]
+        sl = _stack([p[0] for p in parts])
+        center = np.vstack([p[1] for p in parts])
+        mp_seed = None
+        if method == "morgenstern_price":  # describe 와 같은 출발 참고값(같은 면 Spencer)
+            mp_seed = evaluate(sl, center, ("spencer",))["spencer"]["F"]
+            sl.pop("_pre", None)
+        F = evaluate(sl, center, (method,), mp_seed=mp_seed)[method]["F"]
+        return np.where(sl["valid"], F, np.nan)
+
+
+def _stack(parts: list[dict]) -> dict:
+    """한 줄짜리 절편 묶음들 → (S, …) 하나(배열 칸만 잇고 스칼라는 그대로)."""
+    out = {}
+    for key, value in parts[0].items():
+        if key == "_pre":
+            continue
+        if isinstance(value, np.ndarray) and value.ndim >= 1 and value.shape[0] == 1:
+            out[key] = np.concatenate([p[key] for p in parts])
+        else:
+            out[key] = value
+    return out
+
+
+def _chunks(n: int, size: int):
+    for i in range(0, n, size):
+        yield slice(i, min(i + size, n))
+
+
+def sensitivity(
+    fixed: Fixed, variables: list[dict], method: str, points: int = SWEEP_POINTS
+) -> dict:
+    """변수 하나씩 범위 훑기(나머지 평균) · 토네이도 = 끝값 F 폭 큰 순."""
+    means = np.array([float(v["mean"]) for v in variables])
+    rows, spans = [], []
+    for j, v in enumerate(variables):
+        lo, hi = sweep_range(v)
+        xs = np.linspace(lo, hi, points) if hi > lo else np.full(points, means[j])
+        block = np.tile(means, (points, 1))
+        block[:, j] = xs
+        rows.append(block)
+        spans.append(xs)
+    F0 = float(fixed.F(variables, means[None, :], method)[0])
+    if not rows:
+        return {"F_mean": F0, "curves": [], "tornado": []}
+    F = fixed.F(variables, np.vstack(rows), method).reshape(len(variables), points)
+    num = lambda a: [float(x) if np.isfinite(x) else None for x in a]  # noqa: E731
+    curves = []
+    for v, xs, f in zip(variables, spans, F):
+        ends = [f[0], f[-1]]
+        swing = abs(ends[1] - ends[0]) if all(np.isfinite(ends)) else None
+        curves.append(
+            {
+                "name": v["name"],
+                "x": num(xs),
+                "F": num(f),
+                "F_lo": num([ends[0]])[0],
+                "F_hi": num([ends[1]])[0],
+                "swing": swing,
+            }
+        )
+    tornado = sorted(
+        (c for c in curves if c["swing"] is not None), key=lambda c: c["swing"], reverse=True
+    )
+    return {"F_mean": F0 if math.isfinite(F0) else None, "curves": curves, "tornado": [
+        {"name": c["name"], "F_lo": c["F_lo"], "F_hi": c["F_hi"], "swing": c["swing"]}
+        for c in tornado
+    ]}  # fmt: skip
+
+
+# ---- 매회 탐색(search) ------------------------------------------------------------
+_CTX: dict = {}
+
+
+def _init(ctx: dict) -> None:
+    _CTX.clear()
+    _CTX.update(ctx)
+
+
+def _search_one(row) -> float:
+    c = _CTX
+    m = c["method"]
+    methods = tuple(dict.fromkeys(("bishop", m)))
+    res = analyze(
+        perturb(c["model"], c["variables"], row),
+        direction=c["direction"],
+        limits=c["limits"],
+        methods=methods,
+        noncircular=m not in ("ordinary", "bishop"),  # 비원호 Fellenius · Bishop = 참고 칸
+        min_depth=MIN_DEPTH_M,
+    )
+    best = res["methods"][m]["min"]
+    return float(best["F"]) if best and best.get("F") is not None else math.nan
+
+
+def search_F(ctx: dict, X: np.ndarray, workers: int = 1, emit=None) -> np.ndarray:
+    """표본마다 임계면 다시 찾음 — workers > 1 = 프로세스 풀(서버 밖 일꾼에서만)."""
+    out = np.full(len(X), np.nan)
+    if workers <= 1:
+        _init(ctx)
+        for i, row in enumerate(X):
+            out[i] = _search_one(row)
+            emit and emit(i + 1)
+        return out
+    from concurrent.futures import ProcessPoolExecutor, as_completed
+
+    with ProcessPoolExecutor(max_workers=workers, initializer=_init, initargs=(ctx,)) as pool:
+        futures = {pool.submit(_search_one, row): i for i, row in enumerate(X)}
+        for done, future in enumerate(as_completed(futures), 1):
+            out[futures[future]] = future.result()
+            emit and emit(done)
+    return out
+
+
+# ---- 한 번 해석 ---------------------------------------------------------------
+def run(
+    model: dict,
+    *,
+    surface: dict,
+    direction: int,
+    limits=None,
+    variables: list[dict],
+    correlations: list | None = None,
+    n: int = 1000,
+    sampling: str = "lhs",
+    seed: int = 0,
+    mode: str = "fixed",
+    method: str = "bishop",
+    required: float | None = None,
+    workers: int = 1,
+    emit=None,
+) -> dict:
+    """확률 해석 + 민감도 한 벌. emit(끝난 표본 수) = 진행률."""
+    t0 = time.perf_counter()
+    X = sample(variables, n, sampling=sampling, seed=seed, correlations=correlations)
+    fixed = Fixed(model, surface, direction)
+    t1 = time.perf_counter()
+    if mode == "search":
+        ctx = {
+            "model": model,
+            "variables": variables,
+            "direction": direction,
+            "limits": limits,
+            "method": method,
+        }
+        F = search_F(ctx, X, workers, emit)
+    else:
+        F = np.full(n, np.nan)
+        for part in _chunks(n, 250):
+            F[part] = fixed.F(variables, X[part], method)
+            emit and emit(part.stop)
+    t2 = time.perf_counter()
+    sens = sensitivity(fixed, variables, method)
+    t3 = time.perf_counter()
+    realized = [
+        {"name": v["name"], "mean": float(X[:, j].mean()), "std": float(X[:, j].std(ddof=1))}
+        for j, v in enumerate(variables)
+    ]
+    return {
+        "mode": mode,
+        "method": method,
+        "n": n,
+        "sampling": sampling,
+        "seed": seed,
+        "variables": variables,
+        "correlations": correlations or [],
+        "realized": realized,
+        "stats": statistics(F, required),
+        "F": [round(float(f), 4) if np.isfinite(f) else None for f in F],
+        "sensitivity": sens,
+        "elapsed_s": {
+            "prepare": round(t1 - t0, 3),
+            "samples": round(t2 - t1, 3),
+            "sensitivity": round(t3 - t2, 3),
+            "total": round(t3 - t0, 3),
+        },
+    }
+
+
+def main() -> None:
+    """일꾼 — job = run() 인자 + out(결과 파일) · 한 줄마다 진행."""
+    path = Path(sys.argv[1])
+    job = json.loads(path.read_text(encoding="utf-8"))
+    out = Path(job.pop("out"))
+    result = run(**job, emit=lambda k: print(json.dumps({"done": k}), flush=True))
+    out.write_text(json.dumps(result, ensure_ascii=False), encoding="utf-8")
+
+
+if __name__ == "__main__":
+    main()
diff --git a/resources/tester/test_537_1_e03_prob.py b/resources/tester/test_537_1_e03_prob.py
new file mode 100644
index 000000000..af806ae21
--- /dev/null
+++ b/resources/tester/test_537_1_e03_prob.py
@@ -0,0 +1,227 @@
+# -*- coding: utf-8 -*-
+"""537-1 E03 2D 확률 해석 · 민감도 — 엔진(`E03_StabSection_Engine_Prob`).
+
+- 무한사면 해석해: F = c/(γ z sinβ cosβ) + (1 − m γw/γ) tanφ/tanβ · c · φ 정규(MIT OCW 1.151 App. 17 ·
+  z 5 · m 0.5 · β 35° · c N(25, 5) · φ N(30°, 7.5°)) — Pf · μ · σ 정확값 = φ 위 수치 적분 · LHS 1만 ±2 %
+  · c · tanφ 정규면 F 도 정규 → Pf = Φ(−β) · 절편 경로(임계면 고정 · 건조 무한사면) = 닫힌 꼴.
+- Slide2 Tutorial 08(Probabilistic Analysis · LHS 1000 · Global Minimum · Bishop): FS 평균 1.145 · σ 0.1206
+  · PF 10.2 % · RI 정규 1.204 · 대수정규 1.238.
+- Slide2 Verification #36(Li & Lumb 1987 · Hassan & Wolff 1999 · r_u): 결정론 1.340(H&W 1.334) · 평균 1.355 ·
+  PF 1.5 % · RI 대수정규 2.482(H&W 2.336).
+- 고정 시드 재현 · LHS 층 · 절단 · 대수정규 평균 · 민감도(평균 = 결정론 F · 토네이도 순) · 매회 탐색 갈래.
+
+    ./venv/Scripts/python.exe -m pytest resources/tester/test_537_1_e03_prob.py -q
+"""
+
+from __future__ import annotations
+
+import math
+import sys
+from pathlib import Path
+
+import numpy as np
+import pytest
+from scipy import integrate
+from scipy.stats import norm
+
+ROOT = Path(__file__).resolve().parents[2]
+sys.path.insert(0, str(ROOT))
+
+from E03_StabSection import E03_StabSection_Engine_Prob as prob  # noqa: E402
+from E03_StabSection.E03_StabSection_Engine_Search import analyze  # noqa: E402
+
+
+def var(name, prop, mean, std, rel=None, mat="s", dist="normal"):
+    return {
+        "name": name, "material": mat, "prop": prop, "dist": dist,
+        "mean": mean, "std": std, "lo": rel, "hi": rel,
+    }  # fmt: skip
+
+
+def near(v, ref, tol):
+    assert v is not None and abs(v - ref) <= tol * abs(ref), (v, ref)
+
+
+# ---- 무한사면 해석해 ---------------------------------------------------------------
+Z, BETA, G, GW = 5.0, math.radians(35.0), 20.0, 9.81
+C_PHI = [var("c", "c", 25.0, 5.0), var("phi", "phi", 30.0, 7.5)]
+
+
+def _ab(m):
+    return 1 / (G * Z * math.sin(BETA) * math.cos(BETA)), (1 - m * GW / G) / math.tan(BETA)
+
+
+def _exact(m):
+    """φ 위 수치 적분 — Pf = ∫ Φ((1 − b tanφ − a μc)/(a σc)) f(φ) dφ · 평균 · 표준편차."""
+    a, b = _ab(m)
+    f = lambda p: norm.pdf(p, 30.0, 7.5)  # noqa: E731
+    F = lambda p: a * 25.0 + b * math.tan(math.radians(p))  # noqa: E731
+    pf = integrate.quad(lambda p: norm.cdf((1 - F(p)) / (a * 5.0)) * f(p), 0, 89, limit=200)[0]
+    mu = integrate.quad(lambda p: F(p) * f(p), 0, 89)[0]
+    m2 = integrate.quad(lambda p: (F(p) ** 2 + (a * 5.0) ** 2) * f(p), 0, 89)[0]
+    return pf, mu, math.sqrt(m2 - mu * mu)
+
+
+def test_infinite_slope_closed_form_pf_beta():
+    """MIT App. 17 무한사면 — LHS 1만 Pf · β ±2 % · 평균 · 표준편차 ±1 %(φ 정규 · 각도)."""
+    a, b = _ab(0.5)
+    X = prob.sample(C_PHI, 10000, seed=0)
+    s = prob.statistics(a * X[:, 0] + b * np.tan(np.radians(X[:, 1])))
+    pf, mu, sd = _exact(0.5)
+    near(pf, 0.2252, 0.01)  # 같은 문헌 MC 2×10⁶
+    near(s["pf"], pf, 0.02)
+    near(s["beta_normal"], (mu - 1) / sd, 0.02)
+    near(s["mean"], mu, 0.01)
+    near(s["std"], sd, 0.01)
+    V = sd / mu
+    near(
+        s["beta_lognormal"],
+        math.log(mu / math.sqrt(1 + V * V)) / math.sqrt(math.log1p(V * V)),
+        0.02,
+    )
+
+
+def test_infinite_slope_normal_tanphi_exact():
+    """c · tanφ 정규 · 독립 → F 정규 → Pf = Φ(−β) 정확(LHS 2만 ±2 %)."""
+    a, b = _ab(0.5)
+    vs = [var("c", "c", 25.0, 5.0), {**var("t", "c", 0.57735, 0.11547), "prop": "x"}]
+    X = prob.sample(vs, 20000, seed=0)
+    s = prob.statistics(a * X[:, 0] + b * X[:, 1])
+    mu = a * 25.0 + b * 0.57735
+    sd = math.hypot(a * 5.0, b * 0.11547)
+    near(mu, 1.154412, 1e-5)
+    near(sd, 0.163756, 1e-4)
+    near(s["beta_normal"], (mu - 1) / sd, 0.02)
+    near(s["pf"], norm.cdf(-(mu - 1) / sd), 0.02)
+
+
+def _infinite_model():
+    """건조 무한사면 β 35° · 지표 아래 5 m 나란한 면 · 모형 양 끝까지(끝 효과 없음)."""
+    t, L = math.tan(BETA), 100.0
+    model = {
+        "materials": [{"id": "s", "name": "s", "gamma_t": G, "c": 25.0, "phi": 30.0}],
+        "regions": [{"material_id": "s", "polygon": [[0, 0], [L, -t * L], [L, -t * L - 20], [0, -20]]}],
+    }  # fmt: skip
+    return model, {"kind": "polyline", "points": [[0, -Z], [L, -t * L - Z]]}
+
+
+@pytest.mark.parametrize("method", ["ordinary", "bishop", "janbu"])
+def test_fixed_surface_slices_equal_closed_form(method):
+    """임계면 고정 절편 경로 = 무한사면 닫힌 꼴(물성 바꾼 줄마다) · LHS 1만 Pf ±2 %."""
+    model, surface = _infinite_model()
+    fixed = prob.Fixed(model, surface, 1)
+    a, b = _ab(0.0)
+    X = np.array([[25.0, 30.0], [10.0, 20.0], [40.0, 35.0]])
+    got = fixed.F(C_PHI, X, method)
+    assert got == pytest.approx(a * X[:, 0] + b * np.tan(np.radians(X[:, 1])), rel=1e-9)
+    out = prob.run(model, surface=surface, direction=1, variables=C_PHI, n=10000, method=method)
+    near(out["stats"]["pf"], _exact(0.0)[0], 0.02)
+
+
+# ---- Slide2 대조 -------------------------------------------------------------------
+T08 = {
+    "materials": [{"id": "s", "name": "s", "gamma_t": 19.0, "c": 5.0, "phi": 30.0}],
+    "regions": [{"material_id": "s", "polygon": [
+        [0, 0], [130, 0], [130, 50], [80, 50], [50, 30], [0, 30]]}],
+}  # fmt: skip
+T08_VARS = [var("c", "c", 5, 1, 3), var("phi", "phi", 30, 3, 9), var("g", "gamma", 19, 0.5, 1.5)]
+
+
+def _surface(model, method="bishop"):
+    return analyze(model, direction=-1, noncircular=False, methods=(method,))["methods"][method][
+        "min"
+    ]
+
+
+def test_slide2_tutorial08_global_minimum():
+    """Tutorial 08 — LHS 1000 · 고정 시드 셋 평균: FS 평균 ±1 % · σ ±3 % · PF ±2 %p · RI ±3 %."""
+    surface = _surface(T08)
+    near(surface["F"], 1.142, 0.01)
+    runs = [
+        prob.run(T08, surface=surface, direction=-1, variables=T08_VARS, n=1000, seed=k)["stats"]
+        for k in range(3)
+    ]
+    avg = {
+        k: float(np.mean([r[k] for r in runs]))
+        for k in ("mean", "std", "pf", "beta_normal", "beta_lognormal")
+    }
+    near(avg["mean"], 1.145, 0.01)
+    near(avg["std"], 0.1206, 0.03)
+    assert abs(avg["pf"] - 0.102) <= 0.02
+    near(avg["beta_normal"], 1.204, 0.03)
+    near(avg["beta_lognormal"], 1.238, 0.03)
+
+
+def test_slide2_36_li_lumb_ru():
+    """#36 r_u 확률 변수 · MC 1만 — 결정론 ±1 % · 평균 ±1.5 % · PF ±1 %p · RI 대수정규 = Hassan & Wolff ±3 %."""
+    model = {
+        "materials": [{"id": "s", "name": "s", "gamma_t": 18.0, "c": 18.0, "phi": 30.0}],
+        "regions": [{"material_id": "s", "polygon": [[0, 0], [20, 0], [20, 15], [15, 15], [5, 5], [0, 5]]}],
+        "groundwater": {"mode": "ru", "ru": 0.2},
+    }  # fmt: skip
+    surface = _surface(model)
+    near(surface["F"], 1.340, 0.01)
+    vs = [
+        var("c", "c", 18, 3.6), var("phi", "phi", 30, 3), var("g", "gamma", 18, 0.9),
+        var("ru", "ru", 0.2, 0.02, mat=None),
+    ]  # fmt: skip
+    s = prob.run(model, surface=surface, direction=-1, variables=vs, n=10000, sampling="mc")[
+        "stats"
+    ]
+    near(s["mean"], 1.355, 0.015)
+    assert abs(s["pf"] - 0.015) <= 0.01
+    near(s["beta_lognormal"], 2.336, 0.03)
+
+
+# ---- 표본 · 재현 · 민감도 · 갈래 ---------------------------------------------------------
+def test_seed_reproducible_and_lhs_strata():
+    """같은 시드 = 같은 표본 · 다른 시드 = 다름 · LHS = 변수마다 n 칸에 하나씩 · 절단 안 · 상관."""
+    a = prob.sample(T08_VARS, 500, seed=7)
+    assert np.array_equal(a, prob.sample(T08_VARS, 500, seed=7))
+    assert not np.array_equal(a, prob.sample(T08_VARS, 500, seed=8))
+    for j, v in enumerate(T08_VARS):
+        lo, hi = v["mean"] - v["lo"], v["mean"] + v["hi"]
+        assert a[:, j].min() >= lo and a[:, j].max() <= hi
+        pa, pb = norm.cdf(-3), norm.cdf(3)
+        u = (norm.cdf((a[:, j] - v["mean"]) / v["std"]) - pa) / (pb - pa)
+        assert sorted(np.floor(u * 500).astype(int)) == list(range(500))
+    rho = prob.sample(C_PHI, 5000, seed=1, correlations=[(0, 1, -0.5)])
+    assert np.corrcoef(rho.T)[0, 1] == pytest.approx(-0.5, abs=0.03)
+
+
+def test_truncation_and_lognormal_moments():
+    """대수정규 = 평균 · σ 그대로 · 음수 없음 · 정규 c 는 물리 한계 0 에서 절단."""
+    ln = prob.sample([var("c", "c", 10.0, 5.0, dist="lognormal")], 20000, seed=2)[:, 0]
+    near(ln.mean(), 10.0, 0.01)
+    near(ln.std(), 5.0, 0.03)
+    assert ln.min() > 0
+    nc = prob.sample([var("c", "c", 5.0, 5.0)], 5000, seed=2)[:, 0]
+    assert nc.min() >= 0
+
+
+def test_sensitivity_and_tornado():
+    """민감도 — 평균 줄 F = 결정론 F · 곡선 51 점 · c ↑ = F ↑ · 토네이도 = 폭 큰 순."""
+    surface = _surface(T08)
+    out = prob.run(T08, surface=surface, direction=-1, variables=T08_VARS, n=200)
+    sens = out["sensitivity"]
+    assert sens["F_mean"] == pytest.approx(surface["F"], rel=1e-9)
+    curves = {c["name"]: c for c in sens["curves"]}
+    assert len(curves["c"]["F"]) == 51 and curves["c"]["x"][0] == pytest.approx(2.0)
+    assert np.all(np.diff(curves["c"]["F"]) > 0) and np.all(np.diff(curves["phi"]["F"]) > 0)
+    swings = [t["swing"] for t in sens["tornado"]]
+    assert swings == sorted(swings, reverse=True) and sens["tornado"][0]["name"] in ("c", "phi")
+    assert out["stats"]["n"] == 200 and len(out["F"]) == 200
+
+
+def test_search_mode_matches_fixed_at_low_variance():
+    """매회 탐색 — 변동이 아주 작으면 고정면과 같은 F(±1 %) · 표본 수 · 시간 칸."""
+    surface = _surface(T08)
+    tiny = [var("c", "c", 5, 0.01), var("phi", "phi", 30, 0.01)]
+    out = prob.run(T08, surface=surface, direction=-1, variables=tiny, n=4, mode="search")
+    for f in out["F"]:
+        near(f, surface["F"], 0.01)
+    assert out["mode"] == "search" and out["elapsed_s"]["samples"] > 0
+
+
+if __name__ == "__main__":
+    sys.exit(pytest.main([__file__, "-q"]))