李江飛，李巖芳，謝冬梅，段興華，張 康，胡 凱，逯國強，李宏磊

（承德石油高等?？茖W校熱能工程系，河北承德067000）

渦量流函數(shù)法模擬方腔內(nèi)粘性不可壓流動

李江飛，李巖芳，謝冬梅，段興華，張 康，胡 凱，逯國強，李宏磊

（承德石油高等專科學校熱能工程系，河北承德067000）

通過對渦量控制方程進行無量綱化推導和離散，用聯(lián)合迭代方法求解二維方腔流動這一不可壓縮黏性典型流動的解析解.基于Matlab編程，采用渦量流函數(shù)法求解二維方腔流動，計算采用有限體積算法，對流項采用QUICK格式，擴散項采用二階中心差分格式，并采用延遲修正技術的離散格式對該問題進行數(shù)值求解，得到流動達到穩(wěn)定狀態(tài)時各物理量的分布.

數(shù)值模擬；方腔流動；渦量流函數(shù)；有限容積法

Peng等[1]針對二維不可壓縮黏性流體方腔流動邊界條件：流動速度u、v采用無滑移邊界條件，利用動量方程推導了壓力p的邊界條件，本文課題組曾采用交錯網(wǎng)格MAC算法求解二維方腔流動[2]，本文在此基礎上，擬用渦量流函數(shù)法模擬方腔內(nèi)的粘性不可壓流動，以得到流動達到穩(wěn)定狀態(tài)時各物理量的分布.

1 問題分析

1.1 渦量控制方程無量綱化

守恒型渦量控制方程為：

流函數(shù)控制方程為：

用方腔高度H、流體密度ρ,和拖動速度utop作為無量綱標尺，將控制方程無量綱化.無量綱化后的控制方程：

上式中：

其中：H是方腔長度和高度；utop為頂蓋拖動速度；u指變量?在水平x方向的流速；v為變量?在垂直y方向的流速；μ為黏度；ρ為流體密度；ω是渦量；ψ為流函數(shù).

1.2 邊界條件

邊界條件為：流動速度采用無滑移邊界條件，壁面處法向速度恒為0，切向速度也為0.頂蓋u=1, v=0;其余u=v=0.

流函數(shù)的邊界條件為：

渦量的邊界條件為：

2 方程離散與算法分析

2.1 離散控制方程

渦量控制方程為[3-4]：

對流項采用QUICK格式[5]：

對流項渦量內(nèi)邊界點采用一階迎風格式：

擴散項二階中心差分格式為[6]：

離散代數(shù)方程及系數(shù)表達式為：

流函數(shù)的離散采用中心差分格式：

2.2 求解步驟及程序框圖

離散方程可采用聯(lián)合迭代方法求解[7-8]：

（1）根據(jù)ψ0i,j計算u、v,確定離散方程的系數(shù)；

（2）求解渦量代數(shù)方程，得ω0i,j；

（3）由ω0i,j求解流函數(shù)方程，得ψ1i,j；

（4）根據(jù)ψ1i,j和邊界條件確定渦量邊界值ω1B；

（5）根據(jù)ψ1i,j和ω0i,j,迭代計算直至收斂；

（6）利用壓力的Poisson方程計算壓力.

具體流程如圖1所示.

渦量方程式與流函數(shù)方程是互相耦合的，渦量以源項形式出現(xiàn)在流函數(shù)方程中，而渦量方程中對流項的系數(shù)則由流函數(shù)所決定.假設流函數(shù)的一個分布，利用該分布計算u和v，從而確定離散方程的系數(shù)；接下來計算渦量和流函數(shù).在獲得了流函數(shù)的收斂解后，可用各節(jié)點上的流函數(shù)離散處理壓力Poisson方程，從而計算壓力[3].在獲得了流函數(shù)的收斂解后，壓力方程的源項即可用各節(jié)點上的ψ值進行離散處理.

3 結果分析

通過比較不同雷諾數(shù)（100和1000）分別在網(wǎng)格數(shù)為20、40、80、100時的流線圖，可以看出網(wǎng)格數(shù)對結果的影響并不是很大，在精度要求不高時可以采用網(wǎng)格數(shù)較少的計算，節(jié)省時間.穩(wěn)定狀態(tài)時，Re=1000，網(wǎng)格數(shù)=20×20、40×40、80×80、100×100的渦量圖和流函數(shù)圖，計算結果如圖2所示.

圖1 計算程序框圖

穩(wěn)定狀態(tài)時，Re=10、100、1000的流線圖，計算結果如圖3所示.

穩(wěn)定狀態(tài)時，Re=1000，網(wǎng)格數(shù)=20×20，100× 100的壓力、u、v及其中心線速度分布圖，計算結果如圖4所示.

對圖2、圖3、圖4的計算結果進行分析，可知：

（1）二維不可壓縮黏性流體方腔流動的中心大渦并不在中心位置，方腔內(nèi)流動也并不對稱，表明是帶奇性的方腔流動.

（2）本文處理壓力邊界條件時，分別采用動量方程和NeumannBoundary兩種方法，計算結果比較理想.

（3）隨著雷諾數(shù)的增大，在方腔的兩個底角開始出現(xiàn)二次渦流，說明有限體積算法QUICK格式計算精度較高.

（4）通過分析，對流擴散方程中當對流項占主要作用的時候，很容易得到與網(wǎng)格數(shù)無關的解.

（5）本文計算分別采用延遲修正處理和不采用延遲修正處理兩種方法.延遲修正處理能保證方程對角占優(yōu)，減少求解時間.

圖2 不同網(wǎng)格數(shù)下的渦量圖和流函數(shù)圖（Re=1000）

圖3 不同雷諾數(shù)下的流線圖（自左至右：Re=10，100，1000）

圖4 壓力、u、v及其中心線速度分布圖（Re=1000,網(wǎng)格數(shù)=20×20,100×100）

[1]PENGYF,SHIAUYH,HWANGRR.Transitionina2-Dliddrivencavityflow[J].Computer&Fluids,2003,32(3):337-352.

[2]李江飛，石兆東，段興華，等.MAC算法計算二維方腔頂蓋流動[J].宜賓學院學報，2015，15（6）：28-31.

[3]陶文銓.數(shù)值傳熱學[M].第二版.西安:西安交通大學出版社，2001：301-305.

[4]LIUCH,LEUNGDYC.Developmentofafiniteelementsolution fortheunsteadyNavier-Stokesequationsusingprojectionmethod andfractionalθ-scheme[J].Computermethodsinappliedmechan?icsandengineering.2001,190(32-33):4301-4317.

[5]AYDINM,FENNERRT.Boundaryelementanalysisofdrivencav? ityflowforlowandmoderateReynoldsnumbers[J].IntJNumer MethFluids,2001,37(1):45-64.

[6]BRUNEAUCH,SAADM.The2Dlid-drivencavityproblemrevis?ited[J].Computers&Fluids2006,35(3):326-348.

[7]MERCANH,ATALIKK.Vortexformationinlid-drivenarc-shape cavityflowsathighReynoldsnumbers[J].EuropeanJournalofMe?chanicsB/Fluids,2009,28(1):61-71.

[8]WANGJ,LIJF,CHENGWX,etal.ComparisonofFiniteDifference andFiniteVolumeMethodforNumericalSimulationoftheIncom?pressibleViscousDrivenCavityFlow[J].AdvancedMaterialsRe?search,2013.doi:10.4028/www.scientific.net/AMR.732-733.413.

【編校：許潔】

NumericalSimulationofTwo-DimensionalDrivenSquareCavityFlowBasedonStreamFunction-Vor?ticityMethod

LIJiangfei,LIYanfang,XIEDongmei,DUANXinghua,ZHANGKang,HUKai,LUGuoqiang,LIHonglei
(DepartmentofThermalEngineering,ChengdePetroleumCollege,Chengde,Hebei067000,China)

Combinediterationmethodwasusedtogaintheanalyticalsolutionoftwo-dimensionalsquareacvityflowby nondimensionalizedderivationanddiscreteofvorticityequation.BasedonMatlabprogramming,vortexstreamfunction methodsolvedtwo-dimensionalsquarecavityflow.Finitevolumemethodsolvedtheconservationequations.QUICK schemeandtwo-order-central-differenceschemewereappliedintoconvectionanddiffusiontermsseparately.Thedelay correctiontechniqueofdiscreteformatwasappliedintofinitevolumemethodforcalculation.Thedistributionofeach physicalquantitywasobtainedwhentheflowisinasteadystate.

numericalsimulation;cavityflow;vortexstreamfunction;finitevolumemethod

TB126

A

1671-5365(2015)12-0006-04

李江飛，李巖芳，謝冬梅，等.渦量流函數(shù)法模擬方腔內(nèi)粘性不可壓流動[J].宜賓學院學報,2015,15(12):6-9. LIJF,LIYF,XIEDM,etal.NumericalSimulationofTwo-DimensionalDrivenSquareCavityFlowBasedonStreamFunction-VorticityMethod[J].JournalofYibinUniversity,2015,15(12):6-9.

2015-07-06修回：2015-07-31

河北省科技支撐項目“采用電氣液壓調節(jié)的蒸汽透平機運行特性和節(jié)能分析”（13211609）

李江飛（1988-），男，講師，碩士，研究方向為油氣儲運

時間：2015-07-3120:38

http://www.cnki.net/kcms/detail/51.1630.z.20150731.2038.001.html