Ⅰ. 서 론
Ⅱ. 재료 및 방법
1. DAMBRK-HEC-RAS 연계 모형의 개요
2. 연구 대상지
3. DAMBRK를 이용한 저수지 붕괴 모의
4. HEC-RAS 2D를 이용한 저수지 하류부 홍수파 해석
Ⅲ. 결과 및 고찰
1. 저수지 붕괴 모의 결과
2. 하류부 2차원 홍수파 해석 결과
Ⅳ. 요약 및 결론
Ⅰ. 서 론
우리나라 농업용 저수지 17,047개소 중 준공년도 50년 이상인 저수지는 약 87.8%로 나타났으며 (MAFRA, 2024), 저수지의 내구연한이 50년이라는 점에서 농업용저수지의 노후화 문제는 심각한 상황이다 (Jung and Lee, 2023). 또한, 극한강우 발생 시 농업용저수지의 홍수배제능력 부족으로 제방 붕괴, 사면유실, 시설파손과 같은 시설피해가 발생할 수 있으며, 2020년 산양저수지, 북좌저수지, 2024년 법곡저수지, 2025년 송원저수지 등 저수지 붕괴 사례가 꾸준히 발생하고 있다. 농업용저수지의 노후화가 심해지고 기후변화로 인한 극한강우의 빈도 및 강도가 증가함에 따라 농업용저수지 붕괴 위험은 점차 커질 것이다.
저수지 붕괴 시 짧은 시간 내에 저류된 물이 하류부로 빠르게 전파됨에 따라 인적, 물적 피해가 크게 발생할 수 있으므로 관련된 연구가 지속적으로 수행되고 있다. Choi et al. (2011)은 DAMBRK 모형을 활용하여 저수지 붕괴방류 조건에 따른 하류 홍수분석을 수행하였으며, Kim and Han (2016)은 DAMBRK 모형을 활용하여 장현저수지와 동막저수지의 붕괴를 재현하였다. Won et al. (2010)과 Song et al. (2012)은 농업용 필댐 붕괴 시 BREACH 모형의 적용성을 검토한 바 있다.
그러나 흙댐 붕괴 시 발생하는 붕괴류는 제체를 구성하는 토사와 물이 섞여 유출되는 다상 유동이며, 일반적인 뉴턴유체와는 다른 비선형 점성 유동을 보이는 비뉴턴유체이다. 흙댐 붕괴 시 하류부는 토석류 유동, 토사 매몰 등의 피해가 발생하나 기존의 수리해석 모형을 활용하여 비뉴턴유체 유동, 다상 흐름을 고려한 홍수파 연구 사례는 드물었다.
HEC-RAS (Hydrologic Engineering Center – River Analysis System)는 2차원 홍수파 해석에 있어 높은 수리적 안정성과 많은 적용사례를 보이고 있다. Ashraf (2021)은 HEC-RAS 2D가 격자 크기가 커짐에도 지형 정보를 효과적으로 반영하여 FLO-2D보다 정확한 홍수파 해석이 가능함을 보였다. Jeon et al. (2018)은 HEC-RAS 1D와 2D 해석 모듈을 연계하여 하천 인근의 농경지 침수해석을 수행하였다. Yu et al. (2019)은 HEC-RAS 2D를 활용하여 댐 붕괴 매개변수에 따른 붕괴 유출수문곡선을 산정하고, 이에 따른 불갑저수지 하류 홍수범람도를 작성하였다. Hong et al. (2022)은 HEC-HMS와 HEC-RAS를 연계하여 김천부항댐을 대상으로 댐 붕괴 및 하류부 홍수범람 모의를 수행하였다. Kim et al. (2023)은 HEC-RAS 2D를 활용하여 어은저수지 하류부의 댐 붕괴 및 하류부 홍수범람 모의를 수행하였으며, 모의 결과를 바탕으로 홍수피해액을 산정한 바 있다.
최근 HEC-RAS의 비뉴턴유체 해석 기능 개발로 광미댐, 산사태로 발생하는 토석류에서의 적용 가능성이 확인되고 있다. Gibson et al. (2020)은 미국과 브라질의 실제 토석류 사례에 대해 HEC-RAS의 비뉴턴 흐름 해석 기능을 적용하여 도달 시간, 확산 범위, 침수심 등의 주요 변수 예측을 실측 자료와 비교하였고, 모형의 재현성과 현장 적용 가능성을 평가하였다. Gibson et al. (2022)는 Santa Barbara 및 Brumadinho 광미댐 붕괴 사례를 대상으로 항복응력과 점성계수 변화에 따른 민감도 분석을 수행하며 HEC-RAS의 비뉴턴 흐름 해석 기능을 기술적으로 검토하였다. Chang et al. (2024)은 산악지형에서 발생한 토석류의 하천 및 해안 유입 영향을 분석한 바 있다.
그러나 HEC-RAS 2D의 비뉴턴유체 해석 기능을 활용한 선행연구는 대부분 산사태 토석류나 광미댐 붕괴 등에서 발생한 토석류를 대상으로 하였으며, 저수지 또는 댐 붕괴로 인한 비뉴턴유체 기반 홍수파 해석에 HEC-RAS 2D를 적용한 연구는 제한적이다. 이에 따라 HEC-RAS 2D의 비뉴턴 유동 해석 기능을 댐 붕괴 상황에 적용하여 적절성과 재현 가능성을 검토하는 연구가 필요하다.
그러나 HEC-RAS는 강우에 의한 유입량 산정을 모형 내에서 할 수 없다는 한계가 있다 (Hong et al., 2009). 이를 개선하기 위하여 HEC-HMS로 산정한 저수지 유입량, DAMBRK로 산정한 저수지 붕괴 방류량을 경계조건으로 입력하는 형태의 수리수문모형과의 연계가 필요하다. 그중 DAMBRK 모형은 강우-유출 해석으로 산정된 유입수문곡선과 저수지 붕괴 매개변수로부터 댐 붕괴 방류수문곡선을 산정할 수 있다. DAMBRK는 댐 붕괴 모의 및 국내외 비상대처계획 수립에 적용되어 활용되고 있으며, HEC-RAS, FLO-2D 등 다른 홍수파 해석 모형과의 연계성이 입증되었다 (Go et al., 2015; Lee et al., 2017).
따라서 본 연구에서는 농업용 필댐 붕괴 시 발생하는 홍수파의 유변학적 특성을 고려하기 위하여 DAMBRK-HEC-RAS 2D 연계 모형을 구축하고, 적용성을 평가하고자 한다.
Ⅱ. 재료 및 방법
1. DAMBRK-HEC-RAS 연계 모형의 개요
DAMBRK-HEC-RAS 연계 모형을 적용하여 토사흐름을 고려한 저수지 붕괴 홍수파 해석을 위한 모의 흐름도는 Fig. 1과 같다. DAMBRK 모형에서는 강우 발생에 따른 저수지 홍수추적 및 붕괴 모의를 수행한다. 저수지 유입량, 저수지 제원, 붕괴 매개변수 등의 입력자료를 구축하여 저수지 붕괴부 방류량 및 제방 토사 유실량을 산정한다. HEC-RAS 모형은 저수지 붕괴 홍수파의 하류부 전파를 모의한다. 저수지 하류부 지형자료를 구축한 후, DAMBRK 모의 결과인 붕괴부 방류량, 제방 토사 유실량과 비뉴턴유체 매개변수 설정을 통해 흙댐 붕괴류 자료를 생성한다. HEC-RAS 2D의 홍수파 해석을 통해 도출된 침수면적, 침수심, 유속, 홍수파 도달시간 등의 결과자료를 바탕으로 저수지 하류부의 홍수 피해를 분석한다.
2. 연구 대상지
본 연구의 대상지는 충청북도 충주시 엄정면 추평리에 위치한 직동저수지 유역으로 위치는 Fig. 2(a)와 같다. Fig. 2(b), (c)는 각각 추평리의 하천도와 토지피복지도이다.
직동저수지는 1966년 준공되어 2020년 8월 집중호우로 인해 제방 붕괴가 발생하였다. 붕괴 이전의 저수지 제원은 유역면적 86.0 ha, 수혜면적 8.3 ha, 제당 높이 12.0 m, 제당 길이 88.0 m, 유효저수량 7,000 m3이다. 2020년 8월 2일 오전 3시경 집중호우에 의한 월류로 인해 길이 18 m, 높이 12 m의 제방 붕괴가 발생하였으며, 하류 농경지 97,482.0 m2가 침수되었다 (Chungju-si, 2021; Baik et al., 2022).
3. DAMBRK를 이용한 저수지 붕괴 모의
가. DAMBRK 모형
저수지 붕괴에 따른 수문곡선을 산정하기 위하여 DAMBRK 모형을 활용하였다. DAMBRK 모형은 댐의 붕괴유출수문곡선의 유도와 댐 하류부의 홍수추적을 수리학적으로 해석하기 위하여 개발되었다 (Won et al., 2010). 미국 연방재난관리청 (Federal Emergency Management Agency, FEMA), 국제대댐회 (International Commission on Large Dams, ICOLD), 유럽연합 (European Union, EU) 등에서는 DAMBRK를 댐 붕괴 모의 모형으로 추천하고 있으며, 미국 및 우리나라의 댐, 저수지 붕괴로 인한 비상대처계획 (EAP, Emergency Action Plan) 수립에 있어 많은 적용 사례를 보이고 있다 (Kim and Han, 2016; Lee et al., 2017).
DAMBRK 모형은 댐 붕괴로 인한 저수지의 유출수문곡선을 산정한 후, 연속방정식과 운동방정식으로 이루어진 Saint-Venant 방정식을 지배방정식으로 하여 비선형 유한차분법 (Finite Difference Method, FDM)을 통해 하류부 1차원 홍수추적을 수행한다 (Fread, 1988). 본 연구에서는 DAMBRK 모형을 통해 저수지의 붕괴 유출수문곡선을 산정하고, 저수지 붕괴 단면 변화를 계산하여 제방 토사 유실량을 산정하였다. 식 (1)은 넓은 마루웨어 공식으로 월류 붕괴 시의 붕괴부 유출량을 계산하고, 식 (2)는 오리피스 유량공식으로 파이핑 붕괴 시의 붕괴부 유출량을 계산한다.
여기서, C1은 넓은 마루 직사각형 웨어의 유량계수, C2는 넓은 마루 삼각형 웨어의 유량계수, C3는 오리피스 흐름에 대한 유량계수, h는 시간 t에서의 저수지 수위 (EL. m), hb는 시간 t에서의 붕괴부 바닥표고 (EL. m), 는 파이핑 시작 시의 중심부 수두 (EL. m), β는 잠수보정계수를 나타낸다.
나. 입력자료 구축
저수지 붕괴 모의를 위해 필요한 DAMBRK의 입력자료로는 저수지 유입량 자료, 저수지 제원, 저수지 붕괴에 대한 매개변수가 있다. 연구 대상지인 직동저수지는 지자체에서 관리하는 소규모 농업용 저수지로 실제 붕괴 시점에서의 저수지 수위, 저수지 유입량 등의 관측 자료가 존재하지 않으므로, 수문해석기법, 붕괴 후 조사 결과를 바탕으로 입력자료를 구축하였다.
1) 저수지 유입량 산정
저수지 붕괴 당시의 저수지 유입량 자료는 환경부의 홍수량 산정 표준지침 (2019)에 따라 강우관측자료와 강우-유출관계모형을 활용하여 산정하였다. 강우자료는 직동저수지 유역의 대표관측소인 방재기상관측 (AWS, Automatic Weather System) 엄정 관측소의 분단위 강우자료를 활용하였다. 붕괴 당시 2020년 8월 1일 오후 9시 30분부터 8월 2일 오전 10시 30분까지 13시간 동안 총 346.5 mm의 강우가 발생하였으며, 지역빈도해석 결과 엄정 관측소의 500년 빈도 지속시간 13시간 확률강우량 (304.23 mm)을 상회하는 극한강우로 나타났다.
강우 발생에 따른 홍수량 산정을 위하여 유효우량 산정 방법으로 NRCS 방법을, 합성단위도 방법으로 Clark 단위도법을 채택하였다. 유출곡선지수 (Curve Number, CN) 는 MOE (2019)에 따라 산정하였으며, 환경부의 2020년 세분류 토지피복지도와 흙토람 토양도를 활용하였다. Table 1은 CN 산정을 위해 직동저수지 상류 유역의 토지이용과 수문학적 토양군에 따른 CN과 면적 비를 정리한 것으로 CN II는 65로 산정되었다. 2020년 7월 29일부터 7월 30일까지 163.5 mm의 강우가 발생하여 선행토양함수조건에 따라 CN III를 적용하였으며, 81로 나타났다.
Table 1.
Area ratio of land cover classes and hydrologic soil groups for CN estimation
유출곡선지수 외 Clark 단위도법의 매개변수로는 도달시간과 저류상수가 있으며, 본 연구에서는 MOE (2019)에서 제안한 서경대 공식을 적용하였다. 식 (3)과 (4)는 도달시간과 저류상수를 산정하는 서경대 공식으로, 유로연장은 1.33 km, 유역의 표고차는 315 m로 나타났으며, 이를 통해 산정한 도달시간과 저류상수는 각각 0.12 hr, 0.18 hr이다.
여기서, Tc는 도달시간 (hr), K는 저류상수 (hr), L은 유로연장 (km), H는 유역 최원점과 홍수량 산정지점의 표고차 (m), A는 유역면적 (km2), α는 유역계수이다.
Fig. 3은 Clark 단위도법을 적용하여 산정한 붕괴 당시의 직동저수지 유입수문곡선이다. 2020년 8월 2일 오전 3시 40분과 오전 5시 40분에 각각 15.59 m3/s, 19.48 m3/s의 첨두값을 나타냈다.
2) 저수지 제원
직동저수지는 지자체가 관리하는 노후화된 농업용 저수지로 수위-내용적 곡선과 수위-방류량 곡선 등의 상세한 저수지 제원 수집에는 한계가 있었다. 따라서 본 연구에서는 국토지리정보원의 2020년 수치지형도 자료와 충주시 정보공개청구를 통해 제당고는 EL. 175.0 m, 제방 높이는 12.0 m로 설정하였다. 직동저수지의 여수로는 월류식 여수로로 붕괴 이전 관측된 방류량과 상세 여수로 제원의 취득이 어려우므로 여수로를 직사각형 월류부로 가정하고 월류식 여수로의 월류량 산정식인 식 (5)를 통해 수위-방류량 관계를 도출하였다. 여수로 폭은 위성영상을 통해 8 m, 여수로 표고 및 만수위는 유효저수량에 해당하는 EL. 174.0 m로 설정하였다. 이를 통해 추정한 직동저수지의 수위-내용적 곡선과 수위-방류량 곡선은 Fig. 4와 같다.
여기서, Q는 월류식 여수로의 방류량 (m3/s), C는 유량계수, L은 여수로 폭 (m), He는 접근속도수두를 포함한 총수두 (m)이다.
3) 저수지 붕괴 매개변수 설정
DAMBRK 모의 시 필요한 붕괴 매개변수로는 붕괴원인, 붕괴 시점의 저수지 수위, 저수지 붕괴형성시간, 최종 붕괴부 형상 등이 있다. Fig. 5는 직동저수지 제체 붕괴부 형상을 나타낸 것이고, Table 2는 DAMBRK 모의를 위해 설정한 매개변수를 정리한 것이다.
Table 2.
Input dam breach parameter for DAMBRK simulation
직동저수지 붕괴 원인은 극한강우로 인한 제방 월류로 설정하였으며, 2020년 7월 29일부터 7월 30일까지 발생한 강우 (163.5 mm)를 고려하여 붕괴 발생 직전 직동저수지의 저수위는 만수위로 설정하였다. 제체의 최종 붕괴부 형상은 BAI (2021)과 Baik et al. (2022)의 조사 결과를 바탕으로 붕괴 단면을 결정하였다. 저수지 붕괴형성시간은 Froehlich (2008)가 기존의 댐 붕괴 사고 사례를 바탕으로 회귀분석을 통해 개발한 식 (6)을 적용하였다. 해당 경험식은 저수지 붕괴로 인한 하류부 침수 예측을 위한 실무 가이드라인 및 연구에서 범용적으로 사용되며 높은 추정력을 가지고 있다 (Sammen et al., 2017).
여기서, tf는 붕괴형성시간, Vw는 붕괴 시점의 저수용량, hb는 붕괴 높이, g는 중력가속도를 의미한다. 수위-내용적 곡선을 활용하여 붕괴 시점의 저수위 EL. 175.0 m에 해당하는 저수용량 11,332.5 m3를 Vw로 결정하였으며, hb는 12.0 m로 설정하였다. 이에 따라 산정된 붕괴형성시간 tf는 3분으로 산정되었다.
4. HEC-RAS 2D를 이용한 저수지 하류부 홍수파 해석
가. HEC-RAS 2D 모형
HEC-RAS는 미 육군 공병단 수문공학센터 (HEC)에서 개발한 수리수문해석 모델로 현재 6.6버전까지 정식 배포되었으며, 크게 지형 (Geometry), 정상류 (Steady flow), 부정류 (Unsteady flow), 유사이동 (Sediment) 모의, 그리고 공간정보 해석 및 시각화 (RAS Mapper) 모듈로 구성되어 있다. 이외에도 수온 및 수질, 하상변동, 수리구조물의 영향, 토석류 흐름, 하천의 홍수터 침수 및 제방 범람 등을 모의할 수 있어 국내외로 홍수 모의에 널리 이용되고 있다 (Urzică et al., 2020; Jeong et al., 2024).
본 연구에서는 붕괴 초기의 급격한 유속과 수면변화를 해석하기 위하여 2차원 천수방정식 (SWE, Shallow Water Equation)을 지배방정식으로 하였다. 천수방정식은 연속방정식과 운동량방정식으로 구성되며, 물리적으로는 중력에 의해 유동이 발생하고, 바닥경사, 마찰, 압력항, 외부유입에 의해 변화하는 지표수 흐름의 거동을 나타낸다. 천수방정식은 확산파방정식과는 다르게 가속항과 관성항을 고려하기 때문에 연산시간과 모의 불안정성의 한계에도 댐 붕괴 홍수파, 급류, 토석류 등 비선형성을 띠는 유체의 유동 해석에 적합하다 (Scholtz and Chetty, 2021; Pandey et al., 2024). 2차원 천수방정식을 구성하는 연속방정식과 운동량방정식은 식 (7), (8), (9)와 같다.
여기서, t는 시간, h는 수심, u, v는 각각 x, y 방향의 유속, q는 단위면적당 유입 또는 유출량, S0x, S0y는 x, y 방향의 하상경사, Sfx, Sfy는 x, y 방향의 마찰경사, g는 중력가속도를 나타낸다.
나. 저수지 하류부 지형자료 구축
본 연구에서는 홍수파 해석 시 직동저수지 하류의 건축물, 하천, 제방 등 지형지물의 영향을 고려하기 위하여 LiDAR (Light Detection and Ranging) 기법을 이용하여 지형자료를 구축하였다. 본 연구에서 항공사진 취득을 위해 이용한 UAV (Unmanned Aerial Vehicle) 항공촬영방법은 Table 3과 같다. 항공촬영을 통해 취득한 3차원 점군모델 (point cloud)은 DJI Terra 프로그램을 통해 후처리 과정을 거쳐 0.5 m 해상도의 수치표고모델 (Digital Elevation Model, DEM)로 변환하였다. 2020년과 2025년의 수치지형도와 지적도, 위성사진을 검토한 결과 필지, 하천 형상, 구조물 등 홍수파 전파에 있어 큰 영향을 끼치는 지형적 변화는 존재하지 않는다고 판단하였으며, RAS-Mapper와 GIS를 활용하여 농경지, 하천 식생, 교량 등 일부 구조물로 인한 표고 차이를 조정하였다. 농경지의 경우 촬영 당시 재배된 농작물의 높이를 고려하여 필지 경계를 기준으로 DEM 표고를 0.5 m 낮추었으며, 하천 식생, 교량의 경우 하천 흐름이 인위적으로 차단되지 않도록 RAS-Mapper의 Terrain Modification 기능을 활용하여 표고를 조정하였다.
Table 3.
UAV-based LiDAR survey for constructing terrain data
| Flight parameter | Specification |
| Type of UAV | Matrice 350 RTK |
| Type of LiDAR | Zenmuse L2 |
| Flight date | 2025.09.09 |
| Flight time (min) | 15 |
| Flight speed (m/sec) | 10 |
| Ground sampling distance (m) | 0.1 |
Fig. 6은 항공촬영을 통해 구축한 직동저수지 하류부 DEM에 대하여 HEC-RAS 2D에 입력한 지형자료를 나타낸 것으로, 침수해석 영역은 직동저수지 제방을 상류로 하여 직동저수지 하류의 직동천을 따라 추평천 합류부 및 인근 농경지를 포함하도록 설정하였다. 해석 영역은 2 m 크기의 110,473개의 격자로 이루어졌다. 저수지 하류부에 위치한 제내지의 홍수 피해를 분석하기 위하여 4개의 분석지점 (RP#1, RP#2, RP#3, RP#4)을 선정하였다. 각 지점은 저수지 하류부의 농경지로 RP#1은 밭으로 하천 제방으로부터 31.7 m 떨어져있다. RP#2, RP#3, RP#4는 논으로 하천 제방으로부터 각각 3.6, 55.6, 44.9 m 떨어져 있다.
홍수파 흐름은 지표면의 토지피복 특성에 따라 다르게 저항이 발생하여 변화한다. 조도계수는 토지이용상태를 대변하는 수리학적 매개변수로 하도나 지표면 흐름에서 발생되는 유출상황을 예측하는 데 사용된다 (Park and Lee, 2008). Park and Lee (2008)는 Hjelmfelt (1986)의 지표면 조도계수를 우리나라 38개 소분류체계에 맞게 분류하여 적용 시 실측치와 가장 근사한 결과를 보임을 확인한 바 있다. 해석 영역의 조도계수는 환경부의 2020년 세분류 토지피복지도를 Park and Lee (2008)의 분류체계에 적용한 후, 미분류 항목인 하천에 대하여 하천설계기준 (MLTMA, 2009)에 따라 0.04로 결정하였다.
다. 토사 흐름 모의를 위한 비뉴턴유체 홍수파 해석
본 연구에서는 흙댐 붕괴 시 발생하는 물-토사 혼합 홍수파의 하류부 전파 과정을 분석하기 위하여 HEC-RAS 6.6의 비뉴턴유체 해석기능을 활용한 2차원 홍수파 해석을 수행하였다. 입력 유량 조건으로는 DAMBRK를 통해 산정한 저수지 붕괴 홍수파, 여수로 방류량, Clark 단위도법을 통해 산정한 하류 발생 홍수량을 설정하였으며, 계산 간격은 1초로 설정하였다. 홍수파의 유변학적 가정에 따라 저수지 하류부 4개 지점에서의 시간에 따른 침수심, 유속, 홍수파 도달시간 및 전체 영역에서의 침수면적, 최대 침수심을 비교하여 홍수파 전파 특성의 차이를 분석하였다.
홍수파의 비뉴턴유체로서의 거동을 모의하기 위하여 저수지 붕괴로 인한 홍수파 해석 (Gibson et al., 2020; Tang et al., 2022)에 활용된 Herschel-Bulkley 모델을 적용하였다. Herschel-Bulkley 모델에서의 비뉴턴유체의 응력-변형 관계는 식 (10)과 같다. 식 (10)의 전단응력은 식 (8)과 (9)에서의 마찰경사 (Sfx, Sfy)를 구성하는 성분이 된다.
여기서, τ는 전단응력, τy는 항복응력, K는 점성계수, 는 전단변형률, n은 유동지수를 나타낸다.
비뉴턴 유동 해석에서 Herschel–Bulkley 모델의 매개변수는 토사 체적농도 및 입도에 따라 변화한다. 따라서 제체를 구성하는 토사의 입도 특성과 붕괴 홍수파의 토사 체적농도를 기반으로 경험식을 통해 매개변수를 추정하는 절차가 필요하다. 그러나 농업용 필댐 제체 재료는 점토·실트 함량과 함수비에 따라 유변학적 특성이 크게 달라지며, 현장 조건을 대표하는 매개변수의 적용에는 불확실성이 존재한다. Tang et al. (2022)은 dam-break 상황에서의 이류 흐름을 대상으로 Herschel–Bulkley 모델의 매개변수를 실험적으로 도출하고 수리실험 결과와 검증한 바 있다. 본 연구에서는 Tang et al. (2022)에서 도출한 Herschel-Bulkley 모델의 문헌 기반 매개변수 (τy = 5.01 Pa, K = 2.04 Pa·sn, n = 0.37)를 적용하여 비뉴턴유체 가정에 따른 홍수파 해석의 적용성을 검토하고자 하였다.
또한, 실제 제체 재료의 유변학적 특성에 따른 불확실성을 고려하기 위하여, Table 4와 같이 비뉴턴유체 매개변수 조합을 구분하고 항복응력(τy), 점성계수(K), 유동지수(n)에 대한 민감도 분석을 수행하였다. Table 4의 Base case는 민감도 분석에서 각 매개변수 변화에 따른 침수심 및 유속의 상대적 변화를 비교하기 위한 기준 조건으로 설정하였다. 이를 통해 각 매개변수 변화에 따른 홍수파 거동 특성을 추가적으로 검토하였다.
Table 4.
Sensitivity analysis cases for Herschel-Bulkley model’s parameters
| Case | Base | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | ||||||
| Yield stress (τy) | 10 | 0 | 5 | 50 | 100 | 10 | |||||||||||||
| Consistency index (K) | 10 | 0.5 | 1 | 50 | 100 | 10 | |||||||||||||
| Flow index (n) | 0.5 | 0.25 | 0.75 | 1 | 1.5 | ||||||||||||||
한편, 현재 HEC-RAS 6.6은 비뉴턴유체 해석 시 Herschel–Bulkley 모델의 유변학적 매개변수를 계산영역 내 모든 흐름에 동일하게 적용한다. 따라서 본 연구에서 입력한 저수지 붕괴 홍수파뿐만 아니라 토사 혼입이 적을 것으로 예상되는 여수로 방류량 및 하류 유역 홍수량에도 동일한 비뉴턴유체 매개변수가 적용된다. 이에 따라 하류 전파 과정에서 서로 다른 입력조건으로 유입되는 물 흐름에 따른 희석 효과와 토사 농도의 공간적 변화를 반영하는 데 한계가 있다.
라. 홍수 피해면적에 대한 재현성 평가
DAMBRK-HEC-RAS 연계 모형의 홍수 피해의 재현성은 LSSI (Lee Sallee Shape Index)를 활용하여 평가하였다. LSSI는 참값과 추정값에 대하여 참값과 추정값의 합집합에 대한 교집합의 비를 계산하여 공간적인 위치 부합정도를 측정하는 지수로 식 (11)을 통해 산정된다 (Lee and Sallee, 1970; Cho and Park, 2004).
여기서, A는 실제 피해면적, B는 홍수파 해석을 통해 산정된 피해면적을 나타낸다.
LSSI를 계산하기 위하여 침수흔적도를 참값으로, HEC-RAS 2D 모의를 통해 산정된 피해범위를 추정값으로 정하였다. 참값이 되는 침수흔적도는 2020년 9월 11일 한국국토정보공사에서 침수흔적 조사 및 측량을 통해 작성하였으며 (Chungju-si, 2021), 행정안전부의 재난안전데이터공유플랫폼을 통해 연구 대상지의 침수흔적도 자료를 취득하였다.
LSSI 값은 0에서 1까지의 범위로 나타나며, LSSI 값이 1에 가까울수록 관측된 피해범위와 HEC-RAS 2D로 모의된 피해면적의 일치도가 높음을 의미한다. Lee et al. (2019)에서 제시한 LSSI의 평가기준은 Table 5와 같다. Lee et al. (2019)는 대한민국 10개 행정구역에 대하여 침수흔적도와 홍수피해예측지도 간의 LSSI 산정을 수행하였으며, 10개 지역의 평균인 0.252를 중간단계로 하여 단계별 적정범위를 산정하였다.
Ⅲ. 결과 및 고찰
1. 저수지 붕괴 모의 결과
Fig. 7은 직동저수지 붕괴 과정에서의 저수지 유입량, 붕괴부 방류량, 여수로 방류량, 토사 유실량을 나타낸 것이다. DAMBRK 모의 결과 월류에 의한 직동저수지의 붕괴는 붕괴 추정 시각인 2020년 8월 2일 오전 3시부터 오전 5시 사이에 해당하는 오전 3시 46분에 발생하였다. 저수지 붕괴 직후 최대 189.6 m3/s의 유량이 방류되어 붕괴형성시간동안 총 13,501 m3의 물이 유출되었으며, 저수지 제체로부터 유실된 토사량은 총 3,288 m3로 나타났다.
2. 하류부 2차원 홍수파 해석 결과
가. 흙댐 붕괴 홍수파의 유변학적 가정에 따른 홍수파 해석 결과
1) 유변학적 가정에 따른 홍수파 모의 결과 비교
뉴턴유체와 비뉴턴유체 홍수파 해석 시 댐 붕괴 이후 시간에 따른 침수면적 및 침수심 변화는 Fig. 8과 같다. 모의 결과 댐 붕괴 이전에도 여수로 배수량, 하류 유역의 홍수량으로 인해 농경지 저류로 인한 침수가 발생하였다. 이로 인해 저수지 붕괴 홍수파는 저류된 농경지의 침수심을 상승시키는 역할을 하며, 뉴턴유체와 비뉴턴유체 모두 유사한 시기에 각 지점별 붕괴 홍수파가 통과하는 것을 확인할 수 있다.
Fig. 9는 유체의 유변학적 가정에 따른 홍수파 전파 시 분석지점별 수심과 유속 변화를 도시한 것으로 모든 지점에서 뉴턴유체와 비뉴턴유체의 침수심 및 유속 변화는 전반적으로 유사한 경향을 보이며, 첨두 시점 및 첨두값의 차이는 크지 않은 것으로 나타났다.
RP#1의 경우 뉴턴유체에서의 최고 침수심, 최대 유속은 저수지 붕괴 2분 후 발생하였으며, 각각 2.01 m, 1.47 m/s로 산정되었다. 비뉴턴유체에서도 최고 침수심, 최대 유속은 뉴턴유체와 동일한 시점에서 발생하였으며, 각각 2.03 m, 1.46 m/s로 산정되었다. RP#1이 붕괴부 인근의 상류 지점으로서 홍수파가 빠르게 전파됨에 따라 유체의 유변학적 가정에 따른 홍수파의 첨두 시점 및 첨두 값의 차이는 미비하였다. 반면 첨두 이후 감쇠 구간에서는 비뉴턴유체의 경우 유변학적 저항으로 인해 유속이 뉴턴유체에 비해 낮게 나타났다.
RP#2의 경우 뉴턴유체에서의 최고 침수심, 최대 유속은 저수지 붕괴 3분 후 발생하였으며, 각각 2.69 m, 2.44 m/s로 산정되었다. 비뉴턴유체에서도 최고 침수심, 최대 유속은 뉴턴유체와 동일한 시점에서 발생하였으며, 각각 2.69 m, 2.39 m/s로 산정되었다. RP#3의 경우 뉴턴유체에서의 최고 침수심, 최대 유속은 저수지 붕괴 6분 후 발생하였으며, 각각 2.06 m, 1.98 m/s로 산정되었다. 비뉴턴유체에서도 최고 침수심, 최대 유속은 뉴턴유체와 동일한 시점에서 발생하였으며, 각각 2.07 m, 1.91 m/s로 산정되었다. RP#2와 RP#3의 경우 유체 가정에 따른 침수심과 유속 분포의 큰 차이가 나타나지 않았다.
RP#4의 경우 뉴턴유체에서의 최고 침수심과 최대 유속은 9분 후 발생하였으며, 각각 약 2.05 m, 0.63 m/s로 나타났다. 비뉴턴유체에서는 동일 시점에서 각각 약 2.07 m, 0.61 m/s로 나타났다. RP#4에서는 상류 지점에 비해 홍수파 도달이 지연됨에 따라 유변학적 영향이 존재하였고, 유속이 비뉴턴유체에서 더 빠르게 감소하는 것을 확인할 수 있다.
모든 지점의 침수심, 유속의 변화를 분석한 결과 뉴턴유체와 비뉴턴유체의 차이는 첨두값의 발생 시점보다는, 수심, 유속에 대한 첨두 값의 차이, 감쇠부에서 차이가 나타나는 것으로 보인다. 다만 분석지점이 모두 농경지에 분포하고 있어 저류로 인해 홍수파 전파의 특징이 뚜렷하게 나타나지 않은 것으로 보인다.
2) 침수흔적도와의 재현성 평가
행정안전부에서 제공하는 2020년 발생한 직동저수지 하류부의 침수흔적도와 HEC-RAS 2D를 통해 도출된 홍수 피해면적을 Fig. 10에 도시하였다. 뉴턴유체와 비뉴턴유체 해석의 홍수 피해면적은 각각 169,792.6 m2와 169,964.8 m2로 산정되었으며, 저수지 붕괴로 인한 홍수파가 침수흔적도에 기록된 침수영역보다 하류부로 전파되었다. 이는 실제 피해 당시 위성사진으로 기록된 침수흔적도 하류 방향으로도 피해가 발생하였으나 침수흔적도가 지적도를 기반으로 하류부 피해를 집계하지 않고 기록된 것으로 파악된다. 따라서 침수흔적도 작성 범위의 최하류 필지를 기준으로, 직동교부터 추평천 합류부까지의 하류부 초과 전파 침수면적을 제외한 후 LSSI를 산정하였다. 뉴턴유체와 비뉴턴유체 해석에서의 침수흔적도에 대한 LSSI는 각각 0.781, 0.784로 Excellent에 해당하였다.
피해면적에 대한 재현성은 우수한 것으로 나타났으나, 침수심의 경우 모의치와 침수흔적도의 차이가 나타났다. 침수흔적도 상 RP#1에 해당하는 저수지 붕괴부 인근 지역의 침수심은 3.05 m, 하류부로 전파될수록 침수심은 2.62, 2.26 m로 작아지는 것으로 기록되었다. 각 지점별 최고 침수심은 RP#1, RP#2, RP#3, RP#4가 각각 뉴턴유체 기준 2.01, 2.69, 2.06, 2.05로 나타났다. 붕괴부 직하류의 경우 RP#1의 경우 1 m 가량의 침수심 차이가 발생하였으나, 나머지 지점의 경우 침수흔적도와 유사하게 최고 침수심이 모의되었다.
나. 비뉴턴유체 매개변수 민감도 분석 결과
1) 항복응력 (τy)
Fig. 11은 항복응력에 따른 지점별 수심 및 유속 변화를 나타낸 것이다. 항복응력에 따른 민감도 분석 결과, 항복응력의 증가는 모든 관측지점에서 홍수파의 전파를 방해하고, 최대 유속을 감소시키는 것으로 나타났다. 이러한 경향은 상류보다 중하류로 갈수록 뚜렷하게 나타났다.
상류에 위치한 RP#1에서는 항복응력이 0 Pa일 때 최대 유속이 약 1.47 m/s로 나타났으며, 항복응력을 100 Pa까지 증가시켜도 최대 유속은 1.40 m/s로 약간 감소하였다. 또한, 홍수파 도달시간의 차이는 거의 없었다. 이는 붕괴 직후 홍수파가 빠르게 전파되기 때문에 항복응력이 미치는 영향이 제한적인 것으로 보인다.
또한 경사지에 위치한 RP#2에서도 RP#1과 유사한 경향이 나타났다. 항복응력이 증가함에 따라 최대 유속은 2.38 m/s에서 2.3 m/s로 감소하였으며 홍수파 도달시간 역시 붕괴 후 3분을 기준으로 큰 차이가 없었다. 이는 경사가 큰 구간에서 중력에 의한 홍수파의 전파가 지배적으로 작용하여 항복응력의 영향이 작게 나타났기 때문으로 해석된다.
반면 완만한 하류부에 위치한 RP#3와 RP#4에서는 항복응력 증가에 따른 유속 저감과 홍수파의 지연이 뚜렷하게 나타났다. RP#3의 경우, 항복응력이 0 Pa일 때 최대 유속은 1.98 m/s였으나, 항복응력이 5 Pa로 증가하면 1.91 m/s, 50 Pa에서는 1.64 m/s, 100 Pa에서는 1.35 m/s까지 단계적으로 감소하였다. 홍수파의 도달시간은 50, 100 Pa에서 0 Pa 대비 약 0.5분 지연되는 것으로 나타났다.
RP#4에서 항복응력 0 Pa에서 최대 유속은 0.63 m/s로 나타났으나, 항복응력이 10 Pa일 때 0.61 m/s, 100 Pa에서는 0.50 m/s로 감소하여 RP#3와 유사하게 첨두값이 점차 감소하는 경향을 보였다. 또한, 항복응력이 커짐에 따라 홍수파 도달시간이 단계적으로 지연되는 것으로 나타났다. 이는 항복응력 증가에 따라 유체가 흐르기 위해 필요한 최소 전단응력이 커지면서 홍수파 전파가 억제되었기 때문이다.
수심의 경우, 항복응력 증가에 따라 하류부에서 뚜렷한 증가 경향과 감쇠 지연 현상이 확인되었다. RP#4에서는 항복응력이 0 Pa일 때 최대 수심이 2.052 m였으나, 항복응력을 50 Pa로 증가시키면 2.068 m, 100 Pa에서는 2.074 m까지 증가하였다. RP#3에서도 항복응력 증가에 따라 최대 수심이 증가하고, 최대 수심 이후 수심 감소가 완만해지는 경향이 나타났다. 이는 항복응력이 증가할수록 홍수파의 전파 속도가 저하되고, 하류부에서 유량이 정체되기 때문이다.
2) 점성계수 (K)
Fig. 12는 점성계수에 따른 지점별 수심 및 유속 변화를 나타낸 것이다. 점성계수에 따른 민감도 분석 결과, 점성계수도 항복응력과 유사하게 홍수파의 전파와 최대 유속을 변화시키는 매개변수로 나타났다. 점성계수가 증가함에 따라 모든 관측지점에서 유속이 전반적으로 감소하는 경향이 확인되었으며, 이러한 경향은 중·하류 지점인 RP#3와 RP#4에서 특히 두드러졌다.
RP#1에서는 점성계수가 0.5 Pa·sn일 때 최대 유속이 1.46 m/s로 나타났으며, 점성계수를 100 Pa·sn으로 증가시켰을 때 최대 유속은 1.45 m/s로 감소하여 변화가 거의 없었다. 항복응력과 마찬가지로 저수지 붕괴부 인근에서 발생하는 빠른 붕괴류로 인해 점성계수의 영향력이 작은 것으로 사료된다.
RP#2에서도 점성계수 증가에 따른 영향은 제한적으로 나타났다. 점성계수가 증가함에 따라 최대 유속은 약 2.39 m/s에서 2.30 m/s 수준으로 감소하였으며, 도달시간 역시 약 3분 전후로 큰 변화는 나타나지 않았다. 항복응력과 마찬가지로 중력에 의한 흐름 가속이 우세하여 점성의 영향이 상대적으로 작게 나타나기 때문으로 판단된다.
RP#3과 RP#4에서는 점성계수 증가 시 최대 유속이 크게 저감되는 것으로 나타났다. RP#3에서는 점성계수가 0.5 Pa·sn일 때 최대 유속이 1.89 m/s였으나, 점성계수가 10 Pa·sn일 때 1.81 m/s, 50 Pa·sn일 때 1.55 m/s, 100일 때 1.33 m/s로 크게 감소하였다. RP#4에서도 점성계수 0.5 Pa·sn에서 최대 유속은 0.61 m/s였으나, 100 Pa·sn으로 증가 시 0.47 m/s까지 낮아졌다. 이는 점성계수 증가로 인해 전단 변형 과정에서의 에너지 소산이 증대되어 홍수파의 전파와 가속이 억제된 것으로 보인다.
수심의 경우, 점성계수 증가에 따라 하류부에서 최대 수심이 증가하고, 높은 수심이 장시간 유지되었다. RP#4에서 점성계수 0.5 Pa·sn일 때 최대 수심이 2.08 m였으나, 점성계수 50 Pa·sn에서 2.21 m, 100 Pa·sn에서 2.29 m로 증가하였다. 하류 지점으로 갈수록 점성계수 증가 시 홍수파의 정체가 이루어지며 침수심이 상승하는 것으로 나타났다.
3) 유동지수 (n)
Fig. 13은 유동지수에 따른 지점별 수심 및 유속 변화를 나타낸 것이다. 유동지수에 따른 민감도 분석 결과, 유동지수의 영향은 항복응력과 점성계수에 비해 상대적으로 완만하게 나타났으며, 매개변수에 따른 침수심의 차이도 크게 확인되지 않았다.
RP#1부터 RP#3에서는 유동지수가 증가함에 따라 첨두 이후 유속이 더 빠르게 감쇠하는 경향이 나타났다. RP#1의 경우 유동지수가 0.25일 때 최대 유속은 약 1.47 m/s였으며, n이 증가함에 따라 최대 유속은 소폭 감소하였다. RP#2와 RP#3에서도 유사하게 n 증가 시 첨두 이후 유속 저감이 확인되었으며, 특히 RP#3에서는 n이 0.25일 때 최대 유속이 약 1.98 m/s였으나, n이 증가함에 따라 최대 유속 및 감쇠 후 유속이 점진적으로 감소하였다.
반면 RP#4에서는 유동지수 변화에 따른 유속 차이가 상대적으로 작게 나타났다. RP#4의 최대 유속은 n 변화에도 약 0.6 m/s 내외로 유사하게 나타났으며, 첨두 이후 감쇠 곡선 또한 대부분 중첩되는 양상을 보였다. 이는 RP#4는 농경지 침수가 계속 진행되어 홍수파가 이미 상당 부분 감쇠된 이후 도달하는 지점이기 때문에, 유동지수 변화에 따른 추가적인 유속 변화가 제한적으로 나타난 것으로 판단된다.
수심의 경우 모든 지점에서 유동지수 변화에 따른 차이는 크지 않았으며, 최고수심과 도달시간 또한 거의 유사하게 나타났으며, 유동지수 변화에 따른 홍수파 도달 시점은 모든 지점에서 차이가 발생하지 않았다. 따라서 유동지수는 유체의 전파 시점 자체보다는 첨두 이후 홍수파 전파에 대한 감쇠현상에 영향을 주는 매개변수인 것으로 판단된다.
Ⅳ. 요약 및 결론
본 연구는 농업용 필댐의 붕괴 시 발생하는 토사 유출을 고려하여 하류부의 홍수파 특성과 피해 양상을 분석하였다. 이를 위해 충청북도 충주시 직동저수지를 대상으로, DAMBRK-HEC-RAS 2D 연계 모형을 구축하여 저수지 붕괴 시의 유출수문곡선과 제체 토사 유실량을 산정하고, 유변학적 특성을 고려한 2차원 홍수파 해석을 수행하였다.
DAMBRK 모형을 활용하여 붕괴 시나리오별 방류량 및 제체 유실량을 산정한 결과, 최대 방류량은 약 189.6 m3/s, 제체 토사 유실량은 약 3,288 m3로 산정되었다. HEC-RAS 2D 해석에서는 유체의 유변학적 거동을 모의하기 위해 Herschel-Bulkley 비뉴턴유체 모델을 적용하여 하류부 홍수파 해석을 수행하였다. 홍수파 해석 결과, 비뉴턴유체 적용 시 첨두 이후 감쇠부와 침수 지속성에 영향을 주는 것으로 확인되었다. 첨두 홍수파 이후 감쇠 구간에서 비뉴턴유체는 항복응력 이하의 전단응력 조건에서 흐름의 체류가 발생하여 침수 지속시간이 증가하였다. 침수흔적도와의 LSSI는 모두 0.78 이상으로 Excellent 등급을 나타내어 침수영역의 재현은 우수했으나, 전반적인 모의결과 침수심이 일부구간에서 침수흔적도 대비 낮게 모의되는 한계가 있었다.
비뉴턴유체 모델의 매개변수 민감도 분석 결과 항복응력과 점성계수가 증가할수록 하류에서 최대 유속은 크게 감소하고, 홍수파 도달시간은 지연되며, 최대 침수심과 침수 지속시간은 증가하는 것으로 나타났다. 이는 토사 농도가 증가할수록 유체의 흐름에 필요한 전단응력이 커져 홍수파의 전파가 방해되기 때문으로 해석된다. 반면 유동지수는 홍수파의 도달 시점이나 최대 침수심에는 상대적으로 작은 영향을 미치며, 주로 홍수파의 유속 감쇠부에 영향을 주는 것으로 나타났다.
본 연구에서 구성한 DAMBRK–HEC-RAS 2D 연계 모형은 농업용 필댐 붕괴류의 유변학적 가정을 고려하여 토사 혼합 홍수파의 전파와 침수 지속성을 정성적·정량적으로 분석하는 데 활용되었다. 분석 결과, 비뉴턴유체 가정은 홍수파의 전파 속도, 유속 감쇠 및 침수 지속시간에 영향을 미칠 수 있음을 확인하였다. 다만 본 연구에서는 현장 제체 재료의 유변학적 특성 및 붕괴 당시 홍수파 관측자료가 확보되지 않아, 모형의 적용성 평가는 침수흔적도와의 공간적 재현성 및 문헌 매개변수 기반의 민감도 분석을 중심으로 수행하였다. 또한 현재 HEC-RAS 6.6의 비뉴턴유체 해석에서 Herschel–Bulkley 모델의 유변학적 매개변수가 계산영역 전체에 동일하게 적용됨에 따라, 하류 전파 과정에서 발생할 수 있는 토사농도의 공간적 변화와 물 흐름 유입에 따른 희석효과를 직접 반영하는 데 한계가 있다. 향후 붕괴 당시의 수문·수리 조건, 토사농도 변화 및 제체 재료의 유변학적 특성을 고려한 매개변수 설정 및 적용 방안을 검토함으로써, 현장 조건을 보다 적절히 반영한 비뉴턴유체 해석이 가능할 것으로 판단된다. 본 연구 결과는 향후 농업용 저수지 붕괴에 따른 홍수–토사 복합재해 위험도 평가, 비상대처계획(EAP) 수립, 그리고 홍수 피해 저감 대책 마련을 위한 기초 자료로 활용될 수 있을 것으로 기대된다.
