最小二乘最优化反演方法 -D数据和3-D高密度电法勘探反演

www.gaoxiao88.net

1.对称四极装置电阻率测深一维线性反演的基本思路

·第1步：建立目标函数

电阻率测深数据的反演方法可归结为使下面目标函数趋于极小：

地球物理数据处理基础

这里，NS是供电极距数；ρ_ai是第i个极距的实测视电阻率；M（M_j，j＝1，2…，NM）是模型参数（地层的电阻率或厚度），NM是预测模型参数个数；ρ_ci是由预测模型正演计算所得的第i个极距的理论视电阻率；α是范数，当α＝2时即为最小二乘法。

由于地球物理反问题存在多解性。因此在求解时，为减少反演的多解性，常在目标函数式（9－48）中对模型加入先验信息，如加入模型先验信息的目标函数变为

地球物理数据处理基础

式中：M_b为背景模型；C为光滑度矩阵；λ为拉格朗日系数。即要求预测模型正演结果与实测数据最接近，而且要求预测模型与给定的背景模型最接近并且光滑。显然拉格朗日

系数λ可以调控对预测模型的要求偏重于哪边。

·第2步：线性化

对式（9－48）利用泰勒展开，将模型的解在初始模型M处展开，并忽略二次及以上高次项，得

地球物理数据处理基础

要上式趋于极小，则对于各供电极距i（i＝1，2，…，NS）要满足下面的线性方程：

地球物理数据处理基础

·第3步：计算视电阻率对模型参数的偏导数得到偏导数矩阵。

·第4步：解线性方程组（9－50），得到预测模型M的修改量ΔM，建立新的预测模型。

·第5步：计算新模型的正演理论视电阻率，与实测视电阻率进行对比，如果精度满足要求，则新的预测模型即为反演结果；否则回到第3步重新计算偏导数矩阵，重新计算模型修改量，直到精度满足要求。

上述反演过程中有三个问题需要解决：①偏导数的求解；②模型参数的处理，模型参数中有不同量纲的层电阻率和厚度，尤其电阻率参数变化范围很大，如果直接由上述方法求解，不但会导致方程（9－50）严重病态，而且电阻率和厚度参数的修改也不会正确，从而会导致反演方法不收敛；③线性方程组（9－50）的求解，线性方程组（9－50）往往是不相容方程组，求解方法不同，收敛速度和收敛性也不同。下面将对这三个问题的处理方法逐一进行介绍。

2.偏导数的计算

对称四极装置电阻率测深的视电阻率对模型参数的偏导可以从数值滤波正演过程中直接给出。但一般反演中，偏导数通常用差分方法来计算。这里介绍用差分法来求对称四极装置电阻率测深的视电阻率对模型参数的偏导

取ΔM_j＝0.1M_j，则

地球物理数据处理基础

3.模型参数的无量纲化处理

为了解决视电阻率数据和电阻率模型参数变化范围大的问题，不少曾用它们的对数值来进行计算，但是这并没有解决层厚度参数与电阻率参数量纲不同的问题，因此，在反演过程中，出现了视电阻率测深曲线拟合很好，但预测模型与真实模型相差很大的情况，所以必须对方程（9－50）无量纲化。

将式（9－51）的偏导数代入方程（9－50）中可得

地球物理数据处理基础

其中：i＝1，2，…，NS。

上式两端同时除以ρ_ci（M₁，M₂，…，M_j，…，M_NM）可得

地球物理数据处理基础

这样方程组（9－52）可写为矩阵形式

地球物理数据处理基础

上式方程组中，左端系数矩阵A中各系数A_ij、未知向量x中各变量x_j以及右端向量B中各系数B_i都无量纲，从而解决了参数变量量纲不同的问题。

解线性方程组（9－54），则新模型参数Μ^*｛M^*_j，j＝1，2，…，NM｝为

地球物理数据处理基础

另外为了防止模型参数修改过量，实际过程中又作如下规定：x_j>1.5时取x_j＝1.5，x_j<－0.5时取x_j＝－0.5，即每次修改量不超过原有模型参数值的一半，以保证收敛稳定。

4.线性方程组求解

线性方程组（9－54）的求解问题，由于A为NS×NM阶矩阵，且对地球物理反演问题而言A多半是接近奇异的，所以不能用常规方法求解。解这个非方阵线性方程组最常用的方法有：

（1）奇异值分解法：它的基本思想是建立在奇异值分解定理上，即任意NS×NM阶矩阵A均可分解为A＝UWV^T，这里U为NS×NS阶正交阵和V为NM×NM阶正交阵，

地球物理数据处理基础

其中：δ₁，…，δ_r为矩阵A的奇异值，r是矩阵的秩。当A非奇异时，奇异值较大，方程组（9－54）有广义逆解X＝VW^－1U^T，这里

地球物理数据处理基础

当A接近奇异时，有的奇异值就较小，此时由于W^－1中系数过大，上述解的误差就较大。为了解决这个问题，维根斯（Wiggins）提出用最接近的矩阵R来代替A，R＝UW_eV^T，而

地球物理数据处理基础

W中小的奇异值在这里便被零代替了。因此有较精确的广义逆解X＝VW^－1_eU^T。

（2）最小二乘法：在方程组（9－54）两端左乘A^T，则方程变为A^TAx＝A^TB，当A非奇异时，有唯一解X＝（A^TA）^－1A^TB；当A接近奇异时，用上式求解误差就会很大，导致反演不收敛。为了解决这个问题，马奎特提出给系数矩阵A^TA的对角元素上加上一个常数，即A^TAx＋λI＝A^TB，以改善方程的条件，这就是著名的马奎特法又称阻尼最小二乘法。

分析奇异值（维根斯）分解方法和阻尼最小二乘（马奎特）方法，可以认为奇异值分解（维根斯）方法是将导致小的奇异值的方程取消来改善奇异值，而马奎特方法则是通过在A^TA的对角元素上加常数来增大A的奇异值，两者的思想刚好相反。

非线性解析反演与遗传算法的结合反演方法

周辉
（青岛海洋大学海洋地球科学学院，青岛　266003）
何樵登
（长春地质学院地球物理系，长春　130026）
摘要　各向异性介质参数反演通常为非线性优化问题。非线性反演方法可以分为两大类：随机搜索方法，如Monte Carlo法、模拟退火和遗传算法及基于非线性最小平方理论的非线性解析反演方法。遗传算法能寻找到全局最优解，但它为一种较费时的方法。非线性解析反演方法能给出一个与初始模型有关的局部最优解。然而，这种方法具有较快的收敛速度。遗传算法与非线性解析反演方法相结合的反演方法利用这两种反演方法的优点而克服其缺点。因此，结合反演方法既能快速收敛，又能寻找到全局最优解。如何合理地将遗传算法和非线性解析反演方法结合是十分重要的。本文提出一种结合方案，即在连续若干次遗传算法迭代后作一次非线性解析反演。理论算例表明结合反演方法具有上述特点。
关键词　遗传算法　非线性解析反演　非线性结合反演　各向异性介质
1　引言
遗传算法为随机搜索类方法之一，它以概率论为理论基础，用于求解多极值复杂优化问题[9]。遗传算法不要求已知模型空间中后验概率密度的形状并能广泛搜索模型空间。遗传算法模拟自然选择和遗传规律，并遵循适者生存的原则。
遗传算法由Holland在1975年提出[4]。Berg首先将遗传算法应用于地球物理优化问题[1]。Stoffa等系统地研究了种群大小、交叉概率、选择概率和变异概率对多参数优化问题收敛性和收敛速度的影响[11]。Sen等讨论了在选择概率中引入温度参数的作用并提出一些退火方案[10]。周辉等则研究了目标函数与收敛速度和解的精度的关系[16]。
基于最小平方优化理论的非线性反演方法是两大类反演方法之一。当给定的初始模型位于目标函数全局最优解所在的峰谷附近时，这种下降类方法能给出正确解而与初始模型位置无关。下降类算法研究得较深入，应用较广。
Tarantola提出一种基于广义最小二乘法的多维多偏移距声波地震波形解释的一般性非线性地震波形反演方法[12]。随后，Tarantola将该理论推广于各向同性介质的弹性波反演[13]。Gauthier等用理论数据验证了Tarantola提出的方法的正确性[2]。稍后，Tarantola研究非线性解析法反射波弹性反演的策略，指出以纵横波的波阻抗和密度作为反演参数，才尽可能使反演参数之间相互独立[14]。Pan用τ—P变换研究层状声学介质中平面波地震记录非线性解析反演的理论和可行性[6]。为了更多地利用地震数据中的信息，包括VSP资料中反射和转换信息，Mora作了一些工作[5]。当仅用反射数据时反演主要解决引起反射的P波和S波的波阻抗突变。当利用转换数据时，则能分辨大尺度的P波和S波速度变化。Sambridge等改进了修改模型的方法[8]。在子空间中，可同时得到P波、S波波阻抗和密度。周辉等将非线性梯度反演方法推广于多维、多道、多分量任意弹性各向异性介质参数的反演[17]。
非线性解析反演方法和遗传算法结合的反演方法利用非线性解析反演和遗传算法的优点，克服它们的缺点。因此，结合反演方法不仅能搜索到全局最优解，而且能较快地收敛。Porsani等在遗传算法和广义线性反演方法相结合方面作了一些研究[7]。
本文讨论各向异性介质的非线性解析反演方法和遗传算法与非线性解析反演方法相结合的结合反演方法[17]。对于遗传算法读者可参考遗传算法的相关文献[3,9～11]。
2　各向异性介质参数非线性解析反演方法
2.1　共轭梯度法
反演的目的是利用地面或井中测得的位移场ui（xr,t）求取地下介质密度分布ρ（x）和弹性参数分布Cijkl（x）。ρ（x）、Cijkl（x）称为模型参数。x为研究介质中或边界上任一点，x=（x1,x2,x3）,xr为接收点。反演的目标是使目标函数

岩石圈构造和深部作用

取极小值。其中Cd、Cm分别为数据（波场）和模型参数的协方差算子。m0为先验模型参数，m为反演过程中求得的模型参数。由于模型参数有多个，故用向量表示。ucal为给定m的波动方程正演记录，uobs为观测波场，上角标t表示转置。地震记录u和模型参数m之间的函数关系为

岩石圈构造和深部作用

g为非线性算子，（2）式为波动方程的算子形式。记第n次迭代时的模型参数为mn，则有

岩石圈构造和深部作用

及共轭梯度法的迭代公式[15]

岩石圈构造和深部作用

其中Gn为g对mn的Frechet导数，ηn为一常数，可由多种方法计算[5,8]。
梯度为模型空间的对偶空间中的一个元素。模型空间和其对偶空间以模型参数的协方差算子Cm=Diag（Cp，Cc）由式（4d）相联系。在后面将给出和的表达式。
式（4）为梯度反演方法的基本公式。当该公式中的每一量都已知时，迭代就可进行。在这些变量中，最关键的是梯度向量。
2.2　目标函数
在最小二乘理论中，权函数是协方差算子逆的核。假设数据集中的误差是不相关的，它仅取决于时间或源和接收器的位置，那么有[14]

岩石圈构造和深部作用

其中σ为数据的均方差。
2.3　各向异性介质中的弹性波动方程
令fi（x,t;xs）是第s次激发的内体力密度，Ti（x,t;xs）是地球表面S的应力矢量分量，ni（x）是表面的单位法向分量。那么与第s次激发相应的位移由以下微分方程组给出[15]

岩石圈构造和深部作用

2.4　梯度向量
式（4）中梯度向量的分量为[17]

岩石圈构造和深部作用

其中，T为地震记录的长度，为反向传播场，满足

岩石圈构造和深部作用

其中，t∈[T,0]，满足终了时间条件。
3　结合反演方法
3.1　遗传算法和非线性解析反演方法的优缺点
遗传算法是利用概率论来求解多极值复杂优化问题的一种随机搜索方法，由一组随机选取的模型开始，不需要更多的先验信息，广泛而有效地对模型空间的最优部分采样。尽管遗传算法是基于自然选择、遗传规律，搜索模型空间的最优部分而求得最优解，但它是一种计算量很大的方法。由于地震模型空间大，用全局最优化方法估计各向异性介质参数的地震波形反演十分费时。
目标函数的梯度信息是非线性解析反演方法修改模型参数的依据，它能给出一个接近初始模型的一个局部最优解。如果初始模型选择得合适，即当初始模型处在全局最优解所在的目标函数低谷时，非线性解析反演方法能收敛于全局最优解。然而，恰好给出一个接近全局最优解的初始模型的概率是非常小的，尤其对没有模型参数的任何先验信息的情况。但应强调的是，非线性解析反演方法具有较快的收敛速度。
发挥非线性解析反演方法快速收敛和遗传算法能搜索到全局最优解的优点，而克服前者仅能寻找到局部最优解和后者运算量大的缺点是很有意义的。非线性解析反演方法和遗传算法相结合的反演方法可达到上述目的。在结合反演方法中，遗传算法的作用是提供接近全局最优解的模型，非线性解析反演的作用是尽快求出全局最优解。因此，结合反演方法具有搜索到全局最优解的能力和比遗传算法收敛速度快的特点。
3.2　结合方案
遗传算法在优化过程中连续不断地搜索整个模型空间。在每次迭代结束后，得到一个本代的最优模型。根据遗传算法的数学原理[3]，最优模型的数量在下一代中得以增加，同时经交叉和变异作用又有新的模型产生。在下一代种群中，最优模型可能与前一代的相同，也有可能劣于前一代的最优模型。所有这些最优模型可能在目标函数的同一低谷处，也有可能在其它低谷处。遗传算法寻找最优模型要经过多次迭代才能确定一个极值。遗传算法的随机性导致遗传算法是一种费时的方法。然而正是遗传算法的这种随机性保证了它能搜索到全局最优解。
如果将每次遗传算法迭代的最优解作为非线性解析反演的初始模型，非线性解析反演可以找出与初始模型毗邻的局部最优解。由于非线性解析反演是一种确定性的方法，它按目标函数的梯度方向修改模型，所以非线性解析反演方法只需几次迭代即可收敛。非线性解析反演求得的解是否为全局最优解，非线性解析反演方法本身是无法得以保证的。只有当遗传算法提供接近全局最优解的初始模型时，非线性解析方法反演才能收敛到全局最优解。
结合反演方法中遗传算法和非线性解析反演方法的匹配方式是十分重要的。非线性解析反演方法得到接近遗传算法提供的初始模型的局部最优解后，在以后若干代中因遗传算法的随机性而使其最优解与该局部最优解相同。如果每次遗传算法迭代后作非线性解析反演，那么结合反演的结果在几代内都是相同的。显然其中的一些非线性解析反演是没有必要的。因此，结合方式应为在连续多次遗传算法迭代后作一次非线性解析反演，然后将非线性解析反演的结果作为下一代种群中的一个母本模型。图1为结合反演的框图。

图1　结合反演框图

4　算例
为了验证结合反演方法的优越性，对一维多层横向各向同性介质参数的反演理论实例作了分析。
图2是目标函数值与迭代次数的关系图。在该结合反演算例中每次遗传算法迭代后就作一次非线性解析反演迭代。结合反演的误差在开始几次迭代中下降很快，尤其在前3次。结合反演方法在第10次迭代达到的较小误差，遗传算法在第42次迭代才达到。结合反演的误差比遗传算法的跳跃得严重。这是因为非线性解析反演得到的模型在遗传算法中作为母代参加繁衍。这个模型因遗传算法的随机性常常被新的模型替代。这两个模型可能位于目标函数两个不同的低谷中，因此非线性解析反演的结果不同。
尽管结合反演的目标函数有些振荡，但也存在连续几次迭代目标函数几乎不变的现象。这意味着这几次迭代的最优模型是很接近的。在这种情况下非线性解析反演不能提供较大的改进。所以，此时的非线性解析反演是没有必要的，否则只能增加计算量。

图2　结合反演（实线）和遗传算法（虚线）的误差与迭代次数的关系

结合反演中每次遗传算法迭代后作一次非线性解析反演迭代
图3是另一个例子。在该结合反演例子中，每五次遗传算法迭代作一次非线性解析反演。在这里遗传算法占主要地位。此时结合反演的误差函数明显比遗传算法的小。结合反演的误差在第5次迭代末突然下降，并在第10次迭代时的小误差，遗传算法在42代才达到。遗传算法始终没有到达结合反演的最小误差。结合反演的误差在后期迭代过程中平稳下降，这是遗传算法占主导地位的原因。
从该例可知，若遗传算法与非线性解析反演方法比较合理地结合，结合反演方法比遗传算法具有快得多的收敛速度。
5　结论
非线性结合反演方法扬遗传算法和非线性解析反演方法之长，抑其之短，它是一种具有较快收敛速度的全局反演方法。
在结合反演中遗传算法和非线性解析反演方法的结合方式是重要的。从算例可得出，五次遗传算法迭代后作一次非线性解析反演的结合反演的效果明显优于每次遗传算法迭代后都作非线性解析反演的结合反演的效果。但是在结合反演中连续作多少次遗传算法迭代及连续迭代次数在整个迭代过程中的可变性还有待于进一步研究。

图3　结合反演（实线）和遗传算法（虚线）的误差与迭代次数的关系

结合反演中每五次遗传算法迭代后作一次非线性解析反演迭代
在结合反演中遗传算法的作用是提供接近全局最优解的初始模型。结合反演的运算速度主要取决于遗传算法的运算速度。均匀设计理论可以应用于遗传算法以加快随机搜索的速度。
与遗传算法相同，其它随机搜索方法也可用来与非线性解析反演方法形成结合反演方法。
参考文献
[1]E.Berg.Simple convergent genetic algorithm for inversion of multiparameter data.SEG60 Expanded Abstracts,1990，Ⅱ，1126～1128.
[2]O.Gauthier,J.Virieux and A.Tarantola.Two-dimensional nonlinear inversion of seismic waveforms:Numerical results.Geophysics,1986,51,1387～1403.
[3]D.E.Goldberg.Genetic Algorithms in Search,Optimiztion,and Machine Learning.Addison-Wesley,Reading,MA,1989.
[4]J.H.Holland.Adaptation in Natural and Artifical Systems.The University of Michigan Press,Ann Arbor,1975.
[5]P.Mora.2D elastic inversion of multi-offset seismic data.Geophysics,1988,52,2031～2050.
[6]G.S.Pan,R.A.Phinney,and R.I.Odom.Full-waveform inversion of plane-wave seismograms in stratified acoustic media:Theory and feasibility.Geophysics,1988,53,21～31.
[7]M.J.Porsani,P.L.Stoffa,M.K.Sen,et al..A combined Genetic and linear inversion algorithm for seismic wave-form inversion.SEG63 Expanded Abstracts,1993,692～695.
[8]M.S.Sambridge,A.Tatantola and Kennet.An alternative strategy for nonlinear inversion of seismic waveforms.Geophysical Prospecting,1991,39,723～736.
[9]M.Sambridge,and G.Drijkoningen.Genetic algorithms in seismic waveform inversion.Geophys.J.Int.,1992,109,323～342.
[10]M.K.Sen,P.L.Stoffa.Rapid sampling of model space using genetic algorithms:examples from seismic waveform inversion.Geophys.J.Int.,1992,109,323～342.
[11]P.L.Stoffa,M.K.Sen.Nonlinear multiparametre optimization using genetic algorithms:Inversion of plane-wave seismograms.Geophysics,1991,56,1794～1810.
[12]A.Tarantola.Inversion of seismic reflection data in the acoustic approximation.Geophysics,1984（a）,49,1259～1266.
[13]A.Tarantola.The seismic reflection inverse problem.In:F.Santosa,Y.-H.Pao,W.W.System,and C.Holland Eds.Inverse problems of acoustic and elastic waves.Soc.Industr.Appl.Math.,1984（b）,104～181.
[14]A.Tarantola.A strategy for nonlinear elastic inversion of seismic reflection data.Geophysics,1986,51,1893～1903.
[15]A.Tarantola.Inverse problem theory:Methods for data fitting and model parameter estimation.Elsevier Science Publ.Co.Inc.,1987.
[16]周辉，何樵登.遗传算法在各向异性介质参数反演中的应用.长春地质学院学报，1995,25，增刊1,62～67.
[17]周辉.各向异性介质波动方程正演及其非线性反演方法研究.长春地质学院博士论文，1995.

5.3.1 3-D高密度电法勘探数据
在非学术圈，有时会出现两个有趣的问题:定义一个3-D反演方法和3-D数据由什么构成，第一个问题可以很容易回答，但第二个问题就不一定能说得很清楚了。
一个3-D反演方法的特点是允许模型电阻率值在所有3个方向（即x，y，z方向）变化。而2-D反演假设地下电阻率值只在x和z方向变化，y方向保持不变；1-D反演电阻率值仅在z方向变化。3-D反演模型使用独立的矩形单元（图5.18），尤其允许模型值在所有3个方向不同，因此，采用了真3-D反演方法。请注意，一系列平行的2-D测线反演出的模型结构不是真正的3-D反演模型。

图5.18 3-D反演模型a—标准模型，矩形单元格在x和y方向上的宽度等于单位电极距；b—模型前几层在纵、横向上对半细分，目的是提供更好的分辨率；c—模型单元在水平方向上精细划分，但在垂直方向不细分

另一种3-D反演方法的定义特征是用一个3-D正演子程序，如3-D有限差分和有限元（Dey et al.，1979b；Silvester et al.，1990）计算模型视电阻率和雅克比（Jacobian）矩阵值。棘手的问题是“3-D” 数据由什么构成，目前还没有一个普遍能接受的定义，虽然反演数据的方法是3-D，但是，数据是否包含有意义的3-D信息则是另一回事，有人说，为了对地下获得一个好的2-D模型，2-D覆盖的数据必须质量很高，对于2-D勘探来说，可以通过不同电极距观测和不同水平位置观测来获得这样一个2-D覆盖面，但是，前一套数据所需的数据覆盖程度可能考虑 “3 -D” 情况缺少清晰度，以单极-单极（pole-pole）装置勘探为例，下边提出了3-D数据分类系统，按3-D信息内容减少的顺序列出。
第1类:理想的3-D勘探是电极布设成矩形测网，在所有可能的方向进行观测（见图5.2a），即沿不同角度的网格线。
第2类:电极布设成矩形网，所有的观测沿测网线（即在x和y方向），仅在测网某个角度方向进行有限次观测（如45°对角线的方格），如图5.2b所示。
第3类:观测仅沿测网线的两个方向进行，即x和y方向，在测网线之间的夹角不进行观测。这种观测方法通常用在某一时间多通道电极系统没有足够的通道来覆盖整个勘探区域，一种可能的观测序列如图5.13所示。
第4类:观测只在一个方向（如x方向）沿着一系列平行的2-D测线进行，尤其金属矿勘探，这种情况在以往的勘探中是比较常见的。对于这种类型的数据，通常先进行2-D反演，然后组合所有测线的数据进行3-D反演，其目的是从已有的数据中获取新的信息，看看3-D效果是否更显著（即2-D反演结果是否有效）。这种3-D反演效果是否成功，在某种程度上依赖于线距和装置类型，一般来说，线距不得小于2倍测线内的单位电极距。
理想的数据应该是第一类或者至少是一些角度数据可用的第二类，从第3类、第4类数据得到的3-D模型精度会比第1类、第2类数据得到的模型低，在很大程度上依赖于线距和装置类型。但是，即使是类型4的数据，3-D反演结果对有意义的3-D效果也提供了有用的指示，尤其可对不同测线获得的单独2-D反演结果提供有效的检查。
5.3.2 3-D高密度电法数据反演
3-D数据反演可以用类似于2-D反演的光滑约束最小平方反演方法（deGroot-Hedlinet al.，1990；Sasaki，1992），该技术比传统最小平方反演方法快10倍，对于大型数据来说，可以有效节省内存。该方法的另一个优点在于，阻尼因子和光滑滤波系数可适应不同的数据类型，可根据具体情况来调整。光滑约束最小平方法基于下面的方程:

高密度电法勘探方法与技术

式中: 为水平光滑滤波系数矩阵；fz为垂直光滑滤波系数矩阵；J为偏导数矩阵；JT为J的转置矩阵；u为阻尼因子；D为模型参数修改矢量；G为残差矢量。
解该方程可以用拟牛顿法（quasi-Newton）最小二乘法优化技术（Loke et al.，1996），该技术对于大型数据来说，运算速度比传统最小二乘法快10倍，同时需要的内存更少；也可以用高斯-牛顿（Gauss-Newton）法，它的运算速度远慢于拟牛顿法（quasi-Newton）法，但在电阻率差异比超过10:1的地区，该方法可能得到一个比较好的结果；第三种选择方法是先利用高斯-牛顿（Gauss-Newton）法进行2～3次迭代，然后再用拟牛顿法（quasi-Newton）法，在许多情况下，这是一个最好的解决方案（Loke et al.，2002）。除了光滑约束方法外，还有其他一些反演方法也可能是可靠的，反演方法的选择应遵循勘探区域已知的地质情况来决定。
反演时，将地下细分成若干个小矩形棱柱体，并试图确定该棱柱体的电阻值，以减少计算出的视电阻率与观测值之间的差异。Loke和Barker（1996）给出了一种可能的布设，如图5.18a所示，顶层每块的角点处有一根电极。除了这种基本的布设方式外，其他两种布设方式也遵从这种方案，一种方式是前几层在纵、横向上对半细分（图5.18b）；另一种模型单元仅在水平方向上对半细分（图5.18c）。因为电阻率分辨率随深度快速下降，通过研究发现，细分地下模型块仅对顶部两层有利，在许多情况下，仅细分顶层就足够了。通过细分单元格，模型参数的数量和计算机反演数据所需的时间将大大增加。
最优化方法试图通过调整模型块的电阻率来减小计算视电阻率值与测量视电阻率值之间的差异，这一差异通过均方根误差（RMS）来衡量。但是，在模型中，有时最小的模型RMS可能会出现很大和不切实际的模型电阻率值变化，但是，从地质角度来看也不是 “最佳” 的模型。一般来说，最谨慎的做法是选择迭代之后均方根误差不发生明显变化的模型，这种模型通常在4～5次迭代后出现。
5.3.3 电阻率固定
在某些情况下，地下剖面的真正电阻率可能会知道，如从钻孔电阻率测量获得，反演时，允许地下一些区域固定电阻率值，这些区域的形状最好固定为正方体或长方体，区域的形状也应该确定，对于矩形区域，给出顶部左后角x，y，z坐标和底部右前角，y，z坐标，如图5.19所示。然后是固定区域电阻率的阻尼因子权重，该参数可控制反演子程序改变区域电阻率，通常情况，该区域的一些电阻率存在不确定性，钻孔测量只能给出井孔附近有限区域的电阻率，因此，建议可以（有限地）改变区域电阻率。如果使用阻尼因子权重为1.0，该区域的电阻率允许与地下模型其他区域的电阻率改变相同；如果采用较大阻尼因子权重，则固定区域电阻率允许改变较小，如10.0，在反演过程中，该区域的电阻率变化是非常小的，这样大的值仅用电阻率和形状准确知道的区域。

图5.19 一个矩形区域电阻率固定的反演模型

搞笑吧

最小二乘最优化反演方法 -D数据和3-D高密度电法勘探反演

相关参考：

相关评论

相关话题