戴晴琴

(北京師范大學(xué) 物理學(xué)系,北京 100875)

對(duì)于第Ⅱ類超導(dǎo)體,當(dāng)外加磁場(chǎng)大小介于上下臨界場(chǎng)之間時(shí),其正常態(tài)(N)與超導(dǎo)態(tài)(S)的界面能為負(fù)值,從而允許外磁場(chǎng)以磁通渦旋(vortex)的形式穿入超導(dǎo)體內(nèi)部,形成超導(dǎo)與正常態(tài)渦旋共存的混合態(tài)(mixed state, MS)[1]. 此時(shí)如果在垂直于磁場(chǎng)的方向外加電流,則單根磁通在單位長(zhǎng)度上會(huì)受到洛倫茲力f=Φ0J,從而發(fā)生運(yùn)動(dòng)[2].

非理想超導(dǎo)體中總是會(huì)存在一些缺陷,例如晶界、點(diǎn)或面缺陷、雜質(zhì)相等. 為降低整個(gè)體系的自由能,渦旋會(huì)傾向于被釘扎在這些位置[1]. 這一類能夠阻礙渦旋運(yùn)動(dòng)的區(qū)域被稱為釘扎中心,其與渦旋的相互作用被稱為釘扎力.

近年來,利用不對(duì)稱釘扎控制渦旋動(dòng)力學(xué)過程的研究引起了人們極大的興趣,常被用于設(shè)計(jì)磁通泵、整流器或二極管等[3-7]. 渦旋的相互作用及其動(dòng)力學(xué)行為決定了超導(dǎo)體幾乎所有的電、磁學(xué)特性,通過對(duì)渦旋的動(dòng)力學(xué)研究,可以深入了解超導(dǎo)材料的物性和超導(dǎo)機(jī)制,為其工業(yè)化應(yīng)用提供理論支撐.

本文基于此目的開展研究,首先對(duì)渦旋運(yùn)動(dòng)基本模型進(jìn)行簡(jiǎn)單介紹,后加入對(duì)渦渦相互作用力的考慮,利用MATLAB軟件進(jìn)行數(shù)值模擬比較渦旋數(shù)目、洛倫茲力大小、交流電周期對(duì)渦旋運(yùn)動(dòng)的影響,最后加入熱噪聲擾動(dòng)項(xiàng),從而實(shí)現(xiàn)了對(duì)超導(dǎo)體中渦旋動(dòng)力學(xué)過程的簡(jiǎn)單模擬.

1 渦旋運(yùn)動(dòng)基礎(chǔ)模型[8]

1.1 不對(duì)稱周期勢(shì)

考慮如圖1所示幾何形狀的II型超導(dǎo)體薄膜,其釘扎勢(shì)分布滿足U(x,y)=U(x)=U(x+l)即在y方向具有平移不變性,而在x方向呈現(xiàn)以l為周期的變化,包括l1上升沿與l2下降沿. 周期中最大勢(shì)能為ΔU.由于假定l1大l2,且釘扎勢(shì)[8]滿足

圖1 具有不對(duì)稱勢(shì)的超導(dǎo)體示意圖

圖2 不對(duì)稱勢(shì)參數(shù)(x方向)

1.2 洛倫茲力

如圖1所示,給樣品施加y方向交流電,渦旋將受到洛倫茲力[8]：

(2)

1.3 運(yùn)動(dòng)方程與解析解

考慮上述兩種受力,渦旋滿足動(dòng)力學(xué)方程[8]

ηv=fL+fu

(3)

其中η為黏度,fL為洛倫茲力,fu為釘扎力,v為某時(shí)刻受力后產(chǎn)生瞬時(shí)速度.

當(dāng)直流電沿正y方向流動(dòng)時(shí),洛倫茲力使渦流沿正x方向以速度v+移動(dòng);改變電流方向后,渦旋速度反向,記作v-.由于勢(shì)的不對(duì)稱,v+與v-大小并不相等,從而產(chǎn)生沿某一方向的凈速度[8]：

根據(jù)f1、f2與洛倫茲力fL的相對(duì)大小,考慮交流電周期趨于無窮大,可分為3種運(yùn)動(dòng)情況.

當(dāng)fL

v1=0

(5)

當(dāng)f1

(6)

同時(shí)v2-=0,故有[8]：

(7)

當(dāng)f2

(8)

其凈速度[8]為

(9)

利用MATLAB對(duì)以上速度解進(jìn)行繪圖.值得注意的是,本文重點(diǎn)研究各力、速度、位移之間的相對(duì)大小,對(duì)于絕對(duì)值不作嚴(yán)格討論,故全文中所繪圖像均以a0為位移單位,f0為力常數(shù)單位,v0為速度單位,vy0為y方向分速度單位. 在實(shí)際計(jì)算以及生產(chǎn)運(yùn)用中,只需要將各單位代入具體數(shù)值計(jì)算即可.

圖3 v隨fL的變化關(guān)系

渦旋速度在fL=f2時(shí)達(dá)到最大值,若洛倫茲力繼續(xù)增大,則會(huì)出現(xiàn)過驅(qū)動(dòng)現(xiàn)象,即速度隨洛倫茲力增大而減小.

2 渦旋運(yùn)動(dòng)數(shù)值模擬

2.1 不對(duì)稱周期勢(shì)模擬

如圖4所示,通過在樣品上設(shè)置部分重疊的強(qiáng)弱釘扎位點(diǎn),模擬產(chǎn)生不對(duì)稱周期勢(shì).

圖4 單組強(qiáng)弱釘扎示意圖

釘扎中心采用衰減長(zhǎng)度為Rp高斯勢(shì)阱模型,則單個(gè)釘扎的釘扎力[9]可表示為

(10)

其中fp0為釘扎力常數(shù),ri為渦旋位置,R為釘扎中心所在位置,強(qiáng)弱釘扎的不同之處體現(xiàn)在力常數(shù)與釘扎位置上.

如圖5所示,在二維樣品的x與y方向等間距a0設(shè)置一系列強(qiáng)弱釘扎,各組強(qiáng)弱釘扎間距離為d.

圖5 樣品中強(qiáng)弱釘扎位置示意圖

固定x的值,利用MATLAB模擬y方向各位點(diǎn)所受釘扎力fp的受力情況如圖6所示.

圖6 釘扎力隨y方向變化情況

2.2 渦渦相互作用模擬

當(dāng)樣品中的渦旋個(gè)數(shù)大于1時(shí),渦旋不光會(huì)受到釘扎力與洛倫茲力,渦旋之間的排斥相互作用也應(yīng)當(dāng)考慮在內(nèi). 可假設(shè)該力大小與兩渦旋間距離成反比,方向指向兩渦旋連線方向即[9]

(12)

其中ri為第i個(gè)渦旋的位移,rj為第j個(gè)渦旋的位移,fvv0為渦渦相互作用力常數(shù),并設(shè)置截止作用距離為5a0．

2.3 數(shù)值模擬方法

加入對(duì)渦渦相互作用的考慮后,渦旋動(dòng)力學(xué)方程修正為[9]

ηvi=fL+fvv(ri)+fp(ri)

(14)

利用生成隨機(jī)數(shù)的方法獲得渦旋的初始位移,設(shè)置迭代步長(zhǎng)τ0,利用牛頓法的思想進(jìn)行迭代. 首先利用初始位置求得受力,得到初始速度,再利用

x←x+vτ0

(15)

獲得更新后的位置,重復(fù)上述操作,并繪制一個(gè)渦旋相應(yīng)的時(shí)間-y方向位移曲線與時(shí)間-y方向速度曲線.

2.4 單個(gè)渦旋模擬結(jié)果

渦旋運(yùn)動(dòng)情況與洛倫茲力的大小、周期密切相關(guān). 如圖7所示,控制洛倫茲力大小不變,且f1

圖7 單渦旋y方向位移(左)/速度(右)- 時(shí)間曲線

運(yùn)動(dòng)結(jié)果可借助圖8進(jìn)行分析解釋,設(shè)初始時(shí)刻渦旋位置為O,洛倫茲力沿y軸正方向．

圖8 渦旋位置示意圖

當(dāng)周期較小為400τ0時(shí),渦旋首先沿y軸正方向運(yùn)動(dòng),爬升至OA段中間時(shí),洛倫茲力反向,又回到初始位置．

當(dāng)周期增大至500τ0時(shí),渦旋可沿正方向從O點(diǎn)運(yùn)動(dòng)至BC段,后洛倫茲力反向,渦旋會(huì)到B點(diǎn),對(duì)應(yīng)圖中的沿y軸負(fù)方向位移,但由于f1

當(dāng)周期繼續(xù)增大至700τ0時(shí),渦旋可運(yùn)動(dòng)至BC段偏C處段,故對(duì)應(yīng)的位移后退會(huì)更加明顯,且仍然存在平臺(tái)期．

增大洛倫茲力使得f1

圖9 增大洛倫茲力后的渦旋運(yùn)動(dòng)曲線

該模擬結(jié)果對(duì)應(yīng)模型中的過驅(qū)動(dòng)情況,即正負(fù)方向均有速度,但由于v+>|v-|,故仍存在y方向的凈位移．

2.5 渦旋的速度-周期曲線

為進(jìn)一步探究在固定洛倫茲力下,渦旋運(yùn)動(dòng)速度隨周期的變化情況,可引入一渦旋平均速度Vdc,其定義為各渦旋在一個(gè)時(shí)間周期內(nèi)的速度平均即[9]

(18)

控制其他參數(shù)不變,且維持f1

圖10 渦旋的速度-周期曲線

3 熱噪聲下的渦旋運(yùn)動(dòng)

為更好模擬渦旋運(yùn)動(dòng)的真實(shí)情況,對(duì)熱噪聲的考慮是不可避免的,可看作許多比渦旋更輕的粒子對(duì)其碰撞產(chǎn)生了熱噪聲力[10],運(yùn)動(dòng)方程因此變?yōu)閇9]

ηvi=fL+fvv(ri)+fp(ri)+fT(t)

(21)

利用拆分法[11](Splitting method)對(duì)上述朗之萬方程進(jìn)行模擬,將原方程拆分為ABO三項(xiàng)．

B項(xiàng)：

v←v+F/2

(22)

其中F為除熱噪聲外的合力,處理思路與2.3節(jié)類似,通過(12)式來計(jì)算某處受力產(chǎn)生的速度．

A項(xiàng)：

x←x+vτ0/2

(23)

由于此處使用拆分法,故更新位置時(shí)僅計(jì)算半個(gè)迭代步長(zhǎng)對(duì)應(yīng)的位移．

而關(guān)于噪聲項(xiàng)有

(24)

由Ornstein-Uhlenbeck過程表示噪聲擾動(dòng)項(xiàng),其中α為均值回歸速率,回歸均值設(shè)為0,β為波動(dòng)率,dW表示遵從高斯分布的布朗運(yùn)動(dòng)．

解該方程即可得到O項(xiàng)：

(24)

其中R為均值為0,方差為1的高斯隨機(jī)數(shù).

完成上述BAOA過程對(duì)應(yīng)一個(gè)完整的迭代周期,將加入熱噪聲的單個(gè)渦旋周期-速度曲線與為考慮熱噪聲時(shí)的圖像對(duì)比,如圖11所示.

圖11 單渦旋加入噪聲前(左)后(右)的速度-周期曲線

當(dāng)噪聲較小時(shí),圖像趨勢(shì)不變,僅光滑程度受到影響.

4 結(jié)論

本文首先對(duì)渦旋在超導(dǎo)體中的受力進(jìn)行了簡(jiǎn)單分析,包括由不對(duì)稱周期勢(shì)產(chǎn)生的釘扎力、隨時(shí)間周期性變化的洛倫茲力,并由此列出渦旋的運(yùn)動(dòng)學(xué)方程. 進(jìn)一步就各力的相對(duì)大小關(guān)系進(jìn)行了討論,可分為三種情況：當(dāng)洛倫茲力小于釘扎力時(shí),渦旋被困在釘扎中心無法運(yùn)動(dòng);當(dāng)洛倫茲力大于某一方向的釘扎力時(shí),渦旋可受洛倫茲力驅(qū)動(dòng)克服釘扎力進(jìn)行單方向運(yùn)動(dòng);當(dāng)洛倫茲力大于兩個(gè)方向的釘扎力時(shí),渦旋將處于過驅(qū)動(dòng)狀態(tài),向具有凈速度的一方移動(dòng).由此求得了渦旋運(yùn)動(dòng)方程的解析解,并用MATLAB完成圖像繪制,發(fā)現(xiàn)在洛倫茲力大小與較大釘扎力相等時(shí)取得渦旋運(yùn)動(dòng)的最大速度.

其次在上述模型的基礎(chǔ)上,加入對(duì)渦旋間相互作用力的考慮,利用MATLAB進(jìn)行數(shù)值模擬,繪制了渦旋運(yùn)動(dòng)的位置-時(shí)間圖以及速度-時(shí)間圖,并探究了洛倫茲力大小與交流電周期對(duì)其產(chǎn)生的影響,結(jié)合基礎(chǔ)模型對(duì)其分析解釋. 進(jìn)一步地,繪制了不同渦數(shù)所對(duì)應(yīng)的平均速度-交流電周期圖,在T→∞與解析解相對(duì)應(yīng)．

最后用Ornstein-Uhlenbeck過程表示噪聲擾動(dòng)項(xiàng),利用Splitting method對(duì)朗之萬方程進(jìn)行數(shù)值模擬,繪制了噪聲擾動(dòng)下的速度-周期曲線,并與未考慮噪聲時(shí)的曲線進(jìn)行對(duì)比．

電流驅(qū)動(dòng)渦流的動(dòng)力學(xué)研究對(duì)于理解強(qiáng)相互作用渦流的基本集體行為以及在超導(dǎo)體中獲得高非耗散電流具有重要意義[12],且在許多應(yīng)用中發(fā)揮著關(guān)鍵作用,例如高場(chǎng)磁體[13]、超導(dǎo)數(shù)字存儲(chǔ)器和量子位[14]、太赫茲輻射源[15]或用于粒子加速器的諧振腔[16]等．本文基于前人研究,對(duì)超導(dǎo)體中渦旋運(yùn)動(dòng)進(jìn)行了簡(jiǎn)單分析與模擬,并在模型中加入對(duì)熱噪聲干擾項(xiàng)的考慮,使其更貼近真實(shí)應(yīng)用情境,具有一定的現(xiàn)實(shí)意義．

致謝：

感謝北京師范大學(xué)劉海文教授對(duì)本文提出的修改意見。