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

中心圓孔圓盤試件三維最大無量綱應力強度因子的標定

2012-09-26 12:45:04梁亞磊
中國有色金屬學報 2012年8期
關鍵詞:裂紋

張 盛,梁亞磊,

(1. 河南理工大學 能源科學與工程學院,焦作 454003;2. 中國平煤神馬集團十礦,平頂山 467000)

中心圓孔圓盤試件三維最大無量綱應力強度因子的標定

張 盛1,梁亞磊1,2

(1. 河南理工大學 能源科學與工程學院,焦作 454003;2. 中國平煤神馬集團十礦,平頂山 467000)

基于平面假設標定的圓盤試件的無量綱應力強度因子對巖石斷裂韌度的準確測定有一定的影響。使用三維斷裂分析軟件FRANC3D對不同中心孔尺寸的中心圓孔平臺巴西圓盤HFBD試件的三維最大無量綱應力強度因子進行準確標定。結果表明:沿著圓盤厚度方向虛擬三維裂紋前緣各點的最大無量綱應力強度因子并不相等,其中心點與二維標定值最大誤差接近3%,推薦了一種不受三維尺寸影響可簡化為二維分析的中心圓孔平臺圓盤試件尺寸。

中心圓孔圓盤試件;三維應力強度因子;巖石斷裂韌度

在礦井巷道支護、水壓致裂法測定地應力以及爆破施工等工程中都涉及到巖石斷裂問題,巖石斷裂韌度在這些領域已經得到一定的應用,但其準確測定仍然存在很多問題。如裂縫難以預制、無量綱應力強度因子標定不準確、不同測試方法測定的巖石斷裂韌度存在較大的差異等[1]。針對這些問題,1995年國際巖石力學協會推薦了一種人字形切槽巴西圓盤試件(Cracked chevron notched Brazilian disc, CCNBD)測試巖石Ⅰ型斷裂韌度的方法[2];由于該類圓盤試樣具有體積小,加載方便,對試驗設備要求不高等優點,得到了研究者的推廣。但是CCNBD試件構形復雜,不能按平面問題進行分析,其斷裂韌度公式中的無量綱應力強度因子的準確標定就是一個非常復雜的工作,賈學明和王啟智[3]、樊鴻等[4]先后對其進行了大量的標定和修正工作?;趫A盤試件的優點,同時為了避免CCNBD復雜的三維構形,王啟智和賈學明[5]提出采用平臺巴西圓盤(Flattened Brazilian disc, FBD)進行巖石斷裂韌度的測試;隨后,張盛和王啟智[6]進一步比較了5種不同類型圓盤的測試方法。相對于需要提前預制裂縫的圓盤試件,中心圓孔平臺圓盤試件(Holed flattened Brazilian disc, HFBD)試件更容易制作,而且比FBD試件更能引發中心裂紋產生并避免加載端效應[6?7]。但是,由于HFBD試件具有三維構形,采用平面分析的方法對HFBD試件進行無量綱應力強度因子的標定,可能會由于不同的中心孔徑尺寸造成標定值的較大誤差,從而影響巖石斷裂韌度值的準確測定。

長期以來,一些研究者使用不同的方法來獲得裂紋體裂紋前緣的應力強度因子。張玉軍和李治國[8]采用有限元位移法分析了帶裂紋隧道二次襯砌裂紋尖端的應力強度因子。王承強和鄭長良[9]利用平面扇形域哈密頓體系方程,通過分離變量法及共軛辛本征函數向量展開法,推導出一個圓形奇異解析單元列陣用于求解平面裂紋的應力強度因子。黃海燕和劉小健[10]分別使用節點位移解析法和有限元法計算了含中心裂紋有限加筋板的應力強度因子,并通過對比和分析,提出了一個用于評估應力強度因子計算精度的參數。但是,這些研究都基于平面假設,其計算結果并沒有考慮三維裂紋的影響。事實上,受載裂紋體裂紋前緣應力強度因子沿其前緣分布并不一致,樊鴻和張盛[11]采用應力初始單元法進行三維裂紋的分析,但該方法構造裂紋非常復雜,計算工作量較大。鑒于此,使用專業的三維斷裂分析軟件FRANC3D(Fracture analysis code in three dimensions)計算HFBD試件受載起裂時的三維裂紋的應力強度因子[12],并研究圓盤試件中心孔徑尺寸對最大無量綱應力強度因子的影響,對斷裂韌度公式中的最大無量綱應力強度因子系數進行了修正。本文作者對三維裂紋進行無量綱應力強度因子的標定將有助于建立更加合理的斷裂韌度公式,為準確獲得巖石的動態斷裂韌度值奠定基礎。

1 HFBD試件測試巖石斷裂韌度的原理

在試件中預制裂縫的寬度和形狀會給巖石斷裂韌度的測試值帶來很大的誤差[13],而對于HFBD和FBD試件都不需要預先制作裂縫。它們的測試原理是從普通巴西圓盤試件測試巖石斷裂韌度的方法推廣而來的[14?15]。試件受力以后,會產生裂紋,通過臨界裂紋對應的最大無量綱應力強度因子和二次倒拐的局部最小載荷來確定斷裂韌度值[6]。

如圖1(a)所示,HFBD試件受均布壓力載荷作用后,沿加載直徑AB方向中心孔邊引起高度的應力集中。隨著加載載荷的增大,裂紋首先在孔邊應力集中處萌生,進而沿加載直徑AB方向擴展[6?7],HFBD試件利用孔邊高度的應力集中來幫助試件在中心圓孔邊發生起裂,產生裂紋,形成中心圓孔裂紋平臺圓盤試件(Hole-cracked flattened Brazilian disc, HCFBD),如圖1(b)所示。根據裂紋擴展過程中臨界裂紋長度與之對應的局部最小載荷與最大無量綱應力強度因子的一致關系確定斷裂韌度值。

圖1 圓盤試件加載及尺寸參數示意圖Fig. 1 Loading and dimension diagram of disc sample: (a) Holed flatterned Brazilian disc; (b) Holed-cracked flattened Brazilian disc

如圖1所示,HFBD試件受均布壓力合力為p,加載角為2β,半徑為R,裂紋擴展長度為2a,厚度為t,中心孔半徑為r。根據文獻[6],由該類圓盤試件確定的應力強度因子如式(1)所示,斷裂韌度值表達式如式(2)所示:

式中:KI為試樣受載起裂后裂紋前緣應力強度因子;KIC為試驗測試的巖石的斷裂韌度值;Y為無量綱應力強度因子;ρ=r/R為圓盤試件中心孔半徑與圓盤半徑之比;α為無量綱裂紋長度,是裂紋擴展長度與圓盤半徑之比a/R;pmin為加載曲線上的局部最小載荷;Ymax為臨界裂紋長度與之對應的最大無量綱應力強度因子。

當裂紋前緣應力強度因子達到某個臨界值時,裂紋發生擴展,此時這個應力強度因子臨界值稱為斷裂韌度。根據斷裂原理,在受載材料裂紋擴展過程中的任何載荷p以及其對應的裂紋長度可由式(1)計算得出巖石斷裂韌度。但是,脆性材料的裂紋擴展速度非??欤茈y精確確定在某一載荷時的裂紋長度。所以,無法使用式(1)直接計算巖石材料的斷裂韌度。然而,中心圓孔平臺巴西圓盤試件在加載時載荷隨裂紋長度的變化情況對應于無量綱應力強度因子的變化情況,即隨著裂紋的擴展無量綱應力強度因子經歷了上升、最大值、下降的過程,其最大無量綱應力強度因子Ymax正好對應于加載曲線上局部最小載荷pmin(此時,臨界裂紋長度與之對應)[6]。因此,可由pmin和Ymax代入式(2)計算巖石斷裂韌度KIC。

2 三維裂紋模擬分析

通過上面的分析可知,HFBD試件在標定其最大無量綱應力強度因子時,需要在HFBD的基礎上構造虛擬裂紋從而形成HCFBD試件,根據虛擬裂紋計算其受載后裂紋擴展過程中沿厚度方向不同點處的應力強度因子,從而確定其最大無量綱應力強度因子。采用三維斷裂分析軟件FRANC3D可進行三維虛擬裂紋的模擬分析。FRANC3D[12]是美國康奈爾大學斷裂工作組(Cornell fracture group, CFG)自主研發的一套具有建模、應力分析、應力強度因子計算、裂紋自動擴展模擬等功能的軟件。模擬裂紋能力是FRANC3D的特色之一,用戶可以通過內置裂紋庫或者以文件的方式生成所需要的裂紋。對于裂紋的構形,可以設置任意的裂紋前緣。對模型計算結束后,用戶可以儲存應力強度因子計算結果,以便調用。

基于本文作者的前期工作,采用HFBD圓盤試件能夠測試巖石的斷裂韌度[6],但該方法是建立在二維平面分析的基礎上的,在此擬進行三維分析。根據圓盤試件的對稱性,取其一半進行建模,如圖2所示。模型幾何參數為:R=40mm,ρ=0.2,2β=20°,t=30 mm。大理巖材料參數:彈性模量E=1.63 GPa,泊松比μ=0.3。

考慮到裂紋沿模型厚度Z方向的對稱性,主要取沿模型標準化一半厚度方向3個點處的應力強度因子,即取無量綱厚度λ=2Z/t為0、0.25、0.5處的應力強度因子,λ=0代表圓盤沿厚度方向外部表面處應力強度因子,λ=0.25代表圓盤沿厚度方向1/4處點的應力強度因子,λ=0.5代表圓盤沿厚度方向中心處點應力強度因子。

根據FRANC3D的計算,得到中心孔徑ρ=0.2,在不同的裂紋長度條件下,無量綱裂紋應力強度因子Ymax與無量綱裂紋長度a/R的關系如圖3所示,其中,橫坐標代表無量綱裂紋長度即裂紋長度a和圓盤半徑R之比,縱坐標代表無量綱應力強度因子Y。

圖2 圓盤試件三維模型(一半)Fig. 2 Three dimensions model of disc specimen (half)

圖3 無量綱應力強度因子Y和無量綱裂紋長度a/R的關系Fig. 3 Relation between dimensionless stress intensity factorYand dimensionless crack lengtha/R

從圖3可以看出,圓盤試件裂紋前緣中心點,1/4點以及表面點的無量綱應力強度因子隨著裂紋長度的增加,無量綱應力強度因子隨著裂紋長度的增加先增加后減小;沿厚度方向的分布有較大差異;而最大無量綱應力強度因子所對應的裂紋長度并不一致,其值相差也較大。當λ=0時,Ymax=1.214;當λ=0.25時,Ymax=1.048;當λ=0.5 時,Ymax=1.017;裂紋前緣圓盤表面點處的最大無量綱應力強度因子Ymax比中心點處的Ymax大了16.2%。另外,使用有限元軟件ANSYS對該類圓盤試件進行平面問題處理,得到二維計算的Ymax結果為1.038,這比三維方法計算的中心點處Ymax大2.1%。其原因可能是由于三維裂紋考慮了靠近圓盤兩個端面到中心整個裂紋長度的情況,由于靠近兩端面的自由表面上沿裂紋切線方向的位移不受限制,有明顯的泊松比效應,背離了平面應變情況,所以與平面問題處理產生一定的誤差。

3 HFBD試件斷裂韌度公式的三維修正

對圓盤半徑R=40 mm、加載角2β=20°、厚度t=30 mm,9種不同孔徑比(ρ)HFBD試件(ρ為0.1、0.15、0.2、0.25、0.3、0.35、0.4、0.45、0.5),進行的大量三維應力強度因子和二維應力強度因子計算后,得到最大無量綱應力強度因子Ymax和孔徑比ρ的關系如圖4所示。

從圖4可以看出,隨著圓盤試件中心孔徑ρ的增大,裂紋前緣表面點、中心點和采用二維應力強度因子計算方法確定的最大無量綱應力強度因子Ymax都在增大;由二維計算方法確定的Ymax,與三維計算方法圓盤中心確定的Ymax最大誤差為3.0%;與三維計算方法圓盤表面裂紋前緣點確定的Ymax最大誤差為27.8%。

圖4 最大無量綱應力強度因子Ymax和孔徑比ρ的關系Fig. 4 Relationship between biggest dimensional stress intensity factorYmaxand dimensional hole sizeρ

如果采用圓盤中心處點的Ymax代入式(2)計算材料的斷裂韌度,經過曲線擬合如圖5所示,得到三維方法計算的圓盤中心點標定的Ymax與中心孔徑比ρ的函數表達式如式(3)所示:

經驗證,在HFBD試件中心孔徑ρ為0.1~0.5范圍內,由式(3)得出的最大無量綱應力強度因子Ymax與有限元計算標定的值最大誤差小于2.0%,如表1所列。

圖5 圓盤中心點最大無量綱應力強度因子Ymax與中心孔徑比ρ的關系擬合曲線Fig. 5 Fitting curve for biggest dimensional stress intensity factorYmaxof center point with diemnsional hole sizeρ

表1 有限元與式(3)中Ymax計算值的比較Table 1 Comparison ofYmaxcalculation results by FEM and formula (3)

通過比較圓盤試件裂紋前緣表面點和中心處點的Ymax發現:在中心孔徑比ρ為0.1~0.4范圍內,裂紋前緣表面點與裂紋前緣中心點、二維方法確定的Ymax的差值隨ρ的增大而由大變??;在ρ=0.4附近,裂紋前緣表面點、中心點處Ymax與二維方法確定的Ymax基本上相等,二維Ymax比裂紋前緣中心點Ymax大0.41%、比裂紋前緣表面點Ymax小0.46%;在ρ為0.4~0.5范圍內,裂紋前緣表面點、中心處點Ymax的差值隨ρ的增大而由小變大。

從上面的分析可以得出,采用中心孔徑ρ=0.4的HFBD試件,最大無量綱應力強度因子Ymax受三維尺寸的影響可以忽略。因此,作者推薦一種用來測試巖石斷裂韌度的較為理想的中心圓孔圓盤試件,即圓盤的半徑R=40 mm,圓盤厚度t=30 mm,平臺加載角2β=20°,中心孔徑比ρ=0.4,此時由三維方法計算的圓盤表面點、圓盤中心點及二維方法計算的Ymax分別為2.20、2.18、2.19,所采用斷裂韌度的計算公式如下。

由該尺寸和公式得到的斷裂韌度值基本不受三維尺寸的影響,能夠采用平面分析得到準確的測試值。

4 結論

1) HFBD試件虛擬裂紋前緣不同點的最大無量綱應力強度因子值各不相同,隨著圓盤試件中心孔徑比ρ的增大,各點的最大無量綱應力強度因子Ymax都在增大,由二維計算方法標定的值相比于三維方法計算的圓盤中心點確定的Ymax最大誤差接近3%,因此,不能簡單地采用平面問題進行HFBD試件最大無量綱應力強度因子的標定。

2) 推薦了一種能夠基于平面問題進行分析測試巖石斷裂韌度的HFBD試件,該圓盤直徑R=40 mm,厚度R=30 mm,平臺加載角 2β=20°,中心孔徑比ρ=0.4,其斷裂韌度的計算公式如式(4)所示。

REFERENCES

[1]崔振東, 劉大安, 安光明, 周 苗. 巖石Ⅰ型斷裂韌度測試方法研究進展[J]. 測試技術學報, 2009, 3(2): 189?196. CUI Zhen-dong, LIU Da-an, AN Guang-ming, ZHOU Miao. Research progress in mode-Ⅰfracture toughness testing methods for rocks[J]. Journal of Test and Measurement Technology, 2009, 3(2): 189?196.

[2]ISRM Testing Commission. Suggested method for determining mode (Ⅰ): Fracture toughness using cracked chevron notched Brazilian disc (CCNBD) specimens[J]. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, 1995, 32: 57?64.

[3]賈學明, 王啟智. 斷裂韌度試樣CCNBD寬范圍應力強度因子標定[J]. 巖土力學, 2003, 21(6): 907?912. JIA Xue-ming, WANG Qi-zhi. Wide range calibration of the stress intensity factor for the fracture toughness specimen CCNBD[J]. Rock and Soil Mechanics, 2004, 21(6): 764?768.

[4]樊 鴻, 張 盛, 王啟智. 巖石斷裂韌度試樣CCNBD臨界應力強度因子的全新數值標定[C]//中國計算力學大會’2010(CCCM2010)暨第八屆南方計算力學學術會議(SCCM8)論文集. 綿陽: 西南科技大學, 2010?08?20. FAN Hong, ZHANG Sheng, WANG Qi-zhi. A completely new numerical calibration of the critical stress intensity factor for the rock fracture toughness specimen CCNBD[C]//The 2010th China Conference on Computational Mechanics and the 8th South Conference on Computational Mechanics. Mianyang: Southwest University of Science and Technology, 2010?08?20.

[5]王啟智, 賈學明. 用平臺巴西圓盤試樣確定脆性巖石的彈性模量、拉伸強度和斷裂韌度—第一部分: 解析和數值結果[J].巖石力學與工程學報, 2002, 21(9): 1285?1289. WANG Qi-zhi, JIA Xue-ming. Determination of elastic modulus, tensile strength and fracture toughness of brittle rocks by using flattened Brazilian disk specimen—Part Ⅰ: Analytical and numerical results[J]. Chinese Journal of Rock Mechanics and Engineering, 2002, 21(9): 1285?1289.

[6]張 盛, 王啟智. 用五種圓盤試件的劈裂試驗確定巖石斷裂韌度[J]. 巖土力學, 2009, 30(1): 12?18. ZHANG Sheng, WANG Qi-zhi. Determination of rock fracturetoughness by split test using five types of disc specimen[J]. Rock and Soil Mechanics, 2009, 30(1): 12?18.

[7]FISCHER M P, ELSWORTH D, ALLEY R B, ENGELDER T. Finite element analysis of the modified ring rest for determining mode (Ⅰ): Fracture toughness[J]. International Journal of Rock Mechanics & Mining Sciences, 1996, 33(1): 1?15.

[8]張玉軍, 李治國. 帶裂紋隧道二次襯砌承載能力的平面有限元計算分析[J]. 巖土力學, 2005, 26(8): 1201?1206. ZHANG Yu-jun, LI Zhi-guo. 2D finite element analysis of bearing capacity of tunnel lining with cracks[J]. Rock and Soil Mechanics, 2005, 26(8): 1201?1206.

[9]王承強, 鄭長良. 平面裂紋應力強度因子的半解析半有限元法[J]. 工程力學, 2005, 22(1): 33?37. WANG Chen-qiang, ZHENG Chang-liang. Semi-analytical finite element method for plane crack stress intensity factor[J]. Engineering Mechanics, 2005, 22(1): 33?37.

[10]黃海燕, 劉小健. 含中心裂紋有限加筋板應力強度因子的計算精度評估[J]. 船舶力學, 2008, 12(5): 765?769. HUANG Hai-yan, LIU Xiao-jian. Accuracy evaluation of stress intensity factor about stiffened plate with a central crack[J]. Journal of Ship Mechanics, 2008, 12(5): 765?769.

[11]樊 鴻, 張 盛. 求解應力強度因子的單元初始應力法[J].武漢理工大學學報, 2008, 30(11): 90?93. FAN Hong, ZHANG Sheng. An unit initial stress method to calculate stress intensity factor[J]. Journal of Wuhan University of Technology, 2008, 30(11): 90?93.

[12]賈學明, 王啟智. 三維斷裂分析軟件FRANC3D[J]. 計算力學學報, 2004, 21(6): 764?768. JIA Xue-ming, WANG Qi-zhi. Three-dimensional fracture analysis software[J]. Chinese Journal of Computational Mechanics, 2004, 21(6): 764?768.

[13]CHEN M, ZHANG G Q. Laboratory measurement and interpretation of the fracture toughness of formation rocks at great depth[J]. Journal of Petroleum Science and Engineering, 2004, 41: 221?231.

[14]GUO H, AZIZ N I, SCHMIDT L C. Rock fracture toughness determination by the Brazilian test[J]. Engineer Geology, 1993, 33(2): 177?188.

[15]YANG S C, TANG T X, ZOLLINGER D, GURJAR A. Splitting tension tests to determine concrete fracture parameters by peak-load method[J]. Advn Cem Bas Mat, 1997, 5: 18?28.

(編輯 龍懷中)

Calibration of three-dimensional maximum dimensionless stress intensity factor of holed flattened Brazilion disc specimen

ZHANG Sheng1, LIANG Ya-lei1,2
(1. School of Energy Science and Engineering, Henan Polytechnic University, Jiaozuo 454003, China; 2. China Pingmei Shenma Group No.10 Coal Mine, Pingdingshan 467000, China)

There was a certain effect to measure accurately the rock fracture toughness with holed disc specimens, of which the maximum dimensionless of stress intensity factor was calibrated with 2D method. By three-dimensional fracture analysis software FRANC3D, the maximum dimensionless of the stress intensity factor of disc samples was calculated, which has different radius sizes of center hole. The results show that the factors on the crack the front edge along the disc thickness direction are varied, even the maximum relative error of the center point is about 3%, compared to the value made with 2D method. A type of holed flattened Brazilian disc specimen is recommend, which can be simplified to a two-dimensional problem without three-dimensional size effect.

holed flattened Brazilian disc specimen; three-dimensional stress intensity factor; rock fracture toughness

TD313

A

國家自然科學基金青年基金資助項目(51004043);河南理工大學博士基金資助項目(B2008-66648270);河南省教育廳自然科學基礎研究計劃項目(021- 2010A440005)

2011-07-21;

2011-12-10

張 盛,副教授;電話:18603918628;E-mail: zhangs@hpu.edu.cn

1004-0609(2012)08-2347-06

猜你喜歡
裂紋
基于擴展有限元的疲勞裂紋擴展分析
裂紋長度對焊接接頭裂紋擴展驅動力的影響
裂紋圓管彎曲承載能力研究
裂紋敏感性鋼鑄坯表面質量控制
山東冶金(2019年6期)2020-01-06 07:45:58
Epidermal growth factor receptor rs17337023 polymorphism in hypertensive gestational diabetic women: A pilot study
42CrMo托輥裂紋的堆焊修復
山東冶金(2019年3期)2019-07-10 00:54:06
心生裂紋
揚子江(2019年1期)2019-03-08 02:52:34
Overcoming scarring in the urethra:Challenges for tissue engineering
微裂紋區對主裂紋擴展的影響
A7NO1鋁合金退火處理后焊接接頭疲勞裂紋擴展特性
焊接(2015年2期)2015-07-18 11:02:38
主站蜘蛛池模板: 一级毛片在线免费视频| 亚洲天堂色色人体| 精品天海翼一区二区| 日韩国产欧美精品在线| 在线国产综合一区二区三区| 国产手机在线小视频免费观看| 91精品日韩人妻无码久久| 久久无码av一区二区三区| 精品国产三级在线观看| 久久窝窝国产精品午夜看片| 国产高清在线精品一区二区三区 | 久久视精品| 国产成人综合网| 天天躁狠狠躁| 国产亚洲欧美另类一区二区| 国产成人禁片在线观看| 怡红院美国分院一区二区| 久久久久国产精品熟女影院| 伊人久久综在合线亚洲91| 久久精品无码一区二区国产区| 色综合激情网| 国产精品网址你懂的| 亚洲国产中文在线二区三区免| 巨熟乳波霸若妻中文观看免费 | 真人免费一级毛片一区二区| 国产成人无码久久久久毛片| 亚洲综合一区国产精品| 免费看黄片一区二区三区| 一边摸一边做爽的视频17国产| 91丨九色丨首页在线播放 | 婷婷伊人久久| 亚洲第一黄色网| 色香蕉网站| a级毛片免费网站| 亚洲欧美另类日本| 自拍偷拍欧美| 在线人成精品免费视频| www.91中文字幕| 国产99视频精品免费观看9e| 国产精品女同一区三区五区| 国产精品久久久久无码网站| 中文字幕亚洲无线码一区女同| 国产精品男人的天堂| 青青草欧美| 全色黄大色大片免费久久老太| 毛片在线看网站| www.av男人.com| 国产成人三级| 欧美日韩国产在线播放| 精品无码一区二区三区在线视频| 亚洲A∨无码精品午夜在线观看| 鲁鲁鲁爽爽爽在线视频观看| 99久久精品国产自免费| 亚洲综合精品第一页| 成人日韩欧美| 亚洲国产精品日韩欧美一区| 久久青草精品一区二区三区| 亚洲欧美综合精品久久成人网| 日韩精品成人在线| 亚洲动漫h| 国产剧情国内精品原创| www.91在线播放| 国产日韩欧美精品区性色| 91福利片| 多人乱p欧美在线观看| 小蝌蚪亚洲精品国产| 国产高清在线观看91精品| 日韩毛片基地| 欧美精品二区| 亚洲h视频在线| 国产精品女同一区三区五区| 日韩成人午夜| 久久国产成人精品国产成人亚洲 | 欧美在线精品怡红院| 天天色天天综合| 国产真实乱子伦精品视手机观看 | 婷婷亚洲最大| 久久免费精品琪琪| 色视频国产| 亚洲精品麻豆| 久久先锋资源| 亚洲av成人无码网站在线观看|