电法报告-直流电测深正演曲线.doc

关于多层水平地层地电断面电测深曲线的正演的读书报告 姓名：林俊 班级：061084-27 学号：20081003195 指导老师：师学明 日期：二一一年五月目录前言2 目的2 任务要求2 工作过程2 成果2原理3 1-多层水平地层上的对称四极电测深视电阻率表示式3 1.多层水平地层地面点电流源的电场3 2.多层水平地层上电测深的s 表示式和电阻率转换函数5 3.电阻率转换函数的递推公式6 2-水平地层上视电阻率的滤波算法6 3-多层水平地层的电测深曲线类型9 A 二层情况9 B 三层情况9 C 四层及多层情况9编程10 感想18关于多层水平地层地电断面电测深曲线的正演的读书报告前言目的：熟悉并掌握多层水平地层地电断面直流电测深曲线的正演任务要求：编制适用于n层地电断面的正演电测深程序（编程环境不限制，可用C语言，C+，VC，VB，matlab，推荐用matlab）。每个同学计算两个标准地电模型的正演计算第一个模型：二层G型地电模型 第一层地层电阻率 10欧姆米，第一层厚度10米； 第二层地层电阻率100欧姆米第二个模型：三层H型地电模型 第一层地层电阻率：班号（4） 100欧姆米，第一层厚度15米； 第二层地层电阻率：序号（27） 1 欧姆米，第二层厚度20米； 第三层地层电阻率：1000欧姆米AB/2为13个：2, 3, 4.5, 6, 9, 12, 15, 20, 30, 45, 60, 90, 120 （米）。工作过程： 先进行原理分析，再用matlab进行编程，最后小结。成果： 用matlab实现了n层地电断面的直流电测深正演。原理-1多层水平地层上的对称四极电测深视电阻率表示式 1.多层水平地层地面点电流源的电场 如图所示，水平地面下有n层水平地层，各层电阻率分别为1、2 n； 各层厚度分别为h1、h2hn-1； 各层底面到地表的距离分别为H1、H2Hn-1，Hn。求解思路：U E s Ti(r) 用分离变量法求方程 边界条件 通解 各层的电位表达式 利用边界条件求待定系数 A2、A3，,An B1、B2，,Bn-1 设地面点电流源A的强度为I。为求各层 中的电位表达式，将柱坐标系的原点设在A点，Z轴垂直向下。 在所设条件下，电位与角无关，满足如下形式的拉普拉斯方程： 边界条件： （1） （2） （3） （） （4） （5）由分离变量得到零阶贝赛尔方程，其通解为： 式中：A (m), B(m) 为待定的积分变量m的函数；J0(mr) 为零阶贝赛尔函数。利用边界条件，可以得到第一层电位公式： 第二层以下至n-1层，第 i 层的电位为： 第n层内的电位表达式，由 得： 电测深只在地面工作，即z0，故只需求出B1， 式中：J(mr)为零阶第一类贝赛尔函数；B(m)为积分变量m的函数。 对于层数确定的水平地层，根据地层界面上电位和电流密度法向分量连续的边界条件，可具体求出B1(m)的表示式。例如，最简单的二层水平地层，利用1 和2 岩层分界面的相应边界条件可具体求出 n=3时， 其中2.多层水平地层上电测深的s 表示式和电阻率转换函数令 则地面上电位公式为： 若采用MN0的装置测量，相应的s 表达式为： 令 则多层水平地层上的电测深s 公式简写成： 式中，T(m)定义为电阻率转换函数又称核函数。可见，电阻率转换函数与各层的层参数(厚度和电阻率)及积分量m有关。3.电阻率转换函数的递推公式 对于二层水平地层情况，若将式先后代入式和式，便得到二层水平地层的电阻率转换函数:归纳每一层的电阻率转换函数，就可导出电阻率转换函数的递推公式： 式中 为双曲函数 由此可知：当给定n，n-1,n和h1，h2，hn-1。可递推出上层或下层的电阻率转换函数。 电阻率转换函数递推公式的导出，免去应用边界条件解方程组求系数B(m)的计算。-2 水平地层上视电阻率的滤波算法 在电测深视电阻率的表达式中的被积函数可以分为两部分：一是包含地下各层所有信息（厚度及电阻率）的电阻率转换函数T1（m）；此外是与地层参数无关的贝塞尔函数。其虽没有解析计算结果，但可用线性滤波方法进行计算。根据电阻率转换函数的变化规律，对m的抽样取对数等间隔比较合适，因此，首先令ex=mr，则电阻率的表达式变为 根据采样定理，一个函数可以用它在等间距离散抽样点上的抽样值表达： 将电阻率装换函数用它在x数轴上的离散抽样值表达为： 记 则 将Ck=预先计算出来，实际上取有限项求和就可以达到足够的精度了。从频谱分析的观点看，当电阻率装换函数用它在x数轴上的离散抽样值表达式时，相当于滤取了频率高于1/2x的谐波成分，因此这种计算视电阻率的方法称为滤波计算方法，Ck称为滤波系数。为了提高线性滤波计算的精度，减少滤波系数的数目，需要对x的抽样点进行位移，实际使用的线性滤波计算视电阻率的公式为 式中，s（r）为供电极距为r时的视电阻率；T1(ekx+s/r)为mk=ekx+s/r时的电阻率转换函数；Ck为第k个滤波系数；x为抽样间隔；s为位移系数。 实际线性滤波计算常用的抽样间隔有两种。一种为六点是抽样间隔，即x=（ln10）/6=101/6，直直流电阻率测深曲线一般都采用这种抽样间隔进行计算。另一种为十点式抽样间隔，即 x=（ln10）/10=101/10，电磁测深曲线一般都采用这种抽样间隔进行计算。下表为常用的一套六点式抽样间隔的滤波系数，共有20个系数：k=120，其位移系数s=-2.1719，e-2.1719=0.11396。用上述公式就可计算某个供电极距r的视电阻率，只要计算20个对应于这个供电极距r的不同m值的电阻率转换函数T1（mk），将其与下表中对应的20个滤波系数相乘再求和就可以了。 因为要计算一条视电阻率测深曲线就需要计算多个供电极距r的视电阻率，为了减少计算工作量，取ri=eix=10k/6，对应的mk=ekx+s/ri=e(k-i)x+s，这样，计算不同供电极距的视电阻率所需要的电阻率转换函数大多数可以共用。下面是用数值滤波法计算视电阻率测深曲线的计算机流程：（1） 输入层参数，包括层数n、各层的层厚度和电阻率；（2） 输入要计算的供电极距范围，并由此得到r=eix中i的变化范围：imin imax；（3） 计算k-i的变化范围：（k-i）min （k-i）max ；（4） 用电阻率转换函数递推公式，循环计算mj=ejx+s，j=（k-i）min，（k-i）max时的T1（mj）；（5） 用滤波公式，循环计算ri=ei，i=imin，imax是的s（ri）。-3多层水平地层的电测深曲线类型 A 二层情况 二层地点断面上层电阻率为，下层电阻率为。有如下两种电测深曲线类型：（如右图） （1）D型：，基底为低阻； （2）G型：，基底为高阻。 B 三层情况 三层断面共有以下4种类型：（如右图） （1）A型：123，电阻率递增； （2）K型：13，中间高阻层； （3）H型：1223，电阻率递减。 C 四层及多层情况 四层及多层断面电阻率测深曲线类型用三层断面类型的组合表示。四层断面及其电阻率测深曲线类型共有8种，例如1234的四层断面及其电阻率测深曲线类型为HA型。五层断面及其电阻率测深曲线类型共有16中，例如135的五层断面及其电阻率测深曲线类型为KHK型。每多一层，曲线类型增加一倍。N层地电断面的电阻率测深曲线类型数为2n-1。编程要求计算两个标准地电模型的正演：第一个模型：二层G型地电模型 1=10 m h1=10 m 2=100 m第二个模型：三层H型地电模型 1=班号（4）100=400 m h1=15 m 2=序号（27）1 =27 m h2=20 m 3=1000 mAB/2为13个：2,3,4.5,6,9,12,15,20,30,45,60,90,120程序设计流程图程序如下：function ps=dcszhengyan(ab2,n,p,h)c=0.003042 -0.001198 0.01284 0.0235 0.08688 0.2374 0.6194 1.1817 0.4248 -3.4507 2.7044 -1.1324 0.393 -0.1436 0.05812 -0.0252 0.01125 -0.004978 0.002072 -0.000318;%滤波系数C，共20个rn=length(ab2);%极距的个数d=log(10)/6;%取样间距s=-2.1719;%位移系数ps=zeros(rn,1);%视电阻率for i=1:rn for k=1:20 m=exp(k*d+s)/ab2(i); T2=p(n);%电阻率转换函数 for j=n:-1:2 T1=p(j-1)*(p(j-1)*tanh(m*h(j-1)+T2)/(p(j-1)+T2*tanh(m*h(j-1);%电阻率转换函数 T2=T1; end ps(i)=ps(i)+T1*c(k); endend loglog(ab2,ps,-o)axis(1 1000 10 1000)grid on% xx=2:120; % yy=spline(ab2,ps,xx);%三次样条插值，使曲线变光滑% loglog(ab2,ps,-o,xx,yy);% axis(1 1000 10 1000)% grid on xlabel(AB/2 (m);ylabel(s (m);legend(正演结果);title(电测深曲线);return在Command Window中输入：clearclose allclc ab2=2,3,4.5,6,9,12,15,20,30,45,60,90,120;%极距AB/2 n=2;%两层地电断面模型p=10 100;%每一层的电阻率h=10 0;%每一层的厚度 % n=3;%三层地电断面模型% p=400 27 1000;%每一层的电阻率% h=15 20 0;%每一层的厚度ps=dcs_zhengyan(ab2,n,p,h) fid=fopen(psdate.txt,w);% open the file with write permissionfprintf(fid,AB/2 (m) s (m)rn);for i=1:length(ab2) fprintf(fid,%4g %12.4frn,ab2(i),ps(i);endfclose(fid)二层G型地电模型正演电测深曲线图如下： 输出结果为： 三层H型地电模型正演电测深曲线图如下： 输出结果为： 下面是经过三次样条插值光滑处理后的图形： 二层G型地电模型正演电测深曲线图 三层H型地电模型正演电测深曲线图下面是增加极距AB/2个数后的处理结果：ab2=1 2 3 4.5 6 9 12 15 20 30 45 60 90 120 150 200 300 450 600 900 1200 1500 2000 3000 4500 6000 9000 12000 15000 20000 30000 45000 60000 90000 120000 150000;%极距AB/2 二层G型地电模型正演电测深曲线图三层H型地电模型正演电测深曲线图三层电测深曲线具有等值性：对于H型曲线：当1、h1和3相同时，若在一定范围内按比例改变h2和2，保持S2=h2/2值不变，则地中电场分布变化甚小，因而测点附近的电流密度jMN也无明显改变。这将导致不同的地电断面对应形状几乎相同的s电测深曲线。这就是H型曲线的S等值性。 2=271.5=40.5 m h2=20 1.5=30m 2=272=54 m h2=20 2=40m2=273=81 m h2=20 3=60m2=274=108 m h2=20 4=80m感想 本次作为电法勘探资料处理这门课的作业，要求我们随便用一种编程工具来做多层水平地层地电断面电测深曲线的正演。我所编写的程序是基于matlab这样的一个高层次平台，所要实现的功能是多层水平地层地电断面电测深曲线的正演模拟，绘制出电测深曲线。因为在大二时曾经有一定的编程基础，是基于VC+平台的，后来还获得了全国计算机等级考试（二级C）优秀。再后来到大三上，我们又学过数值分析及信号分析与处理等课程，并多次涉及到matlab的相关知识，最后我还报名参加了数理学院开设的“计算方法数学实验”自主学习课程，主要就是针对matlab的各种应用。对比以前编程的经历，有些感想和收获。最深的印象是要编好专业软件需要先研究好专业知识，相对来说编程实现反而比较容易。比如我们在构建电测深曲线时，必须做到生成的电测深曲线符合地质规律，利用电阻率转换函数得到的视电阻率序列也要符合地质规律。为了使得到的电测深曲线更近完美，就必须做好原理分析，这里我先通过阅读程志平编写的电法勘探教程中水平层状大地对称四极电阻率测深曲线这一节，了解到理论上多层水平地层的视电阻率表达式的推导思路，即先求得电势U，再求得电场E，最后就可得到视电阻率s，而首先要求得电势U的话，就得结合边界条件用分离变量的办法解出U所满足的拉普拉斯方程，又由于电测深工作在地面上进行，只需求出地面上（z=0）的U1，因此只要能求得B1（m），便可得到U1，而不一定要求所有的通解中的待定系数。由此引出电阻率转换函数T及其递推公式的推导。