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

大跨度懸索橋施工階段三維顫振精細化分析

2020-04-09 04:08:06
浙江工業大學學報 2020年2期
關鍵詞:風速結構影響

(浙江工業大學 建筑工程學院,浙江 杭州 310023)

懸索橋結構受力性能好,跨越能力強,是千米級主跨橋梁工程的首要選擇。隨著橋梁工程建設由跨越大江大河向近海連島工程及跨越海峽和海洋等更廣闊的水域發展,懸索橋的跨徑將進一步增大,潛在需求在2 000~5 000 m之內[1]。懸索橋結構跨度大,剛度小,風作用下的結構穩定性(主要指顫振)已成為控制懸索橋設計和施工的重要因素。較之于成橋狀態,施工狀態的懸索橋尚未形成最終的結構體系,結構邊界約束更少,結構的整體剛度尤其是扭轉剛度明顯降低,導致結構抗風穩定性更差[2-3]。雖然施工階段結構顫振檢驗風速可以降低,但是一般情況下懸索橋施工期難以避開強風天氣。隨著懸索橋跨徑的進一步增大和施工狀態懸索橋剛度的顯著降低,以下兩個因素對施工期懸索橋顫振穩定性的影響將更加突出:1)靜風作用下結構大變形導致的結構剛度和施加在結構上以結構變形為函數的風荷載的非線性變化及其三維效應,簡稱靜風效應;2)風速空間分布的非一致性。實測資料表明風速沿著豎直高度和水平方向是變化的,但已有大跨度橋梁顫振分析中通常將橋址區域內的風速考慮為均勻分布。懸索橋的主纜矢高和橋塔高度都比較大,橋面主梁采用豎曲線布置,依據風的空間分布特性,橋面主梁、主纜和吊桿以及橋塔等構件上風速的差異性將更加明顯,形成風速的空間非均勻分布。已有分析表明這些因素對成橋狀態懸索橋的顫振存在著重要影響[4],但對施工狀態大跨度懸索橋顫振穩定性的影響如何,則需要進一步分析和明確。

大跨度懸索橋施工階段的顫振穩定性作為一個重要的工程問題已經引起了國內外學者的關注與研究。Tanaka等[5]、Ge等[6]和Diego等[7]結合Hoga Kusten橋的風洞試驗和理論分析研究了不同施工方法對施工期懸索橋顫振穩定性的影響,發現主梁非對稱架設比對稱架設具有更好的穩定性能。劉高等[8]采用多模態顫振分析方法對某懸索橋設計方案架設階段顫振臨界風速的發展趨勢進行分析,并針對其最危險的施工階段分析了主纜與主梁之間設置交叉索對提高顫振穩定性的效果。劉竹釗等[9]用半解析半試驗的方法,分析了主梁非對稱拼裝和主纜間設置臨時水平交叉等措施對提高懸索橋施工過程顫振穩定性的效果。魏志剛等[10]結合西堠門大橋施工過程氣彈模型風洞試驗發現非對稱施工在某些施工階段能明顯提高顫振穩定性,但卻在另一些施工階段中表現出比對稱架設更差的顫振穩定性。李永樂等[11]采用狀態空間法,對主跨3 500 m碳纖維主纜懸索橋施工階段的顫振穩定性進行了分析,討論了加勁梁剛度折減和阻尼比等參數對施工階段顫振穩定性的影響,并研究了交叉索對提高懸索橋施工階段顫振穩定性的效果。Chobsilprakob等[12]將階躍函數應用于丹麥大帶東橋成橋和架設階段的顫振時域分析,并探討了結構幾何非線性的影響。Han等[13]采用ANSYS有限元分析軟件對中都長江大橋進行了施工過程的顫振穩定性數值分析,研究了初始風攻角、主梁安裝順序、加勁剛度折減對顫振穩定性的影響。張新軍等[14]以世界首座跨度超千米的三塔懸索橋-泰州長江公路大橋為工程背景,分別模擬主梁從主跨跨中向兩側橋塔、從兩側橋塔向主跨跨中以及從兩側橋塔和主跨跨中同時向主跨四分點處對稱拼裝的施工順序,在風速空間均布分布情況下采用三維非線性空氣靜力和動力穩定性分析方法,分析主梁拼裝過程結構的空氣靜力和動力穩定性的演變規律,并從抗風穩定性角度提出三塔懸索橋適宜的主梁拼裝施工順序。可以看出,上述研究主要集中在探明施工過程懸索橋顫振穩定性的變化趨勢和改善措施等方面,顫振分析采用基于平板顫振理論的估算公式或基于狀態空間法的線性多模態方法,并假定風速空間一致分布,因此都沒有綜合考慮靜風效應和橋址處風速非均勻分布等因素對施工過程懸索橋顫振穩定性的影響。為此,考慮結構非線性、靜風效應和風速空間非均勻分布等因素,筆者建立了精細化的大跨度橋梁三維非線性顫振分析方法,并編制了計算程序。采用該程序,結合潤揚長江大橋南汊懸索橋,模擬兩種主梁架設順序,分析大跨度懸索橋施工全過程顫振穩定性的演變規律,探索適宜的橋面主梁架設方法,同時探明靜風效應和風速空間非均勻分布等因素對施工狀態懸索橋顫振穩定性的影響,為確保大跨度懸索橋安全施工提供理論依據。

1 三維非線性精細化顫振分析方法及計算程序

1.1 風速空間非均勻分布模型

橋址空間范圍內平均風速可表示為

U=μU0

(1)

式中:U0為參考點處的風速值,一般可以取為中跨主梁跨中處的風速值;μ為風速空間分布系數,依據風場的實測資料[4]可近似地表示為

μ=μH·μV

(2)

式中:μH為風速水平變化系數;μV為風速豎向變化系數;L為橋梁總長;L1為風場分布寬度;e為風場分布非對稱性系數,0≤e≤1,e=0表示風速相對于中跨跨中水平對稱分布;x為風速計算點至橋跨左端的距離;y為風速計算點處的離地高度;y0為參考點處的離地高度;α為地面粗糙度指數。

1.2 風荷載計算模型

風對橋梁結構的作用可以分解為平均風和脈動風作用。在平均風作用下,橋面主梁單位長度上受到的靜風荷載可以分解為如圖1所示的靜力三分力即順風向阻力Fz、橫風向升力Fy和升力矩Mx。由于橋面主梁在靜風作用下產生的變形會反過來改變來流風與橋面主梁間的相對攻角,使得作用在其上的靜力風荷載的非線性變化及三維效應。考慮風速空間分布后單位長度橋面主梁所受到的靜風荷載可表達為

(3)

式中:ρ為空氣密度;D和B分別為橋面主梁在豎直和水平方向的投影高度和寬度;Cz(αe),Cy(αe),CM(αe)分別為體軸下的靜力三分力系數,可由節段模型風洞試驗測得;αe為有效風攻角,等于來流風初始攻角θ0與靜風作用產生的主梁扭轉角θ之和。

圖1 作用在主梁上的靜風荷載

在脈動風作用下,橋面主梁將在靜風作用下的平衡狀態上振動,作用在單位長度橋梁斷面上的自激氣動力可以采用Scanlan提出的用18 個顫振導數表達的自激氣動力計算公式。采用有限元方法將主梁單元受到的均布自激氣動力轉化為單元兩端結點的等效集中荷載,則作用于主梁單元e兩端結點的等效自激氣動力可表示為

(4)

(5)

1.3 三維非線性顫振精細化分析程序

基于線性顫振分析理論,考慮風速空間非均勻分布、靜風作用引起的結構動力特性和作用在結構上的自激氣動力隨結構變形的非線性變化效應,筆者建立了大跨度橋梁三維非線性顫振精細化有限元分析方法,并編制了計算分析程序Nflutter,程序計算流程如圖2所示。

圖2 三維非線性顫振精細化分析流程

程序采用外部風速和內部振動頻率雙參數搜索方法尋求結構系統的顫振臨界狀態,具體分析思路如下:

1)外部風速循環的目的是獲得當前計算風速下結構的靜風平衡狀態,以此作為結構在自激氣動力作用下振動的基準態。靜風作用下,結構的幾何和內力狀態將發生變化,導致結構剛度以及結構單元與來流風之間形成的有效攻角αe發生變化,并將最終影響到結構系統的顫振性能。

2)由于單元自激氣動力是有效風攻角和折算頻率的函數,基于試驗測得的氣動導數,將每一級風速下由步驟1)得到的各構件單元的有效風攻角進行插值計算得到相應的氣動導數,并重新形成式(5)所示的各單元氣動剛度和氣動阻尼矩陣,以此考慮靜風效應對單元自激氣動力的影響。風速的空間非均勻分布則同時考慮風速沿水平向和豎向的變化,結構單元的風速值則取單元兩端結點風速的平均值。

3)由于結構體系的顫振一般發生在以主梁振動主導的低階頻率特別是扭轉頻率上,為提高計算效率,計算中主要跟蹤主梁振動為主的低階復特征值隨風速的變化。因此,內部振動頻率的搜索主要在以主梁振動為主的低階模態頻率范圍內進行。在外部風速循環下,隨著風速的變化,結構的靜風平衡狀態隨之改變,導致結構的剛度及其動力特性發生變化。為加快頻率迭代收斂,在每一級搜索風速下都進行靜風平衡狀態下的結構動力特性分析,由程序自動獲取變化后的以主梁振動為主的頻率和模態。

4)結構系統的顫振特征方程是一個2m階的廣義特征方程(m為參與顫振分析的結構模態數),由于結構系統的氣動阻尼和氣動剛度矩陣為不對稱矩陣,考慮結構本身的阻尼和剛度矩陣后得到的結構總體阻尼和剛度矩陣均為不對稱矩陣。采用多模態分析方法后,特征值方程的階數大大減少,采用Hessenberg矩陣變換和QR變換相結合的方法求解全部特征值及其特征向量。

2 橋梁及主梁架設方案簡介

潤揚長江大橋南汊懸索橋是一座單跨懸索橋[14],中跨1 490 m,兩側邊跨各470 m,見圖3。中跨主纜矢跨比為1/10,主纜橫橋向中心距為34.3 m;吊桿縱橋向間距為16.1 m,共設91 對吊桿;橋面主梁采用全焊扁平流線型鋼箱梁,總寬38.7 m,梁高3 m;橋塔為雙柱三橫梁混凝土門式框架結構,塔高約210 m。

圖3 潤揚長江大橋南汊懸索橋總體布置圖

如圖4所示,懸索橋的主梁架設方法按其推進方式主要有兩種:一種是從跨中位置開始向兩側橋塔對稱拼裝,在施工過程中只有一個架設梁段,梁段兩端是自由的;另一種是從兩側橋塔位置開始向跨中對稱拼裝,在施工過程中有兩個獨立的架設梁段,梁段的一端是自由的,另一端擱置在橋塔上。筆者分別采用這兩種主梁架設方法對潤揚長江大橋施工過程的顫振穩定性進行分析,架設方案一和架設方案二施工階段劃分和梁段拼裝情況分別如表1,2所示。

圖4 主梁架設方案

表1 主梁架設方案一施工階段和梁段拼裝率

Table 1 Construction stages and deck erection ratios under deck erection sequence Ⅰ

施工階段12345678910拼裝率/%10.921.732.643.554.365.276.186.997.8100

表2 主梁架設方案二施工階段和梁段拼裝率

Table 2 Construction stages and deck erection ratios under deck erection sequence Ⅱ

施工階段123456789拼裝率/%13.123.934.845.756.567.478.389.1100

3 大跨度懸索橋施工階段顫振穩定性精細化分析

在0°和±3°初始風攻角下,基于上述兩種主梁架設方案,采用Flutter(三維線性顫振分析)和Nflutter(三維非線性顫振分析)程序,取用該橋主梁節段模型風洞試驗測得的靜力三分力系數和顫振導數[15],結構阻尼比為0.5%,進行施工全過程的顫振穩定性分析,揭示其顫振穩定性的演變規律,同時探明靜風效應和風速空間非均勻分布對大跨度懸索橋施工過程顫振穩定性的影響。分析前,采用大跨度懸索橋施工狀態分析程序(IASB)確定兩種主梁架設方案各施工階段的結構幾何和內力狀態,以此作為顫振分析的基準態。由于鉸接法施工的大跨度懸索橋施工階段加勁梁的豎向剛度比成橋狀態的豎向剛度小很多,而扭轉剛度則與成橋狀態的扭轉剛度比較接近,結合已有研究成果,在施工階段分析中加勁梁的豎向剛度、橫向剛度和扭轉剛度分別折減為成橋狀態的30%,40%,80%[16]。橋位地表粗糙度橫橋向按A類場地考慮,地表粗糙度系數α為0.12。

3.1 線性顫振分析

采用線性顫振分析程序Flutter,在0°和±3°初始風攻角和風速空間均勻分布情況下,對該橋上述兩種主梁架設方案施工期的顫振穩定性進行分析,施工全過程顫振臨界風速的變化趨勢如圖5所示。

圖5 兩種主梁架設方案施工過程顫振臨界風速的變化趨勢

在0°風攻角下,該橋模擬主梁節段從跨中向兩側橋塔對稱拼裝施工順序開展了施工階段的全橋模型風洞試驗,測得的顫振臨界風速如圖5(b)所示。可以看出:筆者計算結果與試驗結果在變化趨勢和數值上均較一致,說明筆者計算方法和程序是可靠的。

主梁架設方案一施工過程的顫振臨界風速變化非常平緩,總體上前三階段先增大,而后減小,當拼裝率達到40%后緩慢地增大直至合攏階段,且不同風攻角下的變化趨勢均相似。施工過程結構的顫振臨界風速以-3°攻角最大,0°攻角次之,+3°攻角最小,這與該主梁斷面的氣動性能一致。

主梁架設方案二施工過程的顫振臨界風速變化情況則較架設方案一大為不同,呈現出單一的下降趨勢,尤以前三階段降幅最為明顯,此后下降較平緩。在主梁架設初期,結構的顫振臨界風速最高,接近200 m/s,此后逐漸下降到合攏時的60 m/s左右。與架設方案一類似,施工過程結構的顫振臨界風速也以-3°攻角最大,0°攻角次之,+3°攻角最小。

從圖5兩種主梁架設方案的結果比較可以看出:架設方案二的絕大部分施工階段的顫振臨界風速均顯著大于架設方案一,尤其當主梁拼裝率在40%以下的早期階段。因此,懸索橋加勁梁采用從兩側橋塔向跨中對稱架設方案可以獲得比較好的顫振穩定性,可以避免施工初期的顫振不穩定狀態,從顫振穩定性角度而言,加勁梁采用架設方案二是比較適宜的。

3.2 靜風效應影響分析

在風速空間分布均勻(μ=1.0)以及0°和±3°初始風攻角下,采用Nflutter程序,進行了兩種主梁架設方案施工過程的顫振穩定性分析,靜風效應對施工過程顫振臨界風速的影響分別如圖6,7所示。

圖6 靜風效應對架設方案一施工過程顫振臨界風速的影響

圖7 靜風效應對架設方案二施工過程顫振臨界風速的影響

從圖6,7兩種主梁架設方案下顫振臨界風速隨施工過程變化情況的分析和比較可以得出以下幾點結論:

1)兩種主梁架設方案的非線性顫振分析得到的施工過程顫振臨界風速的演變規律與前述的線性顫振分析基本一致。

2)靜風作用下結構的大變形及剛度變化對施工過程的顫振臨界風速影響較明顯,特別是對于結構剛度較小的施工初期。總體上看,靜風效應對主梁拼裝率在40%以下的各施工階段影響比較顯著,此后影響逐漸減弱。隨著主梁拼裝節段的增多,結構剛度逐漸增強,靜風作用下的結構變形和剛度變化減小,結構顫振穩定性因此受靜風效應的影響減弱。

3)各風攻角下,靜風效應對顫振臨界風速的影響不一,存在著增大或減小現象,這主要和主梁斷面各風攻角下的氣動性能有關。該橋主梁斷面在負攻角下的氣動性能比正攻角好,因此在0°和-3°風攻角下靜風效應使得施工過程結構的顫振穩定性增強,而在+3°風攻角下則有所降低。

3.3 風速豎向變化影響分析

為了揭示風速沿豎向高度方向變化對施工期大跨度懸索橋顫振穩定性的影響,在風速水平均勻分布(e=0)情形下,采用Flutter和Nflutter程序進行風速沿豎向高度變化的施工過程線性和非線性顫振分析,風速豎向變化對施工過程顫振臨界風速的影響分別如圖8,9所示。

圖8 風速豎向變化對架設方案一施工過程顫振臨界風速的影響

圖9 風速豎向變化對架設方案二施工過程顫振臨界風速的影響

從圖8,9的結果分析和比較可以看出:

1)風速沿豎向高度變化對線性顫振臨界風速基本沒有影響,但對非線性顫振臨界風速略有影響,但影響甚微。這主要是顫振分析中僅考慮主梁自激氣動力的作用,沒有考慮主纜自激氣動力的影響,而主梁雖然有豎曲線,但豎向高差不大,對平均風速的影響非常小,故線性顫振分析時風速沿豎向高度的變化產生的影響不明顯,可以忽略不計。在非線性分析時,則同時考慮了主梁和主纜的靜風荷載作用,主纜的豎向高度大,風速沿豎向高度變化較大,主纜的靜風作用對結構的變形、剛度及氣動性能會產生影響,因此考慮風速沿豎向高度的變化因素后顫振臨界風速受到了一定程度的影響。

2)風速沿豎向高度變化因素對施工過程顫振臨界風速的變化趨勢基本沒有影響,但靜風效應對施工過程顫振臨界風速影響依然非常明顯,特別是主梁架設初期,各風攻角下靜風效應對結構顫振臨界風速的影響規律同3.2節所述。

3.4 風場分布寬度影響分析

為了揭示風場分布寬度對施工期大跨度懸索橋顫振穩定性的影響,在風速空間對稱分布情形下(e=0),采用Flutter和Nflutter程序進行不同風場分布寬度的施工過程結構線性和非線性顫振分析,計算得到的各風場分布寬度下的結構顫振臨界風速隨施工過程的變化趨勢分別如圖10,11所示。

圖10 風場分布寬度對架設方案一施工過程顫振臨界風速的影響

圖11 風場分布寬度對架設方案二施工過程顫振臨界風速的影響

從圖10,11的結果分析和比較可以看出:

1)在風場水平向對稱分布的情況下,隨著分布范圍L1的增大,兩種施工方案下各施工階段的顫振臨界風速逐漸減小,逐漸向風速空間均勻分布情況逼近,但變化幅度非常小。從式(2)可知:隨著風場分布范圍L1的增大,風速分布系數逐漸增大,趨向于均勻流場,致使顫振臨界風速向風速空間均勻分布情況逼近。總體而言,風場分布寬度對施工過程結構的顫振臨界風速影響不大。

2)在各風速分布寬度下,靜風效應對施工過程結構的顫振臨界風速影響依舊比較明顯,特別是主梁架設初期,各風攻角下靜風效應對結構顫振臨界風速的影響規律同3.2節所述。

3.5 風速空間非對稱分布影響分析

為了揭示風速空間非對稱分布對施工期大跨度懸索橋顫振穩定性的影響,在風場分布寬度為10L以及風速空間分布非對稱系數e取0.1,0.2,0.3,0.4的情況下,采用Flutter和Nflutter程序進行施工過程結構線性和非線性顫振分析,計算得到的各風速空間非對稱分布系數下的結構顫振臨界風速隨施工過程的變化情況分別如圖12,13所示。

圖12 風速空間非對稱分布對架設方案一施工過程顫振臨界風速的影響

圖13 風速空間非對稱分布對架設方案二施工過程顫振臨界風速的影響

從圖12,13的結果分析和比較可以得到以下幾點結論:

1)在相同風攻角下,各風速空間非對稱分布情況的顫振臨界風速隨施工過程的演變規律基本一致,風速空間非對稱分布因素不影響懸索橋施工期顫振性能的演變規律。

2)隨著風速空間非對稱分布系數e的增大,即風速對稱中心越偏離跨中,顫振臨界風速逐步提高。從式(2)可知:隨著風速空間非對稱分布系數e的增大,橋址范圍內風速分布的不均勻程度加劇,同時各個位置上的實際風速值都在降低,結構實際受到的風荷載減小,結構的顫振臨界風速因而提高。

3)靜風效應對各風速空間非對稱分布情況的顫振臨界風速影響也比較顯著,尤其是在主梁架設初期(主梁拼裝率在40%以下),而后影響逐漸減弱。各風攻角下,靜風效應的影響規律也同3.2節所述。

4)同風場分布寬度因素相比,風速空間的非對稱分布對大跨度懸索橋施工期顫振臨界風速的影響更為顯著。

4 結 論

考慮結構非線性、靜風效應和風速空間非均勻分布等因素,建立了精細化的大跨度橋梁三維非線性顫振分析方法,并編制了計算程序。采用該程序,結合潤揚長江大橋南汊懸索橋,模擬兩種主梁架設順序進行大跨度懸索橋施工全過程顫振穩定性分析,并得出結論:1)加勁梁采用從兩側橋塔向跨中對稱架設方案結構的顫振穩定性要比從跨中向兩側橋塔對稱架設方案好,因此從抗風穩定性角度而言,懸索橋加勁梁宜采用從兩側橋塔向跨中對稱拼裝施工順序。2)靜風作用下結構的大變形及剛度變化對施工過程的顫振臨界風速影響較明顯,特別是對于結構剛度較小的施工初期。總體上看,靜風效應對主梁拼裝率在40%以下的各施工階段影響比較顯著,而后影響逐漸減弱。因此,懸索橋施工過程的顫振分析必須考慮靜風效應的影響。3)風場分布寬度對懸索橋施工期結構的顫振穩定性影響不大,但風速空間的非對稱分布因素對大跨度懸索橋施工期顫振臨界風速的影響非常顯著,值得重視和考慮。

猜你喜歡
風速結構影響
是什么影響了滑動摩擦力的大小
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
哪些顧慮影響擔當?
當代陜西(2021年2期)2021-03-29 07:41:24
基于Kmeans-VMD-LSTM的短期風速預測
基于最優TS評分和頻率匹配的江蘇近海風速訂正
海洋通報(2020年5期)2021-01-14 09:26:54
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
論《日出》的結構
擴鏈劑聯用對PETG擴鏈反應與流變性能的影響
中國塑料(2016年3期)2016-06-15 20:30:00
基于GARCH的短時風速預測方法
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
主站蜘蛛池模板: 被公侵犯人妻少妇一区二区三区| 国模视频一区二区| 性视频一区| 日韩国产无码一区| 五月激激激综合网色播免费| 亚洲精品国产综合99| 一级高清毛片免费a级高清毛片| 久久精品亚洲中文字幕乱码| 九色综合视频网| 啦啦啦网站在线观看a毛片| 国语少妇高潮| 免费女人18毛片a级毛片视频| 中文字幕有乳无码| 亚洲人成网18禁| 国产农村1级毛片| 亚洲Av激情网五月天| 老司机精品一区在线视频| 2021国产精品自产拍在线| 亚洲精品午夜天堂网页| 欧美一级99在线观看国产| 97在线公开视频| 中文字幕伦视频| 在线观看网站国产| 午夜电影在线观看国产1区| 在线a视频免费观看| 国产一级毛片网站| 啪啪永久免费av| 成人国产精品网站在线看| 国产后式a一视频| 欧美国产日韩在线观看| 2019年国产精品自拍不卡| 国产视频只有无码精品| 天天躁狠狠躁| 欧美日韩导航| 不卡无码网| 国产女同自拍视频| 国产屁屁影院| 亚洲中文在线看视频一区| 婷婷综合亚洲| 国产乱人伦AV在线A| 国产在线视频福利资源站| 欧亚日韩Av| 真实国产乱子伦视频| 亚洲欧美不卡视频| 91美女视频在线| 欧美久久网| 精品自拍视频在线观看| 一区二区影院| 国产欧美日韩在线一区| 久久综合五月| 国产91高跟丝袜| 欧美日本中文| 超薄丝袜足j国产在线视频| 亚洲水蜜桃久久综合网站| 午夜色综合| 久久成人免费| 香蕉99国内自产自拍视频| 国产喷水视频| 久久综合九九亚洲一区| 欧美另类视频一区二区三区| 福利片91| 欧美日韩一区二区三区四区在线观看| 欧美第一页在线| 中文字幕 日韩 欧美| 日韩视频精品在线| 亚洲人精品亚洲人成在线| 国产成人亚洲综合A∨在线播放| 99色亚洲国产精品11p| 亚洲成人在线免费观看| 中国国产A一级毛片| 57pao国产成视频免费播放| 日韩少妇激情一区二区| 国产午夜人做人免费视频| 国产91在线免费视频| 亚洲日本中文字幕天堂网| 色悠久久综合| 91最新精品视频发布页| a级毛片免费播放| 国产永久免费视频m3u8| 亚洲人成色77777在线观看| 久久久精品国产亚洲AV日韩| 婷婷综合亚洲|