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

高Nb-TiAl合金高溫力學行為本構建模仿真技術

2018-02-05 08:09:28董成利于慧臣焦澤輝孔凡濤陳玉勇
航空材料學報 2018年1期
關鍵詞:有限元模型

董成利, 于慧臣, 焦澤輝, 孔凡濤, 陳玉勇

(1.中國航發北京航空材料研究院 先進高溫結構材料重點實驗室,北京 100095; 2.航空材料檢測與評價北京市重點實驗室,北京 100095; 3.材料檢測與評價航空科技重點實驗室,北京 100095; 4.哈爾濱工業大學 材料科學與工程學院,哈爾濱 150001)

TiAl合金具有低密度、高比強度、優異的高溫性能等綜合優點,因此廣泛應用于現代先進航空發動機渦輪葉片的制造,成為代替傳統鎳基高溫合金制造渦輪葉片的唯一候選材料,應用前景廣闊[1-4]。據權威統計,超過40000件TiAl合金用于制造波音787飛機的GEnx 1B發動機低壓渦輪葉片和波音747-8飛機的GEnx 2B低壓渦輪葉片[5]。為了研發能夠承受更高溫度的TiAl合金材料,國內外相關科研人員在傳統TiAl合金材料成分的基礎上添加了一定量的Nb元素,從而研制了具有更高高溫強度、更好抗氧化性的高Nb-TiAl合金,高Nb-TiAl合金的成功研制吸引了國內外業內科研人員的廣泛、持續關注。針對不同成分和組織的高Nb-TiAl合金材料,國內外相關科研人員開展了不同試驗條件下的力學性能測試與表征工作,高Nb-TiAl合金材料的性能測試與表征工作主要體現在室溫和高溫條件下的成分-組織-性能方面,力學性能主要涉及高Nb-TiAl合金材料的拉伸[6-8]、低周疲勞[9-11]、蠕變/持久[12-13]、裂紋擴展[14]和斷裂特性[15-16]等。

盡管目前國內外相關科研人員針對高Nb-TiAl合金材料開展了許多的力學性能測試,但缺乏對高Nb-TiAl合金材料疲勞-蠕變交互力學性能的認知和了解。一方面,渦輪葉片在航空發動機服役工況下不可避免地要承受疲勞-蠕變交互載荷,高Nb-TiAl合金要取代鎳基高溫合金,就必須研究其疲勞-蠕變交互力學行為。另一方面,從渦輪葉片結構強度設計層面上講,若要實現高Nb-TiAl合金的工程化應用,必須建立表征其綜合力學行為的本構模型,從而實現渦輪葉片在服役工況下力學行為的數值仿真模擬和分析。宏觀唯象Chaboche黏塑性統一本構模型常被用于鎳基高溫合金力學行為的有限元數值計算和分析,眾多計算案例表明,宏觀唯象Chaboche黏塑性統一本構模型能夠準確表征和預測鎳基高溫合金的力學行為,計算精度較高。實際上,Chaboche黏塑性統一本構模型本質上是層次結構,可以進行合理的裁剪,即可以任意減少或增加相關內變量及其演化方程來實現對材料或結構從簡單到復雜力學行為的數值模擬。

本研究針對高Nb-TiAl合金開展750 ℃條件下不同應變速率的單軸拉伸、低周疲勞、疲勞-蠕變交互和蠕變試驗,建立考慮Ohno-Wang修正并耦合Kachanov損傷的Chaboche黏塑性統一本構模型,采用自適應的顯式Euler法將微分形式的本構模型離散為差分方程組,并植入到有限元軟件ABAQUS的用戶材料子程序UMAT 中,對高Nb-TiAl合金的不同應變速率的單軸拉伸、低周疲勞、疲勞-蠕變交互和蠕變行為進行了有限元數值仿真模擬。

1 試驗材料與方法

試驗驗材料為高Nb-TiAl合金,合金名義成分為Ti-45Al-8Nb-0.2W-0.2B-0.02Y (原子分數/%),該合金具有均勻的近片層組織,且包含較細小的片層團(平均尺寸為50 μm左右)、分布于片層團周圍的γ相、硼化物以及YAl2等。該高Nb-TiAl合金較之傳統的近片層高Nb-TiAl合金擁有更多分布于片層團周圍的γ相。另外,鑄造過程造成了合金片層團的扭曲變形,具體原始顯微組織見圖1。

圖1 高Nb-TiAl合金原始顯微組織Fig.1 Microstructure of high-Nb TiAl alloy

為了獲得高Nb-TiAl合金的高溫綜合力學性能,并校驗后續本構模型對其綜合力學行為的表征預測能力,針對高Nb-TiAl合金開展了750 ℃條件下的不同應變速率的單軸拉伸、低周疲勞、疲勞-蠕變交互以及蠕變試驗。高Nb-TiAl合金的單軸拉伸、低周疲勞、蠕變-疲勞交互以及蠕變試樣實物如圖2所示。單軸拉伸試驗在Instron-5982試驗機上進行,試驗方法參照GB/T 4338—2006 《金屬材料 高溫拉伸試驗方法》執行,應變速率選擇1×10-3s-1和5×10-3s-1兩種形式;低周疲勞及疲勞-蠕變交互試驗在MTS-810試驗機上進行,試驗方法參照ASTM E606/E606M—2012《應變控制的疲勞試驗方法》執行,低周疲勞采用應變控制,加載速率為1×10-3s-1,疲勞-蠕變交互采用應變控制,加載速率為1×10-3s-1,保載形式分別為上保載60 s、上下各保載30 s和下保載60s三種類型,載荷波形見圖3,其中tT表示上保載時間,tC表示下保載時間;蠕變試驗在機械杠桿式高溫蠕變持久試驗機RJ—30試驗機上進行,試驗方法參照HB 5151—1996《金屬高溫拉伸蠕變試驗方法》執行。

2 Chaboche黏塑性統一本構理論

2.1 彈性本構關系

一般地,金屬材料的彈性變形部分服從廣義虎克定律(Generalized Hook’s Law)。令{ε}=[ε11,ε22,ε33,γ12,γ23,γ31]和{σ}=[σ11,σ22,σ33,τ12,τ23,τ31]分別為應變向量和應力向量,則金屬材料的彈性本構關系可以表示為:

{ε}=[C]{σ}

(1)

式中:[C]為柔度矩陣。

2.2 屈服函數

Chaboche黏塑性本構理論是與屈服準則相關聯的。它總是假設存在某種理想的屈服狀態,屈服函數F作為一個標量出現在本構方程中。當應力水平低于屈服應力時,非彈性變形并不出現,即F≤0,此時應力應變關系服從廣義虎克定量。只有當應力狀態達到足以引起屈服時才會有非彈性變形產生,此時F>0。另外,屈服函數F是與應力、溫度和內變量相關的函數,其一般表達式如下[17]:

圖2 高Nb-TiAl合金試樣實物圖 (a)單軸拉伸試樣;(b)蠕變試樣;(c)低周疲勞及疲勞-蠕變交互試樣Fig.2 Solid specimens of high-Nb TiAl alloy (a)uniaxial tensile specimen; (b)creep specimen;(c)low cycle fatigue and creep-fatigue specimens

圖3 高Nb-TiAl合金疲勞及疲勞-蠕變加載波形示意圖 (a)低周疲勞;(b)上保載60 s; (c)下保載60 s;(d)上下各保載30 sFig.3 Loading waveforms of LCF and creep-fatigue interaction (a)LCF;(b) tensile dwell time 60 s; (c)compressive dwell time 60 s;(d)balanced dwell time 30 s

(2)

式中:T是溫度;Vk(k=1,2,…,N)是N個內變量;σy是與加載速率相關的初始屈服應力;R是一個表示各向同性硬化的標量,表示由于各向同性硬化引起的應力(drag stress);Xij稱為背應力(back stress),表示在黏塑性流動過程中,應力空間內屈服面中心(或等勢面)的移動,描述了與方向有關的硬化。

有效應力偏量的第二不變量為

(3)

2.3 流動法則

在Chaboche黏塑性本構理論中,流動法則是與Von Mises屈服函數相關聯的,即:

(4)

式中:Λ=Λ(F)是一個非負的標量函數。

將屈服函數表達式代入式(4),經整理可得:

(5)

流動法則給出了非彈性應變率和應力、內變量之間所具有的與溫度相關的函數關系。在Chaboche黏塑性本構理論中,將式(4)中的Λ表達為屈服函數F的一個冪函數:

(6)

式中:K和n都是與溫度有關的材料參數;=uH(u);H(u)是Heaviside函數(u≤0,H(u)=0;u>0,H(u)=1)。關于Φ(F)的其他具體表達形式在Perzyna的文章[18-19]中也有所論述,這里不再贅述。

從式(4)~(6)可以看出,在Chaboche黏塑性本構理論中,黏塑性勢函數可以定義為以下形式:

(7)

黏塑性(非彈性)應變率是黏塑性勢關于應力的偏導數:

(8)

2.4 各向同性硬化/軟化模型

材料的各向同性硬化/軟化需要硬化變量R來描述,即通過改變當前彈性域尺寸大小來實現。各向同性硬化/軟化對應于塑性流動過程中應力的緩慢變化,表示屈服面在各個方向以相同的量在擴大或縮小。即通常發生在以下兩種情況中:一是在單調塑性流動過程中,屈服面以較低的速率在各個方向以相同的量在擴大或縮??;二是材料經過大量循環而接近達到循環硬化/軟化的穩定狀態[17]。

一般來說,各向同性變量服從下面的演化方程:

(9)

式中:Q,b,γ和m都是材料參數;Qr是各向同性變量R的漸近值;Q+k是屈服面尺寸的漸近值。注意,Q具有應力的量綱,其他幾個是無量綱的常數。

2.5 隨動硬化模型

材料在循環變形過程中體現出的平均應力松弛和棘輪行為主要是一種各向異性硬化行為的體現,為了更好地描述平均應力松弛和棘輪行為,已有研究主要集中在隨動硬化模型的發展方面。

Chaboche等采用非彈性應變率及其累積作為度量來描述非線性的隨動硬化,具體形式如下:

(10)

式(10)右邊第一項表達了Prager線性運動硬化準則,但是僅限于單調加載;第二項為動態恢復,它使得硬化可以考慮到載荷方向的改變,如Bauschinger效應,該項的引入改善了硬化法則描述遲滯環的能力。第三項是后來單獨加入的,描述了時間硬化恢復效應。

為了擴大模型的適用范圍,可以將隨動硬化變量分為幾個分量,每個分量都遵從相同的非線性演化規律,具體形式如下:

(11)

這樣,總的隨動硬化變量為:

(12)

2.6 蠕變損傷演化方程

蠕變損傷演化律采用Kachanov損傷方程[20-21]表達如下:

(13)

其中,pa,pr和pk是與溫度有關的材料蠕變損傷參數。

對于各向同性材料,其損傷有效應力張量為:

(14)

為了考慮損傷的作用,須在以上所述的本構方程組采用損傷有效張量形式。Chaboche黏塑性統一本構模型的應用價值不僅體現在能準確模擬材料的復雜變形,而且也體現在可很好地與現有大型通用有限元程序相結合,為實際復雜結構和部件的強度分析提供良好的工具。由于ABAQUS有限元軟件在非線性分析領域所具有的獨特優勢,特別是具有很強的數值穩定性和環境友好界面,因此本工作將Chaboche黏塑性統一模型通過ABAQUS提供的用戶材料子程序UMAT,實現Chaboche黏塑性統一本構方程與ABAQUS程序的有機結合,為利用Chaboche黏塑性統一本構模型完成航空發動機結構強度評價和分析奠定基礎。

3 結果與分析

圖4給出了750 ℃不同應變速率條件下(1×10-3s-1和5×10-3s-1)的高Nb-TiAl合金單軸拉伸及有限元模擬曲線。從圖4可以看出,施加的不同應變速率對材料的彈性階段幾乎沒有影響,兩組不同應變速率條件下的彈性拉伸數據幾乎重合,而施加的不同應變速率對材料的屈服階段影響明顯,即隨著施加應變速率的增加,材料發生屈服后的應力增加。以上試驗表明,該高Nb-TiAl合金在750 ℃下具有明顯的拉伸率相關性。另外,采用Chaboche黏塑性統一本構模型對高Nb-TiAl合金在750 ℃下的拉伸行為及率相關性進行了數值模擬計算,可以看出,Chaboche黏塑性統一本構模型能夠準確模擬該材料的拉伸行為,同時也能對材料的拉伸率相關性進行準確表征。

圖4 不同應變速率的單軸拉伸及有限元模擬曲線Fig.4 Experimental and simulated tensile data at different strain rates

圖5給出了750 ℃不同應力水平條件下(360 MPa,380 MPa和420 MPa)的高Nb-TiAl合金蠕變及有限元模擬曲線。從圖5可以看出,隨著施加應力水平的增加,材料蠕變斷裂壽命逐漸減小,材料發生蠕變斷裂時的蠕變應變也逐漸減小,而穩態蠕變速率逐漸增加。對于應力水平為380 MPa和420 MPa的載荷狀況,高Nb-TiAl合金的蠕變曲線僅包含初始蠕變和穩態蠕變兩個典型階段,而對于應力水平為360 MPa的載荷狀況,高Nb-TiAl合金的蠕變曲線包含初始蠕變、穩態蠕變和加速蠕變三個完整階段。采用Chaboche黏塑性統一本構模型對高Nb-TiAl合金在750 ℃下的蠕變行為進行了數值模擬計算,可以看出,Chaboche黏塑性統一本構模型能夠準確模擬該材料的蠕變行為,能夠對材料的初始蠕變、穩態蠕變和加速蠕變三個完整階段進行準確表征。

圖5 不同應力下的蠕變及有限元模擬曲線Fig.5 Experimental and simulated creep data at different stresses

圖6給出了750 ℃不同應變水平條件下(0.30%和0.35%)的高Nb-TiAl合金循環應力響應規律,主要考察了保載條件對材料循環應力響應的影響。從圖中可以看出,不論是低周疲勞還是疲勞-蠕變交互試驗條件,隨著施加應變水平的增加,材料的循環應力增加,而斷裂壽命減小。高Nb-TiAl合金循環應力響應規律因試驗條件的不同而有所不同,對于高Nb-TiAl合金的低周疲勞試驗而言,不同應變水平下的循環應力隨著循環數的增加而增加,表現出明顯的循環硬化現象,且硬化現象隨著施加應變水平的增加而愈發明顯。對于高Nb-TiAl合金的疲勞-蠕變交互試驗而言,不論是上保載、下保載和是上下對稱保載狀態,不同應變水平下的循環應力隨著循環數的增加而減小,表現出明顯的循環軟化現象。這是因為在疲勞試驗應變峰谷處保持一定時間的應變恒定,可產生明顯的應力松弛效應,從而產生循環軟化現象。采用Chaboche黏塑性統一本構模型對高Nb-TiAl合金在750 ℃下的循環應力變化行為進行了數值模擬計算,可以看出,Chaboche黏塑性統一本構模型能夠準確模擬該材料的循環應力響應,能夠對材料的循環硬化和/或軟化行為進行合理表征,雖然擬合曲線與試驗數據存有一定的微小差距,但已完全滿足材料工程化應用的需求。

圖6 不同應變條件下的應力范圍隨循環數變化的試驗及有限元模擬曲線 (a)上保載60 s;(b)上下各保載30 s;(c)下保載60 s;(d)不保載Fig.6 Experimental and simulated stress range-cyclic life data at different strains (a)tensile dwell time 60 s; (b)balanced dwell time 30 s;(c)compressive dwell time 60 s;(d)no dwell time

圖7給出了750 ℃不同應變水平條件下(即0.30%和0.35%)的高Nb-TiAl合金穩態遲滯回線響應規律,主要考察了保載條件對材料穩態遲滯回線響應的影響。從圖7可以看出,不論是低周疲勞還是疲勞-蠕變交互試驗條件,隨著施加應變水平的增加,穩態遲滯回線的面積(即循環應變遲滯能)增加,且穩態遲滯回線的峰谷應力值也增加。對于無保載的低周疲勞和上下各保載30 s的疲勞-蠕變交互試驗,穩態遲滯回線關于應力軸或應變軸基本上呈反對稱結構,而對于上保載60 s和下保載60 s疲勞-蠕變交互試驗,穩態遲滯回線分別沿著應力軸下移和上移,從而產生非零的平均應力。從疲勞-蠕變交互試驗曲線可以看出,應力松弛發生在應變峰谷處的保持階段。采用Chaboche黏塑性統一本構模型對高Nb-TiAl合金在750 ℃下的穩態遲滯回線進行了數值模擬計算,可以看出,Chaboche黏塑性統一本構模型能夠準確模擬該材料的穩態遲滯回線,特別是能夠準確模擬疲勞-蠕變交互試驗下的穩態遲滯回線所反映出來的應力松弛現象。

圖7 不同應變條件下的穩態遲滯回線的試驗及有限元模擬曲線 (a)上保載60 s; (b)上下各保載30 s;(c)下保載60 s;(d)不保載Fig.7 Experimental and simulated stable hysteresis loops at different strains (a)tensile dwell time 60 s; (b)balanced dwell time 30 s;(c)compressive dwell time 60 s;(d)no dwell time

綜上所述,本工作采用的Chaboche黏塑性統一本構模型能夠準確模擬高Nb-TiAl合金的拉伸行為(率相關性)、蠕變行為(完整蠕變曲線三個階段)、低周疲勞行為(循環硬化/軟化現象)和疲勞-蠕變行為(應力松弛現象),且模擬精度較高,為Chaboche黏塑性統一本構模型成功應用于航空發動機TiAl合金葉片的結構優化設計和強度校核提供了技術支撐。

4 結論

(1)通過對高Nb-TiAl合金750 ℃條件下的綜合力學性能試驗,獲得了材料不同應變速率的單軸拉伸、低周疲勞、疲勞-蠕變交互和蠕變試驗數據和曲線。

(2)采用自適應的顯式Euler法將微分形式的Chaboche黏塑性統一本構模型離散為差分方程組,并植入到有限元軟件ABAQUS的用戶材料子程序UMAT,成功實現對高Nb-TiAl合金高溫力學行為的有限元數值模擬仿真。

(3)宏觀唯象Chaboche黏塑性統一本構模型耦合了蠕變損傷演化率并考慮了Ohno-Wang修正,使得該模型能夠更加準確地表征和預測高Nb-TiAl合金的拉伸行為(率相關性)、蠕變行為(完整蠕變曲線三個階段)、低周疲勞行為(循環硬化/軟化現象)和疲勞-蠕變行為(應力松弛現象),且模擬精度較高。

[1] NIEH T G, HSIUNG L M, WADSWORTH J. Superplastic behavior of a powder metallurgy TiAl alloy with a metastable microstructure [J]. Intermetallics, 1999, 7: 163-170.

[2] SONG L, LIN J P, LI J S. Effects of trace alloying elements on the phase transformation behaviors of ordered ω phases in high Nb-TiAl alloys [J]. Materials and Design, 2017, 113: 47-53.

[3] FANG L, LIN J P, LIANG Y F,etal. Microstructural instability in surface layer of a high Nb-TiAl alloy processed by shot peening following high temperature exposure [J]. Intermetallics, 2016, 78: 8-16.

[4] PATRIARCA L, FILIPPINI M, BERETTA S. Digital image correlation-based analysis of strain accumulation on a duplex γ-TiAl[J]. Intermetallics, 2016, 75: 42-50.

[5] CHEN G, PENG Y B, ZHENG G,etal. Polysynthetic twinned TiAl single crystals for high-temperature applications [J]. Nature Materials, 2016, 15: 876-881.

[6] ZHOU H T, KONG F T, WANG X P,etal. High strength in high Nb containing TiAl alloy sheet with fine duplex microstructure produced by hot pack rolling [J]. Journal of Alloys and Compounds, 2017, 695: 3495-3502.

[7] KARTAVYKH A V, ASNIS E A, PISKUN N V,etal. Room-temperature tensile properties of float-zone processed β-stabilized γ-TiAl(Nb,Cr,Zr) intermetallic [J]. Materials Letters, 2017,188: 88-91.

[8] ZHANG S Z, ZHANG C J, DU Z X,etal. Microstructure and tensile properties of hot fogred high Nb containing TiAl based alloy with initial near lamellar microstructure [J]. Materials Science and Engineering: A, 2015, 642:16-21.

[9] DING J, LIANG Y F, XU X J,etal. Cyclic deformation and microstructure evolution of high Nb containing TiAl alloy during high temperature low cycle fatigue [J]. International Journal of Fatigue, 2017, 99: 68-77.

[10] YU L, SONG X P, ZHANG M,etal. Effect of stress ratio on fatigue lifetime of high Nb containing TiAl alloy at elevated temperature [J]. Materials & Design, 2015, 84: 378-384.

[11] KRUML T, OBRTLK K. Microstructure degradation in high temperature fatigue of TiAl alloy [J]. International Journal of Fatigue, 2014, 65: 28-32.

[12] TIAN S G, LV X X, YU H C,etal. Creep behavior and deformation feature of TiAl-Nb alloy with various states at high temperature [J]. Materials Science and Engineering: A, 2016, 651: 490-498.

[13] TIAN S G, WANG Q, YU H C,etal. Microstructure and creep behaviors of a high Nb-TiAl intermetallic compound based alloy [J]. Materials Science and Engineering: A, 2014, 614: 338-346.

[14] ZHANG M, SONG X P, YU L,etal. In situ observation of fatigue crack initiation and propagation behavior of a high-Nb TiAl alloy at 750°C [J]. Materials Science and Engineering: A, 2015, 622: 30-36.

[15] CUI W F, LIU C M. Fracture characteristics of γ-TiAl alloy with high Nb content under cyclic loading [J]. Journal of Alloys and Compounds, 2009, 477: 596-601.

[16] WU Z, HU R, ZHANG T B,etal. Microstructure determined fracture behavior of a high Nb containing TiAl alloy [J] . Materials Science and Engineering: A, 2016, 666: 297-304.

[17] CHABOCHE J L. Development of continuum damage mechanics for elastic solids sustaining anisotropic and unilateral damage [J]. International Journal of Damage Mechanics, 1993, 2: 311-329.

[18] PERZYNA P. The constitutive equations for rate sensitive plastic materials [J]. Quarterly Journal of Applied Mathematics, 1963, 20: 321-332.

[19] PERZYNA P. Fundamental problems in viscoplasticity[J]. Advances in Applied Mechanics, 1966, 9: 243-377.

[20] LEMAITRE J, CHABOCHE J L. Mechanics of solid materials [M]. New York: Cambridge University Press, 1990.

[21] RABOTHOV Y M. Creep problems in structural members [M]. Amsterdam: North-Holland, 1969.

猜你喜歡
有限元模型
一半模型
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
新型有機玻璃在站臺門的應用及有限元分析
上海節能(2020年3期)2020-04-13 13:16:16
基于有限元的深孔鏜削仿真及分析
基于有限元模型對踝模擬扭傷機制的探討
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
磨削淬硬殘余應力的有限元分析
基于SolidWorks的吸嘴支撐臂有限元分析
主站蜘蛛池模板: 91精品国产综合久久香蕉922| 毛片免费在线视频| 亚洲精品成人片在线观看| 日韩中文无码av超清| 毛片久久久| 亚洲视频在线观看免费视频| 精品一区二区三区波多野结衣| 久久亚洲国产最新网站| 亚洲一级毛片在线观播放| 国产精品自在线天天看片| 亚洲一区二区在线无码| 国产乱子伦手机在线| 青青青国产精品国产精品美女| 亚洲欧洲日韩久久狠狠爱| 欧美日韩动态图| 在线无码九区| v天堂中文在线| 国产精品成人一区二区| 成人在线欧美| 在线观看91香蕉国产免费| 国产精品夜夜嗨视频免费视频| 在线观看91香蕉国产免费| 久久免费观看视频| 香蕉网久久| 岛国精品一区免费视频在线观看 | 最新国产午夜精品视频成人| 国产在线观看一区二区三区| 人妻精品久久无码区| 黄色在线不卡| 日韩无码黄色| 亚洲高清免费在线观看| 亚洲欧美人成人让影院| 国产精品自在拍首页视频8| 成人无码一区二区三区视频在线观看 | AⅤ色综合久久天堂AV色综合| 日韩精品成人网页视频在线 | 日韩欧美国产成人| 综合色在线| 亚洲第一成年人网站| 日韩A∨精品日韩精品无码| 免费人成视网站在线不卡| 99久视频| 欧美精品亚洲精品日韩专区| 亚洲精品第一页不卡| 久久精品人人做人人爽电影蜜月 | 无码粉嫩虎白一线天在线观看| 日韩成人免费网站| 国产福利小视频在线播放观看| 88av在线| 国产超碰一区二区三区| 婷婷色一二三区波多野衣| 亚洲日韩国产精品无码专区| 亚洲美女久久| 免费啪啪网址| 亚洲中文字幕在线一区播放| 天堂va亚洲va欧美va国产| 波多野结衣无码AV在线| 亚洲第一色网站| 日韩欧美色综合| 亚洲乱码精品久久久久..| 四虎国产精品永久在线网址| 人妻丝袜无码视频| 污污网站在线观看| 全免费a级毛片免费看不卡| 国产资源站| 午夜a视频| 成人免费一区二区三区| 亚洲一本大道在线| AV不卡在线永久免费观看| 一本大道香蕉高清久久| 色综合热无码热国产| 中文字幕 91| 亚洲日韩国产精品综合在线观看| 99热在线只有精品| 老司机精品99在线播放| 日韩无码视频播放| 国产欧美日韩视频怡春院| 午夜日b视频| 99久视频| 激情综合激情| 婷婷六月色| 无码专区第一页|