李鑫郡，張靖周，譚曉茗

南京航空航天大學(xué) 能源與動力學(xué)院，南京 210016

對于壓電風(fēng)扇激勵流動和傳熱性能的研究大多以單個壓電風(fēng)扇作為研究對象。Kim等[4-5]采用PIV(Particle Image Velocimetry)鎖相測試和煙跡顯示方法對自由空間壓電風(fēng)扇激勵的瞬時流場特征進行了研究，揭示出在每一個振動周期，振動膜片葉尖附近形成一對相互逆轉(zhuǎn)的大尺度渦結(jié)構(gòu)，并且在這兩個相互逆轉(zhuǎn)的渦之間形成高速流動；Kimber等[6]對壓電風(fēng)扇誘導(dǎo)的壓力場和拍動射流流動速率進行了實驗測試，研究表明可達流動速度與膜片葉尖速度近似呈現(xiàn)拋物線分布關(guān)系，當(dāng)膜片經(jīng)歷其上下最大振幅位置時上游泵吸能力加強；譚蕾[7-8]和孔岳[9]等采用數(shù)值模擬方法研究了壓電風(fēng)扇誘導(dǎo)的非定常流場特征。Acikalin[10-12]和Kimber[13-14]等開展了一系列的壓電風(fēng)扇傳熱研究，較為系統(tǒng)地揭示了激勵參數(shù)對于換熱特性的影響；Liu等[15]實驗研究了水平和垂直于壁面放置的單個壓電風(fēng)扇傳熱特性；Lin[16]基于振動測試獲得了諧振風(fēng)扇位移函數(shù)，運用動網(wǎng)格技術(shù)研究了壓電風(fēng)扇激勵的非定常射流沖擊壁面的流動和換熱瞬時特征，之后采用同樣的方法研究了強迫流動環(huán)境中圓柱曲面上壓電風(fēng)扇的作用效果[17]，研究發(fā)現(xiàn)在強迫對流通道中，壓電風(fēng)扇激勵射流與橫流外掠圓柱形成的尾跡流相互作用，能夠形成改善圓柱表面對流換熱的效果，但在某些情形下也可能削弱圓柱表面對流換熱，與主流和壓電風(fēng)扇射流特征速度密切關(guān)聯(lián)；Jeng和Liu[18]實驗研究了軸向流通道中熱沉上游安置壓電風(fēng)扇的傳熱特性，研究表明，在通道低雷諾數(shù)流動條件下，壓電風(fēng)扇的諧振可以增強主流的湍流強度，起到改善熱沉傳熱能力的作用。

本文關(guān)注雙壓電風(fēng)扇系統(tǒng)在存在橫流環(huán)境中的傳熱特性及其影響機制?；跈M流環(huán)境中壓電風(fēng)扇振動測試獲得風(fēng)扇位移規(guī)律，結(jié)合動網(wǎng)格技術(shù)數(shù)值模擬橫流環(huán)境中壓電風(fēng)扇激勵的三維非定常流場。本文研究的雙壓電風(fēng)扇采用并列和串列兩種方式，橫流速度與壓電風(fēng)扇的特征速度相當(dāng)，以期揭示雙風(fēng)扇系統(tǒng)在橫流環(huán)境中的流動傳熱機理。

1 計算物理模型

1.1 壓電風(fēng)扇

本文研究以某個特定結(jié)構(gòu)的壓電風(fēng)扇激勵器作為研究對象，該型壓電風(fēng)扇主要由壓電陶瓷片(PZT)、不銹鋼膜片以及底座組成，如圖1所示。壓電陶瓷片單側(cè)粘附于不銹鋼膜片，其長×寬×厚(Lp×W×tp)為28.0 mm×25.0 mm×0.5 mm。柔性膜片伸出長度×厚度(Lb×tb)為38.0 mm×0.1 mm。記風(fēng)扇葉尖沿振動方向處于前后極限位置時的最大位移為App，該位移是葉尖振幅Ap的2倍。

圖1 壓電風(fēng)扇結(jié)構(gòu)示意圖Fig.1 Geometry schematic of piezoelectric fan

1.2 計算域和邊界條件

圖2 數(shù)值計算模型示意圖Fig.2 Schematic of numerical computational model

為了研究雙壓電風(fēng)扇系統(tǒng)在橫流環(huán)境中，排布方式和振動相位差φ對于流場和加熱壁面換熱能力的影響，分別建立如圖2(a)和圖2(b)所示的雙風(fēng)扇并列、雙風(fēng)扇串列布置的物理模型，對于并列布置模式，兩壓電風(fēng)扇之間的無量綱間距設(shè)為P/W=1.2 (其中P為兩風(fēng)扇中心線之間的距離)；對于串列布置模式，兩壓電風(fēng)扇之間的無量綱間距也設(shè)為P/W=1.2。同時為了對比分析，也建立相同尺寸的單風(fēng)扇物理模型，如圖2(c)所示。根據(jù)已有的壓電風(fēng)扇傳熱研究[13-14]，壓電風(fēng)扇具有局部傳熱強化的特征，壓電風(fēng)扇葉尖距離加熱壁面的間距對于流動和傳熱狀況存在重要的影響，過小的間距由于會引發(fā)風(fēng)扇流和壁面的相互作用而導(dǎo)致風(fēng)扇振幅衰減，因此本文選擇的風(fēng)扇葉尖距離加熱壁面的間距G=3 mm。在三維非定常流場數(shù)值計算中，將壓電風(fēng)扇簡化處理為均勻的振動膜片，其運動按照特定的位移函數(shù)設(shè)定。加熱表面(Heat foil)設(shè)置為無滑移速度邊界和恒熱流密度(q=1 700 W/m2)熱邊界條件，外掠加熱表面橫流入口設(shè)為均勻速度進口，橫流出口設(shè)為壓力出口，其余邊界均設(shè)為壓力邊界條件。根據(jù)已有的針對單個壓電風(fēng)扇振動特性的測試研究[19]，該型壓電風(fēng)扇在一階諧振頻率下振動的葉片位移最大速度或特征速度約為2.3 m/s，因此，本文計算模型中的橫流入口速度設(shè)定為uCF=2.3 m/s。

2 計算方法

2.1 壓電風(fēng)扇位移函數(shù)

壓電風(fēng)扇的類懸臂梁振動過程所激發(fā)的流動過程是典型的非穩(wěn)態(tài)運動，且射流速度較低，因此采用非定常不可壓縮Navier-Stokes方程模擬。將壓電風(fēng)扇柔性膜片的運動軌跡y′(z,t)使用如下振型函數(shù)描述：

y′(z,t)=Y(z)sin(2πft)

(1)

式中：Y(z)為風(fēng)扇膜片位移；f為振動頻率；t為模擬運動時間。壓電風(fēng)扇的位移函數(shù)采用文獻[19]的測試擬合多項式，即

Y(z)=p1z4+p2z3+p3z2+p4z+p5

(2)

式(2)的單位為mm。多項式系數(shù)依次為p1=-1.886×10-6,p2=2.447×10-4,p3=-7.31×10-3,p4=6.078×10-2,p5=-5.092×10-2。

2.2 計算網(wǎng)格及獨立性驗證

進行壓電風(fēng)扇三維非定常流場數(shù)值模擬時，將壓電風(fēng)扇位移函數(shù)通過用戶自定義函數(shù)(UDF)輸入，作為壓電風(fēng)扇周期性運動條件，每個周期劃分為200個時間步長。計算過程中使用動網(wǎng)格技術(shù)描述壓電風(fēng)扇的運動，利用局部重構(gòu)方法實現(xiàn)計算過程中網(wǎng)格再生成[19]。計算區(qū)域采用四面體非結(jié)構(gòu)網(wǎng)格，在壓電風(fēng)扇葉尖附近振動位移最劇烈的區(qū)域以及換熱表面附近進行適當(dāng)?shù)木W(wǎng)格加密及自適應(yīng)處理設(shè)置。圖3所示為壓電風(fēng)扇一個周期內(nèi)最大位移時刻的局部網(wǎng)格示意圖。

計算網(wǎng)格數(shù)經(jīng)過網(wǎng)格獨立性驗證確定。圖4是以單風(fēng)扇模型為例，通過數(shù)值模擬得到的計算網(wǎng)格數(shù)對加熱表面時均對流換熱系數(shù)的影響情況。計算過程中通過監(jiān)測每個循環(huán)中的壁面平均溫度隨時間的變化過程，可以發(fā)現(xiàn)在經(jīng)歷約20個周期以后，壁面平均溫度隨時間呈周期性變化。選擇一個周期進行分析，時均對流換熱系數(shù)采用瞬時值積分平均確定:

(3)

式中：Δt為每個周期的時間；dt為時間步長。

可以看出，當(dāng)計算網(wǎng)格數(shù)在100萬和160萬之間時，時均對流換熱系數(shù)的差異僅在3%以內(nèi)，因此最終單風(fēng)扇模型的計算網(wǎng)格數(shù)選取為132萬。運用同樣的方式，確定雙風(fēng)扇模型的計算網(wǎng)格數(shù)為160萬。

圖3 不同相位時局部網(wǎng)格重構(gòu)Fig.3 Re-meshing of local grids at different phases

圖4 網(wǎng)格獨立性驗證Fig.4 Verification of grid independency

2.3 計算方法和驗證

應(yīng)用三維非定常雷諾時均控制方程(RANS)進行數(shù)值模擬。湍流模型選擇SST(Shear Stress Transport)k-ω兩方程湍流模型，該湍流模型被眾多的研究者應(yīng)用于壓電風(fēng)扇流場計算[16-17, 19]。

數(shù)值模擬采用Fluent軟件，在計算過程中需要采用欠松弛加以控制以確保非定常計算過程的穩(wěn)定性，當(dāng)3個連續(xù)的風(fēng)扇振動周期內(nèi)壁面溫度的相對誤差小于10-4時，認(rèn)為數(shù)值解達到準(zhǔn)穩(wěn)定的周期性狀態(tài)[16-17]。

圖5 無橫流單個風(fēng)扇時均局部對流換熱系數(shù)分布Fig.5 Distributions of time-averaged local convective heat transfer coefficients for single fan without cross flow

為驗證本文的數(shù)值計算方法，選用Kimber等[13]的單個壓電風(fēng)扇在無橫流影響下的表面對流換熱特性實驗進行驗證。圖5為不同風(fēng)扇葉尖距離加熱壁面的間距下數(shù)值計算得到的時均局部對流換熱系數(shù)分布情況，圖中虛線框代表風(fēng)扇振動時所掃過的區(qū)域，實線代表風(fēng)扇停振時的位置?？梢钥吹綁弘婏L(fēng)扇作用下的加熱表面對流換熱系數(shù)在葉尖包絡(luò)區(qū)兩側(cè)呈現(xiàn)啞鈴狀分布，與Kimber等[13]實驗揭示的分布特征相符。

圖6為在不同風(fēng)扇葉尖距離加熱壁面的間距下數(shù)值計算得到的駐點努塞爾數(shù)與Kimber等[13]的實驗數(shù)據(jù)對比。此處的壓電風(fēng)扇雷諾數(shù)和努塞爾數(shù)分別定義為

(4)

式中：uPF為壓電風(fēng)扇的特征速度；DPF為風(fēng)扇振動葉尖包絡(luò)區(qū)的特征長度；k為導(dǎo)熱系數(shù)；ν為運動黏度。uPF、DPF的表達式分別為

(5)

比較表明，本文數(shù)值計算結(jié)果與Kimber等[13]的實驗數(shù)據(jù)相吻合。

圖6 駐點努塞爾數(shù)隨無量綱距離G/Ap的變化Fig.6 Stagnation Nusselt number vs non-dimensional gap G/Ap

3 結(jié)果與分析

3.1 瞬態(tài)流場特征

壓電風(fēng)扇振動誘導(dǎo)一系列非定常的拍動渦環(huán)并不斷地聚合，從而形成風(fēng)力相對集中的“準(zhǔn)連續(xù)”射流。為了揭示壓電風(fēng)扇射流與橫流相干的渦系結(jié)構(gòu)，借鑒Jeong和Hussain[20]提出的λ2渦識別判據(jù)進行分析。該判據(jù)是將流場的速度梯度張量J分解為對稱部分S(應(yīng)變率張量)和非對稱部分Ω(旋轉(zhuǎn)張量)，并通過計算組合張量S2+Ω2的3個特征值(λ1≥λ2≥λ3)，認(rèn)為壓力達到截面最小的充要條件為λ2<0，其中λ2為負(fù)值的點即屬于渦核空間位置。J、S、Ω三者之間的關(guān)系式為

(6)

圖7 單風(fēng)扇不同振動相位瞬時溫度和λ2=-3×104 瞬態(tài)等值面Fig.7 Instantaneous temperature and iso-surface for λ2=-3×104 at different vibrating phases for single fan

圖7為單個壓電風(fēng)扇在一個運動周期內(nèi)4種典型振動相位下λ2=-3×104的瞬態(tài)等值面發(fā)展演變情況。當(dāng)壓電風(fēng)扇運動至其平衡位置時，如圖7(b)和圖7(d)所示，其葉尖位移速度達到最大值，依靠風(fēng)扇葉片對氣流的剪切機制在葉尖以及葉片兩側(cè)緣形成較為強烈的馬蹄渦，而當(dāng)風(fēng)扇運動至其最大振幅位置時，由于風(fēng)扇位移速度逐漸減小，渦結(jié)構(gòu)顯著減弱，此時的橫流效應(yīng)相應(yīng)增強，挾帶葉尖脫落的渦向下游遷移。從葉尖脫落的渦對壁面形成沖擊作用，從而導(dǎo)致壓電風(fēng)扇振動的包絡(luò)區(qū)域附近的局部對流換熱強化；由于周期性的風(fēng)扇振動，葉尖脫落渦在包絡(luò)區(qū)下游形成周期性掃掠，同時流經(jīng)葉尖包絡(luò)區(qū)的橫流也受到壓電風(fēng)扇的脈動激勵，因而在壓電風(fēng)扇下游形成溫度相對較低的條帶狀區(qū)域。

圖8 雙風(fēng)扇并列反相振動的瞬時溫度和λ2=-3×104 瞬態(tài)等值面Fig.8 Instantaneous temperature and iso-surface for λ2=-3×104 for dual out-of-phase fans in face-to-face mode

圖8和圖9分別為雙風(fēng)扇并列模式下反相和同相振動時λ2=-3×104瞬態(tài)等值面以及加熱壁面瞬時溫度分布情況。與單個風(fēng)扇相比，由于雙風(fēng)扇的作用區(qū)域增大，因此加熱表面的冷卻效果得到增強。同時相鄰風(fēng)扇間的振動對于渦系的演變也具有一定的影響，對于并列布置的雙壓電風(fēng)扇，振動相位引起的相鄰風(fēng)扇之間渦結(jié)構(gòu)演化特征具有較大的差異。并列雙風(fēng)扇反相振動時，如圖8所示，相鄰兩個風(fēng)扇相向運動至各自的平衡位置時，葉尖沖擊脫落渦之間形成融合和發(fā)展，并對其他時刻的渦系向下游的發(fā)展也產(chǎn)生影響，從而在壓電風(fēng)扇下游形成溫度更低、范圍更寬的影響區(qū)域。而對于并列雙風(fēng)扇同相振動，雖然存在著相鄰風(fēng)扇間的相互作用，但每個風(fēng)扇誘導(dǎo)的渦結(jié)構(gòu)與單個風(fēng)扇較為接近。

圖10和圖11分別為雙風(fēng)扇串列模式下反相和同相振動時的渦核瞬態(tài)等值面以及加熱壁面的瞬時溫度分布情況。由于兩個風(fēng)扇串列于橫流中，上游的壓電風(fēng)扇誘導(dǎo)的渦在橫流作用下向下游遷移，與下游壓電風(fēng)扇激勵作用形成耦合，因此，相對于單一風(fēng)扇，串列雙壓電風(fēng)扇葉尖脫落的渦對下游表面對流換熱影響作用的區(qū)域更長。對比雙風(fēng)扇串列模式下反相和同相振動時的渦結(jié)構(gòu)演變，可以看出，同相振動情形下相鄰風(fēng)扇的脫落渦之間的融合和發(fā)展更具優(yōu)勢。

圖9 雙風(fēng)扇并列同相振動的瞬時溫度和λ2=-3×104瞬態(tài)等值面Fig.9 Instantaneous temperature and iso-surface of λ2=-3×104 for dual in-phase fans in face-to-face mode

圖10 雙風(fēng)扇串列反相振動的瞬時溫度和λ2=-3×104瞬態(tài)等值面Fig.10 Instantaneous temperature and iso-surface for λ2=-3×104 for dual out-of-phase fans in edge-to-edge mode

圖11 雙風(fēng)扇串列同相振動的瞬時溫度和λ2=-3×104瞬態(tài)等值面Fig.11 Instantaneous temperature and iso-surface for λ2=-3×104 for dual in-phase fans in edge-to-edge mode

3.2 傳熱時均結(jié)果分析

圖12為不同壓電風(fēng)扇作用下的一個周期內(nèi)恒熱流壁面時均局部對流換熱系數(shù)分布情況。對于單個風(fēng)扇，如圖12(a)所示，受橫流的影響，風(fēng)扇葉尖掃掠的包絡(luò)區(qū)的對流換熱系數(shù)分布呈現(xiàn)水滴狀，與單個壓電風(fēng)扇在無橫流環(huán)境中的啞鈴形分布不同(圖5(a))，顯然這是橫流對風(fēng)扇激勵射流沖擊的影響所致。對于并列雙風(fēng)扇，如圖12(b)和圖12(c)所示，風(fēng)扇葉尖掃掠的包絡(luò)區(qū)附近局部對流換熱系數(shù)依然呈現(xiàn)水滴狀分布，在包絡(luò)區(qū)的時均換熱能力與單風(fēng)扇幾乎相同，但相對于單風(fēng)扇，水滴狀在流向上有所延伸，這是相鄰風(fēng)扇脫落渦之間的融合和發(fā)展所致，注意到風(fēng)扇之間的干涉作用主要體現(xiàn)在下游區(qū)域，反相振動時風(fēng)扇間隙下游區(qū)域的換熱能力明顯得到強化。對于串列雙風(fēng)扇排布，如圖12(d)和圖12(e)所示，反相振動時上游風(fēng)扇包絡(luò)區(qū)的換熱能力與單風(fēng)扇情況基本一致，下游風(fēng)扇包絡(luò)區(qū)的換熱能力有所弱化，但下游風(fēng)扇尾跡區(qū)域的換熱能力明顯高于單風(fēng)扇的情況；同相振動時，上游風(fēng)扇與下游風(fēng)扇包絡(luò)區(qū)內(nèi)的換熱能力均高于單風(fēng)扇情況，兩風(fēng)扇間隙區(qū)域的換熱能力也得到明顯加強，綜合而言，在橫流中串列雙風(fēng)扇同相振動時，局部換熱能力可以得到更大程度的強化。

為了進一步對比分析各種排布方式和不同振動相位下壁面平均換熱能力的差異，沿橫流方向定義側(cè)向積分平均對流換熱系數(shù)為

(7)

式中：Ly為側(cè)向積分的參考長度。對于單風(fēng)扇和串列雙風(fēng)扇，Ly分別取App和2App；對于并列雙風(fēng)扇，Ly取為雙風(fēng)扇間距P。

圖13為沿橫流方向的側(cè)向積分平均對流換熱系數(shù)。圖中豎實線、豎虛線分別表示壓電風(fēng)扇中心線位置、風(fēng)扇側(cè)緣位置。對于單個風(fēng)扇，如圖13(a)所示， 在單風(fēng)扇振動包絡(luò)區(qū)， 按照App和2App平均的對流換熱系數(shù)存在較大的差異，而在單風(fēng)扇振動包絡(luò)區(qū)下游1.5倍風(fēng)扇寬度(1.5W)之后，按照App和2App平均的對流換熱系數(shù)相等，說明壓電風(fēng)扇的影響趨于消失。對于并列雙風(fēng)扇，如圖13(b)所示，同相和反相振動對上游以及包絡(luò)區(qū)的對流換熱影響并無明顯差異，但是在下游1W之后反相振動的平均對流換熱系數(shù)明顯高于同相振動；對于串列雙風(fēng)扇，如圖13(c)和圖13(d)所示，同相振動始終優(yōu)于反相振動情況，在反相振動時，相鄰兩風(fēng)扇之間的對流換熱顯著低于同相振動的方式。

圖12 有橫流時加熱表面時均局部對流換熱系數(shù)分布Fig.12 Distributions of time-averaged local convective heat transfer coefficients on heat foil in presence of cross flow

圖13 側(cè)向積分平均對流換熱系數(shù)沿橫流方向的分布Fig.13 Distribution of laterally-averaged convective heat transfer coefficient along cross flow direction

4 結(jié) 論

本文采用動網(wǎng)格技術(shù)對橫流環(huán)境中垂直加熱壁面的單個壓電風(fēng)扇、并列雙風(fēng)扇以及串列雙風(fēng)扇三維非定常流動和傳熱特性進行了數(shù)值模擬，結(jié)果表明：

1) 存在橫流時，單個風(fēng)扇葉尖掃掠包絡(luò)區(qū)的對流換熱系數(shù)分布呈現(xiàn)水滴狀，與無橫流影響下的啞鈴形分布不同。

2) 對于并列布置的雙壓電風(fēng)扇，反相振動有利于相鄰風(fēng)扇脫落渦之間的融合和發(fā)展，雙風(fēng)扇振動包絡(luò)區(qū)下游的換熱能力優(yōu)于同相振動情形。

3) 對于串列布置的雙壓電風(fēng)扇，上游的壓電風(fēng)扇誘導(dǎo)的渦動射流在橫流作用下向下游遷移，與后一個壓電風(fēng)扇激勵作用形成耦合。反相振動時相鄰兩風(fēng)扇之間局部區(qū)域的對流換熱顯著低于同相振動的情形。

[1] YOO J H, HONG J I, CAO W. Piezoelectric ceramic bimorph coupled to thin metal plate as cooling fan for electronic devices[J]. Sensors Actuators A, 2000, 79(1): 8-12.

[2] CHUNG H C, KUMMARI K L, CROUCHER S J, et al. Development of piezoelectric fans for flapping wing application[J]. Sensors and Actuators A: Physical, 2009, 149(1): 136-142.

[3] GILSON G M, PICKERING S J, HANN D B, et al. Piezoelectric fan cooling: A novel high reliability electric machine thermal management solution[J]. IEEE Transactions on Industrial Electronics, 2013, 60(11): 4841-4851.

[4] KIM Y H, WERELEY S T, CHUN C H. Phase-resolved flow field produced by a vibrating cantilever plate between two endplates[J]. Physics in Fluids, 2004, 16(1): 145-162.

[5] KIM Y H, CHUN C H, WERELEY S. Flow field around a vibrating cantilever: Coherent structure education by continuous wavelet transform and proper orthogonal decomposition[J]. Journal of Fluid Mechanics, 2011, 669: 584-606.

[6] KIMBER M, SUZUKI K, KITSUNAI N, et al. Pressure and flow rate performance of piezoelectric fans[J]. IEEE Transactions on Components and Packaging Technologies, 2009, 32(4): 766-775.

[7] 譚蕾， 譚曉茗， 張靖周. 壓電風(fēng)扇激勵非定常流動和換熱特性數(shù)值研究[J]. 航空學(xué)報， 2013, 34(6): 1277-1284.

TAN L, TAN X M, ZHANG J Z. Numerical investigation on unsteady flow and heat transfer characteristics of piezoelectric fan[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(6): 1277-1284 (in chinese).

[8] TAN L, ZHANG J Z, TAN X M. Numerical investigation of convective heat transfer on a vertical surface due to resonating cantilever beam[J]. International Journal of Thermal Sciences, 2014, 80: 93-107.

[9] 孔岳， 李敏， 吳蒙蒙. 壓電風(fēng)扇非定常流場速度分布的數(shù)值研究[J]. 工程力學(xué)， 2016, 33(1): 209-216.

KONG Y, LI M, WU M M. Numerical investigation on the velocity of unsteady flow field induced by piezoelectric fan[J]. Engineering Mechanics, 2016, 33(1): 209-216 (in Chinese).

[10] ACIKALIN T, WAIT S M, GARIMELLA S V, et al. Experimental investigation of the thermal performance of piezoelectric fans[J]. Heat Transfer Engineering, 2004, 25(1): 4-14.

[11] ACIKALIN T, GARIMELLA S V, RAMAN A, et al. Characterization and optimization of the thermal performance of miniature piezoelectric fans[J]. International Journal of Heat and Fluid Flow, 2007, 28(4): 806-820.

[12] ACIKALIN T, GARIMELLA S V. Analysis and prediction of the thermal performance of piezoelectrically actuated fans[J]. Heat Transfer Engineer, 2009, 30(6): 487-498.

[13] KIMBER M, GARIMELLA S V, RAMAN A. Local heat transfer coefficients induced by piezoelectrically actuated vibrating cantilevers[J]. ASME Journal of Heat Transfer, 2007, 129(9): 1168-1176.

[14] KIMBER M, GARIMELLA S V. Measurement and prediction of the cooling characteristics of a generalized vibrating piezoelectric fan[J]. International Journal of Heat and Mass Transfer, 2009, 52(19-20): 4470-4478.

[15] LIU S F, HUANG R T, SHEU W J, et al. Heat transfer by a piezoelectric fan on a flat surface subject to the influence of horizontal/vertical arrangement[J]. International Journal of Heat and Mass Transfer, 2009, 52(11-12): 2565-2570.

[16] LIN C N. Analysis of three-dimensional heat and flow induced by piezoelectric fan[J]. International Journal of Heat and Mass Transfer, 2012, 55(11-12): 3043-3053.

[17] LIN C N. Enhanced heat transfer performance of cylindrical surface by piezoelectric fan under forced convection conditions[J]. International Journal of Heat and Mass Transfer, 2013, 60: 296-308.

[18] JENG T M, LIU C H. Moving-orientation and position effects of the piezoelectric fan on thermal characteristics of the heat sink partially filled in a channel with axial flow[J]. International Journal of Heat and Mass Transfer, 2015, 85: 950-964.

[19] 李鑫郡， 張靖周， 譚曉茗. 單個壓電風(fēng)扇傳熱特性研究[J]. 航空學(xué)報， 2017, 38(7): 120982.

LI X J, ZHANG J Z, TAN X M. Characteristics of heat transfer with a single piezoelectric fan[J]. Acta Aeronautica et Astronautica Sinica, 2017, 38(7): 120982 (in Chinese).

[20] JEONG J, HUSSAIN F. On the definition of a vortex[J]. Journal of Fluid Mechanics, 1995, 285: 69-94.