模拟退火法在地球物理反演中的应用-张明浩

模拟退火法在地球物理反演中的应用简述

王占文（中国矿业大学北京，地球科学与测绘工程学院） 摘要：模拟退火法源于统计热力物理学，它模拟熔融状态下物体缓慢冷却达到结晶状态的物理过程。模拟退火反演算法的基本思想是：生成一系列参数向量模拟粒子的热运动，通过缓慢地减小一个模拟温度的控制参数，使模拟的系统最终冷却结晶达到系统能量最小值的过程。模拟退火反演算法实质是利用了地球物理反演问题求解过程与熔化固体退火过程的相似性，开辟了地球物理反演的新途径，是非线性反演算法中一种最常用的算法。 关键词：模拟退火法；地球物理反演；蒙特卡洛法；

一、 前言

地球物理学中,模拟退火法的使用与高度多解性的资料解释有关。模拟退火不仅能提供很好的拟合模型,而 且也能提供有地质意义的额外信息。例如:可以解释资料的那些模型能有什么样的差异？那些可用的模型是什么样子的？它们的相对概率是多少？模拟退火法（Simulated Annealing Method，简称SA），源于统计热力物理学，是模拟物质退火的物理过程，即统计试验的组合优化过程。模拟退火算法的基本思想是: 生成一系列参数向量模拟粒子的热运动，通过缓慢地减小一个模拟温度的控制参数， 使模拟的热系统最终冷却结晶达到系统能量最小值。模拟退火算法的思想最早由Metropolis等人于1953年提出。Kirkpatrik等人研究证实，冷却后物质变成晶体时的状态，恰好是系统处于整体能量极小的状态，并于1983年成功把它应用于组合优化问题。基于这种类比，Rothman在1985年将退火原理引入到地震资料处理中的剩余静校正问题中，并称之为模拟退火法。紧接着Sen和Stoffa等将其用于一维地震波形反演、一维电测深反演。随后，模拟退火法在大地电磁测深、重磁资料等方法中。

模拟退火反演算法与传统的线性反演方法相比该方法具有：不依赖初始模型的选择、能寻找全局最小点而不陷入局部极小、在反演过程中不用计算雅克比偏导数矩阵等优点。所以近年来，模拟退火法在地球物理资料的反演中得到广泛应用与研究，并且取得了一些令人瞩目的研究成果，它已经成为了目前非线性反演的重要方法之一。

本文主要讨论了自模拟退火法提出以来的一些研究进展和成果，以及现阶段模拟退火法在地球物理反演中的应用情况，所参考文献有26篇，起止年份自1953年至2012年。

二、 模拟退火法的原理

模拟退火算法( Simulated Annealing)源于统计物理学，据统计热力学，物体中的每个分子的状态服从Gibbs分布， 即：

式中: E(ri)为第i个分子的能量函数；；ri为第i个分子所处的状态； k为波耳兹曼常数， T表示温度；(ri)为第i个分子的概率密度，为了方便起见令k= 1。

算法源于对实际固体退火过程的模拟， 即先将固体加温至充分高， 再逐渐冷却。加温时， 固体内部粒子变为无序状态， 内能增大; 而逐渐降温时， 粒子趋于有序， 在每个温度都达到平衡态， 最后在常温时达到基态， 内能减为最小。因此， 算法实际上是将优化问题类比为退火过程中能量的最低状态， 也就是温度达到最低点时， 概率分布中具有最大概率的状态。

模拟退火算法的具体步骤如下:

1) 给定模型每一个参数变化范围, 在这个范围内随机选择一个初始模型

m0 ，并计算相应的目标函数值 E(m0)；

2) 对当前模型进行扰动产生一个新模型 m,计算相应的目标函数值

E(m) , 得到E= E(m) - E(m0)；

3) 若 E< 0, 则新模型 m 被接受; 若 E > 0,则新模型 m 按概率 P = ex

p( - E/ T ) 进行接受,T 为温度。 当模型被接受时, 置m0 = m；

4) 在温度 T 下,重复一定次数的扰动和接受过程, 即重复步骤 2)、3)； 5) 缓慢降低温度 T；

6) 重复步骤 2) 、5) , 直至收敛条件满足为止。为了避免最好的解在优化过程中被忽略掉,可以稍做改进, 即在整个搜索过程中随时记下最好的值, 因为该法的特点决定了最后的优化值并不一定就是最优的值。

三、 模拟退火法的发展和应用

在模拟退火法的应用方面，上世纪90年代出现了一批新技术新应用，但是均处于起步阶段，仍有很多地方并不成熟和完善。和一般地球物理反演问题一样，Rothman 首先将静校正问题转换为目标函数的极小化问题，即:minE ( x) 式中: E 是共中心点道集的叠加能谱，x 为任意炮点到接收点的静校量所组成的矢量。应该看到，由于炮点和接收点非常多，求解一个全局最优的静校正矢量并不是一件非常容易的事情。Rothman引入了Kirkpatrick 等提出的模拟退火法对上述问题进行求解。

Sen和Stoffa尝试将模拟退火方法用于一维地震波形反演。一维地震波形反演问题是要求一个最优的模型，使该模型正演的地震波形与实际地震波形的拟合度最大，即求这样一个模型矢量m，使观测波形和理论波形的拟合度E最大。为了加快对模型空间的搜索效率，Sen 和Stoffa采用了热浴法( Heat - bath) ，即在某个温度上，只变化一个参数而固定其它参数，使该模型参数从最小值依次搜索到最大值，计算不同向量对应的Gibbs条

件概率

根据该概率值对该模型参数进行寻优。然后变化另一个参数而固定其它参数，使另一个参数也从最小值依次搜索到最大值，根据上式计算 Gibbs概率进行寻优，依次类推。整个参数寻优一次后，将温度降温，然后再重复刚才模型寻优的过程。随着温度的降低，反演参数将逐渐向最优值演化。

Paulo M. Carvalho和Arthur B. Weglein在1994年提出可以使用模拟退火算法为界面多次波衰减提供子波估计，且可提供良好的结果。该方法对于求取资料的比例因子，甚至当它变化很大时特别有效。这里，它将涉及早期在这种情况下试图使用梯度最佳化法的失败。应用与模拟退火法有关的梯度尚需进一步研究，有可能改善最终估计。在叠后层速度计算上，也有学者提出了可以利用模拟退火法，并且在早期的PE机上得到了实现。但是这种方法所得到的反演解的分辨率不高，直接解决地质问题是不现实的。不过模拟退火法得到的反演解是正确的，且比初始地质模型更加真实细致的反应底下地质的平均效应。所以可将模拟退火反演解作为共扼梯度法的初始解，获得了可提供更丰富信息并能解决地质问题的最终解。

虽然在地球物理非线性反演中，模拟退火法是一种较为先进的启发式蒙特卡洛法，但是实际地质资料处理中，该方法存在着计算效率不高，甚至有时会失效的弊端。因此，诸多学者提出了对模拟退火法的改进。Rothman，Sen和Stoffa提出的用试验方法估求Tmin。需对系统作若干次迭代运算，然后才能得到合适的Tmin的估计值。姚姚提出由部分局部势能估求Tc*（及Tmin）的方法不仅避免了盲目性，计算效率也大为提高。因为计算局部势能只是系统能量计算中的很少一部分，即使全部局部势能计算