洪德成，王海軍

(吉林大學(xué) 物理學(xué)院，吉林 長春 130012)

隨著計算機模擬技術(shù)的發(fā)展，在科學(xué)研究領(lǐng)域，計算物理學(xué)與理論物理學(xué)、實驗物理學(xué)成為當(dāng)代物理學(xué)發(fā)展的三駕馬車.但在目前的大學(xué)物理教育體系中仍以理論物理為主，實驗物理次之，計算物理再次之.理論物理和實驗物理都已建立起成熟而完備的教學(xué)體系，而計算物理起步晚，發(fā)展緩慢，教學(xué)案例較為有限[1-3].

本文針對磁偶極子源，從麥克斯韋方程組出發(fā)，通過引入赫茲勢推導(dǎo)無限大均勻介質(zhì)中電磁場代數(shù)表達(dá)式和積分表達(dá)式，并給出相應(yīng)的數(shù)值計算方法.代數(shù)解形式簡單、參數(shù)關(guān)系明確，但不能處理含有邊界的非均勻介質(zhì)電磁場問題；積分解適用性較廣，可進(jìn)一步擴展到水平層狀分層、圓柱狀分層及球狀分層介質(zhì)中電磁場分布問題的求解[4-9].根據(jù)不同的邊界類型(水平邊界、圓柱邊界、球形邊界)，積分解的具體表達(dá)形式也不同.

為簡明起見，本文給出適用于水平分層介質(zhì)的電磁場積分解.積分表達(dá)式為包含貝塞爾函數(shù)的高階震蕩無窮積分(也稱Sommerfeld積分).采用三次樣條擬合、高斯-勒讓德數(shù)值積分[7]、以及兩種算法相結(jié)合的方式分別求解Sommerfeld積分.

1 磁偶極子源的電場和磁場

1.1 場的代數(shù)解

無限小電流環(huán)可等價為磁偶極子源，設(shè)電流隨時間變化為e-iωt，麥克斯韋方程組滿足如下形式：

(1)

(2)

介質(zhì)的本構(gòu)關(guān)系：

B=μH+μMδ(r)，D=εE，J=σE

(3)

其中μ、σ和ε分為介質(zhì)磁導(dǎo)率、電導(dǎo)率和介電常量.M代表磁偶極子源強度.將式(3)代入式(1)、式(2)得

(4)

(5)

其中σ*=σ-iωε為復(fù)電導(dǎo)率.對式(5)取散度，并考慮到·(×A)=0，得

(6)

由此引入赫茲勢π[8]：

E=iωμ×π

(7)

同時，令赫茲勢π滿足如下散度方程：

(8)

將式(7)代入式(5)得磁場表達(dá)式：

H=iωμσ*π+φ

(9)

將式(7)、式(8)代入式(4)，得到如下赫茲勢方程：

(10)

(11)

式(10)進(jìn)一步展開為

(12)

可以看出，由于空間對稱關(guān)系，每個分量滿足形式相同的偏微分方程.這里我們只考慮z方向磁偶極子源Mz.即令Mx=My=0，則式(12)的解為[8]

(13)

(14)

將式(13)代入式(7)，得到電場分量表達(dá)式：

Ez=0

(15)

將式(13)和式(14)代入式(9)，得到磁場分量表達(dá)式：

(16)

1.2 場的積分表達(dá)式

在圓柱坐標(biāo)系(ρ，φ，z)中，有如下恒等式[8]：

(17)

(18)

進(jìn)一步將式(18)代入式(8)得標(biāo)勢φ積分表達(dá)式：

(19)

將上面二式分別代入式(7)和式(9)，并應(yīng)用柱坐標(biāo)系下偏微分算子得電場分量：

Eρ=0，

Ez=0

(20)

及磁場分量：

Hφ=0，

(21)

(22)

2 數(shù)值積分方法

式(20)和式(21)中的Sommerfeld積分可統(tǒng)一寫為如下形式：

(23)

本文分別采用三次樣條插值、高斯-勒讓德積分及這2種方法結(jié)合的方式求解式(23)的數(shù)值積分.

2.1 三次樣條擬合積分方法

(24)

其中，多項式系數(shù)由下式給出：

(25)

式中hj=kρ，j+1-kρ，j為離散點間隔，參量Mj由下列矩陣求解：

圖1 被積函數(shù)隨自變量變化關(guān)系

(26)

方程兩側(cè)矩陣元素分別為：

d1=d2，dN+1=dN

(27)

在數(shù)值計算中，通常將式(23)的無窮積分截斷為有限積分.取kρ，N+1為積分上限，kρ，N+1以后對應(yīng)的被積函數(shù)足夠小，在滿足積分精度要求下可以忽略不計.由式(25)—式(27)確定所有區(qū)間內(nèi)的三次函數(shù)系數(shù)后，每一區(qū)間積分結(jié)果為

(28)

將所有區(qū)間積分結(jié)果相加即為式(23)的積分結(jié)果.

2.2 高斯-勒讓德積分方法

高斯-勒讓德積分是一種將積分運算轉(zhuǎn)換為乘積運算的成熟數(shù)學(xué)積分方法，其運算規(guī)則如下：

(29)

其中，方程左側(cè)g(t)是任意形式的被積函數(shù)，積分限為[-1，1]；方程右側(cè)為采樣點對應(yīng)函數(shù)值與加權(quán)系數(shù)的乘積的和，其中tl為自變量分布在[-1，1]區(qū)間的固定采樣點，Wl為對應(yīng)點的加權(quán)系數(shù).

例如：當(dāng)L=5時，有

tl=±0.906 179 845 9，Wl=0.236 926 885 1

tl=±0.538 469 310 1，Wl=0.478 628 670 5

tl=0，Wl=0.568 888 888 9

一般來講，階數(shù)L越大數(shù)值結(jié)果精度越高，但耗時也越多.

在求解式(23)積分時，首先將被積函數(shù)劃分為若干區(qū)間，然后對每一區(qū)間積分上、下限做變量替換.令

(30)

則每一區(qū)間的積分變形為

(31)

與式(29)對照，應(yīng)用高斯-勒讓德積分方法即可求出該區(qū)間的積分結(jié)果，然后將所有區(qū)間結(jié)果相加即為式(23)的積分結(jié)果.

2.3 三次樣條擬合結(jié)合高斯-勒讓德積分方法

三次樣條擬合方法和高斯-勒讓德積分方法都可以獨立求解式(23)的數(shù)值積分問題，為增加對積分方法的深入理解與靈活應(yīng)用，也可將兩種方法結(jié)合起來.采用三次樣條函數(shù)擬合式(23)中被積函數(shù)隨自變量變化相對平滑的部分：

(32)

該區(qū)間內(nèi)的函數(shù)積分可進(jìn)一步表達(dá)為

(33)

式(33)右側(cè)積分可采用上述的高斯-勒讓德積分方法求出.

上述3種積分方法是嚴(yán)格完備的，理論上，計算結(jié)果也可以無限逼近真實值.但計算機編程中總會存在由多種因素引起的計算誤差，其中被積分函數(shù)形態(tài)及采樣點分布規(guī)則是常見的主要因素.對同一積分問題采用不同的積分方法，一方面可以相互驗證數(shù)值積分結(jié)果的正確性，另一方面也可通過對計算精度和計算效率的評估優(yōu)選出最適合該問題的最佳積分方法.

3 數(shù)值結(jié)果

本小節(jié)針對一個給定的物理模型，給出電磁場分量的數(shù)值計算結(jié)果.假定介質(zhì)磁導(dǎo)率為真空磁導(dǎo)率μ=μ0=4π×10-7H/m;介電常量為真空介電常量ε=ε0=8.854×10-12F/m，電導(dǎo)率σ=0.1 S/m.磁偶極子強度Mz=1 A·m2，頻率f=20 kHz.發(fā)射點位置(0，0，0)，場點位置(0.5，0，1).表1給出了代數(shù)解與積分解對比結(jié)果.可以看出3種積分?jǐn)?shù)值結(jié)果與代數(shù)解吻合，證明了理論方法的正確性.數(shù)值結(jié)果有效數(shù)字后幾位的計算誤差主要來自于對積分上下限的截斷及積分區(qū)間分割的位置、數(shù)量多少等.

4 結(jié)論

本文從麥克斯韋方程出發(fā)推導(dǎo)了磁偶極子源在均勻介質(zhì)中激發(fā)的電磁場代數(shù)解與積分解.并采用三次樣條擬合方法、高斯-勒讓德積分方法及它們的結(jié)合方式分別求解Sommerfeld數(shù)值積分.為學(xué)生提供了一個以電動力學(xué)為背景的計算物理教學(xué)案例.

表1 電磁場代數(shù)解與積分解數(shù)值對比