• <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ī)律研究
    在线观看免费视频网站a站| 岛国毛片在线播放| 在线看a的网站| 大香蕉久久成人网| 欧美另类一区| 丝袜脚勾引网站| tube8黄色片| 女性生殖器流出的白浆| 日本av免费视频播放| 2018国产大陆天天弄谢| 亚洲欧美一区二区三区国产| 亚洲人成网站在线观看播放| 视频在线观看一区二区三区| 久久影院123| 老司机亚洲免费影院| 亚洲少妇的诱惑av| 插逼视频在线观看| 日韩不卡一区二区三区视频在线| 18禁在线无遮挡免费观看视频| 国产精品一二三区在线看| 考比视频在线观看| a级毛片在线看网站| 七月丁香在线播放| 日本欧美国产在线视频| 精品视频人人做人人爽| 亚洲精品aⅴ在线观看| 免费不卡的大黄色大毛片视频在线观看| 久久99热这里只频精品6学生| a级毛片黄视频| 日本欧美视频一区| 亚洲精品日韩在线中文字幕| 有码 亚洲区| 五月天丁香电影| av在线老鸭窝| 交换朋友夫妻互换小说| av一本久久久久| 亚洲欧美一区二区三区国产| 在线播放无遮挡| 久久99一区二区三区| 熟女人妻精品中文字幕| av专区在线播放| 午夜福利网站1000一区二区三区| 97在线视频观看| 亚洲欧美中文字幕日韩二区| 蜜桃在线观看..| 亚洲精品乱码久久久久久按摩| 日本猛色少妇xxxxx猛交久久| 日本欧美国产在线视频| 国产精品国产av在线观看| 91午夜精品亚洲一区二区三区| 免费久久久久久久精品成人欧美视频 | 亚洲欧美清纯卡通| 免费少妇av软件| 亚洲欧美一区二区三区黑人 | av专区在线播放| 日日啪夜夜爽| 99九九线精品视频在线观看视频| 亚洲av男天堂| 一本色道久久久久久精品综合| 亚洲,一卡二卡三卡| 香蕉精品网在线| 人妻一区二区av| 汤姆久久久久久久影院中文字幕| 黑人欧美特级aaaaaa片| 欧美亚洲日本最大视频资源| 久久女婷五月综合色啪小说| 如日韩欧美国产精品一区二区三区 | 欧美日韩视频高清一区二区三区二| 日本黄色片子视频| 久久精品人人爽人人爽视色| 黑人高潮一二区| 免费日韩欧美在线观看| 日韩在线高清观看一区二区三区| 2022亚洲国产成人精品| 精品一区在线观看国产| 欧美丝袜亚洲另类| 性高湖久久久久久久久免费观看| 久久久久久久久大av| 国产免费一区二区三区四区乱码| 在线观看一区二区三区激情| 高清午夜精品一区二区三区| 久久99蜜桃精品久久| 亚洲成人一二三区av| av线在线观看网站| 中国美白少妇内射xxxbb| 精品熟女少妇av免费看| 26uuu在线亚洲综合色| 日本vs欧美在线观看视频| 国产国拍精品亚洲av在线观看| 精品少妇久久久久久888优播| 亚洲综合精品二区| 嘟嘟电影网在线观看| 考比视频在线观看| 两个人的视频大全免费| 亚洲色图 男人天堂 中文字幕 | 晚上一个人看的免费电影| 少妇精品久久久久久久| 日韩熟女老妇一区二区性免费视频| 你懂的网址亚洲精品在线观看| 在线观看免费日韩欧美大片 | 十八禁高潮呻吟视频| 97超碰精品成人国产| 日本黄色日本黄色录像| av网站免费在线观看视频| 成人国产麻豆网| 午夜激情久久久久久久| 丰满乱子伦码专区| av免费在线看不卡| 国产精品免费大片| 18禁裸乳无遮挡动漫免费视频| 免费播放大片免费观看视频在线观看| 如日韩欧美国产精品一区二区三区 | 午夜老司机福利剧场| 五月天丁香电影| 极品少妇高潮喷水抽搐| 一个人看视频在线观看www免费| 精品熟女少妇av免费看| 亚洲国产av新网站| 久久99热这里只频精品6学生| 春色校园在线视频观看| 免费黄色在线免费观看| 精品亚洲成a人片在线观看| 黄色欧美视频在线观看| 精品国产露脸久久av麻豆| 看非洲黑人一级黄片| 日韩成人av中文字幕在线观看| 麻豆精品久久久久久蜜桃| 日韩电影二区| 国产男女超爽视频在线观看| 久久精品国产自在天天线| 欧美另类一区| 一区二区三区免费毛片| 亚洲少妇的诱惑av| 99视频精品全部免费 在线| 看非洲黑人一级黄片| 久久影院123| 亚洲高清免费不卡视频| 亚洲av不卡在线观看| 精品一区二区三卡| 在线精品无人区一区二区三| 男女免费视频国产| 日本与韩国留学比较| 国产极品粉嫩免费观看在线 | 亚州av有码| 亚洲色图综合在线观看| 中文字幕人妻熟人妻熟丝袜美| 视频区图区小说| 亚洲,一卡二卡三卡| 精品人妻一区二区三区麻豆| 2022亚洲国产成人精品| 中文字幕av电影在线播放| videos熟女内射| 久久免费观看电影| 91在线精品国自产拍蜜月| 精品亚洲成国产av| 午夜日本视频在线| av卡一久久| 免费黄频网站在线观看国产| 日本与韩国留学比较| 嫩草影院入口| 国产精品国产av在线观看| 国产一区二区在线观看日韩| 曰老女人黄片| 久久久久久久久久久免费av| 日韩强制内射视频| 午夜福利视频在线观看免费| 在线观看免费日韩欧美大片 | 国产精品一区二区在线观看99| 欧美精品人与动牲交sv欧美| 国产成人精品久久久久久| 在线 av 中文字幕| 在线免费观看不下载黄p国产| 亚洲av免费高清在线观看| 国产亚洲精品第一综合不卡 | 秋霞在线观看毛片| 中文字幕久久专区| 欧美xxⅹ黑人| 人体艺术视频欧美日本| 人妻夜夜爽99麻豆av| 亚洲精品aⅴ在线观看| 看非洲黑人一级黄片| 只有这里有精品99| 一级a做视频免费观看| 人体艺术视频欧美日本| 亚洲欧美成人综合另类久久久| 亚洲精品乱久久久久久| 一级片'在线观看视频| 看非洲黑人一级黄片| av国产精品久久久久影院| 午夜免费鲁丝| 国产视频首页在线观看| 少妇猛男粗大的猛烈进出视频| 在线免费观看不下载黄p国产| 狂野欧美白嫩少妇大欣赏| av有码第一页| 国产精品一区二区在线不卡| 午夜福利视频在线观看免费| 看非洲黑人一级黄片| 中文乱码字字幕精品一区二区三区| 麻豆成人av视频| 女人精品久久久久毛片| 欧美成人精品欧美一级黄| 久久人妻熟女aⅴ| 18禁在线播放成人免费| 国产色爽女视频免费观看| 黄片无遮挡物在线观看| 国产高清有码在线观看视频| 赤兔流量卡办理| 国产极品粉嫩免费观看在线 | 肉色欧美久久久久久久蜜桃| 国精品久久久久久国模美| 少妇精品久久久久久久| 日本色播在线视频| 日韩亚洲欧美综合| 亚洲成色77777| 欧美国产精品一级二级三级| 中文字幕人妻熟人妻熟丝袜美| 久久99一区二区三区| 视频中文字幕在线观看| 热99久久久久精品小说推荐| 亚洲精品色激情综合| 婷婷成人精品国产| 这个男人来自地球电影免费观看 | 伊人亚洲综合成人网| 国精品久久久久久国模美| 在线观看免费视频网站a站| 亚洲欧洲国产日韩| 少妇人妻精品综合一区二区| 国产男女超爽视频在线观看| 久久99热这里只频精品6学生| 日韩电影二区| 亚洲精品乱久久久久久| 国产精品久久久久久久久免| 视频在线观看一区二区三区| 五月天丁香电影| 性高湖久久久久久久久免费观看| 天堂8中文在线网| 欧美97在线视频| 少妇被粗大猛烈的视频| 久久99热6这里只有精品| 免费大片黄手机在线观看| 一本色道久久久久久精品综合| 婷婷成人精品国产| av在线观看视频网站免费| 婷婷色综合www| 高清视频免费观看一区二区| 亚洲欧洲精品一区二区精品久久久 | 91精品三级在线观看| 国产亚洲欧美精品永久| 丁香六月天网| 亚洲欧洲国产日韩| 久久久久久久大尺度免费视频| 边亲边吃奶的免费视频| 亚洲四区av| 国产欧美亚洲国产| 校园人妻丝袜中文字幕| 久久亚洲国产成人精品v| 一本—道久久a久久精品蜜桃钙片| 国产精品久久久久久久久免| 街头女战士在线观看网站| 国产老妇伦熟女老妇高清| 久久精品人人爽人人爽视色| 亚洲精品av麻豆狂野| 久久久精品免费免费高清| 精品一区二区免费观看| 亚洲不卡免费看| 国产精品成人在线| 美女国产高潮福利片在线看| 精品熟女少妇av免费看| 国语对白做爰xxxⅹ性视频网站| 99re6热这里在线精品视频| 高清av免费在线| 国产乱人偷精品视频| 日本午夜av视频| 亚洲人与动物交配视频| 肉色欧美久久久久久久蜜桃| 人妻人人澡人人爽人人| 色网站视频免费| 日本av手机在线免费观看| 一区二区日韩欧美中文字幕 | 99九九线精品视频在线观看视频| 欧美97在线视频| 国产成人91sexporn| 亚洲av综合色区一区| 18禁观看日本| 不卡视频在线观看欧美| 精品一区在线观看国产| 另类亚洲欧美激情| 国国产精品蜜臀av免费| 免费人成在线观看视频色| 七月丁香在线播放| 精品亚洲成a人片在线观看| 男人操女人黄网站| 欧美亚洲日本最大视频资源| 午夜av观看不卡| 秋霞在线观看毛片| 久久久午夜欧美精品| 亚洲国产精品成人久久小说| 人体艺术视频欧美日本| 精品久久久精品久久久| 国产精品一区二区在线观看99| 久久99热这里只频精品6学生| 两个人免费观看高清视频| 一级毛片 在线播放| 国产一级毛片在线| 日本色播在线视频| 欧美日韩亚洲高清精品| 日本色播在线视频| 欧美日韩亚洲高清精品| 久久精品国产自在天天线| 在线观看www视频免费| 久久热精品热| 美女福利国产在线| 青春草视频在线免费观看| 大陆偷拍与自拍| 老熟女久久久| 大香蕉久久网| 亚洲怡红院男人天堂| 性色avwww在线观看| 纯流量卡能插随身wifi吗| 国产精品国产av在线观看| 久久热精品热| 日韩不卡一区二区三区视频在线| 大片电影免费在线观看免费| 美女福利国产在线| 91精品国产九色| 人妻一区二区av| 亚洲精品自拍成人| 亚洲内射少妇av| 亚洲精品第二区| 高清视频免费观看一区二区| 少妇人妻 视频| 久久久久久久大尺度免费视频| 亚洲美女黄色视频免费看| 午夜激情福利司机影院| 亚洲精品亚洲一区二区| 国产精品国产三级国产专区5o| 日韩av不卡免费在线播放| 丝袜脚勾引网站| 午夜福利网站1000一区二区三区| 黄色视频在线播放观看不卡| 成人毛片a级毛片在线播放| 亚洲精品日本国产第一区| 高清午夜精品一区二区三区| 亚洲欧美清纯卡通| 满18在线观看网站| 国产精品久久久久久精品电影小说| 欧美三级亚洲精品| 女性被躁到高潮视频| 亚洲欧美精品自产自拍| 男人爽女人下面视频在线观看| 欧美三级亚洲精品| 久久久久网色| 满18在线观看网站| 欧美精品亚洲一区二区| 婷婷色麻豆天堂久久| 水蜜桃什么品种好| 国产黄色免费在线视频| 国产69精品久久久久777片| 欧美亚洲 丝袜 人妻 在线| 国产黄色视频一区二区在线观看| 亚洲精品乱码久久久久久按摩| 国产黄色视频一区二区在线观看| 2022亚洲国产成人精品| 超碰97精品在线观看| 2022亚洲国产成人精品| 高清在线视频一区二区三区| 天天躁夜夜躁狠狠久久av| 成人无遮挡网站| 另类亚洲欧美激情| 自拍欧美九色日韩亚洲蝌蚪91| 欧美xxxx性猛交bbbb| 99精国产麻豆久久婷婷| 免费观看性生交大片5| 日韩一区二区视频免费看| 黄片无遮挡物在线观看| 少妇 在线观看| 99视频精品全部免费 在线| 中文字幕最新亚洲高清| av在线app专区| 中文欧美无线码| 免费日韩欧美在线观看| 欧美日本中文国产一区发布| 国产片内射在线| 91精品国产九色| 亚洲中文av在线| 99热网站在线观看| 午夜老司机福利剧场| 久久国产精品男人的天堂亚洲 | 国产成人精品福利久久| 日本91视频免费播放| 成人毛片a级毛片在线播放| 久久综合国产亚洲精品| 美女视频免费永久观看网站| 免费看不卡的av| 亚洲av男天堂| 另类精品久久| www.av在线官网国产| 2018国产大陆天天弄谢| av在线app专区| 午夜av观看不卡| 国产精品久久久久久精品电影小说| 亚洲,欧美,日韩| 制服诱惑二区| 欧美3d第一页| 少妇高潮的动态图| 欧美国产精品一级二级三级| 日本欧美视频一区| 最新的欧美精品一区二区| 亚洲欧美清纯卡通| 九草在线视频观看| 亚洲精品国产色婷婷电影| 在线观看人妻少妇| 成年美女黄网站色视频大全免费 | 久久久精品区二区三区| 麻豆成人av视频| 国产精品一区二区在线不卡| 亚洲美女搞黄在线观看| 日韩av在线免费看完整版不卡| 男男h啪啪无遮挡| 国产日韩欧美在线精品| 亚州av有码| 在线亚洲精品国产二区图片欧美 | 最近2019中文字幕mv第一页| 99视频精品全部免费 在线| 在线观看美女被高潮喷水网站| 国产成人精品一,二区| 伊人亚洲综合成人网| 亚洲国产最新在线播放| 久久av网站| 搡老乐熟女国产| 在线观看一区二区三区激情| 十分钟在线观看高清视频www| 丝袜在线中文字幕| 亚洲激情五月婷婷啪啪| 国产亚洲精品久久久com| 免费大片黄手机在线观看| 黑人欧美特级aaaaaa片| 天天影视国产精品| 亚洲熟女精品中文字幕| 精品卡一卡二卡四卡免费| 欧美日本中文国产一区发布| 丰满乱子伦码专区| 在线观看免费视频网站a站| 亚洲伊人久久精品综合| 国产在线视频一区二区| 久久久久人妻精品一区果冻| 成年av动漫网址| 超色免费av| 亚洲人与动物交配视频| 你懂的网址亚洲精品在线观看| 亚洲高清免费不卡视频| 亚洲性久久影院| 欧美少妇被猛烈插入视频| 赤兔流量卡办理| 国产男女内射视频| 97在线人人人人妻| 欧美日韩视频精品一区| 久久久精品区二区三区| 婷婷色麻豆天堂久久| 2021少妇久久久久久久久久久| 国产黄片视频在线免费观看| 老女人水多毛片| 中文天堂在线官网| 国产黄片视频在线免费观看| 丝袜喷水一区| 免费高清在线观看视频在线观看| 亚洲精品av麻豆狂野| 黄色欧美视频在线观看| 国产黄频视频在线观看| 日韩欧美精品免费久久| 考比视频在线观看| 国产国语露脸激情在线看| 99久国产av精品国产电影| 亚洲性久久影院| 色婷婷av一区二区三区视频| 欧美精品亚洲一区二区| 国产成人精品福利久久| 一边亲一边摸免费视频| 晚上一个人看的免费电影| 在线看a的网站| 91在线精品国自产拍蜜月| 日韩视频在线欧美| 一区在线观看完整版| 午夜激情福利司机影院| 国产精品一区二区在线观看99| 国产精品欧美亚洲77777| 亚洲精品国产av蜜桃| 交换朋友夫妻互换小说| 日日摸夜夜添夜夜爱| 午夜免费观看性视频| 欧美bdsm另类| a级毛片黄视频| 成人免费观看视频高清| 午夜久久久在线观看| 九九爱精品视频在线观看| 我的老师免费观看完整版| 黄色毛片三级朝国网站| 日本欧美视频一区| 3wmmmm亚洲av在线观看| 国产片特级美女逼逼视频| 男女床上黄色一级片免费看| 亚洲国产欧美在线一区| 午夜视频精品福利| a级毛片在线看网站| 99国产精品99久久久久| 精品人妻在线不人妻| 菩萨蛮人人尽说江南好唐韦庄| 在线观看66精品国产| 国产有黄有色有爽视频| 成人影院久久| 亚洲第一av免费看| 国产欧美日韩一区二区精品| 美国免费a级毛片| 99国产精品免费福利视频| 大码成人一级视频| 精品国产超薄肉色丝袜足j| 精品国产乱码久久久久久男人| 成人18禁在线播放| 亚洲精品在线观看二区| 国产色视频综合| 9热在线视频观看99| 露出奶头的视频| 50天的宝宝边吃奶边哭怎么回事| 精品国产一区二区三区四区第35| 视频区图区小说| 在线观看www视频免费| 亚洲黑人精品在线| 香蕉久久夜色| 亚洲精品国产精品久久久不卡| 成人18禁在线播放| 国产亚洲精品第一综合不卡| 久久午夜亚洲精品久久| 欧美日韩av久久| 香蕉丝袜av| 国产av国产精品国产| 国产淫语在线视频| 午夜老司机福利片| 黄色成人免费大全| 久久中文字幕一级| 久久 成人 亚洲| 天堂动漫精品| 搡老乐熟女国产| 在线观看免费视频网站a站| 亚洲免费av在线视频| 狠狠精品人妻久久久久久综合| 欧美日韩成人在线一区二区| 激情在线观看视频在线高清 | 亚洲,欧美精品.| svipshipincom国产片| 国产日韩欧美在线精品| 一进一出抽搐动态| 妹子高潮喷水视频| 丝袜人妻中文字幕| 亚洲一区中文字幕在线| 精品国产乱码久久久久久男人| 精品高清国产在线一区| 亚洲精品美女久久av网站| 丰满迷人的少妇在线观看| 国产高清激情床上av| 精品午夜福利视频在线观看一区 | www.999成人在线观看| 国产aⅴ精品一区二区三区波| av在线播放免费不卡| e午夜精品久久久久久久| 久久人人97超碰香蕉20202| 久9热在线精品视频| 国产亚洲一区二区精品| 中国美女看黄片| 9热在线视频观看99| 久久99一区二区三区| 欧美 亚洲 国产 日韩一| 久久中文字幕一级| 在线观看舔阴道视频| 90打野战视频偷拍视频| 啦啦啦视频在线资源免费观看| 欧美黄色片欧美黄色片| 最新美女视频免费是黄的| 热99久久久久精品小说推荐| √禁漫天堂资源中文www| 国产精品香港三级国产av潘金莲| 精品欧美一区二区三区在线| 亚洲国产毛片av蜜桃av| 久久久久久亚洲精品国产蜜桃av| 男女午夜视频在线观看| 99热网站在线观看| 日本欧美视频一区| 亚洲国产欧美在线一区| 天天影视国产精品| 91麻豆av在线| 亚洲国产看品久久| 国产男女内射视频| 久久av网站| 久久 成人 亚洲| 午夜福利免费观看在线| 少妇猛男粗大的猛烈进出视频| 久久久久精品人妻al黑| xxxhd国产人妻xxx| 一二三四社区在线视频社区8| 50天的宝宝边吃奶边哭怎么回事| 深夜精品福利| 黄片播放在线免费| 午夜福利,免费看| 久久人人爽av亚洲精品天堂| 久久中文看片网| 精品一区二区三区av网在线观看 | 中亚洲国语对白在线视频| 国产精品1区2区在线观看. | 国产xxxxx性猛交| 国产单亲对白刺激| 国产日韩欧美亚洲二区| 18禁裸乳无遮挡动漫免费视频| 日本欧美视频一区| 法律面前人人平等表现在哪些方面| 日韩欧美一区视频在线观看|