feat(536-1): E06 암반 평사투영 — Markland 판정 · Hoek-Bray 평면 · 쐐기 안전율 · 스테레오넷 · 절리 입력 · 화면

- 엔진(자체 · NumPy · SciPy): 평면 · 쐐기(교선) · 전도 Markland(측방 한계 20° 기본) · 평면 Hoek-Bray(인장균열 없음 · 임계 · 깊이 · 수압 · kh · 수동 볼트) · 쐐기 SWedge 꼴(미끄럼 방향 8 경우 · 두 면 · 한 면 · 수압 중력식 · kh · 볼트)
- 스테레오넷 좌표(등각 · 등면적 하반구 · 극점 · 대원 · 교선 · 마찰원 · 일광 엔벨로프 · 위험 영역 = 판정 부등식)
- 길: 절리 작업본 E06_StabRock/joints.json(초기값 · 저장 · 초기화) · 해석 = 암 측점 깎기 쪽만 · 파생 저장본 stab_rock · 결과 · 측점 스테레오넷
- 화면: 좌 절리군 표 · 설정 / 가운데 위험 띠 · 스테레오넷 / 우 측점 결과 · 안전율 표 · E 단계 막대에 E06(끝에 덧붙임)
- 시험 test_536_1(Hoek Sau Mau Ping 0.8184 · 1.0655 · Sharma 윗면 무게 · Hoek-Bray 쐐기 계수 · Kumsar · Priest 두 예제 · Markland 경계 · 투영 왕복 · 위험 영역 = 판정) · test_528 E 단계 수 4 → 5

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This commit is contained in:
2026-10-25 21:26:37 +09:00
co-authored by Claude Opus 5.5
parent 6029108bdb
commit f2e8bd6d5f
21 changed files with 3420 additions and 7 deletions
+1 -1
View File
@@ -257,7 +257,7 @@ def test_register_ignores_sent_company():
assert "company_name" not in data
@pytest.mark.parametrize("kind,count", [("ROAD", 7), ("PLAN", 3), ("STABILITY", 4)])
@pytest.mark.parametrize("kind,count", [("ROAD", 7), ("PLAN", 3), ("STABILITY", 5)])
def test_stage_seed_by_kind(kind, count):
log = []
asyncio.run(initialize_project_stages(FakeCursor({}, log), "p1", kind))
+426
View File
@@ -0,0 +1,426 @@
"""PLAN 536-1 — E06 암반 평사투영 · Markland · Hoek-Bray 평면 · 쐐기 (요구 345 · 351).
기준 예제(문헌 값 · Rocscience 검증서가 옮긴 값):
- 평면 = Hoek(Sau Mau Ping, 『Practical Rock Engineering』 7장 · RocPlane 검증 #1) 균열 없음 0.8184254 ·
임계 균열 + 수위 0.9z 1.0654738 / Sharma 외(1995) 윗면 0~30° 무게(RocPlane 검증 #2)
- 쐐기 = Hoek & Bray 쐐기 계수(대칭 쐐기 마찰만 · 해석해) · Kovári & Fritz · Kumsar 외(2000) 실험 ·
Priest(1993) 두 예제(SWedge 검증 #1~#7 값)
- Markland 경계 · 투영 왕복 · 위험 영역 = 판정 부등식
"""
import math
import numpy as np
import pytest
from E06_StabRock import E06_StabRock_Engine as eng
from E06_StabRock import E06_StabRock_Model as model
from E06_StabRock import E06_StabRock_Stereonet as net
TOL = 0.02 # 완료 기준 ±2 %
def close(value, expect, tol=TOL):
assert value == pytest.approx(expect, rel=tol), (value, expect)
# ---- 평면 (Hoek-Bray) -------------------------------------------------------------
SAU_MAU = dict(
H=60, slope_dip=50, plane_dip=35, c=100, phi=35, gamma=27, kh=0.08, gamma_w=10
) # MN → kN · 0.1 MPa = 100 kPa
def test_plane_hoek_no_crack():
r = eng.plane_fs(**SAU_MAU, water_fill=1.0)
close(r["F"], 0.8184254, 1e-4)
close(r["W"], 28627.8, 1e-4)
close(r["A"], 104.6068, 1e-5)
close(r["U"], 15691.0, 1e-4)
def test_plane_hoek_tension_crack():
z = eng.plane_critical_depth_hoek(60, 50, 35)
close(z, 14.0092, 1e-5)
r = eng.plane_fs(**SAU_MAU, water_fill=0.9, crack={"mode": "depth", "z": z})
close(r["F"], 1.0654738, 1e-4)
close(r["W"], 24843.9, 1e-4)
close(r["A"], 80.1826, 1e-5)
close(r["crack_b"], 15.3357, 1e-4) # 마루 뒤 거리(RocPlane 검증 #1)
def test_plane_critical_crack_search_matches_hoek():
"""건조 · 윗면 수평 — 수치 탐색 임계 자리 = Hoek-Bray 식 z · b."""
r = eng.plane_fs(
H=60, slope_dip=50, plane_dip=35, c=20, phi=35, gamma=27, crack={"mode": "critical"}
)
close(r["z"], eng.plane_critical_depth_hoek(60, 50, 35), 1e-4)
close(r["crack_b"], 15.3357, 1e-4)
@pytest.mark.parametrize(
"upper,weight,fs,fill",
[
(0, 2392.38, 1.58612, 1.00),
(10, 2446.29, 1.58148, 0.84),
(20, 2503.66, 1.58382, 0.71),
(30, 2568.90, 1.57849, 0.61),
],
)
def test_plane_upper_slope_sharma(upper, weight, fs, fill):
"""윗면 기운 사면 · 연직 균열 마루 뒤 15.336 m · 균열 수위(RocPlane 정수 %)."""
kw = dict(
H=60, slope_dip=50, plane_dip=35, c=12, phi=45, gamma=2.6, gamma_w=1.0, upper_dip=upper
)
r = eng.plane_fs(**kw, crack={"mode": "x", "b": 15.33576}, water_fill=fill)
close(r["W"], weight, 1e-4)
close(r["F"], fs)
def test_plane_bolt_vector():
"""수동 볼트 — 저항 += T cos(θ+ψp) + T sin(θ+ψp) tanφ (힘 벡터 풀이)."""
base = dict(H=10, slope_dip=60, plane_dip=35, c=10, phi=30, gamma=25)
r0 = eng.plane_fs(**base)
r1 = eng.plane_fs(**base, bolt_t=100, bolt_plunge=15)
a = math.radians(50)
expect = (
r0["resist"] + 100 * math.cos(a) + 100 * math.sin(a) * math.tan(math.radians(30))
) / r0["drive"]
close(r1["F"], expect, 1e-9)
def test_plane_not_formed():
assert eng.plane_fs(H=10, slope_dip=40, plane_dip=45, c=0, phi=30, gamma=25) is None
# 윗면이 활동면보다 가파르면 덩이 안 생김(균열 없음)
assert (
eng.plane_fs(H=10, slope_dip=60, plane_dip=30, c=0, phi=20, gamma=25, upper_dip=35) is None
)
# ---- 쐐기 (Hoek-Bray · SWedge 꼴) ------------------------------------------------
FLAT = {"dip_direction": 180, "dip": 0}
def _wedge(j1, j2, slope, phi, H=10, upper=FLAT, **kw):
return eng.wedge_fs(
j1={"dip_direction": j1[0], "dip": j1[1], "c": 0, "phi": phi},
j2={"dip_direction": j2[0], "dip": j2[1], "c": 0, "phi": phi},
slope={"dip_direction": slope[0], "dip": slope[1]},
upper=upper,
H=H,
gamma=26,
**kw,
)
@pytest.mark.parametrize(
"j1,j2", [((141, 45), (219, 45)), ((119, 50), (241, 50)), ((150, 60), (230, 55))]
)
def test_wedge_hoek_bray_wedge_factor(j1, j2):
"""마찰만 · 같은 φ — Hoek & Bray FS = K tanφ / tanψi · K = sinβ / sin(ξ/2)
(ξ = 두 면 사이 각(교선에 수직 단면) · β = 그 단면에서 쐐기 이등분선의 기울기)."""
phi = 35
r = _wedge(j1, j2, (180, 70), phi)
n1 = eng.plane_normal(*j1)
n2 = eng.plane_normal(*j2)
i = eng.line_vector(*r["line"])
# 교선에 수직인 면 안에서 두 면 안쪽 법선(위쪽) 사이 각 → 쐐기 사이 각 ξ = 180 − 그 각
xi = math.pi - math.acos(float(n1 @ n2))
bis = (n1 + n2) / np.linalg.norm(n1 + n2) # 쐐기 이등분선(교선에 수직 단면 안)
h = np.cross(i, [0.0, 0.0, 1.0])
h /= np.linalg.norm(h) # 그 단면의 수평
beta = math.acos(abs(float(bis @ h)))
K = math.sin(beta) / math.sin(xi / 2)
expect = K * math.tan(math.radians(phi)) / math.tan(math.radians(r["line"][1]))
close(r["F"], expect, 1e-3)
@pytest.mark.parametrize(
"j1,j2,slope,phi,plunge,fs",
[
((141, 45), (219, 45), (180, 70), 35, 37.85, 1.0), # Kovári & Fritz · SWedge #1
((119, 50), (241, 50), (180, 70), 35, 30.02, 1.6325), # SWedge #2
((195, 45), (105, 70), (160, 70), 30, None, 0.73), # Kumsar Case 3(Ankara) · SWedge 0.7123
((33.5, 75), (248, 75), (337.5, 75), 40.8, None, 2.02), # Kumsar Case 4(Dinar)
],
)
def test_wedge_kumsar_kovari(j1, j2, slope, phi, plunge, fs):
r = _wedge(j1, j2, slope, phi)
assert r["mode"] == "12"
close(r["F"], fs, 0.03 if fs == 0.73 else 0.01) # 실험값 0.73 = 판독(SWedge 0.7123)
if plunge:
close(r["line"][1], plunge, 1e-3)
def test_wedge_kumsar_seismic():
r = _wedge((33.5, 75), (248, 75), (337.5, 75), 40.8, kh=0.3303)
close(r["F"], 0.99) # Kumsar 0.99 · SWedge 0.9872
def test_wedge_priest_double_plane():
"""Priest(1993) 예제 1 — 두 면 미끄럼 · 수압 5 · 15 kPa · 철탑 180 kN(168/70)."""
r = eng.wedge_fs(
j1={"dip_direction": 203, "dip": 47, "c": 10, "phi": 40},
j2={"dip_direction": 287, "dip": 52, "c": 20, "phi": 35},
slope={"dip_direction": 230, "dip": 60},
upper={"dip_direction": 225, "dip": 5},
H=6.7978,
gamma=26,
water={"pressure": [5, 15]},
external=[{"force": 180, "trend": 168, "plunge": 70}],
)
assert r["mode"] == "12"
close(r["F"], 1.4966, 1e-3) # Priest ≈ 1.5
close(r["volume"], 45.201, 1e-3)
close(r["areas"][0], 41.147, 1e-3)
close(r["areas"][1], 20.428, 1e-3)
close(r["face_area"], 38.955, 1e-3)
close(r["N"][0], 407, 3e-3)
close(r["N"][1], 251, 3e-3)
def test_wedge_priest_single_plane():
"""Priest(1993) 예제 2 — 면 2 하나로 미끄럼(면 1 떨어짐) · 수압 25 · 15 kPa."""
r = eng.wedge_fs(
j1={"dip_direction": 65, "dip": 74, "c": 15, "phi": 32},
j2={"dip_direction": 186, "dip": 41, "c": 5, "phi": 40},
slope={"dip_direction": 134, "dip": 65},
upper={"dip_direction": 122, "dip": 11},
H=6.8471,
gamma=25,
water={"pressure": [25, 15]},
)
assert r["mode"] == "2"
close(r["F"], 0.8493, 1e-3) # Priest ≈ 0.86 · SWedge 0.8493
close(r["volume"], 81.741, 1e-3)
assert r["N"][0] == 0.0
def test_wedge_gravity_water_and_bolt():
"""중력 수압 U = p³ γw Hw a / 6(SWedge · Hoek-Bray) · 볼트는 안전율을 낮추지 않음."""
dry = _wedge((141, 45), (219, 45), (180, 70), 35)
wet = _wedge((141, 45), (219, 45), (180, 70), 35, water={"fill": 1.0})
expect = eng.GAMMA_W * wet["Hw"] / 6 * wet["areas"][0]
close(wet["water_force"][0], expect, 1e-9)
assert wet["F"] < dry["F"]
bolted = _wedge(
(141, 45), (219, 45), (180, 70), 35, bolt={"force": 50, "trend": 0, "plunge": 15}
)
assert bolted["F"] > dry["F"]
def test_wedge_not_formed():
# 교선이 사면보다 가파름 → 사면으로 안 나옴
assert _wedge((141, 75), (219, 75), (180, 50), 35) is None
# ---- Markland 경계 -------------------------------------------------------------
SLOPE = {"dip_direction": 90.0, "dip": 70.0}
def _j(dd, dip, phi=30.0, jid="J"):
return {"id": jid, "dip_direction": dd, "dip": dip, "phi": phi}
@pytest.mark.parametrize(
"dd,dip,hit",
[
(109.9, 45, True),
(110.1, 45, False),
(70.1, 45, True),
(90, 30.1, True),
(90, 29.9, False),
(90, 69.9, True),
(90, 70.1, False),
],
)
def test_markland_planar(dd, dip, hit):
assert (eng.kinematic(SLOPE, [_j(dd, dip)])["planar"] == ["J"]) is hit
def test_markland_planar_apparent_daylight():
"""측방 15° — 겉보기 사면 경사(tanψf cos15°)를 넘는 절리는 안 나옴."""
app = eng.apparent_dip(70, 90, 105)
assert eng.kinematic(SLOPE, [_j(105, app - 0.1)])["planar"] == ["J"]
assert eng.kinematic(SLOPE, [_j(105, app + 0.1)])["planar"] == []
@pytest.mark.parametrize(
"dd,dip,hit", [(270, 50.1, True), (270, 49.9, False), (289.9, 80, True), (290.1, 80, False)]
)
def test_markland_toppling(dd, dip, hit):
"""ψp > 90 − ψf + φ = 50 · 역방향 ±20."""
assert (eng.kinematic(SLOPE, [_j(dd, dip)])["toppling"] == ["J"]) is hit
def test_markland_wedge_plunge_bounds():
a, b = _j(50, 60, jid="A"), _j(130, 60, jid="B")
trend, plunge = eng.intersection(a, b)
assert trend == pytest.approx(90.0, abs=1e-6)
hit = eng.kinematic(SLOPE, [a, b])["wedge"]
assert hit and hit[0]["ids"] == ["A", "B"]
# φ 를 선주각 위로 올리면 사라짐
hi = [dict(a, phi=plunge + 0.1), dict(b, phi=plunge + 0.1)]
assert eng.kinematic(SLOPE, hi)["wedge"] == []
# 사면을 교선보다 완만하게 하면 안 나옴
assert eng.kinematic({"dip_direction": 90, "dip": plunge - 0.1}, [a, b])["wedge"] == []
# ---- 투영 ----------------------------------------------------------------------
@pytest.mark.parametrize("kind", net.PROJECTIONS)
def test_projection_round_trip(kind):
rng = np.random.default_rng(536)
for trend, plunge in zip(rng.uniform(0, 360, 200), rng.uniform(0, 89.9, 200)):
x, y = net.project(trend, plunge, kind)
t2, p2 = net.inverse(x, y, kind)
assert p2 == pytest.approx(plunge, abs=2e-3)
# 좌표는 소수 다섯째 자리 — 중심 가까이는 추향 오차가 1e-5 / r rad 까지
assert eng.ang_diff(t2, trend) <= math.degrees(1e-5 / math.hypot(x, y)) + 1e-6
def test_projection_known_radii():
assert net.radius(0, "equal_angle") == pytest.approx(1.0)
assert net.radius(0, "equal_area") == pytest.approx(1.0)
assert net.radius(90, "equal_angle") == pytest.approx(0.0)
assert net.radius(45, "equal_angle") == pytest.approx(math.tan(math.radians(22.5)))
assert net.radius(45, "equal_area") == pytest.approx(
math.sqrt(2) * math.sin(math.radians(22.5))
)
# 동쪽 수평선 = (1, 0) · 북쪽 = (0, 1)
assert net.project(90, 0, "equal_angle") == pytest.approx([1.0, 0.0], abs=1e-5)
assert net.project(0, 0, "equal_area") == pytest.approx([0.0, 1.0], abs=1e-5)
def _inside(pt, poly):
x, y = pt
hit = False
for (x1, y1), (x2, y2) in zip(poly, poly[1:] + poly[:1]):
if (y1 > y) != (y2 > y) and x < x1 + (y - y1) * (x2 - x1) / (y2 - y1):
hit = not hit
return hit
@pytest.mark.parametrize("kind", net.PROJECTIONS)
def test_zones_match_kinematic(kind):
"""위험 영역 다각형 안 극점 · 교선 = Markland 판정(경계 1° 띠 밖 표본)."""
rng = np.random.default_rng(7)
phi = 30.0
data = net.stereonet(SLOPE, [], phi, kind)
checked = 0
for dd, dip in zip(rng.uniform(0, 360, 1500), rng.uniform(1, 89, 1500)):
kin = eng.kinematic(SLOPE, [_j(dd, dip, phi)])
near = (
abs(eng.ang_diff(dd, 90) - 20) < 1
or abs(eng.ang_diff(dd, 270) - 20) < 1
or abs(dip - phi) < 1
or abs(dip - eng.apparent_dip(70, 90, dd)) < 1
or abs(dip - 50) < 1
)
if near:
continue
pole = net.project(*eng.pole(dd, dip), kind)
assert _inside(pole, data["zones"]["planar"]) == bool(kin["planar"]), (dd, dip)
assert _inside(pole, data["zones"]["toppling"]) == bool(kin["toppling"]), (dd, dip)
checked += 1
assert checked > 1000
# ---- 정본 · 측점 사면 ------------------------------------------------------------
def _doc():
doc = model.build_initial()
doc["sets"] = [
{
"id": "J1",
"name": "J1",
"from_m": 0,
"to_m": 100,
"dip_direction": 90,
"dip": 45,
"phi": 30,
"c": 5,
"persistence_m": None,
}
]
return doc
def test_check_and_clean():
assert model.check_joints(_doc()) == []
bad = _doc()
bad["sets"][0].update(dip=95, to_m=-1, name=" ")
bad["sets"].append(dict(bad["sets"][0]))
bad["settings"]["crack"] = {"mode": "depth", "z": None}
paths = {e["path"] for e in model.check_joints(bad)}
assert {"sets.0.dip", "sets.0.to_m", "sets.0.name", "sets.1.id", "settings.crack.z"} <= paths
clean = model.clean(dict(_doc(), extra=1))
assert set(clean) == {"version", "sets", "settings"}
assert clean["sets"][0]["dip"] == 45.0
def _cross(side_cut: str):
"""좌(+offset) 또는 우 깎기 — 노선 북쪽 진행(left = 서쪽 −x) · 원지반 깎기 쪽으로 0.6 오름 ·
노면 99 · 깎기 1:0.4 · 쌓기 1:1.5."""
sign = 1.0 if side_cut == "left" else -1.0
offs = np.arange(-20, 20.5, 0.5)
ground = [[float(o), 100 + 0.6 * sign * o] for o in offs]
design = []
for o, z in ground:
d = sign * o
if abs(o) <= 2:
z = 99.0
elif d > 2:
z = min(z, 99 + 2.5 * (d - 2))
else:
z = max(z, 99 - (abs(d) - 2) / 1.5)
design.append([o, z])
return (
{
"frame": {"left_xy": [-1.0, 0.0]},
"samples": [{"offset_m": o, "elevation_m": z} for o, z in ground],
},
{
"section_mode": f"{side_cut}_cut",
"cut_slope_ratio": 0.4,
"road_edges": {"left": {"offset_m": 2.0}, "right": {"offset_m": -2.0}},
"design_line": [{"offset_m": o, "elevation_m": z} for o, z in design],
},
)
@pytest.mark.parametrize("side,dd", [("left", 90.0), ("right", 270.0)])
def test_station_slope_orientation(side, dd):
"""좌 = 서쪽 산 · 사면이 동(노면)으로 기욺 = 경사방향 90 · 우는 270 · 경사 = atan(1/0.4)."""
cross, design = _cross(side)
slopes = model.station_slopes(cross, design)
assert len(slopes) == 1 and slopes[0]["side"] == side
assert slopes[0]["dip_direction"] == pytest.approx(dd)
assert slopes[0]["dip"] == pytest.approx(math.degrees(math.atan(2.5)), abs=0.01)
assert slopes[0]["height_m"] > 1
assert slopes[0]["upper_dip"] == pytest.approx(math.degrees(math.atan(0.6)), abs=0.5)
def test_analyze_station(monkeypatch):
cross, design = _cross("left")
monkeypatch.setattr(
model,
"rock_stations",
lambda _s: [{"chainage_m": 20.0, "label": "STA", "ground_type": "ripping_rock"}],
)
monkeypatch.setattr(model.src, "cross_section", lambda _s, _c: {**cross, "design": design})
ground = {
"materials": [],
"stations": [],
"safety_factors": {"dry": {"cut": 1.5}, "wet": {"cut": 1.2}},
"seismic": {"enabled": False},
}
out = model.analyze(None, _doc(), ground)
assert len(out) == 1
e = out[0]
assert e["kinematic"]["planar"] == ["J1"]
assert {c["condition"] for c in e["cases"]} == {"dry", "wet"}
dry = next(c for c in e["cases"] if c["condition"] == "dry")
assert dry["required"] == 1.5 and dry["ok"] == (dry["F"] >= 1.5)
assert e["verdict"] in ("OK", "NG")
assert e["gamma_r"] == model.DEFAULT_GAMMA_R
nets = model.station_net(e, _doc()["settings"])
assert set(nets) == {"equal_angle", "equal_area", "phi"}
assert nets["equal_area"]["zones"]["planar"]