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

非貫通裂隙巖體損傷演化率相關性及變形特征*

2019-09-25 03:23:52鄧正定周尖榮王觀石王月梅
爆炸與沖擊 2019年8期
關鍵詞:裂紋變形模型

鄧正定,向 帥,周尖榮,王觀石,王月梅

(1. 江西理工大學建筑與測繪工程學院,江西 贛州 341000;2. 江西理工大學江西省環境巖土與工程災害控制重點實驗室,江西 贛州 341000;3. 深圳勘察測繪研究院有限公司,廣東 深圳 518000;4. 江西理工大學應用科學學院,江西 贛州 341000)

天然巖體和工程巖體中均普遍存在成因不同的節理、裂隙等軟弱結構面,其結構面幾何特征對巖體的變形破壞機制具有較大影響[1]。同時,應變率效應對巖體的力學特性也影響顯著[2]。因此,建立同時考慮應變率效應及裂隙幾何特征影響的力學模型,對相關巖體工程的安全與穩定具有重要意義。

針對非貫通裂隙巖體,已有較多研究:張平等[3]對斷續裂隙砂巖試件進行單軸靜、動載實驗,發現動載下分支裂紋擴展及貫通具有慣性效應,裂尖翼裂紋易在兩預制裂隙內端部產生直接貫通;劉紅巖等[4]、李地元等[5]分別通過對含裂隙巖體試件進行SHPB 動態沖擊實驗,發現裂隙的幾何特征及應變率對巖體的強度及破壞模式均影響顯著。在實驗基礎上,針對裂隙巖體已構建了動態本構模型:李夕兵等[6]在模型實驗基礎上,將冪函數模型描述節理面動態閉合變形,又在此基礎上建立了考慮應變率效應的裂隙巖體本構模型;張力民等[7]、劉紅巖[8]對已有TCK 模型進行改進,綜合考慮宏細觀缺陷的復合損傷變量(張量),建立了節理巖體動態損傷本構模型;李杰等[9]基于斷裂力學,建立了微裂紋動態擴展與體積擴容的方程,較好地描述了加卸載過程巖體的動態破壞過程。

已有裂隙巖體本構模型雖考慮了宏細觀初始缺陷對巖體的動態力學特性的影響,但對裂隙閉合且裂紋擴展過程中的損傷演化研究仍不盡完善。主要表現為:(1)針對裂隙對巖體的影響大多考慮初始裂隙對巖體的宏觀損傷,而對裂隙翼裂紋動態擴展引起的巖體損傷考慮不足;(2)目前基于斷裂力學與能量理論所構建的損傷模型,其附加應變能大多按裂隙面有效滑動驅動力作用下復合斷裂模式計算[7-8],計算顯然結果會偏大,但按純Ⅰ型或Ⅱ型斷裂擴展模式計算結果會偏小;(3)通常裂隙傾角較緩時,裂隙面有效滑動驅動力為零,已有模型中裂隙存在對巖體損傷也為零,但實際上初始裂隙即使沒有擴展,巖體強度仍較完整巖石強度更低[4]。因此,針對上述問題,引入合適的率效應動力演化方程表征巖體宏細觀缺陷的發展對巖體力學特性的影響,具有較好的學術和工程應用價值。

1 裂隙巖體動態損傷本構模型的建立

采用分離式Hopkinson 壓桿裝置對含非貫通裂隙巖體進行動態沖擊的實驗結果表明[4-5]:裂隙巖體宏觀損傷對巖體力學性能存在明顯軟化效應,應變率對巖體力學性能存在明顯硬化效應。因此,將裂隙巖體動態破壞過程視為具有復合損傷、靜態彈性特性、動態黏滯特性的非均質點組成。對黏彈性響應的Maxwell 體進行改進,將細觀損傷體D1與裂隙宏觀損傷體D2根據應變等價原理并聯組成宏細觀復合損傷體D12,在此基礎上再并聯Maxwell 體構造一個非貫通裂隙巖體動態損傷模型,如圖1 所示。

圖1 中,Maxwell 體具有瞬時變形、等速蠕變和松弛的性質,其主加載方向上的本構方程為:

式中:σM為Maxwell 應力,EM為Maxwell 彈性模量,為 應變率,為加載率,η 為實驗參數。

結合邊界條件,ε(0)=0,σM(0)=0,式(2)變換為:

對式(3)進行Laplace 逆變換,則有:

現假定巖體中未擴展的初始裂隙引起的巖體損傷定義為D0,巖體中微孔洞、微裂紋等細觀損傷定義為損傷變量D1,初始裂隙及其翼裂紋擴展后引起的巖體損傷定義為宏觀損傷變量D2,宏細觀綜合引起巖體的復合損傷變量為D12,根據Lemaitre 應變等價原理[7-8,10]:

式中:σD為損傷應力,E0為無損傷巖體的初始彈性模量。

通過對式(5)化簡,得裂隙巖體宏細觀復合損傷變量表達式為:

假定損傷體元件服從胡克定律,則損傷體元件本構方程為:

根據模型中各元件的串并聯關系,可知:

式(4)、(7)~(9)聯解,可得裂隙巖體動態損傷本構方程:

對于巖體細觀損傷,結合袁小清等[10]的研究成果,細觀損傷D1服從Weibull 分布:

式中:ε(t)為應變,ε0、m 為擬合常數。

下面,主要討論巖體初始裂隙及翼裂紋動態擴展所引起的宏觀損傷變量D2的計算。

2 裂隙巖體宏觀動態損傷演化

2.1 含裂隙巖體翼裂紋起裂與擴展

由實驗結果[11-12]可知:在壓剪狀態下,對于含不同傾角的非貫通裂隙巖體,當裂隙傾角較小時,裂隙面上的切向應力不足以克服法向應力產生的摩擦阻力,此時裂隙面主要發生閉合變形,不會產生沿裂隙面的滑移,整個裂隙巖體試件主要發生與完整巖石相似的張拉破壞模式,但巖體彈性模量及峰值強度均有一定程度的弱化;當裂隙傾角較大時,裂隙面除發生閉合變形外,主要發生沿裂隙面的剪切變形,裂隙尖端產生應力集中,進而導致裂紋擴展。

圖 2 裂隙擴展簡化模型Fig. 2 Simplified model for crack propagation

如圖2 所示,設雙軸應力σ3=kσ1=kσ,在遠場應力作用下,作用于傾角為α 的裂隙表面法向應力σα和切向應力τα分別為:

切向應力τα迫使裂隙上下表面相對滑動,但由于裂隙面是接觸或含充填物的,會產生一個摩擦力來抵抗這個滑動,此有效滑動驅動應力 τe應為:

式中: μ為裂隙面摩擦因數。

式(14)表示在遠場應力作用下的裂隙面上的滑動驅動力。當裂隙面上的切應力不足以克服正應力產生的摩擦力時,有效滑動驅動力為零。當有效滑動驅動力達到一定值時,在主裂紋尖端沿最大周向應力方向產生張拉翼裂紋,最大周向應力產生的應力強度因子為[13]:

式中:a 為原裂隙半長,θ 為與原裂隙呈一定角度的起裂擴展角。

由式(15)可知,KⅠ(θ)隨τe增大而增大,達到斷裂韌度KⅠC時,裂紋開始起裂,則有:

由式(14)、(16),則可得起裂強度:

結合Lee 等[14]的研究,裂隙擴展后翼裂紋尖端應力強度因子為:

式中:l 為翼裂紋擴展長度。

式(18)表示,在滑動驅動應力作用下,翼裂紋擴展后尖端的應力強度因子。式(18)右邊的前半部分表示有效滑動驅動力促使裂隙翼裂紋擴展,而后半部分表示翼裂紋面的法向壓應力阻止翼裂紋擴展。

在壓剪作用下,巖體裂隙在初始階段沿最大周應力方向擴展,但沿一個彎曲路徑后很快趨于最大主應力方向[15]。因此,可將翼裂紋擴展路徑簡化為與主加載方向平行的直線,則式(18)可簡化修正為:

當翼裂紋尖端應力強度因子KⅠ(l)減小到KⅠC時停止擴展,據此可得到翼裂紋最大擴展長度lm:

將式(20)中的前半部分冪級數展開并取前兩項,可解得:

2.2 裂隙巖體準靜態損傷演化能量機制

如圖2 所示,主裂隙擴展前,閉合裂隙的Ⅰ型應力強度因子為零,對于平面應力問題,單裂隙引起的附加應變能為:

式(22)表示翼裂紋擴展前主裂隙尖端Ⅱ型應力強度因子引起的附加應變能,其隱含的條件為加載應力小于或等于起裂強度σ0,因此式(22)可改寫為:

當σ>σ0時,翼裂紋開始擴展,主要為Ⅰ型擴展,其由翼裂紋擴展所引起的附加應變能U2為:

聯解式(15)、(17)、(24),并考慮翼裂紋起裂應力,可得翼裂紋擴展引起的附加應變能:

由式(23)、(25)所得的應變能,均考慮了有效滑動驅動力。如裂隙傾角較小時,裂隙面的有效滑動驅動力為零,按以上方法計算所得的附加應變能顯然也為零,但主裂隙尖端即使不能擴展,裂隙的存在仍對巖體的變形及強度有一定程度的影響,說明它仍儲存一定應變能。許多文獻采用Ⅰ型應力強度因子來計算應變能,在拉剪狀態下,此方法較為適合,但在壓剪狀態下,顯然計算結果會偏大[15]。因此,本文中采用閉合變形的方法計算裂隙受壓應力作用下引起的附加應變能。

裂隙的閉合變形主要由裂隙間充填物法向壓縮變形所致,其法向變形規律可遵循衰減法則[16],即:

式中:dεα為偏應變,dσα為偏應力,Eα為充填物彈性模量,εm為最大法向應變,εα為法向應變。

對式(26)進行積分,可得:

將式(12)代入式(27)中,可得:

則在法向壓應力作用下儲存在裂隙中的彈性應變能U3為:

式中:H 為試件高度。

2.3 裂隙巖體動態損傷演化能量機制

式中: V為翼裂紋擴展速度。

對于裂速函數,Freund 等[18]提出了簡化的近似公式:

其中 cR為瑞利波速度,表示為[19]:

式中:E0、 ν、ρ 分別為巖石彈性模量、泊松比、密度。

寧建國等[20]認為翼裂紋擴展速度與翼裂紋擴展長度、應變及應變率有如下關系:

式中:εcr為臨界應變。

根據應變等價原理,翼裂紋擴展長度在形成的同時伴隨巖體應變,通過對式(33)觀察,翼裂紋擴展長度l 與應變ε 實際上是呈正相關的,這也表明裂紋擴展速度與應變率也呈正相關關系。潘紅宇等[21]通過應變率和裂紋擴展速度的實驗數據擬合也發現了相似規律,因此裂紋擴展速度可定義為:

式中:β 為實驗擬合參數。

聯解式(30)~(31)、(34),可得翼裂紋擴展尖端的動態強度因子:

根據式(35)所得動態強度因子,式(21)所示的靜態擴展長度lm修正為翼裂紋動態擴展長度:

2.4 裂隙巖體宏觀動態損傷變量

假定完整巖石的彈性模量為E0,則雙軸平面應力下儲存在巖體內的彈性應變能U0為:式中:L 為試件直徑,H 為試件高度。

根據能量互易定理,巖體因裂隙存在而引起附加應變能,同時導致巖體的等效彈性模量的弱化,假定裂紋擴展后的宏觀損傷變量為D2,則有:

式中: ?U為初始裂隙及裂隙翼裂紋擴展后的附加應變能。有:

則根據式(39)~(40),可得主加載方向上宏觀損傷變量:

結合式(6)、(10)~(11)、(41),可得到裂隙巖體動態損傷模型:

3 裂隙巖體動態變形特征

3.1 模型驗證

為驗證模型的合理性,討論裂隙幾何特征及應變率效應對巖體力學特性的影響。以文獻[22-23]實驗資料為例,采用的巖石相似材料[4],采用單軸加載的方式(即k=0),試件尺寸高50 mm,直徑50 mm,應變率近似取恒應變率75 s?1,裂紋擴展速度V=0.3cR[24],非貫通裂隙傾角α 分別為0°、30°、45°、60°,裂隙長度均為11 mm。巖石的彈性模量E0=8.7 GPa,泊松比ν=0.28,斷裂韌度KⅠC=0.25 MPa·m1/2;裂隙結構面的彈性模量Eα=198 MPa,最大應變εm=0.02,摩擦因數 μ=0.3。

從巖體動態損傷模型可知,模型中的參數ε0、m、η、EM需通過實驗數據擬合才能得到,其取值對巖石力學性能具有重要影響。采用文獻[22-23]中完整試件的數據資料,將巖體力學參數代入本構模型中,可得應力應變曲線,以此分析各參量對巖體動態本構關系的影響。

圖3 為參量m、η、EM均為常量,改變參量ε0時的巖石本構關系。由圖3 可知,巖石的峰值強度及其所對應的應變均隨參量ε0的增大而增大,且ε0對本構關系的影響主要表現在非彈性的微元體破裂發展階段,說明參量ε0是巖石內部微元強度平均大小的反映[25]。

圖4 為參量ε0、η、EM均為常量,改變參量m 時的巖石本構關系。由圖4 可知,參量m 對巖石的峰值強度的影響較小,且對巖石在彈性變形階段的影響也較小,而對巖石峰前非彈性變形階段及峰后變形階段的曲率影響顯著。隨著m 的增大,巖石本構關系曲線在峰前及峰后階段的曲率均相應增大,即更為陡峭,說明參量m 是巖石內部微元強度集中程度的反映。

圖5 為參量ε0、m、EM均為常量,改變參量η 時的巖石本構關系。由圖5 可知,參量η 對巖石的峰值強度及峰前非彈性變形階段影響顯著,而對巖石的峰值應變及彈性變形階段的曲率影響較小。隨著η 的增大,巖石峰值強度及非彈性變形階段的曲率均增大,巖石的動態黏滯效應也相應增強,且η 越大,粘滯效應的增長幅度也相應更大,說明參量η 是巖石應變率效應相關度的反映。

圖6 為參量ε0、m、η 均為常量,改變Maxwell 彈性模量EM時的巖石本構關系。由圖6 可知,參量EM對巖石的峰值強度及本構關系曲線的曲率均有一定影響,而峰值應變基本未受影響,表明EM也是巖體黏滯效應的重要參數,但與參量η 的影響不同,巖石峰值強度增長幅度隨EM的增大而減小。同時,隨著EM的增大,本構關系曲線峰前、峰后階段的曲率均相應更大,表明參量EM還反映了黏滯效應的集中度,以上分析說明參量EM是巖石應變率效應相關度及集中度的綜合反映。

圖 3 參量ε0 對本構關系的影響Fig. 3 Influence of parameters ε0 on constitutive relation

圖 4 參量m 對本構關系的影響Fig. 4 Influence of parameters m on constitutive relation

圖 5 參量η 對本構關系的影響Fig. 5 Influence of parameters η on constitutive relation

圖 6 參量EM 對本構關系影響Fig. 6 Influence of parameters EM on constitutive relation

圖 7 裂隙巖體的實驗結果與理論計算結果比較Fig. 7 Comparison of experimental and theoretical results of fractured rock mass

綜合以上理論分析可知,各擬合參量對本構關系的影響規律。將文獻[22-23]中的完整巖石實驗數據及力學參數代入式(42),則可擬合得到參量:ε0=0.029 4,m=1.28,η=0.23,EM=3.2 GPa。再將其他裂隙巖體物理力學參數代入式(42),即可得到不同裂隙巖體的理論模型曲線。圖7 為單軸動荷載下用本文理論模型得到的巖體應力應變曲線與文獻[22-23]中的實驗數據對比。如圖7 所示,完整巖石、含30°、45°、60°非貫通裂隙巖體的實驗結果與模型計算結果基本吻合。模型計算所得的巖體動態應力曲線較好地反映了巖體變形的全過程,說明本文模型能較好地描述非貫通裂隙巖體的沖擊力學特性。

3.2 裂隙幾何特征對巖體動態力學特性的影響

計算模型和參數與第3.1 節相同,應變率均為75 s?1,通過改變非貫通裂隙傾角α 及貫通度,研究單軸沖擊荷載作用下非貫通裂隙幾何特征對巖體動態力學特性的影響。圖8~9 為貫通度分別為50%、80%,不同傾角裂隙的存在對巖體初始宏觀損傷變量D0及巖體起裂強度σ0的影響規律。圖10 為非貫通裂隙傾角為45°時,不同裂隙貫通度情況下巖體的動態應力應變曲線。

圖 8 初始損傷變量隨裂隙傾角的變化Fig. 8 Variation of initial damage with slit angle

圖 9 起裂強度隨裂隙傾角的變化Fig. 9 Variation of fracture strength with slit angle

由圖8 可知,當裂隙傾角較小時,巖體初始損傷變量D0隨裂隙傾角增大而增大;當裂隙傾角較大時,巖體初始損傷變量D0隨裂隙傾角α 增大而減小。值得注意的是,由于計算模型考慮了裂隙結構面閉合產生的附加應變能,因而裂隙傾角α 即使為0°,仍會產生初始損傷,巖體強度較完整巖石強度會降低,這也跟諸多文獻實驗結果較一致。當裂隙傾角α 達到90°時,由于裂隙表面法向正應力和切向應力都為零,巖體初始損傷變量D0也為零,因而含90°裂隙巖體與完整巖石力學特性無明顯區別。

圖 10 不同裂隙貫通度的巖體動態應力應變曲線Fig. 10 Dynamic stress-strain curves of rock mass with different fracture penetrability

由圖9 可知,巖體裂隙起裂強度σ0隨裂隙傾角呈“凹”型分布。當裂隙傾角小于裂隙面摩擦角時,作用于裂隙表面的切向應力不足以克服正應力產生的摩擦力,因而裂隙面有效滑動驅動力為零。因此裂隙傾角較小時,巖體裂隙起裂強度σ0無窮大(圖中虛線部分),巖體主要產生張拉破壞及微裂隙、微孔洞裂紋擴展。當裂隙傾角大于裂隙面摩擦角時,有效滑動驅動力大于零,隨著應力的增加,裂隙尖端應力強度因子達到斷裂韌度。由圖9 可知:裂隙傾角45°時,起裂強度最小,這也與實驗結果較為吻合;當裂隙傾角繼續增大時,裂隙有效滑動驅動力又逐漸較小,直至宏觀裂隙起裂強度σ0又達到無窮大。

由圖10 可知,裂隙貫通度對巖體的動態力學特性有顯著影響,巖體峰值強度隨裂隙貫通度的增加而顯著減小。通過應力應變曲線可發現,不同貫通度的裂隙存在對巖體彈性變形階段幾乎沒有影響,而在屈服階段影響明顯。這是由于在彈性變形階段,裂隙尖端應力強度因子尚未達到巖石的斷裂韌度,裂隙的存在對巖體的損傷主要為裂隙面的閉合變形。當裂隙尖端應力強度因子達到斷裂韌度后,翼裂紋開始擴展,根據應變等價原理,裂隙擴展后引起巖體的宏細觀耦合損傷顯著增加,因此應力應變曲線明顯變緩直至破壞。

3.3 應變率對巖體動態力學特性的影響

計算模型和參數與第3.1 節相同,裂隙傾角分別取45°、60°,荷載應變率分別取準靜載0 s?1、近似恒應變率100、200、300、400 s?1,裂隙翼裂紋動態擴展長度及巖體應力應變曲線計算結果分別如圖11~12 所示。

圖 11 應變率對翼裂紋擴展長度的影響Fig. 11 Influence of strain rate on length of wing crack

圖 12 應變率對巖體力學特性的影響Fig. 12 Influence of strain rate on mechanicalproperties of rock

如圖11 所示,應變率對裂隙翼裂紋擴展長度有顯著影響。當裂隙傾角為45°時,由于在其他相同條件下,其裂隙表面有效滑動驅動力較裂隙傾角為60°時更大,因此翼裂紋擴展長度更長。翼裂紋擴展長度隨應變率呈線性分布規律,應變率越大,其翼裂紋最大擴展長度越小。

如圖12 所示,在不同荷載應變率下,巖體在彈性變形階段應變率效應較小,在翼裂紋擴展階段應變率效應顯著。隨著載荷應變率的增大,巖體試件的峰值強度也逐漸變大,5 種不同應變率荷載作用下峰值強度分別為31、38、46、57 和66 MPa,且對應的峰值應變及總應變也隨之增加,這與文獻[24,26]中巖體在不同載荷應變率下實驗結果得到的變形規律一致。

3.4 巖石彈性模量及斷裂韌度對巖體動態力學特性的影響

計算模型和參數與第3.1 節相同,通過改變巖石的彈性模量E0及巖石的斷裂韌度KⅠC,研究巖石性質對巖體動態力學特性的影響規律,分別如圖13~14 所示。

圖 13 E0 對巖體動態力學特性的影響Fig. 13 Influence of E0 on mechanical properties

圖 14 KⅠC 對巖體動態力學特性的影響Fig. 14 Influence of KⅠC on mechanical properties

如圖13 所示,巖石的彈性模量E0對巖體的動態力學特性具有顯著影響。隨著彈性模量E0的增大,應力應變曲線的曲率及峰值強度均顯著增大,巖石彈性模量E0分別為2、4、6、8 和10 GPa 時,對應的巖體峰值強度分別為26.5、49.1、71.8、94.5 和117.2 MPa。值得注意的是,通過圖13 中的巖石彈性模量E0與圖6 中的Maxwell 彈性模量EM對巖體動態力學特性的影響的分析比較,在兩者增加同樣幅度時,巖石的彈性模量對巖體的峰值強度及曲線曲率影響更為明顯。

如圖14 所示,巖石的斷裂韌度KⅠC對巖體的動態力學特性也有較大影響。隨著斷裂韌度KⅠC的增大,巖石的峰值強度也相應增大。值得指出的是,此處斷裂韌度對巖體強度的影響是在裂隙面有效滑動驅動力大于零的前提下,如裂隙傾角趨近0°或90°,巖體強度將不受斷裂韌度KⅠC的影響。由不同斷裂韌度KⅠC下巖體的應力應變曲線可發現,斷裂韌度對巖體的影響主要表現在非彈性階段,而在彈性階段幾乎不受影響,這是由于,斷裂韌度主要影響巖體的起裂強度σ0。當加載應力σ<σ0時,翼裂紋沒有擴展,裂隙的存在對巖體的損傷主要表現為裂隙的閉合變形,而加載應力σ 超過起裂強度σ0后,翼裂紋開始擴展,裂隙對巖體的損傷顯著增加,因而巖體應力應變曲線也迅速變緩。

4 結 論

(1)針對自然界大量存在的非貫通裂隙巖體,基于變形理論、斷裂力學與能量理論相結合,運用模型變形元件的方法,綜合考慮巖體宏細觀缺陷及損傷演化的率相關性,構建了裂隙巖體動態損傷本構模型。

(2)根據Lemaitre 應變等價原理,將巖體的細觀損傷和宏觀損傷耦合,再結合Maxwell 體以描述動態加載過程中巖體的應變率效應,模型中的宏觀損傷變量同時考慮了裂隙翼裂紋的動態擴展以及裂隙面的閉合變形,克服了現有本構模型的不足。

(3)利用模型討論了裂隙幾何特征、應變率及巖石性質對巖體動態變形特征的影響,巖體強度及變形均隨裂隙傾角的增大先變小后增大,應變率效應對巖體具有明顯的硬化效應,且在裂隙擴展階段更為明顯,同時巖石本身的彈性模量及斷裂韌度也對巖體的強度及變形特征有重要影響。

猜你喜歡
裂紋變形模型
一半模型
裂紋長度對焊接接頭裂紋擴展驅動力的影響
重要模型『一線三等角』
談詩的變形
中華詩詞(2020年1期)2020-09-21 09:24:52
重尾非線性自回歸模型自加權M-估計的漸近分布
Epidermal growth factor receptor rs17337023 polymorphism in hypertensive gestational diabetic women: A pilot study
“我”的變形計
例談拼圖與整式變形
會變形的餅
3D打印中的模型分割與打包
主站蜘蛛池模板: 亚洲成人福利网站| 四虎成人精品| 99精品视频九九精品| 免费观看三级毛片| 色天天综合| 日本成人精品视频| 国产免费一级精品视频| 日本高清在线看免费观看| 欧美精品成人一区二区视频一| 在线观看无码av免费不卡网站 | 国产精品毛片一区视频播| 成年人国产网站| 国产黑丝一区| 精品三级在线| 天天综合网站| 五月综合色婷婷| 亚洲欧美h| 制服无码网站| 熟女成人国产精品视频| 2020国产精品视频| 色爽网免费视频| 欧美一级片在线| 青青青伊人色综合久久| 欧美日韩成人在线观看| 国产欧美视频在线观看| 国产精品一线天| 久草中文网| 国产尤物在线播放| 老司国产精品视频91| jizz国产视频| 欧美午夜理伦三级在线观看| 免费a级毛片18以上观看精品| 日韩久久精品无码aV| 国产一二视频| 欧美成人免费一区在线播放| 九九这里只有精品视频| 成人午夜福利视频| 亚洲精品天堂自在久久77| 欧美午夜性视频| 亚洲 成人国产| 日韩高清在线观看不卡一区二区| 最新日本中文字幕| 久久精品无码国产一区二区三区| 伦精品一区二区三区视频| AV不卡国产在线观看| 精品福利网| 国产99视频精品免费视频7| 久久精品人妻中文系列| 亚洲欧美日韩综合二区三区| 国产专区综合另类日韩一区| 四虎成人免费毛片| 四虎永久在线精品影院| 亚洲综合色婷婷中文字幕| 中国国产一级毛片| 久久99精品久久久久久不卡| 久无码久无码av无码| 久久久久免费精品国产| 午夜精品福利影院| 亚洲免费福利视频| 亚洲高清日韩heyzo| 欧美精品啪啪一区二区三区| 久久网综合| 六月婷婷精品视频在线观看| 成人免费视频一区| 久久综合久久鬼| 五月天福利视频| 免费高清毛片| 日韩在线影院| 91亚洲国产视频| 天天摸夜夜操| 香蕉色综合| 午夜福利在线观看入口| 91啪在线| 亚洲人成网站色7777| 无码福利日韩神码福利片| 色综合激情网| 久草热视频在线| 久久这里只有精品66| 波多野结衣视频一区二区| 久久久久久国产精品mv| 国产精品人莉莉成在线播放| 狠狠色噜噜狠狠狠狠色综合久|