auto: 2026-10-26 03:18 (ESD_LAPTOP)

E04_Stab3D/E04_Stab3D_Engine_Search.py

@@ -167,30 +167,43 @@ def trial_columns(
     hard = strata.hard_top[cy, cx]
 
     # 2) 부분 기둥 — sub×sub 점
-    n_act = np.zeros(len(ci))
-    acc = {k: np.zeros(len(ci)) for k in ("depth", "zs", "gx", "gy", "px", "py")}
     fr = ((np.arange(sub) + 0.5) / sub - 0.5) if sub > 1 else np.zeros(1)
-    for fy in fr:
-        for fx in fr:
+    nf = nw * sub  # 점 격자 한 변
+    win_y, win_x = (ci // nw) % nw, ci % nw
+    pts: dict[str, list] = {k: [] for k in ("c", "f", "depth", "zs", "gx", "gy", "px", "py")}
+    for ka, fy in enumerate(fr):
+        for kb, fx in enumerate(fr):
             xs = strata.x[0] + (cx + fx) * dx
             ys = strata.y[0] + (cy + fy) * dx
             rem, z_s, gx, gy = surface(tt, xs, ys)
             zg = strata.ground_at(cy, cx, fy, fx)
             z_b = np.maximum(z_s, hard)
-            act = (rem > 0) & (z_b < zg)
-            n_act += act
-            acc["depth"] += np.where(act, zg - z_b, 0.0)
-            for key, val in (("zs", z_s), ("gx", gx), ("gy", gy), ("px", xs), ("py", ys)):
-                acc[key] += np.where(act, val, 0.0)
-    keep = n_act > 0
-    ci, tt, cy, cx, n_act = ci[keep], tt[keep], cy[keep], cx[keep], n_act[keep]
-    acc = {k: v[keep] for k, v in acc.items()}
-    active = np.zeros(shape, bool)
-    active.ravel()[ci] = True
-    labels, n_lab = ndimage.label(active, structure=_STRUCT)
-    if n_lab == 0:
+            k = np.flatnonzero((rem > 0) & (z_b < zg))
+            pts["c"].append(k)
+            pts["f"].append((tt[k] * nf + win_y[k] * sub + ka) * nf + win_x[k] * sub + kb)
+            for key, val in (
+                ("depth", zg - z_b),
+                ("zs", z_s),
+                ("gx", gx),
+                ("gy", gy),
+                ("px", xs),
+                ("py", ys),
+            ):
+                pts[key].append(val[k])
+    pt = {key: np.concatenate(v) for key, v in pts.items()}
+    n_pt = len(pt["c"])
+    if not n_pt:
         return empty
-    lab = labels.ravel()[ci]
+    lab_pt, n_lab = _split_pieces(strata, pt, tt, cy, cx, hard, n_t * nf * nf, nf, sub, surface)
+    # 기둥 = (셀, 덩이) 짝 — 두 덩이가 한 셀을 나누면 기둥 둘
+    key, inv = np.unique(pt["c"] * n_lab + lab_pt, return_inverse=True)
+    cell_k, lab = key // n_lab, key % n_lab + 1  # 덩이 번호 1.. (0 = 없음)
+    n_act = np.bincount(inv, minlength=len(key)).astype(float)
+    acc = {
+        k: np.bincount(inv, weights=pt[k], minlength=len(key))
+        for k in ("depth", "zs", "gx", "gy", "px", "py")
+    }
+    tt, cy, cx = tt[cell_k], cy[cell_k], cx[cell_k]
     a_h = dx * dx
     depth_cell = acc["depth"] / (sub * sub)
     vol = np.bincount(lab, weights=depth_cell, minlength=n_lab + 1) * a_h
@@ -262,6 +275,51 @@ def _part(job: tuple, k: dict | None = None) -> tuple:
     return light, trial, big, res
 
 
+def _split_pieces(strata, pt, tt, cy, cx, hard, n_fine, nf, sub, surface):
+    """활동면이 지표를 자른 점들을 덩이로 — 이웃 점(4 방향) 사이 가운데 점도 땅속일 때만 이음.
+
+    한 셀 안에서도 덩이가 갈림(발끝을 살짝 비켜 지나는 구 = 사면 몸통 · 발끝 앞 접시 따로).
+    반환 (점별 덩이 번호 0.., 덩이 수).
+    """
+    from scipy.sparse import coo_matrix
+    from scipy.sparse.csgraph import connected_components
+
+    n_pt = len(pt["c"])
+    where = np.full(n_fine, -1, np.int32)  # 점 격자 칸 → 점 번호
+    where[pt["f"]] = np.arange(n_pt, dtype=np.int32)
+    fr = ((np.arange(sub) + 0.5) / sub - 0.5) if sub > 1 else np.zeros(1)
+    half = 0.5 / sub
+    rows, cols_ = [], []
+    for step, axis in ((1, "x"), (nf, "y")):
+        pos = pt["f"] % nf if axis == "x" else (pt["f"] // nf) % nf
+        i = np.flatnonzero(pos < nf - 1)  # 줄 · 창 끝 넘김 막음
+        q = where[pt["f"][i] + step]
+        i, q = i[q >= 0], q[q >= 0]
+        c = pt["c"][i]
+        # 가운데 점 — 점 i 의 셀에서 (fy, fx) + half
+        ka = ((pt["f"][i] // nf) % nf) % sub
+        kb = (pt["f"][i] % nf) % sub
+        mx = pt["px"][i] + (half * strata.dx if axis == "x" else 0.0)
+        my = pt["py"][i] + (half * strata.dx if axis == "y" else 0.0)
+        rem, z_s, _, _ = surface(tt[c], mx, my)
+        zg = np.empty(len(i))
+        for va in range(sub):
+            for vb in range(sub):
+                m = (ka == va) & (kb == vb)
+                if m.any():
+                    fy = fr[va] + (half if axis == "y" else 0.0)
+                    fx = fr[vb] + (half if axis == "x" else 0.0)
+                    zg[m] = strata.ground_at(cy[c[m]], cx[c[m]], fy, fx)
+        ok = (rem > 0) & (np.maximum(z_s, hard[c]) < zg)
+        rows.append(i[ok])
+        cols_.append(q[ok])
+    r_ = np.concatenate(rows)
+    c_ = np.concatenate(cols_)
+    graph = coo_matrix((np.ones(len(r_), np.int8), (r_, c_)), shape=(n_pt, n_pt))
+    n_lab, lab = connected_components(graph, directed=False)
+    return lab, n_lab
+
+
 def governing(res: dict) -> tuple[np.ndarray, np.ndarray]:
     """다스리는 안전율 · 방법(0 Bishop · 1 Janbu 보정) — 머리 주석 규칙."""
     jc = np.where(np.isfinite(res["janbu_corrected"]), res["janbu_corrected"], np.inf)