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

基于改進模態推覆法的導管架平臺彈塑性抗震性能

2021-04-07 08:58:04劉紅兵孫麗萍艾尚茂閆發鎖陳國明
上海交通大學學報 2021年3期
關鍵詞:模態結構能力

劉紅兵,孫麗萍,艾尚茂,閆發鎖,陳國明

(1.哈爾濱工程大學 船舶工程學院,哈爾濱 150001;2.中國石油大學(華東) 海洋油氣裝備與安全技術研究中心,山東 青島 266580)

海底地震作為一種突發性海洋環境災害,具有持續時間短、震級大、破壞性強等特點[1],如1969年渤海灣7.4級地震、2004年印度洋9.4級地震以及2011年日本Tohok-Oki 9.0級地震等,嚴重威脅海洋導管架平臺的服役安全[2].現階段關于導管架平臺的抗震設防尚缺乏統一的專業規范,主要參考美國石油學會(API)以及建筑結構工程等領域相關的設計規范.其中,API規范要求設計地震下,海洋導管架平臺結構需滿足強度要求,即結構和基礎均無破壞;罕遇地震下平臺需滿足韌性要求,即允許結構和基礎有一些損傷,但結構不能倒塌[3].而建筑工程規范則要求結構滿足“兩階段,三水準”原則,其中“兩階段”分別指彈性抗震和彈塑性抗震階段,“三水準”則指“小震不壞、中震可修、大震不倒”的定性準則[4].這些規范準則對于海洋導管架平臺的線彈性抗震設防具有較好的設計指導作用,但對于強震下平臺結構的非線性彈塑性抗震性能以及平臺結構地震失效機制等問題尚缺乏深入研究.而對于海洋導管架平臺,通常結構冗余度較大,傳統的線彈性抗震評估方法無法準確評估這種大冗余度結構的彈塑性抗震性能,因而準確地評估強震下平臺結構的彈塑性地震響應,對于保證海洋導管架平臺的抗震設計及其安全運行具有重要意義.

目前,海洋導管架平臺的抗震設計分析方法主要有靜力法、反應譜法、動力時程法、彈塑性推覆法等.靜力分析法通常將地震作用等效為一種側向慣性載荷作用,并采用地震系數修正地震過程中的動態效應,但由于平臺結構和地震動之間耦合作用的復雜性,導致地震系數難以準確確定,從而容易產生較大的誤差[5-6].反應譜法又稱為等效靜力分析法,通常通過振型分解將平臺結構簡化為單質點體系,然后采用單質點體系最大地震響應描述結構地震響應,該方法可以較好地獲得簡單結構的彈性地震響應,但卻無法考慮強地震下平臺結構的非線性彈塑性地震響應[7-8].動力時程法通過對平臺結構輸入時變地震動并進行動力積分,獲得結構彈塑性地震響應,該方法可較好地反映結構動態地震響應,但由于地震動的隨機性和差異性,即使具有相同控制參數的不同輸入地震動時程,也可能導致結構地震響應相差巨大[9-10].彈塑性推覆法(能力譜法)通過將設計反應譜引入結構抗震性能評估,可較好地識別結構彈塑性抗震性能以及抗震薄弱環節和破壞機制,已在結構工程領域獲得廣泛認可,并被美國、日本、韓國等規范納入,但該方法主要基于結構基本振型和形狀向量保持不變的假設,忽略了高階振型以及塑性抗震階段結構振型向量發生變化的影響,從而在平臺結構實際彈塑性抗震設計中具有一定的局限性[11-13].為此本文以渤海某導管架平臺為研究對象,建立一種適于海洋導管架平臺彈塑性抗震性能評估的改進模態推覆(IMPA)法,對比分析不同分析方法(能力譜、動力時程和改進模態推覆法)下,平臺結構彈塑性抗震性能的差異性,探討高階振型和振型向量對平臺結構彈塑性抗震性能的影響,識別平臺結構抗震薄弱環節,相關結論和研究成果可為海洋石油導管架的平臺抗震設計提供理論依據和工程參考.

1 模態推覆法(MPA)理論基礎

海洋導管架平臺服役環境惡劣,地震作用下除了受到海床震動對結構產生的慣性載荷外,還將受到地震引發海水對結構產生的附加載荷作用.由于海水附加質量產生的慣性載荷對平臺結構的地震響應影響較小,故可將其忽略.根據非線性結構動力學原理可建立地震載荷作用下,海洋平臺結構彈塑性抗震運動的微分方程為

(1)

假設平臺結構地震響應主要由結構前n階振型控制,則根據振型分解的基本思想,對式(1)進行解耦即可得到平臺結構運動的微分方程為

(2)

由式(2)可知,地震載荷作用下的平臺結構恢復力函數依賴于各模態坐標之間的耦合作用,但由于此時的平臺結構響應仍以n階模態為主導,可忽略各模態坐標之間的耦合作用[14-15],所以式(2)可采用類似于彈性體系的強迫解耦算法進行解耦.假定平臺結構地震作用下n階模態單自由度體系(SDOF)的位移響應為

(3)

式中:Dn為平臺結構n階模態單自由度體系位移響應向量;I為單位矩陣;qn為n階振型模態位移.將式(3)代入式(2),整理可得平臺結構單自由度體系n階振型運動方程為

(4)

(5)

式中:Vb,n為平臺結構n階模態剪力.平臺結構等效單自由度體系屈服剪力Vb,ny和屈服位移Dny為

(6)

式中:Fs,ny為n階模態屈服剪力所對應的恢復力;ur,ny為n階模態結構頂點屈服位移;φr,n為n階模態頂點振型值.根據式(6)即可獲得平臺結構n階模態單自由度體系的等效能力曲線,如圖1所示.

圖1 平臺結構等效單自由度體系能力曲線Fig.1 Capacity curve of equivalent single-degree-of-freedom system of platform

根據平臺結構n階模態等效單自由度體系,即可求得n階模態的有效周期為

(7)

進一步將各階振型下平臺結構的地震響應進行疊加組合,即可得到地震載荷作用下的平臺結構響應,表達式為

(8)

2 改進模態推覆法

2.1 改進模態推覆法的基本原理及分析步驟

模態推覆法通過振型組合的方式對平臺結構進行抗震性能評估,該方法可以考慮高階振型對平臺結構抗震性能的影響,但是在各階振型分析過程中采用統一的位移形狀向量,忽略了平臺結構由彈性抗震階段轉變為彈塑性抗震階段形狀向量發生變化的影響,從而導致模態推覆法的評估結果與平臺結構實際抗震性能不相符.對于導管架平臺結構,在彈性抗震階段時,結構形狀向量保持不變,但當平臺結構發生屈服進入彈塑性抗震階段時,其結構剛度和位移形狀向量均發生改變,理想的計算方法是將上一步地震載荷作用下的結構位移向量作為下一步地震載荷作用下的形狀向量,這樣形成的彈塑性抗震分析方法與平臺結構地震載荷作用下的實際受力情況最吻合,但同時會導致計算非常復雜且計算量巨大.根據MPA法獲得的平臺結構抗震能力曲線可知(見圖1),平臺結構進入塑性抗震階段后,抗震能力曲線的切線剛度雖然一直在變化,但其變化量卻較小,故可近似采用塑性階段的等效折線斜率代替平臺結構的塑性剛度.因而IMPA法的基本思想為在彈性抗震階段,采用平臺結構彈性自振形狀向量進行抗震分析;在彈塑性抗震階段,則采用結構塑性階段的等效折線斜率所對應的平臺結構位移向量進行抗震分析.

IMPA法的具體流程如圖2所示,主要分析步驟如下:① 建立平臺結構的有限元模型,計算平臺結構各階的自振頻率ωn和振型φn;② 計算平臺結構各階振型質量參與系數,確定有效參振振型;③ 以平臺結構的n階振型φn作為初始形狀向量,對平臺結構施加側向載荷Fn1=mφn,開展推覆分析,獲得平臺結構n階振型的推覆能力曲線;④ 根據圖1,將平臺結構的n階振型推覆能力曲線簡化為等效雙折線模型(推覆能力曲線和雙折線與x軸圍成的面積相等),并進一步建立平臺結構線彈性n階振型等效單自由度體系能力曲線((Fs,n/Ln)-Dn關系曲線);⑤ 若平臺結構發生屈服,進入塑性抗震階段,可在平臺結構推覆能力曲線上選擇一點,使得該點切線的斜率和簡化后的等效折線斜率相等,然后將該點處所對應的位移向量作為平臺結構塑性抗震階段的形狀向量φn1,對平臺結構施加側向載荷Fn2=mφn1,進行第2次推覆分析,獲得平臺結構塑性抗震階段推覆能力曲線;⑥ 組合平臺結構彈性和塑性抗震階段推覆能力曲線,并建立平臺結構彈塑性n階振型等效單自由度體系能力曲線;⑦ 組合各階振型下平臺結構的地震響應,計算平臺結構的總地震響應,并評估平臺結構的抗震性能.

圖2 IMPA法的基本流程圖Fig.2 Basic flowchart of IMPA method

2.2 改進模態推覆法的有效振型

地震載荷作用下,海洋平臺結構各階振型對平臺結構總體反應的貢獻大小差異較大,如何選取合適的參振振型對于提高平臺結構抗震性能評估的準確度和減小計算量具有重要的意義.目前,關于海洋導管架平臺結構抗震設計過程中有效參振振型的選取沒有具體規范,主要參考建筑結構抗震設計規范,即不考慮扭轉時,若結構基本周期小于1.5 s,可取2~3個振型,若基本周期大于1.5 s或高寬比大于5時,振型數可適當增加;考慮扭轉時,可取前9~15個振型[4],然而海洋導管架平臺結構抗震性能與陸地建筑結構的差異較大,參考陸地建筑規范選取有效參振振型可能會導致較大的誤差.因而,本文采用振型質量參與系數法確定有效參振振型,即保證所有參振振型的參與質量不小于平臺結構總質量的90%[16].

平臺結構n階振型廣義質量為

(9)

(10)

根據式(10),定義平臺結構n階振型參與質量為平臺結構n階振型廣義質量與n階振型參與系數平方的積,具體表達式為

(11)

式中:me,n為平臺結構n階振型的參與質量.平臺結構n階振型質量參與系數an為

(12)

將式(12)按照x、y和θ方向分解,即可得到平臺結構不同自由度下各階振型質量參與系數,表達式為

(13)

式中:ME,nx和ME,ny分別為平臺結構n階振型參與質量在x方向和y方向上的分量;ME,nθ為n階振型參與質量θ方向的轉動慣量;Jn為平臺結構n階振型的轉動慣量.

3 算例分析

3.1 平臺結構有效參振振型

以文獻[1]中渤海某導管架平臺為研究對象進行分析.該平臺為4樁腿導管架平臺,主要由樁腿、導管架和上部組塊等構成,其中平臺底部標高為 -15.0 m,上部組塊頂部標高為20.0 m.為充分了解該平臺結構自振特性及其扭轉效應,系統地分析了平臺結構前20階自振頻率和模態振型特征,其中前6階自振頻率和模態振型如圖3所示,其中f為頻率.由圖3可知,平臺結構基本自振頻率為1.84 Hz,1階振型為平臺沿y軸平動,2階振型為平臺沿x軸平動,3~5階振型為平臺繞3個坐標軸轉動,6階振型為平臺沿z軸平動.

圖3 平臺結構前6階模態振型圖Fig.3 First 6 modal shapes of platform

根據式(13)分別計算平臺結構6個自由度方向各階振型質量參與系數,前10階振型質量參與系數如表1所示,其中anz為n階振型沿z方向的參與系數;ar,nx、ar,ny和ar,nz分別為n階振型沿x、y、z軸旋轉方向參與系數.由表1可知,平臺結構x方向振型質量參與系數在前5階累計達到98.89%,y方向振型質量參與系數在前4階累計達到98.18%,z方向振型質量參與系數在前6階累計達到95.78%,x、y、z軸旋轉方向振型質量參與系數在前9階累計均達到90%以上,因而采用改進模態推覆法對該平臺進行彈塑性抗震性能評估時,可選取前9階或9階以上振型作為有效參與振型.同時由于平臺結構地震響應主要集中于側向x和y方向,其余方向地震響應相對較小,為減小計算量,僅考慮平臺結構x和y兩個方向的地震響應.進一步分析x和y方向振型的質量參與系數可發現,平臺結構x方向振型質量參與系數主要集中于第2、5、9階振型(累計振型質量參與系數為99.20%),y方向振型質量參與系數主要集中于1、4、8階振型(累計振型質量參與系數為99.22%),因而可忽略其余振型對平臺結構x和y方向地震響應的影響,即計算以上6階振型對應方向的地震響應,然后通過振型組合方法獲得平臺結構整體的地震響應.

表1 平臺結構各階模態振型質量參與系數表Tab.1 Mode mass participation factors of platform

3.2 設防烈度8度地震響應

選取設防烈度為8度的地震(50年超越概率為10%,峰值加速度0.2g,g為重力加速度),分別針對上述有效參振振型進行模態推覆分析,獲得各階模態振型下平臺結構能力曲線和等效單自由度雙折線.通過對平臺結構各階模態振型下頂點位移和基底剪力與等效屈服點(二折線交點)所對應的頂點位移和基底剪力之間的比較,判斷不同模態振型下平臺結構是否發生屈服,進入塑性抗震階段,結果如表2所示.其中:ur為平臺結構頂點位移;Vb為基底剪力.

表2 設防烈度8度地震下平臺結構各階模態振型響應Tab.2 Modal responses of platform at 8 degree seismic fortification intensity

由表2可以看出,6個模態振型下的平臺結構頂點位移和基底剪力均小于對應的等效屈服點頂點位移和基底應力響應,平臺結構沒有發生屈服,因而在后續模態推覆分析過程中可采用統一的形狀向量.

上述各階模態振型下的平臺結構x和y方向位移響應圖如圖4所示,其中:h為平臺結構垂向高度;u為側向位移.由圖4可看出,平臺結構2階模態振型下的x方向位移和1階振型下的y方向位移顯著大于其余各階模態振型下的位移響應,平臺結構地震響應以第1和2階模態振型為主.2階模態振型下平臺結構x方向位移在導管架和上部組塊區間段隨著高度的增大而逐漸增大,而在導管架頂層和下甲板之間的斜撐區段,位移值隨著高度的增大而減小,這主要是導管架和上部組塊間質量和剛度矩陣的差異性所致;5階模態振型下,x方向位移值為負,表明此時結構沿著x軸反方向振動;9階模態振型下,x方向位移非常小,說明結構x方向的振動以第2和5階模態振型為主,如圖4(a)所示.與x方向位移類似,1階模態振型下,平臺結構y方向位移隨著高度的增大而逐漸增大,且上部組塊區間段的增大幅值小于導管架區間段;第4和8階模態振型下,y方向的位移值較小,結構y方向振動以第1階模態振型為主,如圖4(b)所示.

圖4 各階模態振型下平臺結構x和y方向位移圖Fig.4 Displacement curves in x and y directions of platform in different modal shapes

圖5 P1和P2的地震加速度譜和時程Fig.5 Acceleration spectrums and time histories of P1 and P2 seismic

分別采用上述3種不同方法獲得的平臺結構地震位移和位移角(η)響應如圖6所示.由圖6可知,平臺結構的最大位移均位于平臺結構頂層,最大位移角均位于導管架頂端位置處.對比MPA法和能力譜法的位移響應可看出,在導管架下部區域,MPA法位移響應小于能力譜法;在導管架中間區段,MPA法和能力譜法的位移響應相差不大;在平臺上部組塊區域,MPA法的位移響應則大于能力譜法,其中MPA法的位移最大值為4.49 cm,能力譜法的位移最大值為4.21 cm,前者比后者增大了約6.65%,說明高階振型對平臺結構底部和上部組塊結構的地震位移響應具有較大的影響,不可忽略.對比P1和P2地震時程分析法的位移響應可看出,不同地震時程作用下,平臺結構的位移響應相差較大,其中P2地震位移響應在導管架區段和MPA法以及能力譜法的位移響應較為接近,在上部組塊區段卻小于MPA法和能力譜法的位移響應,而P1地震位移響應則遠大于P2、MPA法以及能力譜法的位移響應.P1地震最大位移響應為6.27 cm,P2地震最大位移響應為3.97 cm,P1最大位移比P2最大位移約增大了58.12%,表明具有相同峰值加速度的不同輸入地震時程,平臺結構的地震位移響應表現出較大的離散性和差異性,如圖6(a)所示.

圖6 設防8度地震下的平臺結構地震響應Fig.6 Seismic responses of platform at 8 degree seismic fortification intensity

對比3種方法下的平臺結構位移角響應可看出,平臺結構的最大位移角均位于導管架頂端位置處,說明該位置處為平臺結構的抗震薄弱環節,進行彈塑性抗震設計時需重點關注;MPA法和能力譜法下,平臺結構的位移角整體形狀類似,且MPA法的位移角整體略大于能力譜法的位移角,其中MPA法最大位移角為0.24%,能力譜法最大位移角為0.22%,MPA法最大位移角約增大了9.31%,進一步表明高階振型對平臺結構抗震性能的影響;P1和P2地震時程作用下,平臺結構位移角隨著高度的變化表現出更為復雜的變化趨勢,其中P1地震位移角響應整體上大于MPA法和能力譜法位移角響應,P2地震位移角響應在導管架下部區域略大于MPA法和能力譜法響應,在導管架中部區域小于MPA法和能力譜法響應,在上部組塊區域則又大于MPA法和能力譜法響應,如圖6(b)所示.

3.3 罕遇烈度8度地震響應

由3.2分析可知,設防烈度8度地震作用下,平臺結構沒有發生屈服,處于彈性抗震階段,為進一步研究平臺結構彈塑性抗震性能,選取罕遇烈度8度的地震(50年超越概率為2%~3%,峰值加速度為0.4g,水平地震最大影響系數為1.2)進行抗震性能分析.根據前文所述的平臺結構地震載荷作用下有效參振振型,分別開展模態推覆分析,獲得各階有效參振振型下平臺結構能力曲線和等效單自由度雙折線,并通過與等效屈服點的比較判斷平臺結構各階振型下是否發生屈服,進入塑性抗震階段,結果如表3所示.

由表3可看出,x方向的第2階模態振型和y方向的第1階模態振型下平臺結構發生屈服,進入塑性抗震階段,其抗震形狀向量發生改變,故需采用IMPA法進行彈塑性抗震性能分析.

表3 罕遇烈度8度地震下平臺結構各階模態振型響應Tab.3 Modal responses of platform at 8 degree seismic rare intensity

分別采用IMPA法和MPA法,獲得的平臺結構在第1階振型下y方向位移和第2階模態振型下x方向位移的比較如圖7所示.由圖7可看出,兩種方法下的平臺結構x和y方向位移形狀相似,其中在導管架部分,IMPA法和MPA法的位移響應相差不大.但從導管架頂部到上部組塊區域,IMPA法的位移結果顯著大于MPA法的位移結果.這主要是由于導管架的頂端位置處為抗震薄弱環節,結構在該位置處發生屈服進入塑性抗震階段,使得導管架頂端上部組塊位移顯著增大所致.通過改變形狀向量進行二次推覆分析的IMPA法,可較好地體現平臺結構的彈塑性抗震性能特征,與平臺結構的實際抗震狀態較為吻合.

圖7 IMPA法和MPA法分析結果比較圖Fig.7 Analysis results of IMPA and MPA methods

采用SRSS振型組合法組合獲得的平臺結構整體地震響應,并分別同MPA法、能力譜法和時程分析法進行比較,其中時程分析法輸入地震波仍采用調幅至峰值加速度為0.4g的P1和P2地震波,結果如圖8所示.由圖8(a)可看出,4種方法下平臺結構地震位移響應隨高度變化趨勢基本類似,位移最大值均發生在平臺頂點處;對比能力譜、MPA法和IMPA法平臺結構位移可知,在導管架部分三者位移響應基本相等,但在上部組塊區域,3種方法的位移響應表現出明顯的差異,其中IMPA法位移響應最大、MPA法次之,能力譜法位移響應最小,這主要是由于隨著高度增大,高階振型影響逐漸增大以及平臺結構塑性振動特性變化越來越明顯所致.P1和P2地震波時程分析法結果相差較大,其中P1地震波作用下平臺結構的位移響應均大于能力譜、MPA法和IMPA法位移響應,而P2地震波作用下平臺結構位移卻均小于其余3種方法的位移響應,進一步體現了不同的輸入地震動特性對于平臺結構抗震性能評估產生的巨大影響.對比4種方法平臺結構的最大位移可知,IMPA法平臺結構頂點的最大位移值為0.121 m,比MPA法的最大位移增大了7.64%,比能力譜法的最大位移增大了11.42%,比P1地震時程法的最大位移減小了8.01%,比P2地震時程法的最大位移增大了51.01%,充分說明對于罕遇強震,在無法獲取準確的強震加速度時程情況下,采用能力譜法和MPA法評估平臺結構的彈塑性抗震性能均會產生較大誤差,建議采用 IMPA 法.

圖8 罕遇8度地震下平臺結構的地震響應Fig.8 Seismic responses of platform at 8 degree seismic rare intensity

由圖8(b)可知,罕遇烈度8度地震作用下,平臺結構位移角表現出更為復雜的變化特征,在導管架下部區域,P1地震時程法位移角最大,P2地震時程法位移角最小,IMPA法、MPA法和能力譜法基本相同,在導管架上部區域則變為IMPA法位移角最大,MPA法、P1地震時程法、能力譜法和P2地震時程法依次次之,在平臺上部組塊區域,位移角表現出明顯的波動性,這主要是由于進入塑性抗震階段后,平臺結構的振動特性變得更為復雜所致.與設防烈度8度地震類似,平臺結構地震最大位移角也均位于導管架頂端位置處,其中IMPA法最大位移角值為0.59%,MPA法為0.54%,能力譜法為0.45%,P1地震時程法為0.50%,P2地震時程法0.30%,均小于1%,但此時平臺結構已經發生屈服進入塑性抗震階段,與建筑結構工程中1%判據具有一定的差別,這主要是由于平臺結構冗余度大,具有更高的側向位移剛度,所以其側向位移角相對較小.

4 結論

(1) 高階模態振型對于平臺結構彈塑性抗震性能影響較大,地震載荷作用下,平臺結構x方向振型質量參與系數主要集中于第2、5、9階振型,累計振型質量參與系數達99.20%,y方向振型質量參與系數主要集中于第1、4、8階振型,累計振型質量參與系數達99.22%,因而為保證平臺結構彈塑性抗震性能精度,需考慮平臺結構前9階或9階以上的模態振型影響.

(2) 設防烈度8度地震下,平臺結構處于彈性抗震階段,最大位移和最大位移角分別位于平臺結構頂點和導管架頂部,表明導管架頂部為平臺結構抗震薄弱環節,需重點關注;不同分析方法下平臺結構的最大地震響應相差較大,MPA法最大位移和最大位移角分別為4.49 cm和0.24%,比能力譜法最大位移和最大位移角分別增大了6.65%和9.31%,表明高階振型對平臺結構彈性抗震性能的影響.

(3) 罕遇烈度8度地震下,平臺結構發生屈服,進入塑性抗震階段,結構最大位移和位移角仍分別位于平臺頂點和導管架頂部,其中IMPA法平臺頂點最大位移為12 cm,比MPA法和能力譜法分別增大了7.64%和11.42%,IMPA法導管架頂部最大位移角為0.59%,比MPA法和能力譜法分別增大了9.26%和31.11%,進一步表明高階振型和位移形狀向量對平臺結構塑性抗震性能的巨大影響.

(4) 具有相同峰值加速度的不同地震時程下平臺結構地震響應表現出明顯的離散性和差異性,當地震峰值加速度為0.2g時,P1和P2地震時程作用下平臺結構頂點最大位移相差約58.12%;當地震峰值加速度為0.4g時,兩者相差約60.47%.因而在無法準確獲取地震時程情況下,時程分析法、能力譜法和MPA法均會產生較大誤差,建議采用IMPA法評估平臺結構彈塑性抗震性能.

猜你喜歡
模態結構能力
消防安全四個能力
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
大興學習之風 提升履職能力
人大建設(2018年6期)2018-08-16 07:23:10
你的換位思考能力如何
論《日出》的結構
國內多模態教學研究回顧與展望
抄能力
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
基于HHT和Prony算法的電力系統低頻振蕩模態識別
主站蜘蛛池模板: 欧美视频在线不卡| 91福利在线看| 日韩国产黄色网站| 中文字幕自拍偷拍| 91视频青青草| 亚洲精品无码日韩国产不卡| 美女内射视频WWW网站午夜| 永久免费AⅤ无码网站在线观看| 亚洲成人高清无码| 久久国产V一级毛多内射| 中文字幕va| 不卡国产视频第一页| 亚洲国产亚综合在线区| 自拍亚洲欧美精品| 乱人伦中文视频在线观看免费| 国产精品观看视频免费完整版| 亚洲欧美在线综合一区二区三区| 中国一级特黄视频| 亚洲精品男人天堂| 波多野结衣一区二区三区四区| 亚洲人网站| 伊人久久大香线蕉aⅴ色| 欧美精品三级在线| 草逼视频国产| 国产在线第二页| 欧美国产日韩另类| 日韩国产欧美精品在线| 亚洲中文字幕手机在线第一页| av尤物免费在线观看| 日韩精品毛片人妻AV不卡| 国产微拍一区| 国产成年无码AⅤ片在线| 国产经典免费播放视频| 亚洲精品久综合蜜| 欧美国产在线看| 91精品国产丝袜| 欧美亚洲一区二区三区导航| 国产黄色视频综合| 成人自拍视频在线观看| 亚洲成a人片在线观看88| 亚洲一级无毛片无码在线免费视频 | 天堂岛国av无码免费无禁网站| 美女扒开下面流白浆在线试听| 91视频首页| 亚洲天堂在线免费| 99热国产这里只有精品9九 | 亚洲综合国产一区二区三区| 少妇人妻无码首页| 国产亚洲精品91| 东京热高清无码精品| 亚洲中文字幕在线观看| 99热国产这里只有精品无卡顿" | 亚洲精品国产日韩无码AV永久免费网| 亚洲一级毛片免费看| 男人天堂亚洲天堂| 丝袜亚洲综合| 日韩在线欧美在线| 亚洲AV无码乱码在线观看代蜜桃| 亚洲AⅤ永久无码精品毛片| 美女被操黄色视频网站| 亚洲无码电影| 青青草a国产免费观看| 欧美黄网站免费观看| 成年人视频一区二区| 国产综合精品日本亚洲777| av在线手机播放| 一本久道热中字伊人| 久久亚洲美女精品国产精品| 成人国产小视频| 国产日韩精品欧美一区灰| 日韩小视频网站hq| 国产午夜福利亚洲第一| 在线精品自拍| 成人亚洲国产| 亚洲欧美日韩中文字幕在线一区| 天堂va亚洲va欧美va国产| 国产午夜精品鲁丝片| 国产最新无码专区在线| 国内视频精品| 最新国语自产精品视频在| 欧美色99| 91视频精品|