1. 서 론
최근 해안 및 항만구조물은 파랑 및 표사제어와 같은 구조물의 1차적인 기능성에 더해 경제성, 경관성, 생산성, 환경성 등과 같은 다양한 추가적인 기능이 요구되고 있다(Lee et al., 2011). 항내의 정온도 확보를 위해 건설되는 신설 방파제의 경우 제체의 일부구간에 해수소통구를 설치하여 항내 수질개선을 도모하거나 케이슨의 구조를 다양화하여 파랑제어성능과 경관을 동시에 충족시키기 위한 시도들이 지속적으로 이루어지고 있는 것은 이와 같은 사회적인 요구를 반영하는 것으로 판단된다. 또한, 최근 신재생에너지(renewable energy)에 대한 시대적·사회적 관심이 급증함에 따라 연안구조물을 활용하여 파랑 및 조류에너지를 회수하기 위한 다양한 연구와 실증실험 등이 수행되고 있다(Hong et al., 2007; Park et al., 2014; Kim et al., 2016; Lee et al., 2018; Park et al., 2018; Cho and Kim, 2020). 현재 우리나라 연안역에서 시도되고 있는 파력발전의 경우 진동수주형과 같이 주로 중력식 구조물을 기초구조물로 이용하고 있지만, 연안역에 비해 파랑 에너지가 상대적으로 높은 대수심에 설치가 자유롭고 대단위 에너지단지를 용이하게 조성할 수 있는 부유식 구조물의 활용이 매우 기대된다. 더욱이 울산시에서 추진하고 있는 부유식 해상풍력단지 조성사업이 구체화 됨에 따라 향후에는 연안 및 항만에도 가동형 구조물 이용이 증가할 것으로 예상된다.
한편, 최근 계산기 성능의 비약적인 발전과 더불어 다양한 수치모델이 개발됨에 따라 수치해석이 해안 및 항만공학적인 문제를 해결함에 있어서 중요한 해석 툴로 인식되고 있다. 특히, 조파수조를 이용하는 수리모형실험과 유사한 조건을 수치적으로 구현하여 구조물과 파랑의 상호작용을 수치해석으로 모의할 수 있는 수치파동수로(Numerical Wave Tank) 모델은 수치해석의 이용을 가속화 시키고 있다.
수치파동수조는 주로 데카르트좌표계에 의한 유한차분법에 기초하고 있기 때문에 이동구조물의 해석에 있어서 질량보전이나 이동경계면에 대한 경계조건의 취급 등에 있어서 많은 문제점을 포함하고 있고 더욱이 곡면과 같은 복잡한 형상의 구조물로의 적용성이 현저하게 감소한다. 이와 같은 직각격자계의 단점을 보완하기 위해 Hirt and Sicilian(1985)가 제안한 격자 내에 구조물의 점유면적을 고려하여 이산화방정식을 구성할 수 있는 Fig. 1과 같은 FAVOR(Fractional Area Volume Obstacle Representation) 법의 적용이 증가하고 있다(Ataur Rahman et al., 2006; Kermanpur et al., 2008; Kuroiwa et al., 2010; Najafi-Jilani et al., 2014; Abbasi et al., 2018; Tajabadi et al., 2018; Ghanbari and Heidarnejad, 2020).
하지만 FAVOR 법은 직각좌표계를 기초로 이산화가 수행되므로 가동형 구조물에 대한 질량보전이나 이동경계면에 대한 경계조건의 취급 등에 있어서 많은 문제점을 포함하고 있고, 곡면과 같은 복잡한 형상의 구조물로의 적용성이 떨어지는 단점을 지닌다(Lee et al., 2011). 한편, 유한차분법과 대별되는 유한체적법의 경우 구조물의 형상표현이 비교적 자유로운 반면, 높은 격자생성비용과 자유수면의 추적 및 수치조파 등이 문제점으로 지적되고 있다. 이와 같은 전통적인 수치모델에 더해 최근 들어 기존의 데카르트좌표계상에서 임의 형상의 불투과 경계를 갖는 물체와 유체와 연성해석이 가능한 IB(Immersed Boundary) 법이 주목을 받고 있다. 해안공학 분야에 있어서도 Lee and Mizutani(2009)는 IB 법과 VOF 법을 결합한 새로운 수치모델을 수평원주 주변의 파동장에 적용하여 그 유효성을 검증하였다. Lee and Mizutani (2009)가 제안한 IB 법은 Fig. 2에 보인 바와 같이 구조물의 경계면에 위치한 Lagrangian point를 위치시켜 이로부터 구조물의 경계면에서 산정되는 유체력으로부터 수치적으로 경계조건을 만족시키는 계산기법이다.
Lee et al.(2011)과 Peng et al.(2013, 2015)은 Lee and Mizutani(2009)의 모델을 이동구조물로 확장하여 부유체의 파랑응답을 해석하였다. 하지만, Lee et al.(2011)이 제안한 IB법의 경우, 구조물 표면에서의 정식화 가정이 복잡하여 3차원 문제로의 확장이 곤란한 단점이 있다. 따라서, 본 연구에서는 기존의 FAVOR 법에 IB 법의 개념을 도입한 간략한 IB (Simplified IB, 이하 SIB) 법을 이용하여 새로운 수치파동수조모델을 제안한다.
2. SIB 법에 기초한 수치모델 구축
본 연구에서 새롭게 제안하는 SIB 법은 Lee and Mizutani (2009)에 의해 제안된 수치파동수조 개념을 이용한다. Lee and Mizutani(2009)에 의해 제한된 수치파동수조모델은 쇄파 같이 자유수면이 극심하게 변형되는 경계면을 추적하기 위한 VOF 법과 유체와 구조물의 상호작용((Fluid-Structure Interaction, 이하 FSI)을 모의할 수 있는 IB 법으로 구성된다. 그러나, 본 연구에서는 VOF법을 대신하여 이동경계면의 추적에 대한 알고리즘이 상대적으로 간단하고 경계면의 추적성능이 양호한 Takewaki et al.(1985)이 제안한 CIP(Constrained Interpolation Profile or Cubic Interpolated Pseudo-particle)법을 적용하였다. 이하에는 본 연구에서 제안하는 수치모델의 계산기법에 대하여 기술한다.
2.1 지배방정식
본 연구에서 제안하는 SIB법은 Lee and Mizutani(2009)와 유사하게 비압축성 점성유체를 대상으로 식(1)의 연속방정식과 식(2)의 운동방정식, 그리고 자유수면의 추적을 위해 도입하는 유체율 밀도함수 ϕ에 대한 식(3)의 이류방정식으로 구성된다.
여기서, xi 는 데카르트 좌표에서 i = 1, 2는 x, z에 대응하며, ui 는 i 방향의 유속성분(u, w), p는 압력, ρ는 유체의 밀도, ν는 동점성계수, Q는 소스항, t는 시간, g는 중력가속도, Dij는 변형속도텐서, γ는 계산영역의 측면경계에 위치한 감쇠영역에서 양의 값을 갖는 감쇠계수, Fi는 SIB 법에서 고려하는 계산영역내의 구조물에 의한 i 방향의 외력항이다. 식(2)의 운동방정식은 소스항 Q를 통해 목표하는 파랑과 흐름의 발생이 가능하지만 본 연구의 목적은 가동물체형 구조물 거동해석을 위한 SIB 법의 개발이므로 Q = 0를 적용하였다.
한편, Lee and Mizutani(2009)의 IB 법의 경우 액체만을 해석대상으로 하는 단상류의 모델인 반면, 본 연구에서는 자유수면 부근에 설치되는 구조물을 정도 높게 해석하고 자유수면을 포함하는 이동경계의 효과적인 추적을 위해 Lee et al.(2011)에 의해 제안된 혼상류 모델인 단일유체모델(one-field model for immiscible two-phase fluid)을 동시에 적용하였다. 따라서, 식(2)의 유체밀도 ρ ^ ν ^
여기서, 첨자 s, w 및 a는 각각 고체, 액체, 기체의 상태를 나타내며, 유체밀도함수 ϕi는 ϕ1 + ϕ2 + ϕ3 = 1의 관계를 갖는다. 식(1)~(3)의 지배방정식은 직교 교호격자를 적용한 유한차분법에 의해 이산방정식을 구성하였으며, 압력 연성에는 Amsden and Harlow(1970)에 의해 개발된 SMAC(simplified marker and cell) 법에 기초한 다단계(Fractional step) 법을 적용하였다.
2.2 CIP 법에 의한 수치확산의 최소화
이산화방정식은 유한한 격자점에서 구성되므로 식(6)과 같은 이류항의 계산 시 수치진동이나 확산이 발생하는 것은 잘 알려져 있다.
본 연구에서는 이류항에 대한 안정된 계산을 위해 Kim et al.(2020)의 연구에 기초하여 CIP 법을 사용하였다. Fig. 3(a)에 보인 바와 같이 각 격자점에서 초기값이 이류속도 u에 의해 ∆t시간 후 이류된 물리량의 분포형태는 Fig. 3(a)의 점선과 같은 원형 그대로 이류되어야 함에도 불구하고 Fig. 3(b)와 같이 제한된 격자점에서의 정보로 인해 이류대상 함수의 형태가 붕괴되는 것(그림 중의 파란색 실선)을 확인할 수 있다.
CIP 법에서는 이와 같은 물리량의 형태를 유지하면서 이류시키기 위해 식(7)과 같이 이류 물리량의 기울기(gradient)를 추가로 고려하여 이류시키는 방법으로 Fig. 3(c)와 같이 이류대상 물리량의 원형과 근사한 형태를 유지할 수 있다.
여기서, gx, gz는 각각 ∂f/∂x와 ∂f/∂z로 정의된다. 특히, CIP법은 인접한 격자점의 정보만으로 3차 정도의 계산이 가능하여 쌍곡선방정식의 수치해석에 있어서 수치확산이 작고 안정적인 알고리즘으로 알려져 있다(Kim et al., 2020). 인접한 격자점에 f와 g 값이 주어지면 이 격자점들 사이의 물리량 분포형태는 다음과 같은 3차 방정식으로 보간될 수 있다.
여기서, X = (x - xi), Z = (z - zk)이며, Cl, m은 다항식의 계수로 격자점 (i, k)와 인접한 격자점에서의 함수 f(x, z)와 함수의 공간기울기 g(x, z)로부터 결정할 수 있다.
2.3 구조물(고체상태)의 취급
Lee and Mizutani(2009)의 IB 법은 앞서 제시한 Fig. 2와 같이 Lagrangian point를 이용하여 유체영역을 점유하고 있는 구조물의 형상을 표현하고 구조물 경계에서 유체력을 평가하여 유체내의 구조물 영향을 Euler 격자에서 계산되는 유체장으로 피드백하는 과정을 통해 구조물과 유체의 상호작용을 모의하였다. 이와 같은 IB 법의 가장 큰 장점은 Lagrangian point를 통해 구조물의 경계를 표현하므로 임의 형상을 갖는 구조물과 유체와의 상호작용을 모의할 수 있다는 점이다. 하지만, 하지만 Lagrangian 격자와 Euler 격자간의 유체력 교환에 있어서 복잡한 보간이 필요하므로 3차원으로의 확장이 곤란하며 수치해석코드도 매우 복잡하다. 따라서, 본 연구에서는 Xiao et al.(1997)에 의해 제안된 방법인 고체상태를 강체구조물로 취급하는 방법을 적용하였다. 따라서, 물체의 운동은 이동과 회전의 조합으로 구성된다는 가정에 기초하여 해석을 수행한다.
우선 물체의 병진운동은 고체의 밀도함수 ϕ1와 고체 밀도 ρs를 적분하여 계산영역 내를 점유하고 있는 물질의 총질량과 고체의 중심위치를 산정하면 다음과 같이 고체의 중심위치에서의 속도( u i ¯ ) s
여기서, ( x i ¯ ) s
여기서, 고체에 작용하는 총 회전력(torque) Γ s t o t a l
2.4 SIB 법의 개요
본 연구에서 제안한 수치모델은 유체 중에 이동하는 고체입자 해석을 위해 Kashijima et al.(2001)에 의해 제안된 속도장에 대한 체적가중평균법을 Lee et al.(2011)에 의해 제안된 혼상류 모델과 결합하였다. 이하에는 본 연구에서 적용한 SIB 법을 간략히 기술한다. 기존의 IB 법은 유체입자의 운동방정식에 외력항을 도입하여 구조물의 경계면에서 산정되는 유체력으로부터 수치적으로 경계조건을 만족시키는 directingforcing 법이나 구조물 내부에 가상셀(ghost cell)을 위치시켜 보간을 통해 경계조건을 만족시키는 계산기법들이 대부분이다. 즉, IB 법은 이와 같이 구조물의 영향을 고려하기 위해 도입되는 외력항의 평가방법에 의해 구분되며 최근에는 구조물의 소성변형까지 고려할 수 있는 다양한 방법들이 제안되고 있다. 본 연구에서는 Kashijima et al.(2001)가 제안한 고체의 밀도함수로 속도장을 체적가중평균(Volume-weighted average)하는 식(13)을 적용하였다.
여기서, ui는 흐름장에서의 속도, uf는 국부유체속도, us는 국부고체속도이다.
식(13)에 의해 표현되는 흐름장의 속도 ui는 Navier-Stokes 운동방정식을 만족하는 것으로 가정하고 계산영역 내에 존재하는 고체영역에 대해서도 고밀도의 유체로 가정하여 식(2)의 운동방정식으로부터 다음식과 같이 시간발전에 대한 가유속 u ~ i
식(14)의 계산에서는 식(2)의 외력항인 Fi를 고려하지 않고 산정되며 계산된 가유속에 대해 Fi를 이용하여 구조물 내의 유속을 조정한다. 즉, ϕ1 = 1이 되는 구조물셀일 경우 F i = ( u s - u ~ i ) / ∆ t F i = ϕ i ( u s - u ~ i ) / ∆ t
식(15)에 보이는 바와 같이 본 연구에서 적용하고 있는 SIB법에 기초한 최종적인 수치결과는 FAVOR 법과 유사하게 격자내의 고체상태 점유율의 계산정도에 의존하는 것을 알 수 있다. 하지만, Lagrangian point와 Euler 격자를 혼용하는 Lee and Mizutani(2009)의 IB 법에 비해 단순한 알고리즘으로 인해 계산비용이 적고 3차원으로의 확장성이 뛰어난 장점을 지닌다.
3. 수치해석 결과
3.1 아르키메데스 모델
본 논문에서는 새롭게 제안된 SIB 법을 수치파동수조의 적용성을 검토하기 Fig. 2에 보인바와 같이 수면상에 놓인 각주의 액체혼입깊이에 대해 수치해석결과와 아르키메데스 원리에 따른 이론 결과를 비교하였다. 아르키메데스 원리의 이론결과는 개발된 수치해석 코드의 벤치마킹으로 자주 활용되는 문제 중의 하나이다.
Fig. 4에서 보인 기체를 포함한 계산영역은 해상도 ∆x = ∆z = 1의 단위 격자를 갖는 100 × 100의 격자시스템에서 해석되었다. 각주의 높이는 H = 10∆z와 H = 15∆z를 각각 적용하였으며, 각주의 초기위치는 각주의 하단부의 위치가 정수면과 동일하게 위치하여 각주의 초기속도는 0이다. 이와 같은 계산조건 하에서 각주의 질량을 변화시켜 최대혼입깊이를 이론치와 비교한 결과를 Fig. 5에 보인다.
Fig. 3에서 확인되는 바와 같이 수치계산결과는 각주의 밀도에 따른 액체혼입깊이에 대한 이론치의 전체적인 경향을 성공적으로 모의하고 있음을 알 수 있다. 각주의 밀도 변화에 따른 이론값과의 오차는 수치계산에서 각주의 이동에 따라 압력과 속도의 변화로 인한 흐름장의 변화가 발생하는 반면에 이론에서는 구조물 주변의 어떠한 흐름장도 고려되지 않는 이상적인 결과에 따른 것으로 판단된다.
3.2 자유수면으로 낙하는 구조물에 대한 추적성능
본 연구의 목적은 유체 중에 운동하는 가동형 구조물의 경계면 추적을 수치적으로 모의할 수 있는 수치모델을 구축하는 것이므로 밀도가 다른 혼상의 유체(기체와 액체)에서 가동형 구조물의 추적성능을 검토하기 위해 공기중에서 자유수면으로 자유낙하하는 물체(각주)에 대한 수치모의를 진행하였다. 자유낙하 실험은 투명 아크릴로 제작된 수조에 크기 8 × 7.2 c m, 질량 m = 68 g의 물체를 수면에서 5 cm 이격 된 높이에서 자유낙하 시켰으며, 자유낙하하는 과정은 초고속 카메라(MGL 10000, TSI Inc)로 촬영하였다. 자유낙하 실험 시 각주가 공기중에 회전하지 않고 수중에 연직으로 낙하하도록 반복 실험을 수행하였다.
한편, 수치계산에서는 각 방향의 격자해상도를 ∆x = 0.4 와 ∆z = 0.2를 적용하였으며, 본 연구에서는 이동하는 물체의 경계면 추적이 목적이므로 자유수면 부근으로 계산영역을 제한하였다. Fig. 6에서부터 Fig. 8에 자유낙하하는 각주의 실험 결과와 계산된 결과를 시간경과에 따라 비교하였다. Fig. 6은 자유낙하하는 물체가 수면에 부딪히기 직전까지의 결과로 본 연구에서 제안한 SIB 법은 공기 중에서 이동하는 경계면을 잘 추적하고 있음이 확인된다.
Fig. 7과 Fig. 8은 자유낙하하는 각주가 수면과 접촉하여 수면아래로 진입되는 과정을 보여준다. t = 0.1103 s 이후에 각주가 수면으로 진입되면서 낙하속도가 감소하게 되고, 이러한 각주의 속도 변화는 자유표면의 파동을 형성시키는 일련의 과정을 재현하고 있음을 확인할 수 있다. 하지만, 각주의 수면진입에 따라 발생하는 자유수면의 형성패턴은 실험결과와 다소 차이가 존재한다. 이는 수치계산의 경우 단면 2차원에 대한 계산결과로 수리모형실험의 3차원적인 효과가 반영되지 못한 결과에 기인하는 것으로 판단되지만 보다 정량적인 검토가 요구된다. Fig. 8에서는 각주가 수중으로 2/3 이상이 진입하면서부터 수면에서 발생한 파동이 계산영역의 측면으로 이동하는 현상이 관찰된다. 또한, Fig. 7과 유사하게 실험에서 관찰된 각주 주변의 수면파형을 과대평가하고 있음이 확인되지만 전체적으로 각주의 수면 진입에 따라 각주의 형상과 발생되는 수면파형을 양호하게 재현하고 있다고 판단된다. 이상으로부터 본 연구에서 제안한 SIB 법은 다상의 유체가 혼합된 흐름장에서 고체의 밀도함수로 속도장을 체적가중 평균(Volume-weighted average)함으로써 유체중에 놓인 가동형 물체의 해석이 가능함을 확인하였다.
5. 결 론
본 연구에서는 고정된 격자시스템에서 가동물체형 구조물 해석이 가능함과 더불어 3차원으로의 확장도 용이한 SIB (Simplified Immersed Boundary) 법을 제안하였다. 본 연구에서 제안한 SIB 법은 기존의 IB 법이 갖는 있는 복잡한 알고리즘을 개선하기 위해 각 상(phase)의 밀도함수가 국소질량의 중심과 함께 이동하는 것으로 가정한 단일유체모델(one-field model for immiscible two-phase fluid)을 기초로 고체의 밀도함수를 이용한 체적가중평균법을 적용하여 고체상태의 구조물을 취급하였다. 또한, 각 상의 이류계산 시 발생 가능한 수치확산을 최소화하기 위해 인접한 격자만으로 3차 정도의 수치계산이 가능한 CIP 법을 적용하였다. 제안된 SIB법의 해석성능을 검토하기 위해 자유수면으로 낙하하는 물체에 대한 모형실험과 수치모의를 동시에 수행하였다. 그 결과, 본 연구에서 제안한 SIB은 단순한 알고리즘에도 불구하고 자유수면으로 낙하하는 물체와 유체와의 상호작용을 양호하게 재현함을 확인하였다. 하지만, 본 연구에서는 파랑작용과 같이 수면 진동에 따른 구조물 응답은 재현되지 않았으며, 유체 중에서 고체상태의 회전운동에 대해서는 충분히 검토되지 못한 한계점이 있다. 향후 파랑작용 중에 위치한 구조물의 회전운동을 포함한 응답성능에 대한 추가적인 검토를 진행할 계획에 있다.
PDF Links
PubReader
ePub Link
Full text via DOI
Download Citation 
Print 
