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

基于Matlab的三維四步編織仿真

2016-08-06 07:28:49章宇界閻建華
關(guān)鍵詞:仿真

章宇界,閻建華

(東華大學(xué) a. 紡織學(xué)院;b. 研究院,上海 201620)

?

基于Matlab的三維四步編織仿真

章宇界a,閻建華b

(東華大學(xué) a. 紡織學(xué)院;b. 研究院,上海 201620)

摘要:根據(jù)編織參數(shù)對四步法編織進(jìn)行算法演繹,得到編織體的具體工藝參數(shù),然后利用Matlab 編寫腳本追蹤攜紗器的運動,確定紗線運動的空間位置趨勢并存儲到矩陣中,編寫B(tài)均勻樣條曲線腳本并在CATIA軟件中導(dǎo)入空間紗線位置點并生成樣條曲線,實現(xiàn)三維四步法編織的仿真以及矩形編織和圓形編織建模.以4×4矩形編織、3×18圓形編織為例,對具體編織參數(shù)展開計算,對紗線軌跡算法進(jìn)行優(yōu)化,實現(xiàn)編織預(yù)制件的仿真.

關(guān)鍵詞:三維編紗; 紗線軌跡; 編織參數(shù); Matlab; 仿真

三維編織是20世紀(jì)80年代興起的新型紡織技術(shù),不同于其他三維織物,三維編織中的紗線運動較為復(fù)雜,所有紗線按一定規(guī)律相互交織,相互之間交纏緊密,構(gòu)成三維一體成型織物.直接觀察三維編織物結(jié)構(gòu)存在一定難度,而仿真模擬不僅提供了可視化的模型,可根據(jù)編織工藝參數(shù)變化觀察織物的表觀形態(tài),還能進(jìn)行有限元分析,從而得到最優(yōu)參數(shù),節(jié)約大量的時間和成本.充分利用計算機(jī)強(qiáng)大的圖形處理功能,建立數(shù)字化的幾何模型,對三維復(fù)合材料的力學(xué)結(jié)構(gòu)研究具有重大意義[1].目前,常用的三維編織方法為二步法和四步法,其中四步法可以編織許多不同截面的結(jié)構(gòu).當(dāng)前用四步法對矩形、圓形編織的空間紗線軌跡模擬及預(yù)制件仿真還不夠全面,大部分還停留在模擬矩形編織紗線軌跡上,對可視化三維建模特別是圓形編織的研究還不夠深入[1-6].本文通過Matlab軟件與CATIA軟件聯(lián)合建模,采用B均勻樣條原理實現(xiàn)矩形編織和圓形編織的空間紗線軌跡模擬及其預(yù)制件仿真.

1矩形編織

1.1編織參數(shù)

調(diào)整編織參數(shù)會影響織物緊密程度、紗線體積分?jǐn)?shù),進(jìn)而影響織物性能,而編織參數(shù)也給產(chǎn)品的可設(shè)計性提供契機(jī).一般的編織參數(shù)有表面編織角、預(yù)制件尺寸和編織紗線數(shù)等.

參照文獻(xiàn)[2-3],對于矩形編織, 假設(shè)載紗器橫向運動一步, 編織紗線運動的距離為X, 縱向運動一步, 編織紗線運動的距離為Y, 每縱向或橫向運動一步后, 制件成型方向前進(jìn)的節(jié)距為Z.假設(shè)編織樣品寬度為W,厚度為T,行列數(shù)分別為m和n,表面編織角為α,水平角為δ,編織花節(jié)長度為h,則單胞尺寸可表示為

(1)

給定編織樣品寬度W與厚度T,可得

(2)

(3)

此處的W為左邊緣紗線與右邊緣紗線的中心距,編織體實際寬度應(yīng)該為W+d,其中d為紗線直徑,如圖1所示.

圖1 編織體寬度Fig.1 Width of braided preform

1.2算法設(shè)計

攜紗器底盤上的每個攜紗器都有各自的平面位置,且攜紗器的運動軌跡是連續(xù)的,故紗線的軌跡也是可追蹤的.基本原理:以i、j分別表示當(dāng)前某攜紗器處于i排j列并記錄每次運動一步后的平面位置,為便于區(qū)分,每個攜紗器以數(shù)字標(biāo)記序號.只要追蹤每個數(shù)字所代表的攜紗器,就可記錄其運動位置的變化,同時設(shè)置編織速度,得到紗線的空間位置變化即紗線軌跡.以4×4矩形編織為例,需要的紗線數(shù)S為

S=m×n+m+n

(4)

故此例中紗線數(shù)S=24,其紗線排布如圖2所示.

圖2 4×4矩形編織攜紗器排布Fig.2 Yarn carriers align of 4×4 rectangular braiding

以矩陣形式表示該攜紗器初始位置如圖3所示,其中0表示攜紗器上無紗線,1~24分別為24個有效攜紗器.

圖3 攜紗器初始位置矩陣

攜紗器按一定規(guī)律運動,相鄰兩列或兩排的攜紗器運動方向是相反的.如圖4所示,攜紗器進(jìn)行第一步運動時,第三、五排深色攜紗器向左移動一步,第二、四排淺色攜紗器向右移動一步.

圖4 攜紗器第一步排布Fig.4 Align of the first step of yarn carriers

可用矩陣表示攜紗器第一步位置如圖5所示.

圖5 攜紗器第一步位置矩陣

第二步時,第二、四列攜紗器向下移動一步,第三、五列攜紗器向上移動一步.攜紗器運動直到第四步完成一個循環(huán),第二到第四步攜紗器排布矩陣變化過程分別如圖6~8所示.

圖6 攜紗器第二步位置矩陣

圖7 攜紗器第三步位置矩陣

圖8 最終攜紗器排布矩陣

紗線的運動軌跡數(shù)即組數(shù)G由m和n確定,計算式為

(5)

其中:LCM為兩個數(shù)的最小公倍數(shù)[4].

相同運動軌跡紗線數(shù)量N可由S/G得到,以“3”攜紗器為例,一個循環(huán)中其X-Y矩陣位置變化為(2,1)→(2,2)→(3,2)→(3,3)→(4,3).與“3”攜紗器有相同運動軌跡的攜紗器共有4個,與此運動軌跡不同的有3種,共4種運動軌跡.對紗線的運動軌跡進(jìn)行歸類,可簡化對編織體的分析.

根據(jù)X,Y,Z的賦值,“3”攜紗器上的紗線空間位置(2X,Y,0)→(2X,2Y,h/4)→(3X,2Y,h/2)→(3X,3Y,3h/4) →(4X,3Y,h).對于紗線空間位置的變化,可以用矩陣來記錄.X{m}、Y{m}可表示攜紗器在平面內(nèi)的運動位置變化,可轉(zhuǎn)化為紗線的空間位置變化矩陣M.實際上,經(jīng)過上述步驟得到的是構(gòu)成軌跡曲線的點的集合,需要導(dǎo)入建模軟件中生成軌跡曲線和仿真建模.

在Matlab中創(chuàng)建腳本,記錄各攜紗器的平面位置走向,腳本如下:

form=1:S

[i,j]=find(A==m);

X{m}=[X{m}i];

Y{m}=[Y{m}j];

end

其中:A為攜紗器平面位置矩陣;X{m}、Y{m}記錄第m個攜紗器的橫向與縱向位置.將X{m}、Y{m}轉(zhuǎn)化為紗線的空間位置變化,同時存儲到矩陣集合中,得到24根紗線的空間位置矩陣,利用循環(huán)語句,可計算出M,用Matlab腳本可表示為

M=[];

fori=1:S

M=[M;[x1 *X{i}′,y1*Y′,z1*z′]];

end

其中:M為紗線空間運動位置點的集合矩陣;x1和y1分別代表紗線橫向、縱向運動一步的步長;z1表示節(jié)距.

1.3編織體模擬

軟件運行環(huán)境為Windows 7(32位),計算機(jī)CPU為雙核2 GHz,計算機(jī)內(nèi)存為2 GB.

令腳本運行兩個編織循環(huán),記錄紗線每運動一步所處的位置并導(dǎo)入Excel表格中,運行時間為4.07 s, 使用CATIA的宏功能,將每條紗線的軌跡點進(jìn)行插值連接,得到的紗線軌跡如圖9所示.賦予這些紗線一定的截面形狀,生成仿真編織體.若截面面積過大,則致使紗線之間產(chǎn)生穿插,如圖10所示,取紗線截面為圓形,直徑d為0.7mm,觀察可發(fā)現(xiàn),紗線的屈曲較多,編織體較為松弛且纖維體積分?jǐn)?shù)較低,不太符合實際的交織運動軌跡.

圖9 4×4矩形編織體紗線軌跡Fig.9 Yarn traces of 4×4 rectangular braided preform

圖10 4×4矩形編織體渲染圖Fig.10 Rendering figure of 4×4 rectangular braided preform

1.4算法優(yōu)化

在實際編織過程中,紗線是有張力的,不會產(chǎn)生圖9所展現(xiàn)出的松弛結(jié)構(gòu),由于拉力的作用及打緊作用,紗線會自動地沿著最短的路徑展開,折線及松弛曲線形式的路徑是不存在的,表面的紗線緊繞包纏預(yù)制體,內(nèi)部紗線沿最短路線交織,故需對算法進(jìn)行優(yōu)化.文獻(xiàn)[2,5]采用Bézier擬合紗線軌跡,Bézier曲線可以從多邊形的形狀來預(yù)測所要產(chǎn)生曲線的形狀,可通過對多邊形定點的控制來改變曲線的形狀,但在Bézier曲線中任何一個控制頂點發(fā)生改變,會使得整條Bézier曲線形狀發(fā)生變化.基于此,文獻(xiàn)[7]提出了B樣條曲線方法,繼承了Bézier曲線優(yōu)點的同時克服了上述缺點[8],其優(yōu)勢在于曲線可以分段,任意改變一個控制頂點只會影響相鄰曲線的形狀而不影響整條曲線的走向,有效解決了局部控制和連接問題.文獻(xiàn)[6,9-10]使用B樣條曲線進(jìn)行擬合. 因為矩陣M中控制點的分布較為均勻,故采用三次B均勻樣條曲線擬合紗線軌跡.三次B均勻樣條曲線需要4個控制頂點確定三次多項式曲線,計算式為

(6)

其中:pi為控制頂點即空間坐標(biāo)值點;Ni,3為節(jié)點向量.用矩陣形式可表示為

(7)

其中:u∈[0,1].三次B均勻樣條曲線并不經(jīng)過控制頂點,這與實際結(jié)構(gòu)中紗線因相互作用而產(chǎn)生偏移的情況符合,故用B均勻樣條曲線擬合紗線運動軌跡.如圖11所示,仍取紗線直徑d=0.7mm, 編織預(yù)制件經(jīng)優(yōu)化后變得緊湊,由其中某根纖維束的運動趨勢可見紗線屈曲得到改善,仿真預(yù)制件纖維體積分?jǐn)?shù)得到提升,與真實預(yù)制件更為接近.

圖11 優(yōu)化模型Fig.11 Optimized model

在三維編織過程中,預(yù)制件的纖維分?jǐn)?shù)、纖維束的走向等主要由編織角的大小決定,編織角增加,花節(jié)長度減小,仿真模型變得更加緊湊,同時紗線張力也相應(yīng)增加,使得預(yù)制件在高度方向被壓縮.在實際編織中,編織角主要由打緊作用控制,打緊后編織體內(nèi)部紗線的屈曲減少,使得紗線張力增大,打緊力越大,紗線間空隙減少,結(jié)構(gòu)也越為緊湊,編織體不易松弛和垮塌.取W=10mm,T=6mm,m=6,n=10, α=40°,d=0.6mm,進(jìn)行3個循環(huán)后其仿真模型如圖12所示,運行時間為4.46s.

圖12 6×10矩形編織體Fig.12 6×10 rectangular braided preform

在生成空間編織體時,需考慮編織體的整體緊密性.換言之,編織參數(shù)例如編織角的選取需要參考實際編織體結(jié)構(gòu)參數(shù),編織體應(yīng)該是紗線緊密排布,編織過程中紗線間相互擠壓但不穿插以保持編織體不松弛、不垮塌.同時紗線截面的選擇在一定程度上也影響緊密性及編織復(fù)合材料的纖維體積分?jǐn)?shù).

2圓形編織

2.1編織參數(shù)

三維圓形編織的優(yōu)勢在于能夠?qū)A形管件和異形件進(jìn)行一體成型,制成的復(fù)合材料具有優(yōu)異的力學(xué)性能. 當(dāng)前對于圓形編織建模的研究較少,圓形編織紗線軌跡模擬和編織體建模對分析圓形編織結(jié)構(gòu)、評價三維圓形編織復(fù)合材料和結(jié)構(gòu)優(yōu)化設(shè)計具有重要意義.

編織圓形管件的原理和矩形預(yù)制件是相通的,同樣是四步法成型,而不同之處在于縱向運動改為徑向運動,橫向運動變?yōu)橹芟蜻\動,紗線的排列方式與交織結(jié)構(gòu)因此發(fā)生變化.紗線在編織體外圍和內(nèi)圍間往復(fù),故編織角也是動態(tài)變化的.攜紗器周向和徑向間歇運動帶動紗線交織,管件因此成型.底盤上首列與末列相鄰,故周向不再有附加紗線,徑向附加紗線仍存在且只能進(jìn)行徑向運動.經(jīng)過四步法運動后,一個機(jī)器循環(huán)完成,重復(fù)上述編織運動,紗線將相互交織而形成一定長度的管狀編織體.管件的編織需要使用芯模以支撐編織體,芯模的外形也會對編織體的交織結(jié)構(gòu)產(chǎn)生一定的影響,這里只討論常規(guī)圓形編織.

圓形編織中,紗線數(shù)量因周向不再有附加紗線,故紗線數(shù)計算方式變?yōu)椋篠=(m+1)n,其中n為列數(shù)且必須為偶數(shù),否則必有相鄰兩列攜紗器運動軌跡相同.紗線運動軌跡數(shù)即組數(shù)G[11]為

(8)

編織角α是一個重要的編織工藝參數(shù),預(yù)制件表面形成的交織紋路角稱為編織角,其正切等于預(yù)制體表面的花節(jié)寬度與花節(jié)長度之比.假設(shè)某管件的內(nèi)外徑為din、dout,第i列最外部紗線在一個循環(huán)之后會移動到相鄰列,所走過的曲線平面投影為一圓弧,圓弧長度為2πdout/n,同理,最內(nèi)部某紗線在一個循環(huán)后所走過的平面圓弧長度為2πdin/n.故由內(nèi)外編織角的定義可知:

(9)

2.2算法設(shè)計

以3×18管件編織為例,假定din=12 mm,dout=24 mm,αout=40°,則需要的紗線數(shù)量為S=72,紗線組數(shù)G=108/18=6,相同運動軌跡攜紗器數(shù)量N=S/G=12,其攜紗器排布如圖13所示,共有6種圖案分別代表6組攜紗器,相同圖案代表具有相同運動軌跡的攜紗器.相鄰列攜紗器的夾角定義為θ0=360°/n,以[i,j]表示某攜紗器運動若干步后處于第i層和第j列,以r、θ表示該攜紗器上紗線中心與管件中心連線距離以及該連線與垂直方向的夾角,利用三角函數(shù)運算得到其面內(nèi)位置.

圖13 3×18圓形編織攜紗器分布Fig.13 Yarn carriers align of 3×18 tubular braiding

在Matlab中運行創(chuàng)建的腳本文件,讀取攜紗器初始位置和每一步之后的平面位置,并轉(zhuǎn)化為每根紗線的實際空間位置變化并記錄存儲,Matlab腳本如下:

form=1:S

[i,j]=find(A==m);

θ=(j-1)*θ0;

r=(a+2-i)*D+din/2;

x=r*sinθ;

y=r*cosθ;

X{m}=[X{m}x];

Y{m}=[Y{m}y];

end

其中:a為排數(shù);相鄰層同列紗線的中心距D= (dout-din)/2(m+1),X{m}、Y{m}記錄第m條紗線的周向與徑向位置及其變化.

在Matlab腳本程序中輸入編織參數(shù)m、n、din、dout、αout的值,得到72根紗線的空間軌跡變化集合矩陣.對該矩陣同樣進(jìn)行B均勻樣條曲線擬合后的 紗線軌跡如圖14所示,取紗線直徑d=1.5 mm,仿真效果如圖15所示,運行時間為4.24 s.

圖14 3×18管件編織體紗線軌跡Fig.14 Yarn traces of 3×18 tubular braided preform

圖15 3×18管件編織體仿真Fig.15 Simulation of 3×18 tubular braided preform

取6×36管件編織進(jìn)行仿真,假定din=12mm,dout=24mm,αout=30°,紗線共計252根,取紗線直徑d=0.8mm,運行3個循環(huán),仿真編織預(yù)制件如圖16所示,運行時間為6.63s.

圖16 6×36管件編織體仿真Fig.16 Simulation of 6×36 tubular braided preform

實際上,紗線之間存在相互擠壓,擠壓處截面會產(chǎn)生變形,而非規(guī)則的圓形,而且預(yù)制件在經(jīng)過固化成型工藝成為復(fù)合材料后,由于與樹脂的相互作用,紗線變形進(jìn)一步加大,這與仿真結(jié)果是有差異的.另外為了避免紗線的穿插,所取紗線直徑均較小,以致仿真預(yù)制件的纖維體積分?jǐn)?shù)較低,希望在今后的算法研究中得到改善.

3結(jié)語

本文以三維矩形編織和圓形編織為例,利用Matlab記錄紗線的空間運動變化,并用B均勻樣條曲線擬合紗線軌跡重新生成軌跡點,導(dǎo)入CATIA軟件以實現(xiàn)四步法矩形和圓形編織預(yù)制件中紗線軌跡的模擬與預(yù)制件的仿真,并對關(guān)鍵算法進(jìn)行說明.使用B均勻樣條曲線進(jìn)行紗線軌跡擬合與編織預(yù)制件仿真,基本接近真實編織效果,同時矩形編織和圓形編織中紗線間的穿插效果也得到改善.

通過輸入編織工藝參數(shù),實現(xiàn)計算機(jī)編織模擬,用戶可以進(jìn)行編織參數(shù)的調(diào)整,直觀地了解產(chǎn)品的基本外觀以及紗線的空間走向并實時調(diào)整,在一定程度上實現(xiàn)可視化.

參考文獻(xiàn)

[1] 肖來元,汪博峰.三維編織復(fù)合材料幾何成型算法及仿真研究[J].計算機(jī)工程與科學(xué),2008,30(1):82-85.

[2] 邵將,溫衛(wèi)東,崔海濤.三維四步法編織復(fù)合材料結(jié)構(gòu)的計算機(jī)仿真[J].南京航天航空大學(xué)學(xué)報,2009,41(1):36-40.

[3] CHEN L, TAO X M, CHOY C L. On the micro structure of three-dimensional braided preforms [J]. Compos Sci Technol, 1999,59(3):2383-2391.

[4] 陳利,李嘉祿,李學(xué)明.三維編織中紗線的運動規(guī)律分析[J].復(fù)合材料學(xué)報,2002,19(2):71-74.

[5] 周新貴,張長瑞,龔國輝,等.三維編織物的真三維仿真[C]//2000年中國材料研討會論文集.2000:4.

[6] 何紅闖,楊連賀,陳利.矩陣組合截面四步法二次三維編織及其空間模型可視化[J].復(fù)合材料學(xué)報,2010,27(4):160-167.

[7] 張明星.廣義Bézier曲線與B樣條曲線的研究[D].長沙:中南大學(xué)數(shù)學(xué)與統(tǒng)計學(xué)院,2013.

[8] 肖來元,寇曉菲. 三維編織復(fù)合材料編織工藝過程仿真研究[J]. 計算機(jī)工程與科學(xué),2014,36(4):719-724.

[9] KANG T J, KIM S J, JUNG K. Analysis of geometrical parameters using a CAD system for a 3-D braided preform [J]. Textile Research Journal, 2008,78(10): 922-935.

[10] BARSKY B A, BEAUTY J C. Local control of bias and tension in beta-splines [J]. ACM Trans on Graphics, 1983,2(2):109-134.

[11] 曹紅蓓,王君澤,瞿暢,等.管狀三維編織物三維動畫仿真探索[J].紡織學(xué)報,2004,25(5): 71-73.

文章編號:1671-0444(2016)03-0363-07

收稿日期:2015-04-28

基金項目:上海市科委資助項目 (12521102400)

作者簡介:章宇界(1990—),男,浙江嘉興人,碩士研究生,研究方向為紡織復(fù)合材料. E-mail:zhangyujie181@sina.com 閻建華(聯(lián)系人),男,研究員,E-mail:jh_yan@dhu.edu.cn

中圖分類號:TB 332

文獻(xiàn)標(biāo)志碼:A

Simulation of Four-Step Three-Dimensional Braiding Based on Matlab

ZHANGYu-jiea,YANJian-huab

(a. College of Textiles; b. Research Institute, Donghua University, Shanghai 201620, China)

Abstract:The parameters of braiding preform is obtained by an algorithm deduction according to four-step three-dimensional (3-D) braiding parameters. The traces of yarn carriers can be acquired by a Matlab script. Through such a way, the space position trends of yarn movements are gained and can be stored into a matrix for four-step braiding. A script of uniform rational B-spline is compiled to generate points data of yarns and then plot into splines in CATIA software. By means of that, the simulation and modelling of four-step 3-D braiding are realized. 4×4 rectangular braiding and 3×18 tubular braiding are quoted as examples to simulate braiding preforms through calculation of specific braiding parameters and optimization of the algorithm for yarn traces.

Key words:three-dimensional braiding; yarn trace; braiding parameters; Matlab; simulation

猜你喜歡
仿真
Proteus仿真軟件在單片機(jī)原理及應(yīng)用課程教學(xué)中的應(yīng)用
工業(yè)機(jī)器人模擬仿真技術(shù)在職業(yè)教育中的應(yīng)用淺析
一種幫助幼兒車內(nèi)脫險應(yīng)急裝置的仿真分析
科技資訊(2016年18期)2016-11-15 20:09:22
論虛擬仿真實訓(xùn)系統(tǒng)在口腔實驗教學(xué)中的應(yīng)用
科技資訊(2016年18期)2016-11-15 07:55:28
基于機(jī)電設(shè)備電氣控制線路排故的仿真系統(tǒng)設(shè)計
價值工程(2016年29期)2016-11-14 02:01:16
航空電氣系統(tǒng)中故障電弧研究
Buck開關(guān)變換器的基本參數(shù)設(shè)計及仿真分析
試析PLC控制下的自動化立體倉庫仿真情況分析
基于MADYMO的航空座椅約束系統(tǒng)優(yōu)化設(shè)計
科技視界(2016年18期)2016-11-03 21:44:44
中國體態(tài)假人模型與FAA Hybrid Ⅲ 型假人模型沖擊差異性分析
科技視界(2016年18期)2016-11-03 20:31:49
主站蜘蛛池模板: 毛片一区二区在线看| 天天视频在线91频| 国产色爱av资源综合区| 精品国产电影久久九九| 91亚洲视频下载| 亚洲精品国偷自产在线91正片| 中文纯内无码H| 国产一在线观看| 国产www网站| 亚洲国产精品日韩专区AV| 午夜精品一区二区蜜桃| 99在线观看精品视频| 久久综合亚洲鲁鲁九月天| 男人天堂亚洲天堂| 久无码久无码av无码| 拍国产真实乱人偷精品| 亚洲天堂视频网站| 亚洲三级网站| 东京热一区二区三区无码视频| 亚洲无码91视频| 亚洲色无码专线精品观看| 精品国产91爱| 成人中文在线| 亚洲精品麻豆| 亚洲成人一区二区三区| 女人一级毛片| 影音先锋丝袜制服| 超碰91免费人妻| 一本大道视频精品人妻 | 婷婷激情亚洲| 日韩av资源在线| 午夜视频在线观看区二区| 成年人国产网站| 一级毛片基地| 99免费视频观看| 国产男人的天堂| 一区二区影院| 无码福利视频| 国产大片黄在线观看| 97国产在线观看| 中文成人在线视频| 色婷婷亚洲十月十月色天| 少妇精品网站| 漂亮人妻被中出中文字幕久久 | 69综合网| 国产真实乱子伦视频播放| 91啦中文字幕| 国产一级在线观看www色| 国产va在线观看| 97色伦色在线综合视频| 国产真实乱了在线播放| 亚洲欧洲自拍拍偷午夜色无码| 国产精品免费露脸视频| 欧美一区二区三区欧美日韩亚洲| 久久午夜夜伦鲁鲁片无码免费| 国产精品久久久久久影院| 国产成人精品优优av| 狼友视频一区二区三区| 中文字幕乱妇无码AV在线| 亚洲经典在线中文字幕| 免费高清毛片| 美女免费黄网站| 8090午夜无码专区| 久久久精品国产亚洲AV日韩| 亚洲美女久久| 亚洲成人黄色网址| 黄片一区二区三区| 欧美综合在线观看| 日韩不卡高清视频| 日韩欧美综合在线制服| 香蕉国产精品视频| 广东一级毛片| 国产色婷婷| 精品一区二区三区波多野结衣| 91丨九色丨首页在线播放| 久久动漫精品| 精品视频一区在线观看| 韩日无码在线不卡| 国产91久久久久久| 日韩在线播放中文字幕| 女人18毛片久久| 国产一区在线观看无码|