1. 서 론

최근 빈발하고 있는 대규모 홍수와 자연재해는 하천 흐름에 대한 정확도가 높은 수치해석 모델에 대한 관심의 증가로 이어지고 있다. 현재 하천에서 발생하는 일반적인 흐름은 기존에 개발된 여러 형태의 수치기법에 의해 해석되고 있다. 그러나, 연속적이지 않은 형태의 흐름을 해석하거나, 정확한 해석을 필요로 하는 경우에는 기존의 수치해석기법은 많은 한계를 보여 주고 있다. 보통 흐름을 위한 수치해석기법은 공간차분형태에 따라 크게 중앙차분법과 상류이송 기법(upwind scheme)으로 분류될 수 있다. 중앙차분기법은 각 계산지점의 특성만을 반영하여 계산하는 기법으로 차분과정이 비교적 간단하여 지금까지 많이 사용되어 왔으나, 도수, 댐 붕괴류 등 불연속이 나타나는 흐름의 경우 비물리적인 진동이 발생하여 해에 수렴할 수 없게 된다. 반면 상류이송기법은 각 지점의 특성선 특성에 따라 공간차분을 변화시켜 주기 때문에 상류-사류-상류 등으로 변화하는 천이류 흐름을 계산할 수 있는 장점을 가지고 있다. 이로 인하여 댐 붕괴류, 보 월류흐름 등 불연속 흐름의 계산을 위해서 상류이송기법이 많이 이용되고 있으나, 하상경사, 마찰경사, 단면변화 등을 나타내는 생성항(source term)에 대한 보존 특성을 만족하는 수치적 처리의 어려움으로 인하여, 이를 극복할 수 있는 하천 흐름해석 모델의 개발이 요구되고 있다(Kim and Han, 2009).

개수로 흐름의 해석에는 유한차분법(Finite Difference Method), 유한체적법(Finite Volume Method), 유한요소법(Finite Element Method) 등 다양한 수치기법이 사용되고 있다. 유한차분법은 차분과정이 비교적 간단하고 빠르며 고차다항식 적용이 가능하다는 장점이 있으나, 2차원 이상의 복잡한 지형에 대한 해석에 어려움이 있으며, 고차 다항식을 적용할 경우 많은 격자점(grid point)을 필요로 하여 경계면의 처리에 불리하다. 유한체적법은 복잡한 지형에 대한 해석이 가능하나, 고차 정확도의 적용이 어렵다. 유한요소법의 경우에는, 표준적인 갤러킨 유한요소법이 이송(advection)-지배적인 흐름에 대한 방정식의 해석에서 불필요한 진동(wiggle)으로 적용에 한계가 있다는 것이 알려진 이후로, SUPG (Hughes and Brooks, 1982), CBS (Characteristic Based Split, Zienkiewicz and Codina, 1995) 등의 알고리즘이 널리 사용되었다.

불연속 갤러킨 유한요소법은 Reed and Hill (1973)에 의하여 중성자 이송방정식(neutron transport equation)의 모의에 처음 적용된 이래로, 1990년대 이후 Runge-Kutta 불연속 갤러킨(RKDG) 유한요소법으로 적용이 증가하고 있다. 불연속 갤러킨 유한요소법에서는 개별 요소에 대한 약형식(weak form)의 결과로서 나타나게 되는 요소 경계에서 수직으로 표현되는 흐름률을 Riemann 근사해법에 의하여 구하는 것이 일반적이다. 기존의 안정화(stabilized) 유한요소법(예를 들면, SUPG, 특성선 기반 유한요소법(characteristic Galerkin Method) 등)과는 달리, Fig. 1에 보인 바와 같이 요소와 요소사이에 불연속성을 허용하므로 (연속)유한요소법에 비하여 충격파(shock wave)나, 불연속성을 포함하는 문제 해석에 있어서 매우 유리하다. 또한, 이송방정식의 처리에 있어서, 상류이송기법에 의하여 산정되는 흐름률 자체가 필요한 안정성을 제공하기 때문에 안정화 기법에 사용되는 추가적인 안정항(stabilization term)이 필요하지 않다는 점을 장점이라고 할 수 있다. 특히, 개별 요소를 중심으로 계산이 이루어지며, 인접한 요소는 흐름률 만으로 연결되므로, 유한요소법의 대표적인 특징인 전체 강성행렬(global stiffness matrix)과 전체 하중벡터(global force vector)를 구성하는 조립과정(global assembly)이 필요하지 않으며, 2차원, 3차원의 다차원 문제의 경우에는 영역분할(domain decomposition) 방식의 병렬기법의 적용이 용이한 것으로 알려져 있다(Cockburn, 1999). 다만, 요소의 경계에서 둘 이상의 변수가 지정되어 변수의 수가 많아지며, CFL (Courant-Friedrichs-Lewy) 제약조건이 전통적인 (연속)유한요소법에 비하여 엄격히 적용되는 등의 단점이 있다.

천수방정식에 대한 적용은 Chavent and Salzano (1982)에 의하여 시도되었으나, 시간적분에 사용한 Euler 기법으로 인하여 안정적인 계산결과를 제시하지 못하였다. 그후 적분의 정확도를 향상시키기 위하여 공간차분에 고차의 불연속 갤러킨 기법을, 시간적분에 압축성 유체해석 분야의 성과를 반영한 Runge-Kutta기법을 결합하여, Runge-Kutta 불연속 갤러킨(RKDG)기법으로 발전하게 되었다. 최근 Schwanenberg and Harms (2004)에 의하여 2차원 천수방정식에 적용되었으며, Lee and Lee (2013)에 의하여 2차원 사류(supercritical flow), 천이류의 모의에 적용된 바 있다.

요소 경계에서의 흐름률을 계산하기 위한 리만(Riemann) 근사해법으로는 전통적으로 Lax-Friedrichs, Roe, HLL (Harten-Lax- van Leer) 기법이 사용되어 왔으며, 본 연구에서도 이를 적용하였다. 또한, 대부분의 고차 수치기법은 수치해에 수치진동(numerical oscillations)이 발생할 수 있다. 특히, 수심 및 유속 변화의 기울기가 큰 불연속 흐름주위에는 비현실적인 수치진동이 발생될 수 있기 때문에 이를 제거하는 기법을 찾아야 되는 경우가 많다. 수치진동을 제거하는 가장 널리 사용되는 기법은 흐름률 제한자(flux limiter), 기울기 제한자(slope limiter), ENO 차분기법이 있으며, 본 연구에서는 기울기 제한자를 사용하였다.

2차원/3차원 모델에 대한 기술적 요구와 더불어, 1차원 하천 모델은 방재 및 예보, 장기간의 변화 예측 등 향후에도 지속적인 연구과 개발이 요구된다. 본 연구에서는 1차원 천수방정식에 대한 TVD (Total Variation Diminishing) RKDG 기법의 수치모형을 수립하고, 이를 댐 붕괴류의 해석에 적용하였다. 수치진동을 제어하기 위하여 기울기 제한자로 MUSCL 기울기 제한자를 적용하였다. 먼저 본론에서 지배방정식과 수치 기법에 대하여 기술하고, 적용사례로서 삼각형, 사각형 단면에서 댐 붕괴류 흐름과 지형변화에 따른 천이류 흐름을 모의하였으며, 적용 결과와 향후 연구에 대하여 기술하였다.

2. 본 론

2.1 지배방정식

흐름 현상 중에서 많은 경우에는, 수심에 비하여 수평방향 스케일이 크며, 이 경우에 천수방정식(shallow water equation)을 지배방정식으로 적용하여 왔다. 또한, 폭에 비하여 길이가 긴 하천 흐름에 대해서는 1차원 천수방정식인 Saint-Venant 방정식을 이용하는 경우가 많다. Saint-Venant 방정식은 연속방정식과 운동량 방정식으로 다음과 같이 표현할 수 있다(Cunge et al., 1980).

 (1)

여기서, 벡터 U, F, S 는 다음과 같이 정의한다.

, ,   (2)

A는 흐름 면적을, Q는 유량()을, 는 중력 가속도를 나타내며, 그 이외의 기호는 Fig. 2에 정의된 바와 같다. 는 기준면(datum)에 대한 수면의 경사, 는 마찰경사로서 다음과 같이 Manning 공식으로 정의한다.

 (3)

여기서 은 Manning계수이며, 로 동수반경, 는 윤변(wetted perimeter)이다.

지배방정식인 Eq. (1)에 대한 고유치(eigenvalue)와 고유벡터(eigenvector)는 다음과 같다.

, (4)

여기서 는 파속으로 정의된다(Fig. 2). 이미 기술한 바와 같이, 본 연구에서 흐름률 항은 리만 근사해법으로 산정하며, 리만 근사해법은 고유치로 주어지는 특성곡선 방향으로 풍상효과(upwinding effect)를 반영한다. 사실상 고유치가 특성속도, 즉, 교란이 하도를 따라 진행하는 속도를 나타낸다.

2.2 불연속 갤러킨 유한요소법

불연속 갤러킨 유한요소법에서는 (연속)유한요소법에서와 같이 전체영역 를 개의 요소와 개 절점(node)으로, Fig. 3과 같이 겹치지(non-overlapping)않도록 구성한다. 즉, 전체영역 가 일 때, 절점의 분할은 다음과 같다.

 (5)

각 요소의 국지(local) 좌표 에서, 는 요소의 시작 좌표, 는 요소의 끝 좌표를 의미하며 수로의 형상과 조도계수(roughness coefficient)는 각각의 절점에서 지정된다. 요소 절점 값은 형상함수(shape function)로 알려진 선형 보간함수를 이용하여, 단면적 , 유량 는 다음과 같이 표현할 수 있다.

  ;   (6)

여기서 과 는 선형 형상함수(linear shape function)로서, 이고 이다. 은 요소를 구성하는 양단의 절점 기호이다.

천수방정식 Eq. (1)에 대한 수치근사해(numerical approximate solution)를 아래첨자 ‘h’로 표기하면, 다음과 같다.

 (7)

유한요소법에서와 같이 각 요소에서 가중함수(weight function) 를 곱하여 적분하면 다음 Eq. (8)을 얻게 된다.

 (8)

Eq. (8)의 에 대하여 부분적분(integration by

parts)을 하면 다음과 같은 약형식(weak form)을 얻을 수 있다.

 

 (9)

여기서 은 각각의 요소의 경계(boundary)이다. Fig. 3과 같이 각 요소의 경계인 절점(node)에서는 모든 변수가 두 개의 값이 정의되므로, 어떤 값을 사용할 것인가, 또는, 어떤 조합을 사용할 것인가의 문제가 대두된다. 이 조합을 ‘수치 흐름률(, numerical flux)’로 정의하고, 유한체적법에서와 같이 Riemann 근사 해법을 이용하여 산정한다. 위 Eq. (9)의 약형식에서 좌변 두 번째 항을 다시 부분적분하면 다음과 같은 강형식(strong form)을 얻을 수 있다(Hesthaven and Warburton, 2007).

 (10)

 

여기서 는 요소 내부 변수를 이용하여 정의한 흐름률 항이다. (연속)유한요소법에서는 원래 방정식을 ‘strong form(강형식)’이라 하고, 이에 대한 가중잔차형(weighted residual form)을 약형식이라 한다. 이와 달리 불연속 갤러킨 기법에서는, 가중적분(weighted integral) 내부에 원래 방정식을 포함하고 있기 때문에, weak form과 구분하기 위하여 strong form이라 한다. 대체로 weak form을 사용하든, strong form을 사용하든 결과에 차이가 없으나, strong form에는 가중함수(또는 test function)의 미분함수가 포함되어 있지 않으므로, 가중함수에 대한 평활조건(smoothness condition)이 제약조건으로서 필요하지 않다. 따라서, 보다 일반적인 collocation 기법의 이론 전개에도 유리하다고 알려져 있으며(Hesthaven and Warburtin, 2007), 이러한 이유로 선호되기도 한다. 본 연구에서는 Eq. (10)을 기본 지배방정식으로 사용하였다. 갤러킨 기법에 의하여 가중함수에 형상함수(shape function)를 사용하면, 다음과 같은 식을 얻을 수 있다.

 

(11)

여기서 는 질량행렬(mass matrix),

로 이류행렬(advection matrix)이다. 과

은 각각 요소의 오른쪽, 왼쪽 경계이다.

2.3 Riemann 근사 해법, 시간적분, 기울기 제한자

1차원 천수방정식에 대한 불연속 갤러킨 기법의 적용에서는, 요소 경계의 값들이 불연속적이며, 수치 흐름률(Eqs. (9)~(11)의 ) 값을 경계에서 추정해야 하므로, Riemann 문제로 간주할 수 있다(Lai and Khan, 2012). 이와 같이 Riemann 문제에 대한 수치 흐름률을 계산하기 위해 대표적인 Riemann 근사 해법으로 Lax (1954)은 Lax-Friedrichs 흐름률 기법을, Roe (1981)는 Roe 흐름률 기법을 제시하였고, Harten et al. (1983)은 HLL 흐름률 기법을 제시하였다. 본 논문에서는 각각의 흐름률 기법을 사용하여 요소 경계에서의 수치 흐름률을 계산하였으며, 이에 대하여 개략적으로 설명한다. 또한, 시간적분 기법과 기울기 제한자에 대하여 소개한다.

2.3.1 Lax-Friedrichs 흐름률 기법

Lax-Friedrichs 흐름률 기법을 이용한 수치 흐름률의 산정은 다음과 같다(Hesthaven and Warburton, 2007).

 (12)

 (13)

여기서 와 는 Fig. 3에 보인 바와 같이 각각 절점(node)을 중심으로 왼쪽 요소의 오른쪽 값과 오른쪽 요소의 왼쪽 값으로서, 다음과 같이 정의한다.

 (14)

, 로 각각 정의한다.

2.3.2 Roe 흐름률 기법

Roe 흐름률 기법에서는 원래의 비선형 보존법칙(nonlinear conservation law)이 Roe 자코비안 행렬이라는 일정 계수를 가진 선형관계식으로 대체된다. 근사적 흐름률은 다음과 같이 주어진다(Lai and Khan, 1012).

 (15)

여기서 와 , 은 다음과 같다.

,   (16)

,       (17)

,         (18)

 (19)

 (20)

,      (21)

여기서 , , , 각각 경계에서 단면적과 유량으로서, Eq. (14)의 정의에 따른다. Roe 흐름률 기법에서 요소의 경계 왼쪽(L) 하상, 또는 오른쪽(R) 하상이 마른(dry)하상일 경우 파의 전파속도는 각각 다음 Eqs. (22) and (23)과 같다.

,     (22)

,     (23)

2.3.3 HLL (Harten - Lax - van Leer) 흐름률 기법

HLL 수치흐름률 산정은 다음과 같다.

(24)

여기서 , 은 각각 요소 경계의 왼쪽과 오른쪽에서의 파의 전파속도로서 다음과 같이 정의한다.

 (25)

 (26)

 (27)

 (28)

HLL 흐름률 기법에서 요소의 경계 왼쪽(L) 하상, 또는 오른쪽(R) 하상이 마른(dry)하상일 경우 파의 전파속도는 각각 다음과 같다.

,    (29)

,    (30)

2.4 시간 적분

불연속 갤러킨 공간 차분의 정확도에 상응하는 정확도의 시간적분을 적용하기 위해 3차 TVD Runge-Kutta 시간 적분 기법을 적용하였다. 먼저, 시간적분을 위한 Eq. (11)은 다음과 같이 쓸 수 있다.

 (31)

여기서 은 질량행렬(mass matrix)이며 는 절점 값(nodal value)이다. Eq. (31)의 시간 적분에 3차 Runge-Kutta 기법을 적용하면 다음과 같다.

  (32)

 (33)

 (34)

또한, 불연속 갤러킨 기법에서 사용되는 양해법의 특성으로 인하여, 시간간격 는 다음과 같이 CFL 조건을 만족하도록 산정하였다(Cockburn, 1999).

 (35)

여기서 형상함수의 차수는 로서, 본 연구에서는 1차 다항식을 사용하였다.

2.5 기울기 제한자(slope limiter)

고차 정확도 수치기법인 불연속 갤러킨 유한요소법은 수심 및 유속 등의 기울기가 큰 불연속 흐름 주위에서 비현실적인 진동(oscillations)이 발생하기 때문에, TVD 특성을 유지하기 위해서는 이를 억제할 필요가 있다. 수치진동을 제거하는 기법으로 흐름률 제한자(flux limiter), 기울기 제한자(slope limiter)등이 있으며 본 논문에서는 기울기 제한자를 Runge-Kutta 기법의 각 단계에 적용하였다.

전체 1차원 요소중 번째 요소(Fig. 3)에서, 변수 에 대한 기울기 제한로서, 다음과 같은 MUSCL (Monotonic Upstream- centered Scheme for Conservation Laws) 제한자(LeVeque, 2002)가 적용된 바 있다.

 (36)

여기서, 는 요소내 변수의 평균값, 는 요소의 중앙 좌표, 와 는 각각 요소의 시작과 끝점의 좌표이다. minmod 함수는 다음과 같이 정의한다.

(37)

여기서 이다. 본 연구에서는 Cockburn (1999)의 일반화된 기울기 제한자(generalized slope limiter)를 면적 A와 유량 Q에 대하여 적용하였다. 순서는 먼저 진동(oscillation)이 관찰되는 요소에 대하여 다음과 같이

 (38)

를 계산하고, 이에 대하여 다시 MUSCL 기울기 제한자를 적용하는 방식으로 계산을 진행하였다.

3. 수치 모의

댐 붕괴류(dam-break flow)는 개수로(open channel)에서 발생하는 천이류(transcritical flow) 흐름의 대표적인 사례로 알려져 있다. 수로의 단면이 일정한 일양단면 하도(prismatic channel)에서는 임의의 단면에 대하여, 바닥 저면의 경사가 0이고, 마찰이 없는 경우, 특성곡선법에 의한 해석해(exact solution)가 존재한다(Jain, 2001). Fig. 4는 댐 붕괴류 모의를 위한 초기조건으로서, 두께가 0인 가상적인 댐으로부터 상류수심 또는 상류하도 단면적을 각각 , , 하류수심 또는 하류하도 단면적을 , 로 부여한다. 본 연구에서는 삼각형 단면 수로와 사각형 단면 수로에서의 댐 붕괴류 흐름을 모의하여, 이를 해석해와 비교하였다. 또한, 단면의 폭이 변화하는 사각형 단면 수로에서 댐 붕괴류의 예로서 Townson and Al-Salihi (1989)의 수리실험 데이터와 비교하였다.

3.1 일양단면 하도(prismatic channel)에서 댐 붕괴류 해석

3.1.1 삼각형 단면에서의 댐 붕괴류

하상 경사와 마찰이 0인 전체 [-500 m, 500 m]인 총 1,000 m인 삼각형 일양단면(prismatic) 수로에서 댐 붕괴류 흐름을 모의 하였다. 단면은 측면경사(side slope)가 1.0인 삼각형 단면이며 댐은 0 m에 위치하고 있다. 초기조건으로 상류하도 단면적은 1.0 m² (수심 1.0 m)이며 하류하도에는 젖은 하상조건 0.01 m² (수심 0.1 m)에 대하여 댐 붕괴류를 모의하였다. 5 m의 일정 크기 격자를 사용하였으며, 계산 시간간격 는 CFL 조건(Eq. (35))에 의한 상한값의 90%값을 사용하였다. Riemann 근사해법으로 Lax-Friedrichs 흐름률과 Roe, HLL 흐름률을 사용하였다. Fig. 5(a)는 댐 붕괴후 80 sec 일 때 수위를 모의한 것으로서, 세 흐름률 모두 해석해와 잘 일치하는 것을 볼 수 있으며, Fig. 5(b)는 유속을 비교한 것이다. 세 흐름률 모두 해석해와 잘 일치하며, 흐름률 사이의 차이는 없는 것으로 보인다. 다만, x = 200m 부근에서 차이(undershoot)을 보이고 있다. 이는 사각형 단면의 흐름에서는 보이지 않는 현상으로서(Lee and Lee, 2013), 향후 고차(higher-order) 근사 등 다른 방법과 비교 연구가 필요할 것으로 보인다.

3.1.2 사각형 단면에서의 댐 붕괴류

삼각형 수로에서와 같은 조건으로, 하상 경사와 마찰이 0인 전체길이 2,000 m [-1000 m, 1000 m]인 사각형 단면 수로에서의 댐 붕괴류 흐름을 모의하였다. 단면의 폭은 이며 댐은 0 m에 위치하고 있다. 초기조건으로 상류하도 단면적은 10.0 m² (수심 10 m)이며 하류하도에는 단면적 0.5 m² (수심 0.5 m)과 단면적 0.005 m² (수심 0.005 m)를 적용하였다. 5 m의 일정 크기 격자를 사용하였으며, 계산 시간간격 는 삼각형 수로와 같이 CFL 조건(Eq. (35))에 의한 상한값의 90%값을 사용하였다. Riemann 근사해법으로 Lax-Friedrichs 흐름률, Roe 흐름률과 HLL 흐름률을 사용하였다. Figs. 6 and 7은 하류 수심이 각각 0.5 m, 0.005 m 일 때의 모의 결과이다. 댐이 붕괴된 후 시간 50 sec 일 때, 하류 수심 0.5 m인 경우에는 역시 x = 200 m 부근에서 다소 수심을 과소평가하지만, 세 흐름률 모두 수위와 유속에서 해석해와 대체로 잘 일치하고 있다. 그러나, 수심 0.005 m 인 경우에는 잘 일치한다고 보기 어렵다. 보다 정확한 결과를 위해서는 고차 근사를 사용하거나, 마른 하도(dry bed)에 적합한 알고리즘이 필요할 것으로 보인다.

3.2비일양 단면(nonprismatic channel) 하도에서 댐 붕괴류 해석

Townson and Al-Salihi (1989)는 Fig. 8에 보인 바와 같이 사각형 단면의 수축 및 확장 수로에서의 댐 붕괴류 흐름을 모의하였다. 단면의 폭이 수축하여 평행하게 흐르는 수로(Fig. 8의 파선)와 평행한 단면을 유지하다가 확장하는 수로(Fig. 8의 실선)에서 실험을 수행하였으며, 상류경계에서 부터 0.5 m, 1.0 m, 1.8 m, 2.5 m, 3.0 m, 3.5 m에 위치한 관측점에서 수위를 측정하였다. 단면 폭의 수축 및 확장 각도는 5°이며 댐이 위치한 곳의 폭은 두 수로에서 모두 이다. 전체 수로의 길이는 4 m이며, 댐(x = 1.8 m)을 기준으로 왼쪽이 상류, 오른쪽이 하류이며, 상류부와 하류부의 길이는 각각 1.8 m, 2.2 m이다. Townson and Al-Salihi (1989)의 실험 조건에 따라, 초기 상류 수심()은 0.1 m를, 하류 수심()으로 마른 하상으로 10^-5 m, 젖은 하상으로는 0.0176 m를 수심으로 부여하였다. 수로는 경사가 없는 수평()이고, Townson and Al-Salihi (1989)는 수로의 마찰계수를 명시하지 않았으며, 본 연구에서는 Mohapatra and Bhallamudi (1996)가 재현한 바와 같이 마찰이 없는 수로(0)를 가정하였다. Riemann 근사해법으로 Lax-Friedrichs 흐름률, Roe 흐름률과 HLL 흐름률을 사용하였다. 모두 0.05 m의 격자를 사용하였으며, 시간 간격 는 CFL 조건(Eq. (35))에 의한 상한값의 90%값을 사용하였다.

3.2.1 사각형 수축(contracting) 단면에서의 댐 붕괴류

단면의 폭은 0 m에서 최대 폭을 가지며 5°의 각도로 흐름방향에 따라 댐(dam)이 위치한 1.8 m까지 폭이 축소하고 댐 이후로 폭이 0.1 m로 일정한 수로에서 댐 붕괴류 흐름을 모의하였다. Fig. 9(a)은 초기 상류수심 0.1 m, 마른 하상 조건의 하류수심으로 10^-5 m를 부여하여, 댐이 붕괴된 후 1.5 sec일 때의 수심을 보인 것이다. Fig. 9(b)은 초기 상류수심 0.1 m, 하류수심 0.0176 m 일 때, 댐 붕괴 후 2.1 sec에서 수심을 실험 측정 값과 비교한 결과이며, 대체로 잘 일치함을 알 수 있다. 사실상 본 연구에서 적용한 세가지 수치 흐름률 사이의 차이점은 크지 않다.

3.2.2 사각형 확대(expanding) 단면에서의 댐 붕괴류

Fig. 8에 보인 바와 같이 단면의 폭은 0 m에서 댐이 위치한 1.8 m까지 0.1 m로 일정하다. 댐 이후 4 m까지 5°의 각도로 확대되는 사각형 확장 단면 수로에서 댐 붕괴류 흐름을 모의하였다. Fig. 10(a)은 초기 상류수심 0.1 m, 마른 하상 조건을 만족하기 위해 하류수심 10^-5 m를 적용하여, 댐이 붕괴된 후 1.5 sec일 때의 해석결과이다. Fig. 10(b)은 초기 상류수심 0.1 m, 하류수심 0.0176 m를 적용하여, 댐이 붕괴된 후 2.1 sec일 때의 해석결과이다. 수치해석 결과를 Townson and Al-Salihi (1989)의 실험 데이터와 비교하여 잘 일치함을 확인하였다.

사각형 단면 축소 및 확대, 두 경우 모두에서 평행 수로와 접합부인 댐 지점에서 수치모형과 차이를 보이고 있다. 이는 연결부의 방사선상 흐름(radial flow)을 1차원 모의로 재현하기 어렵기 때문인 것으로 보인다. 향후 2차원 모형을 적용한다면 보다 실제 흐름과 근접한 결과를 얻을 것으로 기대된다.

3.3 지형 변화에 따른 1차원 천이류

폭이 1 m인 직사각형 단면으로 전체 구간 [0m,1000m]에 대하여, 바닥 저면의 하상 변화가 다음과 같은 경우에 대하여 흐름을 모의하였다.

 (39)

상류 유입 경계조건으로서 0 m에 수심, (0)10 m, 유량 (0)20 m³/s를, 유출 경계조건으로 (1000)7 m를 적용하였다. 각 요소의 크기는 25 m로 모두 균일하며, 시간 간격 는 CFL 조건(Eq. (35))에 의한 상한값의 90%값을 사용하여 정상상태(steady state)에 도달할 때까지 계산을 수행하였다. 각각의 수치적 흐름률로 Lax-Friedrichs, Roe, HLL 흐름률을 사용한 결과를 해석해, Meselhe et al. (1997)의 수치해와 비교하였다. Fig. 11는 수위/수심 변화를 보인 것으로서, 세 수치 흐름률에서 모두 Meselhe et al. (1997)의 해에 비하여 정확한 거동을 보여준다. 다만, 세 흐름률 사이의 우수성을 비교하기는 어렵다. Lai and Khan (2012)이 지적한 바와 같이 Roe 흐름률과 HLL 흐름률 사이는 물론, Lax-Friedrichs 흐름률에서도 뚜렷한 차이를 발견할 수 없다.

4. 결론 및 향후 연구

본 연구에서는 1차원 천수방정식을 해석하기 위한 수치기법으로 불연속 갤러킨 유한요소법을 적용하였고 천이흐름의 대표적인 사례인 댐 붕괴류와 저면 지형변화에 의한 천이류를 모의하였다. 수치기법을 검증하기 위해 해석해가 존재하는, 하상경사와 마찰, 단면 폭의 변화가 없는 삼각형 단면과 사각형 단면 수로에 대하여 댐 붕괴류 적용하였다. 또한, 해석해가 없는 수로에서의 댐 붕괴류를 모의하기 위해 Townson and Al-Salihi (1989)의 실험에서와 같이 단면 폭이 감소하는 사각형 수축 단면 수로와 단면 폭이 증가하는 사각형 확장 단면 수로의 댐 붕괴류를 모의하였으며, 생성항이 있는 경우의 사례로서 지형 변화가 있는 천이류에 적용하였다.

Townson and Al-Salihi (1989)의 수축/확대 수로에서는 하류수심을 달리하여 수치흐름을 모의하였고, 그 결과를 실험값과 비교하여, 평행수로와 접합부를 제외하면, 비교적 잘 일치함을 확인하였다. 수리실험 결과에서와 같이 단면이 축소하는 하도에서 댐 붕괴 후 깊은 수심이 형성되는 것을 알 수 있었다. 지형 변화가 있는 1차원 천이류 모의에서는 본 연구에서 개발한 기법이 천이류 흐름을 잘 모의함을 확인하였다. 다만, 본 연구의 결과로만 본 다면, 흐름률 사이의 결과의 차이는 크지 않은 것으로 보인다. 추후 단면 폭의 변화뿐만 아니라 하상고와 저면 마찰을 포함한 수로에서의 흐름에 대해서도 본 연구에서 개발된 기법을 적용할 예정이다.