二维电阻率成像的有限元解法
有限单元法在电阻率法地形改正中的应用
［ ５ ］ 徐世浙， 赵生凯， 地形对大地电磁勘探的影 响 ［ Ｊ ］ ． 西北地震学 １ ９ ８ ５ ， ７ （ ４ ） ： ６ ９－ ７ ８ ． 报， ［ ６ ］ 徐世浙 ． 点源二维各 向 异 性 地 电 断 面 直 流 电 场 的 有 限 元 解 法 ［ Ｊ ］ ． 山东海洋学院学报， １ ９ ８ ８ ， １ ８ （ １ ） ： ８ １－ ９ ０ ． ［ ７ ］ 周熙襄， 等． 点源二维电法正演的有限单元 法 ［ Ｊ ］ ． 物化探电子 １ ９ ８ ３ ， （ ３ ） ： １ ９－ ４ ０ ． 计算技术， ［ ８ ］ 黄兰珍， 田 宪 谟， 等． 点源场电阻率法二维地形改正边界元法 ［ Ｊ ］ ． 物探化探计算技术， １ ９ ８ ６ ， （ ３ ） ： ２ ０ １－ ２ ０ ７ ． ［ ９ ］ 梁芳敏， 张小路， 匡应平， 等 ．边界元法在中梯装置视电阻率二 维地形改正中的应用［ Ｊ ］ ． 矿产勘查， ２ ０ １ ０ ， １ （ ４ ） ： ３ ９ ５－ ３ ９ ８ ．
图１ 山谷、 山脊模型联剖实验结果
和３ ０ ° 山谷模型实验得 到 的 联 合 剖 面 纯 地 形 视 电 阻 率曲线。由图 １可 见， 在正地形上（ ） ， 顶部由 φ＜π
Ａ Ｂ 于电流密度减小而表现为低阻异常， ρ ｓ与 ρ ｓ 异常形
成反交点（ 山脊顶部和山谷的两侧） ； 在 负地形 上 （ φ ＞π ） ， 顶部则 形 成 高 阻 正 交 点 异 常 （ 山谷底部和山
２ 有限单元法在地形改正中的应用
２ ． １ 模拟纯地形视电阻率计算
如图 １ （ ａ ） （ ｂ ） 引用于文献［ ２ ］ ， 分别为 ３ ０ ° 山脊
１ 方法原理
［ １ ］
根 据 点 源 二 维 电 场 的 基 本 特 征 以 及 边 值 问 题、 变分问题对相关公 式 进 行 推 导， 并且用傅立叶变换
高密度电阻率法二维和三维有限差分正演计算
高密度电阻率法二维和三维有限差分正演计算高密度电阻率法（High-Resistivity Method，简称HRM）是一种用于精确描述电流在较大空间内散射和复杂反射情况下的分布和阻抗、电磁场和电动势场强度的方法。
它可以精确地计算出沿一个给定路径上电流和电场强度的分布，从而获得某一个时刻电阻地层的特点。
有别于PT和MT方法，HRM有明显优势，可以处理电磁参数和地质参数的复杂变化，适应复杂的反射现象，以及更好地解析各种地表以下界面及其结构形成的脉动反射特征。
HRM二维和三维正演计算都基于相同的数学模型，有一个主要的区别，即边界条件的处理不同。
二维正演计算是基于一维电磁方程组和单边带边界条件，而三维正演计算基于3维电磁方程组，具有等张场及大气层的沿着深度的变化。
HRM的数值正演模拟具有以下特点：（1）快速收敛：基于电磁方程组，通过精确描述电极间耦合时，就可以获得快速收敛。
（2）反射损失小：以解耦形式处理反射，它们的损失更小，从而可以更好地获得精确的结果。
（3）地质参数可调：由于地表以下空间有不同的地质参数，如电导率、磁导率等，只需简单地修改这些参数，就可以很好地模拟不同的深部地质结构。
HRM正演模拟分为2步完成：第一步建立模型，根据真实的地质情况设置正演路径，并确定正演模型的地表，地表以下的地质学参数及边界条件；第二步进行数值正演模拟，根据建立的模型，使用有限差分法的方法解决电磁方程，获得强度场及电场、磁场，并与水文地质诸要素参数建立联系以便反演出水文地质参数。
HRM正演模拟应用非常广泛，能够获得较精确的电流分布和阻抗、电磁场和电动势数据集，对地质勘查提供了有效的技术手段，例如辅助高精度的地质资料的获取和评价，预测浅层以及深部复杂地质结构非常有用，可帮助开发者识别和定位储层分布，有助于钻井现场安全和高效地进行。
二维热传导方程有限元体积法
二维热传导方程有限元体积法（实用版）目录一、二维热传导方程的概述二、有限元体积法的基本原理三、二维热传导方程的有限元体积法求解过程四、有限元体积法的优缺点及应用范围五、结论正文一、二维热传导方程的概述二维热传导方程是描述物体在二维空间中热传导过程的偏微分方程，其基本形式为：T = λT/t其中，T 表示温度分布，λ表示热扩散系数，t 表示时间。
二维热传导方程在工程领域有广泛的应用，例如电子器件冷却、建筑节能等。
二、有限元体积法的基本原理有限元体积法是一种求解偏微分方程的数值方法，其基本思想是将求解域划分为有限个小的子区域（单元），然后在每个子区域中选取一个基函数，用基函数的线性组合表示子区域中的解，最后将子区域的解拼接起来得到整个求解域的解。
有限元体积法的基本步骤包括：划分网格、选择基函数、组装矩阵、求解线性或非线性方程组、计算结果等。
三、二维热传导方程的有限元体积法求解过程1.划分网格：将求解域划分为二维网格，每个网格单元的尺寸可以根据实际问题和求解精度要求进行调整。
2.选择基函数：在每个网格单元内选择一个基函数，如拉格朗日基函数、塞尔达基函数等。
基函数的选择会影响到求解精度和计算复杂度。
3.组装矩阵：根据基函数的选取，组装离散化矩阵和系数矩阵。
离散化矩阵用于表示求解域中的解，系数矩阵用于表示偏微分方程的离散形式。
4.求解线性或非线性方程组：根据组装的矩阵，求解线性或非线性方程组，得到网格点上的温度分布。
5.计算结果：将每个网格点上的温度分布拼接起来，得到整个求解域的温度分布。
四、有限元体积法的优缺点及应用范围有限元体积法具有以下优点：1.适用于大型、复杂的求解域；2.可以灵活选择基函数和网格划分，以适应不同的求解精度和计算资源要求；3.能够处理非线性方程组；4.求解过程具有较好的稳定性和收敛性。
有限元体积法的缺点主要包括：1.计算复杂度较高，尤其是对于大规模问题；2.需要预先划分网格，可能导致计算过程中出现网格畸变。
偏微(20)二维有限元法
ai bi
aj bj
am （2.15） bm
三、单元上的线性插值多项式 Ni ( x, y)ui N j ( x, y)u j N m ( x, y)um （2.9） uh ( x, y )
uh uh ( x, y) Nue
（2.13）
（2.9）式和（2.13）式就是e上根据三角形的
Bue
（2.14）
三、单元上的线性插值多项式
uh uh uh ( x, y) Nue x Bue uh uh uh ( x, y )的梯度向量可以表示为 y （2.14）
N i x B N i y N j x N j y N m x 1 N m 2 e y
1 (ai ui a j u j am um ), a 2 e 1 (bi ui b j u j bm um ), b 2 e 1 (ci ui c j u j cm um ). c 2 e
，
axi byi c ui , ax j by j c u j , axm bym c um .
写成了第一个式子，例如，(2.10)式中N i ( x , y ) 则其他的式子（N j ( x, y)和N m ( x, y) 式）可以 即将第一式的脚标 通过 i , j , m 的脚标轮换得到。
， i换为 j，j 换为 m m 换为 i 就可得到第二式。
再换一次得第三式。
三、单元上的线性插值多项式 uh ( x, y ) Ni ( x, y)ui N j ( x, y)u j Nm ( x, y)um
设单元en为Pi Pj Pm 在其三个顶点上，
二维高密度电阻率法数值模拟与应用_金聪
第5 期
金聪等：二维高密度电阻率法数值模拟与应用
图1 Fig. 1
高密度电阻率法常用排列示意图
Sketch of commonly used arrangements for the highdensity resistivity method
B － 供电电极；M、 N － 测量电极；a － 电极距；n － 隔离系数； （ a） － 温纳 － 斯伦贝谢装置；（ b） － 温纳装置； A、 （ a） － WennerSchlumberger device； （ b） － Wenner device； A，B － current electrode； M，N － measuring electrode； a － electrode spacing； n － isolation coefficient
一定方式相互联接的单元。 根据电场所满足的微分方程以及边界条件， 相 应的变分问题为： F（ u） = σ（ u ） ∫[ 1 2
Ω 2
+ σ’ u0 ·u d Ω +
2
]
n） σcos（ r， u ∫ [1 r 2
Γ∞
+
cos（ r， n） σ’ u0 u d Γ， r
]
（2） δF（ u） = 0 式中：u0 是正常电位； u 是异常电位； σ 为介质 电导率； σ ’ 为异常电导率； n 为边界外法线方向； r 为电源点到边界上点的向径， Γ 为边界 （ 刘斌 等， 2010 ） 。 用矩形单元对区域 Ω 进行剖分， 在各单元上近 似地将变分方程离散化， 并通过单元分析和总体合 成， 导出以各节点电位值为变量的高阶线性方程组 ； 最后解此方程组， 再通过反傅里叶变换算出各节点 的电位值， 得到地下半空间场的分布， 以表征稳定电 。 流场的空间分布 1. 3 光滑约束最小二乘反演 最小二乘反演作为地球物理反演最常用的方 法， 它对地球物理资料反演具有很好的稳定性和适 应性。最小二乘反演问题可以表示为求使目标函数 Ψ 极小时的最优模型参数改正值矢量 Δm。 为了解 决反演问题常常表现出的混定问题， 将光滑约束引 入反演方程， 有效地提高了反演计算效率。 施加光 2012 ） ： 滑约束后的目标函数 Ψ 为（ 刘斌等，
二维均匀电场的直流电阻率法的边界单元法
二维均匀电场的直流电阻率法的边界单元法一、边界元法计算均匀场中水平地形条件下二维不均匀体的异常设在均匀场0E 中。
有一电阻率为2ρ的二维不均匀体，由于地面水平，故可由镜像法原理求空间任一点的电场（如图(1.1)所示）上述电场的边值问题满足：１ Ω∈=∇02u (1.1) １ Ω∈=∇02v （1.2)ΓΓ=v u (1.3)ΓΓ∂∂=∂∂nv nu2111ρρ(1.4)x E u 0-=∞Γ (1.5)（一）积分方程的建立由格林公式 图.1.1区域与边界Γ∂∂-∂∂=Ω∇-∇⎰⎰∞Γ+Γ+ΓΩd nun ud u u )()('122φφφφ (1.6) 设r1ln 21πφ=为二维拉普拉斯方程的基本解 (1.7) 式中r 为Ω中任意点至电场计算点P 的距离则 )(2p δφ-=∇ (1.8) )(p δ为以p 为中心的δ函数，(1.1)、（1.3)式代入(1.6)式左边 ⎪⎩⎪⎨⎧Γ-Ω--=Ω∇-∇⎰⎰ΩΩ上在 内在 p u p p u p u d u u p p πωδφφ2)()(1221 (1.3) 式中p u 是p 点的电位，p ω为p 对区域1Ω的张角。
格林公式的右边为Γ∂∂-∂∂+Γ∂∂-∂∂=Γ∂∂-∂∂⎰⎰⎰∞∞ΓΓ+ΓΓ+Γ+Γd n un u d n u n u d n u n u)(()(''φφφφφφ而 r1ln 21πφ=),cos(21n r rn r r n πφφ-=∂∂∂∂=∂∂由于边界∞Γ在无穷远处，我们可近似将p 点视作∞Γ的圆心，因而在∞Γ上，下列关系式成立⎪⎪⎪⎩⎪⎪⎪⎨⎧-=∂∂-=-=∂∂+-=-=r n y r E x n E nuy r r x E x E u p πφ21),sin(),cos()],sin([0000 (1.10) 式中p x 为p 点的x 坐标 (r,y)为失径r 与坐标y 的夹角 (n,x)为法向与坐标x 坐标的夹角 (1.10)式代入Γ++=Γ∂∂-∂∂⎰⎰∞∞ΓΓd y r E rr y r r x E d n u n up )},sin(1ln 2121),sin([{)(00ππφφp u x E d x E rp p 00021-==Γ=⎰∞Γπ 式中p u 0为点p 点的初始电位Γ∂∂-Γ-=Γ∂∂-∂∂⎰⎰⎰Γ+ΓΓ+ΓΓ+Γd rn u d r n r d n u n u '''1ln 212),cos()(ππφφ于是当p 点在Ω内时，则(1.6)式变为Γ∂∂+Γ+=⎰⎰Γ+ΓΓ+Γd rn u d r n r u u u op p ''1ln 212),cos(ππ (1.11)当P 在地面时Γ∂∂+Γ+=⎰⎰ΓΓd rn u d r n r uu u op p 1ln 2122),cos(2ππ (1.12)由此可见，只要知道边界Γ上的u 和nu∂∂，便可由积分方程(3.1.12)求出地表的u.。
python有限元差分求解二维偏微分方程
python有限元差分求解二维偏微分方程摘要：1.引言2.有限元方法简介3.差分方法简介4.二维偏微分方程的求解方法5.Python 有限元差分求解二维偏微分方程的实现6.总结正文：【引言】本文旨在介绍如何使用Python 有限元差分法求解二维偏微分方程。
有限元方法和差分方法是数值计算领域中常用的方法，它们可以用来解决偏微分方程问题。
在本文中，我们将使用Python 编程语言来实现这些方法，以解决二维偏微分方程。
【有限元方法简介】有限元方法是一种数值计算方法，它将连续的求解区域离散化为有限个单元，通过在每个单元内求解局部问题，最后将各单元的解合并得到原问题的解。
有限元方法广泛应用于固体力学、流体力学、热传导等领域。
【差分方法简介】差分方法是一种数值计算方法，它通过将连续的函数值用离散的点表示，然后通过差分公式将函数在某点的导数表示为有限差分形式。
差分方法主要包括前向差分、后向差分和中心差分等。
差分方法在数值计算中具有重要地位，例如在求解常微分方程的初值问题、边值问题以及偏微分方程等方面都有广泛应用。
【二维偏微分方程的求解方法】二维偏微分方程的求解方法主要包括有限元方法和差分方法。
有限元方法将求解区域离散化为有限个单元，通过在每个单元内求解局部问题，最后将各单元的解合并得到原问题的解。
而差分方法通过将连续的函数值用离散的点表示，然后通过差分公式将函数在某点的导数表示为有限差分形式。
这两种方法在求解二维偏微分方程时，可以相互结合，提高求解精度。
【Python 有限元差分求解二维偏微分方程的实现】Python 作为一门功能强大的编程语言，提供了丰富的数值计算库，如NumPy 和SciPy。
利用这些库，我们可以方便地实现有限元差分法求解二维偏微分方程。
下面是一个简单的示例，使用Python 求解二维Laplace 方程：```pythonimport numpy as npfrom scipy.spatial import griddatadef Laplace(x, y):return np.sin(x) + np.cos(y)def Laplace_dirichlet(x, y, x0, y0):return np.sin(x) + np.cos(y) - np.sin(x0) - np.cos(y0) # 网格划分x = np.linspace(0, 2 * np.pi, 100)y = np.linspace(0, 2 * np.pi, 100)X, Y = np.meshgrid(x, y)# 边界条件x0, y0 = 0, np.pi / 2# 求解Z = griddata(X, Y, Laplace(X, Y), (x, y), method="cubic")Z_dirichlet = griddata(X, Y, Laplace_dirichlet(X, Y, x0, y0), (x, y), method="cubic")# 绘制结果import matplotlib.pyplot as pltplt.figure()plt.contourf(x, y, Z)plt.colorbar()plt.title("二维Laplace 方程的数值解")plt.show()plt.figure()plt.contourf(x, y, Z_dirichlet)plt.colorbar()plt.title("二维Laplace 方程的有限元差分解")plt.show()```【总结】本文介绍了Python 有限元差分求解二维偏微分方程的方法。
海洋可控源电磁法二维有限元正演及反演
海洋可控源电磁法二维有限元正演及反演海洋可控源电磁法（Marine Controlled-Source Electromagnetic，MCSEM）是一种新兴的海洋地球物理勘探技术。
它通过在海洋中设置控制电源，利用电磁场对海底地质进行探测，可以获取地底结构的信息。
本文将从MCSEM的原理入手，介绍MCSEM的二维有限元正演和反演方法，并讨论其应用及优缺点。
一、MCSEM原理MCSEM利用控制电源发出高频电信号，该信号在海洋中传播时，会激发海底地质物体中的电流。
这些电流在海水中会产生电磁场，通过检测电磁场的变化，可以解析出地底物质的电导率、磁导率等物理参数，从而获取地底结构信息。
MCSEM主要有两种控制电源：直流电源和交流电源。
直流电源具有较大的侵入深度和较好的低频响应，适用于大区域的浅层勘探；交流电源具有较好的高频响应，适用于小区域的深部勘探。
控制电源的设置可以根据勘探的需求进行调整。
二、二维有限元正演MCSEM二维有限元正演是指将地下介质分布描述为二维平面的电导率分布，采用有限元理论计算出该模型产生的电磁响应。
MCSEM的二维有限元正演主要包括以下步骤：1.建立数学模型建立地下介质的二维平面电导率模型，将所要研究的地质构造物体分为多个区域。
2.建立有限元网格将地下介质划分为若干小块，每个小块中的电导率均为常数。
在每个小块内部建立一个节点，并通过有限元网格连接所有的节点。
3.设置边界条件和时间离散在模型的边缘处设置边界条件，确定控制电源和检测电极的位置。
对时间进行离散化，并设置时间步长。
4.求解定态矩阵根据有限元法原理，求解定态矩阵，包括系统矩阵、电势向量和电流向量。
5.计算电场和磁场响应根据电场和磁场的计算公式，使用有限元法计算出电场和磁场的响应。
三、二维有限元反演MCSEM二维有限元反演是指利用已知的电磁响应数据，反推出地下介质的电导率分布。
MCSEM的二维有限元反演主要包括以下步骤：1.误差函数的定义确定反演的误差函数，通常采用观测数据与模拟数据之间的二次差值作为误差函数。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
990318岩石力学与工程学报 CHINESE JOURNAL OF ROCK MECHANICS AND ENGINEERING 1999年 第18卷 第3期 Vol.18 No.3 1999
二维电阻率成像的有限元解法*
底青云 王妙月摘要 为了研究电流线追踪电阻率成像的可靠性， 首先开展了二维电阻率成像的有限元方法研究， 在传统有限元成像的基础上，做了三方面的研究工作： 一， 推导了有限元反演方程， 从理论上对比了有限元方法和电流线追踪方法获得反演系数矩阵元素的优缺点； 二， 采用分块压缩的方式求反演系数， 节省计算机内存和计算机时间； 三， 用改进的塞德尔方法快速求解反演方程。 从试算的数值模型实例来看， 结果比较可靠。 关键词 有限元， 电阻率成像 分类号 O242.1
TWO DIMENSION RESISTIVITY TOMOGRAPHY USING THE FINITE ELEMENT METHOD
Di Qingyun Wang Miaoyue(Institute of Geophysics, The Chinese Academy of Sciences, Beijing 100101)Abstract In order to verify the reliability of resistivity tomography by tracing the electric current ray, the research of 2D resistivity tomography using the finite element method is completed. The work is done from three aspects. First, the reverse FEM equations with artificial boundary terms are derived, and their advantages and shortages are pointed out compared with the method of tracing electric current ray. Second, through solving the reversed coefficient using the method of section by section and using the sparse and symetrial property of the coefficient, the computing time and computer memory are saved. At last, the reverse equations are solved by the improved Gauss-Sadal iteration method. The results of model simulation are much reliable. Key words FEM, resistivity tomography
1 引 言
直流电法勘探一直是资源勘探、水文工程、地质勘察以及地下埋设物调查的一种
有效方法［1,2］。反演求解中最常用的方法有一维自动反演；对于二维地电结构主要是
采用有限元法、有限差分法和边界元法等技术进行模拟分析；对三维地电结构的求解
file:///E|/qk/yslxygcxb/yslx99/yslx9903/990318.htm（第 1／11 页）2010-3-23 13:16:27990318方法则较少见［3,4］。这些传统的方法主要是对地表的观测资料进行处理求解。近年来，随着高密度电法采集系统的发展，采用阵列式布极，需要发展能够对跨孔、孔～地、以及地表的观测结果进行反演求解的方法。因此，借鉴地震波和电磁波的层析成像技术，不少学者开始了对电阻率的重建成像研究，并取得了可喜的成果［5,6］。
电阻率成像兼顾了地震波和电磁波成像的优点，野外作业不但方便，而且穿透深度也较大，可产生比常规电阻率勘探方法更高的分辨率，电阻率成像对于解决识别断层、破碎带、油气层、地下水及其污染等问题非常有用。 随着高密度电阻率成像仪的引进，高密度电阻率野外观测可以顺着地形面，涉及的深度可以从米的量级直到1km的深度，采用二维电阻率有限元成像方法可以获得这个深度范围内电性细结构图像，这对解决岩土工程中的岩土力学问题及工程设计问题非常有用。 国外二维和三维电阻率成像技术［7～9］，在广义环境勘探以及地质调查的众多领
域已经取得重要进展。发表的方法主要是有限元方法，电流线方法极为少见。在国内，利用直流电场的观测结果对地下电阻率结构进行成像的研究则刚刚起步。本文借鉴前人的工作经验，试图就电阻率成像问题的原理、方法及实现过程作一探讨。文中首先从方程的推导出发，给出了电阻率成像的有限元反演公式和电流线反演公式。通过对有限元成像与电流追踪成像方法的理论对比，表明两种方法各有优缺点。考虑到有限元成像国外已有比较多的文章发表，其可靠性已经得到了论证，本文重点讨论2D有限元成像方法的实现过程。通过几个数值模拟的实例，进一步论证有限元方法的有效性及可行性，为电流线追踪等其他电阻率成像方法的可靠性研究提供了资料基础和方法对比的研究基础。
2 电阻率反演方程的推导
设2D研究区域可以分割成N个矩形网格块，每个网格块内电阻率均匀，则待求的
2D研究区域内的电阻率分布是一个N维矢量ρ
设ρ0为ρ的初值，则稳定电流场中任一点R的电位Uk(ρ)可以表示成ρ0点的泰勒级数，取一次项有
file:///E|/qk/yslxygcxb/yslx99/yslx9903/990318.htm（第 2／11 页）2010-3-23 13:16:27990318写成矩阵形式是［A］δρ=δU (1)式中：δU为一系列观测点的电位矢量与这些点初始模型为ρ0时的理论电位矢量的差
值；A是系数矩阵，它的元素是；δρ为待求的电阻率分布的修正值。一旦δρ求
得，则
ρ=ρ0+δρ (2)
式(2)中的ρ即为所求的解答。 可见为了求得ρ，求得系数矩阵A是关键，电阻率反演方程的推导归结为求得可
靠的。下面从二个方面阐述求解的方法，并比较它们的优缺点。
2.1 电流线方法 稳定电流场中电位分布满足U=-ρJ 即 U=-E (3)式中：U为电位，ρ为电阻率，J为稳定电流场的电流密度。取n0为J的方向，亦为等位
面U的法线方向。沿着电流线则有
(4)S为电流线弧长。于是在标量形式下结合式(3)，(4)成为
(5)对于第i个网格，电阻率为ρi，δSi为第k条电流线在第i个网格内的长度，则由式(5)，
file:///E|/qk/yslxygcxb/yslx99/yslx9903/990318.htm（第 3／11 页）2010-3-23 13:16:27990318δSi二端之间的电位差δUi为
(6) (7)Aki即为我们所求的系数矩阵A的元素，其中F0为量纲因子，Ji为第i个网格块处的电流密
度，为该处电流密度随电阻率的变化值。
2.2 有限元方法 2D线源稳定电流场的微分方程的边值问题的有限元解法可归结为有限元方程［10］：
KU+F边=I (8)
(9)式中：矩阵K称为刚度矩阵，U为节点的电位矢量，I为节点供电电流矢量，U边为边界
上的节点的电位矢量，n为边界法线方向坐标，F边为人为截断边界的贡献。对于自由
面边界，F边恒等于零。σ为电导率即电阻率的倒数。对于截断边界，节点i对节点i自
身的贡献为
(10)节点i对其相邻节点m的贡献为file:///E|/qk/yslxygcxb/yslx99/yslx9903/990318.htm（第 4／11 页）2010-3-23 13:16:27990318 (11)其中σ右，σ左分别为节点i右、左单元的电导率，当相邻节点m在节点i右边时σ=σ
右，在左边时，σ=σ左。X，Y为网格块的边长。
对式(8)等式两边求ρi的偏导数，则有
(12)式中：
(13)式(12)中即为待求的系数矩阵A，在式(12)中由式(13)可知β已知，K为刚度矩阵，当给定ρi的初值后，亦为已知，因此类似于式(8)的有限元方程可以求得所需的系数矩
阵A。
2.3 电流线方法和有限元方法的比较 电流线方法求系数矩阵A归结为求解式(7)。式(7)表明对于系数矩阵A的每个元素Aki和第i个网格块的参数有关，包括第i个网格块的电阻率数ρi，电流密度Ji，电流密度
变化率以及电流线越过第i个网格块的长度δSi。通过电流线追踪可以获得δSi。，Ji这两个量显然和电阻率的总体分布有关，和几何扩散有关，因此比较难
于确定，作为一级近似可用ρi的初值ρi0的理论电流密度替代。
有限元方法求系数矩阵A归结为求解式(12)，式(12)表明对于每一个网格块的ρi，
i=1,2,…，N，都需求解一次类似于式(8)的有限元方程。众所周知，求解有限元方程是比较费时间的。由于式(8)的右端只有供电点有值，因此只有少数分量不为零，一些快速有效的迭代方法，例如赛德尔迭代方法很难应用，需要用比较费时的QR分解、LLT分解方法。当模型比较大时(N比较大)，计算量是非常巨大的。 式(5)，(6)表明对于电流线方法，直接求电阻率比较方便，而对于有限元方程，由
file:///E|/qk/yslxygcxb/yslx99/yslx9903/990318.htm（第 5／11 页）2010-3-23 13:16:27