（检测技术与自动化装置专业论文）基于matlab的电阻层析成像（ert）与仿真平台设计.pdf

基于m a t l a b 的电阻层析成像( e r t ) 与仿真平台设计 张亮( 检测技术与自动化装置) 指导教师：石天明( 教授) 摘要 电阻层析成像( e r t ) 技术是近年来发展起来的一种基于电阻传感机 理的过程层析成像( 盯) 技术，适用于两相流中以导电性介质为连续相 的工业过程，可提供封闭的管道或过程容器设备内部多相组分物质参数的 二维三维可视化信息，具有非侵入、响应速度快、成本低、安全性能好、 适应范围广等优点。针对我校对电阻层析成像技术的研究尚处于初步阶段 的情况，在总结了电阻层析成像技术的研究现状、技术特点，分析了其数 学物理模型的基础上，基于m a t l a b 编写了一个电阻层析成像仿真软件包。 该软件包完成了通过有限元方法( f e m ) 对e r tj 下问题的求解；实现了 利用线性反投影算法与修正的牛顿拉夫逊类算法等e r t 反问题算法 进行图像重建，同时就反问题中的电极数、正则化因子、迭代初值以及噪 声的影响做了深入探讨：并提出了相关系数与方差作为重建图像的质量评 价指标，使重建图像的优劣有了在数值上的比较。另外，将遗传算法与粒 子群优化算法这两种人工智能优化算法引入到图像重建算法中，极大地提 高了两相流重建算法的收敛性和准确性。 关键词：电阻层析成像，图像重建，图像评价，人工智能算法 d e s i g no fam a t l a bb a s e ds i m u l a t i o np l a t f o r mf o r e l e c t r i c a lr e s i s t a n c et o m o g r a p h ym e a s u r e m e n t z h a n gl i a n g ( m e a s u r e m e n tt e c h n o l o g ya n d a u t o m a t i o ne q u i p m e n t ) d i r e c t e db yp r o f e s s o rt i a n m i n gs h i a b s t r a c t e l e c t r i c a lr e s i s t a n c et o m o g r a p h y ( e r t ) i so n eo fp r o c e s st o m o g r a p h y ( p dt e c h n i q u e sb a s e do nt h ep r i n c i p l eo fi m p e d a n c em e 勰u r e i n e 咄a n dh a s b e e nd e v e l o p e dv e r yf a s to v e rt h el a s td e c a d e e r ti sp a r t i c u l a r l ys u i t a b l ef o r t h e t w o p h a s e ( o rm u l t i - p h a s e ) f l o wa p p l i c a t i o n sh a v i n g ac o n d u c t i v e c o n t i n u o u sp h a s e ，a n di sa b l et op r o v i d e2 d 3 dv i s u a l i z a t i o ni n f o r m a t i o no f t h ep a r a m e t e r so f am u l t i p h a s em i x t u r ef l o w i n gi n s i d eac l o s e dp i p eo rv e s s e l e r to f f e r ss o m ee x c e l l e n tc u t t i n ge d g e si nm e a s l l r e m e ms u c ha sc a p a b l eo f y i e l d i n gp r o c e s si m a g e s ，n o n - i n t r u s i v et ot h ef l o wm e d i a , f a s t - r e s p o n d i n g ， l o w - c o s t ，n o n - r a d i a t i o n ，f l e x i b l em e a s u r e m e mo p t i 0 1 1 8e t c t h ew o r kd e s c r i b e d i nt h i sp a p e ra d d r e s s o st h er e s e a r c hn e e do fe r ti n i t i a t e di nc h i n au n i v e r s i t y o fp e t r o l e u m ，b yd e s i g n i n gam a t l a bb a s e ds i m u l a t i o ns o f t w a r ep a c k a g e c a p a b l eo fp e r f o r m i n gs o m eb a s i c c a l c u i a t i o n ss u c ha se r tp r o b l e m s m o d e l i n g ，s o l v i n g f o r w a r dp r o b l e m sa n dp r o v i d i n gs o m er e c o n s t r u c t i o n a l g o r i t h m s t h ep a c k a g ei n c l u d e saf o r w a r ds o l v e rd e v e l o p e db yu s i n gt h e f i m t ee l e m e n tm e t h o d ( f e m ) ，a n ds e v e r a le r ti n v e r s e r e c o n s t r u c t i o n m o d u l e ss u c ha st h el i n e a rb a c kp r o j e c t i o n ( l b p ) m e t h o da n dt h e m o d i f i e d - n e w t o n r a p h s o n ( m n a ) m e t h o d m o r e o v e r i nt h i sp a p e r , t h e i n f l u e n c e so ft h em e a s u r i n ge l e c t r o d en u m b e r s ，t h em a g n i t u d e o f i i i r e g u l a r i z a t i o nf a c t o r , t h ei n i t i a lv a l u e sf o ri t e r a t i o n , a n dt h en o i s ei m p a c t s o i l t h ei n v e r s ep r o b l e ma r ca l s od i s c u s s e di nd e t a i l s a l s oi nt h i sp a p e r , w e p r o p o s e dam e t h o dt ou s ce i t h e rac r o s s - c o r r e l a t i o nc o e f f i c i e n tm e t h o do ra v a r i a n c em e t h o dt oq u a n t i t a t i v e l ya s s e s st h ei n v e r s e di m a g eq u a l i v yo na c c o u n t o ft h el i k e n e s sb e t w e e nt h em o d e l e di m a g ea n dt h er e c o n s t r u c t e di m a g e t h i s w o r ka l l o w su st ob ea b l et oe v a l u a t et h ei m p r o v e m e n t sn u m e r i c a l l yi nt h e m o d i f i c a t i o n sm a d et ob o t hf o r w a r da n di n v e r s em o d e l i n g s i na d d i t i o n , t h e w o r ka l s oi n t r o d u c e st w oa r t i f i c i a li n t e l l i g e n c eo p t i m i z a t i o nm e t h o d s ，g e n e t i c a l g o r i t h ma n dp a r t i c l es w a l t no p t i m i z a t i o na l g o r i t h m t ot h ed e s i r e dp a c k a g e t oe x p a n dt h ei m a g er e c o n s t r u c t i o np o w e r i ti m p r o v e ss i g n i f i c a n t l yt h e c o n v e r g e n c ec a p a b i l i t ya n dt h er e c o n s t r u c t i o na c c u r a c yc o m p a r e dw i t hr e s u l t s o b t a i n e df r o mt h ea l g o r i t h m sw i t h o u tu s i n ga r t i f i c i a li n t e l l i g e n c em e t h o d k e yw o r d s ：e l e c t r i c a lr e s i s t a n c et o m o g r a p h y , i m a g er e c o n s t r u c t i o n , i m a g e q u a l i t ye v a l u a t i o n ，a r t i f i c i a li n t e l l i g e n c eo p t i m i z a t i o nm e t h o d 独创性声明 本人声明所呈交的论文是我个人在导师指导下进行的研究工作及取 得的研究成果。尽我所知，除了文中特别加以标注和致谢的地方外，论文 中不包含其他人已经发表或撰写过的研究成果，也不包含为获得中国石油 大学或其它教育机构的学位或证书而使用过的材料。与我一同工作的同志 对本研究所做的任何贡献均已在论文中作了明确的说明并表示了谢意。 签名： 关于论文使用授权的说明 本人完全了解中国石油大学有关保留、使用学位论文的规定，即：学 校有权保留送交论文的复印件及电子版，允许论文被查阅和借阅；学校可 以公布论文的全部或部分内容，可以采用影印、缩印或其他复制手段保存 论文。 ( 保密论文在解密后应遵守此规定) 学生签名： 导师签名： ；长彩 诵年s 目z b 明年岁月z 日 中国石油大学( 华东) 硕士论文第1 章绪论 第1 章绪论 两相流多相流体系广泛存在于自然界、工业生产、科学研究、环保 及日常生活等各方面。电阻层析成像技术是一种可实现两相流参数在线检 测的高新检测技术，可应用于连续相具有导电性的两相流参数可视化监 测，本论文将探索新颖的电阻层析成像技术及其在两相流参数检测中的应 用。 随着现代工业和科学技术的发展对两相流参数测量的要求日趋迫切， 分析两相流参数测量技术的发展现状，研究电阻层析成像技术的基本原理 及应用技术是非常有意义的。 1 1 两相流检测方法的背景及现状 1 1 1 两相流的概念及特点 “相”的概念通常是指某一系统中具有相同成分及相同物理、化学性 质的均匀物质部分，各相之间有明显可分的界面。自然界的物质宏观上可 分成三种相，即气相、液相和固相。单相物质的流动称为单相流。两相流 ( t w o p h a s ef l o w ) 或多相流( m u l t i p h a s ef l o w ) ，是指同时存在两种或多 种不同相的物质的流动。常见的有气固、气液、液固、液液两相流， 以及气液液、气固液多相流等i l 】。 在两相流动研究中，根据流动介质的连续与否，可把两相流动介质分 为连续相和离散相，其中连续不断地构成一个整体的相称为连续相( 如气 体、液体等) ，而以不连续介质如固体颗粒、液滴或水泡等形态存在于连 续相中的物质称为离散相( 或非连续相) 。 两相流比单相流具有更广泛的普遍性和实用性，在工业设备、锅炉以 及核反应堆等许多过程或设备内都有两相流存在。其主要有如下特点：流 型复杂多变、两相界面有相间作用力、相间存在着相对速度、物性变化较 大、能耗有增有减、界面有扰动、声速与临界速度不同、有松弛现象、有 电磁效应、数学描述难度大等。由于两相流流动的复杂多样性，使得两相 中国石油大学( 华东) 硕士论文 第1 章绪论 流需要测量的参数很多，包括流动参数( 浓度、流量、压力等) ，传质、传 热情况，流型以及特征参数( 相的粒度分布、相速度和相含率，液膜厚度 分布等) ，并且还常常要求能够给出许多参数的在线实时测量结果，考虑 到测量对象处于复杂的运动状态，又常常伴有传热、传质及化学反应同时 并存的情况。因此，可以认为两相流参数测量具有高难度、多门类交叉的 特点【2 l 。 1 1 2 两相流多相流检测技术的现状 迄今为止，已有的检测技术可分为以下几类【2 】： ( 1 ) 采用传统的单相流仪表和两相流测试模型组合的测量方法。 把成熟的单相流仪表如差压式流量计、涡轮流量计、靶式流量计、容 积式流量计、涡街流量计、电磁流量计、超声波流量计和科里奥利力流量 计等单相流量计应用于两相流的测试，已取得了较大的进展。单相流中已 有的传统的光学、电学、热学等探头和传感器，经改造广泛地应用到两相 流测试系统中来了。采用双( 多) 传感器组合，进行双( 多) 参数组合测量确 定流量或干度等也获得了较多的成功应用。据已发表的研究成果表明这类 组合测量多为组分浓度仪表、速度仪表和流体流动的动量通量仪表的交叉 排列组合应用。 ( 2 ) 使用传统的单相流检测仪表，采用新的信号检测和处理方法( 如噪 声检测理论、模糊评判理论、谱分析技术) 进行两相流流动状态的辨识及 参数估计。 ( 3 ) 应用近代新技术检测两相流参数。 ( 4 ) 在成熟的硬件基础上，以计算机技术为支撑平台，应用基于软测量 技术的软测量方法( 如状态估计、过程参数辨识、入工神经两络和模式识 别等) 。 具体来说，检测两相流参数的主要方法有：相关法，容积法，节流法， 涡轮流量计法，激光多普勒法，p i v 法( p a r t i c l ei m a g ev e l o e i m e t e r ，粒子 成像测速法) ，热线、热膜风速仪，过程层析成像技术( p r o c e s st o m o g r a p h y ， 简称p t ) ，核磁共振法( m a g n e t i cr e s o n a n c el m 硒n g ，简称m m ) ，直接使用 质量流量计测量等i j j 。 2 中国石油大学( 华东) 硕士论文第1 章绪论 1 1 3 两相流检测的意义 随着现代工业的发展及两相流在工业中应用越来越广泛，对生产过程 中两相流的计量、控制要求越来越高，两相流参数的准确测量变得至关重 要。它可揭示两相流流动机理，建立两相流模型并对流动过程预测、设计 和优化控制提供必要的依据，对大型系统生产过程的安全高效运行具有指 导意义。 然而，目前的两相流参数检测技术还远远不能达到工业应用的要求。 迄今为止，现有的测试技术和方法大多处在实验室应用研究阶段，已商品 化可在线检测的工业型仪表为数还很少，因此，开发出一种精度高，实时 性好、成本相对低廉的两相流参数检测仪表是两相流检测中的重要课题。 1 2 过程层析成像技术 1 2 1 过程层析成像技术的概念 过程层析成像( p r o c e s st o m o g r a p h y ，简称p t ) 技术的产生和发展与 科学研究和工程实践中对两相流多相流过程内部信息的迫切需求密切相 关，是2 0 世纪8 0 年代初期，随着计算机技术和检测技术的进步正式形成和 发展起来的一种新型的过程参数检测技术，是以两相流或多相流为主要研 究对象的过程参数二维或三维分布状况的在线实时检测技术，是医学诊断 中的计算机层析成像( c o m p u t e r i z e dt o m o g r a p h y ，简称c t ) 技术在工业 过程中的应用1 4 ”。 1 2 2 过程层析成像技术的组成与特点 p t 系统的构成如图卜l 所示，主要由以下三部分构成【4 】： 、获取被测物场信息的敏感器空间阵列，即传感器； 、传感器控制及数据采集系统； 、图像重建与物场参数提取单元，即图象重建及处理计算机。 中国石油大学( 华东) 硕士论文 第1 章绪论 图卜1 过程层析成像技术的系统构成 虽然p t 技术是医学工程中c t 技术在工业领域中的发展，但由于测量对 象、测量目的以及运行环境的不同，使得p t 技术无论是在信息的获取方式 和信息的处理方法上，还是在测量结果的解释和应用上与c t 技术有着显著 的不同： 面向工业过程，具有更快的数据采集速度和数据处理能力； 电学p t 是一种“软场”( s o f tf i e l d ) ； p t 系统获得的投影数据少； p t 技术的传感器等要装在现场，要有更好的适应能力； p t 技术能根据重建的图像提取有关被测物场的特征信息，并将这些 特征信息反馈给被测过程，实现被测过程的调节和控制。 1 3 电阻层析成像( e r t ) 技术概述 电阻层析成像( e l e c t r i c a lr e s i s t a n c et o m o g r a p h y ，简称e r t ) 技术是 电阻抗层析成像( e l e c t r i c a li m p e d a n c et o m o g r a p h y ，简称e i t ) 技术的一 种简化情况，即只利用了电阻的实部信息。 1 3 1e r t 技术的产生与发展 二十世纪8 0 年代中期形成和发展起来的电阻层析成像技术属于过程 层析技术的一种，是以两相流或多相流为主要对象的过程参数二维或三维 4 中国石油大学( 华东) 硕士论文第1 章绪论 分布状况的在线实时监测技术。所测得的边界电压数据包含了丰富的两相 流场的信息，可以用来提取流型、相含率等特征参数1 6 - 9 1 。 在国外，自上世纪8 0 年代中期兴起的p t 技术因采用非侵入或非接触方 式测量，能在线连续地提供两相流体截面状况的二维可视化信息，并可经 过进一步处理提取若干有关被测两相流体的特征参数，受到了科技工作者 的普遍关注。早在1 9 8 8 年英国曼彻斯特大学( u m i s t ) 理工学院研制成的 8 电极电容层析成像系统，1 9 9 0 年该系统发展成为1 2 电极，配备t r a n s p u t e r 系列高速并行处理器件的电容层析成像系统，在油气水混合流体试验装置 上稳定运行，在线图像重建速率达至l j 4 0 帧s ，并可以1 0 0 次s 的速率更新所 存储的图像象素灰度数据。1 9 9 8 年，u m i s t 大学研制出m k 2 a 系统，其数据 采集速率可达1 0 0 幅秒，数据传输率达于i m b s ，并己在实验室的石油开 采用泥浆测量装置上进行了实验，其结果较为令人满意。2 0 0 1 年，又开 发出利用双极性脉冲电流源作为激励源的电阻层析成像系统，这一技术克 服了常规直流激励带来的介质电极化效应，为解决系统实时性不高提供了 新的途径【o j 。 与此同时，一种安全性好、价格低廉的用于成像人体部位的电阻抗层 析成像技术( e i t ) 引起了医学界科学家的广泛重视。1 9 8 8 年u m i s t 也 开始开发研制对含有导电流体的容器进行成像的e i t 系统。由于医学e i t 和工业过程e r t 存在许多共同之处，因而医学e i t 的技术积累使得工业 过程e i t 技术的研究取得了迅速进展，e r t 技术作为e i t 技术的一种简 化形式也得到了广泛研究。另外，e r t 技术在环境监测中的应用也取得了 成功。 在国内，对e r t 技术的研究开始于8 0 年代后期徐苓安教授带领的天 津大学p t 研究小组。目前国内许多领域都对e r t 技术进行深入的研究： 北京航空航天大学和东北大学的e r t 小组在算法上做了大量仿真实验， 取得了很好的成绩；浙江工学院开发的电阻层析成像系统已用于土壤测 量；天津大学的e r t 科研小组，目前已开发出t e r t - i ，2 ，3 ，4 型系统 样机，其采样速率可达5 0 幅( 4 0 k h z , 1 0 4 点幅，模拟解调) 。其中，t e r t - 4 型成像系统已完成动态试验。 我校目前对e r t 技术的研究还处在初始阶段，尚无成型的科研设备及 5 中国石油大学( 华东) 硕士论文 第l 章绪论 技术支持。在研究生学习期间，我们首先开发了一个基于m a t l a b 的e r t 仿真软件包，为利用e r t 技术进行多相流检测方面提供一个界面友好e r t 仿真平台，并为软硬件的实现奠定了理论数据基础，代表了信息与控制工 程学院检测技术与自动化装置专业在多相流检测方面又开拓出一个新的 方向。 1 3 2e r t 技术实现的系统结构与特点 如图1 2 所示：e r t 系统包括用于激励测量的电极阵列、数据采集与 处理单元、图像重建与分析显示单元三部分 7 1 。 图1 - 2 e r t 系统的构成 e r t 系统的电极阵列由特殊设计的电极等间隔排布，控制单元( 计算 机) 向数据采集单元发出指令，给某一对电极施加激励电流，在过程对象 内部建立起敏感场。测量边界上的电压信号，将得到的测量数据送图像重 建单元，以适当的算法重建出对象内部的电阻率分布，从而得到媒质分布 图像( - - 维或三维) 。最后送图像分析单元，对图像的物理意义加以解释， 提取有关的特征参数，为过程控制或实验研究提供必要的依据。 e r t 技术主要有以下几个特点： 被检测物场的连续相必须是导电的，一般必须是含有水的生产过程； 敏感阵列为非侵入式，由一系列等间隔排布的电极构成。敏感阵列的 设计对于敏感场的性能有直接的影响。电极的设计、排布以及电极的 性能是整个系统至关重要的一部分； 设计敏感电极阵列时，必须考虑电极阵列所形成的空间敏感场的非均 6 中国石油大学( 华东) 硕上论文第1 章绪论 匀性的影响，应尽可能使其灵敏度均匀性好； 敏感场的激励信号为低频的恒流源。频率不能过高或者过低，过低容 易引起电极的电化学反应；过高，电磁感应及分布阻抗带来很强的测 量噪声； 敏感场是“软场。，敏感场的分布受到场内介质分布的影响，敏感场与 物场的相互作用为非线性，导致图像重建与分析困难： 检测信号为微弱的电压信号或其微小变化，要求测量电路必须具有高 的灵敏度和信噪比。 1 3 3e r t 技术存在的主要问题和发展趋势 e r t “软场”问题的克服，通过对敏感场分布特性的认识与改进，尽 可能提高敏感场灵敏度分布的均匀性【s 】； 提高数据采集精度和增大数据采集量； 研究有效降低不适定性、快速、准确的图像重建算法； 在实现二维图像的基础上向三维图像转化，由于e r t 系统的敏感场 是三维非均匀分布的，简单地用二维场近似地分析必然引起误差，这 就要求必须研究三维场内的敏感场分布、图像重建算法等等吧 研究过程参数的提取方法和系统性能评价准则。 1 4e r t 技术的工作原理 由于两相流体各相组分具有不同 的电阻率( 由于电导率与电阻率互为倒 数，所以这两个概念在研究e r t 技术时 可以互换使用) ，因此电阻率分布与相 的组分分布一一对应。如图1 3 所示， 在两相流工艺装置或管道某一截面的 边界等间隔布置多个电极( 这些电极在 工作时同时分别作为激励电极和测量 电极) ，轮流给某一相邻的电极对施加 7 电 激励电极依次激励 图1 - 3 e r t 测量原理图 v i 中国石油大学( 华东) 硕士论文第1 章绪论 激励电流i s 作为输入信号，在过程对象内部建立起敏感场，再测量其它电 极对上的输出信号v j ( 电压大小) 。当场内电阻率分布变化时，电流 场的分布会随之变化，导致场域内电势分布变化，从而场域边界上的测量 电压也要发生变化。利用边界上的测量电压，通过一定的成像算法，可以 重建出场内的电阻率分布或反映电阻率分布情况的灰度分布，实现可视化 测量。 1 5 本课题的主要研究任务 本研究课题的任务是用m a t l a b 编程建立一个e r t 仿真软件包，本课题 的软件包针对二维问题，主要研究层流、核心流、环形流等流型。课题内 容包括： ( 1 ) 实现两相流场域的有限元剖分、计算边界电压值等e r t 正问题的求 解，并通过与a n s y s 的仿真结果比较，验证基于m a t l a b 的仿真结果的正确性， 再根据正问题数据画出场域内的电势分布图； ( 2 ) 实现e r t 反问题的应用，包括线性反投影算法( l b p ) 、修正的牛 顿拉夫逊类算法( m n r ) 等基本的图像重建算法，并详细分析各自 的优缺点； ( 3 ) 在图像重建算法中引入遗传算法与粒子群算法等人工智能算法，并 将其与经典的图像重建算法结合，形成一种成像精度很高的组合型两相流 成像算法； ( 4 ) 讨论信息量( 测量电压数量) 、噪声等对图像重建算法仿真结果的 影响。并提出一种合理的量化指标对各种成像结果进行评判； ( 5 ) 通 过_ m a t l a b 的图形用户界面( g u i ) 建立一个简单的e r t 仿真平台。 8 中国石油大学( 华东) 硕士论文第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 2 1e r t 技术的数学物理模型 通常的，如果给定电阻分布p ( x ，j ，) 和边界电势分布妒( x ，力，求物体内 部的电势分布u ( x ，力，这就是d r i c h l e t 边值问题： ( d ) v 伊砚) 2 0 ，( 薯y ) 刨 ( 2 _ 1 ) 、7 【甜= 妒( x ，j ，) ，( x ，y ) 硷 如果给定电阻分布p ( x ，力和边界电流密度j ( x ，力，求物体内部的电势 分布u ( x ，y ) ，这就是n e u m a n n 边值问题： iv ( p 。1 v u ) = o ，( z ，y ) q p 0 u ( m y g 小施 q 。2 电磁场的基本方程是麦克斯韦( m a x w e l l ) 方程，在电阻层析成像系 统中也同样满足。由于求解一般形式的麦克斯韦方程是相当困难的，所以 常常根据实际情况来简化方程，对实际的电阻层析成像系统作两点假设： 、电阻层析成像系统的场域遵循敏感场( 静电场) 的规律，即不考虑电 磁场在空间上的时间滞后。当测量管道的直径小于六分之一电磁波波长 时，这种假设可以成立。一般情况下，激励源的频率选择小于1 0 0 k h z ， 电磁波长为旯v f = 3 x 1 0 3 m ，实际e r t 系统直径小于5 米，远小于 2 6 = 5 0 0 ，因此敏感场条件成立，可以将e r t 传感器的敏感场看成是静 电场。 、研究的敏感场域内没有激励源频率范围内的电流源，从而敏感场的各 处电流散度为零。 根据以上的假设，由m a x w e l l 方程，该静电场域归结为第二类边界求 解问题，即n e u m a n n 问题，可描述为如下方程i l 】： 9 中国石油大学( 华东) 硕士论文第2 章电阻层析成像( e 盯) 的数学模型及其正问题求解 v c o v ) = 0 盯娑i 。：， 功 妒( ，y 。) = 0 在q 内 ( x o ，y o ) q ( 2 3 ) ( 2 - 4 ) ( 2 5 ) 其中，矶和v 分别为散度和梯度算子，q 为物体所在的空间区域，仃 为待求的电导率分布，a q 为其边界，妒为未知电势分布函数，j 为边界电 流密度分布函数，取q 内一点4 9 0 ，y o ) 作为参考节点，其电势妒( ，) 为 零。 2 2e r t 技术的实现路线 e r t 技术的实现分为正问题( f o r w a r dp r o b l e m ) 的求解与反问题 ( i n v e r s ep r o b l e m ) 的实现两部分，如图2 1 所示。 图2 - 1e r t 正问题与反问题 e r t 正问题即是由已知电阻率分布p ( 或电导率分布o ) 及敏感场的 初始条件( 输入激励电流i o ) 求解敏感场内任意一点的电势分布由的 过程。 e r t 反问题则是由已知边界的测量电压v 及初始条件( 输入激励电流 1 0 ) 求解敏感场内的电阻率分布p ( 或电导率分布o ) ，即是获得物场 介质分布图像的图像重建算法。 本章将详细论述e r t 正问题的求解方法。 1 0 中国石油大学( 华东) 硕士论文 第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 2 3e r t 的正问题求解 正问题的解法大体分为两种，一种是直接解法，即解析解；另一种是 目前求解正问题最常用的有限元法( f i n i t e e l e m e n t m e t h o d ，简称f e m ) 。 2 3 1 正问题的解析解 采用分离变量法可求取两种媒介分布情况下e r t 正问题解析解1 1 。 在均匀介质中，电导率o 在整个场域内为常数，敏感场满足式( 2 3 ) ， 为表述方便，采用极坐标系。设圆形场域半径为r ，激励电流强度i ，盯 为场内介质电导率，口为极角，p 为极径，为电极所对圆心角( 采取弧 度形式) ，易、西分别为注入电流电极、流出电流电极与坐标轴正实轴的 夹角，如图2 - 2 所示： 图2 2 二维均匀敏感场不葸图 则场内任意一点( p ，0 ) 处的电势 妒( 胪) = 型o z 篙尹k r ) s i n 警脚( o - q ) - - c o s 一( 口一岛) ( 2 6 ) 若已知激励电流i ，则电极上电流密度_ ，= i ，于是电位分布可以进一 a 【 步写成 蚍咖羔萋( 盯s i n 拿c o s ”( o - b ) - - c o s n ( 0 - 岛) 沿z ， 很显然，以上仅是给出了一种简单情况下e r t 正问题的解析解，对 于电阻率分布复杂的情况，则很难求其解析解，而只能利用数值方法，例 如下面介绍的有限元法求解。 中国石油大学( 华东) 硕士论文 第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 2 3 2 正问题的有限元解 有限元法( f e m ) 的基本思想是将连续的 求解域离散为有限个单元的组合体，这样的组 合体能较好地模拟或逼近求解域，随着单元数 目的增加，亦即单元尺寸的缩小，解的近似程 度将不断改进。如果单元满足收敛要求，近似 解最后将收敛于精确f e m 。分析步骤一般包 括： 区域剖分； 选取单元内的场变量插值函数； 单元特性矩阵及单元节点方程形成； 总矩阵的组装及整体方程的形成； 解线性方程组求出场变量值。 有限元法求解e r t 正问题的流程图如图 2 3 所示，下面将详细介绍应用于e r t 中的有 限元法 1 1 - 1 3 1 。 在二维e r t 成像区域，一般将区域剖分 为m 个三角单元，m 的大小取决于阻抗图像 i 有限元剖分 k 一。；。 障元矩阵计算r i 总刚矩阵计算y y 谚 i 谚= y - 1a qi 边界电压值v ” + c j c j ： 图2 3 有限元法流程图 空间分辨率，受限于外围边界的电极个数e ，为了给e r t 反问题提供足 够的信息量，一般应满足m s n = e ( e 一1 ) 2 ，n 为独立测量的次数。 图2 4 有限元剖分及其三角单元 中国石油大学( 华东) 硕上论文第2 章电阻层析成像( e rr ) 的数学模型及其芷问题求解 如图2 - 4 所不，取出剖分区域的一个单兀，对于三结点三角形单兀， 共有- - 个自由度，其阻抗为p e ，当剖分单元数n l 足够多时，可选用简单 的线性插值函数来近似求解场变量电位的分布，插值多项式设为1 1 1 1 矿( 并，y ) = p o + 届x + 岛y ( 2 8 ) 将三节点处的坐标值代入( 2 8 ) 式，得三节点的电位值办、办、九， 旃2 妒“，y 1 ) = p o + 届而+ , a 2 y , 以2 矿也，儿) = p o + 届而+ 厦y 2 ( 2 - 9 ) 办2 似，乃) = p o + 届而+ 岛y 3 联立上面三式，得： p o2 寺如磊+ 口2 ：+ q 珐) 届2 西1 ( 6 1 破+ b 2 z + b 3 庐3 ) ( 2 - 1 0 ) 展2 寺“办z + c 3 九) 其中，a 为三角形的面积，a 、b 、c 可通过轮换下标1 、2 、3 求得： 彳：圭5 屹x l 羔l ：圭( q + 口：+ 口，) c 2 t ， 2 工2 乃一而儿日2 2 x 3 m 一而乃 h = y 2 一y 3b z = y 3 一y l c l2 x ，一x 2c 2 。x l x 3 妒k y ，= 石g ，y 溉+ 五k y 耽+ 石k j ，玩= 【，；k y x 五g ，y x 工g ，y 4 墨 j s ， ，g ，y ) 2 寺o ，+ 螺y ) | 扛1 、2 3 ( 2 )屹( y 屯胜而 一 一 一 儿以砭 而= i i = 以岛 吩 中国石油大学( 华东) 硕士论文第2 章电阻层析成像( e r r ) 的数学模型及其正问题求解 = 去纠鼍鼍+ 蒡号卜 计算得单元导纳矩阵元素： 。去以o + q 勺l _ 1 ，2 ，3 融潮= 目 其中，办( f = 1 ，2 ，3 ) 为节点电位，q ( 扛1 ，2 ，3 ) 为流入节点电流。 y 为n x n 维矩阵，、c 为n x l 维矩 阵，n 为网格中的总节点个数。 叠加总体刚度矩阵的方法可以简 单描述如下，设某一区域被剖分为4 个三角形单元，共有5 个结点，如图 2 5 所示，其中每个单元的单元矩阵局 部结点编号用( 1 ) 、( 2 ) 、( 3 ) 表示。 对于、单元，在按式( 2 1 1 ) 计算单元面积时，各单元节点按逆时针 编号顺序代入，保证单元面积计算结果 为正( 如果采用顺时针编号，则在计算 ( 2 1 5 ) ( 2 1 6 ) ( 2 1 7 ) ( 2 1 8 ) 图2 - 5 总刚矩阵整合示意图 面积时加绝对值，确保面积计算结果为正值) ；而在计算单位刚度矩阵时， 由于矩阵本身的性质决定，无论按哪种顺序编号得到的矩阵，总可以等效 的转换成编号按从小到大顺序得到的单位刚度矩阵，因此，4 个单元的单 位刚度矩阵分别为 1 4 中国石油大学( 华东) 硕士论文 第2 章电阻层析成像( e r r ) 的数学模型及其正问题求解 离 慝y110) 3 1 2 ( i ) 款2 ( i ) ) 1 3 2 ( 1 ) y 1 4 ( 3 ) 几4 ( 3 ) 如( 3 ) y 1 3 ( 1 ) y 2 3 0 ) y 3 3 ( i ) y 1 5 ( 3 ) y 4 5 0 ) y 5 5 ( 3 ) 医 医 y 1 3 ( 2 ) 乃3 ( 2 ) y 4 3 ( 2 ) y j s ( , o y 5 爿4 ) 咒5 ( 4 ) ( 2 1 9 ) 其中矩阵每个元素的下标括号内数字表示第几个单元，例如下标“1 2 ( 4 ) ” 表示第4 个单元中的对应于总节点的编号“l ”、“2 ”。 将第四个单元矩阵等效变换，得 y l l ( 1 ) 呦) 嘲) l 咒l ( 1 ) 耽2 0 ) ，奴1 ) l 乃i 【1 ) 乃2 ( 1 ) 乃划 k 3 ) y l 。( 3 ) 确 1 只1 ( 3 ) ) ，叫3 ) 儿 l 乃t ( 3 ) 乃4 ( 3 ) 儿印) k ：) 咒 l 乃h 2 ) 乃3 ( 2 ) 【儿1 ( 2 ) y 叫2 ) k 。) 峋。) i 乃1 ( 4 ) 耽2 ( 4 ) 【乃1 ( 4 ) 弘2 ( 4 ) ( 2 - 2 0 ) 由于该区域共有5 个结点，因此总体刚度矩阵y 是5 5 矩阵，分别 将各个单元刚度矩阵扩展成5 5 矩 咒1 1 3 ) 0 0 y 4 1 ( 3 ) 欺1 ( 3 ) y 1 2 ( i ) 咒2 ( t ) 乃2 ( 1 ) 0 0 00 00 00 00 o 0 0 0 00 00 00 00 y 1 5 ( 3 ) 0 0 以h 3 ) 乃m 3 ) y 1 1 ( 2 ) 0 y 3 1 ( 2 ) y 4 1 ( 2 ) 0 y l i ( 4 ) y 2 1 ( 4 ) 0 0 乃1 ( 4 ) 阵，得 0 咒m 2 ) oo 0 乃x 2 ) 0 y 4 3 ( 2 ) 00 乃2 【4 ) 0 y 2 2 ( 4 ) 0 oo 0o 乃2 ( 4 ) 0 1 5 咒“2 ) 0 00 y 3 4 ( 2 ) 0 2 ) 0 00 0 咒爿4 ) o 咒5 ( 4 ) 00 o0 0 乃5 ( 4 ) ( 2 2 1 ) 对 0 町 q m 脚眦 肋埘脚 1j_1j 1j 1，j o o o 砷 对 对 m 脚脚o o m o o m 脚 刖脚聊o o 中国石油大学( 华东) 硕士论文第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 y = 再将这4 个扩展后的5 5 矩阵相加，即得到总体刚度矩阵y 。 y j l o ) + 咒l ( 2 ) + 7 h o ) + 一l ( 砷y , 2 ( o + 7 1 2 ( 4 )7 1 3 0 ) + m 玎4咒q 习+ m 4 ( 3 )y 1 5 0 ) + m x ql y 2 1 ( o + 咒l ( 砷y 叫1 ) + ，瓴4 ) 如( 1 ) 0 ，葡4 )l 7 3 1 ( 1 ) + y 3 l ( 4乃砷)7 3 3 ( 1 ) + y 3 3 ( 2 )y34(2)01 只峭+ 乃聃0 ) 瓴2 ) ) 钆2 ) + 只y 4 5 0 )i 咒1 p ) + 乃l f 4 ，y 5 2 ( 4 ) 0 钮97 5 s o ) + y 5 5 0 ) j ( 2 - 2 2 ) 在仿真计算过程中，可以将上述方法适当简化，先设总体刚度矩阵为 n x n 的零阵，然后按照各单元刚度矩阵内的元素编号直接写入总刚矩阵 即可。 最后，在得到总刚矩阵y 以后，代入己知条件注入电流c ( 电流流入 的节点处为+ 1 ，电流流出的节点处为1 ，其余节点的电流均为o ) ，规定坐 标原点处的节点为参考零点( 巾1 = 0 ) ，即去掉y 矩阵的第一行第一列、 巾矩阵的第一行以及c 矩阵的第一行，得 y ( n - i ) x ( ) 砧一l m2 - 1 阳 ( 2 2 3 ) 进而得到正问题中除去原点其余各节点的电势值 砧一 2 ：t m ) q 一- ( 2 - 2 4 ) 其中，誓：。卜卜。) 的求解一般有直接法或者迭代法两种，在仿真计算中，根 据誓m ) 是对称正定矩阵的特性，利用c h o l e s k y 分解或l u 分解，得 m 。) = 三u ，畋。卜。= u 、( 三c 0 一，阳) ，再叠加上原点的电势值巾i = o ， 求得所有节点的电势 伽怯。 ( 2 2 5 ) m 。1 即为一次激励下的节点电势，为了求得边界电极间的电压v ，给定 ( e x n ) 的转换矩阵t e 。，该矩阵为每一行只有两个元素分别为1 和一1 ， 分别对应当前用于测量的两电极节点编号，其余均为0 的稀疏矩阵，有 e m 0 ) = 7 赢唬。l ( e 一一电极数，n m - 节点数) ( 2 - 2 6 ) 即求得第一次激励下边界测量电压矩阵唾2 ，将每次激励下计算得到的测 1 6 中国石油大学( 华东) 硕士论文 第2 章电阻层析成像( e r t ) 的数学模型及其正问题求解 量电压矩阵叠加，得 。= l 唾2 ，吃： ，疃0i ( 2 - 2 7 ) 根据对称原理，l 、2 号电极激励时测量3 、4 号电极间的电压，与3 、 4 号电极激励时测量1 、2 号电极间的电压相等，该电压矩阵。是一个 对称矩阵。再除去1 、2 号电极激励时l 、2 号电极间的电压，有效测量电 压值即为矩阵k 。的下三角矩阵( 或上三角矩阵) ，有效测量电压数为 e x ( e i ) 2 。但由于矩阵。的次对角线上的电压值在测量时离激励电 极很近，导致其数值比最小电压值大几十倍以上，一旦在测量中引入噪声， 会对算法的收敛性产生过大影响，详细论述参见2 4 1 节讨论( 2 ) ，因此 有时也舍弃不用，此时的有效测量电压数= e x ( e - 3 ) 2 。 2 4e r t 正问题的仿真结果与讨论 e r t 正问题的仿真结果是否正确将直接关系到后面求