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

    展向旋轉(zhuǎn)平板Couette湍流中的多態(tài)現(xiàn)象

    2020-03-13 01:29:36夏振華
    關(guān)鍵詞:流態(tài)雷諾數(shù)流向

    夏振華

    (浙江大學(xué) 航空航天學(xué)院工程力學(xué)系, 浙江 杭州 310027)

    0 引 言

    旋轉(zhuǎn)條件下的壁湍流問(wèn)題是自然界和工業(yè)工程中常見(jiàn)的湍流現(xiàn)象之一??紤]地球自轉(zhuǎn)作用,那么地球上的壁湍流問(wèn)題幾乎都可以視為旋轉(zhuǎn)壁湍流問(wèn)題。而航空發(fā)動(dòng)機(jī)的轉(zhuǎn)子、渦輪間的流動(dòng)也是典型的旋轉(zhuǎn)壁湍流問(wèn)題。為了設(shè)計(jì)性能更好的航空發(fā)動(dòng)機(jī)等旋轉(zhuǎn)機(jī)械,我們急需了解旋轉(zhuǎn)條件下的壁湍流的統(tǒng)計(jì)規(guī)律,以及對(duì)應(yīng)的流動(dòng)結(jié)構(gòu)隨旋轉(zhuǎn)強(qiáng)度的變化規(guī)律。另一方面,航空發(fā)動(dòng)機(jī)等旋轉(zhuǎn)機(jī)械中關(guān)鍵湍流問(wèn)題有很多,也非常復(fù)雜。為了更好地獨(dú)立出旋轉(zhuǎn)因素的影響,人們提出了一系列簡(jiǎn)化的標(biāo)準(zhǔn)模型問(wèn)題,其中展向旋轉(zhuǎn)平板Couette湍流(Spanwise rotating plane Couette flow, RPCF)就是其中之一[1-2]。

    在旋轉(zhuǎn)系統(tǒng)下,由于科氏力的存在,基本的穩(wěn)定性分析告訴我們,旋轉(zhuǎn)會(huì)根據(jù)其方向的選擇性使得剪切流動(dòng)變得更穩(wěn)定或者更不穩(wěn)定。如果旋轉(zhuǎn)方向和平均流渦量同號(hào)(cyclonic rotation),那么旋轉(zhuǎn)會(huì)使得流動(dòng)更穩(wěn)定;如果兩者異號(hào)(anticyclonic roation),那么旋轉(zhuǎn)會(huì)讓流動(dòng)變得更不穩(wěn)定。在展向旋轉(zhuǎn)的槽道湍流中,由于其平均速度存在極值,平均流渦量會(huì)在槽道中變號(hào),因此旋轉(zhuǎn)會(huì)在一部分槽道中使得流動(dòng)變得不穩(wěn)定;而在另一部分槽道中,其會(huì)使流動(dòng)變得更穩(wěn)定。而在RPCF中,其平均流渦量基本是相同符號(hào)的,因此旋轉(zhuǎn)在整個(gè)槽道保持相同的效應(yīng),即或全場(chǎng)讓流動(dòng)更不穩(wěn)定;或全場(chǎng)讓流動(dòng)更穩(wěn)定[3-5]。

    圖1為RPCF流動(dòng)的示意圖。相距為2h的兩塊平板以相對(duì)速度2Uw沿x方向運(yùn)動(dòng),而整個(gè)系統(tǒng)沿展向(z方向)以常角速度Ωz旋轉(zhuǎn)。該系統(tǒng)有兩個(gè)控制參數(shù):雷諾數(shù)Rew=Uwh/ν和旋轉(zhuǎn)數(shù)Ro=2Ωzh/Uw。

    圖1 RPCF簡(jiǎn)圖[5]Fig.1 Sketch of the flow geometry in RPCF[5]

    在過(guò)去的幾十年中,很多研究組針對(duì)RPCF開(kāi)展了系列的實(shí)驗(yàn)和數(shù)值模擬研究。在實(shí)驗(yàn)方面,Alfredsson教授和其團(tuán)隊(duì)做出了一系列工作[6-9]。特別需要指出的是,Tsukahara、Tillmark和Alfredsson[7]在2010年時(shí)總結(jié)并擴(kuò)充了之前的實(shí)驗(yàn)工作,系統(tǒng)地研究了參數(shù)范圍在0

    在數(shù)值研究上,Andersson教授及其合作者首先針對(duì)Rew=1300時(shí)的RPCF開(kāi)展了直接數(shù)值模擬(Direct numerical simulation, DNS)研究[13-17]。他們研究了弱旋轉(zhuǎn)時(shí)的流動(dòng)統(tǒng)計(jì)以及流動(dòng)結(jié)構(gòu),中高旋轉(zhuǎn)數(shù)時(shí)的旋轉(zhuǎn)效應(yīng)以及極強(qiáng)旋轉(zhuǎn)時(shí)的流動(dòng)統(tǒng)計(jì)。他們發(fā)現(xiàn)在一定旋轉(zhuǎn)數(shù)范圍里,湍流和二次流隨著Ro增加而增強(qiáng),但是超過(guò)一定的旋轉(zhuǎn)數(shù)范圍,其湍流和二次流則會(huì)隨著Ro增加而減弱。在極高旋轉(zhuǎn)數(shù)條件下,二次流消失,這個(gè)時(shí)候的湍流特征與無(wú)旋轉(zhuǎn)時(shí)的流動(dòng)特征差別很大。Salewski和Eckhardt[18]用DNS研究了四個(gè)不同雷諾數(shù)(Rew= 650、1300、2600和5200)下,以壁面剪切為表征的動(dòng)量流隨Ro的變化行為。他們考察了-0.04≤Ro≤0.6之間的算例,發(fā)現(xiàn)變化曲線在兩個(gè)較低雷諾數(shù)情況下只有一個(gè)極值點(diǎn);而在兩個(gè)較高雷諾數(shù)下,變化曲線有兩個(gè)極值點(diǎn),其中一個(gè)極值點(diǎn)位于Ro=0.02附近(較小旋轉(zhuǎn)數(shù))。Gai等[5]也針對(duì)Rew=1300時(shí)的RPCF流動(dòng)開(kāi)展了一系列DNS模擬研究,他們分析了中心線處的平均剪切率、壁面摩擦雷諾數(shù)、體平均的總湍動(dòng)能、體平均的二次流湍動(dòng)能、體平均的剩余場(chǎng)湍動(dòng)能等隨著旋轉(zhuǎn)數(shù)的變化規(guī)律。他們還分析了不同旋轉(zhuǎn)數(shù)下的流動(dòng)結(jié)構(gòu)的異同。

    近些年,本文作者與其合作者還發(fā)現(xiàn)了RPCF中存在多湍態(tài)現(xiàn)象[19-22]。在經(jīng)典的湍流理論中,人們都假設(shè)雷諾數(shù)足夠高的時(shí)候,湍流是各態(tài)遍歷的,也就是說(shuō)所有的湍流態(tài)之間都是等價(jià)的,或者說(shuō)湍流態(tài)是唯一的。在這個(gè)假設(shè)下,雖然瞬時(shí)場(chǎng)跟初始條件有關(guān),但是湍流的統(tǒng)計(jì)規(guī)律是與其初始流場(chǎng)無(wú)關(guān)的[23]。湍流的唯一統(tǒng)計(jì)態(tài)假設(shè)是經(jīng)典湍流理論的基石,也是很多湍流結(jié)論和假設(shè)能夠從低雷諾數(shù)外推到更高雷諾數(shù)的原因之一[24]。然而,近些年來(lái),一些實(shí)驗(yàn)和數(shù)值模擬的結(jié)果顯示,部分湍流中可能存在不唯一的湍流態(tài)[25]。在多態(tài)條件下,相同的控制參數(shù),可能會(huì)出現(xiàn)不同的湍流統(tǒng)計(jì)和流動(dòng)結(jié)構(gòu)。例如,在Taylor-Couette流動(dòng)中,Huismann等[24]就通過(guò)實(shí)驗(yàn)發(fā)現(xiàn)在Re為106量級(jí)下,通過(guò)不同途徑達(dá)到相同的內(nèi)外筒轉(zhuǎn)速比a時(shí),流動(dòng)的Taylor渦對(duì)數(shù)和扭矩等結(jié)果在一定的參數(shù)范圍里(0.17

    本文將就RPCF中的多態(tài)現(xiàn)象的研究現(xiàn)狀進(jìn)行綜述介紹。

    1 RPCF的控制方程和數(shù)值解法

    1.1 控制方程

    RPCF的幾何構(gòu)型如圖1所示。在通常的設(shè)置中,我們假設(shè)流向?yàn)閤,法向?yàn)閥,展向?yàn)閦(或者x1、x2、x3),其對(duì)應(yīng)的速度為u、v、w(或者u1、u2、u3)。用Uw,h無(wú)量綱化之后,系統(tǒng)的控制方程為:

    (1)

    (2)

    這里p是包含了離心力的有效壓力(effective pressure)。其中雷諾數(shù)Rew和旋轉(zhuǎn)數(shù)Ro的定義前文已給出。

    1.2 數(shù)值解法介紹

    在實(shí)際計(jì)算中,為了保證程序的精度,我們采用高精度Fourier-Chebyshev偽譜法對(duì)空間項(xiàng)進(jìn)行離散。具體而言,在流向和展向周期方向,我們采用Fourier級(jí)數(shù)展開(kāi);而在法向采用Chebyshev多項(xiàng)式展開(kāi)。在時(shí)間推進(jìn)上,我們采用了半隱式離散格式,即線性項(xiàng)采用2階精度的Crank-Nicolson根式而非線性項(xiàng)采用的是2階的Adams-Bashforth格式。此外,我們用3/2法則消除混淆誤差。值得說(shuō)明的是,我們的方法比Bech和Andersson[14-15]和Barri和Andersson[17]的數(shù)值解法精度更高。Bech和Andersson用的是2階精度的有限差分程序ECCLES[14-15],而B(niǎo)arri和Andersson則用的是2階精度的有限體積程序MGLET[17]。

    在實(shí)際計(jì)算中,我們?cè)诖笮長(zhǎng)x×Ly×Lz的計(jì)算域上開(kāi)展計(jì)算,離散的網(wǎng)格為Nx×Ny×Nz。具體的計(jì)算網(wǎng)格,在后續(xù)會(huì)做介紹。

    1.3 雷諾分解和三重分解

    (3)

    由于在RPCF中存在二次流,因此我們?cè)诜治龅臅r(shí)候,還可以開(kāi)展三重分解[5,14,19]:

    (4)

    (5)

    (6)

    即湍動(dòng)能和雷諾應(yīng)力都可以分解為由二次流貢獻(xiàn)的部分和剩余場(chǎng)貢獻(xiàn)的部分。

    2 結(jié)果與討論

    本部分將著重介紹相關(guān)的結(jié)果與討論部分。

    2.1 旋轉(zhuǎn)數(shù)對(duì)湍流統(tǒng)計(jì)和流動(dòng)結(jié)構(gòu)的影響

    Gai等人[5]針對(duì)RPCF開(kāi)展了一系列的DNS模擬,基本控制參數(shù)如表1所示。其主要研究目標(biāo)是旋轉(zhuǎn)數(shù)對(duì)湍流統(tǒng)計(jì)和流動(dòng)結(jié)構(gòu)的影響。在分析該問(wèn)題時(shí),研究者們一般重點(diǎn)關(guān)注兩個(gè)量:平均速度在中心線處的法向梯度Ψ和壁面摩擦雷諾數(shù)Reτ。

    表1 計(jì)算參數(shù)Table 1 Computational parameters

    圖2(a、b)分別給出了計(jì)算中20個(gè)算例的中心線速度的法向梯度的平均值及其標(biāo)準(zhǔn)差以及壁面摩擦雷諾數(shù)Reτ隨旋轉(zhuǎn)數(shù)Ro的變化規(guī)律??梢钥吹?,在Ro從0增加到0.02附近時(shí),Ψ逐漸減?。欢?dāng)Ro繼續(xù)從0.02增加到0.9時(shí),Ψ則隨著Ro增大而增加。非常有意思的是,在Ro=0.02時(shí),Ψ是負(fù)值,即在該旋轉(zhuǎn)數(shù)條件下,平均流動(dòng)存在局部回流現(xiàn)象。Salewski和Eckhardt[18]在Ro=0.04時(shí)也發(fā)現(xiàn)了平均速度的局部回流現(xiàn)象。不過(guò)需要說(shuō)明的是,他們所采用的流向計(jì)算域過(guò)小,對(duì)流動(dòng)結(jié)構(gòu)和統(tǒng)計(jì)有一定的影響(在下文,我們將討論這種差別)。而Kawata和Alfredsson[10]他們也通過(guò)實(shí)驗(yàn)測(cè)到過(guò)這種平均速度局部回流現(xiàn)象。相反地,壁面摩擦雷諾數(shù)Reτ隨旋轉(zhuǎn)數(shù)Ro的變化則呈現(xiàn)先增加后減小的趨勢(shì)。其最大值約在Ro=0.15附近。在該問(wèn)題中,Rew保持固定,因此Reτ可以表征壁面摩擦的大小。由圖可以看出,在00.611時(shí),旋轉(zhuǎn)可以起到減阻的作用。

    (b)

    圖3給出了體平均的湍動(dòng)能[k]y以及二次流和剩余場(chǎng)貢獻(xiàn)的部分[ks]y和[k″]y隨旋轉(zhuǎn)數(shù)Ro的變化曲線。其中,算子“[]y”定義如下

    (7)

    圖3 體平均的湍動(dòng)能及其分量隨旋轉(zhuǎn)數(shù)Ro的變化規(guī)律。相關(guān)量都用Ro=0時(shí)的無(wú)量綱化[5]Fig.3 Volume-averaged kinetic energy and its two shears (the contributions from the secondary flows and the residual field) as a function of Ro. All quantities are normalized by at Ro=0[5]

    由于湍動(dòng)能等項(xiàng)已經(jīng)是時(shí)間、流向和展向平均之后的量,因此其法向平均可以認(rèn)為是體平均的結(jié)果。由圖中可以看出,總湍動(dòng)能隨著Ro是先增加后減小的,在Ro=0.25附近達(dá)到最大值。即當(dāng)Ro<0.25時(shí),旋轉(zhuǎn)會(huì)增強(qiáng)湍流;而當(dāng)Ro>0.25時(shí),旋轉(zhuǎn)相反會(huì)抑制湍流。二次流部分的湍動(dòng)能的結(jié)果同總湍動(dòng)能基本保持類似的趨勢(shì),在Ro>0.6時(shí),二次流基本消失。剩余場(chǎng)部分的湍動(dòng)能則體現(xiàn)出同總動(dòng)能相反的趨勢(shì)。它是先減小,后增大,最后又減小的過(guò)程。在Ro=0.1附近其達(dá)到極小值,而在Ro=0.5時(shí),其已經(jīng)非常接近總湍動(dòng)能的值。

    圖4是六個(gè)不同旋轉(zhuǎn)數(shù)下二次流的統(tǒng)計(jì)圖。該圖中,橫截面速度向量由(vs,ws)表征;流向速度分量us則由等值面表征。從圖中可以看出,在Ro<0.25時(shí)隨著Ro增加,二次流逐漸增強(qiáng);而當(dāng)Ro>0.25時(shí),二次流則隨著Ro增加而減弱,在Ro=0.7時(shí)已經(jīng)幾乎無(wú)法看到明顯的二次流了。這和前面的圖3的結(jié)論基本一致。值得注意的是,當(dāng)Ro=0.25和0.32時(shí),二次流的大渦兩側(cè)存在次生的二次流。其中Ro=0.25時(shí),次生二次流強(qiáng)于Ro=0.32時(shí)的結(jié)果。

    (a) Ro=0

    (b) Ro=0.01

    (c) Ro=0.02

    (d) Ro=0.25

    (e) Ro=0.32

    (f) Ro=0.7

    2.2 弱旋轉(zhuǎn)條件下的流動(dòng)結(jié)構(gòu)和統(tǒng)計(jì)

    在上一節(jié)中,我們提到Rew=1300、Ro=0.02時(shí)平均流在槽道中心位置有局部回流現(xiàn)象。本部分將介紹在弱旋轉(zhuǎn)Ro=0.02時(shí)的流動(dòng)結(jié)構(gòu)和統(tǒng)計(jì)分析。在Huang等人[20]的文章中,他們開(kāi)展了一系列DNS研究,分析了網(wǎng)格分辨率、流向計(jì)算域、展向計(jì)算域和計(jì)算時(shí)間的影響。在計(jì)算結(jié)果中,他們發(fā)現(xiàn),在Ro=0.02時(shí),流動(dòng)結(jié)構(gòu)可能會(huì)隨著計(jì)算時(shí)長(zhǎng)而改變,不過(guò)這種計(jì)算結(jié)果還跟流向計(jì)算域有關(guān)。當(dāng)流向計(jì)算域較小,例如Lx=4π時(shí),展向Lz=4π的槽道中二次流的渦對(duì)(Roll-cells)數(shù)恒定為3對(duì);而當(dāng)Lx≥8π時(shí),展向Lz=4π的槽道中會(huì)先出現(xiàn)3對(duì)渦對(duì)而后流動(dòng)會(huì)演變?yōu)?對(duì)渦對(duì)。在流動(dòng)處在3對(duì)渦對(duì)時(shí),平均速度在中心區(qū)存在局部回流;而當(dāng)流動(dòng)處于2對(duì)渦對(duì)時(shí),該局部回流消失。如果Lx=8π,Lz=6π時(shí),該流態(tài)轉(zhuǎn)換依舊存在,流動(dòng)會(huì)從4對(duì)渦對(duì)轉(zhuǎn)換成3對(duì)渦對(duì)。圖5給出了該計(jì)算域下x=0位置處的中心線瞬時(shí)流向速度在展向和時(shí)間上的變化規(guī)律??梢钥闯觯陂L(zhǎng)達(dá)5000個(gè)無(wú)量綱時(shí)間中,流動(dòng)從層流加擾動(dòng)初始解首先演化出4對(duì)流向渦對(duì),該4對(duì)渦對(duì)的狀態(tài)持續(xù)了大約1800無(wú)量綱時(shí)間(200~2000)。之后該流動(dòng)的渦對(duì)開(kāi)始出現(xiàn)合并,最后在3000個(gè)無(wú)量綱時(shí)間的時(shí)候,流動(dòng)變成3對(duì)渦對(duì)結(jié)構(gòu)。值得說(shuō)明的是,我們嘗試過(guò)計(jì)算超過(guò)10 000個(gè)無(wú)量綱時(shí)間,流動(dòng)后期依舊保持在3對(duì)渦對(duì)的狀態(tài)。

    圖5 在x=0處中心線流向瞬時(shí)速度信號(hào)u(0,0,z,t)隨時(shí)間的變化關(guān)系[20]Fig.5 Contours of u(0,0,z,t) as time evolves[20]

    圖6給出了4對(duì)渦和3對(duì)渦結(jié)構(gòu)下,流場(chǎng)的湍流統(tǒng)計(jì)規(guī)律。可以看出,兩者之間是存在明顯差別的。特別的,在4對(duì)渦時(shí),平均速度在中心線處存在局部回流(從圖6(d)中雷諾剪應(yīng)力大于1可以推出來(lái));而當(dāng)流動(dòng)變成3對(duì)渦的時(shí)候,該平均速度的局部回流消失,雷諾剪應(yīng)力的最大值小于1。

    (a) 平均速度

    (d) -〈u′v′〉+

    需要說(shuō)明的是,在前人的工作中,Lz=4π的槽道寬度中只發(fā)現(xiàn)了3對(duì)渦的狀態(tài)[5,14-15]。我們認(rèn)為這主要原因在于過(guò)去的計(jì)算時(shí)間不夠長(zhǎng)。在過(guò)去的計(jì)算中,最多只算了幾百無(wú)量綱時(shí)間[5],圖5顯示了計(jì)算時(shí)間如果不夠長(zhǎng),容易忽略掉后面穩(wěn)定的流態(tài)。因此,我們建議算湍流問(wèn)題時(shí),要用更長(zhǎng)的計(jì)算時(shí)間。

    2.3 多態(tài)現(xiàn)象及其統(tǒng)計(jì)結(jié)果、魯棒性分析

    事實(shí)上,RPCF中的多湍流態(tài)現(xiàn)象在文獻(xiàn)[19]中已有論述。在該文章中,本文作者和合作者針對(duì)Rew=1300,Ro=0.2的RPCF流動(dòng)開(kāi)展了DNS模擬。由于在文獻(xiàn)[5]的基礎(chǔ)上,我們知道隨著Ro從0.02增加到0.2時(shí),二次流是會(huì)變得更強(qiáng)更穩(wěn)定的;而在前期研究Ro=0.02時(shí),我們也發(fā)現(xiàn)了Ro=0.02時(shí)存在流態(tài)轉(zhuǎn)換過(guò)程。因此,在設(shè)計(jì)計(jì)算的時(shí)候,我們選取了兩個(gè)不同的初始場(chǎng),一個(gè)初始場(chǎng)是Ro=0.02時(shí)較前期的瞬時(shí)場(chǎng),此時(shí)該流場(chǎng)擁有3對(duì)渦對(duì)結(jié)構(gòu)(Lz=4π);另一個(gè)初始場(chǎng)則選取Ro=0.02時(shí)后期的瞬時(shí)場(chǎng),該流場(chǎng)擁有2對(duì)渦對(duì)結(jié)構(gòu)。在計(jì)算中,兩個(gè)初始場(chǎng)分別按照相同的程序,在相同的控制參數(shù),相同的計(jì)算域和相同的網(wǎng)格分辨率條件下進(jìn)行時(shí)間推進(jìn)。經(jīng)過(guò)超過(guò)3000個(gè)無(wú)量綱時(shí)間的計(jì)算,我們發(fā)現(xiàn)流動(dòng)的最終統(tǒng)計(jì)結(jié)果和流動(dòng)結(jié)構(gòu)不一樣。另外,我們還在流向和展向進(jìn)行了網(wǎng)格加密,結(jié)果依舊出現(xiàn)多態(tài)現(xiàn)象。

    圖7給出了穩(wěn)態(tài)條件下,兩種不同流態(tài)下的二次流的橫流速度矢量圖和流向分量的等值線圖??梢钥闯觯瑑山M不同的初始條件會(huì)等到不同的流動(dòng)結(jié)構(gòu):一個(gè)會(huì)得到2對(duì)渦對(duì)結(jié)構(gòu);另一個(gè)則會(huì)保持3對(duì)渦對(duì)的狀態(tài)。

    (a) 2對(duì)渦對(duì)的結(jié)果

    (b) 3對(duì)渦對(duì)的結(jié)果

    兩種不同流態(tài)下,其平均速度存在差異。為了說(shuō)明問(wèn)題,表2中給出了兩組不同網(wǎng)格下,2對(duì)渦和3對(duì)渦的一些統(tǒng)計(jì)結(jié)果。從表中的數(shù)據(jù)可以看出:相同流態(tài)下,兩組不同的網(wǎng)格下壁面摩擦雷諾數(shù)Reτ相差很小,誤差都小于1‰。而不同流態(tài)之間的Reτ之間的差別則明顯增大,差別約為6%。類似地,平均速度在中心線處的法向梯度也呈現(xiàn)類似的行為,只是誤差相對(duì)較大。在3對(duì)渦流態(tài)下,兩組不同網(wǎng)格的計(jì)算給出的值約為0.2;而在2對(duì)渦流態(tài)下對(duì)應(yīng)的值減小到約0.14。從穩(wěn)定性分析的角度看,前者在中心區(qū)趨向于中性穩(wěn)定;而后者則處于穩(wěn)定區(qū)域,這也可以解釋為何2對(duì)渦條件下湍動(dòng)能在中心區(qū)會(huì)相對(duì)較小。

    表2 不同流態(tài)下的湍流特征統(tǒng)計(jì)[19]Table 2 Turbulent statistics at different states[19]

    圖8給出了2對(duì)渦和3對(duì)渦時(shí)的湍流二階統(tǒng)計(jì)結(jié)果??梢钥闯觯瑑煞N不同流態(tài)之間的統(tǒng)計(jì)存在非常明顯的差異,而同一種流態(tài)在不同網(wǎng)格下的統(tǒng)計(jì)差異幾乎可以忽略。在靠近壁面的區(qū)域,3對(duì)渦條件下的流向和展向速度脈動(dòng)更強(qiáng);同時(shí),其對(duì)應(yīng)的法向速度脈動(dòng)和雷諾剪應(yīng)力則在槽道中心區(qū)更強(qiáng)。

    為了驗(yàn)證多態(tài)現(xiàn)象的存在是物理的,而不是計(jì)算域造成的虛假結(jié)果,我們?cè)谖墨I(xiàn)[19]中擴(kuò)大了流向和展向的計(jì)算域,計(jì)算結(jié)果基本還能保持多態(tài)的結(jié)果。為了進(jìn)一步驗(yàn)證多態(tài)現(xiàn)象的魯棒性,我們還考察了多態(tài)現(xiàn)象對(duì)初始流場(chǎng)的依賴性。我們構(gòu)造了一個(gè)新的流場(chǎng)ui:

    (8)

    更進(jìn)一步,我們還模擬實(shí)驗(yàn)的方式[25],在計(jì)算中逐漸增加旋轉(zhuǎn)數(shù)或者逐漸減小旋轉(zhuǎn)數(shù),進(jìn)而來(lái)研究多態(tài)現(xiàn)象的存在區(qū)間。在文獻(xiàn)[21]中,Huang等人開(kāi)展了兩組DNS模擬。在相同的計(jì)算域和計(jì)算網(wǎng)格下,兩組計(jì)算的初始場(chǎng)都是層流解加擾動(dòng)的無(wú)散場(chǎng)。其中一組的旋轉(zhuǎn)數(shù)Ro逐漸從0.01增加到0.02,0.03,0.04,0.05,0.1,0.2,0.3,0.5,0.6;而另一組則按照相同的Ro從0.6逐漸減小到0.01。在每個(gè)Ro上,都計(jì)算1250個(gè)無(wú)量綱時(shí)間,后一個(gè)算例(Ro)是以前一個(gè)算例的最后一個(gè)時(shí)刻的瞬時(shí)流場(chǎng)作為初始場(chǎng)。

    (a) 〈u′u′〉

    (b) 〈v′v′〉

    (c) 〈w′w′〉

    (d) -〈u′v′〉

    圖9 不同參數(shù)α下初場(chǎng)的最終流態(tài)示意圖[22]Fig.9 A sketch of final states with different initial conditions defined by α[22]

    圖10給出了法向瞬時(shí)速度的前四個(gè)模態(tài)的幅值隨計(jì)算時(shí)間的變化曲線。其中,模態(tài)幅值定義如下:

    (9)

    其中,F(xiàn)z(f(z))表征的是函數(shù)f(z)的離散傅里葉變換,“〈〉x”表征流向平均。通常情況下,幅值最大的kz值可以視為渦對(duì)數(shù)。由圖10可以看出,在0.03≤Ro≤0.3范圍里,隨著旋轉(zhuǎn)數(shù)Ro逐漸增加,kz=2的模態(tài)逐漸占主導(dǎo);而當(dāng)Ro逐漸減小時(shí),kz=3的模態(tài)逐漸占主導(dǎo)。

    (a) 隨著Ro增加的結(jié)果

    (b) 隨著Ro減小的結(jié)果

    圖11給出了在Ro在0.01到0.5之間的流向瞬時(shí)速度在x=0、y=0處的時(shí)間歷程。其中,圖11(a)給出的是隨著Ro增加的情況,圖11(b)給出的是隨著Ro減小時(shí)的情況。可以看出在0.03≤Ro≤0.3范圍里,隨著Ro增加時(shí),流場(chǎng)呈現(xiàn)2對(duì)渦;而隨著Ro減小時(shí),流場(chǎng)則有3對(duì)渦。圖11的結(jié)果與Huisman等人[25]在Taylor-Couette流動(dòng)的實(shí)驗(yàn)中的結(jié)果非常相似。

    (a) Ro逐漸增加

    (b) Ro逐漸減小

    我們還分析了壁面摩擦雷諾數(shù)Reτ、體平均的湍動(dòng)能以及體平均的二次流湍動(dòng)能隨著Ro增加和減小的結(jié)果。圖12中給出了體平均的湍動(dòng)能隨著Ro增加和減小時(shí)的變化規(guī)律。由圖可以看出,在0.03≤Ro≤0.3范圍里,兩組計(jì)算出來(lái)的體平均湍動(dòng)能在相同的Ro下存在明顯差別,3對(duì)渦時(shí)(Ro減小方向)的湍動(dòng)能明顯大于2對(duì)渦時(shí)的湍動(dòng)能。在文章中,我們還發(fā)現(xiàn)Reτ和體平均的二次流的湍動(dòng)能都存在類似于體平均湍動(dòng)能一樣的遲滯環(huán)[21]。

    需要說(shuō)明的是,在這兩組計(jì)算中,10個(gè)Ro的具體數(shù)值并非預(yù)先有意設(shè)計(jì)的,而是在計(jì)算中隨機(jī)選取的。計(jì)算時(shí)間間隔1250個(gè)無(wú)量綱時(shí)間也不是特意選取的。在所有計(jì)算結(jié)果出來(lái)之前,我們也不能保證我們能觀測(cè)到遲滯環(huán)現(xiàn)象。最終結(jié)果我們發(fā)現(xiàn)了遲滯環(huán),這至少說(shuō)明了多態(tài)現(xiàn)象可以在很大的Ro范圍里存在。我們承認(rèn)如果改變Ro的具體值選取或者改變時(shí)間間隔,結(jié)果可能不一樣。事實(shí)上,在Huisman等人的實(shí)驗(yàn)研究中,他們也是在一些情況下遇到了多態(tài)現(xiàn)象[25]。

    圖12 體平均的湍動(dòng)能隨著Ro增加和減小的變化規(guī)律[21]Fig.12 Changes of the volume-averaged turbulent kinetic energy as Ro increases or decreases[21]

    2.4 多態(tài)條件下的能量流動(dòng)分析

    在Taylor-Couette湍流中多態(tài)現(xiàn)象的工作中,Huisman等人曾經(jīng)猜測(cè)大尺度相干結(jié)構(gòu)的選擇性作用在多態(tài)現(xiàn)象中占據(jù)重要地位[24]。在RPCF中,二次流就是這種大尺度的流動(dòng)結(jié)構(gòu)。為了驗(yàn)證Huisman等人的假設(shè),我們研究了二次流和剩余場(chǎng)之間的能量傳輸過(guò)程。參考Bech和Andersson[14]和Cai等人[26]的工作,我們將速度場(chǎng)分成四部分:

    (u,v,w)=(〈u〉,〈v〉,〈w〉)+(us,0,0)+

    (0,vs,ws)+(u″,v″,w″)

    (10)

    (11)

    (12)

    (13)

    (a) 3對(duì)渦的流態(tài)

    (b) 2對(duì)渦的流態(tài)

    3 結(jié)論與展望

    我們通過(guò)一系列直接數(shù)值模擬在Rew=1300的RPCF中發(fā)現(xiàn)了多湍流態(tài)現(xiàn)象。在不同的湍流態(tài)下其流動(dòng)結(jié)構(gòu)和湍流統(tǒng)計(jì)差別明顯。該多態(tài)現(xiàn)象比較魯棒,可以在非常大的Ro范圍里被觀察到。

    我們還分析了不同湍流態(tài)下,能量在平均場(chǎng)、二次流場(chǎng)和剩余場(chǎng)之間的傳遞關(guān)系。通過(guò)數(shù)據(jù)分析我們發(fā)現(xiàn),3對(duì)渦時(shí)(Lz=4π),二次流較強(qiáng),此時(shí)二次流相當(dāng)于能量傳輸?shù)囊粋€(gè)橋梁:平均場(chǎng)將大部分能量傳輸給二次流,再由二次流傳給剩余場(chǎng);而當(dāng)流動(dòng)處在2對(duì)渦的流態(tài)時(shí),此時(shí)二次流相對(duì)較弱,平均流雖然傳輸了很多能量給二次流,但是大部分都用于維持二次流本身,只有少部分能量能傳輸給剩余場(chǎng)。而剩余場(chǎng)的能量主要是通過(guò)平均場(chǎng)直接傳輸過(guò)來(lái)的。這種能量傳輸?shù)牟煌憩F(xiàn),側(cè)面驗(yàn)證了大尺度結(jié)構(gòu)在多態(tài)現(xiàn)象中占據(jù)重要地位。

    經(jīng)典的湍流理論假設(shè)湍流是各態(tài)遍歷的,是跟初始場(chǎng)無(wú)關(guān)的。這是湍流建模和模擬的基礎(chǔ)之一。多態(tài)現(xiàn)象的出現(xiàn)將對(duì)湍流建模和模擬提出新的挑戰(zhàn)。如果在模型的構(gòu)建(生成項(xiàng)的大小都不一樣)和驗(yàn)證中,將多態(tài)現(xiàn)象考慮進(jìn)去是值得研究的課題。此外,多湍態(tài)現(xiàn)象在工程問(wèn)題中是否存在?如果存在,它會(huì)對(duì)工程設(shè)計(jì)和安全有什么影響?這些也是非常值得研究的問(wèn)題。

    猜你喜歡
    流態(tài)雷諾數(shù)流向
    側(cè)邊機(jī)組故障對(duì)泵站前池流態(tài)的影響
    小溪啊!流向遠(yuǎn)方
    井岡教育(2020年6期)2020-12-14 03:04:42
    基于Transition SST模型的高雷諾數(shù)圓柱繞流數(shù)值研究
    改進(jìn)邊界條件的非恒定流模型在城市河流橡膠壩流態(tài)模擬中的應(yīng)用
    十大漲幅、換手、振副、資金流向
    失穩(wěn)初期的低雷諾數(shù)圓柱繞流POD-Galerkin 建模方法研究
    基于轉(zhuǎn)捩模型的低雷諾數(shù)翼型優(yōu)化設(shè)計(jì)研究
    民機(jī)高速風(fēng)洞試驗(yàn)的阻力雷諾數(shù)效應(yīng)修正
    流向逆轉(zhuǎn)的啟示
    動(dòng)態(tài)流態(tài)冰蓄冷系統(tǒng)在千級(jí)凈化廠房的應(yīng)用
    国产精品久久久av美女十八| 亚洲av日韩在线播放| 嫩草影视91久久| 精品人妻在线不人妻| 欧美亚洲日本最大视频资源| av电影中文网址| 日日夜夜操网爽| 操出白浆在线播放| 成在线人永久免费视频| 亚洲精品美女久久久久99蜜臀| 国产成人系列免费观看| 日韩制服丝袜自拍偷拍| 欧美成人午夜精品| av电影中文网址| 久久精品国产综合久久久| 精品国产乱子伦一区二区三区| 色尼玛亚洲综合影院| 99热网站在线观看| 69精品国产乱码久久久| 欧美乱妇无乱码| 成年人黄色毛片网站| 亚洲人成伊人成综合网2020| 99久久99久久久精品蜜桃| 999久久久精品免费观看国产| 国产精品香港三级国产av潘金莲| 欧美在线一区亚洲| 极品人妻少妇av视频| 50天的宝宝边吃奶边哭怎么回事| 超碰97精品在线观看| 两个人免费观看高清视频| 亚洲九九香蕉| 热99re8久久精品国产| 欧美久久黑人一区二区| 激情视频va一区二区三区| 在线观看免费视频网站a站| 91成年电影在线观看| 欧美色视频一区免费| 老司机福利观看| 色综合婷婷激情| 亚洲精品国产区一区二| 精品人妻在线不人妻| 久久精品国产亚洲av高清一级| 操美女的视频在线观看| 国产精品电影一区二区三区 | 亚洲精品中文字幕一二三四区| bbb黄色大片| 亚洲在线自拍视频| 国产精品一区二区精品视频观看| 欧美精品高潮呻吟av久久| 亚洲av日韩精品久久久久久密| 如日韩欧美国产精品一区二区三区| 熟女少妇亚洲综合色aaa.| 久久人人爽av亚洲精品天堂| 中文字幕av电影在线播放| 亚洲熟女毛片儿| 久久久精品免费免费高清| 亚洲欧美一区二区三区黑人| bbb黄色大片| av国产精品久久久久影院| 婷婷成人精品国产| 亚洲男人天堂网一区| 欧美日韩成人在线一区二区| √禁漫天堂资源中文www| 757午夜福利合集在线观看| 人人妻,人人澡人人爽秒播| 高清视频免费观看一区二区| 日韩欧美免费精品| 国产国语露脸激情在线看| svipshipincom国产片| 国产成+人综合+亚洲专区| 飞空精品影院首页| 精品人妻在线不人妻| 狠狠婷婷综合久久久久久88av| 丝袜人妻中文字幕| 婷婷丁香在线五月| 亚洲 欧美一区二区三区| 久久国产精品影院| 国产男靠女视频免费网站| 欧美激情极品国产一区二区三区| cao死你这个sao货| 女警被强在线播放| 国内久久婷婷六月综合欲色啪| 王馨瑶露胸无遮挡在线观看| 亚洲专区国产一区二区| 亚洲美女黄片视频| 日本黄色日本黄色录像| 中亚洲国语对白在线视频| 国产成人欧美在线观看 | 1024香蕉在线观看| 亚洲av日韩在线播放| e午夜精品久久久久久久| 人人妻人人添人人爽欧美一区卜| 亚洲欧美精品综合一区二区三区| 国产野战对白在线观看| 露出奶头的视频| 别揉我奶头~嗯~啊~动态视频| 欧美国产精品一级二级三级| 欧美 亚洲 国产 日韩一| 成人国产一区最新在线观看| 亚洲av成人av| 老鸭窝网址在线观看| 久99久视频精品免费| 国产单亲对白刺激| 久久久国产成人精品二区 | 亚洲专区国产一区二区| 国产亚洲欧美精品永久| 亚洲精品乱久久久久久| 女人精品久久久久毛片| 国产精品偷伦视频观看了| 亚洲国产欧美一区二区综合| 免费高清在线观看日韩| 国产成人精品在线电影| videosex国产| 露出奶头的视频| 啦啦啦 在线观看视频| 亚洲成人国产一区在线观看| 亚洲熟妇熟女久久| 亚洲一码二码三码区别大吗| 十分钟在线观看高清视频www| 亚洲中文日韩欧美视频| 婷婷成人精品国产| 日本黄色日本黄色录像| 免费女性裸体啪啪无遮挡网站| 一区在线观看完整版| 制服诱惑二区| 夜夜夜夜夜久久久久| 最近最新中文字幕大全电影3 | 国产真人三级小视频在线观看| 国产欧美亚洲国产| 久久精品人人爽人人爽视色| 在线观看免费午夜福利视频| 久久人妻av系列| 两性午夜刺激爽爽歪歪视频在线观看 | 亚洲专区中文字幕在线| 深夜精品福利| 精品免费久久久久久久清纯 | 自线自在国产av| 天堂动漫精品| 亚洲精华国产精华精| 国产成人一区二区三区免费视频网站| av免费在线观看网站| 精品久久久久久电影网| 中出人妻视频一区二区| 精品国产乱子伦一区二区三区| 91九色精品人成在线观看| 国产精品98久久久久久宅男小说| 纯流量卡能插随身wifi吗| 久久人人97超碰香蕉20202| 亚洲精品久久成人aⅴ小说| 久久久久精品人妻al黑| 天天操日日干夜夜撸| 一区二区三区精品91| 国产极品粉嫩免费观看在线| 亚洲精品国产一区二区精华液| 高潮久久久久久久久久久不卡| 成年版毛片免费区| 波多野结衣av一区二区av| 欧美老熟妇乱子伦牲交| 国产精品久久久久久人妻精品电影| 国产麻豆69| 黄色丝袜av网址大全| 变态另类成人亚洲欧美熟女 | 久久精品国产99精品国产亚洲性色 | 日本wwww免费看| 操美女的视频在线观看| 99精品在免费线老司机午夜| av免费在线观看网站| 女人高潮潮喷娇喘18禁视频| 色在线成人网| tube8黄色片| 亚洲第一青青草原| 51午夜福利影视在线观看| 欧美成狂野欧美在线观看| av片东京热男人的天堂| 国产成人av激情在线播放| 男男h啪啪无遮挡| 亚洲午夜理论影院| 丝袜美足系列| 一区二区三区激情视频| 国产高清激情床上av| 日本精品一区二区三区蜜桃| 久久久久久久精品吃奶| 黑人巨大精品欧美一区二区蜜桃| 两性午夜刺激爽爽歪歪视频在线观看 | 黄色成人免费大全| 99热网站在线观看| 国产成人av激情在线播放| 天堂俺去俺来也www色官网| 亚洲成人免费电影在线观看| 麻豆乱淫一区二区| 亚洲精品国产区一区二| 日韩欧美一区视频在线观看| 91成年电影在线观看| 免费在线观看黄色视频的| 国产男女内射视频| 亚洲专区字幕在线| 久久精品国产a三级三级三级| 日韩三级视频一区二区三区| 桃红色精品国产亚洲av| 国产极品粉嫩免费观看在线| 一级毛片女人18水好多| 欧美成人午夜精品| 巨乳人妻的诱惑在线观看| 日韩视频一区二区在线观看| 久久中文字幕一级| 精品久久久久久久毛片微露脸| 成年人黄色毛片网站| 国产成人av激情在线播放| 下体分泌物呈黄色| 亚洲精品成人av观看孕妇| 午夜精品国产一区二区电影| 亚洲精品久久成人aⅴ小说| 国产不卡av网站在线观看| 亚洲精品粉嫩美女一区| 亚洲成av片中文字幕在线观看| tube8黄色片| 自线自在国产av| 亚洲 国产 在线| 三上悠亚av全集在线观看| 波多野结衣av一区二区av| 午夜免费成人在线视频| 韩国精品一区二区三区| 欧美午夜高清在线| 国产精品一区二区在线观看99| 亚洲人成电影观看| 99久久99久久久精品蜜桃| 久久狼人影院| a在线观看视频网站| 91在线观看av| 在线永久观看黄色视频| 久99久视频精品免费| 久久久久国产精品人妻aⅴ院 | 伦理电影免费视频| 国产精品影院久久| 国产黄色免费在线视频| 久久国产亚洲av麻豆专区| 妹子高潮喷水视频| 国产欧美亚洲国产| 18禁裸乳无遮挡动漫免费视频| 亚洲欧美精品综合一区二区三区| 国产亚洲欧美在线一区二区| 美女扒开内裤让男人捅视频| 成人av一区二区三区在线看| 久久99一区二区三区| 国产成人av激情在线播放| 欧美av亚洲av综合av国产av| 亚洲五月婷婷丁香| 国产亚洲精品一区二区www | 男男h啪啪无遮挡| 久热这里只有精品99| 免费在线观看黄色视频的| 777久久人妻少妇嫩草av网站| 成人免费观看视频高清| 多毛熟女@视频| 乱人伦中国视频| 1024视频免费在线观看| 亚洲av成人一区二区三| 欧美午夜高清在线| 欧美丝袜亚洲另类 | 两个人看的免费小视频| 欧美精品高潮呻吟av久久| 欧美精品一区二区免费开放| 99国产精品一区二区三区| 亚洲精品美女久久av网站| 婷婷丁香在线五月| 亚洲全国av大片| 久久精品成人免费网站| 色播在线永久视频| 美女高潮喷水抽搐中文字幕| 国产成人欧美在线观看 | 在线av久久热| 国产一区二区三区视频了| 亚洲一码二码三码区别大吗| 美女福利国产在线| 少妇被粗大的猛进出69影院| 十八禁人妻一区二区| av福利片在线| 可以免费在线观看a视频的电影网站| 亚洲午夜理论影院| 涩涩av久久男人的天堂| 久久国产精品大桥未久av| 女人精品久久久久毛片| 国产成+人综合+亚洲专区| 在线视频色国产色| 免费女性裸体啪啪无遮挡网站| 老司机午夜福利在线观看视频| 中文字幕精品免费在线观看视频| 在线观看一区二区三区激情| 亚洲aⅴ乱码一区二区在线播放 | 亚洲av第一区精品v没综合| 老司机午夜福利在线观看视频| 欧美 亚洲 国产 日韩一| 男女午夜视频在线观看| videosex国产| x7x7x7水蜜桃| 一级毛片精品| 亚洲专区中文字幕在线| 色94色欧美一区二区| 日本vs欧美在线观看视频| 王馨瑶露胸无遮挡在线观看| 村上凉子中文字幕在线| 亚洲七黄色美女视频| 成熟少妇高潮喷水视频| 成年版毛片免费区| 国产蜜桃级精品一区二区三区 | 中文字幕制服av| 中国美女看黄片| 成人永久免费在线观看视频| 国产午夜精品久久久久久| 黑人猛操日本美女一级片| 亚洲少妇的诱惑av| 这个男人来自地球电影免费观看| 欧美性长视频在线观看| 国产免费av片在线观看野外av| 69精品国产乱码久久久| 亚洲五月色婷婷综合| 亚洲七黄色美女视频| 丝袜美足系列| 亚洲成av片中文字幕在线观看| www.精华液| 午夜免费成人在线视频| 一级作爱视频免费观看| 国产亚洲精品久久久久久毛片 | 一区二区三区激情视频| 亚洲av欧美aⅴ国产| 一本综合久久免费| 女人久久www免费人成看片| www.精华液| 波多野结衣av一区二区av| 人人妻人人澡人人看| 热re99久久精品国产66热6| 亚洲av欧美aⅴ国产| 国产成+人综合+亚洲专区| 国产亚洲欧美在线一区二区| 黄色a级毛片大全视频| 日韩欧美一区视频在线观看| 成年女人毛片免费观看观看9 | 大香蕉久久成人网| 国产乱人伦免费视频| 午夜亚洲福利在线播放| 一区二区三区国产精品乱码| 国产一区二区三区综合在线观看| 久久国产精品人妻蜜桃| 国产免费av片在线观看野外av| 中文字幕色久视频| av天堂在线播放| 欧美日韩福利视频一区二区| 丰满饥渴人妻一区二区三| 大型av网站在线播放| 色播在线永久视频| 岛国毛片在线播放| 精品国产乱子伦一区二区三区| 精品国产美女av久久久久小说| 久久草成人影院| 新久久久久国产一级毛片| 国产精品美女特级片免费视频播放器 | 久久性视频一级片| 又黄又爽又免费观看的视频| 久久人妻福利社区极品人妻图片| 国产成人免费无遮挡视频| 天堂动漫精品| 欧美久久黑人一区二区| 一进一出好大好爽视频| 多毛熟女@视频| 亚洲国产看品久久| 亚洲成国产人片在线观看| 俄罗斯特黄特色一大片| 法律面前人人平等表现在哪些方面| 精品国内亚洲2022精品成人 | 国产精品亚洲av一区麻豆| 免费在线观看日本一区| 国产xxxxx性猛交| 在线观看免费午夜福利视频| 午夜91福利影院| 亚洲精品一二三| 精品久久久久久电影网| 亚洲欧洲精品一区二区精品久久久| av一本久久久久| 日韩欧美国产一区二区入口| 久久精品成人免费网站| 亚洲欧美日韩高清在线视频| 中文字幕制服av| 最近最新中文字幕大全免费视频| 欧美日韩福利视频一区二区| 91在线观看av| 精品视频人人做人人爽| 麻豆国产av国片精品| 国产精品亚洲一级av第二区| 国产精品久久视频播放| 亚洲专区国产一区二区| 国产精品 国内视频| 国产精品乱码一区二三区的特点 | 国产精品影院久久| 国产成人一区二区三区免费视频网站| 男人舔女人的私密视频| 欧美人与性动交α欧美精品济南到| 建设人人有责人人尽责人人享有的| av超薄肉色丝袜交足视频| 老熟妇乱子伦视频在线观看| 亚洲av日韩在线播放| 一本大道久久a久久精品| 俄罗斯特黄特色一大片| 女人爽到高潮嗷嗷叫在线视频| 91成年电影在线观看| 成人精品一区二区免费| 精品人妻熟女毛片av久久网站| 可以免费在线观看a视频的电影网站| 丝瓜视频免费看黄片| 人人妻人人澡人人看| 欧美+亚洲+日韩+国产| 99久久国产精品久久久| 亚洲 国产 在线| 一级作爱视频免费观看| 精品电影一区二区在线| 午夜精品久久久久久毛片777| 少妇粗大呻吟视频| 久久九九热精品免费| 老司机靠b影院| 久9热在线精品视频| 国产高清视频在线播放一区| 99精品久久久久人妻精品| 男人舔女人的私密视频| 99热网站在线观看| 搡老岳熟女国产| 亚洲综合色网址| 国产又色又爽无遮挡免费看| 国产成人啪精品午夜网站| 久久精品国产亚洲av香蕉五月 | 国产精品.久久久| 精品国产国语对白av| 91av网站免费观看| 手机成人av网站| 国产精品久久视频播放| 两性午夜刺激爽爽歪歪视频在线观看 | 午夜福利免费观看在线| 中文字幕精品免费在线观看视频| 亚洲黑人精品在线| 一区二区三区激情视频| 免费在线观看黄色视频的| 久久精品国产a三级三级三级| 最新在线观看一区二区三区| 国产又色又爽无遮挡免费看| 成人18禁高潮啪啪吃奶动态图| 99re6热这里在线精品视频| 精品少妇久久久久久888优播| 岛国毛片在线播放| 亚洲专区字幕在线| 久久亚洲真实| 久久久久久免费高清国产稀缺| 精品熟女少妇八av免费久了| 亚洲成人免费av在线播放| 成年女人毛片免费观看观看9 | e午夜精品久久久久久久| 一级a爱片免费观看的视频| √禁漫天堂资源中文www| 亚洲,欧美精品.| 免费在线观看日本一区| 99国产综合亚洲精品| 高清欧美精品videossex| 啪啪无遮挡十八禁网站| 动漫黄色视频在线观看| 亚洲人成伊人成综合网2020| 日本欧美视频一区| 久久影院123| 黄色视频不卡| 又黄又粗又硬又大视频| 母亲3免费完整高清在线观看| 午夜福利免费观看在线| 欧美国产精品va在线观看不卡| 啪啪无遮挡十八禁网站| 不卡一级毛片| av片东京热男人的天堂| 亚洲精品久久成人aⅴ小说| 9色porny在线观看| 国产一区二区三区在线臀色熟女 | 精品久久久久久,| 久久久国产一区二区| 国产淫语在线视频| 两个人看的免费小视频| 免费高清在线观看日韩| 亚洲国产精品合色在线| 成人手机av| 飞空精品影院首页| 久久久国产成人精品二区 | 久久久久久久久久久久大奶| 亚洲精品国产区一区二| 国产亚洲欧美98| 午夜福利欧美成人| 9色porny在线观看| 天天添夜夜摸| 成人永久免费在线观看视频| 亚洲精品自拍成人| 一进一出抽搐gif免费好疼 | 国产成+人综合+亚洲专区| 国产精品 国内视频| 欧美亚洲日本最大视频资源| 成人国语在线视频| 亚洲视频免费观看视频| 久久国产亚洲av麻豆专区| 51午夜福利影视在线观看| 在线观看免费视频日本深夜| 美女 人体艺术 gogo| 精品国产一区二区三区四区第35| 淫妇啪啪啪对白视频| 99精国产麻豆久久婷婷| 一级黄色大片毛片| 热99re8久久精品国产| 午夜福利欧美成人| 一级毛片女人18水好多| 一夜夜www| 建设人人有责人人尽责人人享有的| 国产国语露脸激情在线看| 国产一区在线观看成人免费| 一级片'在线观看视频| 国产xxxxx性猛交| 久久草成人影院| ponron亚洲| 777久久人妻少妇嫩草av网站| 高清黄色对白视频在线免费看| 757午夜福利合集在线观看| 天堂中文最新版在线下载| 亚洲欧美一区二区三区久久| av有码第一页| 精品久久久精品久久久| 51午夜福利影视在线观看| 亚洲色图av天堂| 啪啪无遮挡十八禁网站| 亚洲欧美一区二区三区久久| 精品人妻熟女毛片av久久网站| 久久午夜综合久久蜜桃| 久久天堂一区二区三区四区| 大码成人一级视频| av天堂在线播放| 日本撒尿小便嘘嘘汇集6| 在线视频色国产色| 国产精品 欧美亚洲| 宅男免费午夜| 成人18禁高潮啪啪吃奶动态图| 精品一区二区三卡| 99香蕉大伊视频| 悠悠久久av| 两性夫妻黄色片| 国产精品亚洲一级av第二区| 这个男人来自地球电影免费观看| 亚洲色图av天堂| 中文字幕另类日韩欧美亚洲嫩草| 亚洲专区字幕在线| 狠狠婷婷综合久久久久久88av| 亚洲欧美日韩高清在线视频| 18禁国产床啪视频网站| 最新在线观看一区二区三区| 精品人妻在线不人妻| 我的亚洲天堂| 天天添夜夜摸| 久久精品亚洲精品国产色婷小说| 亚洲精品美女久久av网站| 看黄色毛片网站| 国产欧美日韩精品亚洲av| 国产aⅴ精品一区二区三区波| 中文字幕人妻丝袜一区二区| 黄色怎么调成土黄色| 精品久久久久久久毛片微露脸| 国产精品国产高清国产av | 欧洲精品卡2卡3卡4卡5卡区| 亚洲国产看品久久| 成人三级做爰电影| 久久狼人影院| 免费久久久久久久精品成人欧美视频| 午夜福利欧美成人| svipshipincom国产片| 在线永久观看黄色视频| 国产高清国产精品国产三级| 男女床上黄色一级片免费看| 两个人看的免费小视频| 99在线人妻在线中文字幕 | 国产一区二区三区综合在线观看| 国产激情欧美一区二区| av有码第一页| 咕卡用的链子| 日本a在线网址| 激情在线观看视频在线高清 | 亚洲精品一卡2卡三卡4卡5卡| 亚洲人成电影观看| 国产成人欧美| bbb黄色大片| 久久久国产欧美日韩av| 欧美另类亚洲清纯唯美| 18禁国产床啪视频网站| 国产99久久九九免费精品| 精品久久久久久久毛片微露脸| 一本综合久久免费| 亚洲aⅴ乱码一区二区在线播放 | 国产av一区二区精品久久| 午夜日韩欧美国产| 国产精华一区二区三区| 欧美在线黄色| 中文字幕高清在线视频| 久久天躁狠狠躁夜夜2o2o| av天堂在线播放| 91老司机精品| 高清黄色对白视频在线免费看| 亚洲 国产 在线| 91大片在线观看| 欧美日韩乱码在线| ponron亚洲| 国产人伦9x9x在线观看| 狠狠婷婷综合久久久久久88av| 女人被躁到高潮嗷嗷叫费观| 一级片'在线观看视频| 一个人免费在线观看的高清视频| 日本黄色视频三级网站网址 | 男女床上黄色一级片免费看| 亚洲午夜精品一区,二区,三区| 国产精品秋霞免费鲁丝片| 国产一卡二卡三卡精品|