基于非负最小二乘法的一维MT正则化反演研究
周绍民;柳建新;孙欢乐
【摘 要】In order to solve theill-posed,instability and uniqueness in magnetotelluric magnetotelluric regularized inversion algorithmbased onnon-negative least squares.The ill-posed problem has been effectively improved by introducing regularization thought and the model constraint object function u the smooth constrained regularization matrix.In order to avoid the logarithm of resistivity,the fast non negative least squares method is introduced to simplify the calculation of partial derivative matrix.Finally stable inversion result the different geo models it shows that this inversion algorithm is effective.%为了有效地解决在大地电磁反演过程中出现的不适定性、解的不稳定性和非唯一性等问题,提出了基于非负最小二乘法的MT正则化反演算法,该算法的思想是通过引入正则化思想,模型约束目标函数被采用最光滑模型约束正则化矩阵,有效地改善了不适定问题得到有效改善.为了避免常规方法需要对电阻率取对数的情况,引入了快速非负最小二乘法(这里把反演参数设为电阻率,因而必须为正值),从而简化了反演过程中偏导数矩阵的计算.最后试算不同的地电模型,得到稳定的反演结果,表明了此反演算法是有效的.
【期刊名称】《工程地球物理学报》
【年(卷),期】2017(014)003
【总页数】9页(P253-261)
【关键词】大地电磁测深法;非负最小二乘法;正则化因子;反演
【作 者】正则化最小二乘问题周绍民;柳建新;孙欢乐
【作者单位】中南大学 地球科学与信息物理学院,湖南 长沙 410083;中南大学 有资源与地质灾害探查湖南省重点实验室,湖南 长沙 410083;中南大学 地球科学与信息物理学院,湖南 长沙 410083;中南大学 有资源与地质灾害探查湖南省重点实验室,湖南 长沙 410083;中南大学 地球科学与信息物理学院,湖南 长沙 410083;中南大学 有资源与地质灾害探查湖南省重点实验室,湖南 长沙 410083
【正文语种】中 文
【中图分类】P631.3
大地电磁测深法(Magnetotelluric,简称MT)是由前苏联学者Tikhonov[1](1950)和法国学者Cagniard[2](1953)提出来的,通过观测和研究天然的交变电磁场,从而得到地球内部电性结构信息的一种地球物理勘探方法。由于大地电磁测深法是以天然场源为研究基础,不需要外加场源,因此野外勘探工作容易开展,同时能够降低野外勘探的成本;天然场源信号不会被高阻层所屏蔽,对于低阻层具有较高的分辨率[3];根据趋肤效应准则,其趋肤深度仅仅与地下介质的电阻率和电磁场的频率有关。在岩石导电性一定的条件下,频率越高,对应的勘探深度越浅,频率越低,对应的勘探深度越大,这就构成了大地电磁测深的工作基础[4];其探测深度的能力是其他地球物理方法所不能比拟的,因此得到了广泛的应用。MT正演数值模拟主要采用有限单元法[5]、有限差分法[6-7]和积分方程法[8],这三种数值模拟方法都有自己独特的优势与不足,有限单元法应用较广泛。MT反演方面,早期采用的是一维渐近线近似法和手工量板法,目前应用较广泛的有OCCAM一维反演[9]、一维连续介质反演的曲线对比法[10,11]和ARIA一维反演[12]。
大地电磁测深反演问题都是不适定问题,因此在反演过程中,得到的解具有不稳定性和非唯一性,地球物理工作者通常在反演过程中引入正则化思想,能够有效地解决这一类问题。在大地电磁测深反演计算过程中,对模型施加必要的先验约束条件,解的稳定性得到很好的改
善,同时可以减少反演结果的非唯一性,这是目前MT反演研究者最常采用的,也是最为有效的方法。根据不同的模型约束条件类型,比较常用的模型约束目标函数有三种,分别是最小模型约束、最平缓模型约束和最光滑模型约束。在反演求解方程组过程中,为了消除在反演过程中电阻率出现负值,常规的处理思路是对模型参数取以10为底的对数,而在本论文中引入了数学理论方法非负最小二乘法,反演模型参数直接取为电阻率参数,提供了一种新的处理方法思路,对于不满足要求的负值解,直接将其略去,从而加快了方程组的求解速度,反演效率和精度得到了有效的提高。
在本文中,对MT一维问题分别给出了解析解和数值模拟结果,并进行了对比,证明了数值模拟方法在地球物理中的可行性;采用Bostick[13]反演结果作为反演的初始模型,反演结果是直接由观测数据计算得到的,没有进行数据的优化过程,其反演结果的精度虽然不高,但反演计算过程简便,其解释具有唯一性,不存在人为的干扰,能较好地反映地电断面的基本特性,因此能为大地电磁测深资料的精确反演提供可靠的初始模型;针对不同的一维层状介质模型,得到满足精度的反演结果,最后对三种不同的模型约束目标函数的反演结果进行对比研究。
2.1 正演理论
在大地电磁测深方法研究中,传导电流占支配地位,因此忽略了位移电流的影响,即研究的是准静态极限问题[14],假设时间因子为eiωt。于是,可以得到Maxwell方程组可表示为:
式中,E和H分别为电场(V/m)和磁场分量(A/m);ω为角频率(rad/s);μ(μ0=4π×10-7H/m)均匀半空间磁导率(B/H);σ为介质的电导率(μs/m)。
将上述方程组展开成分量形式,对于一维层状模型,可得到视电阻率ρ和相位φ的公式表示为:
φ=arctan
2.2 有限单元法分析
电场所满足的偏微分方程为
+iωμσEx=0
给出边界条件,其边值问题所对应的变分问题为
δF(Ex)=0
对研究区域进行剖分时,为了实现无穷远边界条件的近似处理,在保证计算精度的前提下,采用非均匀网格剖分,同时节约了计算时间。
利用矩形单元双二次差值构造相应的形函数[15],求解上述变分问题,最后将得到一个复系数方程组:
KEx=0
其中,K为总体刚度矩阵。将电场所满足的上边界条件和下边界条件代入到方程组当中,解方程组可得节点处的电场值,再利用数值方法求出电场值沿垂向的偏导数值,从而可以计算出模型响应的视电阻率和阻抗相位。
为了验证有限单元法模拟大地电磁响应的可行性,选用了一个电阻率为20 Ω·m的均匀半空间,分别得到了解析解和FEM数值解。通过图1分析比较得出,两者拟合一致,因此有限单元法在模拟计算大地电磁响应理论上是可行的。
大地电磁反演问题可以概括为已知观测数据求对应模型参数的过程[16],并且是不适定问题,得到不稳定的反演结果,同时反演解释具有非唯一性,也就是说同一组观测数据对不同
的地电模型具有一致的拟合效果,也就是出现多解性问题。因此为了改善和解决这一问题,通常是引入正则化思想
Pα(m)=φ(m)+αs(m)
其中,Pα(m)为总体目标函数;φ(m)为数据目标函数;s(m)为模型约束目标函数;α为正则化因子。
因此,反演问题的总目标函数可表示为
Pα(m)=[do bs-F(m)]2+αWm(m-mref)2
式中:do bs为观测值;F(m)为正演响应函数;Wm为光滑度矩阵,也称为模型权系数矩阵;mref 为先验模型。
设mk是模型的第k次迭代结果,那么基于先验模型的正则化反演迭代更新公式可写为mk+1=mk+Δm。
同时设Jk为Jacob矩阵,写为将上述总目标函数对Δm求导并令其等于零,通过解方程组得到
模型修正量Δm,将其加到预测模型参数矢量中,得到新的模型参数矢量;如此重复这一过程,直到总体目标函数满足要求为止,得到迭代形式为

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。