《電阻抗斷層成像的求解方法分析》4700字

電阻抗斷層成像的求解方法分析綜述目錄TOC\o"1-2"\h\u10363電阻抗斷層成像的求解方法分析綜述 1312731.1EIT技術(shù)的正問題 15666Σ?φ?n=J (2-3) 1306861.1.1EIT正問題的基本實(shí)現(xiàn)原理 2142391.2.2有限元法求解EIT 4306352.2EIT技術(shù)的逆問題 5284282.2-1EIT逆問題 55251.2.2EIT的圖像重建算法 5160141.2.3牛頓-拉夫遜算法 6178361.2.4正則化修正技術(shù) 81.1EIT技術(shù)的正問題EIT研究具有特殊邊界條件的電場(chǎng)，其數(shù)學(xué)模型可以等效于Laplace方程： ??σ??=0 (2-1)加強(qiáng)邊界條件：?=?0 自然邊界條件： Σ?φ?n=J 其中，?是物體內(nèi)部的電壓分布，σ為電導(dǎo)率，?0和J分別為邊界區(qū)域的電壓和電流密度REF_Ref1122394767\n\h[13].根據(jù)驅(qū)動(dòng)信號(hào)和測(cè)量模型的阻抗分布計(jì)算內(nèi)部電流和電壓分布，這在EIT的研究中統(tǒng)稱為正問題(ForwardProblem)，為了解決正問題，可以通過求解拉普拉斯方程得到該區(qū)域的節(jié)點(diǎn)電壓REF_Ref1119739261\n\h[12]，然后利用已知的條件計(jì)算出內(nèi)部電流密度，從而達(dá)到綜合分析電場(chǎng)的目的.當(dāng)前科學(xué)研究團(tuán)隊(duì)中，用于解決正問題的數(shù)值方法是有限元法，有限體積元法，邊界元法等.有限元法(FiniteElementMethod，F(xiàn)EM)，通過將離散和連續(xù)求解區(qū)域簡(jiǎn)化為有限數(shù)量的元素節(jié)點(diǎn)，該方法適用于分割復(fù)雜字段，該方法更全面，更靈活.邊界元法(BoundaryElementMethod，BEM)適用于開域問題.用高斯定理降階復(fù)雜高階問題，使得待求解方程被轉(zhuǎn)化成邊界積分方程，從而將邊界細(xì)分成有限個(gè)單元，這種方法由于Gauss定理，保證了精確度，減少了計(jì)算工作量REF_Ref745565020\n\h[9].有限體元法(FiniteVolumeMethod，F(xiàn)VM)，它將場(chǎng)域劃分為大量的主單元元素，并圍繞單元主元素的每個(gè)節(jié)點(diǎn)構(gòu)建一個(gè)封閉的控制體，稱為次單元元素.電流通量由主單元中每個(gè)節(jié)點(diǎn)的電勢(shì)來計(jì)算和表示.這樣就可以形成一個(gè)矩陣方程來求解每個(gè)節(jié)點(diǎn)的電勢(shì).該方法精度高，計(jì)算量少，可以容易地解決電導(dǎo)率突然變化的精確模擬問題.圖1.2.1-1EIT正問題EIT問題無法通過解析方法解決.這與許多其他非線性逆問題相同.必須借助迭代計(jì)算.首先，假設(shè)電導(dǎo)率分布為已知的感應(yīng)區(qū)域，解決正問題以獲取場(chǎng)邊界處的電勢(shì)分布.根通過基于已知邊界條件迭代地校正假定的電導(dǎo)率分布，并進(jìn)行無限次，直到正問題的結(jié)果的誤差小于某個(gè)值，也就是說，此時(shí)的正問題解已經(jīng)被認(rèn)定為收斂到滿足已知邊界條件，從而獲得導(dǎo)電率的分布情況，然后在利用現(xiàn)有的圖像重建算法對(duì)求解出的結(jié)果進(jìn)行成像處理.雖然已經(jīng)能將這個(gè)過程通過有限元法變形，將式(2-1)轉(zhuǎn)化成求解線性方程組，但是，正問題求解過程迭代一直是最大的問題，它會(huì)花費(fèi)很長(zhǎng)的時(shí)間.有限元法求解過程：開始，場(chǎng)域剖分，取得各種信息(節(jié)點(diǎn)總數(shù)、強(qiáng)加邊界條件節(jié)點(diǎn)編號(hào)、單元節(jié)點(diǎn)編號(hào)、三角元總數(shù)以及單元電導(dǎo)率)，形成系數(shù)矩陣，構(gòu)成有限元方程，邊界條件、激勵(lì)模式，求解方程，結(jié)束REF_Ref1128932690\n\h[14].1.1.1EIT正問題的基本實(shí)現(xiàn)原理EIT技術(shù)的主要研究問題是EIT正問題計(jì)算，EIT逆問題計(jì)算和硬件系統(tǒng)設(shè)計(jì)REF_Ref1151807017\n\h[15].EIT正向問題的數(shù)學(xué)模型主要從麥克斯韋的電磁場(chǎng)方程開始，推導(dǎo)它，然后使用有限元方法求解.現(xiàn)有研究已經(jīng)證實(shí)，生物組織中的載體材料具有非均勻分布的各向異性，并且生物組織的電學(xué)性質(zhì)不僅僅是電阻，而是電阻和電容，這證實(shí)了生物組織中存在著虛部信息.電磁場(chǎng)中的電流和電導(dǎo)率分布之間基礎(chǔ)約束關(guān)系可通過麥克斯韋方程組進(jìn)行求解表達(dá)： ??L=ρ (2-4) ?×E=?B?t ??B=0 (2-6) ?×M=J+?D?t其中： B=μM ( D=εE (2-9) I=σE (2-10)電位移——L.電場(chǎng)強(qiáng)度——E.磁感應(yīng)強(qiáng)度——B，磁場(chǎng)強(qiáng)度——M.電流密度——I.電荷密度——ρ。介質(zhì)的磁導(dǎo)率——μ.介質(zhì)的介電常數(shù)——?.介質(zhì)電導(dǎo)率——σ.由于電阻抗成像問題中敏感場(chǎng)特性包含實(shí)部信息和虛部信息，當(dāng)向場(chǎng)域內(nèi)加入頻率為??的電流激勵(lì)時(shí)，其電場(chǎng)強(qiáng)度為： E=E可得： ?×M=σE+ε?E?t又因?yàn)椋篍=?φ是場(chǎng)內(nèi)電勢(shì)的分布，故上式可以變形為： ?×M=?γ?Φ (2-13)將上式兩端取散度得到的Laplace方程如下： ??(γ?Φ=0) (2-14)其中：γ=σ(x,y)+jωε(x,y)，γ是敏感長(zhǎng)內(nèi)復(fù)電導(dǎo)率.一旦確定了邊界條件，它也可以由上式唯一地確定：對(duì)于EIT問題，由于激勵(lì)源是電流源，因此滿足Neumann邊界條件，即，給定邊界. γ??Φ?其中，n為邊界外法線向量，S為場(chǎng)域邊界，I為電流密度.1.2.2有限元法求解EIT當(dāng)給定了σ和在實(shí)際的測(cè)量中，因?yàn)榇骐娢粊頊y(cè)量物體上各點(diǎn)的電壓，所以以O(shè)點(diǎn)為基準(zhǔn)點(diǎn)，其電勢(shì)為0.如果將其他邊界點(diǎn)和點(diǎn)O之間的電位差設(shè)定為電位值.將點(diǎn)O的坐標(biāo)設(shè)為(x0，y0)，則測(cè)得邊界電勢(shì)分布函數(shù)φ(x,y)(首先將問題化為變分形式.記Ω=Ω∪?Ω，假設(shè)C’C`(Ω)={u|u及uU=u(x,y)是方程(N’)的解，則對(duì)任意的v=v(x,y)∈ ν?(??ρ?1即??(v在Ω上積分后再利用高斯公式可得：Ωρ?1?ν??udxdy=Ω??(νρ?1?u)dxdy記D(u,ν)=Ωρ?1?ν??udxdy，G(ν)=(G)求u∈在原始問題中，需要u的二次偏導(dǎo)數(shù)的存在，但是在積分形式中，積分和格林公式的應(yīng)用大大減少了解的平滑度，ux面對(duì)一般的u∈C2.2EIT技術(shù)的逆問題1.2.1EIT逆問題的基本實(shí)現(xiàn)原理2.2-1EIT逆問題EIT的逆問題是一個(gè)病態(tài)的非線性問題，由圖2.2-1以及2.1節(jié)公式2-1，2-2，2-3可知，逆問題就是給出邊界的電壓測(cè)量數(shù)值和激勵(lì)電流反解模型內(nèi)的電阻抗分布，即由?0，J求解ρ1.2.2EIT的圖像重建算法EIT最關(guān)鍵的一環(huán)，就是最終的成像，成像的好壞往往決定著它的使用價(jià)值.EIT圖像重建算法主要是根據(jù)邊界測(cè)量值來確認(rèn)測(cè)量體內(nèi)部電導(dǎo)率分布情況，隨著國(guó)際上對(duì)EIT技術(shù)的關(guān)注，EIT圖像重建算法逐漸成為EIT技術(shù)中最頻繁的研究方向.當(dāng)前，智能算法和非智能算法是EIT圖像重建算法的兩個(gè)主要類別.非智能算法通常通過使用線性模型模擬逆EIT問題來解決，該問題可大致分為迭代算法和非迭代算法.非迭代類算法不需要經(jīng)歷迭代來優(yōu)化目標(biāo)函數(shù)就可以獲得EIT逆問題的解.1985年，Muria等人提出了一種基于Geselowitz的靈敏度算法，以通過靈敏度矩陣對(duì)線性方程進(jìn)行求解，從而簡(jiǎn)化了問題.1990年，Cheny等人創(chuàng)新出基于最小二乘法的單步牛頓殘差重構(gòu)法.1991年，Somersalo等人提出層剝法.1998年，Vauhkonen提出了對(duì)應(yīng)阻抗分布先驗(yàn)假設(shè)的Tikhonov正則化矩陣的構(gòu)造法，等等.迭代算法是通過更新參數(shù)的連續(xù)迭代來使目標(biāo)函數(shù)達(dá)到最小值的算法的一種REF_Ref1161487849\n\h[16].1987年，Yorkey提出了改進(jìn)Newton-Raphson法，并提出用標(biāo)準(zhǔn)方法和補(bǔ)償定理計(jì)算雅可比矩陣.1991年，Huda率領(lǐng)眾人基于改進(jìn)的Newton-Raphson法，提出一種先驗(yàn)信息融入圖像重建的正則化法，等等.而在具體的實(shí)現(xiàn)里，眾多科學(xué)家都在致力于研究出一個(gè)有效的求解器.目前，EIT的大多數(shù)研究小組使用的是在MATLAB中電阻抗層析成像和漫反射光學(xué)層析成像重建軟件EIDORS.EIDORS提供了一組有用的性質(zhì)，例如二維和三維正演模擬，以及一組廣泛的重建算法、可視化功能等等.Horesh等人使用不同的預(yù)處理?xiàng)l件和更有效的例程對(duì)EIDORS進(jìn)行了調(diào)整，產(chǎn)生了一個(gè)稱為SuperSolver的版本，目前仍在UCL的小組中使用.然而，對(duì)于大網(wǎng)格，MATLAB缺乏高效的并行編程能力，這使得計(jì)算正解成為一項(xiàng)冗長(zhǎng)的任務(wù).而后，Borsic等人將Jacobian計(jì)算(但不是系統(tǒng)矩陣的集合)轉(zhuǎn)移到稀疏并行直接求解器庫(kù)PARDISO中，以超越這些限制.與Horesh等人相比，他們能夠?qū)⒄蚰M的速度提高約5.3倍.他們把它用在含有大約50萬個(gè)元素的網(wǎng)格上.在較大的網(wǎng)格上，直接求解需要大量?jī)?nèi)存，這些內(nèi)存通常會(huì)限制可計(jì)算的網(wǎng)格大小.此外，在本文中，直接求解器的裝配速度要比一個(gè)好的預(yù)處理器慢得多，這使得迭代方法的執(zhí)行速度更快，而這取決于唯一電流注入模式的數(shù)目.特別地，代數(shù)多重網(wǎng)格預(yù)處理程序已經(jīng)被證明可以顯著地提高求解時(shí)間.基于圖形處理單元(GPU)的計(jì)算已經(jīng)成功地應(yīng)用于對(duì)內(nèi)存的快速訪問至關(guān)重要的Jacobian矩陣的計(jì)算中.1.2.3牛頓-拉夫遜算法當(dāng)使用牛頓迭代法作為EIT圖像重建算法時(shí)，它基于最優(yōu)化理論的思想，以重建模型邊界的電壓的計(jì)算值與測(cè)量得到的邊緣電壓的值兩者間差的二范數(shù)最小，即得到一個(gè)合適的靜態(tài)電導(dǎo)率分布REF_Ref1164244197\n\h[17]σ?，使得(2-18)最小化 f(σ)=1其中σ作為電導(dǎo)率分布是一個(gè)非線性函數(shù).F(σ):Cn→假設(shè)T(σ)=F(σ)?V，由于電阻抗是非線性的，有牛頓迭代法，在點(diǎn)( T(σ+h)=T(σ)+T`采用線性最小二乘法，找到適合的h，讓目標(biāo)函數(shù)達(dá)到最小，即T(σ+h)≈0， T(σ+h)=T(σ)+T`(σ)h≈0求解： h=?(T`(σ))?1T同時(shí)： h=σk+1?σk=[F即： σk+1=σk+[F`基于之前的理論，更異于使趨近于0，也就是可以使得目標(biāo)函數(shù)達(dá)到最小值.而通常就是Jacobian矩陣，定義為[F`由上述推導(dǎo)可知，最終目標(biāo)是最小化的F(σ)?V，但是，在展開之后，不計(jì)算高階導(dǎo)數(shù)項(xiàng)，從而引起較大的誤差，導(dǎo)致最終結(jié)果的質(zhì)量.所以公式(2-18)是最小化的研究對(duì)象，令?f(σ)= ?f(σ)=對(duì)f(σ) f``(σ)=F`(σ)(2-25)中的F`(σ)TF將公式(2-18)在(σ+h)處展開成泰勒級(jí)數(shù)，可得： f(σ+h)=f(σ)+f`求上式的梯度并將其設(shè)為0，就可以使其最小化，可得： h??f=求解得： h=(F`(σ)T綜合上述推導(dǎo)，可得： h=(F在上述的推導(dǎo)中，本文還是做了簡(jiǎn)化計(jì)算，例如不計(jì)算高階變量，在實(shí)際應(yīng)用中，由于Hessian陣的條件數(shù)一般都是比較大的，所以這樣的處理會(huì)使得誤差較大，對(duì)最終結(jié)果影響也很大.因此，還要引入正則化法來修正.1.2.4正則化修正技術(shù)Hessian矩陣的廣義逆條件的數(shù)目較大，導(dǎo)致該解的病態(tài)需要一些額外的先驗(yàn)和解的附加約束條件，為了將不適當(dāng)?shù)膯栴}轉(zhuǎn)變?yōu)檫m當(dāng)?shù)膯栴}，可以通過歸一化方法獲得穩(wěn)定的近似解并進(jìn)行校正.常見的歸一化算法包括Tikhonov歸一化和Marquatdt方法.Marquatdt方法是通過引入和歸一化單位矩陣來實(shí)現(xiàn)的.此參數(shù)減少了Hessian矩陣中的廣義逆條件的數(shù)量，以減少不良條件.本文主要利用Tikhonov正則化規(guī)則來修改Newton-Raphson算法，這也是目前使用的最成熟，最常用的修改方法，Tikhonov正則化的基本原理是最小二乘問題.在實(shí)際中，求解方法是在目標(biāo)函數(shù)中加入一定的約束元素，在一定程度上提高了解的穩(wěn)定性.把最小化函數(shù)看作： Z(a,L)=argmin||h(Z)?d|(2-30)中，Z0為參數(shù)Z的初始值，α是修正參數(shù)，F(xiàn)(Z)=||L(Z?Z0 Z(a,L)(k+1)=Z(a,L)通過上述公式，原來求Jaccobian矩陣Jk的逆的問題轉(zhuǎn)化成了求JkTJk所以將這種修正方法使用到逆問題求解中，對(duì)f(σ) f(σ)=||F(σ)?V||其中λ為正則化參數(shù)，R是一個(gè)調(diào)制矩陣，它包含一些關(guān)于電