• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    湍流系統(tǒng)的能量最小多尺度模型研究進展

    2022-07-06 08:14:14王利民郭舒宇向星付少童
    化工學(xué)報 2022年6期
    關(guān)鍵詞:大渦層流黏性

    王利民,郭舒宇,向星,付少童

    (1 中國科學(xué)院過程工程研究所多相復(fù)雜系統(tǒng)國家重點實驗室,北京 100190;2中國科學(xué)院大學(xué)化工學(xué)院,北京 100049)

    引 言

    過程工程中存在著不同尺度的湍流,其對物質(zhì)傳遞與反應(yīng)效率發(fā)揮著重要作用。自雷諾實驗以來,湍流科學(xué)研究已有一個多世紀的歷史[1],但其仍為經(jīng)典物理學(xué)中尚未解決的主要難題之一[2]。Reynolds將充分發(fā)展的湍流運動分解成時均運動和脈動運動,在對Navier-Stokes 方程進行雷諾分解時,產(chǎn)生了未知的雷諾應(yīng)力項,造成了雷諾時均運動方程不封閉的根本性困難[3]。為了封閉雷諾時均運動方程,需要對雷諾應(yīng)力項進行?;;姆绞娇梢远喾N多樣,構(gòu)造了各種湍流模型。

    當(dāng)前文獻報道的湍流模型已達數(shù)百種。普遍認為,湍流模型的關(guān)鍵在于模型封閉的合理性。傳統(tǒng)的湍流模型可分為兩類[4]。第一類是渦黏法,假設(shè)湍流脈動對平均流有耗散作用,Boussinesq 假設(shè)下的雷諾應(yīng)力在數(shù)學(xué)上類似于牛頓流體的應(yīng)力-應(yīng)變關(guān)系,并與平均應(yīng)變率成正比。典型的渦黏模型主要有零方程模型[5]、一方程模型[6]和兩方程模型[7-9]。渦黏性模型中,與分子理論中的平均自由程類似,Prandtl率先提出了混合長度模型[10],通過引入湍流特征長度來模化雷諾應(yīng)力,但該模型存在很多不足,在模型中特征長度的取值根據(jù)簡單工況的實驗確定,在應(yīng)用到比較復(fù)雜工況時其適用性依然存疑;湍流場的尺度范圍很大,用單一長度尺度來描述雷諾應(yīng)力不合適。隨后Taylor 的渦量輸運理論、von Karman的相似性理論也通過假設(shè)雷諾應(yīng)力和時均速度梯度之間的關(guān)系來封閉雷諾時均(Reynoldsaveraged Navier-Stokes,RANS)方程[11]。隨著湍流建模和模擬技術(shù)的發(fā)展,這些唯象湍流模型逐漸被更復(fù)雜的模型所取代。由于湍流過程中有序和無序共存的復(fù)雜性,目前還沒有一個公認的普適湍流模型。另一類是雷諾應(yīng)力模型(Reynolds stress model,RSM)[12]或二階矩閉合模型,該模型考慮了雷諾應(yīng)力的各向異性和復(fù)雜的湍流相互作用,比渦黏模型更詳細、更通用。Chou[13]和Rotta[14]在雷諾應(yīng)力模型做出了開創(chuàng)性工作,直接從雷諾應(yīng)力輸運方程或其簡化代數(shù)方程計算Navier-Stokes 方程中的雷諾應(yīng)力張量,從而摒棄了湍流黏性系數(shù)的概念。然而這種建模方法中雷諾應(yīng)力存在擴散輸運、湍流壓力-應(yīng)變相互作用,以及產(chǎn)生項和耗散項等,仍存在封閉問題[15]。此外,由于雷諾應(yīng)力的計算需要另外求解六個輸運方程,使得RSM 的計算成本非常高。因此,與直接數(shù)值模擬(direct numerical simulation,DNS)[16]和大渦模擬(large eddy simulation, LES)[17]一樣,RSM 產(chǎn)生的信息量遠超工程中計算平均流動的所需,不具備成本效益[18]。

    各種層次的湍流模擬方法如圖1 所示,除DNS外均需要構(gòu)建湍流模型。湍流模型是影響工程湍流模擬準(zhǔn)確性的關(guān)鍵因素之一,在計算流體力學(xué)(computational fluid dynamics,CFD)模擬中起著關(guān)鍵作用。當(dāng)前各種湍流模型形式千差萬別,但封閉雷諾時均方程的方式仍是經(jīng)驗性的,缺乏對復(fù)雜流動的各種相互作用的競爭中協(xié)調(diào)及控制機制的深入研究。而介科學(xué)理論致力于闡明復(fù)雜系統(tǒng)中各種控制機制的相互作用,本文基于介科學(xué)框架,以能量最小多尺度(energy-minimization multi-scale,EMMS)的思想為核心,介紹了湍流中介尺度行為的共性原理,包括湍流中的黏性控制機制與慣性控制機制,以及二者在流動中的競爭中協(xié)調(diào)以實現(xiàn)穩(wěn)定性條件。基于該思想構(gòu)建了具備工程應(yīng)用價值的EMMS 湍流模型,顯著改進了RANS 模擬預(yù)測的精度,在層湍轉(zhuǎn)捩預(yù)測以及全球氣候仿真等具備重要應(yīng)用價值,為介科學(xué)理論作為復(fù)雜系統(tǒng)的普適理論奠定了基礎(chǔ)。

    圖1 湍流模擬的通用多尺度框架Fig.1 A generic multiscale framework for modelling turbulent flows

    1 EMMS建模思路

    20 世紀80~90 年代,李靜海等[19-21]從氣固流態(tài)化中的顆粒聚團現(xiàn)象出發(fā),認為:氣體和顆粒都擁有自己獨立的運動趨勢,即氣體趨向于選擇阻力最小的途徑流動,而顆??偸潜M可能處于最小位能的位置,介尺度聚團的形成源于氣體和顆粒的各自運動趨勢在競爭中的協(xié)調(diào)。李靜海等原創(chuàng)性地提出了基于顆粒尺度、顆粒聚團尺度和設(shè)備尺度的多尺度分析方法和介尺度聚團結(jié)構(gòu)應(yīng)滿足的穩(wěn)定性條件,將該穩(wěn)定性條件與氣體、固體的守恒方程和團聚物尺寸方程同時求解,從而得到相應(yīng)的結(jié)構(gòu)參數(shù),建立了具有普遍意義的EMMS 模型[19-21]。EMMS 模型成功實現(xiàn)了非均勻氣固系統(tǒng)的準(zhǔn)確模擬,顯著提升了預(yù)測的準(zhǔn)確性和精度[22-23]。邱小平等[24]耦合了EMMS 曳力與簡化雙流體模型,實現(xiàn)了工業(yè)規(guī)模氣固反應(yīng)器快速模擬的潛力。陳愷成等[25]開展了基于EMMS 的循環(huán)流化床流域研究,確定了快速床的操作邊界。佟穎等[26]基于雙分散顆粒EMMS 曳力模型,模擬了兩種不同的雙分散鼓泡流化床,探討了模型的適用性。

    EMMS 原理成功推廣到其他復(fù)雜系統(tǒng)[27-28]。在氣液體系中,應(yīng)用EMMS模型能量分解的基本思路,Ge等[28-29]對氣液體系微觀和介觀能耗過程進行了區(qū)分和定義,提出了氣液體系的穩(wěn)定性條件。在此基礎(chǔ)上,Yang 等[30]假設(shè)鼓泡塔體系中存在兩種氣泡特征直徑,以反映體系內(nèi)的非均勻結(jié)構(gòu),提出了雙氣泡(dual-bubble-size, DBS)模型,從理論上預(yù)測了鼓泡塔中的流型轉(zhuǎn)變,在結(jié)構(gòu)參數(shù)范圍內(nèi),全局微尺度能耗最小值對應(yīng)的解會發(fā)生跳躍,從物理上解釋了流動結(jié)構(gòu)的宏觀演化。Han 等[31]根據(jù)介區(qū)域分析了主導(dǎo)機制,并基于介尺度框架研究了鼓泡塔中的流域轉(zhuǎn)變。陳衛(wèi)等[32]將流態(tài)化過程中鼓泡、湍動和快速流化過程類比物質(zhì)的固、液和氣三態(tài),嘗試基于EMMS 原理構(gòu)建相變理論。馬永麗等[33]對氣液固流化床內(nèi)復(fù)雜三相流動結(jié)構(gòu)的介尺度機理模型進行了分析和總結(jié)。汪帆等[34]基于介尺度研究范式及理論,探討了利用EMMS 模型對現(xiàn)有成核數(shù)學(xué)模型進行修正和優(yōu)化的思路。

    上述研究的非平衡復(fù)雜系統(tǒng)都遵循控制機制之間的競爭中協(xié)調(diào),文獻中的大量模擬和實驗結(jié)果都證實了EMMS原理對所研究非平衡復(fù)雜系統(tǒng)的適用性[35-37],表明了EMMS 模型背后蘊含的原理具有一定的普遍性。EMMS建模思路從理解復(fù)雜系統(tǒng)的多尺度結(jié)構(gòu)和不同機制之間的相互作用和控制機制開始,通過分析各控制機制的競爭中協(xié)調(diào)獲得系統(tǒng)的穩(wěn)定性條件,在數(shù)學(xué)上可表述為多目標(biāo)變分問題。因此,通過相應(yīng)的穩(wěn)定性條件和復(fù)雜系統(tǒng)的封閉模型,可以實現(xiàn)不同尺度下動力學(xué)方程的關(guān)聯(lián)[38-39]。經(jīng)過三十多年的努力,EMMS思路已經(jīng)成為解決復(fù)雜系統(tǒng)問題的一個非常有前途的工具,這一原創(chuàng)性學(xué)術(shù)思想形成了創(chuàng)新性的EMMS原理與介科學(xué)理論體系[40]。

    2 湍流控制機制與穩(wěn)定性條件

    雷諾實驗通過臨界Reynolds數(shù)來區(qū)分流動的層流和湍流狀態(tài),率先開展了流動穩(wěn)定性的研究。Heisenberg[41]將Orr-Sommerfeld 方程用于分析平面Poiseuille 流的穩(wěn)定性,發(fā)現(xiàn)在接近湍流起始流速的有界剪切流中會出現(xiàn)Tollmien-Schlichting 波,該方法也被用于解決流動穩(wěn)定性問題[42]。Lorenz[43]發(fā)現(xiàn),流體確定性模型的所有非周期解均有界但存在不規(guī)則的波動,這與控制方程非線性導(dǎo)致的無窮多序列分叉和擾動有關(guān)。Helmholtz[44]提出了恒定外力作用下流體緩慢穩(wěn)定流動的最小能量耗散率原理。Korteweg[45]發(fā)現(xiàn)層流狀態(tài)下流體的流動分布使得黏性的能量損失最小。Reyleigh[46]將最小能量耗散率原理擴展到許多其他流體流動系統(tǒng)。根據(jù)非平衡熱力學(xué)理論,Bejan[47]提出了層流中最小熵產(chǎn)生原理,這與最小能量耗散率原理一致。

    一般認為,在層流中,系統(tǒng)完全由流體黏性決定,速度分布受黏性耗散率最小的條件控制。而在單相湍流中,黏性控制機制和慣性控制機制共同作用,比層流中黏性作用主導(dǎo)的情況要復(fù)雜很多,因此Helmholtz[44]提出的最小能量耗散率原理不再適用于湍流。

    目前,基于各種最大化流動耗散函數(shù),提出了幾種變分準(zhǔn)則來預(yù)測單相湍流中的平均速度場。Malkus[48]提出了最大總黏性耗散率D的變分準(zhǔn)則,以預(yù)測湍流管流中的軸向速度分布。隨后,該準(zhǔn)則被修正為最大“效率”E=DDf/Dm(這里Dm是平均流的黏性耗散率,Df=D-Dm是脈動的耗散率),由最大“效率”E預(yù)測的最佳速度場比最大總黏性耗散率D得出的結(jié)果更接近實際流動[49]。Bertram[50]提出了一種新的最小湍動能變分判據(jù),并試圖闡明不同變分判據(jù)在不同約束條件下的聯(lián)系。自Malkus 的最大總黏性耗散率D變分準(zhǔn)則出現(xiàn)以來,湍流系統(tǒng)中出現(xiàn)了各種變分原理,但均無法解釋為什么一些變分原理可以產(chǎn)生實際結(jié)果。此外,對于給定最大/最小化的特定流動量缺乏明確的物理基礎(chǔ)。

    由于單相湍流系統(tǒng)中存在黏性和慣性兩種不同的控制機制,變分原理不能僅用黏性控制機制的極值趨勢或慣性控制機制的極值趨勢來表示;相反,這兩種控制機制的競爭中協(xié)調(diào)對系統(tǒng)穩(wěn)定起著重要作用。如果認為單相湍流由許多湍流渦組成,且渦與周圍環(huán)境的相互作用類似于顆粒和氣體之間的相互作用,則湍流系統(tǒng)和氣固系統(tǒng)的耗散行為類似[27]。黏性和慣性控制機制都不能完全支配一個系統(tǒng),因此它們不能完全實現(xiàn)各自的趨勢;相反,它們通過競爭中協(xié)調(diào)來共同實現(xiàn)系統(tǒng)穩(wěn)定。Li等[27]發(fā)現(xiàn)了黏性控制機制和慣性控制機制之間的競爭中協(xié)調(diào),類似于氣固系統(tǒng)中顆粒主導(dǎo)機制和氣體主導(dǎo)機制之間的競爭中協(xié)調(diào),提出黏性耗散率Wv最小和總耗散率WT最大,即

    利用關(guān)系式(1)作為穩(wěn)定性條件成功預(yù)測了圓管湍流的徑向速度分布。Wang 等[28,51-52]發(fā)現(xiàn),湍流中黏性力和慣性力分別反映了速度場在空間和時間上的非均勻性。慣性力和黏性力同時作用于流體微元,流體微元提供載體,黏性和慣性在流體微元內(nèi)競爭中協(xié)調(diào),其中慣性機制和黏性機制均與耗散有關(guān)。將總能量耗散分解為平均剪切耗散和湍流脈動耗散:

    式(2)右側(cè)第一項為平均剪切耗散,定義為Wv;第二項為湍流脈動耗散,定義為Wte。假設(shè)Wv與完全層流中能量耗散的本質(zhì)相同,由黏性占優(yōu),則根據(jù)最小黏性耗散原理,與時均速度場相關(guān)的平均剪切耗散Wv應(yīng)趨于最小。

    當(dāng)流動處于穩(wěn)定狀態(tài)時,流體微元受到的慣性力為:

    即表明流體的脈動決定其慣性運動,湍流脈動越劇烈則慣性作用越強。湍流慣性機制與湍流脈動的耗散存在密切關(guān)聯(lián)。根據(jù)Li等[27]認為總耗散率WT趨于最大,由于黏性作用使黏性平均剪切耗散率Wv趨于最小,所以湍流脈動耗散Wte=WT-Wv趨于最大,這樣可推導(dǎo)出如下關(guān)系:

    建立合理湍流模型的關(guān)鍵在于如何更深刻地理解計算網(wǎng)格中的湍流結(jié)構(gòu),而這與穩(wěn)定性條件密切相關(guān)。如果不考慮這種穩(wěn)定性條件,很難開發(fā)出更好的模型。EMMS 湍流模型考慮穩(wěn)定性條件,利用式(4)進行封閉湍流動力學(xué)方程。

    3 EMMS湍流模型

    工程中的湍流主要表現(xiàn)為層流和湍流共存。然而傳統(tǒng)的湍流模型假設(shè)流體在計算網(wǎng)格中處于均勻和充分湍流狀態(tài),從而忽略了流動中的層流部分,導(dǎo)致工程中實際流動的湍流模擬并不精確。湍流模擬中,面對計算網(wǎng)格內(nèi)的非均勻性,通常采取的方法為平均化處理、通過濾波只解析大渦和解析所有尺度的渦結(jié)構(gòu),而忽視了湍流中的介尺度渦團及黏性和慣性之間競爭中協(xié)調(diào)的主導(dǎo)機制,因而難以同時兼顧模擬的準(zhǔn)確性和計算量。

    3.1 EMMS湍流模型的建立

    在由不同大小“旋渦”組成的單相湍流中,存在多尺度結(jié)構(gòu)和復(fù)雜的相互作用[11]。湍流中渦的運動在空間尺度和時間尺度上是動態(tài)耦合[53],由于大渦并不穩(wěn)定,最終會分裂成許多小渦,伴隨著湍動能由大渦向小渦傳遞。眾多小渦繼承了前一個旋渦的能量,經(jīng)歷同樣的過程分裂成更小的渦,再將能量向下傳遞。通過這種方式,能量從大的運動尺度一直傳遞到足夠小尺度。

    微尺度:范圍從分子尺度到Kolmogorov 尺度η,該尺度下黏性和慣性交替變化。當(dāng)黏性占主導(dǎo),一組分子作為一個整體進行相同的運動;而當(dāng)慣性占主導(dǎo)時,考慮兩個旋渦之間的界面高速剪切作用。

    介尺度:一定范圍的渦代表湍流的介尺度結(jié)構(gòu)。湍流慣性子區(qū)尺度(η?l?L)遠大于Kolmogorov 尺度,但遠小于宏觀的流動尺度。介尺度渦從大尺度渦的波動中獲取動能,并將其傳遞給更小的渦。介尺度渦可看作是將大尺度含能渦與小尺度耗散渦聯(lián)系起來的橋梁。

    宏尺度:包含大尺度流體的流動,L。大尺度流動受邊界特定的幾何特征控制,因而是各向異性的。大多數(shù)大尺度波動從平均流接收能量,并將其轉(zhuǎn)換為小尺度波動。宏尺度包含了尺度范圍內(nèi)渦的大部分動能。宏觀尺度渦的作用是從主流獲取能量以維持湍流。因此,宏尺度行為體現(xiàn)在與壁面和邊界效應(yīng)相關(guān)的全局波動。

    隨著EMMS模型和DBS模型在氣固和氣液兩種系統(tǒng)中的成功發(fā)展,EMMS 建模思路擴展到了單相湍流體系[52,54],即EMMS湍流模型(圖2),其使用了穩(wěn)定條件來封閉湍流。具體來說,視流動為由湍流流體相和非湍流流體相組成的兩相流問題,引入描述該兩相流系統(tǒng)的介尺度結(jié)構(gòu)參數(shù),特別是湍流流體成分所占的體積分數(shù),使湍流有效黏性系數(shù)進一步增強所含的湍流介尺度結(jié)構(gòu)信息。隨后,根據(jù)湍流渦級串理論,將湍流中包含的能量耗散分解為不同部分,然后進行量化,與氣固系統(tǒng)EMMS模型的原理類似,通過慣性和黏性協(xié)調(diào)控制機制形成的湍流穩(wěn)定性條件封閉湍流模型,優(yōu)化出湍流的介尺度結(jié)構(gòu)參數(shù),改進湍流模式理論中的渦黏性系數(shù)。

    圖2 湍流系統(tǒng)的介尺度框架[52]Fig.2 A mesoscale framework for turbulence system[52]

    3.1.1 流動分相與對應(yīng)的約束方程 最新發(fā)展的雙渦EMMS 湍流模型[55]借鑒兩相相互滲透的雙流體模型,假設(shè)在流動中具有層流流體成分、大渦流體成分和小渦流體成分組分(圖3)。則該系統(tǒng)可通過下述結(jié)構(gòu)參數(shù)來描述:層流、大渦和小渦的體積分數(shù)(fl,fL和fS);大渦和小渦對應(yīng)的等效直徑(dL和dS);層流、大渦和小渦的表觀速度(Ul,UL和US)。

    圖3 雙渦EMMS湍流模型結(jié)構(gòu)示意圖[55]Fig.3 Structure diagram of dual-eddy EMMS-based turbulence model[55]

    當(dāng)流動處于穩(wěn)定狀態(tài)時,渦團所受的曳力和有效浮力相平衡,大、小渦團的力平衡方程分別如下:

    其中,g是重力加速度。大、小渦團和周邊層流流體之間的滑移速度分別為:

    本模型通過將渦團類比于氣泡,從而得到其曳力系數(shù)。大渦的曳力系數(shù)計算公式如下:

    小渦的曳力系數(shù)可以根據(jù)式(9)以此類推計算而得。

    假設(shè)求解范圍內(nèi)的所有湍動能k均由大渦貢獻,其湍流強度為終端湍流強度I0,則大渦的特征方程表示為:

    由耗散渦Reynolds數(shù)約為1可得:

    此外,流體速度和相體積分數(shù)的守恒方程分別為:

    3.1.2 能耗分解 假設(shè)流體流動處于穩(wěn)態(tài),則總能量耗散率WT等于能量輸入率,可分解為四個部分:湍流渦團與周邊流體表面互相振蕩而引起的能量耗散率Ws,大尺度渦團破碎而產(chǎn)生的能量耗散率Wbk,耗散渦尺度的小渦能量耗散率WKolmogorov,以及發(fā)生在層流(非湍流)成分流體內(nèi)的分子黏性耗散Wv。根據(jù)上述描述得到方程:

    WT可近似為初始大渦的能量耗散率:

    Ws與湍流渦團表面的劇烈振蕩而引起的阻力FD,surf有關(guān):

    Wbk可寫成如下形式:

    其中,ωL,λ(dL,λ)為抵達頻率;PL?(dL|fBV,λ)為碰撞破碎的概率密度分布函數(shù);dL為被撞擊的渦團尺度;λ為撞擊的渦團尺度;EDIS為渦團破碎時的黏性耗散。

    WKolmogorov根據(jù)均勻各向同性湍流的耗散譜計算得出:

    通過輸入U、k、σ、μl、ρl、ρL、ρS求解以上代數(shù)方程,結(jié)合式(4)即可優(yōu)化得到該工況下的層流體積分數(shù)、層流速度、大渦體積分數(shù)、大渦速度等一系列結(jié)構(gòu)參數(shù)。

    3.2 EMMS湍流模型與CFD耦合

    雙渦EMMS 湍流模型假設(shè)在一般湍流中有層流、大渦和小渦三部分共存,并相互滲透。在CFD求解RANS 方程中,首先獲得計算網(wǎng)格中的流體速度U和湍動能k,通過求解EMMS湍流模型即可獲得該網(wǎng)格單元下的結(jié)構(gòu)參數(shù)fL,fS,fl?;诖鬁u體積分數(shù)fL對RANS 模型中的湍流黏性系數(shù)進行修正。以耦合RANS方程中的k-ωSST模型為例,湍流黏性系數(shù)μt被修正為:

    即將控制方程中的μt乘以由fL歸一化得到的湍流控制因子f*,同時,k方程中的湍動能生成項也需要乘以f*來抑制湍動能在局部區(qū)域的生成。基于EMMS模型的計算結(jié)果,f*已擬合為代數(shù)形式:

    3.3 EMMS模型顯著改進了預(yù)測精度

    基于EMMS 原理的湍流模型顯著提高了RANS在湍流建模和仿真中的精度[56]。例如,在頂蓋方腔流中,EMMS 湍流模型可以成功捕捉到左下角和右下角的三次角渦,而標(biāo)準(zhǔn)的k-ε模型無法預(yù)測(圖4)。此外,EMMS 湍流模型預(yù)測的流線相比于標(biāo)準(zhǔn)k-ε模型預(yù)測的流線更符合DNS結(jié)果[52,54]。

    圖4 基于EMMS的湍流模型預(yù)測的Re=10000時頂蓋驅(qū)動方腔流的流線圖與DNS和k-ε的流線比較[52]Fig.4 Streamline patterns for lid-driven cavity flow at Re=10000 predicted by the EMMS-based turbulence model in a comparison with those from DNS and k-ε model[52]

    邊界層轉(zhuǎn)捩預(yù)測在高速飛行器的設(shè)計中起著重要作用[57]。與CFD 發(fā)展的其他領(lǐng)域相比,轉(zhuǎn)捩區(qū)的CFD 預(yù)測仍然是一個難點。RANS 方法用于預(yù)測轉(zhuǎn)捩問題具有計算成本低的優(yōu)勢[58]。文獻[59-61]中報道了RANS 下三種計算層流到湍流轉(zhuǎn)捩的方法:(1)在實驗確定的轉(zhuǎn)捩位置打開湍流模型或使用湍流黏度;(2)利用低Reynolds 數(shù)湍流模型;(3)使用間歇的概念來混合從層流到湍流區(qū)域的流動。這些方法大多缺乏真實的物理基礎(chǔ),適用范圍相當(dāng)有限。

    為了耦合CFD 進行計算,基于OpenFOAM 進行了二次開發(fā),構(gòu)建EMMS湍流模型庫、不可壓縮湍流模型動態(tài)庫和可壓縮湍流模型動態(tài)庫?;陔p渦EMMS 湍流模型優(yōu)化湍流體積分數(shù)f,并與k–ωSST湍流模型耦合[9],實現(xiàn)了平板及NACA0012翼型邊界層轉(zhuǎn)捩的預(yù)測,結(jié)果如圖5 和圖6 所示。算例計算網(wǎng)格與邊界條件如圖7 和圖8 所示。f*由0 逐漸增大到1表明流體微元由完全層流控制逐步轉(zhuǎn)變?yōu)橥耆牧骺刂啤D9 和圖10 分別為平板和翼型邊界層的表面摩擦系數(shù),Cf突增的點即為轉(zhuǎn)捩位置,與圖5 和圖6 中f*突變的位置相對應(yīng),且均與實驗數(shù)據(jù)[62-63]吻合較好,驗證了雙渦EMMS 湍流模型預(yù)測轉(zhuǎn)捩問題的可行性。

    圖5 平板邊界層轉(zhuǎn)捩算例的湍流控制因子f*分布[55]Fig.5 Turbulence control factor distribution for flat plate boundary transition case[55]

    圖6 NACA0012翼型上表面轉(zhuǎn)捩點附近湍流控制因子f*分布Fig.6 Turbulence control factor distribution near transition point on the NACA0012 top-surface

    圖7 T3A 算例計算網(wǎng)格與邊界條件Fig.7 Computational mesh and boundary conditions for T3A case

    圖8 NACA0012翼型算例計算網(wǎng)格Fig.8 Computational mesh for NACA0012 airfoil case

    圖9 雙渦EMMS湍流模型和k-ω SST湍流模型預(yù)測平板邊界層表面摩擦系數(shù)Cf及與實驗數(shù)據(jù)的比較Fig.9 Skin friction Cf predicted by two-eddy EMMS-based turbulence model and k-ω SST turbulence model in a comparison with experimental data on the flat plate boundary

    圖10 雙渦EMMS湍流模型和k-ω SST湍流模型預(yù)測NACA0012翼型表面摩擦系數(shù)Cf[55]Fig.10 Skin friction Cf predicted by two-eddy EMMS-based turbulence model and k-ω SST turbulence model in comparison with experimental data on the NACA0012 airfoil top-surface[55]

    復(fù)雜的氣候模型由描述大氣、海洋、陸地表面和冰的數(shù)學(xué)方程組成[64]。氣象大氣和海洋模式涉及典型的工程湍流模擬,計算網(wǎng)格通常達到幾十公里量級,在粗網(wǎng)格里湍動和非湍動(即層流)區(qū)域共存,而傳統(tǒng)模式都是假設(shè)網(wǎng)格內(nèi)是均勻狀態(tài),從而導(dǎo)致模型預(yù)測不準(zhǔn)確。通過應(yīng)用EMMS湍流模型對全球海洋環(huán)流進行了模擬,得到了初步結(jié)果(圖11)。湍流結(jié)構(gòu)只是網(wǎng)格內(nèi)非均勻性的一部分,未來的工作可以將介尺度理論內(nèi)涵擴展到更一般的情況,如計算網(wǎng)格內(nèi)除了湍流,還可以考慮云、植被生長、雨雪、地形高程等在粗網(wǎng)格內(nèi)的不均勻分布等,這些都屬于介尺度理論的研究范疇?;诮榭茖W(xué)概念考慮計算網(wǎng)格內(nèi)的非均勻性,避免平均化方法的均勻處理,將會對天氣預(yù)報和氣候模擬的計算精度有很大改進。

    圖11 全球海洋環(huán)流的EMMS模擬Fig.11 EMMS modeling in ocean currents

    3.4 EMMS湍流模型驗證了介科學(xué)理論的普適性

    介科學(xué)認為系統(tǒng)的復(fù)雜性源于共存的兩種(或更多)主導(dǎo)機制之間競爭中協(xié)調(diào)[65]。不同機制的相對主導(dǎo)作用隨著系統(tǒng)狀態(tài)的變化而變化。以兩種主導(dǎo)機制為例,系統(tǒng)將出現(xiàn)三個區(qū)域,即機制A 主導(dǎo)、機制A 和機制B 之間競爭中協(xié)調(diào)以及機制B 主導(dǎo)[66]。復(fù)雜性出現(xiàn)于機制A 和機制B 協(xié)調(diào)(共存)的介區(qū)域中。當(dāng)系統(tǒng)狀態(tài)處于兩個極端區(qū)域之間的介區(qū)域時,不同共存主導(dǎo)機制之間競爭中協(xié)調(diào)導(dǎo)致了介尺度的復(fù)雜性[67]。因此,介科學(xué)很好地闡述了這種復(fù)雜性,對介尺度和介區(qū)域而言具有重要意義。

    基于雙渦EMMS 湍流模型同樣計算出了Li等[40,66-67]提出的復(fù)雜系統(tǒng)中不同區(qū)域的行為。如圖12 所示,隨著流速的增加,依次出現(xiàn)了完全黏性控制、黏性和慣性之間競爭中協(xié)調(diào)以及完全慣性控制三種區(qū)域,大渦等效直徑dL在中間位置出現(xiàn)了劇烈波動。說明流動的復(fù)雜性同樣出現(xiàn)在黏性控制機制與慣性控制機制競爭中協(xié)調(diào)(層流成分和湍流成分流體共存)的介區(qū)域。這表明,介科學(xué)的概念在湍流系統(tǒng)中得到進一步證實,同時也為介科學(xué)理論作為復(fù)雜系統(tǒng)普適理論提供證據(jù)。

    圖12 湍流系統(tǒng)的三種區(qū)域預(yù)測Fig.12 Three-regime prediction identified in turbulence system

    4 結(jié)論與展望

    介科學(xué)聚焦于單元尺度和系統(tǒng)尺度之間存在的控制機制之間的競爭中協(xié)調(diào)導(dǎo)致的結(jié)構(gòu),并建立其穩(wěn)定性條件,從而實現(xiàn)對介尺度結(jié)構(gòu)及其性能的定量描述[65]。在揭示競爭中協(xié)調(diào)的主導(dǎo)機制及其關(guān)系時,應(yīng)特別注意介尺度問題的層次特征,主導(dǎo)機制之間的競爭中協(xié)調(diào)導(dǎo)致了不同尺度上的結(jié)構(gòu)差異和時空多尺度行為,每種主導(dǎo)機制自身的規(guī)律及其極值趨勢。

    湍流一直被認為是經(jīng)典物理中最復(fù)雜的問題之一。一個多世紀以來,吸引了眾多科學(xué)家的持續(xù)努力,但人們對它的認識還只是局部的,目前還沒有完全令人滿意的湍流理論。而在湍流模擬中,湍流模型對于模擬結(jié)果的準(zhǔn)確性至關(guān)重要,而當(dāng)前湍流模型大都是基于實驗數(shù)據(jù)或經(jīng)驗關(guān)聯(lián)式封閉湍流, 對物理模型的內(nèi)部機理缺乏進一步探索,因而難以找到一種普適的湍流模型。介尺度湍流模型的發(fā)展將新興的介科學(xué)思想與CFD 相結(jié)合,對模擬過程工程中的湍流問題和介尺度結(jié)構(gòu)具有理論和實際意義。EMMS 建模思路成功地推廣到湍流系統(tǒng),通過建立工程湍流非均勻結(jié)構(gòu)的物理模型,闡明其控制機制及其相互關(guān)系,解釋了介尺度行為——非均勻渦團結(jié)構(gòu)的形成機制及其對動量傳遞過程的影響。目前基于介科學(xué)視角下的EMMS湍流模型主要有如下進展:(1)在頂蓋方腔流中捕捉到標(biāo)準(zhǔn)的k-ε模型無法預(yù)測到的三次角渦,與DNS結(jié)果更加符合;(2)基于雙渦EMMS 湍流模型成功預(yù)測了平板及NACA0012 邊界層轉(zhuǎn)捩;(3)湍流系統(tǒng)中成功復(fù)現(xiàn)了介科學(xué)理論預(yù)測的完全黏性控制、黏性和慣性之間競爭中協(xié)調(diào)以及完全慣性控制三種區(qū)域,進一步驗證了介科學(xué)理論在復(fù)雜系統(tǒng)中的普適性。

    介科學(xué)有助于解決非均勻湍流系統(tǒng)的定量模擬,EMMS原理改進了湍流建模和模擬的精度,使其能夠解決工程實際問題。湍流系統(tǒng)的EMMS模型為介科學(xué)作為共性原理提供了證據(jù)。EMMS原理對介尺度問題具有普遍性,EMMS 思路框架有助于交叉學(xué)科的融合發(fā)展,介科學(xué)的發(fā)展將極大提高不同學(xué)科解決復(fù)雜問題的能力。

    致謝:感謝中國科學(xué)院過程工程研究所李靜海院士長期以來對本研究給予的指導(dǎo)及鼓勵!

    猜你喜歡
    大渦層流黏性
    層流輥道電機IP56防護等級結(jié)構(gòu)設(shè)計
    防爆電機(2022年5期)2022-11-18 07:40:18
    摻氫對二甲醚層流燃燒特性的影響
    層流切應(yīng)力誘導(dǎo)microRNA-101下調(diào)EZH2抑制血管新生
    富硒產(chǎn)業(yè)需要強化“黏性”——安康能否玩轉(zhuǎn)“硒+”
    如何運用播音主持技巧增強受眾黏性
    傳媒評論(2019年4期)2019-07-13 05:49:28
    基于壁面射流的下?lián)舯┝鞣欠€(wěn)態(tài)風(fēng)場大渦模擬
    軸流風(fēng)機葉尖泄漏流動的大渦模擬
    玩油灰黏性物成網(wǎng)紅
    華人時刊(2017年17期)2017-11-09 03:12:03
    基層農(nóng)行提高客戶黏性淺析
    基于大渦模擬的旋風(fēng)分離器錐體結(jié)構(gòu)影響研究
    国产亚洲av嫩草精品影院| 午夜免费激情av| 精品国产三级普通话版| 久久国产乱子免费精品| av在线观看视频网站免费| 久久久久久九九精品二区国产| 欧美日韩瑟瑟在线播放| 中文资源天堂在线| 99热精品在线国产| 亚洲精品乱码久久久v下载方式| 18禁黄网站禁片午夜丰满| 美女cb高潮喷水在线观看| 午夜福利在线观看吧| 超碰av人人做人人爽久久| 琪琪午夜伦伦电影理论片6080| 久久久久久久久大av| 久久久成人免费电影| 搡老岳熟女国产| 国产精品福利在线免费观看| 男人舔奶头视频| 在线观看av片永久免费下载| 国产免费一级a男人的天堂| 亚洲 国产 在线| 久久欧美精品欧美久久欧美| av中文乱码字幕在线| 男女啪啪激烈高潮av片| 午夜免费成人在线视频| 五月伊人婷婷丁香| 国产精品,欧美在线| 18禁黄网站禁片免费观看直播| 久久国产乱子免费精品| 美女大奶头视频| 成人国产麻豆网| 久久久成人免费电影| 久久亚洲精品不卡| 亚洲美女视频黄频| 亚洲精品成人久久久久久| 一级a爱片免费观看的视频| 亚洲精品成人久久久久久| 欧美bdsm另类| 男女下面进入的视频免费午夜| 少妇被粗大猛烈的视频| 精品国产三级普通话版| 久久人人精品亚洲av| 精品人妻1区二区| 亚洲乱码一区二区免费版| 久久国产乱子免费精品| 男女做爰动态图高潮gif福利片| 日本色播在线视频| 久久精品夜夜夜夜夜久久蜜豆| 日韩欧美三级三区| 亚洲18禁久久av| 亚州av有码| 亚洲国产欧洲综合997久久,| 国产成年人精品一区二区| 国产成年人精品一区二区| 亚洲熟妇熟女久久| 麻豆国产97在线/欧美| 最近视频中文字幕2019在线8| 国产精品自产拍在线观看55亚洲| 国内毛片毛片毛片毛片毛片| 欧美成人性av电影在线观看| 小蜜桃在线观看免费完整版高清| 91精品国产九色| 成人鲁丝片一二三区免费| 性欧美人与动物交配| 欧美成人免费av一区二区三区| 精品久久久久久成人av| av在线老鸭窝| 免费观看精品视频网站| 免费av观看视频| 国内精品一区二区在线观看| 伦理电影大哥的女人| 精品不卡国产一区二区三区| 不卡一级毛片| 精品一区二区三区视频在线观看免费| 亚洲电影在线观看av| 日韩欧美 国产精品| 国产男靠女视频免费网站| 无遮挡黄片免费观看| 国产一区二区在线观看日韩| 在线观看一区二区三区| 久久亚洲精品不卡| 我的女老师完整版在线观看| 亚洲熟妇熟女久久| 俄罗斯特黄特色一大片| 俄罗斯特黄特色一大片| 一区二区三区免费毛片| 欧美xxxx黑人xx丫x性爽| 亚洲人与动物交配视频| 女的被弄到高潮叫床怎么办 | 窝窝影院91人妻| 中国美白少妇内射xxxbb| 一个人看视频在线观看www免费| 成人国产综合亚洲| 亚洲第一区二区三区不卡| 在线观看av片永久免费下载| 国产男靠女视频免费网站| 久久这里只有精品中国| 久久精品影院6| 国产高清有码在线观看视频| 毛片女人毛片| 亚洲一区高清亚洲精品| av黄色大香蕉| 亚洲成av人片在线播放无| 黄色日韩在线| 两个人的视频大全免费| av.在线天堂| 此物有八面人人有两片| 亚洲成人精品中文字幕电影| 国产主播在线观看一区二区| 国产精品国产三级国产av玫瑰| 亚洲av日韩精品久久久久久密| 男女之事视频高清在线观看| 免费人成视频x8x8入口观看| bbb黄色大片| 亚洲精品乱码久久久v下载方式| 天堂av国产一区二区熟女人妻| 精品福利观看| 三级国产精品欧美在线观看| 深夜a级毛片| 午夜福利欧美成人| 国产精品久久久久久精品电影| 国产高清视频在线观看网站| 亚洲性夜色夜夜综合| 桃红色精品国产亚洲av| 热99re8久久精品国产| 午夜福利成人在线免费观看| 成人综合一区亚洲| 日韩在线高清观看一区二区三区 | 欧美三级亚洲精品| 亚洲精品一区av在线观看| 亚洲美女视频黄频| 亚洲最大成人av| 国产精品久久久久久av不卡| 在现免费观看毛片| 少妇人妻精品综合一区二区 | 色综合婷婷激情| 国产精品一区二区三区四区免费观看 | 日韩欧美免费精品| 日日干狠狠操夜夜爽| 高清毛片免费观看视频网站| 国产欧美日韩一区二区精品| 一级av片app| 国产又黄又爽又无遮挡在线| 麻豆一二三区av精品| 婷婷亚洲欧美| 亚洲 国产 在线| 在线观看一区二区三区| 国产一区二区激情短视频| 成人无遮挡网站| 无人区码免费观看不卡| 日韩欧美国产在线观看| 国产精品不卡视频一区二区| 黄色日韩在线| 欧美日本亚洲视频在线播放| 高清日韩中文字幕在线| 偷拍熟女少妇极品色| 俺也久久电影网| 一卡2卡三卡四卡精品乱码亚洲| 男女那种视频在线观看| 国产免费男女视频| 成年女人看的毛片在线观看| 日本与韩国留学比较| 免费高清视频大片| 免费人成视频x8x8入口观看| 热99re8久久精品国产| 97人妻精品一区二区三区麻豆| 国产伦精品一区二区三区视频9| 成人特级av手机在线观看| 一边摸一边抽搐一进一小说| 中文字幕精品亚洲无线码一区| 日本在线视频免费播放| 国产精品人妻久久久久久| 少妇高潮的动态图| 成人毛片a级毛片在线播放| 久久精品夜夜夜夜夜久久蜜豆| 九九热线精品视视频播放| 男女之事视频高清在线观看| 免费观看的影片在线观看| 干丝袜人妻中文字幕| 国产色爽女视频免费观看| 白带黄色成豆腐渣| 日韩精品青青久久久久久| 国产蜜桃级精品一区二区三区| 黄色视频,在线免费观看| 亚洲最大成人手机在线| 亚洲欧美精品综合久久99| 亚洲三级黄色毛片| 成人高潮视频无遮挡免费网站| 国产高潮美女av| 精品一区二区免费观看| 国产精品女同一区二区软件 | 国产乱人伦免费视频| 久久欧美精品欧美久久欧美| 亚洲一级一片aⅴ在线观看| 搡老妇女老女人老熟妇| 国产麻豆成人av免费视频| 又黄又爽又刺激的免费视频.| 欧美绝顶高潮抽搐喷水| 99国产极品粉嫩在线观看| 免费不卡的大黄色大毛片视频在线观看 | 两人在一起打扑克的视频| 老师上课跳d突然被开到最大视频| 成人永久免费在线观看视频| 22中文网久久字幕| 亚州av有码| 国产精品久久久久久av不卡| 国产亚洲欧美98| 啦啦啦韩国在线观看视频| 日韩中文字幕欧美一区二区| 国产精品av视频在线免费观看| 久久久久久国产a免费观看| av视频在线观看入口| 夜夜夜夜夜久久久久| 99精品在免费线老司机午夜| 国产极品精品免费视频能看的| 少妇丰满av| 极品教师在线免费播放| 久久久久性生活片| 国产av麻豆久久久久久久| 国产极品精品免费视频能看的| 欧美精品国产亚洲| 一本久久中文字幕| 久久久久精品国产欧美久久久| 日韩国内少妇激情av| 久久人人精品亚洲av| 色综合站精品国产| 桃红色精品国产亚洲av| 人妻久久中文字幕网| 日日摸夜夜添夜夜添av毛片 | 高清毛片免费观看视频网站| 中文字幕免费在线视频6| 国产精品人妻久久久久久| 99久久久亚洲精品蜜臀av| 久久久久久久久久黄片| 免费观看的影片在线观看| 亚洲男人的天堂狠狠| 欧美人与善性xxx| 日韩欧美精品v在线| 舔av片在线| 人妻少妇偷人精品九色| 啦啦啦观看免费观看视频高清| 九九在线视频观看精品| 午夜a级毛片| 欧美激情在线99| 蜜桃久久精品国产亚洲av| 一个人看的www免费观看视频| 国内揄拍国产精品人妻在线| 国产高清有码在线观看视频| 欧美成人免费av一区二区三区| 看免费成人av毛片| 国产真实伦视频高清在线观看 | 国产精品一区二区免费欧美| 精品福利观看| 男女下面进入的视频免费午夜| 欧美zozozo另类| 欧美日韩亚洲国产一区二区在线观看| 少妇熟女aⅴ在线视频| 深夜a级毛片| 91久久精品国产一区二区成人| 女人被狂操c到高潮| 亚洲无线观看免费| 亚洲乱码一区二区免费版| 黄色视频,在线免费观看| 精品人妻熟女av久视频| 99热这里只有是精品50| 一级黄色大片毛片| 草草在线视频免费看| 91在线精品国自产拍蜜月| 成人二区视频| 日韩欧美在线乱码| 身体一侧抽搐| 亚洲国产欧洲综合997久久,| 国产亚洲av嫩草精品影院| 精品久久久久久久久久免费视频| 女同久久另类99精品国产91| 99久久精品热视频| 欧美xxxx黑人xx丫x性爽| 日本在线视频免费播放| 国语自产精品视频在线第100页| 久久这里只有精品中国| av视频在线观看入口| 校园春色视频在线观看| 少妇高潮的动态图| 亚洲欧美日韩无卡精品| 国国产精品蜜臀av免费| 1024手机看黄色片| 成人精品一区二区免费| 中文字幕av成人在线电影| 国产国拍精品亚洲av在线观看| 国产精华一区二区三区| 久久国产精品人妻蜜桃| 成年女人看的毛片在线观看| АⅤ资源中文在线天堂| 国产真实乱freesex| 免费一级毛片在线播放高清视频| 亚洲一区二区三区色噜噜| 久久精品国产清高在天天线| 欧美zozozo另类| 欧洲精品卡2卡3卡4卡5卡区| 欧美日韩综合久久久久久 | 成人午夜高清在线视频| 特大巨黑吊av在线直播| 国产黄片美女视频| 他把我摸到了高潮在线观看| 中亚洲国语对白在线视频| 高清在线国产一区| 国产视频一区二区在线看| 高清日韩中文字幕在线| 一本久久中文字幕| 国产在线精品亚洲第一网站| 1000部很黄的大片| 亚洲精品国产成人久久av| 免费看日本二区| 少妇熟女aⅴ在线视频| 精品一区二区免费观看| 午夜免费激情av| 国产伦精品一区二区三区四那| 此物有八面人人有两片| 狂野欧美激情性xxxx在线观看| 干丝袜人妻中文字幕| 又黄又爽又免费观看的视频| 国产免费av片在线观看野外av| 欧美色视频一区免费| 久久久久久久久中文| 色在线成人网| 校园人妻丝袜中文字幕| 日本在线视频免费播放| 全区人妻精品视频| 国产亚洲91精品色在线| 午夜激情欧美在线| 国产精品国产高清国产av| 婷婷六月久久综合丁香| 国产精品嫩草影院av在线观看 | 中文字幕av成人在线电影| 两个人视频免费观看高清| 久99久视频精品免费| 免费一级毛片在线播放高清视频| 日韩欧美精品v在线| 国产精品久久视频播放| 黄色欧美视频在线观看| 男女啪啪激烈高潮av片| 99九九线精品视频在线观看视频| 成人亚洲精品av一区二区| 国产亚洲精品久久久久久毛片| 亚洲av免费在线观看| 午夜精品久久久久久毛片777| 国产精品不卡视频一区二区| 久久精品夜夜夜夜夜久久蜜豆| 国产精华一区二区三区| 欧美潮喷喷水| 亚洲欧美日韩无卡精品| 最后的刺客免费高清国语| 大又大粗又爽又黄少妇毛片口| 97碰自拍视频| 狠狠狠狠99中文字幕| 直男gayav资源| 国产精品综合久久久久久久免费| av在线天堂中文字幕| 内地一区二区视频在线| 免费观看在线日韩| 国产一区二区在线av高清观看| 久久久国产成人精品二区| 国产白丝娇喘喷水9色精品| 色吧在线观看| 热99re8久久精品国产| av在线亚洲专区| 久久精品久久久久久噜噜老黄 | 精品无人区乱码1区二区| 成人国产综合亚洲| 亚洲av美国av| 国产中年淑女户外野战色| 国产伦在线观看视频一区| 亚洲综合色惰| 少妇人妻精品综合一区二区 | 国产精品久久久久久精品电影| 午夜精品久久久久久毛片777| 搡老熟女国产l中国老女人| 男女做爰动态图高潮gif福利片| 久久久色成人| 特级一级黄色大片| 99热网站在线观看| 国内精品久久久久久久电影| 欧美最新免费一区二区三区| 久久久久久国产a免费观看| 女的被弄到高潮叫床怎么办 | 久9热在线精品视频| 日韩国内少妇激情av| 丰满人妻一区二区三区视频av| 亚洲欧美激情综合另类| .国产精品久久| 日韩欧美免费精品| 老女人水多毛片| 人人妻人人澡欧美一区二区| 久久午夜福利片| 午夜日韩欧美国产| 国内揄拍国产精品人妻在线| 午夜精品久久久久久毛片777| 不卡视频在线观看欧美| 亚洲av美国av| xxxwww97欧美| 国产高清不卡午夜福利| 亚洲精品国产成人久久av| 亚洲美女搞黄在线观看 | 999久久久精品免费观看国产| 精品免费久久久久久久清纯| 日韩一区二区视频免费看| avwww免费| 国产美女午夜福利| 又黄又爽又免费观看的视频| 亚洲av美国av| 亚洲精品456在线播放app | 99在线视频只有这里精品首页| 一个人看视频在线观看www免费| 亚洲一级一片aⅴ在线观看| 欧美日韩瑟瑟在线播放| 天堂动漫精品| 男人的好看免费观看在线视频| 精品一区二区三区av网在线观看| 国产一区二区三区av在线 | 中文字幕高清在线视频| 国产不卡一卡二| 18禁黄网站禁片午夜丰满| 麻豆成人av在线观看| 99九九线精品视频在线观看视频| 日本三级黄在线观看| 久久欧美精品欧美久久欧美| 午夜福利在线观看吧| 日日摸夜夜添夜夜添av毛片 | 国产综合懂色| 3wmmmm亚洲av在线观看| 久99久视频精品免费| 波野结衣二区三区在线| 最近视频中文字幕2019在线8| 国产又黄又爽又无遮挡在线| 久久人人精品亚洲av| 免费观看的影片在线观看| 老司机深夜福利视频在线观看| 麻豆国产97在线/欧美| 欧美性猛交╳xxx乱大交人| 国产高清三级在线| 村上凉子中文字幕在线| 国产亚洲精品久久久com| 日韩精品中文字幕看吧| 国国产精品蜜臀av免费| 日韩欧美国产在线观看| 久久国产乱子免费精品| 亚洲欧美激情综合另类| 1000部很黄的大片| 俺也久久电影网| 白带黄色成豆腐渣| 欧美一区二区亚洲| 国产伦精品一区二区三区四那| 国产免费av片在线观看野外av| 国产视频内射| 国产精品久久久久久久电影| 精品免费久久久久久久清纯| 成人欧美大片| 亚洲av.av天堂| 特大巨黑吊av在线直播| 日日啪夜夜撸| 嫩草影院入口| 日韩在线高清观看一区二区三区 | 亚洲av第一区精品v没综合| 国产伦精品一区二区三区四那| 亚洲综合色惰| 老司机福利观看| 午夜福利高清视频| 久久久久国内视频| 欧美精品啪啪一区二区三区| 国产黄a三级三级三级人| 亚洲天堂国产精品一区在线| www日本黄色视频网| 99精品久久久久人妻精品| 九九在线视频观看精品| 亚洲av第一区精品v没综合| 中文字幕av在线有码专区| 99热网站在线观看| 亚洲av不卡在线观看| 国产av在哪里看| 别揉我奶头~嗯~啊~动态视频| 国产免费男女视频| 丰满乱子伦码专区| 国产在线精品亚洲第一网站| avwww免费| 老司机福利观看| 国产av在哪里看| 麻豆成人午夜福利视频| 国产精品亚洲一级av第二区| 91在线精品国自产拍蜜月| 国产精品人妻久久久久久| 国产亚洲91精品色在线| 中文在线观看免费www的网站| 一a级毛片在线观看| 久久久国产成人精品二区| 国产av在哪里看| 很黄的视频免费| 男人舔奶头视频| 国产精品美女特级片免费视频播放器| av.在线天堂| 男人的好看免费观看在线视频| 在线a可以看的网站| 一区二区三区四区激情视频 | av视频在线观看入口| 伦精品一区二区三区| 丰满乱子伦码专区| 啦啦啦观看免费观看视频高清| 91麻豆精品激情在线观看国产| 国产久久久一区二区三区| 黄片wwwwww| 99久国产av精品| 深夜a级毛片| 性欧美人与动物交配| 人人妻人人澡欧美一区二区| 韩国av在线不卡| 日韩精品青青久久久久久| 麻豆成人午夜福利视频| 在线观看66精品国产| 欧美zozozo另类| 又黄又爽又免费观看的视频| 一级毛片久久久久久久久女| 亚洲精华国产精华液的使用体验 | 成人精品一区二区免费| 欧美日韩中文字幕国产精品一区二区三区| 久久久精品大字幕| 美女高潮喷水抽搐中文字幕| 乱码一卡2卡4卡精品| 日日撸夜夜添| 一本精品99久久精品77| 黄色女人牲交| 九九爱精品视频在线观看| avwww免费| 免费黄网站久久成人精品| 国产不卡一卡二| aaaaa片日本免费| 身体一侧抽搐| 国产精品嫩草影院av在线观看 | 美女高潮喷水抽搐中文字幕| 十八禁国产超污无遮挡网站| 在线免费观看的www视频| 亚洲人成网站高清观看| 色吧在线观看| 中文在线观看免费www的网站| 色综合亚洲欧美另类图片| av黄色大香蕉| 88av欧美| 国产探花在线观看一区二区| 国产精品三级大全| 日日摸夜夜添夜夜添小说| 成年女人永久免费观看视频| av专区在线播放| 有码 亚洲区| 午夜视频国产福利| 狂野欧美激情性xxxx在线观看| 少妇的逼水好多| 欧美性猛交╳xxx乱大交人| 在线观看舔阴道视频| 免费大片18禁| 99久久九九国产精品国产免费| 男女视频在线观看网站免费| 亚洲经典国产精华液单| 国产淫片久久久久久久久| 免费在线观看影片大全网站| 亚洲四区av| 国产一级毛片七仙女欲春2| 亚洲va在线va天堂va国产| 亚洲第一区二区三区不卡| 搡老妇女老女人老熟妇| 日韩大尺度精品在线看网址| 精品久久久久久久久av| 麻豆一二三区av精品| 久久久久久国产a免费观看| 看黄色毛片网站| 深爱激情五月婷婷| 日本熟妇午夜| 欧美zozozo另类| 亚洲专区国产一区二区| 性色avwww在线观看| 亚洲第一区二区三区不卡| 麻豆国产av国片精品| 91久久精品电影网| 精品无人区乱码1区二区| 神马国产精品三级电影在线观看| 久久久久性生活片| 人人妻人人澡欧美一区二区| 欧美丝袜亚洲另类 | 国产精品一区www在线观看 | 1024手机看黄色片| 中文字幕av成人在线电影| 国产亚洲欧美98| 日韩中字成人| 亚洲美女搞黄在线观看 | 亚洲av美国av| 国产精品电影一区二区三区| 丰满乱子伦码专区| 熟妇人妻久久中文字幕3abv| 欧美激情在线99| av福利片在线观看| 给我免费播放毛片高清在线观看| 91久久精品电影网| 日本撒尿小便嘘嘘汇集6| 亚洲精品日韩av片在线观看| 日韩欧美在线二视频| 97超视频在线观看视频| 日韩大尺度精品在线看网址| a级毛片a级免费在线| 久久精品国产亚洲av涩爱 | 色哟哟哟哟哟哟| 国产精品美女特级片免费视频播放器| 亚洲精品在线观看二区| 精品日产1卡2卡| 色综合亚洲欧美另类图片| 免费观看精品视频网站| 极品教师在线免费播放|