夏舸 楊立 寇蔚 杜永成

(海軍工程大學(xué)動(dòng)力工程學(xué)院,武漢 430033)

基于變換熱力學(xué)的三維任意形狀熱斗篷設(shè)計(jì)?

夏舸 楊立?寇蔚 杜永成

(海軍工程大學(xué)動(dòng)力工程學(xué)院,武漢 430033)

(2016年11月22日收到;2017年3月9日收到修改稿)

在變換熱力學(xué)的基礎(chǔ)上,通過(guò)坐標(biāo)變換的方法,推導(dǎo)出三維任意形狀熱斗篷導(dǎo)熱系數(shù)的通解表達(dá)式,并進(jìn)行了全波仿真驗(yàn)證.結(jié)果表明∶熱流均能繞過(guò)保護(hù)區(qū)域流出,保護(hù)區(qū)域的溫度保持不變,而且熱斗篷外的溫度場(chǎng)并沒(méi)有破壞,具有很好的熱保護(hù)和熱隱身的效果.這一方法把變換熱力學(xué)從二維拓展到三維,具有普遍的適用性.同時(shí),這種技術(shù)為熱流流動(dòng)路徑和目標(biāo)溫度場(chǎng)的控制奠定了理論基礎(chǔ),在微芯片、電動(dòng)機(jī)的保護(hù)以及目標(biāo)熱隱身上有潛在應(yīng)用.

∶變換熱力學(xué),熱斗篷,三維任意形狀,溫度場(chǎng)

PACS∶44.10.+I,05.70.—a,81.05.Xj,07.05.TpDOI∶10.7498/aps.66.104401

1 引 言

2006年,Pendry等[1]和Leonhardt[2]從麥克斯韋方程的形式不變性出發(fā),通過(guò)坐標(biāo)變換的數(shù)學(xué)方法推導(dǎo)出具有各向異性特性的人工超材料的參數(shù),使得電磁波按照預(yù)先設(shè)定的路徑傳播,從而達(dá)到電磁隱形的目的.從此,這種控制電磁波傳播的理論方法引起了廣大學(xué)者極大的關(guān)注,并逐漸發(fā)展成為變換光學(xué)的一般理論[3?5].變換光學(xué)理論不但可以應(yīng)用于制作隱形器件[6,7],還能用來(lái)設(shè)計(jì)波束集中器[8,9]、完美透鏡[10]、波束旋轉(zhuǎn)器[11]以及指向性天線等[12].后來(lái),人們又把光學(xué)變換理論的思想應(yīng)用到聲波[13,14]、物質(zhì)波[15]、彈性波[16]和熱流[17]的傳播路徑設(shè)計(jì)上.

與電磁波、聲波以及彈性波等不同,熱流并不是一種波,僅僅是從熱端到冷端的流動(dòng)[18].基于導(dǎo)熱方程的不變性,通過(guò)構(gòu)造出特定的熱物性參數(shù),可以實(shí)現(xiàn)對(duì)熱流傳遞方向的控制與利用,滿足工程實(shí)際應(yīng)用.2008年,Fan等[19]通過(guò)利用分級(jí)材料,設(shè)計(jì)出球形和橢球型熱隱身衣,并對(duì)其進(jìn)行了仿真驗(yàn)證,這是變換熱力學(xué)理論第一次被提出并驗(yàn)證.后來(lái),Guenneau等[20]提出一種通過(guò)20層同心的各向同性材料構(gòu)建出二維非穩(wěn)態(tài)熱斗篷和熱聚集器的方法.在此基礎(chǔ)上,Schittny等[21]利用銅片和聚二甲基硅氧烷兩種材料交替排列構(gòu)建10個(gè)等厚度的環(huán),從實(shí)驗(yàn)上驗(yàn)證了二維熱斗篷的隱身效果.2014年,毛福春等[22]在變換熱力學(xué)的基礎(chǔ)上,設(shè)計(jì)出具有非共形任意橫截面的柱形熱斗篷熱導(dǎo)率的表達(dá)式,并通過(guò)COMSOL Multiphysics仿真驗(yàn)證了熱斗篷對(duì)隱身區(qū)域具有熱保護(hù)的功能.同時(shí),秦春雷等[23]將連續(xù)介質(zhì)力學(xué)的方法引入熱力學(xué)中,數(shù)值求解連續(xù)介質(zhì)位移的拉普拉斯方程獲得封閉曲線內(nèi)場(chǎng)的形變信息,進(jìn)而可以計(jì)算任意形狀的熱斗篷的材料參數(shù).

熱斗篷由于其特殊的導(dǎo)熱系數(shù)分布使得熱流繞過(guò)保護(hù)區(qū)域,這一特質(zhì)使其能夠應(yīng)用于計(jì)算機(jī)芯片、高速電動(dòng)機(jī)的熱保護(hù)以及軍事目標(biāo)的偽裝中.然而,目前對(duì)熱斗篷導(dǎo)熱系數(shù)研究大多數(shù)局限在二維空間里,而對(duì)于三維形狀熱斗篷也基本局限于球體、橢球體以及柱體等規(guī)則共形熱斗篷的研究.而在實(shí)際應(yīng)用過(guò)程中,這些簡(jiǎn)單規(guī)則的三維模型往往不能滿足要求.而對(duì)基于拉普拉斯方程的熱斗篷研究方法由于其介質(zhì)參數(shù)全部是數(shù)值結(jié)果,在精確提取參數(shù)進(jìn)行材料實(shí)現(xiàn)時(shí)困難較大,不利于在工程上實(shí)現(xiàn).因此,本文基于變換熱力學(xué)的方法,推導(dǎo)出三維任意形狀熱斗篷的導(dǎo)熱系數(shù)的通用表達(dá)式,并通過(guò)全波仿真分析球體、橢球體、共形任意形狀熱斗篷和非共形任意形狀熱斗篷,驗(yàn)證了表達(dá)式的正確性以及熱斗篷的隱身和保護(hù)雙重功能,這些結(jié)果為三維熱斗篷熱物性參數(shù)的隱身設(shè)計(jì)和改進(jìn)奠定了理論基礎(chǔ).

2 理論和方法

熱傳導(dǎo)是物體只依靠微觀粒子(分子、原子或者自由電子)的熱運(yùn)動(dòng)產(chǎn)生的熱量傳遞的現(xiàn)象.在無(wú)內(nèi)熱源情況下,穩(wěn)態(tài)熱傳導(dǎo)方程的一般形式可以寫(xiě)成[24]

式中,T表示區(qū)域溫度,λ表示導(dǎo)熱系數(shù),?表示梯度算符.

與基于麥克斯韋方程的變換光學(xué)類(lèi)似,根據(jù)變換熱力學(xué)理論[20],由于方程(1)具有形式不變性,導(dǎo)熱方程在變換空間可表示為

式中,T′和λ′分別表示物體在變換空間中的溫度和導(dǎo)熱系數(shù),?′表示變換空間中的梯度算符.變換空間與原空間介質(zhì)的導(dǎo)熱系數(shù)的關(guān)系為[25]

式中,A為變換空間與原空間之間的雅克比矩陣,ATAT是A的轉(zhuǎn)置矩陣,det(A)為矩陣的行列式的值.不同坐標(biāo)系下的雅克比矩陣是不同的,在球坐標(biāo)系下的雅克比矩陣∶

對(duì)于三維熱斗篷,假設(shè)內(nèi)外邊界分別為r=R1(θ,φ)和r=R2(θ,φ), 為了實(shí)現(xiàn)對(duì)保護(hù)區(qū)域的完美隱形,把0＜ r＜ R2(θ,φ)區(qū)域壓縮到R1(θ′,φ′) ＜ r ＜ R2(θ′,φ′)區(qū)域, 得到保護(hù)區(qū)域r ＜ R1(θ′,φ′). 區(qū)域r ＞ R2(θ,φ) 保持不變. 因此變換方程如下[22]∶

把(5)式代入(4)和(3)式,經(jīng)過(guò)推導(dǎo)可以得到所需變換空間的導(dǎo)熱系數(shù)為

式中

將球坐標(biāo)系下的材料導(dǎo)熱系數(shù)轉(zhuǎn)換為材料在笛卡爾坐標(biāo)系下的導(dǎo)熱系數(shù),具體轉(zhuǎn)換關(guān)系如下∶

化簡(jiǎn)得到

根據(jù)(6)和(8)式即可得到內(nèi)外邊界為r=R1(θ′,φ′)和r=R2(θ′,φ′)的三維任意形狀熱斗篷導(dǎo)熱系數(shù)的一般表達(dá)式.當(dāng)消去θ′,簡(jiǎn)化為二維任意形狀熱斗篷時(shí),與參考文獻(xiàn)[22]一致.

3 仿真驗(yàn)證與分析

為了驗(yàn)證熱斗篷導(dǎo)熱系數(shù)表達(dá)式的正確性,本文采用有限元仿真軟件COMSOL Multiphysics進(jìn)行全波仿真.假設(shè)整個(gè)區(qū)域?yàn)橐粋€(gè)1 m×1 m×1 m的空間,前后左右四個(gè)面保持絕熱,上下兩個(gè)面分別保持恒溫,上面為高溫保持在400 K,下面為低溫保持在300 K;在被保護(hù)的隱形區(qū)域和熱斗篷外的區(qū)域的材料為銅.在常溫下銅的各項(xiàng)基本參數(shù)如下∶密度ρ=8960 kg/m3,導(dǎo)熱系數(shù)λ=400 W/(m·K),常壓熱容Cp=385 J/(kg·K).在熱斗篷區(qū)域,材料的導(dǎo)熱系數(shù)λ′由(6)式和(8)式計(jì)算得出.

3.1 共形熱斗篷

熱斗篷區(qū)域的形狀主要取決于內(nèi)外邊界函數(shù)的選取, 當(dāng)外邊界函數(shù)R2(θ′,φ′)與內(nèi)邊界R1(θ′,φ′)的比值為常數(shù), 且r=R1(θ′,φ′)全部區(qū)域都在r=R2(θ′,φ′) 內(nèi)部時(shí), 則此為共形熱斗篷.先取典型案例球體熱斗篷作為分析對(duì)象,假設(shè)內(nèi)球半徑R1(θ′,φ′)=0.2 m, 外球半徑為R2(θ′,φ′)=2R1(θ′,φ′)=0.4 m.

圖1(a)—(c)為穩(wěn)態(tài)球體熱斗篷的模擬結(jié)果.從圖1(a)和圖1(b)可以看出,未加熱斗篷時(shí),所有等溫面保持平行,與熱流方向垂直.在加入熱斗篷后,進(jìn)入熱斗篷之前所有的等溫面保持平行,與未加入斗篷的情況一樣.在進(jìn)入熱斗篷區(qū)域后,等溫面開(kāi)始發(fā)生變化從而影響整個(gè)溫度場(chǎng)的分布,當(dāng)?shù)葴孛孢M(jìn)入上半球殼熱斗篷時(shí),由于球殼內(nèi)不同位置導(dǎo)熱系數(shù)不同,越接近內(nèi)球,徑向?qū)嵯禂?shù)越小,所以等溫面像一頂帽子扣在內(nèi)球之上,而下半球殼熱斗篷情況與上半球殼正好對(duì)稱(chēng).內(nèi)球并沒(méi)有任何熱流穿過(guò),說(shuō)明內(nèi)球保護(hù)區(qū)域的溫度保持一致.當(dāng)?shù)葴孛娉隽藷岫放窈?又恢復(fù)了平行狀態(tài).圖1(c)為yoz平面上的熱流與等溫線分布圖.熱流從上到下均勻流動(dòng),在斗篷外側(cè),等溫線間距大小一致,熱流密度處處相等.熱流能夠平滑地繞過(guò)內(nèi)部的保護(hù)區(qū)域而對(duì)斗篷外部區(qū)域不產(chǎn)生任何影響,這與文獻(xiàn)[26]結(jié)論一致.由于內(nèi)球保護(hù)區(qū)域目標(biāo)的存在不會(huì)破壞斗篷外部的溫度場(chǎng),同時(shí)外部熱流無(wú)法進(jìn)入目標(biāo)的區(qū)域,從上表面方向上無(wú)法發(fā)現(xiàn)目標(biāo)的存在.因此,此球體熱斗篷具有較好的熱保護(hù)和熱隱身效果.

圖1 (網(wǎng)刊彩色)球體熱斗篷(a),(b)等溫面分布圖;(c)yoz平面的等溫線和熱通量分布Fig.1.(color online)Spherical thermal cloak:(a),(b)Distribution of isothermal surface;(c)isotherms and heatflux profile in the yoz plane.

圖2和圖3為三維橢球和共形三維任意形狀熱斗篷的仿真結(jié)果.橢球的球坐標(biāo)方程為

共形三維任意形狀熱斗篷的仿真參數(shù)方程為

仿真時(shí),R1(θ′,φ′)= R(θ′,φ′),R2(θ′,φ′)=1.5R1(θ′,φ′), 由圖3 的仿真結(jié)果可以看出, 該熱斗篷同樣具有熱保護(hù)和熱隱身雙重功能.

圖2 (網(wǎng)刊彩色)橢球熱斗篷 (a)三維空間溫度分布;(b)yoz平面的等溫線和熱通量分布Fig.2. (color online)Ellipsoidal thermal cloak:(a)Three-dimensional temperature distribution;(b)isotherms and heatflux profile in the yoz plane.

圖3 (網(wǎng)刊彩色)共形三維任意形狀熱斗篷 (a)三維空間溫度分布;(b)yoz平面的等溫線和熱通量分布Fig.3. (color online)Three-dimensional conformal thermal cloak with arbitrary shapes:(a)Threedimensional temperature distribution;(b)isotherms and heatflux profile in the yoz plane.

3.2 非共形熱斗篷

為了驗(yàn)證方便,本文首先選取外邊界函數(shù)為球體,內(nèi)邊界保護(hù)區(qū)域?yàn)闄E球體.仿真時(shí),球體半徑R2(θ′,φ′)=0.3 m, 橢圓的具體參數(shù)a=0.15 m,b=0.1 m,c=0.05 m.具體仿真結(jié)果如圖4所示.

圖4(a)和圖4(b)是球加橢球非共形熱斗篷的整個(gè)三維空間溫度分布以及yoz平面的熱擴(kuò)散圖.球與橢球是實(shí)際工程中常見(jiàn)的幾何體,由圖可以看出,當(dāng)熱流經(jīng)過(guò)斗篷時(shí)發(fā)生彎曲,從而繞過(guò)了內(nèi)橢球,在斗篷外恢復(fù)原傳遞方向.因此,該斗篷也同樣具有熱保護(hù)和熱隱身雙重功能.

圖5給出了非共形任意三維形狀熱斗篷的仿真結(jié)果. 當(dāng)外邊界函數(shù)R2(θ′,φ′)與內(nèi)邊界函數(shù)R1(θ′,φ′)比值不為常數(shù),且內(nèi)邊界所圍成的區(qū)域全部位于外邊界之內(nèi)時(shí),則此熱斗篷為非共形熱斗篷.仿真時(shí)將內(nèi)外邊界參數(shù)設(shè)置如下∶

由圖5可以看出,內(nèi)部保護(hù)區(qū)域沒(méi)有熱流流入,溫度保持不變,與外部熱流強(qiáng)度大小無(wú)關(guān).

圖4 (網(wǎng)刊彩色)球加橢球非共形熱斗篷 (a)三維空間溫度分布;(b)yoz平面的等溫線和熱通量分布Fig.4. (color online)Non-conformal thermal cloak with the sphere outside the ellipsoid:(a)Threedimensional temperature distribution;(b)isotherms and heatflux profile in the yoz plane.

圖5 (網(wǎng)刊彩色)非共形三維任意形狀熱斗篷 (a)三維空間溫度分布;(b)yoz平面的等溫線和熱通量分布Fig.5.(color online)Three-dimensional non-conformal thermal cloak with arbitrary shapes:(a)Threedimensional temperature distribution;(b)isotherms and heatflux profile in the yoz plane.

4 結(jié) 論

根據(jù)變換熱力學(xué)的基礎(chǔ)理論,推導(dǎo)出了三維任意形狀熱斗篷的導(dǎo)熱系數(shù)通用表達(dá)式,在此基礎(chǔ)上設(shè)計(jì)了球體、橢球體、共形三維任意形狀熱斗篷、球加橢球以及非共形三維任意形狀熱斗篷,并通過(guò)有限元軟件Comsol Multiphysics進(jìn)行全波仿真,結(jié)果驗(yàn)證了本文所推導(dǎo)的三維任意形狀熱斗篷導(dǎo)熱系數(shù)的正確性.在這五個(gè)算例中,熱流均能在熱斗篷的引導(dǎo)下平滑地繞過(guò)保護(hù)區(qū)域,使得保護(hù)區(qū)域不受外部熱流的影響,從而起到了保護(hù)作用;同時(shí),在熱斗篷之外的區(qū)域與原先的溫度場(chǎng)一樣,并沒(méi)有被保護(hù)區(qū)域的目標(biāo)破壞,保持“隱形”的狀態(tài).

相比于先前對(duì)二維任意形狀熱斗篷以及球體、橢球體等規(guī)則三維熱斗篷的研究,本文將變換熱力學(xué)應(yīng)用到了三維任意形狀的熱斗篷中去,具有更好的普適性,更加貼近工程實(shí)際應(yīng)用,同時(shí)也為目標(biāo)的熱隱身和熱保護(hù)提供了一種有效途徑.

[1]Pendry J B,Schurig D,Smith D R 2006 Science 312 1780

[2]Leonhardt U 2006 Science 312 1777

[3]Pendry J B,Schurig D,Smith D R 2006 Opt.Express 14 9794

[4]Pendry J B,Schurig D,Smith D R 2007 Opt.Express 15 14772

[5]Liu Y,Zentgraf T,Bartal G,Zhang X 2010 Nano Lett.10 1991

[6]Schurig D,Mock J J,Justice B J,Cummer S A,Pendry J B,Starr A F,Smith D R 2006 Science 314 977

[7]Cummer S A,Popa B,Schurig D,Smith D R,Pendry J B 2006 Phys. Rev. E 74 036621

[8]Rahm M, Schurig D, Roberts D A, CummerSA,SmithDR,PendryJB2008Photon. Nanostruct. Fundam. Appl.6 87

[9]Aubry A,Lei D Y,Fernándezdomínguez A I,Sonnefraud Y,Maier S A,Pendry J B 2010 Nano Lett.10 2574

[10]Xu H X,Wang G M,Qi M Q,Li L,Cui T J 2013 Adv.Opt.Mater.1 495

[11]Chen H,Chan C T 2007 Appl.Phys.Lett.90 241105

[12]Enoch S,Tayeb G,Sabouroux P,Guérin N,Vincent P 2002 Phys.Rev.Lett.89 213902

[13]Chen Y,Yang F,Xu J Y,Liu X J 2008 Appl.Phys.Lett.92 151913

[14]Wei Q,Chen Y,Liu X J 2012 Appl.Phys.A 109 913

[15]Zhang S,Genov D A,Sun C,Zhang X 2008 Phys.Rev.Lett.100 123002

[16]Farhat M,Guenneau S,Enoch S 2009 Phys.Rev.Lett.103 024301

[17]Hu R,Wei X L,Hu J Y,Luo X B 2014 Sci.Rep.4 3600

[18]Li T H,Zhu D L,MaofC,Huang M,Yang J J,Li S B 2016 Front.Phys.11 1

[19]Fan C Z,Gao Y,Huang J P 2008 Appl.Phys.Lett.92 251907

[20]Guenneau S,Amra C,Veynante D 2012 Opt.Express 20 8207

[21]Schittny R,Kadic M,Guenneau S,Wegener M 2013 Phys.Rev.Lett.110 195901

[22]MaofC,Li T H,Huang M,Yang J J,Chen J C 2014 Acta Phys.Sin.63 014401(in Chinese)[毛春福,李廷華,黃銘,楊晶晶,陳俊昌2014物理學(xué)報(bào)63 014401]

[23]Qin C L,Yang J J,Huang M 2014 Acta Phys.Sin.63 194402(in Chinese)[秦春雷,楊晶晶,黃銘2014物理學(xué)報(bào)63 194402]

[24]Yang S M,Tao W Q 2006 Heat Transfer(4th Ed.)(Beijing:Higher Education Press)p43(in Chinese)[楊世銘,陶文銓 2006傳熱學(xué)(第四版)(北京:高等教育出版社)第43頁(yè)]

[25]Yang T Z,Huang L J,Chen F,Xu W K 2013 J. Phys. D: Appl. Phys.46 305102

[26]Chen T,Weng C N,Tsai Y L 2015 J.Appl.Phys.117 054904

[27]Wu Q,Zhang K,MengfY,Li L W 2010 Acta Phys.Sin.59 6071(in Chinese)[吳群,張狂,孟繁義,李樂(lè)偉2010物理學(xué)報(bào)59 6071]

PACS∶44.10.+I,05.70.—a,81.05.Xj,07.05.TpDOI∶10.7498/aps.66.104401

*Project supported by the National Natural Science Foundation of China(Grant No.11504426)and the National Defense Foundation of China(Grant No.1010502020202).

?Corresponding author.E-mail:yangli123123@126.com

Design and research of three-dimensional thermal cloak with arbitrary shape based on the transformation thermodynamics?

Xia GeYang Li?Kou WeiDu Yong-Cheng
(College of Power Engineering,Naval University of Engineering,Wuhan 430033,China)

22 November 2016;revised manuscript

9 March 2017)

Based on the form-invariance of the thermal conduction equation different from wave equation,transformation thermodynamics has opened up a new area for the arbitrarily manipulating of heatfluxes at discretion by using thermal metamaterials.Moreover,it can help researchers to design different kinds of thermal devices with many unique properties that cannot be simply realized by natural materials,such as thermal cloaking,thermal concentrating,thermal rotating and thermal illusion.Among these devices,the conventional thermal cloak enabling heatfluxes to travel around the inner region,has attracted the most significant attention so far.At the present time,the studies of the thermal cloak mainly focus on two-dimensional space with arbitrary shape and three-dimensional space with regular shape,which appear to be far from enough to meet the engineering requirements.In this paper,we derive the general expression of the thermal conductivity for three-dimensional thermal cloak with arbitrary shape according to the transformation thermodynamics.In this paper,the thermal conductivity in the polar coordinate system is transformed into that in the Cartesian coordinate system by means of coordinate transformation.On the basis of the expression of the thermal conductivity,we adopt full-wave simulation by using the software COMSOL Multiphysics to analyze the cloaking performances offive designed thermal cloaks,i.e.,spherical thermal cloak,ellipsoidal thermal cloak,three-dimensional conformal thermal cloak with arbitrary shapes,non-conformal thermal cloak with the sphere outside the ellipsoid,and three-dimensional non-conformal thermal cloak with arbitrary shapes.The results show that the heatfluxes travel around the protection area,and eventually return to their original paths.The temperature profile inside the thermal cloak keeps unchanged,and the temperature field outside the thermal cloak is not distorted,which proves that the cloak has a perfect thermal invisible effect.Both the conformal and non-conformal thermal cloak have perfect thermal protection and invisible function.In this paper,the transformation thermodynamics is extended from two-dimensional thermal cloak to three-dimensional thermal cloak with better universality.At the same time,this technology provides morefl exibility in controlling heat flow and target temperature field,which will have potential applications in designing microchip,motor protection and target thermal stealth.It is believed that the method presented here can be applied to other branches of physics,such as acoustics,matter waves and elastic waves.

∶transformation thermodynamics,thermal cloak,three-dimensional arbitrary shape,temperature field

?國(guó)家自然科學(xué)基金(批準(zhǔn)號(hào):11504426)和國(guó)防預(yù)研基金(批準(zhǔn)號(hào):1010502020202)資助的課題.

?通信作者.E-mail:yangli123123@126.com

?2017中國(guó)物理學(xué)會(huì)Chinese Physical Society