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

多種邊界條件下自旋圓柱殼振動特性研究

2024-10-24 00:00:00董有恒李映輝李翔宇茅曉晨
振動工程學報 2024年10期

摘要: 自旋圓柱殼作為實際工程結構中的重要部件,其殼體端部邊界條件復雜多樣,對部件的振動特性有著重要影響。為研究多種邊界條件下自旋圓柱殼的振動特性,基于Lagrange方程和Novozhilov殼體理論建立自旋圓柱殼的結構動力學模型;根據圓柱殼邊界條件的數學描述,采用Chebyshev多項式構建滿足邊界條件且與殼體結構參數無關的位移場離散函數;通過求解與振動特性相關的特征方程得到靜止圓柱殼的固有頻率和模態,分析殼體轉動慣性對固有頻率的影響,討論殼體理論在不同殼體幾何尺寸下的適用性;此外,發現了環向波數相關的模態函數,并將其用于圓柱殼零階環向模態固有頻率的求解及多種邊界條件下自旋圓柱殼行波固有頻率的求解,討論了自旋圓柱殼行波固有頻率隨結構參數變化的情況。

關鍵詞: 自旋圓柱殼; 多種邊界條件; 轉動慣性; 位移場離散

中圖分類號: O321; TH113.1 文獻標志碼: A 文章編號: 1004-4523(2024)10-1731-08

DOI:10.16385/j.cnki.issn.1004-4523.2024.10.011

引 言

憑借其中心對稱特性和較高的比剛度,圓柱殼結構作為重要零部件在機械領域和航空航天領域得到了廣泛應用。高速離心機的內膽和外殼、風電機組聯軸器的合金鋼中間管、航空發動機的葉片輪轂等零部件均可簡化為自旋圓柱殼力學模型。不同于靜止圓柱殼的駐波特性,由自旋運動引起的科氏效應會導致自旋圓柱殼結構出現行波振動[1?2]。同時,實際工程中的圓柱殼結構邊界條件種類較多,常見的簡支邊界條件難以滿足振動特性求解精度的要求。因此,針對自旋圓柱殼結構,發展除簡支邊界條件外的其他多種邊界條件描述方法,并準確分析自旋圓柱殼的行波振動特性,能夠為實際工程結構中多種邊界自旋圓柱殼的動力學分析提供研究基礎。

由于單個三角函數能準確描述簡支邊界條件下圓柱殼關于軸向坐標的模態構型,現有研究中關于圓柱殼振動特性的研究大多在兩端簡支邊界條件下開展[3?9]。徐港輝等[3]通過正、余弦模態函數研究了兩端簡支圓柱殼動力學響應中的模態參與問題。Bich等[4]研究了簡支邊界功能梯度材料圓柱殼在外加載荷作用下的非線性振動,首先將殼體位移分量在軸向和環向離散為單模態下的三角函數形式,然后在非線性運動控制方程的基礎上忽略殼體面內慣性項,得到關于殼體徑向位移的非線性微分方程,最后得到了圓柱殼的頻響函數和非線性動力學響應的數值解。Wang[5]運用Donnell非線性殼體理論建立了復合材料自旋圓柱殼關于徑向位移的非線性運動控制方程,將殼體徑向位移用三角函數表示,運用諧波平衡法研究了簡支邊界條件下圓柱殼的強迫振動,研究發現較高的自旋速度會導致復雜的頻響曲線并產生新的穩態。基于一階剪切變形殼體理論,Sofiyev等[6]研究了彈性地基上功能梯度材料圓柱殼的非線性振動,考慮殼體端部簡支邊界,假設殼體單模態振動,將殼體位移場離散為三角函數形式,通過數值算例討論了系統參數對圓柱殼頻響曲線的影響規律。

為研究多種邊界條件下圓柱殼的振動特性,學者們用三角函數和雙曲函數組合形式的梁模態函數來描述圓柱殼關于軸向坐標的振型[10?12]。孫述鵬等[10]采用固支邊界梁的模態函數描述圓柱殼模態函數在軸向的構型,分別研究了橫向簡諧力和恒力作用下兩端固支自旋圓柱殼的行波振動響應。此外,Jin等[13]和Lin等[14]通過在圓柱殼端部添加虛擬彈簧來模擬彈性邊界,其中不同的殼體邊界條件由不同形式的彈簧剛度系數組合決定,此特定的剛度系數組合隨殼體幾何尺寸和材料性質而發生變化。石先杰等[15]采用邊界虛擬彈簧模擬功能梯度圓柱殼的任意經典邊界,研究了殼體自由振動和瞬態振動特性。龐福振等[16]研究了邊界虛擬彈簧約束下半球殼的受迫振動響應。

為建立多種邊界條件下與圓柱殼結構參數無關的殼體模態函數,Kurylov等[17]將圓柱殼位移場離散為正交多項式或三角函數的多項求和形式,在此基礎上對比分析了簡支和固支邊界圓柱殼的固有頻率。Dong等[18]使用Chebyshev多項式分別給出了滿足圓柱殼8種邊界條件的位移場離散函數,得到了較高精度的殼體固有頻率。

綜上所述,學者們對多種邊界條件下圓柱殼的振動特性分析已開展廣泛研究,本文將進一步研究多種邊界條件下自旋圓柱殼的行波振動特性。給出含有自旋運動的殼體結構模態分析流程,同時分析殼體轉動慣性對固有頻率求解精度的影響,并通過數值算例對比不同殼體理論的適用范圍,研究不同邊界條件對行波振動特性的影響。

1 自旋圓柱殼動力學模型

自旋圓柱殼力學模型及其橫截面微元如圖1所示,其中紅色虛線圓所在的環面為圓柱殼中面,圓柱殼以自旋速度(r/min)繞其中心軸線作自旋運動。圓柱殼長度、壁厚、中面半徑分別為L,h,r;圓柱殼材料的楊氏模量、泊松比、密度分別為E,,。為建立自旋圓柱殼動力學模型,并使用Chebyshev多項式表示位移場離散函數,需將圓柱殼軸向坐標轉移到多項式的定義區間[-1,1]。將慣性坐標系(,,)的軸置于圓柱殼中心軸線,柱坐標系(,,z)固結于殼體中面,其中,[-1,1],[0,2],[-0.5h,0.5h]。在固結于中面的柱坐標系下,圓柱殼中面任意點沿,,z三個方向的位移分別為,,。

考慮圓柱殼橫截面轉動,則圓柱殼任意點沿,,z三個方向的位移,,可表示為關于圓柱殼中面位移場的函數:

(1)

式中 變量下標中的“,”表示對下標“,”后面變量的微分。

本文使用薄壁殼體理論中較為精確的Novozhilov殼體理論[19],其應變場表示為:

(2)

其中,殼體中面應變和曲率變化分量分別表示為:

(3)

圓柱殼應力?應變關系則表示為:

(4)

式中 (i,j=1,2,4)為與楊氏模量和泊松比有關的彈性常數[2]。

為方便對圓柱殼邊界條件進行描述,接著給出圓柱殼橫截面處軸力和彎矩的定義式:

(5)

根據圖1所示的殼體橫截面微元,圓柱殼應變能可表示為:

(6)

同時,自旋圓柱殼的動能K包含結構的振動動能和自旋運動引起的動能,具體表示為:

(7)

式中 變量上方加‘·’表示對時間變量t的微分。

由離心慣性力產生的環向拉伸應變能為:

(8)

結合式(1)和(5),可給出圓柱殼在端部可移動簡支(Sm)、端部不可移動簡支(Si)、滑移(Sd)、固支(C)、自由(F)等多種邊界條件下,關于殼體邊界處(=-1或1)的位移場、軸力和彎矩的數學描述:

(9)

2 振動特性

為了研究多種邊界條件下自旋圓柱殼的行波振動特性,需要對位移場進行離散,再通過求解特征值問題來得到殼體結構固有頻率。

由于圓柱殼位移場關于環向坐標的模態構型呈花瓣狀,且相互正交,故可通過單一的正、余弦函數準確描述。而圓柱殼位移場關于軸向坐標的離散函數與殼體邊界條件有關,文獻[18]給出的數值結果證明,若使用梁模態函數(三角函數和雙曲函數的組合形式)近似描述多種邊界條件下殼位移場關于軸向坐標的離散函數,得到的固有頻率存在較大誤差。

基于此,本文選取Chebyshev多項式描述圓柱殼位移場關于軸向坐標的離散函數。當存在自旋運動時,多項式描述的陀螺矩陣具有奇異性和稀疏性,無法直接求解相對應的特征值問題。因此,首先通過對靜止圓柱殼的模態分析,得到多種邊界條件下各階模態的離散函數描述形式,在此基礎上求解多種邊界條件下自旋圓柱殼的行波振動特性,具體流程圖如圖2所示。

首先給出多種邊界條件下圓柱殼單模態振動的殼體中面位移場離散形式:

(10)

式中 U(),V()和W()表示位移場關于軸向和環向坐標的離散函數;f(t)為與時間有關的簡諧函數。

考慮第一類Chebyshev多項式,其前j階多項式依次表示為:,,。式(10)中的離散函數可表示為的線性組合形式:

(11)

式中 n表示圓柱殼各階模態的環向波數;表示Chebyshev多項式的項數;,,表示Chebyshev多項式的系數。

將式(10)和(11)代入式(9),通過運算即可得到滿足圓柱殼不同邊界條件的Chebyshev多項式組合形式。此時,滿足邊界條件的圓柱殼離散函數為:

(12)

式中 fun表示滿足邊界條件的第j個關于Chebyshev多項式的線性組合函數;M表示U(),V()和W()關于軸向坐標的離散函數的求和項數;,,表示Chebyshev多項式的修正系數。

本文使用如下Lagrange方程求解圓柱殼固有振動特性:

(13)

式中 q為向量={,…,,,…,,,…,}中的元素。

將式(10)和(12)代入能量表達式(6)~(8),再結合式(13),并假設圓柱殼的簡諧振動f(t)=,其中表示圓柱殼的固有頻率。當自旋速度=0時,可得到靜止圓柱殼的關于固有頻率和固有振型的特征方程:

(14)

式中 ,分別表示靜止圓柱殼的質量矩陣和剛度矩陣。

當自旋速度0時,自旋圓柱殼出現行波振動,固有頻率分為前行波頻率和后行波頻率。對于環向波數為n的行波模態,此時圓柱殼關于軸向坐標和環向坐標的離散函數為:

(15)

式中 ,和為滿足殼體邊界條件的模態函數,與對應的靜止圓柱殼模態函數一致,由式(12)和(14)得到,如圖2流程圖所示。此時,q為向量q={,,,,,}中的廣義變量,根據Lagrange方程得到:

(16)

式中 M,G和K分別表示自旋圓柱殼的質量矩陣、陀螺矩陣和剛度矩陣。值得注意的是,由自旋運動引起的科式效應和離心效應,分別在矩陣G和K中體現。

假設簡諧振動形式q=,其中和分別為待定標量和向量,將此式代入式(16),即可得到二次特征值問題。為方便求解行波固有頻率,使用狀態空間描述,令,廣義坐標q在狀態空間中的運動變為,得到關于自旋圓柱殼行波固有頻率的特征值方程:

(17)

式(17)與(16)具有相同的特征值。求解式(17)可得到6個頻率值,其中最小的兩個頻率對應該邊界條件下自旋圓柱殼的前行波頻率和后行波頻率。

2.1 數值驗證與殼體轉動慣性的影響

為驗證本文中基于Chebyshev多項式的模態離散函數的有效性,及多種邊界條件下圓柱殼振動特性的理論推導與數值模擬的正確性,首先將本文中基于Novozhilov殼體理論的圓柱殼固有頻率與有限元結果、文獻[2]中基于Amabili?Reddy三階剪切變形殼體的理論結果進行對比,如表1所示。通過數值解對比發現,當式(12)中多項式項數M=20時,固有頻率值收斂,故在本文數值解求解過程中M均取為20。

值得注意的是,在C?C和Si?Si邊界條件下,表1中針對圓柱殼零階環向模態(n=0)的固有頻率計算,在C?C和Si?Si邊界條件下,使用式(12)和(14)無法直接得到其固有頻率。此時,采用(1,1)模態下求得的離散函數作為(1,0)模態頻率計算的離散函數,可得到較高精度的固有頻率。由表1可知,本文解與有限元解、文獻[18]的解有較好的對照,驗證了Chebyshev多項式在模擬圓柱殼多種邊界條件時的有效性。

現有文獻在關于圓柱殼振動特性和非線性振動響應的分析中,為了便于計算分析,一般忽略殼體轉動慣性[4,6,10?11,17],在式(7)中僅考慮中面位移場相關的動能。基于Novozhilov殼體理論和Donnell殼體理論,表2分別給出了考慮殼體轉動慣性和忽略殼體轉動慣性情況下圓柱殼的固有頻率,同時對比了圓柱殼在不同半徑?厚度比(r/h)情況下,上述兩種殼體理論的適用性。為了更直觀地進行分析,將表2中基于上述兩種殼體理論的固有頻率計算結果與有限元結果對比的相對誤差展示在圖3中。

從表2和圖3中可以看出,針對兩種不同的殼體理論,忽略殼體轉動慣性造成的固有頻率計算誤差較小,這是因為殼體厚度與殼體長度相比較小,且轉動慣性是比殼體中面慣性更高階的量,故在較大半徑?厚度比情況下,忽略殼體轉動慣性的假設是合理的。此外,對比Novozhilov殼體理論和Donnell殼體理論的適用性,Novozhilov殼體理論比Donnell殼體理論具有更高的固有頻率計算精度。在薄壁殼情況下(r/h200),與有限元結果對比發現,基于Donnell殼體理論的計算結果最大誤差不超過1.6%。在中厚殼情況下(),基于Novozhilov殼體理論的計算結果仍具有較高的精度,與Amabili?Reddy三階剪切變形殼體理論的結果相差較小,而基于Donnell殼體理論的計算結果存在較大誤差。故Novozhilov殼體理論在中厚殼描述中也是一個較好的選擇。

2.2 環向波數相關的模態函數

依據圖2所示的流程求解特征方程式(14),可得到圓柱殼在不同邊界條件下的固有頻率和對應的模態構型。圖4給出了C?C邊界條件下,圓柱殼的(1,2)、(1,3)、(1,4)、(1,5)、(1,6)、(1,7)模態構型沿軸向的變化情況。

在相同軸向半波數條件下,分別對比圖4中,和在不同環向波數情況下隨圓柱殼軸向坐標的變化,發現圓柱殼模態構型在不同環向波數中存在明顯差異,稱此現象為環向波數相關的模態函數。

為進一步驗證圓柱殼環向波數相關的模態函數,表3分別給出了采用(1,2)、(1,3)、(1,4)、(1,5)模態下已知的關于軸向坐標的離散函數作為模態構型時,求解得到的固有頻率。由表3的計算結果及其與準確值間的相對誤差可知,通過自身模態去求解固有頻率可以得到準確值,當用于求解固有頻率的模態函數所對應的環向波數與待求模態波數相差較大時,會造成固有頻率求解誤差較大。

值得注意的是,與圓柱殼的環向波數相關的模態函數在除簡支?簡支邊界外的其他邊界條件中均存在,此現象將在多種邊界條件下,圓柱殼零階環向模態的固有頻率求解、自旋圓柱殼行波振動特性求解、圓柱殼的多模態耦合振動分析中有重要應用。

圖5對比了C?C,Sm?Sm邊界條件下,圓柱殼固有頻率隨環向波數的變化。從圖5中可以看出,不同于梁、板結構,圓柱殼的最低固有頻率一般不在最小模態波數處取得。相同模態波數情況下,固支邊界的固有頻率遠高于端部可移動簡支邊界的固有頻率。

2.3 自旋圓柱殼行波固有頻率

依據圖2所示的流程,考慮Sm?Sm,Sd?Sd邊界條件,表4給出了自旋圓柱殼行波振動的固有頻率,并與現有文獻[2]的結果進行了對比。從表4中可以看出,本文采用Chebyshev多項式模擬圓柱殼多種邊界條件下的位移場離散函數,并使用環向波數相關的模態函數求解多種邊界條件下的自旋圓柱殼行波振動特性是可行的,具有較高的精度。圖6則給出了Sm?Sm,C?C邊界條件下,自旋圓柱殼(1,5)行波模態的固有頻率隨自旋速度的變化。

3 結 論

本文基于Novozhilov殼體理論建立了圓柱殼動力學模型,運用Chebyshev多項式組合形式的位移場離散函數研究了多種邊界條件下自旋圓柱殼的行波振動特性。主要結論如下:

(1)Novozhilov殼體理論可用于中厚殼的理論分析,在圓柱殼振動特性分析中能達到與高階剪切變形殼體理論相近的準確度,且具有較簡便的形式(包括非線性項),在圓柱殼多模態耦合振動分析中,Novozhilov殼體理論是一個較好的選擇。

(2)在較大半徑?厚度比情況下,殼體轉動慣性可忽略不計,以簡化振動特性求解過程。圓柱殼的環向波數相關的模態函數在除簡支?簡支邊界外的其他邊界條件中均存在。

(3)在多種邊界條件下,由Chebyshev多項式構成的模態離散函數可以精確滿足圓柱殼的多種邊界條件。環向波數相關的模態函數在圓柱殼零階環向模態的固有頻率求解、自旋圓柱殼行波振動特性求解、圓柱殼的多模態耦合振動分析具有重要作用。

在后續的研究中,進一步將環向波數相關的模態函數應用于自旋圓柱殼的非線性多模態耦合振動,研究多種邊界條件下自旋圓柱殼的非線性動力學響應。

參考文獻:

[1]曹航,朱梓根. 轉動殼體行波振動的有限元分析方法[J]. 航空動力學報,2002,17(2): 222-225.

GAO Hang,ZHU Zigen. Travelling-wave vibration of rotating shells by finite element method[J]. Journal of Aerospace Power,2002,17(2): 222-225.

[2]Dong Y H,Liu H,Hu H Y,et al. Semi-analytical and experimental studies on travelling wave vibrations of a moderately thick cylindrical shell subject to a spinning motion[J]. Journal of Sound and Vibration,2022,535: 117095.

[3]徐港輝,祝長生. 圓柱殼體動力響應中的模態參與問題研究[J]. 振動工程學報,2024,37(1): 83-94.

XU Ganghui,ZHU Changsheng. Study on modal participation in dynamic responses of cylindrical shells[J]. Journal of Vibration Engineering,2024,37(1): 83-94.

[4]Bich D H,Nguyen N X. Nonlinear vibration of functionally graded circular cylindrical shells based on improved Donnell equations[J]. Journal of Sound and Vibration,2012,331(25): 5488-5501.

[5]Wang Y Q. Nonlinear vibration of a rotating laminated composite circular cylindrical shell: traveling wave vibration[J]. Nonlinear Dynamics,2014,77(4): 1693-1707.

[6]Sofiyev A H,Hui D,Haciyev V C,et al. The nonlinear vibration of orthotropic functionally graded cylindrical shells surrounded by an elastic foundation within first order shear deformation theory[J]. Composites Part B: Engineering,2017,116: 170-185.

[7]Amabili M,Reddy J N. The nonlinear,third-order thickness and shear deformation theory for statics and dynamics of laminated composite shells[J]. Composite Structures,2020,244: 112265.

[8]Yadav A,Amabili M,Panda S K,et al. Forced nonlinear vibrations of circular cylindrical sandwich shells with cellular core using higher-order shear and thickness deformation theory[J]. Journal of Sound and Vibration,2021,510: 116283.

[9]Teng M W,Wang Y Q. Spin-induced internal resonance in circular cylindrical shells[J]. International Journal of Non-Linear Mechanics,2022,147: 193-206.

[10]孫述鵬,曹登慶,初世明. 轉動薄壁圓柱殼行波振動響應分析[J]. 振動工程學報,2013,26(3): 459-466.

SUN Shupeng,CAO Dengqing,CHU Shiming. Analysis of travelling wave vibration response for thin rotating cylindrical shell[J]. Journal of Vibration Engineering,2013,26(3): 459-466.

[11]Dong Y H,Zhu B,Wang Y,et al. Nonlinear free vibration of graded graphene reinforced cylindrical shells: effects of spinning motion and axial load[J]. Journal of Sound and Vibration,2018,437: 79-96.

[12]Song X Y,Ren Y P,Han Q K. Nonlinear vibration of rotating cylindrical shell due to unilateral contact induced tip rubbing impact: theoretical and experimental verification[J]. Mechanical Systems and Signal Processing,2022,164: 108244.

[13]Jin G Y,Ye T G,Ma X L,et al. A unified approach for the vibration analysis of moderately thick composite laminated cylindrical shells with arbitrary boundary conditions[J]. International Journal of Mechanical Sciences,2013,75: 357-376.

[14]Lin H G,Cao D Q,Shao C H. An admissible function for vibration and flutter studies of FG cylindrical shells with arbitrary edge conditions using characteristic orthogonal polynomials[J]. Composite Structures,2018,185: 748-763.

[15]石先杰,左朋. 熱環境下功能梯度圓柱殼振動特性分析[J]. 振動工程學報,2023,36(2): 526-533.

SHI Xianjie,ZUO Peng. Vibration analysis of functionally graded cylindrical shell under thermal environment[J]. Journal of Vibration Engineering,2023,36(2): 526-533.

[16]龐福振,張明,高聰,等. 復雜邊界條件下半球殼受迫振動響應分析[J]. 振動工程學報,2024,37(3): 374-383.

PANG Fuzhen,ZHANG Ming,GAO Cong,et al. Forced vibration response analysis of hemispherical shell under complex boundary conditions[J]. Journal of Vibration Engineering,2024,37(3): 374-383.

[17]Kurylov Y,Amabili M. Polynomial versus trigonometric expansions for nonlinear vibrations of circular cylindrical shells with different boundary conditions[J]. Journal of Sound and Vibration,2010,329: 1435-1449.

[18]Dong Y H,Hu H Y,Wang L F. A comprehensive study on the coupled multi-mode vibrations of cylindrical shells[J]. Mechanical Systems and Signal Processing,2022,169: 108730.

[19]Amabili M. Nonlinear vibrations of circular cylindrical shells with different boundary conditions[J]. AIAA Journal,2003,41(6): 1119-1130.

Vibration characteristics of spinning cylindrical shells under various boundary conditions

DONG You-heng1,LI Ying-hui2,LI Xiang-yu2,MAO Xiao-chen1

(1.College of Mechanics and Engineering Science,Hohai University,Nanjing 211100,China; 2.School of Mechanics and Aerospace Engineering,Southwest Jiaotong University,Chengdu 611756,China)

Abstract: Spinning cylindrical shells are critical components in practical engineering structures. The boundary conditions at the shell ends are diverse and significantly influence the vibration characteristics of the shell. To study these characteristics under various boundary conditions,a dynamic model of the spinning cylindrical shell is established using Lagrange equations and Novozhilov’s shell theory. The mathematical description of the boundary conditions for the cylindrical shell is combined with the discretized displacement functions,which are constructed based on a linear combination of Chebyshev polynomials. These functions satisfy the boundary conditions and are independent of the cylindrical shell's parameters. The vibration characteristics of stationary cylindrical shells are determined by solving the eigenvalue problems,revealing the influence of rotary inertia on the vibration characteristics. The applicability of different shell theories with respect to various geometrical parameters of the shell is discussed. Additionally,circumferential wave-dependent mode functions are identified and used to compute the natural frequencies of shell modes with the zero circumferential waves,as well as the travelling waves of the spinning cylindrical shell under different boundary conditions. The impact of structural parameters on the natural frequencies of the travelling waves is also analyzed.

Key words: spinning cylindrical shell;various boundary conditions;rotary inertia;displacement field discretization

作者簡介: 董有恒(1992―),男,博士,副研究員。E-mail:dyhe@hhu.edu.cn。

通訊作者: 李映輝(1964―),男,博士,教授。E-mail:yhli2007@sina.com。

主站蜘蛛池模板: 成年人国产视频| 伊人久久影视| 国产精品女熟高潮视频| 亚洲国产成人精品无码区性色| 成人国内精品久久久久影院| 亚洲综合欧美在线一区在线播放| 中国国语毛片免费观看视频| 亚洲精品天堂在线观看| 免费一级成人毛片| 亚洲一区二区三区国产精华液| 丰满人妻一区二区三区视频| 国产制服丝袜91在线| 极品性荡少妇一区二区色欲 | 亚洲视频一区在线| 91 九色视频丝袜| 久久国产高潮流白浆免费观看| 久久无码av三级| 国产香蕉在线视频| 国产精品亚洲综合久久小说| 日韩不卡高清视频| 欧美国产成人在线| 91精品国产情侣高潮露脸| AV天堂资源福利在线观看| 国产成人久久综合777777麻豆| 任我操在线视频| aa级毛片毛片免费观看久| 国内精品九九久久久精品| 亚洲成人精品在线| 好吊日免费视频| 国产黄在线观看| 亚洲乱码在线视频| 日韩无码视频播放| 在线欧美国产| 日本日韩欧美| 久久精品91麻豆| 国产网友愉拍精品| 国产精品香蕉| 日韩午夜福利在线观看| 中文字幕第4页| 黄网站欧美内射| 色悠久久综合| 试看120秒男女啪啪免费| 福利在线免费视频| 亚洲无码视频图片| 亚洲精品视频免费| 国产h视频在线观看视频| 亚洲精品无码AⅤ片青青在线观看| 亚洲精品成人片在线观看| 午夜日本永久乱码免费播放片| 日韩第九页| 一本综合久久| 日韩av手机在线| 国产在线一区视频| 国产自产视频一区二区三区| 国产精品永久免费嫩草研究院| 在线国产毛片| 夜色爽爽影院18禁妓女影院| 青青青草国产| jizz在线免费播放| 美女黄网十八禁免费看| 国产三级毛片| 毛片基地美国正在播放亚洲| 国产丝袜无码精品| 国产女人喷水视频| 亚洲综合色吧| 久久99精品国产麻豆宅宅| 成人伊人色一区二区三区| 99久久精品免费看国产免费软件 | 国产尤物jk自慰制服喷水| 999精品色在线观看| 国产精品真实对白精彩久久| 一级黄色片网| 成人午夜网址| 欧美午夜视频在线| 精品久久久久久久久久久| 亚洲成人网在线观看| 欧美性爱精品一区二区三区| 久久青草视频| 九色最新网址| 亚洲综合欧美在线一区在线播放| 亚洲AⅤ永久无码精品毛片| 色精品视频|