999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

突變流道內空化準周期特性研究

2022-07-19 03:49:34王江云侯琳倩
石油學報(石油加工) 2022年4期

王 壯, 王江云, 薛 凱, 王 娟, 侯琳倩

(1.中國石油大學(北京) 重質油國家重點實驗室,北京 102249;2.過程流體過濾與分離技術北京市重點實驗室,北京 102249; 3.徐州徐工基礎工程機械有限公司,江蘇 徐州 221001;4.中國石油大學(北京) 克拉瑪依校區,新疆 克拉瑪依 834000; 5.西安益翔航電科技有限公司,陜西 西安 710065)

當高速液流流經節流孔口等突變流道造成局部壓力低于該溫度下飽和蒸汽壓時會產生汽化空泡,這種現象稱為空化。空化是一種涉及相變、非定常、湍流、可壓縮等多種復雜影響因素的流動現象[1],且具有顯著的準周期性[2-3]。空化過程中強烈的非定常特性以及大尺度空泡的發展、斷裂和脫落等會引起空化噪聲、振動和空蝕等[4-5],但合理利用空泡潰滅時產生的能量亦可為液體殺菌、水力開采等提供便利[6-7]。

為詳細探究空化過程中空泡演變的行為特征和準周期特性,學者們進行了大量研究。曹彥濤等[8]、曹友銓等[9]、劉登成等[10]和Ji等[11]先后采用實驗和數值模擬方法探究了三維扭曲水翼空化的演變過程和發展機制,分析了空泡脫落與漩渦結構間的緊密關系,得到了空泡脫落的周期性演化行為和特征;彭熾等[12]和李偉等[13]利用高速攝像機對空化射流過程進行分析,發現云狀空泡演化過程具有明顯的空泡產生、發展、脫落、潰滅4階段周期性。楊龍等[14]指出不同的三維水翼翼型下空泡演化的準周期特性不同;趙靜等[15]采用DES方法模擬發現空泡演化的周期特性會使翼型表面的升阻力表現出同步的周期性波動變化;Sato等[16]的研究表明云狀空泡演化頻譜中存在明顯的低頻和高頻部分,推測分別是由固有載荷和空泡脫落的振動所造成的;陳瑛等[17]也在繞翼型空泡的周期性流動現象中發現翼型振動頻率與空泡脫落頻率一致。在對空泡脫落演化機制的研究中,Soyama等[18]認為空化射流剪切層中的復雜壓力梯度可能是造成空泡云周期性脫落的原因;Watanabe等[19]和Hutli等[20]認為空泡周期性脫落的原因在于周期性的回射流。此外,在對空泡脫落特征的微觀研究中,Kubota等[21]闡述了空泡的微觀結構,認為在空泡的形成和發展過程中包含了許多微小的空泡,一定程度上會促進回射流發展進而引發空泡脫落。

上述研究雖然較為詳細地分析討論了空泡演化特征及演化機理,并指出空化過程具有準周期特性和非定常特性。但對準周期內的空泡初生、發展、匯聚和潰滅等具體行為特征的剖析還不夠深入,未能全面揭示不同空化狀態下準周期演化特征的變化趨勢和內在機理。因此,筆者以Krella等[22]設計的空化室內動態演化過程的實驗和數值模擬分析為基礎,深入探究不同空化狀態下空泡演化的具體行為特征和準周期特性,以期獲悉空化準周期特性的變化規律和內在機理。

1 空化室幾何模型及網格劃分

圖1(a)為200 mm×30 mm×30 mm方腔空化室結構示意圖(展向厚度30 mm)[22],內部設置有上下2個半徑為9 mm的半圓形凸臺,并可通過移動下凸臺改變節流縫寬,在研究中,調節縫寬為3 mm。圖1(b)和圖1(c)為空化室中心軸面上的監測點位置及三維計算網格。監測點位置分別為靶區內:b1、b2、b3、b4、b5、b6;下游流場內:d1、d2、d3、d4、d5、d6。利用ANSYS ICEM軟件進行完全結構化的六面體網格劃分。為保證湍流的充分發展和數值迭代計算的穩定性,將入口段延長為3倍當量直徑長度,出口段延長為5倍當量直徑長度,對兩凸臺處進行邊界層設置和網格加密處理,越靠近凸臺邊壁,其法向網格越密。

R—Radius圖1 空化室結構及網格劃分Fig.1 Structure and grid of cavitation chamber(a) 2D structure; (b) Monitoring points; (c) Grid generation

2 實驗裝置

為探究空化室內空泡形態和演化特征,搭建如圖2所示的實驗裝置,包括:電磁閥、儲水罐、泵、電磁流量計、壓力表、空化室,并在空化室同水平高度安裝數碼相機用于拍攝空化過程。當自來水流經空化室喉部狹縫處時節流降壓產生空化,可通過調節流量大小和兩凸臺間的狹縫間距改變空化強弱,實現對不同空化狀態的研究分析。

1—Electric valve; 2—Water storage tank; 3—Pump; 4—Electromagnetic meter; 5—Pressure gauge; 6—Cavitation test section; 7—Camera圖2 空化實驗裝置示意圖Fig.2 Schematic diagram of cavitation test devices

3 計算模型與邊界條件

3.1 計算模型

3.1.1 湍流模型

為捕捉三維空化流動中復雜的非定常脈動特性,可以采用近壁處耦合雷諾時均(RANS)、湍流充分發展區運用大渦模擬(LES)的分離渦模型(DES)[23]。LES模型屬于過濾小渦后對N-S方程的直接求解,相比其他基于雷諾時均的湍流模型可以更好地捕捉由于空化現象所引起的流場內的非定常特性。延遲分離渦模型(DDES)屬于DES的一種,其中在近壁處耦合k-ωSST模型的方法可以較好地阻止由RANS向LES轉換過程中“灰色區域”的出現[24]。雖然k-ωSST模型在近壁邊界層的網格要求較高,近壁處第一層網格的無量綱距離y+通常在1~5 范圍內,但恰好滿足了近壁處對空化初生現象捕捉的精細化網格要求。

3.1.2 多相流模型

空化發展演變的過程復雜,大尺度空泡團內氣、液兩相混合共存,且二者在隨主流移動過程中存在速度滑移。Mixture模型屬于Euler均相流模型,適用于相間充分混合的多相流動的模擬計算,可以較好地實現對空泡群的捕捉和空泡均勻度的刻畫,在對空化問題的研究分析中已被廣泛應用[25-26]。

3.1.3 空化模型

基于Rayleigh-Plesset方程推導來的Zwart-Gerber-Belamri(ZGB)模型[27]和Schneer-Sauer(SS)模型[28]是較為常用的空化模型,2種模型對空泡內蒸發凝結的傳質過程和不可凝結氣體的計算方式不同。其中,ZGB模型充分考慮了空泡密度對蒸發項的影響,且假設空泡的初始半徑相同;SS模型雖然也將空泡密度考慮在內,但為定值,且并未考慮不可凝結氣體對空化流動的影響。在研究中,筆者重點考察不同空化狀態下空化演變的準周期特性,對空泡密度的變化較為敏感,因此,SS模型并不適用,故選擇ZGB空化模型。

3.2 邊界條件

為捕捉空泡形態和演化行為,模擬計算過程中采用隱式非穩態求解,并參照Courant & Friedrichs & Lewy(CFL)計算收斂準則設置時間步長為1×10-6s[29-30]。設置液相為25 ℃水,密度為999.19 kg/m3,動力黏度為0.001139 Pa·s,氣相為理想狀態水蒸氣。兩凸臺間縫寬為3 mm,無滑移固體壁面。離散項采用控制容積積分法,壓力差分格式為PRESTO!格式,動量、體積分數、湍動能、湍流耗散率都采用高階QUICK格式。湍流強度為0.5%,氣相體積分數初始值為0,飽和蒸汽壓為3169 Pa[2]。

3.3 網格無關性驗證

在對三維非定常空化特性的模擬計算時網格數量會影響計算精度,筆者采用了3套網格進行了網格無關性驗證。從兩凸臺處開始加密布置網格,遠離凸臺網格逐漸變疏,為保證數值傳遞的穩定性,相鄰網格漸變率不超過1.15,可根據空泡初生半徑的大小[31]設置第一層網格厚度為1×10-5m,試算并保證近壁網絡第一層厚度與邊界層厚度之比(y+)≤5。空化室展向方向初始布置30個節點,得到第1套網格數量為200860個,參照Ji等[32]的研究結果通過增加展向節點數加密所得第2套和第3套網格數量分別為412570個和630650個。

設置入口壓力0.43 MPa、出口壓力0.1 MPa,圖3為非穩態求解下同一時刻3套網格的Q判據[33-34]渦結構(閾值Q為12000 s-2)、速度矢量和空泡氣相體積分數(αv,表征空泡內氣體的體積分數)的云圖。由圖3可知,3套網格下速度分布相似,但隨著網格數量增加,對空泡捕捉和渦結構刻畫能力增強,其中圖3(b)和(c)相比圖3(a)都較為完整地刻畫出了空泡形態和渦旋結構。圖4為同一時刻下兩凸臺中心軸線(如圖1(b)中所示)上的壓力分布。由圖4可以看出,第2套網格和第3套網格壓力分布相似,而第1套網格的結果偏差較大,考慮到計算的經濟性,后續求解中都采用第2套網格,即412570個。

ν—Mixture velocity; αv—Gas volume fraction圖3 不同網格數下流場云圖Fig.3 Cloud chart of flow field under different grid numbersGrid number: (a) 200860; (b) 412570; (c) 630650

pw—Near wall pressure; x—Axis length of boss圖4 不同網格數下兩凸臺軸線上壓力分布Fig.4 Pressure distribution on the axis of two bosses under different grid numbers(a) Upper boss; (b) Lower boss

4 空泡演化特征

以空泡演化準周期特征最為突出的云狀空化為例,詳細分析空泡演化細節和流場波動的頻譜特征,并以此為基礎分析探究不同空化狀態下的準周期演化規律。在此云狀空化狀態下的進口壓力設置為0.6 MPa,出口壓力設置為0.1 MPa。

4.1 云狀空化過程

圖5為空化室一個周期內云狀空化演變過程,可觀察到脫落空泡呈現橢球狀,融合空泡團呈現云狀,且空泡以下凸臺為主要初生依附位置。圖5(a)~(d)對應從空化初生到最劇烈空化狀態,而后圖5(e)~(h)對應空化程度逐漸減弱。整個過程呈現出顯著的準周期特性,并測得周期約為74 ms。

圖5 模擬和實驗空化過程Fig.5 Simulations and experiments of cavitation processCavitation moments and characteristics: (a) t=t0, Cavitation inception; (b) t=(t0+10) ms, Growth; (c) t=(t0+21) ms, Shedding and fusion; (d) t=(t0+32) ms, Internal faults; (e) t=(t0+43) ms, Separation; (f) t=(t0+54) ms; Collapse; (g) t=(t0+64) ms, Extinction; (h) t=(t0+74) ms, Re-inception

具體演化過程為由圖5(a)~圖5(b)出現空泡的初生依附并逐漸長大。由圖5(c)開始,空泡不斷生長,空泡的體積和氣相體積分數持續增大,并逐漸發展為主體大空泡。當主體空泡達到足夠體積后開始脫落、融合,并逐漸向下游移動進而匯聚成大尺度的空泡團,同時空泡團內還含氣相體積分數不均一的大量小尺度空泡。從圖5(e)~(h)空化程度逐漸變弱,在圖5(e)中融合空泡團內部開始出現明顯的分離斷層,隨后斷層逐漸擴大,出現空泡的斷裂分離。隨著空化強度減弱,空泡的初生依附處向下游提供的脫落空泡數量減少,下游空泡團體積變小,且內部持續出現較小尺度空泡潰滅與融合的交互作用,小尺度空泡數量減少,氣相體積分數降低,空泡團變得透明。在最后階段內空泡初生依附處已完全不能向下游提供脫落空泡,空泡團內氣相體積分數繼續降低,內部空泡持續潰滅、受壓變小,并最終消失殆盡。

4.2 準周期特性分析

4.2.1 空泡氣相體積分數

圖6為靶區和下游區內模擬計算所得空泡氣相體積分數隨時間變化曲線。由圖6可知,在靶區和下游區內的氣相體積分數呈現出團簇狀的準周期性變化,這與云團狀的大尺度空泡團的演化形態相吻合,且兩區域內的氣相體積分數變化周期基本相等約74 ms,但二者的變化趨勢相差較大。在靶區內從b1到b6氣相體積分數的變化規律相同,都呈現出團簇狀,且其值也緩慢增大,表明在靶區空間內空泡向下游移動過程中出現匯聚,逐漸形成氣泡群密集的空泡團,但并未出現明顯地空泡融合的現象。在下游區內,與d1位置相比,其他各處的氣相體積分數都稍有減小,且變化曲線由規整變得混亂再變得規整。表明在下游空間內空泡由脫落時的橢球形狀在向下游移動過程中逐漸發生空泡間的依附、黏連,最終使得融合空泡團變得畸形,空泡團內部的氣相體積分數也分布不均勻,而當空泡團內部小尺度空泡完全融合后又逐漸均一。

αv—Gas volume fraction in monitoring point;b1—6—Monitoring points in target area; Tb—Long period of cavitation evolution in target area;d1—6—Monitoring points in downstream area; Td—Long period of cavitation evolution in downstream area圖6 靶區及下游區內氣相體積分數波動變化曲線(長周期)Fig.6 Fluctuation curves of gas volume fraction in target and downstream area (long period)(a) Target area; (b) Downstream area

假定在圖6中所呈現出的大尺度空泡團演變時氣相體積分數變化的周期為長周期即74 ms。圖7所示即為兩監測區域內一個空泡團演變長周期內的氣相體積分數隨時間變化曲線,各點的氣相體積分數變化規律依舊相同。在圖7(a)中,氣相體積分數的變化周期性明顯且擁有幾乎相同的周期約為0.98 ms,如Tb2~Tb6,將此類周期定義為短周期。分析可知,該短周期間接表明了微小空泡的數量密度和體積大小,每個波峰則代表一個微小空泡。結合長周期和短周期可推斷一個大尺度空泡或空泡團內約有近百個微小空泡,且小泡的形態幾乎相同并隨著空泡團向下游移動。此外,靶區各監測點處氣相體積分數波形變化還出現了明顯的延遲,相鄰等距的兩監測點間延遲間隔時間也基本相同,約為0.39 ms,如Δt2~Δt5。而且氣相體積分數的峰值和波寬在向下游移動過程中逐漸變大,表明在大尺度空泡或空泡團內的眾多微小空泡在向下游移動、匯聚過程中逐漸長大。相比其他監測位置,b6監測點處的小空泡內氣相體積分數始終大于零,表明已經出現了微小空泡之間的均勻化匯聚和融合。

圖7(b)為下游區內大尺度空泡團的內部和外部的空泡氣相體積分數變化曲線,可知,在空泡團內部的氣相體積分數值較高且微小空泡匯聚緊密,并擁有與靶區內微小空泡相同的短周期約1 ms,表明該空化狀態下的微小空泡再生短周期相同,而在空泡團外部則較為混亂,表明相鄰空泡團黏附明顯,但靠近下游這種現象又會逐漸消失。

Δt2—5—Evolution delay time of microbubbles between two adjacent monitoring points in target area;Tb2—6—Short period of cavitation evolution in target area; Td4—Short period of cavitation evolution at d4圖7 靶區及下游區內氣相體積分數波動變化曲線(短周期)Fig.7 Fluctuation curves of gas volume fraction (αv) in target and downstream area (short period)(a) Target area; (b) Downstream area

4.2.2 流場內壓力波動

圖8所示為靶區和下游區內各監測點的壓力波動曲線,與圖6中空泡的氣相體積分數變化規律相吻合,壓力的波動變化也呈現出團簇狀,且間隔長周期也約為74 ms。圖9為兩監測區域內氣相體積分數、壓力波動和出口質量流量隨時間波動變化曲線。由圖8 和圖9可以得出:各點處壓力團簇狀波動時刻正好對應空泡氣相體積分數為零,且此時出口流量為最大;當壓力基本保持不變時對應空泡氣相體積分數最大,此時出口流量較低。表明了空泡的潰滅消亡會引起局部高壓脈動,并使得壓力波動呈現出周期特性,而空泡團的脫落演化也導致出口流量呈現周期性變化。

p—Pressure in monitoring point圖8 靶區和下游區內壓力波動變化曲線Fig.8 Pressure fluctuation curves in target and downstream area(a) Target area; (b) Downstream area

圖9 氣相體積分數(αv)、壓力(p)和波動(Qm,out)與出口質量流量對應規律曲線Fig.9 Corresponding law curves of gas volume fraction (αv),pressure (p) fluctuation and outlet mass flow rate (Qm,out)

圖10為靶區和下游區內各監測點處的壓力波動進行快速傅里葉變換后求取的功率譜密度(PSD)。由圖10可知,在兩區域內都擁有第一主頻率14 Hz,并且在靶區內還存在第二主頻率1020 Hz。由上文分析可知:空泡的演化行為與流場壓力波動相對應,則第一主頻率14 Hz對應大尺度空泡團的融合、脫落演化過程,與其演化周期74 ms相對應;而第二主頻率1020 Hz對應微小空泡的再生演化過程,周期為0.98 ms,因靶區內微小空泡數量較多,黏連依附較少,演化規律性強,所以同時呈現出了第二主頻率。

PSD—Power spectral density; f—Pressure fluctuation frequency圖10 壓力波動的功率譜密度Fig.10 Power spectral density of pressure fluctuation(a) Target area; (b) Downstream area

4.3 不同空化狀態下頻譜特征及規律

為研究分析空化室內各類空化狀態的準周期特性,參照上文云狀空化分析方法,實驗過程中通過改變入口流量調節空化狀態為無空化、泡狀空化、云狀空化和超空化,并設置對應各自空化狀態下的進、出口壓力,模擬計算不同空化狀態下靶區內的空泡團脫落頻率(對應長周期)和微小空泡再生頻率(對應短周期)。圖11為四類典型空化狀態下靶區內氣相體積分數和壓力波動變化曲線及實驗和模擬所得空泡形態。由圖11可知,同一空化狀態下的氣相體積分數和壓力波動變化特征相對應,但不同空化狀態之間差異明顯。在無空化時,即空化室內無空泡存在,氣相體積分數為零,且壓力波動僅受流道結構影響,波動較緩;在泡狀空化時,空泡可見且孤立存在而并不匯聚,高峰值氣相體積分數的空泡隨機出現且壓力波動出現極峰,波動程度也加劇;在云狀空化時,空泡尺度較大且可見并大面積匯聚成團,氣相體積分數和壓力波動呈團簇狀,且波動變化都較劇烈;在超空化時,空泡尺度很小呈現出霧狀,幾乎不可見,但微小空泡充滿流場,氣相體積分數和壓力波動變化十分劇烈且已無明顯規律。分析可知,不同空化狀態下氣相體積分數和壓力波動的分布變化差異的原因是空泡的形態及周期性演化行為特征的不同所造成的。

(1)Monitoring fluctuation curves; (2)Experimental results; (3)Simulation results圖11 各空化狀態下氣相體積分數(αv)和壓力(p)波動變化曲線及空泡形態Fig.11 Gas volume fraction (αv), pressure (ρ) fluctuation curve and cavitation shape under different cavitation statesCavitation states: (a) Non cavitation; (b) Bubble cavitation; (c) Cloud cavitation; (d) Super cavitation

圖12為不同空化狀態下的兩類空泡演化頻率變化曲線。依據對空泡團和微小空泡的演化特征分析,并結合圖11發現在不同的空化狀態內兩類頻率差異較大。在無空化狀態下無高頻存在,其中的空泡團脫落頻率指代流道結構的固有頻率約25 Hz。進入泡狀空化后再生頻率基本穩定在110 Hz(因該空化狀態下只有較大尺度孤立空泡的脫落演化,比微小空泡大而又不至于匯聚成團,則該狀態下空泡演化的頻率既可稱脫落頻率也可稱再生頻率)。如前文分析,在云狀空化內呈現出兩類頻率,且脫落頻率穩定在14 Hz,而再生頻率波動較大但穩定在900 Hz。當進入超空化狀態后,脫落頻率消失,再生頻率波動變小仍穩定在900 Hz。

f1—Vacuole shedding frequency;f2—Microbubble regeneration frequency;Q—Inlet flow of cavitation chamber圖12 不同空化狀態下兩類空泡演化頻率的變化規律曲線Fig.12 Evolution frequency curves of two kinds of cavitation under different cavitation states

仔細分析兩類頻率的變化規律發現,當空化從無到有,脫落頻率消失是因為游離空泡還未匯聚成團,無法監測大尺度空泡團的脫落演化,而在無空化狀態時是由于流道固有頻率使得流場內壓力出現低頻波動。當空化發生后,流態轉變,空化帶來的劇烈波動因素掩蓋了流道固有結構所引起的壓力波動。在泡狀空化內,小尺度空泡陸續出現掩蓋了類似云狀空化或超空化狀態內產生的微小空泡,該階段空泡尺度稍大,這也是在泡狀空化內所呈現出的和超空化內微小空泡再生頻率不同的主要原因。在云狀空化狀態內,較小尺度空泡或微小空泡匯聚成團,導致了脫落頻率的出現。而當進入超空化時,只有呈霧狀的微小空泡存在,且小尺度空泡相互混合趨于均勻化,使得在超空化狀態下隨著流量增大,微小空泡的再生頻率逐漸趨于穩定。

5 結 論

結合動態實驗拍攝,以云狀空化過程為例,采用數值模擬方法詳細分析了以空泡演變、氣相體積分數、壓力波動和頻譜特征為依據的空化準周期特性,并據此對比分析了不同空化狀態下的規律特性,得出如下結論:

(1)云狀空化狀態下,空泡自初生處長大脫落匯聚成團,氣相體積分數和壓力波動變化曲線呈簇狀且周期性顯著;空化室內呈現出空泡團的脫落和內部微小空泡的再生、匯聚現象,并擁有各自的演化頻率。

(2)不同空化狀態內氣相體積分數、壓力波動變化和空泡形態各有特點,并與各自的空泡演化行為相關聯,且都符合空泡消亡時壓力突增的規律特性。不同空化狀態內兩類空泡演化頻率的變化規律揭示了準周期演化過程是以微小空泡的再生發展為基礎,空化程度越高,微小空泡的主導地位越明顯,且大尺度空泡團也是由微小空泡的高頻初生、脫落、致密匯聚、長大融合所形成的。

主站蜘蛛池模板: 操国产美女| 国产精品专区第一页在线观看| 亚洲区视频在线观看| 亚洲第一综合天堂另类专| 中文字幕啪啪| 久久 午夜福利 张柏芝| 99久久人妻精品免费二区| 国产欧美高清| 亚洲一区毛片| 午夜福利亚洲精品| 亚洲第一网站男人都懂| 国产高清色视频免费看的网址| 久久久久久尹人网香蕉| 国产91在线|中文| 国产免费a级片| 一级一级一片免费| 新SSS无码手机在线观看| 亚洲欧洲日韩综合| 妇女自拍偷自拍亚洲精品| 亚洲床戏一区| 妇女自拍偷自拍亚洲精品| 日本精品视频一区二区| 免费中文字幕一级毛片| 性做久久久久久久免费看| 精品国产三级在线观看| 欧美成人精品在线| 日韩精品免费一线在线观看| 青青热久免费精品视频6| 玖玖精品视频在线观看| 亚洲成人手机在线| 成年人久久黄色网站| 国产精品毛片一区| 国产美女丝袜高潮| аⅴ资源中文在线天堂| 老司机久久99久久精品播放| 国产乱视频网站| 在线一级毛片| 国产欧美日韩在线一区| 亚洲精品桃花岛av在线| 激情影院内射美女| 97在线视频免费观看| 在线看免费无码av天堂的| 波多野结衣一区二区三区四区| 在线播放真实国产乱子伦| 高清色本在线www| 老色鬼久久亚洲AV综合| 国产乱人伦精品一区二区| 超清无码熟妇人妻AV在线绿巨人 | 国产无码制服丝袜| 欧美亚洲国产精品第一页| 色综合网址| 久久特级毛片| 亚洲欧州色色免费AV| 自拍中文字幕| 国产成人综合日韩精品无码首页| 五月天香蕉视频国产亚| 成人一区在线| 最新精品国偷自产在线| 亚亚洲乱码一二三四区| 久久免费观看视频| 日韩成人免费网站| 亚洲香蕉在线| 久久精品国产亚洲麻豆| 国产亚洲精品资源在线26u| 99视频全部免费| 国产亚洲精品无码专| 亚洲三级色| 日韩成人在线网站| 久操线在视频在线观看| 国产十八禁在线观看免费| 在线播放真实国产乱子伦| 国产一级二级在线观看| 露脸国产精品自产在线播| 日韩在线1| 国产高清国内精品福利| 欧美精品啪啪| 国产精品熟女亚洲AV麻豆| 欧美精品综合视频一区二区| 无码中文字幕乱码免费2| 91人妻日韩人妻无码专区精品| 全部免费特黄特色大片视频| 精品国产一区91在线|