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

基于流體-固體耦合的液壓滑閥開啟過程仿真分析與優化

2017-04-11 06:47:02張勝王強何曉暉高亞明王新文
兵工學報 2017年3期
關鍵詞:優化

張勝, 王強, 何曉暉, 高亞明, 王新文

(解放軍理工大學 野戰工程學院, 江蘇 南京 210007)

基于流體-固體耦合的液壓滑閥開啟過程仿真分析與優化

張勝, 王強, 何曉暉, 高亞明, 王新文

(解放軍理工大學 野戰工程學院, 江蘇 南京 210007)

為改善液壓滑閥開啟過程中閥芯的工作性能,基于流體- 固體耦合理論和動網格技術,建立滑閥開啟過程三維流體- 固體耦合數值分析模型,對閥芯在驅動力、彈簧力及液動力作用下的動作過程進行了模擬分析。針對閥芯受力分析結果,改進滑閥內部流道結構;以改進后閥芯槽口和凸臺的4個結構參數為優化對象,利用最小二乘擬合和反向傳播神經網絡構建最大馮米塞斯應力與液動力峰值的目標函數,借助遺傳算法確定了槽口和凸臺的4個結構參數最優值。研究結果表明:優化后閥芯所受液動力和最大馮米塞斯應力的峰值分別減小了16.3%和22.0%;優化設計閥芯的結構參數可明顯提高滑閥開啟性能。

流體力學; 液壓滑閥; 流體- 固體耦合; 動網格; 遺傳算法

0 引言

滑閥結構被液壓系統的流量、方向和壓力控制閥廣泛采用,其開啟過程中流場不穩定,會導致閥芯控制精度不高。閥芯開啟過程的動力學特征直接影響液壓滑閥的控制性能,對于液壓系統有著較大的影響。考慮流體- 固體耦合的研究可追溯到19世紀初的機翼氣動彈性問題[1]。目前,流體- 固體耦合的研究主要集中于水輪機、冷卻器蓋[2-4]等,針對液壓閥的研究并不多。而對于液壓閥閥芯的動態仿真大多是直接給定閥芯的開口度或閥芯的運動速度,并未考慮閥芯開啟過程的瞬態受力情況和流體- 固體耦合的影響[5-6]。因此,有必要對液壓滑閥的閥芯及其開啟過程進行重點研究,分析其瞬態受力情況并加以優化。針對上述問題,借助流體- 固體耦合理論和動網格技術,以流體力學軟件CFD-ACE+為平臺,對液壓滑閥開啟過程進行流體- 固體耦合特性分析;針對分析結果,改進液壓滑閥的內部流道結構并進行優化;建立優化后的模型,對優化結果進行仿真驗證,并比較優化前后液壓滑閥的開啟性能。

1 幾何模型和網格劃分

1.1 幾何模型

液壓滑閥通常由多個閥腔組成,基于相似性原理,用其中的一個閥腔就能說明滑閥內部流體的流動狀態。借鑒文獻[7-9]的簡化方法,歸納出具有共性的液壓滑閥作為研究對象,其結構如圖1所示。

圖1 液壓滑閥結構Fig.1 Structure diagram of hydraulic spool valve

在該滑閥閥芯的開啟過程中,液壓油的流動具有對稱性,為提高計算效率,取流動區域和閥芯的一半作為研究對象,建立用于流體- 固體耦合分析的三維幾何模型,如圖2所示。

圖2 液壓滑閥幾何模型Fig.2 Geometric model of hydraulic spool valve

1.2 網格劃分

流體- 固體耦合會導致網格變形較大,出現三維瞬態流場解析發散的問題,單一的網格變形方法很難同時具有較大的網格變形能力和得到較高的變形質量。為此,本文將彈簧法和超限插值法相結合,解決結構網格變形問題。基于網格分塊的思想,將計算區域分成16個網格塊;采用彈簧連接不同網格塊的角點,通過彈簧驅動各角點的變形;對于網格塊內部節點,采用超限插值法將角點位移通過插值函數分布到網格塊內部節點[10];網格劃分采用前處理軟件GEOM實現;為更好地捕捉閥口處的流場變化,對閥口區域網格進行細化。在大量試算的基礎上,最終確定網格總數為95 252個。其中,流域網格51 029個,閥芯網格44 223個,如圖3所示。

圖3 網格劃分Fig.3 Grid division

2 控制方程

2.1 動網格條件下的流場控制方程

在網格運動過程中,流場內控制體將發生改變,一般采用積分形式來表示諸守恒方程:

(1)

2.2 流體- 固體耦合控制方程

在流體- 固體耦合交界面處應力:

(2)

式中:ΔFi為在i方向上的力(N);ΔAi為在i方向的受力面積(m2)。

力平衡方程

(3)

式中:σxx、σyy、σzz為應力分量;τxy、τxz、τyx、τyz、τzx、τzy為應變分量。

幾何方程

(4)

式中:ξxx、ξyy、ξzz為應變分量;u、v、w為x、y、z方向的速度分量;γxy、γyx、γzy、γyz、γxz、γzx為切應變分量。

本構方程

(5)

式中:λ為泊松比;E為楊氏模量(GPa)。

3 邊界條件與參數選取

3.1 流域參數

入口壓力取6.00 MPa,出口壓力取0.15 MPa;液壓油型號VG46,密度ρ=885 kg/m3,運動黏度ν=46 mm2/s,由于液壓油壓縮性很小,視為不可壓縮流體;彈簧剛度80 N/mm;滑閥閥芯兩端驅動壓差1 MPa;采用k-ε湍流模型。

3.2 閥芯參數

閥芯材料為45號鋼。密度ρs=7 850 kg/m3;楊氏模量Es=210 GPa;泊松比λs=0.3.

4 仿真結果分析與優化

4.1 閥芯運動情況分析

閥芯的運動狀態由驅動力、彈簧力及液動力共同決定,仿真得到的閥芯運動情況如圖4所示。圖中點1為閥芯開啟初期(t=0.42 ms),點2為閥芯位移最大時(t=1.80 ms),點3為閥芯位移達到穩態時(t=47.00 ms)。

圖4 閥芯開啟過程位移變化Fig.4 Displacement change of spool during opening

從圖4可知,液壓滑閥開啟過程中,閥芯位移存在一定的波動,最大位移達到1.49 mm,穩態位移為0.78 mm. 閥芯位移的波動會引起流量的波動,從而影響執行元件的運動狀態;位移的波動也會產生振動和噪聲,影響流場及閥芯瞬態受力情況。

液壓滑閥閥芯穩態位移和穩態液動力的理論計算公式[11]為

(6)

式中:D為流道直徑,D=16 mm;xv為閥芯穩態位移(mm);W為閥口的面積梯度(mm2/mm);kt為彈簧剛度,kt=80 N/mm;F為滑閥閥芯所受的驅動力(N);Fs為閥芯所受的穩態液動力(N);Δp為滑閥進出口壓差,Δp=5.85 MPa;Δpq為驅動壓差,Δpq=1.00 MPa.

聯立(6)式,可以解得閥芯穩態位移xv=0.71 mm,穩態液動力Fs=-43.73 N.

仿真得到的閥芯穩態位移xvf=0.78 mm,穩態液動力Fsf=-47.96 N. 比較理論計算與仿真計算結果,閥芯穩態位移和穩態液動力的誤差分別為9.9%和9.7%,可以認為仿真計算結果與理論計算結果基本吻合,證明了仿真的可靠性。誤差的產生主要是因為在理論計算中穩態液動力采用經驗公式求解,沒有考慮內部結構不同及流場微觀流動帶來的影響。

4.2 流體- 固體耦合分析

圖5 不同時刻流體- 固體耦合分析Fig.5 Fluid-structure interaction analysis at different moments

針對閥芯開啟初期(t=0.42 ms,圖4中點1),最大位移(t=1.8 ms,圖4中點2),穩態位移(t=47 ms,圖4中點3),三維流體- 固體耦合仿真結果如圖5所示,圖中給出了液壓滑閥流域的壓力場、壓力等值線分布以及液壓滑閥閥芯所受的馮米塞斯應力分布。其中,馮米塞斯應力是根據第四強度理論得到的一種當量應力,考慮了第1、第2和第3主應力。馮米塞斯應力值越大的地方材料變形就越大,往往也是越容易受到破壞的地方,因而可以用來對疲勞、破壞等進行評價[12]。

從圖5可以看出:

1)在流場的A區域和B區域會產生負壓,負壓區域面積隨著閥芯位移的增大而變大。流場負壓區是引起氣蝕的主要原因,氣蝕會導致閥芯產生氣蝕破壞;

2)在液壓滑閥開啟初期,閥芯應力集中于閥芯的1、2、3和4區域且四周對稱分布。當閥口位移較大后,應力集中于3和4區域。

閥芯在開啟過程中所受最大馮米塞斯應力的大小和位置隨時間變化情況如圖6所示。由圖6可看出,最大馮米塞斯應力主要集中在1、2、3、4區域,在閥芯開啟瞬間達到峰值。

圖6 不同時刻最大馮米塞斯應力分析Fig.6 Analysis of maximum von Mises stresses at different moments

4.3 液壓滑閥結構優化分析

為改善滑閥開啟過程中閥芯受力情況,需減小閥芯所受的液動力(包括穩態液動力與瞬態液動力)和因流體- 固體耦合作用所產生的最大馮米塞斯應力。現改進滑閥內部流道結構,在閥芯上增設凸臺和槽口,新結構如圖7所示。

圖7 新結構滑閥Fig.7 New structure of spool valve

圖8分別給出了原結構、僅有凸臺、僅有槽口和新結構閥芯在開啟過程中最大馮米塞斯應力和液動力的變化情況。其中,液動力是通過求解閥芯各壁面軸向力的合力得到[13],為瞬態液動力和穩態液動力之和,本文將其作為其中一個優化指標。從圖8中可看出,閥芯凸臺結構能有效地減小液動力,槽口則對最大馮米塞斯應力減小較為明顯,而同時增設槽口和凸臺后,最大馮米塞斯應力和液動力均減小。

圖8 結構變化對閥芯受力影響Fig.8 Effect of structure change on force of spool

經大量試算,選取如圖7所示的L1、L2、L3、L4的尺寸為優化對象,以滑閥開啟過程中不同瞬時閥芯所受最大馮米塞斯應力和液動力的峰值為目標函數,建立33個滑閥試驗模型,進行三維流體- 固體耦合仿真分析,結果如表1所示。為了便于表達,實驗因子L1、L2、L3、L4分別用x1、x2、x3、x4表示,評價指標σmax、Fsmax分別用y1、y2表示。其中,1.0 mm≤x1≤5.0 mm,0.1 mm≤x2≤2.0 mm,0.1 mm≤x3≤2.0 mm,0.1 mm≤x4≤2.0 mm.

表1 液壓滑閥開啟過程實驗設計表Tab.1 Test design of opening process of spool valve

優化閥芯結構要求最大馮米塞斯應力和液動力均最小,取(7)式作為優化目標函數:

y=wvfv+wtft,

(7)

式中:fv、ft分別為最大馮米塞斯應力和液動力峰值的隸屬度(無因次量);wv、wt為對應權重,均取0.5. 隸屬度函數[14]按(8)式計算:

(8)

式中:si、so和sw分別為表1中某一樣本的指標值、該指標的最優值和最差值,i=1,2,…,33.

為取得每個指標的最優值和最差值,采用2階函數分別對最大馮米塞斯應力和液動力峰值的仿真結果進行最小二乘擬合,擬合值和期望值如圖9所示。擬合方程分別為

圖9 最大馮米塞斯應力和液動力峰值最小二乘擬合情況Fig.9 Least square fitting of peak value of maximum von Mises stress and transient flow force

擬合方程y1和y2評價標準為誤差決定系數[15],檢驗誤差決定系數為

(9)

利用遺傳算法[16]對擬合出的兩個方程求解最值。在處理約束時,采用罰函數的方法,將其轉化為無約束優化問題。罰函數具體操作是:在約束條件內,適應度值高;約束條件外,適應度值低,容易被淘汰。得出最大馮米塞斯應力峰值的最優值和最差值分別為68 MPa、130 MPa,液動力峰值的最優值和最差值分別為-68 N、-101 N. 根據(7)式計算得到目標函數(見表1)。

由于最小二乘擬合各設計因子與目標函數之間的關系時效果不太理想,采用反向傳播(BP)神經網絡構建其響應面模型。選取表1中27組數據作為訓練數據,其余6組與預測結果進行對比分析誤差,擬合情況如圖10所示。

圖10 目標函數BP神經網絡模型預測結果Fig.10 Prediction of BP neural network for objective function

從圖10可看出,預測誤差絕對值在8%以內,因此構建的BP神經網絡預報精度完全可以接受。再次借助遺傳算法尋求目標函數的最優值,尋優適應度收斂過程曲線如圖11所示。最終得到優化結果如表2所示。

圖11 遺傳算法適應度收斂曲線Fig.11 Fitness convergence curves of GA

表2 閥芯優化結果Tab.2 Optimimized results of spool

根據尋優結果,選取L1、L2、L3、L4的尺寸分別為3.1 mm、0.9 mm、1.6 mm、0.9 mm,建立優化后的三維流體- 固體耦合模型,得出相同設定條件下閥芯所受的最大馮米塞斯應力和液動力變化情況,與優化前的結果對比,如圖12所示。原結構、優化設計及仿真結果對比如表3所示。

圖12 滑閥開啟過程優化前后對比情況Fig.12 Comparison of results during opening before and after optimization

表3 閥芯優化設計驗證對比Tab.3 Verification and comparison of optimization design for spool

由圖12和表3可知,優化后的閥芯在開啟過程中不同瞬時受到最大馮米塞斯應力的峰值為94.7 MPa,減小了22.0%,液動力峰值為83.2 N,減小了16.3%,與優化設計結果較為接近,但是均比優化前的值要小,且在整個開啟過程中最大馮米塞斯應力變化較為平穩,明顯改善了滑閥開啟性能。

5 結論

1)液壓滑閥開啟過程中流域會產生負壓區,負壓區主要集中在閥口處靠近閥芯和閥座的兩個區域,且隨閥芯位移的增加而擴散;

2)液壓滑閥開啟過程中不同位置閥芯所受應力是不斷變化的,且在閥芯開啟初期最大馮米塞斯應力達到峰值;

3)液壓滑閥開啟過程中,閥芯所受最大馮米塞斯應力值主要集中于4個截面突變區域。改進閥芯結構,可減小這些區域的最大馮米塞斯應力,同時也能減小液動力;

4)內部流道結構的優化設計值與仿真結果基本吻合,證明仿真分析方法和優化方法具有可行性。

References)

[1] 邢景棠, 周盛, 崔爾杰. 流固耦合力學概述[J]. 力學進展, 1997, 27(1):19-37. XING Jing-tang, ZHOU Sheng, CUI Er-jie. Overview of fluid solid coupling mechanics[J]. Advances in Mechanics, 1997, 27(1):19-37.(in Chinese)

[2] Saeed R A, Galybin A N, Popov V. 3D fluid-structure modelling and vibration analysis for fault diagnosis of Francis turbine using multiple ANN and multiple ANFIS[J]. Mechanical Systems and Signal Processing, 2013, 34(1): 259-276.

[3] 胡丹梅, 張志超, 孫凱, 等. 風力機葉片流體- 固體耦合計算分析[J]. 中國電機工程學報, 2013, 33(17):98-104. HU Dan-mei, ZHANG Zhi-chao, SUN Kai, et al. Computational analysis of wind turbine blades based on fluid-structure interaction[J]. Proceedings of the CSEE, 2013, 33(17):98-104.(in Chinese)

[4] 張俊紅, 郭遷, 王健, 等. 基于流固耦合與多目標拓撲優化的噪聲塑料機油冷卻器蓋優化設計[J]. 振動與沖擊, 2016, 35(7):186-191. ZHANG Jun-hong, GUO Qian, WANG Jian, et al. Low noise optimization design of plastic oil cooler cover based on liquid-solid coupling and multi-objective topological optimization[J]. Journal of Vibration and Shock, 2016, 35(7):186-191.(in Chinese)

[5] 郭熛, 解寧, 郭津津, 等. 滑閥液動力研究及結構分析[J]. 液壓氣動與密封, 2012, 32(4):11-15. GUO Piao, XIE Ning, GUO Jin-jin, et al. Research on flow force of the spool valve and structure analysis[J]. Hydraulics Pneumatics & Seals, 2012, 32(4):11-15. (in Chinese)

[6] 郭昌盛. 基于動網格的調壓閥流場模擬及閥芯運動分析[D].重慶: 重慶大學,2013. GUO Chang-sheng. Moving mesh based on flow field simulation of pressure regulating valve and spool motion analysis[D]. Chongqing: Chongqing University, 2013.(in Chinese)

[7] 趙蕾, 陳青, 權龍. 閥芯運動狀態滑閥內部流場的可視化分析[J]. 農業機械學報, 2008, 39(11):142-145. ZHAO Lei, CHEN Qing, QUAN Long. Visualization analysis of internal flow field of a moving spool valve[J]. Transactions of the Chinese Society for Agricultural Machinery, 2008, 39(11):142-145.(in Chinese)

[8] 王安麟, 吳小鋒, 周成林, 等. 基于CFD的液壓滑閥多學科優化設計[J]. 上海交通大學學報, 2010, 44(12):1767-1772. WANG An-lin, WU Xiao-feng, ZHOU Cheng-lin, et al. Multidisciplinary optimization of a hydraulic slide valve based on CFD[J]. Journal of Shanghai Jiao Tong University, 2010, 44(12):1767-1772. (in Chinese)

[9] Song X G, Wang L, Park Y C. Analysis and optimization of butterfly valve disc[J]. Journal of Process Mechanical Engineering, 2009, 223(2):81-89.

[10] 張來平, 鄧小剛, 張涵信. 動網格生成技術及非定常計算方法進展綜述[J]. 力學進展, 2010, 40(4):424-447. ZHANG Lai-ping, DENG Xiao-gang, ZHANG Han-xin. Review on the development of the dynamic mesh generation technique and the unsteady computation method[J]. Advances in Mechanics, 2010, 40(4):424-447. (in Chinese)

[11] 房小立. 基于CFD有限元法的滑閥穩態液動力研究[D]. 武漢:武漢科技大學, 2012. FANG Xiao-li. Research on steady flow force of slide valve based on CFD finite element method[D]. Wuhan: Wuhan University of Science and Technology, 2012. (in Chinese)

[12] 金奕山, 李琳. 隨機振動結構von Mises應力過程峰值概率密度函數的研究[J]. 應用力學學報, 2006, 23(4):645-648, 701. JIN Yi-shan, LI Lin. Research on the peak probability density function of von Mises stress process in random vibration structure[J]. Chinese Journal of Applied Mechanics, 2006, 23(4):645-648, 701.(in Chinese)

[13] 張曉俊, 權龍, 趙斌. 內流式滑閥壁面壓力分布可視化計算及試驗驗證[J]. 機械工程學報, 2016, 52(14):196-203. ZHANG Xiao-jun, QUAN Long, ZHAO Bin. Visualization calculation and experimental verification of steady axial flow force on the converged flow sliding valve[J]. Journal of Mechanical Engineering, 2016, 52(14):196-203. (in Chinese)

[14] 張俊紅, 郭遷, 王健, 等. 基于流固耦合與多目標拓撲優化的噪聲塑料機油冷卻器蓋優化設計[J]. 振動與沖擊, 2016, 35(7):186-191. ZHANG Jun-hong, GUO Qian, WANG Jian, et al. Low noise optimization design of plastic oil cooler cover based on liquid-solid coupling and multi-objective topological optimization[J]. Journal of Vibration and Shock, 2016, 35(7):186-191. (in Chinese)

[15] 吳小鋒, 干為民, 劉春節, 等. 液壓換向滑閥內部結構的健壯性設計[J]. 中國機械工程, 2015, 26(15):2030-2035, 2040. WU Xiao-feng, GAN Wei-min, LIU Chun-jie, et al. Robust design of hydraulic slide valve internal structure[J]. China Mechanical Engineering, 2015, 26(15):2030-2035, 2040.(in Chinese)

[16] 陳崢, 潘二東, 申江衛, 等. 基于遺傳算法的PHEV動力參數優化與仿真[J]. 昆明理工大學學報:自然科學版, 2016, 41(1):46-51. CHEN Zheng, PAN Er-dong, SHEN Jiang-wei, et al. Optimization and simulation of PHEV dynamic parameters based on genetic algorithm[J]. Journal of Kunming University of Science and Technology: Natural Science Edition, 2016, 41(1):46-51. (in Chinese)

Simulation and Optimization of Hydraulic Spool Valve Opening Process Based on Fluid-structure Interaction

ZHANG Sheng, WANG Qiang, HE Xiao-hui, GAO Ya-ming, WANG Xin-wen

(College of Field Engineering, the PLA University of Science and Technology, Nanjing 210007, Jiangsu, China)

To improve the working capability of the spool during the opening process of hydraulic spool valve, a three-dimensional numerical model is established based on fluid-structure interaction theory and dynamic mesh technology, and the movement of spool under the action of driving force, spring force and flow force is simulated. The internal flow channel structure is improved in view of stress on spool. The objective functions of peak value of the maximum von Mises stress and flow force are set up through least squares fit and BP neural network, with four structural parameters of improved spool notches and convex being taken as the optimization object. The optimal values of four structure parameters of notches and convex platform are determined by using genetic algorithm. After optimization, the maximum von Mises stress is decreased by 22.0%, and the maximum flow force is decreased by 16.3%. Results show that the spool valve opening performance is obviously improved through the study of the optimization design of structural parameters of spool.

fluid mechanics; hydraulic spool valve; fluid-structure interaction; dynamic mesh; genetic algorithm

2016-06-15

國家自然科學基金項目(51175511)

張勝(1992—), 男, 碩士研究生。E-mail:18205152450@163.com

王強(1964—),男,教授,碩士生導師。E-mail:gcbhxh@163.com

TH137.52+1

A

1000-1093(2017)03-0608-08

10.3969/j.issn.1000-1093.2017.03.026

猜你喜歡
優化
超限高層建筑結構設計與優化思考
房地產導刊(2022年5期)2022-06-01 06:20:14
PEMFC流道的多目標優化
能源工程(2022年1期)2022-03-29 01:06:28
民用建筑防煙排煙設計優化探討
關于優化消防安全告知承諾的一些思考
一道優化題的幾何解法
由“形”啟“數”優化運算——以2021年解析幾何高考題為例
圍繞“地、業、人”優化產業扶貧
今日農業(2020年16期)2020-12-14 15:04:59
事業單位中固定資產會計處理的優化
消費導刊(2018年8期)2018-05-25 13:20:08
4K HDR性能大幅度優化 JVC DLA-X8 18 BC
幾種常見的負載均衡算法的優化
電子制作(2017年20期)2017-04-26 06:57:45
主站蜘蛛池模板: 久久成人免费| 91国内外精品自在线播放| 超薄丝袜足j国产在线视频| 成人免费午间影院在线观看| 国产网站一区二区三区| 国产一区成人| 91探花在线观看国产最新| 日本免费a视频| 国产极品嫩模在线观看91| 中国精品自拍| 国产精品免费露脸视频| 亚洲国产精品成人久久综合影院| 一级黄色欧美| 国产日韩av在线播放| 97精品久久久大香线焦| 欧美a在线看| 免费A∨中文乱码专区| 美女高潮全身流白浆福利区| 欧美日韩北条麻妃一区二区| 国产精品无码久久久久久| 久久黄色毛片| 国产黄网站在线观看| 国产一区在线视频观看| 福利姬国产精品一区在线| 超碰免费91| 国产精品片在线观看手机版| 青青草国产免费国产| 国产在线观看人成激情视频| 91午夜福利在线观看| 国产毛片基地| 国产精品九九视频| 波多野结衣一区二区三区AV| 中文字幕在线看| 全部无卡免费的毛片在线看| 亚洲成a人片7777| av色爱 天堂网| 欧美成a人片在线观看| 国产午夜一级毛片| 手机在线看片不卡中文字幕| 国产美女精品一区二区| 久久五月视频| 九九视频免费看| 99re在线免费视频| 成年人福利视频| 国产在线视频二区| 日韩精品一区二区深田咏美| 国产麻豆精品久久一二三| 日韩欧美综合在线制服| 久久国产黑丝袜视频| 91精品国产91欠久久久久| 欧美成在线视频| 亚洲欧洲日产国码无码av喷潮| 欧美第九页| 一级做a爰片久久毛片毛片| 国语少妇高潮| 麻豆精品在线播放| 国产www网站| 成人在线亚洲| 久久精品亚洲热综合一区二区| 高清不卡毛片| 成人免费一级片| 国产黄色免费看| 999在线免费视频| 91精品专区国产盗摄| 在线亚洲小视频| 亚洲天堂网在线播放| 亚洲三级视频在线观看| 欧美在线国产| av无码久久精品| 她的性爱视频| 蜜臀av性久久久久蜜臀aⅴ麻豆| 国产精品偷伦视频免费观看国产 | 国产精品v欧美| 尤物特级无码毛片免费| 国产又色又爽又黄| 色综合综合网| 国产精品无码久久久久久| 久久精品无码国产一区二区三区| 97久久免费视频| 欧美成人手机在线观看网址| 2024av在线无码中文最新| 国产成人精品一区二区免费看京|