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

雙閉室復合材料薄壁梁的結構阻尼細觀分析

2012-09-15 10:17:32任勇生宋玉壁孫雙雙
振動與沖擊 2012年23期
關鍵詞:模態復合材料結構

任勇生,宋玉壁,孫雙雙

(1.山東科技大學 機械電子工程學院,青島 266510;2.青島科技大學,青島 266061)

在過去的幾年間,風力發電機的尺寸持續增加,直徑超過120 m的風輪樣機正在樣機試驗階段。隨著風輪葉片尺寸的增加以及碳纖維復合材料的使用,在周圍風場的作用下,不可避免地存在氣動力、慣性力和彈性力間的耦合,有可能誘發葉片模態之間的耦合振動,導致葉片發生顫振。近來,關于風力機葉片顫振抑制,已經成為現代大型風力機結構動力學研究的一個重要的方面,日益受到風能工程領域的關注與重視[1-3]。

葉片顫振抑制可以采用主動與被動的方法。目前由于主動顫振抑制存在附加控制系統重量與能耗問題,特別是主動控制系統花費過大,勢必增加風力發電的成本,因此,被動控制方法開始受到關注[2-3]。研究表明,被動阻尼在改善柔性結構的動力學特性,包括實施有效振動控制、提高疲勞壽命和改善氣動彈性穩定性等方面,已經顯示出較強的實用性。氣彈阻尼作為氣動阻尼和結構阻尼的總和,是衡量與控制葉片是否發生顫振的一個基本指標。樹脂基纖維復合材料結構由于高分子材料內部晶體的內摩擦以及纖維/基體間的中間相阻尼,展示出高結構阻尼特性。此外,復合材料阻尼還能夠通過提高組分阻尼、纖維混雜以及設置阻尼層進行調整,以獲得理想的被動阻尼效果。復合材料結構阻尼的預測,離不開先進復合材料阻尼理論的指導。近年來,宏/細觀阻尼分析理論在復合材料殼、板和實心梁的阻尼研究和設計以及復合材料結構的被動振動控制中,得到了廣泛的應用[4],然而對風力機葉片這種復合材料結構阻尼特性的研究,至今尚不多見。風力機葉片具有閉合截面復合材料薄壁梁的結構特點,它不僅在風力機葉片上采用,在先進飛機固定機翼和直升機旋翼上也被使用。因此,建立閉合截面復合材料薄壁梁的結構阻尼分析模型,對于評價此類典型復合材料結構的阻尼耗散能力并且據此進行精確結構阻尼設計,以便于建立相應的被動顫振抑制有效方法,具有理論研究價值和工程實用性。

Suresh等[5]在不考慮復合材料薄壁箱形梁的結構特點以及截面翹曲變形情況下,采用經典的層合薄板有限單元法進行結構離散,復合材料的阻尼機理采用彈性-黏彈性對應原理(復模量法)進行描述,分析了邊界條件和鋪層角對固有頻率和結構阻尼的影響。Saravanos等[6]從空心截面復合材料薄壁梁的幾何和變形特征出發,基于Rehfied薄殼位移場方程[7],提出了一個管狀層合復合材料梁的三維有限元阻尼分析模型,其中,復合材料阻尼通過耗散能進行度量,不計薄壁梁拉-剪耦合、拉-扭耦合和彎-扭彈性耦合的作用,得到三種規則截面梁阻尼分析結果。Chortis等[8]在Saravanos等[6]工作的基礎上進一步研究了復合材料彈性耦合的影響,采用有限元方法計算了單/雙閉室兩種翼型截面的復合材料薄壁梁的模態阻尼,并對照實驗結果進行了驗證。Ramkumar等[9]分別采用三維梁和層合薄板兩種有限元模型,研究薄壁復合材料箱形梁的阻尼特性。結構建模同樣也是基于Rehfied薄壁梁理論,而復合材料阻尼的描述則采用復模量法并且不考慮薄壁梁扭轉翹曲的影響。由于上述阻尼模型所依據的Rehfield薄壁梁理論非漸進準確,難以保證橫截面剛度特性計算結果一致精確性[10]。

目前關于閉合截面復合材料薄壁梁的阻尼預測模型,僅限于有限元分析模型。雖然有限元離散在復雜結構的力學建模已經得到廣泛的應用,但是由于計算量大、內存需求高,對于此類結構的阻尼預測,缺乏實用性。任勇生等[11]基于復合材料薄壁梁理論并且結合最大應變能理論,建立了彎扭耦合復合材料單閉室薄壁梁的結構宏觀阻尼預測的分布參數模型,采用Galerkin法求解振動分析模型,分析了纖維鋪層角、薄壁梁長寬比和截面寬高比對阻尼性能的影響。計算過程表明,基于上述分布參數模型,具有計算速度快和內存需求小的優點,并且所得到的計算結果與有限元相比可以滿足對計算精度的要求。

單閉室與雙閉室截面梁的主要差別表現在對扭轉的分析,由于雙閉室截面梁屬于靜不定問題,存在無法由平衡方程確定的獨立的剪力流,特別是對于各向異性復合材料梁而言,由于扭轉與其它變形之間的耦合,上述差別變得更為明顯,因此,雙閉室復合材料薄壁梁的橫截面剪切與扭轉特性分析相對具有較大的難度。

本文在前期研究的基礎上,提出一個彎扭耦合復合材料雙閉室薄壁梁的結構阻尼連續分布模型。基于單層混雜材料的細觀力學阻尼計算方法[12]和多胞模型[13],分別獲得單層復合材料的等效阻尼特性和等效彈性特性,薄壁梁采用VAM[10]進行結構力學建模,基于Hamilton原理導出復合材料薄壁梁的自由振動模型,將位移按廣義坐標進行模態展開,采用Galerkin法求解振動分析模型。在導出復合材料薄壁梁的應變能和耗散能表達式的基礎上,根據最大應變能理論,對雙閉室薄壁復合材料梁模態阻尼性能進行理論分析,雙閉室截面剛度系數采用文獻[14]給出的公式進行計算,通過數值近似計算,揭示了纖維含量、纖維鋪層角對阻尼性能的影響。

1 分析模型

1.1 位移場和應變場

本節采用VAM[10]導出各向異性雙閉室薄壁截面梁的位移場和二維截面分析模型。圖1表示一細長的雙閉室復合材料薄壁梁,長度為L,厚度h,中面曲率半徑 r,d表示最大橫截面尺寸。假定 d?L,h?d,h?r,坐標ξ沿中面法線方向度量,滿足 -h/2≤ξ≤h/2。(x,y,z)表示薄壁梁整體坐標系;(x,s,ξ)表示薄壁梁局部坐標系,s沿截面中線(即,薄壁梁的中面與橫截面的交線)切向,ξ沿截面中線外法線方向。

薄壁梁上的任意點沿坐標系(x,y,z)的三個坐標軸方向的位移分量為[10]:

其中:截面翹曲函數g(s,x)由位移場的連續性、環向力為零以及定常剪力的條件決定。

與位移場(1)相對應的二階近似應變場為:

圖1 雙閉室坐標系統及運動學變量圖Fig.1 Two-cell coordinate systems and kinematic variables

1.2 復合材料薄壁梁的應變能密度

復合材料薄壁梁的應變能密度:

對式(3)沿厚度積分,可得:

其中:

hk、hk-1分別為第k層的上、下表面坐標,N 為總層數。

假設薄壁梁的環向應力很小,可以忽略不計,則有:

由式(4)、(6),得:

將式(7)代入(4)消去 γ22,得:

其中:

復合材料薄壁梁的應變能為:

利用應變-位移關系(2),得位移/扭轉角表示應變能為:

1.3 復合材料薄壁梁的耗散能密度

復合材料薄壁梁的耗散(阻尼)應變能密度:

其中:阻尼(損耗因子)矩陣[ψ]為下列對角矩陣:

其中:ψ1,ψ2,ψ6分別表示纖維復合材料單層的軸向、橫向拉壓和面內剪切損耗因子。

假設:① 基體/纖維界面相互不脫離,為理想粘結;② 復合材料各組分的應變沿纖維方向是相等的,材料應力在與纖維垂直方向是均勻的;③ 纖維與基體具有粘彈性性質,且纖維單向埋設。于是按照復合材料細觀力學阻尼分析理論,由纖維和基體構成的復合材料介質,其阻尼性能可直接由下式計算[12]:

其中:下標n,s分別表示拉壓和剪切方向,m表示基體,Vf表示纖維含量。

復合材料介質的宏觀彈性常數采用多胞模型計算如下[13]:

一般來說,矩陣(13)中的元素ψi是材料內部應力幅值的函數,但是在理想的線性阻尼分析的框架下,可以近似地認為ψi與應力幅值無關,于是遲滯環的形狀為橢圓形,ψi具有定常的特性。

由式(12)積分得:

其中:

類似地,利用式(7)對式(19)進行簡化,得:

其中:

復合材料薄壁梁的耗散能為:

其中:C為復合材料薄壁梁橫截面的4×4阻尼矩陣,其矩陣元素表達式類似于剛度矩陣K,只需要將A(s),B(s),C(s)替換為 Ad(s),Bd(s),Cd(s)即可得到。

1.4 復合材料薄壁梁的自由振動方程及其求解

為了導出復合材料薄壁梁的運動方程,利用Hamilton原理:

其中:動能T可由下式確定:

ρ為材料密度,V表示變形后的梁上任意一點的速度矢量,它與變形梁任意一點的位置矢量:r=(x+u1)i之間滿足關系:表示對時間t求偏導。

由Hamilton方程(21),可導出復合材料薄壁梁的無阻尼自由振動偏微分方程組[15]。采用Galerkin法可進一步導出復合材料薄壁梁的自由振動問題的特征方程:

其中:

{Xm}是復合材料薄壁梁的模態矢量。

1.5 結構阻尼

復合材料薄壁梁的模態阻尼比可以定義為每個振動周期內的耗散(阻尼)能與最大應變能之比:

與式(24)中的第一式類似,復合材料薄壁梁的阻尼矩陣C—具有下列形式:

2 結果及討論

首先,為了驗證雙閉室復合材料箱形薄壁梁剛度系數結果的正確性,表1給出了如圖3所示構型,即沿薄壁梁截面的周線的鋪層方室分別為[+45]2和[-45]2,雙閉室復合材料箱形薄壁梁的剛度系數計算結果,并且與文獻[16]及[17]的結果進行了比較,模型的幾何參數和材料參數均取自兩篇文章。可以看出,三者的計算結果較吻合。

圖2 雙閉室箱型薄壁梁截面Fig.2 Two-cell thin-walled box beam cross section

表1 雙閉室復合材料箱形薄壁梁的剛度系數計算結果Tab.1 Stiffiness results for two-cell beam.

其次,為了檢驗本文建立的復合材料薄壁梁阻尼模型及其近似計算方法的正確性,表2、表3給出在兩種CUS構型,即沿薄壁梁截面的周線的鋪層方式分別為[0]16和[90]16,復合材料箱型截面薄壁懸臂梁的模態阻尼和固有頻率計算結果,并且與文獻[6]的3D剪切梁阻尼有限元結果進行了比較,結構的幾何參數和材料參數均取自文獻[6]。可以看出,二者符合得很好。

表2 復合材料箱形梁的模態頻率和阻尼,L/a=14.36,a/b=5,[0]16Tab.2 Modal frequency and damping of cantilever composite

表3 復合材料箱形梁的模態頻率和阻尼,L/a=14.36,a/b=5,[90]16.Tab.3 Modal frequency and damping of cantilever composite

下面分別計算圖3所示四種構型下雙閉室復合材料薄壁梁的結構阻尼。模型的幾何尺寸為:外截面寬度 a=0.025 m,外截面高度 b=0.025 m,長度 L=1.524 m,單層厚度0.000 127 m,鋪層數為6 層,復合材料基體與纖維性能參數如表4所示。雙閉室箱形復合材料薄壁梁結構如圖3所示。

表4 復合材料性能參數Tab.4 Mechanical properties of composite material

圖3 雙閉室箱型梁模型Fig.3 The diagram of double-cell thin-walled box beam

圖4 表示具有圖3(a)構形的雙閉室復合材料薄壁懸臂梁的第一階彎曲為主模態阻尼隨鋪層角的變化曲線,其中也顯示出纖維含量的影響。結果表明,彎曲為主模態阻尼隨著鋪層角的增加呈現先增加后減少的變化趨勢,最大阻尼發生在鋪層角θ為50°附近,在鋪層角90°的阻尼大約是0°阻尼值的5倍左右。增加纖維含量,可以在一定的鋪層角范圍內提高或者降低薄壁梁的阻尼,但是這種作用的效果并非很明顯。圖5表示具有圖3(a)構形的雙閉室復合材料薄壁懸臂梁的第一階扭轉為主模態阻尼隨鋪層角的變化曲線。顯然,扭轉為主模態阻尼隨著鋪層角的變化趨勢,與彎曲為主模態阻尼是相反的,最小阻尼發生在35°~55°之間,鋪層角90°的阻尼大和0°阻尼值是相同的。鋪層角在0°~90°之間,增加纖維含量可以明顯提高阻尼值。

圖6表示具有圖3(b)構形的雙閉室復合材料薄壁懸臂梁的第一階揮舞彎曲模態阻尼隨鋪層角的變化曲線,結果表明,阻尼隨鋪層角的變化規律與具有圖3(a)構形的雙閉室的結果(見圖4)是類似的。圖7表示具有圖3(b)構形的雙閉室復合材料薄壁懸臂梁的扭轉-揮舞彎曲耦合模態阻尼隨鋪層角的變化曲線,阻尼隨鋪層角的變化以及纖維含量的影響,與具有圖3(a)構形的雙閉室情形下的結果(見圖5)是類似的。

圖8和圖10分別表示對應于圖3(c)和3(d)構形的雙閉室復合材料薄壁懸臂梁的第一階揮舞彎曲模態阻尼隨鋪層角的變化曲線,結果表明,阻尼隨鋪層角的變化規律與具有圖3(a)構形的雙閉室的結果(見圖4)是類似的。圖9和圖11分別表示對應于圖3(c)和3(d)構形的雙閉室復合材料薄壁懸臂梁的扭轉-揮舞彎曲耦合模態阻尼隨鋪層角的變化曲線,阻尼隨鋪層角的變化以及纖維含量的影響,與具有圖3(a)構形的雙閉室情形下的結果(見圖5)是類似的。

圖4 不同纖維含量的薄壁梁(圖3(a))的第一階彎曲為主模態阻尼隨鋪層角變化曲線Fig.4 The first bending modal damping for the beam(Fig.3 (a))vs.ply angle for various fibre volume fraction

圖5 不同纖維含量的薄壁梁(圖3(a))薄壁梁的第一階扭轉為主模態阻尼隨鋪層角變化曲線Fig.5 The first twisting modal damping for the beam(Fig.3 (a))vs.ply angle for various fibre volume fraction

圖6 不同纖維含量的薄壁梁(圖3(b))的第一階彎曲為主模態阻尼隨鋪層角變化曲線Fig.6 The first bending modal damping for the beam(Fig.3 (b))vs.ply angle for various fibre volume fraction

圖7 不同纖維含量的薄壁梁(圖3(b))薄壁梁的第一階扭轉為主模態阻尼隨鋪層角變化曲線Fig.7 The first twisting modal damping for the beam(Fig.3 (b))vs.ply angle for various fibre volume fraction

圖8 不同纖維含量的薄壁梁(圖3(c))的第一階彎曲為主模態阻尼隨鋪層角變化曲線Fig.8 The first bending modal damping for the beam(Fig.3 (c))vs.ply angle for various fibre volume fraction

圖9 不同纖維含量的薄壁梁(圖3(c))薄壁梁的第一階扭轉為主模態阻尼隨鋪層角變化曲線Fig.9 The first twisting modal damping for the beam(Fig.3 (c))vs.ply angle for various fibre volume fraction

圖10 不同纖維含量的薄壁梁(圖3(d))的第一階彎曲為主模態阻尼隨鋪層角變化曲線Fig.10 The first bending modal damping for the beam(Fig.3 (d))vs.ply angle for various fibre volume fraction

圖11 不同纖維含量的薄壁梁(圖3(d))薄壁梁的第一階扭轉為主模態阻尼隨鋪層角變化曲線Fig.11 The first twisting modal damping for the beam(Fig.3 (d))vs.ply angle for various fibre volume fraction

圖12 雙閉室風力機葉片截面翼型輪廓線Fig.12 Two cell airfoil profile composite beam

圖12 給出一個帶有腹板的雙閉室風力機葉片翼型輪廓線,它由文獻[19]翼型曲線的計算公式畫出。翼型輪廓線以及腹板的鋪層方式取[θ]6。梁長為0.762 m,復合材料性能參數取自表4。圖13和14分別表示纖維含量以及鋪層角對第一階揮舞彎曲模態阻尼和第一階扭轉-揮舞耦合模態阻尼的影響規律,可以看出,在上述鋪層方式下,鋪層角的影響與雙閉室箱型截面梁的結果是相似的,但是纖維含量的影響與雙閉室箱型截面梁的結果有所不同。

3 結論

圖13 雙閉室翼型截面復合材料薄壁梁的第一階彎曲為主模態阻尼隨鋪層角變化曲線Fig.13 The first bending modal damping for the beam vs.ply angle for various fibre volume fraction

圖14 雙閉室翼型截面復合材料薄壁梁的第一階扭轉為主模態阻尼隨鋪層角變化曲線Fig.14 The first twisting modal damping for the beam vs.ply angle for various fibre volume fraction

本文基于單層混雜材料的細觀力學阻尼計算方法和多胞模型,綜合采用VAM法和Hamilton原理并且借助于Galerkin法,建立了彎扭耦合復合材料雙閉室薄壁梁的結構阻尼分析模型。研究結果表明,采用本文建立的解析模型及其近似計算方法所獲得的薄壁箱形懸臂梁的模態阻尼與固有振動頻率數值近似解,與已有文獻的有限元解相當一致。雙閉室薄壁復合材料梁的結果阻尼的分析結果表明:

(1)根據本文模型能夠對雙閉室薄壁復合材料梁的結構阻尼性能進行參數分析,以確定影響薄壁梁阻尼性能的主要因素,為進一步改善或提高薄壁梁的阻尼耗散性能,提供有用的信息。

(2)纖維鋪層角對薄壁梁的結構阻尼能夠產生明顯的影響。對于彎曲為主模態來說,相應的阻尼最大值出現在50°鋪層角附近;對于扭轉為主模態,相應的阻尼最小值出現在35°~55°鋪層角附近,并且在上述兩種情形下阻尼隨纖維鋪層角的變化規律是相反的。

(3)隨著纖維含量的增加,可以在一定的鋪層角范圍內提高或者降低薄壁梁彎曲為主和扭轉為主振動模態的結構阻尼,纖維含量的影響效果將取決于雙閉室薄壁復合材料梁的截面的鋪層方式以及構型。

[1]Barlas T K,Van Kuik G A M.Review of state of the art in smart rotor control research for wind turbines[J].Progress in Aerospace Sciences,2010,46(1):1 -27.

[2]Hansen M H Aeroelastic instability problems for wind turbines[J].Wind Energy,2007,10:551 -577.

[3]Chaviaropoulos P K,Politis E S,Lekou D J,et al.Enhancing the damping of wind turbine rotor blades,the DAMBLADE Project[J].Wind Energy,2006,9:163 -177.

[4]任勇生,劉立厚.纖維增強復合材料結構阻尼研究進展[J].力學與實踐,2004,26(1):9 -16.

[5]Suresh R,Malhotra S K.Vibration and damping analysis of thin-walled box beams[J].Journal of Sound and Vibration,1998,215(2):201-210.

[6]Saravanos D A,Varelis D,Plagianakos T S,et al.A shear beam finite element for the damping analysis of tubular laminated composite beams[J].Journal of Sound and Vibration,2006,291:802 823.

[7]Rehfield R W.Design analysis methodology for composite rotor blades[C].Proceedings of the Seventh DoD/NASA Conference on Fibrous Composites in Structural Design.Denver,CO,Grant NAG-2-238,1985.

[8]Chortis D I,Chrysochoidis N A,Saravanos D A.Damped structural dynamics models of large wind-turbine blades including material and structural damping[J].Journal of Physics:Conference Series,2007,75:1 -11.(doi:10.1088/1742-6596/75/1/012076).

[9]Ramkumar K,Ganesan N,Kannan R.Global and local behaviour based composite damping studies on thin-walled box structure[J].European Journal of Mechanics A/Solids,2010,29:253-265.

[10]Berdichevsky V L,Armanios E,Badir A M.Theory of anisotropic thin-walled closed-cross-section beams[J].Composites Engineering,1992,2(5 -7):411 -432.

[11]任勇生,杜向紅,孫雙雙,等.單閉室復合材料薄壁梁的結構阻尼[J].振動與沖擊,2012,21(3):141-146.

[12]Saravanos D A,Chamis C C.Unified micromechanics of damping for unidirection and off-axis fiber composites[J].J CompositesTechnology & Research. JCTRER, 1990,12(1):31-40.

[13]Murthy P L N. Second generation integrated composite analyzer(ICAN)computer code[R].NASA TP 3290,1993.

[14]Badir A M.Analysis of two-cell composite beams[C].Proceedings of the 36th AIAA Structures,Structural Dynamics and Materials Conferences,New Orleans,AIAA 95-1208 -CP,1995.

[15]任勇生,杜向紅,楊樹蓮.風力機復合材料柔性葉片的顫振分析[J].振動與沖擊,2011,30(9):64-69,128.

[16]Cesnik C E S,Shin S J.On the modeling of integrally actuated helicopter blades[J].International Journal of Solids and Structures,2001,38:1765-1789.

[17]Cesnik C E S,Hodges D H.VABS:A new concept for composite rotor blade cross-sectional modeling[J].Journal of American Helicopter Society,1997,42(1):27 -38.

[18]Kaliske M, Rothert H. Damping characterization of unidirectional fiberreinforced polymercomposites[J].Composites Engineering,1995,5(5):551 -567.

[19]王旭東,陳 進,Shen W Z,等.風力機葉片翼型型線集成設計理論研究[J].中國機械工程,2009,20(2):211-213,228.

附錄:剛度系數的求解[14]

圖15 雙閉式截面Fig.15 Schematic of two-cell section

如圖15所示,雙閉式截面被分成了三部分。位移場中提到的一些系數可以通過圖15每個部分求解:

其中:

下標Ⅰ和Ⅱ表示圖中的左右兩個閉室。分別用AeⅠ和AeⅡ表示下標為Ⅰ和Ⅱ的左右兩個閉室的面積。剛度矩陣[K]是4×4對稱陣,其中的元素組成Kij( i,j=1,2,3,4)由橫截面的幾何變化和材料的性質決定,表達式為:

猜你喜歡
模態復合材料結構
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
民機復合材料的適航鑒定
復合材料無損檢測探討
電子測試(2017年11期)2017-12-15 08:57:13
論《日出》的結構
國內多模態教學研究回顧與展望
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
基于HHT和Prony算法的電力系統低頻振蕩模態識別
TiO2/ACF復合材料的制備及表征
應用化工(2014年10期)2014-08-16 13:11:29
由單個模態構造對稱簡支梁的抗彎剛度
計算物理(2014年2期)2014-03-11 17:01:39
主站蜘蛛池模板: 国产精品午夜电影| 99无码中文字幕视频| 99精品久久精品| 国产在线拍偷自揄拍精品| 国内精品九九久久久精品| 2020国产在线视精品在| 97人人做人人爽香蕉精品| 国产极品粉嫩小泬免费看| 亚洲av无码久久无遮挡| 久久精品亚洲热综合一区二区| 99精品热视频这里只有精品7| 国产真实乱子伦视频播放| 人妻熟妇日韩AV在线播放| 欧洲熟妇精品视频| 国产av色站网站| 国产福利2021最新在线观看| 色久综合在线| 欧美综合成人| 伊人久久福利中文字幕| 欧美日韩导航| 无码网站免费观看| 日本a级免费| 久久精品日日躁夜夜躁欧美| A级毛片高清免费视频就| 少妇被粗大的猛烈进出免费视频| 99久久人妻精品免费二区| 3D动漫精品啪啪一区二区下载| 欧美爱爱网| 日本久久久久久免费网络| 亚洲欧美日本国产综合在线| 欧美在线黄| 亚洲中文字幕23页在线| 欧美中文字幕无线码视频| 国内精品视频在线| 成人福利一区二区视频在线| 欧美日韩v| 色婷婷电影网| 六月婷婷综合| 国产无遮挡猛进猛出免费软件| 国产性猛交XXXX免费看| 在线观看精品自拍视频| 青青草一区| 国产18在线播放| 欧美亚洲国产精品第一页| 国产自在线播放| 久久特级毛片| 中文字幕色站| 亚洲精品图区| 国产精品永久不卡免费视频| 99热这里只有成人精品国产| 久久毛片基地| 精品久久国产综合精麻豆| 亚洲欧州色色免费AV| 久久这里只精品热免费99| 国产一级特黄aa级特黄裸毛片| 欧美a在线| 中文无码日韩精品| 91精品福利自产拍在线观看| A级全黄试看30分钟小视频| 日韩视频免费| 亚洲福利视频一区二区| 女人一级毛片| 69av在线| 久久性妇女精品免费| 狠狠干综合| 日韩欧美国产中文| 国产美女91视频| 91啪在线| 99免费在线观看视频| 无码高清专区| 91精品网站| 亚洲人成网站色7777| 日韩经典精品无码一区二区| 亚洲国产成人综合精品2020 | 国产精品自拍露脸视频| 国产精品免费入口视频| 亚洲欧美另类色图| 国产www网站| 漂亮人妻被中出中文字幕久久| 欧美啪啪网| 日韩AV无码一区| 日韩在线视频网|