Research Paper

Journal of Korean Tunnelling and Underground Space Association. 30 September 2026. 555-570
https://doi.org/10.9711/KTAJ.2026.28.5.555

ABSTRACT


MAIN

  • 1. 서 론

  • 2. 물리적 시나리오와 해석모델

  •   2.1 관 누수에 의한 공동 발달과 수치모델

  •   2.2 수리이력 전달과 포화도에 따른 지반 물성 저하

  •   2.3 침투·침식 기준과 공동 확대의 판정

  •   2.4 비교조건과 수치해석의 적정성 평가

  • 3. 해석결과

  •   3.1 누수 유무에 따른 공동 주변 응답

  •   3.2 침투·침식 이후 파괴 과정의 전이

  •   3.3 침투·침식 임계유속과 누수 경과시간별 수리상태

  •   3.4 공동 주변 파괴영역의 공간적 발달

  •   3.5 메쉬 크기에 따른 결과 변화

  •   3.6 선행연구와의 정성적 비교

  •   3.7 물성 저하 가정의 민감도와 적용범위

  • 4. 결 론

1. 서 론

도심지 매설관의 결함부에서 누수가 지속되면 침투력에 의해 토립자가 이탈·이송되어, 관 주변 지반에 침식공동이 형성될 수 있다(Sato and Kuwano, 2015). 공동 확대로 상부 지반의 하중전달 경로가 약화되면 지표침하 또는 도로함몰로 이어질 수 있다. 모형실험에 따르면 공동의 형상과 지반 변형범위는 토사의 입도분포와 상대밀도에 크게 좌우된다. 비슷한 크기의 입자로 이루어져 세립분이 거의 없는 토사는 입자 크기의 분포가 넓고 세립분을 포함한 토사보다 공동 형성에 취약하였으며, 느슨한 지반에서는 변형이 지표까지 확산된 반면 조밀한 지반에서는 공동 주변에 집중되었다(Kwak et al., 2019). Indiketiya et al. (2019)의 반복 누수실험에서도 관 주변 되메움재의 입도분포와 결함 크기가 침식저항과 공동 발달을 지배하였다. 이러한 결과는 누수 유발 공동이 이미 형성된 빈 공간의 단순 붕괴가 아니라 토립자 이탈·이송, 지반 물성 변화 및 응력재분배가 연속적으로 작용하여 공동이 확대되는 현상임을 나타낸다.

기존 수치연구는 내부침식을 다공질 연속체의 질량교환으로 해석하는 방법과 개별 입자의 이동을 추적하는 방법으로 구분된다. 다공질 연속체 모델은 유체 흐름에 따른 세립분의 이탈·이송과 이에 수반되는 수리·역학적 물성 변화를 연계할 수 있으나, 공동 규모가 커지거나 불연속면이 형성된 이후에는 연속체 가정의 적용범위가 제한된다(Zhang et al., 2013). 토체를 지지하는 입자 구조, 침식 가능한 세립분, 이동 중인 세립분 및 간극수 사이의 질량교환을 고려한 수리–역학 연계모델도 제안되었다(Yang et al., 2019). Zhang et al. (2020)은 손상관 주변의 침식률을 사용자 요소에 도입하여 수리수두가 결함형상보다 공동 형성과 토사유실에 더 큰 영향을 미친다는 실험에서 관찰된 경향을 재현하였다. 다만 이 연구에서는 세립분 유실에 따른 지반의 강성 및 강도 변화를 고려하지 않았다.

입자기반 해석은 토립자의 이동경로와 공동 상부의 하중이 주변 지반으로 전달되는 아칭 현상을 상세히 추적하는 데 유리하다. CFD–DEM 해석에서는 지하수위 상승이 일정 범위까지 침식을 가속하고 관 결함의 방향이 최종 침식범위를 변화시키는 것으로 나타났다(Ibrahim and Meguid, 2023). 반복 누수로 형성되는 공동은 급속 확대, 완만한 성장 및 안정화 단계를 거치며, 성장이 멈춘 공동도 강우 침투로 불포화토의 모관흡수력이 감소하면 붕괴할 수 있음이 모형실험과 FDM–DEM 해석에서 보고되었다(Wang et al., 2026). 그러나 입자기반 방법은 다수의 입자와 유체–입자 상호작용을 계산해야 하므로 해석영역이 넓고 변형이 큰 문제에는 상당한 계산비용이 요구된다. Eulerian 해석은 입자 크기에 따른 선택적 유실을 직접 재현하지 못하지만, 고정된 Eulerian 격자에서 재료가 이동하므로 심한 요소왜곡 없이 공동이 크게 확대되는 과정을 해석할 수 있다.

본 연구에서는 작은 초기 공동에서 침투·침식 기준에 따라 토체 이탈이 시작되고, 공동 주변의 응력분포 변화와 인장응력 손상 및 전단변형 집중으로 공동이 확대되는 과정을 Eulerian 대변형 해석으로 구현하였다. 최대 이력 유효포화도에 따라 지반의 강성, 강도 및 파괴에너지를 감소시키고, 각 요소를 처음 도달한 이탈기준에 따라 침투·침식형, 인장응력 손상형 및 전단변형 집중형으로 분류하였다. 공동 확대의 시간적·공간적 전개는 에너지비, 유량수지 및 메쉬 크기에 따른 결과 변화와 함께 평가하였다. 연구의 목적은 입자별 침식량이나 실제 공동크기를 예측하는 것이 아니라, 전달된 수리이력에 따른 침투·침식 개시와 후속 역학적 확대를 구분하고 수리조건별 상대경향을 평가하는 데 있다.

2. 물리적 시나리오와 해석모델

2.1 관 누수에 의한 공동 발달과 수치모델

누수에 따른 지중공동의 대변형을 해석하기 위해 CEL 해석체계의 Eulerian 정식화를 적용하였다. 본 모델에는 관로 등의 Lagrangian 구조부재를 포함하지 않았으며, 토체와 물의 대변형 및 이동은 고정된 Eulerian 격자에서 계산하였다. Fig. 1(a)는 관 결함부의 국부 누수에 따른 유효포화도와 Darcy 유속의 증가, 포화도 의존 지반 물성 저하 및 임계유속에 따른 침투·침식이 토체 이탈로 이어지는 과정을 개념적으로 나타낸다. 초기 이탈로 자유면이 형성되면 공동 주변의 응력이 재분배되고, 천단과 측벽에서 인장응력에 의한 손상과 전단변형 집중이 발생하면서 공동이 상부와 측방으로 확대될 수 있다. Fig. 1(a)의 관로는 현장 누수원을 설명하기 위한 개념도이며, 수치모델에서는 관 자체를 구조부재로 포함하지 않고 누수 유입부와 하류 배출부로 단순화하였다.

Fig. 1(b)의 계산영역은 길이 20.0 m, 높이 9.5 m, 두께 0.25 m의 얇은 3차원 Eulerian 모델이다. 전 영역의 두께방향 속도를 구속하여 평면변형률에 가까운 거동을 모사하였으며, 토체와 물은 동일한 Eulerian 격자에서 재료 체적분율로 구분하였다. 누수 위치에는 반경 0.15 m의 원형 초기 공동만을 설정하여 후속 이탈·손상영역이 해석 과정에서 발달하도록 하였다. 좌우 경계와 저류·배출영역의 측벽은 수평방향 속도를, 상·하단 경계는 연직방향 속도를 구속하였다. 하류 배출영역에서는 지정 출구에만 하향 유속을 부여하고 나머지 경계의 유출을 제한하였다. 상단 경계의 연직변위를 구속하였으므로 지표침하와 함몰은 해석대상에서 제외하였으며, 평가대상은 초기 공동 주변의 토체 이탈과 파괴영역 확대이다.

Fig. 1(c)는 기준 메쉬와 세분화된 메쉬의 구성을 나타낸다. 공동 주변의 요소크기는 각각 0.125 m와 0.0625 m이며, 두 메쉬의 결과를 비교하여 최초 이탈시점, 파괴영역의 발달방향 및 정량지표에 대한 메쉬 크기의 영향을 평가하였다.

물의 압축성은 충격파 속도와 입자속도의 선형 관계를 사용하는 Us–Up 상태방정식으로 나타냈다. 토체에는 Drucker–Prager 항복기준과 비연관 유동법칙(non-associated flow rule)을 적용하였으며, 공동 주변의 인장균열과 전단변형에 따른 강도 저하는 손상모델로 표현하였다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F1.jpg
Fig. 1.

Leakage-induced cavity-development mechanism, numerical idealization, and meshes

2.2 수리이력 전달과 포화도에 따른 지반 물성 저하

독립된 유한체적 침투해석에서 누수 개시 후 0, 1, 3, 7, 14 및 30-day의 유효포화도 Se와 Darcy 유속 q를 산정하였다. 각 시점의 공간분포를 역학모델에 전달하고, 그 사이의 수리상태는 시간에 따라 선형 보간하였다. 이하에서는 이 시간별 수리해석 결과를 수리이력이라 하고, 수리해석의 누수 경과시간을 th로 표기한다. 30-day 시점의 수리상태를 적용한 뒤에는 해당 유효포화도와 Darcy 유속 분포를 고정하여 후속 역학응답을 계산하였다. 따라서 th는 역학해석의 계산시간이나 현장 공동 발생시간이 아니라, 역학모델에 전달된 수리이력이 대표하는 누수 경과시간이다.

역학해석에서 발생한 공동 확대, 간극률 및 투수경로 변화는 수리해석에 되돌려 반영하지 않았으므로, 본 해석은 수리응답이 역학응답에만 영향을 미치는 단방향 수리–역학 연계(one-way coupling)에 해당한다. 수리이력 전달, 포화도에 따른 물성 저하, 요소 이탈기준, 고정격자에서의 응답지표 산정 및 제외된 상호작용을 Fig. 2에 정리하였다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F2.jpg
Fig. 2.

One-way transfer of time-dependent hydraulic fields and calculation of fixed-grid response measures

침투에 따른 모관흡수력 감소가 불포화토의 강도 및 강성을 낮출 수 있다는 실험적 관찰(Abeyrathne et al., 2023)을 참고하여, 본 해석에서는 포화도 이력에 따른 지반 물성 저하를 모델 입력관계로 설정하였다. 다만 해당 실험은 본 해석의 토질이나 저하계수를 검증한 자료가 아니다. Se는 유효포화도이며, Se,max는 해당 재료점이 해석 중 경험한 최대 유효포화도이다. 포화도가 감소한 뒤 물성이 즉시 회복되는 것을 방지하기 위해 Se,max를 이력변수로 사용하고, 무차원 물성 잔존계수 η를 식 (1)과 같이 정의하였다.

(1)
η=1Se,max≤Se,0η=exp-kdSe,maxndSe,max>Se,0

Se,0 = 0.3은 초기 배경 유효포화도 0.2보다 큰 값으로 정한 모델의 물성 저하 활성화 기준이며, 실험으로 확인한 임계포화도가 아니다. 실제 해석 구현에서는 이 값을 넘으면 최대 이력 유효포화도 자체를 지수함수에 대입한다. 계수 kd = 1.0과 지수 nd = 1.5는 활성화된 저하 관계의 크기와 비선형성을 조절하는 입력값이다. 최대 유효포화도를 이력값으로 저장하므로 건조과정에서도 저하된 물성이 회복되지 않는다. 아래첨자 0과 h는 각각 초기값과 저하 후 값을 나타내며, 탄성계수 E, 점착력 c 및 내부마찰각 φ의 잔존비는 식 (2), (3), (4)로 계산하였다.

(2)
EhE0=maxrE,1-aE(1-η)
(3)
chc0=maxrc,1-ac(1-η)
(4)
tanφhtanφ0=maxrφ,1-aφ(1-η)

식 (2), (3), (4)에서 aE, ac 및 aφ는 유효포화도 증가에 대한 각 물성의 저하 민감도이며, rE, rc 및 rφ는 최소 잔존비이다. 본 해석에는 aE = 0.6875, ac = 0.9, aφ = 0.2, rE = 0.1, rc = 0.15 및 rφ = 0.7을 적용하였다. 초기 물성은 탄성계수 E0 = 50 MPa, 포아송비 ν = 0.3, 점착력 c0 = 25 kPa, 내부마찰각 φ0 = 30°, 팽창각 ψ = 5° 및 인장강도 ft0 = 2 kPa이다. 유효포화도 증가에 따라 인장강도와 인장 파괴에너지는 각각 ft = ft0η와 Gft = Gft0ηmG로 감소시켰으며, Gft0 = 15 J/m2와 mG = 1.0을 사용하였다. 초기 전단 파괴에너지 Gfs0는 125 J/m2이다. 이들 값은 현장 보정값이 아니라 파괴 과정과 상대경향을 분석하기 위한 모델 입력값이다. rE = 0.1과 rc = 0.15는 각각 잔존비의 하한으로, 포화도 저하만으로 탄성계수와 점착력이 반드시 초기값의 10%와 15%까지 감소한다는 뜻은 아니다. 본 연구의 포화도 범위(0–1)와 입력계수에서 포화도 저하 성분만 고려한 최소 잔존비는 각각 약 0.565와 0.431이다. 이 잔존비는 제시한 재료와 물성 저하식의 입력조건에서만 성립한다.

물성 저하식의 형태에 따른 결과 변화를 확인하기 위해 기존의 불연속 활성화식과 이를 연속화한 대안식을 비교하였다. 대안식에서는 활성화 기준 이하일 때 물성 잔존계수를 1로 두고, 그 기준을 넘으면 최대 이력 유효포화도의 증가분을 활성화 기준에서 완전포화까지의 범위로 나누어 식 (1)의 지수항에 대입하였다. 활성화 기준은 0.2와 0.4로, 식 (1)의 저하계수는 0.75와 1.25로 각각 한 번에 하나씩 변경하였다. 활성화 기준을 바꿀 때에는 물성 저하뿐 아니라 요소 이탈기준을 평가하는 포화도 조건에도 동일한 값을 적용하였다. 모든 사례에서 기준 메쉬, 경계조건, 수리이력 및 해석시간은 같게 유지하였다. 이 값들은 시험으로 확인한 재료변수의 범위가 아니라 모델 입력값의 변화에 대한 민감도 검토범위이다.

2.3 침투·침식 기준과 공동 확대의 판정

이탈기준은 초기 공동과 연결되고 포화도가 증가한 요소에서만 평가하였다. Darcy 유속 q가 임계유속 qcrit에 도달하면 침투·침식 기준에 도달하며, 이후 인장응력에 의한 손상 또는 전단변형 집중 기준에 먼저 도달한 요소는 각각 인장응력 손상형과 전단변형 집중형으로 분류하였다. 기준 임계유속은 qcrit = 4.0×10-7 m/s (400 nm/s)이다. 이는 실내 침식시험으로 얻은 값이 아니라 200, 400, 600 nm/s의 비교조건에서 가운데에 놓인 해석 입력값이다. 이 분류는 손상 정도의 단계를 의미하지 않으며, 각 요소가 처음 도달한 이탈기준을 나타낸다.

일반적인 내부침식은 세립분의 이탈·이송과 이에 따른 간극률 및 투수특성 변화를 포함한다(Zhang et al., 2013; Yang et al., 2019). 본 모델은 토립자의 실제 이동을 직접 해석하지 않으므로, 본 논문에서는 침투·침식을 초기 공동과 연결된 요소에서 Darcy 유속이 설정 임계값에 도달할 때 하중전달 능력을 상실하도록 하는 해석상의 이탈기준으로 정의하였다. 파괴영역의 확대는 각 요소가 처음 도달한 이탈기준을 기록하고 해당 요소의 평면 투영면적을 누적하여 정량화하였다. 이 면적은 하중전달 능력을 상실한 영역의 상대적 크기를 나타내며, 실제 개방 공동면적이나 유실토량과 동일하지 않다.

2.4 비교조건과 수치해석의 적정성 평가

주요 비교조건은 수리상태 적용 유무, 침투·침식 임계유속, 누수 경과시간 및 메쉬 크기로 구성하였다(Table 1). 물성 저하식의 전환형태와 입력값은 별도의 민감도 해석으로 검토하였다. 각 비교에서는 요소 이탈 여부, 파괴영역의 확대 및 최대 상향 확장거리를 평가하였다.

수치결과의 적정성은 에너지비(KE/IE ≤ 5%, AE/IE ≤ 15%), 유량수지 오차(5% 이하) 및 비지정 경계 유출량(공급수량의 1% 이하)으로 확인하였다. 수리상태를 적용하지 않은 비교해석은 최종 IE가 매우 작으므로 해석의 정상 종료, 수리해석 상태변수의 범위, 요소 이탈 여부와 질량·물 체적 변화로 평가하였다.

Table 1.

Principal comparison cases and their physical purpose

Comparison case Variable Purpose
Hydraulic field 
excluded
Initial effective saturation Se = 0.2; 
Darcy velocity = 0 m/s
Determine whether the detachment criterion is reached solely 
because of the initial cavity and self-weight
Critical Darcy 
velocity
qcrit = 200 (numerical diagnostic case), 400, and 600 nm/s Evaluate the response to changes in the seepage-erosion 
threshold; exclude the 200 nm/s case from quantitative 
comparison because it exceeds the energy-acceptance criterion
Hydraulic fields at 
different elapsed 
leakage times
qcrit = 400 nm/s at 
th = 8.6, 17.1, and 30 days
Assess changes in the extent of tensile damage and 
shear localization with elapsed leakage time
Mesh-size effect h = 0.125 and 0.0625 m Evaluate mesh-size effects on the onset time, growth direction, 
and quantitative response metrics

3. 해석결과

3.1 누수 유무에 따른 공동 주변 응답

기준 누수해석과 수리상태를 적용하지 않은 비교해석은 수리상태 적용 여부만 달리하고 형상, 재료, 초기 공동 및 경계조건은 동일하게 유지하였다. 이를 통해 초기 공동과 자중만으로 발생하는 변형을 전달된 수리상태에 따른 요소 이탈과 구분하였다.

기준 누수해석에서는 누수 경과시간 th ≈ 1.0 day의 수리상태에서 첫 요소가 이탈기준에 도달하였다. 30-day 수리상태를 고정한 최종 역학응답에서 파괴영역 면적은 1.875 m2, 초기 공동 천단 기준 최대 상향 확장거리는 0.638 m였다(Fig. 3). 기준 누수해석의 KE/IE, AE/IE 및 유량수지 오차는 각각 0.0517%, 10.18% 및 0.0983%로 허용기준을 만족하였다. 비교해석에서도 요소 이탈이 없고 질량과 물 체적 변화가 각각 -6.15×10-7과 -0.0127%로 작아 수치적으로 안정된 결과를 확인하였다. 따라서 기준 누수해석에서 확인된 요소 이탈은 초기 공동과 자중만으로 발생한 현상이 아니라 전달된 수리상태에 의해 유발된 응답으로 해석하였다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F3.jpg
Fig. 3.

Response metrics with and without transferred hydraulic fields

3.2 침투·침식 이후 파괴 과정의 전이

Fig. 4(a)는 요소별로 처음 도달한 이탈기준에 따른 누적면적을 나타내며, 공동 확대의 주된 과정이 침투·침식에서 역학적 손상으로 전환되는 양상을 보여준다. 침투·침식 기준은 th ≈ 1.0 day에 처음 도달하였고, th ≈ 1.8 day에는 인장응력 손상 기준에 도달한 요소가 발생하였다. 침투·침식 기준 요소의 누적면적 증가는 th ≈ 5.0 day 이후 제한적이었으며, th ≈ 5.7 day부터 인장응력 손상 기준 요소의 누적면적이 침투·침식 기준 요소의 누적면적을 초과하였다. 전단변형 집중 기준은 th ≈ 23.2 day에 처음 도달하였다.

Fig. 4(b)의 초기 공동 천단 기준 최대 상향 확장거리는 연속적으로 증가하지 않고 th ≈ 5.6, 14.6 및 29.1 day 부근에서 계단형으로 증가하였다. 이는 초기 침투·침식으로 자유면이 형성된 뒤 공동 주변의 응력이 재분배되고, 포화도 증가로 물성이 저하된 천단부에 인장응력 손상이 누적되고, 누수 경과시간 약 23.2 day부터 전단변형이 좁은 영역에 집중된 결과로 해석된다.

침투·침식으로 자유면이 확대되면 공동 주변의 하중전달 경로가 변하고, 포화도 증가로 강성과 강도가 감소한 천단·측벽부에 인장응력과 전단변형이 집중될 수 있다. 누수 경과시간 약 23.2 day에 전단변형 집중 기준에 도달한 요소가 나타난 것은 공동 주변에서 전단변형이 국부적으로 집중되는 영역이 발달한 것으로 해석된다. 이러한 공동 확대 순서는 누수에 의한 토립자 유실과 공동 주변에 변형이 집중되는 과정이 연속적으로 나타난다는 모형실험의 관찰(Indiketiya et al., 2019; Kwak et al., 2019)과 정성적으로 부합한다. 다만 선행실험과 본 해석은 토질, 누수조건 및 경계조건이 다르므로, 비교는 현상 전개순서가 정성적으로 일치하는지를 확인하는 수준으로 한정하였다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F4.jpg
Fig. 4.

Evolution of cavity-growth indicators with elapsed leakage time

3.3 침투·침식 임계유속과 누수 경과시간별 수리상태

에너지 허용기준을 만족한 조건에서 qcrit가 400에서 600 nm/s로 증가하면 최초 이탈기준 도달시점은 th ≈ 1.0 day에서 1.3 day로 지연되었다. th = 8.6 day에서 파괴영역 면적은 0.625 m2에서 0.219 m2로 감소하였고, 초기 공동 천단 상부에서 추가로 이탈기준에 도달한 요소는 없었다(Fig. 5(a)). qcrit = 200 nm/s 조건에서는 파괴영역 면적 1.000 m2와 상향 확장거리 0.638 m가 계산되었으나 최종 인공에너지/내부에너지비가 53.59%로 허용기준을 초과하였다. 이 결과는 정량 비교에서 제외하고 참고값으로만 제시하였다. 기준값 400 nm/s는 현장이나 실내시험에서 보정한 침식 개시 유속이 아니라 설정값에 따른 모델 응답을 비교하기 위한 기준조건이다. 200 nm/s 조건은 에너지 허용기준을 벗어나므로, 400 nm/s를 검증된 침식 임계값으로 해석해서는 안 된다.

qcrit = 400 nm/s에서 th = 8.6, 17.1 및 30-day의 수리상태를 비교하면 파괴영역 면적은 각각 0.625, 1.000 및 1.875 m2로 증가하였다(Fig. 5(b)). 최대 상향 확장거리도 0.388, 0.513 및 0.638 m로 증가하였다. 면적과 길이의 단위가 다른 점을 고려하여 Fig. 5에서는 파괴유형별 면적과 상향 확장거리를 상·하단 그래프로 구분하였다. 누수 경과시간이 증가해도 침투·침식 기준 요소의 면적 증가는 제한적이었으나, 인장응력 손상 기준 요소의 면적은 증가하였다. 결과적으로 본 모델에서는 임계유속을 높일수록 침투·침식 개시가 지연된 반면, 누수 경과시간이 길어질수록 자유면 형성 이후의 인장응력 손상과 전단변형 집중 범위가 확대되었다. 30-day 조건의 값은 해당 수리상태를 고정한 뒤 후속 역학응답이 종료될 때까지 계산한 결과이다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F5.jpg
Fig. 5.

Effects of critical Darcy velocity and elapsed leakage time

침투·침식 임계유속에 따라 이탈 개시시점과 파괴영역이 달라진 결과는 수리조건이 토사유실과 공동 발달에 중요하다는 Zhang et al. (2020)의 관찰과 정성적으로 부합한다. 다만 선행연구는 구동력인 수리수두를 변화시킨 반면, 본 연구는 침투·침식 기준의 임계조건인 qcrit를 변화시켰다. 두 비교변수의 물리적 의미가 다르므로 본 결과를 선행연구의 정량 검증이나 직접적인 경향 일치로 해석할 수 없다. 또한 qcrit는 토질시험으로 산정한 재료 고유의 임계유속이 아니므로, 400 nm/s는 설계값이 아니라 동일한 수치모델에서 민감도를 비교하기 위한 기준조건이다.

3.4 공동 주변 파괴영역의 공간적 발달

최초 이탈기준의 공간분포에서는 침투·침식 기준 요소가 초기 공동 부근에서, 인장응력 손상 기준 요소가 측벽부에서 먼저 나타났다. 전단변형 집중 기준 요소는 누수 경과시간 약 23.2 day에 천단·측벽에서 확인되었다(Fig. 6). 각 결과에는 대응하는 누수 경과시간 th ≈ 1.0, 1.8, 5.0, 5.7 및 23.2 day를 표시하였으며, 최종 상태는 30-day 수리상태를 고정한 후의 역학응답이다. 모든 시점은 동일한 메쉬와 범례를 적용하여 비교하였다. 회색 영역은 이탈기준에 도달하지 않은 토체이며, 흰색 영역은 토체가 존재하지 않는 초기·확대 공동 또는 배출공간이다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F6.jpg
Fig. 6.

Spatial distribution of the detachment criteria first met with hydraulic elapsed time

범례의 정수 분류값 0, 1, 2 및 3은 각각 이탈기준 미도달, 인장응력 손상 기준, 전단변형 집중 기준 및 침투·침식 기준을 의미한다. 중간색은 인접 절점의 값을 평균하여 표시하는 75% 절점평균화 과정에서 생긴 색상 보간이며 별도의 파괴유형이나 손상수준을 의미하지 않는다. 파괴기준별 요소 수와 파괴영역 면적은 등고선 면적이 아니라 요소별로 처음 도달한 이탈기준의 원자료로 산정하였다(Table 2).

Table 2.

Failure indicators at the reference critical Darcy velocity

Elapsed leakage 
time, th
(day)
Elements meeting 
seepage–erosion 
criterion
Elements meeting 
tensile-damage 
criterion
Elements meeting 
shear-localization 
criterion
Failed area
(m2)
Maximum 
upward growth
(m)
8.6 12 28 0 0.625 0.388
17.1 14 50 0 1.000 0.513
30 (final) 14 98 8 1.875 0.638

The failed area was calculated by multiplying the total number of elements meeting the three criteria by the area of one element in the reference mesh, (0.125 m)2. The 30-day row reports the final mechanical response after fixing the hydraulic state.

Fig. 6(f)의 상향 확대 결과는 공동이 커진 뒤의 지하수 흐름을 재산정하지 않은 단방향 연계와 상단 연직구속을 함께 고려하여 해석해야 한다. 상단 경계의 연직변위가 허용되지 않으므로 지표침하나 함몰의 발생 여부를 이 결과에서 판정할 수 없다. 따라서 최대 상향 확장거리는 본 모델 내부 파괴영역의 상대적 지표로 한정하였다.

3.5 메쉬 크기에 따른 결과 변화

공동 주변 요소크기를 0.125 m에서 0.0625 m로 줄인 기준 메쉬와 세분화된 메쉬에서 최초 토체 이탈은 누수 경과시간 th ≈ 1.0 day에 발생하였고 최대 수평 확장폭은 3.0 m로 같았다. 30-day 수리상태를 고정한 최종 응답에서 파괴영역 면적과 최대 상향 확장거리는 세분화된 메쉬에서 각각 11.67%와 14.71% 감소하였다(Fig. 7). 최초 이탈시점과 주된 발달방향은 유지되었으나 정량지표의 변화가 10%를 초과하였다. 따라서 두 메쉬의 비교만으로 메쉬 수렴성이 입증되었다고 볼 수 없으며, 본 결과는 메쉬 크기의 영향을 확인한 예비 검토에 해당한다. 세분화된 메쉬의 최대 운동에너지/최종 내부에너지, 최종 인공에너지/내부에너지 및 유량수지 오차는 각각 0.0248%, 8.79% 및 0.1024%로 허용기준을 만족하였다.

https://cdn.apub.kr/journalsite/sites/kta/2026-028-05/N0550280512/images/kta_2026_285_555_F7.jpg
Fig. 7.

Effect of mesh refinement on the final failure region for the 30-day hydraulic state

정량적 메쉬 수렴성을 판단하려면 동일한 세분비를 유지한 h = 0.03125 m 메쉬를 추가하고, 파괴영역 면적, 상향 확장거리, 최대 수평 확장폭 및 개시시점의 변화율을 평가해야 한다. 따라서 현재 결과는 공동 확대 과정과 조건별 상대경향을 해석하는 데 활용할 수 있으나, 현장 공동크기나 함몰시점을 직접 예측하는 근거로는 충분하지 않다.

3.6 선행연구와의 정성적 비교

본 연구와 동일한 토질, 누수 및 경계조건의 실험자료가 없으므로 선행연구와의 비교는 정량 검증이 아니라 현상 전개의 정성적 비교로 한정하였다(Table 3). 비교항목은 외부 수리작용에 따른 공동 개시, 침투·침식 이후 역학적 손상이 이어지는 순서 및 파괴영역의 발달방향이다. 모형척도, 토질, 유량 정의 및 관 결함조건이 서로 다르므로 공동 발생시간과 최종 크기는 직접 비교하지 않았다.

Guo et al. (2024)의 대형 모형실험에서는 침투력이 증가하면서 토립자 배열이 흐트러지고 공동과 붕괴영역이 단계적으로 발달하였으며, 공동 발달이 진행된 뒤에는 관 결함부에서 상부로 넓어지는 원추형 붕괴영역이 관찰되었다. Tan et al. (2022)은 특정 원형 천공조건에서 폭이 높이보다 큰 타원형 공동을 보고했으나, 관 파손형태와 조건에 따라 공동형상이 달라질 수 있음을 함께 제시하였다. 따라서 선행연구에서 보고된 공동 형상은 보편적인 단일 형상이 아니라 각 실험의 결함 형상과 수리·지반 조건에 따른 결과로 해석해야 한다.

본 해석의 등고선 영역은 실제 공동 경계가 아니라 최초 이탈기준에 도달한 고정 Eulerian 요소의 집합이다. 최종 파괴영역은 초기 공동에서 측방과 상부로 발달하였으며, 이는 선행연구에서 관찰된 결함부에서 시작되는 확대 방향과 정성적으로 부합한다. 다만 공동 확대에 따른 지하수 흐름 변화와 토립자 이송을 직접 해석하지 않았으므로, 문헌의 공동 형상과 본 해석의 파괴영역을 정량적으로 대응시킬 수 없다.

Table 3.

Qualitative phenomenon comparison with published leakage-induced cavity studies

Comparison item Numerical response in this study Observations from 
previous studies
Comparison
Effect of transferred 
hydraulic fields
No element reached a detachment 
criterion when the hydraulic field 
was excluded, whereas detachment 
initiated at th ≈ 1 day under 
leakage conditions
Groundwater flow and pipe 
defects provide the driving and 
discharge conditions for 
soil-particle migration 
(Guo et al., 2024; 2025)
Qualitatively consistent with 
cavity initiation induced by 
the transferred hydraulic field
Sequence of 
cavity growth
The seepage-erosion criterion was 
followed by the expansion of 
tensile damage and subsequent 
shear localization
Erosion and the collapse zone 
develop progressively after 
hydraulic disturbance of the soil 
skeleton (Guo et al., 2024; 2025)
Qualitatively consistent with 
the reported sequence of 
phenomenon development
Growth direction and 
geometry
The failure region represented by 
fixed Eulerian elements expanded 
laterally and upward from the 
initial cavity
Laterally dominant elliptical 
cavities and subsequent conical 
collapse zones under specific 
circular-orifice conditions (Tan et al., 2022; Guo et al., 2024)
Only the growth direction is 
compared qualitatively 
because the cavity boundary 
and failure region represent 
different physical quantities
Direct comparison of 
time and size
Model response metrics for the 
30-day hydraulic field
Model scale, soil properties, 
discharge definition, and 
pipe-defect conditions differ 
among the studies
Occurrence time and final 
size are excluded from 
quantitative comparison

3.7 물성 저하 가정의 민감도와 적용범위

불연속 활성화식을 연속식으로 바꾸면 기준 메쉬의 최종 이탈요소 수는 120개에서 118개로 줄었고, 파괴영역 면적은 1.875 m2에서 1.844 m2로 약 1.7% 감소하였다. 첫 침투·침식 기준 도달은 두 식 모두 누수 경과시간 약 1.0 day였고 최대 상향 확장거리는 0.638 m로 같았다. 연속식의 활성화 기준을 0.2와 0.4로 바꾼 경우 파괴영역 면적은 각각 1.813과 1.844 m2였으며, 첫 침투·침식 기준 도달은 각각 약 1.0과 1.2 day였다. 식 (1)의 저하계수를 0.75와 1.25로 바꾼 경우의 파괴영역 면적은 모두 1.844 m2였다. 모든 민감도 해석에서 침투·침식, 인장응력 손상, 전단변형 집중 기준에 도달한 순서는 유지됐고 최대 상향 확장거리도 같았으며, 에너지·유량수지 허용기준을 만족하였다. 다만 활성화 기준 0.2와 0.4에서 전단변형 집중 기준에 처음 도달한 누수 경과시간은 각각 약 20.6과 28.4 day로 달랐다.

본 연구의 누수 경과시간 th는 독립 수리해석의 입력상태를 구분하는 변수이며, 현장 공동의 발생시간 또는 성장속도와 동일하지 않다. 또한 본 모델은 세립분의 선택적 이탈·이송, 침식에 따른 투수특성 변화, 공동 확대에 따른 지하수 흐름의 재산정, 실제 관로의 구조거동 및 지표함몰을 해석에 포함하지 않았다. 물성 저하 활성화 기준과 저하계수는 해당 토질의 시험으로 보정하지 않았으며, 활성화 기준을 넘을 때 식 (1)의 물성 잔존계수가 불연속적으로 변하는 단순화도 포함된다. 침투·침식 임계유속 역시 현장시험으로 보정하지 않았고, 탄성계수·점착력·내부마찰각의 개별 저하계수에 대한 독립 민감도 검토도 수행하지 않았다. 이번 추가해석은 공통 물성 저하계수와 활성화 기준을 바꾼 수치적 범위검토이므로 물성식의 현장 타당성을 입증하지는 않는다. 따라서 결과의 활용범위는 제시한 경계조건과 입력범위에서 공동 확대 과정을 해석하고 수리조건별 상대응답을 비교하는 데 한정된다.

4. 결 론

1. 수리상태를 적용하지 않은 비교해석에서는 이탈기준에 도달한 요소가 없었으나, 기준 누수해석에서는 th ≈ 1.0 day의 수리상태에서 첫 요소가 이탈기준에 도달하였다. 30-day 수리상태를 고정한 최종 역학응답에서 파괴영역 면적과 초기 공동 천단 기준 최대 상향 확장거리는 각각 1.875 m2와 0.638 m였다.

2. 기준 임계유속(qcrit = 400 nm/s)에서 침투·침식 기준 요소의 누적면적 증가는 th ≈ 5.0 day 이후 제한적이었던 반면, 인장응력 손상 기준 요소의 누적면적은 누수 경과시간이 긴 조건에서도 증가하였다. 이는 침투·침식으로 자유면이 형성된 뒤 공동 주변의 응력분포가 변하고, 천단·측벽의 인장응력 손상과 누수 경과시간 약 23.2 day의 전단변형 집중이 차례로 발생하는 확대 과정과 일관된다.

3. qcrit가 400에서 600 nm/s로 증가하면 th = 8.6 day의 파괴영역 면적은 0.625 m2에서 0.219 m2로 감소하였고, 초기 공동 천단 상부에서 추가로 이탈기준에 도달한 요소는 없었다. 200 nm/s 조건은 최종 인공에너지/내부에너지비가 53.59%로 허용기준을 초과하여 정량 비교에서 제외하였다.

4. 요소크기 0.0625 m의 세분화된 메쉬에서도 첫 이탈기준 도달시점과 최대 수평 확장폭은 동일하였으나, 파괴영역 면적과 최대 상향 확장거리는 각각 11.67%와 14.71% 감소하였다. 주된 확대방향은 유지되었으나 정량적 메쉬 수렴성은 확인되지 않았다.

5. 본 모델은 토립자의 실제 이동, 공동 확대에 따른 지하수 흐름의 재산정 및 지표함몰을 해석에 포함하지 않았으며, 물성 저하계수와 침투·침식 임계유속도 현장시험으로 보정하지 않았다. 물성 저하식의 연속화와 입력값 변동에도 검토한 범위에서 이탈기준에 도달한 순서는 유지되었다. 따라서 결과의 활용범위는 현장 공동크기나 발생시간의 직접 예측이 아니라 수리이력, 지반 물성 저하, 침투·침식과 공동 주변에 집중되는 역학적 손상의 연계 과정을 해석하고 해석조건별 응답을 비교하는 데 한정된다.

Acknowledgements

본 연구는 국토교통부/국토교통과학기술진흥원의 지원으로 수행되었음(스마트 지하안전 조사·분석·평가 기술개발, 과제번호: RS-2026-25528114).

저자 기여도

반호기는 원고 작성, 연구 개념 및 설계, 데이터 분석을 하였고, 강재모는 연구내용 및 원고 검토를 하였다.

References

1

Abeyrathne, A., Jayasundara, C., Costa, S., Kodikara, J. (2023), “Deviatoric behaviour of compacted unsaturated soils within (q, vw, p) space—An experimental study”, Acta Geotechnica, Vol. 18, pp. 1651-1669. https://doi.org/10.1007/s11440-022-01608-3

10.1007/s11440-022-01608-3
2

Guo, J., Zhang, Y., Cheng, Y., Zhang, X., Shi, H., Zheng, J., Ma, Y. (2025), “Study on urban ground collapse induced by defective pipelines based on physical model experiments and numerical simulation”, Scientific Reports, Vol. 15, Art. 6085. https://doi.org/10.1038/s41598-025-90146-5

10.1038/s41598-025-90146-539971999PMC11840151
3

Guo, J., Zhang, Y., Li, Y., Zhang, X., Zheng, J., Shi, H., Zhang, Q., Chen, Z., Ma, Y. (2024), “Model experimental study on the mechanism of collapse induced by leakage of underground pipeline”, Scientific Reports, Vol. 14, Art. 17717. https://doi.org/10.1038/s41598-024-68824-7

10.1038/s41598-024-68824-739085627PMC11291715
4

Ibrahim, A., Meguid, M.A. (2023), “CFD-DEM simulation of sand erosion into defective gravity pipes under constant groundwater table”, Tunnelling and Underground Space Technology, Vol. 131, Art. 104823. https://doi.org/10.1016/j.tust.2022.104823

10.1016/j.tust.2022.104823
5

Indiketiya, S., Jegatheesan, P., Rajeev, P., Kuwano, R. (2019), “The influence of pipe embedment material on sinkhole formation due to erosion around defective sewers”, Transportation Geotechnics, Vol. 19, pp. 110-125. https://doi.org/10.1016/j.trgeo.2019.03.001

10.1016/j.trgeo.2019.03.001
6

Kwak, T.Y., Woo, S.I., Kim, J., Chung, C.K. (2019), “Model test assessment of the generation of underground cavities and ground cave-ins by damaged sewer pipes”, Soils and Foundations, Vol. 59, No. 3, pp. 586-600. https://doi.org/10.1016/j.sandf.2018.12.011

10.1016/j.sandf.2018.12.011
7

Sato, M., Kuwano, R. (2015), “Influence of location of subsurface structures on development of underground cavities induced by internal erosion”, Soils and Foundations, Vol. 55, No. 4, pp. 829-840. https://doi.org/10.1016/j.sandf.2015.06.014

10.1016/j.sandf.2015.06.014
8

Tan, F., Tan, W., Yan, F., Qi, X., Li, Q., Hong, Z. (2022), “Model test analysis of subsurface cavity and ground collapse due to broken pipe leakage”, Applied Sciences, Vol. 12, No. 24, Art. 13017. https://doi.org/10.3390/app122413017

10.3390/app122413017
9

Wang, Z.Y., Tan, Y., Liu, J.C. (2026), “Collapse of underground cavities: Initiated by cyclic pipeline leakage and triggered by rainfall infiltration”, Tunnelling and Underground Space Technology, Vol. 173, Art. 107612. https://doi.org/10.1016/j.tust.2026.107612

10.1016/j.tust.2026.107612
10

Yang, J., Yin, Z.Y., Laouafa, F., Hicher, P.Y. (2019), “Internal erosion in dike-on-foundation modeled by a coupled hydromechanical approach”, International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 43, No. 3, pp. 663-683. https://doi.org/10.1002/nag.2877

10.1002/nag.2877
11

Zhang, D.M., Du, W.W., Peng, M.Z., Feng, S.J., Li, Z.L. (2020), “Experimental and numerical study of internal erosion around submerged defective pipe”, Tunnelling and Underground Space Technology, Vol. 97, Art. 103256. https://doi.org/10.1016/j.tust.2019.103256

10.1016/j.tust.2019.103256
12

Zhang, X.S., Wong, H., Leo, C.J., Bui, T.A., Wang, J.X., Sun, W.H., Huang, Z.Q. (2013), “A thermodynamics-based model on the internal erosion of earth structures”, Geotechnical and Geological Engineering, Vol. 31, No. 2, pp. 479-492. https://doi.org/10.1007/s10706-012-9600-8

10.1007/s10706-012-9600-8
페이지 상단으로 이동하기