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

    基于球元重建與大渦模擬耦合并行算法的數(shù)值模擬

    2017-04-12 03:31:22郭飛宏王澤宇仲兆平王肖祎
    東南大學學報(自然科學版) 2017年2期
    關鍵詞:曳力大渦柱形

    郭飛宏 王澤宇 仲兆平 王肖祎

    (1東南大學能源與環(huán)境學院, 南京 210096)(2東南大學建筑設計院有限公司, 南京 210096)

    基于球元重建與大渦模擬耦合并行算法的數(shù)值模擬

    郭飛宏1王澤宇1仲兆平1王肖祎2

    (1東南大學能源與環(huán)境學院, 南京 210096)(2東南大學建筑設計院有限公司, 南京 210096)

    為了更好地模擬柱形顆粒與石英砂的流動,利用球元重建方法對柱形顆粒進行球元重建.采用大渦模擬中的SGS(sub-grid scale)亞格子模型,對柱形顆粒造成的湍流進行解析.模擬中耦合并行算法,縮短了模型計算時間.基于上述方法,模擬了表觀氣速為1.5 m/s的柱形顆粒與石英砂的混合流動過程. 通過流動瞬時圖像和壓力脈動分析表明,模擬結果與實驗結果相似,表明耦合并行算法能較好地模擬流化床中柱形顆粒與石英砂的流動.當柱形顆粒尺寸較大時,加劇了湍流效應對混合流動的影響.模擬的z方向氣體平均速度呈現(xiàn)中間高兩邊低的對稱分布狀態(tài),這與床內的氣泡和柱形顆粒分布有關,符合實際流化規(guī)律.

    柱形顆粒;石英砂;模擬;球元重建;大渦模擬;并行算法

    傳統(tǒng)的離散元和計算流體力學方法(DEM-CFD)只能對球形顆粒的流動進行模擬,然而在實際的硫化床氣固反應過程中,有很多柱形顆粒[1-3].柱形顆粒形狀改變了顆粒與湍流的相互作用,不同形狀顆粒在流場中的碰撞和受力有顯著的差異[4-5].在研究方法上,缺少對柱形顆粒系統(tǒng)的實驗研究,針對球形顆粒的研究方法應用于柱形顆粒系統(tǒng)時,容易產生很大的誤差.對流態(tài)化下柱形顆粒的運動開展數(shù)值模擬也具有一定的難度,如何構建柱形顆粒,采用較為統(tǒng)一的數(shù)理模型描述柱形顆粒的受力和運動,是至今仍未完全解決的難題.國內外很多研究者對DEM-CFD方法進行了改進,Hilton等[6]采用超曲面法建立非球形顆粒模型,模擬了非球形顆粒在流化床中的流動.Zhong等[7]對氣固流化床中的柱形顆粒流動進行了模擬,將柱形顆粒抽象為球元的組合體.Ren等[8-9]對流化床中玉米狀顆粒流動進行了模擬,與Zhong等[7]的方法一樣,仍然是將非球形顆粒抽象為球元的結合體,顆粒所受氣相曳力近似等于球元所受曳力之和.盧洲等[10]采用DEM-CFD方法模擬了彎管中輸送過程的柱狀顆粒.Oschmann等[11]采用球元聚集法和多面體法構造柱形顆粒并描述顆粒之間的運動.當顆粒尺寸較大時,無法忽視流化床湍流對氣固流動的影響.大渦數(shù)值模擬對于尺寸大于網格的湍流,采用直接數(shù)值求解[12-14];對于尺寸小于網格的小尺度的湍流脈動,采用SGS模型求解.由于大渦模擬能在較短的計算時間內獲得比雷諾平均模擬更多更準確的湍流信息.因此,大渦模擬在湍流解析上得到愈來愈廣泛的運用.

    本文采用球元重建方法,將柱形顆粒表面所受的曳力轉化為虛擬球元所受的曳力.采用大渦模型,并引入并行算法,從而實現(xiàn)柱形顆粒與石英砂流動的DEM-CFD并行模擬.

    1 基本理論

    1.1 球元重建

    計算CFD氣固曳力時,采用直徑為0.8 mm虛擬球元重建一個長10 mm、直徑10 mm的柱形顆粒模擬真實柱形顆粒(見圖1).本文采用0.8 mm的虛擬球元,一方面虛擬球元的直徑與石英砂顆粒直徑大小相同,利于計算;另一方面,球元之間成正四面體緊密堆積,此時網格的空隙率約為0.2.當CFD計算網格完全處于柱形顆粒內部時,既能保證計算過程有很好的收斂性,又能借助虛擬球元實現(xiàn)在柱形顆粒內部對流體形成較大的阻力,使得柱形顆粒內部流體速度接近于零.在計算柱形顆粒與其他顆粒及壁面碰撞時,由直徑為1.5 mm的虛擬球元重建形成一個中空的柱形殼體,如圖1(c)所示.選用1.5 mm的虛擬球元有利于計算,同時保證了柱形顆粒重構的準確性.采用殼體布置方式的原因是:DEM計算過程中,計算顆粒碰撞力所用的時間很長,且計算用時隨著顆粒數(shù)量的增加呈指數(shù)級增加.研究發(fā)現(xiàn),在雙組分DEM-CFD模擬中,可以忽略柱形顆粒內部球元的受力狀況,這樣可以縮短計算時間.為了盡可能地接近柱形顆粒的表面,球元與球元之間存在重疊.

    (a) 柱形顆粒

    (b) CFD球元重建

    (c) DEM球元重建

    為了驗證球元重建的有效性,對柱形顆粒曳力系數(shù)Cd進行定量分析.通過改變虛擬球元大小,進而改變柱形顆粒內部的空隙率,模擬柱形顆粒曳力系數(shù)隨雷諾數(shù)的變化.如圖2所示,柱形顆粒曳力系數(shù)Cd隨著雷諾數(shù)的增加而逐漸減小,并且最后趨于穩(wěn)定;柱形顆粒曳力系數(shù)Cd隨著虛擬球元所占體積分數(shù)a的增加而增加.當Re<3 000,α=0.7時,曳力系數(shù)與實驗曲線吻合良好;而當Re>3 000,α=0.7時,曳力系數(shù)與實驗值有輕微偏離,而α=0.5時,與實驗曲線吻合良好.本次模擬的截面氣速在1.0~2.0 m/s之間,因此對應的雷諾數(shù)一般小于3 000,所以控制虛擬球元大小使α=0.7,可以獲得較好的流動模擬結果.

    圖2 柱形顆粒曳力系數(shù)Cd隨雷諾數(shù)Re的變化

    1.2 大渦模擬

    對于湍流的解析,采用大渦模型中普遍采用的SGS(sub-grid scale)亞格子模型[15].在對傳統(tǒng)N-S方程進行過濾后,大渦模型的連續(xù)性方程和動量方程為

    (1)

    (2)

    (3)

    (4)

    (5)

    式中,Cs為SGS常數(shù),在本次模擬研究中取0.13;Δ為特征長度,根據孔隙率的不同,其表達式為

    (6)

    1.3 氣固作用力

    (7)

    (8)

    本研究中,根據空隙率的不同分別采用Ergun和en&Yu模型來表示氣固曳力系數(shù)β,表達式為

    (9)

    式中,Cd為流場內的顆粒曳力系數(shù),由Shiller和Naumann方程以及牛頓曳力法則計算,即

    (10)

    式中, 雷諾數(shù)Re的表達式為

    (11)

    (12)

    式中,αp和αc分別為單個計算網格內石英砂和生物質的體積分數(shù);αfic為生物質內虛擬顆粒的體積分數(shù);dp和dfic為石英砂顆粒和虛擬顆粒的直徑.

    1.4 離散元方法

    Tsuji等[5]提出了DEM-CFD方法,將DEM模型引入到稠密氣固兩相流動的模擬研究中.其基本思想是假設顆粒在碰撞過程中可以有輕微的重疊現(xiàn)象并可以持續(xù)一定的時間.如圖3所示,采用彈簧、阻尼器和滑移器來表示顆粒與顆粒間或顆粒與壁面間的碰撞.

    (a) 法線方向 (b) 切線方向

    1—彈簧;2—滑移器;3—阻尼器

    圖3 DEM顆粒碰撞示意圖

    顆粒發(fā)生碰撞時其作用力在碰撞點上,這個力稱為顆粒接觸力Fc.顆粒接觸力被分為法向力Fcn和切向力Fct,法向力使顆粒發(fā)生平動,切向力由于作用方向線與指點不在同一直線上,會產生力矩使顆粒發(fā)生轉動.法向力Fcn和切向力Fct由下式計算:

    Fcn=-kndn-ηnvn

    (13)

    vn=(vr·n)n

    (14)

    Fct=-ktdt-ηtvt

    (15)

    vt=vr-(vr·n)n+rwp×n

    (16)

    式中,dn和dt分別表示法向和切向的相對位移;vr為顆粒間相對速度;vt表示碰撞點的相對滑移速度;kn和kt分別表示法向和切向方向彈簧的彈性系數(shù);ηn和ηt分別為法向和切向的阻尼系數(shù);n為單位法向量;wp為顆粒間相對角速度.阻尼系數(shù)η可由恢復系數(shù)e確定,即

    (17)

    (18)

    α=-(1/π)lne

    (19)

    式中,m為質量;k為彈性系數(shù).在計算柱形顆粒接觸力時,只需計算柱形顆粒表面的球元與其他顆?;虮诿娴呐鲎?柱形顆粒表面球元碰撞接觸力計算方法與球形顆粒是相同的.在柱形顆粒與其他顆粒碰撞時,可能只有一個球元受力,也可能有多個球元同時受力.球形顆粒與柱形顆粒的碰撞視為球形顆粒與球元的碰撞.柱形顆粒總的受力為柱形顆粒內所有球元的受力.柱形顆粒的運動可以認為是由球元的平動和轉動疊加產生的.根據牛頓第二定律,球形顆粒i的平動和轉動方程為

    (20)

    (21)

    式中,mp為顆粒質量;vpi為顆粒i的速度矢量;fpi為顆粒i所受氣相曳力;∑Fci為顆粒i所受接觸力總和;Ip為顆粒的轉動慣量;ωpi為顆粒i的角速度矢量;∑Mpi為顆粒i所受力矩總和.

    1.5 并行算法

    對柱形顆粒與石英砂流動進行模擬,采用了區(qū)域分解技術和顆粒分配技術相結合的方法[16-17].基于信息傳遞的并行編程技術(message passing interface,MPI)對CFD求解區(qū)域和顆粒碰撞計算進行多線程并行算法優(yōu)化.在CFD計算中,連續(xù)相的速度、壓力等重要參數(shù)需要耦合迭代計算,一個時間步長需要進行上千次的迭代計算,因此連續(xù)相的速度、壓力的計算是一個十分復雜且耗時的過程,對整個計算區(qū)域按照x方向進行區(qū)域劃分.程序運行后,在同一時刻存在N個線程對不同計算區(qū)域進行CFD求解,這樣單位時間內進行的迭代次數(shù)就會明顯增加.而每個計算域都比原先劃分多了左右2列虛擬邊界作為數(shù)據緩沖區(qū),這些區(qū)域存儲了與當前計算域相鄰的邊界信息.在一個迭代計算步驟中,分區(qū)的邊界信息通過通信句柄MPI_send將計算域的邊界信息傳輸至虛擬邊界,通過通信句柄MPI_recv從相鄰區(qū)域的緩沖區(qū)獲取邊界信息以保證整個流場迭代計算的連續(xù)性.在DEM計算中,由于顆粒數(shù)量較多,通常都有幾十萬個,因此顆粒碰撞力占用了絕大部分的DEM計算時間,需要對這部分計算進行并行算法優(yōu)化.系統(tǒng)統(tǒng)一對柱形顆粒中的虛擬球元和石英砂顆粒分配ID號,然后按照ID號的大小將不同顆粒的計算任務分配給不同的線程處理,從而提高計算效率.

    2 模型的建立

    2.1 基本流程圖

    本文模擬的計算流程如圖4所示,采用Fortran語言編寫程序.在一個時間步長內,模擬程序主要進行了兩大部分計算:① 連續(xù)相求解(CFD).基于正交結構化交錯網格和有限體積離散方法求解連續(xù)相的流場.將動量方程離散后采用SIMPLE算法進行計算.首先假設一個壓力場已知,求解流場的速度,然后修正壓力和速度,最后進行迭代計算判斷收斂性.② 離散相求解(DEM).首先檢索網格內的顆粒,判斷顆粒是柱形顆粒中的球元顆粒還是石英砂顆粒.采用顆粒碰撞軟球模型和球元重建理論,計算顆粒間的碰撞力,得到柱形顆粒(虛擬球元)和石英砂的運動軌跡.

    圖4 DEM-CFD模擬計算流程

    2.2 計算區(qū)域及定解條件

    構建80 mm×40 mm×500 mm的計算區(qū)域.如果網格尺寸過大,會使模擬計算精度大大降低;如果網格尺寸過小,則氣固曳力公式不適用.在CFD計算中,網格尺寸一般為顆粒尺寸的2.5~4.0倍.在DEM計算中,網格尺寸的要求相對較低,一般為顆粒尺寸的1.0~1.5倍.模擬的計算區(qū)域與實驗中的流化床幾何尺寸一致,在實驗中觀察到流化床上半部分幾乎沒有顆粒的存在,為節(jié)省模擬時間,在模擬中只計算床內高度小于500 mm的流動情況.流化床的幾何尺寸見圖5.流化床模擬計算中,壁面采用無滑移邊界條件.入口邊界采用速度入口邊界條件,布風板上氣孔內設置為速度入口邊界條件,其余部分為壁面的無滑移邊界條件.出口邊界采用常用的壓力出口邊界條件.對于冷態(tài)模擬,床內壓強和氣相速度的初始條件對模擬結果的影響很小,本研究設置初始壓強為0,初始氣相速度為0.

    (a) 主視圖

    (b) 側視圖

    2.3 基本參數(shù)

    流化床柱形顆粒-石英砂雙組分流動模擬主要參數(shù)見表1.

    表1 柱形顆粒-石英砂雙組分流動模擬參數(shù)

    3 結果與討論

    3.1 流動圖像分析

    模擬結果采用Tecplot360圖像后處理軟件輸出.圖6為柱形顆粒質量分數(shù)為10%、表觀氣速為1.5 m/s時,柱形顆粒與石英砂混合流動模擬圖.圖7為該工況下的流動瞬時圖.對比模擬圖和流動瞬時圖可以發(fā)現(xiàn),模擬結果與實驗結果得到的圖像較為相似.當t=0 s時,柱形顆粒分布在石英砂上面,流化床處于靜止狀態(tài),床高約為60 mm.隨著時間的增加,顆粒開始流化:當t=0.3 s時,床層下部形成比較明顯的大氣泡,柱形顆粒和石英砂在氣泡的作用下被吹起,床層“膨脹”的高度為150 mm左右,約為靜止床高的2.5倍,此時柱形顆粒主要集中在石英砂上部,雖然被氣泡帶到一定高度,但與石英砂并未較好地混合.當t=0.5 s時,氣泡逐漸增大和上升,在床層中部形成較大的氣泡區(qū)間,此時柱形顆粒已經到達一定的高度,最大位置變化不大,但是由于中間氣泡的上升,部分石英砂和柱形顆粒受到中間力的作用,開始向床兩側的壁面移動.當t=1.2 s時,氣泡破裂成較小氣泡,柱形顆粒和石英砂回落,并再次與上升的氣泡相聚,造成氣泡重復破裂(見圖6),整個過程有利于柱形顆粒和石英砂的混合.在0~1.2 s的過程中,氣泡的形成—生長—破裂、柱形顆粒的上升和下落、柱形顆粒與石英砂的混合程度在逐漸增加.且圖6的氣固流動行為與圖7中實驗拍攝的圖像行為十分相似,模擬效果較好.

    (a) t=0 s

    (b) t=0.3 s

    (c) t=0.5 s

    (d) t=1.2 s

    (a) t=0 s

    (b) t=0.3 s

    (c) t=0.5 s

    (d) t=1.2 s

    當t=1.0 s時,柱形顆粒與石英砂的流動與模擬如圖8所示.此時,較大柱形顆粒的回落會對區(qū)域A內的大氣泡造成破裂,形成較小氣泡.回落的顆粒再次與上升的氣泡相聚,造成氣泡重復破裂,促進柱形顆粒和石英砂的混合.圖8中的流動圖像和模擬圖像均較好地反映出了這一現(xiàn)象.同時,由于柱形顆粒體積較大,四周的石英砂與其碰撞后容易向外側運動,在柱形顆粒周圍形成較大空隙.這樣,在計算氣固曳力時就需要考慮湍流對氣固流動的影響.

    (a) 模擬圖 (b) 流動圖

    3.2 壓力脈動分析

    圖9為不同床層高度h的壓降模擬值與實驗值對比.由于最小流化速度為0.5 m/s,因此模擬沒有包含固定床工況(表觀氣速v<0.5 m/s).當床高h=50 mm時,實驗值的床層壓降隨著氣速的增加,先上升后下降并趨于穩(wěn)定,流化后的平均壓降ΔP為1.3 kPa.模擬值主要處于流化階段,床層平均壓降約為1.2 kPa.當床高h=85 mm時,壓降

    (a) h=50 mm

    (b) h=85 mm

    隨著氣速的增加,先上升后下降并趨于穩(wěn)定,流化后的平均壓降約為1.7 kPa,模擬的平均壓降約為1.6 kPa.由于在模擬過程中布風條件、床料形狀和尺寸等方面做了一系列簡化假設,因而模擬結果與實驗結果存在一定的偏差,相同氣速條件下床層壓差的偏差量在10%以內.當表觀氣速大于臨界流化風速時,模擬方法能較為準確地模擬分析真實柱形顆粒與石英砂混合流動時的壓力脈動.

    3.3 湍流分析

    在傳統(tǒng)的流化床中,氣速分布均勻,雷諾數(shù)比較小,湍流對氣固流動的影響不很明顯.但是在柱形顆粒和石英砂的混合流動過程中,柱形顆粒尺寸比較大,加劇了湍流效應對流動的影響,因此有必要對湍流進行分析.大渦模擬的基本思想是對大尺度湍流用直接數(shù)值模擬方法求解,對小尺度湍流在添加亞格子尺度應力后,采用SGS模型求解.小尺度湍流對大尺度湍流的影響主要表現(xiàn)在耗散效應上,也就是說大渦模擬中的能量耗散是從大渦向小渦串級傳輸?shù)男问匠霈F(xiàn).SGS湍流黏度的大小在一定程度上反映了小渦對大渦的影響.利用SGS湍流黏度來定量分析柱形顆粒對氣體湍流耗散產生的影響.圖10為不同氣速下SGS湍流黏度的分布.圖中紅色為較大的湍流黏度,藍色為較小的湍流黏度,顏色越深對應的SGS湍流黏度越大.由于噴口附近進氣速度較大,雷諾數(shù)較高,在噴口上方很小區(qū)域會形成湍流黏度最大的區(qū)域.當表觀氣速v=1.0 m/s時,氣體速度較低,床內柱形顆粒和石英砂的流動不很充分.此時,氣體呈氣拴式由床層底部向上緩慢運動,顆粒上下混合不均勻.柱形顆粒聚集在石英砂上部,氣體穿過此區(qū)域后在顆粒層上部形成了較強的尾渦,從而增加了氣體的黏性耗散,SGS湍流黏度較大.當表觀氣速v=1.5 m/s時,流化床內產生的氣泡變大,混合流動加劇,柱形顆粒和石英砂混合較為充分,顆粒在整個床層內分布較為分散,此時SGS湍流黏度分布相對于表觀氣速v=1.0 m/s時較為均勻.

    (a) v=1.0 m/s

    (b) v=1.5 m/s

    3.4 氣相速度

    圖11為不同表觀氣速下,不同床層高度z方向氣體平均速度vz的分布.流化床x方向的長度為80 mm.當表觀氣速v=1.0 m/s時,低床層(h=20,40 mm)z方向氣體平均速度受噴口進氣影響較大,此時的氣相速度范圍為1.4~1.8 m/s,平均速度分別為1.6和1.5 m/s.當h=80 mm時,z方向氣體平均速度明顯下降,呈現(xiàn)中間低兩邊高的波動分布狀態(tài),說明氣體穿過底部床料層后向壁面偏移,并沿著壁面向上運動.這主要是由于柱形顆粒的布置和表觀氣速較低造成的.如圖7(a)所示,柱形顆粒開始集中在x方向中心位置,顆粒本身的重力會對底部中心造成附加壓力,這種附加的阻力導致氣體向壓力較低的兩側壁面運動.當表觀氣速較小(1.0 m/s)時,氣體呈氣拴式由床層底部向上緩慢運動,顆粒上下混合不均勻,這種附加壓力會持續(xù)一段時間.z方向氣體平均速度的波動分布還與布風板的設計有關.

    (a) v=1.0 m/s

    (b) v=1.5 m/s

    當表觀氣速v=1.5 m/s時,床內z方向氣體平均速度呈現(xiàn)中間高兩邊低的對稱分布狀態(tài).當床層高度從40 mm增加到200 mm時,z方向氣體最大速度分別從3.0 m/s下降到1.5 m/s;當h=40,80 mm時,中間區(qū)域氣體縱向速度是左、右兩側靠近壁面處的2倍.這是由于在中間床層產生了較大的氣泡(見圖7(b)),床料顆粒受到氣泡的抬升和排擠作用向兩側回落,因此氣體在中間受到床料的阻力比靠近壁面要小,造成了氣體平均速度呈現(xiàn)中間高兩邊低的對稱分布狀態(tài).表觀氣速增加后,h=20 mm的低床層同時受到噴口進氣和氣固混合的影響.當表觀氣速為1.5 m/s時,z方向氣體平均速度波動性明顯加大,靠近壁面x=70~80 mm位置處,氣相速度較低.這主要是因為表觀氣速的增加,加劇了柱形顆粒與石英砂的混合,而這種混合并不均勻.運動過程中受到多種因素的影響,柱形顆粒和石英砂的分布不均勻,體現(xiàn)在z方向氣體平均速度的波動(見圖7(d)).當h=200 mm時氣體已經穿過床料層,縱向速度分布較為均勻,接近表觀氣速1.5 m/s.

    4 結論

    1) 采用球元重建方法,構建了柱形顆粒,利用大渦模擬中的SGS亞格子模型,對湍流進行解析,同時耦合并行算法,縮短了計算時間.

    2) 基于球元重建與大渦模擬耦合并行算法對柱形顆粒和石英砂的混合流動進行了模擬,模擬結果與實驗結果相似.

    3) 模擬結果表明:較大柱形顆粒的上升和回落會造成氣泡的破裂,促進柱形顆粒和石英砂的混合;柱形顆粒尺寸較大,加劇了湍流效應對流動的影響,有必要對湍流進行解析;當表觀氣速v=1.5 m/s時,床內z方向氣體平均速度呈現(xiàn)中間高兩邊低的對稱分布狀態(tài),這主要與床內的氣泡和柱形顆粒分布有關,符合實際流化規(guī)律.

    References)

    [1]Khan A A, de Jong W, Jansens P J, et al. Biomass combustion in fluidized bed boilers: Potential problems and remedies[J].FuelProcessingTechnology, 2009, 90(1): 21-50. DOI: 10.1016/j.fuproc.2008.07.012.

    [2]Cheng H, Hu Y. Municipal solid waste (MSW) as a renewable source of energy: Current and future practices in China[J].BioresourceTechnology, 2010, 101(11): 3816-3824. DOI:10.1016/j.biortech.2010.01.040.

    [3]Hilton J E, Cleary P W. The influence of particle shape on flow modes in pneumatic conveying[J].ChemicalEngineeringScience, 2011, 66(3): 231-240. DOI:10.1016/j.ces.2010.09.034.

    [4]Wang C, Zhong Z, Wang X. Microscopic flow characteristics in fluidized bed of cylinder-shaped particles[J].KoreanJournalofChemicalEngineering, 2015, 32(12): 2384-2393. DOI:10.1007/s11814-015-0033-y.

    [5]Tsuji T, Higashida K, Okuyama Y, et al. Fictitious particle method: A numerical model for flows including dense solids with large size difference[J].AIChEJournal, 2014, 60(5): 1606-1620. DOI:10.1002/aic.14355.

    [6]Hilton J E, Mason L R, Cleary P W. Dynamics of gas-solid fluidised beds with non-spherical particle geometry[J].ChemicalEngineeringScience, 2010, 65(5): 1584-1596. DOI:10.1016/j.ces.2009.10.028.

    [7]Zhong W Q, Zhang Y, Jin B S, et al. Discrete element method simulation of cylinder-shaped particle flow in a gas-solid fluidized bed[J].ChemicalEngineering&Technology, 2009, 32(3): 386-391. DOI: 10.1002/ceat.200800516.

    [8]Ren B, Zhong W, Chen Y, et al. CFD-DEM simulation of spouting of corn-shaped particles[J].Particuology, 2012, 10(5): 562-572. DOI:10.1016/j.partic.2012.03.011.

    [9]Ren B, Zhong W, Jiang X, et al. Numerical simulation of spouting of cylindroid particles in a spouted bed[J].TheCanadianJournalofChemicalEngineering, 2014, 92(5): 928-934. DOI:10.1002/cjce.21900.

    [10]盧洲,劉雪東,潘兵.基于CFD-DEM方法的柱狀顆粒在彎管中輸送過程的數(shù)值模擬[J].中國粉體技術,2011,17(5):65-69.DOI:10.3969/j.issn.1008-5548.2011.05.017. Lu Zhou, Liu Xuedong, Pan Bing. Numerical simulation of the process of cylindrical particles in elbow pipe based on CFD-DEM method [J].ChinaPowderTechnology, 2011, 17(5): 65-69. DOI:10.3969/j.issn.1008-5548.2011.05.017. (in Chinese)

    [11]Oschmann T, Hold J, Kruggel-Emden H. Numerical investigation of mixing and orientation of non-spherical particles in a model type fluidized bed[J].PowderTechnology, 2014, 258: 304-323. DOI: 10.1016/j.powtec.2014.03.046.

    [12]Yamamoto Y, PotthoffM, Tanaka T, et al. Large-eddy simulation of turbulent gas-particle flow in a vertical channel: Effect of considering inter-particle collisions[J].JournalofFluidMechanics, 2001, 442: 303-334.DOI:10.1017/s0022112001005092.

    [13]Koukouvinis P, Gavaises M, Li J, et al. Large eddy simulation of diesel injector including cavitation effects and correlation to erosion damage[J].Fuel, 2016, 175: 26-39. DOI: 10.1016/j.fuel.2016.02.037.

    [14]Yin B, Yu S, Jia H, et al. Numerical research of diesel spray and atomization coupled cavitation by large eddy simulation (LES) under high injection pressure[J].InternationalJournalofHeatandFluidFlow, 2016, 59: 1-9. DOI: 10.1016/j.ijheatfluidflow.2016.01.005.

    [15]Hoover K A, Foley M G, Heasler P G, et al. Sub-grid-scale characterization of channel lengths for use in catchment modeling[J].WaterResourcesResearch, 1991, 27(11): 2865-2873. DOI: 10.1029/91WR01551.

    [16]Li T, Garg R, Galvin J, et al. Open-source MFIX-DEM software for gas-solids flows: Part Ⅱ—Validation studies[J].PowderTechnology, 2012, 220(SI): 138-150. DOI: 10.1016/j.powtec.2011.09.020.

    [17]Garg R, Galvin J, Li T, et al. Open-source MFIX-DEM software for gas-solids flows: Part Ⅰ—Verification studies[J].PowderTechnology, 2012, 220(SI): 122-137. DOI: 10.1016/j.powtec.2011.09.019.

    Numerical simulation based on sphere reconstruction and large eddy simulation coupled parallel algorithm

    Guo Feihong1Wang Zeyu1Zhong Zhaoping1Wang Xiaoyi2

    (1School of Energy and Environment, Southeast University, Nanjing 210096, China)(2Architectural Design and Research Institute Co., Ltd., Southeast University, Nanjing 210096, China)

    To simulate the mixing flow of cylindrical particles and quartz sand, the spherical element was used to reconstruct cylindrical particles. A sub-grid scale (SGS) model was used to analyze the turbulent flow. The coupled parallel algorithm in simulation reduced the computation time. Based on the above method, the mixing flow process of cylindrical particles and the quartz sand at 1.5 m/s was simulated. The simulation results were in good agreement with the experimental results on instantaneous flow and pressure fluctuation. The result shows that the proposed method can simulate the mixing flow of cylindrical particles and the quartz sand in fluidized bed. When the size of cylindrical particles is larger, the influence of the turbulence on mixing flow is aggravated. The average velocity of gas in thezdirection is simulated, presenting symmetrical distribution of high in the middle and low on both sides. This is related to the distribution of bubbles and cylindrical particles in fluidized bed, thus it is in line with the actual fluidization law.

    cylindrical particles; quarts sand; simulation; spherical element reconstruction; large eddy simulation; parallel algorithm

    10.3969/j.issn.1001-0505.2017.02.015

    2016-07-25. 作者簡介: 郭飛宏(1986—),男,博士生;仲兆平(聯(lián)系人),男,博士,教授,博士生導師,zzhong@seu.edu.cn.

    國家自然科學基金資助項目(U13161115,51276040)、東南大學優(yōu)秀博士學位論文培育基金資助項目(YBJJ1644).

    郭飛宏,王澤宇,仲兆平,等.基于球元重建與大渦模擬耦合并行算法的數(shù)值模擬[J].東南大學學報(自然科學版),2017,47(2):283-290.

    10.3969/j.issn.1001-0505.2017.02.015.

    TQ051

    A

    1001-0505(2017)02-0283-08

    猜你喜歡
    曳力大渦柱形
    預測天然氣斜井臨界攜液流量新方法
    循環(huán)流化床鍋爐爐膛流動特性數(shù)值模擬進展
    煤炭學報(2021年10期)2021-11-10 03:25:22
    基于壁面射流的下?lián)舯┝鞣欠€(wěn)態(tài)風場大渦模擬
    非柱形容器的壓力和壓強
    軸流風機葉尖泄漏流動的大渦模擬
    風機技術(2018年1期)2018-03-21 09:15:10
    從“柱形凸透鏡成像”實驗感受體驗教學的魅力
    基于EMMS模型的攪拌釜內氣液兩相流數(shù)值模擬
    化工學報(2016年7期)2016-08-06 07:11:55
    基于大渦模擬的旋風分離器錐體結構影響研究
    溪洛渡水電站GIL柱形絕緣子局部放電原因分析
    “柱形”鋁內襯纖維纏繞復合材料氣瓶自緊分析
    火箭推進(2014年3期)2014-03-16 05:49:00
    又黄又粗又硬又大视频| 美女视频免费永久观看网站| 精品国内亚洲2022精品成人 | 菩萨蛮人人尽说江南好唐韦庄| 国产三级黄色录像| 老司机靠b影院| 黄色怎么调成土黄色| 欧美日韩一级在线毛片| 午夜福利免费观看在线| 亚洲国产毛片av蜜桃av| 97精品久久久久久久久久精品| 亚洲欧洲日产国产| 亚洲av国产av综合av卡| 日韩欧美免费精品| 久久天堂一区二区三区四区| 性色av一级| 久久国产精品人妻蜜桃| 嫁个100分男人电影在线观看| 国产一区二区三区综合在线观看| 亚洲欧美日韩另类电影网站| 青草久久国产| kizo精华| 久久人人爽人人片av| 黑人欧美特级aaaaaa片| 久久天躁狠狠躁夜夜2o2o| 国产成人精品在线电影| 午夜福利视频在线观看免费| 精品少妇一区二区三区视频日本电影| 精品国产乱码久久久久久小说| 日韩视频一区二区在线观看| h视频一区二区三区| 中文字幕制服av| 丰满少妇做爰视频| 一级片'在线观看视频| cao死你这个sao货| 欧美黄色片欧美黄色片| 飞空精品影院首页| 欧美人与性动交α欧美软件| 爱豆传媒免费全集在线观看| 黄色视频不卡| 午夜影院在线不卡| 欧美日韩av久久| 亚洲欧美日韩高清在线视频 | 欧美黑人精品巨大| 日韩 亚洲 欧美在线| 成人三级做爰电影| 99国产综合亚洲精品| 国产野战对白在线观看| 亚洲专区国产一区二区| 亚洲av片天天在线观看| 男女国产视频网站| 日本一区二区免费在线视频| 亚洲一码二码三码区别大吗| 精品人妻在线不人妻| 少妇粗大呻吟视频| 超碰成人久久| 50天的宝宝边吃奶边哭怎么回事| 国产麻豆69| a在线观看视频网站| 日韩熟女老妇一区二区性免费视频| 久久精品aⅴ一区二区三区四区| 999久久久精品免费观看国产| 一本综合久久免费| 欧美中文综合在线视频| 精品无人区乱码1区二区| 美女午夜性视频免费| 欧美+亚洲+日韩+国产| 国产黄片美女视频| 淫妇啪啪啪对白视频| 国产激情久久老熟女| 成人手机av| 中文字幕人妻丝袜一区二区| 久久人人精品亚洲av| 小说图片视频综合网站| 日韩欧美在线二视频| 国产亚洲欧美98| 国产一区二区在线观看日韩 | 欧美日韩乱码在线| 亚洲熟女毛片儿| 国产区一区二久久| 丰满的人妻完整版| 美女大奶头视频| 欧美乱色亚洲激情| 91老司机精品| www.精华液| 波多野结衣高清作品| av中文乱码字幕在线| 欧美三级亚洲精品| 搡老熟女国产l中国老女人| 成年人黄色毛片网站| 天堂动漫精品| 每晚都被弄得嗷嗷叫到高潮| 日日爽夜夜爽网站| 校园春色视频在线观看| x7x7x7水蜜桃| 99国产精品一区二区蜜桃av| 99久久精品热视频| av国产免费在线观看| 人人妻人人澡欧美一区二区| 麻豆久久精品国产亚洲av| 一进一出抽搐gif免费好疼| 99在线视频只有这里精品首页| 亚洲avbb在线观看| 久久久水蜜桃国产精品网| 国产欧美日韩一区二区三| 在线观看www视频免费| 日韩国内少妇激情av| 久久 成人 亚洲| 夜夜躁狠狠躁天天躁| 亚洲在线自拍视频| 一本精品99久久精品77| 伦理电影免费视频| 国产v大片淫在线免费观看| 狂野欧美白嫩少妇大欣赏| 黄片大片在线免费观看| 岛国在线免费视频观看| 丝袜美腿诱惑在线| 久久久水蜜桃国产精品网| 国产精品九九99| 天堂动漫精品| 天堂√8在线中文| 亚洲精品久久成人aⅴ小说| 国产精品香港三级国产av潘金莲| 精品久久久久久,| 日本五十路高清| 精品人妻1区二区| 日韩欧美 国产精品| 日韩欧美国产一区二区入口| 精品不卡国产一区二区三区| 亚洲欧美精品综合一区二区三区| 国产激情偷乱视频一区二区| 亚洲av中文字字幕乱码综合| 看免费av毛片| 久久伊人香网站| 欧美日韩亚洲国产一区二区在线观看| 国产区一区二久久| 可以在线观看的亚洲视频| 精品欧美国产一区二区三| 日韩有码中文字幕| 午夜精品久久久久久毛片777| 亚洲成av人片免费观看| 成人手机av| a在线观看视频网站| 男人的好看免费观看在线视频 | 男男h啪啪无遮挡| 免费在线观看日本一区| 亚洲av成人精品一区久久| www.熟女人妻精品国产| 亚洲第一电影网av| 黑人欧美特级aaaaaa片| 人妻丰满熟妇av一区二区三区| 淫妇啪啪啪对白视频| 天堂影院成人在线观看| 久久久久精品国产欧美久久久| 最近最新免费中文字幕在线| 一个人免费在线观看的高清视频| 国产午夜精品久久久久久| 91麻豆精品激情在线观看国产| 久久这里只有精品中国| 日韩精品免费视频一区二区三区| 日韩成人在线观看一区二区三区| 日韩欧美国产一区二区入口| av天堂在线播放| 午夜影院日韩av| 精品久久久久久,| 97人妻精品一区二区三区麻豆| 欧美在线黄色| 亚洲精品久久成人aⅴ小说| 一区福利在线观看| 国产主播在线观看一区二区| 色噜噜av男人的天堂激情| 国产精品av视频在线免费观看| 99久久精品热视频| 中文字幕人妻丝袜一区二区| av有码第一页| 久久中文字幕一级| 久久久久久久午夜电影| 老汉色av国产亚洲站长工具| 亚洲av日韩精品久久久久久密| 变态另类丝袜制服| 免费看美女性在线毛片视频| 亚洲精品中文字幕一二三四区| 久久久久久久精品吃奶| 高清毛片免费观看视频网站| a在线观看视频网站| 后天国语完整版免费观看| 日韩大码丰满熟妇| 一级作爱视频免费观看| 两人在一起打扑克的视频| 亚洲男人的天堂狠狠| 91老司机精品| 国产精品亚洲美女久久久| 手机成人av网站| 久久香蕉精品热| 又粗又爽又猛毛片免费看| 最好的美女福利视频网| 巨乳人妻的诱惑在线观看| 久久精品aⅴ一区二区三区四区| 黄色女人牲交| 国产精品免费一区二区三区在线| 哪里可以看免费的av片| 长腿黑丝高跟| 欧美色视频一区免费| 九色国产91popny在线| 免费在线观看日本一区| 看黄色毛片网站| 成人永久免费在线观看视频| 日本成人三级电影网站| 免费看a级黄色片| 欧美午夜高清在线| 嫩草影院精品99| 男女之事视频高清在线观看| 午夜免费激情av| 亚洲免费av在线视频| 日韩欧美在线二视频| 久久 成人 亚洲| 手机成人av网站| 国产精品综合久久久久久久免费| 亚洲乱码一区二区免费版| 午夜免费观看网址| or卡值多少钱| 日韩高清综合在线| 欧美日韩福利视频一区二区| 可以免费在线观看a视频的电影网站| 久久久久久久久免费视频了| 搡老熟女国产l中国老女人| 熟女电影av网| 2021天堂中文幕一二区在线观| 亚洲av电影在线进入| 久久久久性生活片| 亚洲成av人片在线播放无| 搡老岳熟女国产| 一个人免费在线观看电影 | 搡老妇女老女人老熟妇| 久久久久国内视频| 精品一区二区三区四区五区乱码| 国产精品 欧美亚洲| 黄色丝袜av网址大全| 国产av麻豆久久久久久久| 国产精品98久久久久久宅男小说| 久久久久久人人人人人| 在线观看一区二区三区| 亚洲人成77777在线视频| 日本 欧美在线| 欧美成人免费av一区二区三区| 亚洲欧洲精品一区二区精品久久久| 久久人妻av系列| 麻豆一二三区av精品| 一本综合久久免费| 日韩国内少妇激情av| 中文资源天堂在线| 国产视频一区二区在线看| 在线a可以看的网站| 两个人的视频大全免费| 伊人久久大香线蕉亚洲五| 日韩欧美精品v在线| 久久99热这里只有精品18| 亚洲欧美日韩无卡精品| 在线永久观看黄色视频| 欧美午夜高清在线| 长腿黑丝高跟| av福利片在线观看| 欧美 亚洲 国产 日韩一| 亚洲av片天天在线观看| 无限看片的www在线观看| 欧美精品亚洲一区二区| 老熟妇乱子伦视频在线观看| 成年版毛片免费区| 动漫黄色视频在线观看| 精品国产亚洲在线| 国产三级在线视频| 亚洲国产精品999在线| 午夜免费观看网址| 51午夜福利影视在线观看| 午夜精品在线福利| 91麻豆av在线| 高潮久久久久久久久久久不卡| 亚洲欧美日韩东京热| 又粗又爽又猛毛片免费看| 久久婷婷人人爽人人干人人爱| 亚洲国产欧洲综合997久久,| 精品免费久久久久久久清纯| 国产精品美女特级片免费视频播放器 | 麻豆国产av国片精品| 精品日产1卡2卡| 一级毛片女人18水好多| 老熟妇乱子伦视频在线观看| 又粗又爽又猛毛片免费看| 香蕉av资源在线| 亚洲人成77777在线视频| 在线观看免费日韩欧美大片| 国产97色在线日韩免费| 怎么达到女性高潮| 欧美性长视频在线观看| 中文亚洲av片在线观看爽| 亚洲欧美日韩无卡精品| 丰满人妻一区二区三区视频av | 三级毛片av免费| 黑人欧美特级aaaaaa片| 久久 成人 亚洲| 黄色成人免费大全| 神马国产精品三级电影在线观看 | e午夜精品久久久久久久| 中文在线观看免费www的网站 | 久久精品国产清高在天天线| 淫秽高清视频在线观看| 亚洲五月婷婷丁香| 中文在线观看免费www的网站 | 午夜影院日韩av| 国产精品九九99| 制服诱惑二区| 欧美大码av| 亚洲欧美精品综合一区二区三区| 巨乳人妻的诱惑在线观看| 国产精品爽爽va在线观看网站| 一本一本综合久久| 久久伊人香网站| 又黄又爽又免费观看的视频| 亚洲精品一区av在线观看| 床上黄色一级片| 青草久久国产| 成人亚洲精品av一区二区| netflix在线观看网站| 91字幕亚洲| 亚洲熟妇熟女久久| 在线观看免费日韩欧美大片| 最新美女视频免费是黄的| 国产视频一区二区在线看| avwww免费| 亚洲黑人精品在线| 欧美日韩一级在线毛片| 免费看十八禁软件| 欧美日韩黄片免| 黄色视频,在线免费观看| 免费在线观看黄色视频的| 亚洲中文日韩欧美视频| 人妻丰满熟妇av一区二区三区| 国产伦在线观看视频一区| 此物有八面人人有两片| 久久香蕉国产精品| 一本久久中文字幕| 正在播放国产对白刺激| 国产精品1区2区在线观看.| 五月伊人婷婷丁香| 久久久精品欧美日韩精品| 中文字幕精品亚洲无线码一区| 成人一区二区视频在线观看| 免费在线观看影片大全网站| 亚洲激情在线av| 母亲3免费完整高清在线观看| 成人三级做爰电影| 亚洲自偷自拍图片 自拍| 日韩精品青青久久久久久| 欧美黑人欧美精品刺激| 丝袜人妻中文字幕| 国内久久婷婷六月综合欲色啪| 老司机在亚洲福利影院| 丁香欧美五月| 性色av乱码一区二区三区2| 久久这里只有精品19| 亚洲 欧美一区二区三区| 身体一侧抽搐| 亚洲精品色激情综合| 亚洲精华国产精华精| 国内久久婷婷六月综合欲色啪| 日韩国内少妇激情av| 18禁黄网站禁片免费观看直播| 国产伦人伦偷精品视频| 女人被狂操c到高潮| netflix在线观看网站| 国产精品影院久久| АⅤ资源中文在线天堂| 在线观看日韩欧美| 啦啦啦韩国在线观看视频| 亚洲全国av大片| 久久香蕉国产精品| 国产区一区二久久| 日日夜夜操网爽| 女人高潮潮喷娇喘18禁视频| 最近视频中文字幕2019在线8| bbb黄色大片| 一本综合久久免费| 男女下面进入的视频免费午夜| 免费人成视频x8x8入口观看| 欧美中文综合在线视频| 久久香蕉精品热| 久久久水蜜桃国产精品网| 久久久国产成人免费| 久久久精品大字幕| 老熟妇仑乱视频hdxx| 午夜福利免费观看在线| 国产精品亚洲av一区麻豆| 日韩国内少妇激情av| 日韩 欧美 亚洲 中文字幕| 日韩高清综合在线| 看免费av毛片| 大型黄色视频在线免费观看| 欧美日韩中文字幕国产精品一区二区三区| 亚洲成a人片在线一区二区| 欧美乱色亚洲激情| 日韩精品免费视频一区二区三区| 可以在线观看毛片的网站| 91字幕亚洲| 欧美性猛交╳xxx乱大交人| 日本撒尿小便嘘嘘汇集6| 国产免费男女视频| 99re在线观看精品视频| 老鸭窝网址在线观看| 国产av一区在线观看免费| 精品一区二区三区视频在线观看免费| 国产一区二区激情短视频| 首页视频小说图片口味搜索| 精品福利观看| 国产精品久久久久久精品电影| 禁无遮挡网站| 亚洲av片天天在线观看| 久久国产乱子伦精品免费另类| 欧美在线一区亚洲| 中文在线观看免费www的网站 | 少妇裸体淫交视频免费看高清 | 日韩欧美在线乱码| 日韩成人在线观看一区二区三区| 亚洲一码二码三码区别大吗| 精品久久久久久久久久久久久| 欧美成人午夜精品| 久久久水蜜桃国产精品网| 国产成人一区二区三区免费视频网站| 搡老熟女国产l中国老女人| 中文字幕人成人乱码亚洲影| 在线观看午夜福利视频| 国产精品电影一区二区三区| 国产欧美日韩一区二区三| e午夜精品久久久久久久| av免费在线观看网站| 久久亚洲精品不卡| 欧美成人性av电影在线观看| 国产精品永久免费网站| 999精品在线视频| 国产精品精品国产色婷婷| 亚洲国产精品久久男人天堂| 色哟哟哟哟哟哟| 久久人妻福利社区极品人妻图片| www.熟女人妻精品国产| 久久久久久九九精品二区国产 | 级片在线观看| ponron亚洲| cao死你这个sao货| 成人18禁在线播放| 在线观看www视频免费| 90打野战视频偷拍视频| 中出人妻视频一区二区| 欧美成人午夜精品| 欧美一区二区国产精品久久精品 | 精品久久久久久久久久免费视频| 91字幕亚洲| 免费人成视频x8x8入口观看| 中文字幕久久专区| 一级a爱片免费观看的视频| 这个男人来自地球电影免费观看| 久久久久久久久中文| 国产精品久久久av美女十八| 亚洲精品久久国产高清桃花| 亚洲,欧美精品.| 美女午夜性视频免费| 青草久久国产| 搡老熟女国产l中国老女人| 免费av毛片视频| 一级作爱视频免费观看| 国产精品亚洲美女久久久| 日本精品一区二区三区蜜桃| 久久草成人影院| 亚洲成人久久性| cao死你这个sao货| 中文字幕久久专区| 18禁黄网站禁片免费观看直播| 亚洲成av人片免费观看| 特级一级黄色大片| 99在线人妻在线中文字幕| 免费看日本二区| 黄色 视频免费看| 欧美日韩乱码在线| 一级毛片女人18水好多| 99re在线观看精品视频| av福利片在线| 欧美日韩瑟瑟在线播放| 亚洲18禁久久av| 91麻豆av在线| 精品乱码久久久久久99久播| 中文字幕精品亚洲无线码一区| 国产黄片美女视频| or卡值多少钱| 欧美久久黑人一区二区| 久久国产乱子伦精品免费另类| 久久午夜亚洲精品久久| 亚洲欧美日韩东京热| 亚洲无线在线观看| 亚洲人成网站高清观看| 久久精品夜夜夜夜夜久久蜜豆 | 亚洲在线自拍视频| 日本 av在线| 亚洲欧美激情综合另类| 欧美三级亚洲精品| 亚洲aⅴ乱码一区二区在线播放 | 国产午夜精品久久久久久| 国产精品永久免费网站| www日本黄色视频网| 女警被强在线播放| 国产午夜福利久久久久久| 岛国在线免费视频观看| 91九色精品人成在线观看| 床上黄色一级片| 91大片在线观看| 国产精品香港三级国产av潘金莲| 99久久国产精品久久久| 不卡一级毛片| 国模一区二区三区四区视频 | 又黄又粗又硬又大视频| 亚洲中文av在线| 中文字幕久久专区| 美女黄网站色视频| 国产精品久久电影中文字幕| 免费看a级黄色片| 欧美日韩福利视频一区二区| 亚洲电影在线观看av| 中亚洲国语对白在线视频| 不卡av一区二区三区| 国产精品精品国产色婷婷| 黄色毛片三级朝国网站| 88av欧美| 欧美大码av| 香蕉av资源在线| 精品久久久久久成人av| 人人妻人人澡欧美一区二区| 18禁黄网站禁片午夜丰满| 我的老师免费观看完整版| 午夜福利视频1000在线观看| 亚洲精品一卡2卡三卡4卡5卡| 精品免费久久久久久久清纯| 久久精品国产综合久久久| 亚洲av电影在线进入| 日日夜夜操网爽| 在线免费观看的www视频| 欧洲精品卡2卡3卡4卡5卡区| 巨乳人妻的诱惑在线观看| e午夜精品久久久久久久| 日本免费一区二区三区高清不卡| 国产亚洲av嫩草精品影院| 亚洲欧美精品综合一区二区三区| 亚洲精品中文字幕一二三四区| 91九色精品人成在线观看| 九色成人免费人妻av| 91老司机精品| 免费电影在线观看免费观看| 美女扒开内裤让男人捅视频| 亚洲国产精品sss在线观看| 99久久精品国产亚洲精品| 日本免费一区二区三区高清不卡| 两性夫妻黄色片| 欧美成人午夜精品| 亚洲五月天丁香| 天天一区二区日本电影三级| 免费看美女性在线毛片视频| 亚洲午夜理论影院| 欧美一区二区国产精品久久精品 | 不卡一级毛片| 禁无遮挡网站| 三级国产精品欧美在线观看 | 国产高清有码在线观看视频 | 黄色 视频免费看| 国产亚洲av嫩草精品影院| 欧美成人午夜精品| 免费人成视频x8x8入口观看| 亚洲最大成人中文| 亚洲av第一区精品v没综合| 午夜激情av网站| 久久性视频一级片| 男人舔女人的私密视频| 91麻豆av在线| 国产成年人精品一区二区| 亚洲国产欧美网| 最近视频中文字幕2019在线8| 欧美高清成人免费视频www| 午夜福利免费观看在线| 色综合婷婷激情| 亚洲av电影在线进入| 一级毛片高清免费大全| 久久久久久久午夜电影| 久久 成人 亚洲| 色综合站精品国产| 亚洲国产欧美网| 丝袜美腿诱惑在线| 亚洲自偷自拍图片 自拍| 黄频高清免费视频| 亚洲欧美日韩东京热| 国产精品爽爽va在线观看网站| 亚洲av成人一区二区三| 亚洲性夜色夜夜综合| 一本大道久久a久久精品| 黄频高清免费视频| 欧美激情久久久久久爽电影| 夜夜爽天天搞| 欧美成人午夜精品| 毛片女人毛片| 免费在线观看影片大全网站| av福利片在线| 神马国产精品三级电影在线观看 | 久久婷婷成人综合色麻豆| 欧美日韩一级在线毛片| 97碰自拍视频| 首页视频小说图片口味搜索| 美女午夜性视频免费| 法律面前人人平等表现在哪些方面| 日本 av在线| 九色国产91popny在线| 亚洲精品美女久久av网站|