水平管氣液兩相分層流界面剪切預(yù)測(cè)研究進(jìn)展

鄭平，趙梁（遼寧石油化工大學(xué)石油天然氣工程學(xué)院，遼寧 撫順 113001）

?

綜述與專論

水平管氣液兩相分層流界面剪切預(yù)測(cè)研究進(jìn)展

鄭平，趙梁
（遼寧石油化工大學(xué)石油天然氣工程學(xué)院，遼寧 撫順 113001）

摘要：水平管氣液兩相分層流雖流型簡(jiǎn)單，但由于界面存在復(fù)雜的動(dòng)量和能量傳遞，分層流的界面剪切預(yù)測(cè)至今沒有一致的結(jié)論。本文從理論模型、實(shí)驗(yàn)?zāi)Ｐ汀?shù)值計(jì)算 3個(gè)角度出發(fā)，詳細(xì)闡述水平管氣液兩相分層流界面剪切預(yù)測(cè)的研究現(xiàn)狀，得出不同研究方法的優(yōu)勢(shì)和缺陷。針對(duì) 3種研究方法，指出理論模型通過模型簡(jiǎn)化和經(jīng)驗(yàn)關(guān)聯(lián)式來建立封閉模型，實(shí)驗(yàn)?zāi)Ｐ蛣t在封閉關(guān)系上修正經(jīng)驗(yàn)關(guān)聯(lián)式，但由于簡(jiǎn)化假設(shè)和實(shí)驗(yàn)條件的限制，使得這兩種研究方法對(duì)界面剪切應(yīng)力的預(yù)測(cè)具有一定的局限性；數(shù)值計(jì)算能夠彌補(bǔ)機(jī)理模型在流場(chǎng)細(xì)節(jié)等方面的不足，但能夠提供界面剪切預(yù)測(cè)或封閉關(guān)系的工作很少。此外，對(duì)比了 5種不同形式的已有模型對(duì)氣液兩相分層流持液率和壓降預(yù)測(cè)的結(jié)果。最后展望了水平管氣液兩相分層流界面剪切預(yù)測(cè)的研究趨勢(shì)，提出理論和實(shí)驗(yàn)研究需要提出更詳細(xì)的局部模型，并考慮工程實(shí)際工況進(jìn)行研究，發(fā)展針對(duì)氣液界面計(jì)算的新方法，并為分層流提供封閉關(guān)系則是數(shù)值計(jì)算研究面臨的挑戰(zhàn)。

關(guān)鍵詞：水平管；氣液兩相流；分層流；界面剪切應(yīng)力；持液率；壓降

第一作者及聯(lián)系人：鄭平（1973—），女，博士，副教授，主要從事多相流體力學(xué)應(yīng)用研究。E-mail zhengping0413@163.com。

水平管氣液兩相分層流是指存在密度差的兩相流體，在質(zhì)量力（如重力或其他慣性力）作用下，密度較大的相在下層流動(dòng)，密度較小的相在上層流動(dòng)，相界面存在清晰可辨的薄膜［1-2］。水平管氣液兩相分層流是化工工程和石油天然氣工程中最常見的兩相流流型之一，分層流的特點(diǎn)是流場(chǎng)除受到氣壁和液壁的剪切作用外，還受到氣液界面相互作用的影響。由于氣液界面存在復(fù)雜的動(dòng)量和能量傳遞，因此目前對(duì)分層流界面剪切應(yīng)力的預(yù)測(cè)至今沒有一致的結(jié)論。持液率和壓降是分層流研究中最重要的兩個(gè)參數(shù)，而氣液界面剪切應(yīng)力的計(jì)算與分層流持液率和壓降的預(yù)測(cè)密切相關(guān)，對(duì)工程實(shí)際應(yīng)用具有重要意義。

理論、實(shí)驗(yàn)和數(shù)值計(jì)算是研究水平管分層流界面剪切預(yù)測(cè)的3類基本方法。理論建模和實(shí)驗(yàn)的間接測(cè)量都是在嘗試分析界面相互作用機(jī)理基礎(chǔ)上提出不同的計(jì)算模型，因此也被稱為機(jī)理性模型［3］，機(jī)理性模型的共同特點(diǎn)是需要人為提供與經(jīng)驗(yàn)關(guān)系同樣數(shù)量的本構(gòu)關(guān)系式。數(shù)值計(jì)算即計(jì)算流體力學(xué)（computational fluid dynamics，CFD）方法，通常依賴?yán)碚摵蛯?shí)驗(yàn)的結(jié)果，在機(jī)理性模型的支持下深入觀察流場(chǎng)細(xì)節(jié)特征，但討論氣液界面剪切應(yīng)力本構(gòu)關(guān)系的工作并不多。

本文從理論、實(shí)驗(yàn)、數(shù)值計(jì)算3類基本方法入手，對(duì)氣液兩相分層流界面剪切預(yù)測(cè)問題進(jìn)行綜述，通過對(duì)比分析5種已有模型得到的持液率和壓降結(jié)果，找出各自的優(yōu)勢(shì)與缺陷，提出此問題的研究方向和展望。

1 氣液兩相分層流界面剪切預(yù)測(cè)理論

氣液兩相分層流的理論模型從數(shù)學(xué)模型和幾何模型兩方面論述。

1.1 數(shù)學(xué)模型

目前，已經(jīng)建立了兩類數(shù)學(xué)模型，即整體模型和與之互補(bǔ)的局部模型［4］。

早期的兩相流模型對(duì)兩相流動(dòng)作了很大程度的簡(jiǎn)化假設(shè)，如最早提出水平管分層流壓降和持液率的LOCKHART等［5］修正公式，此模型通過單相流模型計(jì)算兩相流流動(dòng)，適應(yīng)性差，所以預(yù)測(cè)準(zhǔn)確率偏低。

TAITEL等［6-7］建立了一維兩流體模型，該模型假設(shè)氣液界面為水平的光滑壁面，從而將氣相流動(dòng)假設(shè)為單相流，氣壁和氣液界面具有相同的摩擦因子，同時(shí)假設(shè)液相為明渠流動(dòng)，相界面為一水平面，為此對(duì)氣、液兩相分別建立動(dòng)量方程，這是整體模型的典型代表，許多理論模型及實(shí)驗(yàn)都建立在這一模型之上。但它忽略了流場(chǎng)的結(jié)構(gòu)細(xì)節(jié)，對(duì)于大雷諾數(shù)湍流，通過平均速度前提下的經(jīng)驗(yàn)關(guān)系式計(jì)算剪切應(yīng)力往往偏小。因此該模型適用于層流和雷諾數(shù)較小的湍流分層流［8］。

HANSEN等［9］建立了二維移動(dòng)壁面數(shù)學(xué)模型，通過雙對(duì)數(shù)律速度分布假設(shè)和界面的無滑移自由移動(dòng)壁面假設(shè)，按照近壁速度分布規(guī)律分析，針對(duì)界面光滑和波動(dòng)兩種情形，得到了平穩(wěn)湍流的氣液兩相流流場(chǎng)分布規(guī)律。計(jì)算結(jié)果與實(shí)驗(yàn)相對(duì)比表明，二維移動(dòng)壁面模型可以較好地預(yù)測(cè)光滑分層流和波狀分層流，但是對(duì)后者的氣液界面速度分布預(yù)測(cè)效果較差。

ULLMANN等［10-11］對(duì)層流分層流和湍流分層流（光滑與波狀分層流）建立的封閉關(guān)系，是在TAITEL& DUKLER的動(dòng)量方程基礎(chǔ)上，結(jié)合已有的模型結(jié)果導(dǎo)出的半理論半經(jīng)驗(yàn)分層流封閉關(guān)聯(lián)式。但是，由于已有的不同模型對(duì)界面和液壁剪切應(yīng)力的預(yù)測(cè)偏差很大，必然導(dǎo)致關(guān)聯(lián)式在一定條件下（如當(dāng)界面和液壁剪切應(yīng)力預(yù)測(cè)同時(shí)偏高時(shí)）不能準(zhǔn)確預(yù)測(cè)持液率。

兩流體模型的缺陷在于將分層流看作疊置的兩相分離的流動(dòng)，其界面波動(dòng)的幅度小于界面厚度，減弱了預(yù)測(cè)效果。為此，KANG等［12］將氣液界面剪切應(yīng)力視為附著于光滑界面的剪切和由波浪移動(dòng)形成的剪切的總和，提出了剪切應(yīng)力預(yù)測(cè)關(guān)聯(lián)式。為了簡(jiǎn)化模型，此模型只考慮界面波附近流動(dòng)分離較弱的情形，因此只適用于波動(dòng)現(xiàn)象相對(duì)穩(wěn)定的分層流。此外，對(duì)界面形狀和界面速度分布形式，此模型并沒有考慮。

LINé等［4］發(fā)展了一種局部?jī)闪黧w模型來計(jì)算氣液界面的波動(dòng)區(qū)域，他們并不預(yù)測(cè)界面形狀，而是將界面假設(shè)為一已知的先驗(yàn)情形，該模型用于計(jì)算界面由壓強(qiáng)和黏性剪切分布造成的動(dòng)量傳遞。波動(dòng)理論在分層流流動(dòng)行為研究中的應(yīng)用也是實(shí)驗(yàn)的重要課題之一，但對(duì)波動(dòng)產(chǎn)生、變形、傳播機(jī)理以及波動(dòng)上下的湍流區(qū)域并沒有進(jìn)行分析。

劉夷平等［13］分析了液壁、氣液間剪切應(yīng)力預(yù)測(cè)的不確定度，并利用最大不確定度原理，對(duì)影響其準(zhǔn)確度的各因素進(jìn)行誤差分析，修正了ANDRITSOS等［14］的實(shí)驗(yàn)?zāi)Ｐ汀＝Y(jié)果表明，氣相摩擦因子可以用單相流的指數(shù)關(guān)系式表示，而液相摩擦因子受氣液相界面剪切應(yīng)力的影響較大。

1.2 幾何模型

TAITEL & DUKLER模型是基于水平界面幾何模型建立的數(shù)學(xué)模型，對(duì)于液相表觀速度不大或彎管處的低持液率問題，由于BOND數(shù)減小，一些學(xué)者提出了非水平面幾何模型，著重考察橫向界面形狀對(duì)分層流封閉模型的可能性。

HART等［15］將液相橫截面區(qū)域處理為均勻厚度、緊貼管底的凹面，而 CHEN等［16］以及之后的VLACHOS等［17-18］將液相橫截面區(qū)域處理為雙圓交叉形成的月牙形，見圖 1。這些非水平界面模型很大程度上彌補(bǔ)了早期理論模型在低持液率分層流的預(yù)測(cè)誤差。

圖1 分層流幾何模型示意圖

水平界面幾何模型應(yīng)用范圍比較廣，等厚度和凹形界面幾何模型通常用于低持液率問題，凹形界面問題在討論半環(huán)狀流時(shí)非常有效，對(duì)于非牛頓流體分層流，凹形界面和與之對(duì)應(yīng)的凸形界面幾何模型都有廣泛的研究［19］。

BRAUNER等［20］在預(yù)測(cè)界面形狀時(shí)考慮了能量，并分析了液體的物性、持液率、管道尺寸、壁面附著力和重力對(duì)界面形狀的影響。界面曲率預(yù)測(cè)提供了完美解決彎曲界面分層流所需的封閉關(guān)系式。后來，結(jié)合BRAUNER等［21］的流動(dòng)關(guān)系式，又得到了一系列層流兩流體系統(tǒng)界面形狀和相關(guān)的局部和整體流動(dòng)特性［22-23］，但只準(zhǔn)確地分析了二維層流分層流的流速情況。1998年，BRAUNER等［24］發(fā)展了一種基于彎曲氣液界面的兩流體模型，確定界面曲率與給定操作條件下持液率的封閉關(guān)系式。此模型能夠?qū)Τ忠郝屎蛪航堤峁┮环N可靠的預(yù)測(cè)，對(duì)于湍流流動(dòng)，此兩流體模型的結(jié)果更準(zhǔn)確。但當(dāng)兩流體模型在反常的高黏度環(huán)境中應(yīng)用時(shí)，會(huì)發(fā)生最大的誤差。

ROVINSKY等［25］對(duì)完全偏心的核心環(huán)形層流兩相流進(jìn)行了研究，用解析法求出持液率和壓降。結(jié)果表明，當(dāng)黏性相與管壁之間的連接區(qū)域增大時(shí)，潤(rùn)滑相面的形狀就可能是凸形，這時(shí)壓降減小的潛能會(huì)進(jìn)一步減小。

SOLEIMANI［26］通 過 數(shù) 值 方 法 計(jì) 算 YOUNG-LAPLACE方程，確定液-液系統(tǒng)界面形狀。結(jié)果表明，界面形狀取決于接觸角、BOND數(shù)和兩流體當(dāng)?shù)伢w積分?jǐn)?shù)的比率。

此外，非線性理論氣液兩相流界面特性（如界面波與不穩(wěn)定性）的研究成果很有價(jià)值［16，27-29］，這些方法不僅為建立更接近實(shí)際情況的理論模型提供依據(jù)和思路，也為實(shí)驗(yàn)數(shù)據(jù)處理提供了方便有效的方法。

2 氣液兩相分層流界面剪切預(yù)測(cè)實(shí)驗(yàn)

理論與實(shí)驗(yàn)相結(jié)合一直是研究水平管氣液兩相分層流界面剪切預(yù)測(cè)的主要方法，測(cè)量方法和數(shù)據(jù)處理是實(shí)驗(yàn)的基本內(nèi)容，理論模型通常由實(shí)驗(yàn)對(duì)比結(jié)果得到合理解釋。

KOWALSKI［30］測(cè)量了空氣-水和氟利昂 12-水系統(tǒng)的氣壁和液壁剪切應(yīng)力，得出新的壁面摩擦因子關(guān)聯(lián)式（稱為KOWALSKI形式公式）。結(jié)果表明，BLASIUS方程預(yù)測(cè)的氣壁摩擦因子與實(shí)驗(yàn)結(jié)果比較吻合；在得到氣壁與液壁剪切應(yīng)力后，通過外推界面附近的雷諾數(shù)剪切分布，分別得到了光滑和波動(dòng)界面的剪切應(yīng)力，建立了KOWALSKI形式表達(dá)的分層流封閉關(guān)系。

ANDRITSOS等［14］、HAGIWARA等［31］通過實(shí)驗(yàn)分別分析了界面行為（如界面波）對(duì)氣壁和界面剪切應(yīng)力的影響，并討論了氣壁和界面摩擦因子的本構(gòu)關(guān)系。ANDRITSOS等［14］利用動(dòng)量守恒預(yù)測(cè)了界面剪切應(yīng)力，同時(shí)提出過渡氣相表觀速度概念，以此作為光滑與波狀分層流的界限。結(jié)果表明，對(duì)于光滑分層流，界面摩擦因子等于氣壁摩擦因子；當(dāng)出現(xiàn)界面波，即氣相表觀速度大于過渡值，界面與氣壁的摩擦因子之比為液高（實(shí)際液高/管徑）和氣相表觀速度的本構(gòu)關(guān)系式。HAGIWARA 等［31］用熱膜電導(dǎo)探針和熱線測(cè)定了充分發(fā)展的氣液兩相流的壁面剪切應(yīng)力、液膜厚度和氣體速度，并觀察了氣壁剪切與界面波的關(guān)系。實(shí)驗(yàn)結(jié)果表明，氣液界面上的波幅增大時(shí)，氣壁剪切應(yīng)力也隨之增大，為界面摩擦因子與氣壁摩擦因子之比作為分層流封閉關(guān)系提供了合理依據(jù)。

TZOTZI等［32］在內(nèi)徑為0.025m的水平管內(nèi)，以空氣、二氧化碳、氦氣作為氣相，水、丁醇-水溶液作為液相，進(jìn)行氣液兩相分層流動(dòng)實(shí)驗(yàn)，并在ANDRITSOS & HANRATTY模型的基礎(chǔ)上，著重考慮氣體密度、表面張力兩個(gè)因素，完善從光滑流到2D波、從2D波到大幅度波轉(zhuǎn)型的半理論相關(guān)性和兩類波狀流界面摩擦因子的經(jīng)驗(yàn)關(guān)系式。結(jié)果表明，轉(zhuǎn)型相關(guān)性與實(shí)驗(yàn)數(shù)據(jù)符合良好。

STRAND［33］的實(shí)驗(yàn)在長(zhǎng)度為 35m、內(nèi)徑為10cm的水平管道中進(jìn)行，利用探針測(cè)量了空氣-水分層流流場(chǎng)分布，得到了氣壁、液壁和界面剪切應(yīng)力，并通過摩擦因子與壁面粗糙度關(guān)系［34］計(jì)算了分層流界面等效的粗糙度，其結(jié)果為許多新的CFD方法提供了數(shù)據(jù)基礎(chǔ)。

LI等［35］用兩種方法測(cè)定界面剪切應(yīng)力，對(duì)氣相采用動(dòng)量平衡法，對(duì)氣-液界面采用雷諾剪切應(yīng)力外推法，得到一種新的預(yù)測(cè)界面剪切應(yīng)力公式。結(jié)果表明，界面摩擦因子與空隙率和氣體表觀雷諾數(shù)有關(guān)，預(yù)測(cè)結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好。

3 氣液兩相分層流界面剪切預(yù)測(cè)數(shù)值計(jì)算

目前，數(shù)值計(jì)算方法（CFD）在預(yù)測(cè)水平管氣液界面剪切應(yīng)力的研究應(yīng)用已經(jīng)越來越普遍［36］。DOMINIQUE［37］在討論氣液兩相流所有流型的CFD建模的困難時(shí)指出，由于氣液兩相流界面的離散性質(zhì)，CFD方法受限于空間和時(shí)間的精確計(jì)算。

較早的應(yīng)用 CFD研究氣液兩相分層流的工作主要集中于利用MARKER & CELL或Volume of Fluid（VOF）方法模擬分析2D/3D自由表面流動(dòng)或氣泡的生長(zhǎng)。AKAI等［38］提出了局部一維數(shù)值分層流兩方程湍流模型，預(yù)測(cè)了氣液兩相的壓力梯度、持液率和速度，與實(shí)驗(yàn)數(shù)據(jù)吻合較好，但未能給出界面波動(dòng)于氣相雷諾數(shù)的經(jīng)驗(yàn)關(guān)系式。

SHOHAM等［8］最早用 CFD建立了二維分層湍-湍氣液兩流體模型，應(yīng)用渦黏理論求解液相N-S動(dòng)量方程，用零方程模型計(jì)算湍流黏度，在修正界面應(yīng)力時(shí)考慮了界面結(jié)構(gòu)，預(yù)測(cè)出液相的速度場(chǎng)、壓降和持液率，得到了比 LOCKHART & MARTINELLI和TAITEL & DUKLER更好的結(jié)果。

ISSA［39］對(duì)充分發(fā)展的氣液兩相分層流建立了二維湍流-湍流k-ε兩方程模型，采用曲線正交網(wǎng)格和雙極坐標(biāo)系，得到了接近實(shí)驗(yàn)數(shù)據(jù)的氣液兩相壓力梯度、速度場(chǎng)和界面高度。但如何指定氣液界面（特別是波狀分層流界面）的湍動(dòng)能和湍動(dòng)能耗散率邊界條件，尚未得到解決。模型中，包括氣液界面的所有邊界都要使用壁面函數(shù)，并需要借助額外的經(jīng)驗(yàn)公式，且只適用于內(nèi)徑小于或等于25mm的管道。

MEKNASSI等［40］的工作很有代表性，將氣液界面處理為具有等效粗糙度的粗糙壁面，壁面函數(shù)采用低雷諾數(shù)k-ε兩方程湍流模型，通過有限差分法成功地模擬了STRAND［33］和LOPEZ［41］的實(shí)驗(yàn)。結(jié)果表明，當(dāng)考慮二次流時(shí)，模擬結(jié)果與實(shí)驗(yàn)結(jié)果符合良好，氣相的模擬結(jié)果證實(shí)了界面粗糙度能夠反映氣液界面的相互作用。

NEWTON等［42］通過動(dòng)量方程改進(jìn)了50mm和80mm內(nèi)徑氣液水平管流KOWALSKI形式的界面摩擦因子關(guān)聯(lián)式，其思路是盡可能少地使用經(jīng)驗(yàn)關(guān)系式，以提高結(jié)論的普適性。后來，NEWTON等［43-44］又分別對(duì)湍流氣液兩相和波動(dòng)氣液兩相分層流進(jìn)行了數(shù)值模擬。前者氣液界面按照光滑處理，湍流模型只引用了單相流阻尼函數(shù)一個(gè)經(jīng)驗(yàn)公式，二維軸向動(dòng)量方程式在氣液各相中分別求解，湍流黏度由低雷諾數(shù)k-ε模型求得；后者將界面處理為具有等效粗糙度的壁面，除阻尼函數(shù)外，還引用了粗糙度與摩擦因子經(jīng)驗(yàn)關(guān)系式，給出了界面剪切應(yīng)力經(jīng)驗(yàn)分布函數(shù)，由k-ε湍流模型進(jìn)行求解。

MOUZA等［3］利用CFX軟件模擬了VLACHOS等［45］的實(shí)驗(yàn)，將氣液兩相分開進(jìn)行模擬，并計(jì)算了氣壁、液壁和氣液界面剪切應(yīng)力分布情況，發(fā)現(xiàn)模擬計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好。

BERTHELSEN等［46-47］借助 NEWTON等［44］的界面剪切分布函數(shù)，二維穩(wěn)態(tài)動(dòng)量方程采用有限差分格式復(fù)合求解，在管壁和氣液界面附近使用重疊網(wǎng)格技術(shù)進(jìn)行局部加密，近壁處采用 k-l一方程模型，遠(yuǎn)離壁面時(shí)采用 k-ε兩方程模型，并在氣液界面處采用 Level-Set界面重構(gòu)技術(shù)對(duì)水平管和微傾管分層流進(jìn)行了數(shù)值模擬，討論了剪切應(yīng)力值及其分布情況，模擬結(jié)果與湍流流動(dòng)的實(shí)驗(yàn)結(jié)果吻合較好。

GHORAI等［48］利用FLUENT軟件建立VOF模型，模擬了STRAND［33］和LOPEZ［41］的實(shí)驗(yàn)，使用STRAND得到的界面等效粗糙度對(duì)氣液界面進(jìn)行了修正。采用控制體積有限元法（CVFEM）求解控制方程，CVFEM的優(yōu)點(diǎn)在于既保留了有限單元法（FEM）網(wǎng)格劃分的靈活性，又未犧牲控制容積有限差分法（CVFDM）的魯棒性和經(jīng)濟(jì)性，網(wǎng)格劃分采用非結(jié)構(gòu)網(wǎng)格。通過數(shù)值計(jì)算得到了滿意的模擬效果，并分析了氣液界面摩擦因子和液壁摩擦因子的本構(gòu)關(guān)系，給出了界面摩擦因子與氣壁摩擦因子之比與表觀流速的關(guān)聯(lián)式。

DE SAMPAIO等［49］利用有限單元法和k-ω湍流模型對(duì)層-層分層流進(jìn)行了數(shù)值模擬，迭代過程中，將液位高度和壓降的估計(jì)值直接帶入控制方程，并用NEWTON-RAPHSON法封閉方程，即不斷循環(huán)修正液位高度和壓降值，直到滿足收斂精度。最后將計(jì)算的壁面剪切應(yīng)力、界面剪切應(yīng)力和液高預(yù)測(cè)與已有實(shí)驗(yàn)進(jìn)行了對(duì)比，認(rèn)為k-ω湍流模型能很好地模擬高含液率的光滑分層流。

SIDI-ALI等［36］的工作值得關(guān)注，利用FLUENT軟件，在HANSEN等［8］二維移動(dòng)壁面模型基礎(chǔ)上，建立了三維氣相的移動(dòng)壁面模型，即將氣液界面定義為一個(gè)無滑移的勻速移動(dòng)壁面，通過對(duì)比STRAND［33］實(shí)驗(yàn)的氣速分布，模擬得到了移動(dòng)壁面的移動(dòng)速度，利用粗糙速度通過有限容積法計(jì)算出界面摩擦因子和界面剪切應(yīng)力。不過他們的模型存在兩個(gè)明顯的缺陷；第一是移動(dòng)壁面的移動(dòng)速度是液相速度平均速度，這與實(shí)際情況不同；二是他們計(jì)算的氣壁剪切應(yīng)力明顯大于理論和實(shí)驗(yàn)結(jié)果（平均誤差大于 90%），因此他們并沒有進(jìn)一步討論氣相封閉關(guān)系。

4 5種現(xiàn)有模型的壓降與持液率預(yù)測(cè)對(duì)比分析

在對(duì)氣壁、液壁和氣液界面剪切應(yīng)力的預(yù)測(cè)基礎(chǔ)上，封閉氣液兩相動(dòng)量方程，從而達(dá)到通過已知條件（如管徑、氣液兩相物性值和氣液兩相表觀速度等）實(shí)現(xiàn)對(duì)持液率和壓降的預(yù)測(cè)，這就是對(duì)分層流封閉關(guān)系研究的主要意義，也是工程實(shí)踐最關(guān)心的問題。

4.1 氣壁摩擦因子關(guān)聯(lián)式

分層流氣壁剪切應(yīng)力的計(jì)算基礎(chǔ)是單相流氣壁剪切應(yīng)力的計(jì)算，早期氣液兩相流都采用單相流計(jì)算方法，由于分相流模型在分層流模型預(yù)測(cè)中的優(yōu)勢(shì)，這里只對(duì)兩相流的氣壁剪切作簡(jiǎn)單介紹。TAITEL等［6-7］分層流氣壁剪切應(yīng)力由BLASIUS形式公式計(jì)算仍然有很高精度，說明分層流氣壁剪切與單相流壁面剪切非常相似。表1給出了4個(gè)具有代表性的分層流氣壁摩擦因子計(jì)算公式。

4.2 液壁摩擦因子關(guān)聯(lián)式

BLASIUS公式和TAITEL & DUKLER公式在較低的雷諾數(shù)或者光滑分層流下，預(yù)測(cè)的液壁剪切應(yīng)力誤差不大，但隨著雷諾數(shù)的增加，上述兩個(gè)關(guān)聯(lián)式誤差較大，不能滿足預(yù)測(cè)精度。KOWALSKI的液壁摩擦因子關(guān)聯(lián)式考慮了持液率關(guān)系，被更多研究者所使用。4個(gè)湍-湍分層流液壁摩擦因子關(guān)聯(lián)式見表2。

表1 分層流氣壁摩擦因子關(guān)聯(lián)式

表2 分層流液壁摩擦因子關(guān)聯(lián)式

4.3 界面摩擦因子關(guān)聯(lián)式

早期TAITEL等認(rèn)為界面摩擦因子等于氣壁摩擦因子，后來的研究發(fā)現(xiàn)這一假設(shè)只在雷諾數(shù)較小的湍流和層流中才能較好地滿足實(shí)際結(jié)果［8，52］。對(duì)于界面剪切摩擦因子關(guān)聯(lián)式，大致分以下幾種典型形式:

（1）常數(shù)，如SHOHAM & TAITEL［8］；

（2）BLASIUS形式，如 TAITEL 等［6-7］，SIDI-ALI等［36］；

（3）KOWALSKI形式，如KOWALSKI［30］，NEWTON等［53］；

（4）fi/fG經(jīng)驗(yàn)公式，如ANDRITSOS等［14］，TZOTZI等［32］。

表3給出了7個(gè)典型的湍-湍分層流界面摩擦因子關(guān)聯(lián)式。TAITEL & DUKLER模型認(rèn)為氣壁摩擦因子等于界面摩擦因子；SHOHAM等［8］沿用了COHEN & HANRATTY矩形管（0.305m×0.0254 m）分層流實(shí)驗(yàn)結(jié)果，介質(zhì)為空氣和水，在波狀分層流下，平均界面摩擦因子為0.0142；KOWALSKI對(duì)空氣-水和氟利昂12-水系統(tǒng)進(jìn)行實(shí)驗(yàn)測(cè)量，管徑0.051m，管口壓力分別為2.25kPa和4.2kPa，在氣壁和液壁的摩擦因子關(guān)聯(lián)式與動(dòng)力方程的本構(gòu)關(guān)系基礎(chǔ)上，提出了KOWALSKI形式的界面摩擦因子關(guān)聯(lián)式；HART等使用ECK［54］單相流壁面剪切應(yīng)力公式，將氣液界面假設(shè)為粗糙壁面；CROWLEY等［55］在研究微傾管流型轉(zhuǎn)變時(shí)發(fā)現(xiàn)，當(dāng)界面摩擦因子使用5倍或10倍由TAITEL & DUKLER公式計(jì)算的氣壁摩擦因子時(shí)，預(yù)測(cè)的流型與實(shí)驗(yàn)結(jié)果能很好地吻合；SIDI-ALI等［36］通過三維移動(dòng)壁面的CFD模型得到的移動(dòng)壁面摩擦因子關(guān)聯(lián)式；TZOTZI等以空氣、二氧化碳、氦氣作為氣相，水、丁醇-水溶液作為液相，通過波動(dòng)理論完善從光滑流到2D波、從2D波到大幅度波的轉(zhuǎn)型。

表3 分層流界面摩擦因子關(guān)聯(lián)式

4.4 不同模型預(yù)測(cè)結(jié)果對(duì)比

選擇5個(gè)具有代表意義的模型進(jìn)行對(duì)比分析，它們是TAITEL & DUKLER［6-7］（T-D）、SHOHAM & TAITEL［8］（S-T）、KOWALSKI［30］（K）、SIDI-ALI & GATIGNOL［36］（ SA-G） 和TZOTZI & ANDRITSOS［32］（T-A），其中，SA-G模型為半封閉模型，且氣壁摩擦因子計(jì)算誤差很大，沒有達(dá)到改進(jìn)T-D公式的目的，而并沒有考慮液壁摩擦因子，因此，這里對(duì)比計(jì)算中，SA-G模型的氣壁和液壁摩擦因子采用TAITEL & DUKLER公式；T-A的氣壁摩擦因子直接使用TAITEL & DUKLER公式。實(shí)驗(yàn)數(shù)據(jù)取自ANDRITSOS［14］（D=0.09525m）、ANDRITSOS［56］（D=0.08m）、STRAND［33］（D=0.1m）和AYATI［57］（D=0.1m）等，介質(zhì)為空氣和水，表觀雷諾數(shù)范圍滿足 8000≤ReSG≤90000和 5000≤ReSL≤17000。圖2為5種不同模型壓降和持液率預(yù)測(cè)結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比。

從圖2（a）可以看出，T-D模型和T-A模型預(yù)測(cè)的持液率偏大，前者的誤差主要來自界面剪切的過小預(yù)測(cè)，后者的誤差來自對(duì)液壁剪切的過大預(yù)測(cè)；K模型預(yù)測(cè)效果較好，但預(yù)測(cè)值分布的總體趨勢(shì)與實(shí)驗(yàn)測(cè)量結(jié)果不吻合，適用范圍有限；S-T模型和由T-D壁面摩擦因子關(guān)聯(lián)式封閉的SA-G模型在很大程度上修正了T-D模型的假設(shè)造成的誤差，這一結(jié)果從界面剪切應(yīng)力計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合得很好可以得到解釋。

從圖2（b）可以看出，T-D模型和T-A模型預(yù)測(cè)的壓降誤差最大，前者的誤差主要來自界面剪切的過小預(yù)測(cè)，后者的誤差來自液壁剪切的過大預(yù)測(cè)，這與持液率預(yù)測(cè)結(jié)果相對(duì)應(yīng)；S-T模型預(yù)測(cè)的壓降誤差較小，說明fi=0.0142在所討論的雷諾數(shù)流動(dòng)中能保證預(yù)測(cè)結(jié)果誤差在20%以內(nèi)，但由T-D壁面摩擦因子關(guān)聯(lián)式封閉的SA-G模型誤差卻很大，這與SA-G界面摩擦因子預(yù)測(cè)公式在較大雷諾數(shù)下預(yù)測(cè)值明顯偏低有關(guān)；K模型預(yù)測(cè)效果很好，可見基于KOWALSKI形式的封閉關(guān)聯(lián)式在壓降預(yù)測(cè)方面具有較高精度。

5 結(jié)語(yǔ)與展望

（1）理論模型、實(shí)驗(yàn)測(cè)量、CFD數(shù)值計(jì)算 3類不同研究方法都可以較好地預(yù)測(cè)分層流界面剪切性質(zhì)，進(jìn)而預(yù)測(cè)分層流的持液率和壓降。其中理論模型為實(shí)驗(yàn)和CFD方法提供了理論依據(jù)，實(shí)驗(yàn)結(jié)果驗(yàn)證并修正了理論模型，而CFD數(shù)值計(jì)算可以通過合理假設(shè)，在較少引入已有本構(gòu)關(guān)聯(lián)式的前提下預(yù)測(cè)分層流，并提供更詳細(xì)的流場(chǎng)細(xì)節(jié)。但3類研究方法都具有一定的局限性，理論模型受限于模型的簡(jiǎn)化條件和方程的封閉性，實(shí)驗(yàn)間接測(cè)量受限于實(shí)驗(yàn)測(cè)量技術(shù)與數(shù)據(jù)分析，CFD則受限于計(jì)算技術(shù)和方法。

圖2 5種模型預(yù)測(cè)結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比

（2）理論模型取得了一定的成果，并在實(shí)驗(yàn)和CFD仿真計(jì)算中起到了關(guān)鍵的指導(dǎo)作用。然而理論模型受到兩方面限制:一是理論模型的簡(jiǎn)化假設(shè)條件，這些簡(jiǎn)化導(dǎo)致在實(shí)際應(yīng)用中適用范圍較窄；二是理論模型的封閉性。對(duì)于絕大多數(shù)的理論模型，都需要額外提供本構(gòu)關(guān)系式來封閉求解方程，這些本構(gòu)關(guān)系通常是基于實(shí)驗(yàn)的經(jīng)驗(yàn)關(guān)系式。如何建立合理的分層流理論模型，是分層流界面剪切預(yù)測(cè)的關(guān)鍵。應(yīng)根據(jù)工程普遍問題提出合理假設(shè)，并針對(duì)界面建立更精確的局部模型，使模型預(yù)測(cè)具有更大的適應(yīng)性和準(zhǔn)確性；同時(shí)大力開展有效的實(shí)驗(yàn)工作，得到更好的符合實(shí)際問題的經(jīng)驗(yàn)關(guān)系式，使理論模型的求解更接近實(shí)際情況。

（3）由于實(shí)驗(yàn)條件在很大程度上影響了氣液兩相分層流界面剪切預(yù)測(cè)的結(jié)論，因此很有必要在實(shí)驗(yàn)過程中考慮更多的影響因素，但這樣也會(huì)使研究難度無形加大。此外，實(shí)驗(yàn)室中的大部分實(shí)驗(yàn)與實(shí)際工況環(huán)境相差較大，因而會(huì)產(chǎn)生實(shí)驗(yàn)數(shù)據(jù)與實(shí)際工況不同程度的誤差。因此，如何在實(shí)驗(yàn)室中得到與實(shí)際工況（特別是高溫高壓等特殊工況）相符的研究結(jié)論，仍需要做更多的基礎(chǔ)工作。

（4）分層流流場(chǎng)細(xì)節(jié)和界面行為直接影響理論模型建立和CFD數(shù)值計(jì)算結(jié)果的正確性。應(yīng)完善更加細(xì)致的分層流流場(chǎng)細(xì)節(jié)和界面行為的實(shí)驗(yàn)觀察，特別是對(duì)界面波動(dòng)與氣泡/液滴的觀察，為理論模型的合理建立和 CFD數(shù)值計(jì)算的正確求解提供更詳細(xì)的科學(xué)數(shù)據(jù)。

（5）CFD數(shù)值計(jì)算仿真技術(shù)可以考慮被理論模型簡(jiǎn)化的條件，也可以觀測(cè)到實(shí)驗(yàn)難以觀測(cè)的流程細(xì)節(jié)特征。但是，CFD模型的建立離不開理論和實(shí)驗(yàn)得到的已知結(jié)論，并且，計(jì)算技術(shù)直接影響到計(jì)算結(jié)果，對(duì)于界面行為的模擬一直是CFD研究的難點(diǎn)和重點(diǎn)。隨著CFD技術(shù)和方法的不斷發(fā)展，在更精確的流場(chǎng)特征前提下，氣液兩相分層流界面行為模擬和界面剪切應(yīng)力的預(yù)測(cè)是未來 CFD數(shù)值計(jì)算最主要的挑戰(zhàn)，也是大勢(shì)所趨。

參 考 文 獻(xiàn)

［1］ ULLMANN A，BRAUNER N. Closure relations for two-fluid models for two-phase stratified smooth and stratified wavy flows［J］. International Journal of Multiphase Flow，2006，32（1）:82-105.

［2］ 郭烈錦. 兩相與多相流動(dòng)力學(xué)［M］. 西安:西安交通大學(xué)出版社，2002.

［3］ MOUZA A A，PARAS S V，KARABELAS A J. CFD code application to wavy stratified gas-liquid flow［J］. Chemical Engineering Research and Design，2001，79（5）:561-568.

［4］ LINé A，LOPEZ D. Two-fluid model of wavy separated two-phase flow［J］. International Journal of Multiphase Flow，1997，23（6）: 1131-1146.

［5］ LOCKHART R W，MARTINELLI R C. Proposed correlation of data for isothermal two phase，two component flow in pipes［J］. Chemical Engineering Progress，1949，45（1）:38-48

［6］ TAITEL Y，DUKLER A E. A model for predicting flow regime transitions in horizontal and near horizontal gas-liquid flow［J］. AIChE Journal，1976，22（1）:47-55.

［7］ TAITEL Y，DUKLER A E. A theoretical approach to the Lockhart-Martinelli correlation for stratified flow［J］. International Journal of Multiphase Flow，1976，2（4）:591-595.

［8］ SHOHAM O，TAITEL Y. Stratified turbulent-turbulent gas-liquid flow in horizontal and inclined pipes［J］. AlChE Journal，1984，303 （3）:377-385.

［9］ HANSEN E A，VESTED H J. Liquid hold-up，pressure drop，and velocity profiles in steady uniform stratified flow［J］. Journal of Energy Resources Technology，1991，113（2）:87-93.

［10］ ULLMANN A，ZAMIR M，GAT S，et al. Multi-holdups in co-current stratified flow in inclined tubes［J］. International Journal of Multiphase Flow，2003，29（10）:1565-1581.

［11］ ULLMANN A，GOLDSTEIN A，ZAMIR M，et al. Closure relations for the shear stresses in two-fluid models for laminar stratified flow［J］. International Journal of Multiphase Flow，2004，30（7/8）:877-900.

［12］ KANG H C，KIM M H. The relation between the interfacial shear stress and the wave motion in a stratified flow［J］. International Journal of Multiphase Flow，1993，19（1）:35- 49.

［13］ 劉夷平，王茹金，陳超，等.氣液兩相分層流剪切應(yīng)力的不確定度分析［J］. 實(shí)驗(yàn)流體力學(xué)，2012，26（4）:43-47.

［14］ ANDRITSOS N，HANRATTY T J. Influence of interfacial waves in stratified gas-liquid flows［J］. International Journal of Multiphase Flow，1987，13（5）:583-603.

［15］ HART J，HAMERSMA P J，F(xiàn)ORTUIN J M H. Correlations predicting frictional pressure drop and liquid holdup during horizontal gas-liquid pipe flow with a small liquid holdup［J］. International Journal of Multiphase Flow，1989，15（6）:947-964.

［16］ CHEN X T，CAI X D，BRILL J P. Gas-liquid stratified-wavy flow in horizontal pipelines［J］. Journal of Energy Resources Technology，1997，119:209-216.

［17］ VLACHOS N A，PARAS S V，KARABELAS A J. Prediction of holdup，axial pressure gradient and wall shear stress in wavy stratified and stratified/atomization gas/liquid flow［J］. International Journal of Multiphase Flow，1999，25（2）:365-376.

［18］ VLACHOS N A. Studies of stratified/atomization two-phase flow in a horizontal pipe［D］. Thessaloniki:Aristotle University of Thessaloniki，1997

［19］ NG T S，LAWRENCE C J，HEWITT G F. Interface shapes for two-phase laminar stratified flow in a circular pipe［J］. International Journal of Multiphase Flow，2000，27（7）:1301-1311.

［20］ BRAUNER N，ROVINSKY J，MOALEM M D. Determination of the interface curvature in stratified two-phase systems by energy considerations［J］. International Journal of Multiphase Flow，1996，22（6）:1167-1185.

［21］ BRAUNER N，ROVINSKY J，MOALEM M D. Analytical solution for Laminar-Laminar two-phase flow in circular conduits［J］. Chemical Engineering Communication，1996，141/142（1）:103-143.

［22］ BRAUNER N，ROVINSKY J，MOALEM M D. Analytical solution of laminar stratified flow with curved interfaces［C］// Proceedings of the NURETH-7 Meeting，ANS 1，1995:192-211.

［23］ BRAUNER N，ROVINSKY J，MOALEM M D. Characteristics of annular and stratified two-phase flow in the limit of fully eccentric core annular configuration［C］// Proc. of the ExHFT 4，Brussels，1997.

［24］ BRAUNER N，MOALEM M D，ROVINSKY J. A two-fluid model for stratified flows with curved interfaces［J］. International Journal of Multiphase Flow，1998，24（6）:975-1004.

［25］ ROVINSKY J，BRAUNER N，MOALEM M D. Analytical solution for laminar two-phase flow in a fully eccentric core-annular configuration［J］. International Journal of Multiphase Flow，1997，23 （3）:523-543.

［26］ SOLEIMANI A. Phase distribution and associated phenomena in oil-water flows in horizontal tubes［D］. London:Imperial College，University of London，1999.

［27］ BIRVALSKI M，TUMMERS M J，DELFOS R，et al. PIV measurements of waves and turbulence in stratified horizontal two-phase pipe flow［J］. International Journal of Multiphase Flow，2014，62:161-173.

［28］ FULLMER W D，RANSOM V H，LOPEZ de Bertodano M A. Linear and nonlinear analysis of an unstable， but well-posed，one-dimensional two-fluid model for two-phase flow based on the inviscid Kelvin-Helmholtz instability［J］. Nuclear Engineering and Design，2014，268:173-184.

［29］ LIU G，WANG Y S，ZANG G J，et al. Viscous Kelvin- Helmholtz instability analysis of liquid-vapor two-phase stratified flow for condensation in horizontal tubes［J］. International Journal of Heat and Mass Transfer，2015，84:592-599.

［30］ KOWALSKI J E. Wall and interfacial shear stress in stratified flow in a horizontal pipe［J］. AIChE Journal，1987，32（2）:274-281.

［31］ HAGIWARA Y，ESMAEILZADEH E，TSUTSUI H，et al. Simultaneous measurements of liquid film thickness，wall shear stress and gas flow turbulence of horizontal wavy two-phase flow ［J］. International Journal of Multiphase Flow，1989，15（3）:421- 431.

［32］ TZOTZI C，ANDRITSOS N. Interfacial shear stress in wavy stratifed gas-liquid flow in horizontal pipes［J］. International Journal of Multiphase Flow，2013，54:43-54.

［33］ STRAND ?. An experimental investigation of stratified two-phase flow in horizontal pipes［D］. Oslo:University of Oslo，1993.

［34］ HAALAND S E. Simple and explicit formulas for the friction factor in turbulent pipe flow［J］. Journal of Fluids Engineering，1983，105: 89-90.

［35］ LI W D，SUN K X，ZHOU F D. Interfacial shear stress of stratified flow in a horizontal pipe［J］. Chinese Journal of Chemical Engineering，1999，7（3）:263-270.

［36］ SIDI-ALI K，GATIGNOL R. Interfacial friction factor determination using CFD simulations in a horizontal stratified two-phase flow［J］. Chemical Engineering Science，2010，65（18）:5160- 5169.

［37］ DOMINIQUE B. The difficult challenge of a two-phase CFD modelling for all flow regimes［J］. Nuclear Engineering and Design，2014，279:116-125.

［38］ AKAI M，INOUE A，AOKI S. The prediction of stratified two-phase flow with a two equation model of turbulence［J］. International Journal of Multiphase Flow，1981，7（1）:21-39.

［39］ ISSA R I. Prediction of turbulent，stratified，two-phase flow in inclined pipes and channels［J］. International Journal of Multiphase Flow，1988，14（2）:141-154.

［40］ MEKNASSI F，BENKIRANE R，LINé A，et al. Numerical modelling of wavy stratified two-phase flow in pipes［J］. Chemical Engineering Science，2000，55（20）:4681-4697.

［41］ LOPEZ D. Ecoulements diphasiques a phases separees a faible contenu de liquide［D］. France:I.N.P. Toulouse，1994.

［42］ NEWTON C H，BEHNIA M. On the use of the stratified momentum balance for the deduction of shear stress in horizontal gas-liquid pipe flow［J］. International Journal of Multiphase Flow，1998，24（8）: 1407-1423.

［43］ NEWTON C H，BEHNIA M. Numerical calculation of turbulent stratified-gas-liquid pipe flows［J］. International Journal of Multiphase Flow，2000，26（2）:327-337.

［44］ NEWTON C H，BEHNIA M. A numerical model of stratified-wavy gas-liquid pipe flow［J］.Chemical Engineering Science，2001，56 （24）:6851-6861.

［45］ VLACHOS N A，PARAS S V，KARABELAS A J. Liquid-to-wall shear stress distribution in stratified/atomization flow［J］. International Journal of Multiphase Flow，1997，23（5）:845-863.

［46］ BERTHELSEN P A，YTREHUS T. Calculations of stratified wavy two-phase flow in pipes［J］. International Journal of Multiphase Flow，2005，31（5）:571-592.

［47］ BERTHELSEN P A，YTREHUS T. Stratified smooth two-phase flowusing the immersed interface method［J］. Computers & Fluids，2007，36（7）:1273-1289.

［48］ GHORAI S，NIGAM K D P. CFD modelling of flow profiles and interfacial phenomena in two-phase flow in pipes［J］. Chemical Engineering and Processing，2006，45（1）:55-65.

［49］ DE SAMPAIO P A B，F(xiàn)ACCINI J L H，SU Jian. Modelling of stratified gas-liquid two-phase flow in horizontal circular pipes［J］. International Journal of Heat and Mass Transfer，2008，51: 2752-2761.

［50］ HAND N P. Gas-liquid co-current flow in a horizontal pipeline［D］. Belfast:The Queen's University of Belfast，1991.

［51］ SRICHAI S. High pressure separated two-phase flow［D］. London: Imperial College，1994.

［52］ BARNEA D，TAITEL Y. Structural and interfacial stability of multiple solutions for stratified flow［J］. International Journal of Multiphase Flow，1992，18（6）:821-830.

［53］ NEWTON C H，BEHNIA M. On the use of the stratified momentum balance for the deduction of shear stress in horizontal gas-liquid pipe flow［J］. International Journal of Multiphase Flow，1998，24（8）: 1407-1423.

［54］ ECK B. Technische Stromunglehre［M］. New York:Springer，1973.

［55］ CROWLEY C J，WALLIS G B，BARRY J J. Validation of a one-dimensional wave model for the stratified-to-slug flow regime transition，with consequences for wave growth and slug frequency［J］. International Journal of Multiphase Flow，1992，18（2）:249-271.

［56］ ANDRITSOS N. Statistical analysis of waves in horizontal stratified gas-liquid flow［J］. International Journal of Multiphase Flow，1992，18（3）:465-473.

［57］ AYATI A J，KOLAAS A，JENSEN G W，et al. A PIV investigation of stratified gas-liquid flow in a horizontal pipe［J］. International Journal of Multiphase Flow，2014，61:129- 143.

Progress of interfacial shear stress in the gas-liquid stratified flow in horizontal pipes

ZHENG Ping，ZHAO Liang
（College of Petroleum Engineering，Liaoning Shihua University，F(xiàn)ushun 113001，Liaoning，China）

Abstract:The gas-liquid two-phase stratified flow has complex momentum and energy transmission phenomena at interface，although its pattern is simple. There are still no unified conclusions for interfacial shear stress in stratified flow. Current progress of interfacial shear stress in horizontal pipes were elaborated by theoretical models，experimental models and numerical simulations. As for theoretical models，closure models were established by models simplification and empirical correlations. Empirical correlations were corrected with closure relations in experiments. Due to simplifying assumptions and experimental conditions，there are some limitations in predicting the interfacial shear stress with theoretical and experimental ways. The detailed flow fields were studied by numerical simulations，but closure relations have been less obtained so far. Five existing models are compared according to liquid holdup and pressure drop. The future research trends of gas-liquid two-phase stratified flow interfacial shear stress in horizontal pipes were further discussed. More detailed local models are needed to be proposed and engineering practice should be taken into account. It is necessary to develop new methods of gas-liquid interface calculations and closure relations for numerical simulations.

Key words:horizontal pipe；gas-liquid flow；stratified flow；interfacial shear stress；liquid holdup；pressure drop

中圖分類號(hào)：O 359.1

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1000-6613（2016）05-1289-09

DOI：10.16085/j.issn.1000-6613.2016.05.004

收稿日期：2015-11-09；修改稿日期：2016-02-01。

基金項(xiàng)目：國(guó)家自然科學(xué)基金項(xiàng)目（51201009）。