二维均匀介质中积分方程的快速精确算法

采用稳定型双共轭梯度快速Fourier变换(BCGS-FFT)算法精确计算二维均匀介质中的积分方程.采用一种新的插值函数作为基函数和试探函数对积分方程进行弱化离散,离散后的积分方程采用稳定型双共轭梯度迭代方法进行求解,从而得到异常体内电场的分布.计算时采用快速Fourier变换技术将积分方程内Green函数与电场的乘积表示成褶积形式以加快计算速度.数值计算举例说明了算法的精确性和有效性.     

基于改进Krylov子空间算法的井中激电反演

井中激电是二次找矿重要的地球物理勘探手段,快速而稳定的正反演算法有助于方法的推广和应用.本文在正演模拟中,给出了考虑井眼影响下的网格剖分方式;用右端项校正技术减小边界效应和源点奇异性引起的模拟误差;并采用循环Krylov子空间算法提高多线性方程组的求解效率.反演用Gauss-Newton法结合Jacobian-free Krylov迭代求解技术,给出了Jacobian矩阵向量积的简化计算方法;用不精确预处理共轭梯度法对模型修正量方程近似求解以减少计算量;采用不同于正演的反演网格剖分降低不适定性.数值算例验证了相关算法的有效性和可靠性.     

二维介质柱的电磁成像研究

从体等效原理出发得到用介质柱介电参数表示的体等效电流,并得到一个积分方程组求解散射问题．利用变分原理,对此方程组求变分,并化为矩阵方程组,消去含介质柱内总场变化的项,得到散射场变化与介质柱介电参数变化之间关系的矩阵方程．即反演方程．用算例检验此法的反演精度、收敛性能、多目标反演能力、抗噪声性能及利用不完备测量散射信息反演目标的能力．结果表明,此法是一种实用性很强的二维介质柱电磁成像方法．     

一种随机粗糙海面电磁散射的极化核函数

从磁场的积分方程出发,应用迭代方法推导了二次迭代条件下粗糙海面的近似电磁散射极化核函数,并对其进行了分析.这种函数的表示是在光学近似即KA近似的基础上的一种扩展形式,因而可称为扩展KA近似(EKA).Elfouhaily采用了与Holliday类似的方法分析了双站配置条件下散射的极化核函数,试图保留所有的高次项,但是由于处理的复杂性他只给出了斜率项在散射项中的表示形式,并称之为KA近似的"补充项".本文在此基础上进一步考虑了Elfouhaily近似的形式,给出了表示双站散射极化核函数的完整表达形式,并以此为基础推算了后向散射条件下极化核函数的理论结果,并与Elfouhaily的结果进行了对比,说明了极化核函数完整表达形式的合理性.近似散射极化核函数的分析对于深入理解粗糙表面(特别是海面)电磁散射特性以及散射近似模型的建立等问题的研究具有重要意义.     

基于Born敏感核函数的VTI介质多参数全波形反演

本文基于VTI介质拟声波方程,利用散射积分原理,在Born近似下导出了速度与各向异性参数的敏感核函数,同时结合作者前期研究提出的矩阵分解算法实现了一种新的VTI介质多参数全波形反演方法.矩阵分解算法通过对核函数-向量乘进行具有明确物理含义的向量-标量乘分解累加运算实现目标函数一阶方向或二阶方向的直接求取,从而避免了庞大核函数矩阵与Hessian矩阵的存储,该方法同时可以大大降低常规全波形反演在计算二阶方向时的庞大计算量.为了克服不同参数对波场影响程度的不同,本文利用作者前期在VTI介质射线走时层析成像研究中提出的分步反演策略实现了多参数联合全波形反演.理论模型实验表明,本文提出的基于Born敏感核函数的各向异性矩阵分解全波形反演方法可以获得较好的多参数反演结果.     

正则化预处理迭代算法在频率域声波模拟中的应用

高精度及高效频率域声波数值模拟的关键在于高效求解声波方程经离散化后得到的大型稀疏线性方程组.该方程组系数矩阵具有很强的稀疏性,非对称性和非正定性等特征,常用的迭代算法难以准确、高效地求解.为了改善数值模拟迭代算法的收敛性与稳定性,在算法基础上添加预条件算子是求解该类方程的常用方案.本文基于以上思路,引入正则化技术来构造合适的预条件算子,提出正则化预条件迭代算法,以加速求解方程组.通过包含有均匀介质和高非均匀度介质(Marmousi)模型的数值模拟实验结果表明:与单独使用迭代算法相比,本文提出的正则化预条件迭代算法在计算量方面仅多了一次矩阵-矢量相乘,内存消耗未增加;同时,基于该算法的数值模拟结果能够满足精度要求,较单独使用迭代法能够有效改善收敛性质,加快收敛速度;而且,在二维模型算例下,与LU分解算法相比,基于该算法的内存消耗大幅下降.     

用变分玻恩迭代方法重建二维非均匀介质结构

提出了用于二维轴对称非均匀介质结构的反演和成像的一种新的反演迭代方法──变分玻恩迭代方法（ＶＢＩＭ）．首先利用玻恩近似将非线性积分方程线性化,然后应用变分方法导出用于反演的电场积分方程．正演数据则利用高效的数值模式匹配方法获得．数值结果表明,ＶＢＩＭ与ＢＩＭ相比,其收敛速度、成像质量等均得到较大的改善。     

位场曲化平积分方程的迭代解

提出了位场曲化平的新方法. 给定观测曲面S上的位场、S对下方水平面P的相对高程,确定P上的位场. 利用由P向上延拓到S的积分式,建立这两个面上位场及相对高程三者所满足的方程,它是第一类Fredholm积分方程. 用Fourier逆变换式把这一空间域积分式化为波数域积分式,再由指数函数的Taylor展开进一步化为级数式. 积分方程的解采用逐次逼近法迭代计算,即用S上的位场观测值作为P上位场的初始迭代值,用导出的级数式求得S上的位场计算值、由S上的位场观测值与计算值之差校正P上的位场,多次迭代,直到满足迭代终止准则. 我们还给出该积分方程的波数域迭代计算方法. 模型算例表明,重力异常曲化平的均方差和磁异常曲化平的均方差分别为0.0008 mGal和0.0019 nT,在主频为2.26 GHz的笔记本电脑运行,2048×2048数据量,计算时间是975 s. 野外磁场实际资料处理也证实这种方法的有效性.     

用变分坡恩迭代方法重建二维非均匀介质结构

提出了用于二维轴对称非均匀介质结构的反演和成像的一种新的反演迭代方法－变分坡恩迭代方法（ＶＢＩＭ）,首先利用玻恩近似将非线民生积分方程线性化,然后应用变分方法导出用于反演的电场积分方程,正演数据则利用高效的数值模式匹配方法获得,数值结果表明,ＶＢＩＭ与ＢＩＭ相比,其收敛速度,成像质量等均得到较大的改善。     

基于高斯牛顿法的频率域可控源电磁三维反演研究

三维反演解释是电磁法勘探发展的重要趋势,而如何提高三维反演的可靠性、稳定性和计算效率是算法开发者们目前的研究重点.本文实现了一种频率域可控源电磁(CSEM)三维反演算法.其中正演基于拟态有限体积法离散化,利用直接矩阵分解技术来求解大型线性系统方程,不仅准确、稳定,而且特别有利于含有大量发射场源位置的CSEM勘探情况;对目标函数的最优化采用高斯牛顿法(GN),具有近似二次的收敛性;使用预条件共轭梯度法(PCG)求解每次GN迭代所得到的法方程,避免了显式求解和存储灵敏度矩阵,减小了计算量.以上这些方法的结合应用,使得本文的三维反演算法准确、稳定且高效.通过陆地和海洋CSEM勘探场景中的典型理论模型的反演测试,验证了本文算法的有效性.     

Inversion of surface gravity data for 3-D density modeling of geologic structures using total variation regularization

We develop an inversion procedure using the total variation (TV) regularization method as a stabilizing function to invert surface gravity data to retrieve 3-D density models of geologic structures with sharp boundaries. The developed inversion procedure combines several effective algorithms to solve the TV regularized problem. First, a matrix form of the gradient vector is designed using the Kronecker product to numerically approximate the 3-D TV function. The piecewise polynomial truncated singular value decomposition (PP-TSVD) algorithm is then used to solve the TV regularized inverse problem. To obtain a density model with depth resolution, we use a sensitivity-based depth weighting function. Finally, we apply the Genetic Algorithm (GA) to select the best combination of the PP-TSVD algorithm and the depth weighting function parameters. 3-D simulations conducted with synthetic data show that this approach produces sub-surface images in which the structures are well separated in terms of sharp boundaries, without the need of a priori detailed density model. The method applied to a real dataset from a micro-gravimetry survey of Gotvand Dam, southwestern Iran, clearly delineates subsurface cavities starting from a depth of 40 m within the area of the dam reservoir.     

用改进的光滑约束最小二乘正交分解法实现电阻率三维反演

对三维电阻率反演问题进行了深入研究,提供了一种利用地表观测数据实现三维反演的实用算法.该方法应用有限差分求正演解,并通过对粗糙度矩阵元素进行适当改进,使之适用于各种情况下粗糙度矩阵的求取,进而建立在模型的总粗糙度极小条件下的反演方程.对反演方程采用收敛速度快且稳定的最小二乘正交分解（LSQR）法进行迭代求解,在迭代求解过程中只需利用偏导数矩阵和其转置矩阵乘以一个向量的结果,回避了直接求偏导数矩阵的繁琐计算,节省了内存,加快了反演的计算速度.不同的计算实例表明上述方法是求解大规模三维电阻率反演问题的有效方法.     

2-D Versus 3-D Magnetotelluric Data Interpretation

In recent years, the number of publications dealing with the mathematical and physical 3-D aspects of the magnetotelluric method has increased drastically. However, field experiments on a grid are often impractical and surveys are frequently restricted to single or widely separated profiles. So, in many cases we find ourselves with the following question: is the applicability of the 2-D hypothesis valid to extract geoelectric and geological information from real 3-D environments? The aim of this paper is to explore a few instructive but general situations to understand the basics of a 2-D interpretation of 3-D magnetotelluric data and to determine which data subset (TE-mode or TM-mode) is best for obtaining the electrical conductivity distribution of the subsurface using 2-D techniques. A review of the mathematical and physical fundamentals of the electromagnetic fields generated by a simple 3-D structure allows us to prioritise the choice of modes in a 2-D interpretation of responses influenced by 3-D structures. This analysis is corroborated by numerical results from synthetic models and by real data acquired by other authors. One important result of this analysis is that the mode most unaffected by 3-D effects depends on the position of the 3-D structure with respect to the regional 2-D strike direction. When the 3-D body is normal to the regional strike, the TE-mode is affected mainly by galvanic effects, while the TM-mode is affected by galvanic and inductive effects. In this case, a 2-D interpretation of the TM-mode is prone to error. When the 3-D body is parallel to the regional 2-D strike the TE-mode is affected by galvanic and inductive effects and the TM-mode is affected mainly by galvanic effects, making it more suitable for 2-D interpretation. In general, a wise 2-D interpretation of 3-D magnetotelluric data can be a guide to a reasonable geological interpretation.     

大地电磁二维对称各向异性介质的有限元数值模拟

假定垂直轴为二维对称各向异性介质主轴之一，构造走向与另一主轴方向成任一夹角，用伽勒金（Ｇａｌｅｒｋｉｎ）有限元法和矩形网格，优化地合成总体刚度矩阵，使二维各向异性介质的基本方程形成有限元代数方程组，求出各节点场值，并利用ＭＯＭ法求出辅助场，进而求出张量阻抗等响应函数．同时还对有关文献的计算模型进行了数值模拟检验     

三维地形直流电阻率有限元法模拟

基于稳定电流场的基本方程、三维区域满足的边值问题以及相应的变分问题,研究了三维起伏地形条件下电阻率的有限单元数值模拟算法. 离散积分区域时,以三棱柱为最小研究单元,推导了含有地形特征信息的三线性插值型函数以及单元刚度矩阵. 采用变带宽、一维数组方式只存储稀疏刚度矩阵中非零元素,能够节约内存. 利用Cholesky分解法只分解一次大型稀疏矩阵,通过回代可以求出方程组的全部解,当求解有多个供电点的测深问题时可以缩短计算时间. 模型计算表明,在水平层状介质模型上,三维计算结果与解析解或二维数值解十分吻合,计算精度满足误差要求. 在二维山脊上的二极剖面或三维山谷上的中间梯度剖面上,其三维计算结果与相应模型的土槽实验结果或边界元法计算结果也非常接近.     

Parallel 3-D Simulation of a Fault Gouge Using the Lattice Solid Model

Despite the insight gained from 2-D particle models, and given that the dynamics of crustal faults occur in 3-D space, the question remains, how do the 3-D fault gouge dynamics differ from those in 2-D? Traditionally, 2-D modeling has been preferred over 3-D simulations because of the computational cost of solving 3-D problems. However, modern high performance computing architectures, combined with a parallel implementation of the Lattice Solid Model (LSM), provide the opportunity to explore 3-D fault micro-mechanics and to advance understanding of effective constitutive relations of fault gouge layers. In this paper, macroscopic friction values from 2-D and 3-D LSM simulations, performed on an SGI Altix 3700 super-cluster, are compared. Two rectangular elastic blocks of bonded particles, with a rough fault plane and separated by a region of randomly sized non-bonded gouge particles, are sheared in opposite directions by normally-loaded driving plates. The results demonstrate that the gouge particles in the 3-D models undergo significant out-of-plane motion during shear. The 3-D models also exhibit a higher mean macroscopic friction than the 2-D models for varying values of interparticle friction. 2-D LSM gouge models have previously been shown to exhibit accelerating energy release in simulated earthquake cycles, supporting the Critical Point hypothesis. The 3-D models are shown to also display accelerating energy release, and good fits of power law time-to-failure functions to the cumulative energy release are obtained.     

3-D Numerical modeling of an electromagnetic field in geoelectrics using the Trefftz method

The algorithm for numerically solving the direct 3-D problem of calculating the electromagnetic field varying harmonically in an arbitrary inhomogeneous 3-D media is developed on the base of the Trefftz method, which has not been used previously in geoelectrics. The corresponding system of algebraic equations has been solved with the use of a modification of the well-known Kaczmarz iterative method. A cyclic method of equalization is used as a procedure of preconditioning the matrix of the system.     

青藏高原东部及其边缘地区的地壳上地幔三维速度结构

利用在青藏高原东部及其邻近地区记录到的１万余条近震到时资料,反演该地区的地壳上地幔三维速度结构。采用网格点模型描述三维速度结构,模型维数为２２２２６,网格点间距水平向为１００ｋｍ,垂直向为２０ｋｍ,网格点之间的速度值通过线性插值给出。采用改进了的快速三维射线追踪方法,确定三维非均匀介质中的地震射线路径和理论走时。反演结果显示,青藏高原南部的上地壳中（３０ｋｍ左右的深度）存在一低速区,这和面波反演的结果一致,羌塘块体下地壳有明显的低速异常带,青藏公路沿线的垂直速度剖面显示出岩石层受挤压增厚的构造特征。     

Grounded-source transientE-field modeling of a shallow resistive layer with 3-D inhomogeneity

There is growing interest in the use of transient electromagnetic (TEM) sounding for shallow geotechnical, environmental and groundwater investigations. Two commonly used transmitterreceiver configurations for TEM sounding are 1) loop-loop or its variation, in-loop configuration and 2) wire-loop configuration. The less common configuration of a horizontal electric dipole (HED) transmitter and receiver is treated in this study and called wire-wire configuration.Two important problems of shallow investigation in hard and soft rocks respectively are, defining 1) a fractured/fissured zone of medium resistivity, sandwiched between an overlying surface weathered rock of low resistivity and an underlying fresh compact rock of high resistivity and 2) a body of resistive sand buried in conductive clay. Lateral change in the middle layer resistivity is modeled by including a 3-D body of anomalous resistivity. The effect of perturbing the resistivity of the 3-D inclusion and the host middle layer for the wire-wire configuration is compared with that of the commonly used loop-loop configuration. The wire-wire configuration is found more sensitive to the model perturbations than the loop-loop configuration.1-D inversion of synthetic 1-D data sets for the wire-wire configuration finds resolution and estimation errors to be less than 10 percent for all the model parameters. For 3-D models, 1-D inversion results give a resolution error of 10 percent or less for the depth to, resistivity and thickness of, the 3-D inclusion. The estimate is within 10 percent of the true value for the first parameter but 40 percent for the other two. Resolution as well as estimation of the basement resistivity is always very poor.Using the wire-wire configuration, it is theoretically possible to define a buried resistive layer and any lateral change in its resistivity, subject to the above limitations of 1-D inversion. However, the basement resistivity cannot be estimated with reasonable accuracy in the presence of a lateral inhomogeneity in the overlying layer.     

利用共轭梯度算法的电阻率三维有限元正演

引入对称超松弛预条件共轭梯度(SSOR PCG)迭代算法求解电阻率三维有限元计算形成的大型线性方程组，并有机结合系数矩阵的稀疏存储模式，使得三维有限元正演计算的速度大大提高而内存需求则大大减少. 该算法可方便地应用于求三维异常电位，在保持快速计算的基础上，正演计算的精度显著提高.