1. 서 론

최근 세계 각국이 기후변화와 지구온난화 위기에 대응하기 위한 탄소중립을 도모하고 있는 가운데 유럽, 일본, 호주 등의 친환경 선도국에서는 기존 에너지 인프라를 친환경에너지 시스템으로 전환하기 위한 다양한 수소경제 활성화 정책들을 발표하고 있다(FCH JU, 2019; CSIRO, 2018). 일례로, 유럽연합(EU)은 2019년 수소 로드맵인 Fuel Cells and Hydrogen Joint Undertaking의 일환으로 2050년까지 수소 및 연료전지의 보급을 확대하는 계획을 발표하였고, 호주에서도 호주 수소전략(COAG Energy Counsil, 2019)과 같은 정책을 통해 그린수소 생산 및 공급인프라를 구축하여 수소경제를 활성화하고 있다(Choi et al., 2022). 이러한 연구 동향 가운데 우리나라도 수소경제 활성화 정책을 추진하는 등 수소경제로의 전환을 준비하고 있으며, 인구밀도가 높은 도심지에서의 수소 연료의 원활한 공급을 위해 중·대규모의 안정적인 지하 수소 저장 방안을 마련하고 있다(Go et al., 2024). 수소 연료의 대중화를 위해서는 수소에너지의 생산 및 충전 기술 못지않게 저장 기술개발이 필수적으로 요구되며, 수소의 생산 및 공급 확대는 필연적으로 중·대규모의 안정적인 수소 저장 기술을 요구하게 될 것이기 때문이다. 즉 수소저장기술은 산업 원료, 교통수단의 연료 및 건물 에너지 공급체로서 수소의 잠재력을 최대한 보장하고 생산과 소비 간의 격차를 메우는 데 필수적인 기술이라 할 수 있다. 수소 저장기술은 다양한 방식이 고려되고 있지만, Material 방식과 Chemical 방식은 효율이 매우 낮아서 공학적 목적에 부합하지 않다고 알려져 있다. Physical 방식으로 중·대규모의 수소를 현행과 같이 지상에 저장하는 경우 각종 재난재해 및 폭발 사고에 의한 안전상의 우려가 존재할 뿐만 아니라 부지선정에 대한 어려움이 있어 최근에는 세계 각국에서는 Geological 방식, 즉 수소의 지중저장 방식에 주목하고 있다(Go et al., 2022). 전 세계적으로 지중수소저장과 관련된 사업이나 선행연구들은 주로 지하 내 자연적으로 형성된 암염층을 이용해 조성된 인공공동을 활용하는 암염돔(salt cavern) 지중저장방식에 대하여 수행되고 있다(Simon et al., 2015; Heinemann et al., 2018). 공동을 둘러싸고 있는 암염층은 화강암이나 점토성 암반보다 상대적으로 낮은 투과성을 가지고 있어서 가스 유출 가능성이 낮고 높은 압력에 대한 구조적 안정성이 높다고 알려져 있다. 하지만, 암염돔을 활용한 수소 저장시설은 극히 일부 국가에서만 운영 중이며, 우리나라와 같이 암염돔을 활용할 수 있는 지질학적 조건을 지니고 있지 못한 국가들이 많다(Park et al., 2022). 따라서 국내의 지질학적인 여건과 상황에 맞는 수소 지하저장시스템 구축이 필요한 실정이라 할 수 있다.

이에 국내에서는 한국지질자원연구원이 수소의 지하저장시스템 기술 개발을 목표로 중대용량 액화수소 지하저장시스템 구축 핵심 설계기술개발을 연구 목표로 설정하였으며, 관련된 연구를 수행하고 있다(Lee et al., 2025). Fig. 1은 한국지질자원연구원에서 제안한 액화수소 지하저장시스템의 개념도를 나타내며, 이에 대한 자세한 설명은 2.3절에서 후술한다(KIGAM, 2025).

Fig. 1

Development of core design technologies for a large-capacity underground liquid hydrogen storage system conducted by the Korea Institute of Geoscience and Mineral Resources (KIGAM, 2025)

액화수소 지하저장시스템 구축에 있어서 장기적 저온 환경에 노출되는 저장공동 인근 암반에서는 공극 또는 불연속면에서의 동결과 유체 이동에 의한 불연속면 팽창 이슈로 인해 액화수소 저장공동시스템을 구성하는 여러 인공 구조물의 안정성에 큰 영향을 미칠 것으로 예상된다. 이 경우 수치해석을 이용한 저온 암반 복합 거동 분석기법은 암반의 장기적 동결 팽창 특성을 이해하는 데 매우 효과적일 것으로 판단된다. 또한, 향후 수소 지하저장시스템 사용의 대중화를 위해서는 다양한 시나리오에 따른 액화수소 지하저장시스템 안정성 평가에 대한 기준자료(reference data)가 요구된다. 따라서, 저장소 인근 저온 암반의 동결 팽창에 따른 액화수소공동의 역학적 안정성을 평가하고 다양한 불연속면 조건에 따른 저장공동의 피해를 저감시킬 수 있는 최적 설계연구가 필요하다. 이에 본 연구에서는 액화수소 지하저장소 인근 저온 암반에서의 공극수 상변화를 고려한 복합해석(열-수리 완전 연계)기법 개발하고 지하수 흐름 조건에서의 저온 암반의 냉각 거동을 분석하고자 한다.

2. 수치해석

2.1 수학적 모형

본 연구에서는 액화수소 지하저장소 인근 암반에서의 극저온 동결에 따른 온도변화를 평가하기 위해 다음과 같은 사항을 가정하였다. 암반은 미세균열만 있고 불연속면이 없는 다공질 건전암(porous intact rock)으로 간주하고 미세공극은 지하수로 완전히 포화된 상태로 가정하였다. 이에 다공질암은 암석입자, 미세공극을 채우는 비동결 공극수와 동결공극수를 포함하며 이들의 체적분율에 따른 열물성 변화와 동결 시 공극수 상변화에 따른 잠열 영향 등이 고려되어야 한다. 국부적인 온도평형상태를 가정할 때 동결에 따른 다공질암의 온도변화는 다음과 같은 지배방정식(governing equation)에 의해 모사될 수 있다(Michalowski and Zhu, 2006; Grenier et al., 2018).

여기서, n은 다공질암의 공극률(-), T는 다공질암의 온도(℃)이며, t는 시간, ρ_w는 비동결 공극수의 밀도(kg/m3), ρ_i는 동결공극수의 밀도(kg/m³), ρ_s는 암석입자의 밀도(kg/m³), C_w는 비동결 공극수의 비열(J/(kg·K)), C_i는 동결 공극수의 비열(J/(kg·K)), C_s는 암석 입자의 비열(J/(kg·K)), S_w는 미동결함수(Unfrozen water content curve), L은 비동결 공극수의 단위질량당 잠재 융해열(kJ/kg)을 의미한다. 한편, 미세균열 내 공극수의 동결 과정에서 공극 내의 일부 수분이 동결되지 않은 상태로 존재하는 것을 부동수분이라 하는데(Wang et al., 2017) 저온 암반의 온도해석에 있어서 중요한 인자로 작용한다(Kim and Go, 2023). 본 연구에서는 대표적인 부동수분곡선 모형인 Mckenzie et al.(2007)의 모델을 수치해석에 적용하였다(식 (2) 참조).

여기서, S_w,res는 잔류포화도, T는 다공질암의 온도(K), T_f는 어는점(K), W는 상변화전이온도구간(K)을 의미한다. 한편, 식 (1)의 λ_eff는 다공질암의 유효 열전도도로서 각 구성요소 j의 열전도도 λ_j를 해당 체적율 θ_j로 가중한 기하평균(geometric mean) 함수로 적용된다.

여기서, θ는 체적분율(volume fraction)을 나타낸다. 또한, 구성요소 j로 표현되는 s, w, i는 각각 암석입자, 미세균열 내 비동결 공극수, 동결 공극수를 의미한다.

2.2 모델 검증

본 연구에서는 구축된 해석모델의 예측 신뢰성 검증을 위해 다공질 재료를 활용한 동결토 토조실험을 수행한 후 이를 해석모델의 결과와 비교하였다. 사용된 동결토조시스템은 Fig. 2와 같이 토조시스템, 동결관과 연결된 냉각칠러, 항온수조와 연결된 칠러, 데이터로거 및 모니터 등으로 구성되어 있다. 다공질암석을 토조 안으로 삽입하기 어려운 한계가 있어 다공질암을 모래(주문진 표준사)로 대체하여 토조를 채운 후 동일 조건에서의 다공질재료의 동결거동을 실험과 해석에서 비교 검증하고자 하였다. 또한, 토조 외벽은 아크릴로 제작하여 시험 도중의 동결벽체의 형상을 관측할 수 있도록 하였다. Fig. 3은 모형 토조동결실험을 모사하기 위한 수치해석의 경계조건을 보여준다. 먼저 수리학적 경계조건의 경우, 실험에서와 동일하게 해석 도메인의 하단에서 상단으로의 지하수 흐름이 발생하도록 토조 상·하단에 △H만큼의 수두차를 부여하였다. 한편, 열적 경계조건으로는 하단에서 상단으로 물의 흐름에 의한 이류(heat advection)가 발생하도록 설정하되 토조 하단 경계에 초기 지반온도(T_ini)를 Dirichlet 경계조건으로 부여하여 지하수 유동이 존재하는 지반의 실질적인 열특성을 모사하였다. Fig. 4에서 알 수 있듯이 동일 동결실험 종료 시점(3시간)에서의 실험과 해석결과에서 동결벽체의 모양 및 형성 패턴이 유사하게 나타나는 것을 확인하였고, 이는 구축된 해석모델이 다공질 재료의 동결 거동에 대해 신뢰성 있는 예측 정확도를 지니고 있음을 확인하였다.

Fig. 2

Artificial freezing chiller system for porous media

Fig. 3

Boundary conditions of finite element model for the indoor freezing test. (a) Hydraulic BCs; (b) Thermal BCs

Fig. 4

Comparison results of frozen zone of porous material after 3 hours. (a) frozen zone in experiment, (b) frozen zone in FEM results, (c) comparison results of frozen zone between experiment and FEM results

2.3 내조시스템 개념모델

2.3.1 내조시스템 개념

본 연구에서 고려하고 있는 액화수소 지하저장시스템은 국내에서 채택가능한 유일한 선택지인 지하공동(rock cavern) 방식이며, 이에 대한 배경과 상세설명은 Lee et al.(2025)에 제시되어 있다. 해당 저장방식은 복공식 암반공동(lined rock cavern, LRC)으로 불린다(Fig. 1). 지하공동을 활용한 저장시스템의 경우 기밀성 및 단열을 위하여 공동 내에 단열재, 멤브레인, 콘크리트 라이너 등을 구성요소로 하는 내조시스템을 적용한다. 이에 따라 각 구성 재료의 두께, 재질, 재료적 특성, 암반 공동의 형상, 암반의 심도 등에 대한 설계가 필수적이며, 설계 과정에서는 효율과 안전성, 경제성을 종합적으로 고려해야한다(Park et al., 2010). 내조시스템에 주입되는 액화수소의 기화온도가 영하 253℃이기 때문에 암반공동 주변은 극저온 환경에 노출된다. 특히 지하수가 존재하는 지하암반 조건에서는 온도가 어는점 이하로 저하될 경우 주변 암반 내에서 얼음이 형성될 수 있다. 물이 얼음으로 상변화하는 과정에서 부피 팽창이 발생하며, 이는 내조시스템의 역학적 안정성에 영향을 미칠 수 있다. 반면, 형성된 얼음은 암반의 균열이나 공극을 메워 차수 역할을 할 수 있으며, 단열 역할을 할 수 있다. 이와 같이 얼음이 형성된 영역을 ice ring이라 하며, ice ring의 형성/영향 범위 평가는 내조시스템의 효율과 안정성을 결정하는 중요한 요소 중 하나이다(Jung et al., 2006).

2.3.2 수치해석 모델 설정

본 연구에서는 해석의 간편화를 위해 복공식 암반공동의 정중앙을 기준으로 2차원 반단면만 해석영역으로 설정하였고(2D axial symmetric), 저장소는 길이 2.7 m(L/2), 높이 9 m(H) 규격의 마제형으로 설정하였다. 또한, 콘크리트 라이너의 두께는 0.35 m, 절연체 두께는 0.2 m, 스틸 라이너 두께는 0.15 m로 설정하였고, 암반은 불연속면(절리면)이 없는 포화된 다공질암으로 간주하였다. 현재 전 세계적으로 액화수소 저장을 위한 복공식 지하저장시스템에 대한 구체적인 설계안이 제안된 바가 없으며, 이에 따라 본 연구에서는 한국지질자원연구원에서 목표로 하는 원천기술 기반 설계안을 적용하였다. 다공질암 및 LRC 시스템 구성요소에 대한 재료물성은 Table 1에 제시하였다. 한편, 지표(y = 0 m) 경계에서의 온도는 10℃로 일정한 조건으로 가정하였고(T_surface = 10℃), 암반의 초기온도는 0.03℃/m의 온도경사로 선형적으로 증가하는 형태로 부여하였다(T_ini = T_surface-0.03 [℃/m] × y). 또한, 개념모형에서, 콘크리트 라이너, 절연체, 스틸 라이너의 두께는 각각 0.3 m, 0.2 m, 0.15 m로 간주하였고, 액화수소와 스틸 라이더가 접촉하는 경계면에 –253℃의 온도를 Dirichlet 경계조건으로 부여하였다.

Fig. 5

Model domain and specification of LRC system used in the simulation

3. 결과 및 고찰

3.1 LRC 모델의 예비 해석 결과

Fig. 6은 액화수소용 복공식 저장소 인근 암반에서의 시간에 따른 온도분포 양상을 보여준다. 온도분포는 COMSOL에 내장된 Filter 기능을 사용하여 영하의 조건에서만 plot이 되도록 설정하여 ice ring의 분포를 시각적으로 도시하였다. 액화수소 저장 후 처음 5년까지는 LRC 인근 암반에는 ice ring이 전혀 만들어지지 않았고, 50년이 지나야 최대 5 m 두께를 지니는 방사형의 ice ring이 형성됨을 확인하였다. 액화수소가 매우 낮은 기화 온도를 지니고 있음에도 불구하고 절연체를 포함한 LRC 구조체의 다양한 구조물의 영향으로 인해 저장소 인근 암반에 적절한 ice ring이 만들어지는 데는 50년 이상의 시간이 소요된다는 것을 예측할 수 있다.

Fig. 6

Temperature distribution over time in the rock mass near LRC for liquefied hydrogen: the red line represents the 0°C isotherm

Fig. 7은 복공식 암반공동의 콘크리트 라이닝 최외각 경계로부터 인근 암반의 수평 거리에 따른 암반 온도 분포 양상을 보여준다. 저장 후 1년이 지나면 라이닝 최외각 경계에서 영하의 온도에 진입하며 50년 후에는 –7℃까지 경계온도가 냉각됨을 알 수 있다. 또한, 0℃ 라인을 기준으로하여 ice ring의 최대 두께를 산정하였는데, 저장 후 30년 후에 약 2.5 m, 50년 후에 약 5 m 정도 두께가 형성되었다.

Fig. 7

Temperature change along the ML1 (horizontal distance from the outer boundary of the LRC, see Fig. 5)

3.2 LRC의 다양한 설계조건에 따른 ice ring 분포 양상

3.2.1 LRC 절연체의 단열성능이 ice ring 반경에 미치는 영향

Fig. 8은 복공식 암반공동 내 절연체의 단열성능에 따른 50년 후의 ice ring의 분포 변화를 보여준다. 3.1절의 예비해석에서 가정된 절연체의 열전도도는 0.0184 W/(mK)이며, 이는 기존 연구를 참고하여 유사한 개념에서 사용되는 절연체의 열전도도 값으로 설정하였다(Jung et al., 2011; KIGAM, 2003). 본 연구에서는 이 값을 기준 값으로 설정하고, 기준 값 대비 80% 수준, 100% 수준, 120% 수준, 150% 수준으로 달리하여 50년 후의 ice ring 분포 양상을 비교하였다. 여기서 기준 값 대비 절연체의 열전도도가 증가하는 경우는 단열 성능이 저하됨을 의미한다. Fig. 8에서 알 수 있듯이 절연체의 열전도도의 단열 성능 차이(80%–150%)에 따라 ice ring의 최대 두께는 3 m에서 15 m까지 달라짐을 알 수 있다. 따라서 복공식 암반공동 주변의 ice ring의 두께를 신뢰성 있게 예측하기 위해서는 절연체의 열성능을 매우 정밀하게 산정할 필요가 있다. 또한, 목표로 하는 ice ring의 두께를 제어하기 위해서는 절연체의 단열 성능을 제어하는 것이 가장 효과적임을 알 수 있다.

Fig. 8

Changes in ice ring distribution after 50 years according to the different thermal conductivity of LRC insulators: the red line represents the 0°C isotherm. (a) 80% of the reference thermal conductivity, (b) 100% of the reference thermal conductivity, (c) 120% of the reference thermal conductivity, (d) Insulation 150% of the reference thermal conductivity

3.2.2 LRC 형상이 ice ring 반경에 미치는 영향

Fig. 9는 복공식 암반공동의 형상에 따른 ice ring 분포 양상을 동일시점(저장 50년 후)에서 비교한 결과이다. 저장소의 형상은 원형터널과 마제형 터널로 구분하였고, 각 저장소 체적(단면적)은 동일하게 설정하였다. Fig. 9에서 알 수 있듯이 LRC의 형상에 따라 인근 암반에서 형성되는 ice ring(0°C 등온선)의 확산 범위에 뚜렷한 차이가 발생하였는데, 마제형 저장소에 비해 원형 저장소에서 ice ring이 더 넓게 발달하였다. 이는 원형 형상이 열전달 측면에서 보다 기하학적으로 균일한 방사형의 경계 조건을 제공하기 때문이다. 원형 저장소의 경우, 저장소–암반 경계에서의 곡률 반경이 일정하여 저온 열유속이 전 방향으로 균일하게 분산되며, 이에 따라 냉기가 방사형으로 효과적으로 전파된다. 그 결과, 암반 내 등온선이 완만하고 연속적으로 형성되면서 동결 전선이 보다 먼 거리까지 안정적으로 확산되는 경향이 나타난 것으로 분석된다. 반면, 마제형 저장소는 곡률 변화가 크고 형상 불연속부가 존재하여 열유속이 국부적으로 집중되는 비균일한 분포를 보인다. 이러한 열전달 집중 현상은 저장소 인근에서 급격한 온도 구배를 형성하고, 결과적으로 저온 에너지의 공간적 확산을 제한하여 ice ring의 확장을 억제하는 요인으로 작용한다.

Fig. 9

Effect of LRC’s shape on ice ring distributions

따라서 본 해석 결과는 ice ring의 분포양상이 복공식 암반공동의 단열 성능 뿐 아니라 저장소 형상에 크게 의존함을 보여주며, 특히 원형 형상이 마제형에 비해 동결 영향 범위를 확대시키는 열역학적으로 보다 효율적인 구조임을 시사한다.

4. 결 론

본 연구에서는 액화수소 지하저장시스템 개념모델을 대상으로 예비 수치해석을 수행하고, 저장소 인근 암반에서의 장기 온도분포를 평가하였다. 본 연구를 통해 도출된 주요 결론은 다음과 같다.

(1) 저온 해석을 위하여 연속체 기반의 열-수리 해석 모델을 구축하였으며, 모델의 검증을 위하여 다공질 재료를 활용한 동결토 토조 실험을 수행하였다. 실험 결과와 수치해석 결과 비교시, 동결벽체 형상 및 형성 패턴이 서로 유사하게 나타났으며, 이를 통해 구축된 해석모델이 다공질 재료의 동결 거동을 신뢰성 있게 예측 가능함을 확인하였다.

(2) 검증된 수치해석 모델을 이용하여 액화수소 지하저장을 가정한 예비 해석을 수행하였다. 본 연구에서 설정한 복공식 지하공동 조건에서 저장소 인근 암반은 액화수소 저장 후 약 1년 경과 시 라이닝 최외곽 경계에서 영하 온도 영역에 진입하는 것으로 나타났다. 이후 냉각이 지속되어 50년 후에는 해당 경계의 온도는 약 –7℃까지 저하되었다. ice ring의 최대 두께는 액화수소 저장 후 약 30년 경과 시 약 2.5 m, 50년 경과 시에는 약 5 m 두께가 형성되었다.

(3) 설계 조건에 따른 해석 결과, ice ring의 분포 양상은 지하공동의 단열 성능 뿐만 아니라 공동 형상에 크게 의존하는 것으로 나타났다. 특히 원형 형상은 마제형에 비해 동결 영향 범위를 확대가능하여 열역학적으로 보다 효율적일 수 있음을 시사한다.

본 연구에서 제시한 해석 결과와 검증된 수치해석 모델은 액화수소 지하저장시스템의 개념 설계를 위한 기초 자료로 활용될 수 있다. 향후 연구에서는 저온 암반에서의 열-수리-역학적 복합거동 평가 모델을 개발하여, 장기 운영 조건에서의 역학적 안정성 평가가 수행되어야 할 것이다.