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

        ?

        運用改進的無單元Galerkin方法分析機場道面斷裂力學(xué)問題?

        2017-08-07 07:59:20鄒詩瑩席偉成彭妙娟程玉民
        物理學(xué)報 2017年12期
        關(guān)鍵詞:道面尖端機場

        鄒詩瑩 席偉成 彭妙娟 程玉民?

        1)(上海大學(xué)上海市應(yīng)用數(shù)學(xué)和力學(xué)研究所,上海 200072)

        2)(上海大學(xué)土木工程系,上海 200444)

        運用改進的無單元Galerkin方法分析機場道面斷裂力學(xué)問題?

        鄒詩瑩1)席偉成2)彭妙娟2)程玉民1)?

        1)(上海大學(xué)上海市應(yīng)用數(shù)學(xué)和力學(xué)研究所,上海 200072)

        2)(上海大學(xué)土木工程系,上海 200444)

        (2017年2月26日收到;2017年3月26日收到修改稿)

        在改進的無單元Galerkin方法的基礎(chǔ)上,將能反映裂紋尖端附近應(yīng)力奇異性的特征項引入改進的移動最小二乘法的基函數(shù)中,將斷裂力學(xué)和改進的無單元Galerkin方法結(jié)合,研究了線彈性斷裂力學(xué)的改進的無單元Galerkin方法,并對含反射裂縫的機場復(fù)合道面層狀體系結(jié)構(gòu)進行了數(shù)值分析.本文的理論為機場復(fù)合道面斷裂力學(xué)分析提供了一種新方法.

        改進的無單元Galerkin方法,機場復(fù)合道面,層狀體系,反射裂縫

        1 引 言

        無網(wǎng)格方法的研究已有二十多年的歷史,國際上目前將基于點的近似構(gòu)造試函數(shù)、不需要考慮節(jié)點與單元間關(guān)聯(lián)條件的各種數(shù)值方法稱為無網(wǎng)格方法[1,2].

        程玉民等[3,4]在移動最小二乘法的基礎(chǔ)上,建立了改進的移動最小二乘法. 引入復(fù)變量理論,程玉民等建立了復(fù)變量移動最小二乘法[5,6]、改進的復(fù)變量移動最小二乘法[7]和復(fù)變量重構(gòu)核粒子法[8].在此基礎(chǔ)上,提出了改進的無單元Galerkin方法[9,10]、插值型無單元Galerkin方法[11?14]、邊界無單元法[15?18]、復(fù)變量無單元Galerkin方法[19?23]、改進的復(fù)變量無單元Galerkin方法[7,24?26]和復(fù)變量重構(gòu)核粒子法[8,27?29].

        由于改進的移動最小二乘法形成的方程組易于求解并且不形成病態(tài)方程組,因而改進的無單元Galerkin方法可以提高無單元Galerkin方法的計算精度和計算效率.

        機場復(fù)合道面兼顧了剛性路面和柔性路面的優(yōu)點.水泥混凝土板提供了一個穩(wěn)定、堅實的基礎(chǔ),瀝青混凝土則作為表面功能特性優(yōu)良的磨耗層.然而,在水泥混凝土板間各種形式的接縫或裂縫使得加鋪的瀝青面層易產(chǎn)生反射裂縫.在瀝青混凝土加鋪層出現(xiàn)反射裂縫后,會使地表水進入瀝青混合料內(nèi)部,從而降低集料與瀝青之間的黏結(jié)性能,形成松散、剝落等病害.特別是在飛機輪載的作用下容易導(dǎo)致應(yīng)力集中現(xiàn)象,反射裂縫周邊道面又會形成新的交叉次生裂縫,發(fā)展到一定程度后將使道面出現(xiàn)大面積剝落,嚴重影響飛機的運行安全.反射裂縫是機場復(fù)合道面破壞的主要表現(xiàn)形式之一.

        關(guān)于機場復(fù)合道面反射裂縫,Garzon等[30,31]利用三維廣義有限元法對荷載作用下機場道面的反射裂縫進行了研究;李淑明等[32]采用有限元法分析了土工布、玻璃纖維網(wǎng)對舊水泥混凝土路面加鋪瀝青混凝土層復(fù)合路面結(jié)構(gòu)的內(nèi)力影響;凌建明等基于彈性層狀體系理論,利用有限元法對飛機荷載作用下機場道面地基附加應(yīng)力進行了分析[33],提出了機場復(fù)合道面剩余壽命預(yù)估方法[34];馬翔等以B777-200飛機荷載為計算荷載,利用有限元法分析了復(fù)合機場道面荷載應(yīng)力[35],提出了以荷載與溫度綜合疲勞彎拉應(yīng)力和瀝青面層反射裂縫疲勞壽命為設(shè)計指標的復(fù)合機場道面結(jié)構(gòu)設(shè)計方法[36].

        綜上研究進展,目前關(guān)于機場復(fù)合道面反射裂縫的計算均采用有限元法,尚未見到采用無網(wǎng)格方法分析機場復(fù)合道面反射裂縫,而無網(wǎng)格方法在求解斷裂問題方面優(yōu)于有限元法.因而,對機場復(fù)合道面來說,研究應(yīng)用于反射裂縫分析的斷裂力學(xué)的無網(wǎng)格方法是非常必要的.

        2 改進的移動最小二乘法

        在移動最小二乘法中,取函數(shù)u(x)的逼近函數(shù)為

        其中m是基函數(shù)的個數(shù),pi(x)是基函數(shù),ai(x)是相應(yīng)的系數(shù).

        對應(yīng)于(1)式的整體逼近,在點x鄰域內(nèi)的局部逼近定義為

        其中xI為影響域覆蓋點x的節(jié)點,n是影響域覆蓋點x的節(jié)點數(shù),uI=u(xI).

        對J取極值,即得

        可得

        對于點集{xi}和權(quán)函數(shù){wi},若一組函數(shù)φ1(x),φ2(x),···,φm(x)滿足

        則稱φ1(x),φ2(x),···,φm(x)是關(guān)于點集{xi}帶權(quán){wi}的正交函數(shù)族(k,j=1,2,···,m).

        若{pi(x)},i=1,2,···,m,為Hilbert空間span(p)上的關(guān)于點集{xi}的帶權(quán)的正交基函數(shù)族,即

        則方程(6)可寫成

        這樣,我們可以直接得到系數(shù)ai(x),即

        寫成矩陣形式

        其中

        將(10)式代入(2)式可得

        其中Φ?(x)為形函數(shù),

        這樣,系數(shù)ai(x)可以簡單、直接地得到,不需要求矩陣A(x)的逆,避免了求解病態(tài)或奇異的方程組,既提高了效率,又提高了精度.

        對如下基函數(shù)

        利用Schmidt正交化方法得到的正交基函數(shù)為

        3 裂紋尖端試函數(shù)

        基函數(shù)擴展法是通過在基函數(shù)中加入擴展函數(shù)項,來達到對試函數(shù)擴展的目的.對于斷裂力學(xué)的改進的無單元Galerkin方法,基函數(shù)中必須加入裂紋尖端位移場表達式中的一些重要項和一些重要的梯度項.基函數(shù)擴展的項數(shù)決定了數(shù)值解的精度,當基函數(shù)中加入了能反映包含裂紋尖端位移場中的所有項,即完全基函數(shù)擴展時,取如下形式:

        類似于改進的移動最小二乘法的推導(dǎo),可得

        類似(13)式可得裂紋尖端試函數(shù)為

        由此可以看出基函數(shù)擴展法沒有增加新的節(jié)點未知變量.但是由于基函數(shù)的項數(shù)增加了,形函數(shù)的計算量也會相應(yīng)地增加.為了提高計算效率,本文采用局部基函數(shù)擴展法,即在基函數(shù)中加入即

        這樣,可以得到裂紋尖端的試函數(shù)為

        可得裂紋尖端節(jié)點的形函數(shù)

        4 斷裂力學(xué)的改進的無單元Galerkin方法

        假設(shè)二維彈性體的求解域為?,邊界為Γ,彈性力學(xué)問題的平衡方程為

        其中σij為應(yīng)力分量,bi是單位體積上的體力分量,i,j=1,2.

        幾何方程為

        其中εij是應(yīng)變分量,ui是位移分量.

        本構(gòu)方程為

        式中D=(Dijkl)是彈性矩陣,對于平面應(yīng)力問題,

        對于平面應(yīng)變問題,

        其中E為材料的彈性模量.

        對應(yīng)的邊界條件為

        本節(jié)采用罰函數(shù)法施加本質(zhì)邊界條件,則彈性力學(xué)問題對應(yīng)的泛函為

        如果位移邊界條件是沿xi方向的,則si為0,否則為1.

        將求解域離散為有限個節(jié)點,節(jié)點總數(shù)為M.利用改進的移動最小二乘法建立逼近函數(shù),可得位移的逼近函數(shù)為

        對非裂紋尖端節(jié)點,形函數(shù)如(14)式所示;裂紋尖端節(jié)點的形函數(shù)見(24)式.

        將(37)式代入(32)式,可得

        矩陣Kα是根據(jù)節(jié)點信息得到的全局罰函數(shù)矩陣,

        向量Fα由本質(zhì)邊界條件得到,

        裂縫尖端的應(yīng)力強度因子采用以下方法進行計算.對裂紋尖端選取適當?shù)膢OC|=L得點C,|OB|=L/4得點B.由點O,B,C即構(gòu)成二維問題裂紋尖端的虛擬奇異二次元,裂紋尖端局部坐標如圖1所示,應(yīng)力強度因子表達式為

        其中a為裂紋長度的一半,vB和vC分別為點B和點C豎直方向的位移.

        圖1 裂紋尖端虛擬奇異元Fig.1.The virtual singular element at the tip of a crack.

        5 機場道面反射裂縫的數(shù)值分析

        5.1 機場復(fù)合道面反射裂縫的計算模型

        為了進行反射裂縫影響因素分析,將復(fù)合機場道面結(jié)構(gòu)簡化為平面應(yīng)變問題的層狀體系,采用二維模型進行數(shù)值計算.將荷載簡化為二維平面應(yīng)變問題下的條形均布荷載,采用B777-200B,輪壓為1.45 MPa,輪距為1.40 m.此外假定各結(jié)構(gòu)層都由線彈性的各向同性、均質(zhì)材料組成,地基由彈性半空間地基假設(shè),在模型中采用有限尺寸.模型的邊界條件為:地基底部完全約束,各結(jié)構(gòu)層兩側(cè)鉸接,具體幾何參數(shù)如圖2,各層的物理參數(shù)見表1,反射裂縫起裂位置如圖3,裂縫寬度取0.5 cm.

        圖2 機場復(fù)合道面模型幾何參數(shù)(mm)Fig.2.The geometric parameters of the airport composite pavement(mm).

        表1 結(jié)構(gòu)的物理參數(shù)Table 1.The physical parameters of the structure.

        圖3 反射裂縫的位置Fig.3.The position of the re fl ective crack.

        5.2 機場復(fù)合道面反射裂縫力學(xué)分析

        首先考慮起裂位置處于接縫處的裂縫計算模型.無網(wǎng)格方法的節(jié)點布置采用5370個節(jié)點,輪載下方區(qū)域節(jié)點進行加密,并對起裂位置進行局部加密,沿水泥板接縫位置布置2×27個加密節(jié)點,具體節(jié)點分布和局部加密如圖4和圖5.

        圖4 模型節(jié)點布置Fig.4.The node distribution of the model.

        圖5 裂縫附近節(jié)點加密分布(2×27)Fig.5.The node re fi nement distribution at the tip of the crack(2×27).

        計算結(jié)束后,選取裂紋尖端適當計算點,采用(47)式計算,得到此處應(yīng)力強度因子KI為9.42 MPa·mm1/2,裂縫尖端位于混凝土板層,回彈模量為30000 MPa.

        在此基礎(chǔ)上,進一步考慮裂縫反射到瀝青加鋪層的反射裂縫模型,如圖6取反射至瀝青面層的裂縫長度為4 cm,并在裂紋區(qū)域進行節(jié)點加密(3×31個加密節(jié)點),其節(jié)點布置局部如圖7所示.

        通過計算得到的KI為0.36 MPa·mm1/2,裂縫尖端位于瀝青層,回彈模量為2000 MPa.

        圖6 反射至瀝青層的裂縫Fig.6.The re fl ective crack in the asphalt layer.

        圖7 反射裂縫附近節(jié)點加密分布(3×31)Fig.7.The node re fi nement distribution at the tip of the re fl ective crack(3×31).

        圖8 瀝青層計算點豎向位移Fig.8.The vertical displacement at the points in the asphalt layer.

        選取計算點后對模型的位移和應(yīng)力進行計算,計算結(jié)果如圖8—圖10所示.從圖8可以看出,在荷載作用位置附近,豎向位移較大,與路面實際情況符合.從圖9可以看出,當混泥土板開裂后,在荷載作用下裂縫正上方面層的正應(yīng)力具有奇異性,從而使得路面反射裂縫繼續(xù)向上發(fā)展,以致最后貫通整個路面.可見,當基層出現(xiàn)裂縫后,由于裂縫尖端存在一定的應(yīng)力集中,這種應(yīng)力集中對于路面正常工作是不利的.所以,路面工程施工中要盡可能減少混凝土板的收縮裂縫,以防止反射裂縫發(fā)生.

        從圖9和圖10可以看出,在對稱荷載作用下,裂縫尖端正應(yīng)力遠大于剪應(yīng)力,裂縫為張開型,即裂縫的擴展受正應(yīng)力或KI的控制.

        圖9 輪載下各層計算點正應(yīng)力σ11Fig.9.The normal stress σ11at the points in layers under the wheel load.

        圖10 輪載下各層計算點剪應(yīng)力τ21Fig.10.The shear stress τ21at the points in layers under the wheel load.

        6 結(jié) 論

        本文針對機場復(fù)合道面的反射裂縫問題,建立了斷裂力學(xué)改進的無單元Galerkin方法.數(shù)值算例表明,該方法通過對試函數(shù)的擴展能有效地反映裂紋尖端應(yīng)力場的奇異性,提高了無網(wǎng)格方法的求解精度.

        在對試函數(shù)進行改進時,采用基函數(shù)的局部擴展法,不僅能有效地提高對裂紋尖端應(yīng)力場的計算精度,也比試函數(shù)外部擴展法以及基函數(shù)完全擴展法具有更好的求解效率.

        使用改進無單元Galerkin方法對機場復(fù)合道面模型進行計算,在布置節(jié)點時,在裂紋局部進行加密布置.對于計算裂紋擴展的模型時,可以在擴展路徑的區(qū)域內(nèi)增加布置相應(yīng)的節(jié)點,而無需進行有限元法的網(wǎng)格重構(gòu).

        本文的理論為機場復(fù)合道面反射裂縫的分析提供了一種新的方法.

        [1]Cheng Y M 2015 Meshless Methods(Beijing:Science Press)pp1–13(in Chinese)[程玉民2015無網(wǎng)格方法 (北京:科學(xué)出版社)第1—13頁]

        [2]Cheng Y M,Ji X,He P F 2004 Acta Mech.Sin.36 43(in Chinese)[程玉民,嵇醒,賀鵬飛 2004力學(xué)學(xué)報 36 43]

        [3]Cheng Y M,Chen M J 2003 Acta Mech.Sin.35 181(in Chinese)[程玉民,陳美娟 2003力學(xué)學(xué)報 35 181]

        [4]Cheng Y M,Liew K M,Kitipornchai S 2009 Int.J.Numer.Meth.Eng.78 1258

        [5]Cheng Y M,Peng M J,Li J H 2005 Chin.J.Theor.Appl.Mech.37 719(in Chinese)[程玉民,彭妙娟,李九紅2005應(yīng)用力學(xué)學(xué)報37 719]

        [6]Cheng Y M,Li J H 2005 Acta Phys.Sin.54 4463(in Chinese)[程玉民,李九紅 2005物理學(xué)報 54 4463]

        [7]Bai F N,Li D M,Wang J F,Cheng Y M 2012 Chin.Phys.B 21 020204

        [8]Chen L,Cheng Y M 2008 Acta Phys.Sin.57 1(in Chinese)[陳麗,程玉民 2008物理學(xué)報 57 1]

        [9]Peng M J,Li R X,Cheng Y M 2014 Eng.Anal.Bound.Elem.40 104

        [10]Cheng R J,Cheng Y M 2016 Chin.Phys.B 25 020203

        [11]Sun F X,Wang J F,Cheng Y M 2013 Chin.Phys.B 22 120203

        [12]Cheng Y M,Bai F N,Peng M J 2014 Appl.Math.Modell.38 5187

        [13]Cheng Y M,Bai F N,Liu C,Peng M J 2016 Int.J.Comput.Mater.Sci.Eng.5 1650023

        [14]Sun F X,Wang J F,Cheng Y M 2016 Int.J.Appl.Mech.8 1650096

        [15]Cheng Y M,Peng M J 2005 Sci.China Ser.G 48 641

        [16]Peng M J,Cheng Y M 2009 Eng.Anal.Bound.Elem.33 77

        [17]Ren H P,Cheng Y M,Zhang W 2009 Chin.Phys.B 18 4065

        [18]Wang J F,Wang J F,Sun F X,Cheng Y M 2013 Int.J.Comput.Meth.10 1350043

        [19]Cheng Y M,Li J H 2006 Sci.China Ser.G 49 46

        [20]Peng M J,Li D M,Cheng Y M 2011 Eng.Struct.33 127

        [21]Li D M,Peng M J,Cheng Y M 2011 Sci.Sin.:Phys.Mech.Astron.41 1003(in Chinese)[李冬明,彭妙娟,程玉民2011中國科學(xué):物理學(xué)力學(xué) 天文學(xué)41 1003]

        [22]Cheng Y M,Li R X,Peng M J 2012 Chin.Phys.B 21 090205

        [23]Cheng Y M,Wang J F,Li R X 2012 Int.J.Appl.Mech.4 1250042

        [24]Cheng Y M,Wang J F,Bai F N 2012 Chin.Phys.B 21 090203

        [25]Deng Y J,Liu C,Peng M J,Cheng Y M 2015 Int.J.Appl.Mech.7 1550017

        [26]Cheng Y M,Liu C,Bai F N,Peng M J 2015 Chin.Phys.B 24 100202

        [27]Chen L,Ma H P,Cheng Y M 2013 Chin.Phys.B 22 050202

        [28]Weng Y J,Cheng Y M 2013 Chin.Phys.B 22 090204

        [29]Chen L,Cheng Y M,Ma H P 2015 Comput.Mech.55 591

        [30]Garzon J,Duarte C A,Buttlar W 2010 Road Mater.Pavement Design 11 459

        [31]Garzon J,Kim D,Duarte C A 2013 Int.J.Comput.Meth.10 1350045

        [32]Li S M,Cai X M,Xu Z H 2005 J.Tongji Univ.(Nat.Sci.)33 1616(in Chinese)[李淑明,蔡喜棉,許志鴻 2005同濟大學(xué)學(xué)報(自然科學(xué)版)33 1616]

        [33]Wo R H,Ling J M 2001 J.Tongji Univ.29 288(in Chinese)[咼潤華,凌建明2001同濟大學(xué)學(xué)報29 288]

        [34]Zhou Z F,Ling J M,Yuan J 2007 J.Traffic Transport.Eng.7 50(in Chinese)[周正峰,凌建明,袁捷 2007交通運輸工程學(xué)報7 50]

        [35]Ma X,Ni F J,Chen R S 2010 J.Changan Univ.(Nat.Sci.)30 23(in Chinese)[馬翔,倪富健,陳榮生2010長安大學(xué)學(xué)報(自然科學(xué)版)30 23]

        [36]Ma X,Ni F J,Gu X Y 2010 J.Traffic Transport.Eng.10 36(in Chinese)[馬翔,倪富健,顧興宇2010交通運輸工程學(xué)報10 36]

        PACS:02.60.Cb,02.60.Lj,46.50.+aDOI:10.7498/aps.66.120204

        Analysis of fracture problems of airport pavement by improved element-free Galerkin method?

        Zou Shi-Ying1)Xi Wei-Cheng2)Peng Miao-Juan2)Cheng Yu-Min1)?

        1)(Shanghai Institute of Applied Mathematics and Mechanics,Shanghai University,Shanghai 200072,China)
        2)(Department of Civil Engineering,Shanghai University,Shanghai 200444,China)

        26 February 2017;revised manuscript

        26 March 2017)

        Using the improved element-free Galerkin(IEFG)method,in this paper we introduce the characteristic parameterwhich can re fl ect the singular stress near the crack tip into the basic function of the improved moving least-squares(IMLS)approximation.Combining fracture theory with the IEFG method,we present an IEFG method of treating the elastic fracture problems,and analyze a numerical example of two-dimensional layered system of airport composite pavement with re fl ective crack.

        In the IEFG method,the IMLS approximation is used to form the shape function.The IMLS approximation is presented from the moving least-squares(MLS)approximation,which is the basis of the element-free Galerkin(EFG)method.Compared with the MLS approximation,the IMLS approximation uses the orthonormal basis functions to obtain the shape function,which leads to the fact that the matrices for obtaining the undetermined coefficients are diagonal.Then the IMLS approximation can obtain the solutions of the undetermined coefficients directly without the inverse matrices.The IMLS approximation can overcome the disadvantages of the MLS approximation,in which the ill-conditional or singular matrices are formed sometimes.And it can also improve the computational efficiency of the MLS approximation.

        Because of the advantages of the IMLS approximation,the IEFG method has greater computational efficiency than the EFG method which is based on the MLS approximation,and can obtain the solution for arbitrary node distribution,even though the EFG method cannot obtain the solution due to the ill-conditional or singular matrices in the MLS approximation.

        Paving the asphalt concrete layer on the cement concrete pavement is an e ff ective approach to improving the structure and service performance of an airport pavement,which is called airport composite pavement.The airport composite pavement has the advantages of rigid pavement and fl exible pavement,but there are various forms of joints or cracks of cement concrete slab,which makes the crack re fl ect into the asphalt overlay easily under the plane load and environmental factors.Re fl ective crack is one of the main failure forms of the airport composite pavement.Therefore,it is of great theoretical signi fi cance and engineering application to study the generation and development mechanism of re fl ective crack of the airport composite pavement.

        For the numerical methods of solving the fracture problems,introducing the characteristic parameterwhich can re fl ect the singular stress near the crack tip into the basic function is a general approach.In this paper,we use this approach to obtain the IEFG method for fracture problems,and the layered system of airport composite pavement with re fl ective crack is considered.The numerical results of the displacements and stresses in the airport composite pavement are given.And at the tip of the crack,the stress is singular,which makes the crack of the airport composite pavement grow.

        This paper provides a new method for solving the re fl ective crack problem of airport composite pavement.

        improved element-free Galerkin method,airport composite pavement,layered system,re fl ective crack

        10.7498/aps.66.120204

        ?國家自然科學(xué)基金委員會-中國民航局民航聯(lián)合研究基金(批準號:U1433104)資助的課題.

        ?通信作者.E-mail:ymcheng@shu.edu.cn

        ?2017中國物理學(xué)會Chinese Physical Society

        http://wulixb.iphy.ac.cn

        *Project supported by the National Natural Science Foundation of China(Grant No.U1433104).

        ?Corresponding author.E-mail:ymcheng@shu.edu.cn

        猜你喜歡
        道面尖端機場
        機場罷工
        淺談機場跑道道面性能變化規(guī)律
        通航機場跑道道面承載能力分析
        如何避免GSM-R無線通信系統(tǒng)對機場電磁干擾
        郭紹?。核枷肱鲎苍炀图舛巳瞬?/a>
        面部識別使機場安檢提速
        濕滑跑道飛機著陸輪胎-水膜-道面相互作用
        最有創(chuàng)意的機場
        鏡頭看展
        基于位移相關(guān)法的重復(fù)壓裂裂縫尖端應(yīng)力場研究
        斷塊油氣田(2014年5期)2014-03-11 15:33:49
        少妇av免费在线播放| 免费a级毛片在线播放不收费| 日韩av一区二区网址| 亚洲欧美中文字幕5发布| 中文字幕日本特黄aa毛片| 久久国产品野战| 国产精品一区二区久久精品蜜臀| 国内精品嫩模av私拍在线观看| 国产精品国产三级农村妇女| 青草久久婷婷亚洲精品| 在线亚洲妇色中文色综合| 国产精品国产三级国产专播下| 无码人妻久久久一区二区三区 | 久久亚洲中文字幕精品一区 | 日本韩国亚洲三级在线| 亚洲gay片在线gv网站| 极品av麻豆国产在线观看| 久久狠狠高潮亚洲精品暴力打 | 麻豆久久久9性大片| 无码日韩人妻AV一区免费| 91亚洲最新国语中文字幕| 国产一区二区三区在线影院| 中文字幕一区二区三区视频| 内射白浆一区二区在线观看| 少妇aaa级久久久无码精品片| www国产无套内射com| 亚洲aⅴ久久久噜噜噜噜| 日产一区一区三区区别| 日日碰狠狠添天天爽超碰97久久| 国产人妻熟女高跟丝袜图片| 亚洲 暴爽 av人人爽日日碰| 亚洲日日噜噜噜夜夜爽爽| 毛片成人18毛片免费看| 日本边添边摸边做边爱喷水| 特级婬片国产高清视频| 精品久久亚洲一级α| 久久久精品国产老熟女| 伊人久久大香线蕉午夜av | 99久久99久久精品免费看蜜桃| 久久精品国产亚洲av四虎| 四虎影视在线观看2413|