Practical Numerical Model for Wave Propagation and Fluid-Structure Interaction in Infinite Fluid

Jeong-Rae Cho; Seong-Wook Han; Jin Ho Lee

doi:10.7734/COSEIK.2021.34.6.427

Preview

Research Paper

Journal of the Computational Structural Engineering Institute of Korea. 31 December 2021. 427-435
https://doi.org/10.7734/COSEIK.2021.34.6.427

Practical Numerical Model for Wave Propagation and Fluid-Structure Interaction in Infinite Fluid

무한 유체 영역에서의 파전파 해석 및 유체-구조물 상호작용 해석을 위한 실용적 수치 모형

Jeong-Rae Cho¹^*

Seong-Wook Han²

Jin Ho Lee³

조 정래¹^*

한 성욱²

이 진호³

¹Research Fellow, Department of Structural Engineering Research, Korea Institute of Civil Engineering and Building Technology, Goyang, 10223, Korea

²Post-Doctoral Researcher, Department of Structural Engineering Research, Korea Institute of Civil Engineering and Building Technology, Goyang, 10223, Korea

³Associate Professor, Department of Ocean Engineering, Pukyong National University, Busan, 48513, Korea

¹한국건설기술연구원 구조연구본부 연구위원

²한국건설기술연구원 구조연구본부 박사후연구원

³부경대학교 해양공학과 부교수

^{*Corresponding Author}

ABSTRACT

An analysis considering the fluid-structure interaction is required to strictly evaluate the seismic behavior of facilities such as, environmental facilities and dams, that store fluids. Specifically, in the case of an infinite domain in the upstream direction, such as a dam-reservoir system, this should be carefully considered. In this study, we proposed a practical numerical model for both wave propagation and fluid-structure interaction analyses of an infinite domain, for a system with a semi-infinite domain such as a dam-reservoir system. This method was applicable to the time domain, and enabled accurate boundary analysis. For an infinite fluid domain, a small number of mid-point integrated acoustic finite elements were applied instead of a general acoustic finite element, and a viscous boundary was imposed on the outermost boundary. The validity and accuracy of the proposed method were secured by comparing analytic solutions of a reservoir having infinite domain, with the parametric analysis results, for the number of elements and the size of the modeling region. Furthermore, the proposed method was compared with other fluid-structure interaction methods using additional mass.

Keywords

environmental facility

dam

infinite domain

fluid-structure interaction

wave propagation analysis

seismic response

환경시설물, 댐과 같은 유체를 저장하는 시설물을 대상으로 엄밀하게 지진 거동을 평가하기 위해서는 유체-구조물 상호작용을 고려한 해석이 필요하다. 특히, 댐-호소 계와 같이 상류 방향으로 무한 영역을 가지는 경우에는 이를 적절히 고려해야 할 필요가 있다. 본 연구에서는 댐-호소 계와 같은 반무한 유체 영역을 갖는 시스템을 대상으로 무한 영역의 파전파 해석 및 유체-구조물 상호작용 해석을 위한 실용적인 수치 모형을 제시하였다. 시간영역에 적용가능한 방법으로 정확하면서도 경계적인 해석이 가능하다. 무한 유체 영역에 대해서는 일반 acoustic finite element 대신 작은 개수의 mid-point integrated acoustic finite element를 적용하고 최종 경계에는 점성경계를 부과한다. 제안하는 방법의 유효성과 정확성을 검증하기 위해 강체 댐체를 가정한 반무한 호소계를 대상으로 적용하는 요소의 개수, 모델링 영역 크기 등을 매개변수로 해석해와 비교·검증하였다. 제안된 방법을 적용하여 댐-호소 계의 유체-구조물 상호작용을 부가질량을 사용하는 경우와 비교하였다.

키워드

환경 시설물

댐

무한 영역

유체-구조물 상호작용

파전파 해석

지진응답

MAIN

1. 서 론
2. 유체 및 유체-구조물 상호작용 정식화
3. 유체 무한 영역 처리를 위한 실용적 모델
4. 검증 및 적용
4.1 강체 댐체 가정시 지진응답
4.2 제안 기법의 유체-구조물 상호작용 해석 적용예
5. 결 론

1. 서 론

최근 경험한 두 개의 큰 지진인 경주 및 포항 지진은 더 이상 우리나라가 지진에 대한 안전지대가 아니라는 점을 시사하고 있다. 특히 2017년 11월 15일 포항지진(규모 5.4)에서는 양덕·유강·약성 등 3개 정수장, 장량·흥해 등 2개 하수처리장, 청하농공단지폐수처리장 등의 유체 저장 시설을 갖춘 환경시설물에 많은 피해가 발생하였다. 유체 저장 시설의 경우 유체를 합리적으로 고려한 해석법을 적용할 필요가 있다. 유사하게 댐과 같은 국가적인 중요 시설물의 경우 유체-구조물 상호작용을 고려한 해석이 필요하다.

지진 해석시 유한요소법을 적용하는 경우 유체에는 압력에 대한 파동방정식으로 유도한 acoustic finite element를 적용하고, 구조물에는 일반 고체 유한요소를 적용한다. 유체와 구조물의 경계에 적절한 연결요소로 적용하면 일반 상용 유한요소 프로그램에서 큰 어려움없이 유체-구조물 상호작용해석이 가능하다(ABAQUS, 2019; Cook et al., 2002). 하지만, 댐-호소 계와 같이 상류 방향으로 무한 영역을 가지는 경우에는 이를 적절히 고려해야 할 필요가 있다. 유체의 무한 영역은 (1) 시간 영역에서 충분히 큰 영역의 유체를 모델링하고, 경계에 감쇠기를 부과하는 방법(Cook et al., 2002), (2) 주파수 영역에서 해석해에 기반해 유도된 전달경계(transmiting boundary)를 적용하는 방법을 고려해 볼 수 있다(Cho, 1998). 첫 번째 방법은 시간영역해석이므로 댐체의 비선형성을 고려할 수 있다는 장점이 있으나 상당히 큰 영역의 유체를 모델링해야 하므로 계산 시간 등에서 비효율적이다. 두 번째 방법은 해석해에 기반하여 무한 영역을 처리하므로 정확한 방법이지만, 댐체의 비선형성을 고려할 수 없고, 일반 상용 유한요소 프로그램에서 사용할 수 없는 한계를 가진다.

시간영역에서 반무한 영역에서 발생하는 파전파 문제 및 유체-구조물 상호작용 문제의 해를 효율적으로 구하기 위해 고차 흡수경계조건(Givoli, 2004)과 perfectly matched layer(PML) (Basu and Chopra, 2004)이 개발되어 다양하게 사용되고 있다. 하지만, 고차 흡수경계조건과 PML을 구현하기 위해서는 추가 내부변수를 고려하기 때문에 구현이 어렵거나 변위 기반 유한요소기법 적용 시에는 변위를 시간에 대하여 적분한 항을 고려해야 하는 단점이 있다. 또한, 이 방법들의 수치 안정성에 대해서는 아직 명확한 결론이 나지 않았다. 이 연구에서는 무한 매질에서의 시간영역 파전파 해석 및 유체-구조물 상호작용 해석에 활용할 수 있는 정확하고 효율적이면서도 구현하기 쉬운 수치적 방법을 제안하고자 한다. 이를 위해 반무한 매질에 대한 mid-point integrated finite element(Guddati et al., 2016)와 점성경계에 근거한 방법을 제안할 것이다. Mid-point integrated finite element는 파가 전파하는 방향으로 선형의 형상함수를 사용하고 그 방향으로 mid-point integration을 사용하여 시스템 행렬을 계산하는 방법이다. 이 방법을 사용하면 통상적인 유한요소와 비슷한 형태의 시스템 행렬을 얻을 수 있어, 기존 해석코드에 사용자 요소로 쉽게 구현이 가능하다. 그러므로 제안하는 모델은 무한 매질에서의 파전파 문제에 대하여 아주 “실용적”인 접근법이라고 할 수 있다. 본 연구에서는 댐-호소 계와 같은 반무한 유체 영역을 갖는 시스템을 대상으로 무한 영역의 파전파 해석 및 유체-구조물 상호작용 해석을 위한 실용적인 수치 모형을 제시하고자 한다. 제안된 방법은 시간영역에 적용가능한 방법으로 정확하면서도 경제적인 해석이 가능하다. 무한 유체 영역에 대해서는 일반 acoustic finite element 대신 작은 개수의 mid-point integrated acoustic finite element를 적용하고 최종 경계에는 점성경계를 부과한다. 제안하는 방법의 유효성과 정확성을 검증하기 위해 적용하는 요소의 개수, 모델링 영역 크기 등을 매개변수로 해석해와 비교·검증하였다. 기존의 해석에서 유체의 영향을 단순히 부가질량으로 모델링하는 경우가 많은데, 제안된 방법을 적용하여 댐-호소 계의 유체-구조물 상호작용 해석을 수행하고 기존의 부가질량을 사용하는 경우와 비교하여 그 우수함을 보였다.

2. 유체 및 유체-구조물 상호작용 정식화

이 장에서는 Fig. 1과 같이 댐-호소 계를 시간 영역에서 충분히 큰 영역의 유체를 모델링하고, 경계에 감쇠기를 부과하는 통상적인 방법(Cook et al., 2002)의 이론을 소개한다. 유체는 비점성 압축성 이상유체로 가정한다. 이때 유체의 지배방정식은 다음과 같다.

(1)

\nabla^{2} p = \frac{1}{c^{2}} \dot{p}

여기서, $p (x, y, z, t)$ 는 유체의 동수압이고, $c$ 는 유체내 압축파의 속도로 $\sqrt{K} / ρ, K$ 는 체적탄성계수, $ρ$ 는 밀도이다. 물의 경우 $K = 2.1904 \times 10^{9} Pa, ρ = 1, 000 kg / m^{3}, c = 1, 480 m / s$ 이다. 경계조건으로는 다음과 같은 유형이 가능하다.

(2a)

p = p_{0}

(2b)

\frac{\partial p}{\partial n} = {(\frac{\partial p}{\partial n})}_{0}

(2c)

\frac{\partial p}{\partial n} = - \frac{1}{g} \ddot{p}

(2d)

\frac{\partial p}{\partial n} = - \frac{β}{c} \dot{p}

식 (2a), (2b)는 경계에서 필드 또는 그 플럭스에 대해 기지의 값을 지정하는 일반적인 경계조건이다. 식 (2c)는 유체 자유 표면에서의 표면 유동(sloshing)을 고려한 조건이다. 식 (2d)는 무한 영역의 비반사(non-reflecting) 경계를 표현한다. $n$ 은 유체 경계면의 외향 단위 법선 벡터, $g$ 는 중력가속도이고, $β$ 는 실험에서 결정되는 무차원 계수이며, 1이면 경계면에서 압축파에 대한 점성경계(viscous boundary) 역할을 한다(Cook et al., 2002). 한편 고체와의 경계면에서는 다음의 조건을 만족해야 한다.

(3)

\frac{\partial p}{\partial n} = - ρ {\ddot{u}}_{n}^{t}

여기서, ${\ddot{u}}_{n}^{t}$ 은 구조물의 운동으로 인한 유체 경계면의 전체 가속도(total acceleration)이다.

식 (1), (2)를 대상으로 동수압 $p (x, y, z, t) = Np$ 와 같이( $N$ 과 $p$ 는 각각 유체의 형상함수와 절점 압력) 이산화하는 유한요소 정식화를 수행하면 다음과 같은 유체 유한요소 모형의 지배방정식을 얻을 수 있다(Cook et al., 2002).

(4)

(M_{p} + M_{i}) \ddot{p} + C_{i} \dot{p} + K_{p} p = F_{p}

여기서, $M_{p}$ 와 $K_{p}$ 는 유체 영역내에서 발생하는 질량 및 강성행렬이고, $M_{i}$ 는 식 (2c)의 자유표면 출렁임 경계조건에 의해 발생하는 부가질량행렬, $C_{i}$ 는 식 (2d)의 비반사(non-reflecting) 경계조건에서 발생하는 감쇠행렬, $F_{p}$ 는 절점력이며 다음과 같다.

(5a)

M_{p} = \frac{1}{c^{2}} \int N^{T} N d V

(5b)

M_{i} = \frac{1}{g} \int N^{T} N d S

(5c)

C_{i} = \frac{β}{c} \int N^{T} N d S

(5d)

K_{p} = \int [\frac{\partial N^{T}}{\partial x} \frac{\partial N}{\partial x} + \frac{\partial N^{T}}{\partial y} \frac{\partial N}{\partial y} + \frac{\partial N^{T}}{\partial z} \frac{\partial N}{\partial z}] d V

(5e)

F_{p} = \int N^{T} \frac{\partial p}{\partial n} d V

위 식에서 $V$ 는 유한요소 영역에 대한 적분을, $S$ 는 경계면에 대한 적분을 의미한다.

한편 고체인 구조물에 대해 표준 유한요소모델을 적용하면 다음과 같다.

(6)

M {\ddot{u}}^{t} + C {\dot{u}}^{t} + {Ku}^{t} = F

여기에서 $M, C, K$ 는 구조물의 질량, 감쇠, 강성행렬이고, $u^{t}$ 는 절점의 전체 변위(total displacement), $F$ 는 절점 하중이다. 식 (4)와 (6)에 정의된 절점력을 외부에서 가력하는 하중벡터와 상호작용에 의한 하중벡터와 구분하여 표기하면 $F_{p} = F_{p}^{o} + F_{p}^{i}, F = F^{o} + F^{i}$ 이다. 위첨자 $o$ 는 외부 가력 하중을, $i$ 는 상호 작용력을 의미한다. 유체와 고체를 대상으로 서로에게 전달되는 상호작용력 $F_{p}^{i}, F^{i}$ 는 식 (3)을 도입하여 유도하면 다음과 같다.

(7)

F_{p}^{i} = - ρ S \ddot{u}

(8)

F_{}^{i} = S^{T} p

$S$ 는 유체와 구조물의 동적 거동을 연계시켜 주는 연결 요소(interface element)의 거동을 나타내는 행렬이다.

(9)

S = \int N^{T} n^{T} N_{s} d S

여기에서 $n = [n_{x} n_{y} n_{z}]^{T}$ 로 연결면에서의 유체영역의 외향 법선 벡터이고, $N_{s}$ 는 구조물의 형상함수이다. 상호작용력은 서로 평형을 이루므로 $F_{p}^{i} + F^{i} = 0$ 이다. 따라서 식 (4), (6)에 식 (7), (8)을 적용하면 유체-구조물 상호작용을 고려한 전체 운동방정식을 유도할 수 있다.

(10)

[\begin{matrix} M & 0 \\ ρ S & M_{p} + M_{i} \end{matrix}] \{\begin{matrix} {\ddot{u}}^{t} \\ \ddot{p} \end{matrix}\} + [\begin{matrix} C & 0 \\ 0 & C_{i} \end{matrix}] \{\begin{matrix} {\dot{u}}^{t} \\ \dot{p} \end{matrix}\} + [\begin{matrix} K & - S^{T} \\ 0 & K_{p} \end{matrix}] \{\begin{matrix} u^{t} \\ p \end{matrix}\} = \{\begin{matrix} F^{o} \\ F_{p}^{o} \end{matrix}\}

만약 외부 가력 하중이 없는( $F^{o} = 0, F_{p}^{o} = 0$ ) 지진 해석을 가정하고, 응답을 지반운동 ${\ddot{u}}_{g}$ 대한 상대운동으로 표현한다면 식 (10)은 다음과 같이 정리할 수 있다.

(11)

[\begin{matrix} M & 0 \\ ρ S & M_{p} + M_{i} \end{matrix}] \{\begin{matrix} \ddot{u} \\ \ddot{p} \end{matrix}\} + [\begin{matrix} C & 0 \\ 0 & C_{i} \end{matrix}] \{\begin{matrix} \dot{u} \\ \dot{p} \end{matrix}\} + [\begin{matrix} K & - S^{T} \\ 0 & K_{p} \end{matrix}] \{\begin{matrix} u \\ p \end{matrix}\} = \{\begin{matrix} - M \\ - ρ S \end{matrix}\} r {\ddot{u}}_{g}

식 (10), (11)에서 유체-구조물 상호작용계는 비대칭 시스템 행렬을 가지는 것을 알 수 있다. 한편, 구조물을 강체 고체를 가정하고 유체만을 모델링하여 지진해석하는 경우, 강체 고체로 취급하는 구조물과의 경계에서 지진하중이 작용하게 된다. 따라서 지배방정식은 다음과 같다.

(12)

[M_{p} + M_{i}] \ddot{p} + C_{i} \dot{p} + K_{p} p = - ρ Sr {\ddot{u}}_{g}

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F1.jpg

Fig. 1.

Fluid-structure interaction analysis of dam-reservoir system

3. 유체 무한 영역 처리를 위한 실용적 모델

반무한 유체 영역을 갖는 시스템을 대상으로 시간 영역에서 파전파 해석 및 유체-구조물 상호작용 해석을 수행하는 일반적인 방법은 시간 영역에서 충분히 큰 영역의 유체를 acoustic element로 모델링하고, 경계에는 식 (2d)의 비반사경계조건을 적용하는 방법(Cook et al., 2002)이다. 식 (2d)에 근거한 비반사 경계조건은 식 (5c)의 감쇠행렬로 모사되며, 결국 경계면에 분포하는 감쇠기(distributed damper)를 설치하는 것과 동등하므로 점성경계(viscous boundary)라고도 불린다. 이 방법은 많은 수의 요소가 필요함에 따라 경제적인 해석이 어려운 단점이 있다. 본 연구에서는 무한 유체 영역에 대해서는 일반 acoustic finite element 대신 작은 수의 mid-point integrated acoustic finite element (Guddati et al., 2016)를 적용하고 최종 경계에는 기존 방법과 같이 점성경계를 부과한다. Mid-point integrated acoustic finite element는 일반 acoustic finite element에서 파의 진행 방향으로 선형의 형상함수를 사용하고 1개 적분점을 적용한 것을 의미한다. 이 요소를 사용하면 파의 진행 방향으로의 분산 특성(dispersion property)을 Padé 근사를 이용하여 근사할 수 있다.

Fig. 2가 본 연구에서 제안한 실용적 모델을 도식화하여 나타낸 것이다. 유체 영역은 근역(near field)과 원역(far field)으로 구분한다. 근역은 일반 acoustic finite element를 적용하는데, 유체내 동수압의 관찰이 필요하거나 불규칙한 영역을 모델링할 때 적용한다. 원역은 파의 전파방향으로 적분점을 1개 적용한 mid-point integrated acoustic finite element를 적용한다. 그림과 같이 4절점 선형 요소라면 강성행렬을 적분할 때 일반 acoustic finite element는 2×2 적분을, mid-point integreated acoustic finite element는 1×2를 적용한다. 질량행렬에 대해서는 3×3과 1×3을 각각 적용하면 된다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F2.jpg

Fig. 2.

Fluid-structure interaction analysis of dam-reservoir system

4. 검증 및 적용

4.1 강체 댐체 가정시 지진응답

제안하는 방법의 유효성과 정확성을 검증하기 위해 적용하는 Fig. 3과 같이 높이 100m의 반무한 영역을 갖는 호소계를 대상으로 강체 댐체를 가정하여 해석해와 비교하였다. 수치 모델에서 AC2D4 요소는 4절점의 일반 acoustic finite element를 의미하며, 1:1의 형상비를 적용하였다. AC2D4Mx 요소는 4절점의 mid-point integrated element를 의미하며, 1:5의 형상비를 적용하였다. AC2D4 요소를 적용한 수치 모델이 기존 방법을 의미하며, AC2D4Mx 요소를 적용한 수치 모델이 제안 방법을 의미한다. 수치 모델의 지배방정식은 강체 댐체를 가정하므로 식 (12)에 해당한다. 정해와의 비교를 위해 자유수면에서 동수압을 0을 적용하다. 높이 방향으로는 10개 요소로 이산화하였으며 길이 $L$ 과 유체의 재료감쇠비 $ξ$ 를 매개변수로 정해와 비교한다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F3.jpg

Fig. 3.

Fluid-structure interaction analysis of dam-reservoir system

강체 댐체에 작용하는 반무한 유체의 동수압 및 그 합력에 대한 정해는 자유수면에서의 동수압을 0을 적용할 경우 주파수 영역에서 다음과 같다(Cho, 1998).

(13)

p (y, ω) = \frac{4 ρ {\ddot{u}}_{g}}{π} \sum_{n = 1}^{\infty} \frac{(- 1)^{n - 1}}{(2 n - 1) \sqrt{λ_{n}^{2} - k^{2}}} \cos λ_{n} y

(14)

F (ω) = \frac{8 H ρ {\ddot{u}}_{g} (ω)}{π^{2}} \sum_{n = 1}^{\infty} \frac{1}{(2 n - 1)^{2} \sqrt{λ_{n}^{2} - k^{2}}}

여기에서 $p (y, ω)$ 는 높이 $y$ 에서의 동수압이고, $F (ω)$ 는 동수압을 높이 방향으로 적분한 합력이다. $H$ 는 댐체의 높이, $k = ω / c$ 는 파수(wave number)이고, $λ_{n}$ 은 다음과 같다.

(15)

λ_{n} = \frac{(2 n - 1) π}{2 H}

한편, 고유진동수는 다음 식으로 계산할 수 있으며, 높이 100m인 경우 $f_{1} = 3.7 Hz, f_{2} = 11.1 Hz$ 이다.

(16)

f_{n} = \frac{(2 n - 1) c}{4 H}

해석은 2개의 수치모델(AC2D4 및 AC2D4Mx 요소 적용), 길이 30개(L=1,000, 2,000, ... 30,000m), 3종의 재료감쇠비 $ξ$ (0.5%, 1%, 2%)에 수행하였다. 수치모델에 사용된 입력지진은 Fig. 4에서 제시한 것과 같이 1989년 Loma Prieta 지진 기록 중 RSN765_LOMAP_G01090를 사용하였다(PEER Ground Motion Database, 2021). 유체의 재료 감쇠는 복소 영역의 정해에 적용할 때는 체적탄성계수 $K$ 대신 $K (1 + 2 ξ i)$ 를 적용하는 것으로 쉽게 고려가능하다. 수치모델에서는 다음과 같이 stiffness-proportional Rayleigh 감쇠로 근사적으로 고려할 수 있다(Chopra, 2017).

(17)

C^{s} = \frac{2 ξ}{ω_{0}} K^{p}

여기에서 $C^{s}$ 는 수치모델에 추가되는 재료 감쇠행렬, $ω_{0}$ 는 재료 감쇠비가 일치되는 각진동수로 여기에서는 유체의 1차고유진동 3.7Hz에 대응하는 각진동수를 적용하였다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F4.jpg

Fig. 4.

Input ground motion

Fig. 5와 Fig. 6은 감쇠비 0.01을 적용한 경우 합력의 시간이력을 해석해에서 구한 시간이력과 비교한 것이다. 해석해의 시간이력은 식 (14)에서 구한 주파수 응답을 푸리에 역변환을 통해 계산한 것이다. 그림의 범례에서 해석 경우를 model-L-n- $ξ$ 형태로 표시하였으며, model은 적용 요소, L은 길이, n은 길이방향 요소수이다. 예를 들어 AC2D4Mx-5000-50-0.01는 AC2D4Mx를 적용하고, L=5000m, L방향 50개의 요소, 재료감쇠비는 1%임을 의미한다. 그림에서 더 큰 영역을 모델링할수록(L이 클수록) 점차 정해와 근접함을 보여주고 있다. 또한 AC2D4Mx를 적용하는 경우 일반 유한요소인 AC2D4를 사용하는 경우와 비교할 때 같은 길이 L을 적용하지만 5배 작은 요소 개수를 적용하더라도 충분히 정해에 접근합을 보여준다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F5.jpg

Fig. 5.

Hydrodynamic resultant force – AC2D4 model

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F6.jpg

Fig. 6.

Hydrodynamic resultant force – AC2D4Mx model

Fig. 7과 Fig. 8은 감쇠비 0.01을 적용한 경우 합력의 전달함수를 정해와 비교한 것이다. 이 연구에서 동수압 및 그 합력의 전달함수는 다음과 같이 입력 지반 가속도의 비로서 정의하였다.

(18a)

T^{p} (y, ω) = p (y, ω) / {\ddot{u}}_{g} (ω)

(18b)

T^{F} (ω) = F (ω) / {\ddot{u}}_{g} (ω)

수치모델에 대한 시간이력결과를 계산한 후 주파수영역으로 푸리에변환하고, 그 결과를 입력 지반 가속도를 푸리에변환한 결과로 나누어 전달함수를 계산하였다. 그림에서 L이 작을 때는 해석해를 중심으로 변동하지만 L이 클수록 정해에 근접합을 알수 있다. 또한 AC2D4 대비 AC2D4Mx를 동일한 L을 적용하고 5배 작은 요소개수를 적용하여도 유사한 정확도의 결과를 보임을 알 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F7.jpg

Fig. 7.

Transfer function of resultant force – AC2D4 model

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F8.jpg

Fig. 8.

Transfer function of resultant force – AC2D4Mx model

수치모델 적용시 모델링 영역의 크기인 L에 따라 해석해에 수렴하는 정도를 정량적으로 표시하기 위해 주파수 영역에서 입력 지반 운동과 무관하게 정의할 수 있는 전달함수를 대상으로 다음과 같이 오차를 정의하였다.

(19)

e (ω_{1}, ω_{2}) = \frac{\int_{ω_{1}}^{ω_{2}} | | T^{F} (ω) | - | T_{0}^{F} (ω) | | d ω}{\int_{ω_{1}}^{ω_{2}} | T_{0}^{F} (ω) | d ω} \times 100 [%]

여기에서 아래첨자 0는 해석해를 의미하며, 오차는 관심있는 주파수 영역 $[ω_{1}, ω_{2}]$ 을 대상으로 해석해의 전달함수와 수치해의 전달함수의 차이를 면적비로 정의한 것이다. Fig. 9와 Fig. 10과 Table 1은 AC2D4를 적용한 경우와 AC2D4Mx를 적용한 경우 감쇠비와 길이 L에 따른 오차를 나타낸 것이다. 오차를 계산할 때 관심 영역은 0~25Hz로 설정하였다. 그림에서 유체의 재료감쇠비가 클수록 짧은 L에 대해 수치모델링을 하여도 빠르게 해석해에 접근함을 알수 있다. 예를 들어 약 2% 오차를 허용한다면 감쇠비가 0.5%인 경우 9000m를, 감쇠비 2%인 경우 3000m를 모델링할 필요가 있다. 한편 AC2D4Mx를 적용한 경우 AC2D4와 같은 길이 L을 5배 작은 수의 요소로 모델링하여도 유사한 오차를 보임을 알 수 있다. 이는 AC2D4Mx를 이용하여 보다 경제적인 해석이 가능함을 의미한다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F9.jpg

Fig. 9.

Error of transfer function – AC2D4 model

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F10.jpg

Fig. 10.

Error of transfer function – AC2D4Mx model

Table 1.

Error of transfer fuction of resultant force with therange of 0~25Hz [%]

L(m)	AC2D4			AC2D4Mx
	$ξ$			$ξ$
	0.5%	1%	2%	0.5%	1%	2%
1000	18.4	12.0	6.7	24.4	17.5	9.7
2000	11.7	6.3	2.9	17.8	10.1	4.3
3000	8.2	3.8	1.7	13.0	6.1	2.6
4000	6.1	2.6	1.3	9.8	4.0	1.9
5000	4.6	1.8	1.2	7.6	2.9	1.7
6000	3.7	1.4	1.1	5.9	2.3	1.6
7000	2.9	1.2	1.1	4.7	1.9	1.6
8000	2.5	1.0	1.1	3.9	1.6	1.6
9000	2.0	0.9	1.1	3.2	1.5	1.6
10000	1.8	0.9	1.1	2.8	1.4	1.6
11000	1.5	0.8	1.1	2.5	1.3	1.6
12000	1.3	0.8	1.1	2.2	1.3	1.6
13000	1.2	0.8	1.1	1.9	1.3	1.6
14000	1.1	0.8	1.1	1.8	1.3	1.6
15000	1.0	0.8	1.1	1.6	1.2	1.6
16000	0.9	0.8	1.1	1.5	1.3	1.6
17000	0.9	0.8	1.1	1.4	1.2	1.6
18000	0.8	0.8	1.1	1.4	1.3	1.6
19000	0.8	0.8	1.1	1.3	1.3	1.6
20000	0.8	0.8	1.1	1.3	1.3	1.6
21000	0.7	0.8	1.1	1.2	1.3	1.6
22000	0.7	0.8	1.1	1.2	1.3	1.6
23000	0.7	0.8	1.1	1.2	1.3	1.6
24000	0.7	0.8	1.1	1.2	1.3	1.6
25000	0.7	0.8	1.1	1.2	1.3	1.6
26000	0.7	0.8	1.1	1.2	1.3	1.6
27000	0.7	0.8	1.1	1.2	1.3	1.6
28000	0.7	0.8	1.1	1.2	1.3	1.6
29000	0.7	0.8	1.1	1.2	1.3	1.6
30000	0.7	0.8	1.1	1.2	1.3	1.6

4.2 제안 기법의 유체-구조물 상호작용 해석 적용예

제안된 방법을 적용하여 댐-호소 계의 유체-구조물 상호작용 해석을 수행하고, 그 결과를 단순히 부가질량을 사용하는 경우와 비교하였다. 해석 대상은 높이 120m, 폭 96m의 Pine Flat 댐으로 모델링 편의를 위해 유체와의 접촉면은 수직하게 가정하였다. 유체의 높이는 100m로 가정하였다. 댐체는 선형탄성을 가정하였고, 탄성계수는 34.45GPa, 밀도는 2,480kg/m³이다. 유체-구조물 상호작용을 적용하는 경우 Fig. 11과 같이 유체 무한 영역의 길이 L=5000m를 적용하였고. AC2D4Mx 요소로 길이 방향으로 100개의 요소를 사용였다(10m×50m 요소 적용). 자유수면에는 sloshing boundary를 적용하였다. 부가질량은 Westergaard(1933)가 제시한 다음 식을 적용하였으며, 댐체 만을 모델링한 후 유체의 접촉면에 질량을 추가한다.

(20)

m (y) = \frac{7}{8} ρ \sqrt{H (H - y)}

여기에서 $m (y)$ 는 높이 $y$ 에서의 댐체의 접촉면에 작용하는 단위면적당 부가질량이다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F11.jpg

Fig. 11.

Flud-structure interation of dam-reservoir system

Fig. 12는 각 경우에 대해 댐 마루(Fig. 11의 A점)에서의 변위시간이력을 나타낸 것이다. 이 연구에서 제시한 방법에 비해 부가질량을 적용한 경우 일부 오차가 발생함을 알 수 있으며, 이는 Chopra(1967)에서 지적한 바와 같이 Westergaard(1933)가 제시한 부가질량 모델은 유체의 압축성을 고려하지 않기 때문이다. Fig. 13은 댐 마루 변위에 대한 전달함수를 나타낸 것이다. 전달함수는 지반가속도가 아닌 지반변위에 대해 나타낸 것이다. 부가 질량을 사용하는 경우 첫 번째 지배주파수는 어느정도 정확하게 모사하지만, 이후 진동수 성분에서 오차가 발생함을 알 수 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F12.jpg

Fig. 12.

Displacement time history at dam crest

https://cdn.apub.kr/journalsite/sites/jcoseik/2021-034-06/N0040340610/images/Figure_jcoseik_34_06_10_F13.jpg

Fig. 13.

Transfer function of displacement at dam crest

5. 결 론

본 연구에서는 댐-호소 계와 같은 반무한 유체 영역을 갖는 시스템을 대상으로 무한 영역의 파전파 해석 및 유체-구조물 상호작용 해석을 위한 실용적인 수치 모형을 제시하였다. 무한 유체 영역에 대해서는 일반 acoustic finite element 대신 작은 개수의 mid-point integrated acoustic finite element를 적용하고 최종 경계에는 점성경계를 부과한다. 제안하는 방법의 유효성과 정확성을 검증하기 위해 적용하는 요소의 개수, 모델링 영역 크기 등을 매개변수로 해석해와 비교·검증하였다. 검증 결과 유체의 감쇠비에 따라 모델링 영역의 크기는 변경되지만 mid-point integrated acoustic finite element를 사용하면 일반 acoustic finite element를 적용하는 영역을 보다 경제적으로 계산할 수 있음을 확인하였다. 또한, 댐-호소 계의 유체-구조물 상호작용을 부가질량을 사용하는 경우와 비교하여, 기존의 유체-구조물 상호작용 해석에 사용되고 있는 부가질량 모델보다 더 정확한 결과를 얻을 수 있음을 확인하였다. 이 연구에서 제시된 mid-point integrated acoustic finite element와 점성경계를 사용하는 실용적인 수치 모형을 사용하면 무한 매질에서의 acoustic wave 전파 해석 및 유체-구조물 상호작용 해석을 시간영역에서 정확하게 수행할 수 있고, 향후 비선형 유체-구조물 상호작용 해석 등에 적용할 수 있을 것으로 기대된다.

Acknowledgements

본 결과물은 환경부의 재원으로 한국환경산업기술원의 환경시설 재난재해 대응기술개발사업의 지원을 받아 연구되었습니다(2019002850003).

References

ABAQUS (2019) ABAQUS Documentation, Dassault Systèmes, Providence, RI, USA.

Basu, U., Chopra, A.K. (2004) Perfectly Matched Layers for Transient Elastodynamics of Unbounded Domains, Int. J. Numer. Methods Eng., 59, pp.1039~1074. 10.1002/nme.896

Cho, J.-R. (1998) Seismic Response Analysis of Dam-Reservoir System Considering the Interaction between the Flexible Dam and the Compressible Impound Water, Master's Thesis, Seoul National University.

Chopra, A.K. (1967) Hydrodynamic Pressures on Dams During Earthquakes, J. Eng. Mech. Div., 93(6), pp.205~223. 10.1061/JMCEA3.0000915

Chopra, A.K. (2017) Dynamics of Structures: Theroy and Applications to Earthquake Engineering, 5th Edition, Pearson Education Inc.

Cook, R.D., Malkus, D.S., Plesha, M.E., Witt, R.J. (2002) Concepts and Applications of Finite Element Analysis, 4th Edition, John Wiley & Sons. Inc.

Givoli, D. (2004) High-Order Local Non-Reflecting Boundary Conditions: A Review, Wave Motion, 39, pp.319~326. 10.1016/j.wavemoti.2003.12.004

Guddati, M.N., Druskin, V., Astaneha, A.V. (2016) Exponential Convergence through Linear Finite Element Discretization of Stratified Subdomains, J. Comput. Phys., 322, pp.429~447. 10.1016/j.jcp.2016.06.045

PEER Ground Motion Database, accessed Nov 30, 2021, https://ngawest2.berkeley.edu/.

Westergaard, H.M. (1933) Water Pressures on Dams During Earthquakes, Trans. Am. Soc. Civil Eng., 98(2), pp.418~433. 10.1061/TACEAT.0004496

Journal of the Computational Structural Engineering Institute of Korea ISSN:1229-3059(Print) 2287-2302(Online) 한국전산구조공학회 논문집

Preview

Practical Numerical Model for Wave Propagation and Fluid-Structure Interaction in Infinite Fluid

ABSTRACT

MAIN

(1)

(2a)

(2b)

(2c)

(2d)

(3)

(4)

(5a)

(5b)

(5c)

(5d)

(5e)

(6)

(7)

(8)

(9)

(10)

(11)

(12)

Fig. 1.

Fluid-structure interaction analysis of dam-reservoir system

Fig. 2.

Fluid-structure interaction analysis of dam-reservoir system

Fig. 3.

Fluid-structure interaction analysis of dam-reservoir system

(13)

(14)

(15)

(16)

(17)

Fig. 4.

Input ground motion

Fig. 5.

Hydrodynamic resultant force – AC2D4 model

Fig. 6.

Hydrodynamic resultant force – AC2D4Mx model

(18a)

(18b)

Fig. 7.

Transfer function of resultant force – AC2D4 model

Fig. 8.

Transfer function of resultant force – AC2D4Mx model

(19)

Fig. 9.

Error of transfer function – AC2D4 model

Fig. 10.

Error of transfer function – AC2D4Mx model

Table 1.

Error of transfer fuction of resultant force with therange of 0~25Hz [%]

(20)

Fig. 11.

Flud-structure interation of dam-reservoir system

Fig. 12.

Displacement time history at dam crest

Fig. 13.

Transfer function of displacement at dam crest

Acknowledgements

References