基于Matlab的非均勻介質靜電特性分析

章三妹　唐正明

摘 要：針對傳統的靜電場分析主要局限于同種介質，在介紹有限差分法求電位分布原理的基礎上，導出了非均勻介質有限差分法的計算公式，并運用Matlab編程分析了所設定的非均勻介質區域靜電特性。仿真結果與靜電場理論符合較好，有助于加深對電介質靜電特性的理解。

關鍵詞：有限差分法；非均勻介質；靜電場；Matlab

中圖分類號：TP301 文獻標識碼：A

Analyzing of Electrostatic Field of Inhomogeneous Medium based on Matlab

ZHANG Sanmei1，TANG Zhengming2

（1 Experiment Center， China West Normal University， Nanchong Sichuan 637009， China；

2 School of Electronics and Information Engineering， China West Normal University， Nanchong Sichuan 637009， China）

Abstract： Aiming at the situation that traditional analysis of electrostatic field is limited to single medium， through the introduction of the principle of finite difference method， the calculation formula for inhomogeneous medium is presented， simulation for the electrostatic field distribution of inhomogeneous medium is carried on Matlab as well. The results are consistent with electrostatic field theory， and this study is good for deeply understanding the characteristic of electrostatic field.

Keywords： Finite Difference Method； Inhomogeneous Medium；Electrostatic Field； Matlab

0引 言

隨著工程問題復雜度的提升及計算機處理能力的顯著提高，應用數值法分析電磁場問題日趨廣泛。Matlab具有程序設計簡單，圖像顯示和處理功能強大等優點，比較適合應用于電磁場問題的數值分析中。二維靜態電磁場的邊值問題是求解電磁場的基礎，然而，現有文獻所涉及的常為同種媒質的靜電場分析[1-4]。由于，有限差分法求解靜電場邊值問題較為有效，同時，其也是廣泛應用于時諧場分析的時域有限差分法的基礎[5-7]。本文即運用Matlab結合有限差分法分析了非均勻介質填充區域的靜電場特性，并進一步給出了處理該類問題的程序框圖及程序實例。

1有限差分法的原理[8]

設函數 的獨立變量 存在微小增量 ，則該函數的增量可表示為：

（1）

其一階差商為：

（2）

由于在增量 較小的情況下，微分可近似于差商，即：

（3）

此式常稱為前向差商，類似地可以得到后向差商、中心差商。不同差商方式近似精度也不同，可通過泰勒級數來加以比較。由泰勒公式有：

（4）

（5）

（4）、（5）兩式相減可得：

（6）

可見，中心差商略去的是3階以上的無窮小，其它差商形式的截斷誤差可按類似原理進行分析。相對前向差商和后向差商來說，中心差商截斷誤差最小，在問題無特殊限定情況下較為常用，故本文以中心差分為例展開分析。按照類似方法，可得到二階導數的差商形式：

（7）

由于偏微分方程也可以做類似的差分表示，而Maxwell方程由偏微分方程所描述。因此，有限差分法可用于靜電場問題分析中。

2非均勻介質中的有限差分法形式

將求解區域劃分成網格，然后將區域內連續的場分布，按照有限差分法理論，用網格節點上的離散的數值解代替。差分網格除常規的矩形、正方形網格外，還有三角形等網格[9]。如圖1所示，為了描述方便，將區域采用邊長為 的正方形網格離散（ 即為離散步長，因其直接影響到數值解的精度，要求 充分?。?，并將相應節點放大顯示，其中，坐標為（i，j）的節點電位的電位為 ，其周圍4個節點的電位分別為 、 、 、 。

圖1 節點順序圖

Fig.1 Node sequence diagram

各個節點的電位可表示成以 為基點的泰勒級數形式[10]。

（8）

（9）

（10）

（11）

將以上4式相加后，由于 充分小，略去4階以上的項已可保證相當高的計算精度，可得

（12）

由于區域中的任何點均滿足泊松方程：

（13）

其中， 為場源。為簡化分析，假定問題域為無源區域，根據場強 與電位的關系 及電位移矢量所滿足的方程 ，利用恒等式 ，得：

（14）

（15）

聯立（12）、（13）和（15）得非均勻介質中電位的有限差分形式：

（16）

運用Matlab對（16）編程，可以得出求解區域的電位分布，同時，利用其特有的場分析命令，可以較為便捷地得到其電場分布，并予以顯示。

3非均勻介質填充區域的程序實例

加載介質的軸對稱電磁場問題是現今電工設備設計分析中的典型問題。如圖2所示，無限長矩形腔，四周為理想接地導體，橫截面為3m×2m，其中央有一截面方向1m長，直流電壓為1V的無限長平板，由相對介電常數為29的介質區域所支撐，腔的其余部分為空氣填充，相對介電常數為1。

圖2 求解區域

Fig.2 Solving region

將計算區域劃分為30*20個單元格，非均勻介質中電位的有限差分形式用賽德爾迭代法表示，在程序設計上參照文獻[11]采用Matlab編程。程序流程圖如圖3所示，主要程序代碼如下：

clear all；

hx=31；hy=21；

cx1=6；cx2=26；cy1=11； %設置網格節點數和結構參數

dx=cx2-cx1+1；

v1=zeros（hy，hx）； %初始化輸出結果

mark=ones（hy，hx）*2； %設置標志初值

v1（cy1，cx1：cx2）=ones（）； %初始化邊值（電位，標志及介電常數）

mark（cy1，cx1：cx2）=zeros（）；

epsl=ones（hy，hx）*8.85e-12；

epsl（1：cy1，cx1：cx2）=ones（cy1，dx）*19*8.85e-12；

epsl（1：cy1，cx1）=10*8.85e-12；

epsl（1：cy1，cx2）=10*8.85e-12；

v2=v1；maxt=1； %初始化

k=0；

while （maxt>1e-12） %設置由v1迭代算出v2的迭代誤差

k=k+1； %改變迭代次數

maxt=0；

for i=2：hy-1

for j=2：hx-1

m=mark（i，j）；

if（m>=2） %計算滿足標志條件的區域 v2（i，j）=（v1（i，j+1）+v1（i+1，j）+v2（i-1，j）+v2（i，j-1））/4+（ epsl（i+1，j）-epsl（i-1，j））*（v1（i+1，j）-v2（i-1，j））/16/epsl（i，j） +（ epsl（i，j+1）-epsl（i，j-1））*（v1（i，j+1）-v2（i，j-1））/16/epsl（i，j）； %非均勻介質填充區的拉普拉斯方程差分公式

maxt=abs（v2（i，j）-v1（i，j））；

end

end

end

v1=v2；

end

subplot（1，2，1），mesh（v1） %畫三維曲面圖

axis（[0，31，0，21，0，1]）

subplot（1，2，2），contour（v1，15） %畫等位線

hold on

x=1：1：hx；y=1：1：hy

[xx，yy]= meshgrid（x，y）；

[Gx，Gy]=gradient（v1，0.6，0.6） %計算梯度

quiver（xx，yy，-Gx，-Gy，'black'） %畫電場線圖

圖3 程序流程圖

Fig.3 Program flow chart

為了做到直觀顯示，程序選擇輸出了電位分布三維曲面圖和等位線、電場線分布圖，如圖4所示。

圖4 電位、電場分布圖

Fig.4 Potential， electric field distribution

從結果來看，等位線在介電常數較高的區域相對密集，場強相對較小，因為電介質的極化使得總場有所減弱；而在帶電平板兩端場強極大，其為尖端放電所致，這些均與靜電場基礎理論相吻合。

4結束語

在對有限差分法做簡要介紹的基礎上，詳細推導了非均勻介質填充區域拉普拉斯方程的差分形式。同時，基于Matlab在分析電磁場問題方面的獨特優勢，本文給出了用Matlab編程計算該類問題的程序流程圖及程序設計實例。所做分析有助于加深將Matlab和有限差分法應用于電磁場數值計算的理解和掌握，亦使得相對抽象的靜電場理論具體化且易于理解。

參考文獻：

[1]趙德奎，劉勇.MATLAB 在有限差分法數值計算中的應用[J].四川理工學院學報（自然科學版），2005，18（4）：61-64.

[2]熊彬，阮百堯.MATLAB 在有限差分法中的應用[J].桂林工學院學報，2001，21（2）：104-109.

[3]王潔，陳超波.基于MATLAB的靜態場邊值問題有限差分法的研究[J].微計算機應用，2010，31（3）：1-5.

[4]宋燎原，王平，張海峰，等.靜態電磁場邊值問題計算方法[J].大學物理，2007，26（8）：23-26.

[5] OHTANI T， KANAI Y， COLE J B. A stability improvement technique using PML condition for the three-dimensional nonuniform mesh nonstandard FDTD method[J]. Magnetics， IEEE Transactions on， 2013， 49（5）： 1569-1572.

[6]王為，覃宇建，劉培國，等.基于高階時域有限差分法與改進節點分析法混合求解復雜傳輸線網絡瞬態響應[J].電子與信息學報，2012，34（12）：2999-3005.

[7]朱小敏，任新成，郭立新.指數型粗糙地面與上方矩形截面柱寬帶電磁散射的時域有限差分法研究[J].物理學報， 2014，63（5）：1-7.

[8] Richard C. Booton. Computational Methods for Electromagnetics and Microwaves[M].New York： Wiley， 1992.

[9]孫小東，李振春，王小六.三角網格有限差分法疊前逆時偏移方法研究[J].地球物理學進展，2012，27（5）：2077-2083.

[10]馮慈璋.電磁場（第二版）[M].北京： 高等教育出版社，1983.

[11]何紅雨.電磁場數值計算法與 MATLAB 實現[M].武漢： 華中科技大學出版社，2004.