1. 서 론
2. 유한요소 수식화
2.1 메탄 하이드레이트 해리
2.2 지배방정식(Governing equations)
2.3 주변수(main variables)
3. 메탄가스 생산 수치해석
4. 결 론
1. 서 론
전세계적으로 해양 퇴적지와 영구 동토층, 특히 울릉도 부근해역에 매장된 메탄가스 하이드레이트(methane gas hydrate)는 방대한 저장량으로 차세대 에너지원으로 주목 받고 있다. 이에 차세대 청정 에너지 자원으로 현재 미국, 일본 등 선진국에서 관련 연구가 활발히 진행하고 있다. 하지만 가스 생산 과정의 하이드레이트의 해리(dissociation)는 필연적으로 주위 지반의 변형을 유발하게 되며, 이로 인하여 대규모 해저사면파괴에 대한 우려가 끊임없이 제기되고 있다. 또한 가스 생산과정은 급격한 주위 환경변화를 초래할 뿐만 아니라, 메탄가스가 이산화탄소보다 100배 이상 온실효과를 유발한다는 점에서 가스의 대기유출은 심각한 지구온난화를 초래할수도 있다.
따라서, 안전하고 지속 가능한 가스 생산을 위해서는 생산 절차를 체적화하기 위한 실질적인 평가수단이 필요이다. 그리고 자연상태에서 매우 불균질하게 분포하는 메탄 하이드레이트를 생산하기 위해서는 각각의 지층에 가장 부합되는 경제적인 생산정의 설계와 생산절차를 제시하고 또한 장기적인 생산량을 예측할 수 있어야 한다. 아직까지 수치해석이 가스 잠재 생산량을 예측하고 지반안정성을 평가할 수 있는 유일한 수단이다 (Moridis et al., 2009). 메탄 하이드레이트의 해리 해석에 대한 신뢰성을 높이기 위해서는 지반의 Thermo-Hydro- Mechanical(열-수리-역학적, THM) 물성치에 대한 적절한 평가와 개발된 프로그램의 현장 또는 실내실험에 대한 검증이 매우 중요하다.
국외에서 가스 하이드레이트를 포함한 지반의 거동 해석을 위한 연구용 simulator로는 TOUGH+HYDRATE, MH21, STOMP-HYD, CMG-STARS 등이 있다. TH해석에 관한 이론해가 존재하지 않아 국제공동연구를 통하여 프로그램들 간의 상호 비교검증을 진행하고 있으나, 예측 값들에 상당한 차이를 보이고 있다(Fig. 1). 또한, 수치해석 결과에 대한 수렴성이 현저히 떨어지는 문제도 있다. 하이드레이트를 포함한 지층의 역학적 거동은 대부분 탄성으로 가정하고 수치해석을 수행하고 있으며, 주 변수들 사이의 결합된(coupled) 현상에 대하의 모사가 매우 미흡한 실정이다. 역학적 거동과 이와 결합된 현상에 대해서는 추가적인 연구가 필요한 실정이다.
메탄가스 생산에 대한 국내의 수치해석적 연구수준은 초보적인 단계이며, 하이드레이트를 포함한 지층의 역학적 거동 특성을 모사할 수 있는 모델은 개발되지 있지 않다. 다만, HM 현상을 해석하기 위한 상용 및 연구용 프로그램들이 사용되고 있으며, 최근에 방폐장의 완충재 및 동토의 거동해석을 위한 THM 해석 도구의 개발이 진행중이다(Shin et al., 2012). 그리고, 가스 하이드레이트 생성과 해리, 치환 등 메탄･이산화탄소 가스 하이드레이트에 대한 기초 연구는 상당히 진행되어 있으며, 가스 하이드레이트의 거동을 화학적, 열/수리학적 관점에서 분석･파악하기 위한 실험 및 수치해석적 연구가 각기 진행 중에 있다(Cha et al., 2008). Moridis et al.(2009)는 TOUGH+HYDRATE를 이용하여 동해 울릉유역에 잠재하는 하이드레이트 가스의 생산성을 규모, 유형에 따라 해석을 수행하였다.
본 연구에서는 메탄 하이드레이트의 해리 현상에 관한 fully coupled THM 유한요소 프로그램을 개발하고자 한다. 메탄 하이드레이트, 흙입자, 유체, 및 메탄가스의 질량보존의 법칙, 에너지 보존의 법칙, 그리고 힘평형 방정식으로부터 지배방정식을 유도하고, 지배방정식의 Galerkin 수식화와 시간적분을 통하여 Newton의 반복과정을 적용한다. 그리고 개발된 프로그램의 안정성 평가를 위하여 실내모형실험 결과에 대한 비교검증을 수행한다.
2. 유한요소 수식화
2.1 메탄 하이드레이트 해리
하이드레이트의 형성 조건에 적합한 낮은 온도와 높은 압력조건에서 생성되는 천연 메탄가스 하이드레이트는 가스 저류층(reservoir)에 많은 량의 메탄가스를 포함하고 있다. 물분자간의 수소결합은 공동을 형성하고, 유입된 메탄가스는 물과 가스 분자간의 van der Waals force을 통하여 공동을 안정화 시키게 된다(Sloan & Koh, 2008). 이러한 하이드레이트의 생성(formation)과 해리(dissociation)는 다음과 같은 화학 방정식으로 표현할 수 있다.
(1)
여기서 
는 가스 하이드레이트를 형성하기 위한 1개의 가스 분자당 필요한 물분자의 개수로 일반적으로 6.0이다.
하이드레이트와 얼음은 국부적인 온도와 압력조건에 따라 가스 생산과정에서 생성될수도 있고, 사라질수도 있다. 다공질 재료 공극내의 메탄 하이드레이트의 해리는 다양한 상(phase) 상호간에 질량과 열 교환, 그리고 운동학적 해리(kinetic hydrate dissociation)현상을 유발한다. 하이드레이트의 해리를 예측하는 모델은 Equilibrium model과 Kinetic model이 있다(Moridis et al., 2008; Sloan and Koh, 2008).
Equilibrium model에서는 얼음과 하이드레이트의 생성과 해리가 매우 빠른 속도로 진행된다고 가정한다. 주
어진 온도에서 국부적인 가스압력이 평형압 
보다 낮
으지면, 하이드레이트는 즉시 가스와 물로 해리된다고 보는 것이다. 따라서 이 model에서 메탄은 공극내에서 단일상(single phase)으로만 존재하게 된다.
Kinetic model에서 하이드레이트의 생성과 해리의 속
도는 가스의 도산성(fugacity) 
과 평형 상태 하이드레이트의 도산성 
의 차이에 비례한다고 가정한다(Gamwo & Liu, 2010). 즉, 
와 
의 차이가 하이드레이트 해리
의 driving force가 되며 속도를 결정한다(Fig. 2). Kim- Bishnoi의 제안식(Clarke & Bishnoi, 2001; Kim et al., 1987)은 hydrate 해리에 의한 가스의 몰 변화율을 다음과 같이 나타내고 있다.
(2)
여기서, 메탄 하이드레이트의 단위체적당 몰수는
, Kinetic constant 
, 반응에 참여하는 하이드레이트의 전체 면적비는 
이다. 메탄 하이드레이트의 몰 반응상수는 
이다(Moridis et al., 2008).
가스-하이드레이트-유체의 평형 도산성인 
는 메탄가스의 P-T diagram에서 평형압력 
과 근사하며, 가스의 도산성 
은 다공질 재료의 가스압 
으로 대신할 수 있다(Fig. 2, Yousif et al., 1991; Sun et al., 2005).
2.2 지배방정식(Governing equations)
지배방정식을 유도하기 위해 거시적 접근법(macro-scopic approach)을 사용하였으며, 다공질 재료는 국부적인 양을 체적평균으로 대체하는 등가 연속체로 가정하였다. 일반적인 다공질 재료는 기체(gas), 액체(liquid), 그리고 고체(solid)의 3상으로 이루어져 있으며, 흙입자(solid) 이 외의 공극은 물(water)과 공기(gas)로 채워져 있다. 그리고 흙입자는 고체(solid)상태로만 존재하며, 물(water)은 액체(liquid water) 상태로 존재하게 된다.
본 논문의 다공질 재료 공극내의 메탄은 기체 상태와 물에 용해된 상태, 그리고 하이드레이트의 고체상태로 존재할 수 있다. 기체상태와 액체상태는 이상기체 방정식과 Henry’s law를 이용하여 정의할 수 있다. 그리고 물의 동결에 의한 상변화를 고려할 수 있도록 하였다 (Fig. 3).
|
|
(a) | (b) |
Fig. 3. THM interaction in hydrate bearing sediment, (b) Phase diagram of hydrate bearing sediment | |
다공질 재료에서 흙입자를 제외한 공극을 차지하는 각각에 대한 포화도는 다음과 같다.
(3)
여기서, 공극내의 hydrate saturation은 
이고, 온도에 의존하는 ice saturation은 
이다.
고체상태인 흙입자와 얼음, 그리고 메탄 하이드레이트는 입자이동이 없는 것으로 가정하였다. 그리고 기체와 액체 상의 흐름은 darcy법칙과 상대투수곡선을 적용할 수 있으며, 간극을 차지하고 있는 얼음과 메탄 하이드레이트는 투수성을 감소시킨다고 가정하였다(Ruan et al., 2012).
위와 같은 가정들을 바탕으로 5가지의 지배방정식들과 일련의 구성관계식을 통하여 메탄을 포함한 다공질 재료의 THM 거동특성을 수식화할 수 있다.
다양한 종류의 물질이동에 대한 거시적인 평형방정식은 다음과 같은 일반화된 편미분 방정식으로 정리할 수 있다.
(4)
~ ~
여기서, 
는 다공질 재료의 종류를 나타내며, 
는 그
종류의 단위체적당 질량을 나타낸다. 그리고 
는 해당
~
종류의 질량 흐름(mass flux)으로 이류(advective)형태와 비-이류(non-advective)형태 2가지를 모두 포함한다.
a. Soil mass balance
고체상태로만 존재하는 흙입자에 대한 질량보존의 법칙은 다음과 같이 정리할 수 있다.
(5)
~ ~
여기서, 
는 흙입자의 단위 질량이며, 
는 흙의 간극
률(porosity), 그리고 
는 흙입자의 변위 속도이다.
~
b. Mass balance of methane hydrate (Solid)
고체상태인 메탄 하이드레이트에 대한 질량보존의 법
칙은 단위체적당 무게와 질량흐름 
, 
~ ~
을 이용하여 지배방정식을 유도하였다. 이때 source 항목에 추가되는 단위체적당 무게는 
이다.
(6)
~ ~
c. Water mass balance
다공성 재료내의 유체는 액체 상태의 물로 존재하며, 메탄 하이드레이트의 생성과 해리에 의해서 변화하게 된다. 단위체적당 유체의 무게와 질량흐름은 
,
이고, 하이드레이트의 생성과 해리에 의하
~ ~
여 source 항목에 추가되는 항목은 
이다.
~ ~
(7)
~ ~
여기서 단위면적당 유체의 유량 
은 유체압력(
)를 이용하여 다음과 같이 
표현할 수 있다.
d. Mass balance of methane (gas & liquid state)
다공질내 메탄의 기체 및 액체 상태의 질량에 대한 연속방정식은 다음같이 정리할 수 있다. 단위체적당 무게와 질
량흐름은 
, 
~ ~ ~
이며, source 항목은 하이드레이트의 생성과 소산에 의한 질량이다.
~ ~
(8)
~ ~ ~
여기서 단위면적당 기체의 유량 
은 가스압 (
)를 이용하여 다음과 같이 
표현할 수 있다.
e. Energy balance equation
에너지 보존의 법칙은 내부 에너지(internal energy)를 이용하여 정리하였으며, 다공질 재료내 전체 종(species)의 고체, 액체, 그리고 기체에 대하여 적용하였다. 단위체적당 무게와 질량흐름은 식 (13)과 같다.
(9)
~ ~ ~ ~ ~
여기서, 재료의 열전도에 의한 열흐름은 Fick의 법칙
(
)을 이용하였다.
~ ~
연속방정식(식 (4))에 적용하여 다음과 같은 에너지 보존에 관한 방정식을 유도하였다.
~ ~
(10)
~ ~ ~ ~
여기서, 
는 물의 단위질량당 잠재 융해열(specific latent heat of fusion)이며, 
는 methane의 단위질량당 잠재 융해열이다.
에너지방정식에 대한 수치해석은 다음의 잠재 융해열과 이류흐름에 의하여 쉽게 발산하는 경하는 있다.
1) 고체 재료의 열전도에 의한 에너지 연속방정식은 재료의
열전도도와 열용량에 의하여 정의되며(
),
~
기본 방정식은 토질역학의 탄성압밀해석과 동일하다. 따라서, Newton의 반복과정 없이 수렴하게 되며, 열적 물성의 범위가 제한적이므로 해석결과의 변화폭이 상대적으로 적은 편이다(예로, 압밀해석에서는 투수계수가 범위가 넓어 입력물성의 오차에 의한 해석결과의 차이가 크게 나타난다). 하지만, 다공질내 물질의 상변화에 의하여 잠재 융해열이 발생하게 되면, 열용량보다 큰 잠재열로 인하여 심각한
수치적 오류의 누적이 발생할 수 있다. 예로, 물의 잠재열과 열용량의 비는
으
로, 얼음을 녹이는 잠재열로 같은 양의 물을 80℃ 정도 데울 수 있다. 메탄의 해리 과정 또한 잠재열과 열용량비가 매우 높은 편으로(메탄의 해리에 필요한 잠재열=446.12kJ/kg), 충분히 작은 시간간격에서 Newton의 반복과정을 수행하지 않는 경우 수치해석 결과는 발산하기 쉽다.
2) TH coupled 문제에서 유속이 빠른 경우 유체의 흐름(advection)에 의한 열전달이 전도(conduction)에 의한 열전달보다 크게 되면 기존의 Galerkin수식화에 의한 결과는 발산하게 된다(Brooks & Hughes, 1982). 안정적인 해의 수렴을 위해서는 Petrov-Galerkin 수식화를 통하여 해의 안정성을 확보하여야 한다.
2.3 주변수(main variables)
기존의 하이드레이트를 포함한 다공질 재료의 거동을 해석하기 위한 국외의 전산 프로그램들(TOUGH+ HYDRATE, STOMP-HYD, CMG-STARS 등)은 서로 다른 주변수를 사용하고 있다. 또한 Fig. 2의 상(phase)의 영역에 따라 주변수를 바꾸는 PVSM(Primary Variable Switch Method) 기법을 사용하고 있다. TOUGH-HYDRATE는 gas(water)영역에서 가스압(Pg, gas pressure), 유체 포화도(Pl, liquid saturation), 온도(T)를 사용하고, hydrate(water)영역에서는 압력(P, pressure), 유체 포화도(Pl, liquid saturation), 온도(T)를 사용한다(Moridis et al., 2009). 그리고 STOMP-HYD는 gas(water)영역에서 가스압(Pg, gas pressure), 유체압(Pl, liquid saturation), 온도(T), 그리고 hydrate(water)영역에서는 hydrate 포화도(Sh, hydrate saturation), 유체압(Pl, liquid saturation), 온도(T)를 사용한다(White et al., 2011). 상변화 영역에서 이러한 주변수의 변화는 수치해석 프로그램의 수렴성을 크게 저하시키다(White et al. 2011). 그리고, 물의 동결점이하 영역에 대해서는 TOUGH-HYDRATE만이 얼음 포화도(Si, ice saturation)를 사용하고 있다고 명시하고 있다.
본 연구에서는 다양한 주변수들의 조합에 대한 수치적 안정성 검토를 수행하였다. 이를 통하여 주변수를 u, Sg, Pl, T, Sh으로 채택하여, Fig. 2의 상변화 전영역에서 해석이 가능한 상호의존적 지배방정식을 유도하였다. 상호의존적 지배방정식들을 임의의 시간 
과 
사이에서 사다리꼴 법칙으로 시간적분을 수행하고, Newton의 반복과정을 이용할 수 있도록 근사화시켰다. 시간 적분된 지배방정식의 수식화에 보간함수를 도입하여, 주 변수들에 대한 Galerkin 수식화를 수행하였다. 그리고 에너지 방정식의 convective term에 관해서는 Petrov-Galerkin 수식화를 통하여 해의 안정성을 확보하도록 하였다.
3. 메탄가스 생산 수치해석
기존에 수행된 실내모형 실험결과에 대한 수치해석을 통하여 개발된 프로그램의 비교 검증 해석을 수행하였다. 메탄 하이드레이트를 포함한 Berea sandstone 코아 샘플에 대하여 Masuda et al.(1999)는 하이드레이트의 해리 과정에서의 메탄가스 생산량과 온도변화를 계측하는 실험을 수행하였다. 해석 단면(4,000요소, 12,441절점)은 높이 30cm, 반경 5cm로 2차원 축대칭 요소를 사용하였다. 8절점 변위 요소(u)와 Sg, Pl, T, Sh는 4절점을 사용하였다.
코아의 초기 가스압-온도조건은 382.39t/m2-275.45 K이며, 공극내의 포화도는 Sh = 0.443, Sl = 0.206, Sg = 0.351이다. 가스의 생산은 감압법을 사용하였으며, 코아 출구부의 압력을 289.6t/m2 (△P=92.79t/m2)까지 감압한 후 일정하게 유지시켰다(Case A, Fig. 4). 감압에 의한 가스 생산과정에서 감압 출구(“a”), 중간높이(“b”), 그리고, 바닥(“c”)에서 주요 변수를 시간경과에 따라 모니터링하였다.
수치해석에 이용한 메탄 하이드레이트 해리모델은 kinetic 모델을 사용하였다. Equilibrium model의 경우 P-T 값이 phase boundary 하부로 진입하게 되면 공극내의 모든 하이드레이트가 순간적으로 해리하게 된다. 이로 인하여 에너지 방정식이 매우 불안정한 상태로 되어 시간간격을 최대한 작은 값을 사용해야 하며, 메쉬 의존적인 해석결과가 도출되어 본 논문에서는 비교를 수행하지 않았다.
Sandstone 코아의 감압실험에 대한 수치해석을 수행하고 그 결과를 Fig. 4에 도식화 하였다. 실내실험에서 출구부(“a”)의 가스 누적량은 200분 경과후에 0.405mol(CH4 6.48g)에 이르고, 수치해석 결과는 실험결과와 매우 유사한 가스 생산 경향을 보이면서 최종 생산량은 0.387mol(CH4 6.19g)을 예측하였다(Fig. 4b). 온도변화는 감압에 의하여 해리가 시작되는 초반부(감압시간 △t=16.7min)에 해리의 흡열반응(endothermic reaction)에 의하여 2℃ 이상 급격히 떨어지게 된다. 코아의 외부에서 대류(convection)에 의한 열전달이 이루어지면서 온도가 지속적으로 상승한 후, 250분 경과후에 외기 온도에 수렴하는 것을 알 수 있다(Fig. 4c). 출구부(“a”)에서 해리에 필요한 에너지는 외부와 인접 시료에서 공급되기 때문에 중앙부(“b”)의 온도가 더 낮아지며, 바닥부(“c”) 역시 외부에서 열원이 공급되므로 중앙부(“b”)보다 온도 감소가 작게 나타났다. Fig. 4d에서 하이드레이트 포화도는 감압에 의하여 해리되면서 출구부부터 순차적으로 포화도가 감소됨을 알 수 있다. 개발된 프로그램에 의한 수치결과는 실내 모형실험 결과와 비교적 유사하였으며, 해의 수렴성과 안정성을 확인할 수 있었다.
코아 실내실험에 대하여 감압속도와 감압량의 영향을 분석하기 위하여 감압속도를 10배 증가시킨 경우(Case B, 감압시간=10×△t)와 감압량을 20% 증가시킨 경우(Case C, 감압량=1.2×△P)에 대하여 수치해석을 수행하였다. Fig. 5는 감압속도와 감압량의 변화에 따른 주 변수들의 시간적 변화를 나타내고 있다. 동일한 감압량에 대하여 감압속도의 증가(Case B)는 가스의 생산량과 나머지 변수들에 미치는 영향이 미비한 것으로 나타났다. 다만, 감압량의 증가(Case C)는 가스 생산을 빠른 시간내에 완료시켰으며, 이는 core 내부의 온도변화 및 하이드레이트 포화도의 변화에 의해서도 확인할 수 있다. 이처럼 감압량의 증가에 의한 메탄 생산의 촉진은 높은 감압량으로 메탄가스의 해리경계에 먼저 도달하였을 뿐만 아니라, 내부 온도가 0℃ 이하로 내려가면서 얼음이 형성되면서 나오는 잠재열(exothermic reaction)에 의한 것으로 판단된다.
하지만 실제현장에서 감압속도나 감압량을 지나치게 크게 하면, 케이싱에 과다한 하중이 재하될 뿐만 아니라, 주위지반의 붕괴될 위험이 있으므로 효율적인 감압 계획을 수립하여야 한다. 다양한 감압 크기와 감압 시간에 따른 메탄 생산량과 주위 지반의 안정성 해석을 통하여 최적의 가스 생산계획을 수립해야 한다.
4. 결 론
방대한 저장량으로 차세대 에너지원으로 주목받고 있는 메탄가스 하이드레이트는 가스 생산과정에서 유발될 수 있는 문제를 최소화하면서 최적의 생산조건을 선정하기 위한 하이드레이트 포함한 다공질 재료의 THM 현상에 대한 프로그램의 개발이 절실하다. 국외의 해석 프로그램들은 국제공동연구를 통하여 프로그램들 간의 상호 비교검증을 진행하고 있으나, 예측값의 불일치와 수렴성에 문제가 있는 것으로 나타났다. 그리고 역학적 거동해석에 대한 추가적인 연구가 필요한 실정이다.
본 연구에서는 다공질 재료내의 메탄 하이드레이트의 해리 현상을 해석할 수 있는 fully coupled THM 유한요소 프로그램을 개발하였다. Methane hydrate, solid, water, 및 methane gas의 질량보존의 법칙, 에너지 보존의 법칙, 그리고 힘평형 방정식으로부터 지배방정식을 유도하였다. 다양한 주변수들의 조합을 통하여 주변수 u, Sg, Pl, T, Sh으로 채택하였으며, 상변화 전영역에서 해석이 가능하도록 하였다. 하이드레이트의 해리를 예측하는 모델은 kinetic model을 이용하였다. 개발된 THM 유한요소 프로그램을 이용한 메탄가스 생산 수치해석결과는 실내 모형실험 결과와 비교적 유사하였으며, 해의 수렴성과 안정성을 확인할 수 있었다. Masuda의 실내모형실험 조건하에서 감압속도의 증가는 가스의 생산 속도에 미치는 영향이 미비한 것으로 나타났다.
하이드레이트 포함한 지반의 THM 거동 실험과 다양한 모형생산실험을 통하여 THM 구성모델의 개선이 필요하며, 기존에 개발된 수치해석 프로그램들과의 상호 비교검증을 통한 기술력 공유 필요하다. 추후 울릉분지에 계획중인 시험시공의 안정성 예측-분석에 활용하고, 가스 하이드레이트 상업생산에 계측시스템과 연계하여 위험성 및 생산성 평가할 수 있는 도구로 활용할 수 있을 것으로 판단된다.
