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

感應(yīng)加熱工藝參數(shù)對Q345鋼彎曲角度和曲率半徑的影響

2017-11-14 05:55:33張繼祥
關(guān)鍵詞:有限元

張繼祥,周 偉,徐 昱,鐘 厲,2

(1.重慶交通大學(xué) 機(jī)電與汽車工程學(xué)院,重慶 400074; 2.重慶交通大學(xué) 重慶市特種船舶數(shù)字化設(shè)計與制造工程技術(shù)研究中心,重慶 400074)

感應(yīng)加熱工藝參數(shù)對Q345鋼彎曲角度和曲率半徑的影響

張繼祥1,2,周 偉1,徐 昱1,鐘 厲1,2

(1.重慶交通大學(xué) 機(jī)電與汽車工程學(xué)院,重慶 400074; 2.重慶交通大學(xué) 重慶市特種船舶數(shù)字化設(shè)計與制造工程技術(shù)研究中心,重慶 400074)

基于ANSYS有限元分析軟件,采用APDL語言實現(xiàn)高頻感應(yīng)加熱成形中電-磁-熱-力等多物理場耦合,建立了高頻感應(yīng)加熱彎板的三維移動式有限元模型,并系統(tǒng)研究了高頻感應(yīng)加熱彎板成形目標(biāo)與工藝參數(shù)的關(guān)系。結(jié)果表明:加熱功率較大時,鋼板彎曲角度隨著熱源移動速度增加先增大再減小,而加熱功率較小時鋼板彎曲角度隨著熱源移動速度增加直線下降;鋼板彎曲曲率半徑隨熱源移動速度增加而增大,隨加熱功率增大而減小。最后分別建立了彎曲角度、曲率半徑與加熱功率和熱源移動速度的函數(shù)關(guān)系式。

船舶工程;有限元模擬;高頻感應(yīng);彎板成形;函數(shù)關(guān)系

水火彎板成形在船板制造領(lǐng)域效率低,形狀達(dá)不到要求,而高頻感應(yīng)加熱彎板成形具有溫度控制精準(zhǔn)、氧化少、污染少、加熱效率極高、工作環(huán)境清潔安全、有利于實現(xiàn)自動化加工等優(yōu)勢[1-3],已成為一種廣泛應(yīng)用的板材熱應(yīng)力成形方法,國內(nèi)外相關(guān)學(xué)者均在開展這方面的研究。

在國外,K.Y.BAE等[4]運用統(tǒng)計方法對感應(yīng)加熱變形與加熱參數(shù)之間的關(guān)系進(jìn)行了研究,并證實了其方程的合理性;J.H.LEE等[5]建立有限元模型,分析了鋼板高頻感應(yīng)加熱和可燃?xì)怏w加熱下鋼板的殘余應(yīng)力和變形;C.D.JANG等[6]則提出了鋼板變形預(yù)測可通過仿真模擬技術(shù)來實現(xiàn)的思想。

在國內(nèi),相關(guān)學(xué)者也開展了相關(guān)的研究,并取得了許多成果。殷筱依等[7]、劉芳平等[8]研究了薄壁件的翹曲和屈曲;范平等[9]、張繼祥等[10]和安國銀等[11]分別建立了靜止的電、磁、熱耦合熱源模型,但所建立的感應(yīng)熱源是通過溫度場與磁場間接耦合得到,該模型只能實現(xiàn)二維靜止加熱模擬,無法移動;張雪彪等[12-13]建立了三維靜止式磁熱耦合熱源模型,實現(xiàn)了三維靜態(tài)的耦合模擬,但沒有進(jìn)行移動加熱模擬;周宏等[14]、張繼祥等[15-16]和徐昱[17]分別進(jìn)行了三維移動式有限元模型的研究。

筆者在相關(guān)學(xué)者提出的電-磁-熱-力等多物理場耦合實現(xiàn)方法的基礎(chǔ)上[15-16],基于ANSYS二次開發(fā)APDL語言,建立了三維移動式的高頻感應(yīng)加熱彎板有限元模型;系統(tǒng)研究了成形工藝參數(shù)對鋼板彎曲角度、曲率半徑等的影響;并建立成型參數(shù)與成形目標(biāo)之間的關(guān)系模型。

1 有限元建模

1.1 限元模型

模擬實驗采用400 mm×200 mm×10 mm船用鋼板,材質(zhì)為Q345。實際圓形感應(yīng)線圈是使用直徑6 mm、壁厚1 mm的銅管纏繞而成。為簡化有限元建模模型,將圓形線圈改為34 mm×34 mm×6 mm面積相等的矩形線圈,線圈與鋼板間距2 mm,且截面電流加載方向和圓形線圈的方向一致。圖1為高頻感應(yīng)成形實的有限元模型。

圖1 鋼板高頻感應(yīng)加熱成形模型Fig.1 High frequency induction heating formation model of steel plate

1.2 網(wǎng)格劃分

在計算諧波磁場時,采用SOLID97磁矢量單元,溫度場分析時采用8節(jié)點SOLID70熱實體單元,在進(jìn)行結(jié)構(gòu)應(yīng)力場計算時采用SOLID185單元等效替換SOLID70單元。劃分后的模型如圖2。

圖2 有限元網(wǎng)格劃分Fig.2 Finite element meshing

1.3 模擬計算流程

筆者基于ANSYS軟件,耦合溫度場和磁場熱源,使用APDL語言循環(huán)語句來完成對熱源的移動,并對移動式三維感應(yīng)加熱電磁-熱-應(yīng)力多物理場進(jìn)行模擬,耦合模擬計算流程如圖3。

圖3 模型編程算法流程Fig.3 Algorithm flow chart of model programming

1.4 邊界及約束條件(圖4)

如圖4,在鋼板左右兩側(cè)面和底部兩端交線上的所有點施加UY=0 約束。為更加準(zhǔn)確的模擬出實際狀況下磁力線情況,對所有邊界上的點施加磁力約束。在Z軸方向上,為了模擬X-Y平面上的磁場,對Z=0和Z=M的平面上分別施加與其平行的AZ=0約束;為了模擬X平面上的磁場,在X=-L和X=L的平面上,分別施加與其平行的AZ=0約束。

散熱邊界條件。在鋼板上下、前后表面以及右側(cè)面分別施加換熱系數(shù)為15、10 W/(m2·℃)的散熱邊界條件[13]。

1.5 材料參數(shù)

Q345鋼熱物理參數(shù)中導(dǎo)熱系數(shù)K、屈服強度σs、切變模量G、彈性模量E、磁導(dǎo)率μ等參數(shù)[14]如表1,假定銅線圈和空氣的磁導(dǎo)率μ=1 H/m。

圖4 模型約束條件Fig.4 Model constraints

溫度/℃導(dǎo)熱系數(shù)K/[w·(k·m)-1]溫度/℃屈服強度σs/Pa切變模量G/Pa彈性模量E/Pa溫度/℃磁導(dǎo)率μ/(H·m-1)050.0003.45E81.37E92.06E11030012547.502003.00E81.32E91.92E11160.029025043.404001.50E81.22E71.72E11291.528237540.006008.00E71.04E31.50E11477.626150036.708004.00E701.22E11635.023562533.3010003.00E707.68E10698.020475030.00709.018487525.96720.3135100026.80742.017761.011000.01

1.6 研究方案

電源頻率為20 kHz,設(shè)定的功率分別為20、25、30、35 kW,采用10~30 mm/s 的移動速度對彎板工藝進(jìn)行了數(shù)值模擬。

2 數(shù)值模擬結(jié)果及分析

2.1 工藝參數(shù)對彎曲角度影響(圖5)

由圖5可看出:加熱功率較大時(25~35 kW),鋼板彎曲角度隨著熱源移動速度增加先增大再減小;加熱功率較小時(20 kW),鋼板彎曲角度隨著熱源移動速度增加直線下降,沒有出現(xiàn)峰值;不同功率下彎曲角度峰值出現(xiàn)位置不同,功率越大,峰值對應(yīng)的加熱速度越大。

熱源移動速度較小時,加熱功率對彎曲角度影響不大,但隨著熱源移動速度的增加,功率對彎曲角度的影響變大,功率越大,成形角度越大。熱源的移動速度較小時,鋼板的透熱性較好,厚度方向上溫度梯度變化小,熱應(yīng)力都較小,所以彎曲角度在各功率下變化不大;隨著熱源移動速度不斷增加。加熱功率較低時,鋼板透熱性變差,厚度方向上溫度梯度變化較大;熱應(yīng)力增加,鋼板的彎曲應(yīng)力也隨之增加,但鋼板底部的溫度較低,強度較大,底部材料不能夠發(fā)生塑形拉伸。所以在低功率加熱時,彎曲角度隨著熱源移動速度的增加而減小。

隨著熱源移動速度的不斷增加,在加熱功率較高時,鋼板透熱性變好,厚度方向上溫度梯度較小,熱應(yīng)力增加,鋼板總的彎曲應(yīng)力也隨之增加,此時底部鋼板溫度較高,屈服應(yīng)力較低,發(fā)生了塑性伸長變形,所以彎曲角度不斷增加。由于線能量的原因,鋼板底部加熱溫度隨熱源移動速度的增加而不斷降低,在彎曲應(yīng)力與塑形變形抗力的綜合作用下,彎曲角度會達(dá)到一個峰值,然后逐漸下降。

圖5 彎曲角度與熱源移動速度關(guān)系Fig.5 Relationship between bending angle and heat source moving speed

2.2 工藝參數(shù)對曲率半徑影響(圖6)

由圖6可看出:鋼板彎曲曲率半徑隨著熱源移動速度增加而不斷增大;加熱功率越大,上升幅度越小,即加熱功率越大曲率半徑越小。

在熱源的移動速度較小時,鋼板的透熱性較好,厚度方向上溫度梯度變化小,熱應(yīng)力都較小,所以彎曲角度在各功率下變化不大。各功率下鋼板的彎曲曲率半徑都隨熱源移動速度的增加而不斷增大,同一熱源移動速度下功率越大鋼板彎曲曲率半徑越小。這是因為鋼板的彎曲曲率半徑是由彎曲部分的弧線和彎曲角度共同決定的。在熱源移動速度相同,弧長一定時,功率越小使得鋼板彎曲角度越小,鋼板曲率半徑則越大;當(dāng)加熱功率一定時,鋼板彎曲角度隨熱源移動速度增大而越小,鋼板的曲率則半徑越大。

圖6 曲率半徑與熱源移動速度關(guān)系Fig.6 Relationship between curvature radius and heat source moving speed

2.3 彎曲角度與熱源移動速度、加熱功率關(guān)系

彎曲角度與熱源移動速度、加熱功率之間的函數(shù)關(guān)系如式(1):

θ=f(x,y)=f(v,P)

(1)

式中:θ為目標(biāo)值彎曲角度,(°);v為熱源移動速度,mm/s;P為加熱功率,kW。

筆者建立上述函數(shù)關(guān)系如式(1),并采用MATLAB曲面擬合功能進(jìn)行3次多項式函數(shù)擬合,擬合結(jié)果如圖7。

圖7 彎曲角度-加熱功率-熱源移動速度擬合曲面Fig.7 Bending angle-heating power-heat source moving speed fitting surface

圖7中,x軸表示熱源移動速度,mm/s;y軸表示加熱功率,kW;z軸表示鋼板冷卻后的彎曲角度,(°)。曲面頂點位置即最大彎曲角度以及對應(yīng)的源移動速度、加熱功率工藝參數(shù)值。

在熱源移動速度為v=22 mm/s、加熱功率為P=35 kW的熱加工參數(shù)下,鋼板得到最大彎曲角度。運用MATLAB軟件擬合結(jié)果得到的彎曲角度與熱源移動速度、加熱功率的函數(shù)關(guān)系如式(2):

f(v,P)=p00+p10v+p01P+p20v2+p11vP+p02P2+p30v3+p21v2P+p12vP2+p03P3

(2)

式中:p00=0.344 9,p10=-0.045 57,p-01=0.030 59,p20=0,p11=0.003 26,p02=-0.898 6×10-3,p30=9.306×10-5,p21=0,p12=0.139 7×10-3,p03=-3.378×10-5。

上述擬合的置信度為95%,方差為0.009 236。

2.4 曲率半徑與熱源移動速度、加熱功率關(guān)系

函數(shù)曲率半徑與熱源移動速度、加熱功率之間的關(guān)系函數(shù)如式(3):

R=g(x,y)=g(v,P)

(3)

式中:R為目標(biāo)值曲率半徑,m。

筆者建立上述函數(shù)關(guān)系如式(3),采用MATLAB曲面擬合功能進(jìn)行3次多項式函數(shù)擬合,擬合結(jié)果如圖8。

圖8 曲率半徑-加熱功率-熱源移動速度擬合曲面Fig.8 Curvature radius-heating power-heat source moving speed fitting surface

由圖8可知:x軸表示熱源移動速度,mm/s;y軸表示加熱功率,kW;z軸表示鋼板冷卻后的彎曲曲率半徑,m。

經(jīng)MATLAB擬合處理結(jié)果得到的曲率半徑與熱源移動速度、加熱功率的函數(shù)關(guān)系如式(4):

g(v,P)=p00+p10v+p01P+p20v2+p11vP+p02P2+p30v3+p21v2P+p12vP2+p03P3+p40v4+p31v3P+p22v2P2+p13vP3+p50v5+p41v4P+p32v3P2+p23v2P3

(4)

式中:p00=-146.4,p10=26.7,p01=11.88,p20=-1.193,p11=-2.129,p02=-0.319 1,p30=0.014 97,p21=0.092 06,p12=0.052 5,p03=0.003 385,p40=0.000 134 3,p31=-0.001 389,p22=-0.001 739,p13=-0.000 510 8,p50=-2.799×10-6,p41=-6.955×10-6,p32=1.17×10-6,p23=1.335×10-5。

上述擬合的置信度為95%,方差為4.085。

3 結(jié) 論

1)熱源移動速度較低時,加熱功率對彎曲角度影響較小。隨著熱源移動速度的增加,功率越大,成形角度越大。加熱功率較大時,鋼板彎曲角度隨著熱源移動速度增加先增大再減小,而加熱功率較小時鋼板彎曲角度隨著熱源移動速度增加直線下降,沒有出現(xiàn)峰值。

2)熱源移動速度較低時,加熱功率對彎曲曲率半徑影響不大;隨著熱源移動速度的增加,鋼板彎曲曲率半徑增大;加熱功率越小,曲率半徑越大。

[1] 周宏,羅宇,蔣志勇.船舶板材邊界條件對高頻感應(yīng)彎板成形的影響[J].船海工程,2009,38(5):19-23.

ZHOU Hong, LUO Yu, JIANG Zhiyong. Influence of boundary condition on high frequency inducting plate bending [J].Ship&OceanEngineering, 2009, 38(5): 19-23.

[2] 紀(jì)卓尚,劉玉君,孫振烈,等.船體曲面鋼板水火加工成型工藝的理論與應(yīng)用研究[J].中國造船,1998(增刊):118-124.

JI Zhuoshang, LIU Yujun, SUN Zhenlie, et al. Study on theory and application for hull steel plate bending by line heating[J].ShipbuildingofChina, 1998(Sup): 118-124.

[3] 張雪彪.船體曲面鋼板完全線加熱成形研究[D].大連:大連理工大學(xué),2006.

ZHANG Xuebiao.ResearchonShip-hullCurvedPlateFormingbyPureLineHeating[D].Dalian:Dalian University of Technology,2006.

[4] BAE K Y, YANG Y S, HYUN C M, et al. Derivation of simplified formulas to predict deformations of plate in steel forming process with induction heating [J].InternationalJournalofMachineToolsandManufacture, 2008, 48(15): 1646-1652.

[5] LEE J H, LEE K H, YUN J S. An electromagnetic and thermo-mechanical analysis of high frequency induction heating for steel plate bending[J].KeyEngineeringMaterials,2006,326/328:1283-1286.

[6] JANG C D, KIM H K, HA Y S. Prediction of plate bending by high-frequency induction heating [J].JournalofShipProduction, 2002, 18(4): 226-236.

[7] 張繼祥,殷筱依,劉鳳芝,等.纖維增強薄壁塑件翹曲變形分析[J].重慶交通大學(xué)學(xué)報(自然科學(xué)版),2016,35(4):173-178.

ZHANG Jixiang, YIN Xiaoyi, LIU Fengzhi et al. The analysis of the warpage of thin-wall plastic part for fiber reinforced plastic [J].JournalofChongqingJiaotongUniversity(NaturalScience), 2016, 35(4): 173- 178.

[8] 劉芳平,周建庭.變截面開口薄壁鋼箱穩(wěn)定性分析[J].重慶交通大學(xué)學(xué)報(自然科學(xué)版),2012,31(6):1105-1108.

LIU Fangping, ZHOU Jianting. Stability analysis of open thin-walled steel box with variable cross-section[J].JournalofChongqingJiaotongUniversity(NaturalScience), 2012, 31(6): 1105-1108.

[9] 范平,陳明和,周兆峰,等.板料無模成形技術(shù)—熱應(yīng)力成形[J].新技術(shù)新工藝,2007(11):62-64.

FAN Ping, CHEN Minghe, ZHOU Zhaofeng, et al. Dieless thermal stress forming of plate [J].NewTechnology&NewProcess, 2007(11): 62-64.

[10] 張繼祥,安國銀,李政君,等.不同厚度下船舶板材高頻感應(yīng)自由彎曲成形[J].重慶交通大學(xué)學(xué)報(自然科學(xué)版),2012,31(4):900-904.

ZHANG Jixiang, AN Guoyin, LI Zhengjun et al. Free bending of ship plate of different thickness on high-frequency induction[J].JournalofChongqingJiaotongUniversity(NaturalScience),2012,31(4): 900- 904.

[11] 安國銀,張繼祥,李政君,等.船板高頻感應(yīng)加熱關(guān)鍵點變形規(guī)律研究[J].船舶工程,2013,35(5):87-89.

AN Guoyin, ZHANG Jixiang, LI Zhengjun, et al. Study of deformation law of key points on high-frequency induction heating for ship plate [J].ShipEngineering, 2013, 35(5): 87-89.

[12] 張雪彪,楊玉龍,劉玉君,等.鋼板彎曲工藝中的高頻感應(yīng)加熱過程數(shù)值分析[J].中國造船,2011,52(2):108-116.

ZHANG Xuebiao, YANG Yulong, LIU Yujun, et al. Process applied to steel plate bending technology[J].ShipbuildingofChina,2011, 52(2): 108-116.

[13] 張雪彪,楊玉龍,劉玉君.鋼板高頻感應(yīng)加熱過程電磁-熱耦合場分析[J].大連理工大學(xué)學(xué)報,2012,52(5):676-682.

ZHANG Xuebiao, YANG Yulong, LIU Yujun. Analysis of electromagnetic - thermal coupling field for high frequency induction heating process of steel plate [J].JournalofDalianUniversityofTechnology, 2012, 52(5): 676-682.

[14] 周宏,蔣志勇.表面換熱系數(shù)對高頻感應(yīng)彎板成型的影響[J].材料熱處理學(xué)報,2013,34(增刊Ⅱ):262-266.

ZHOU Hong, JIANG Zhiyong. Effect of high frequency induction technology-based model surface heat exchange coefficient on bending plate forming of ship hull [J].TransactionsofMaterialsandHeatTreatment, 2013, 34(SupⅡ): 262-266.

[15] 張繼祥,徐昱,王智祥,等.鋼板高頻感應(yīng)移動加熱成形有限元建模及模擬研究[J].鍛壓技術(shù),2014,39(7):40-46.

ZHANG Jixiang, XU Yu, WANG Zhixiang, et al. Finite element modeling and simulation research of high frequency induction heating forming of steel plate based on moving heat source [J].Forging&StampingTechnology, 2014, 39(7): 40-46.

[16] 張繼祥,楊泮,劉鳳芝,等.跟蹤水冷下鋼板高頻感應(yīng)三維移動加熱成形有限元模型研究[J].鍛壓技術(shù),2015,40(7):113-119.

ZHANG Jixiang, YANG Pan, LIU Fengzhi, et al. Research on finite element model of 3D mobile tracking water-cooling steel plate with high frequency induction in the heating forming process with high frequency induction [J].Forging&StampingTechnology, 2015, 40(7): 113-119.

[17] 徐昱.三維移動式感應(yīng)加熱彎板成形有限元模擬研究[D].重慶:重慶交通大學(xué),2015.

XU Yu.StudyonFiniteElementSimulationofThree-DimensionMovingInductiveHeatingPlateBendingForming[D].Chongqing:Chongqing Jiaotong University, 2015.

EffectsofParametersofInductionHeatingProcessonBendingAngleandRadiusofCurvatureofQ345Steel

ZHANG Jixiang1, 2, ZHOU Wei1, XU Yu1, ZHONG Li1, 2

(1.School of Mechatronics & Vehicle Engineering, Chongqing Jiaotong University, Chongqing 400074, P. R. China; 2.Chongqing Engineering Research Center for Special Ship Digital Design and Manufacturing, Chongqing Jiaotong University, Chongqing 400074, P. R. China)

Based on the finite element analysis software ANSYS and its APDL language, an electric-magnetic-heat-stress physic field coupling in the HFV (high frequency induction) heating formation was realized and a finite element model of three-dimensional mobile high frequency induction heating bending plate was established. The relationship between HFV bending plate forming targets and their process parameters was systematically studied. The results show that the steel plate bending angle increases first and then decreases with the increase of the heat source moving speed when the heating power is large; while, the bending angle decreases sharply with the increase of heat source moving speed when the heating power is relatively small. The radius of curvature of steel plate increases with the increase of heat source moving speed, and decreases with the increase of heating power. Finally, the function formula of the bending angle with the heating power and the heat source moving speed, and that of the curvature radius with the heating power and the heat source moving speed were established respectively.

ship engineering; FEM simulation; high frequency induction; bending plate formation; function formula

10.3969/j.issn.1674-0696.2017.11.20

2016-04-13;

2016-06-29

重慶市基礎(chǔ)前沿研究計劃項目(cstc2013jcyjA70015);重慶市教育委員會科學(xué)技術(shù)研究項目(KJ080407)

張繼祥(1971—),男,山東單縣人,教授,博士,主要從事材料加工和機(jī)械工程方面的研究。E-mail: jixiangzhang@163.com。

周 偉(1990—),男,安徽碭山人,碩士研究生,主要從事機(jī)械工程方面的研究。E-mail: 42528256@qq.com。

U671.6;TG 161

A

1674-0696(2017)11-106-05

(責(zé)任編輯:劉韜)

猜你喜歡
有限元
基于擴(kuò)展有限元的疲勞裂紋擴(kuò)展分析
非線性感應(yīng)加熱問題的全離散有限元方法
TDDH型停車器制動過程有限元分析
新型有機(jī)玻璃在站臺門的應(yīng)用及有限元分析
基于I-DEAS的履帶起重機(jī)主機(jī)有限元計算
基于有限元模型對踝模擬扭傷機(jī)制的探討
10MN快鍛液壓機(jī)有限元分析
磨削淬硬殘余應(yīng)力的有限元分析
基于SolidWorks的吸嘴支撐臂有限元分析
箱形孔軋制的有限元模擬
上海金屬(2013年4期)2013-12-20 07:57:18
主站蜘蛛池模板: 国产日本一线在线观看免费| 亚洲无线视频| 日韩少妇激情一区二区| 亚洲高清无码精品| 精久久久久无码区中文字幕| 国产爽歪歪免费视频在线观看| 亚洲综合色婷婷| 九九热精品视频在线| 国内熟女少妇一线天| 黄色网址免费在线| 亚洲一区毛片| 最新日本中文字幕| 成人国产一区二区三区| 久久99精品久久久久纯品| 91精品国产自产在线老师啪l| 久草热视频在线| 久草国产在线观看| 一级毛片网| 亚洲国产在一区二区三区| 91成人免费观看| 999国内精品视频免费| 高清欧美性猛交XXXX黑人猛交 | 麻豆精品国产自产在线| 国产成人av一区二区三区| 日韩成人午夜| 成年女人a毛片免费视频| 亚洲av片在线免费观看| 91精品日韩人妻无码久久| 怡春院欧美一区二区三区免费| 99re视频在线| 中文字幕精品一区二区三区视频 | 露脸真实国语乱在线观看| 香蕉eeww99国产在线观看| 在线欧美a| 国产全黄a一级毛片| 国产精品第一区| 伊人久综合| 在线免费看黄的网站| 日韩精品亚洲人旧成在线| 亚洲女同一区二区| 久久久久免费看成人影片 | 久操中文在线| 999国内精品久久免费视频| 亚洲国产第一区二区香蕉| 一级一级一片免费| 亚洲伦理一区二区| 人人爱天天做夜夜爽| 三上悠亚在线精品二区| 欧美日韩中文字幕在线| 五月天久久综合| 全部免费特黄特色大片视频| 久久综合九九亚洲一区| 伊人久久大线影院首页| 日本精品中文字幕在线不卡| 久久无码高潮喷水| 亚洲天堂精品在线| 视频一本大道香蕉久在线播放| 国产又粗又猛又爽视频| 国产精品专区第一页在线观看| 四虎AV麻豆| 22sihu国产精品视频影视资讯| 91精品国产自产在线观看| 国产成人综合网| 欧美va亚洲va香蕉在线| 五月天福利视频| 色AV色 综合网站| 在线国产资源| 欧美日韩一区二区在线免费观看| 亚洲天堂在线免费| 中文字幕欧美日韩| 无码一区中文字幕| 人妻21p大胆| 国产在线精品香蕉麻豆| 99在线视频免费| 1024国产在线| 在线看国产精品| 国产久草视频| 亚洲另类色| 亚洲日韩精品无码专区97| 欧美日韩一区二区三区在线视频| 国产麻豆精品在线观看| 伊人五月丁香综合AⅤ|