• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    高黏度樹脂在纖維織物中浸漬過程的CFD模擬

    2016-11-21 02:48:21周云飛沈春銀戴干策
    關鍵詞:氣泡滲透率織物

    方 荀, 周云飛, 沈春銀, 戴干策

    (華東理工大學化學工程聯(lián)合國家重點實驗室,上海 200237)

    ?

    高黏度樹脂在纖維織物中浸漬過程的CFD模擬

    方 荀, 周云飛, 沈春銀, 戴干策

    (華東理工大學化學工程聯(lián)合國家重點實驗室,上海 200237)

    采用VOF(Volume of Fluid)方法結合多孔介質(zhì)模型對樹脂在纖維織物中的流動前沿進行追蹤,并對截留在纖維織物中的氣泡進行捕獲。纖維織物中的氣泡主要由結構不均勻性造成,結構的不均勻性導致流體速度差異并形成多個流動前沿,氣泡在動壓的作用下會發(fā)生遷移。研究了壓力、黏度以及雙尺度孔隙率對浸漬速率的影響,發(fā)現(xiàn)浸漬時間和壓力成反比,和黏度以及纖維束孔隙率與橫向滲透率的比值成正比;纖維織物的浸漬速率主要取決于微觀浸漬速率,而與纖維束的間距無關。

    數(shù)值模擬; 氣泡; VOF; 多孔介質(zhì); 浸漬時間

    在纖維增強復合材料中,氣泡的存在對復合材料的力學性能有著很大影響,當氣泡體積分數(shù)小于4%時,氣泡體積分數(shù)每增加1%,復合材料的層間剪切強度降低7%[1],因此減少甚至消除復合材料中的氣泡是復合材料成型的重要目的。排除氣泡的過程又被稱為浸漬。氣泡產(chǎn)生的因素如下:一種為揮發(fā)物所導致的氣泡[2],主要由環(huán)境濕度引起,研究認為高成型溫度下,吸濕對熱塑性復合材料的氣泡形成有很大影響。另一種為纖維氈的不均勻性所導致的氣泡形成[3],對于纖維織物,存在兩個尺度的流動,即纖維束內(nèi)和纖維束間的流動,其流速的差異會導致氣泡裹挾。

    對熱固性復合材料的浸漬過程,有大量學者進行了相關研究,基本是針對RTM(樹脂傳遞模塑)等成型工藝。Chen等[4]基于實驗的基礎提出了雙尺度孔隙率的空氣截留模型,在研究樹脂通過兩束纖維的流動時,將束間通道看成收縮-擴張的喉管,而纖維束表面則可以看作是多孔壁面,進而求解,對纖維束內(nèi)的流動則通過滲透速率來計算,并將氣泡壓縮性考慮其中。Patel等[5]提出了雙通道多孔介質(zhì)模型,并用兩相Darcy定律來求解。Chang等[6]構建了2-D模型研究單向排列的纖維氈的橫向流動,他們考慮了兩種形式的流動,繞著纖維束的宏觀流動和繞著纖維的微觀流動,使用單相Darcy定律進行計算,流體前沿壓力可以為真空、一個大氣壓或者毛細壓力,當宏觀流動的滲透率和微觀流動的滲透率之比超過20,則會在纖維束內(nèi)形成橢圓形的氣泡。Kang等[7]提出了一種數(shù)學模型,研究RTM成型過程中的氣泡形成,認識到不均勻的速度場是由纖維預制體結構的不均勻性所導致的,其中速度和毛細壓力的作用可用毛細數(shù)來表示,通過校正,該模型可以預測孔隙含量和尺寸。

    隨著計算機技術的發(fā)展,使用數(shù)值模擬的方法研究熱固性復合材料成型過程中氣泡的形成和排除越來越普遍。其中有限元(FE)/控制體積(CV)用來解決充填問題已經(jīng)非常通用和高效,求解區(qū)域使用固定的有限元網(wǎng)格,控制體積由網(wǎng)格節(jié)點決定,每個控制體積包含一個與之相關的填充因子。Simacek等[8]使用該方法求解了雙尺度多孔介質(zhì)的流動問題,其中大孔道可以使用一般的Darcy定律和質(zhì)量守恒來求解,而小孔道的求解則使用匯項。Lawrence等[9]在此基礎上添加了毛細作用力的影響,Tan等[10-11]進一步研究非等溫過程和反應過程。Lee等[12]則提出了一種真實時間模擬的方法來對空氣泡的形成和運動進行控制,使用光學傳感器檢測流動前沿,實時計算流動前沿的可能發(fā)展,并根據(jù)發(fā)展調(diào)整壓力,從而降低氣泡的含量。VOF方法同樣有著比較廣泛的應用,常用來模擬模腔中的樹脂流動前沿。Yang等[13]針對單向排列的纖維預制體,在纖維束內(nèi)使用Brinkman方程,纖維束間使用Stokes方程求解,使用VOF方法追蹤氣液界面。Devalve等[14]使用VOF方法結合多孔介質(zhì)模型模擬了纖維織物中氣體含量的變化,得到最小氣含量的工藝條件。Lundstr?m等[15]則使用Level set方法對無褶皺編織纖維中的樹脂流動前沿進行模擬。

    熱塑性復合材料的浸漬研究也有不少,但大多數(shù)是模型計算,文獻[16-19]提出了關于粉體浸漬或者混合紗的浸漬模型。Steggall等[20]提出了纖維織物熔體浸漬模型,Jespersen等[21]提出了薄膜疊層浸漬模型,考慮了纖維的運動。在熱塑性復合材料雙鋼帶壓機連續(xù)成型過程中,一定壓力下樹脂熔體在纖維織物的厚度方向上滲透,在這一過程中,熱塑性樹脂在纖維束內(nèi)和纖維束間的流動速度不同,使氣泡截留。氣泡的形成和遷移通過模型的計算不能直觀顯示出來,通過CFD數(shù)值計算不僅可以計算浸漬時間,同時可以觀察氣泡的形成和遷移。本文采用VOF方法對氣液兩相界面進行追蹤,對宏觀的壓力、熔體的物性以及纖維織物的結構等因素的影響進行探索。

    1 模擬對象及數(shù)學模型

    1.1 物理模型

    本文主要基于無褶皺雙軸向纖維織物的幾何結構,其代表單元中包含相互垂直的兩束纖維,且為多層結構。纖維束的邊緣假設為拋物線形,如圖1所示。其中纖維束內(nèi)包含成百上千根纖維單絲,將其看成具有一定孔隙率的多孔介質(zhì),在纖維軸向和橫向上具有不同的滲透率,該滲透率為孔隙率和纖維半徑的函數(shù)。樹脂在壓力的作用下從織物上方滲流,樹脂在該纖維織物代表單元中的流動前沿則使用VOF方法[22]進行追蹤。VOF是一種基于固定網(wǎng)格的方法,它定義了一個流體體積函數(shù),即網(wǎng)格單元中目標流體體積和網(wǎng)格體積之比(α),通過該流體體積函數(shù)來確定自由面,追蹤流體的變化。很顯然α的范圍可以從0到1變化,其中α=0表示該網(wǎng)格單元中全部為空氣,α=1表示網(wǎng)格單元中全部為流體,0<α<1時則為交界面單元。

    圖1 無褶皺纖維織物的幾何模型Fig.1 Geometric model of non-crimp fabrics

    1.2 纖維束間的流動

    纖維束間的通道比較大,可以看成宏觀通道,因此使用常規(guī)的層流VOF模型即可。流體體積分數(shù)守恒方程為

    (1)

    忽略重力的作用,總體質(zhì)量守恒方程和動量守恒方程為

    (2)

    (3)

    其中:ρ為流體和空氣平均密度,ρ=αρr+(1-α)ρa;μ為黏度,樹脂熔體假定為牛頓流體,μ=αμr+(1-α)μa。

    1.3 纖維束內(nèi)的流動

    由于纖維束是由許多根單向排列的纖維構成的,纖維之間構成了它的流動通道,該通道的尺寸比較小,一般為微米尺度,因此可看成是微觀通道。該微觀通道中的流動可以使用多孔介質(zhì)模型來計算。Darcy定律常常用來計算多孔介質(zhì)中的流體流動,將Darcy定律中的多孔介質(zhì)黏性阻力添加到動量源項,忽略慣性阻力,動量守恒方程可以寫成:

    (4)

    1.4 纖維束滲透率

    對于單向排列的纖維束,很早就有學者對滲透率進行了研究,Happel[23]取圓形的代表單元,將圓柱邊界指定為零剪切應力,獲得了纖維軸向和橫向流動的滲透率。后來的學者開始用非圓形的代表單元,如正方形排列或六邊形排列,求解Stokes方程,獲得一系列滲透率與纖維半徑a、孔隙率φ的關系。Gebart[24]得到:

    (5)

    (6)

    1.5 邊界條件

    樹脂熔體從上表面以一定的壓力流入纖維織物,入口設為壓力入口,壓力值為p0,纖維織物的下表面為出口,壓力始終保持一個大氣壓,表壓為0,因此出口邊界條件為壓力出口,壓力值為0。四周由于具有對稱性,因此周圍4個面為對稱平面,設置為對稱邊界條件。當t=0時,纖維織物的孔隙中無樹脂熔體,因此樹脂體積分數(shù)α=0,而入口處一直保持α=1。

    1.6 數(shù)值模擬環(huán)境

    采用商用軟件Fluent14.0,用有限體積法模擬,其中多相流模型選擇VOF兩相模型,纖維束區(qū)域設置為多孔介質(zhì)。對VOF瞬態(tài)計算,壓力-速度耦合方法采用PISO,壓力離散格式為PRESTO!,速度離散采用二階迎風格式,相體積分數(shù)離散采用Geo-Reconstruct方法,時間離散采用一階隱式。計算方程殘差設為10-3。時間步長采用可變時間步長,最大時間步長為1.0×10-3s,最小時間步長為1.0×10-9s,CourantNumber保持1.0。

    1.7 網(wǎng)格劃分及其無關化

    纖維束區(qū)域一共由6個體構成,其中包含2個1/2纖維束和4個1/4纖維束,其網(wǎng)格劃分方法為Cooper方法,即在纖維束的端面上劃分三角形網(wǎng)格,并沿著纖維束軸掃描,構成楔形網(wǎng)格。纖維束間的區(qū)域為了便于劃分Cooper網(wǎng)格,將其分割成3個體,其中2個沿著x方向掃描,1個沿著z方向掃描。劃分Cooper網(wǎng)格的優(yōu)點主要為網(wǎng)格生成比較規(guī)則,網(wǎng)格質(zhì)量較高。

    圖2 無褶皺纖維織物代表單元網(wǎng)格劃分Fig.2 Mesh of representative unit of non-crimp fabrics

    圖3 網(wǎng)格無關化檢驗Fig.3 Mesh independence

    2 結果與討論

    2.1 氣泡形成與運動

    本文采用的無褶皺雙軸向纖維織物代表單元長度L=3 mm,厚度h=1 mm,纖維束的寬度w=2.6 mm,纖維束重疊部分寬度w′=2 mm,代表單元上方有熔體在p0=1.0 MPa的壓力作用下向下流動,熔體假定為牛頓流體,黏度μ=2.5 Pa·s。纖維束的孔隙率φ=0.5,纖維半徑a=20 μm,此時可以通過關系式求得纖維束軸向的滲透率K∥=3.02×10-11m2,橫向的滲透率K⊥=6.54×10-12m2。由于纖維束間的通道比纖維束內(nèi)的大很多,因此樹脂熔體會優(yōu)先通過束間的通道,并同時緩慢地向纖維束內(nèi)滲透。正是由于兩者之間的速率差異顯著,纖維束間快速流動的樹脂會對纖維束內(nèi)的空氣進行圍追堵截,從而形成了氣泡,如圖4所示。實際上,纖維束內(nèi)部會產(chǎn)生毛細作用力,如果毛細作用力夠大的話,也可能導致纖維束內(nèi)的流動快于纖維束間的流動,從而使得氣泡被截留在纖維束間,本文中由于樹脂的黏度比較大,因此忽略毛細作用力的影響,氣泡一般形成于纖維束內(nèi)。

    圖4 熱塑性復合材料中氣泡的形成Fig.4 Bubble formation in thermoplastic composites

    液體浸漬纖維束間宏觀通道的速度非???在t=0.000 3 s時就已經(jīng)將上層的宏觀通道完全填充,而下層的宏觀通道浸漬速度并不快,其原因可以通過流線進行解釋。圖5所示為z=0平面上不同時間的相分布情況。圖6所示為從宏觀通道流入流體的流線,該宏觀通道類似于一個喉管,樹脂熔體收縮進入喉道,從喉道出來的流體只能向兩邊到達一定的距離,因此t=0.000 3 s時下層宏觀通道并不能完全填充。下層的宏觀通道如果需要繼續(xù)被填充,需要通過滲入纖維束中的流體來實現(xiàn),隨著滲入纖維束的熔體越來越多,下層宏觀通道能到達的區(qū)域則越遠,圖5中0.01~0.05 s反映了這一現(xiàn)象。微觀通道中的流體通過宏觀通道緩慢向其滲透,沿著纖維束周圍滲透的樹脂逐漸排除其中的空氣,因此宏觀通道也必須留出通道使空氣排出。

    圖7所示為z=1 mm平面上不同時間的相分布情況。當t=0.000 3 s時,樹脂已經(jīng)通過宏觀通道將中心處填充,而上方纖維束中的樹脂由于多孔介質(zhì)黏性阻力的作用滲透較慢;當t=0.002 s時,宏觀通道中的樹脂向周圍的纖維束橫向滲透,當兩股流體前沿相遇時,則裹挾了一段氣泡。接著氣泡受到流體的擠壓由扁平狀收縮為三角狀,并從纖維束內(nèi)向纖維束間遷移,接著被驅(qū)趕進入下層的纖維束內(nèi),同時又被擠壓成扁平狀,在t=0.05 s時,該氣泡被排除。兩側纖維束中的空氣則不斷被滲入的樹脂所排擠,從邊緣處被排出。

    圖5 z=0處樹脂相分布圖Fig.5 Resin phase distribution at z=0

    圖6 宏觀通道的流線和x-y平面投影圖Fig.6 Streamline of macro channels and projection at x-y plane

    圖7 z=1 mm處樹脂相分布圖Fig.7 Resin phase distribution at z=1 mm

    考慮到纖維床中的氣泡周圍被樹脂包圍,其外部和內(nèi)部的壓力是平衡的,氣泡的自由面具有相同的曲率。當樹脂受到驅(qū)動時,會產(chǎn)生額外的動壓和黏性應力場,并導致自由面曲率發(fā)生變化,下游方向的曲率增大,上游方向的曲率降低,從而使壓差適應外界的壓力梯度,然而自由面對壓力的調(diào)整是有限的,一旦外界的壓力梯度超過了自由面的可調(diào)整范圍,那么氣泡就開始運動。因此,當氣泡兩端的毛細壓差小于動壓,則氣泡就會移動,反之,氣泡則不會被排出。毛細壓差(Δpcap)為

    (7)

    其中σ為表面張力。而黏性流動導致的壓差(Δp)為

    (8)

    其中,θ為氣泡長度與來流方向的夾角。如果Δp/Δpcap>1,則氣泡可以移動。對復合材料中的氣泡遷移過程,Gangloff等[25]通過數(shù)值計算發(fā)現(xiàn),對于給定尺寸的氣泡,一定的操作條件可以通過毛細管數(shù)來決定氣泡能否遷移,很顯然和上面的分析相符合。在熱塑性復合材料中,由于黏性流動所導致的壓差一般大于毛細壓差,因此氣泡能夠發(fā)生遷移。

    圖8所示為x=0平面上不同時間的相分布圖。結合圖7和圖8可以看出,氣泡在纖維織物中的流動方向和樹脂的流動方向基本一致。氣泡剛從纖維束內(nèi)進入纖維束間時,氣泡在向下遷移的同時向兩側運動,當氣泡從纖維束間又進入下層的纖維束內(nèi)時,氣泡在向下遷移的同時又向中心運動,可以看出氣泡的遷移主要受流體黏性壓差的作用驅(qū)動。同時可以發(fā)現(xiàn)氣泡在遷移的過程中會發(fā)生破裂,特別是由纖維束間到纖維束內(nèi)時,大氣泡分裂成小氣泡再排出。

    圖8 x=0處樹脂相分布圖Fig.8 Resin phase distribution at x=0

    圖9所示為纖維束間、纖維束內(nèi)和整個纖維織物中的樹脂體積分數(shù)隨時間的變化。很顯然,纖維束間的浸漬時間很短,而纖維束內(nèi)的浸漬時間相對較長,而整個纖維織物的浸漬時間主要取決于纖維束的浸漬時間。纖維束間在0~0.000 3 s之間樹脂體積分數(shù)變化很快,而在0.000 3~0.05 s之間仍有一段緩慢增長,這是由于下層的宏觀通道并不能迅速填充完全,和圖4所顯示的情況相符。

    圖9 纖維束間、纖維束內(nèi)、預制體中 樹脂體積分數(shù)與時間的關系Fig.9 Volume fraction of resin in intra-tow and inter-tow and preform vs. time

    熱塑性復合材料中的氣泡形成主要是由于纖維織物中局部阻力不同導致流速差異,從而產(chǎn)生多個流動前沿,由于宏觀通道的流速大于微觀通道,因此流動前沿一般在纖維束內(nèi)相遇,從而截留氣泡。由于熱塑性復合材料中的黏性流動產(chǎn)生的壓差大于毛細壓差,因此氣泡能夠發(fā)生遷移,只要時間足夠長,氣泡就可以被排除,這意味著保壓時間對于氣泡的排除十分重要。需要注意的是,延長保壓時間雖然能夠?qū)馀莺拷档偷阶畹椭?但對工業(yè)上需要的高效率生產(chǎn)方式相違背,此時可能需要權衡制品質(zhì)量和保壓時間的關系。就本例而言,當樹脂體積分數(shù)達到0.90時需要的時間為0.03 s,達到0.95時所需要的時間延長50%,達到0.99時所需的時間延長100%,而達到0.995時所需的時間延長200%,由此可見浸漬過程越往后時間越長,如果制品的氣含量能夠允許0.5%,浸漬時間只需要0.09 s即可。

    2.2 壓力對浸漬速度的影響

    隨著壓力的增加,樹脂在纖維織物中的流速也增加,因此排除空氣的速度會加快。圖10示出了不同的壓力對樹脂充填速度的影響。壓力對宏觀通道的浸漬影響不太顯著,如圖10(b)所示,因為宏觀通道中的流速很快,在3×10-3s之前樹脂體積分數(shù)就超過0.9,最終的差別僅僅是下層的宏觀通道的浸漬速率。對比圖10(a)和圖10(c),壓力對微觀通道浸漬速率和纖維織物整體的浸漬速率基本一致,可見整個纖維織物的浸漬速率取決于纖維束的浸漬速率。纖維束內(nèi)的浸漬速率又基本取決于橫向的浸漬速率,這是因為橫向的滲透率要比軸向滲透率低。

    圖10 壓力對樹脂浸漬速度的影響Fig.10 Effects of pressure on impregnation velocity

    纖維束模型常常使用矩形或者圓形。假設纖維束為矩形,如圖11所示,纖維束橫向滲透速率和壓力的關系可用Darcy定律來表示。

    (9)

    圖11 矩形纖維束模型和圓形纖維束模型Fig.11 Rectangular and circle models of fiber tows

    因此,可求得浸漬時間分析解為

    (10)

    其中,b為樹脂滲入纖維束的距離。如果假設纖維束為圓形,流體從四周向中心滲透,Darcy定律寫成柱坐標形式為

    (11)

    其中,徑向流動的壓力梯度為

    (12)

    通過積分求得浸漬時間隨流動前沿變化的分析解為

    (13)

    對一般形狀的纖維束,可以得到推廣的浸漬時間:

    (14)

    對于給定的纖維織物結構,由于纖維束的幾何形狀、尺寸和樹脂流動方向基本不變,可認為f(φ)僅是總的樹脂體積分數(shù)φ的函數(shù)。圖10(d)中示出了不同樹脂填充體積分數(shù)下浸漬時間和壓力的關系,也可以發(fā)現(xiàn)其變化規(guī)律基本和上述的分析結果相符,即在相同的樹脂體積分數(shù)下,浸漬時間和壓力的乘積是定值。對該纖維織物而言,壓力增大到一定程度時,浸漬時間的變化則不太明顯,壓力從0.5 MPa增加到1.0 MPa,樹脂體積分數(shù)達到0.99時的浸漬時間降低了50%,而壓力從0.5 MPa增加到1.5 MPa時,浸漬時間降低66%,和50%相比變化不大,因此在壓力和浸漬時間中權衡的話,1.0 MPa的壓力比較合適,繼續(xù)增加壓力雖能降低浸漬時間,但意義不顯著。

    2.3 樹脂黏度對浸漬速率的影響

    樹脂黏度對浸漬時間的影響很大,黏度越大意味著黏性阻力越大,從之前的分析可知,浸漬時間主要取決于纖維束內(nèi)的浸漬時間,在本文中為纖維橫向的浸漬時間。不管是矩形還是圓形的纖維束模型,在獲得同樣的樹脂體積分數(shù)條件下,浸漬時間都與黏度成正比,如圖12所示。在本文的物理模型中,并沒有忽略慣性項的計算,而實際的結果表明,在該黏度范圍慣性項的作用可忽略,只考慮黏性項時,對宏觀通道的流動僅用Stokes方程,而對微觀的纖維束內(nèi)的流動則使用Darcy定律即可。此時,從方程上也可以得到在相同的壓力下,流速與黏度成反比,因而浸漬時間和黏度成正比。這一結果僅針對黏度較高的情形,對于熱固性樹脂,因為不可忽略慣性項,浸漬時間也就與黏度不成正比。當然,針對實際的聚合物熔體,簡單的牛頓黏性定律不能充分描述流體的性質(zhì),常常使用剪切稀化或者黏彈性本構方程描述,因此實際的情況將更加復雜。需要說明的是,熱塑性樹脂的真實黏度比本文設定的黏度還要高一到兩個數(shù)量級以上,本文設定的黏度雖遠高于一般熱固性樹脂的黏度,但仍存在缺憾。從圖12(b)看出,當樹脂體積分數(shù)達到0.99,黏度從1.5 Pa·s增加到4.5 Pa·s時,浸漬時間從0.04 s增加至0.12 s,可以預計到當黏度提高一到兩個數(shù)量級后,浸漬時間會增加幾十甚至上百倍,此時可以通過提高壓力的方法來縮短浸漬時間。

    圖12 樹脂黏度對浸漬速率的影響Fig.12 Effects of viscosity on impregnation velocity

    2.4 孔隙率對浸漬速率的影響

    在纖維織物中包含了宏觀通道和微觀通道兩種尺度的通道,因此它包含了兩種尺度的孔隙率:一種為大孔孔隙率,另一種為微孔孔隙率,它們都會對樹脂的浸漬以及氣泡的形成產(chǎn)生影響。本文對纖維束中的微孔采用多孔介質(zhì)方法,大孔采用真實的通道進行模擬。纖維束的孔隙率越大,橫向和軸向的滲透率越大,當纖維束完全緊密堆積時,橫向的滲透率則減為零,而軸向仍可以大于零。圖13所示為通過經(jīng)驗式得出的橫向與軸向滲透率隨孔隙率的變化關系。當纖維束孔隙率在0.4~0.7變化時,橫向滲透率與軸向滲透率之比在5左右。當孔隙率減小時,橫向滲透率會產(chǎn)生突降,而軸向滲透率降低相對較慢,兩者之間差別可達多個數(shù)量級。如果纖維束排列特別緊密,橫向滲透將非常困難,尤其本文中樹脂的主流動方向即為纖維束橫向流動,纖維束排列不能過于緊密,否則浸漬時間會增加多個數(shù)量級。

    圖 13 纖維束軸向和橫向滲透率與孔隙率的關系Fig.13 Relationship between axial,transverse permeability and porosity

    圖14所示為纖維束孔隙率對浸漬速率的影響。很顯然,隨著孔隙率的增大,纖維束的滲透率增大,黏性阻力降低,有利于微觀浸漬。由于纖維織物的浸漬速率主要取決于纖維束的浸漬速率,因此孔隙率增大,微觀浸漬速率加快,同時整個纖維織物的浸漬速率也隨之加快。為了更好地預測纖維束孔隙率對浸漬時間的影響,結合式(14)將浸漬時間和孔隙率對橫向滲透率的比值(φ/K⊥)作圖,可以得出浸漬時間和孔隙率對橫向滲透率的比值基本成正比。因此,幾何結構相同的條件下,可以估算出不同的纖維束孔隙率浸漬所需要的時間。當纖維束孔隙率從0.7降低至0.4時,相同樹脂體積分數(shù)下的浸漬時間提高了10倍左右。如果想要獲得高性能的復合材料,必須要提高纖維含量,這也意味著纖維束的孔隙率降低,從而浸漬速率降低,加工難度提高,此時需要通過提高溫度、降低黏度或者提高壓力的方法來提高浸漬速率。

    圖14 纖維束孔隙率對浸漬速率的影響Fig.14 Effects of porosity of fiber tows on impregnation velocity

    除了微孔孔隙率,纖維織物還存在纖維束間的大孔孔隙率。在纖維束尺寸不變的情況下,調(diào)整纖維束間的距離,大孔的孔隙率將發(fā)生變化。圖15顯示了不同的纖維束間距離對樹脂浸漬速度的影響。很顯然,大孔的孔隙率對樹脂浸漬速率影響很小,尤其間距d從0.4 mm到0.6 mm變化時,幾乎不變,這是因為當d=0.2 mm時,宏觀通道的尺寸相對較小,流動阻力較大,因此宏觀通道的浸漬速率仍然會對整個纖維織物的浸漬速率產(chǎn)生微量的影響,而隨著纖維束的間距的增大,宏觀通道的浸漬速率變快,整個纖維織物的浸漬速率幾乎只取決于微觀浸漬速率,因此它們的浸漬速率曲線幾乎重合。但總體來說,纖維織物的浸漬速率還是取決于微觀浸漬速率,與之前的分析一致。

    圖15 纖維束間距對浸漬速率的影響Fig.15 Effects of tow gaps on impregnation velocity

    此外,隨著纖維束間距的減小,受到截留的氣泡運動速度反而加快,圖16對比了3種纖維束間距的受裹挾氣泡的遷移情況。在t=0.03 s,間距d=0.2 mm時,纖維束間和纖維束內(nèi)被截留的氣泡幾乎都被排除;當d=0.4 mm時,被截留的氣泡仍然存在,但已經(jīng)向下遷移到纖維束中,并處于排除狀態(tài);d=0.6 mm時,被截留的氣泡還停在纖維束間的通道中。可以看出d=0.2 mm截留氣泡遷移的速度大于d=0.4 mm,氣泡遷移速度最慢的是d=0.6 mm。在宏觀通道中,雖然垂直通道中的速度規(guī)律是隨著纖維束間距的減小而降低,但是在水平通道中,速度卻隨著纖維束間距的減小而有所上升。

    圖16 t=0.03 s時的相分布圖 (z=0)Fig.16 Phase distribution at t=0.03 s (z=0)

    圖17的速度分布圖可以說明這一情況。從總體速度云圖來看,顯然間距越大,在垂直通道中的流速越快,從局部速度云圖可以看出,間距越小在水平通道中的流速越快,因此流體驅(qū)趕氣泡的速度也會變大,氣泡反而更容易排出。

    圖17 總體速度和局部速度(0~0.2 m/s)云圖(z=0)Fig.17 Global and local (0~0.2 m/s) velocity contour (z=0)

    2.5 模型驗證

    本文主要采用VOF結合多孔介質(zhì)模型的方法對纖維多孔介質(zhì)中的氣-液兩相流動進行模擬,為了驗證CFD模型的正確性,對多孔介質(zhì)方腔,分別采用數(shù)值計算和模型公式求解時間和填充分數(shù)的關系,如圖18所示。其中多孔介質(zhì)滲透率為1.0×10-9,孔隙率為0.5,樹脂黏度為5 Pa·s,方腔長度為10 mm,寬度為5 mm,樹脂沿著方腔長度方向流動,壓力為105Pa。和式(10)類似,推導求得時間和填充分數(shù)的關系為

    (15)

    對比用式(15)和CFD方法得到的時間與樹脂填充分數(shù)的關系,可以發(fā)現(xiàn),采用CFD方法求解的結果和模型公式所得結果吻合很好,可以證明采用VOF追蹤氣-液界面的方法在多孔介質(zhì)中是適用的。

    圖18 模型法和CFD結果對比Fig.18 Comparison between CFD and model results

    Wang等[26]對雙鋼帶壓機的運行速度和浸漬效果進行了定量測試,如果將運行速度轉(zhuǎn)化為浸漬時間,即運行速度正比于浸漬時間的倒數(shù),Wang得到的孔隙率隨浸漬時間的變化趨勢和本文模擬得到的變化趨勢基本一致。

    3 結 論

    本文采用數(shù)值模擬方法對無褶皺雙軸纖維織物代表單元中的樹脂浸漬過程進行模擬,其中VOF方法能夠追蹤氣液兩相的界面,而纖維束內(nèi)則采用多孔介質(zhì)模型進行模擬,纖維束間為一般流動通道。采用CFD方法能夠?qū)椢锏慕n進行模擬分析,得到的規(guī)律能夠?qū)Τ尚推鸬揭欢ǖ闹笇ё饔谩?/p>

    (1) 纖維束間的宏觀浸漬速率比微觀浸漬速率大兩個數(shù)量級,纖維織物的浸漬速率主要取決于纖維束內(nèi)的浸漬速率。

    (2) 樹脂對纖維織物的浸漬速率開始時很快,然后不斷降低,達到100%浸漬程度的時間相對較長,可以根據(jù)氣泡含量要求選擇浸漬時間,過長的時間會降低生產(chǎn)效率。

    (3) 壓力越大,浸漬速率越大,獲得同樣的樹脂體積分數(shù)所需時間和壓力成反比。壓力的選擇也需要權衡考慮,隨著壓力的不斷增大,它對浸漬速率的影響越來越小,然而它對設備的要求越來越高。

    (4) 黏度越大,所需的浸漬時間越長,浸漬所需時間基本和黏度成正比,對于更高黏度的樹脂熔體,浸漬時間可能達到幾十、上百秒,需要增加浸漬壓力來減少浸漬時間??紫堵蕪?.7降低至0.4,浸漬時間提高10倍,可考慮適當提高成型溫度使黏度降低,或者增加壓力來提高浸漬速率。

    符號說明:

    a——纖維半徑,μm

    b——樹脂滲入纖維束的距離,mm

    d——纖維束間距,mm

    h——纖維束厚度,mm

    K——滲透率張量,m2

    K⊥——橫向滲透率,m2

    K∥——軸向滲透率,m2

    L——代表單元長度,mm

    l——氣泡的長度,mm

    p——壓力,MPa

    p0——入口壓力,MPa

    pcap——毛細作用力,MPa

    r——纖維束浸漬前沿的半徑,mm

    r0——纖維束的半徑,mm

    rb——氣泡半徑,mm

    rt——孔隙收縮處半徑,mm

    t——浸漬時間,s

    tcir——圓形纖維束浸漬時間,s

    trec——矩形纖維束浸漬時間,s

    u——流動速度,m/s

    ub——矩形纖維束模型的滲透速率,m/s

    ur——圓形纖維束模型的滲透速率,m/s

    w ——纖維束寬度,mm

    w′——纖維束重疊部分寬度,mm

    α——樹脂相體積分數(shù)

    ρ ——平均密度,kg/m3

    ρa——空氣密度,kg/m3

    ρr——樹脂密度,kg/m3

    μ ——平均黏度,Pa·s

    μa——空氣黏度,Pa·s

    μr——樹脂黏度,Pa·s

    φ ——纖維束孔隙率

    φc——緊密堆積時的孔隙率

    φ——總的樹脂體積分數(shù)

    [1] JUDD N,WRIGHT W.Voids and their effects on the mechanical properties of composites:An appraisal[J].SAMPE Journal,1978,14:10-14.

    [2] ROYCHOWDHURY S,GILLESPIE Jr J,ADVANI S.Computer Aided Design in Composite Material Technology III[M].Netherlands:Springer,1992:89-107.

    [3] YU M C,MIDDLEMAN S.Air entrapment during liquid infiltration of porous media[J].Chemical Engineering Communications,1993,123(1):61-69.

    [4] CHEN Y T,MACOSKO C W,DAVIS H T.Wetting of fiber mats for composites manufacturing:II.Air entrapment model[J].AIChE Journal,1995,41(10):2274-2281.

    [5] PATEL N,LEE L J.Modeling of void formation and removal in liquid composite molding:Part II.Model development and implementation[J].Polymer Composites,1996,17(1):104-114.

    [6] CHANG C Y, HOURNG L W.Study on void formation in resin transfer molding[J].Polymer Engineering & Science,1998,38(5):809-818.

    [7] KANG M K,LEE W I,HAHN H T.Formation of microvoids during resin-transfer molding process[J].Composites Science and Technology,2000,60(12):2427-2434.

    [8] SIMACEK P, ADVANI S G.A numerical model to predict fiber tow saturation during liquid composite molding[J].Composites Science and Technology,2003,63(12):1725-1736.

    [9] LAWRENCE J M,NEACSU V,ADVANI S G.Modeling the impact of capillary pressure and air entrapment on fiber tow saturation during resin infusion in LCM[J].Composites:Part A.Applied Science and Manufacturing,2009,40(8):1053-1064.

    [10] TAN H,PILLAI K M.Multiscale modeling of unsaturated flow of dual-scale fiber preform in liquid composite molding:II.Non-isothermal flows[J].Composites:Part A.Applied Science and Manufacturing,2012,43(1):14-28.

    [11] TAN H,PILLAI K M.Multiscale modeling of unsaturated flow in dual-scale fiber preforms of liquid composite molding:III.Reactive flows[J].Composites:Part A.Applied Science and Manufacturing,2012,43(1):29-44.

    [12] LEE D H,LEE W I,KANG M K.Analysis and minimization of void formation during resin transfer molding process[J].Composites Science and Technology,2006,66(16):3281-3289.

    [13] YANG J,JIA Y,SUN S,etal.Mesoscopic simulation of the impregnating process of unidirectional fibrous preform in resin transfer molding[J].Materials Science and Engineering:A,2006,435:515-520.

    [14] DEVALVE C,PITCHUMANI R.Simulation of void formation in liquid composite molding processes[J].Composites:Part A.Applied Science and Manufacturing,2013,51:22-32.

    [15] LUNDSTROM T S,FRISHFELDS V,JAKOVICS A.Bubble formation and motion in non-crimp fabrics with perturbed bundle geometry[J].Composites:Part A.Applied Science and Manufacturing,2010,41(1):83-92.[16] LEE W I,SPRINGER G S.A model of the manufacturing process of thermoplastic matrix composites[J].Journal of Composite Materials,1987,21(11):1017-1055.

    [17] BERNET N,MICHAUD V,BOURBAN P,etal.An impregnation model for the consolidation of thermoplastic composites made from commingled yarns[J].Journal of Composite Materials,1999,33(8):751-772.

    [18] SALA G,CUTOLO D.The pultrusion of powder-impregnated thermoplastic composites[J].Composites:Part A.Applied Science and Manufacturing,1997,28(7):637-646.

    [19] VAN West B,PIPES R B,ADVANI S.The consolidation of commingled thermoplastic fabrics[J].Polymer Composites,1991,12(6):417-427.

    [20] STEGGALL-MURPHY C,SIMACEK P,ADVANI S G,etal.A model for thermoplastic melt impregnation of fiber bundles during consolidation of powder-impregnated continuous fiber composites[J].Composites:Part A.Applied Science and Manufacturing,2010,41(1):93-100.

    [21] JESPERSEN S,WAKEMAN M,MICHAUD V,etal.Film stacking impregnation model for a novel net shape thermoplastic composite preforming process[J].Composites Science and Technology,2008,68(7):1822-1830.

    [22] HIRT C W,NICHOLS B D.Volume of fluid (VOF) method for the dynamics of free boundaries[J].Journal of Computational Physics,1981,39(1):201-225.

    [23] HAPPEL J.Viscous flow relative to arrays of cylinders[J].AIChE Journal,1959,5(2):174-177.

    [24] GEBART B.Permeability of unidirectional reinforcements for RTM[J].Journal of Composite Materials,1992,26(8):1100-1133.

    [25] GANGLOFF J J,DANIEL C,ADVANI S G.A model of two-phase resin and void flow during composites processing[J].International Journal of Multiphase Flow,2014,65:51-60.

    [26] WANG X,MAYER C,NEITZEL M.Some issues on impregnation in manufacturing of thermoplastic composites by using a double belt press[J].Polymer Composites,1997,18(6):701-710.

    CFD Simulation of Impregnation Process with Highly Viscous Resin in Fiber Fabrics

    FANG Xun, ZHOU Yun-fei, SHEN Chun-yin, DAI Gan-ce

    (State Key Laboratory of Chemical Engineering,East China University of Science and Technology,Shanghai 200237,China)

    Flow front in fiber fabrics was tracked by VOF (volume of fluid) method combined with a porous medium model,and the entrapped air bubbles were captured by simulation.The formation of bubbles is due to the heterogeneous structure which causes the velocity difference with several flow front.The bubbles will migrate under the dynamic pressure.Moreover,the effects of pressure,viscosity and dual-scale porosity on impregnation velocity were studied,and the impregnation time was in propotional to the reciprocal of pressure,viscosity,or transverse permeability divided by tow porosity.The impregnation velocity was determined by the impregnation of micro channels while had nothing to do with the distance between two fiber tows.

    numerical simulation; bubbles; VOF; porous medium; impregnation time

    1006-3080(2016)05-0615-10

    10.14135/j.cnki.1006-3080.2016.05.005

    2016-05-04

    國家自然科學基金(21376086)

    方 荀(1989-),男,江蘇泰州人,博士生,研究方向為熱塑性復合材料加工成型。E-mail:shawnfang@msn.cn

    沈春銀,E-mail:ichuny@ecust.edu.cn

    TB332; TQ021.1

    A

    猜你喜歡
    氣泡滲透率織物
    檸檬氣泡水
    欣漾(2024年2期)2024-04-27 15:19:49
    無Sn-Pd活化法制備PANI/Cu導電織物
    《紡織品織物折痕回復角的測定》正式發(fā)布
    SIAU詩杭便攜式氣泡水杯
    新潮電子(2021年7期)2021-08-14 15:53:12
    浮法玻璃氣泡的預防和控制對策
    竹纖維織物抗菌研究進展
    冰凍氣泡
    中煤階煤層氣井排采階段劃分及滲透率變化
    不同滲透率巖芯孔徑分布與可動流體研究
    SAGD井微壓裂儲層滲透率變化規(guī)律研究
    国内少妇人妻偷人精品xxx网站 | 免费人成视频x8x8入口观看| 操出白浆在线播放| 久久久久久久久免费视频了| 曰老女人黄片| 亚洲美女黄片视频| 国产精品自产拍在线观看55亚洲| 一进一出抽搐动态| 欧美大码av| 久久久久免费精品人妻一区二区| 999久久久国产精品视频| 国产精品一及| 男女视频在线观看网站免费 | 国产不卡一卡二| 亚洲狠狠婷婷综合久久图片| 免费在线观看成人毛片| 欧美最黄视频在线播放免费| 一二三四在线观看免费中文在| 国产av在哪里看| 哪里可以看免费的av片| 制服诱惑二区| 欧美久久黑人一区二区| 亚洲国产高清在线一区二区三| 99re在线观看精品视频| 88av欧美| 88av欧美| 日本 av在线| 人人妻,人人澡人人爽秒播| 老汉色∧v一级毛片| 三级国产精品欧美在线观看 | 一进一出抽搐gif免费好疼| 国产精品一区二区精品视频观看| 欧美国产日韩亚洲一区| 成人亚洲精品av一区二区| 国产69精品久久久久777片 | 久久久久国产一级毛片高清牌| 97超级碰碰碰精品色视频在线观看| 很黄的视频免费| bbb黄色大片| 国产欧美日韩一区二区精品| 午夜福利视频1000在线观看| 欧美在线一区亚洲| 日韩有码中文字幕| 69av精品久久久久久| 极品教师在线免费播放| 欧美黑人巨大hd| 一级作爱视频免费观看| 欧美成人免费av一区二区三区| 一本精品99久久精品77| 免费观看精品视频网站| 一级毛片高清免费大全| 99久久99久久久精品蜜桃| 国产伦在线观看视频一区| 一级黄色大片毛片| 成人特级黄色片久久久久久久| ponron亚洲| 免费在线观看影片大全网站| 99久久久亚洲精品蜜臀av| 午夜免费成人在线视频| 欧美+亚洲+日韩+国产| 一本大道久久a久久精品| 午夜日韩欧美国产| 国内精品一区二区在线观看| 两性午夜刺激爽爽歪歪视频在线观看 | 18禁黄网站禁片午夜丰满| 亚洲五月婷婷丁香| 制服丝袜大香蕉在线| 波多野结衣高清作品| 国产91精品成人一区二区三区| 美女午夜性视频免费| 亚洲国产日韩欧美精品在线观看 | 精品乱码久久久久久99久播| 青草久久国产| 精品久久久久久久久久久久久| 国内精品一区二区在线观看| 国产私拍福利视频在线观看| 91国产中文字幕| 黑人欧美特级aaaaaa片| 精品一区二区三区视频在线观看免费| 精品电影一区二区在线| 国产亚洲精品综合一区在线观看 | 国产亚洲精品一区二区www| 成人18禁高潮啪啪吃奶动态图| 黑人欧美特级aaaaaa片| 97超级碰碰碰精品色视频在线观看| 欧美大码av| 女同久久另类99精品国产91| 日韩欧美三级三区| 国产av又大| 99热这里只有是精品50| 91字幕亚洲| 九色成人免费人妻av| 国产精华一区二区三区| 成人精品一区二区免费| 日韩成人在线观看一区二区三区| 亚洲狠狠婷婷综合久久图片| 很黄的视频免费| 搡老熟女国产l中国老女人| 男人舔奶头视频| 女生性感内裤真人,穿戴方法视频| 99精品久久久久人妻精品| 成人一区二区视频在线观看| 久久欧美精品欧美久久欧美| 午夜亚洲福利在线播放| 人妻久久中文字幕网| 色老头精品视频在线观看| 18禁黄网站禁片免费观看直播| 男女之事视频高清在线观看| 国产高清videossex| 女同久久另类99精品国产91| 五月玫瑰六月丁香| 操出白浆在线播放| 一本大道久久a久久精品| 久久国产乱子伦精品免费另类| 91成年电影在线观看| 国产成人av教育| 国产激情偷乱视频一区二区| 女人高潮潮喷娇喘18禁视频| 老鸭窝网址在线观看| 在线播放国产精品三级| 老汉色∧v一级毛片| 中文字幕人妻丝袜一区二区| 午夜精品久久久久久毛片777| 夜夜躁狠狠躁天天躁| 亚洲熟女毛片儿| 美女黄网站色视频| 在线观看免费午夜福利视频| 女人爽到高潮嗷嗷叫在线视频| 亚洲自拍偷在线| 51午夜福利影视在线观看| 精品久久蜜臀av无| 又紧又爽又黄一区二区| 亚洲男人天堂网一区| tocl精华| 国产亚洲精品一区二区www| 在线视频色国产色| 成人手机av| 国产探花在线观看一区二区| 色在线成人网| 老鸭窝网址在线观看| 中文字幕最新亚洲高清| 午夜福利高清视频| 少妇裸体淫交视频免费看高清 | 国产一级毛片七仙女欲春2| 国产欧美日韩精品亚洲av| 一级作爱视频免费观看| 香蕉久久夜色| 久久精品夜夜夜夜夜久久蜜豆 | netflix在线观看网站| 亚洲自拍偷在线| 免费搜索国产男女视频| 亚洲第一电影网av| 日韩欧美 国产精品| 亚洲无线在线观看| 一卡2卡三卡四卡精品乱码亚洲| 国产亚洲av高清不卡| 亚洲一区高清亚洲精品| 两性午夜刺激爽爽歪歪视频在线观看 | a在线观看视频网站| 日本在线视频免费播放| 国产爱豆传媒在线观看 | 日韩精品中文字幕看吧| 波多野结衣高清作品| 别揉我奶头~嗯~啊~动态视频| 麻豆一二三区av精品| 亚洲精品在线观看二区| 日韩大码丰满熟妇| 一级毛片精品| 国产在线观看jvid| 午夜福利视频1000在线观看| 99久久99久久久精品蜜桃| av天堂在线播放| 在线观看美女被高潮喷水网站 | 麻豆久久精品国产亚洲av| 亚洲无线在线观看| www.熟女人妻精品国产| 日韩欧美 国产精品| av在线天堂中文字幕| 亚洲精品国产一区二区精华液| 国产亚洲av高清不卡| 国产av不卡久久| 美女 人体艺术 gogo| 免费在线观看视频国产中文字幕亚洲| 久久天堂一区二区三区四区| 久久久精品国产亚洲av高清涩受| 两性午夜刺激爽爽歪歪视频在线观看 | 国产日本99.免费观看| 超碰成人久久| 国产成年人精品一区二区| 99热6这里只有精品| 搞女人的毛片| 床上黄色一级片| 91国产中文字幕| 可以免费在线观看a视频的电影网站| tocl精华| 久久久久国产一级毛片高清牌| 久久精品国产99精品国产亚洲性色| 无限看片的www在线观看| 一进一出好大好爽视频| 三级毛片av免费| 亚洲va日本ⅴa欧美va伊人久久| 在线永久观看黄色视频| 欧美中文日本在线观看视频| 一二三四在线观看免费中文在| 国产午夜精品久久久久久| 性欧美人与动物交配| 三级男女做爰猛烈吃奶摸视频| 国产亚洲欧美在线一区二区| 色av中文字幕| 国产一级毛片七仙女欲春2| av在线天堂中文字幕| 精品久久久久久久久久久久久| 国产三级中文精品| a在线观看视频网站| 99久久99久久久精品蜜桃| 久久精品国产亚洲av香蕉五月| 久久久国产成人免费| 少妇熟女aⅴ在线视频| 99久久精品热视频| e午夜精品久久久久久久| 国产精品亚洲一级av第二区| 最近视频中文字幕2019在线8| 欧美中文综合在线视频| 国产亚洲精品第一综合不卡| 亚洲色图 男人天堂 中文字幕| 欧美乱码精品一区二区三区| 亚洲18禁久久av| 少妇的丰满在线观看| 中文亚洲av片在线观看爽| 淫秽高清视频在线观看| 亚洲男人的天堂狠狠| 成人av一区二区三区在线看| 一边摸一边抽搐一进一小说| 欧美性长视频在线观看| 国产一区二区激情短视频| 色av中文字幕| 国产亚洲精品第一综合不卡| 亚洲成av人片免费观看| 成年版毛片免费区| 国产三级黄色录像| 亚洲精品国产一区二区精华液| 999精品在线视频| 搡老妇女老女人老熟妇| 国产1区2区3区精品| 国产精品香港三级国产av潘金莲| 日本一区二区免费在线视频| 变态另类成人亚洲欧美熟女| 少妇裸体淫交视频免费看高清 | 国产不卡一卡二| 在线观看日韩欧美| 久久国产精品影院| 国产成人精品久久二区二区91| 国产成人av激情在线播放| 国产精品影院久久| 欧美日韩亚洲综合一区二区三区_| 99精品久久久久人妻精品| 精品一区二区三区av网在线观看| 国产熟女xx| 成人三级黄色视频| bbb黄色大片| 黑人巨大精品欧美一区二区mp4| 午夜福利免费观看在线| 淫妇啪啪啪对白视频| 国产精品久久久av美女十八| 女警被强在线播放| 我的老师免费观看完整版| 亚洲欧美一区二区三区黑人| 看黄色毛片网站| 欧美午夜高清在线| 成年免费大片在线观看| 亚洲中文字幕一区二区三区有码在线看 | 色播亚洲综合网| 一级毛片高清免费大全| 天堂√8在线中文| 亚洲avbb在线观看| 伊人久久大香线蕉亚洲五| 色精品久久人妻99蜜桃| 亚洲午夜精品一区,二区,三区| 国产三级黄色录像| 看免费av毛片| 国产精品,欧美在线| 久久久国产成人精品二区| 一区二区三区高清视频在线| 国产成人影院久久av| 欧美性长视频在线观看| 国产精品av久久久久免费| 精品一区二区三区av网在线观看| 一a级毛片在线观看| 久久久久性生活片| 嫩草影视91久久| 高清毛片免费观看视频网站| 欧美日韩精品网址| 香蕉丝袜av| 久久久精品大字幕| 国产精品久久久人人做人人爽| 国产视频内射| 老熟妇仑乱视频hdxx| 99国产精品一区二区蜜桃av| 国产亚洲精品av在线| 99精品欧美一区二区三区四区| 老司机午夜十八禁免费视频| 欧美一区二区精品小视频在线| 欧美成人一区二区免费高清观看 | 波多野结衣高清无吗| 久久天躁狠狠躁夜夜2o2o| 午夜日韩欧美国产| 欧美黑人欧美精品刺激| 国产av一区在线观看免费| 色在线成人网| x7x7x7水蜜桃| 亚洲av成人不卡在线观看播放网| 两个人视频免费观看高清| 麻豆一二三区av精品| 18禁国产床啪视频网站| 亚洲男人的天堂狠狠| 亚洲真实伦在线观看| 给我免费播放毛片高清在线观看| 老司机在亚洲福利影院| 岛国在线免费视频观看| 特级一级黄色大片| 亚洲av电影在线进入| 久久人妻av系列| 久久久国产欧美日韩av| 99久久国产精品久久久| 男人舔女人的私密视频| 1024视频免费在线观看| 久久久久免费精品人妻一区二区| 久久久久久大精品| 精品第一国产精品| 久久精品国产清高在天天线| 人妻久久中文字幕网| 精品人妻1区二区| 国产亚洲av嫩草精品影院| 久久久久国产精品人妻aⅴ院| 亚洲av日韩精品久久久久久密| 亚洲天堂国产精品一区在线| 国产精品久久久人人做人人爽| 狂野欧美激情性xxxx| 18禁观看日本| 嫁个100分男人电影在线观看| 国产精品国产高清国产av| 午夜福利免费观看在线| 日日爽夜夜爽网站| 亚洲国产欧洲综合997久久,| 国产aⅴ精品一区二区三区波| 国产精品香港三级国产av潘金莲| 美女扒开内裤让男人捅视频| 美女大奶头视频| 亚洲免费av在线视频| 18禁观看日本| 可以在线观看的亚洲视频| 999久久久精品免费观看国产| 一边摸一边抽搐一进一小说| 国产精品亚洲一级av第二区| 嫁个100分男人电影在线观看| 久热爱精品视频在线9| 麻豆国产97在线/欧美 | 亚洲国产中文字幕在线视频| 香蕉久久夜色| 女生性感内裤真人,穿戴方法视频| 成人18禁在线播放| 欧美午夜高清在线| 久久人妻av系列| 亚洲激情在线av| 在线播放国产精品三级| 在线永久观看黄色视频| 九色国产91popny在线| 美女黄网站色视频| 脱女人内裤的视频| 亚洲无线在线观看| 一本综合久久免费| 99国产精品一区二区蜜桃av| 久热爱精品视频在线9| 亚洲五月天丁香| 在线观看66精品国产| 亚洲国产欧美一区二区综合| 午夜两性在线视频| 99在线视频只有这里精品首页| 我的老师免费观看完整版| 欧美色视频一区免费| 身体一侧抽搐| 三级毛片av免费| 午夜福利在线在线| 欧美黄色淫秽网站| av欧美777| 欧美日韩亚洲国产一区二区在线观看| 国产免费av片在线观看野外av| 欧美一级毛片孕妇| 国内揄拍国产精品人妻在线| www.熟女人妻精品国产| 深夜精品福利| 亚洲精品av麻豆狂野| 色综合站精品国产| 欧美乱色亚洲激情| 久久久久久久精品吃奶| 国产精品免费视频内射| 在线观看免费日韩欧美大片| 免费在线观看影片大全网站| 中文字幕精品亚洲无线码一区| 此物有八面人人有两片| 免费看日本二区| 午夜福利欧美成人| 国产aⅴ精品一区二区三区波| 岛国视频午夜一区免费看| 啦啦啦观看免费观看视频高清| 麻豆国产97在线/欧美 | 国产亚洲欧美在线一区二区| 国产亚洲精品一区二区www| 搡老妇女老女人老熟妇| 激情在线观看视频在线高清| 欧美黑人巨大hd| 亚洲成人免费电影在线观看| 亚洲欧美日韩高清专用| 麻豆成人av在线观看| 叶爱在线成人免费视频播放| 欧美zozozo另类| 人人妻人人看人人澡| 久久久精品欧美日韩精品| 两个人看的免费小视频| 国产一区二区在线av高清观看| 中文在线观看免费www的网站 | 久久久国产成人精品二区| 一区二区三区激情视频| 亚洲人与动物交配视频| 不卡一级毛片| 国产日本99.免费观看| 免费在线观看影片大全网站| 免费在线观看成人毛片| 老司机靠b影院| 黄色片一级片一级黄色片| 又紧又爽又黄一区二区| 狠狠狠狠99中文字幕| 亚洲色图av天堂| 欧美日韩乱码在线| 在线a可以看的网站| 亚洲狠狠婷婷综合久久图片| 日本熟妇午夜| 黄色视频,在线免费观看| 搞女人的毛片| 国产一区二区在线av高清观看| 观看免费一级毛片| av有码第一页| 国产在线精品亚洲第一网站| 国产成人精品久久二区二区91| 99国产精品99久久久久| 一二三四社区在线视频社区8| 欧美日韩亚洲国产一区二区在线观看| 美女黄网站色视频| 国产高清videossex| 长腿黑丝高跟| 午夜久久久久精精品| 波多野结衣高清作品| 天堂影院成人在线观看| 在线a可以看的网站| 久久亚洲真实| 国产aⅴ精品一区二区三区波| 黄色视频不卡| 国产又黄又爽又无遮挡在线| 久久久水蜜桃国产精品网| 黄色丝袜av网址大全| 久久人人精品亚洲av| 精品无人区乱码1区二区| 精品少妇一区二区三区视频日本电影| 在线观看免费午夜福利视频| 亚洲熟女毛片儿| 成人av在线播放网站| 国产99白浆流出| 亚洲五月天丁香| 老司机靠b影院| 国产成人影院久久av| 999久久久国产精品视频| 嫁个100分男人电影在线观看| 日本免费一区二区三区高清不卡| 真人一进一出gif抽搐免费| 国产三级中文精品| 亚洲成人久久爱视频| 大型黄色视频在线免费观看| 在线观看免费午夜福利视频| 午夜视频精品福利| 婷婷丁香在线五月| 精华霜和精华液先用哪个| 黄色丝袜av网址大全| 欧美黑人精品巨大| 啦啦啦韩国在线观看视频| 他把我摸到了高潮在线观看| 国产精品久久电影中文字幕| 国产亚洲精品av在线| 欧美黑人巨大hd| 不卡av一区二区三区| 男插女下体视频免费在线播放| 日韩精品中文字幕看吧| 男人的好看免费观看在线视频 | 欧美黑人巨大hd| 一本大道久久a久久精品| 99riav亚洲国产免费| 亚洲一区中文字幕在线| 免费看美女性在线毛片视频| 午夜两性在线视频| 国模一区二区三区四区视频 | 成年女人毛片免费观看观看9| 欧美黑人精品巨大| 99久久综合精品五月天人人| 91麻豆精品激情在线观看国产| 欧美日韩亚洲国产一区二区在线观看| www.999成人在线观看| 国产人伦9x9x在线观看| 舔av片在线| 韩国av一区二区三区四区| 一级片免费观看大全| 国内精品一区二区在线观看| 久久久精品欧美日韩精品| 国产精品免费视频内射| 一个人观看的视频www高清免费观看 | 免费搜索国产男女视频| 两个人的视频大全免费| 日韩有码中文字幕| 午夜a级毛片| 人妻久久中文字幕网| 老汉色av国产亚洲站长工具| 人妻丰满熟妇av一区二区三区| 欧美日韩乱码在线| 熟妇人妻久久中文字幕3abv| 亚洲五月婷婷丁香| 亚洲精品在线美女| 女人被狂操c到高潮| 午夜影院日韩av| 老司机午夜十八禁免费视频| 亚洲激情在线av| 在线看三级毛片| 99久久精品热视频| 99热这里只有精品一区 | 一二三四在线观看免费中文在| 成人av在线播放网站| 国产视频一区二区在线看| 国产成人一区二区三区免费视频网站| 白带黄色成豆腐渣| 国产精品久久久久久人妻精品电影| 成人亚洲精品av一区二区| 亚洲成人精品中文字幕电影| 在线观看免费日韩欧美大片| 国产乱人伦免费视频| 久久精品综合一区二区三区| 看黄色毛片网站| 美女午夜性视频免费| 伦理电影免费视频| 欧美一区二区精品小视频在线| 亚洲精品国产一区二区精华液| 精品一区二区三区视频在线观看免费| 狂野欧美激情性xxxx| 国产欧美日韩精品亚洲av| 亚洲 欧美一区二区三区| 69av精品久久久久久| 亚洲最大成人中文| 日韩有码中文字幕| 色综合亚洲欧美另类图片| 伦理电影免费视频| 久久人妻福利社区极品人妻图片| 国产精品乱码一区二三区的特点| 后天国语完整版免费观看| 成人三级黄色视频| 亚洲国产中文字幕在线视频| 国产精品爽爽va在线观看网站| 一级毛片精品| 1024视频免费在线观看| 国产亚洲精品第一综合不卡| 久久精品国产99精品国产亚洲性色| 一二三四社区在线视频社区8| 亚洲欧美日韩高清在线视频| 亚洲中文av在线| 黑人操中国人逼视频| 在线观看美女被高潮喷水网站 | 熟女电影av网| 给我免费播放毛片高清在线观看| 99re在线观看精品视频| 女人被狂操c到高潮| 国内揄拍国产精品人妻在线| 欧美高清成人免费视频www| 欧美日韩乱码在线| ponron亚洲| 国产成人aa在线观看| 婷婷六月久久综合丁香| 亚洲成av人片在线播放无| 最近最新中文字幕大全免费视频| 久久久久久人人人人人| 久久久精品国产亚洲av高清涩受| 精品乱码久久久久久99久播| 岛国在线观看网站| 亚洲国产欧美一区二区综合| 午夜影院日韩av| 啦啦啦韩国在线观看视频| 男女视频在线观看网站免费 | 女警被强在线播放| 一区二区三区激情视频| 午夜福利成人在线免费观看| 中文字幕人妻丝袜一区二区| 黑人操中国人逼视频| 真人一进一出gif抽搐免费| 久久这里只有精品中国| 99热这里只有是精品50| 久久亚洲真实| 久久国产精品影院| 老司机午夜福利在线观看视频| 久久精品综合一区二区三区| 一本大道久久a久久精品| 国产精品美女特级片免费视频播放器 | 日韩精品青青久久久久久| 少妇的丰满在线观看| 欧美日韩中文字幕国产精品一区二区三区| 国产亚洲精品久久久久久毛片| 精品欧美国产一区二区三| 午夜a级毛片| 国产亚洲精品av在线| av超薄肉色丝袜交足视频| 免费看a级黄色片| 看免费av毛片|