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

拉桿松弛導致預緊失諧組合轉子性能退化研究

2016-01-15 06:09:56蘇永雷,王艾倫,曾海楠
振動與沖擊 2015年20期

第一作者蘇永雷男,碩士生,1989年生

通信作者王艾倫男,教授,博士生導師,1959年12月生

拉桿松弛導致預緊失諧組合轉子性能退化研究

蘇永雷1,2,王艾倫1,曾海楠3

(1.中南大學高性能復雜制造國家重點實驗室,長沙410083; 2.重慶長安汽車股份有限公司(研究總院),重慶401120;3.天津大學機械工程學院,天津300072)

摘要:由結構損傷導致性能退化角度出發,分析拉桿松弛引起的預緊失諧組合轉子性能退化特性。利用跨尺度計算方法獲得輪盤界面接觸剛度并對接觸界面等效處理,提出考慮粗糙界面及預緊的組合轉子動力學建模方法。在此基礎上進行組合轉子動力學分析,獲得組合轉子固有頻率隨時間變化規律;建立以組合轉子頻率相對降低量表征其性能退化量的定量評估方法,對退化數據進行處理并描述組合轉子退化軌跡,實現拉桿松弛引起的組合轉子性能退化定量評估;對初始預緊失諧程度相同的組合轉子進行不平衡響應分析,獲得預緊失諧組合轉子退化特征。結果表明,組合轉子性能退化量參數服從指數函數分布,負失諧對組合轉子不平衡響應影響較大。

關鍵詞:拉桿松弛;組合轉子;退化量;失諧;固有頻率

基金項目:國家973計劃資助項目(2013CB035706)

收稿日期:2014-03-25修改稿收到日期:2014-06-24

中圖分類號:TH113.1

文獻標志碼:A

DOI:10.13465/j.cnki.jvs.2015.20.005

Abstract:A combined rotor is the key part of an engine, and it is important to reveal its degradation mechanism and assess the performance degradation degree for making maintenance during operation. Performance degradation characteristics of the combined rotor were discussed, from the perspective of performance degradation caused by structural damage. The dynamic features of the combined rotor were analyzed, considering the rod relaxation and interface contact effect and the changing of natural frequency along with the time was observed. Then, by defining the degradation degree of combined rotor, the degradation path was described through processing the degradation data and the degradation assessment was completed as well. At the end, unbalance responses of combined rotors with the same mistuned degree of preload were investigated and the degradation characteristics of combined rotors were obtained. The results show that: degradation parameters of combined rotors have exponential distribution, and the negative mistuning has great effect on the unbalance response amplitude of combined rotors.

Performance degradation of mistuned combined rotor considering rod relaxation

SUYong-lei1,2,WANGAi-lun1,ZENGHai-nan3(1. Key Laboratory of High performance and Complex Manufacturing, Changsha 410083, China;2. Chongqing chang’an automobile Co., Ltd., Chongqing 401120, China;3. School of Mechanical Engineering, Tianjin University, Tianjin 300072, China)

Key words:rod relaxation; combined rotor; degradation degree; mistuning; natural frequency

運行中由于蠕變、裂紋、磨損等原因設備性能會逐步退化,不但降低其可靠性更增加其發生故障的可能性。作為長期嚴酷工況運行造成組合轉子性能退化的有力佐證,GE公司10余臺9FA機組連續發生嚴重振動超標,且故障大多出現于機組運行10000小時以上,主要原因為組合轉子間預緊力在長期運行中逐漸缺失導致預緊力不足,在輪盤間發生滑移,改變轉子的應有剛度;該公司另一主打機組LM5000燃氣輪機運行中也出現壓氣機隨機組負荷增加、振動加劇現象[1]。而某臺燃機中間軸及壓氣機24個連接螺栓因其中1個螺栓破壞導致轉子振動過大,導致重大事故發生[2]。

工業實踐中因制造誤差、材料缺陷等因素使轉子系統循環對稱性遭受破壞導致系統失諧。轉子最典型的失諧形式為葉盤系統失諧,主要反映為模態局部化及振動傳遞局部化。結構動力學中此兩種局部化稱為失諧周期結構振動局部化問題,且廣受關注[3]。而對周向拉桿組合轉子則存在另種形式失諧即各拉桿預緊不均勻導致的失諧。若一根拉桿預緊力與其它拉桿不同,則該拉桿在轉子彎曲時的受力亦不同于其它拉桿,因而導致轉子彎曲剛度改變,并產生各向異性。該失諧不同于葉盤的弱耦合結構失諧,而類似裂紋轉子產生的剛度各向異性[4]。

燃氣輪機長期服役于高溫、高轉速工況,對其進行性能退化特性研究、掌握性能退化程度及趨勢對保障燃氣輪機安全運行具有重要意義。Diakunchak[5]認為燃氣輪機性能退化不僅與運行時間相關,亦受工作模式影響。Brook[6]認為燃氣輪機大部分能量損失由壓氣機部分性能退化所致。Kurz等[7]通過對燃氣輪機性能退化原因詳細闡述,研究不同退化成因導致的燃氣輪機系統性能參數變化。據燃氣輪機性能退化特征,Brothertom等[8]提出用浴盆曲線近似描述燃氣輪機性能退化規律。Li等[9]利用氣路分析方法對退化數據進行擬合,獲得燃氣輪機性能退化曲線,并用于預測燃氣輪機剩余壽命。Venturini等[10]基于蒙特卡洛法,提出可靠診斷與預測燃氣輪機退化方法。

現有研究存三方面問題:①燃氣輪機退化研究主要集中在不同退化成因對整機性能參數影響及診斷方式,而對其組合轉子性能退化研究較少。組合轉子作為燃氣輪機能量轉換的核心部件,其退化會對整機性能產生重大影響,但組合轉子性能退化機理并不明確。②其退化量參數所選整機性能參數(功率、燃耗、壓比)不適合作為組合轉子退化量指標。③雖通過建立具有界面開閉效應的集中參數模型研究組合轉子動力學問題,但對性能隨時間退化研究尚少。為此,本文由結構損傷導致性能退化角度出發,研究拉桿松弛導致預緊失諧組合轉子退化特征;定義退化量指標并建立組合轉子退化量與時間的函數關系,揭示拉桿松弛引起的組合轉子性能退化機理。

1組合轉子計算模型

燃氣輪機壓氣機段周向均布拉桿組合轉子簡化模型(簡稱組合轉子)見圖1,該轉子由兩端軸頭、4個中間輪盤及8根拉桿組成。拉桿周向均布于直徑Drod=80 mm圓上,將六級輪盤預緊組成整體,輪盤間具有多個接觸界面,組合轉子材料密度7800 kg/m3,輪盤、轉軸彈性模量210 GPa。組合轉子各物理參數見表1。

表1 組合轉子物理參數

圖1 組合轉子結構示意圖 Fig.1 Structure of compressor combined rotor

2組合轉子界面剛度等效

接觸界面的存在導致組合結構剛度發生變化,結構通過接觸界面進行動力傳遞,不同界面狀態直接影響結構的振動特性[11]。為準確進行組合轉子動力學分析,在獲得粗糙表面接觸剛度基礎上計算接觸剛度與預緊力關系,進行接觸剛度等效。

2.1界面剛度確定

直接在全尺寸模型中研究微米級粗糙度對界面剛度影響尺寸效應會較明顯從而導致分析困難。本文采用具有相同界面粗糙度的微元體模型分析與宏觀尺寸相結合的跨尺度計算方法獲得輪盤間接觸剛度,即①計算輪盤間(圖 1)A接觸界面面積S1并進行受力分析,獲得預緊后界面壓力P;②據文獻[12]方法建立界面面積S2、粗糙度0.4 μm、壓力P作用下微元體模型進行應力應變分析并提取該微元體模型法向接觸剛度kn;③借鑒文獻[13],在獲得微觀有限元模型接觸剛度基礎上通過面積擴展,獲得組合轉子在預緊力F作用下界面法向接觸剛K,見圖2。面積擴展表達式為

F=PS1

(1)

(2)

圖2 組合轉子界面接觸剛度與預緊力的關系圖 Fig.2 Stiffness of interface versus pretension for combined rotor

2.2界面剛度等效

(3)

式中:h為虛擬材料層厚度。

3考慮界面的預緊失諧組合轉子動力學建模

為進行組合轉子諧調、失諧的退化研究,本文提出考慮粗糙界面和預緊情況的組合轉子建模方法建立p根拉桿預緊(其中q根拉桿失諧)的失諧組合轉子動力學模型,見圖3。

圖3 預緊失諧的組合轉子建模 Fig.3 The model of preload mistuning combined rotor

據失諧的一般表達式,定義預緊失諧量為

ηi=(Fdi-Fh)/Fh

(4)

式中:Fh為諧調預緊力,Fih為組合轉子第i根拉桿預緊力。則ηi>0為預緊正失諧,ηi<0為預緊負失諧,ηi=0為預緊諧調。

將第i根拉桿預緊力Fi等效為所有拉桿中最小預緊力F0與附加預緊力ΔFi兩部分。即

Fi=F0+ΔFi

(5)

本文以8根拉桿預緊、1根拉桿發生失諧(P=8,q=1)為例,對失諧組合轉子進行分析。

4拉桿螺栓松弛應力計算

聯接件受蠕變、疲勞等多種作用,但在穩定載荷狀態下最基本損傷機制為恒應變下應力松弛[14]。

松弛中,設ε0為聯接件初始彈性應變且為常數。在絕對溫度為T的環境中經時間t后彈性應變為εe,蠕變應變為εc,則有

ε0=εe+εc=σ/E+εc

(6)

式中:E為彈性模量。

對式(5)求導得

(7)

針對長期工作的拉桿松弛,用 Norton 模型表達蠕變應變速度與應力關系,即

(8)

式中:C1,C2,C3為材料常數。

拉桿初應力σ0符合

(9)

式中:σ0為拉桿初始應力;Δl為預緊量;lrod為總長;Arod為橫截面積;F0為初始預緊力。

結合式(4)、式(7)~式(9)得拉桿松弛剩余應力為

σ={[(1+η)ε0E]1-C2-

(10)

拉桿預緊力F在t時刻符合

F=σArod

(11)

結合式(4)、(5)、(11)得預緊失諧產生的附加預緊力為

(12)

式中:σh為預緊諧調(η=0)時據式(10)所得應力值。

本文燃氣輪機組合轉子拉桿材料采用GH4169,利用蠕變實驗數據確定Norton蠕變模型參數值[15]。其中C1=2.2E-9,C2=10.2,C3=50825。由材料特性[16]可得彈性模量E在不同溫度T的數值。

參照工程中重型燃氣輪機拉桿預緊情況,將初始預緊量定為拉桿總長的千分之一,即ε0=1/1000。據式(9)所得失諧拉桿在873 K溫度下應力松弛曲線見圖4,據式(11)所得預緊失諧量隨時間變化見圖5。由兩圖看出,材料為GH4169的拉桿在同一溫度環境下失諧量隨時間增加逐步減小,即各拉桿預緊力隨時間增加趨于諧調;初始失諧量η*為10%(5%)的拉桿較初始失諧為-10%(-5%)的拉桿失諧量隨時間增加更接近0值,即初始失諧量相當情況下隨時間增加初始正失諧拉桿預緊力失諧程度較負失諧小。

圖4 失諧拉桿在873K溫度下的應力松弛曲線 Fig.4 Stress relaxation curve of mistuning rod

圖5 失諧拉桿在873K溫度下失諧量隨時間的變化 Fig.5 Mistuning degree change with relaxation time at 873K

5組合轉子固有特性分析

拉桿初始預緊不均產生的廣義附加力矩影響轉子剛度。為進一步研究拉桿松弛引起預緊失諧組合轉子退化特征,對初始預緊失諧的組合轉子在簡單彈性支承下進行動力學分析。拉桿諧調預緊值由拉桿伸長量(拉桿長度的1/1000)確定。本文只分析組合轉子彎振特性,不考慮界面切向剛度,虛擬材料層泊松比為0,厚度0.5 mm;轉子工作溫度873 K,轉速3000 r/min。進行動力學分析時將拉桿螺母與輪盤間及虛擬材料層兩端與輪盤接觸界面粘合。

在同一時間段內組合轉子受預緊降低、預緊失諧兩因素耦合作用;通過調節虛擬材料層彈性模量表達組合轉子預緊程度,調節拉桿附加預緊力表達預緊失諧量變化。以初始預緊失諧量-10%、0、10% 為例,在溫度為873 K、轉速3000 r/min 條件下對8根拉桿、6級輪盤組合轉子進行含預應力的模態分析。附加預緊力據式(12)求得,計算所得組合轉子固有頻率見表2。由表2看出,預緊失諧對組合轉子固有頻率影響不大,且隨松弛時間延長而減弱。此因失諧量為10%時預緊失諧絕對量F占拉桿總預緊力的1/80,對界面剛度影響較小;另外,隨時間延長預緊失諧量逐步降低,失諧對組合轉子固有頻率影響降低。組合轉子前三階固有頻率相對降低量在不同時間分布值相近,且頻率相對降低量隨時間呈遞增趨勢。

表2 組合轉子固有頻率隨時間變化關系

6松弛引起失諧組合轉子性能退化評估

為判斷設備退化情況,本文選組合轉子頻率相對降低量作為拉桿松弛引起的組合轉子性能退化量參數,退化量D定義式如下:

D=(f0-ft)/f0

(13)

式中:f0,ft分別為0時刻、t時時刻組合轉子固有頻率。

各退化量參數值即為退化數據,對退化數據進行曲線擬合所得退化參數隨時間變化的函數即為退化軌跡[17],見圖6,獲得組合轉子(工作溫度873 K)退化軌跡服從指數函數分布,即

D=D1e(-t/D2) +D3

(14)

式中:D1,D2,D3為常數。

圖6 組合轉子在不同失諧量下的退化軌跡 Fig.6 Degradation path at different mistuning degree

由圖6看出,預緊失諧對組合轉子退化軌跡影響不大。由式(14)知,隨時間增加退化量趨近退化參數D3值;在1×105h,D3在負失諧、諧調、正失諧時分別為0.00656、0.00662、0.00666。

7預緊失諧組合轉子性能退化特征分析

拉桿初始預緊不均產生的廣義附加力矩相當于給轉子施加附加力矩,將造成轉子初始彎曲[18]。為研究拉桿松弛所致預緊失諧組合轉子退化特征,對其在簡單彈性支承下進行不平衡響應分析。軸承支撐剛度為108N/m。結合組合轉子幾何特性與振型曲線,不平衡響應觀測點位置見圖7。

圖7 激振點與觀測點位置 Fig.7 Position of exciting points and observation points

軸向力8F0作用下的輪盤界面剛度K通過虛擬材料層彈性模量表達,將大小為F的附加預緊力施加于拉桿上,對失諧量在-10%~10%內的組合轉子進行含預應力模態分析。在此基礎上用模態綜合法對失諧組合轉子進行諧響應分析。獲得3000 r/min轉速、預緊失諧下組合轉子不平衡響應。失諧量為10%的不平衡響應見圖8。由圖8看出,預緊失諧組合轉子一階振動幅值最大,其中沿組合轉子軸向振幅由中間向軸端逐步減小,觀測點3處幅值最大。負失諧下組合轉子不平衡響應規律相同,取觀測點3處振幅進行分析。

圖8 失諧的組合轉子不平衡響應 Fig.8 Unbalance response about mistuning rotor

圖9 組合轉子振幅隨時間的變化 Fig.9 The amplitude of combined rotor change with time

同前,在溫度873K、轉速3000 r/min條件下用模態綜合法對組合轉子進行含預應力的諧響應分析,獲得組合轉子最大幅值隨退化時間變化見圖9。由圖9看出,隨時間增加預緊失諧組合轉子響應幅值減小,但遠大于預緊諧調下振幅。此因在失諧組合轉子退化過程中受預緊降低及預緊失諧兩因素耦合作用中預緊失諧起主導作用,失諧量隨時間減小對組合轉子振幅影響較大。初始預緊失諧程度相同時,隨時間增加負失諧對組合轉子振動影響較大,在104~105h內負失諧所致不平衡響應幅值為正失諧的1.38~1.66倍。

8拉桿松弛所致組合轉子性能退化試驗

圖10 組合轉子實驗裝置 Fig.10 The experiment equipment of combined rotor

本文建立有限元模型分析組合轉子動力學特性,計算中存在一定簡化及假設,建立與圖1同尺寸實驗模型轉子,見圖10。在不同預緊情況下實測振動模態,驗證有限元方法的分析結果。

實驗采用錘擊法,即用力錘敲擊組合轉子上某點施加激勵,由LMS模態分析儀器記錄試件測點響應,經數據處理獲得組合轉子固有頻率。通過測量拉桿軸向拉伸應變獲得不同擰緊力矩下拉桿預緊力,采用120Ω電阻式應變片、HX3811靜態應變儀,測試分辨率為1 με,誤差±2 με。

組裝、預緊組合轉子過程中用立式組裝、對稱交叉預緊及分步拉緊三種方式實現拉桿預緊諧調。在預緊諧調基礎上擰緊或擰松某根拉桿獲得預緊失諧組合轉子。本實驗取873 K溫度下拉桿的等效預緊力,通過調整預緊力值表征拉桿松弛程度,測試獲得組合轉子彎曲振動頻率隨時間變化規律,見表3。由表3看出,實驗組合轉子固有頻率與有限元分析結果相差不大,取組合轉子頻率相對降低量作為退化數據,實驗組合轉子的退化軌跡與有限元分析結果相差不大。一定程度上可證實本文轉子模型的正確性及組合轉子性能退化分析結果的可靠性。

9結論

(1)通過所提考慮粗糙界面的預緊失諧組合轉子動力學建模方法發現,失諧拉桿隨時間增加有趨于諧調趨勢;初始失諧量相當情況下隨時間增加初始正失諧拉桿預緊力失諧程度小于負失諧。

(2)建立以組合轉子頻率相對降低量表征性能退化量的定量評估方法,以組合轉子頻率相對降低量作為組合轉子性能退化量參數發現,拉桿松弛引起的組合轉子退化量參數符從指數函數分布;預緊失諧對退化量分布影響不大。

表3 實驗組合轉子在不同預緊力作用下的彎振頻率

(3)由預緊失諧組合轉子性能退化特征發現,時間增加預緊失諧組合轉子響應幅值減小。初始預緊失諧程度相同時,隨時間增加負失諧對組合轉子振動影響較大,在104~105h范圍內負失諧所致不平衡響應幅值為正失諧的1.38~1.66倍。

參考文獻

[1]朱啟寶, 汪鴻振. LM5000 燃機振動研究[J]. 噪聲與振動控制,2006, 26(2): 16-19.

ZHU Qi-bao,WANG Hong-zhen. Vibration analysis of LM5000 gas turbine[J]. Noise and Vibration Control 2006, 26(2): 16-19.

[2]Mohammadi M, Salimi H R. Failure analysis of a gas turbine marriage bolt[J]. Journal of Failure Analysis and Prevention, 2007, 7(2): 81-86.

[3]Castanier M P, Pierre C. Modeling and analysis of mistuned bladed disk vibration: current status and emerging directions [J]. Journal of Propulsion and Power, 2006, 22(2): 384-396.

[4]袁淑霞,張優云,蔣翔俊,等.拉桿失諧模型及其對端面弧齒應力分布的影響[J].哈爾濱工業大學學報,2013,45(5):64-69.

YUAN Shu-xia, ZHANG You-yun, JIANG Xiang-jun,et al. Analysis of bolt preload mistuned model and its impact on stress distribution of curvic couplings[J]. Journal of Harbin Institute of Technology, 2013,45(5):64-69.

[5]Diakunchak I S. Performance degradation in industrial gas turbines[J]. ASME Journal of Engineering for Gas Turbines and Power, 1992, 114:161-168.

[6]Brooks F J. GE gas turbine performance characteristics[R]. GE Power Systems GER-3567H,2003.

[7]Kurz R, Brun K. Degradation in gas turbine systems[J]. Journal of Engineering for Gas Turbines and Power, 2001,123(1):70-77.

[8]Brotherton T, Jahns G, Jacobs J, et al. Prognosis of fault in gas turbine engines[J].IEEE Aerospace Conference Proceedings, 2000(6):163-171.

[9]Li Y G, Nilkitsaranont P. Gas turbine performance prognostic for condition-based maintenance[J]. Applied Energy, 2009(86): 573-584.

[10]Venturini M, Puggina N. Prediction reliability of a statistical methodology for gas turbine prognostics[J].Journal of Engineering for Gas Turbines and Power,2012,134(10):101601.

[11]肖會芳,卲毅敏,徐金梧. 粗糙界面法向接觸振動響應與能量耗散特性研究[J]. 振動與沖擊, 2014,33(4) : 149-155.

XIAO Hui-fang, SHAO Yi-min, XU Jin-wu. Dynamic response and energy dissipation characteristics for normal contact vibration of a rough interface[J]. Journal of Vibration and Shock,2014,33(4) : 149-155.

[12]李輝光,劉恒,虞烈. 粗糙機械結合面的接觸剛度研究[J]. 西安交通大學學報,2011, 45(6): 69-74.

LI Hui-guang, LIU Heng, YU Lie. Study on contact stiffness of rough mechanical joint surface [J]. Journal of Xi’an Jiao Tong University,2011,45(6): 69-74.

[13]廖伯瑜,周新民,尹志宏. 現代機械動力學及其工程應用[M]. 北京:機械工業出版社,2004:243-245.

[14]郭進全,軒福貞,王正東,等. 基于蠕變的高溫構件應力松弛損傷模型[J].核動力工程,2009,30(4):9-12.

GUO Jin-quan, XUAN Fu-zhen, WANG Zheng-dong,et al. Creep based stress relaxation damage model for high temperature components[J]. NuclearPower Engineering, 2009, 30(4):9-12.

[15]李玲. 渦軸發動機轉子連接螺栓的強度與蠕變壽命研究[D]. 南京:南京航空航天大學, 2008.

[16]于慧臣,吳學仁. 航空發動機設計用材料數據手冊[M]. 北京: 航空工業出版社,2010:85-107.

[17]Chao M T. Degradation analysis and related topics: some thoughts and a review[J]. Proc, Natl, Sci. Counc. ROC(A), 1999,23(5):555-566.

[18]Lu Y J, Zhang Y F, Yu Y B, et al. Nonlinear dynamics of flexible rotor system supported on fixed-tilting pad combination journal bearing[J]. Journal of Central South University of Technology, 2011, 18(3): 22-29.

主站蜘蛛池模板: 九色视频最新网址| 欧美不卡二区| 色首页AV在线| 亚洲Aⅴ无码专区在线观看q| 玩两个丰满老熟女久久网| 国产AV毛片| 欧美成人a∨视频免费观看| 香蕉视频国产精品人| 欧美成人精品高清在线下载| 亚洲狼网站狼狼鲁亚洲下载| 四虎影视8848永久精品| 国产精品视频第一专区| 国产成人精品男人的天堂下载| 2020精品极品国产色在线观看 | 久久国产V一级毛多内射| 婷婷午夜天| 三级视频中文字幕| 在线免费观看AV| 一级香蕉视频在线观看| 在线a网站| 久久频这里精品99香蕉久网址| 国产精品福利在线观看无码卡| 欧美不卡视频一区发布| 色哟哟国产精品一区二区| 日韩欧美高清视频| 亚洲成网站| 亚洲AV人人澡人人双人| 亚洲综合片| 中文字幕av无码不卡免费| 国产精品永久免费嫩草研究院| 亚洲欧洲日韩久久狠狠爱| 2021国产精品自产拍在线| a级毛片免费网站| 国产在线自揄拍揄视频网站| 日韩最新中文字幕| 欧美精品亚洲二区| 黄片在线永久| 57pao国产成视频免费播放| 久久精品视频一| 国产精品视频系列专区| 在线观看无码a∨| 久久综合九色综合97网| 性69交片免费看| 日韩精品专区免费无码aⅴ| 成人国产精品网站在线看| 2021国产精品自拍| 精品综合久久久久久97| 永久免费AⅤ无码网站在线观看| 亚洲成av人无码综合在线观看| 1024你懂的国产精品| 国产精品露脸视频| 亚洲天堂视频在线免费观看| 人妻精品久久无码区| 欧美福利在线| 欧美在线黄| 欧美三级视频在线播放| 欧美激情福利| 无码丝袜人妻| 欧美高清国产| 国产网站免费观看| 国产毛片高清一级国语| 亚洲最大在线观看| 亚洲欧洲日韩国产综合在线二区| 亚洲黄色视频在线观看一区| 夜夜拍夜夜爽| 久久成人免费| 欧美成一级| 啊嗯不日本网站| 国产爽爽视频| www.亚洲一区| 日a本亚洲中文在线观看| 2021国产精品自产拍在线| 伊人久久精品无码麻豆精品| 国产91麻豆视频| 性网站在线观看| 久久精品人人做人人爽| 免费一级成人毛片| 国产真实乱子伦精品视手机观看| 性做久久久久久久免费看| 日韩AV无码免费一二三区| 视频一本大道香蕉久在线播放| 国产精品高清国产三级囯产AV|