海洋电磁法三维正反演研究进展

非线性共轭梯度法（NLCG）
搜索方向：
p0 C0g0 pk Ck g k k pk 1
T pk ( gk gk 1 ) 0
需满足的正交条件：
gT k Ck (g k g k 1 ) Polak-Ribiere 公式： k gT k 1Ck 1g k 1
韩波
• 2014 年我国进口原油 3.1 亿吨，同比增长 9.5% ，对外依存度为
59.6%！ ——中国海关总署2015年1月
• 海洋被认为蕴藏了全球超过70%的石油，据估计，全球深海区潜
在石油储量可能超过1000亿桶（1桶≈160升）。 - “海底的油气如同埋在地里的马铃薯一样等待我们去发掘” —— Michael Economides, 2007 • 海洋CSEM近些年成为海洋油气勘探中重要性仅次于地震勘探的
韩波
本研究的矢量FE法与基于A-Φ分解的矢量FE法 (Ansari et al., 2014)正演 结果对比（李建慧等，2015）
韩波
三维反演算法研究进展
目标函数
• 数据拟合： d Wd d F (m) • 模型约束： m Wm m m r
2
+
d i Fi (m) i i 1
目前被广泛应用于海洋板块构造研究、海底油气资源
勘探等领域。
天然场源：海洋大地电磁法（MMT） 人工场源：海洋可控源电磁法（MCSEM）
韩波
MCSEM
MMT
（Constable, 2010）
韩波
MMT工作原理 • • • • 与陆地MT基本相同 垂直入射，平面波 测量正交电磁场 计算卡尼亚响应
m k J C J k C
T k 1 d

1 1 m
J C d F(m ) C m
T k 1 d k 1 m
k
mr

CG迭代求解
特点： • • • • 利用了一阶导数和二阶导数（Hessian）信息 每次外层GN迭代包含多次内层CG迭代 每次CG迭代只需计算灵敏度矩阵或其转置与向量的乘积， 无需显式存储矩阵，计算快，节约内存 单次迭代收敛程度高，需要的迭代次数少，正演可采用直 接解法 韩波
韩波
数值离散化方法
电磁法中最常用的数值模拟方法包括积分方程法，有限
差分法，有限体积法和有限单元法。我们目前实现了三
种不同的数值模拟方法：
（1）基于二次场方程的交错网格有限体积法（韩波等，2015）；
（2）基于总场方程的交错网格有限体积法（彭荣华等，2015）； （3）基于总场方程的矢量有限单元法（李建慧等，2015）。
韩波
以交错网格FV法为例：
5×5×5
韩波
以交错网格FV法为例：
10×10×10
韩波
以交错网格FV法为例：
20×20×20
韩波
以交错网格FV法为例：
50×50×50
韩波
迭代解法
Krylov子空间方法（GMRES, BiCG, QMR等等） e0 , e1 , e2 , ...
e n e 0 K n (r0 , A),
1123
1119
10-8
10-11 10-5
7.30×10-21
7.32×10-21 2.25×10-22 2.21×10-22 2.23×10-22
4.63×10-13
4.65×10-13 3.52×10-14 3.46×10-14 3.49×10-14
韩波
1138
1125 1111 1118 1134
10
10-8 10-11
对于多场源问题，基于直接解法的FV三维正演计算时间（韩波等，2015）
韩波
基于直接解法的FV三维正演中，直 接解法的内存需求与解方程时间随 未知数个数的变化（韩波等，2015）
韩波
MUMPS求解器与迭代解法的并行可伸缩性（Scalability）比较 （ Puzyrev et al.，2015）
韩波
3D模型正演算例
模型一
• 模拟三维油气储集层 • 高阻层在x-y平面内的截面为正方形 • 发射频率1Hz，inline观测
韩波
FV1法正演计算结果：不含空气层（上排）和包含空气层（下排） （韩波等，2015） 韩波
模型二
(a) Sea 3.2 S/m Wire 100 m 100 m (0,0,0) 950 m Hydrocarbon 0.01 S/m 100 m 2000 m Seabed 1 S/m
• • • • 最速下降法（Steepest descent） 牛顿法（Newton） 高斯牛顿法（Gauss-Newton） 拟牛顿法（Quasi-Newton）
•
非线性共轭梯度法（Nonlinear Conjugate gradient）
目前，我们实现了非线性共轭梯度（韩波，2015）和高斯牛 顿-共轭梯度（Peng et al., 2015）两种反演方法 韩波
韩波
MCSEM三维正演算法研究进展
控制方程 总场方程
• •
直接数值法求解总场 场源处具有奇异性（J→∞），场源附近需加密网格
韩波

二次场方程
其中：
• •
数值法求解二次场，解析法求解一次场，二者相加 选取恰当的背景模型来消除场源的奇异性，无需加密网格
•
背景模型的选取不方便，尤其对于多个场源位置的问题
•
频率测深
Ex Z xx E Z y yx
Z xy H x Z yy H y
H z Tzx
Hx Tzy H y
韩波
MCSEM工作原理
•
• • •
与典型的陆地CSEM如CSAMT相差巨大
特点：
•
•
仅利用了一阶导数（梯度）信息
每次迭代只需计算灵敏度矩阵或其转置与向量的乘积，无 需显式存储矩阵，计算快，节约内存
•
与利用二阶导数信息的方法相比，单次迭代收敛程度低， 因此迭代次数多，模型更新频繁，正演不宜采用直接解法
韩波
高斯牛顿-共轭梯度法（GN-CG）
GN法模型更新方式：
1 m k 1 m k H k gk
2015年中国地球科学联合学术年会
海洋电磁法三维正反演 研究进展
主要内容
1 2 3 4 5
海洋电磁法简介 研究的背景、目的与意义
三维正演算法研究进展 三维反演算法研究进展 结语
韩波
海洋电磁法简介
• 海洋电磁法（marine electromagnetic methods, MEM） 通过在海底布设电磁场接收机来记录天然存在的或人 工激发的交变电磁信号，以推断海底地层的电性情况，
韩波
矢量有限单元法
1
(a)
1 2 2 4 3 5 α
6
-1000 0 x/m 1000
• • •
基于电场总场Helmholtz方程 非结构四面体网格 容易模拟复杂的地电结构
韩波
1000
4
z/m
0
3
大型线性方程组的求解
控制偏微分方程经FV或FE离散化后都得到线性系统方程：
Ae b
A：大型、对称、稀疏复矩阵。
Nd
2
2
d m Wd d F(m) Wm m m r
2
2
韩波
最优化算法——线搜索方法
1
确定一个搜索方向 pk 确定沿 pk方向“行走”的距离，即步长 k 更新旧模型参数，获得新模型： mk 1 mk k pk
基本过程：
2 3
根据搜索方向的不同，或对同一类搜索方向的近似方式与近 似程度不同，基本的线搜索方法可分为：
地球物理技术，其原因在于地震勘探虽然能准确圈定构造，但难
以确定构造中是否含有油气，而海洋CSEM技术能识别高阻油气 储集层，从而提高钻井成功率，达到直接找油气的目的。
韩波
•
为了通过EM观测数据尽可能真实地还原出介质的电性结构以提升 勘探效果，三维反演解释是必需的。
Zhdanov等（2014）对挪威Troll海洋油气田进行了三维CSEM反演， 得到了与电阻率模型与地震数据反演模型很好地吻合
电偶极子或磁偶极子发射 收发距0至几十公里，非平面波 对高阻层反应灵敏
韩波
MSCEM中电磁波的传播路径
韩波
研究的背景、目的与意义
海洋电磁法的发展历程
● 起步阶段（ 1960s-1980s ） （Edwards et al., 1985; Cox et al., 1986; Edwards et al., 1988.） 试验阶段（ 1980s-2000s ） （胡文宝, 1991; Constable et al., 1996, 1998; MacGregor et al., 1998, 2000, 何继善, 1999）

主要缺点
• • • 内存需求巨大 计算速度慢 不易并行化，可并行程度有限
韩波
大范围改变空气电导率和频率时MUMPS直接求解器的精度和矩阵分解时间 （韩波等，2015）
频率（Hz）
空气层电导 率（S/m）
10-5 10-8 10-11 10-5
残差 9.10×10-21 9.09×10-21 8.99×10-21 9.36×10-21 9.58×10-21
实用化、快速发展阶段（ 2000s至今 ） （魏文博, 2001; Baba, 2005; Edwards, 2005; 何展翔, 2006, 2008; 沈金松, 2009; Constable, 2010; Key, 2012; 李予国, 2012, 2014）
●
பைடு நூலகம்
●
韩波
与MEM有关的文献出版情况. 引自Key 2012
韩波
Zhdanov等（2009）对墨西哥湾 Gemini勘探区的海洋MT数据进行
了三维反演，还原出了海底高阻盐
体结构的三维图像，与地震深度偏 移结果基本吻合。
韩波
电磁法的3D反演在最近10年取得了长足的进展，尤其是在MT方面， 发展出了多种性能优良的3D反演算法，已有学者公开代码供科学研 究使用。相较而言，CSEM的3D反演研究成果略少，国内的研究工 作更少，并且目前国内外没有任何公开的CSEM三维反演代码。再 者，由于地球物理反演所固有的非唯一性问题，没有任何一种反演 方法在所有情况下是都是准确的，不同反演方法的相互验证是削弱 非唯一性问题的有效途径。因此，在目前阶段，研究高精度、高效 率的海洋EM的3D反演算法与程序是一项非常有意义的工作。
J1 J J NT

J1 v Jv J v NT
韩波
单一广义场源的灵敏度的计算
韩波
交错网格有限体积法 I
• •

基于二次电场Helmholtz方程 规则六面体网格

s s p p E d i E d i ( ) E d 0 0
韩波
交错网格有限体积法 II
•
基于电场总场Helmholtz方程
•
规则六面体网格
r=1000 m 100 m (0,0,0) (0,0)
x z
(b) y
x
Hydrocarbon 2000 m
• 高阻层在x-y平面内的截面为圆形
韩波
(a)
(b)
-1000
0
1000
-1000
0 x/m
1000
1000
y/m 0
z/m
-1000
0 x/m
1000
• 矢量FE法 • 网格剖分：1,342,046个四面体， 1,562,622条棱边
多频率、多场源的灵敏度计算
• 定义“广义场源”：将频率视为场源，即将同一位置和方位的 发射源发射的不同频率视为不同的广义场源（Egbert, 2012） • 实际发射源个数Nt，第i个发射源发射的频率数Nfi，则广义场 源个数为：
NT Nf i
i 0
Nt
•
总体灵敏度由NT个广义场源的灵敏度合成：
相对残差 6.12×10-11 6.12×10-11 6.06×10-11 5.57×10-12 5.70×10-12
矩阵分解耗 时（s）
1122 1138 1133 1123 1146
0.01
0.1
10-8
10-11
10-5 1
9.37×10-21
7.20×10-21
5.57×10-12
4.58×10-13
主要优点
• • •
n 1, 2, ...
无需显式存储矩阵，内存需求小 计算速度快 容易并行化

主要缺点
• • 收敛性不稳定，取决于系数矩阵的性质 不利于求解多场源（右端项）问题
韩波
直接解法

主要优点
• • • 稳定、精度高，与系数矩阵的谱性质基本无关 多右端项问题，矩阵分解因子可重复利用 多个免费软件包：MUMPS，PARDISO，SuperLU等等