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

水熱變化對季節(jié)凍土區(qū)牽引式滑坡變形的作用機制
——以2020年8月11日西和縣大型滑坡為例

2023-01-15 06:13:02梁小鵬陳虹舉燕強珍楊佳亮馬敏
科學技術(shù)與工程 2022年33期
關(guān)鍵詞:變形模型

梁小鵬, 陳虹舉, 燕強珍, 楊佳亮, 馬敏

(1. 甘肅工程地質(zhì)研究院, 蘭州 730030; 2. 中國科學院西北生態(tài)環(huán)境資源研究院冰凍圈科學國家重點實驗室, 蘭州 730000; 3. 東北林業(yè)大學寒區(qū)科學與工程研究院, 哈爾濱 150040)

凍融引發(fā)的滑坡是季節(jié)凍土區(qū)滑坡的一種主要形式,該類型滑坡的主要特點是滑動速率較小且變形具有季節(jié)特性[1-3]。其主要變形機制是季節(jié)性凍融引發(fā)的水分遷移導致邊坡下部產(chǎn)生了水分集聚,從而在集聚區(qū)附近形成軟弱層[4-5]。在外部荷載或者土體自重作用下,軟弱層可能會演化成滑動面,并最終產(chǎn)生滑坡[6-7],這是目前季節(jié)凍土區(qū)滑坡發(fā)生機制的主要理論成果。然而,該理論僅僅來自定性分析[8],并無定量的解釋和說明,因此難以指導工程實踐。

要對該類型滑坡進行定量分析,必須要建立相應的水熱力耦合數(shù)值模型[9-11]。目前,飽和狀態(tài)下凍融土的水熱耦合、水熱力耦合模型已經(jīng)取得了一定的成果,并且主要用于各種工程建設的工后評價,如青藏公路和青藏鐵路多年凍土融化后的路基沉降問題等[12-16]。其中,Zhang等[8]運用飽和狀態(tài)下凍土的水熱模型計算了滑坡發(fā)生過程中的水熱變化,然后運用單獨的本構(gòu)模型對該滑坡變形進行了模擬,進而對該類型滑坡的變形機制做了初步解釋。然而,該模型假定土體飽和,并且運用兩個單獨的模型分別模擬了水熱和變形過程,因此無法反映水熱變化對變形的作用機制。即便如此,也極少有研究運用該飽和土模型對凍融引發(fā)的滑坡進行數(shù)值分析[4,6,17-18]。主要原因在于飽和狀態(tài)的假設太過于絕對,不符合凍土中水分遷移的實際情況[19-22]。另外,目前對各類型的滑坡來說,數(shù)值研究的主要作用都是對滑坡災害的評價,極少研究其誘發(fā)機制,這也源于理論模型的欠缺[9,23-24]。

由于缺少完整的適用于滑坡變形的非飽和凍融土的水熱力耦合數(shù)值模型,結(jié)果導致該類型滑坡的變形機制只能停留在定性分析層面。比如,一種解釋認為季節(jié)凍融作用會導致邊坡產(chǎn)生季節(jié)性凍結(jié)滯水促滑效應,這種作用會從3個方面影響邊坡穩(wěn)定性:一是凍結(jié)會導致邊坡下部水分集聚;二是水分增加會降低土體強度;三是水分集聚還會增大邊坡的靜、動水壓力[25]。上述分析從概念模型的角度解釋了季節(jié)凍土區(qū)滑坡產(chǎn)生的原因,具有一定的合理性。但是,水分如何作用于滑體、水分的來源和消散過程并不清楚[26-27]。此外,水分富集的位置及其原因,水分遷移如何影響變形等原因尚不清楚,這也導致現(xiàn)有的邊坡處置措施一直未能起到良好的效果[26,28-31]。如果能夠建立一種符合季節(jié)凍土區(qū)邊坡水熱和變形相互作用過程的理論模型,并結(jié)合上述定性分析結(jié)果,對凍融作用引發(fā)滑坡的過程、變形特點和機制進行定量研究,那么可以對定性分析的結(jié)果進行驗證和補充,從而為該類型滑坡處置措施的探索提供更加充分的理論指導。

為此,現(xiàn)從微觀角度首次基于熱力學觀點,建立季節(jié)凍土區(qū)非飽和土水分、溫度和變形相互耦合的數(shù)值模型,闡述水分對邊坡變形的作用機制。該模型的優(yōu)點是能夠反映凍土邊坡中凍融作用引發(fā)的水分遷移過程和水分不均勻分布對土體強度的影響,以此來分析邊坡的變形特性。并以2020年8月11日甘肅隴南市白馮村發(fā)生的大型滑坡H3為例,分析水熱變化對季節(jié)凍土區(qū)牽引式滑坡的作用機制。對滑坡處置前、后邊坡的溫度場、水分場和變形場進行計算,并對滑坡產(chǎn)生的原因和處置前后邊坡的變形特點進行分析。最后討論目前滑坡治理存在的問題,并結(jié)合研究結(jié)果給出該類型滑坡可能的治理方式。

1 模型與方法

水分作為滑坡產(chǎn)生的主要原因,水分的來源和去向至關(guān)重要。季節(jié)凍結(jié)會導致水分向凍結(jié)鋒面遷移從而形成水分聚集。當冰透鏡體融化后水分液態(tài)水增加必然導致活動層土體強度大幅降低,從而產(chǎn)生蠕動。由此,基于質(zhì)能平衡的非等溫水熱氣耦合模型,建立了季節(jié)凍土區(qū)邊坡的水分遷移方程如下。

(1)

式(1)中:θw、θv和θi分別為液態(tài)水、氣態(tài)水和冰體積含量,如無特殊標注,下標w、i、v和s分別指的是液態(tài)水、冰、氣態(tài)水和土的性質(zhì);ρw、ρv和ρi分別為液態(tài)水、氣態(tài)水和冰的密度;h為水頭高度;Kwh和KwT分別為基質(zhì)勢和溫度梯度引起的等溫和非等溫水力傳導率;Kvh和KvT分別為等溫和非等溫氣態(tài)水傳導率。

根據(jù)Van Genuchten模型,非飽和土的土水特征曲線表達式為

Θ=[1+(-αh)n]-m

(2)

式(2)中:Θ為有效飽和度;α、n和m為擬合參數(shù)。

等溫條件下,非飽和土滲透系數(shù)[32]為

Kwh=10-ΩθiKsΘl[1-(1-Θ1/m)m]2

(3)

式(3)中:Ω為量綱為1的經(jīng)驗參數(shù);l為Mualem模型調(diào)整參數(shù);Ks為飽和滲透系數(shù)。

非等溫液態(tài)水滲透系數(shù)、溫度梯度引起的氣態(tài)水滲透系數(shù)和基質(zhì)勢梯度引起的氣態(tài)水滲透系數(shù)分別為KwT、KvT和Kvh,表達式分別為

(4)

式(4)中:GwT為評價土水特征曲線受溫度影響的經(jīng)驗參數(shù);γ為表面張力,γ0為25℃時的表面張力;M和g分別為水的摩爾質(zhì)量和重力加速度;R為氣體常數(shù);Hr為相對濕度;η為水氣擴散增強因子;D為水氣在土中的擴散度。

溫度梯度作為水分遷移的另一個主要驅(qū)動力,凍融產(chǎn)生的溫度梯度會導致水分向活動層遷移,從而產(chǎn)生水分聚集。根據(jù)能量守恒,傳熱方程[32-33]為

(5)

式(5)中:C為土體的等效比熱;λ為土體等效導熱系數(shù);Lv為水的蒸發(fā)潛熱。其余模型參數(shù)表達式見表1。

表1 模型相關(guān)參數(shù)表達式Table 1 Expression of relevant parameters of the model

牽引式滑坡的變形主要為蠕變變形,本構(gòu)方程必須考慮水分變化的影響,這也是本模型的一個主要創(chuàng)新點。基于黏塑性理論,建立季節(jié)凍土區(qū)邊坡的本構(gòu)方程為

Δσ=DT(Δε-Δεvp-Δεv)

(6)

(7)

式中:DT為與溫度和含水量相關(guān)的彈性矩陣;Δε為總應變增量矢量;Δεvp為黏塑性應變增量矢量;Δεv為水分相變引起的應變增量矢量;ET為彈性模量;vT為泊松比。其他主要力學參數(shù)表達式如下。

(8)

(9)

(10)

(11)

式中:Tp為土體凍結(jié)溫度,273.15 K;其余參數(shù)均為試驗擬合結(jié)果。

式(1)~式(11)組成了非飽和凍土水-熱-力耦合數(shù)值模型。

2 案例分析——以2020年8月11日西和縣大型滑坡為例

2.1 研究區(qū)地質(zhì)概況

研究區(qū)位于隴南市西和縣,該地區(qū)處于兩組大斷裂帶的交匯處,其中一組斷裂帶經(jīng)過成縣-西和-禮縣,橫穿西禮盆地南部,呈北西走向,另一組斷裂帶走向北東,即為西禮斷陷的東緣主斷裂,該斷裂帶向南西方向延伸,與武都化馬一帶的北東向斷裂相接。從這兩組大斷裂對新生代地層的控制和對現(xiàn)代地形地貌的影響可以推斷,它們均屬于活動性斷裂構(gòu)造。主要地貌有構(gòu)造侵蝕低中山地貌和侵蝕堆積河谷地貌。該區(qū)域受溫帶大陸性季風氣候影響,年平均氣溫為8.4 ℃,標準凍深0.42 m,最大凍深約為1.0 m。

2.2 滑坡介紹

滑坡發(fā)生后,甘肅工程地質(zhì)研究院立即對該滑坡進行全面勘測,查明了新老滑坡周界(圖1),鉆探確定了滑動面深度,對巖土體結(jié)構(gòu)特性進行了實驗分析。經(jīng)研究發(fā)現(xiàn),H3為老滑坡,該滑坡自2019年8月2日以來一直處于極慢速滑動狀態(tài)。本次因連續(xù)強降雨影響,于2020年8月11日在H3滑坡體上又發(fā)育1處次級滑坡H3-1,該滑坡截止到2020年8月26日一直處于蠕變滑動狀態(tài)(表2)。

表2 模型相關(guān)參數(shù)表達式Table 2 Expression of relevant parameters of the model

圖1 白馮村滑坡全貌圖Fig.1 Overall view of Baifeng Village landslide

H3老滑坡(東經(jīng)105°16′50″~105°17′44″,北緯33°59′54″~34°00′21″)長約390 m,寬約320 m,滑坡前緣高程為1 604 m,滑坡后緣高程為1 729 m,相對高差為125 m,滑向為140°,鉆孔顯示滑體平均厚度為11.5 m,滑坡體積143.52×104m3,屬于大型滑坡。

H3-1滑坡長約132 m,寬約125 m,滑坡前緣高程為1 603 m,滑坡后緣高程為1 646 m,相對高差為43 m,滑向為139°,根據(jù)鉆孔顯示滑體平均厚度為9.2 m,滑坡體積15.18×104m3,屬于小型滑坡。

2019年8月2日—2020年8月11日期間該邊坡一直處于蠕變滑動狀態(tài),H3-1滑坡中部發(fā)育多條橫向張拉裂縫,截止到2020年8月26日,H3滑坡導致白馮村三社9戶居民75間房屋倒塌,136間房屋出現(xiàn)不同程度豎向裂縫,通村硬化道路中斷損壞約240 m,經(jīng)濟損失約460萬元,災情等級為中型。H3滑坡威脅西和縣白村民三社69戶村民361人的生命財產(chǎn)安全,同時威脅十天高速西和隧道。

2020年10月28日,甘肅省自然資源廳立項對該滑坡進行治理,具體措施包括:在H3滑坡前部布設29根2 m×3 m的C30砼抗滑樁,中心間距5 m,抗滑樁走向垂直于滑坡滑動方向,局部根據(jù)地形調(diào)整,同時輔已截排水和局部擋墻加固措施。

2.3 幾何模型

為了分析滑坡治理前后邊坡的穩(wěn)定性,選擇1-1斷面進行分析(圖2),該斷面的幾何尺寸如圖3所示。鉆孔數(shù)據(jù)顯示,邊坡淺層為厚度變化較大的馬蘭黃土層,土質(zhì)疏松、顆粒較均勻,以粉砂為主,呈塊狀,大孔隙顯著,垂直節(jié)理發(fā)育,一般具有濕陷性。馬蘭黃土下部為弱風化泥巖層,巖芯呈短柱狀,采取率高,干鉆困難,遇水易軟化,風干易崩解。根據(jù)現(xiàn)場勘查和鉆孔取樣分析結(jié)果,勘察單位給出了滑移線的大體位置(圖3)。

圖2 滑坡平面布置圖Fig.2 Landslide layout plan

AB和BC邊界為絕熱、不可滲透邊界;AED為水分可滲透邊界圖3 1-1斷面幾何模型及計算域劃分Fig.3 Geometric model of the section 1-1 and division of calculation domain

2.4 土層參數(shù)

根據(jù)原位檢測和室內(nèi)試驗結(jié)果,土層的物理參數(shù)見表3[34-36],力學參數(shù)見表4[37-39]。

表3 土層物理參數(shù)Table 3 Physical parameters of the soil layers

表4 土層力學參數(shù)Table 4 Mechanical parameters of the soil layers

2.5 邊界條件和初始條件

邊坡溫度邊界條件表達式為

T=T0+Asin(2πt+α0)

(12)

式(12)中:α0為初始相位角,由滑坡開始滑動時間決定。年平均地溫T0為10.5 ℃,年平均氣溫變化幅度為14.7 ℃。AED邊界的水平位移和豎向位移均不受限制,AB和CD邊界的豎向位移不受限制,其他邊界及其位移均受限制。

假設2019年8月2日之前邊坡的水、熱和變形狀況穩(wěn)定,根據(jù)邊坡初始邊界條件對邊坡進行穩(wěn)定性計算,當水熱和變形不再發(fā)生變化,此時邊坡的水熱狀況將作為模型起算的初始條件。

3 結(jié)果和討論

由于新滑坡是在老滑坡的基礎上繼續(xù)滑動的結(jié)果,因此必須對兩者同時進行穩(wěn)定性分析。根據(jù)分析結(jié)果選擇在滑坡前緣布設抗滑樁的方式來控制邊坡的變形。

3.1 滑坡處理前邊坡穩(wěn)定性分析

2019年8月2日—2020年8月11日期間該邊坡一直處于極慢速滑動狀態(tài)。圖4(a)~圖4(d)為2020年8月11日該邊坡的溫度、水分和變形狀況。由圖4(a)可見,推測滑移線以上溫度梯度較大,推測滑移線以下溫度梯度快速減小,溫度梯度突變位置與推測滑移線位置接近。這是因為溫度梯度作為水分遷移的主要驅(qū)動力,其突變位置的水分含量必然會產(chǎn)生明顯變化。圖4(b)中推測滑移線附近的含水率最大值約為0.37,大于周圍土體的含水率,這說明水分在此處發(fā)生了集聚(滯水層)。由于含水率越大,泥巖的強度越低,越容易發(fā)生滑移破壞。因此,含水率最大的位置與推測滑移線的位置相吻合。圖4(c)、圖4(d)為邊坡垂直變形和水平變形狀況,其中最大垂直變形為72.41 cm,最大水平變形為68.57 cm。其中,藍色線代表0-位移線與地表的距離,垂直變形和水平變形的0-位移線距離地表的距離分別為11.11 m和11.02 m,其中推測滑移線距離地表的距離為10.67 m。另外,根據(jù)現(xiàn)場勘查資料顯示,滑坡前緣垂直位移大約為0.4 m,圖4(c)中滑坡前緣最大垂直位移為0.42 m,與實測結(jié)果接近。

圖4 2020年8月11日(滑坡處理前)邊坡的溫度、水分和變形狀況Fig.4 Temperature, moisture and deformation of the slope on August 11, 2020

上述結(jié)果表明水分和變形計算結(jié)果與實際情況相吻合,這也驗證了該數(shù)值模型的準確性。

需要強調(diào)的是,該邊坡一直處于蠕變狀態(tài),主要是由于水分不均勻分布引起的。這與水分在不同季節(jié)的遷移方式不同有關(guān),主要過程包括兩個方面:一方面,在冷季,地表溫度降低導致下層滯水在溫度梯度作用下向活動層內(nèi)遷移并凍結(jié)成冰,造成地表附近體積含冰量和含水量增大[40];另一方面,在暖季,地表溫度上升,活動層內(nèi)冰晶融化導致體積液態(tài)水含量迅速增大。其中小部分液態(tài)水通過地氣交換排放到大氣中,大部分則在溫度梯度和基質(zhì)勢作用下向下遷移并在水分遷移驅(qū)動力最小的地方發(fā)生集聚[35],導致該處土體的彈性模量、黏聚力和內(nèi)摩擦角等力學指標大幅降低,加之在上層土體的重力勢作用下,勢必導致斜坡沿力學強度薄弱位置產(chǎn)生滑移破壞。

3.2 設置抗滑樁后邊坡穩(wěn)定性分析

針對該滑坡的特點,選擇在滑坡前緣設置抗滑樁的方式來減小邊坡變形。為此,對設置抗滑樁后邊坡的變形進行穩(wěn)定性計算,當邊坡達到穩(wěn)定狀態(tài)時,其水熱和變形狀況如圖5所示。

圖5(a)、圖5(b)為設置抗滑樁后邊坡的水熱和變形狀況,與未設置抗滑樁時相比,設置抗滑樁后邊坡的溫度場和水分場略有不同。首先,邊坡整體溫度略微升高,最明顯的是2 ℃等溫線范圍有所增加。這是因為黃土的孔隙度較大,導熱系數(shù)較低;當抗滑樁的阻力有效阻止滑體向下移動時,水分遷移的結(jié)果會使得黃土發(fā)生濕陷下沉而非滑移。此時黃土的密實度增加導致其導熱系數(shù)增大,從而使得進入邊坡內(nèi)的熱量增加。因此,邊坡的溫度出現(xiàn)了一定程度的升高。此外,30%含水率等值線高度略有上升,這是因為黃土層密實度增加致使表層土體的溫度梯度增大,進而使溫度梯度突變點所在高度增加,最終導致暖季水分集聚點位置升高。

圖5(c)、圖5(d)為設置抗滑樁后邊坡的垂直和水平變形,其中,最大垂直變形為2.60 cm,最大水平變形為2.72 cm。與未設置抗滑樁時相比,變形明顯減小。由此可見,設置抗滑樁可以有效提升邊坡穩(wěn)定性。此外,垂直變形和水平變形的0-位移線距離地表的距離分別為9.06 m和8.99 m,和圖4(c)、圖4(d)相比,0-位移線分別提升了2.05 m和2.03 m。由此可見,設置抗滑樁后邊坡的0-位移線更高,這說明抗滑樁還可以提升滑移線高度,從而減小了潛在滑體體積,這對提升邊坡穩(wěn)定性有很大幫助。

圖5 滑坡治理后邊坡的溫度、水分和變形狀況Fig.5 Temperature, moisture and deformation of the slope after landslide treatment

總的來說,抗滑樁主要是通過增大滑體阻力來提升邊坡的力學穩(wěn)定性[41-42]。除此之外,由于水熱和變形之間的相互耦合機制,設置抗滑樁會使得邊坡的地溫升高、最大含水量所在位置提升,最終導致邊坡滑移線升高。這在一定程度上減小了潛在滑體的體積,從而減小了邊坡的垂直和水平變形,進而提高了邊坡穩(wěn)定性。

然而,事實證明通過設置抗滑樁的方式對滑坡的治理效果非常有限。畢竟滑移線始終存在,這也是該類型滑坡復活的主要原因。在外界條件(包括暴雨、地震和人為因素)影響下,這種方式難以杜絕滑坡復活。但是,通過詳細論述水熱變化,尤其是水分不均勻分布對邊坡變形的影響,并對牽引式滑坡產(chǎn)生的機理進行了解釋,發(fā)現(xiàn)可以從減小潛在滑體體積(或者說提升潛在滑移線)的角度提升邊坡穩(wěn)定性,這也為滑坡治理提供了一種新思路,也是下一步要從事的研究。

4 結(jié)論

建立了季節(jié)凍土區(qū)非飽和土的水熱力耦合模型,以2020年8月11日甘肅隴南市白馮村發(fā)生的大型滑坡為研究對象,對滑坡處理前后邊坡中的水熱變化和變形進行了分析,討論了邊坡中水分遷移和水分聚集的原因和來源及其對滑坡變形的影響。得到以下主要結(jié)論。

(1)詳細論述了水熱變化,尤其是水分不均勻分布對邊坡變形的影響,并對牽引式滑坡產(chǎn)生的機理進行了解釋。結(jié)果表明,冷季水分向活動層遷移,在冰透鏡體作用下水分在活動層下部發(fā)生集聚,暖季活動層融化導致液態(tài)水大量增加,導致冰透鏡體所在位置變成潛在滑移面。在外部荷載(暴雨等)作用下滑體會沿著該滑移面產(chǎn)生滑動,這便是季節(jié)凍土區(qū)邊坡中的水分遷移過程及其對滑體滑動的作用機制。

(2)黃土地區(qū)牽引式滑坡復發(fā)的主要原因在于滯水層的存在,而季節(jié)凍融作用加劇了水分不均勻分布和滯水效應。因此,在外界條件(暴雨和地震等)影響下,原本相對穩(wěn)定或者已經(jīng)處于不穩(wěn)定狀態(tài)的老滑坡會出現(xiàn)明顯滑動。

(3)從力學層面來說,抗滑樁能明顯改善滑坡的穩(wěn)定性。滑坡處理前最大垂直變形和最大水平變形分別為72.41 cm和68.57 cm。在滑坡前緣設置抗滑樁后,邊坡達到穩(wěn)定狀態(tài)時的最大垂直變形和最大水平變形分別為2.60 cm 和2.72 cm,穩(wěn)定性顯著提升。

(4)抗滑樁只是從力學角度提升了邊坡穩(wěn)定性,但是潛在滑移線仍然存在,在其他外部外荷載作用下,滑坡可能會復活。發(fā)現(xiàn)可以從減小潛在滑體體積(或者說提升潛在滑移線)的角度提升邊坡穩(wěn)定性,這也為滑坡治理提供了另一種思路。

此外,目前對滑坡的變形機理研究尚未有比較具體的量化研究結(jié)果,尤其是地基中的水分變化非常難確定,從而導致理論研究停滯不前。首次從微觀角度基于熱力學原理分析了水分變化對季節(jié)凍土區(qū)滑坡變形的作用機制,這對滑坡變形的防治和量化研究來說具有理論借鑒意義和實踐指導價值。

猜你喜歡
變形模型
一半模型
重要模型『一線三等角』
談詩的變形
中華詩詞(2020年1期)2020-09-21 09:24:52
重尾非線性自回歸模型自加權(quán)M-估計的漸近分布
“我”的變形計
變形巧算
例談拼圖與整式變形
會變形的餅
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
主站蜘蛛池模板: 国产毛片网站| 午夜日韩久久影院| 在线观看网站国产| 欧美成人第一页| 亚洲男人的天堂久久香蕉| 波多野结衣中文字幕一区二区| 久久五月天综合| 日韩欧美国产综合| 国产精品视频公开费视频| 亚亚洲乱码一二三四区| 一区二区三区四区精品视频 | 91精品国产91久无码网站| a级毛片免费播放| 久久www视频| 日韩高清一区 | 亚洲人成网站日本片| 亚洲成人高清无码| 真实国产精品vr专区| 国产精品午夜福利麻豆| 一区二区三区在线不卡免费| 激情视频综合网| 日韩乱码免费一区二区三区| 亚洲系列中文字幕一区二区| 宅男噜噜噜66国产在线观看| 国产丝袜无码一区二区视频| 亚洲第一色网站| 久久婷婷国产综合尤物精品| 亚洲美女操| 久久成人18免费| 亚洲精品男人天堂| 亚洲精品成人福利在线电影| 国产资源免费观看| 人妻丰满熟妇αv无码| 91在线日韩在线播放| 日韩精品成人网页视频在线| 欧美在线综合视频| 好吊色国产欧美日韩免费观看| 久久99国产综合精品1| 福利国产微拍广场一区视频在线 | 成人年鲁鲁在线观看视频| 无码国产伊人| 97精品伊人久久大香线蕉| 欧美人人干| 国产成人精品在线| 国产精品理论片| 亚洲精品欧美日本中文字幕| 精品国产三级在线观看| 日韩无码视频专区| 亚洲人成影视在线观看| 国产v欧美v日韩v综合精品| 亚洲精品日产精品乱码不卡| 欧美中文一区| 欧美精品在线视频观看| 国产肉感大码AV无码| 8090午夜无码专区| 亚洲成人高清在线观看| 亚国产欧美在线人成| 欧美成在线视频| AV无码无在线观看免费| 东京热av无码电影一区二区| 亚洲男人的天堂久久香蕉网| 国产99视频在线| 福利姬国产精品一区在线| 97国产在线视频| 精品福利网| 婷五月综合| 国产成人高清在线精品| aa级毛片毛片免费观看久| 国产午夜看片| 国产办公室秘书无码精品| 亚洲国产系列| 日本欧美一二三区色视频| 一级毛片免费不卡在线| 伊人久久婷婷五月综合97色| 亚洲精品无码在线播放网站| 日韩亚洲综合在线| 亚洲欧美极品| 欧美日韩一区二区三区在线视频| 啪啪啪亚洲无码| 成·人免费午夜无码视频在线观看 | 5555国产在线观看| 久久毛片网|