基于分形理論的雷電先導(dǎo)建模與仿真畢設(shè)論文終極版

1、<p>　　基于分形理論的雷電先導(dǎo)建模與仿真</p><p><b>　　摘要</b></p><p>　　雷電是一種自然放電現(xiàn)象，它是由雷云(帶電的云層)對(duì)地面建筑物及大地的自然放電引起的。雷電放電的隨機(jī)性以及雷電先導(dǎo)在其產(chǎn)生和發(fā)展過程中易受周圍環(huán)境影響的特點(diǎn)導(dǎo)致雷電活動(dòng)對(duì)電力系統(tǒng)產(chǎn)生巨大危害。由于雷電放電的隨機(jī)性，需要用統(tǒng)計(jì)學(xué)的方法來進(jìn)行防雷性能的評(píng)估。

2、為此采用分形模型模擬了先導(dǎo)的發(fā)展過程，在二維和三維中描述雷電放電通道的分支和扭曲。通過比較自然界中雷電放電和仿真結(jié)果的分形維數(shù)，可以找到最接近實(shí)際的分形系數(shù)。與傳統(tǒng)的電氣幾何模型、改進(jìn)電氣幾何模型以及先導(dǎo)發(fā)展模型對(duì)比，分形模型不僅能夠給出逼真的先導(dǎo)通道圖形，而且能夠從統(tǒng)計(jì)學(xué)意義上對(duì)一些不合理的側(cè)面繞擊現(xiàn)象作出解釋。</p><p>　　本課題基于對(duì)雷電發(fā)展過程準(zhǔn)靜態(tài)場(chǎng)的假設(shè)，采用有限差分法作為電場(chǎng)計(jì)算的方法，將偏

3、微分方程的問題化為線性方程組問題，并用迭代法進(jìn)行線性方程的數(shù)值解。結(jié)合以上計(jì)算出的空間電場(chǎng)，采用分形方法完成了二維和三維雷電先導(dǎo)放電通道的建模與仿真，并實(shí)現(xiàn)了雷電先導(dǎo)發(fā)展過程及相應(yīng)的空間電勢(shì)圖形靜態(tài)和動(dòng)態(tài)顯示。在此基礎(chǔ)上，本課題還討論了避雷針(線)對(duì)閃電軌跡的影響。</p><p>　　關(guān)鍵詞：分形模型；三維建模；雷電先導(dǎo)；有限差分法；MATLAB</p><p>　　Modeling a

4、nd Simulation of Lightning Leader based on the Fractal Theory</p><p><b>　　Abstract</b></p><p>　　Lightning is a natural discharge phenomenon, which is caused by thundercloud’s (electr

5、ically charged clouds) natural discharge to ground-level structures and earth. Lightning causes great damage to power system due to the stochastic nature of lightning discharges and the vulnerability of lightning leaders

6、 to various environments in the process of inception and propagation. As a result of the randomness of lightning discharge, it’s necessary to introduce statistical methods into the evaluation of</p><p>　　Thi

7、s topic is based on the assumption of quasi-static field in lightning development process, using the finite difference method to calculate the electric field, changing partial differential equations into linear equations

8、, using the iterative method to get the numerical solution of linear equations. Combining with the above electric field, this topic uses the fractal theory to model and simulate the two-dimensional and three-dimensional

9、lightning leader discharge channel, achieving development</p><p>　　Keywords: Fractal model; Three-dimensional modeling; Lightning leader; Finite difference method; MATLAB</p><p><b>　　目 錄&l

10、t;/b></p><p><b>　　摘要I</b></p><p>　　AbstractII</p><p><b>　　第1章 緒論1</b></p><p>　　1.1 課題背景1</p><p>　　1.2 國(guó)內(nèi)外研究動(dòng)態(tài)2</p>

11、<p>　　1.3 本文主要內(nèi)容和難點(diǎn)3</p><p>　　1.3.1 工作內(nèi)容4</p><p>　　1.3.2 課題難點(diǎn)4</p><p>　　第2章 有限差分法和格點(diǎn)式隨機(jī)分形模型5</p><p>　　2.1 有限差分法概述5</p><p>　　2.2 DBM概述7</p>

12、;<p>　　2.1.1 先導(dǎo)發(fā)展規(guī)則7</p><p>　　2.1.2 分形先導(dǎo)發(fā)展模型關(guān)鍵參數(shù)的選取8</p><p>　　第3章 雷電先導(dǎo)放電通道的二維模型10</p><p>　　3.1 雷電先導(dǎo)放電的物理過程及參數(shù)10</p><p>　　3.1.1雷電放電概述10</p><p>

13、　　3.3 程序流程13</p><p>　　3.3.1 模塊114</p><p>　　3.3.2 模塊214</p><p>　　3.3.3 模塊316</p><p>　　3.4 動(dòng)畫效果和電勢(shì)的變化19</p><p>　　3.4.1 動(dòng)畫效果19</p><p>　　3.4

14、.2 電勢(shì)的變化21</p><p>　　第4章 雷電先導(dǎo)放電通道的三維模型22</p><p>　　4.1 模型的建立與簡(jiǎn)化22</p><p>　　4.2 程序流程22</p><p>　　4.2.1 模塊1 初始情況22</p><p>　　4.2.2 模塊2 計(jì)算電勢(shì)23</p>&

15、lt;p>　　4.2.3 模塊3 計(jì)算概率，決定下一個(gè)擊穿點(diǎn)24</p><p>　　4.3 動(dòng)畫效果和電勢(shì)的變化24</p><p>　　4.3.1 動(dòng)畫效果24</p><p>　　4.3.2 電勢(shì)的變化25</p><p>　　第5章 避雷針(線)對(duì)先導(dǎo)通道發(fā)展的影響28</p><p>　　5

16、.1 避雷針對(duì)先導(dǎo)通道發(fā)展的影響28</p><p>　　5.1.1 二維模擬避雷針對(duì)雷電先導(dǎo)的影響28</p><p>　　5.1.2 三維模擬避雷針對(duì)雷電先導(dǎo)的影響29</p><p>　　5.2 避雷線對(duì)架空線路的保護(hù)作用32</p><p>　　第6章 結(jié)論與展望36</p><p><b&g

17、t;　　6.1 結(jié)論36</b></p><p><b>　　6.2 展望36</b></p><p><b>　　參考文獻(xiàn)37</b></p><p><b>　　致謝39</b></p><p><b>　　第1章 緒論</b>&

18、lt;/p><p><b>　　1.1 課題背景</b></p><p>　　雷電是伴有閃電和雷鳴的一種雄偉壯觀而又有點(diǎn)令人生畏的放電現(xiàn)象[]，是發(fā)生在大氣層中的聲、光、電物理現(xiàn)象，常見的雷電現(xiàn)象是一部分帶電云層與另一部分帶異種電荷的云層與大地之間的迅猛的放電過程[]。</p><p>　　空中的塵埃、冰晶等物質(zhì)在云層中翻滾運(yùn)動(dòng)的時(shí)候，經(jīng)過一些復(fù)雜

19、過程，使這些物質(zhì)分別帶上了正電荷與負(fù)電荷。經(jīng)過運(yùn)動(dòng)，帶上相同電荷的質(zhì)量較重的物質(zhì)會(huì)到達(dá)云層的下部(一般為負(fù)電荷)，帶上相同電荷的質(zhì)量較輕的物質(zhì)會(huì)到達(dá)云層的上部(一般為正電荷)。這樣，同性電荷的匯集就形成了一些帶電中心，當(dāng)異性帶電中心之間的空氣被其強(qiáng)大的電場(chǎng)擊穿時(shí)，就形成“云間放電”(即閃電)。</p><p>　　帶負(fù)電荷的云層向下靠近地面時(shí)，地面的凸出物、金屬等會(huì)被感應(yīng)出正電荷，隨著電場(chǎng)的逐步增強(qiáng)，雷云向下形成

20、下行先導(dǎo)，地面的物體形成向上閃流，二者相遇即形成對(duì)地放電。這就容易造成雷電災(zāi)害。雷電災(zāi)害是“聯(lián)合國(guó)國(guó)際減災(zāi)十年”公布的最嚴(yán)重的十種自然災(zāi)害之一。最新統(tǒng)計(jì)資料表明，雷電造成的損失已經(jīng)上升到自然災(zāi)害的第三位。全球每年因雷擊造成人員傷亡、財(cái)產(chǎn)損失不計(jì)其數(shù)。據(jù)不完全統(tǒng)計(jì)，我國(guó)每年因雷擊以及雷擊負(fù)效應(yīng)造成的人員傷亡達(dá)3000～4000人，財(cái)產(chǎn)損失在50億元到100億元人民幣。</p><p>　　雷電災(zāi)害所涉及的范圍幾乎遍

21、布各行各業(yè)。現(xiàn)代電子技術(shù)的高速發(fā)展，帶來的負(fù)效應(yīng)之一就是其抗雷擊浪涌能力的降低。以大規(guī)模集成電路為核心組件的測(cè)量、監(jiān)控、保護(hù)、通信、計(jì)算機(jī)網(wǎng)絡(luò)等先進(jìn)電子設(shè)備廣泛運(yùn)用于電力、航空、國(guó)防、通信、廣電、金融、交通、石化、醫(yī)療以及其它現(xiàn)代生活的各個(gè)領(lǐng)域，以大型CMOS集成元件組成的這些電子設(shè)備普遍存在著對(duì)暫態(tài)過電壓、過電流耐受能力較弱的缺點(diǎn)，暫態(tài)過電壓不僅會(huì)造成電子設(shè)備產(chǎn)生誤操作，也會(huì)造成更大的直接經(jīng)濟(jì)損失和廣泛的社會(huì)影響。</p>

22、<p>　　雷電的發(fā)生是不能避免，不能阻止的，因?yàn)檫@是自然現(xiàn)象。但是，我們可以設(shè)法避雷，設(shè)法使雷電天氣對(duì)我們不造成災(zāi)害。在雷雨天氣，我們有時(shí)會(huì)—看到一些些參天大樹被雷電擊倒，而周圍的一些高塔、高樓等高層建筑卻安然無(wú)恙。這是因?yàn)檫@些參天大樹由于受到帶大量電荷的云層的感應(yīng)，也帶上大量電荷，積累的電荷過多時(shí)，樹就被擊倒。而同樣倩況下．高層建筑安然無(wú)恙恐怕要?dú)w功于避雷針了。在許多高塔上都有于個(gè)金屬做的、狀如針頭向上直立的東西，這就

23、是避雷針。原來，避雷針應(yīng)叫“引雷針”。在雷雨天氣，高樓上空出現(xiàn)帶電云層時(shí)，迅雷針和高樓頂部都被感應(yīng)上大量電荷，由于避雷針針頭是尖的，而靜電感應(yīng)時(shí)，導(dǎo)體尖端總是聚集了最多的電荷。這樣，避雷針就聚集了大部分電荷。避雷針又與這些帶電云層形成了一個(gè)電容器，由于它較尖，即這個(gè)電容器的兩極板正對(duì)面積很小，電容也就很小，也就是說它所能容納的電荷很少。而它又聚集了大部分電荷，所以，當(dāng)云層上電荷較多時(shí)，避雷針與云層之間的空氣就很容易被擊穿，成為導(dǎo)體。這樣

24、，帶電云層與避雷針形成通路，而避雷針又是接地的。避雷針就可以把云層上的電荷導(dǎo)人大地，使其不對(duì)高層建筑構(gòu)成危險(xiǎn)，保證了它的安全。</p><p>　　國(guó)內(nèi)外超特高壓輸電線路的運(yùn)行統(tǒng)計(jì)結(jié)果表明，雷擊跳閘已經(jīng)成為特高壓線路故障的主要原因，而繞擊則是造成雷擊跳閘的主要因素。另外，超高壓線路雷擊跳閘故障中很大部分也是由繞擊造成的[]。隨著特高壓電網(wǎng)的建設(shè)，如何防護(hù)輸電線路的繞擊逐漸成為線路防雷分析中需要考慮的主要問題，而如

25、何正確的描述雷電的發(fā)展過程又決定了防雷分析的準(zhǔn)確性和合理性。自然界的雷電放電是一種超長(zhǎng)間隙的放電現(xiàn)象。由于其范圍巨大，其發(fā)展過程中會(huì)受到各種不確定因素的影響，如自然界的背景噪聲、空間電場(chǎng)、氣象條件等，從而呈現(xiàn)出彎曲分叉的隨機(jī)放電路徑，給研究者和工程師研究雷擊過程帶來了很大的困難。目前我國(guó)輸電線路的雷電防護(hù)設(shè)計(jì)及其性能估算以我國(guó)電力行業(yè)標(biāo)準(zhǔn)DL/T620-1997《交流電氣裝置的過電壓保護(hù)和絕緣配合》為依據(jù)。規(guī)程法的繞擊計(jì)算未考慮雷電發(fā)展

26、過程、雷電流大小和地面傾角對(duì)屏蔽效果的影響，而是根據(jù)經(jīng)驗(yàn)和小電流試驗(yàn)?zāi)Ｐ驮囼?yàn)結(jié)果提出的綜合平均法，常不能反映線路的具體特性。因此，非常需要找到一種最合理的雷電分析模型。</p><p>　　1.2 國(guó)內(nèi)外研究動(dòng)態(tài)</p><p>　　目前，輸電線路雷擊屏蔽模型主要有：電氣幾何模型(EGM)[、]、改進(jìn)電氣幾何模型[、]、先導(dǎo)發(fā)展模型(LPM)[、、]。</p><p&g

27、t;　　對(duì)于EGM模型及改進(jìn)EGM模型，將雷電流幅值與雷擊距離聯(lián)系起來，將雷擊輸電線路用幾何作圖的方法來描述，繞擊率與雷電流幅值相關(guān)，能夠考慮線路結(jié)構(gòu)和雷電參數(shù)等對(duì)繞擊率的影響，其通過計(jì)算雷電擊距，以繞擊弧的方式對(duì)輸電線路繞擊性能進(jìn)行評(píng)估，但該方法沒有考慮雷電放電的分散性及線路運(yùn)行電壓對(duì)線路繞擊耐雷性能的影響。雖然在不同線路參數(shù)及環(huán)境下的應(yīng)用中提出了一些修正的方法，但其防雷設(shè)計(jì)與超特高壓輸電線路實(shí)際運(yùn)行結(jié)果仍存在較大的誤差[]。但該方法

28、的理論是由保護(hù)角較大、桿塔高度較低的線路運(yùn)行經(jīng)驗(yàn)總結(jié)而來的，在超高壓線路設(shè)計(jì)時(shí)，理論計(jì)算的繞擊率總是遠(yuǎn)小于實(shí)際運(yùn)行值。如20世紀(jì)90年代中后期，華東地區(qū)連年屢次發(fā)生雷害事故，而遭受雷害事故的桿塔多為ZM1型桿塔，其邊相導(dǎo)線保護(hù)角為7.2°，滿足EGM的完全屏蔽條件。根據(jù)現(xiàn)場(chǎng)調(diào)查情況，經(jīng)分析認(rèn)為：由于放電的分散性，超出EGM繞擊區(qū)定位的雷閃仍可能發(fā)生繞擊[3]，這一現(xiàn)象顯然無(wú)法用EGM來解釋。</p><p&

29、gt;　　LPM模型基于實(shí)驗(yàn)室試驗(yàn)的基礎(chǔ)上，闡述了雷電發(fā)展的整個(gè)過程，是較EGM模型更為先進(jìn)的雷電屏蔽模型，LPM模型提出輸電線路繞擊側(cè)面距離及屏蔽失效寬度這兩個(gè)基本參數(shù)使得模型在對(duì)輸電線路繞擊耐雷性能的評(píng)估有較大的適用性，該模型考慮了上下行先導(dǎo)相對(duì)發(fā)展的實(shí)際過程，是對(duì)雷電發(fā)展過程的近似模擬，但對(duì)雷電先導(dǎo)的發(fā)展僅考慮了雷電放電的確定性因素，認(rèn)為雷電先導(dǎo)一直沿先導(dǎo)發(fā)展空間最大場(chǎng)強(qiáng)方向發(fā)展，忽略了隨機(jī)性因素，其模型中用于表征輸電線路繞擊耐雷

30、水平的屏蔽失效寬度這一參數(shù)也無(wú)法體現(xiàn)出雷云放電的分散性。因此傳統(tǒng)的LPM模型及其線路繞擊耐雷性能評(píng)估方法具有一定的局限性。國(guó)內(nèi)學(xué)者基于雷電放電和長(zhǎng)空氣間隙放電物理過程的相似性，建立了實(shí)驗(yàn)室的等比例縮小輸電線路雷擊模型[、]，通過輸電線路繞擊概率曲線的形式對(duì)輸電線路繞擊耐雷性能進(jìn)行研究。但實(shí)驗(yàn)室試驗(yàn)中長(zhǎng)間隙放電的先導(dǎo)發(fā)展并不充分，如何將實(shí)驗(yàn)室數(shù)據(jù)應(yīng)用到實(shí)際的防雷設(shè)計(jì)還有待深入地研究。</p><p>　　因此，尋求

31、一種能夠充分考慮雷電放電確定性及隨機(jī)性因素的輸電線路繞擊屏蔽模型及其繞擊耐雷性能評(píng)估方法對(duì)特高壓的防雷設(shè)計(jì)及運(yùn)行至關(guān)重要?？紤]到雷電先導(dǎo)的發(fā)展呈現(xiàn)樹狀和狹窄的隨機(jī)放電通道結(jié)構(gòu)[]，該結(jié)構(gòu)有明顯的自相似特征，而分形正是描述該結(jié)構(gòu)的最好數(shù)學(xué)方法。分形理論是當(dāng)今世界十分風(fēng)靡和活躍的新理論、新學(xué)科。分形概念是美籍?dāng)?shù)學(xué)家曼德布羅特(B.B.Mandelbort)首先提出的。分形理論既是非線性科學(xué)的前沿和重要分支，又是一門新興的橫斷學(xué)科。作為一種方

32、法論和認(rèn)識(shí)論，其啟示是多方面的：一是分形整體與局部形態(tài)的相似，啟發(fā)人們通過認(rèn)識(shí)部分來認(rèn)識(shí)整體，從有限中認(rèn)識(shí)無(wú)限；二是分形揭示了介于整體與部分、有序與無(wú)序、復(fù)雜與簡(jiǎn)單之間的新形態(tài)、新秩序；三是分形從一特定層面揭示了世界普遍聯(lián)系和統(tǒng)一的圖景。</p><p>　　分形先導(dǎo)發(fā)展模型是一種基于分形理論的雷電先導(dǎo)模型，它在傳統(tǒng)雷電先導(dǎo)模型的基礎(chǔ)上，充分考慮雷電上下行先導(dǎo)路徑的分形特性以及先導(dǎo)發(fā)展過程的隨機(jī)性和確定性因素，將

33、LPM模型的理論基礎(chǔ)和分形理論相結(jié)合而建立的一種基于分形理論的超特高壓輸電線路耐雷性能評(píng)估方法。Niemeyer,Pientronero和Wielsmann于1983年研究了平行玻璃之間SF6中的電擊現(xiàn)象，他們得到了多枝杈的徑向放電圖形(1µs的30kV電脈沖)，他們對(duì)圖形進(jìn)行了計(jì)算機(jī)模擬，提出了介電擊穿模型(dielectric breakdown model，DBM)。繼Niemeyer等人提出DBM模型后[]，不少研究者

34、又將其用在對(duì)雷電閃道的模擬上[、、、、、]。分形模擬不僅其結(jié)果能夠體現(xiàn)通道的分岔和彎曲，而且其過程也與雷電放電過程基本類似。閃電現(xiàn)象屬于氣體電介質(zhì)擊穿現(xiàn)象，用類似DBM模型的方法進(jìn)行計(jì)算機(jī)模擬是比較自然的想法，因此一些學(xué)者以此為依據(jù)來模擬自然界中的閃電現(xiàn)象。金秀儒的《用計(jì)算機(jī)模擬閃電形成的嘗試》即為一次不錯(cuò)的嘗試[]，但作者最終沒有得出一個(gè)完善的模擬。中國(guó)科學(xué)技術(shù)大學(xué)物理2班的王寧遠(yuǎn)在他的基礎(chǔ)上更近一步[]，用</p>&

35、lt;p>　　1.3 本文主要內(nèi)容和難點(diǎn)</p><p>　　本課題的主要任務(wù)是基于MATLAB軟件編程，采用分形方法實(shí)現(xiàn)二維和三維雷電先導(dǎo)放電通道的建模與仿真，充分了解介質(zhì)電擊穿的具體原理以及性質(zhì)，在此基礎(chǔ)上利用其核心公式進(jìn)行條件的確立與簡(jiǎn)化，找到合適的邊界條件并進(jìn)行建模。模型首先通過在二維空間和三維空間劃分網(wǎng)格把實(shí)際問題離散化，然后通過適當(dāng)?shù)乃惴▽?duì)閃電形成過程進(jìn)行模擬。MATLAB編程時(shí)采用有限差分法不

36、斷更新電位，并結(jié)合隨機(jī)函數(shù)使程序不停地尋找生長(zhǎng)點(diǎn)，直到閃電擊穿到地面。</p><p>　　1.3.1 工作內(nèi)容</p><p>　　(1)充分了解介質(zhì)電擊穿的具體原理以及性質(zhì)[]；</p><p>　　(2)理解有限差分法的具體原理[23]；</p><p>　　(3)學(xué)會(huì)使用MATLAB軟件，尤其是MATLAB程序設(shè)計(jì)與圖形處理[]；&l

37、t;/p><p>　　(4)編寫二維程序，實(shí)現(xiàn)雷電先導(dǎo)打到地面或邊界后程序結(jié)束運(yùn)行；</p><p>　　(5)編寫三維程序，實(shí)現(xiàn)雷電先導(dǎo)打到地面或邊界后程序結(jié)束運(yùn)行；</p><p>　　(6)程序中添加雷電擊中地面的時(shí)間；</p><p>　　(7)實(shí)現(xiàn)雷電的動(dòng)態(tài)效果；</p><p>　　(8)通過圖形動(dòng)態(tài)顯示電勢(shì)的

38、變化；</p><p>　　(9)研究避雷針（線）對(duì)雷電先導(dǎo)的影響；</p><p>　　(10)研究一根輸電線上場(chǎng)強(qiáng)的變化。</p><p>　　1.3.2 課題難點(diǎn)</p><p>　　(1)對(duì)介質(zhì)電擊穿的核心公式進(jìn)行條件的確立與簡(jiǎn)化，找到合適的邊界條件；</p><p>　　(2)如何使程序不停地尋找生長(zhǎng)點(diǎn)；&l

39、t;/p><p>　　(3)編程實(shí)現(xiàn)動(dòng)態(tài)效果；</p><p>　　(4)編程實(shí)現(xiàn)三維電勢(shì)的動(dòng)態(tài)變化。</p><p>　　第2章 有限差分法和格點(diǎn)式隨機(jī)分形模型</p><p>　　2.1 有限差分法概述</p><p>　　經(jīng)典場(chǎng)的邊值問題在數(shù)學(xué)上表達(dá)為泊松方程和拉普拉斯方程，但解偏微分方程往往是比較困難的[23]。

40、但是很多時(shí)候?qū)τ诰唧w問題我們需要的不是解析解，而是數(shù)值解，所以可以考慮用連續(xù)變量離散化的方法求出數(shù)值解，在足夠的精度上進(jìn)行逼近，這就引出了有限差分法。</p><p><b>　　有限差分法：</b></p><p>　　微分： (2-1)</p><p>　　用有限的h

41、代替，使得</p><p><b>　　(2-2)</b></p><p><b>　　差分的種類：</b></p><p><b>　　一階差分：</b></p><p><b>　　或者或者</b></p><p><b&

42、gt;　　二階差分：</b></p><p><b>　　(2-3)</b></p><p>　　設(shè)U(x,y,z)為空間電勢(shì)的函數(shù)。</p><p><b>　　泊松方程：</b></p><p><b>　　(2-4)</b></p><p&

43、gt;　　使用二階差分代替泊松方程中的拉普拉斯算符，有：</p><p><b>　　(2-5)</b></p><p>　　Σ表示分別對(duì)三個(gè)變?cè)蟛罘种停韵孪嗤?lt;/p><p>　　矩陣(數(shù)組)是計(jì)算機(jī)中重要的數(shù)據(jù)結(jié)構(gòu)，為了方便用矩陣去存儲(chǔ)數(shù)據(jù)，我們網(wǎng)格去劃分空間，從而不僅使方程化為簡(jiǎn)單的有限差分形式，而且這樣的數(shù)據(jù)類型在計(jì)算機(jī)中易于儲(chǔ)

44、存和運(yùn)算。那么h=k=l=1，并且令f(x,y,z)=u(x,y,z)就有：</p><p><b>　　(2-6)</b></p><p>　　這就是泊松方程的有限差分形式，以下估計(jì)該方程的精度：</p><p>　　由泰勒公式，易知有以下結(jié)果：</p><p><b>　　(2-7)</b>&l

45、t;/p><p><b>　　若考慮離散的點(diǎn)：</b></p><p><b>　　(2-8)</b></p><p><b>　　上述六式相加</b></p><p><b>　　(2-9)</b></p><p><b>

46、　　代入</b></p><p><b>　　得：</b></p><p><b>　　(2-10)</b></p><p>　　(2-1)式的精度為三階，四階及更高階項(xiàng)被略去。</p><p><b>　　若滿足</b></p><p>&

47、lt;b>　　有</b></p><p><b>　　(2-11)</b></p><p>　　此即拉普拉斯方程的有限差分形式。</p><p>　　這里，我們通過有限差分的方法，把偏微分方程在三階精度下簡(jiǎn)化為形式易于計(jì)算的代數(shù)方程，從而使之易于在計(jì)算機(jī)上實(shí)現(xiàn)。</p><p>　　注：有時(shí)我們需要解二

48、維靜電場(chǎng)，則方程退化為：</p><p><b>　　(2-12)</b></p><p><b>　　2.2 DBM概述</b></p><p>　　2.1.1 先導(dǎo)發(fā)展規(guī)則</p><p>　　為簡(jiǎn)潔起見，以二維平板放電為例，介紹DBM的主要算法和具體步驟[24]。圖2-1中給出了運(yùn)行若干步的模

49、擬放電圖樣，其中黑點(diǎn)表示放電通道中的點(diǎn)，黑點(diǎn)以實(shí)線相連形成通道。</p><p>　　圖2-1 平板放電過程的DBM模擬方法示意圖</p><p>　　模擬區(qū)域內(nèi)的電場(chǎng)在準(zhǔn)靜態(tài)近似下，電位滿足Poisson方程。將上下邊界和放電通道都作為第一類邊界條件，即其中格點(diǎn)的電位是固定的。于是無(wú)需考慮空間電荷，因?yàn)槠湫?yīng)已包含在邊界條件中。圖2-1中每一格點(diǎn)電位可由離散的Laplace方程求取，即

50、：</p><p><b>　　(2-13)</b></p><p>　　每一步的放電過程是鄰近已生成的放電通道的未放電格點(diǎn)(i′,j′)，如果與放電通道中某點(diǎn)(i,j)之間的平均場(chǎng)強(qiáng)超過放電臨界場(chǎng)強(qiáng)Ec，即：</p><p><b>　　(2-14)</b></p><p>　　該點(diǎn)就是可能放電點(diǎn)

51、，(2)式中di′j′,ij為兩點(diǎn)之間的距離。可能放電點(diǎn)在圖2-1中用白圈表示。一個(gè)可能放電點(diǎn)成為新的通道點(diǎn)的概率p與局部的場(chǎng)強(qiáng)有關(guān)，其表達(dá)式如(3)式所示：</p><p><b>　　(2-15)</b></p><p>　　其中(i,j)和(i′,j′)分別為一個(gè)已放電點(diǎn)和一個(gè)可能放電點(diǎn)，η稱為發(fā)展概率指數(shù)，放電通道由(i,j)向(i′,j′)發(fā)展的概率與兩點(diǎn)之

52、間平均場(chǎng)強(qiáng)的η次方成正比。(3)式右端分母是對(duì)所有(i,j)和(i′,j′)的組合求和，需要說明的是，每一個(gè)通道點(diǎn)和每一個(gè)與之鄰近的可能發(fā)展點(diǎn)都構(gòu)成一個(gè)這樣的組合。按計(jì)算出的概率隨機(jī)抽取下一步的放電發(fā)展方向，如果新的放電通道是由(i,j)向(i′,j′)發(fā)展，假設(shè)沿著放電通道電場(chǎng)強(qiáng)度為Ech，則(i′,j′)處更新后的電位為：</p><p><b>　　(2-16)</b></p&g

53、t;<p>　　然后開始下一步的模擬，如此循環(huán)，直到通道聯(lián)通上下邊界，表示此時(shí)間隙擊穿。</p><p>　　2.1.2 分形先導(dǎo)發(fā)展模型關(guān)鍵參數(shù)的選取</p><p>　　由(2-15)式可以看出，η是影響上式的關(guān)鍵參數(shù)[]。通過計(jì)算不同取值下模擬出的雷電先導(dǎo)通道的分形維數(shù)，對(duì)比實(shí)際觀測(cè)到的雷電先導(dǎo)通道的分形維數(shù)，可以確定η的合理取值。按照(2-15)式計(jì)算出的發(fā)展概率隨機(jī)

54、抽取下一步的發(fā)展方向，如此循環(huán)，直到間隙擊穿。</p><p>　　分形維數(shù)反映復(fù)雜圖形占有空間的有效性，是分形圖形有序性與無(wú)序性，確定性與隨機(jī)性的量度。由式(2-15)式可見，分形參數(shù)η其值大小影響了分形先導(dǎo)通道周圍各處發(fā)展概率，決定了雷電先導(dǎo)發(fā)展的隨機(jī)性及確定性。選取合適的分形參數(shù)η值對(duì)整個(gè)雷電路徑模擬的精確性有非常重要的影響。因此，本文對(duì)各種分形參數(shù)η值下雷電放電路徑進(jìn)行模擬，以期得出分形參數(shù)η值對(duì)雷電模擬

55、路徑分形維數(shù)的影響，確定分形先導(dǎo)發(fā)展模型合適的分形參數(shù)η值。</p><p>　　基于盒維法的雷電發(fā)展分形維數(shù)為：</p><p><b>　　(2-17)</b></p><p>　　式中，Db為分形維數(shù)，ε為測(cè)量尺度，N(ε)為分形圖形能夠被邊長(zhǎng)為ε的正方形完全覆蓋的最小個(gè)數(shù)。</p><p>　　基于盒維法，文[2

56、6]對(duì)不同分形參數(shù)η值下雷電模擬路徑的分形維數(shù)進(jìn)行了計(jì)算，如表2-1所示。并選取了3個(gè)不同η值下雷電路徑模擬圖(見圖2-2)。</p><p>　　表2-1 不同η值下的分形維數(shù)Db</p><p>　　當(dāng)η取值較小時(shí)，通過式(2-15)可以看出，先導(dǎo)頭部各方向發(fā)展概率趨于相同。因此，導(dǎo)致先導(dǎo)通道分支過多，如圖2-2(a)所示。此時(shí)，隨機(jī)性在雷電先導(dǎo)的發(fā)展中占主要作用，相應(yīng)的分形維數(shù)較大

57、。反之當(dāng)η值取值較大時(shí)，先導(dǎo)頭部最大電場(chǎng)方向獲得了很大的發(fā)展概率，而其它方向的先導(dǎo)發(fā)展概率較小，先導(dǎo)主要沿最大場(chǎng)強(qiáng)方向發(fā)展，分支發(fā)展的概率也相對(duì)較小，如圖2-2(c)所示。此時(shí)，確定性在雷電先導(dǎo)的發(fā)展中占主要地位，相應(yīng)的分形維數(shù)較小。對(duì)于實(shí)際放電分形維數(shù)，Charles I Richman 運(yùn)用盒維法對(duì)雷云放電圖片進(jìn)行分形維數(shù)的測(cè)量得出，實(shí)際雷云放電路徑的分形維數(shù)在1.05~1.4之間，平均值為1.213。而由表2-1可以看出，當(dāng)η=1

58、時(shí)，雷電發(fā)展分形維數(shù)Db=1.20±0.04，與實(shí)際放電路徑的分形維數(shù)最接近。</p><p>　　圖2-2 雷電先導(dǎo)模擬圖</p><p>　　但是對(duì)于η值得選取，其他文獻(xiàn)還有不同的見解，例如文[23]提出根據(jù)生活經(jīng)驗(yàn)，η大約取2或3之間時(shí)應(yīng)該更符合客觀情況。由于這不是本文研究的重點(diǎn)，所以在下文的編程中均采取文[23]的觀點(diǎn)，取η=3。</p><p>

59、　　第3章 雷電先導(dǎo)放電通道的二維模型</p><p>　　3.1 雷電先導(dǎo)放電的物理過程及參數(shù)</p><p>　　3.1.1雷電放電概述</p><p>　　雷電放電包括雷云對(duì)大地、雷云對(duì)雷云和雷云內(nèi)部的放電現(xiàn)象[]。大多數(shù)雷電放電是在雷云與雷云之間進(jìn)行的，只有少數(shù)是對(duì)地進(jìn)行的。我們主要研究雷云對(duì)大地的放電，因?yàn)檫@是造成雷害事故的主要因素。</p>

60、<p>　　雷云帶有大量電荷，由于靜電感應(yīng)作用，在雷云下方的地面或地面上的物體將感應(yīng)聚集與雷云相反的電荷，雷云與大地間就形成了電場(chǎng)。當(dāng)雷云中的電荷逐漸積累，達(dá)到一定的電荷密度，使其表面空間的電場(chǎng)強(qiáng)度足夠大(25~30kV∕cm)時(shí)，就發(fā)展局部放電。如果此地最大場(chǎng)強(qiáng)方向主要是對(duì)地的，就會(huì)發(fā)展對(duì)地的放電，形成下行雷。</p><p>　　雷電的極性是按照從雷云流入大地的電荷極性決定的。廣泛實(shí)測(cè)表明，90

61、%左右是負(fù)極性雷。</p><p>　　正雷出現(xiàn)的機(jī)會(huì)較少，故對(duì)正雷的研究也較少。正雷的下行先導(dǎo)階段沒有明顯的逐級(jí)發(fā)展的特征，正雷通常只有一個(gè)分量，十分難得有兩個(gè)分量的。正雷一般有很長(zhǎng)的波頭(長(zhǎng)達(dá)幾百微秒)和很長(zhǎng)的波尾(長(zhǎng)達(dá)上千微秒)，這樣，它所傳遞的電荷可能比多分量的負(fù)雷還多得多，而其電流陡度則比相應(yīng)的負(fù)雷小得多。</p><p>　　若低面上存在特別高聳的導(dǎo)電性能良好的物體時(shí)，也可能首

62、先從該物體頂端出發(fā)發(fā)展向上的先導(dǎo)(這種雷稱為上行雷)。但上行先導(dǎo)到達(dá)雷云時(shí)，一般不會(huì)發(fā)生主放電過程，這是因?yàn)槔自频姆烹娦阅鼙却蟮夭畹枚?，難以在極短時(shí)間內(nèi)供應(yīng)為中和先導(dǎo)通道中電荷所需的極大的主放電電流，而只能向雷云深處發(fā)展多分支的云中先導(dǎo)，通過寬廣區(qū)域的電暈流注，從分散的水性質(zhì)點(diǎn)上卸下電荷，匯集起來，以中和上行先導(dǎo)中的部分電荷。這樣的放電過程顯然只能是較緩和的，而不可能具有大沖擊電流的特性，而放電電流一般不足千安，而延續(xù)時(shí)間則較長(zhǎng)，有的可

63、能長(zhǎng)達(dá)10-1s。此外，上行先導(dǎo)從一開始就出現(xiàn)分支的概率較大。</p><p>　　3.1.2 雷電的先導(dǎo)過程及參數(shù)</p><p>　　關(guān)于雷電先導(dǎo)的情況[27]，主要由展開照相得來，它與長(zhǎng)間隙擊穿的先導(dǎo)有定性的相似，也是由先導(dǎo)通道、先導(dǎo)頭部和流注區(qū)三部分組成。先導(dǎo)過程延續(xù)約幾毫秒，以逐級(jí)發(fā)展的、高電導(dǎo)的、高溫的、具有高電位的先導(dǎo)通道將雷云到大地的氣隙擊穿。</p><

64、;p>　　整塊雷云可以有若干個(gè)電荷中心[]，負(fù)電荷中心位于雷云的下部，距地面500~10000m的范圍內(nèi)，直接擊向地面的放電通常從負(fù)電荷中心邊緣開始。</p><p>　　下行負(fù)先導(dǎo)具有明顯的分級(jí)發(fā)展的特點(diǎn)[27]。每級(jí)長(zhǎng)度在10~200m范圍內(nèi)，平均為40m.。相鄰兩級(jí)間歇時(shí)間為30~90μs，平均為60μs。分級(jí)前伸速度約為(1~5)×107m/s，延續(xù)約為1μs。由于存在間歇，其總的平均速度

65、約為(1~8)×105m/s，先導(dǎo)過程歷時(shí)約0.005~ 0.01s[]。下行負(fù)先導(dǎo)向地面推進(jìn)時(shí)還可能出現(xiàn)一些分支，但通常只有其中的一支到達(dá)地面。</p><p>　　梯級(jí)先導(dǎo)形成的原因一般解釋為[29]：在先導(dǎo)通道形成之前有一看不見的引路先導(dǎo)(流注)。先導(dǎo)頭部有很強(qiáng)的電場(chǎng)，引路先導(dǎo)在電場(chǎng)的驅(qū)動(dòng)下向下發(fā)展，由于空間電荷區(qū)的消耗，又因先導(dǎo)頂端聚集了大量正電荷，使得先導(dǎo)通道頂端處的電場(chǎng)大為減弱甚至出現(xiàn)反向。

66、此時(shí)引路先導(dǎo)便停止發(fā)展，而電子會(huì)退到先導(dǎo)頭部的正電荷區(qū)，產(chǎn)生強(qiáng)烈的復(fù)合，形成很強(qiáng)的光電離，先導(dǎo)就沿著引路先導(dǎo)形成的通道向前發(fā)展，完成一次梯級(jí)過程，其頭部強(qiáng)電場(chǎng)又將恢復(fù)。這個(gè)過程重復(fù)下去，就是梯級(jí)先導(dǎo)過程經(jīng)典的物理圖象。從另一個(gè)角度看，下行雷的放電是一種超長(zhǎng)間隙放電，放電過程也由先導(dǎo)-流注系來完成。通常的下行負(fù)雷擊從雷云下部負(fù)電荷中心的某處開始，以梯級(jí)先導(dǎo)放電的形式向地面發(fā)展，各梯級(jí)的方向有明顯的垂直分量，但一般來說，每個(gè)梯級(jí)在空間的取向

67、是隨機(jī)的，最終使放電通道成為大家習(xí)見的曲折形狀。各個(gè)梯級(jí)的方向很可能由放電通道頂端周圍的局部電場(chǎng)分布所確定，所以雷電的第一次閃擊幾乎總是有重重分枝，而這些分枝也一樣按梯級(jí)形式發(fā)展，以幾十米的梯級(jí)向下延伸。因此，雷電的發(fā)展具有明顯的分形特性。</p><p>　　很難直接測(cè)出雷云電荷中心對(duì)地的電勢(shì)[27]，但可以間接推算出：先導(dǎo)根部(在云中)對(duì)地的電勢(shì)約為50~100MV；先導(dǎo)頭部對(duì)地電勢(shì)為20~80MV。<

68、/p><p>　　當(dāng)下行雷電導(dǎo)從雷云向建筑物方向發(fā)展時(shí)，從接地的建筑物上可能產(chǎn)生向上的迎面先導(dǎo)。產(chǎn)生初始迎面先導(dǎo)的條件與上述產(chǎn)生上行雷的條件相似。區(qū)別在于：前者不僅考慮雷云電場(chǎng)的作用，還應(yīng)考慮下行先導(dǎo)電荷造成局部電場(chǎng)的作用。</p><p>　　迎面先導(dǎo)在相當(dāng)大的程度上影響著下行先導(dǎo)的發(fā)展路線并決定雷擊點(diǎn)的所在，所以，它在對(duì)地雷擊的發(fā)展中具有很重要的意義。</p><p&g

69、t;　　直擊到平地的雷電，實(shí)際上不存在從地面向上發(fā)展的迎面先導(dǎo)。實(shí)際觀測(cè)統(tǒng)計(jì)指出，平原地區(qū)高度在200m以上的建筑物上，觀測(cè)到相當(dāng)多的上行雷電。</p><p><b>　　3.2 模型的建立</b></p><p>　　3.2.1 理想的模擬場(chǎng)景、步驟和參數(shù)</p><p>　　在文[24]中只考慮比較常見的下行負(fù)極性雷的情形，給出了一個(gè)理想

70、的雷電先導(dǎo)模型。雷電通道的形成分為下行先導(dǎo)和上行先導(dǎo)兩部分，因而模擬時(shí)在一般氣體放電模擬的基礎(chǔ)上考慮多條上行先導(dǎo)起始、發(fā)展以及上、下行先導(dǎo)頭部相接(末躍)等的實(shí)現(xiàn)。下行先導(dǎo)從雷云底部起始，當(dāng)其逼近大地時(shí)，上行先導(dǎo)會(huì)從地面及輸電線路的某位置起始并向上發(fā)展。模擬中，大地以及地面上構(gòu)筑物(如桿塔、避雷線等)設(shè)為零電位，雷云底部等效成一個(gè)平板電極，電位為?U0。另外要在雷云底部隨機(jī)選擇一點(diǎn)，產(chǎn)生一段初始的下行先導(dǎo)通道。在每一時(shí)間步長(zhǎng)，檢查下行先

71、導(dǎo)頭部與地面之間的平均場(chǎng)強(qiáng)，如果超過擊穿臨界值Eg，則雷電將直擊大地。檢查地面目標(biāo)上可能產(chǎn)生上行先導(dǎo)的點(diǎn)處的場(chǎng)強(qiáng)，如果其中某點(diǎn)與一鄰近點(diǎn)之間的平均場(chǎng)強(qiáng)大于上行先導(dǎo)起始的臨界場(chǎng)強(qiáng)Eu，則此處出現(xiàn)擊穿，上行先導(dǎo)起始。對(duì)已經(jīng)起始的每條上行先導(dǎo)，按照類似于下行先導(dǎo)的方式發(fā)展。計(jì)算每條上行先導(dǎo)的頭部與下行先導(dǎo)的頭部之間的平均場(chǎng)強(qiáng)，設(shè)其最大值E(k)對(duì)應(yīng)于第k條上行先導(dǎo)。如果E(k)不小于Ef(Ef為末躍場(chǎng)強(qiáng)判據(jù))，則發(fā)生末躍。單次模擬結(jié)果顯然具有

72、很大的隨機(jī)性。為使得模擬具有實(shí)際價(jià)值，須多次模擬得到統(tǒng)計(jì)</p><p>　　3.2.2 模型的簡(jiǎn)化</p><p>　　雖然現(xiàn)有雷電放電路徑發(fā)展觀察資料表明，雷電先導(dǎo)發(fā)展過程表現(xiàn)出梯級(jí)向下發(fā)展的特點(diǎn)。但由于在計(jì)算機(jī)先導(dǎo)發(fā)展模擬中，忽略先導(dǎo)發(fā)展梯級(jí)間隔對(duì)先導(dǎo)路徑的選擇及最后雷擊目標(biāo)物選擇等關(guān)鍵問題并沒有明顯影響。因此，本模型將雷電先導(dǎo)梯級(jí)發(fā)展簡(jiǎn)化為連續(xù)向下發(fā)展。上文指出，通過對(duì)實(shí)際雷云向地

73、面放電過程的觀察得出，雷電先導(dǎo)梯級(jí)長(zhǎng)度為10~ 200m。對(duì)于計(jì)算機(jī)模擬計(jì)算，將計(jì)算步長(zhǎng)設(shè)置小于1個(gè)梯級(jí)發(fā)展長(zhǎng)度便可滿足模擬先導(dǎo)路徑精度要求，進(jìn)一步減小計(jì)算步長(zhǎng)可以提高模擬精度，但過小的計(jì)算步長(zhǎng)會(huì)造成較大的計(jì)算量。因此,本模型綜合考慮模擬精度及計(jì)算機(jī)計(jì)算量,將計(jì)算步長(zhǎng)取為雷電梯級(jí)長(zhǎng)度的最小值，即為10m。</p><p>　　圖3-1 二維格點(diǎn)電勢(shì)圖</p><p>　　本模型是一個(gè)從概率

74、角度去描述介質(zhì)電擊穿的模型。如圖3-1所示，每一個(gè)點(diǎn)都有可能擊穿離它最近的八個(gè)點(diǎn)，每擊穿一個(gè)點(diǎn)，擊穿點(diǎn)的電勢(shì)就可以認(rèn)為和擊穿源相同。由于雷雨天氣中云層頂部帶正電，底部帶負(fù)電，所以會(huì)使得地面成為高電勢(shì)，云層成為低電勢(shì)，所以在計(jì)算中我們不妨假設(shè)云層的電勢(shì)為0，地面的電勢(shì)為50MV(中等雷電)。天地間簡(jiǎn)化為一個(gè)平面。閃電的路徑為擊穿的點(diǎn)，由于閃電從云端一直到地面，第一個(gè)被擊穿的點(diǎn)在云端，可以認(rèn)為所有擊穿的點(diǎn)的電勢(shì)都為0，其概率滿足(2-15)

75、式。但是因?yàn)槲覀兲幚淼亩际请妱?shì)，所以改寫(2-15)式可以得到：</p><p><b>　　(3-1)</b></p><p>　　從點(diǎn)(i,j)出發(fā)，如果僅僅考慮對(duì)(i+1,j)，(i-1,j)，(i,j-1)，(i,j+1)點(diǎn)的擊穿，很明顯，當(dāng)網(wǎng)格均勻時(shí)，所有di都相等，(3-1)式就蛻變?yōu)橄率剑?lt;/p><p><b>　　(

76、3-2)</b></p><p>　　再考慮加入(i+1,j-1)，(i-1,j-1)，(i+1,j+1)，(i-,j+1)四個(gè)點(diǎn)，考慮到(3-1)式中距離di 的影響，我們對(duì)這四個(gè)點(diǎn)只要把φ改寫為φ∕就可以直接應(yīng)用式(3-2)求解。</p><p>　　本課題模擬的二維雷電先導(dǎo)從1000米高空開始，分級(jí)前伸速度300000m/s。我們假設(shè)第一個(gè)擊穿的點(diǎn)在云端的中間，每擊穿一個(gè)

77、點(diǎn)，我們把已被擊穿的點(diǎn)計(jì)入邊界，這滿足拉普拉斯方程，所以我們?cè)儆糜邢薏罘址ǖ匦掠?jì)算平面的電勢(shì)分布，然后依照新算出來的電勢(shì)再配合隨機(jī)函數(shù)得到下一個(gè)擊穿的點(diǎn)。每個(gè)已被擊穿的點(diǎn)都可能會(huì)作為生長(zhǎng)點(diǎn)，使最近的八個(gè)點(diǎn)之一被擊穿。最終擊穿點(diǎn)到達(dá)邊界則計(jì)算結(jié)束。因?yàn)橛?jì)算電勢(shì)需要完整的邊界條件，所以我們做一個(gè)近似的假定，我們認(rèn)為在整個(gè)過程中平面左右兩個(gè)邊界上的電勢(shì)一直不會(huì)改變，即滿足最初的線性變化的分布。</p><p>&l

78、t;b>　　3.3 程序流程</b></p><p>　　根據(jù)以上的理論和假設(shè)可以做出程序流程圖，如圖3-2。</p><p>　　顯然，該算法包括以下3大模塊：</p><p><b>　　模塊1：初始情況</b></p><p><b>　　模塊2：計(jì)算電勢(shì)</b></p

79、><p>　　模塊3：計(jì)算概率，決定下一個(gè)擊穿點(diǎn)</p><p>　　圖3-2 程序流程圖</p><p><b>　　3.3.1 模塊1</b></p><p>　　clc %清除命令窗口中顯示的所有輸入和輸出</p><p>　　tic;

80、 %啟動(dòng)秒表計(jì)時(shí)器，輸出tic使用后的流</p><p>　　時(shí)間，與toc結(jié)合記錄程序運(yùn)行時(shí)間</p><p>　　lx=100;ly=100; %定義網(wǎng)格的最大行、列數(shù) </p><p>　　w=2/(1+sin(pi/(lx-1))); %加速收斂

81、因子</p><p>　　v1=zeros(ly,lx); %v1初始化，使其成為全零矩陣</p><p>　　v=v1; %將v1賦給v</p><p>　　r=3; %發(fā)展概率指數(shù)</p><p>

82、;　　c=3e5; %分級(jí)前伸速度</p><p>　　T=0; %時(shí)間矩陣</p><p>　　tt=0; %返回雷電先導(dǎo)歷時(shí)</p><p><b>　　temp=0;</b><

83、/p><p>　　e1=zeros(1,97); %求取一條線路上相鄰兩點(diǎn)之間的場(chǎng)強(qiáng)，放在e1中</p><p>　　%-----------------------------------------使電勢(shì)從上到下線性遞增----------------------------------------</p><p>　　for i

84、=1:ly</p><p>　　for j=1:lx</p><p>　　v1(i,j)=(ly-i)*50/(ly-1); </p><p><b>　　end </b></p><p><b>　　end </b></p><p>　　o=sqrt(2);<

85、/p><p>　　v(1,:)=1;v(:,1)=1;v(ly,:)=1;v(:,lx)=1; %邊界</p><p>　　k=0; %迭代次數(shù)</p><p>　　i=99; %起始點(diǎn)的縱坐標(biāo)</p><p>　　

86、j=fix(lx/2); %起始點(diǎn)的橫坐標(biāo)</p><p>　　第3行中天地間開辟了一個(gè)100×100的空間，由于上節(jié)中指出模擬的雷電從1000米高空開始，時(shí)間步長(zhǎng)設(shè)置為10米，所以每一電位坐標(biāo)代表10米。</p><p><b>　　3.3.2 模塊2</b></p><p>　　

87、命xi為第i點(diǎn)電勢(shì)，Laplace方程的離散形式可以寫成：</p><p><b>　　(3-3)</b></p><p>　　則對(duì)全部n個(gè)點(diǎn)，可以列n個(gè)方程：</p><p><b>　　(3-4)</b></p><p>　　方程組中，第i個(gè)方程中xj前的系數(shù)aj如此給定：當(dāng)j=i時(shí)，aj=-4

88、；當(dāng)j為i鄰近的四個(gè)點(diǎn)之一時(shí)，aj=1；其他情況，aj=0。即有n個(gè)方程，每個(gè)方程中的n個(gè)元素中，除了這個(gè)點(diǎn)和與它最鄰近的4個(gè)點(diǎn)外，都為0元；把給定的m個(gè)邊界條件(包括擊穿點(diǎn))的電勢(shì)值代入方程組，剩余的n-m獨(dú)立方程，恰能給出剩余的n-m個(gè)未知網(wǎng)格的電勢(shì)；這樣就把求各點(diǎn)電勢(shì)問題轉(zhuǎn)化為求解線性方程組的問題了。</p><p>　　由于我們用到的差分公式本身就是近似公式，求出該方程的精確解沒有實(shí)際意義，我們只需要一個(gè)

89、符合精度要求的近似解即可，我們所以我們應(yīng)該考慮一個(gè)效率更高的，更容易實(shí)現(xiàn)的算法，那么可以考慮使用迭代法[23]。</p><p><b>　　迭代程序如下：</b></p><p>　　while(maxt>1e-5) %精度要求，達(dá)到精度跳出循環(huán)</p><p><b>　　k=k+1; </b></p

90、><p><b>　　maxt=0; </b></p><p>　　for m=2:ly-1 </p><p>　　for n=2:lx-1; </p><p>　　if v(m,n)==0 </p><p>　　v2(m,n)=(v1(m,n+1)+v1(m+1,n)+v1(m-1,n)+v1(m,

91、n-1))/4;</p><p>　　t=abs(v2(m,n)-v1(m,n)); </p><p>　　if(t>maxt) maxt=t; </p><p><b>　　end </b></p><p><b>　　end</b></p><p><b

92、>　　end</b></p><p><b>　　end</b></p><p><b>　　v1=v2;</b></p><p><b>　　end</b></p><p>　　為了加快迭代速度，可以提前使用新值進(jìn)行迭代，因?yàn)樵谏鲜龀绦蛑?，我們的掃描賦值方式

93、是一行一行掃描，在對(duì)v2(i，j)賦值時(shí)，它周圍四個(gè)點(diǎn)其中有2個(gè)已經(jīng)被掃描過了，即已經(jīng)獲得了新的數(shù)值，這個(gè)數(shù)值應(yīng)該更優(yōu)，所以我們盡量使用新算出的數(shù)值進(jìn)行迭代，因此，只要把上述代碼中：</p><p>　　v2(i,j)=(v1(i,j+1)+v1(i+1,j)+v1(i-1,j)+v1(i,j-1))/4</p><p><b>　　改為：</b></p>

94、<p>　　v2(i,j)=(v1(i,j+1)+v1(i+1,j)+v2(i-1,j)+v2(i,j-1))/4</p><p>　　就可以將迭代次數(shù)減少。</p><p>　　為了能使迭代次數(shù)進(jìn)一步減少，考慮讓每次迭代獲得更多的增量，可以采用逐次超松弛法[]：</p><p><b>　　那么將：</b></p>

95、<p>　　v2(i,j)=(v1(i,j+1)+v1(i+1,j)+v2(i-1,j)+v2(i,j-1))/4</p><p>　　變形為：v2(i,j)=v1(i,j)+(v1(i,j+1)+v1(i+1,j)+v2(i-1,j)+v2(i,j-1)-4*v1(i,j))/4</p><p>　　再考慮：v2(i,j)=v1(i,j)+w*(v1(i,j+1)+v1(i+

96、1,j)+v2(i-1,j)+v2(i,j-1)-4*v1(i,j))/4</p><p>　　其中w稱為“加速收斂因子”，是一個(gè)供選擇的參數(shù)，其值在1~2之間，該方法的快慢與w有著明顯的關(guān)系。實(shí)踐表明，如果w選得好，可以較快地加快迭代速度。經(jīng)驗(yàn)表明正方形場(chǎng)域由正方形網(wǎng)格劃分，每邊節(jié)點(diǎn)數(shù)為(p+1)時(shí)，最佳收斂因子w為</p><p><b>　　(3-5)</b>&

97、lt;/p><p>　　程序最終的電位迭代程序?yàn)椋?lt;/p><p>　　while(maxt>1e-5) </p><p><b>　　k=k+1; </b></p><p><b>　　maxt=0; </b></p><p>　　for m=2:ly-1 </p&

98、gt;<p>　　for n=2:lx-1; </p><p>　　if v(m,n)==0 </p><p>　　v2(m,n)=v2(m,n)+w*(v1(m,n+1)+v1(m+1,n)+v2(m-1,n)+v2(m,n-1)-4*v2(m,n))/4;</p><p>　　t=abs(v2(m,n)-v1(m,n)); </p>

99、<p>　　if(t>maxt) </p><p><b>　　maxt=t; </b></p><p><b>　　end </b></p><p><b>　　end</b></p><p><b>　　end</b></p

100、><p><b>　　end</b></p><p><b>　　v1=v2;</b></p><p><b>　　end</b></p><p><b>　　3.3.3 模塊3</b></p><p>　　結(jié)合式(3-2)，設(shè)變量X

101、表示可能擊穿的n個(gè)點(diǎn)，則這n個(gè)點(diǎn)的概率可以寫成：</p><p><b>　　(3-6)</b></p><p>　　我們以列表的形式表示出：</p><p>　　表3-1 可能擊穿點(diǎn)的概率分布表</p><p>　　從中顯而易見， 作為概率函數(shù)的pi滿足：</p><p><b>　　

102、(3-7)</b></p><p>　　于是，我們命F(x)作為概率分布函數(shù)有：</p><p><b>　　(3-8)</b></p><p>　　由rand(1,1)生成一個(gè)隨機(jī)數(shù)，這個(gè)數(shù)必然在[0,1]區(qū)間中隨機(jī)選取，而且這個(gè)數(shù)必然對(duì)應(yīng)于唯一一個(gè)可能擊穿的點(diǎn)，由此即可判定x所在的區(qū)間，進(jìn)而判定pi以及所對(duì)應(yīng)的ai，這樣，我們就

103、確定了下一個(gè)擊穿點(diǎn)，命這個(gè)擊穿點(diǎn)Φ=0。</p><p><b>　　相關(guān)程序如下：</b></p><p><b>　　p=0;</b></p><p>　　for i=2:ly-1 </p><p>　　for j=2:lx-1 </p><p>　　if (v(i,j)

104、==2) </p><p>　　if ( v(i-1,j)==2) </p><p><b>　　p1=0; </b></p><p>　　else p1=v1(i-1,j)^r; </p><p><b>　　end; </b></p><p>　　if (v(i+1,j)

105、==2) </p><p><b>　　p2=0; </b></p><p>　　else p2=v1(i+1,j)^r;</p><p><b>　　end; </b></p><p>　　if (v(i,j-1)==2) </p><p><b>　　p3=0;

106、 </b></p><p><b>　　else</b></p><p>　　p3=v1(i,j-1)^r;</p><p><b>　　end; </b></p><p>　　if (v(i,j+1)==2) </p><p><b>　　p4=0;

107、</b></p><p><b>　　else</b></p><p>　　p4=v1(i,j+1)^r;</p><p><b>　　end; </b></p><p>　　if ( v(i-1,j-1)==2) </p><p><b>　　p5=0

108、; </b></p><p><b>　　else</b></p><p>　　p5=(v1(i-1,j-1)/o)^r; </p><p><b>　　end; </b></p><p>　　if ( v(i-1,j+1)==2) p6=0; </p><p>

109、　　else p6=(v1(i-1,j+1)/o)^r; </p><p>　　end; if ( v(i+1,j-1)==2) </p><p><b>　　p7=0; </b></p><p>　　else p7=(v1(i+1,j-1)/o)^r; </p><p><b>　　end; </b&g

110、t;</p><p>　　if ( v(i+1,j+1)==2) p8=0; </p><p>　　else p8=(v1(i+1,j+1)/o)^r; </p><p>　　end; p=p1+p2+p3+p4+p5+p6+p7+p8+p; </p><p><b>　　end</b></p><p

111、><b>　　end;</b></p><p><b>　　end; </b></p><p>　　q=rand(1,1); %取0~1之間的隨機(jī)數(shù)</p><p><b>　　s=0; </b></p><p>　　%---------

112、-------------------------------計(jì)算每一個(gè)可能擊穿點(diǎn)的概率-----------------------------------</p><p>　　for i=2:ly-1 </p><p>　　for j=2:lx-1 </p><p>　　if (v(i,j)==2) </p><p>　　if ( v(i

113、-1,j)==2) p1=0; </p><p>　　else p1=(v1(i-1,j)^r/p); end; </p><p>　　if (v(i+1,j)==2) p2=0;</p><p>　　else p2=(v1(i+1,j)^r/p);end; </p><p>　　if (v(i,j-1)==2) p3=0;</p>

114、;<p>　　else p3=(v1(i,j-1)^r/p);end; </p><p>　　if (v(i,j+1)==2) p4=0; </p><p>　　else p4=(v1(i,j+1)^r/p);end; </p><p>　　if ( v(i-1,j-1)==2) p5=0; </p><p>　　else p5

115、=((v1(i-1,j-1)/o)^r)/p; end; </p><p>　　if ( v(i-1,j+1)==2) p6=0; </p><p>　　else p6=((v1(i-1,j+1)/o)^r)/p; end; </p><p>　　if ( v(i+1,j-1)==2) p7=0; </p><p>　　else p7=((v

116、1(i+1,j-1)/o)^r)/p; end; </p><p>　　if ( v(i+1,j+1)==2) p8=0; </p><p>　　else p8=((v1(i+1,j+1)/o)^r)/p; end; </p><p>　　%--------------------------------------結(jié)合隨機(jī)函數(shù)選擇下一個(gè)擊穿點(diǎn)-----------

117、-----------------------</p><p>　　s=s+p1+p2+p3+p4+p5+p6+p7+p8; </p><p>　　if (q>=s-p1-p2-p3-p4-p5-p6-p7-p8) && (q<s-p2-p3-p4-p5-p6-p7-p8) </p><p>　　m=i-1; n=j; x0=j;y0=i