Research Paper

Journal of the Korean Society of Mineral and Energy Resources Engineers. 30 April 2022. 161-172
https://doi.org/10.32390/ksmer.2022.59.2.161

ABSTRACT


MAIN

  • 서 론

  • 횡등방성을 고려한 4차원 탄성파 모니터링

  •   탄성파 이방성

  •   횡등방성을 고려한 3차원 역시간 구조보정

  •   이방성을 고려한 4차원 탄성파 모니터링

  •   4차원 탄성파 모니터링에서의 횡등방성의 영향

  •   모델링에서의 횡등방성의 영향

  •   시간경과 탄성파 모니터링에서의 횡등방성의 영향

  •   북해 볼브 유전지역의 가상의 4차원 탄성파 모니터링

  • 결 론

서 론

최근 2050 탄소중립 실현을 위한 이산화탄소 저감기술로써, 이산화탄소 지중저장 기술(Carbon Capture and Storage, CCS)이 큰 주목을 받고 있다. 국내 이산화탄소 지중저장 연구는 2010년 이후로 육상과 해양에 적합한 부지를 탐사하기 위한 연구를 중심으로 활발히 수행되어 왔다. 2013년 시작한 포항 분지 해상 중소규모 CO2 저장실증 프로젝트를 통하여 2017년 국내 최초로 이산화탄소 주입 실증에 성공하였으나(Kwon, 2018), 지진 발생 가능성에 대한 사회적 우려로 국내 CCS 프로젝트는 한동안 중지된 상태였다. CCS에 적합한 저류층을 보유한 동해 가스전의 생산이 종료되어 감에 따라 지진 피해 가능성이 적은 먼바다에서의 CCS 프로젝트가 다시 활발하게 논의 중이지만, 국내 이산화탄소 모니터링 기술 실증은 아직 미흡한 상태이다. 본 연구의 목적은 지중 저장된 이산화탄소의 모니터링 역량이 우수한 노르웨이와의 국제 공동연구(2022 한-북유럽 협력 사업)를 위한 기초연구로, 역시간 구조보정 기반의 4차원 이방성 탄성파 모니터링 알고리듬을 확보하는 것이다.

4차원 탄성파 모니터링 기술은 반복적인 3차원 탐사를 통해 획득한 모니터링 자료로부터 3차원 탄성파 영상화를 수행하여 저장층의 물성 변화에 따른 반사도의 변화를 추정하는 기법이다(Lumley, 2010). 노르웨이 슬라이프너 프로젝트에서 성공 가능성이 검증된 4차원 탄성파 모니터링 기술은 최근 이산화탄소 지중저장 사업에서 주입 후 이산화탄소의 거동을 모니터링하기 위한 방법으로 주목을 받고 있다(Chadwick et al., 2004; Arts et al., 2008).

4차원 탄성파 모니터링 기술의 핵심은 3차원 탄성파 탐사자료를 반복 영상화하여, 탐사 자료의 변화를 통해 지하의 변화를 추정하는 것으로, 대량의 3차원 탄성파 탐사자료를 영상화하기 위한 구조보정 기술을 선택하는 것이 중요하다. 4차원 탄성파 모니터링에서는 주입 이산화탄소에 의한 약한 신호의 변화를 통해 이산화탄소의 정확한 공간적인 분포를 추정해야 하기 때문에 가능한 많은 자료를 중합하여 신호 대 잡음비를 향상시킬 수 있면서 공간적으로 정확하게 영상화 할 수 있는 구조보정 기술이 요구된다. 따라서 3차원 구조보정에 따른 계산량, 정확성, 신호 대 잡음비 향상 가능성을 고려하면, 키르히호프 구조보정(Kirchhoff migration) 방법과 역시간 구조보정(Reverse-Time Migration, RTM) 방법이 가장 유망할 것으로 보인다. 대규모 3차원 탐사자료의 영상화 과정에서 계산 효율성과 영상의 분해능을 고려하면 키르히호프 구조보정(Buske, 1999)이 더 선호되는 반면, 영상화의 반복성과 정확성을 고려하면 수치모델링을 기반으로 한 역시간 구조보정(Baysal et al., 1983)이 더 적합할 수 있다. 역시간 구조보정은 왕복 주시 파동 방정식을 사용하기에 계산시간이 오래 걸린다는 단점이 존재하나, 모든 반사파를 구현하기 때문에 복잡한 지질구조를 파악하는데 용이하다는 장점이 존재한다. 최근 고성능 컴퓨팅 기술의 발달과 더불어 3차원 역시간 구조보정이 점점 키르히호프 구조보정을 대체하고 있으며(Chang and McMechan, 1989; Yoon et al., 2003; Dussaud et al., 2008), 역시간 구조보정의 성능을 향상시키기 위한 다양한 기술이 활발하게 개발되고 있다(Perrone et al., 2011; Crawley et al., 2012; Oh et al., 2021). 역시간 구조보정은 파동장 연산자를 계산하여 반복적으로 활용하므로, 정확한 지층 영상 도출이 반복 요구되는 시간경과 탄성파 모니터링에 적합하다. 따라서 앞으로의 발전 가능성을 고려하면 역시간 구조보정을 기반으로 한 4차원 탄성파 모니터링 기술에 대한 지속적인 연구가 필요하다.

역시간 구조보정 결과의 품질은 입력 속도 정보에 의해 크게 영향을 받는다. 지하공간의 방향에 따라 탄성파의 속도가 달라지는 이방성이 강한 지역의 지층을 영상화 할 경우, 탄성파 이방성에 의해 오프셋에 따라 도달시간이 달라지기 때문에(Veeken, 2013), 등방성(isotropy)을 가정한 영상화에서는 먼 오프셋 자료를 배제해야 한다. 하지만 신호 대 잡음비가 낮은 탐사자료의 경우에는 가능한 많은 자료를 이용해야 중합(stacking) 과정에서 신호 대 잡음비가 향상될 수 있다. 또한 4차원 시간경과 모니터링 시, 주입 이산화탄소에 의한 반사파의 진폭 변화는 작기 때문에, 보다 많은 자료를 이용해야만 주입 이산화탄소에 의한 반사도 변화를 증폭시킬 수 있다. 따라서 이방성이 강한 지역의 시간경과 모니터링에서는 이방성을 고려한 영상화 기술을 이용하여 모든 오프셋 자료를 영상화에 활용하는 것이 바람직하다.

역시간 구조보정 기반의 시간경과 모니터링 알고리듬은 기존의 등방성 파동방정식을 이방성을 고려한 파동방정식으로 대체함으로써 구현 가능하다. 이때, 어느 정도 수준의 이방성까지 고려하는지에 따라 역시간 구조보정에 소요되는 계산량에 큰 차이가 있다. 고성능 컴퓨팅 기술의 발전과 함께, 최근에는 회전된 사방정계 형태의 탄성 이방성 매질까지 영상화가 가능하다(Oh and Alkhalifah, 2018). 실제 지구 지층은 거의 대부분이 이방성 매질로 이루어져 있기 때문에 이방성이 강한 지층 하부의 이산화탄소의 거동을 모니터링하기 위해서는 이방성 매질을 고려한 알고리듬을 이용하여 구조보정을 해주어야 한다. 4차원 모니터링에서는 부지의 특성에 따라 다양한 형태의 탄성파 이방성이 나타날 수 있으며, 역시간 구조보정의 계산량을 줄이기 위해서는 해당 CCS 부지의 이방성-탄성 특성과 모니터링 자료에 미치는 영향을 고려하여 적합한 이방성 파동방정식을 이용하는 것이 중요하다.

Alkhalifah(2000)는 유사미분연산자를 이용한 횡등방성 매질에 대한 유사음향파동방정식을 제안하였으며, Duveneck et al.(2008)은 유사음향파동방정식의 공간 4차 미분 항을 피하기 위하여 탄성파동방정식으로부터 횡등방성 매질에 대한 음향파동방정식을 유도하였다. 하지만 Duveneck et al.(2008)의 방법은 상대적으로 계산량이 많고, 비타원형 횡등방성 매질에 송신원이 존재할 경우, S파 잡음이 발생한다는 한계점이 있다. Xu and Zhou(2014)Alkhalifah (2000)의 유사음향파동방정식을 개선하여, 음향파동방정식을 등방성항과 이방성 보정항으로 분리한 구형분리 기반의 음향파동방정식을 제안하였으며, Xu et al.(2015)는 이를 더 개선하여 음향파동방정식을 타원형 이방성항과 비타원형 보정항으로 분리한 타원형분리 기반의 음향파동방정식을 제안하였다.

북해 지역은 횡등방성이 강한 지역으로(Banik, 1984), 본 연구에서는 횡등방성을 가정하기에 용이한 타원형 분리 기반의 이방성 음향파동방정식(Xu et al., 2015)을 기반으로 한 4차원 탄성파 모니터링 알고리듬을 제안한다. 단순한 횡등방성 주입 모델을 가정한 시간경과 모니터링을 통해 이방성이 강한 지역에서 획득된 탐사자료에 대한 이방성 모니터링 알고리듬의 정확성을 검토하고, 타원형 이방성을 가정하는 것이 계산량 측면에서 더 효율적임을 확인한다. 마지막으로 횡등방성이 강한 북해지역의 3차원 주입모델에 대한 시간경과 이방성 탄성파 모니터링의 적용성을 검증한다.

횡등방성을 고려한 4차원 탄성파 모니터링

탄성파 이방성

이방성이란 방향에 따라 매질의 물성이 달라지는 것을 의미한다. 매질의 이방성에 의해 탄성파의 전파 과정에서 방향에 따라 속도가 달라지는 탄성파 이방성 현상이 발생할 수 있다. 기존의 탄성파 탐사는 미임계반사(subcritical reflection) 영역 내의 반사파를 대상으로 수행되었으나, 최근에는 영상 품질 향상 및 역산 적용성 향상을 위해 먼 오프셋까지 수진기를 배치하여 초임계반사(supercritical reflection) 영역의 자료까지 획득하고 있다. 미임계반사는 반사각도가 작아서 거의 수직 방향의 전파 성분이 우세한 반면, 굴절파나 초임계반사파는 수평 방향의 전파 성분이 상대적이 커지게 된다. 탄성파 이방성이 존재하는 지역에서는 파동의 수직 방향과 수평 방향의 속도가 다르기 때문에, 만약 먼 거리의 오프셋에서 기록된 자료를 영상화에 이용하고자 한다면, 등방성 가정을 통한 하나의 P파 속도로는 미임계반사파와 초임계반사파 모두를 정확하게 보정하는 것이 어렵다. 따라서 이방성이 강한 지역에서 유가스 탐사나 이산화탄소 지중저장 모니터링을 위한 탄성파 탐사 자료가 획득된다면, 이방성을 적절하게 고려할 수 있는 영상화 알고리듬이 적용되어야 한다.

유가스 탐사나 이산화탄소 지중저장 부지 탐사의 대상이 되는 퇴적분지에서 이방성이 발생하는 원인은 크게 세 가지가 있다(Oh, 2020). 첫 번째는 셰일처럼 암석(또는 이를 구성하는 광물)이 가지고 있는 고유한 이방성이다. 두 번째는 얇은 지층이 교호하는 퇴적분지에서, 이러한 수직적인 층서 구조의 변화가 탄성파 신호의 파장에 비해 너무 작을 때 발생하는 겉보기 이방성이다(Backus, 1962). 세 번째는 퇴적분지 내에 수직적인 균열대가 일정한 방향성을 가지고 존재할 때 유발되는 이방성이 있다. 퇴적 분지에서 흔하게 나타나는 이방성은 종등방성(Horizontal Transverse Isotropic, HTI)과 횡등방성(Vertical Transverse Isotropic, VTI)이다. 앞서 언급한 세 가지 원인 중 첫 번째와 두 번째 원인에 의해 횡등방성이 흔하게 나타나며, 세 번째 원인은 종등방성을 유발할 수 있다.

실제 퇴적분지는 대부분 퇴적층으로 되어있어 횡등방성이 흔하게 나타난다. 횡등방성 매질에서는 일반적으로 수평 방향의 탄성파 속도가 수직 방향의 속도보다 빠르다. 횡등방성 매질에서의 P파의 위상속도 변화는 이방성 변수 𝜀, 𝛿를 이용하여 쉽게 표현할 수 있다(Thomsen, 1986). 횡등방성 매질에서의 전파 각도에 따른 P파의 위상속도(Vh, phase velocity)는 다음과 같이 표현된다(Fomel and Grechka, 2001):

(1)
Vh2(θ)=Vz2(1+2δsin2θcos2θ+2εsin4θ).

Vz는 P파의 수직 속도로 0도의 각도로 전파하는(즉, 수직 방향) 위상속도 Vh(0)와 같다. 횡등방성의 정도는 이방성 변수 𝜀과 𝛿의 크기에 의해 정의된다. 식 (1)에서, 𝜃는 파동이 입사하는 각도로 이방성 변수 𝜀은 sin4θ의 계수이기 때문에 수직 방향에 대해 90도의 각도(즉, 수평 방향)로 전파하는 파동의 속도를 가장 빨라지게 한다. 이방성 변수 𝛿는 sin2θcos2θ의 계수이기 때문에 45도의 전파 속도를 달라지게 한다. 2차원에서 𝜀과 𝛿의 크기가 같을 경우, 파면은 완전한 타원형이 되기 때문에, 이러한 특수한 경우를 타원형 이방성(elliptical anisotropy)라고 한다. 따라서 횡등방성이 존재하는 지역에서는 오프셋에 따라 이방성의 영향이 달라지기 때문에, 먼 오프셋 자료까지 활용한 탄성파 영상화 및 모니터링에서는 이방성에 대한 고려가 필요하다.

미임계반사 영역의 작은 반사 각도를 가지는 반사파는 주로 수직 방향으로 전파하기 때문에 이방성의 영향이 적다. 하지만 오프셋이 커질수록 이방성 변수 𝛿에 의해 파동의 전파 속도가 빨라지거나 느려질 수 있고, 더 먼 오프셋에서는 일반적으로 양의 값을 가지는 이방성 변수 𝜀에 의해 파동의 전파 속도가 빨라지게 된다. Fig. 1은 여러 가지 이방성 변수에 의한 파면의 변화를 보여준다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F1.jpg
Fig. 1.

The wavefronts of the P-wave depending on different anisotropic parameters (modified from Fomel and Grechka, 2001). The dashed lines indicate the wavefronts from the isotropic P-wave in (a).

횡등방성을 고려한 3차원 역시간 구조보정

탄성파 영상화 과정에서 이방성을 구현하기 위해서는 파동방정식을 기반으로 한 역시간 구조보정 방법이 유리하다. 등방성 파동방정식을 기반으로 한 기존의 3차원 역시간 구조보정 기술은 다음과 같은 3차원 등방성 훅의 법칙(Hooke’s law)에 의한 음향 파동방정식을 기반으로 파동장을 예측한다:

(2)
1v22pt2=2px2+2py2+2pz2+f.

이 때, p는 압력파동장, x,y,z는 공간 변수, t는 시간 변수이며, 𝜈는 매질의 P파 속도, f는 에어건과 같은 음향 송신원의 파형을 나타낸다. 등방성 파동방정식 기반의 역시간 구조보정은 식 (3)과 같이 송신원으로부터 온 송신원 파동장(SISO)와 수신기에서의 기록된 신호를 역시간으로 전파시켜 획득한 수신기 파동장(RISO)의 상호상관(Fig. 2)을 통하여 지하의 반사면 정보인 구조보정 영상(Imaging condition, IC)을 얻을 수 있다(Zhou et al., 2018):

(3)
ICISO(x,y,z)=stSISO(t,x,y,z)RISO(t,x,y,z).

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F2.jpg
Fig. 2.

Schematic diagram illustrating the mechanism of reverse-time migration: (a) source wavefields, (b) receiver wavefields and (c) imaging condition (red point) obtained by the cross-correlation between source and receiver wavefields.

하지만, 앞서 설명한 것처럼 이방성이 강한 지역에서는 수평 방향의 전파 속도가 수직 방향의 속도보다 빠르기 때문에, 오프셋이 커짐에 따라서 이방성의 영향이 커지게 된다. 등방성 파동방정식으로 계산한 구조보정 영상(ICISO)은 모든 오프셋이 자료에 대해서 동일한 수직 속도로 깊이 정보를 추정하기 때문에, 오프셋이 커질수록 실제보다 더 느린 속도가 적용되어 추정된 반사면의 깊이가 더 얕아지게 된다. 이러한 횡등방성에 의한 오차는 오프셋이 커짐에 따라서 반사면의 깊이가 얕아지는 결과를 준다, 따라서 오프셋 영역(offset domain)이나 각도 영역(angle domain)에서 반사면 정보를 확인할 경우, 하키 스틱 효과(hockey stick effect)가 발생하여 실제 반사점 위치에서 모든 신호가 중합되기 어렵다(Veeken, 2013).

파동방정식 기반의 역시간 구조보정에서는 식 (2)로 표현된 파동방정식을 이방성을 고려할 수 있는 파동방정식으로 바꾸어 주면 손쉽게 이방성을 고려할 수 있다. 이를 위한 다양한 유사음향파동방정식(Kosloff et al., 1989; Zhou et al., 2006; Bakker and Duveneck, 2011)이 개발되었는데, 본 연구에서는 Xu et al., (2015)가 제안한 타원형분리 기반의 이방성 음향파동방정식을 사용하였다:

(4)
1v022pt2=(1+2ε)2px2+(1+2ε)2py2+2pz2Se+f,
(5)
Se=121+1-8(ε-δ)nh2nz2((1+2ε)nh2+nz2)2.

식 (2)와 비교하여, ν0는 수직 방향 P파 속도이고, nn=k|k|로 정의되는 파수 영역에서의 위상 방향 벡터이며, nx,nynz는 각 방향에 따른 성분이다. 이때, 수평 방향 성분은 nh2=nx2+ny2이다. 이방성 변수 𝜀과 𝛿는 전파 방향에 따른 P파 속도의 차이를 결정하며, Se는 각 지점에서의 비타원적인 전파 특성을 보정하기 위한 상수이다. 𝜀과 𝛿가 같은 타원형 이방성에서는 Se=1이 되고, 이방성 음향파동방정식 (4)는 등방성 음향파동방정식 (2)와 매우 유사해진다. 따라서 프로그래밍 관점에서 타원형 이방성에 대한 이방성 음향파동방정식을 계산하기 위해 필요한 메모리와 시간은 등방성 음향파동방정식을 계산하는데 필요한 메모리, 시간과 거의 비슷하다. 이는 타원형 이방성으로 횡등방성을 근사할 경우, 3차원 영상화 또는 4차원 모니터링에서 추가되는 계산량이 거의 없어 매우 효율적일 수 있음을 지시한다. 최종적으로 이방성을 고려한 구조보정 영상 ICVTI은 다음과 같이 이방성 음향파동방정식을 이용해 계산한 송신원 파동장과 수신기 파동장의 상호상관으로 구현될 수 있다:

(6)
ICVTI(x,y,z)=stSVTI(t,x,y,z)RVTI(t,x,y,z).

타원형 이방성을 가정한 역시간 구조보정 영상은 식 (4)에서 타원형 이방성을 가정하여(𝜀=𝛿 즉, Se=1) 계산한 송신원 파동장과 수신기 파동장으로 계산될 수 있다. 본 연구에서는 역시간 구조보정 영상을 구축하기 위하여 완전파형역산(Full-waveform inversion, FWI)을 기반으로 한 알고리듬을 이용하였다(Virieux an Operto, 2009).

이방성을 고려한 4차원 탄성파 모니터링

4차원 탄성파 모니터링을 이용하면 지구 내부가 시간에 따라 어떻게 변하는지 추정할 수 있다(Lumley, 2010; Park et al., 2021). 4차원 탄성파 모니터링을 하려면 2가지 가정이 이루어져야 한다. 첫 번째는 지구의 지질학적인 시간 척도에 비하여 이산화탄소 주입 행위가 너무 짧아서 주입이 일어나지 않은 부분은 정적인 상태로 근사 될 수 있다는 것이고, 두 번째는 주입된 이산화탄소가 매질의 물성(특히, P파 속도)을 변화시킨다는 것이다. 즉, 저장층의 주입 전 물성과 이산화탄소의 물성을 고려했을 때, 주입된 이산화탄소에 의해 P파 속도의 변화가 발생해야 하고, 이산화탄소가 주입되지 않은 지역에서는 물성 변화가 없어야 한다. 일반적으로 액체로 포화된 염대수층에 기체인 이산화탄소를 주입하면, P파 속도가 감소하기 때문에, 4차원 시간경과 모니터링 기술을 적용하기에 용이하다.

Fig. 3은 역시간 구조보정을 기반으로 한 4차원 탄성파 모니터링의 적용 과정을 보여준다. 주입 전 기준탐사(baseline)로부터 구조보정 영상(IC0ISO,IC0VTI)을 획득하고, 주입 후 반복적인 모니터링 탐사(monitoring survey)를 통해 구조보정 영상(IC1ISO,IC1VTI 등)을 획득한다. 이방성을 고려한 역시간 구조보정 기반의 4차원 탄성파 모니터링 알고리듬은 기존의 등방성 역시간 구조보정을 통해 획득한 구조보정 영상(ICISO)을 이방성 역시간 구조보정을 통해 획득한 구조보정 영상(ICVTI)으로 바꾸어 주기만 하면 된다. 4차원 탄성파 모니터링에서도 타원형 이방성만을 고려했을 때, 구조보정 영상을 계산하기 위해 필요한 메모리와 시간이 등방성 역시간 구조보정과 거의 같으므로, 매우 효율적이다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F3.jpg
Fig. 3.

Flowchart showing the overall procedure in 4D time-lapse seismic monitoring based on reverse time migration. Blue and green arrows indicate isotropic and anisotropic monitoring procedures, respectively.

4차원 탄성파 모니터링에서의 횡등방성의 영향

이번 장에서는 지하 매질의 횡등방성이 4차원 탄성파 모니터링에서 어떤 영향을 주는지 확인한다. 4차원 탄성파 모니터링 기술은 3차원 구조보정을 통한 지하 반사면의 변화를 추정하는 것이며, 3차원 구조보정은 송신원 파동장과 수신기 파동장의 수치모델링을 기반으로 수행되기 때문에, 결국 이방성의 영향을 수치모델링 측면에서부터 이해하는 것이 중요하다.

모델링에서의 횡등방성의 영향

파동 전파 과정에서의 횡등방성의 영향을 확인하기 위하여, 간단한 3차원 균질 매질을 가정하였다. 수치 모델은 가로, 세로, 깊이가 모두 4 km인 정육면체 형태로 설정하였고, 각 방향으로의 공간격자 간격은 20 m로 설정하였다. 모델의 정중앙에 위치한 송신원의 파형은 최대 주파수가 10 Hz인 1차 미분 가우스함수로 가정하였으며, 매질의 P파 속도는 1.5 km/s로 설정하고, 이방성 변수에 따라 총 세 가지 모델을 만들어 비교하였다. 각 모델에서 이방성 변수 (𝜀, 𝛿)는 각각 (0, 0), (0.5, 0.5), (0.5, 0)이다. 첫 번째 모델은 등방성 모델, 두 번째 모델은 타원형 이방성 모델, 세 번째 모델은 횡등방성 모델에 해당한다.

Fig. 4는 각각의 모델에서 파동 전파 순간포착 영상을 보여준다. 등방성 모델의 경우 모든 방향으로 동일한 속도로 파면이 구형으로 형성되는 것을 확인할 수 있다. 각 순간포착 영상 단면에서 빨간색 원은 등방성 모델링에서 획득된 파면이다. 이방성 모델에서는 𝜀이 0.5로 매우 큰 값을 가지므로 수평 방향으로 훨씬 빠르게 전파함을 확인할 수 있다. 타원형 이방성 모델에서는 x-z, y-z 평면에서 타원형의 파면이 형성되고, VTI의 경우 𝛿가 –0.1로 설정되어 식 (1)에 따라 대각 방향의 속도가 느려지면서 수직 단면에서 파면이 마름모 형태로 형성됨을 확인할 수 있다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F4.jpg
Fig. 4.

Snapshots from (a) isotropic model, (b) elliptical anisotropic model, and (c) VTI model at 2 sec. Epsilon and delta are 0.5 and 0 in the elliptical model (b). In VTI (c) model, epsilon is the same as that in the elliptical model, and the delta is set at -0.1. Red circles indicate the wavefronts obtained by isotropic modeling in (a). To compare wave speeds, we used red circles to indicate the wavefronts in the isotropic media.

이방성 변수 𝜀의 영향에 비해 𝛿의 영향은 상대적으로 크지 않다. 이는 전파 각도에 따른 P파 위상속도에 대한 식 (1)에서도 직관적으로 유추할 수 있다. 이방성 변수 𝛿의 각도에 따른 영향을 반영하는 sin2θcos2θ은 45도 방향에서 최대값이 0.25이다. 따라서 동일한 이방성 변수 𝜀과 𝛿가 주어진다고 하더라도, 𝛿가 P파 속도에 미치는 영향은 𝜀의 1/4배이다. 따라서 등방성 음향파동방정식과 비교하면 거의 같은 계산 비용으로, 타원형 이방성을 가정한 이방성 음향파동방정식만으로도 횡등방성에 의한 P파 속도의 변화를 어느 정도 반영할 수 있으며(Alkhalifah, 2016; Oh and Alkhalifah, 2019), 횡등방성을 구현하기 위해 추가되는 계산량과 기대효과를 고려하여 횡등방성으로의 완전한 근사 여부가 결정되어야 한다.

시간경과 탄성파 모니터링에서의 횡등방성의 영향

이제 역시간 구조보정 기반의 시간경과 탄성파 모니터링에서의 횡등방성의 영향을 확인한다. 횡등방성 매질에서는 모든 수평방향 속도가 동일하기 때문에, 수평층에 대한 3차원 횡등방성 모델은 임의의 방위각에 대한 2차원 횡등방성 모델로 근사할 수 있다. 따라서 효율적인 계산을 위하여 본 수치 예제에서는 2차원 공간 기반의 시간경과 탄성파 모니터링 수치 예제를 수행하였다.

2차원 모델(Fig. 5)은 가로 8 km, 깊이 2 km로, 공간격자 간격은 가로와 깊이 방향 모두 20m, 총 기록 시간은 2초로 설정하였다. 시간 격자의 간격은 0.002초이며, 최대 주파수가 10 Hz인 리커파형(ricker wavelet)을 송신원으로 사용하였다. 총 100개의 송신원이 표면에 80 m 간격으로 설치되고, 400개의 수신기가 모든 표면 격자에 설치되었다고 가정하였다. 수치모델의 배경 속도는 2층 모델을 가정하였으며, 첫 번째 층의 P파 속도는 1.5 km/s, 두 번째 층의 P파 속도는 3 km/s로 설정하였다. 가로 2 km 너비, 두께 0.2 km로 이산화탄소가 주입되었다고 가정하여 두 번째 층 상부에 P파 속도를 감소시켰다. 첫 번째 층의 이방성 변수 (𝜀, 𝛿)는 (0.3, 0.1)이고 두 번째 층은 등방성으로 가정하였다. 주입 유체인 이산화탄소도 등방성을 보이기 때문에 주입 후에도 두 번째 층에서의 이방성은 나타나지 않는다고 가정하였다. 따라서 이 수치 예제를 통해서 강한 이방성이 존재하는 첫 번째 층에 의해 두 번째 층 상부에 주입된 이산화탄소의 모니터링이 적절하게 수행될 수 있는지 확인하고자 한다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F5.jpg
Fig. 5.

2-layers VTI model: P-wave velocities (a) before and (b) after injection and anisotropic parameters (c) 𝜀 and (d) 𝛿.

Fig. 6은 횡등방성 자료에 대한 등방성 기반, 타원형 이방성 기반, 횡등방성 기반의 시간경과 탄성파 모니터링 기술의 오프셋에 따른 정확도 차이를 보여준다. 오프셋이 짧은(0~2 km) 자료만 사용하였을 경우에는 수직 속도 성분이 우세한 미임계반사파 위주로 모니터링이 수행되기 때문에, 모든 방법이 유사한 모니터링 결과를 보여준다. 즉, 탐사자료의 오프셋이 짧을 경우, 등방성 기반의 시간경과 탄성파 모니터링 알고리듬으로도 충분히 좋은 결과를 기대할 수 있다. 하지만 긴(2~4 km) 오프셋 길이의 자료만 사용하였을 경우에는 수평 속도 성분이 우세하기 때문에 등방성 기반의 모니터링 자료에서 반사파의 위치가 왜곡되어 나타난다. 이는 오프셋이 길어질 경우, 등방성 기반의 시간경과 탄성파 모니터링 기술은 하키 스틱 효과에 의해 먼 오프셋 자료에서 추정된 주입 이산화탄소의 깊이가 상부로 올라가기 때문에, 모든 오프셋 자료가 정확한 깊이에서 반사점을 형성하지 않는다. 이러한 이유로 긴 오프셋 자료만을 사용할 경우, 파장이 긴 반사파가 주로 영상화 되어 정확한 지층 내부 구조의 파악이 어렵다. 횡등방성 기반의 시간경과 탄성파 모니터링 기술은 모든 오프셋 자료가 주입 이산화탄소의 깊이를 동일하게 추정할 수 있다. Fig. 7은 획득된 모든 오프셋(0~4 km)을 사용하여 모니터링한 결과로, 횡등방성을 가정했을 때 이산화탄소 주입으로 인한 속도 변화를 반영하는 반사신호들이 잘 중합됨을 관찰할 수 있다. 실제 현장 자료에 대한 적용을 고려해 본다면, 주입 이산화탄소에 의한 진폭 변화는 배경 잡음이나 다른 비반복성 요인에 비해 매우 작을 수도 있다. 따라서 보다 많은 자료를 활용하여 주입 이산화탄소에서 반사된 신호의 신호 대 잡음비를 향상시키는 것이 중요하며, 횡등방성 기반의 시간경과 탄성파 모니터링 기술은 보다 많은 오프셋 자료를 이용함으로써 신호 대 잡음비를 향상시키는 것이 가능하다. 타원형 이방성 기반의 탄성파 모니터링 기술은 횡등방성 기반의 탄성파 모니터링과 거의 유사한 결과를 준다. 이는 앞서 언급한 것처럼 이방성 변수 𝛿의 의한 영향이 그리 크지 않기 때문이다. 하지만 횡등방성 기반의 이방성 음향파동방정식에서 매 시간간격마다 식 (5)Se를 계산하는 데 추가적인 비용이 소요되는 것을 생각하면, 타원형 이방성으로 횡등방성을 근사하는 것이 대부분의 상황에서 효율적일 것이다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F6.jpg
Fig. 6.

Time-lapse seismic monitoring results: isotropic monitoring (a and d), elliptical anisotropic monitoring (b and e), and VTI monitoring (c and f) with only short-offset (0~2 km) data (a, b, and c) and with long-offset (2~4 km) data (d, e, and f).

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F7.jpg
Fig. 7.

Time-lapse seismic monitoring results: (a) isotropic monitoring, (b) elliptical anisotropic monitoring, and (c) VTI monitoring with the offsets from 0 to 4 km.

북해 볼브 유전지역의 가상의 4차원 탄성파 모니터링

현실적인 3차원 지질 모델에 대한 이방성을 고려한 4차원 탄성파 모니터링 기술의 적용성을 검토하기 위하여, 북해 볼브 유전지역의 저류층에 가상으로 이산화탄소가 주입되었다고 가정하였다. 4차원 시간경과 모니터링에 사용한 변수는 Table 1과 같다. 400개의 송신원과 1600개의 수신기가 표면에 일정한 간격으로 설치되어 있다고 가정하였다.

Table 1.

Parameters for 4D time-lapse seismic monitoring of virtual CO2 injection in Volve oilfield

Parameter Value
Survey area Inline 11.0 km
Crossline 6.0 km
Depth 4.5 km
Grid size 50 m
Recording time 5 s
Time sampling 0.002 s
Number of sources 400
Number of receivers 1600
Source wavelet Ricker function
Fmax of source wavelet 10 Hz
Frequency range 5~10 Hz
Frequency interval 0.2 Hz

Fig. 8은 노르웨이 국영석유회사 Equinor에 의해 제공된 볼브 유전지역의 P파 속도, 이방성 변수 모델이다(Szydlik et al., 2006). Szydlik et al.(2006)에 따르면 북해 볼브지역의 저류층인 약 3 km 깊이에 존재하는 Hugin 사암층 상부로 강한 횡등방성이 존재한다. 저류층인 Hugin 사암층에 가상으로 2단계의 이산화탄소 주입이 이루어졌다고 가정하여 P파 속도를 감소시킨 주입 후 모델을 구축하였다. 이전 수치 예제와 마찬가지로, 주입 전후 배경 지층은 정적인 상태이고 이산화탄소의 주입이 저장층의 이방성 변수에 영향을 끼치지 않을 것이라 가정하였기 때문에, 이방성 모델은 모든 주입 단계에서 동일하다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F8.jpg
Fig. 8.

(a) P-wave velocity model before injection, differences in P-wave velocity after (b) 1st injection, and (c) 2nd injection and the anisotropic parameters, (d) 𝜀 and (e) 𝛿.

Fig. 9는 사용한 오프셋에 따른 등방성 기반 시간경과 탄성파 모니터링과 횡등방성 기반 시간경과 탄성파 모니터링의 차이를 보여준다. 오프셋은 각각 0~2 km, 2~4 km, 4~6 km로 총 3가지 경우를 비교하였다. 0~2 km 오프셋은 이방성의 영향이 적은 미임계반사파 위주로 시간경과 탄성파 모니터링이 수행되기 때문에 두 알고리듬이 거의 동일한 결과를 보여준다(Fig. 9a and d). 2~4km 오프셋 거리에서부터 이방성 영향이 나타나기 시작하였으며, 등방성 기반의 시간경과 탄성파 모니터링에서는 주입 이산화탄소의 상부 반사면이 실제 반사면보다 더 상부에서 나타났다(Fig. 9b and e). 이는 실제로 이방성에 의해 빠른 속도로 전파한 반사파를, 등방성 가정에 의해 더 느린 속도로 구조보정했기 때문에 나타난 하키스틱 효과에 의한 것이다. 이러한 오차는 4~6 km 오프셋의 자료를 이용한 결과에서 더 명확하게 구분된다(Fig. 9c and f).

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F9.jpg
Fig. 9.

4D time-lapse seismic monitoring results after 1st injection using isotropic (a, b and c) and VTI (d, e and f) approximations with offset ranges of 0~2 km (a and d), 2~4 km (b and e) and 4~6 km (c and f).

Fig. 10은 0~6 km 의 전체 오프셋 자료를 사용했을 때, 등방성 기반, 타원형 이방성 기반, 횡등방성 기반의 시간경과 탄성파 모니터링 결과를 보여준다. Fig. 9에서 확인했던 먼 오프셋 자료에서의 왜곡에 의해 등방성 기반 시간경과 모니터링 결과는 주입 이산화탄소의 수직적인 분포를 추정하는데 오차가 발생하였다. 하지만 횡등방성 기반의 시간경과 탄성파 모니터링 기술은 보다 정확하게 주입 이산화탄소의 깊이 정보를 추정하였으며, 보다 많은 오프셋 자료가 중합되면서 짧은 오프셋 자료만을 이용한 등방성 기반 시간경과 탄성파 모니터링 결과(Fig. 9a)와 비교했을 때, 반사면이 더 선명하게 드러났다. 따라서 신호 대 잡음비가 더 낮은 상황에서도 더 좋은 결과를 줄 것으로 기대된다. 타원형 이방성 기반의 시간경과 탄성파 모니터링 결과는 횡등방성 기반의 탄성파 모니터링 결과와 매우 유사하다. 즉, 등방성 기반 시간경과 탄성파 모니터링과 거의 유사한 계산량으로 이방성에 의한 왜곡을 줄일 수 있기 때문에 실제 대규모 4차원 탄성파 모니터링 문제에서 보다 효율적인 알고리듬이 될 것이다.

/media/sites/ksmer/2022-059-02/N0330590204/images/ksmer_59_02_04_F10.jpg
Fig. 10.

4D time-lapse seismic monitoring results using isotropic (a and d), elliptical anisotropic (b and e), and VTI (c and f) approximations after 1st injection (a, b, and c) and 2nd injection (d, e, and f).

결 론

본 연구에서는 이방성 매질을 고려한 4차원 탄성파 모니터링 알고리듬을 제안하였다. 이방성이 강한 상부 지층 아래로 이산화탄소가 주입될 때, 기존의 등방성 모니터링 알고리듬은 이방성에 의해 오프셋에 따라 추정되는 반사면의 깊이가 달라지고, 이로 인해 주입 이산화탄소의 깊이 정보의 왜곡이 발생되며, 주입층 반사면의 해상도가 떨어지는 것을 확인하였다. 타원형분리 기반의 이방성 음향파동방정식을 이용한 이방성 모니터링 알고리듬은 모든 오프셋 자료로부터 동일한 심도 정보를 추정할 수 있기 때문에, 주입층 반사면에서 중합이 잘 되어 보다 정확한 깊이 정보와 선명한 영상을 제공할 수 있었다. 또한 타원형 이방성을 가정한 모니터링 알고리듬은 등방성 모니터링과 거의 같은 계산량으로 횡등방성을 근사할 수 있음을 확인하였다. 실제 모니터링 자료에 대해 적용하기 위해서는 지하의 이방성 정보를 추정하기 위한 다변수 역산 기술이 같이 활용되어야 하며, 노르웨이와의 후속 공동연구를 통하여 이를 검토할 계획이다.

Acknowledgements

이 논문은 2022년도 한국연구재단의 국제협력사업의 지원을 받아 연구되었습니다(NRF-2022K2A9A2A23000300). 또한 이 연구는 2021년도 정부(과학기술정보통신부) 재원으로 한국연구재단의 지원을 받아 수행된 연구입니다(No. 2021R1A2C1092820). 볼브 유전지역의 3차원 탐사자료를 제공해주신 Equinor에 감사드립니다.

References

1
Alkhalifah, T., 2000. An acoustic wave equation for anisotropic media, Geophysics, 65(4), p.1239-1250. 10.1190/1.1444815
2
Alkhalifah, T., 2016. Research note: Insights into the data dependency on anisotropy: An inversion prospective, Geophysical Prospecting, 64(2), p.505-513. 10.1111/1365-2478.12345
3
Arts, R., Chadwick, A., Eiken, O., Thibeau, S., and Nooner, S., 2008. Ten years' experience of monitoring CO2 injection in the Utsira Sand at Sleipner, offshore Norway, First Break, 26(1), p.65-72. 10.3997/1365-2397.26.1115.27807
4
Backus, G.E., 1962. Long‐wave elastic anisotropy produced by horizontal layering, Journal of Geophysical Research, 67(11), p.4427-4440. 10.1029/JZ067i011p04427
5
Bakker, P.M. and Duveneck, E., 2011. Stability analysis for acoustic wave propagation in tilted TI media by finite differences, Geophysical Journal International, 185, p.911-921. 10.1111/j.1365-246X.2011.04986.x
6
Banik, N.C., 1984. Velocity anisotropy of shales and depth estimation in the North Sea basin, Geophysics, 49(9), p.1411-1419. 10.1190/1.1441770
7
Baysal, E., Kosloff, D.D., and Sherwood, J.W., 1983. Reverse time migration, Geophysics, 48(11), p.1514-1524. 10.1190/1.1441434
8
Buske, S., 1999. Three-dimensional pre-stack Kirchhoff migration of deep seismic reflection data, Geophysical Journal International, 137(1), p.243-260. 10.1046/j.1365-246x.1999.00789.x
9
Chadwick, R.A., Zweigel, P., Gregersen, U., Kirby, G.A., Holloway, S., and Johannessen, P.N., 2004. Geological reservoir characterization of a CO2 storage site, The Utsira Sand, Sleipner, northern North Sea, Energy, 29(9-10), p.1371-1381. 10.1016/j.energy.2004.03.071
10
Chang, W.F. and McMechan, G.A., 1989. 3D acoustic reverse- time migration 1, Geophysical Prospecting, 37(3), p.243-256. 10.1111/j.1365-2478.1989.tb02205.x
11
Crawley, S., Whitmore, N.D., Sosa, A., and Jones, M., 2012. Improving RTM images with angle gathers, In 82nd SEG Technical Program Expanded Abstract, Society of Expolration Geophysicists, 4609p. 10.1190/segam2012-1370.122840967
12
Dussaud, E., Symes, W.W., Williamson, P., Lemaistre, L., Singer, P., Denel, B., and Cherrett, A., 2008. Computational strategies for reverse-time migration, In 78th SEG Technical Program Expanded Abstracts, Society of Exploration Geophysicists, p.2267-2271. 10.1190/1.3059336
13
Duveneck, E., Milcik, P., Bakker, P.M., and Perkins, C., 2008. Acoustic VTI wave equations and their application for anisotropic reverse-time migration, In 78th SEG Technical Program Expanded Abstracts, Society of Exploration Geophysicists, p.2186-2190. 10.1190/1.3059320
14
Fomel, S. and Grechka, V., 2001. Nonhyperbolic reflection moveout of P waves. An overview and comparison of reasons. CWP-372, Colorado School of Mines.
15
Kosloff, D., Filho, A.Q., Tessmer, E., and Behle, A. 1989. Numerical solution of the acoustic and elastic wave equations by a new rapid expansion method, Geophysical Prospecting, 37, p.383-394. 10.1111/j.1365-2478.1989.tb02212.x
16
Kwon, Y.K., 2018. Demonstration-scale Offshore CO2 Storage Project in the Pohang Basin, Korea, The Journal of Engineering Geology, 28(2), p.133-160.
17
Lumley, D., 2010. 4D seismic monitoring of CO2 sequestration, The Leading Edge, 29(2), p.150-155. 10.1190/1.3304817
18
Oh, J.W. and Alkhalifah, T., 2018. Optimal full-waveform inversion strategy for marine data in azimuthally rotated elastic orthorhombic media, Geophysics, 83(4), p.R307-R320. 10.1190/geo2017-0762.1
19
Oh, J.W. and Alkhalifah, T., 2019. Study on the full‐waveform inversion strategy for 3D elastic orthorhombic anisotropic media, Application to ocean bottom cable data, Geophysical Prospecting, 67(5), p.1219-1242. 10.1111/1365-2478.12768
20
Oh, J.W., 2020. Parallel elastic wave modeling considering the effects of orthorhombic anisotropy, Journal of Korean Society of Mineral and Energy Resources Engineers, 57(1), p.86-96 (in Korean with English abstract). 10.32390/ksmer.2020.57.1.086
21
Oh, J.W., Cheng, J., and Min, D.J., 2021. Diffraction-angle filtering of gradient for acoustic full-waveform inversion, Geophysics, 86(2), p.R173-R185. 10.1190/geo2019-0801.1
22
Park, S.E., Li, X., Kim, B.Y., Oh, J.W., Min, D.J., and Kim, H.S., 2021. Seismic Imaging of Ocean-bottom Seismic Data for Finding a Carbon Capture and Storage Site: Two- dimensional Reverse-time Migration of Ocean-bottom Seismic Data Acquired in the Pohang Basin, South Korea, Geophysics and Geophysical Exploration, 24(3), p.78-88.
23
Perrone, M.P., Lu, L., Liu, L., Fedulova, I., Semenikhin, A., and Gorbik, V., 2011. High performance RTM using massive domain partitioning, In 73rd EAGE Conference and Exhibition incorporating SPE EUROPEC 2011, European Association of Geoscientists & Engineers, p.cp-238-00150.
24
Szydlik, T.J., Way, S., Smith, P., Aamodt, L., and Friedrich, C., 2006. 3D PP/PS prestack depth migration on the Volve field, In 68th EAGE Conference and Exhibition incorporating SPE EUROPEC 2006, European Association of Geoscientists & Engineers, p.cp-2-00185.
25
Thomsen, L., 1986. Weak elastic anisotropy, Geophysics, 51(10), p.1954-1966. 10.1190/1.1442051
26
Veeken, P.P., 2013. Seismic stratigraphy and depositional facies models, Academic Press. 10.3997/9789073834439PMC3869884
27
Virieux, J. and Operto, S. 2009. An overview of full-waveform inversion in exploration geophysics. Geophysics, 74(6), p.WCC1-WCC26. 10.1190/1.3238367
28
Xu, S. and Zhou, H. 2014. Accurate simulations of pure quasi-P-waves in complex anisotropic media, Geophysics, 79(6), p.T341-T348. 10.1190/geo2014-0242.1
29
Xu, S., Tang, B. Mu, J., and Zhou, H., 2015. Elliptic decomposition of quasi-P wave equation, In 77th EAGE Conference and Exhibition incorporating SPE EUROPEC 2015, European Association of Geoscientists & Engineers, p.1-5. 10.3997/2214-4609.201413134
30
Yoon, K., Shin, C., Suh, S., Lines, L.R., and Hong, S., 2003. 3D reverse-time migration using the acoustic wave equation: An experience with the SEG/EAGE data set, The Leading Edge, 22(1), p.38-41. 10.1190/1.1542754
31
Zhou, H., Zhang, G., and Bloor, R., 2006. An anisotropic acoustic wave equation for VTI media, In 68th EAGE Conference and Exhibition incorporating SPE EUROPEC 2006, European Association of Geoscientists & Engineers, p.cp-2-00318. 10.3997/2214-4609.201402310
32
Zhou, H.W., Hu, H., Zou, Z., Wo, Y., and Youn, O., 2018. Reverse time migration: A prospect of seismic imaging methodology, Earth-Science Reviews, 179, p.207-227. 10.1016/j.earscirev.2018.02.008
페이지 상단으로 이동하기