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

輕鋼屋面自攻釘節點在臺風作用下時變可靠度研究

2021-12-01 12:31:22楊娜劉威白凡葛瀚文
湖南大學學報(自然科學版) 2021年11期
關鍵詞:風速

楊娜,劉威,白凡,2,葛瀚文

(1.北京交通大學 土木建筑工程學院,北京 100044;2.結構風工程與城市風環境北京市重點實驗室(北京交通大學),北京 100044)

隨著我國制造業的發展,輕鋼結構以其自重輕、工廠預制化程度高等特點,在廠房建筑中被廣泛應用[1-2].其中部分輕鋼屋面系統是由檁條及壓型鋼板通過自攻釘連接構成.由于屋面系統材質輕、柔度大,在局部屋面板承受較大的風吸力后,其自攻釘節點處往往會發生拔出破壞,從而導致整個輕鋼屋面系統失效[3-4],因此,針對自攻釘節點的災后性能評估方法亟待研究,從而為后期的屋面節點體系設計與服役期限評估提供理論支撐.

目前,學者們從不同層面對自攻釘節點的抗風性能進行了研究.在靜力抗風性能方面,Xu 等[5-6]和Mahendran[7]對澳大利亞常用的輕鋼屋面自攻釘節點進行了試驗研究,驗證了輕鋼屋面系統的失效是由于自攻釘節點處的局部破壞所致,分析了不同類型壓型鋼板的節點破壞過程及其抗風性能的影響因素.在節點的疲勞抗風性能方面,Morgan 等[8]和Beck 等[9]探討了自攻釘連接式輕鋼屋面發生疲勞損傷破壞的可能性,并通過疲勞試驗表明,臺風作用會使自攻釘節點周圍的鋼板局部破壞.在以往試驗中,雙跨足尺的屋面模型是研究自攻釘節點疲勞所常用的試驗構件.然而,雙跨足尺的試驗模型對于建造工藝要求較高,且試驗成本高、耗時長.因此,Mahendran[10]通過研究模型較小且更易于操作的單節點試驗模型替代了傳統雙跨足尺模型并得到廣泛應用.Xu[11]通過對循環荷載范圍與均值分布的數據分析,驗證了在疲勞損傷計算中,參考循環荷載范圍的意義遠大于參考循環荷載均值.伴隨風壓系數循環范圍與均值的上升,風荷載循環次數集中在均值等于1/2 范圍處.在輕鋼系統的災后性能評估方面,宋芳芳等[12]分析了臺風作用期間輕鋼結構的破壞過程,提出了風災破壞全過程的預測分析方法等.黃國慶等[13]、趙明偉等[14]建立了構件抗力模型,通過預先設定構件的破壞準則,基于蒙特卡羅模擬,得到了各類構件及整體的超越破壞概率曲線.Mahendran 等[15]通過恒幅循環測試對自攻釘連接的檁條-腹板系統進行研究,給出了輕鋼屋面自攻釘結構的抗疲勞曲線以及疲勞破壞設計方程.Myuran 等[16]、葛瀚文[17]結合靜力及疲勞試驗,考慮了在風致作用下自攻釘節點的靜力極限承載力受疲勞的影響,得到其靜力-疲勞統一方程,從而為該類節點抗疲勞設計和評估提供參考.

在時變可靠度研究方面,張建仁等[18]在既有混凝土橋梁荷載和抗力的時變性基礎上,建立了在役橋梁的時變可靠度計算模型.王草等[19]考慮非平穩車載過程研究了橋梁時變可靠度評估新方法,并給出了考慮非平穩車載過程的橋梁結構時變可靠度顯式計算公式.劉威等[20]基于大氣銹蝕模型和簡支檁條的屈曲計算理論,結合時變可靠度方法,提出了銹蝕檁條在臺風作用下的可靠度計算模型.以上研究的共同特點是結構抗力的退化均服從連續分布,并均能推導出對應結構的時變可靠度顯式計算公式.而對于本文主要研究的臺風致損傷,輕鋼屋面自攻釘節點抗力退化隨每次臺風荷載的發生而突變.此時基于顯式的時變可靠度公式推導將不再適用,因此選擇蒙特卡羅方法進行研究.

綜上所述,本文首先通過模擬臺風(風向和風速)方法,即基于風壓系數時程數據,通過雨流分析法得到臺風通過期間的風壓系數矩陣;再結合自攻釘節點的承載力隨循環荷載幅值的退化規律,計算出不同風速下損傷值的變化情況,從而得到自攻釘節點在臺風作用下的抗力退化模型;繼而進行基于概率的自攻釘節點抗風性能評估,分析并預測各沿海城市自攻釘節點的時變可靠度指標及服役壽命.

1 時變可靠度計算模型

本文通過蒙特卡羅模擬,設計了屋面自攻釘節點在臺風致疲勞作用下的時變可靠度計算程序.圖1(a)給出了自攻釘節點疲勞破壞的時變可靠度計算模型,即自攻釘節點抗力R 及臺風荷載S(不考慮其他荷載效應)的變化模型.每當輕鋼屋面系統遭遇一次臺風荷載Si,自攻釘節點抗力Ri下降ΔRi.fR(r)、fS(s)分別為抗力、荷載的概率密度函數,μR(t)、μS(t)分別為抗力均值、荷載均值隨時間的變化關系.

計算流程如圖1(b)所示,其中Z=R-S 為結構承載力極限狀態方程,對應的極限狀態是輕鋼屋面系統自攻釘節點的疲勞破壞.具體步驟如下.

圖1 時變可靠度計算方法Fig.1 Time-dependent reliability assessment method

步驟1:輸入抗力信息及臺風概率模型.

步驟2:對于單次臺風事件L(j),隨機生成對應的臺風發生次數n 及對應的n 個臺風荷載Sn.

步驟3:計算由每次臺風荷載Sn引起的抗力退化值ΔRn.

步驟4:每次臺風過后,檢查極限狀態方程Zn,若Zn>0,則繼續計算,否則停止計算,并將事件L(j)記錄為失效事件.

步驟5:使用蒙特卡羅模擬,重復計算100 萬次,使計算結果趨于穩定的失效概率,然后將失效事件總次數除以事件發生的總次數,得出該抗力與荷載分布模型下的失效概率.

1.1 基于概率分布的臺風荷載模型

1.1.1 臺風發生次數

在沿海地區,輕鋼屋面系統經常遭遇風荷載影響.本文采用平穩Poisson 過程來模擬臺風發生,即在任一時間段(0,tL]內,臺風事件發生次數N(t)=x 的概率為:

式中:x 為自然數;λ 為Poisson 強度,即臺風頻率.

1.1.2 臺風風速分布

在臺風致疲勞作用下,不同風速導致的承載力退化程度不同.而Weibull 概率分布在風荷載強度分析中廣泛應用,其雙參數曲線能擬合實際風速分布模型,故臺風風速的累積分布函數可表示為:

式中:v 為10 min 計平均風速,m/s;β 和u 分別為服從Weibull 分布的形狀及尺度參數.

1.1.3 風荷載標準值

參考《建筑結構荷載規范》[21],風荷載標準值wk(單位為kN/m2)為:

式中:βgz為陣風系數;μs1為風荷載局部體型系數;μz為風壓高度變化系數;w0為基本風壓,kN/m2.

1.2 臺風致自攻釘節點抗力退化模型

為研究臺風期間屋面自攻釘節點的疲勞性能,從而得到其損傷值與風速的關系.本文參照Jancauskas 等[22]的模擬臺風方法,首先基于風洞試驗的屋面時程風壓系數數據,結合臺風通過時的風速和風向信息,利用雨流分析法計算臺風過境期間的實際尺寸屋面風壓系數矩陣.結合自攻釘節點抗力退化特性,計算不同風速下自攻釘節點由臺風作用導致的抗力退化情況,從而得到輕鋼屋面中自攻釘節點抗力退化模型.

1.2.1 風洞實驗屋面風壓系數矩陣

本文屋面風壓系數時程數據基于東京工藝大學(Tokyo Polytechnic University,TPU)的空氣動力學數據庫[23],該數據庫基于模擬郊區低矮建筑風洞試驗.其中,試驗縮尺模型比例為1 ∶100,速度縮尺比例為1 ∶3,時間縮尺比例為3 ∶100.模擬風場的風速對應于實際10 m 高度處22 m/s.為考慮風荷載作用下輕鋼屋面的最不利位置,選擇平屋面屋角處為測試點.采樣頻率為15 Hz,采樣時長600 s.包括了0~180°風向內間隔為15°的13 組數據,每組數據由9 000 個時程風壓系數組成.房屋模型的長、寬、高分別為240 mm、160 mm、40 mm.

在估算結構因疲勞損傷而導致的抗力退化時,諸如Minner 損傷準則等的線性損傷理論是基于荷載的循環次數.因此,需要將本文獲取的風壓系數時程數據轉換成循環荷載作用次數的分布形式.利用雨流分析法[24],可將風荷載時程曲線轉換為若干個荷載循環,進而用于抗力退化值的計算.

雨流分析以雙參數法為基礎,基本原理如圖2所示.對于風壓系數時程記錄的每個荷載循環,根據該循環荷載的均值和范圍進行分類.通過雨流分析法得到了風向為0~180°,間隔為15°的13 個風壓系數循環矩陣.矩陣的第一行和第一列分別為由“風壓系數”表示的循環荷載的范圍和均值,每個單元內給出了對應的荷載循環次數.表1 列出了45°風向處的風壓系數矩陣.有些正壓循環單元(屋面向下風壓作用),對屋頂造成的疲勞損傷可忽略[22],故不包括在分析矩陣中.

圖2 雨流分析法原理圖Fig.2 Principle of the rain flow count method

表1 45°風向處風壓系數矩陣Tab.1 Matrix of wind pressure coefficient at 45°wind direction kPa

1.2.2 模擬臺風法

1)模擬臺風概述

通常僅根據最大風速來指定設計風,而疲勞對應的是一段時間內風壓的作用,因此有必要根據風速和風向在其持續時間內的變化來指定設計臺風,假設為5 h[22].在臺風通過期間,通常會在風眼路徑的邊緣生成最大破壞程度的風荷載,其風速根據臺風行進速度矢量進行疊加.

臺風過境時,風眼到達不同位置的風速取決于中心壓力p、臺風風速最大時的半徑R′、風眼前進速度U 以及地形和其他氣象因素.當r′=R′時(如圖3所示),所產生的峰值風速方向與最大峰值壓力方向一致.根據從模型房屋屋面角部得到的不同風向下的風壓系數時程矩陣,當臺風行進2.5 h 后,模型房屋位于風眼路徑的邊緣,此時行進速度與切向速度夾角a 為0°,臺風風向角度θ 為75°,也就是自攻釘節點的最不利風向情況.根據計算房屋與臺風的相對運動,得到風眼距房屋不同位置處的風向及風速.

圖3 模擬臺風通過期間的風速和風向分析模型Fig.3 Analytical model of wind speed and direction during simulated typhoon passage

由于自攻釘節點的疲勞損傷值與對應荷載幅值下循環次數有關,而根據風洞試驗數據得到的風壓系數矩陣中,荷載循環次數取決于分析時間間隔的時長.文獻[22]指出,在臺風過境期間,選取分析時間間隔(Δt)為15 min 最接近真實工況,其中臺風總作用時間5 h 內的第一次、最后一次分析時間間隔選為7.5 min.并假定在該時間間隔內,建筑物所受到的風向角度和平均風速恒定.

2)模擬臺風風速的確定

為了確定模擬臺風的風速,Walker 等[25]給出經驗公式用于計算距臺風中心為r′,10 m 高度處的最大10 min 計平均風速V.

式中:Vp=,p 為臺風中心大氣壓,kPa;R′為峰值風速處的半徑,km;切向速度Vp(R′/r′)k及前進速度Uf為向量相加(見圖3);C 和k 為常數,可由臺風數據擬合得到;K 為地形摩擦損失因子;以下參數假定[11]為當r′=R′時出現的峰值風速參數,適用于郊區一類的地形條件,陣風因子取1.7,3 s 最大峰值風速為70 m/s.其中,臺風中心大氣壓p=93 kPa,風眼前進速度Uf=2 m/s,峰值風速處半徑R′=25 km,C=4.5,k=0.67.

3)模擬臺風風速及風向隨過境時間的關系

基于上述模擬臺風模型及特征:(a)可計算臺風通過期間不同時段的風速大小.其中,由于風洞試驗模擬的是4 m 高度處,所以還需要將經驗公式計算的10 m 高度處10 min 計風速按公式(6)轉化為房屋模型所對應的4 m 高度處的10 min 計平均風速;(b)利用三角幾何關系和向量加法法則可計算不同時段的風向.

以第11 個分析時間間隔,即圖3 中T=2.5 h、切向速度與行進速度夾角為75°、r′=R′=25 km 這一情況為例,計算可得10 m 高度處最大10 min 計平均風速約為42 m/s,即臺風最大切向速度為40 m/s,前進速度Uf為2 m/s.按照《建筑結構荷載規范》,轉化為房屋模型所對應4 m 高度處10 min 計平均風速為42 m/s,即:

式中:Vz為所求高度z m 處的10 min 計平均風速,m/s;V10為10 m 高度10 min 計平均風速,地形條件為郊區,風壓高度變化系數截斷高度為10 m,風速剖面指數α 取0.15.

式(7)和式(8)分別給出了本文模擬臺風方法所得到的風速(4 m 高度處10 min 計)及風向的解析解:

式中:θ 為風向角度;a 為行進速度與切向速度夾角;r為建筑物至風眼距離,m;d 為建筑物與臺風風眼在路徑上的投影距離,m;T 為臺風作用時間,s;n 為分析時間間隔,s;n 為分析次序.且式(7)(8)所涉及的其他距離單位均為m,時間單位均為s.

圖4 給出了模擬臺風通過期間,3 s 計極值風速50~85 m/s 下,陣風因子為1.7 時,4 m 高度處屋面角位置10 min 計平均風速和風向變化.

由于臺風最大切向風速與前進速度Uf存在較大差距,使得風向夾角受極值風速影響較小,故在本次分析中,極值風速對各風向角度影響趨于同一條直線.

1.2.3 基于模擬臺風法轉換的實際尺寸屋面風壓系數矩陣

根據風洞試驗時程風壓系數數據,對每個時間間隔內的風壓系數矩陣(見表1)進行如公式(9)的轉換,從而得到實際尺寸下屋面風壓系數矩陣(見表2),具體轉換方法如下:

1)通過對1.2.1 節13 個不同風向角度的屋面風壓系數矩陣之間選取相近的矩陣進行線性插值,從而確定圖4 中21 個時間間隔內對應風向角度的風壓系數矩陣.

2)使用式(9)將風壓系數矩陣(荷載循環次數)從風洞模型條件轉換為實際尺寸條件.

式中:Nf為實際尺寸屋面風壓系數矩陣的循環次數;Nwt為風洞實驗屋面風壓系數矩陣的循環次數;V1為圖4 中的最大10 min 平均風速,m/s;T1為速度V1的持續時間,s,即分析時間間隔;Vwt為風洞中的平均風速,約為10 m 高度處22 m/s;Twt為風洞作用時間,約為實際條件下20 000 s.

對于15 min 分析時間間隔,通過上述計算方法可以得到21 個實際尺寸建筑的風壓系數矩陣.以3 s 極值風速為70 m/s 的臺風為例,45°風向下的實際尺寸風壓系數矩陣見表2.通過觀察風壓系數矩陣,發現大量循環集中在低循環均值和低循環范圍處.隨著循環均值或循環范圍的增大,循環次數變得更小.

表2 45°實際尺寸風壓系數矩陣Tab.2 Matrix of wind pressure coefficients for actual dimensions at 45° kPa

1.2.4 臺風致自攻釘節點抗力退化模型

根據1.2.3 節可得到一次臺風過境期間21 個時間間隔(對應于21 個風向角度)的實際尺寸建筑風壓系數矩陣.為結合抗力退化值與荷載的關系,需將風壓系數Cp轉化為已知風速下的風荷載值w,kN/m2:

式中:w0=v2/1 600.

考慮常見輕鋼屋面檁條間距為1 m,V125 型壓型鋼板單跨長度為0.25 m,則風荷載作用面積為0.25 m2.為便于計算抗力退化值,將式(10)中風荷載值w 寫為用力表達的形式,即:將表2 風壓系數矩陣中的第一行循環荷載范圍寫成力表達的形式:

基于相關自攻釘節點疲勞和靜力的統一試驗結果[17],在不同荷載范圍作用下,單次循環荷載作用下的抗力退化值ΔF(單位為kN)與荷載范圍w′有如下關系:

對于每一個轉化后的實際尺寸屋面風壓系數矩陣,通過累積該矩陣中不同循環荷載范圍w′下的抗力退化值ΔF·Nf,可以得到n 個循環荷載范圍下的累積抗力退化值[22]:

對于以某一速度過境的臺風,通過計算可得21個時間間隔的實際尺寸屋面風壓系數矩陣的抗力退化值.以3 s 計極值風速70 m/s 為例,給出了其抗力退化的時程分布(如表3 所示).其中為便于計算,表格中風速均為10 min 計風速.通過疊加每個時間間隔內的抗力退化值,進而得到抗力退化值隨臺風作用時間的分布曲線.

表3 70 m/s 臺風致抗力退化值分布Tab.3 Distribution of resistance degradation value caused by typhoon of 70 m/s

基于上述計算流程,可計算出不同臺風風速下自攻釘節點的抗力退化值隨時間的分布,進而得到結構在不同風速臺風侵襲后自攻釘節點承載力的下降情況,通過擬合可確定其退化模型,如圖5 所示.

圖5 臺風致自攻釘節點抗力退化模型Fig.5 Resistance degradation model of self-tapping screw joints caused by typhoon

最后通過擬合得到臺風致自攻釘節點抗力退化數學模型,如式(14)所示,線性相關系數為0.995 8.

該數學模型表示輕鋼屋面自攻釘節點在服役期間,每經歷一次以10 min 計風速為vi(單位為m/s)的臺風,其抗力退化ΔRi(單位為kN).

2 算例分析

根據上述臺風荷載模型以及臺風致節點抗力退化模型,以我國東南沿海11 個城市的輕鋼平屋頂屋面體系為例,分析自攻釘節點在考慮臺風致抗力退化作用下的時變可靠度.選取節點形式為V125 型,場地為郊區,地面粗糙程度為B 類.陣風系數取1.7,局部體形系數取1.65,風壓高度變化系數取1.0.輕鋼屋面房屋尺寸為24 m×16 m×4 m(長×寬×高),檁條間距為1 m,V125 型壓型鋼板單跨長度為0.25 m,由此可使用公式(3)(4)得到臺風風速與風荷載的關系.其中,自攻釘節點初始抗力值為1.5 kN(不考慮初始抗力的隨機性).

2.1 各沿海城市臺風概率模型

通過統計我國東南沿海11 個城市的臺風數據[26],其中包含了各地臺風重現期及相應的臺風風速,通過聯立方程計算出服從Weibull 分布的雙參數u 和β 以及臺風頻率λ,結果整理如表4 所示.

表4 東南沿海城市臺風參數Tab.4 Typhoon parameters in southeast coastal cities

2.2 各沿海城市時變可靠度計算

根據不同地區的臺風荷載概率模型,結合臺風致自攻釘節點抗力退化模型,利用蒙特卡羅模擬,計算輕鋼屋面自攻釘節點在不同服役年限下的時變可靠度指標,即結構在(0,tL]時段內的可靠概率Rp(tL)和可靠指標分別為:

其中,若不考慮疲勞作用,則可進行初始抗力不變情況下的時變可靠度計算,此時:

通過對比抗力退化與不退化兩種情況,計算tL=5,10,15,…,50 年時的結構可靠度.繪制11 個城市的時變可靠度指標(如圖6 所示),其中目標可靠指標基于我國《建筑結構可靠度設計統一標準》(GB 50068—2018)[27]以及AISI-S100-07 規范[28],從而進行中美規范下可靠度指標的對比分析.

2.3 可靠度計算結果與分析

時變可靠度理論反映了臺風致疲勞損傷對自攻釘節點承載力退化的影響及臺風的隨機特性,使得輕鋼圍護系統的抗風設計更加科學合理,分析圖6可得:

1)考慮抗力退化后,輕鋼屋面自攻釘節點時變可靠度明顯下降.隨著服役時間的增加,臺風發生次數的增加,臺風荷載的影響也更加明顯,可靠度下降程度相對增大.在本文的臺風致抗力退化模型中,建筑物與臺風風眼之間的最近距離為25 km,作用時間為5 h.而在11 個城市風荷載模型中[23],臺風荷載的統計數據來自建筑物周圍500 km 范圍內.故在本文時變可靠度分析模型下,在建筑物500 km 半徑范圍內發生的臺風,其風眼均會行進至距離建筑物25 km處,此時代表了該區域為受臺風影響的極端危險區域,說明該假定涵蓋了臺風評估中的最不利情況.

2)本文分別參考中美兩國的目標可靠指標,其中我國《建筑結構可靠度設計統一標準》所提出的可靠度規范,對應的是整個圍護結構在三級延性破壞下的承載極限狀態,不應小于可靠指標2.7,但未單獨針對輕鋼屋面體系做出明確規定.而在美國AISI規范中該可靠指標規定為1.5,且針對輕鋼屋面系統,設定的可靠指標明顯低于我國規范指標.基于我國規范給出的目標可靠指標,針對我國11 個東南沿海城市,即使在不考慮抗力退化的情況下,寧波、福州、廈門、臺北等城市在5~10 年內將達到目標可靠指標,其余7 個城市則會在30 年內達到該指標;而在考慮臺風致自攻釘節點抗力退化的情況下,11 個城市在前5 年服役期內均會達到該指標.因此,對于由自攻釘連接的檁條-壓型鋼板系統而言,我國規范的可靠指標較為保守.

3)基于AISI 規范給出的檁條-壓型鋼板系統可靠指標,分析11 個沿海城市的時變可靠度可得:在不考慮抗力退化的情況下,各城市在50 年服役期內均不會達到目標可靠指標;而在考慮臺風致抗力退化的情況下,臺北將在服役10 年左右達到目標可靠指標,上海、寧波、溫州、福州、廈門、深圳、香港則會在服役20~30 年內達到該指標,而廣州、湛江、海口則有至少35~50 年的服役期.

3 結論

針對輕鋼屋面實現了基于概率的自攻釘節點抗風性能評估.對于抗力,首先利用模擬臺風方法,即基于相關風洞試驗的屋面時程風壓系數數據,結合臺風通過時的風速和風向信息,利用雨流分析法計算臺風過境期間的實際尺寸屋面風壓系數矩陣;結合相關自攻釘節點疲勞和靜力的統一試驗結果,計算臺風致作用引起的抗力退化情況,從而得到臺風致自攻釘節點抗力退化模型.對于荷載,參考臺風荷載統計數據模型.最后,通過蒙特卡羅模擬,結合算例完成了輕鋼屋面自攻釘節點的臺風致時變可靠度研究,結論如下:

1)在評估臺風作用對自攻釘節點的性能削弱方面,基于模擬臺風方法給出了模擬臺風的風速及風向的解析解,同時還給出了臺風致自攻釘節點抗力退化數學模型.

2)通過時變可靠度分析,對于輕鋼屋面系統,我國規范中針對輕鋼屋面自攻釘節點的目標可靠指標相對保守.

3)考慮臺風致抗力退化的情況,輕鋼屋面自攻釘節點時變可靠指標顯著下降.基于算例參數及AISI 規范,僅考慮臺風作用下,臺北的輕鋼屋面自攻釘節點服役年限為10 年左右;上海、寧波、溫州、福州、廈門、深圳、香港則會在20~30 年內到達服役年限;而廣州、湛江、海口則有至少35~50 年的服役年限.

猜你喜歡
風速
邯鄲市近46年風向風速特征分析
基于Kmeans-VMD-LSTM的短期風速預測
基于最優TS評分和頻率匹配的江蘇近海風速訂正
海洋通報(2020年5期)2021-01-14 09:26:54
基于時間相關性的風速威布爾分布優化方法
陜西黃土高原地區日極大風速的統計推算方法
陜西氣象(2020年2期)2020-06-08 00:54:38
基于GARCH的短時風速預測方法
快速評估風電場50年一遇最大風速的算法
風能(2016年11期)2016-03-04 05:24:00
考慮風切和塔影效應的風力機風速模型
電測與儀表(2015年8期)2015-04-09 11:50:06
GE在中國發布2.3-116低風速智能風機
考慮風速分布與日非平穩性的風速數據預處理方法研究
主站蜘蛛池模板: 国产成人精品一区二区秒拍1o| 很黄的网站在线观看| 欧美精品一区在线看| 亚洲国产理论片在线播放| 欧美激情,国产精品| 亚洲国产欧洲精品路线久久| 在线中文字幕日韩| 丰满人妻中出白浆| 日韩高清成人| 精品偷拍一区二区| 日韩欧美中文亚洲高清在线| 欧美中文一区| 毛片视频网址| 欧美a级在线| 在线a网站| 呦女亚洲一区精品| 丝袜美女被出水视频一区| www.亚洲天堂| av大片在线无码免费| 日本在线视频免费| 国产区人妖精品人妖精品视频| 成年午夜精品久久精品| 国产午夜看片| 国产三级a| 亚洲成人精品在线| 亚洲精品国产乱码不卡| 欧美无专区| 国内毛片视频| 青草免费在线观看| 欧美一区二区自偷自拍视频| 色噜噜狠狠色综合网图区| 国产福利在线免费| www.91中文字幕| 成人久久18免费网站| 国产精品美女免费视频大全| 成人亚洲视频| a毛片基地免费大全| 极品国产在线| 99re视频在线| 亚洲五月激情网| 亚洲无码电影| 免费无码网站| 992tv国产人成在线观看| 日本一区二区不卡视频| 亚洲精品麻豆| 亚洲av日韩av制服丝袜| 亚洲成肉网| 国产91丝袜在线播放动漫| 国产乱人伦AV在线A| 91精品情国产情侣高潮对白蜜| 亚洲天堂久久| 在线观看国产网址你懂的| 色婷婷丁香| 国产农村妇女精品一二区| 国产一区在线视频观看| 美女视频黄频a免费高清不卡| 亚洲天堂视频在线播放| 久久综合干| 青青草国产在线视频| 99久久精品免费视频| 亚洲精品午夜天堂网页| 超薄丝袜足j国产在线视频| AV在线天堂进入| 日韩在线视频网| 亚洲人成网址| 欧美激情视频二区| 国产黄网永久免费| 久久久久青草大香线综合精品| 九月婷婷亚洲综合在线| 波多野结衣在线se| 91香蕉视频下载网站| 制服丝袜亚洲| 日韩无码视频网站| 91久久偷偷做嫩草影院精品| 黄色成年视频| 亚洲欧洲日韩久久狠狠爱| 91无码人妻精品一区| 国产又大又粗又猛又爽的视频| 97se亚洲综合| 国产精品久久国产精麻豆99网站| 精品第一国产综合精品Aⅴ| 国产成人福利在线视老湿机|