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

基于缺陷敏感度的加筋短板承載能力

2020-11-06 06:44:08劉存趙冬強
航空學報 2020年10期
關鍵詞:短板有限元變形

劉存,趙冬強

1.航空工業第一飛機設計研究院 強度設計研究所,西安 710089 2.西北工業大學 航空學院,西安 710072

層流機翼設計技術可以使得機翼上表面維持較長的層流區域,從而降低飛機的阻力和使用成本,提高飛機性能。然而,層流技術對翼面結構制造、使用過程中存在的缺陷相當敏感。一旦翼面蒙皮具有初始缺陷,就引起層流的提前轉捩,然而,對于實際翼面蒙皮結構,初始缺陷的存在是難以避免的。而幾何缺陷對薄殼結構后屈曲承載特性的影響是該領域研究的重要方向[1]。一味追求結構減重和承載能力最大所獲得的最優設計往往表現為對缺陷非常敏感,因此基于缺陷敏感度的承載能力評估具有重要意義[2]。幾何缺陷影響結構的極限承載能力,具有幾何缺陷的加筋板結構穩定性問題通常是結構大變形和材料彈塑性的雙重非線性的復雜后屈曲問題,而非簡單的分叉問題。

目前針對幾何缺陷對加筋板穩定性及承載能力的影響問題國內外眾多學者進行了相關理論和試驗研究。Campbell等[3]通過掃描壁板表面實測幾何缺陷,并通過ABAQUS命令*IMPERFECTION將缺陷引入有限元模型,分析了缺陷對蒙皮和筋條厚度均為0.9 mm的薄加筋板“坍塌”載荷的影響。Bernard等[4]研究了3種測量冷軋薄壁鋼板幾何缺陷方法,并評價了各自的優劣。Lanzi[5]實測了2塊壁板的初始幾何缺陷,并對測量數據進行了分析,得到了最大正負缺陷幅值以及幅值絕對值的平均值,通過參考真實位置改變蒙皮長桁節點位置來引入缺陷,同時比較了有無測量缺陷下2塊壁板的試驗和數值模擬的屈曲載荷及破壞載荷。Houston等[6]提出了采用整體次加筋板結構來考慮幾何缺陷對壁板設計的影響。Paulo等[7]測量了不同形狀和幅值的中長柱加筋板初始幾何缺陷,并采用有限元模擬技術預測了具有測量缺陷的加筋板破壞載荷。Xu和Soares[8]使用ANSYS分析了實測缺陷與幾何缺陷公式下窄加筋板的極限強度,結果與試驗吻合高度一致,表明初始缺陷公式對分析極限強度是有效的。加筋板的幾何缺陷形式和幅值未知的情況下,可以用屈曲特征模態來模擬其形態,用制造或者裝配可以接受的缺陷上限來限定其幅值。劉存等[9]用有限元線性屈曲特性值分析獲取特征模態,經過屈曲模態矢量歸一化,乘以缺陷偏移量等計算處理后,用一階屈曲模態代替局部初始缺陷,雖然此方法計算結構的后屈曲極限承載能力時計入初始缺陷的影響,但并不能保證與實際結構的缺陷一致,結果難免存在偏差。Hilburger等[10]基于實測幾何缺陷數據標記典型制造缺陷,引入實測缺陷形狀的有限元分析結果與試驗屈曲載荷符合性良好。文獻[11]指出考慮缺陷的非線性有限元分析可以得到高精度的屈曲和破壞強度,屈曲響應可影響制造的公差量。Rigo[12]和Couto[13]等研究了不同形狀和幅值的缺陷對焊接薄壁結構極限強度的影響。

然而,國內外針對幾何缺陷對加筋板的承載能力的研究仍有所欠缺,以上研究采用單個類型的缺陷公式或者特征模態,而均未對加筋短板的幾何缺陷開展全面的研究和敏感度分析,且對加筋短板的失穩形態缺乏更精細的分析和認識,因此,對幾何缺陷敏感度的深入研究就非常必要了。敏感度分析方法可以研究模型的輸出變化對模型參數的敏感程度。篩選方法常用于篩除對模型輸出影響較弱的輸入,而用于敏感度分析的最初階段。局部敏感度分析方法是指每次分析時只有被研究的輸入變化,而其余輸入固定不變的分析法,也稱OAT(One-factor-at-A-Time) 方法[14],其主要分析方法有:差分法、求導法、格林函數法等。全局靈敏度分析(Global Sensitivity Analysis,GSA)方法檢驗多個輸入的同時變化對模型輸出的影響,并分析每項輸入以及各輸入間的交互效應對模型輸出的影響[15],其主要方法有非參數方法[16-17]、基于方差的分析方法[18]、矩獨立分析法[19]、基于失效概率的矩獨立分析法[20-21]等。

本文采用原理簡單、工程上便于使用的局部靈敏度分析法,結合工程經驗,著重考慮靈敏度系數較大的變量,進行初彎曲、初偏心以及初變形對加筋短板承載能力的影響分析。給出了3種幾何缺陷下加筋短板承載能力的設計許用值以及表征加筋短板初彎曲和初變形缺陷的計算公式,提出了通過加工制造要求的控制來降低幾何缺陷對極限承載能力敏感度的設計建議,對指導加筋短板的結構強度設計具有重要的工程應用價值。

1 加筋短板承載能力分析

1.1 加筋短板模型

軸壓載荷下加筋板后屈曲承載能力高精度的預測需要考慮初始幾何缺陷、材料非線性、幾何非線性以及高精度有限元模型。然而有限元模型的細節,如網格細化、理想化邊界條件、分析方法同樣需要事先進行研究。

選取的模型為長桁兩側各1/2長桁間距寬度的蒙皮的單加筋短板結構。根據研究方案,加筋短板的寬度150 mm,長度100 mm。橫剖面如圖1所示。材料為鋁合金7150-T7751預拉伸厚板,性能參數[22]如下:壓縮彈性模量Ec=73 723 MPa,泊松比μ=0.33,極限強度σb=565 MPa,壓縮屈服強度σ0.2c=530 MPa。

圖1 結構剖面簡圖Fig.1 Cross section of structure

具體剖面參數見表1。其中,t為長桁腹板厚度;t1和b0分別為長桁自由凸緣的厚度和寬度;t2和b1分別為長桁底邊凸緣的厚度和寬度;h為長桁高度;R和R1為倒圓半徑。

表1 剖面參數Table 1 Parameters of cross section

有限元模型建模時采用可以考慮大變形和材料非線性的殼單元。網格劃分時,根據模型收斂性計算,為準確模擬結構因受壓而產生的屈曲波,加筋之間蒙皮沿橫向至少劃分6個單元,加筋腹板沿高度方向至少劃分3個單元,加強筋自由凸緣至少劃分2個單元。為了保證建模與計算的準確性,單元數的設置盡可能使得各單元呈現正方形。

邊界條件表征加筋板的實際受載情況,面內和面外邊界條件都需考慮。選取能夠代表機翼上翼面加筋板的典型結構為計算模型,在壓縮載荷作用下,蒙皮自由邊在與壓縮載荷垂直的方向上的面內變形不受約束,長桁自由凸緣則可認為完全自由的。同時施加的邊界條件應保證加筋板變形后相對兩邊相互平行。由于是整體加筋板,加筋對蒙皮的轉動約束以及支撐作用通過共節點來實現。

有限元模型加載及約束如圖2所示。蒙皮沿橫向劃分140個單元,長桁腹板沿高度劃分60個單元,長桁自由凸緣沿橫向劃分20個單元,模型單元均為殼單元CQUAD4。

圖2 有限元模型加載及約束圖Fig.2 Loading and constrain condition of FEM

采用類型為RBE2的多點約束單元模擬試驗機將載荷直接施加在加筋短板端面,加載點為主動節點,加載點約束XY方向的位移,加筋短板端面為從動節點,從動點約束蒙皮、長桁腹板、長桁凸緣各自的面外位移。約束端約束加筋短板端面三向XYZ位移。考核段約束蒙皮自由邊面外Y向位移。有限元模型在加載點施加強迫位移載荷。

為了準確地預測加筋短板的承載能力,采用材料真實的本構關系。鋁合金7150-T7751板材[14]應力-應變曲線如圖3所示。

圖3 7150-T7751板材應力-應變曲線Fig.3 Stress-strain curve of 7150-T7751 plate

1.2 求解策略

鋁合金加筋板在航空領域飛機結構上使用已久,加筋板屈曲后剛度下降,在達到破壞載荷之前,加筋板進入后屈曲,之后發生連續的失效和坍塌。此領域的設計方法是基于Euler柱屈曲理論和Timoshenko板殼彈性穩定性理論,此類前期方法受邊界簡化和彈性假定等限制,適用范圍有限,與試驗相比,也不總能很好地預測臨界載荷。而采用GMNIA(Geometrically and Material Nonlinear Analysis with Imperfections)分析方法則可高精度的預測,同時具備通用性。

求解加筋板的極限承載能力,可能需要跨越加載過程中載荷—位移曲線的屈曲分叉點,達到結構“坍塌”的極限載荷點,進而獲得較為準確的結構的后屈曲承載能力,因此需要制定求解策略。通常情況下,非線性隱式分析可以較為準確地模擬加筋板結構的后屈曲行為,但當結構出現局部屈曲后,求解步長變小,求解時間激增,收斂性難以保證,需要設置合理的求解參數。Lanzi和Giavotto[23]則用動態顯式分析求解復合材料加筋板的后屈曲特性。顯式求解方法不存在收斂問題,但求解時間受模型復雜度、單元尺寸等參數影響較大。Abramovich等[24]用3種商用有限元求解器ABAQUS/Standard、ABAQUS/Explicit、NASTRAN模擬了帶缺陷的加筋板后屈曲特性,并從計算時間、與試驗數據吻合度等方面做了比較分析。

加筋板的幾何缺陷在未知的情況下,可以用屈曲特征模態來替代。需要采用MSC.Nastran SOL105進行線性屈曲分析,獲取有限元模型的屈曲模態。其次將一階屈曲模態矢量歸一化,乘以缺陷基矢量得到缺陷偏移矢量,采用SPCD模型數據卡施加強迫位移,將缺陷添加到有限元模型。最后基于Newton-Raphson迭代的弧長法(Arc-Length Method)SOL600隱式求解器進行非線性有限元計算。針對加筋短板的有限元模型,采用GMNIA分析方法進行求解。模型使用NLPARM卡定義非線性分析,載荷被分為87個等增量,采用ITER方法控制切線剛度修正,矩陣修正之前的迭代次數為40,每個載荷增量的總迭代限為1 000。采用NLPCI卡定義非線性靜態分析中弧長增量求解策略,弧長法類型為MRIKS,最小弧長比為0.25,最大弧長比為4.0,期望收斂的迭代次數為40,每步最大迭代次數500。

1.3 結果分析

采用上述求解策略,對包含線性屈曲特征值法計算得到的屈曲模態轉化為含初始缺陷的有限元模型進行求解,求解時開啟大變形(LGDISP)考慮幾何非線性,引入材料的彈塑性本構關系,進行非線性迭代分析計算,得到加筋短板的極限承載能力為368.3 kN。

表1中參數的加筋板試驗件進行了5件,破壞載荷的試驗與有限元結果比較見表2,同時給出了試驗載荷平均值,表中序號對應試驗件編號。PTEST為試驗值,PFEA為有限元計算值,ε為有限元計算值相對于試驗值的偏差:

表2 破壞載荷的試驗與有限元結果比較Table 2 Comparison of failure load between test data and calculation results

ε=(PFEA-PTEST)/PTEST×100%

(1)

可見,采用有限元計算的加筋壁短板的破壞載荷與試驗值誤差均在10%以內,與試驗值平均值誤差為-0.46%。

提取載荷及位移的計算結果,繪制載荷-位移曲線如圖4所示。

圖4 加筋短板加載端載荷-位移曲線Fig.4 Load-displacement curve of loading end of stiffened short plate

由加筋短板加載端載荷-位移曲線可見,曲線在初始階段呈線性,加載到368.3 kN時,曲線達到頂點,曲線斜率由正值變為負值,即為結構的極限載荷。隨后位移繼續增加,載荷下降,結構喪失承載能力。進一步分析加筋短板的破壞模式可見,長度方向1/2處蒙皮的面外變形最大且一側蒙皮突起而另一側蒙皮凹陷,兩者組成一個完整的屈曲波;長桁腹板有明顯的屈曲變形并向一側彎曲。從整個破壞過程來看,長桁剖面發生局部失穩,壓縮載荷繼續增加,增加部分由較剛硬的角區承受,直到應力增加到足夠的數值造成破壞。而從蒙皮和長桁腹板的失穩先后次序來看,蒙皮先失穩。有限元模擬失效時刻位移云圖如圖5所示,試驗第1#件的失效形貌如圖6所示。可見兩者的失效模式一致。

圖6 試驗失效形貌Fig.6 Failure morphology of specimen

2 幾何缺陷對加筋短板承載能力敏感度分析

2.1 幾何缺陷的表征

軸壓載荷下的無缺陷擾動的薄板,其屈曲特性只可能發生分叉型失穩(Snap-through),而非此類結構設計要追求的極限承載能力(Collapse),只有施加面外的初始缺陷擾動,才能將分叉失穩轉化為與之接近的極限承載能力進行非線性數值計算。初始缺陷在加筋板的屈曲特性中扮演重要角色[25],選擇合適的初始缺陷,獲得合理的極限承載結果。對于加筋板而言,初始幾何缺陷可以分為整體幾何缺陷和局部幾何缺陷。具體而言,初始幾何缺陷主要包括加筋板初始彎曲、載荷初始偏心和加筋板的初始變形。在有限元模型中,加筋板的整體幾何缺陷往往以初彎曲的方式計入,通過將模型節點按照正弦波形式的偏移來實現。對于局部幾何缺陷的施加,將柱型、板型和側移型缺陷分別賦予蒙皮和長桁,同樣通過施加強迫位移的方式來施加。

對于加筋板來說,在制造、安裝、運輸和服役過程中不可避免地產生某種幾何缺陷,因此在進行有限元仿真分析時,引入幾何缺陷是符合工程實際的。除了用線性屈曲特征值法計算得到的屈曲模態作為初始缺陷外,關于初始缺陷的規范[26-27]多來自國際船級社協會(IACS)和挪威船級社(DNV),而對于專業的航空工程技術人員,航空加筋板結構引入何種模式、多大幅值的幾何缺陷是值得花費時間和精力去深入研究的問題。為此,以加筋短板為例開展初彎曲、初偏心、初變形等缺陷對極限承載能力的敏感度分析。

2.2 初彎曲

基于歐洲鋼結構協會(ECCS)委員會對通常情況下比較平直的鋁合金擠壓型材的初彎曲的測量結果,初彎曲總小于L/1 300,L為型材總長度,ECCS委員會在計算壓桿的穩定性時采用了初彎曲為L/1 000的正弦擾曲線[28]。為研究初彎曲對加筋短板極限承載能力的影響,初彎曲采用正弦擾曲線式(2),考慮加筋短板3種初彎曲方向,對擾曲線的幅值進行敏感度研究。

(2)

式中:wb為初彎曲擾曲線;a為加筋短板長度,mm;x為加筋短板沿長度方向坐標,mm。

根據圖6所示的加筋短板試驗件加工要求:零件A面(端面)與B面(端面)平行度0.05,長桁軸線面與A面、B面垂直度0.08。可計算參考擾曲線幅值分別為a/2 000、a/1 250。為考慮加筋短板初彎曲帶來的影響,假設加筋短板初彎曲的幅值分別為a/1 000、a/1 250、a/2 000、a/2 500、a/3 000、a/3 500、a/4 000和a/4 500。在空間中的加筋短板彎曲方向有3種:彎向蒙皮正面、彎向蒙皮反面(即長桁凸緣方向)以及彎向腹板正面/反面。由于文中計算的加筋短板結構關于長桁腹板軸線對稱,故彎向腹板的正面和反面計算結果相同。

各種幅值下,3種初始彎曲方向的載荷-位移曲線如圖7所示。圖中:skin-代表彎曲圓弧中心在蒙皮側(即彎向蒙皮正面),skin+代表彎曲圓弧中心在長桁側(即彎向蒙皮反面),web代表彎曲圓弧中心在腹板一側(即彎向腹板正面/反面)。

圖7 初彎曲加筋短板加載端載荷-位移曲線圖Fig.7 Load-displacement curves of loading end for initial bending of stiffened short plate

由圖7可見,加筋短板初彎曲的幅值為a/1 000、a/1 250、a/2 000、a/2 500、a/3 000、a/3 500、a/4 000和a/4 500,3種彎曲類型中skin-和skin+線性段斜率一致,web線性段斜率略高于skin-和skin+。說明在相同的初彎曲幅值和壓縮位移下,加筋短板彎向蒙皮或彎向長桁自由凸緣方向時的應變小于彎向腹板方向的應變,在加筋短板設計中,更應關注彎向蒙皮或長桁自由凸緣方向的初彎曲,此方向也正是加筋板在機翼盒段中的主要受力方向。在機翼盒段的上翼面,加筋板受機翼翼型的制約,實際結構中初彎曲彎向長桁自由凸緣一側。因此,研究skin+的承載能力具有直接的工程應用意義。

為了便于分析3種初彎曲幅值對加筋短板極限承載能力敏感度,將計算結果匯總見表3。其中,CDS(Curve Decrease Significantly)用于判斷載荷—位移曲線是否出現拐點,即是否獲得極限載荷,是代表有限元仿真結果得到了極限載荷,否代表未得到極限載荷。

由表3可見,skin+在加筋短板初彎曲的幅值為a/2 500時,極限承載能力達到最大值,為此,建議航空結構鋁合金加筋板結構的初彎曲幅值取a/2 500。此時式(2)修正為式(3):

表3 初彎曲對極限承載能力的影響Table 3 Effect of initial bending on ultimate bearing capacity

(3)

此外,3種彎曲類型中skin+的極限載荷最小,在結構設計階段,初彎曲彎向長桁自由凸緣方向的極限承載能力可用作加筋短板承載能力的設計許用值,而在結構件的服役使用中,應盡量避免。

2.3 初偏心

為研究加筋短板初始偏心對極限承載能力的影響,計算時引入了有限元模型的一階屈曲特征模態形狀。初偏心距取0.025倍、0.050倍、0.075倍、0.100倍、0.125倍的回轉半徑。初偏心的方向定義:向長桁自由凸緣側偏心為正偏心,向蒙皮側偏心為負偏心。形心、偏心后的位置及有限元計算的極限承載能力見表4。其中,回轉半徑為23.29 mm。

表4 初始偏心及極限承載能力計算結果Table 4 Calculation results of initial eccentricity and ultimate bearing capacity

由表4可見,各種匹配偏心距和偏心方向下,有限元仿真得到的極限承載能力均小于試驗破壞載荷平均值370.0 kN,說明偏心是不利因素。在結構設計中應避免截面出現偏心,在強度計算中應計入結構真實的剛度特性,并考慮偏心帶來的附加載荷。

正偏心加載端載荷-位移曲線如圖8所示,負偏心加載端載荷-位移曲線如圖9所示。可見,無論載荷施加正偏心還是負偏心,對加筋短板而言,線性段的斜率是一致的,極限承載能力相差不大,均在2%以內。加筋短板長細比小于20,處于短柱范圍內,有限元模型通過RBE2將載荷施加在端面上,初始偏心值對極限承載能力的影響有限。

圖8 正偏心加載端載荷-位移曲線Fig.8 Load-displacement curves for positive eccentric loading end

圖9 負偏心加載端載荷-位移曲線Fig.9 Load-displacement curves for negative eccentric loading end

2.4 初變形

由于加工、制造等方面的原因,加筋板難以避免地帶有一定的初始變形,而這種初始變形通常對壓縮載荷作用下加筋板的大位移或后屈曲特性有明顯的影響。加筋板的初始變形如圖10所示。加筋板常見的初始變形缺陷包括柱型變形缺陷(圖10(a)),板型變形缺陷(圖10(b)和加強筋側移型變形缺陷(圖10(c))。

圖10 加筋板初始變形Fig.10 Initial deformation diagram of stiffened panel

加筋短板的柱型變形缺陷、板型變形缺陷和加強筋側移型變形缺陷的變形函數表達式分別為式(4)~式(6)。其中板型變形缺陷作用在加筋短板的蒙皮上,側移型變形缺陷作用在加筋短板的加強筋上,柱型變形缺陷則是作用在加筋短板上,與初彎曲是相同的,2.2節已經進行了研究,并給出了推薦值。

板型變形缺陷:

(4)

柱型變形缺陷:

(5)

側移型變形缺陷:

(6)

文獻[18-19]給出了船體結構wp、ws的常用推薦值,而航空結構加筋板變形缺陷的幅值大小是值得深入研究的問題。為此,采用有限元方法計算了加筋短板的板型變形缺陷和加強筋側移型變形缺陷在不同幅值下的極限承載能力,得到了各種變形缺陷下幅值對加筋短板極限承載能力的影響程度。不同幅值的板型變形缺陷下加筋短板極限承載能力如圖11所示,不同幅值的側移型變形缺陷下加筋短板極限承載能力如圖12所示。

圖11 板型變形缺陷下加筋短板載荷-位移曲線Fig.11 Load-displacement curves of thin-horse mode with initial imperfection of stiffened short plate

圖12 側移型變形缺陷下加筋短板載荷-位移曲線Fig.12 Load-displacement curves of side-ways with initial imperfection

由圖11和圖12可見,加筋短板初變形的wp、ws在不同的取值下,線性段斜率一致,而極限承載能力是不同的,說明初變形的wp、ws直接影響加筋短板的承載能力。因此,在加筋短板設計中,更應關注wp、ws的取值。圖6中的加筋短板采用航空領域常用的制造加工公差量:蒙皮、筋條厚度公差為±0.15 mm,筋高公差±0.2,零件平面直線度每100 mm間隙不大于0.1 mm。經計算可得wp=4.28%、ws=5.82%。

為了分析板型變形缺陷和側移型變形缺陷下幅值對加筋短板極限承載能力影響,將計算結果匯總見表5。

從表5中可知,加筋短板的板型變形缺陷的幅值wp為a/3 000~a/5 000時,極限承載能力相對穩定,且接近試驗平均值。而側移型變形缺陷ws的在a/4 000~a/10 000區域變化較小且接近試驗平均值。為此,建議航空飛機鋁合金加筋板結構板型變形缺陷幅值取a/5 000,側移型變形缺陷取值為a/6 000。此時的式(4)修正為式(7),式(6)修正為式(8):

表5 初變形缺陷對極限承載能力的影響Table 5 Effect of initial deflection on ultimate bearing capacity

(7)

(8)

將試驗件長度100 mm代入式(7)、式(8)可得幅值wp=2.0%、ws=1.67%。由此可見,要獲得板型變形缺陷和側移型變形缺陷下加筋短板的高承載能力,必須提高制造加工的精度:蒙皮厚度公差為±0.07 mm,零件平面直線度每100 mm間隙不大于0.028 mm。

3 結 論

通過加筋短板承載能力的分析和幾何缺陷對加筋短板承載能力敏感度研究,得到以下結論:

1) 幾何缺陷對加筋短板承載能力均有不同程度的影響,進行加筋短板極限承載能力分析時,應計入幾何缺陷的影響。

2) 通過初彎曲對加筋短板承載能力敏感度分析,給出了航空結構鋁合金加筋短板初彎曲的推薦公式,提出了結構設計許用值。

3) 通過初偏心對加筋短板承載能力敏感度分析,得到了在均勻受壓載荷下,偏心距和偏心方向對加筋短板極限承載能力的影響有限的結論。

4) 通過初變形對加筋短板承載能力敏感度分析,給出了航空結構鋁合金加筋短板的板型變形缺陷幅值以及側移型變形缺陷的推薦公式,提出了通過控制加工制造公差量來提高結構承載能力的設計建議。

猜你喜歡
短板有限元變形
執行“強制休假”還需“補齊三個短板”
談詩的變形
中華詩詞(2020年1期)2020-09-21 09:24:52
立新標 補齊勞動防護短板
勞動保護(2019年7期)2019-08-27 00:41:04
DCT的優勢與短板并存
汽車觀察(2018年12期)2018-12-26 01:05:40
“我”的變形計
補齊短板 建好“四好農村路”
中國公路(2017年17期)2017-11-09 02:25:00
例談拼圖與整式變形
會變形的餅
磨削淬硬殘余應力的有限元分析
基于SolidWorks的吸嘴支撐臂有限元分析
主站蜘蛛池模板: 免费人成在线观看成人片| 日韩经典精品无码一区二区| 天堂va亚洲va欧美va国产| 国产一区亚洲一区| 91视频精品| 亚洲女同欧美在线| 久久成人免费| 国产精品自在线天天看片| 久久综合婷婷| 中文字幕 91| 波多野结衣视频网站| 国产日本一线在线观看免费| 国产三级国产精品国产普男人| 正在播放久久| 亚洲熟妇AV日韩熟妇在线| 又粗又大又爽又紧免费视频| 久久一日本道色综合久久| 国产成人精品一区二区三区| 国产精品偷伦视频免费观看国产 | 免费无码在线观看| 午夜欧美理论2019理论| 免费毛片网站在线观看| 亚洲日本在线免费观看| 91青青草视频| 日本国产精品| 天天综合网色| 午夜视频在线观看区二区| 国产成人8x视频一区二区| 国产一区二区三区日韩精品| 成年人国产视频| 国产高清自拍视频| 国产日韩欧美视频| 在线观看亚洲人成网站| 久久精品人人做人人爽97| 亚洲精品在线影院| 在线观看免费AV网| 这里只有精品免费视频| 国产精品久久自在自线观看| 精品国产成人a在线观看| 毛片基地美国正在播放亚洲| 青青青视频91在线 | 精品99在线观看| 国产爽妇精品| 97在线碰| 热re99久久精品国99热| 国产精品三区四区| 尤物视频一区| 欧美亚洲日韩中文| 免费看一级毛片波多结衣| 美女免费黄网站| 99视频有精品视频免费观看| 在线免费观看a视频| 尤物精品视频一区二区三区| 真实国产乱子伦视频 | 91精品国产情侣高潮露脸| 欧美在线导航| 欧美一级色视频| 久久综合九色综合97婷婷| 伊人色在线视频| 天天综合网在线| 97在线公开视频| 高清免费毛片| 国产黄视频网站| 看你懂的巨臀中文字幕一区二区| 黄色a一级视频| 网久久综合| 好紧好深好大乳无码中文字幕| www亚洲天堂| 精品国产91爱| 亚洲综合久久成人AV| 色综合中文| 亚洲Aⅴ无码专区在线观看q| 国产91丝袜在线播放动漫| 国产真实乱人视频| 国产91蝌蚪窝| 欧美日韩免费观看| 综合天天色| 超碰91免费人妻| 伊人久久大香线蕉综合影视| 欧美激情成人网| 91无码人妻精品一区二区蜜桃 | 激情综合网址|