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

        ?

        基于微觀圖像及內(nèi)聚力模型的復(fù)合材料裂紋擴(kuò)展模擬①

        2017-09-15 09:14:47秦麗媛孟松鶴李金平
        固體火箭技術(shù) 2017年4期
        關(guān)鍵詞:內(nèi)聚力晶界微觀

        秦麗媛,孟松鶴,李金平,金 華

        (哈爾濱工業(yè)大學(xué) 復(fù)合材料與結(jié)構(gòu)研究所,哈爾濱 150008)

        基于微觀圖像及內(nèi)聚力模型的復(fù)合材料裂紋擴(kuò)展模擬①

        秦麗媛,孟松鶴,李金平,金 華

        (哈爾濱工業(yè)大學(xué) 復(fù)合材料與結(jié)構(gòu)研究所,哈爾濱 150008)

        為研究材料微觀結(jié)構(gòu)及晶界強(qiáng)度對(duì)材料力學(xué)性能的影響,在晶界處引入內(nèi)聚力單元模型,模擬晶間破壞過程。以ZrB2-SiC復(fù)合材料為研究對(duì)象,將其掃描的微觀結(jié)構(gòu)圖片進(jìn)行矢量化處理,并導(dǎo)入ABAQUS有限元軟件中建立模型,同時(shí)在其晶界處,設(shè)置內(nèi)聚力單元模擬晶界破壞過程。通過改變ZrB2與SiC相界面強(qiáng)度,得到了晶界及材料不均勻?qū)Σ牧蠎?yīng)力分布及裂紋擴(kuò)展的影響。結(jié)果表明,由于晶界的存在,材料內(nèi)部出現(xiàn)應(yīng)力分布不均勻現(xiàn)象并產(chǎn)生應(yīng)力集中。隨著晶界強(qiáng)度的改變,裂紋起始位置及擴(kuò)展方向發(fā)生改變,且裂紋沿低強(qiáng)度的界面進(jìn)行擴(kuò)展。隨著ZrB2-SiC界面強(qiáng)度增大,材料的強(qiáng)度提高,拉伸模量不變。

        微觀結(jié)構(gòu);內(nèi)聚力單元法;晶界強(qiáng)度;ZrB2-SiC陶瓷

        0 引言

        復(fù)合材料是兩種或兩種以上的組分組成的多相材料,且各相之間具有明顯的界面。由于復(fù)合材料組成相的物理化學(xué)性質(zhì)不同,當(dāng)材料受到外力作用時(shí),其對(duì)外力的響應(yīng)也存在很大的不同。因此,復(fù)合材料的微觀結(jié)構(gòu)對(duì)其力學(xué)性能,如應(yīng)力分布、裂紋擴(kuò)展以及破壞模式有很大的影響。

        目前,國(guó)內(nèi)外很多研究通過改變材料的微觀結(jié)構(gòu),進(jìn)而得到理想的材料性能。以ZrB2-SiC陶瓷為例,可在ZrB2基體中添加第二相,并改變第二相形貌,來提高材料性能。Monteverde F在ZrB2陶瓷中引入10 % αSiC晶粒,SiC顆粒可控制硼化物晶粒的長(zhǎng)大,并達(dá)到致密化,從而提高其斷裂韌性[1]。Laura Silvestroni[2]和DilettaSciti[3]研究了SiC晶須和SiC纖維對(duì)ZrB2基陶瓷的增韌機(jī)制。發(fā)現(xiàn)在高溫的條件下,SiC晶須發(fā)生破壞,但可提高其室溫強(qiáng)度。而SiC纖維的加入,可增強(qiáng)陶瓷材料的韌性,且其增韌效果優(yōu)于晶須。王明福等研究了SiC晶片(SiCpl)對(duì)ZrB2復(fù)合材料的影響。隨著SiCpl的加入,SiCpl/ ZrB2復(fù)合材料相對(duì)致密度和斷裂韌性都得到了較大的提高[4]。張幸紅等還研究了在ZrB2中加入其他材料[5],如碳短纖維[6]或石墨片[7-8]對(duì)材料力學(xué)性能的影響。因此,材料的微觀結(jié)構(gòu)對(duì)其力學(xué)性能影響很大。然而,這些微結(jié)構(gòu)對(duì)材料力學(xué)性能的影響,大多是通過宏觀力學(xué)試驗(yàn)嘗試對(duì)其破壞機(jī)制進(jìn)行定性分析得到的,而缺少定量的計(jì)算。

        一般情況下,裂紋是引起各種結(jié)構(gòu)及零部件失效根源。因此,發(fā)現(xiàn)各種裂紋現(xiàn)象、了解裂紋擴(kuò)展及失穩(wěn)擴(kuò)展的條件、掌握裂紋擴(kuò)展的規(guī)律及控制裂紋的擴(kuò)展非常有必要。為了研究裂紋萌生及擴(kuò)展的真正物理機(jī)理,必須從介觀或微觀尺度下進(jìn)行數(shù)值計(jì)算。隨著計(jì)算機(jī)科學(xué)的快速發(fā)展,有限元法作為工程分析的有效方法在研究材料微觀結(jié)構(gòu)中得到廣泛的應(yīng)用。在傳統(tǒng)有限元基礎(chǔ)上,又發(fā)展了處理非連續(xù)問題的內(nèi)聚力單元模型。內(nèi)聚力模型最先由Dugdale和Barenblatt在1960年提出[9]。1976年,Hillerborg等又加入了拉伸強(qiáng)度,已經(jīng)允許已有裂紋的增長(zhǎng),并且還允許新裂紋的萌生和演化[10]。Xu和Needleman等又提出了粘接則[11],進(jìn)一步完善了該模型。Horacio D Espinosa等使用內(nèi)聚力單元研究了多晶Al2O3中的裂紋萌生及擴(kuò)展,并研究了內(nèi)聚力單元中各參數(shù)的影響[12]。

        由于陶瓷材料本身強(qiáng)度很大,在燒結(jié)過程中由于雜質(zhì)和晶界的存在,往往晶界處成為材料破壞的弱項(xiàng),而產(chǎn)生晶間斷裂。本文將通過直接對(duì)材料的微觀結(jié)構(gòu)進(jìn)行建模,模擬材料內(nèi)部的應(yīng)力分布情況,并使用內(nèi)聚力單元對(duì)晶界進(jìn)行模擬,得到裂紋的形成及擴(kuò)展。以ZrB2-SiC陶瓷為研究對(duì)象,基于其微觀形貌直接建立有限元模型。由于ZrB2-SiC不同相之間的匹配程度及熱膨脹系數(shù)不同,會(huì)造成ZrB2-ZrB2與ZrB2-SiC的界面強(qiáng)度不同。為研究晶界強(qiáng)度對(duì)材料力學(xué)性能及裂紋擴(kuò)展的影響,在其晶界處設(shè)置內(nèi)聚力單元,并改變ZrB2-SiC晶界的強(qiáng)度進(jìn)行計(jì)算,得到不同的界面強(qiáng)度對(duì)材料的整體強(qiáng)度、拉伸模量及裂紋擴(kuò)展的影響。

        1 基于圖像建模

        基于數(shù)字圖像建模技術(shù),在生物技術(shù)領(lǐng)域由Holister等[13]首次提出用于模擬人類骨骼的多孔微結(jié)構(gòu)。這種建模方法廣泛用于巖土工程材料研究領(lǐng)域[14-16]。數(shù)字圖像處理技術(shù)可直接獲取材料內(nèi)部微觀結(jié)構(gòu),得到其內(nèi)部組織的分布特征。

        1.1 圖像獲取

        將熱壓燒結(jié)ZrB2-SiC陶瓷材料使用不同型號(hào)砂紙進(jìn)行研磨、拋光。將拋光后的試樣放入熔融的強(qiáng)堿中進(jìn)行堿腐蝕。再通過掃描電鏡,得到材料的數(shù)值圖像。

        圖1為ZrB2-SiC陶瓷的微觀結(jié)構(gòu)形貌圖,其中,黑色部分為SiC,灰色部分為ZrB2。

        1.2 晶界提取

        將掃描圖像進(jìn)行矢量化處理。將晶界的明暗分界面的像素點(diǎn)數(shù)據(jù)轉(zhuǎn)換成幾何矢量信息,進(jìn)行有限元計(jì)算。微觀結(jié)構(gòu)圖的不同晶粒,由于其反光性質(zhì)不同,呈現(xiàn)出不同的灰度,且對(duì)進(jìn)行腐蝕后的微觀結(jié)構(gòu)呈現(xiàn)出不同的灰度,可明顯地看出晶界的形貌。根據(jù)其灰度的不同,將材料微觀結(jié)構(gòu)進(jìn)行矢量化處理,從而將晶界提取出來,且對(duì)不同組成相在計(jì)算時(shí),賦予不同的材料屬性。對(duì)于不太清晰的圖片,應(yīng)先進(jìn)行圖片處理,增加其明暗對(duì)比度,以便得到清晰的晶界輪廓。圖片處理晶界處理過程如圖2所示。為方便計(jì)算,選取圖1方塊部分進(jìn)行圖片處理及計(jì)算。由掃描圖片可很清楚地分辨出ZrB2-SiC多晶陶瓷材料中不同相間的顏色差異,但由于ZrB2晶粒間的對(duì)比度較低,晶界不是很明顯。將圖片導(dǎo)入Image-Pro Plus軟件,進(jìn)行增加對(duì)比度處理,如圖2(b)所示。然后,使用AlgoLab photo vector,將增強(qiáng)對(duì)比度后的圖像進(jìn)行矢量化處理,如圖2(c)所示。

        由于材料微觀形貌圖中晶界處存在些陰影,因此矢量化后的圖像晶界處存在一些噪音。需要將矢量化后的圖像,根據(jù)原始圖像在CAD中對(duì)晶界進(jìn)行修整,以便進(jìn)行計(jì)算。得到修整后的矢量圖,便可導(dǎo)入有限元軟件進(jìn)行計(jì)算,如圖2(d)所示。

        2 內(nèi)聚力單元模型

        內(nèi)聚力模型避免了線彈性斷裂力學(xué)中的裂紋尖端應(yīng)力奇異性,而被廣泛用于斷裂力學(xué)模擬中。內(nèi)聚力損傷模型通過界面周圍材料之間內(nèi)聚力關(guān)系來模擬裂紋的開裂和擴(kuò)展。一般陶瓷材料具有晶體結(jié)構(gòu),且多數(shù)的破壞為晶間斷裂。因此,可通過在晶間設(shè)置內(nèi)聚力單元來模擬材料的穿晶破壞過程。

        將矢量化的微觀圖像導(dǎo)入ABAQUS有限元軟件中。根據(jù)材料本身的性能對(duì)實(shí)體材料部分賦予不同材料屬性。在晶界處設(shè)置內(nèi)聚力單元,可將其看作被一個(gè)厚度隔開的兩個(gè)面,這兩個(gè)面分別和其他實(shí)體單元連接,如圖3所示。

        內(nèi)聚力單元根據(jù)力-位移模型,來描述材料失效。當(dāng)內(nèi)聚力單元的應(yīng)力或應(yīng)變狀態(tài)滿足損傷起始準(zhǔn)則后,單元將發(fā)生損傷。在界面應(yīng)力或應(yīng)變達(dá)到最大值以前,兩個(gè)晶粒間的界面是完整的。一旦達(dá)到最大拉力,界面開始惡化,且拉力減小至零,界面被認(rèn)為破壞,晶界產(chǎn)生微裂紋。

        3 計(jì)算

        根據(jù)掃描的微觀結(jié)構(gòu)圖,不同組成相及界面相賦予不同的材料性能,組成相的楊氏模量及泊松比如表1所示。

        為有效模擬界面開裂過程,選用了ABAQUS中二維四節(jié)點(diǎn)內(nèi)聚力單元(COH2D4)來定義界面單元。界面內(nèi)聚力單元采用最大應(yīng)力準(zhǔn)則作為損傷起始準(zhǔn)則。當(dāng)內(nèi)聚力單元的應(yīng)力滿足損傷起始準(zhǔn)則后,單元進(jìn)入損傷演化階段。本文損傷演化的力和位移模型使用雙線性損傷演化模型,如圖4所示。

        表1 ZrB2和SiC的楊氏模量E及泊松比ν[17]

        (1)

        圖4中的三角形面積即為界面單元完全破壞過程中的應(yīng)變能釋放率Gc,即

        (2)

        由式(1)和式(2)得到界面單元最大分離位移δ為0.2 μm。

        為研究不同晶界結(jié)合強(qiáng)度對(duì)材料力學(xué)性能及裂紋擴(kuò)展的影響,改變ZrB2-SiC界面的最大牽引力。由材料的拉伸強(qiáng)度值,選取ZrB2-ZrB2界面最大牽引力為300 MPa,選取ZrB2-SiC相界面的最大牽引力為100、200、300、400、500 MPa。

        對(duì)模型施加單軸拉伸作用,模擬材料整體宏觀力學(xué)行為。模型右端施加均勻的位移邊界條件,模型左端固定,兩側(cè)邊界自由。在位移邊界條件作用下,通過提取下邊界所有節(jié)點(diǎn)沿加載方向的平均反力,求得復(fù)合材料整體宏觀應(yīng)力-應(yīng)變關(guān)系曲線,從而得到其強(qiáng)度及拉伸模量。

        4 模擬結(jié)果及討論

        圖5所示為基于微觀結(jié)構(gòu)的數(shù)字圖像,并在晶界處插入內(nèi)聚力單元的多晶模型,在單軸拉伸載荷下的有限元計(jì)算結(jié)果。由于對(duì)不同的組成相賦予了不同的材料屬性,從圖5中可看出,材料受到載荷作用時(shí),呈現(xiàn)出非均勻受力情況,并在晶界處產(chǎn)生應(yīng)力集中。如圖5中箭頭所示,在ZrB2晶粒中部分應(yīng)力約為383 MPa,并整個(gè)晶粒受力均勻。而SiC晶粒邊緣應(yīng)力為778 MPa,如圖5中箭頭所示。由此可見,由材料不均勻及顆粒形狀不規(guī)則造成的應(yīng)力集中很明顯。

        保持ZrB2-ZrB2相界面強(qiáng)度不變,對(duì)ZrB2-SiC相界面選取不同的強(qiáng)度值進(jìn)行計(jì)算,得到材料的應(yīng)力-應(yīng)變曲線,如圖6所示。材料的強(qiáng)度和模量隨界面強(qiáng)度的變化曲線如圖7所示。材料的拉伸模量不隨強(qiáng)度的改變發(fā)生改變,計(jì)算得到的材料模量為420 GPa,且不隨ZrB2-SiC相界面強(qiáng)度的變化有所改變,并小于組成相及界面的模量。

        整體模型的強(qiáng)度隨ZrB2-SiC界面的強(qiáng)度的增大而增大,如圖7所示。當(dāng)ZrB2-SiC界面的強(qiáng)度為100 MPa時(shí),與材料整體的強(qiáng)度234 MPa相差較大。材料的整體強(qiáng)度更接近于ZrB2-ZrB2界面的強(qiáng)度300 MPa。當(dāng)界面強(qiáng)度為500 MPa時(shí),與材料整體的強(qiáng)度535 MPa較接近。材料整體強(qiáng)度是不同相界面強(qiáng)度共同作用的結(jié)果。

        圖8給出了模型中ZrB2-SiC界面強(qiáng)度不同時(shí),裂紋隨拉伸位移U變化的擴(kuò)展情況。

        當(dāng)模型中ZrB2-SiC界面強(qiáng)度為100 MPa時(shí),裂紋首先在ZrB2-SiC界面產(chǎn)生,且容易沿低強(qiáng)度的ZrB2-SiC界面進(jìn)行擴(kuò)展。當(dāng)模型中ZrB2-SiC界面強(qiáng)度為300 MPa時(shí),裂紋起始位置發(fā)生改變,由于ZrB2-SiC界面強(qiáng)度與ZrB2-ZrB2界面強(qiáng)度相同,裂紋擴(kuò)展傾向于沿小的偏轉(zhuǎn)角方向擴(kuò)展。當(dāng)ZrB2-SiC界面強(qiáng)度為500 MPa時(shí),由于對(duì)不同的組成相賦予不同的材料屬性,且各相形狀不同,容易產(chǎn)生應(yīng)力集中,如圖5頂端SiC顆粒處,裂紋仍在ZrB2-SiC界面產(chǎn)生。裂紋沿低強(qiáng)度ZrB2-ZrB2界面進(jìn)行擴(kuò)展。

        5 結(jié)論

        使用基于數(shù)字圖像技術(shù)進(jìn)行有限元建模,并在晶界處插入內(nèi)聚力單元進(jìn)行計(jì)算,可很好地模擬材料的不均勻性對(duì)材料內(nèi)部應(yīng)力分布情況的影響,得到材料內(nèi)部更真實(shí)的受力情況。由于不同相之間的匹配程度及熱膨脹系數(shù)不同,會(huì)造成ZrB2-ZrB2與ZrB2-SiC的界面強(qiáng)度不同。分別選取不同的ZrB2-SiC的界面強(qiáng)度,得到其對(duì)材料的整體強(qiáng)度、拉伸模量及裂紋擴(kuò)展的影響。

        (1)由于材料組成相的異質(zhì)性及晶界的存在,晶粒邊緣及晶界處出現(xiàn)很嚴(yán)重的應(yīng)力集中。

        (2)材料的拉伸模量不隨ZrB2-SiC晶面強(qiáng)度的變化發(fā)生改變,而整體模型的強(qiáng)度隨晶界強(qiáng)度的增大而增大,且高于晶界強(qiáng)度。

        (3)晶界強(qiáng)度的變化對(duì)裂紋的起始位置及擴(kuò)展方向有重要影響。由于組成相的異質(zhì)性及晶體形狀不同,造成材料內(nèi)部應(yīng)力集中,裂紋的起始位置具有偶然性。裂紋容易沿低強(qiáng)度晶界進(jìn)行擴(kuò)展。

        [1] Monteverde F.Beneficial effects of an ultra-fine α-SiC incorporation on the sinterability and mechanical proper-ties of ZrB2[J].Applied Physics A,2006,82(2):329-337.

        [2] Silvestroni L,Sciti D,Melandri C,et al.Toughened ZrB2-based ceramics through SiC whisker or SiC chopped fiber additions[J].Journal of the European Ceramic Society,2010,30(11):2155-2164.

        [3] Sciti D,Pienti L,Dalle Fabbriche D,et al.Combined effect of SiC chopped fibers and SiC whiskers on the toughening of ZrB2[J].Ceramics International,2014,40(3):4819-4826.

        [4] 王明福,汪長(zhǎng)安,尉磊,等.SiC 晶片增韌ZrB2復(fù)合陶瓷材料的制備與性能[J].稀有金屬材料與工程,2009,38(2):913-915.

        [5] Sun X,Han W,Hu P,et al.Microstructure and mechanical properties of ZrB2-Nb composite[J].International Journal of Refractory Metals and Hard Materials,2010,28(3):472-474.

        [6] 吳世平,楊飛宇,張幸紅,等.碳短纖維對(duì)ZrB2-SiC 基超高溫陶瓷力學(xué)性能的影響[J].材料工程,2007(5):15-18.

        [7] Zhang X,Wang Z,Sun X,et al.Effect of graphite flake on the mechanical properties of hot pressed ZrB2-SiC ceramics[J].Materials Letters,2008,62(28):4360-4362.

        [8] Wang Z,Wang S,Zhang X,et al.Effect of graphite flake on microstructure as well as mechanical properties and thermal shock resistance of ZrB2-SiC matrix ultrahigh temperature ceramics[J].Journal of Alloys and Compounds,2009,484(1):390-394.

        [9] Dugdale D S.Yielding of steel sheets containing slits[J].Journal of the Mechanics and Physics of Solids,1960,8(2):100-104.

        [10] Hillerborg A,Modéer M,Petersson P E.Analysis of crack formation and crack growth in concrete by means of fracture mechanics and finite elements[J].Cement and Concrete Research,1976,6(6):773-781.

        [11] Xu X P,Needleman A.Numerical simulations of fast crack growth in brittle solids[J].Journal of the Mechanics and Physics of Solids,1994,42(9):1397-1434.

        [12] Espinosa H D,Zavattieri P D.A grain level model for the study of failure initiation and evolution in polycrystalline brittle materials.Part II:numerical examples[J].Mechanics of Materials,2003,35(3):365-394.

        [13] Rezaie A,Fahrenholtz W G,Hilmas G E.Effect of hot pressing time and temperature on the microstructure and mechanical properties of ZrB2-SiC[J].Journal of ,ate-rials Science,2007,42(8):2735-2744.

        [14] Pienti L,Sciti D,Silvestroni L,et al.Effect of milling on the mechanical properties of chopped SiC Fi-ber-reinforced ZrB2[J].Materials,2013,6(5):1980-1993.

        [15] Silvestroni L,Sciti D.Microstructure evolution upon annealing of a ZrB2-SiC composite containing lanthana and magnesia[J].Journal of the European Ceramic Society,2013,33(2):403-412.

        [16] 王江洋,錢振東,王亞奇.細(xì)觀尺度下大孔隙環(huán)氧瀝青混合料損傷演化分析[J].北京工業(yè)大學(xué)學(xué)報(bào),2013,39(8):1223-1229.

        [17] Zhi W,Zhanjun W,Guodong S.Effect of annealing treatment on mechanical properties of a ZrB2-SiC-graphite ceramic[J].Materials Science and Engineering:A,2011,528(6):2870-2874.

        [18] Rezaie A,Fahrenholtz W G,Hilmas G E.Effect of hot pressing time and temperature on the microstructure and mechanical properties of ZrB2-SiC[J].Journal of Materials Science,2007,42(8):2735-2744.

        [19] Akin I,Hotta M,Sahin F C,et al.Microstructure and densification of ZrB2-SiC composites prepared by spark plasma sintering[J].Journal of the European Ceramic Society,2009,29(11):2379-2385.

        (編輯:崔賢彬)

        Simulation of crack propagation in composite materials using cohesive zone model and image-based numerical model

        QIN Li-yuan,MENG Song-he,LI Jin-ping,JIN Hua

        (Center for Composite Materials and Structure,Harbin Institute of Technology,Harbin 150008,China)

        The cohesive zone model is used to simulate the intergranular fracture propagation to get the effects of the microstructure and the interfacial strength on the mechanical properties of the composites.The microstructure morphology of ZrB2-SiC ceramics illustrated by SEM was imported into ABAQUS finite element software,with the cohesive elements embedded along grain boundaries,to predict the intergranular fracture propagation.The effects of grain boundary and material inhomogeneity on the stress distribution and crack propagation were investigated by changing the ZrB2-SiC interfacial strength.The result shows an uneven stress distribution and stress concentration are observed due to the presence of grain boundaries.The locations of the crack initiation and propagation direction vary with the interfacial strengths and the crack propagates along the grain boundaries with low strength.The strength of polycrystals increases as ZrB2-SiC interfacial strength increases, and the elastic modulus remains substantially unchanged.

        microstructure;cohesive zone model;interfacial strength;ZrB2-SiC ceramics

        2016-06-01;

        2016-07-06。

        國(guó)家自然科學(xué)基金青年科學(xué)基金(11502058)。

        秦麗媛(1986—),女,博士,研究方向?yàn)槎嘞嗵沾刹牧隙喑叨确治觥-mail:qinliyuan@hit.edu.cn

        TB332

        A

        1006-2793(2017)04-0501-05

        10.7673/j.issn.1006-2793.2017.04.018

        猜你喜歡
        內(nèi)聚力晶界微觀
        CRTS Ⅱ型軌道板/CA 砂漿界面內(nèi)聚力模型研究
        晶界工程對(duì)316L不銹鋼晶界形貌影響的三維研究
        上海金屬(2022年4期)2022-08-03 09:52:00
        基于截?cái)嗲驙钅P偷腇e扭轉(zhuǎn)晶界的能量計(jì)算
        鐵/鎳基奧氏體多晶合金晶界彎曲研究進(jìn)展
        基于內(nèi)聚力模型的輪盤破裂轉(zhuǎn)速預(yù)測(cè)方法研究
        大學(xué)英語教學(xué)中影響閱讀教學(xué)的因素淺析
        一種新的結(jié)合面微觀接觸模型
        Inconel 600 合金的晶界工程工藝及晶界處碳化物的析出形貌
        上海金屬(2015年6期)2015-11-29 01:09:02
        微觀的山水
        詩(shī)選刊(2015年6期)2015-10-26 09:47:10
        微觀中國(guó)
        浙江人大(2014年8期)2014-03-20 16:21:15
        国产精品成人无码a 无码| 亚洲国产成人av二区| 精品国产成人av久久| 国产av无码专区亚洲版综合| 中文字幕久久精品一二三区 | 国产剧情一区二区三区在线| 国产精品久久久久aaaa| 色狠狠av老熟女| 国产亚洲欧美精品一区| 邻居少妇张开腿让我爽视频| 人妻精品视频一区二区三区| 国产精品亚洲欧美大片在线看| 最近中文字幕视频高清| 青春草在线视频精品| 国内精品极品久久免费看| 福利视频一区二区三区| 蜜桃久久精品成人无码av| 久久无码高潮喷水| 欧美破处在线观看| 日韩av在线亚洲女同| 国产在线无码一区二区三区视频 | 成人欧美日韩一区二区三区| а中文在线天堂| 91久久精品人妻一区二区| 色婷婷久久精品一区二区| 亚洲av无码一区二区三区不卡| 日本午夜免费福利视频| av天堂线上| 精品人妻69一区二区三区蜜桃| 国产精品久久久久aaaa| 国产美女在线精品亚洲二区| 蜜臀久久久精品国产亚洲av| 亚洲97成人在线视频| 人妻夜夜爽天天爽| 久久精品岛国av一区二区无码| 午夜天堂精品一区二区| 精品极品视频在线观看| 人妻无码一区二区三区免费| 日韩精品无码区免费专区| 午夜香蕉av一区二区三区| 亚洲国产精品一区二区成人av|