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

汶川地震近斷層地震動峰值加速度衰減關系研究

2024-12-31 00:00:00李平徐建元鞠雅倩陳鈺鑫歐陽剛壘宣雨童沙明卓
振動工程學報 2024年11期

摘要: 本文選取2008年汶川8.0級特大地震斷層距41 km內13個土層臺站的強震記錄,結合以隨機有限斷層法模擬得到的29個計算點的強震記錄,采用三種近斷層地震動衰減模型和第五代地震動參數區劃圖衰減模型進行了非線性最小二乘法擬合,得到了汶川地震近斷層地震動峰值加速度衰減關系,并通過±1倍標準差和相關性系數進行可靠性分析。結果表明,采用不同衰減模型得到的衰減關系有一定的差異:邵廣彪模型得到的峰值加速度偏低;王國權模型因沒有考慮震級項,其預測結果針對其他震級的地震會產生較大偏差;第五代地震動參數區劃圖衰減模型得到的峰值加速度偏低。

關鍵詞: 汶川地震;"近斷層地震動;"隨機有限斷層法;"衰減關系;"標準差

中圖分類號: P315.9 """文獻標志碼: A """文章編號: 1004?4523(2024)11?1875?09 "DOI:10.16385/j.cnki.issn.1004?4523.2024.11.008

引""言

中國學者對近斷層地震動衰減規律的研究因缺乏近斷層強震記錄,以及受限于近斷層地震動的復雜性而不夠深入。隨著近些年破壞性地震的頻發,中國學者獲取了大量近斷層強震記錄,關于近斷層地震動的研究也取得了一些實質性的突破,得到了一系列研究成果1?3。文獻[4?5]以中國臺灣集集地震發震斷層55 km范圍內的130個臺站記錄為數據庫,分析了地震動峰值加速度、速度、位移和持時的特征,提出了適用于集集地震的近場地震動衰減關系;邵廣彪等6基于全球范圍內6.0~7.6級且震源深度小于20 km的地震,通過對強震記錄峰值加速度進行統計分析,提出了近斷層地震動的衰減模型;KOKETSU等7利用汶川地震發震斷層150 km內54個臺站的加速度記錄,擬合回歸得到了汶川地震近斷層地震動的衰減關系。

除了依靠現有強震記錄進行回歸分析確定地震動衰減關系之外,隨機有限斷層法在建立地震動衰減關系中也得到廣泛應用。傅磊等8建立了龍門山地區與高程相關的κ0模型,模擬得到了汶川地震動PGA(峰值加速度)分布;GHASEMI等9基于隨機有限斷層法,模擬了汶川地震近斷層150 km內54個臺站的PGA,并對模擬結果和經驗衰減模型進行了對比分析,驗證了模擬結果的準確性;喻煙10利用隨機有限斷層法,模擬出了汶川地震主震地震動場,但并未得到北川老縣城的地表地震動。關于汶川地震衰減關系的成果中,并未給出明確的近斷層區域劃分和適用于汶川地震近斷層地震動參數衰減關系的衰減模型。因此,基于已有研究成果,本文利用李明等11對近斷層區域的劃分范圍,確定了汶川地震近斷層區域在斷層距41 km內,并選取了13個土層臺站的強震記錄,結合以隨機有限斷層法模擬得到的29個模擬計算點的強震記錄,采用王國權衰減模型4、邵廣彪衰減模型6、本文模型和第五代地震動參數區劃圖衰減模型12?13進行非線性最小二乘法擬合,得到了汶川地震近斷層地震動峰值加速度衰減關系,并通過±1倍的標準差及相關性系數對擬合曲線的離散程度進行可靠性分析。

1 近斷層區域劃分及臺站選取

2008年5月12日,四川省汶川縣發生8.0級特大地震,震中位于阿壩藏族羌族自治州汶川縣映秀鎮,震源深度為14 km。根據地質學研究表明:控制本次地震的是北川—映秀斷裂帶,地表破裂長度約240 km,以逆斷層破裂為主并伴有右旋走滑分量14?15。此次地震共計觸發四川省強震臺站141個,參考李明等1116在2010年給出的近斷層地震動區域的劃分標準,根據USGS場地分類標準、斷層類型和矩震級三個指標,確定汶川地震近斷層區域為斷層距41 km內,選取了近斷層41 km內的13個土層臺站和2個基巖臺站(選取兩個基巖臺站作為后文的隨機有限斷層法模擬地震動的參照臺站)。本文臺站信息來自中國地震局工程力學研究所強震動觀測與工程振動研究中心,使用的強震記錄數據均已進行基線校正和濾波。本文計算衰減關系采用斷層距作為距離項,具體計算方法是以汶川地震地表破裂帶作為發震斷層在地表的投影位置,根據多個地表破裂帶現場調查點的經緯度坐標和臺站的經緯度坐標計算地表破裂帶現場調查點和臺站位置兩點之間的距離,取其中的最小值為斷層距。臺站分布如圖1所示,臺站信息如表1所示。

2 汶川地震隨機有限斷層模擬

本文以2008年汶川特大地震為例,基于文獻[17?18]反演得到的單斷層位錯模型,確定了沿走向280 km和沿傾向40 km的矩形計算區域,計算點共計9060個,如圖2所示。參考文獻[8?10]的研究成果,確定了模擬汶川特大地震地震動場的參數,如表2所示,表中R為斷層距。

由于隨機有限斷層模擬的是水平向基巖地震動,為了和13個土層臺站(Ⅱ類場地)的PGA進行對比分析20,利用《中國地震動參數區劃圖》(GB 18306—2015)21中的場地地震動峰值加速度調整系數,采用線性內插法調整得到臺站地表的PGA模擬計算結果,如表3所示。表3中除51SFB,51JYH,51JYD臺站誤差稍大,其余臺站的強震記錄數據和模擬誤差平均在20%左右,誤差較大的三個臺站可能是由于土層臺站位于山地地區,地形地貌對地震動影響顯著。

為了驗證利用隨機有限斷層法模擬汶川地震主震地震動場的準確性,本文利用《中國地震烈度表》(GB/T 17742—2020)22中的地震烈度對應的地震動峰值(地震動峰值為所對應的儀器測定的地震烈度中值),將汶川地震Ⅸ,"Ⅹ,"Ⅺ烈度等值線轉換為地震動峰值等值線,并與模擬的地震動場等值線進行對比分析,結果如圖3所示。模擬得到的汶川地震近場加速度等值線和烈度等值線基本吻合,其中Ⅹ和Ⅺ烈度區的等值線和峰值加速度分區吻合度較高。

汶川地震中近斷層41 km內有13個土層臺站和2個基巖臺站。由于臺站僅分布在4~6 km,"18~25 km,"30~36 km,強震記錄數據少且分布不均,所以僅用強震記錄數據擬合近斷層地震動衰減關系不夠精確。本文通過汶川地震隨機有限斷層法在斷層距41 km范圍內增加了29個模擬計算點(依據實測臺站的位置、模擬計算點的PGA值和斷層距確定補充計算點的位置,使其均勻分布在近斷層41 km內,便于后文衰減關系的擬合),汶川地震隨機有限斷層法計算區域和模擬計算點如圖4所示。模擬補充基巖處的計算點通過《中國地震動參數區劃圖》(GB 18306—2015)21中的場地地震動峰值加速度調整系數進行調整,采用線性內插法調整得到模擬地表PGA的計算結果如表4所示。

綜上分析,以隨機有限斷層法模擬近斷層地震動是可靠的,所得數據可用于后文中的汶川地震近斷層衰減關系擬合。

3 近斷層地震動衰減關系模型

目前常用的近斷層地震動衰減關系模型有邵廣彪衰減模型6、王國權衰減模型4,本文衰減模型以及第五代地震動參數區劃圖青藏區短軸衰減模型,具體闡述如下:

邵廣彪衰減模型6形式為:

式中""Y為地震動峰值PGA;R為斷層距;M為矩震級;H為震源深度;c0c3為待求的參數;ε為隨機誤差。

由于本文所研究的是汶川地震主震,式(1)中的變量矩震級、震源深度都是常數,所以經過修改后的衰減模型為:

式中""為擬合值的標準差。

王國權等4基于1999年中國臺灣集集(Chi?Chi)地震提出的衰減模型為:

文獻[9]提出的PGA衰減模型為:

式中""a0a3為待求的參數。

依據已有強震記錄提出的近斷層地震動PGA衰減模型為(本文所使用的衰減模型):

第五代地震動參數區劃圖青藏區短軸衰減模型20為:

采用以上三種衰減模型和第五代地震動參數區劃圖衰減模型,對實測臺站的水平向記錄、豎向記錄、實測臺站記錄和模擬記錄進行非線性最小二乘法擬合(最小二乘法應用到曲線擬合時,當自變量和因變量同時存在均值為零,方差相同的隨機誤差時,此方法能給出在統計意義上最好的參數擬合結果),回歸結果、擬合優度R2和標準差σ如表5所示。第五代地震動參數區劃圖青藏區短軸衰減模型參數已確定。

4 回歸結果分析

4.1 本文衰減關系與其他衰減關系對比

本文基于以上三種衰減模型和第五代地震動參數區劃圖衰減模型,分別對實測強震記錄和模擬強震記錄進行非線性最小二乘法擬合(水平向峰值加速度衰減關系擬合使用13個土層臺站的峰值加速度和隨機有限斷層模擬得到的29個補充記錄,豎向地震動峰值加速度進行衰減關系擬合時,僅使用13個土層臺站的峰值加速度進行回歸擬合),如圖5(a),(b)所示,邵廣彪6和王國權模型4的結果對地震動峰值加速度估計過低,而本文模型和第五代地震動參數區劃圖衰減模型地震動峰值加速度較為接近。如圖5(c)所示,使用隨機有限斷層法對模擬計算點強震記錄和實測強震記錄進行擬合,除了邵廣彪模型6對峰值估計過低之外其他三種模型峰值均達到了1.0g左右,所以隨著隨機有限斷層法模擬計算點的補充,近斷層區域內衰減規律較為明顯,而第五代地震動參數區劃圖衰減模型和本文模型較為接近,但預測峰值加速度偏低。根據整個汶川地震的宏觀現象以及臺站記錄到的強震記錄來看,距離震中最近的臥龍臺PGA高達957.7 gal,而遠離震中的北川震害情況非常嚴重,與震中相比甚至更加顯著。所以北川老縣城的地表地震動估計值要大于臥龍臺所記錄到的結果,至少要達到1.0g以上。故根據北川老縣城的實際震害以及地表地震動估計,本文模型的結果更為合理。

4.2 標準差分析

本文將實測強震記錄和模擬強震記錄分別用邵廣彪衰減模型6、王國權衰減模型4和本文模型進行預測值和觀測值的對比分析,結果如圖6~8所示,圖中實線為擬合曲線,虛線為估計值的±1倍標準差值。由圖6可知,邵廣彪模型6的觀測值大部分均處在±1倍標準差內,只有極少數大于預測值范圍,增加隨機有限斷層模擬點強震記錄之后擬合結果有明顯提高但是整體趨勢變化不大23?24。由圖7可知,王國權衰減模型4擬合曲線與其他三種模型擬合曲線相比整體趨勢過低,尤其是豎向觀測值不符合實測記錄以及震害特征,增加了模擬計算點強震記錄之后整體趨勢變化較大,對峰值的估計增大了0.2g左右,但是由于王國權衰減模型4沒有考慮震級項,故而預測結果對于其他震級的地震會產生較大偏差。由圖8可知,本文結果對PGA的估計值較為合理,符合北川老縣城的震害特征,整體衰減趨勢較快,觀測值大部分處在±1倍標準差內,增加了模擬計算點強震記錄之后,整體趨勢相似性較高,且擬合優度接近0.8,標準差也在0.2左右,不僅說明本文模型可以很好地估計汶川地震近斷層區域內的PGA,同時也反映了北川老縣城地表地震動的經驗估計。所以本文模型的估計值整體上較為合理,符合汶川地震近斷層PGA特征。

5 結論與展望

本文選取了汶川地震近斷層內13個土層臺站的強震記錄,結合以隨機有限斷層法模擬得到的29個模擬計算點的強震記錄,使用邵廣彪模型6、王國權模型4、本文模型以及第五代地震動參數區劃圖衰減模型11?12,對水平和豎向峰值加速度進行非線性最小二乘法擬合,利用相關性系數和±1倍標準差進行了可靠性分析,得出以下結論:

(1)"本文衰減模型對實測臺站水平向記錄、豎向記錄的擬合回歸結果較為合理,絕大多數峰值加速度觀測值分布在本文衰減模型±1倍標準差以內。斷層距1 km范圍內峰值加速度達到了1.0g,符合北川老縣城的實際震害情況,基本達到了預測值范圍。

(2)"在實測臺站強震記錄和模擬強震記錄中,邵廣彪模型的估值過低不符合實際震害特征,王國權模型因沒有考慮震級項導致預測結果對于其他震級的地震會產生較大偏差,而本文模型和第五代地震動參數區劃圖衰減模型的結果較為合理,但是第五代地震動參數區劃圖衰減模型對于本文模型來說峰值較低,所以本文模型較適用于汶川特大地震近斷層區域的峰值加速度預測。

本文并未討論模擬地震動計算點參與到衰減關系擬合中可能存在的誤差,本文模型適用性較為單一,僅適用于汶川地震近斷層41 km內的衰減關系擬合。對隨機有限斷層法補充計算點位置的確定并未做深入研究,可能會造成一些誤差。斷層上下盤的震害程度不同,本文因數據點數量有限并未在衰減關系中考慮上下盤效應的影響。未來工作可以從以上問題入手,進行更加深入的討論和研究。

參考文獻:

[1] COUSINS W J,"ZHAO J X,"PERRIN N D. A model for the attenuation of peak ground acceleration in New Zealand earthquakes based on seismograph and accelerograph data[J]. Bulletin of the New Zealand Society for Earthquake Engineering,"1999,"32(4):193?220.

[2] 江輝,"慎丹,"倪永軍,等. 近斷層地震地面運動峰值衰減規律研究[J]. 北京交通大學學報,"2011,"35(4):83?87.

JIANG Hui,"SHEN Dan,"NI Yongjun,"et al. Research on peak attenuation relationship of near?fault earthquake ground motion[J]. Journal of Beijing Jiaotong University,"2011,"35(4):83?87.

[3] 趙金鑫. 汶川地震近斷層效應研究[D]. 哈爾濱:中國地震局工程力學研究所,"2013.

ZHAO Jinxin. Study on near?fault effects of Wenchuan earthquake[D]. Harbin:"Institute of Engineering Mechanics,"China Earthquake Administration,"2013.

[4] WANG G Q,"ZHOU X Y,"ZHANG P Z,"et al. Characteristics of amplitude and duration for near fault strong ground motion from the 1999 Chi?Chi,"Taiwan earthquake[J]. Soil Dynamics and Earthquake Engineering,"2002,"20(1):73?96.

[5] SHABESTARI K T,"YAMAZAKI F. Attenuation relationships of the ground motion parameters considering directivity effects in the 1999 Chi?Chi,"Taiwan earthquake[C]∥Proceedings of Seventh US National Conference on Earthquake Engineering. Boston,"Massachusetts,"2002:"21?25.

[6] 邵廣彪,"馮啟民. 近斷層地震動加速度峰值衰減規律的研究[J]. 地震工程與工程振動,2004,"24(3):30?37.

SHAO Guangbiao,"FENG Qimin. Research on attenuation of near?fault peak strong ground motion acceleration[J]. Earthquake Engineering and Engineering Vibration,"2004,24(3):30?37.

[7] KOKETSU K,"YOKOTA Y,"GHASEMI H,"et al. Source process and ground motions of the 2008 Wenchuan earthquake[C]∥Proceedings of International Conference on Earthquake Engineering—The First Anniversary of Wenchuan Earthquake. 2009:"11?12.

[8] 傅磊,李小軍. 龍門山地區的Kappa(κ0)模型及汶川Ms8.0地震的強地震動模擬[J]. 地球物理學報,"2017,"60(8):2935?2947.

FU Lei,"LI Xiaojun. The Kappa(κ0)"Model of the Longmenshan region and its application to simulation of strong ground?motion by the Wenchuan Ms8.0 earthquake[J]. Chinese Journal of Geophysics,"2017,"60(8):2935?2947.

[9] GHASEMI H,"FUKUSHIMA Y,"KOKETSU K,"et al. Ground?motion simulation for the 2008 Wenchuan,"China,"earthquake using the stochastic finite?fault method[J]. Bulletin of the Seismological Society of America,2010,100(5B):"2476?2490.

[10] 喻煙. 汶川地震區地震動估計經驗模型[D]. 哈爾濱:"中國地震局工程力學研究所,"2012.

YU Yan. Empirical estimate model for ground motion of Wenchuan earthquake zone[D]. Harbin:"Institute of Engineering Mechanics,"China Earthquake Administration,"2012.

[11] 李明,"謝禮立,"翟長海,等. 近斷層地震動區域的劃分[J]. 地震工程與工程振動,"2009,"29(5):20?25.

LI Ming,"XIE Lili,"ZHAI Changhai,"et al. Scope division of near?fault ground motion[J]. Earthquake Engineering and Engineering Vibration,"2009,"29(5):20?25.

[12] 俞言祥. 新一代地震區劃圖地震動參數衰減關系的建立與特點分析[J]. 城市與減災,"2016(3):34?38.

YU Yanxiang. Establishment and characterization of ground motion parameter attenuation relationships in new generation seismic zoning maps[J]. City and Disaster Reduction,"2016(3):34?38.

[13] 俞言祥,李山有,肖亮.為新區劃圖編制所建立的地震動衰減關系[J]. 震災防御技術,"2013,"8(1):24?33.

YU Yanxiang,"LI Shanyou,"XIAO Liang. Development of ground motion attenuation relations for the new seismic hazard map of China[J]. Technology for Earthquake Disaster Prevention,2013,"8(1):24?33.

[14] 袁一凡,"田啟文. 工程地震學[M]. 北京:"地震出版社,"2012.

[15] 李小軍,于愛勤,"甘朋霞,"等. 汶川8.0級地震北川縣城區災害調查與分析[J]. 震災防御技術,"2008,"3(4):352?362.

LI Xiaojun,"YU Aiqin,"GAN Pengxia,"et al. Survey and analysis of the disaster and engineering damage of Beichuan County seat in Ms8.0 Wenchuan earthquake[J]. Technology for Earthquake Disaster Prevention,"2008,"3(4):352?362.

[16] 李明. 近斷層地震動對結構抗震設計的影響研究[D]. 哈爾濱:"中國地震局工程力學研究所,"2010.

LI Ming. Study of the effect of near?fault ground motion on structural seismic design[D]. Harbin:"Institute of Engineering Mechanics,"China Earthquake Administration,"2010.

[17] BJERRUM L W,"SΦRENSEN M B,"ATAKAN K. Strong ground?motion simulation of the 12 May 2008 Mw7.9 Wenchuan earthquake,"using various slip models[J]. Bulletin of the Seismological Society of America,"2010,"100(5B):"2396?2424.

[18] JI C. Preliminary result of the May 12,"2008 Mw7.97 Sichuan earthquake[DB/OL].https://ji.faculty.geol.ucsb.edu/big_earthquakes/2008/05/12/ShiChuan.html,2008.

[19] 鐘菊芳,"袁峰. 震源參數對合成時程持時的影響分析[J]. 防災減災工程學報,"2019,"39(5):"733?747.

ZHONG J F,"YUAN F. Influence of source parameters on the duration of simulated ground motion[J]. Journal of Disaster Prevention and Mitigation Engineering,"2009,"39(5):"733?747.

[20] 于海英,"王棟,"楊永強,"等. 汶川8.0級地震強震動加速度記錄的初步分析[J]. 地震工程與工程振動,"2009,"29(1):"1?13.

YU Haiying,"WANG Dong,"YANG Yongqiang,et al. The preliminary analysis of strong ground motion records from the Ms8.0 Wenchuan earthquake[J]. Earthquake Engineering and Engineering Vibration,"2009,"29(1):1?13.

[21] 中國地震局. 中國地震動參數區劃圖:"GB 18306—2015[S]. 北京:中國建筑工業出版社,"2015.

China Earthquake Administration. Seismic ground motion parameters zonation map of China:"GB 18306—2015[S]. Beijing:"China Architecture amp; Building Press,2015.

[22] 中國地震局. 中國地震烈度表:"GB/T 17742—2020 [S]. 北京:中國建筑工業出版社,"2021.

China Earthquake Administration. The Chinese seismic intensity scale:"GB/T 17742—2020[S]. Beijing:"China Architecture amp; Building Press,"2021.

[23] 張斌,俞言祥,李小軍,"等. 西南地區水平向峰值速度、峰值位移衰減關系研究[J].地球物理學報,"2021,"64(8):"2733?2748.

ZHANG Bin,"YU Yanxiang,"LI Xiaojun,"et al. Ground motion attenuation relationship of horizontal component of PGV and PGD in Southwest China[J]. Chinese Journal of Geophysics,"2021,"64(8):"2733?2748.

[24] 張斌,李小軍,俞言祥,等. 魯甸地震強震動記錄與地震動衰減模型的對比研究[J].地球物理學報,"2020,"63(8):2999?3014.

ZHANG Bin,"LI Xiaojun,"YU Yanxiang,"et al. Comparison of strong ground motion records from Ludian,"China,"earthquake with ground?motion attenuation models[J]. Chinese Journal of Geophysics,"2020,63(8):2999?3014.

Attenuation relationship of peak ground acceleration near the fault of the Wenchuan earthquake

LI"Ping12,"XU"Jian?yuan13,"JU"Ya?qian12,"CHEN"Yu?xin12,"OUYANG"Gang?lei14,"XUAN"Yu?tong12SHA Ming-zhuo5

(1.College of Geological Engineering,Institute of Disaster Prevention,"Sanhe 065201,"China;2.Key Laboratory of Earthquake Disaster Prevention and Risk Assessment of Hebei Province,"Sanhe 065201,"China;"3.Municipal Planning and Engineering Design Institute,"Jiangsu Provincial Planning and Design Group Co.,"Ltd.,"Nanjing 210019,"China;"4.Nuclear Industry Shaft Construction Group Co.,"Ltd.,"Huzhou 313001,"China;"5. Shenyang Water Supply amp; Drainage Prospecting Design Research Institute Co.,"Ltd.,"Shenyang 110023,"China)

Abstract: In this paper,"the strong earthquake records of 13 soil stations within 41 km from the fault of the 2008 Wenchuan 8.0 magnitude earthquake are selected. These records are combined with the strong earthquake records of 29 calculation points obtained by using the stochastic finite fault method simulation. Moreover,"nonlinear least squares fitting is performed by using three near-fault ground motion attenuation models and the fifth-generation ground motion parameter zoning map attenuation model to obtain the peak ground motion acceleration of the Wenchuan earthquake near-fault. The attenuation relationships are obtained and analyzed by considering ±1 times standard deviation and correlation coefficients for reliability. The results show that there are some differences in the attenuation relations obtained by different attenuation models. The peak acceleration obtained by the Shao Guangbiao model is low. The prediction results of Wang Guoquan model has a large deviation for other magnitude earthquakes because the magnitude term is not considered. Besides,"the peak acceleration of the fifth-generation ground motion parameter zoning map attenuation model is low.

Key words: Wenchuan earthquake;"near-fault seismic ground motion;"stochastic finite tomography;"attenuation relationship;"standard deviation

作者簡介: 李""平(1981―),男,博士,教授。"E-mail:"chinaliping1981@126.com。

主站蜘蛛池模板: 福利在线一区| 亚洲第一成年人网站| 久久精品亚洲专区| 亚洲综合色区在线播放2019| 成人精品亚洲| 依依成人精品无v国产| 伊人激情综合| 精品伊人久久久大香线蕉欧美| 激情综合激情| 97精品国产高清久久久久蜜芽| 欧美在线一二区| 国产精品尤物铁牛tv| 97青草最新免费精品视频| 欧美激情伊人| 国产一区二区三区在线无码| 日本成人精品视频| 色亚洲激情综合精品无码视频| 伊人中文网| 亚洲成a人在线观看| 亚洲一区免费看| 日韩精品中文字幕一区三区| 日本免费福利视频| 国产成人一区| 伊人AV天堂| 成人福利一区二区视频在线| 国产精品网曝门免费视频| 97综合久久| 国产精品久久久免费视频| 欧美亚洲国产日韩电影在线| 熟女成人国产精品视频| 久久人搡人人玩人妻精品| 久久一级电影| 国产在线拍偷自揄拍精品| 午夜国产理论| 这里只有精品在线播放| 亚洲午夜国产精品无卡| 九月婷婷亚洲综合在线| 久无码久无码av无码| 综合成人国产| 香蕉视频在线精品| 中文一区二区视频| 99色亚洲国产精品11p| 久久毛片基地| 久久久久久久久亚洲精品| 欧美中日韩在线| 欧美一级黄片一区2区| 国产在线观看第二页| 91福利在线看| 91福利免费| 亚洲国产黄色| 亚洲综合色婷婷| 国内精品91| 91最新精品视频发布页| 免费激情网址| 国产高颜值露脸在线观看| 99精品在线视频观看| 自拍欧美亚洲| 99re在线免费视频| 婷婷久久综合九色综合88| 亚洲精品成人片在线观看| 在线看免费无码av天堂的| 91福利免费视频| 91www在线观看| 又猛又黄又爽无遮挡的视频网站| 狠狠躁天天躁夜夜躁婷婷| 欧美日韩国产在线人| 国产成人91精品免费网址在线| 人妻21p大胆| 日韩国产亚洲一区二区在线观看| 毛片免费观看视频| 在线国产毛片| 美美女高清毛片视频免费观看| 国产高清又黄又嫩的免费视频网站| 免费高清a毛片| 一级片一区| 亚洲综合片| 99这里精品| 2048国产精品原创综合在线| 久久9966精品国产免费| 欧美一级夜夜爽| 国产丰满大乳无码免费播放 | 国模极品一区二区三区|