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

一種基于應變能理論的加速蠕變本構(gòu)模型

2014-06-07 05:55:21沈才華王文武
煤炭學報 2014年11期
關(guān)鍵詞:理論模型

沈才華,張 兵,王文武

(1.河海大學土木與交通學院,江蘇南京 210098;2.中國礦業(yè)大學深部巖土力學與地下工程國家重點實驗室,江蘇徐州 221008)

一種基于應變能理論的加速蠕變本構(gòu)模型

沈才華1,2,張 兵1,王文武1

(1.河海大學土木與交通學院,江蘇南京 210098;2.中國礦業(yè)大學深部巖土力學與地下工程國家重點實驗室,江蘇徐州 221008)

采用連續(xù)介質(zhì)應變能理論分析蠕變宏觀力學行為突變的過程,定義臨界應變能密度值作為預測加速蠕變發(fā)生時刻的控制參數(shù),采用西原正夫元件模型與Perzyna黏塑性理論相結(jié)合,考慮應力狀態(tài)對加速蠕變的影響,用過屈服應力比函數(shù)反映加速階段蠕變應變速率變化,建立了一種能描述蠕變3個階段全過程的加速蠕變本構(gòu)模型。該方法確定的加速蠕變時刻不僅能反映累積蠕變應變量的影響,而且能有效反映應力狀態(tài)對加速蠕變的影響,為實現(xiàn)預測加速蠕變提供了可能,并且確定本構(gòu)模型參數(shù)所需試驗數(shù)量大大減少,為實際運用提供了方便。最后結(jié)合三軸蠕變試驗,采用D-P屈服準則結(jié)合塑性流動法則,分析了不同流動法則對預測加速蠕變發(fā)生時刻的影響規(guī)律,數(shù)值模擬結(jié)果與試驗數(shù)據(jù)基本吻合。

應變能理論;蠕變本構(gòu)模型;應變能密度;長期強度

巖石流變包括蠕變、應力松弛、長期強度等現(xiàn)象。其中,蠕變[1]特指恒定應力作用下應變隨時間持續(xù)增長的現(xiàn)象,因其與工程實踐聯(lián)系緊密而成為流變研究領(lǐng)域中最重要的內(nèi)容。其中加速蠕變階段對工程穩(wěn)定性起主要決定性作用,因此引起了學術(shù)界的重視,Vёiìalov考慮應力與應變速率之間的非線性關(guān)系,認為Bingham模型中黏滯元件的黏滯系數(shù)不是常數(shù),而是變化的,并建立了修正Bingham定律[2];王軍保、劉新榮等[3]基于非線性流變力學理論,提出了一種非線性黏滯體,其黏滯系數(shù)是所加應力水平和蠕變時間的函數(shù),有效的模擬了等速蠕變狀態(tài);宋勇軍等[4]借鑒經(jīng)典元件模型的建模思路,將巖石的初始屈服強度作為蠕變硬化的應力閾值,巖石的長期強度作為損傷軟化的應力閾值,引入能反映巖石硬化效應的硬化函數(shù)和損傷效應的損傷變量,建立能夠全面反映蠕變機制的巖石非線性蠕變模型。

總的說來,目前巖石流變的理論研究,國外研究多集中于冪律本構(gòu)或者內(nèi)變量模型,而國內(nèi)多數(shù)基于元件理論,即按照巖石的彈性、塑性和黏滯性設(shè)置基本元件(Hooke體、St.Venant體、Newton體)[5]建立由其串聯(lián)或并聯(lián)組合而成的模型,目前使用較多的有Kelvin模型、Burgers模型、西原模型、Sterpi和Gioda模型等[6-7]。其中,西原模型[1]因簡單且能夠比較清楚且全面地反映流變中的彈-黏彈-黏塑性而受到廣大研究者的重視。

對于蠕變加速階段如何判別和模擬研究還存在爭議,目前大部分得到的都是經(jīng)驗性公式,理論研究不夠,有的需要很多試驗才能獲得其本構(gòu)公式,結(jié)果很難推廣。Sih.G.C[8]采用應變能理論對裂紋擴展特性進行了研究,很好的解釋了裂紋開展引起的宏觀應力應變過程,成為分析含裂紋體應變場變化的有效方法,因此本文基于應變能理論,研究蠕變斷裂過程,提出了基于Perzyna黏塑性理論結(jié)合西原模型的一種非線性蠕變統(tǒng)一本構(gòu)模型。

1 基于應變能理論的加速蠕變模型

傳統(tǒng)的西原模型由胡克體(H),黏彈性體(N/ H)和黏塑性體(N/St.V)組成,一維的西原本構(gòu)方程可表示為

式中,E0為彈性模量;E1為黏彈性模量;η1和η2為黏滯系數(shù);σ為總應力;ε為總應變;σs為長期強度; t0為初始屈服時刻。

西原模型中黏塑性階段采用應變率為常數(shù)進行模擬,不能反映黏塑性加速蠕變,而Perzyna黏塑性理論[9-10]能很好地模擬材料的非線性黏塑性力學行為,其黏塑性應變率本構(gòu)方程定義如下:

根據(jù)流變學理論,決定蠕變力學行為特征的主要因素應從材料內(nèi)部考慮,根據(jù)應變能理論[11-13],材料的力學行為改變認為與其內(nèi)能變化有關(guān),可以采用應變能密度變化來描述。因此假設(shè)瞬時彈塑性應變與蠕變應變(總應變減去瞬時應變)產(chǎn)生的內(nèi)部機理不同,即蠕變應變能密度變化控制了材料的蠕變行為的變化。結(jié)合巖石蠕變?nèi)^程特點,采用臨界應變能密度值作為巖石蠕變過程中不同階段的控制閥值,因此定義非零等速蠕變臨界應變能密度值,其下限值為,當不同應力狀態(tài)情況下,只要其應變能密度值超過其臨界值,則巖石就將進入非零等速蠕變階段。同理定義加速蠕變臨界應變能密度值,其下限值為,當不同應力狀態(tài)情況下,只要其應變能密度值超過其臨界值,則巖石就將進入加速蠕變階段(圖 1)。因此式(3)中的φ(F)函數(shù)可用蠕變應變能密度函數(shù)[8,14-15]重新定義如下:

式中,F為屈服函數(shù)。

圖1 蠕變臨界應變能密度控制閥值示意Fig.1 Schematic diagram of control threshold creep strain energy density

由圖1所示,巖石宏觀蠕變應變變化過程一般分為3個階段:第1段為衰減蠕變階段,如果蠕變應變能密度值小于非零等速蠕變臨界應變能密度值,巖石將不會進入等速蠕變階段;第3階段為巖石應變能密度值大于等于加速蠕變臨界應變能密度值,巖石進入加速蠕變階段;第2階段處于兩者之間為等速蠕變階段。根據(jù)定義,非零等速蠕變臨界蠕變應變能密度計算公式可以寫為

其中,F?(σ)為非零等速蠕變對應的應力狀態(tài)函數(shù);τ1為非零等速蠕變起始時刻。對于一維恒定荷載下的蠕變過程,式(5)可簡寫為

式中,σi為能進入非零等速蠕變的應力,即大于長期強度荷載時對應的應力;為非零等速蠕變起始時刻對應的蠕變應變值。

因此式(2)改寫為

根據(jù)Perzyna黏塑性理論,{m}可取塑性流動方向,當采用相關(guān)聯(lián)流動法則[16]時,可定義為與塑性流動方向一致,一維狀態(tài)時等于1(為了容易理解,下面仍在公式中保留{m})。其中W為蠕變應變能密度值,一維情況下可簡化為

其中,εs(σ,t)為蠕變應變,在漸進破壞過程中是一個復雜函數(shù),如果只在某一恒定加載條件下(σ為常數(shù))進行蠕變過程研究,則εs(σ,t)可簡化為時間的函數(shù),可通過擬合試驗結(jié)果獲得εs(t)的表達式,假設(shè)采用多項式進行無限逼近蠕變過程曲線可得

如果k=0,則f(t)=a0,方程(12)簡化為Perzyna黏塑性本構(gòu)模型[5]。

因此將式(12)的黏塑性方程與西原模型結(jié)合可得改進的西原模型黏彈性、黏塑性階段的本構(gòu)方程如下:

因此根據(jù)不等式遞推理論,當k=0(即f(t)= a0),采用彈塑性理論中的流動法則[16]確定蠕變流動方向,則方程(13)退化為傳統(tǒng)的西原模型。

等速蠕變與加速蠕變是巖石蠕變過程中截然不同的兩種力學行為,根據(jù)前面加速蠕變臨界應變能密度(圖2)的定義,不同應力狀態(tài)下應變能密度大于時,巖石進入加速蠕變階段,此時式(11)中的f(t)明顯是時間相關(guān)的函數(shù),由于在同一種類型的應力狀態(tài)下(如單軸受拉類型、單軸受壓類型或三軸受壓類型等)巖石進入加速蠕變時的破壞模式很相似,因此假設(shè)同一種巖石在兩種不同應力狀態(tài)下的加速蠕變臨界應變能密度成函數(shù)關(guān)系,即屬于同一種類型的任意進入加速蠕變的應力狀態(tài)的臨界應變能密度可表示成式(15)。

圖2 蠕變過程臨界應變能密度示意Fig.2 Schematic diagram of the strain energy density of accelerated creep

試驗表明巖石在同一類型應力狀態(tài)下(如三軸受壓)的破壞模式基本是一致的,因此可以采用靜力學連續(xù)介質(zhì)應變能理論進行宏觀應力應變分析[10]。應變能密度是標量值,如果把巖石視為連續(xù)均勻介質(zhì),則作為同一類型應力狀態(tài)下的加速蠕變應變能密度的臨界值(標量)不應受應力方向的改變而變化,因此可以視為常數(shù),即式(16)中f(σi, σs-long)可近似取1。

因此可獲得一維狀態(tài)下基于應變能理論的能描述加速蠕變的統(tǒng)一黏彈塑性蠕變本構(gòu)方程如下:

當n=1,k=1時,f(t)=a0+a1(t-τ2),蠕變統(tǒng)一本構(gòu)方程式退化為Qinghui Jiang等[17]采用增加非線性黏滯壺的方式建立的加速蠕變本構(gòu)模型。

根據(jù)實際情況,應力狀態(tài)增加對加速蠕變的快慢也將起到影響,因此式(12)中k的取值(蠕變非線性特征參數(shù),反映了蠕變加速的快慢)應該與過屈服應力比值大小有關(guān)[18],因此有式(20)成立。

因此結(jié)合一維本構(gòu)方程推廣到三維矩陣形式的類比法,得到基于應變能理論的可以描述加速蠕變過程的簡化實用三維加速蠕變本構(gòu)模型如下:

式中,帶上劃線的表示三維情況下對應的屈服強度計算值

本構(gòu)模型的參數(shù)可以通過巖石蠕變試驗獲得。假設(shè)采用 Drucker-Prager屈服準則,則

2 三軸蠕變試驗的加速蠕變本構(gòu)模型模擬分析

對某砂巖(圖3)進行三軸蠕變試驗,軸向壓力61 MPa,圍壓6.0 MPa,采用最小二乘法,結(jié)合方程(21)~(23)進行擬合獲得本構(gòu)模型參數(shù)。采用D-P屈服準則和相關(guān)聯(lián)流動法則建立三維本構(gòu)模型(本構(gòu)模型Ⅰ),根據(jù)實測曲線擬合的參數(shù)見表1,實測數(shù)據(jù)與本構(gòu)模型Ⅰ的擬合曲線如圖4所示。采用D-P屈服準則和非相關(guān)聯(lián)流動法則(α=0)建立三維本構(gòu)模型(本構(gòu)模型Ⅱ)。實測加速蠕變發(fā)生時刻(τ2=36.2 h)的加速蠕變應變能密度值為0.169 J/m3。

圖3 砂巖試樣Fig.3 Sandstone specimens

表1 本構(gòu)模型計算參數(shù)Table 1 Parameters of the triaxial creep tests and the parameters of constitutive model

圖4 實測值與本構(gòu)模型Ⅰ擬合曲線Fig.4 The data of triaxial test and the fitting curve of modified Nishihara model

表1顯示,采用不同的流動法則時對預測加速蠕變有較大影響。采用本構(gòu)模型Ⅰ計算的實際加速蠕變起始時刻(τ2=36.2 h)對應的加速蠕變應變能密度值(0.149 9 J/m3)小于實測值(0.169 J/m3);采用本構(gòu)模型Ⅱ計算的實際加速蠕變起始時刻(τ2= 36.2 h)對應的加速蠕變應變能密度值(0.184 3 J/m3)大于實際值(0.169 J/m3)。圖5顯示采用本構(gòu)模型Ⅰ(α=0.331,相關(guān)聯(lián)流動法則)預測的加速蠕變發(fā)生時刻比實際的晚,模擬結(jié)果偏于危險;采用本構(gòu)模型Ⅱ(α=0)預測的加速蠕變發(fā)生時刻比實際的早,模擬結(jié)果偏于安全,可見非相關(guān)聯(lián)流動法則的膨脹角對預測結(jié)果有很大影響,這與實際情況一致。因此取非相關(guān)聯(lián)流動法則α=0.232時作為本巖樣的本構(gòu)模型(本構(gòu)模型Ⅲ)。圖5顯示本構(gòu)模型Ⅲ模擬計算曲線與實測值基本一致。

為了驗證本構(gòu)模型的合理性,采用本構(gòu)模型Ⅲ對不同軸向應力作用下的蠕變?nèi)^程進行預測,結(jié)果如圖6所示。

圖5 不同本構(gòu)模型模擬值與實測值對比曲線Fig.5 The data of triaxial creep test and the different constitutive model simulation curves

圖6 軸向蠕變應變-時間曲線(本構(gòu)模型Ⅲ,不同軸壓)Fig.6 The time-strain curves of creep(constitutive modelⅢ and different axial pressures)

圖6顯示基于應變能理論建立的蠕變本構(gòu)模型能有效描述蠕變?nèi)A段全過程,預測不同軸向應力下的蠕變過程規(guī)律與實際情況基本相似。

3 結(jié) 論

(1)采用非線性黏塑性理論研究與元件方法相結(jié)合建立巖石蠕變的有限元本構(gòu)模型,能很好地反映巖石蠕變?nèi)^程,為實現(xiàn)巖石蠕變過程數(shù)值模擬分析提供了一種途徑。

(2)基于連續(xù)介質(zhì)應變能理論,采用臨界應變能密度值描述材料發(fā)生蠕變力學行為突變的特征量,該方法確定的加速蠕變時刻不僅能反映累積蠕變應變量的影響,而且能有效反映應力狀態(tài)對加速蠕變的影響,為實現(xiàn)預測加速蠕變提供了可能,并且確定本構(gòu)模型參數(shù)所需試驗數(shù)量大大減少,為實際運用提供了方便。

(3)最后通過數(shù)值模擬與試驗結(jié)果對比分析說明:采用本文方法獲得的巖石蠕變本構(gòu)方程能有效模擬衰減蠕變、等速蠕變和加速蠕變的全過程,采用非相關(guān)聯(lián)流動法則的蠕變本構(gòu)模型預測加速蠕變發(fā)生時刻比實測值早,相對偏于安全。本文提出的用蠕變臨界應變能密度值預測加速蠕變發(fā)生時刻的方法為研究巖石蠕變破壞的宏觀本構(gòu)模型提供了新的思路。

[1] Nishihara M.Creep of shale and sandy-shale[J].The Journal of Geological Magazine of Japan,1952,58:373-377.

[2] 齊亞靜,姜清輝,王志儉,等.改進西原模型的三維蠕變本構(gòu)方程及其參數(shù)辨識[J].巖石力學與工程學報,2012,31(2):347-355.

Qi Yajing,Jiang Qinghui,Wang Zhijian,et al.3D creep constitutive equation of modified Nishihara model and its parameters identification[J].Chinese Journal of Rock Mechanics and Engineering,2012, 31(2):347-355.

[3] 王軍保,劉新榮,郭建強,等.鹽巖蠕變特性及其非線性本構(gòu)模型[J].煤炭學報,2014,39(3):445-451.

Wang Junbao,Liu Xinrong,Guo Jianqiang,et al.Creep properties of salt rock and its nonlinear constitutive model[J].Journal of China Coal Society,2014,39(3):445-451.

[4] 宋勇軍,雷勝友,劉向科.基于硬化和損傷效應的巖石非線性蠕變模型[J].煤炭學報,2012,37(S2):287-292.

Song Yongjun,LeiShengyou,Liu Xiangke.Non-linearrock creep model based on hardening and damage effect[J].Journal of China Coal Society,2012,37(S2):287-292.

[5] De Borst R,Crisfield M A,Remmers J J C,et al.Time-dependent material models[A].Non-Linear Finite Element Analysis of Solids and Structures,Second Edition[C].2012:281-304.

[6] Tang M,Wang Z.Experimental study on rheological deformation and stress properties of limestone[J].Journal of Central South University of Technology,2008,15(1):475-478.

[7] Sterpi D,Gioda G.Visco-plastic behaviour around advancing tunnels in squeezing rock[J].Rock Mechanics and Rock Engineering, 2009,42(2):319-339.

[8] Sih G C.From monoscale to multiscale modeling of fatigue crack growth:Stress and energy density factor[J].Science China Physics,Mechanics and Astronomy,2014,57(1):39-50.

[9] Zhou H W,Wang C P,Han B B,et al.A creep constitutive model for salt rock based on fractional derivatives[J].International Journal of Rock Mechanics and Mining Sciences,2011,48(1):116-121.

[10] Zhou H W,Wang C P,Mishnaevsky Jr L,et al.A fractional derivative approach to full creep regions in salt rock[J].Mechanics of Time-Dependent Materials,2013,17(3):413-425.

[11] Zhu Hehua,Ye Bin,Cai Yongchang,et al.An elasto-viscoplastic model for soft rock around tunnels considering overconsolidation and structure effects[J].Computers and Geotechnics,2013,50:6-16.

[12] Xu W Y,Wang R B,Wang W,et al.Creep properties and permeability evolution in triaxial rheological tests of hard rock in dam foundation[J].Journal of Central South University of Technology, 2012,19(1):252-261.

[13] Shibata K,Tani K,Okada T.Creep behavior of tuffaceous rock at high temperature observed in unconfined compression test [J].Soils and Foundations,2007,47(1):1-10.

[14] Huang Y,Zhang Y,Zhang M,et al.Method for predicting the failure load of masonry wall panels based on generalized strain-energy density[J].Journal of Engineering Mechanics,2014,140(8):1-12.

[15] Ayatollahi M R,Sedighiani K.Mode I fracture initiation in limestone by strain energy density criterion[J].Theoretical and Applied Fracture Mechanics,2012,57(1):14-18.

[16] Lubliner J.Plasticity theory[M].New York:Dover Publications, 2008.

[17] Jiang Qinghui,Qi Yajing,Wang Zhijian,et al.An extended Nishihara model for the description of three stages of sandstone creep [J].Geophysical Journal International,2013,193(2):841-854.

[18] Rutter E H.On the creep testing of rocks at constant stress and constant force[J].International Journal of Rock Mechanics and Mining Sciences&Geomechanics Abstracts,1972,9(2):191-195.

[19] 鄭穎人,沈珠江,龔曉南.巖土塑性力學原理[M].北京:中國建筑工業(yè)出版社,2002.

Zheng Yingren,Shen Zhujiang,Gong Xiaonan.The principles of geotechnical plastic mechanics[M].Beijing:China Architecture& Building Press,2002.

A new accelerated creep constitutive model based on the strain energy theory

SHEN Cai-hua1,2,ZHANG Bing1,WANG Wen-wu1
(1.College of Civil and Transportation Engineering,Hohai University,Nanjing 210098,China;2.State Key Laboratory of Geomechanics and Deep Underground Engineering,China University of Mining and Technology,Xuzhou 221008,China)

The anthors analyzed the mutation process of the macro mechanical creep behavior from the perspective of the continuous strain energy theory,and defined the critical strain energy density to judge the occurrence time of accelerated creep.A accelerated creep constitutive model which can describe the whole process of transient creep,steady state creep,and accelerated creep was established by combining the Perzyna Visco-plasticity theory and Nishihara model,taking into account the influence of stress state on the accelerated creep and using the yield stress ratio to reflect the creep strain rate change of accelerated creep.The determined accelerating creep occurring time using this method not only reflectes the effect of cumulative strain of the creep,but also the influence of stress states on the accelerating creep,which provides possibility of prediction of the accelerating creep.In addition,the number of experiments to determine constitutive model parameters was significantly reduced,which is more convenient in practical applications.Finally,according to the triaxial creep test,the law that different flow rules affect the occurrence time of accelerated creep was analyzed by adopting Drucker-Prager yield criterion and plastic flow rule.As the simulation result is basically consistent with the test data,the research results provide a new method for the study of rock creep fracture process.

strain energy theory;creep constitutive model;strain energy density;long-term strength

TU452

A

0253-9993(2014)11-2195-06

2014-04-03 責任編輯:王婉潔

江蘇省自然科學基金面上資助項目(BK20141419);深部巖土力學與地下工程國家重點實驗室開放基金資助項目(SKLGDUEK1211);中央高校基本科研業(yè)務費專項資金資助項目(2014B04914)

沈才華(1976—),男,江蘇宜興人,講師,博士。E-mail:shencaihua@163.com

沈才華,張 兵,王文武.一種基于應變能理論的加速蠕變本構(gòu)模型[J].煤炭學報,2014,39(11):2195-2200.

10.13225/j.cnki.jccs.2014.0445

Shen Caihua,Zhang Bing,Wang Wenwu.A new accelerated creep constitutive model based on the strain energy theory[J].Journal of China Coal Society,2014,39(11):2195-2200.doi:10.13225/j.cnki.jccs.2014.0445

猜你喜歡
理論模型
一半模型
堅持理論創(chuàng)新
當代陜西(2022年5期)2022-04-19 12:10:18
神秘的混沌理論
理論創(chuàng)新 引領(lǐng)百年
相關(guān)于撓理論的Baer模
重要模型『一線三等角』
重尾非線性自回歸模型自加權(quán)M-估計的漸近分布
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
理論宣講如何答疑解惑
學習月刊(2015年21期)2015-07-11 01:51:44
主站蜘蛛池模板: 99re精彩视频| 精品国产污污免费网站| 激情在线网| 九九九精品成人免费视频7| 久久国产乱子| 亚洲不卡网| 日本伊人色综合网| a国产精品| 亚洲一区二区三区国产精品| vvvv98国产成人综合青青| 国产精品视频a| 毛片最新网址| 欧美在线视频不卡第一页| 欧美午夜精品| 嫩草影院在线观看精品视频| 亚洲日韩欧美在线观看| 伊人久久影视| 在线国产毛片| 看国产毛片| 人人爱天天做夜夜爽| 91成人试看福利体验区| 国产www网站| 日韩免费中文字幕| 草逼视频国产| 久久精品人人做人人爽| 91黄色在线观看| 91色老久久精品偷偷蜜臀| 精品国产网站| 91在线国内在线播放老师| 人人看人人鲁狠狠高清| 中文字幕精品一区二区三区视频| 成人亚洲视频| 国产精品大尺度尺度视频| 精品99在线观看| 亚洲精品视频免费看| 91视频青青草| 国产亚洲欧美日韩在线一区二区三区| 综合色区亚洲熟妇在线| 欧美精品一区二区三区中文字幕| 国产一区二区三区夜色| 欧美福利在线观看| 日韩欧美一区在线观看| 日韩国产综合精选| 18黑白丝水手服自慰喷水网站| 激情无码字幕综合| av在线手机播放| 欧美一级专区免费大片| 国产精品香蕉在线| 亚洲成人一区二区| 亚洲精品成人片在线观看| 在线看免费无码av天堂的| 重口调教一区二区视频| 久久久久久国产精品mv| 国产乱子伦精品视频| 2024av在线无码中文最新| 免费va国产在线观看| 在线观看国产黄色| 欧美一区日韩一区中文字幕页| 国产JIZzJIzz视频全部免费| 熟妇丰满人妻| 在线免费无码视频| 久久国产热| 日韩东京热无码人妻| www.91中文字幕| 91色老久久精品偷偷蜜臀| 狠狠做深爱婷婷久久一区| 国产微拍精品| 久久精品国产精品青草app| 网久久综合| 国产日本一区二区三区| 国产a v无码专区亚洲av| 日韩福利在线观看| 亚洲中文字幕97久久精品少妇| 五月婷婷综合色| 免费播放毛片| 国产青青草视频| 亚洲国产综合自在线另类| 亚洲精品麻豆| 欧美在线网| 亚欧美国产综合| 日韩专区欧美| 久久五月视频|