亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        基于Barton-Bandis非線性破壞準則的邊坡破壞概率分析

        2020-12-21 04:46:58陳歡歡劉子航
        金屬礦山 2020年11期
        關鍵詞:準則巖體邊坡

        陳歡歡 雍 睿 劉子航

        (紹興文理學院土木工程學院,浙江紹興 312000)

        邊坡穩(wěn)定性問題涉及到公路工程、水利工程、建筑工程等諸多領域,采用合理有效的邊坡穩(wěn)定性分析方法準確評價邊坡穩(wěn)定性,對于保證人民生命財產(chǎn)安全,確保工程安全有著重要的現(xiàn)實意義。

        目前邊坡穩(wěn)定性評價的常用方法可大致分為兩類,即定性分析方法和定量分析方法。定性分析方法包括自然(成因)歷史分析法、工程地質(zhì)類比法和圖解法等,定量分析方法包括極限平衡分析法(如Bishop 法[1]、Janbu 法[2]、Morgenstern-Prince 法[3]、Spencer法[4]、剩余推力法[5]等)和數(shù)值分析法(如有限元法[6]、離散元法[7]、塊體理論與不連續(xù)變形分析[8]等)。由于傳統(tǒng)評價方法存在一定的局限性,國內(nèi)眾多學者在以上兩類方法的基礎上,引入了一些新的理論方法,如不確定性分析法,其中包括可靠性評價法[9]、模糊綜合評價法[10]、灰色系統(tǒng)理論[11]等。目前邊坡穩(wěn)定性分析方法中的理論公式及數(shù)值分析多是基于Mohr-Coulomb準則(簡稱“M-C準則”)。經(jīng)研究發(fā)現(xiàn),巖體內(nèi)部廣泛分布著結構面,其抗剪強度受粗糙度和壁巖強度的影響呈現(xiàn)非線性特征,當法向應力σn較小時,巖體結構面的抗剪性能并不能充分滿足M-C準則[12]。而Barton-Bandis準則(簡稱“B-B準則”)能夠描述巖體內(nèi)部結構面的非線性力學破壞行為。在估算節(jié)理抗剪強度時,B-B準則既能考慮到節(jié)理表面的粗糙起伏特征,又能考慮到應力環(huán)境及巖壁強度,計算結果較為可靠實用。

        根據(jù)可靠性設計理論,邊坡穩(wěn)定評價可用概率來表達,而不必用“絕對穩(wěn)定”或“必然破壞”等結論表示,具體可表示為“具有90%的安全概率”,“90%的破壞可能性”等結論[13]。徐洪[14]以傳統(tǒng)邊坡穩(wěn)定性分析中的穩(wěn)定性系數(shù)為基礎,結合可靠性理論對邊坡穩(wěn)定進行了常規(guī)可靠性分析;陳坤杰[15]通過Bishop條分法進行穩(wěn)定性分析,建立了邊坡的可靠度計算模型,并利用MATLAB編程實現(xiàn)了邊坡穩(wěn)定性的蒙特卡羅模擬計算;張穎[16]運用蒙特卡羅法(Monte Carlo)和Rosenblueth法進行了邊坡可靠性分析,有效彌補了傳統(tǒng)的邊坡穩(wěn)定性評價方法沒有考慮設計參數(shù)變異性的不足;張露[17]將瑞典條分法與Monte Carlo法相結合進行了邊坡可靠度分析;周超海等[18]將穩(wěn)定性系數(shù)考慮為具有一定破壞概率的隨機分布函數(shù),通過蒙特卡羅模擬及Janbu模型對邊坡進行了穩(wěn)定性分析;周春梅等[19]采用Monte Carlo法模擬分析影響邊坡穩(wěn)定性系數(shù)的參數(shù)分布,在剩余推力法的基礎上建立了滑坡破壞概率模型,經(jīng)過敏感性分析確定了影響邊坡破壞概率的主要因素。以上研究在進行邊坡穩(wěn)定性分析時多是將M-C準則作為強度準則。

        為更好地適應巖體結構面抗剪強度所呈現(xiàn)的非線性特征和抗剪強度參數(shù)的離散性特征,本研究基于B-B準則的邊坡破壞概率分析方法,建立了巖質(zhì)邊坡穩(wěn)定性分析模型,對最危險滑動面進行破壞概率分析,并將改進的破壞概率分析方法與傳統(tǒng)極限平衡分析法的計算結果進行對比分析,驗證基于B-B準則的邊坡破壞概率分析方法評價邊坡穩(wěn)定性的有效性。

        1 基于B-B準則的邊坡破壞概率分析法

        基于B-B準則的邊坡破壞概率分析方法是在非線性B-B準則的基礎上選用表征巖體結構面特征的3種抗剪強度參數(shù)及其分布來建立巖質(zhì)邊坡穩(wěn)定性分析模型,結合Monte Carlo法和各參數(shù)分布生成多組隨機數(shù),通過Morgenstern-Price法(簡稱“M-P法”)對邊坡巖體內(nèi)部潛在的滑移面進行穩(wěn)定性分析,得到潛在滑移面上的多組穩(wěn)定性系數(shù),采用Monte Carlo法對多組穩(wěn)定性系數(shù)進行統(tǒng)計分析,確定邊坡穩(wěn)定性系數(shù)均值和破壞概率。

        1.1 結構面抗剪強度參數(shù)選取

        BARTON等[20]等在研究多組結構面的直剪特性和試驗結果的基礎上,提出了用于估計不規(guī)則、無充填結構面峰值抗剪強度的非線性JRC-JCS經(jīng)驗公式:

        式中,τ為結構面抗剪強度,MPa;σn為作用于結構面上的法向應力,MPa;?b為基本摩擦角,(°);JRC為結構面粗糙度系數(shù);JCS為結構面壁巖強度,MPa。

        M-C準則與B-B準則對比如圖1所示。圖1中,c為結構面黏聚力,?為結構面內(nèi)摩擦角,ct為結構面等效黏聚力,?t為結構面等效內(nèi)摩擦角,兩者可由BB準則等效線性擬合M-C準則得到??紤]到M-C準則的局限性,本研究選取非線性JRC-JCS經(jīng)驗公式中的JRC、JCS、?b等作為邊坡巖體結構面的特征參數(shù)進行邊坡建模。

        1.2 穩(wěn)定性系數(shù)計算

        M-P法在進行邊坡穩(wěn)定性分析時不僅考慮了條塊之間力的平衡和力矩的平衡,而且考慮了滑動巖體的邊界條件,可對非圓弧滑動面上可能滑動的巖體中微分條塊列出同時滿足力與力矩平衡條件的微分方程式。假定每一條塊切向條間力T與法向條間力E之間存在一個對水平坐標x的函數(shù)關系:

        式中,λ為條間力比例常數(shù);T為每個條塊的切向條間力,kN;E為每一條塊的法向條間力,kN。

        于是,每一條塊兩側的邊界條件為

        邊坡整體力平衡穩(wěn)定性系數(shù)Ff和力矩平衡穩(wěn)定性系數(shù)Fm可分別進行如下計算[22]:

        式中,c為結構面黏聚力,kPa;?為結構面內(nèi)摩擦角,(°);Wi為滑體各條塊的重力,kN;u為孔隙水壓力,kPa;Ni為條塊底部法向力,kN;D為集中荷載,kN;ω為集中荷載與豎直方向的夾角,(°);αi為滑體底面傾角,(°);R為滑弧圓心到滑弧的距離,m;x為水平坐標;f為摩擦系數(shù);d為集中荷載的力矩,kN·m;β為折算系數(shù);i為滑體各條塊編號,1≤i≤n。

        結合邊界條件及微分方程的迭代計算,求得Ff和Fm隨λ變化的分布曲線,兩分布曲線的交點即為穩(wěn)定性系數(shù)Fs。由于微分計算過程較繁瑣,故使用SLIDE軟件中的M-P法模塊進行計算和分析。

        1.3 破壞概率計算

        Monte Carlo法又稱隨機模擬法或統(tǒng)計試驗法,是邊坡可靠性分析時常采用的一種方法。影響穩(wěn)定性系數(shù)的各種參數(shù)(如JRC、JCS、?b等)在獲取過程中會因外界因素或人為因素的作用而產(chǎn)生隨機性誤差,這些參數(shù)從嚴格意義上可以看作隨機變量。穩(wěn)定性系數(shù)Fs可看作隨機變量JRC、JCS、?b的狀態(tài)函數(shù),其函數(shù)表達式為

        若已知狀態(tài)變量JRC、JCS、?b的概率分布,根據(jù)邊坡的極限狀態(tài)條件Fs=0,利用Monte Carlo法產(chǎn)生符合狀態(tài)變量概率分布的3組隨機 數(shù)(JRC1,JRC2,…,JRCn),(JCS1,JCS2,…,JCSn),,并將其代入狀態(tài)函數(shù)(式(6))計算得到狀態(tài)函數(shù)Fs的y個隨機數(shù)。如果在這y個隨機數(shù)中有x個小于或等于1(以穩(wěn)定性系數(shù)表示邊坡狀態(tài)),當y足夠大時,根據(jù)大數(shù)定律,此時的頻率已近似于概率,因而可得到邊坡的破壞概率:

        2 工程實例

        2.1 礦山概況

        柏楊廟礦山位于貴州省遵義市桐梓縣婁山關鎮(zhèn)沙嘴村,中心點地理坐標為東經(jīng)106°49′55.20″,北緯28°6'28.80″。研究區(qū)位于黔北山地與四川盆地的銜接地帶,冰川作用顯著,溶蝕侵蝕并存,長期的地質(zhì)作用使巖體結構多呈現(xiàn)層狀結構。研究區(qū)地處中亞熱帶高原季風濕潤性氣候區(qū),年平均降雨量1 038.8 mm,夏季降水量多,冬季降水量少??辈槠陂g未發(fā)現(xiàn)該礦山附近存在地表水系,礦山開挖面也未見地下水溢出。該區(qū)不處于地震帶,屬于弱震環(huán)境,根據(jù)中國地震烈度區(qū)劃,區(qū)內(nèi)為Ⅴ度區(qū)。礦山邊坡整體輪廓如圖2所示。該邊坡近NE走向,傾向SE??傮w邊坡高度100 m,寬度145 m,進深100 m,坡向137°,坡角75°,無臺階邊坡。邊坡整體主要出露石灰?guī)r,以發(fā)育層面、節(jié)理為主,邊坡范圍內(nèi)沒有斷層發(fā)育,自然容重為27 kN/m3,飽和容重為28 kN/m3。層面B1產(chǎn)狀132°∠66°,規(guī)模120 m,無填充物。邊坡發(fā)育兩組節(jié)理J1、J2,其中節(jié)理J1產(chǎn)狀54°∠38°,規(guī)模0.2~1 m,間距0.1~0.6 m,無充填物;節(jié)理J2產(chǎn)狀225°∠73°,規(guī)模0.8~2 m,間距0.1~0.4 m,無充填物。經(jīng)調(diào)查,沒有發(fā)現(xiàn)該邊坡巖性、構造、工程地質(zhì)和水文地質(zhì)條件存在明顯變化的情況。

        2.2 滑移面確定

        該礦山邊坡潛在滑移面可采用露天礦山邊坡穩(wěn)定性分級分析法確定[23]。通過對層面B1與總體邊坡坡面進行赤平投影(圖3)可知,層面(貫穿性結構面)控制著邊坡的整體穩(wěn)定性??傮w邊坡坡向137°,坡角75°。層面B1傾向132°,傾角66°,層面傾向與邊坡傾向夾角5°,小于30°,為順坡向層面,且層面傾角(66°)小于邊坡坡角(75°),邊坡整體很有可能沿該結構面發(fā)生單平面型滑移破壞,故需對總體邊坡進行穩(wěn)定性評價。

        2.3 結構面統(tǒng)計測量與數(shù)據(jù)分析

        依據(jù)赤平投影分析結果,在層面B1的作用下,該邊坡整體有可能沿該結構面發(fā)生單平面型滑移而產(chǎn)生破壞,層面B1為潛在滑移面。因此,需要對層面B1開展結構面取樣調(diào)查。

        2.3.1JRC統(tǒng)計測量與取值

        在潛在滑移面所對應的結構面上,使用輪廓曲線儀繪制了106條結構面輪廓曲線(圖4)。采用工程高清掃描儀掃描結構面輪廓曲線圖,將圖中曲線轉換為圖片格式文件,通過計算程序對圖中曲線進行異常點排除、形態(tài)學濾波去噪和圖像歸一化,提取掃描圖片的灰度值矩陣。根據(jù)結構面實測長度與數(shù)字化灰度矩陣的大小關系,將輪廓曲線上各像素點坐標轉化為實際坐標數(shù)據(jù),采用全域搜索法在每條輪廓曲線上提取10 cm的試樣,分別讀取并存儲每條輪廓曲線10 cm試樣的坐標數(shù)據(jù)。本研究通過杜時貴[23-24]所提出的Barton直邊法的簡明公式計算每條結構面試樣的JRC:

        式中:RA為相對幅度;Ln為輪廓曲線長度,cm;Ry為結構面表面輪廓曲線齒凸幅度,cm。

        2.3.2JCS統(tǒng)計測量與取值

        本研究采用回彈法確定JCS值,選用指針直讀式回彈儀在潛在滑移面所對應的結構面上測量了48組回彈值(圖5),依據(jù)回彈值與壁巖強度關系和壁巖狀態(tài)(含水、風化情況),確定結構面壁巖強度大?。?3]。為方便起見,基于MATLAB軟件程序計算JCS值。首先輸入巖石容重、回彈值及結構面真傾角;然后設定曲線選型范圍、容重范圍、坐標轉換參數(shù)及強度插值參數(shù);最后輸入回彈參數(shù)進行坐標轉換,進行強度計算、偏差分析,確定JCS值和飽和狀態(tài)回彈值。

        2.3.3 ?b取值

        將回彈值代入基本摩擦角的經(jīng)驗公式,即理查茲提出的結構面回彈值與基本摩擦角的線性關系式[25],得到?b:

        式中,N為回彈值。

        2.4 抗剪強度參數(shù)統(tǒng)計分布檢驗

        2.4.1 概率密度函數(shù)

        在進行邊坡穩(wěn)定性分析之前,需要驗證參數(shù)的概率分布。各參數(shù)的頻率分布如圖6所示。由圖6可知:各參數(shù)分布圖形與正態(tài)分布函數(shù)圖形相似,可對3組數(shù)據(jù)分別進行正態(tài)分布檢驗。

        2.4.2 K-S檢驗

        由圖6可知:樣本數(shù)據(jù)概率分布的擬合程度與正態(tài)分布相近,但精確檢驗還需要通過數(shù)量方式實現(xiàn),可選擇單樣本K-S檢驗法檢驗樣本數(shù)據(jù)是否服從某種理論分布。檢驗結果如表1所示。正態(tài)分布形式的雙尾概率P(JRC)=0.052>0.05,P(JCS)=0.2>0.05,P(?b)=0.2>0.05,表明三者均接受該形式假設。因此,JRC、JCS、?b的樣本數(shù)據(jù)均符合正態(tài)分布特征的假設。

        2.5 邊坡穩(wěn)定性分析

        2.5.1 基于B-B準則的邊坡破壞概率分析

        2.5.1.1 邊坡模型及荷載設定

        按照礦山邊坡所對應的原型幾何尺寸建立的分析模型如圖7所示。該邊坡為自然邊坡,根據(jù)《建筑抗震設計規(guī)范》(GB 50011—2010)[26]要求,可不考慮地震作用,邊坡表面沒有任何外部荷載,周圍沒有明顯的水文地質(zhì)條件,故模型表面不需設置荷載,可不考慮地下水的滲流作用,但需考慮降雨的影響。降雨工況下,巖體內(nèi)部水、氣兩相壓力變化導致巖體結構面抗剪強度參數(shù)發(fā)生變化,故需對降雨工況下的抗剪強度參數(shù)作適當折減處理[27]。

        2.5.1.2 強度模型及參數(shù)設置

        本研究選用B-B準則作為穩(wěn)定性計算的強度準則,所采用的計算參數(shù)為JRC、JCS、?b,將三者在干燥狀態(tài)和飽和狀態(tài)下的平均值、標準差、相對最小值及相對最大值分別導入模型,基本數(shù)據(jù)如表2所示。JCS在飽和狀態(tài)下的樣本數(shù)據(jù)可通過軟化系數(shù)k=0.7對干燥狀態(tài)的樣本數(shù)據(jù)進行折減獲得[28],?b在飽和狀態(tài)下的樣本數(shù)據(jù)可根據(jù)2.3.2節(jié)確定的飽和狀態(tài)回彈值代入式(9)計算獲得。

        2.5.1.3 模擬過程

        注:表中各參數(shù)樣本數(shù)據(jù)均符合正態(tài)分布。

        本研究采用M-P法來計算邊坡穩(wěn)定性系數(shù),模擬方法為Monte Carlo法,滑動方向為從右向左,迭代計算次數(shù)小于50次,滑動面類型為非圓弧型滑動。在模型中預先設置好潛在滑動面,根據(jù)參數(shù)統(tǒng)計數(shù)據(jù),采用Monte Carlo法進行隨機抽樣,在考慮計算精度的前提下,選擇模擬次數(shù)為1 000次[29],選用自動搜索方式對每一組參數(shù)隨機樣本值進行最不利滑裂面搜索,分別求得總體邊坡在自然工況、降雨工況下的多組穩(wěn)定性系數(shù),采用Monte Carlo法統(tǒng)計分析穩(wěn)定性系數(shù)值,確定邊坡破壞概率。

        2.5.1.4 計算結果分析

        經(jīng)過模擬計算,得到兩種工況下各1 000組穩(wěn)定性系數(shù)值,對其進行統(tǒng)計分析,得到穩(wěn)定性系數(shù)的頻率分布(圖8),圖中灰色部分為Fs<1。不同工況下的計算結果如圖9所示,圖中PF為邊坡破壞概率。自然工況下,穩(wěn)定性系數(shù)服從正態(tài)分布,其最大值為1.683,最小值為0.619,平均值為1.029,破壞概率為49.424%。降雨工況下,穩(wěn)定性系數(shù)服從正態(tài)分布,其最大值為1.438,最小值為0.561,平均值為0.875,破壞概率為83.852%。根據(jù)穩(wěn)定性等級劃分標準[30]判斷,自然工況下的邊坡處于中等危險狀態(tài),穩(wěn)定等級為3;降雨工況下的邊坡處于高危險狀態(tài),穩(wěn)定等級為2。

        2.5.2 基于傳統(tǒng)極限平衡分析法的邊坡穩(wěn)定分析

        2.5.2.1 基于B-B準則的極限平衡分析

        基于B-B準則的極限平衡分析法進行邊坡穩(wěn)定性分析時選用JRC、JCS、?b的平均值作為邊坡穩(wěn)定性分析參數(shù)。將兩種工況下各參數(shù)的平均值分別導入基于B-B準則的邊坡模型,采用M-P法計算得到兩種工況下的穩(wěn)定性系數(shù)(圖10)。自然工況下,邊坡的穩(wěn)定性系數(shù)為1.029,大于1小于1.05;降雨工況下,邊坡的穩(wěn)定性系數(shù)為0.875,小于1。根據(jù)《滑坡防治工程勘查規(guī)范》(GBT 32864—2016)[31],該邊坡整體在自然工況下處于欠穩(wěn)定狀態(tài),在降雨工況下處于不穩(wěn)定狀態(tài)。

        2.5.2.2 基于M-C準則的極限平衡分析

        (1)法向應力取值。由野外地質(zhì)調(diào)查可知,桐梓縣地區(qū)地質(zhì)構造穩(wěn)定。因此,可按式(10)來計算邊坡垂直方向主應力和水平方向主應力[32]:

        式中,H為滑體高度,m;γ為結構面上覆巖石容重,(×10-3MN/m3);σV為結構面所處部位垂直主應力,MPa;σH為結構面所處部位水平主應力,MPa。

        在已知垂直主應力和水平主應力的情況下,可通過矢量法計算獲得結構面上的法向應力,選取5個試算點1/3σn、2/3σn、σn、4/3σn和5/3σn。

        (2)穩(wěn)定性系數(shù)計算。將試算的5個法向應力值和JRC、JCS、?b的平均值分別代入式(1)得到多組抗剪強度值,結果如表3所示。根據(jù)記錄結果,以最小二乘法進行線性擬合,結果如圖11所示。由直線的傾角確定巖體結構面自然工況下的峰值內(nèi)摩擦角?=40.59°,降雨工況下的峰值內(nèi)摩擦角?=35.42°;由直線的截距確定結構面自然工況下黏聚力c=0.541 2 MPa,降雨工況下黏聚力c=0.456 8 MPa。將兩種工況下的峰值內(nèi)摩擦角和黏聚力分別導入基于M-C準則的邊坡模型,采用M-P法計算得到兩種工況下的穩(wěn)定性系數(shù)。

        (3)計算結果分析。極限平衡分析法的計算結果如圖12所示。自然工況下,邊坡的穩(wěn)定性系數(shù)為1.037,大于1小于1.05;降雨工況下,邊坡的穩(wěn)定性系數(shù)為0.884,小于1。根據(jù)《滑坡防治工程勘查規(guī)范》(GBT 32864—2016)[31],該邊坡整體在自然工況下處于欠穩(wěn)定狀態(tài),在降雨工況下處于不穩(wěn)定狀態(tài)。

        2.5.3 計算結果對比分析

        通過上述兩類方法對巖質(zhì)邊坡進行穩(wěn)定性分析,發(fā)現(xiàn)其計算結果與所采用的強度準則密切相關,在計算方法相同的前提下,不同的強度準則可以求解出不同的穩(wěn)定性系數(shù)。傳統(tǒng)極限平衡分析法以線性M-C準則作為強度模型進行邊坡穩(wěn)定性分析時,天然節(jié)理的粗糙度和壁巖強度的影響使M-C準則不再適用于巖體結構面非線性破壞的情況,故計算結果存在一定的誤差;以非線性B-B準則作為強度模型時,盡管適用于巖體結構面抗剪強度的非線性特征,但僅以穩(wěn)定性系數(shù)來評價邊坡穩(wěn)定性不夠精確?;贐-B準則的破壞概率分析方法以非線性B-B準則作為強度模型,考慮到了巖體結構面抗剪強度的非線性特征,所采用的抗剪強度參數(shù)獲取方便,以多種隨機變量的統(tǒng)計分布特征為基礎進行建模,避免了抗剪強度參數(shù)離散性的影響,計算結果以穩(wěn)定性系數(shù)和破壞概率形式來表達可靠性較強。

        3 結 論

        本研究將可靠度理論和Barton-Bandis準則相結合,根據(jù)巖質(zhì)邊坡破壞極限狀態(tài)函數(shù),討論了以破壞概率和穩(wěn)定性系數(shù)表征邊坡臨界狀態(tài)的邊坡穩(wěn)定性分析方法,指出了傳統(tǒng)的極限平衡分析法在進行邊坡穩(wěn)定性分析時所存在的不足,結合工程實例分析得出以下結論:

        (1)進行邊坡穩(wěn)定性分析時,基于B-B準則的破壞概率分析方法選用可表征巖體結構面特征的3種抗剪強度參數(shù)JRC、JCS、?b,并充分考慮了它們的統(tǒng)計分布特征,而傳統(tǒng)極限平衡分析法選用抗剪強度參數(shù)均值來表征巖體物理力學指標,忽略了巖體結構面特征參數(shù)的離散性。

        (2)傳統(tǒng)的極限平衡分析法采用“穩(wěn)定性系數(shù)小于1”作為邊坡破壞的臨界狀態(tài)不夠精確;基于B-B準則的邊坡破壞概率分析方法以穩(wěn)定性系數(shù)和破壞概率形式來表達邊坡的穩(wěn)定程度,可為現(xiàn)場施工提供可靠的理論依據(jù)。

        猜你喜歡
        準則巖體邊坡
        這邊山 那邊坡
        具非線性中立項的二階延遲微分方程的Philos型準則
        基于無人機影像的巖體結構面粗糙度獲取
        甘肅科技(2020年20期)2020-04-13 00:30:18
        水利水電工程高邊坡的治理與加固探討
        基于Canny振蕩抑制準則的改進匹配濾波器
        基于SLOPE/W的邊坡穩(wěn)定分析
        平泉縣下營坊雜巖體分異演化及其成巖成礦
        基于不同軟件對高邊坡穩(wěn)定計算
        一圖讀懂《中國共產(chǎn)黨廉潔自律準則》
        單一層狀巖體和軟硬復合巖體單軸壓縮破損特征試驗研究
        国产成人精品蜜芽视频| 天天影视性色香欲综合网| 亚洲欧美综合在线天堂| 国产91在线精品福利| av免费在线播放一区二区| 未满十八18禁止免费无码网站| 人妻丰满熟妇av无码区hd| 无码一区二区三区AV免费换脸| 国产女主播视频一区二区三区| 日韩极品视频免费观看| 巨胸喷奶水www视频网站| 538在线啪在线观看| 日韩av免费在线不卡一区| 有坂深雪中文字幕亚洲中文| 欧美亚洲熟妇一区二区三区| 国产精品一区二区在线观看99| 亚洲精品中文有码字幕| 人成综合视频在线播放| 狠狠做深爱婷婷久久综合一区| 国产精品亚洲A∨天堂| 亚洲精品乱码久久麻豆| 免费av一区二区三区| 亚洲乱亚洲乱少妇无码99p| 免费观看久久精品日本视频| 国内精品毛片av在线播放| 亚洲av成人无遮挡网站在线观看| 法国啄木乌av片在线播放| 久久久精品国产亚洲麻色欲| 国产一区二区三区色哟哟| 人人妻人人爽人人澡欧美一区| 国产亚洲欧美在线观看的| 蜜臀一区二区av天堂| 久久成人国产精品一区二区| 国产一区二区三区在线观看免费 | 日日摸夜夜添夜夜添无码免费视频| 免费人成视频网站在线| 青青草免费手机视频在线观看| 国产精品无码久久久久| 亚洲午夜久久久久中文字幕久| 日韩av在线亚洲女同| 日韩av无码精品一二三区|