刘　庆，徐凤姣，Osborne Kachaje

(长江大学 油气资源与勘探技术教育部重点实验室，湖北 武汉　430100)

20世纪50年代初Tikhonov(1950)[1]和Cagnird(1953)[2]分别提出利用天然电磁场进行勘探的方法。

有限单元法由于精度高，能够模拟复杂地形，因此在国内外得到了广泛应用。徐世浙[3-4]对有限单元法进行了深入研究，为国内学者使用有限单元法进行MT的正演奠定了理论基础。

刘鹏茂[5]、顾关文[6]等基于CPU实现了大地电磁的二维与三维反演，取得了一定的加速效果。本文采用CPU(OpenMP)和GPU(CUDA)异构并行的方式进行MT的二维正演模拟，并采用全稀疏存储的方式来存储有限元中的刚度矩阵，有效地减少了单次MT二维正演所需的时间和对内存的占用。

1　MT二维正演

根据麦克斯韦方程，角频率为ω(时谐因子为e-iωt)的定态电磁场的方程是

▽×E=iwμH

(1)

▽×H=(σ-iωε)E

(2)

假定电性结构是二维的，取走向为y轴，x轴与y轴垂直，保持水平，z轴在垂直方向上。可知∂u/∂y=0，将(1)、(2)展开可得到两个独立的方程组，即E型波和H型波，方程组的具体形式如下。

E型：

(3)

H型：

(4)

本文选用H型波的方程，并基于有限单元法进行大地电磁二维正演。其边值问题可以归纳为：

(5)

边值问题(5)与下列变分问题等价：

(6)

2　有限单元法

利用4节点的等参单元来进行数值模拟，采用双线性插值，将式(6)中的区域积分分解为各单元积分之和：

(7)

单元积分：

(8)

(9)

单元积分：

(10)

(11)

单元积分：

(12)

(13)

(14)

δF(u)=δuTKu=0

(15)

由δu得任意性，可得Ku=0。解线性方程组前，将AB线上得边界值代入。解得线性代数方程组后，得各节点得u。

3　稀疏存储

本文采用三元组的存储方式来储存总体刚度矩阵K。并分别采用全稀疏存储和不做稀疏处理的存储对不同网格大小的二维区域进行MT二维正演，所得的刚度矩阵所占用的内存对比如图1所示。

图1　完全存储与稀疏存储对比

由图1可知，采用完全存储的方式时，刚度矩阵中的元素将随着网格节点的个数以二次方的量级增长，而全稀疏存储的方式则对内存的需求较小。

4　基于CPU和GPU的并行计算

本文采用CPU(OpenMP)和GPU(CUDA)异构并行的方式来进行MT二维正演计算，加速效果如表1所示。

表1　MT二维正演并行效率对比

5　MT二维并行正演结果分析

在二维情况下，设计了如图2所示的二维正演模型，图3为其二维正演视电阻率断面图。

图2　二维正演模型

图3　正演结果

从图3中可以看出正演的视电阻率断面图可以很好的对正演模型做出响应，进一步地验证了MT二维并行正演程序的正确性。

6　结　语

1)采用全稀疏存储的方式避免了由于存储刚度矩阵中的零元素而造成的内存的浪费，为模拟范围更大、网格更多、精度更高的MT二维正演提供了便利。

2)采用CPU和GPU异构并行的方式较明显地提高了单次大地电磁二维正演的速度。