1、<p><b>　　本科畢業(yè)設(shè)計(jì)</b></p><p><b>　?。?0 屆）</b></p><p>　　有限差分法在細(xì)胞介電譜中的應(yīng)用</p><p>　　所在學(xué)院 </p><p>　　專(zhuān)業(yè)班級(jí) 數(shù)學(xué)與應(yīng)用數(shù)學(xué)

2、 </p><p>　　學(xué)生姓名 學(xué)號(hào) </p><p>　　指導(dǎo)教師 職稱(chēng) </p><p>　　完成日期 年 月 </p><p><b>　　摘　要</b></p&g

3、t;<p>　　【摘要】交流電作用于細(xì)胞體系產(chǎn)生介電弛豫效應(yīng)。文中將該系統(tǒng)等效成電容與電阻并聯(lián)而成的模型，根據(jù)復(fù)介電常數(shù)虛部的定義（即表示在電場(chǎng)變化一個(gè)周期中的能量損失），通過(guò)等效電路中的電阻及相關(guān)參數(shù)來(lái)計(jì)算。在等效模型中，電容在電場(chǎng)變化一個(gè)周期所儲(chǔ)存的能量復(fù)介電常數(shù)實(shí)部有關(guān)，建立相應(yīng)函數(shù)關(guān)系，通過(guò)所給的相關(guān)參數(shù)計(jì)算復(fù)介電常數(shù)實(shí)部。利用此電磁學(xué)原理建立三維模型，以有限差分法對(duì)生物細(xì)胞介電弛豫進(jìn)行模擬。最后將模擬結(jié)果與PS方

4、程對(duì)比。對(duì)比結(jié)果的相對(duì)誤差為，的相對(duì)誤差為。</p><p>　　【關(guān)鍵詞】有限差分法；介電弛豫；PS方程；復(fù)介電常數(shù)；復(fù)電導(dǎo)率。</p><p><b>　　Abstract</b></p><p>　　【ABSTRACT】While applying an alternating current electric field on biol

5、ogical cells,the cells produce dielectric relaxation effect.In this paper,this system is equivalent to the model which includes capacitance and resistors in parallel.According to the definition of the imaginary part of c

6、omplex permittivity which means the energy loss in a cycle,we can calculate the imaginary part by resistance in equivalent circuit and parameters.The real part of complex permittivity is related to energy chang</p>

7、<p>　　【KEYWORDS】Finite difference method;Dielectric relaxation;PS equation;Complex permittivity ;Complex conductivity.</p><p><b>　　目　錄</b></p><p><b>　　摘　要I</b><

8、;/p><p>　　AbstractI</p><p><b>　　目　錄II</b></p><p><b>　　1引言1</b></p><p>　　1.1細(xì)胞介電弛豫效應(yīng)1</p><p>　　1.2三維有限差分法1</p><p>

9、　　2介電基本理論2</p><p>　　2.1常用符號(hào)2</p><p>　　2.2宏觀電性質(zhì)3</p><p>　　2.2.1真空介電常數(shù)3</p><p>　　2.2.2電導(dǎo)率3</p><p>　　2.2.3界面極化3</p><p>　　2.2.4電磁學(xué)方程

10、3</p><p>　　3三維有限差分模型與算法5</p><p>　　3.1三維有限差分模型5</p><p>　　3.1.1模型原理5</p><p>　　3.1.2三維有限差分模型6</p><p>　　4三維有限差分模型的算法9</p><p>　　4.1三維差分模

11、型算法9</p><p>　　4.1.1建立位圖和導(dǎo)納元胞數(shù)組矩陣9</p><p>　　4.2超松弛迭代法求解電勢(shì)10</p><p>　　4.2.1Jacobi迭代10</p><p>　　4.2.2Seidel迭代11</p><p>　　4.2.3松弛迭代法11</p>&l

12、t;p>　　4.3超松弛程序框圖12</p><p>　　5單球殼模型計(jì)算12</p><p>　　5.1單球殼模型13</p><p>　　5.1.1單球殼模型13</p><p>　　5.1.2其他模型14</p><p>　　5.2算法改進(jìn)14</p><p>

13、;　　5.2.1改進(jìn)理論14</p><p>　　5.2.2改進(jìn)結(jié)果15</p><p>　　6結(jié)果與討論15</p><p>　　6.1Cole-Cole理論15</p><p>　　6.1.1Cole-Cole理論發(fā)展15</p><p>　　6.2PS理論16</p><

14、;p>　　6.2.1PS理論發(fā)展16</p><p>　　6.2.2基于參數(shù)PS求解16</p><p>　　6.3基于參數(shù)三維有限差分法求解17</p><p>　　6.4三維有限差分法與PS方程理論結(jié)果對(duì)比18</p><p>　　6.5三維有限差分法與PS方程理論三維結(jié)果對(duì)比21</p><

15、p><b>　　7結(jié)束語(yǔ)23</b></p><p><b>　　參考文獻(xiàn)24</b></p><p>　　致謝錯(cuò)誤!未定義書(shū)簽。</p><p><b>　　引言</b></p><p><b>　　細(xì)胞介電弛豫效應(yīng)</b></p&g

16、t;<p>　　細(xì)胞普遍存在著介電弛豫效應(yīng)[1](dielectric dispersion),即細(xì)胞的復(fù)介電常數(shù)和復(fù)電導(dǎo)率隨外加電場(chǎng)頻率變化發(fā)生的介電響應(yīng)。研究細(xì)胞介電弛豫可以了解細(xì)胞界面極化的性質(zhì)和細(xì)胞導(dǎo)電、絕緣的電生理學(xué)特性，有助于人們驗(yàn)證生物細(xì)胞的數(shù)學(xué)模型和物理模型的正確性，深人了解生物細(xì)胞的結(jié)構(gòu)和功能。通常利用計(jì)算機(jī)控制的寬頻帶阻抗分析儀非損傷地測(cè)量細(xì)胞的介電常數(shù)和電導(dǎo)率,然后對(duì)計(jì)算機(jī)采集的數(shù)據(jù)進(jìn)行多種圖形處理和

17、非線(xiàn)性數(shù)值計(jì)算分析，建立生物細(xì)胞的介電參數(shù)。</p><p><b>　　三維有限差分法</b></p><p>　　有限差分方法(FDM：Finite Difference Method)是計(jì)算機(jī)數(shù)值模擬最早采用的方法，至今仍被廣泛運(yùn)用。該方法將求解域劃分為差分網(wǎng)格，用有限個(gè)網(wǎng)格節(jié)點(diǎn)代替連續(xù)的求解域。三維有限差分法以Taylor級(jí)數(shù)展開(kāi)等方法，把控制方程中的導(dǎo)數(shù)用網(wǎng)

18、格節(jié)點(diǎn)上的函數(shù)值的差商代替進(jìn)行三維空間離散，從而建立以網(wǎng)格節(jié)點(diǎn)上的值為未知數(shù)的代數(shù)方程組[2]。該方法是一種直接將微分問(wèn)題變?yōu)榇鷶?shù)問(wèn)題的近似數(shù)值解法，表達(dá)簡(jiǎn)單，數(shù)學(xué)概念直觀，是發(fā)展較早且比較成熟的數(shù)值方法。</p><p><b>　　介電基本理論</b></p><p><b>　　常用符號(hào)</b></p><p>　　

19、介電常數(shù)： 介電常數(shù)：</p><p>　　復(fù)介電常數(shù)實(shí)部： 復(fù)介電常數(shù)虛部：</p><p>　　電導(dǎo)率： 復(fù)電導(dǎo)率：</p><p>　　復(fù)電導(dǎo)率實(shí)部：

20、 復(fù)電導(dǎo)率虛： </p><p>　　真空中介電常數(shù)： 測(cè)定角頻率：</p><p>　　所加交流電場(chǎng)頻率： 虛數(shù)單位：</p><p>　　直流電壓： 交流電壓：</p><p>　　電容器電容：

21、 所帶電荷：</p><p>　　大立方體變成： 細(xì)胞體積分?jǐn)?shù)：</p><p>　　相鄰網(wǎng)格點(diǎn)間距： 細(xì)胞半徑：</p><p><b>　　細(xì)胞膜厚度為：</b></p><p>

22、　　網(wǎng)格點(diǎn)元素坐標(biāo)： 網(wǎng)格點(diǎn)的電勢(shì)： </p><p>　　溶液復(fù)電導(dǎo)率實(shí)部： 細(xì)胞質(zhì)復(fù)電導(dǎo)率實(shí)部：</p><p>　　細(xì)胞膜復(fù)電導(dǎo)率實(shí)部： 溶液復(fù)介電常數(shù)實(shí)部：</p><p>　　細(xì)胞質(zhì)復(fù)介電常數(shù)實(shí)部：

23、 細(xì)胞膜復(fù)介電常數(shù)實(shí)部： </p><p>　　流經(jīng)和邊界電流： 流經(jīng)和邊界電流： 流經(jīng)和邊界電流： 流經(jīng)和邊界電流： 流經(jīng)和邊界電流： 流經(jīng)和邊界電流：</p><p>　　在和之間元素導(dǎo)納： 在和之間元素導(dǎo)納：</p><p>　　在和

24、之間元素導(dǎo)納： 在和之間元素導(dǎo)納：</p><p>　　在和之間元素導(dǎo)納： 在和之間元素導(dǎo)納：</p><p><b>　　宏觀電性質(zhì)</b></p><p>　　一般的，材料的宏觀電性質(zhì)可以由介電常數(shù)和電導(dǎo)率來(lái)表征。</p><p><b>　　真空介電常數(shù)</b><

25、/p><p>　　將直流電壓加在一真空平行板電容器上，在極板上產(chǎn)生一定數(shù)量的電荷 ，此時(shí)真空電容器的電容為：</p><p><b>　　（2.1）</b></p><p>　　而電容大小與電壓無(wú)關(guān)，只與其幾何結(jié)構(gòu)有關(guān)，每個(gè)極板的面積為，極板間距為，則</p><p><b>　　（2.2）</b>&l

26、t;/p><p>　　比例系數(shù)稱(chēng)為真空介電常數(shù)，其值為[3]。</p><p>　　當(dāng)電容器中充滿(mǎn)電介質(zhì)時(shí)，極板上的電荷增加到，電容器電容比真空狀態(tài)下增加了倍，即</p><p><b>　?。?.3）</b></p><p>　　比例系數(shù)稱(chēng)為相對(duì)介電常數(shù)，表征電介質(zhì)儲(chǔ)存電能的能力大小。</p><p&g

27、t;<b>　　電導(dǎo)率</b></p><p>　　表征溶液傳導(dǎo)電流的能力，在直流電導(dǎo)中主要由質(zhì)子、離子和電子傳導(dǎo)機(jī)制引起，電導(dǎo)通過(guò)介電損失來(lái)發(fā)生，介電損失是一種交流電現(xiàn)象。</p><p><b>　　界面極化</b></p><p>　　當(dāng)材料中的自由電荷（正、負(fù)離子或電子）在兩種材料的相界面上或者在一種材料內(nèi)兩個(gè)不同

28、區(qū)域積聚時(shí)，將形成材料中空間電荷的不均勻分布，產(chǎn)生宏觀偶極矩，這種極化稱(chēng)為界面。界面極化在一定程度上可以等效地看成偶極矩取向極化，界面極化主要存在于具有相界面的不均勻材料和具有缺陷、顆粒和雜質(zhì)的材料中錯(cuò)誤!未找到引用源。。界面極化機(jī)制如圖2.1所示。</p><p>　　圖2.1 界面極化機(jī)制示意圖</p><p><b>　　電磁學(xué)方程</b></p>

29、<p>　　細(xì)胞在交變電場(chǎng)作用下產(chǎn)生極化，由于存在介電質(zhì)致使極化過(guò)程產(chǎn)生相位差，而此時(shí)電束密度也產(chǎn)生相應(yīng)的相位差[15]。為方便推導(dǎo)，交流電場(chǎng)用復(fù)數(shù)表示：</p><p><b>　?。?.4）</b></p><p><b>　　電束密度為：</b></p><p><b>　?。?.5）</

30、b></p><p>　　若定義復(fù)介電常數(shù)為：</p><p><b>　?。?.6）</b></p><p><b>　　則，</b></p><p>　　這與通過(guò)直流電時(shí)形式一致。在復(fù)介電常數(shù)中實(shí)部表示在靜態(tài)時(shí)介電常數(shù)。</p><p><b>　　通過(guò)系

31、統(tǒng)的電流</b></p><p><b>　　（2.7）</b></p><p>　　參考靜電場(chǎng)中形式，定義復(fù)電導(dǎo)率：</p><p><b>　?。?.8）</b></p><p>　　所以復(fù)介電常數(shù)和復(fù)電導(dǎo)率關(guān)系：</p><p><b>　?。?.

32、9）</b></p><p>　　三維有限差分模型與算法 </p><p><b>　　三維有限差分模型</b></p><p>　　有限差分法在數(shù)值計(jì)算中應(yīng)用廣泛，本文將其推廣應(yīng)用于三維空間模擬，對(duì)細(xì)胞介電參數(shù)進(jìn)行計(jì)算分析。在之前三維有限差分法已經(jīng)應(yīng)用于計(jì)算混合系統(tǒng)的等效介電常數(shù)[6][7]，在這里將其推廣至覆蓋有原生質(zhì)膜的生物細(xì)

33、胞[8]。</p><p><b>　　模型原理 </b></p><p>　　對(duì)細(xì)胞介電性質(zhì)研究包括復(fù)介電常數(shù)和復(fù)電導(dǎo)率，復(fù)介電常數(shù)和電導(dǎo)率的關(guān)系為：</p><p><b>　　（3.1）</b></p><p>　　而復(fù)電導(dǎo)率的計(jì)算公式為：</p><p><b&

34、gt;　　（3.2）</b></p><p>　　其中為電導(dǎo)率，為相對(duì)介電常數(shù)，在實(shí)驗(yàn)中可有相應(yīng)實(shí)驗(yàn)儀器測(cè)出，對(duì)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行計(jì)算可得相應(yīng)的細(xì)胞介電特性。</p><p>　　將包含細(xì)胞的培養(yǎng)液進(jìn)行極化時(shí)，把系統(tǒng)看成平行板電容器與電阻并聯(lián)的系統(tǒng)。電容器的極板面積為，極板間距為，而流經(jīng)整個(gè)系統(tǒng)的電流：</p><p>　　圖一.系統(tǒng)的等價(jià)電路</p&g

35、t;<p><b>　?。?.3）</b></p><p>　　在該模型中，電容器在交流電場(chǎng)中儲(chǔ)存能量：</p><p><b>　?。?.4）</b></p><p>　　而復(fù)介電常數(shù)的虛部表示電場(chǎng)變化一個(gè)周期的能量損失，主要以電阻產(chǎn)生焦耳熱的損失。</p><p><b>

36、;　?。?.5）</b></p><p>　　根據(jù)計(jì)算公式，在此建立三維有限差分模型對(duì)細(xì)胞進(jìn)行數(shù)字模擬，從而進(jìn)行細(xì)胞介電特性分析，對(duì)細(xì)胞結(jié)構(gòu)分析、細(xì)胞介電弛豫分析提供可靠的依據(jù)。在對(duì)模型進(jìn)行數(shù)字分析，計(jì)算出相應(yīng)的參數(shù)和與之實(shí)驗(yàn)進(jìn)行對(duì)比。</p><p><b>　　三維有限差分模型</b></p><p>　　細(xì)胞置于三維平行板電容

37、器中，將所在區(qū)域三維分割實(shí)現(xiàn)差分網(wǎng)格，在每個(gè)差分網(wǎng)格點(diǎn)都位于立方體元素的中心，其分割圖如圖（3.1）所示。整個(gè)三維系統(tǒng)擁有個(gè)網(wǎng)格點(diǎn)，每個(gè)網(wǎng)格點(diǎn)相距為，每個(gè)網(wǎng)格所在的立方體元素體積為。</p><p>　　圖3.1差分網(wǎng)格圖 圖3.2網(wǎng)格點(diǎn)</p><p>　　對(duì)于每個(gè)網(wǎng)格點(diǎn)來(lái)說(shuō)，其電勢(shì)可由周?chē)鶄€(gè)網(wǎng)格點(diǎn)的電勢(shì)和所在的立方體的復(fù)電導(dǎo)率計(jì)算得到。如圖3.2

38、所示，為網(wǎng)格所在點(diǎn)，為各網(wǎng)格點(diǎn)的電勢(shì)，灰色部分為細(xì)胞內(nèi)部（復(fù)電導(dǎo)率），白色部分為細(xì)胞外部培養(yǎng)液（復(fù)電導(dǎo)率），兩者相交處為細(xì)胞膜（復(fù)電導(dǎo)率，厚度）。</p><p>　　假設(shè)一個(gè)網(wǎng)格點(diǎn)的立方體元素處于均勻平行電場(chǎng)中，電場(chǎng)位于網(wǎng)格點(diǎn)和之間，流過(guò)立方體元素和邊界的電流的計(jì)算公式為[9]：</p><p><b>　?。?.6）</b></p><p>

39、;　　上式中為處于立方體元素和的導(dǎo)納，其計(jì)算公式為：</p><p><b>　?。?.7）</b></p><p>　　同理流過(guò)網(wǎng)格點(diǎn)的其他五個(gè)邊界的電流和相應(yīng)的導(dǎo)納計(jì)算為：</p><p><b>　?。?.8）</b></p><p><b>　?。?.9）</b><

40、;/p><p><b>　　（3.10）</b></p><p><b>　?。?.11）</b></p><p><b>　　（3.12）</b></p><p><b>　?。?.13）</b></p><p><b>　

41、?。?.14）</b></p><p><b>　?。?.12）</b></p><p><b>　　（3.13）</b></p><p><b>　?。?.14）</b></p><p>　　在網(wǎng)格點(diǎn)處電流的流進(jìn)和流出量相同，其矢量和為零，因此得出的一個(gè)差分方程：

42、</p><p><b>　　（3.15）</b></p><p><b>　　整理得：</b></p><p><b>　　（3.16）</b></p><p>　　當(dāng)有細(xì)胞膜處于立方體元素和之間,導(dǎo)納的計(jì)算可有下式表示:</p><p><b&

43、gt;　　(3.17)</b></p><p>　　對(duì)有限差分劃分的每個(gè)網(wǎng)格電流分析，得出相應(yīng)的差分方程，聯(lián)立方程得到關(guān)于電勢(shì)的差分方程組，通過(guò)超松弛法求解。將相應(yīng)的電勢(shì)代入方程（3.4），（3.5）求得相應(yīng)的介電特性參數(shù)。</p><p>　　三維有限差分模型的算法 </p><p><b>　　三維差分模型算法</b></

44、p><p>　　在計(jì)算三維差分模型時(shí)，首先得建立原始位圖以便于導(dǎo)納的計(jì)算，其次計(jì)算各個(gè)網(wǎng)格點(diǎn)電勢(shì)，求出相應(yīng)電勢(shì)后計(jì)算系統(tǒng)總電流，最后對(duì)細(xì)胞復(fù)電導(dǎo)率和復(fù)介電常數(shù)進(jìn)行求解，分析相應(yīng)的特性。</p><p>　　建立位圖和導(dǎo)納元胞數(shù)組矩陣</p><p>　　將細(xì)胞置于立方體內(nèi)進(jìn)行劃分差分網(wǎng)格，網(wǎng)格點(diǎn)整個(gè)系統(tǒng)內(nèi)要么在細(xì)胞內(nèi)部，或在細(xì)胞外部的培養(yǎng)液中。通過(guò)簡(jiǎn)化模型，處于細(xì)胞內(nèi)和

45、培養(yǎng)液之間的立方體元素為細(xì)胞膜，注意計(jì)算時(shí)導(dǎo)納的區(qū)分。</p><p>　　建立原始位圖矩陣，其內(nèi)部元素只有或，判斷標(biāo)準(zhǔn)為網(wǎng)格點(diǎn)處于細(xì)胞內(nèi)部為，否則為，對(duì)位圖矩陣計(jì)算可得出相應(yīng)細(xì)胞占溶液的體積分?jǐn)?shù)。</p><p>　　對(duì)導(dǎo)納計(jì)算公式由式（3.4），（3.10）-（3.14），（3.17）可得，我們由原始位圖矩陣分析可計(jì)算出每個(gè)網(wǎng)格點(diǎn)相應(yīng)的個(gè)導(dǎo)納，建立導(dǎo)納矩陣。</p>&l

46、t;p>　　考慮邊界問(wèn)題，我們以網(wǎng)格點(diǎn)為例進(jìn)行處理以便于統(tǒng)一導(dǎo)納矩陣元素的格式，網(wǎng)格點(diǎn)的示意圖如圖3.3所示。</p><p>　　圖3.3 元素(1,1,1)位置</p><p>　　網(wǎng)格點(diǎn)的上部不存在與之相鄰的網(wǎng)格點(diǎn)元素，我們將其令為，故，</p><p>　　同理，，以此為原則對(duì)系統(tǒng)內(nèi)所有網(wǎng)格點(diǎn)元素導(dǎo)納進(jìn)行計(jì)算，網(wǎng)格點(diǎn)的導(dǎo)納在Matlab中以元胞數(shù)組進(jìn)

47、行存儲(chǔ)：。每個(gè)網(wǎng)格點(diǎn)內(nèi)擁有六個(gè)導(dǎo)納即元胞內(nèi)擁有六個(gè)元素為：，，，，，。所以網(wǎng)格點(diǎn)的導(dǎo)納為：。對(duì)系統(tǒng)的導(dǎo)納以上述規(guī)則計(jì)算并存儲(chǔ)，以備求解電勢(shì)所用。</p><p>　　超松弛迭代法求解電勢(shì)</p><p>　　在均勻介質(zhì)中，電流流過(guò)導(dǎo)體時(shí)其電勢(shì)呈現(xiàn)線(xiàn)性關(guān)系。由于在稀薄溶液中，細(xì)胞對(duì)電勢(shì)影響不大，我們假設(shè)在兩塊電極板之間的電勢(shì)變化為線(xiàn)性梯度，即從電極板正極到負(fù)極呈單調(diào)遞減形，以此建立初始電勢(shì)矩

48、陣。</p><p><b>　　設(shè)原始方程組為：</b></p><p><b>　?。?.1）</b></p><p><b>　　以矩陣形式表示為：</b></p><p><b>　?。?.2）</b></p><p>　　

49、其中 為系數(shù)矩陣,非奇異,是右端項(xiàng),為解向量。</p><p><b>　　假設(shè),令</b></p><p><b>　　(4.3)</b></p><p>　　則方程(3.22)可整理為:</p><p><b>　　(4.4)</b></p><p>

50、;<b>　　Jacobi迭代</b></p><p>　　令 ， ， （4.5）</p><p>　　所以方程式（4.4）可表示為：</p><p><b>　?。?.6）</b></p><p>　　選取初始矩陣，代入式（4.4）右端可得一組新向量，記為：.如此迭代下去，迭代

51、格式為：</p><p>　?。?.7）以上迭代為簡(jiǎn)單迭代，即Jacobi迭代。</p><p><b>　　Seidel迭代</b></p><p>　　基于上述Jacobi迭代法則，假設(shè)先采用代入式（4.4）第一個(gè)方程右端，計(jì)算出，然后用新算出的第一個(gè)分量代替，用代入式（4.4）第二個(gè)方程右端，計(jì)算出，用代入式（4.4）第三式右端計(jì)算，如此

52、類(lèi)推直到所有變量用取代。</p><p><b>　　迭代過(guò)程為：</b></p><p><b>　?。?.8）</b></p><p><b>　　令</b></p><p>　　， （4.9）</p>&l

53、t;p>　　迭代過(guò)程用用矩陣表示：</p><p><b>　　（4.10）</b></p><p><b>　　將上式整理得：</b></p><p><b>　?。?.11）</b></p><p>　　上述迭代稱(chēng)為Seidel迭代法，這種迭代既可節(jié)省存儲(chǔ)空間又給編程

54、帶來(lái)方便[10]。</p><p><b>　　松弛迭代法</b></p><p>　　對(duì)Seidel迭代法關(guān)系處理，令</p><p>　?。?.12）于是我們可得 （4.13） <

55、/p><p>　　對(duì)Seidel迭代法來(lái)說(shuō)，是加上一修正項(xiàng)而得到的。若在修正項(xiàng)前加入一參數(shù)便得到松弛法的計(jì)算公式：</p><p><b>　?。?.14）</b></p><p>　　其中稱(chēng)為松弛因子，當(dāng)時(shí)稱(chēng)為超松弛，當(dāng)時(shí)即為Seidel迭代法[11]。</p><p>　　定理1：對(duì)任何初始向量,和常數(shù)項(xiàng),有迭代格式&l

56、t;/p><p>　　產(chǎn)生向量序列收斂且極限與初值無(wú)關(guān)的充分必要條件為（為矩陣譜半徑.）</p><p>　　收斂最快的松弛因子稱(chēng)為最優(yōu)松弛因子1,記為,其計(jì)算公式為:</p><p><b>　　(4.15)</b></p><p><b>　　超松弛程序框圖</b></p><p

57、><b>　　單球殼模型計(jì)算</b></p><p><b>　　單球殼模型</b></p><p>　　細(xì)胞本身是含有細(xì)胞膜和細(xì)胞質(zhì)的非均勻體系，細(xì)胞膜含有蛋白質(zhì)磷脂雙分子層，它具有很低的離子通透性，為非導(dǎo)電的外殼。而細(xì)胞質(zhì)具有導(dǎo)電性，圖過(guò)細(xì)胞質(zhì)在電性質(zhì)上均一的，那么電模型就是一個(gè)被絕緣性的殼包裹的到點(diǎn)球，稱(chēng)為單球殼模型[12]。<

58、/p><p><b>　　單球殼模型</b></p><p>　　單球殼模型主要用來(lái)解釋細(xì)胞懸浮液中源于細(xì)胞膜的界面極化而產(chǎn)生的色散。細(xì)胞膜內(nèi)部的細(xì)胞質(zhì)同時(shí)包含細(xì)胞器和蛋白質(zhì)等大分子成份，也是個(gè)非均勻體系，因此細(xì)胞質(zhì)的電性質(zhì)也依賴(lài)測(cè)量時(shí)的頻率，應(yīng)該產(chǎn)生弛豫現(xiàn)象。 </p><p>　　細(xì)胞的單球殼模型如圖4.1所示，</p><

59、;p>　　圖4.1 細(xì)胞單球殼模型示意圖</p><p>　　球殼半徑為，細(xì)胞膜厚度為，將復(fù)介電常數(shù)分別為和的細(xì)胞質(zhì)和細(xì)胞膜構(gòu)成的球殼</p><p>　　模型視為均一體，其等效復(fù)介電常數(shù)為：</p><p><b>　?。?.1）</b></p><p><b>　?。?.2）</b><

60、;/p><p>　　根據(jù)Wagner理論，該粒子以體積分?jǐn)?shù)分散到復(fù)介電常數(shù)為的介質(zhì)中所構(gòu)成的懸浮液的復(fù)介電常數(shù)為：</p><p><b>　?。?.3）</b></p><p>　　通常情況下， ，以及的條件滿(mǎn)足，高頻弛豫條件的強(qiáng)度很小，幾乎可忽略，而且 介電色散可以近似地看成單一弛豫時(shí)間表示。</p><p><b

61、>　　其他模型</b></p><p>　　對(duì)于具有細(xì)胞膜外基質(zhì)的細(xì)胞，例如菌類(lèi)、酵母菌等，用單球殼模型不能很好解釋其相應(yīng)的介電性質(zhì)，此時(shí)將細(xì)胞壁和連續(xù)的介質(zhì)視為具有相同電性質(zhì)，提出帶細(xì)胞壁的細(xì)胞模型[13]。</p><p>　　該模型如圖5.2所示：</p><p>　　圖5.2 細(xì)胞雙球殼模型示意圖</p><p>　

62、　為整個(gè)球體的等價(jià)復(fù)介電常數(shù)，為細(xì)胞壁的復(fù)介電常數(shù)，為細(xì)胞壁厚度，為整個(gè)球體的外半徑。整個(gè)球體的復(fù)介電常數(shù)：</p><p><b>　?。?.5）</b></p><p><b>　?。?.6）</b></p><p><b>　　算法改進(jìn)</b></p><p>　　以單球

63、殼形為例進(jìn)行計(jì)算，在進(jìn)行劃分差分網(wǎng)格時(shí)，由于該系統(tǒng)為三維系統(tǒng)，劃分次數(shù)與總網(wǎng)格點(diǎn)數(shù)呈立方級(jí)上升。</p><p><b>　　改進(jìn)理論</b></p><p>　　我們考慮計(jì)算機(jī)儲(chǔ)存能力，對(duì)之前算法中的導(dǎo)納以元胞數(shù)組進(jìn)行儲(chǔ)存。對(duì)迭代方程式（3.16）研究，其系數(shù)矩陣為：</p><p>　　對(duì)迭代系數(shù)矩陣觀察得，共有個(gè)元素，但是每行最多只有六個(gè)

64、非零元素，系數(shù)陣為稀疏陣，為減少存儲(chǔ)量，設(shè)計(jì)壓縮存儲(chǔ)方式進(jìn)行存儲(chǔ)，格式見(jiàn)表圖5.3 。</p><p>　　圖5.3壓縮存儲(chǔ)格式</p><p>　　空格放位置方程的系數(shù)，中放置，以此類(lèi)推到存放；</p><p>　　空格放位置方程中未知數(shù)在行向量中的位置（對(duì)應(yīng)系數(shù)的相應(yīng)位置），中放置在行向量中的位置，以此類(lèi)推到。當(dāng)處在立方體邊界位置時(shí)，對(duì)應(yīng)不存在的等的位置設(shè)定為，

65、對(duì)應(yīng)的系數(shù)值為。</p><p><b>　　改進(jìn)結(jié)果</b></p><p>　　按照?qǐng)D5.3壓縮存儲(chǔ)格式，在計(jì)算中不僅節(jié)省了存儲(chǔ)空間，也提高了運(yùn)算效率。將改進(jìn)存儲(chǔ)算法與原算法進(jìn)行對(duì)比，結(jié)果見(jiàn)表5.4 。</p><p>　　表5.4 壓縮存儲(chǔ)與原算法存儲(chǔ)對(duì)比表</p><p>　　有上表分析得，當(dāng)劃分差分格數(shù)越多，壓

66、縮存儲(chǔ)的優(yōu)勢(shì)越明顯。在理論上當(dāng)時(shí)，原算法僅迭代系數(shù)矩陣所占內(nèi)存為，目前對(duì)于普通的PC來(lái)說(shuō)達(dá)不到計(jì)算的要求。</p><p><b>　　結(jié)果與討論</b></p><p>　　Cole-Cole理論</p><p>　　在對(duì)生物細(xì)胞介電譜研究中，K. S. Cole 和 R. H. Cole于1941年提出的經(jīng)驗(yàn)公式[14]，為細(xì)胞介電弛豫研究

67、提供介電特性的方法。</p><p>　　Cole-Cole理論發(fā)展</p><p>　　如果在交流電場(chǎng)作用下，可得到頻率依存的平均偶極矩：</p><p><b>　?。?.1）</b></p><p>　　進(jìn)而得到復(fù)介電常數(shù)的頻率依存性：</p><p><b>　?。?.2）<

68、;/b></p><p>　　其中為靜態(tài)介電常數(shù)，光學(xué)介電常數(shù)，為：</p><p><b>　　（6.3）</b></p><p>　　滿(mǎn)足式（5.2）的弛豫現(xiàn)象稱(chēng)為Debye型弛豫[15]。</p><p>　　K. S. Cole 和 R. H. Cole在1941年提出的經(jīng)驗(yàn)公式，由于電導(dǎo)率包括直流電導(dǎo)率，

69、損耗因子通過(guò)從中減去低頻限制電導(dǎo)率計(jì)算得到（見(jiàn)方程（6.2））</p><p><b>　　(6.4)</b></p><p><b>　　(6.5)</b></p><p>　　其中是的實(shí)部，為弛豫時(shí)間， 為Cole—Cole參數(shù)，特征頻率定義為。該曲線(xiàn)擬合已經(jīng)通過(guò)最小二乘法算出，該計(jì)算方法最大限度的使誤差的平方和最小，

70、方差的計(jì)算如下式：</p><p><b>　　(6.6)</b></p><p>　　其中是第個(gè)角頻率，下標(biāo)從（6.4）式中估計(jì)得到。</p><p><b>　　PS理論</b></p><p>　　在單球殼模型中，Pauly 和 Schwan將 式（5.1）和式（5.3）聯(lián)立稱(chēng)為PS方程[12

71、]。這也是對(duì)單球殼模型介電性質(zhì)研究的理論方程，對(duì)于其他方法對(duì)單球殼模型介電性質(zhì)研究起著基準(zhǔn)的作用。</p><p><b>　　PS理論發(fā)展</b></p><p>　　P．S理論是針對(duì)稀薄形系統(tǒng)提出，此后Hanai提出將方程（5.1）代入Hanai方程[16]：</p><p><b>　　（6.7）</b></p

72、><p>　　將式（5.1）和式（6.7）聯(lián)立，稱(chēng)為HAK方程。該方程是針對(duì)細(xì)胞高濃度下（），PS方程不適用的情形提出的，在一定程度對(duì)PS的繼承和發(fā)展。</p><p><b>　　基于參數(shù)PS求解</b></p><p>　　在參數(shù)下[7]，對(duì)PS進(jìn)行計(jì)算，可以求出在特定體積分?jǐn)?shù)下細(xì)胞介電弛豫，其復(fù)相對(duì)介電常數(shù)實(shí)部，復(fù)電導(dǎo)率實(shí)部如圖6.1所示。&

73、lt;/p><p>　　圖6.1 PS方程復(fù)相對(duì)介電常數(shù)實(shí)部和復(fù)電導(dǎo)率實(shí)部</p><p>　　有上圖分析可得，對(duì)于復(fù)介電常數(shù)實(shí)部在頻率從時(shí)，其趨勢(shì)為平緩。從理論上來(lái)說(shuō)，在低頻段時(shí)細(xì)胞的內(nèi)部電荷活動(dòng)不大，對(duì)于復(fù)電導(dǎo)率分析同樣可得電導(dǎo)性變化緩慢，在低頻段影響較小。當(dāng)頻率從時(shí)，電性質(zhì)變化不大，此時(shí)在高頻段細(xì)胞壁和膜已經(jīng)擊穿，整個(gè)細(xì)胞的電性質(zhì)和培養(yǎng)液一樣。在頻率從時(shí)，細(xì)胞在電場(chǎng)作用下其電荷發(fā)生定向移

74、動(dòng)，出現(xiàn)極化作用，細(xì)胞兩頭電壓變化明顯致使細(xì)胞在高頻段發(fā)生擊穿，此時(shí)的電特性為過(guò)渡階段隨著頻率的改變而發(fā)生變化。</p><p>　　基于參數(shù)三維有限差分法求解</p><p>　　在參數(shù)下，對(duì)PS進(jìn)行計(jì)算，可以求出在特定體積分?jǐn)?shù)下細(xì)胞介電弛豫，其復(fù)相對(duì)介電常數(shù)實(shí)部，復(fù)電導(dǎo)率實(shí)部如圖6.2所示。</p><p>　　圖6.2 PS方程復(fù)相對(duì)介電常數(shù)實(shí)部和復(fù)電導(dǎo)率實(shí)部

75、</p><p>　　將細(xì)胞進(jìn)行數(shù)字模擬其介電性質(zhì)的變化如圖所示，在頻率從時(shí)細(xì)胞介電性質(zhì)基本不變，其頻率對(duì)細(xì)胞產(chǎn)生影響很小，當(dāng)頻率從時(shí)，同樣電性質(zhì)變化不大，此時(shí)在高頻段細(xì)胞壁和膜已經(jīng)擊穿，整個(gè)細(xì)胞的電性質(zhì)和培養(yǎng)液一樣。</p><p>　　三維有限差分法與PS方程理論結(jié)果對(duì)比</p><p>　　我們分析細(xì)胞在理論方程下和數(shù)字模擬情況系，其整體介電性質(zhì)基本相同，所以

76、將其對(duì)照比較，其圖形見(jiàn)圖6.3。 </p><p>　　圖6.3 PS方程與3D-FDM模擬對(duì)照?qǐng)D</p><p>　　在上圖中，實(shí)線(xiàn)表示PS理論方程圖，圓圈表示3D-FMD數(shù)字模擬的圖。對(duì)細(xì)胞復(fù)介電常數(shù)實(shí)部對(duì)照，其吻合程度很高。在頻率從時(shí)細(xì)胞介電性質(zhì)基本不變，外加電壓頻率對(duì)細(xì)胞產(chǎn)生影響很小，當(dāng)頻率從時(shí)，同樣電性質(zhì)變化不大，此時(shí)在高頻段細(xì)胞壁和膜已經(jīng)擊穿，整個(gè)細(xì)胞的電

77、性質(zhì)和培養(yǎng)液一樣。對(duì)于細(xì)胞的電導(dǎo)率在頻率時(shí)出現(xiàn)一段不吻合情況，3D-FMD模擬細(xì)胞在此時(shí)出現(xiàn)上翹情況，而PS理論曲線(xiàn)則是平緩。對(duì)此結(jié)果分析，當(dāng)細(xì)胞在高頻段時(shí)，細(xì)胞幾乎完全電解，細(xì)胞的結(jié)構(gòu)已經(jīng)打破可以近似看成和培養(yǎng)液的電性質(zhì)相同，這是理論上的分析[17]。但實(shí)際上細(xì)胞在極化后還存在著相應(yīng)的細(xì)胞器，此時(shí)的細(xì)胞培養(yǎng)液的電導(dǎo)性和之前細(xì)胞未極化存在差別，而3D-FMD模擬在這里很好的模擬出來(lái)這樣的情況。</p><p>　

78、　我們分析兩種方法的電性質(zhì)相似性可驗(yàn)證數(shù)字模擬很高準(zhǔn)確和吻合度，下圖6.4為復(fù)相對(duì)介電常數(shù)實(shí)部，復(fù)電導(dǎo)率實(shí)部的誤差圖：</p><p>　　圖6.5 前圖為的誤差圖，后圖為的誤差圖</p><p>　　對(duì)模擬數(shù)據(jù)和理論數(shù)據(jù)的誤差進(jìn)行計(jì)算，的平均誤差為，其中最大誤差點(diǎn)誤差為，的平均誤差為，其中最大誤差點(diǎn)誤差為。從誤差上來(lái)看，數(shù)字模擬3D-FDM的效果和PS方程吻合度很高。</p>

79、<p>　　細(xì)胞介電特性中，在復(fù)平面上復(fù)相對(duì)介電常數(shù)的點(diǎn)即（）如下圖6.6所示。</p><p><b>　　圖6.6 </b></p><p>　　將細(xì)胞復(fù)介電常數(shù)中虛部進(jìn)行對(duì)比，虛部即為極化過(guò)程的介電損耗值，反映電介質(zhì)在交變電場(chǎng)下，將一部分電能轉(zhuǎn)換成熱能損耗。如圖6.7所示。</p><p><b>　　圖6.7

80、 </b></p><p>　　在交變電場(chǎng)作用下，電介質(zhì)在低頻段時(shí)產(chǎn)生的能量損耗3D-FDM模擬比理論上的損失多，在高頻段時(shí)其能量損失值幾乎相同[18]。當(dāng)細(xì)胞處于高頻段時(shí)，細(xì)胞完全極化，其介電性質(zhì)和培養(yǎng)液相同所以?xún)煞N模型所反映的能量損失值相同，而在低頻段時(shí)細(xì)胞由于自身的所帶電荷干擾，極化過(guò)程沒(méi)有理論上來(lái)的迅速和均勻，所以理論模型和3D-FDM模擬有差距。</p><p>　　

81、三維有限差分法與PS方程理論三維結(jié)果對(duì)比</p><p>　　綜合上述考慮細(xì)胞在極化過(guò)程與外加交流電頻率有關(guān)，同時(shí)在不同的細(xì)胞體積分?jǐn)?shù)下，細(xì)胞整體系統(tǒng)的介電弛豫結(jié)果也不相同。在之前利用給定參數(shù)計(jì)算細(xì)胞介電特性參數(shù)，結(jié)果比較準(zhǔn)確，在此我們利用所建立的3D-FDM模型模擬三維參數(shù)，在不同細(xì)胞體積分?jǐn)?shù)和外加交流電頻率計(jì)算相應(yīng)細(xì)胞節(jié)點(diǎn)特性參數(shù)，同時(shí)與理論模型對(duì)照。</p><p>　　將PS理論方

82、程在頻率作用下，不同細(xì)胞體積分?jǐn)?shù)情況下的細(xì)胞介電弛豫，由圖6.8三維曲面所示。黑色網(wǎng)格所示為PS理論方程下的曲面，彩色網(wǎng)格為3D-FDM模型的曲面。從圖上可以看出細(xì)胞處于低頻段時(shí)，不同的作用下，其復(fù)介電常數(shù)實(shí)部吻合效果不是很好，而在其他頻率段時(shí)吻合效果相對(duì)較好。這與之前分析在特定下，低頻段的吻合有點(diǎn)差距類(lèi)似，與細(xì)胞內(nèi)在電荷有關(guān)，在外加電場(chǎng)作用下不能迅速發(fā)生改變，而是通過(guò)內(nèi)在電荷調(diào)整。在同一頻率下，不同的的吻合程度在相對(duì)低時(shí)吻合較好，這與

83、Wagner理論稀疏系統(tǒng)理論吻合[19]。</p><p>　　圖6.8 三維對(duì)照?qǐng)D</p><p>　　同理按上述方法對(duì)系統(tǒng)復(fù)電導(dǎo)率進(jìn)行多變量研究，在頻率作用下，不同細(xì)胞體積分?jǐn)?shù)情況下的細(xì)胞介電弛豫，由圖6.9三維曲面所示。黑色網(wǎng)格所示為PS理論方程下的曲面，彩色網(wǎng)格為3D-FDM模型的曲面。從圖上可以看出細(xì)胞處于高頻段時(shí)，不同的作用下，復(fù)電導(dǎo)率的吻合效果不佳，而在低頻率段時(shí)吻合效果

84、較好，與圖6.3中分析的結(jié)果類(lèi)似。當(dāng)細(xì)胞極化后其內(nèi)部細(xì)胞器對(duì)培養(yǎng)液電導(dǎo)率的影響，導(dǎo)致所測(cè)電導(dǎo)率與理論上有差距。</p><p>　　圖6.9 三維對(duì)照?qǐng)D</p><p><b>　　結(jié)束語(yǔ)</b></p><p>　　本文從介電參數(shù)的定義出發(fā)，利用三維有限差分法模擬細(xì)胞介電弛豫現(xiàn)象，建立模型對(duì)細(xì)胞所在整個(gè)系統(tǒng)的介電特性進(jìn)行模擬計(jì)算，同時(shí)以單

85、球殼模型為例對(duì)細(xì)胞的介電特性和計(jì)算，以PS理論方程為基準(zhǔn)進(jìn)行對(duì)比。在對(duì)模型分析求解過(guò)程中對(duì)原先模型算法提出改進(jìn)，大大提高了計(jì)算速度和減少了內(nèi)存容量。對(duì)結(jié)果分析，三維有限差分法數(shù)字模擬與PS理論方程的結(jié)果很近，誤差在之內(nèi)。最后將以頻率作為單變量的模擬擴(kuò)展值雙變量（頻率，體積分?jǐn)?shù)）模擬對(duì)比，三維曲面吻合結(jié)果較好。</p><p>　　但是鑒于三維有限差分法在計(jì)算時(shí)的誤差比較大，在之后研究中以三維有限元法對(duì)其進(jìn)行完善和

86、補(bǔ)充，提高結(jié)果的吻合度。</p><p><b>　　參考文獻(xiàn)</b></p><p>　　Foster K R, Schwan H P, Dielectric Properties of tissues [M].Boca Ratun FL: CRC Press 1996.</p><p>　　曾攀.有限元分析及應(yīng)用[M].北京：清華大學(xué)出版社

87、，2004，1～5.</p><p>　　CODATA. Electric constant. 2006 CODATA recommended values. NIST [2007-08-08]. </p><p>　　Wagner K W 1914 Erklarung der Dielectrishen Nachwirkungsvorgange auf Grund Maxwel

88、lscher Vorstellungen Arch. Electrotech. 2 371–87.</p><p>　　王立洪,馬青.應(yīng)用MAPLE軟件分析生物細(xì)胞介電弛豫效應(yīng)[J].寧波大學(xué)學(xué)報(bào),2003,16(4):305~307.</p><p>　　Calame J P .Evolution of Davidson-Cole relaxation behavior in rando

89、m conductor-insulator composites J. Appl.Phys. 94 5945–57.</p><p>　　Asami K .Simulation of dielectric relaxation in periodic binary systems of complex geometry J. Colloid Interface Sci. 292 228–35.</p>

90、<p>　　Koji Asami. Dielectric properties of biological tissues in which cells are connected by communicating junctions [J].J. Phys. D: Appl. Phys. 2007, 40, 3718~3727.</p><p>　　Koji Asami. Dielectric d

91、ispersion in biological cells of complex geometry simulated by the three–dimensional finite difference method [J] .J. Phys. D: Appl. Phys. 39 (2006) 492~499.</p><p>　　施吉林，劉淑珍等.計(jì)算機(jī)數(shù)值方法(第二版）[M].北京：高等教育出版社,2008

92、(5)：211～215.</p><p>　　徐翠薇，孫繩武.計(jì)算方法引論[M].北京：高等教育出版社，2001（3）：188～189.</p><p>　　Pauly H,Schwan H P.Z Naturforsch,1959,14b:125～131.</p><p>　　Asami K，Hanai T，Koizumi N.J Membrane Biol，19

93、76，28：169.</p><p>　　Cole K Scale R H.D.dispersion and absorption in dielectrics [J].J Chem.Phys,1941, 9.</p><p>　　趙孔雙.介電譜方法及應(yīng)用[M].北京:化學(xué)工業(yè)出版社，2008（6）：19～20.</p><p>　　Hanai T, Asami K

94、, Koizumi N. Dielectric theory of concentrated suspension of shell-spheres in particular reference to the analysis of biological cell suspensions Bull [J]. Inst. Chem. Res., Kyoto University 1979, 57, 297~305.</p>

95、<p>　　馬青，湯治元，趙偉紅.紅細(xì)胞被動(dòng)電生理特征的數(shù)據(jù)分析[J].中國(guó)醫(yī)學(xué)物理雜志，2001,</p><p>　　18(4):231~233.</p><p>　　馬青，渡邊牧夫，洲崎敏伸.應(yīng)用COLE-COLE公式分析蛙血液細(xì)胞介電譜[J].中國(guó)生物醫(yī)學(xué)工程學(xué)報(bào)，2003,22(4),309~312.</p><p>　　Yoichi Kats

96、umoto, Yoshihito Hayashi, Ikuya Oshige,Shinji Omori,Noriyuki Kishii,Akio Yasuda,and Koji Asami.Dielectric Cytometry with Three-Dimensional Cellular Modeling[J]. Biophysical Journal Volume 95 September (2008) 3043~3047.&l

97、t;/p><p><b>　　附錄</b></p><p>　　matlab計(jì)算電勢(shì)分布</p><p>　　%%建立seidel迭代矩陣中B（B表示方程的系數(shù)矩陣）采用壓縮存儲(chǔ)</p><p>　　Y1=cell(n,n);%把Y第一層賦給Y1的第一層（節(jié)省空間），在之后的迭代中有用，而Y在下面給出Y是將會(huì)改變，</

98、p><p>　　%Y1=Y; % 下面要修正Y，暫時(shí)將Y備份到Y(jié)1中。其實(shí)，這里保留Y的第一層數(shù)據(jù)即可，可壓縮。</p><p><b>　　for j=1:n</b></p><p><b>　　for k=1:n</b></p><p>　　Y1{j,k}=Y{1,j,k};</p>

99、<p><b>　　end</b></p><p><b>　　end</b></p><p>　　for j=1:n % 修正Y第一層，因?yàn)樵贐中對(duì)應(yīng)第一層的phi為已知的，已經(jīng)放在g中故令為0</p><p><b>　　for k=1:n</b></p><p&g

100、t;　　Y{1,j,k}(1)=[0];Y{1,j,k}(2)=[0];Y{1,j,k}(3)=[0];Y{1,j,k}(4)=[0];Y{1,j,k}(5)=[0];Y{1,j,k}(6)=[0];</p><p><b>　　end</b></p><p><b>　　end</b></p><p>　　B=zeros

101、((n^3-2*n^2),12);h=1;S1=0;%B表示方程的系數(shù)矩陣，h控制行數(shù)，S1用來(lái)放sum（Y{i,j,k}）</p><p>　　for i=2:n-1</p><p><b>　　for j=1:n</b></p><p><b>　　for k=1:n</b></p><p>　

102、　S1=sum(Y{i,j,k});</p><p>　　%計(jì)算{i,j,k}，在一行矩陣中對(duì)應(yīng)的位置</p><p>　　i1=(i-2)*n^2+(j-1)*n+k; %phi(i-1,j,k)在n^3維向量中的坐標(biāo)</p><p>　　i2= i*n^2+(j-1)*n+k; %phi(i+1,j,k)在n^3維向量中的坐標(biāo)</p>&

103、lt;p>　　i3=(i-1)*n^2+(j-2)*n+k; %phi(i,j-1,k)在n^3維向量中的坐標(biāo)</p><p>　　i4=(i-1)*n^2+ j*n+k; %phi(i,j+1,k)在n^3維向量中的坐標(biāo)</p><p>　　i5=(i-1)*n^2+(j-1)*n+k-1;%phi(i,j,k-1)在n^3維向量中的坐標(biāo)</p><p

104、>　　i6=(i-1)*n^2+(j-1)*n+k+1;%phi(i,j,k+1)在n^3維向量中的坐標(biāo) </p><p>　　%n^3+1表示系統(tǒng)邊界外的，技巧性設(shè)置</p><p>　　if (j==1) i3=n^3+1; end </p><p>　　if (j==n) i4=n^3+1; end</p><p

105、>　　if (k==1) i5=n^3+1; end</p><p>　　if (k==n) i6=n^3+1; end</p><p>　　%第一到第六列放Y{i,j,k}（1..6）</p><p>　　B(h,1)=(Y{i,j,k}(1))/S1;</p><p>　　B(h,2)=(Y{i,j,k}(2))/S1;</p

106、><p>　　B(h,3)=(Y{i,j,k}(3))/S1;</p><p>　　B(h,4)=(Y{i,j,k}(4))/S1;</p><p>　　B(h,5)=(Y{i,j,k}(5))/S1;</p><p>　　B(h,6)=(Y{i,j,k}(6))/S1;</p><p>　　%7--12列放對(duì)應(yīng)的在B中位

107、置</p><p>　　B(h,7)=i1;</p><p>　　B(h,8)=i2;</p><p>　　B(h,9)=i3;</p><p>　　B(h,10)=i4;</p><p>　　B(h,11)=i5;</p><p>　　B(h,12)=i6; </p><p

108、>　　h=h+1;%行數(shù)增一</p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　end</b></p><p>　　4 %步驟%%進(jìn)行迭代計(jì)算 %計(jì)算w %超松弛迭代%初始化電勢(shì)圖</p&

109、gt;<p>　　A=zeros(1,n^3+1); %電勢(shì)， n^3+1技巧性設(shè)置為0</p><p>　　%A的初始值為線(xiàn)性的，最上面為1，最下面為0</p><p>　　for i=1:n^2 %最上層電勢(shì)初始化為1</p><p><b>　　A(1,i)=1;</b></p><p><

110、;b>　　end</b></p><p>　　for i=n^3-n^2+1:n^3 %最下層電勢(shì)初始化為0</p><p><b>　　A(1,i)=0;</b></p><p><b>　　end</b></p><p>　　%A(1,n^3+1)=0</p>

111、<p>　　for i=2:n-1 % 線(xiàn)性電勢(shì)梯度 </p><p>　　for j=(i-1)*(n^2)+1:i*n^2</p><p>　　A(1,j)=1-i/n;</p><p><b>　　end</b></p><p><b>　　end</b></p>

112、<p>　　A1=zeros(1,n^3+1);A1=A;%建立A1用于迭代存放電勢(shì)，功能和A相同</p><p><b>　　%迭代過(guò)程</b></p><p>　　for j=1:8%迭代次數(shù)</p><p>　　j %打印當(dāng)前迭代次數(shù)</p><p>　　for i=1:n^3-2*n^2

113、 </p><p>　　j1=B(i,7);j2=B(i,8);j3=B(i,9);j4=B(i,10);j5=B(i,11);j6=B(i,12);</p><p>　　A(1,n^2+i)=(1-w)*A(1,n^2+i)+w*(B(i,1)*A(1,j1)+B(i,2)*A(1,j2)+B(i,3)*A(1,j3)+B(i,4)*A(1,j4)+B(i,5)*A(1,j5)

114、+B(i,6)*A(1,j6));</p><p>　　if (i<=n^2) A(1,n^2+i)=A(1,n^2+i)+g(i,1)*w;</p><p><b>　　end</b></p><p><b>　　end</b></p><p>　　if max(abs(A1-A))<

115、;0.0001, break;%判斷條件</p><p>　　else A1=A;</p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　clear A1</b></p><p>　　A1=z

116、eros(n,n,n);%將迭代后的結(jié)果轉(zhuǎn)化成[n,n,n]形式</p><p><b>　　for i=1:n</b></p><p><b>　　for j=1:n</b></p><p><b>　　for k=1:n</b></p><p>　　A1(i,j,k)=A(

117、1,n^2*(i-1)+n*(j-1)+k);</p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　%Y=Y1;</b></p>&l

118、t;p><b>　　for i=1:n</b></p><p><b>　　for j=1:n</b></p><p>　　Y{1,i,j}=Y1{i,j};</p><p><b>　　end</b></p><p><b>　　end</b>&l

119、t;/p><p>　　I1=0 %計(jì)算電流，給定初始0</p><p>　　for i=1:n-1</p><p><b>　　for j=1:n</b></p><p><b>　　for k=1:n</b></p><p>　　I1=I1+Y{i,j,k}(2)*(A1(i,

120、j,k)-A1(i+1,j,k));% 電流方向? 我們這里是以i正方向?yàn)闇?zhǔn)，與文中相反</p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　end</b></p><p><b>　　l1=l1+1;&

121、lt;/b></p><p><b>　　end</b></p><p><b>　　end</b></p><p>　　2.maple 畫(huà)結(jié)果圖</p><p>　　> restart;</p><p>　　with(plots):</p><

122、;p>　　> #A:=((ε[b]-ε[p])/(ε[a]-ε[p]))^3*(ε[a]/ε[b])-(1-P)^3=0;</p><p>　　> ε[p]:=ε[m]*(2*(1-v)*ε[m]+(1+2*v)*ε[c])/((2+v)*ε[m]+(1-v)*ε[c]):</p><p>　　> #ε[p]=p1+q1*I;</p><p>

123、;　　ε[m]:=p2+q2*I:</p><p>　　ε[c]:=p3+q3*I:</p><p>　　ε[a]:=p4+q4*I:</p><p>　　#ε[b]:=p+q*I;</p><p>　　q2:=-k[m]/(2*Pi*f*ε[0]):</p><p>　　q3:=-k[c]/(2*Pi*f*ε[0])

124、:</p><p>　　q4:=-k[a]/(2*Pi*f*ε[0]):</p><p>　　> v:=(1-d/R)^3:</p><p>　　> p2:=5: p3:=80:p4:=80: d:=5*10^(-9): R:=5*10^(-6): P:=0.29: p2 := 5:k[m]:=10^(-7): k[c]:=1: k[a]:=1:ε[0

125、]:=8.854*10^(-12):evalc(Re(ε[p])):</p><p>　　> q1:=evalc(Im(ε[p])):</p><p>　　> ε[b]:=ε[a]*(2*(1-P)*ε[a]+(1+2*P)*ε[p])/((2+P)*ε[a]+(1-P)*ε[p]):</p><p><b>　　> ε[b]:</