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

注水誘發地震的譜振幅相關系數及視應力研究*

2015-12-14 09:31:26張致偉喬慧珍趙小艷
地震研究 2015年1期
關鍵詞:機制區域研究

張致偉,喬慧珍,吳 朋,宮 悅,趙小艷

0 引言

地震是地殼構造應力積累超過了巖層的受力極限,應變能快速釋放,斷裂帶破裂失穩的結果,大地震的發生必然受到構造應力場的制約,而震源機制是研究構造應力場的基礎資料,不同階段的震源機制解在一定程度上反映了大震孕育過程中構造應力場的變化。國內外諸多學者研究認為基于震源機制解反演應力場的時空變化能夠捕捉到強震前的震兆信息,震源機制一致性參數是判斷強震危險性的一個有用判據 (Michael et al.,1990;Hauksson,1994;Wiemer et al.,2002;陳颙,1978;刁桂苓等,1994;程萬正等,2006b;刁桂苓等,2011)。

震源機制一致性參數雖然物理含義明確,但由于受地震臺網布局和小震波形信噪比的制約,往往難以給出比較精確的小震震源機制解。鑒于上述問題,Lund和B?evarsson(2002)提出了采用微震體波譜振幅相關系數來研究地震序列震源機制是否相似的問題。國內學者隨后也將該方法運用于強震震例研究中 (朱航等,2006;崔子健等,2012),并取得了較好的成果。研究表明強震前的小震群序列與一般性小震群序列相比,譜振幅相關系數接近于1.0,說明震源區應力場增強,震源機制相似性更為顯著。且譜振幅相關系數的測定和分析,只需地震序列周圍若干數字地震臺站即可,具有一定的可操作性。

四川東南部存在多處采鹽、采氣及廢水回注井,地震活動普遍表現為強度低、頻次高、持續時間長的特征,且與當地加壓注水及采氣 (鹵)等工業活動具有很好的相關性 (阮祥等,2008;張致偉等,2012)。本文選取四川自貢和長寧兩個注水地區作為研究區域,基于測震臺網記錄的數字地震波形資料和觀測報告,分析注水誘發地震的譜振幅相關系數特征及其與當地加壓注水數據的關系,并結合表征絕對應力水平的視應力綜合分析研究區域的應力高低。

1 資料和方法

1. 1 研究區地震活動概況

本文以自貢和長寧兩個注水地區為例 (圖1a中圓形區域)。2000年以來研究區域ML≥2.0地震時序圖 (圖1b)顯示,自貢地區在2000~2008年地震活動較弱,未發生ML≥3.0地震,2009年1月上旬該區小震活動明顯增強,并于2009年2月16日、5月22日先后發生了ML4.4和 ML4.2地震,且在2011、2012年又發生了兩次比較顯著的4級地震震群活動 (圖1中圓形區域b),張致偉等 (2012)研究表明該區地震活動主要受當地廢水回注井的注水壓力調控。長寧地區在2000~2006年3月地震活動同樣呈現相對較弱狀態,2006年4月以來該區域地震活動頻次和強度明顯升高,并于2008年2月1日發生了ML4.8地震,隨后該區域地震頻次略有下降(圖1c),但ML≥4.0地震依然持續發生,并于2013年4月25日發生了ML5.2地震 (圖1a中圓形區域c),阮祥等 (2008)指出對該區域的地震活動阮祥等 (2008)研究表明該區域的地震活動主要受長寧雙河鹽礦注 (出)水量的影響。

1. 2 譜振幅相關系數測定方法

地震波能量在傳播過程中不但隨著傳播距離出現幾何擴散,還受到傳播路徑上介質的吸收和散射等影響,在到達臺站接收儀器前地表下方淺層的介質也會對地震波產生影響,因此觀測得到的地震數據是震源激發的信息經過上述各種過程的產物。在頻率域,臺站j記錄到的地震i的觀測位移譜Uij(f)可表示為

式中,S'i(f)為地震i的震源譜;φij為地震i的震源輻射圖型因子;Pij(f)為傳播路徑效應,包括地震波的幾何擴散和非彈性衰減;L'j(f)為第j個臺站的局部場地效應;N'j(f)為第j個臺站附近的地面運動噪聲;I'j(f)為第j個臺站儀器響應;Surj為臺站附近地表自由表面效應。

根據Brune模型 (1970),S'i(f)可表示為

式中,Ω'0i為地震i的震源譜低頻漸近線值,即零頻極限值;fc為低頻漸近線與高頻漸近線交點處的頻率,稱為拐角頻率。

定義譜振幅Ω0ij為臺站j記錄到的地震i的包含輻射圖型因子的震源譜零頻極限值,可表示為

式 (1)中的Pij(f)可表示為

式中,Gij為幾何擴散因子;Rij為震源距;Q(f)為介質的品質因子;v為地震波 P波或 S波傳播速度。

從地震記錄中扣除噪聲和儀器響應,并令Lj(f)=L'j(f)Surj)(Lj稱為廣義的臺站場地效應),地震i在臺站j的地面運動位移譜Uij(f)為

可見,若臺站j記錄的兩次地震x、y震源位置足夠接近,即其間距比震源距Rij小得多,且震源機制 (輻射圖型因子)相同,則臺站j記錄的這兩次地震的地面運動位移譜Uxj(f)與Uyj(f)的相對大小只與這兩次地震的Ω'0相對大小相關;若震源機制不同,則還應與這兩次地震震源機制的差異有關。反過來說,當由地震記錄來反演譜振幅時,若這兩次地震震源機制相同,則反演得到的Ω0xj與Ω0yj的相對大小只與Ω'0x和Ω'0y的相對大小有關;如果震源機制不同,則與震源機制有關。因此可由臺站的地震記錄反演地震x、y的振幅譜Ω0xj與Ω0yj,通過計算其相關系數rxy來描述震源機制是否相似的問題,若震源機制相似,譜振幅相關系數rxy應較大,接近1;若震源機制不相似,則相關系數較小。地震x和y的譜振幅相關系數 rxy可表示為(Lund,B?evarsson,2002)

式中,xjl、yjl分別表示臺站j記錄的地震x、y的l分量波形資料所反演得到的譜振幅Ω0xj與Ω0yj的對數值;x、y分別為xjl與tjl的平均值;n為使用的臺站數目;l取值為1到5,依次為垂向和徑向P波,垂向、徑向和切向S波。

本研究考慮到持續地震震源機制是否相似是針對多數地震而言的,因此采用崔子健等 (2012)提出的計算方法,將地震按發生時間順序排列,地震m與前面m-1個地震為一組,計算組內每兩個地震的rxy,得出N=m(m-1)/2個相關系數rxy,對rxy求算數平均值,其結果代表了地震m時刻的組內譜振幅的相關程度,以步長為1進行滑動,計算每個組內的譜振幅相關系數的算數平均值,最后可以得到隨時間變化的譜振幅相關系數。

2 譜振幅相關系數結果分析

2. 1 場地響應和區域Q值反演

地震儀器記錄到的地震波形資料包含震源信息、傳播路徑及臺站場地響應等,要想獲得趨于更真實的震源參數,必須對觀測信號做地震波的幾何擴散、非彈性衰減以及臺站場地響應校正 (劉杰等,2003)。基于2009-1-1~2012-12-31四川東南地區20次ML≥3.5地震波形資料及觀測報告,結合Atkinson和Mereu(1992),Moya等 (2000)提出的方法,反演獲得了川東南地區P、S波Q值和7個測震臺站的場地響應。

Q值是描述地球介質特性的主要參數,與地質構造、地震活動性及地域熱流有密切的關系 (周連慶等,2008)。Atkison方法是考慮頻率依賴的Q模型,通常用頻率的冪函數,即Q(f)=Q0fη的形式來擬合Q值與頻率f的關系,參數η反映了Q值對頻率的依賴程度。圖2給出了川東南地區P、S波Q值隨頻率的變化及其擬合曲線,P波Q=53.8f0.879,S 波 Q=101.9f1.025,與阮祥等 (2008)基于早期資料計算獲得長寧地區的尾波Q值比較一致,而與川西北地區的Q值結果 (喬慧珍等,2006)差別較大,表明Q值具有地區差異特征。

圖3給出了川東南地區7個測震臺站的場地響應結果,圖中藍線是由每個地震記錄得到的臺站場地響應,紅線則是該臺站場地響應的平均值。臺站場地響應總體表現為高頻衰減特征,JLI、MBI、LBO和LZH臺在1~4 Hz的放大倍數在1~2之間,而其余臺站的放大效應都大于2。臺站場地響應形態之所以存在差異,主要受臺基、周邊地形等因素的影響。臺站場地響應的計算不僅對研究區臺站場地選點有一定的參考意義,而且也為下一步計算譜振幅相關系數奠定基礎。

2. 2 注水誘發地震的譜振幅相關系數

在獲得川東南地區P、S波的Q值及臺站場地響應的基礎上,以自貢和長寧兩個注水地區為例,采用上述方法研究注水誘發地震的譜振幅相關系數特征及其與當地加壓注水數據的關系。

選取自貢地區2009-1-1~2013-4-30共32次 ML≥3.0地震,基于馬邊 (MBI)、漢王山(HWS)、花馬寺 (HMS)、五馬坪 (WMP)和瀘州臺 (LZH)記錄的地震波資料和觀測報告,計算獲得該區域地震譜振幅相關系數在0.84~0.92之間波動,其平均值為0.89。結合研究區地震活動及注水井的加壓數據,將研究時段分為兩個階段進行分析 (圖4a),第一階段為2009-1~2010-7,該時段注水壓力持續升高,地震活動顯著,譜振幅相關系數由0.84逐漸上升至0.92,呈現隨加壓數據的升高而增大的特征;第二階段為2010-8~2013-4,該時段注水井停止加壓,地震活動明顯減弱,譜振幅相關系數也有逐漸降低的趨勢。

選取長寧地區2006-1-1~2013-4-30共88次 ML≥3.0地震,基于馬邊 (MBI)、雷波(LBO)、花馬寺 (HMS)、五馬坪 (WMP)和瀘州臺 (LZH)記錄的地震波資料和觀測報告計算獲得了譜振幅相關系數,該區域地震相關系數總體上在0.87~0.97之間波動,其平均值為 0.90。結合研究區地震活動及鹽礦注 (出)水量數據,同樣將研究時段分為兩個階段 (圖4b),第一階段為2006-4~2007-11,相比前期注 (出)水量,該時段內注水量大幅增加,但出水量的增加并不明顯,造成注水量和出水量之差大幅上升,小震活動明顯增強,相應的譜振幅相關系數也在高值(0.90~0.97)波動,而且在注 (出)水量差值變化最大時相關系數出現了明顯的階躍;第二階段為2007-12~2013-4,該時段注水量和出水量相當,且有逐漸下降的趨勢,但地震活動依然持續,譜振幅相關系數隨著注 (出)水量的下降也有所降低,但仍然大于0.87,且比較穩定。

綜上所述,自貢和長寧地區注水誘發地震的譜振幅相關系數均在高值波動,其均值分別為0.89和0.90,且變化過程與加壓注水數據具有很好的相關性。Seeber和Armbruster(2000)認為地震可作為應力變化的標志,雖然小地震對總的變形貢獻很小,但它們在空間、時間和運動上的分布對應力變化非常敏感。影響地殼力學狀態的各種現象,包括諸如斷層蠕動或巖漿注入的自然現象,如水庫蓄水或油田注水的大型工程活動。自貢和長寧地區的介質因加壓注水、采鹽等工業活動已發生了變化,注水對局部區域應力場產生擾動影響,使地下淺層裂隙呈優勢取向排列,引發的中小地震震源機制表現出較好的相似性,體現在較高的譜振幅相關系數上。

3 地震視應力結果分析

地震視應力是在有關震源介質的均勻彈性、地震破裂的脆性性質假定下 (Kanamori,Anderson,1975;Kanamori,1977;Wu,2001),可由地震波信息獲得的、與區域平均應力呈正比的物理量。對一個地區中引起地震滑動的應力水平進行區域平均,則可作為當地絕對應力水平的一個間接估計 (吳忠良等,2002)。

地震視應力σapp被定義為單位地震矩輻射出的地震波能量與震源區介質的剪切模型的積 (Wyss,Brune,1968),即

式中,ES為地震波輻射能量;M0為地震矩;μ為震源區介質剪切模量,通常取 μ=3.5×1011dyn/cm2。

基于Brune(1970)圓盤模型,地震震源譜可由震源譜的低頻水平Ω'oi和拐角頻率fc表示 (式(2))。而地震波輻射能量ES和地震矩M0分別表示為

式中,β為S波波速,取β=3.5 km/s;ρ為介質密度,取ρ=2.71 g/cm3;Rθφ為輻射圖型因子,取S波的均方根輻射圖形因子■2/5≈0.63。

計算過程中,一般選取區域地震臺網震中距≤200 km的臺站記錄的具有較高信噪比的寬頻帶數字化地震波形資料,利用S波頻帶范圍在1.0~20 Hz的波形計算震源譜,在獲得地震波輻射能量和地震矩的基礎上,再由 (7)式逐臺求取地震視應力,最終取各臺計算結果的平均值作為相應地震的視應力。

自貢和長寧地區注水誘發地震的譜振幅相關系數較高,可能反映了中小地震震源機制具有很好的相似性,為了配合譜振幅相關系數分析研究區域的應力水平,基于上述方法計算獲得了2000~2013年川東南地區ML≥3.0地震的視應力結果。圖5分別給出了獲得地震視應力的ML≥3.0地震及其視應力空間分布圖,由圖5a可見,3級地震空間分布比較分散,ML≥4.0地震則主要集中分布在四川自貢、長寧、重慶榮昌和云南東北地區。視應力相對高值區均在發生4級地震的地區,唯獨四川自貢和長寧4級地震集中區 (圓圈)的視應力較低 (圖5b),表明川東南局部地區地震視應力存在差異,自貢和長寧地區的應力水平相對偏低。

為了進一步研究自貢和長寧注水地區地震視應力的特征,圖6給出了自貢、長寧及川東南地區ML≥3.5地震的視應力對比圖。由圖可見,川東南地區地震視應力總體呈現隨震級增大而增大的特征,但與程萬正等 (2006a)給出的四川其它地區視應力相比,川東南地區同等震級的地震視應力相對偏低。自貢地區不同震級檔的地震視應力均在川東南地震視應力擬合線下方,長寧地區不同震級檔的地震視應力分布在擬合線附近,與川東南其它區域地震視應力相比,自貢和長寧地區同等震級的地震視應力相對偏低。地殼介質應力水平有隨深度增加而變大的特征 (丁健民,高莉青,1981;李方全,祁英男,1988),自貢和長寧地區小震精定位均顯示其震源深度較淺 (阮祥等,2008;張致偉等,2013),可作為研究區域應力水平較低的佐證。

結合自貢和長寧地區地震活動及加壓、注水數據,分析研究區域地震視應力在不同階段的變化特征 (圖7)。自貢地區除發生在第Ⅱ時段內的14、15和17日ML≥4.5地震視應力比較顯著外,加壓時段 (Ⅰ)與停止加壓時段 (Ⅱ)內同等震級檔的地震視應力相當 (圖7a);長寧地區與上述現象相似,除第Ⅱ時段內的13和17日ML≥4.5地震外,整個時段內同等震級的地震視應力相差不大 (圖7b)。可見自貢和長寧地區地震視應力同樣依賴于震級大小,加壓注水等工業活動對該地區的地震視應力影響不大。

4 結論與討論

本文以四川自貢和長寧兩個注水地區為例,基于區域臺網記錄的數字地震波形資料和觀測報告,分析注水誘發地震的譜振幅相關系數特征及其與當地加壓注水數據的關系,并結合表征絕對應力水平的視應力綜合分析研究區域的應力高低。

自貢和長寧地區注水誘發地震的譜振幅相關系數均在高值波動,且變化過程與加壓注水數據具有很好的相關性。自貢地區地震譜振幅相關系數在0.84~0.92之間波動,其平均值為0.89,在注水壓力持續升高的背景下,譜振幅相關系數呈現隨加壓數據的升高而增大的特征;長寧地區地震譜振幅相關系數在0.87~0.97之間波動,其平均值為0.90,在注水量和出水量之差大幅上升時段,相應的譜振幅相關系數也維持在高值。而在研究區域注 (出)水量相當或停止加壓時,譜振幅相關系數也有相應降低,但依然在高值波動。分析認為加壓注水對局部區域應力場產生擾動影響,使地下淺層裂隙呈優勢取向排列,引發的中小地震震源機制表現出較好的相似性,從而表現為較高的譜振幅相關系數。

川東南地區地震視應力存在局部區域差異,自貢和長寧地區同等震級的地震視應力與其它區域結果相比相對偏低,可能反映了上述區域的應力水平相對偏低,研究區域的地震視應力具有隨震級增大而增大的特征,加壓注水等工業活動對所在區域地震視應力影響不大。

自貢和長寧兩個注水地區地震的譜振幅相關系數整體偏高,反映了中小地震震源機制具有很好的相似性,而視應力結果顯示研究區域的應力水平并不高,且歷史上該區域的地震活動水平也不高。崔子健等 (2012)研究發現強震前的小震群同樣具有較高的譜振幅相關系數,其反映了大震孕育過程中震源區介質應力水平增強。注水地區及其它強震孕育區雖然譜振幅相關系數具有相同的特征,但其孕震機理不同,局部區域應力水平也就存在差異。

鑒于小震震源機制不易求解,理論上譜振幅相關系數能夠反映震源機制的相似性。但小震活動具有隨機性,其受局部應力場及局部構造等因素的影響,應結合表征應力水平的參數進行綜合判定。

本項工作得到了中國地震局“地震分析預測研究青年工作組”的大力支持;中國地震臺網中心周龍泉博士和四川省地震局程萬正研究員給予了重要指導和有益討論;中國地震局預測研究所崔子健和云南省地震局鄔成棟在程序上給予了幫助和指導,在此一并致謝!

陳颙.1978.用震源機制一致性作為描述地震活動性的新參數[J].地球物理學報,21(2):146-159.

程萬正,陳學忠,喬慧珍.2006a.四川地震輻射能量和視應力的研究[J].地球物理學進展,21(3):692-699.

程萬正,阮祥,張永久.2006b.川滇次級地塊震源機制類型與一致性參數[J].地震學報,28(6):561-573.

崔子健,李志雄,陳章立,等.2012.判別小震群序列類型的新方法研究——譜振幅相關分析法[J].地球物理學報,55(5):1718-1724.

刁桂苓,徐錫偉,陳于高,等.2011.汶川 MW7.9和集集MW7.6地震前應力場轉換現象及其可能的前兆意義[J].地球物理學報,54(1):128-136.

刁桂苓,于利民,李欽祖.1994.強震前后震源區應力場變化一例[J].地震學報,16(1):64-69.

丁健民,高莉青.1981.地殼水平應力與垂直應力隨深度的變化[J].地震,2:46-48.

李方全,祁英男.1988.地殼應力隨深度的變化規律[J].巖石力學與工程學報,7(4):301-309.

劉杰,鄭斯華,黃玉龍,等.2003.利用遺傳算法反演非彈性衰減系數、震源參數和場地響應[J].地震學報,25(2):211-218.

喬慧珍,張永久,程萬正,等.2006.川西北地區介質衰減特征研究[J].地震地磁觀測與研究,27(4):1-7.

阮祥,程萬正,張永久,等.2008.四川長寧鹽礦井注水誘發地震研究[J].中國地震,24(3):226-234.

吳忠良,黃靜,林碧蒼.2002.中國西部地震視應力的分布[J].地震學報,24(3):293-301.

張致偉,程萬正,梁明劍,等.2012.四川自貢—隆昌地區注水誘發地震研究[J].地球物理學報,55(5):1635-1645.

張致偉,程萬正,吳朋,等.2013.自貢—隆昌地區地震重新定位及P波速度結構研究[J].中國地震,29(1):37-47.

周連慶,趙翠萍,修濟剛,等.2008.利用天然地震研究地殼Q值的方法和進展[J].國際地震動態,(2):1-11.

朱航,劉杰,陳天長.2006.采用體波譜振幅相關系數方法研究地震序列的震源機制變化過程[J].地震,26(2):1-11.

Atkinson G.M.,Mereu R.F..1992.The Shape of Ground Motion Attenuation Curves in Southeastern Canada[J].Bull.Seism.Soc.Amer.,82(5):2014-2031.

Brune J.N..1970.Tectonic Stress and the Spectra of Seismic Shear Waves from Earthquakes[J].J.Geophys.Res.,75(26):4997-5009.

Hauksson E..1994.State of Stress from Focal Mechanisms before and after the 1992 Landers Earthquake Sequence[J].Bull.Seism.Soc.Am,84(3):917-934.

Kanamori H.,Anderson D.L..1975.Theoretical Basis of Some Empirical Relations in Seismology[J].Bull.Seism.Soc.Am.,65(5):1073-1095.

Kanamori H..1977.The Energy Release in Great Earthquakes[J].J.Geophys.Res.,82(20):2981-2987.

Lund B.,B?evarsson R..2002.Correlation of Microearthquake Bodywave Spectral Amplitudes[J].Bull.Seism.Soc.Am.,92(6):2419-2433.

Michael A.J.,Ellsworth W.L.,Oppenheimer D.H..1990.Coseismic Stress Changes Induced by the 1989 Loma Prieta,California Earthquake[J].Geophys.Res.Lett.,17(9):1441-1444.

Moya C.A.,Aguirre J.,Irikura K..2000.Inversion of Source Parameters and Site Effects from Strong Ground Motion Records using Genetic Algorithms[J].Bull.Seism.Soc.Amer.,90(4):977-992.

Seeber L.,Armbruster J.G..2000.Earthquakes as Beacons of Stress Change[J].Nature,407(6800):69-72.

Wiemer S.,Gerstenberger M.,Hauksson E..2002.Properties of the Aftershock Sequence of the 1999 MW7.1 Hector Mine earthquake:Implications for Aftershock Hazard[J].Bull.Seism.S.Am.,92(4):1227-1240.

Wu Z.L..2001.Scaling of Apparent Stress from Broadband Radiated Energy Catalogue and Seismic Moment Catalogue and its Focal Mechanism Dependence[J].Earth.Planets.Space.,53:943-948.

Wyss M.,Brune J.N..1968.Seismic Moment,Stress,and Source Dimensions forEarthquakes in theCalifornia-NevadaRegion[J].J.Gephys.Res.,73(14):4681-4694.

猜你喜歡
機制區域研究
FMS與YBT相關性的實證研究
遼代千人邑研究述論
視錯覺在平面設計中的應用與研究
科技傳播(2019年22期)2020-01-14 03:06:54
EMA伺服控制系統研究
自制力是一種很好的篩選機制
文苑(2018年21期)2018-11-09 01:23:06
關于四色猜想
分區域
破除舊機制要分步推進
中國衛生(2015年9期)2015-11-10 03:11:12
基于嚴重區域的多PCC點暫降頻次估計
電測與儀表(2015年5期)2015-04-09 11:30:52
注重機制的相互配合
中國衛生(2014年3期)2014-11-12 13:18:12
主站蜘蛛池模板: 欧美啪啪网| 国产18页| 丁香综合在线| 成人午夜视频网站| 欧美精品二区| 狠狠做深爱婷婷久久一区| 久久无码免费束人妻| 国产伦精品一区二区三区视频优播 | 欧美视频二区| 91视频青青草| 国产av一码二码三码无码| 欧美性天天| 欧美α片免费观看| 国产成人在线无码免费视频| 高清无码手机在线观看 | 国产精品v欧美| 99精品在线视频观看| 黄色网在线免费观看| 欧美日韩va| yjizz视频最新网站在线| 国产成人精品一区二区三在线观看| 91热爆在线| 欧美成人综合在线| 中文字幕无码中文字幕有码在线| 成人在线亚洲| 亚洲国产中文精品va在线播放| 2021国产精品自产拍在线| 亚洲不卡无码av中文字幕| 蝴蝶伊人久久中文娱乐网| 狠狠综合久久| 国产青榴视频| 黄片在线永久| 97精品伊人久久大香线蕉| 国产精品入口麻豆| 中文字幕1区2区| 91丝袜在线观看| 亚洲视频三级| 亚洲欧美自拍中文| 99精品免费在线| 无码AV动漫| 亚洲精品片911| a在线观看免费| 国产精品爽爽va在线无码观看 | 中文字幕久久亚洲一区| 国产精品欧美日本韩免费一区二区三区不卡| 日韩AV无码一区| 最近最新中文字幕免费的一页| 亚洲精品第一在线观看视频| 成人蜜桃网| 久久精品人人做人人爽| 91一级片| 国产福利在线免费| 亚洲全网成人资源在线观看| 国产成人精品综合| 99热这里都是国产精品| 无码国内精品人妻少妇蜜桃视频| 波多野结衣第一页| 免费播放毛片| 国产免费久久精品99re丫丫一| 青青热久免费精品视频6| 在线观看免费黄色网址| 91九色国产在线| 久久久噜噜噜| 欧美精品xx| 国产成人综合网| 中文字幕在线观看日本| 欧美一级高清片久久99| 亚洲精品成人片在线观看| 中文天堂在线视频| 日韩欧美国产精品| 日韩精品成人在线| 亚洲天堂网在线观看视频| 国内精品视频在线| 欧美色视频网站| 韩日无码在线不卡| 日韩在线网址| 国产自在线播放| 国产性生交xxxxx免费| 91探花在线观看国产最新| 国产黄色片在线看| 熟妇无码人妻| 天天做天天爱天天爽综合区|