極坐標系下靜態(tài)場的求解方法

楊琳

(商洛學院電子信息與電氣工程學院，陜西商洛726000)

極坐標系下靜態(tài)場的求解方法

楊琳

(商洛學院電子信息與電氣工程學院，陜西商洛726000)

針對圓形場域在直角坐標系下難以求解這一問題，提出了在極坐標系下應用FDTD對圓形靜態(tài)場求解。介紹了差分法的基本思想，在此基礎(chǔ)上導出極坐標系下的差分方程和算法，并在MATLAB環(huán)境下仿真，驗證了算法的有效性。

FDTD；極坐標；靜態(tài)場；MATLAB

靜態(tài)場的求解方法有很多種，鏡像法、分離變量法、有限差分法等，其中FDTD[1]在求解靜態(tài)場問題方面具有相當強的優(yōu)越性。但這種方法的計算量很大，人工計算費時費力，且精度難以到達實際需求。而Matlab[2]軟件在數(shù)據(jù)運算方面很有優(yōu)勢，可以迅速地編制出有效的數(shù)據(jù)處理程序。結(jié)合FDTD與Matlab的優(yōu)勢，可以在計算機中迅速得到高精度的運算結(jié)果。近些年來，很多學者在利用MATLAB求解靜態(tài)場方面進行了研究。熊彬等[3]研究了矩形場域的FDTD算法。趙德奎等[4]研究了正六邊形場域的FDTD算法。這些算法均在直角坐標系下實現(xiàn)。而實際工程應用當中，場域形狀多種多樣，圓形、環(huán)形、扇形等場域更適合在極坐標系下研究。本文在MATLAB環(huán)境中實現(xiàn)了在極坐標系下應用FDTD對圓環(huán)形靜態(tài)場的求解。一方面，在前人研究的基礎(chǔ)之上導出了極坐標系下靜態(tài)場的差分公式，把復雜的二階偏微分方程問題轉(zhuǎn)化為簡單的線性方程組問題，還用迭代法來求解線性方程組，并且用加入松弛因子的方式加快收斂速度，減少運算時間。另一方面，在MATLAB軟件環(huán)境中仿真求解環(huán)形區(qū)域二維場，不僅減少了運算量，提高了計算結(jié)果的精確程度，還可以模擬出場值分布的圖形，形象直觀地反映出場域性質(zhì)。另外，這種方法具有良好的通用性和普遍性，僅僅需要修改網(wǎng)格的初始參數(shù)和邊界條件，即可輕松地移植到其他類似圓形區(qū)域，從而為這一類場域的靜態(tài)場求解提供了一個參考。

1 FDTD的基本思想

有限差分法[5]既簡單又直觀，其思路是將求解場域劃分成網(wǎng)格，把場域的場值連續(xù)分布轉(zhuǎn)化為網(wǎng)格節(jié)點的場值離散分布。該方法求解出的離散點的場值，不是精確值，而是近似值。而其核心則是將復雜的場域的偏微分問題轉(zhuǎn)化為簡單的相鄰離散點之間場值差值問題[6]。求解結(jié)果精度取決于網(wǎng)格劃分方式。網(wǎng)格劃分得越細，求解結(jié)果精度越高，耗時越長；與之相反網(wǎng)格劃分得越粗，求解結(jié)果精度越低，耗時越少。

2 極坐標系下差分方程的推導

2.1 網(wǎng)格劃分

為了方便地求解，采用極坐標系來將圓形或扇形區(qū)域劃分成圓環(huán)狀網(wǎng)格。如圖1所示，取ρ作徑向變量，θ作角度變量。角度步長為Δθ，徑向步長為Δρ，等ρ線為半徑相等的弧線，等θ線為角度相等的徑向射線，ρ=mΔρ，θ=nΔθ(m,n=0,1,2,3,…)。整張平面被這兩種線劃分成網(wǎng)格。每個部分圓環(huán)格為一個單元，網(wǎng)格線的交叉點是節(jié)點，某節(jié)點的位置坐標是(ρ，θ)。圖1給出了扇形區(qū)域的網(wǎng)格劃分。設(shè)0點電位φ0的序號坐標是(m,n)，那么四周四個點電位φ1，φ2，φ3，φ4的序號坐標分別為(m+1,n)、(m,n+1)、(m-1,n)、(m,n-1)。

圖1 扇形區(qū)域差分網(wǎng)格

2.2 差分方程推導

根據(jù)泰勒級數(shù)展開式，一個函數(shù)可以在某一點展開成(1)式的形式，

假設(shè)φ0=φm，n，忽略比二次更高的項，那么四周四個點電位φ1，φ2，φ3，φ4可以用φ0表示成(2)(3)(4) (5)式。

由(1)+(3)式，可得

由(1)-(3)式，可得

由(2)+(4)式，可得

由(2)-(4)式，可得

由(6)(7)(8)(9)式，可得

由于拉普拉斯方程適用于電位函數(shù)，可得

將(10)(11)(12)式代入到(14)式，可得

約分可得點的差分方程

2.3 差分方程求解

根據(jù)(17)式，φ0的第t+1次結(jié)果可由φ1，φ2，φ3，φ4的第t次結(jié)果得出。進行N次迭代之后就可保證相鄰兩次迭代結(jié)果的差值在允許的范圍之內(nèi)，此時結(jié)束迭代即可。但這種迭代方式收斂慢、迭代次數(shù)多、所需要的計算時間比較長。為了節(jié)省資源，減少計算時間，使用加權(quán)迭代法[7]，即加入松弛因子的方式來加快收斂速度，從而減少計算時間、節(jié)省資源。

該方法可用(18)式來表示

松弛因子的取值范圍是[0,2]，在這個范圍之內(nèi)存在一個最優(yōu)松弛因子[8]使得收斂速度最快。獲得最優(yōu)松弛因子的方法有很多種，常見的有二分法、逐步搜索法、黃金分割法等。根據(jù)相關(guān)資料對于符合一定條件的區(qū)域，最優(yōu)松弛因子可以通過(20)式來確定，式中n1、n2代表兩個網(wǎng)格方向的點數(shù)。

3 場域的求解與仿真

根據(jù)上述分析，可知圓、環(huán)、扇形區(qū)域的靜態(tài)場求解過程如下：

1)劃分場域網(wǎng)格；

2)根據(jù)邊界條件對場域賦初值；

3)確定最優(yōu)松弛因子；

4)根據(jù)(18)式、(19)式以及最佳松弛因子，進行場域的若干次迭代，在相鄰的兩次迭代結(jié)果的差值小于規(guī)定的精度時結(jié)束迭代。

假設(shè)有如圖2所示的圓環(huán)狀導體槽，其內(nèi)環(huán)半徑為50 mm，外環(huán)半徑為100 mm，內(nèi)環(huán)電位為100 V，外環(huán)電位為0 V?？梢杂蒙鲜龇椒ㄇ蠼獬鰣A環(huán)區(qū)域內(nèi)的電位分布。

圖2 圓環(huán)形導體槽

求解區(qū)域是一個圓環(huán)形狀，設(shè)置其徑向步長為5 mm，角度步長為15°。根據(jù)上述求解方法，編制出Matlab程序來求解圓環(huán)區(qū)域內(nèi)靜態(tài)場的電位分布。其Matlab仿真程序如下：

仿真結(jié)果如圖3所示，在內(nèi)環(huán)上電位為100 V，外環(huán)上電位為0 V，中間區(qū)域的電位呈遞減趨勢，仿真結(jié)果與理論分析相一致。

圖3 仿真結(jié)果

在仿真過程中，徑向點數(shù)為11，角度點數(shù)為25，網(wǎng)格總點數(shù)為275。迭代精度為10-6V。在不加松弛因子的情況下，需要329次可完成迭代。在加入松弛因子的情況下，僅需要55次迭代即可得到求解結(jié)果?？梢钥闯觯褂眉訖?quán)迭代法可以很大程度地減小迭代次數(shù)，節(jié)省計算時間。

4 結(jié)語

本文根據(jù)圓環(huán)形場域的特征和FDTD的基本思想（即將偏導數(shù)替換為差商），在極坐標系下對場域進行離散化，推導出相應的差分方程，并在MATLAB環(huán)境下，求解了圓環(huán)形場域靜電場問題，不僅大大減少了計算量，還使計算結(jié)果更加精確?？傊疚乃鏊惴ê头椒?，可擴展應用到更復雜的扇形、環(huán)形、柱形等其他類似圓形區(qū)域中，使得求解此類問題變得更加便捷。

[1]Dennis M.Sullivan,electromagnetic simulation using the fdtd method[M].New York:IEEE Publications,2000:3-6.

[2]Mehmet U,Levent S.Matlab-based virtual wedge scattering tool for the comparison of high frequency asymptotics and FDTD method[J].Applied Computational Electromagnetics Society journal,2012,27(9):697-705.

[3]熊彬,阮百堯.MATLAB在有限差分法中的應用[J].桂林工學院學報,2001,21(2):104-109.

[4]趙德奎,劉勇.MATLAB在有限差分法數(shù)值計算中的應用[J].四川理工學院學報:自然科學版,2005,18(4):61-64.

[5]葛德彪,閆玉波.電磁場時域有限差分方法[M].西安:西安電子科技大學出版社,2005:1-5.

[6]何紅雨.電磁場數(shù)值計算法與MATLAB實現(xiàn)[M].武漢:華中科技大學出版社，2004:93-99.

[7]孫慧玲.電磁場中的自適應超松弛迭代技術(shù)研究[D].南京：南京理工大學,2008:8-10.

[8]胡楓,金遠平.SOR最優(yōu)松弛因子選取方法研究[J].西南師范大學學報:自然科學版,2008,33(5):48-50.

（責任編輯：張國春）

The Study of Solution of Static Field Solver in the Polar Coordinates

YANG Lin
（CollegeofElectroncInformationandElectricalEngineering,ShangluoUinversity,Shangluo726000,Shaanxi）

As for the difficult solution concerning round field in the Cartesian coordinates,a new method to solve round static field by FDTD in the Polar Coordinates is proposed.The principle of the finite difference method is introduced,by which the differential equations and the algorithm in the polar coordinate are obtained.Through simulation by MATLAB,the algorithm was verified.

FDTD;polar coordinates;the static field;MATLAB

O412.3

：A

：1674－0033(2014)02－0024－04

10.13440/j.slxy.1674-0033.2014.02.008

2014-01-09

商洛學院教育教學改革研究項目(13JYJX105)

楊琳，女，陜西山陽人，碩士，助教