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

    兩個非磁性顆粒在磁流體中的沉降現(xiàn)象研究?

    2017-09-07 20:55:16陳木鳳李翔牛小東李游Adnan山口博司
    物理學(xué)報 2017年16期
    關(guān)鍵詞:磁流體磁化比值

    陳木鳳 李翔?牛小東 李游 Adnan 山口博司

    1)(汕頭大學(xué)工學(xué)院,汕頭 515063)

    2)(同志社大學(xué)能源轉(zhuǎn)換與研究中心,京都630-0321,日本)

    兩個非磁性顆粒在磁流體中的沉降現(xiàn)象研究?

    陳木鳳1)李翔1)?牛小東1)李游1)Adnan1)山口博司2)

    1)(汕頭大學(xué)工學(xué)院,汕頭 515063)

    2)(同志社大學(xué)能源轉(zhuǎn)換與研究中心,京都630-0321,日本)

    (2017年3月31日收到;2017年6月2日收到修改稿)

    在磁場作用下,在磁流體里添加非磁性顆粒(non-magnetic particles,NPs),可以使得NPs形成不同的結(jié)構(gòu),操控NPs的運動從而影響磁流體的特性,這種應(yīng)用逐漸受到了研究者的關(guān)注.為了更好地操控磁流體里NPs的運動,本文采用一種多物理模型研究在外加磁場作用下,磁流體中兩個NPs沉降的運動過程.其中,用格子玻爾茲曼方法模擬磁流體的運動,外加磁場對磁流體的影響用一種自修正方法求解泊松方程,這個自修正方法可以使歐姆定律滿足守恒定律.NPs之間的偶極干擾力采用偶極力模型,同時采用一種相對過渡平滑的共軛邊界條件處理NPs與磁流體交界面的流固干擾以避免磁場密度過渡的突變.本文主要探究兩個NPs在磁流體中的沉降,揭示磁場作用下NPs的相互干擾原理;同時,對控制NPs運動時的參數(shù)進行調(diào)節(jié),得到NPs不同的運動軌跡,達到操控顆粒運動的目的.本研究可對NPs在磁流體中的應(yīng)用提供定量的分析結(jié)果,對NPs在工業(yè)上的應(yīng)用提供有力的理論支撐.

    磁流體,非磁性顆粒,自修正方法,共軛邊界條件

    1 引 言

    磁流體是由納米級磁性顆粒(magnetic particles,MPs,直徑約10 nm)基載液和分散劑三者混合而成的一種穩(wěn)定的溶液.與一般流體相比,磁流體不僅具有液體的流動性,而且具有磁化性能,利用磁流體的磁化特性,可以通過外加磁場實現(xiàn)對磁流體運動的操控.迄今磁流體的應(yīng)用越來越廣泛,包括磁流體密封、減震、潤滑、醫(yī)療器械、聲音調(diào)節(jié)、傳熱散熱、磁流體選礦等領(lǐng)域.含非磁性顆粒(NPs)的磁流體被稱為反磁流體(inversemagnetic fluid).這是因為NPs的尺寸遠大于磁流體中的納米級磁性顆粒,NPs與磁流體之間的相互作用可以看作是固相顆粒與牛頓流體之間的流固耦合作用;外加磁場時,NPs被附近磁流體反向磁化并表現(xiàn)出各向異性.在磁流體中放置大量NPs時,NPs之間由于磁矩而具有偶極力,這種各向異性使得NPs在磁場方向組裝成鏈狀.由于磁流體的宏觀特性是由磁流體的微觀結(jié)構(gòu)決定的,在磁流體中加入NPs在物理、化學(xué)、生物工程等領(lǐng)域具有重要的工程應(yīng)用價值.這些應(yīng)用包括轉(zhuǎn)軸的密封、揚聲器線圈的散熱、晶體工程、拋光和研磨、提高導(dǎo)熱性能等[1?7].

    很多研究者通過實驗或者理論模擬研究NPs在磁流體中的微觀結(jié)構(gòu)形態(tài)[5?14].其中Kaiser等[8]在1976年首先提出,在磁化磁流體中,沿磁力線方向放置的NPs會相互吸引,垂直于磁力線方向上放置的NPs會相互排斥.Fu jita和MaMiya[10]在1987年通過實驗測量和偶極力模型計算兩種方法,獲得兩個NPs的相互吸引力,兩種方法的結(jié)果相互符合,再次驗證了Kaiser[8]現(xiàn)象.Skjeltorp[9]在1983年提出磁性液體中加入非磁性膠體混合而成的液體成為反磁性液體.他將磁流體中的NPs稱為磁空穴(magnetic holes),并研究了在含聚苯乙烯的磁流體薄片中外加磁場的作用,發(fā)現(xiàn)摻入到磁性液體中的聚苯乙烯小球可以結(jié)晶成不同的格子.

    Zhu等[3]用Stokesian動力學(xué)方法模擬了在磁流體中加入微米級圓形MPs和圓形NPs,結(jié)果表明NPs可以被MPs驅(qū)趕,并最終分布在磁性顆粒團的空隙里,這樣可以實現(xiàn)對NPs的操控.Ido等[4,5]使用同樣的方法,混合了不同比例的球形和膠囊形MPs和NPs研究顆粒的自組裝形態(tài),發(fā)現(xiàn)NPs在磁力線方向上也會發(fā)生輕微的自組裝現(xiàn)象.通過理論上分析NPs受到的磁感應(yīng)力,Gao等[12]發(fā)現(xiàn)NPs在磁流體中受到的磁感應(yīng)力可以通過選擇合適的格子大小、NPs的體積分數(shù)和幾何形狀等調(diào)節(jié).Peng等[13]采用二維蒙特卡羅法模擬均勻磁場作用下的磁流體中懸浮的MPs和NPs,他們發(fā)現(xiàn)顆粒本身運動受顆粒本身的屬性影響大于外磁場的影響;同時,NPs會阻礙磁性顆粒的運動,并且NPs越多,阻礙作用越明顯.Iwamoto等[7]最近通過實驗發(fā)現(xiàn)在磁流體中的銀納米線可以改變磁流體的熱導(dǎo)率,具體表現(xiàn)為當(dāng)磁場平行于熱梯度方向時熱導(dǎo)率提高,垂直時熱導(dǎo)率減少.

    NPs在磁流體中的反向磁化特性不僅豐富了磁流體的研究領(lǐng)域,而且為實現(xiàn)NPs的磁操控開拓了新的思路,使得磁流體的應(yīng)用更加廣泛.但是,目前為止,對磁流體非磁性顆粒的研究大多數(shù)都是定性的分析,很少得到定量驗證和分析.為了使反磁性液體的特性在工業(yè)上得到更加合適的應(yīng)用,本文主要研究兩個NPs在磁流體中的沉降運動,揭示磁場作用下NPs的相互干擾原理;同時對顆粒本身的特性以及外界磁場運動時的參數(shù)進行調(diào)節(jié),得到顆粒的不同的運動軌跡,達到調(diào)節(jié)顆粒運動的目的.這可為調(diào)節(jié)參數(shù)操控NPs運動提供定量的分析結(jié)果,對NPs在工業(yè)上的應(yīng)用提供有力的理論支撐.

    在本研究中,我們提出一種多物理模型研究非磁性顆粒在均勻磁化磁流體中的自組裝現(xiàn)象.磁流體的運動采用介觀格子玻爾茲曼(LBM)方法模擬[15?21];外加磁場用一種自修正方法[22,23]求解泊松方程,其中磁場泊松方程從麥克斯韋磁場方程推導(dǎo)而來,自修正過程可以使歐姆定律滿足守恒定律.對于固體顆粒和磁流體的流固耦合,我們采用最近提出的分布函數(shù)修正的浸沒邊界LBM方法[21].在磁流體中,非磁性顆粒受到的力除了慣性力、流固干擾力、顆粒之間的接觸力(近距離時)、流場邊界碰撞力(與流場邊界近距離時),還受到其他顆粒的偶極力.偶極力是由于顆粒附近磁流體對NPs的磁化而產(chǎn)生的,采用偶極力模型計算[10].另外,由于磁流體與非磁性顆粒之間的磁透率不一樣,這使得顆粒附近的磁感應(yīng)強度不一樣,我們在邊界附近采用一種相對過渡平滑的共軛邊界條件[24?26]避免磁場密度過渡的突變.這些方法使得處理流固耦合的邊界更加精確.與文獻[24—26]的共軛邊界解法不同,本文把共軛邊界結(jié)合磁勢的泊松方程,并采用直接離散的方法求解共軛邊界條件.

    以下本文的結(jié)構(gòu)分為兩大部分:第2部分是對模擬方法進行闡述;第3部分對兩個NPs的吸引力算例進行數(shù)值驗證[10],首先驗證了不加磁場的情況下兩個NPs沉降[17,27?29]的算例,接著加入磁場,改變磁場與顆粒的磁透率比值和磁場參數(shù),探究顆粒的不同運動軌跡.本研究可為實現(xiàn)對磁流體中的NPs進行精準操控提供有力的理論支撐.

    2 模擬方法

    2.1 基于分布函數(shù)修正浸沒邊界的LBM

    單松弛時間的基于分布函數(shù)修正浸沒邊界的LBM方程[17,21]如下:

    其中fα(r,t)是密度分布函數(shù);eα,r,δt,τf,ωα分別為離散格子速度、空間位置矢量、時間步長、松弛時間和權(quán)重函數(shù);ωα=4/9(α =0),1/9(α=1,2,3,4),1/36(α=5,6,7,8);Fα為外力項,(r,t)是由固體邊界產(chǎn)生的擾動反應(yīng)到流場中[21].平衡態(tài)密度分布函數(shù)fαeq(r,t)依據(jù)D2Q9離散速度模型[17]表示如下;

    其中ρ,u分別是無量綱密度和速度;δx=1/N是格子單元(N為格子數(shù)目),c=δx/δt=1是單位格子速度,且c=cscs為格子聲速).在方程(1)中的離散外力項可以定義如下;

    FM是Kelvin力,來源于外加磁場對磁流體的影響,可以表示為FM(X,t)=M?B,其中M,B分別為磁化強度和磁感應(yīng)強度.它的無量綱形式為

    我們定義磁雷諾數(shù)RM=χμ0H20/(ρc2),它是由磁流體和外加磁場強度共同決定的,其中χ,μ0分別為是磁流體的磁化率和真空磁透率;H為磁化強度矢量,H0和分別為它的模和無量綱矢量.

    在玻爾茲曼方程(1)中,介觀虛擬顆粒在流場中經(jīng)歷了碰撞和遷移兩個演化過程,并且在每個時間步長內(nèi)滿足動量守恒定律.宏觀密度和動量是由密度分布函數(shù)求得的,具體如下:

    其中無量綱松弛時間τf如下:

    其中η為流體的動力黏度.

    2.2 一種磁場標量勢方程的自修正過程

    計算磁流體流場的磁流體動力學(xué)采用一種磁場麥克斯韋方程的自修正過程[22,23].與壓力泊松方程質(zhì)量守恒定律類似,自修正過程可以使得磁場標量勢的泊松方程中的磁場強度滿足守恒定律.根據(jù)文獻[22],基于歐姆定律和磁標量勢泊松方程的自修正過程如下:

    μf和?分別為磁流體的磁透率和磁標量勢,其中H=???.ωB是欠松弛因子(0<ωB<1,本文取ωB=0.04).根據(jù)靜磁場的麥克斯韋方程,有

    J為電流密度.另外,由LBM計算的是整個流場,磁流體和固體的磁透率不一樣會影響流固耦合邊界周圍的磁感應(yīng)強度,所以這里我們采用共軛邊界條件處理邊界.具體如下:

    其中,μnp是NPs的磁透率.對于非磁性顆粒,有μnp<μf,這個共軛邊界的實施可以使得非磁性顆粒與磁流體的邊界過度相對平滑,避免突變引起流場發(fā)散.方程(8)的離散形式為

    磁流體在外加磁場作用下,磁場感應(yīng)強度可以由導(dǎo)磁物質(zhì)內(nèi)的磁場公式求得,即B=B0+B′=B0+μfχH,其中B0是外部磁化場,B′是磁流體中導(dǎo)磁物質(zhì)(納米顆粒)和非磁性顆粒的貢獻.所以對于NPs,磁矩方向與外加磁場方向相反,表現(xiàn)出反磁性.

    2.3 NPs在磁流體中受到的力

    在本研究中,顆粒在流體受到的合力(Ftotal)[17]包括流固干擾力(F(Xl,t)),顆粒之間的接觸力(Fp-p)(采用Lennard-Jones(L-J)勢模型計算顆粒間的相互作用力),顆粒與固體邊界的接觸力(Fp-w)(與邊界近距離時,避免穿透邊界),慣性力(FI)和顆粒由于被磁化而對其他顆粒產(chǎn)生的偶極力(FM).第i個顆粒受到的合力可以通過以下公式求得:

    這里的Mf是體積為一個顆粒大小的流體的質(zhì)量;i和j分別代表第i,j個顆粒;ri和rj是顆粒位置坐標;a是顆粒的半徑,并且ε=a2,ri,j=|ri?rj|,ri,w=|ri?rwj|,rwj是顆粒離流場邊界的距離;l與格子長度相等.

    與文獻[17]相比,NPs在外加磁場作用下,還受到被周圍的磁流體磁化而與其他顆粒產(chǎn)生相互偶極干擾力[4,5,10],偶極力可以由偶極力模型[10]求得.第i個顆粒受到的偶極力是其他NPs對它產(chǎn)生的偶極力總和,如下:

    其中,Mi代表第i個顆粒的磁矩,它是由磁流體和其他非磁性顆粒共同磁化的結(jié)果[11,13,14],可以由以下公式求得:

    3 數(shù)值模擬結(jié)果與分析

    首先數(shù)值驗證兩個在磁場方向平行排列的圓形NPs在磁化磁流體中的最大吸引力并與文獻[10]比較.接著在管道兩個NPs沉降[17,27?29]的算例中加磁場,探究顆粒軌跡及其受力的變化.

    3.1 算例驗證——兩個NPs的吸引力驗證

    為了驗證上述提出的混合LBM、偶極力模型和自修正方法的多物理模型,本文數(shù)值驗證兩個在磁場方向平行排列的圓形NPs在磁化磁流體中的最大吸引力,并與文獻[9]中的實驗測量結(jié)果和數(shù)值結(jié)果進行比較.

    圖1 兩個NPs的吸引力計算模型Fig.1.Calcu lation Model force attracted force of 2 NPs.

    圖2 (網(wǎng)刊彩色)不同磁化強度下兩個NPs的吸引力Fig.2.(color on line)Attracted force between 2 NPs in d iff erent Magnetic intensity.

    流場區(qū)域為Lx×Ly=2×4的方腔,如圖1所示,網(wǎng)格為100×200.非磁性顆粒為銅球,半徑為r=0.555/2,初始位置為 P1(1.0,2.0+1.1r),P2(1.001,2.0?1.1r).松弛時間為τf=0.575.本研究采用與文獻[10]相同條件的量綱形式進行數(shù)值模擬.我們?nèi)〔煌帕黧w密度比重(比重為1.2,1.3和1.4),以及每個比重下外加不同的磁化強度進行模擬,獲得兩個NPs的最大吸引力,并與文獻[10]比較.從圖2可以看出,在不同的磁化強度下,兩個NPs的最大吸引力在文獻[10]的計算結(jié)果范圍內(nèi)并且與實驗結(jié)果相符合,證明本文的混合物理模型有效可行.

    3.2 管道中兩個顆粒(NPs)沉降的算例中加磁場

    為了探究非磁性顆粒在磁流體中的運動行為,本文采用上述的多物理模型模擬了兩個圓形非磁性顆粒在垂直管道中的沉降問題[17,27?29],并加磁場研究在磁化磁流體中NPs沉降的軌跡、速度和受力.與經(jīng)典DKT(D rafting-Kissing-Tumbling)研究一樣,本文選取文獻[17]相同的條件,垂直管道為2 cm(寬)×8 cm(高),流體的密度為ρf=1 g/cm3,黏性為0.001 g/cm.s.非磁性顆粒的密度為ρnp=1.01 g/cm3,半徑為r=0.1 cm,顆粒的初始位置為P1(0.999 cm,7.2 cm,實心),P2(1.0 cm,6.8 cm,空心),速度為零.我們選取200×800的網(wǎng)格,每個顆粒有20個網(wǎng)格單元覆蓋,松弛時間為τf=0.65.顆粒受到的力如2.3部分公式所示.

    在本研究中,我們首先不加磁場,對顆粒下落算例進行數(shù)值驗證;接著加入磁場,改變顆粒與磁流體的磁透率比值,探究磁透率不同的非磁性顆粒下落的現(xiàn)象;最后改變RM值,探究RM值對顆粒下落軌跡的影響.同時,探究磁場對顆粒軌跡、速度和受力的影響.

    圖3展示了不加磁場作用下圓形顆粒下落過程中的幾個不同時刻(t =0,1.0,1.5,2.0,2.5,3.0,4.5 s)的位置變化現(xiàn)象.首先,兩個顆粒同時下落,在下面的顆粒P2作為引導(dǎo),為后面的顆粒P2創(chuàng)造一個低壓的尾跡;P2在這個尾跡的影響下,下落的速度加快;經(jīng)過一段時間的拖曳(drafting),P2趕上了P1并且兩個顆粒接觸碰撞(kissing).碰撞后,兩個顆粒為一個整體下落一段時間后,由于流體動力的作用,使得兩個顆粒分開并且各自翻轉(zhuǎn)下落(tumbling).其位置與文獻[17,27—29]的DTK過程相符合.

    圖3 不同時刻兩個顆粒下落的位置Fig.3. SediMentation of 2 NPs at diff erent tiMe stages.

    圖4、圖5和圖6、圖7分別是圓形顆粒在水平方向和垂直方向的位置和速度隨時間的變化.從圖4和圖5可以看出,在0.8 s前兩個顆粒中心幾乎保持相同的距離,之后顆粒P1在顆粒P2創(chuàng)造的尾跡的影響下追趕P2.當(dāng)兩個顆粒距離較近時,我們采用的是L-J勢模型計算顆粒間的相互作用力,避免顆粒之間相互穿透,與文獻[17,29]的結(jié)果一樣.在碰撞的過程中,兩個顆粒的中心逐漸偏離管道中心,然后做翻轉(zhuǎn)運動.可以清晰地看到本文的結(jié)果與文獻符合.由于選擇的模擬方法和碰撞模型不同,會對顆粒碰撞后的軌跡和速度有不同的影響,這些在位置和速度的變化上都有體現(xiàn)[29],由于速度波動的地方數(shù)據(jù)波動大,這里取的是文獻[29]采用L-J勢模型計算的速度平均值.垂直方向碰撞后速度的不穩(wěn)定與文獻[29]選的L-J勢模型計算結(jié)果符合得比較好.

    圖4 (網(wǎng)刊彩色)兩個顆粒下落水平方向的位置隨時間的變化Fig.4.(color on line)Transverse coordinates of the centers of 2 NPs.

    圖5 (網(wǎng)刊彩色)兩個顆粒下落的垂直方向的位置隨時間的變化Fig.5.(color on line)Longitudinal coordinates of the centers of 2 NPs.

    圖6 (網(wǎng)刊彩色)兩個顆粒下落水平方向的速度隨時間的變化Fig.6.(color on line)TeMporal variations With transverse velocity of 2 NPs.

    圖7 (網(wǎng)刊彩色)兩個顆粒下落垂直方向的速度隨時間的變化Fig.7.(color on line)TeMporal variations With longitud inal velocity of 2 NPs.

    接著我們加入外部磁場,探究外加磁場作用下磁流體中非磁性顆粒下落的軌跡和速度的變化,為非磁性顆粒組裝、結(jié)晶等工程問題探究提供新的視覺向?qū)Ш筒灰粯拥乃悸?首先選取參數(shù)RM=10?4,改變顆粒與磁流體的磁透率比值,探究磁透率不同的非磁性顆粒下落的現(xiàn)象.

    圖8展示了RM=10?4時,在磁透率比值不同的情況下,顆粒下落過程中的7個不同時刻(t=0,1.0,1.5,2.0,2.5,3.0,4.5 s)的位置.從圖8可以發(fā)現(xiàn),當(dāng)磁透率比值小于2.9時,兩個顆粒在拖曳、接觸(DK)后,往右邊翻轉(zhuǎn),當(dāng)比值大于3時,DK后往左邊翻轉(zhuǎn).其詳細的下落位置變化如圖9和10所示,可以看出顆粒DK后的翻轉(zhuǎn)位置.當(dāng)比值小于2.9時,顆粒DK后向右翻轉(zhuǎn),并且比值越小,越接近于無磁場時的運動軌跡;當(dāng)比值大于3時,顆粒DK后向左翻轉(zhuǎn),并且比值越大,DK后翻轉(zhuǎn)的時間越早.

    圖8 兩個顆粒下落不同時刻的位置Fig.8. SediMentation of 2 NPs at diff erent tiMe stages.

    圖9 (網(wǎng)刊彩色)磁場作用下兩個顆粒下落水平方向的位置隨時間的變化Fig.9.(color on line)Transverse coordinates of the centers of 2 NPs in Magnetic field.

    圖11和圖12是顆粒下落速度隨時間的變化情況.外加磁場后,顆粒在DK前后的速度變化波動較大.這是由于當(dāng)顆粒放在磁流體里,外加磁場后,磁流體的流體動力受到由外磁場引起的Kelvin力的影響,這個力與磁流體的磁化率和外加磁場的強弱有關(guān).當(dāng)外加磁場時,流體顆粒由于被磁化而帶磁矩,兩個顆粒相互靠近時,除了顆粒之間的接觸力,還受到偶極力的作用,所以受的影響是這些力的疊加.

    圖10 (網(wǎng)刊彩色)磁場作用下兩個顆粒下落垂直方向的位置隨時間的變化Fig.10.(color on line)Longitud inal coordinates of the centers of 2 NPs in Magnetic field.

    圖11 (網(wǎng)刊彩色)磁場作用下兩個顆粒下落水平方向的速度隨時間的變化Fig.11.(color online)TeMporal variationsWith transverse velocity of 2 NPs in Magnetic field.

    圖12 (網(wǎng)刊彩色)磁場作用下兩個顆粒垂直方向的速度隨時間的變化Fig.12.(color on line)TeMporal variations With longitudinal velocity of 2 NPs in Magnetic field.

    由于磁透率比值不一樣,會引起顆粒發(fā)生不一樣的翻轉(zhuǎn),我們分析顆粒水平方向的受力,探究其對顆粒運動的影響.為了分析磁場對顆粒DK后翻轉(zhuǎn)的軌跡變化,分析了顆粒在不加磁場、磁透率比值為2.9和3這三種情況下,顆粒x方向的受力,結(jié)果如圖13—圖15所示.

    圖13 (網(wǎng)刊彩色)無磁場時,兩個顆粒水平方向受力隨時間的變化Fig.13.(color on line)Transverse force of 2 NPs Without Magnetic field.

    圖14 (網(wǎng)刊彩色)μf/μnp=2.9時,兩個顆粒下落水平方向的受力隨時間的變化Fig.14.(color on line)Transverse force of 2 NPs atμf/μnp=2.9.

    圖15 (網(wǎng)刊彩色)μf/μnp=3時,兩個顆粒下落水平方向的受力隨時間的變化Fig.15.(color on line)Transverse force of 2 NPs atμf/μnp=3.

    當(dāng)不加磁場時,如圖13,顆粒在接觸(約t=1.5 s)之間受到的力逐漸由流體動力主導(dǎo)變?yōu)轭w粒之間的接觸力為主導(dǎo),由于顆粒之間的接觸力為排斥力,大小相等方向相反,并且與顆粒之間的中心距離密切相關(guān).當(dāng)顆粒P1在P2的尾跡引導(dǎo)下追趕上P2并且發(fā)生碰撞,以兩個顆粒為整體的系統(tǒng)動量守恒,所以顆粒間碰撞后動量交換,受力方向與原來的方向相反.碰撞后,距離變大;由于流體動力和慣性力的影響,顆粒又開始逐漸靠近,在下落的過程中發(fā)生了多次碰撞、分離、再次碰撞的運動.由于流體動力和慣性力的影響使得兩個顆粒在水平方向上速度先后變?yōu)榱?t=3 s和t=3.75),兩個顆粒分道揚鑣,繼續(xù)下落.

    當(dāng)磁透率比值較小時(μf/μnp=2.9),如圖14,顆粒受到的流體動力受磁場Kelvin力影響,并且顆粒由于被逆磁化而具有磁矩,所以顆粒間的受力包括顆粒偶極力.但是由于磁透率比值較小,所以這些影響比較微弱.顆粒在下落時,碰撞過程前后會受到微小偶極的相互吸引力,所以顆粒接觸后分離的時間(約t=3.02 s)延后.

    當(dāng)磁透率比值較大時(μf/μnp=3),顆粒更容易受到Kelvin力的影響,兩個顆粒的磁矩對磁場更敏感,所以兩個顆粒距離越來越近時,顆粒受到的偶極力的影響較大.兩個顆粒的偶極力為相互吸引力,部分抵消了由于顆粒接觸而產(chǎn)生的接觸力,見圖15.顆粒在碰撞后并沒有被強烈地排斥,而是分開一小段距離后又被吸引了.這樣的過程也會反復(fù)多次,直到顆粒P1在水平方向上的速度變?yōu)榱?約t=3.4 s),在重力、流體動力和偶極力的作用下各自下落.

    當(dāng)磁透率比值更大時(μf/μnp=10),顆粒受力的敏感源包括由Kelvin力影響的流體動力和顆粒間的偶極力,在二者的共同作用下,顆粒的碰撞過程比較急劇,速度波動較大,如圖11和圖12所示.

    最后改變RM值,探究RM值對顆粒下落軌跡的影響.圖16和圖17展示了不同RM值時顆粒下落的位置變化情況.我們選了四個RM值,分別是RM=10?3,5×10?4,10?4和10?5,發(fā)現(xiàn)顆粒的運動軌跡變化不大.具體表現(xiàn)為:在DK前軌跡幾乎重合,DK后翻轉(zhuǎn)的方向相同,軌跡略有偏差.這說明改變RM,會使顆粒的運動軌跡有細微的差別.

    在實際工程應(yīng)用中,對NPs運動軌跡進行宏觀操控時,可以通過改變改變顆粒的性質(zhì),使得磁透率比值不一樣,實現(xiàn)預(yù)想軌跡粗調(diào);接著改變RM值,對NPs的軌跡實現(xiàn)微調(diào),使得NPs的運動情況更貼近我們的預(yù)想軌跡,更好地調(diào)控顆粒的運動.

    圖16 (網(wǎng)刊彩色)不同RM值下兩個顆粒下落水平方向的位置隨時間的變化Fig.16.(color online)Transverse coordinates of the centers of 2 NPs at diff erent RM.

    圖17 (網(wǎng)刊彩色)不同RM值下兩個顆粒下落垂直方向的位置隨時間的變化Fig.17.(color on line)Longitud inal coordinates of the centers of 2 NPs at diff erent RM.

    4 結(jié) 論

    本文混合了浸沒邊界法的LBM方程、自我修正過程解泊松方程和共軛邊界條件,模擬磁化磁流體中的流固耦合現(xiàn)象.在外加磁場作用下,NPs被磁化,并帶有反向磁矩.研究表明,在探究兩個圓形NPs在磁化磁流體中下落時發(fā)現(xiàn),適當(dāng)改變NPs的磁透率,使得磁流體和NPs的磁透率比值不同,會得到不同的顆粒運動軌跡.這個比值的不同,使得NPs在碰撞后的轉(zhuǎn)向發(fā)生改變.具體表現(xiàn)為:當(dāng)比值較小時,顆粒轉(zhuǎn)向變化不大,隨著比值的增大,顆粒逐漸發(fā)生了相反的轉(zhuǎn)向.而改變RM值,對顆粒的運動軌跡影響較小,可以實現(xiàn)對顆粒的微觀調(diào)節(jié).這些發(fā)現(xiàn)對操控磁流體中NPs的運動,從而改變磁流體的特性有著重要的指導(dǎo)意義.

    [1]Halsey T C,Toor W1990 J.Stat.Phys.61 1257

    [2]IslaMMF,Lin K H,Lacoste D,Lubensky T C,Yodh A G 2003 Phys.Rev.E 67 021402

    [3]Zhu Y,UMehara N,Ido Y,Sato A 2006 J.Magn.Magn.Mater.302 96

    [4]Ido Y,Inagaki T,UMehara N 2008 MagnetohydrodynaMics 44 83

    [5]Ido Y,Inagaki T,YaMaguchiT 2010 J.Phys.:Condens.Matter 22 324103

    [6]Chen Q,Bae S C,G ranick S 2011 Nature 469 381

    [7]IwaMoto Y,Yoshioka A,Naito T,Cuya J,Ido Y,Okawa R,YaMaguchi H 2016 Exp.Therm.Fluid Sci.79 111

    [8]Kaiser R,Mir L,Curtis R A 1976 US Patent 3951785

    [9]Skjeltorp At 1983 Phys.Rev.Lett.51 2306

    [10]Fu jitAt,MaMiya M1987 J.Magn.Magn.Mater.65 207

    [11]Fu rst E M,Gast A P 2000 Phys.Rev.E 61 6732

    [12]Gao Y,Jian Y C,Zhang L F,Huang J P 2007 J.Phys.Chem.C 111 10785

    [13]Peng X,Min Y,MAt,Luo W,Yan M2009 T J.Magn.Magn.Mater.321 1221

    [14]Li H,Peng X 2012 J.CoMpu t.Phys.7 1405

    [15]Peskin C S 1977 J.CoMput.Phys.25 220

    [16]Peskin C S 2002 Acta NuMerica 11 479

    [17]N iu X D,Shu C,CheWY T,PeMg Y 2006 Phys.Lett.A 354 173

    [18]He Y L,Wang Y,Li Q 2008 Lattice BoltzMann Method: Theory and Applications(Beijing:Science Press)p31–55(in Chinese)[何雅玲,王勇,李慶 2008格子BoltzMann方法的理論及應(yīng)用(第一版)(北京:科學(xué)出版社)第31—55頁]

    [19]N iu X D,YaMaguchi H,Yoshikawa K 2009 Phys.Rev.E 79 046713

    [20]Hu P,Zhang X W,Niu X D,Meng H 2014 Acta Mech.Sin.46 673(in Chinese)[胡平,張興偉,牛小東,孟輝2014力學(xué)學(xué)報46 673]

    [21]Chen MF,Niu X D,Ma Y R,YaMaguchi H,IwaMoto Y 2015 Procedia Engineering 126 691

    [22]A raseki H,Kotake S 1994 J.CoMpu t.Phys.110 301

    [23]YaMasakiH,YaMaguchiH 2017 J.Magn.Magn.Mater.431 164

    [24]Li L,Chen C,Mei R,K lausner,J F 2014 Phys.Rev.E 89 043308

    [25]Guo K,Li L,X iao G,Au Yeung N,Mei R 2015 In t.J.Heat Mass Transfer 88 306

    [26]Hu Y,Li D,Shu S,N iu X D 2015 CoMput.Math.Appl.70 2227

    [27]Feng J,Hu H H,Joseph D D 1994 J.F luid Mech.261 95

    [28]Feng Z G,Michaelides E E 2004 J.CoMput.Phys.195 602

    [29]Zhang H,Tan Y,Shu S,Niu X D,Trias F X,Yang D,Sheng Y 2014 CoMput.F luids 94 37

    PACS:47.65.Cb,47.11.–j,47.61.JdDOI:10.7498/aps.66.164703

    *Pro ject supported by the National Natural Science Foundation of China(G rant No.11372168).

    ?Corresponding author.E-Mail:15xli1@stu.edu.cn

    Sed iMen tation of tWo non-Magnetic particles in Magnetic fl u id?

    Chen Mu-Feng1)Li Xiang1)?Niu Xiao-Dong1)Li You1)Adnan1)Hiroshi Yamaguchi2)

    1)(College of Engineering,Shantou University,Shantou 515063,China)
    2)(Energy Conversion Research Center,Doshisha University,Kyoto 630-0321,Japan)

    31 March 2017;revised Manuscrip t

    2 June 2017)

    Magnetic fl uid is a stable suspension of solid phaseMagnetic particles of diaMeter about 10 nMin a nonMagnetic carrier fl uid likewater or alcohol.Nowadays,themagnetic fl uid isWidely used in industry areas such as sealing,daMping,lubricating,sound regu lation,heat dissipation,and MHD benefi ciation.Researchers have paid great attention to the behaviors of non-Magnetic particles(NPs)in theMagnetic field becauseMagnetic fluid containing NPs can forMdiff erent Microstructures,which are easily controlled by app lying a magnetic field.In order to app ropriately use the properties of Magnetic fluid in industry,it is necessary to study the interaction aMong NPs in detail.In this paper,a mu ltiphysical nuMerical Model is eMp loyed to investigate the sediMentation of two NPs in Magnetic fluid sub jected to an app lied magnetic field.Themagnetic fluid fl oWis simulated by lattice Boltzmann method,and magneto hyd rodynaMics is calculated With a self-correcting procedure of a Poisson equation solver,which enables the Ohm’s laWto satisfy its conservation law.A dipole forcemodel is used to obtain the dipole interaction force between particles.In addition,as the perMeability of theMagnetic fl uid is quite diff erent froMthose of the NPs and Magnetic fl uid,correctly estab lishing the conjugate boundary condition of theMagnetic intensity at the interface between the particles and surrounding fluid is a key because it aff ects themagnetic induction in the fluid-structure interaction area.A smooth transition scheme of the conjugate boundary condition for Magnetic intensity at the interface between the particles and surrounding fluid is used in this work.The aiMof this work is to investigate sediMentation of two NPs in Magnetized Magnetic fluid.By changing the ratio ofmagnetic permeability and themagnetic parameter,it is found that altering the ratio ofmagnetic permeability isMore eff ective to change the trajectories of NPs,while changing theMagnetic paraMeter can just give rise to a slight transforMof particle trajectories.This can provide good theoretical support for the app lication ofMagnetic fluid in industry area,because the results in the present simu lation can quantitatively analyze the controlling of the MoveMent of NPs.

    magnetic fluid,non-magnetic particles,self-corrected process,conjugate boundary condition

    10.7498/aps.66.164703

    ?國家自然科學(xué)基金(批準號:11372168)資助的課題.

    ?通信作者.E-Mail:15xli1@stu.edu.cn

    ?2017中國物理學(xué)會C h inese P hysica l Society

    http://Wu lixb.iphy.ac.cn

    猜你喜歡
    磁流體磁化比值
    磁流體·吸引力
    中國寶玉石(2024年1期)2024-03-11 04:06:18
    磁流體音箱
    非均勻磁場下磁流體形態(tài)的研究
    電子制作(2019年9期)2019-05-30 09:42:16
    不可壓縮磁流體方程組在Besov空間中的爆破準則
    東北豐磁化炭基復(fù)合肥
    雙色球磁化炭基復(fù)合肥
    比值遙感蝕變信息提取及閾值確定(插圖)
    河北遙感(2017年2期)2017-08-07 14:49:00
    基于磁化能量的鋰電池串模塊化均衡方法
    不同應(yīng)變率比值計算方法在甲狀腺惡性腫瘤診斷中的應(yīng)用
    超強磁場下簡并電子氣體的磁化
    成人毛片a级毛片在线播放| 1024手机看黄色片| 一a级毛片在线观看| 白带黄色成豆腐渣| 色av中文字幕| 日韩欧美精品免费久久| 久久精品影院6| 一进一出抽搐gif免费好疼| 亚洲国产日韩欧美精品在线观看| 欧美xxxx性猛交bbbb| 免费观看的影片在线观看| 日本精品一区二区三区蜜桃| 欧美一区二区精品小视频在线| 老司机深夜福利视频在线观看| 精品久久久久久久久久免费视频| 亚洲一区高清亚洲精品| 亚洲真实伦在线观看| 中文字幕av在线有码专区| 美女高潮的动态| 亚洲国产精品成人综合色| 天堂√8在线中文| 大型黄色视频在线免费观看| 少妇高潮的动态图| 欧美3d第一页| 嫩草影院新地址| 美女免费视频网站| 老熟妇仑乱视频hdxx| 日韩一区二区视频免费看| 一本一本综合久久| 中文字幕熟女人妻在线| 精品久久久久久久久亚洲 | 欧美日韩综合久久久久久 | 国产亚洲精品综合一区在线观看| 很黄的视频免费| 成年人黄色毛片网站| 天堂√8在线中文| 久久99热这里只有精品18| 成人特级黄色片久久久久久久| 美女被艹到高潮喷水动态| 中出人妻视频一区二区| 国产淫片久久久久久久久| 欧美zozozo另类| 欧美性猛交╳xxx乱大交人| 国产亚洲91精品色在线| 国产精品,欧美在线| 草草在线视频免费看| 亚洲av中文av极速乱 | 99久久精品热视频| 亚洲 国产 在线| h日本视频在线播放| av专区在线播放| 极品教师在线视频| 一边摸一边抽搐一进一小说| 成年女人永久免费观看视频| 一a级毛片在线观看| 一区二区三区高清视频在线| 看黄色毛片网站| 麻豆精品久久久久久蜜桃| 亚洲av一区综合| 国产私拍福利视频在线观看| 1000部很黄的大片| 国产成人a区在线观看| 99在线视频只有这里精品首页| 噜噜噜噜噜久久久久久91| 国产一区二区激情短视频| 搡女人真爽免费视频火全软件 | 老女人水多毛片| 人妻久久中文字幕网| 欧美成人免费av一区二区三区| 午夜激情欧美在线| 亚洲欧美日韩高清专用| 97热精品久久久久久| 制服丝袜大香蕉在线| 国产精品综合久久久久久久免费| 在线播放无遮挡| 日日摸夜夜添夜夜添小说| 国产在线精品亚洲第一网站| 内地一区二区视频在线| 国产日本99.免费观看| 午夜激情福利司机影院| 亚洲精品乱码久久久v下载方式| 蜜桃久久精品国产亚洲av| 久久精品国产亚洲av天美| 如何舔出高潮| 九九爱精品视频在线观看| 国内久久婷婷六月综合欲色啪| 非洲黑人性xxxx精品又粗又长| 国产淫片久久久久久久久| 亚洲四区av| 最近在线观看免费完整版| 性色avwww在线观看| 色播亚洲综合网| 亚洲五月天丁香| 中亚洲国语对白在线视频| 亚洲va在线va天堂va国产| 国产黄a三级三级三级人| 精品国产三级普通话版| av天堂在线播放| 三级国产精品欧美在线观看| 国产精品国产高清国产av| 亚洲国产日韩欧美精品在线观看| 国产日本99.免费观看| 男插女下体视频免费在线播放| 一个人免费在线观看电影| 十八禁网站免费在线| 日韩欧美一区二区三区在线观看| 美女免费视频网站| 神马国产精品三级电影在线观看| 不卡一级毛片| 三级男女做爰猛烈吃奶摸视频| 午夜免费男女啪啪视频观看 | 禁无遮挡网站| 天天躁日日操中文字幕| 一个人免费在线观看电影| 日韩强制内射视频| 少妇的逼好多水| 搡女人真爽免费视频火全软件 | 国产激情偷乱视频一区二区| 中文字幕人妻熟人妻熟丝袜美| 久久草成人影院| 51国产日韩欧美| 给我免费播放毛片高清在线观看| 99国产精品一区二区蜜桃av| 制服丝袜大香蕉在线| 九九热线精品视视频播放| 我的女老师完整版在线观看| 黄片wwwwww| 色在线成人网| 午夜免费成人在线视频| 国产精品综合久久久久久久免费| 中文亚洲av片在线观看爽| 国内精品一区二区在线观看| 中文字幕av在线有码专区| 国产色婷婷99| 联通29元200g的流量卡| 又爽又黄a免费视频| 蜜桃久久精品国产亚洲av| 亚洲欧美精品综合久久99| 韩国av在线不卡| 精品一区二区三区视频在线| 男插女下体视频免费在线播放| 色播亚洲综合网| 一边摸一边抽搐一进一小说| 丝袜美腿在线中文| 深爱激情五月婷婷| 久久久久久国产a免费观看| 精品久久久久久久久av| 国产精品,欧美在线| 日韩欧美一区二区三区在线观看| 一区福利在线观看| eeuss影院久久| 深夜a级毛片| 日韩欧美精品v在线| 精品久久久久久,| 亚洲专区中文字幕在线| 此物有八面人人有两片| 午夜免费成人在线视频| 日韩中字成人| 久久婷婷人人爽人人干人人爱| 九九久久精品国产亚洲av麻豆| 亚洲美女视频黄频| av国产免费在线观看| 少妇丰满av| 久久久久久久久久黄片| 国产高清视频在线观看网站| 亚洲av.av天堂| 免费不卡的大黄色大毛片视频在线观看 | 热99在线观看视频| 免费高清视频大片| 在线观看美女被高潮喷水网站| 麻豆成人午夜福利视频| 麻豆精品久久久久久蜜桃| 国产精品久久久久久久久免| 久久久久久久久久久丰满 | 国产欧美日韩精品一区二区| 成人精品一区二区免费| 国产免费男女视频| 欧美精品啪啪一区二区三区| 少妇高潮的动态图| 国产高潮美女av| 久久欧美精品欧美久久欧美| 永久网站在线| 精品一区二区三区人妻视频| 欧美三级亚洲精品| 色哟哟哟哟哟哟| 精品一区二区三区视频在线观看免费| 亚洲中文字幕日韩| 1000部很黄的大片| 嫩草影视91久久| 精华霜和精华液先用哪个| 欧美成人a在线观看| 女人十人毛片免费观看3o分钟| 少妇猛男粗大的猛烈进出视频 | 露出奶头的视频| 久久久久久大精品| 丰满乱子伦码专区| 日本爱情动作片www.在线观看 | 男女视频在线观看网站免费| 性色avwww在线观看| 男人和女人高潮做爰伦理| 热99re8久久精品国产| 亚洲一区高清亚洲精品| 成人亚洲精品av一区二区| 观看免费一级毛片| 国产成人a区在线观看| 波多野结衣巨乳人妻| 国产精品一区二区三区四区免费观看 | 非洲黑人性xxxx精品又粗又长| 51国产日韩欧美| 干丝袜人妻中文字幕| 欧美日韩亚洲国产一区二区在线观看| 亚洲国产精品久久男人天堂| 国产又黄又爽又无遮挡在线| 亚洲,欧美,日韩| 国语自产精品视频在线第100页| 亚洲最大成人中文| 日本在线视频免费播放| 中文字幕久久专区| 内射极品少妇av片p| 韩国av一区二区三区四区| 国产精品亚洲一级av第二区| 色在线成人网| 亚洲最大成人av| 久久久久久国产a免费观看| 天美传媒精品一区二区| 国模一区二区三区四区视频| 一本一本综合久久| 国产视频内射| 亚洲图色成人| 免费观看精品视频网站| 国产午夜精品论理片| 黄色一级大片看看| 亚洲中文字幕日韩| 大型黄色视频在线免费观看| 亚洲av中文av极速乱 | videossex国产| 美女 人体艺术 gogo| 老司机午夜福利在线观看视频| 在线播放国产精品三级| 联通29元200g的流量卡| 久久久久免费精品人妻一区二区| 日韩,欧美,国产一区二区三区 | 国产亚洲精品av在线| 国产精品一区二区三区四区免费观看 | 午夜精品在线福利| 中文字幕精品亚洲无线码一区| 亚洲久久久久久中文字幕| 久久久精品大字幕| a级一级毛片免费在线观看| 国产探花极品一区二区| 麻豆成人av在线观看| 老熟妇仑乱视频hdxx| 2021天堂中文幕一二区在线观| 一a级毛片在线观看| 精品午夜福利在线看| 舔av片在线| 男女视频在线观看网站免费| 1024手机看黄色片| 干丝袜人妻中文字幕| 级片在线观看| 赤兔流量卡办理| 国产熟女欧美一区二区| 搡女人真爽免费视频火全软件 | 夜夜看夜夜爽夜夜摸| 国产精品久久久久久亚洲av鲁大| 在线观看一区二区三区| 美女高潮喷水抽搐中文字幕| 乱人视频在线观看| 桃红色精品国产亚洲av| 免费在线观看成人毛片| 亚洲国产精品久久男人天堂| 久久久久久久久中文| 亚洲色图av天堂| 3wmmmm亚洲av在线观看| 欧美一区二区国产精品久久精品| 99久久精品国产国产毛片| 老司机午夜福利在线观看视频| netflix在线观看网站| 国产探花在线观看一区二区| 亚洲中文日韩欧美视频| 国产黄片美女视频| 国产爱豆传媒在线观看| 国产视频内射| 夜夜夜夜夜久久久久| 真实男女啪啪啪动态图| 国产精品99久久久久久久久| 在线a可以看的网站| 色哟哟哟哟哟哟| 国产男靠女视频免费网站| 高清日韩中文字幕在线| 免费观看人在逋| 国产精品一区www在线观看 | 最好的美女福利视频网| 丰满人妻一区二区三区视频av| 淫秽高清视频在线观看| 亚洲五月天丁香| 精品无人区乱码1区二区| 天堂影院成人在线观看| 熟妇人妻久久中文字幕3abv| 国产免费一级a男人的天堂| 亚洲成人精品中文字幕电影| 日本精品一区二区三区蜜桃| а√天堂www在线а√下载| 简卡轻食公司| 99久久成人亚洲精品观看| 波多野结衣巨乳人妻| 欧美+亚洲+日韩+国产| 久久婷婷人人爽人人干人人爱| 99九九线精品视频在线观看视频| 久9热在线精品视频| 亚洲国产精品久久男人天堂| 久久久久国产精品人妻aⅴ院| 国产精品久久久久久亚洲av鲁大| 国产精品一区二区三区四区久久| 亚洲五月天丁香| 男人狂女人下面高潮的视频| 国产精品电影一区二区三区| 中文亚洲av片在线观看爽| 欧美日韩精品成人综合77777| 床上黄色一级片| 午夜老司机福利剧场| 99久久成人亚洲精品观看| 看十八女毛片水多多多| 国产女主播在线喷水免费视频网站 | 亚洲国产欧美人成| 99riav亚洲国产免费| 亚洲第一电影网av| 国产精品亚洲一级av第二区| 丰满人妻一区二区三区视频av| 天天躁日日操中文字幕| 亚洲不卡免费看| 国产伦一二天堂av在线观看| 男女啪啪激烈高潮av片| 免费av毛片视频| 国产高清有码在线观看视频| 欧美xxxx性猛交bbbb| 国产成人a区在线观看| 国产黄片美女视频| 一级a爱片免费观看的视频| 国产一区二区亚洲精品在线观看| 久久久久国内视频| 欧美另类亚洲清纯唯美| 国产伦在线观看视频一区| 99视频精品全部免费 在线| 亚洲国产高清在线一区二区三| 亚洲精品亚洲一区二区| 人人妻人人澡欧美一区二区| 国产成人福利小说| 亚洲最大成人手机在线| 级片在线观看| 亚洲av五月六月丁香网| 国产精品无大码| 露出奶头的视频| 真人一进一出gif抽搐免费| 久久久久精品国产欧美久久久| 欧美不卡视频在线免费观看| 成人三级黄色视频| 极品教师在线免费播放| 午夜爱爱视频在线播放| 精品久久久噜噜| 久久精品91蜜桃| 我的女老师完整版在线观看| 最近最新免费中文字幕在线| 99热这里只有精品一区| 亚洲欧美日韩东京热| av天堂在线播放| 亚洲无线在线观看| 日韩 亚洲 欧美在线| 淫妇啪啪啪对白视频| 国产亚洲精品久久久久久毛片| 色噜噜av男人的天堂激情| av视频在线观看入口| 欧美3d第一页| 精品福利观看| 国产一区二区亚洲精品在线观看| 99久国产av精品| 亚洲国产高清在线一区二区三| 一区福利在线观看| 日本 欧美在线| 精品国内亚洲2022精品成人| 久久欧美精品欧美久久欧美| 黄色视频,在线免费观看| 成人二区视频| 日本熟妇午夜| 村上凉子中文字幕在线| 我要看日韩黄色一级片| 51国产日韩欧美| 国产亚洲精品综合一区在线观看| 国产精品久久久久久亚洲av鲁大| 亚洲精品国产成人久久av| 国产精品一区www在线观看 | 在线观看av片永久免费下载| 国产国拍精品亚洲av在线观看| 久久久色成人| 99热精品在线国产| 99久久久亚洲精品蜜臀av| 日本精品一区二区三区蜜桃| 男人和女人高潮做爰伦理| 亚洲在线自拍视频| 国产极品精品免费视频能看的| 免费不卡的大黄色大毛片视频在线观看 | 国产午夜福利久久久久久| 国产精品一及| 亚洲最大成人手机在线| 亚洲在线自拍视频| 久久久久久九九精品二区国产| 一级毛片久久久久久久久女| aaaaa片日本免费| 成人永久免费在线观看视频| 国产伦精品一区二区三区视频9| 亚洲四区av| 欧美日韩国产亚洲二区| 久久精品国产自在天天线| 精品久久久久久久末码| 久久人妻av系列| 一级a爱片免费观看的视频| 性插视频无遮挡在线免费观看| 高清在线国产一区| 两个人的视频大全免费| 丰满乱子伦码专区| 如何舔出高潮| 国产爱豆传媒在线观看| 国产黄色小视频在线观看| 国产大屁股一区二区在线视频| 日本免费一区二区三区高清不卡| 精品免费久久久久久久清纯| 真实男女啪啪啪动态图| 国产亚洲精品久久久久久毛片| 久久久久久久午夜电影| 国产69精品久久久久777片| 全区人妻精品视频| 国产黄色小视频在线观看| 乱码一卡2卡4卡精品| 亚洲精品亚洲一区二区| 欧美成人a在线观看| 午夜精品一区二区三区免费看| 国产精品福利在线免费观看| 在线观看午夜福利视频| 午夜免费男女啪啪视频观看 | 精品人妻偷拍中文字幕| 91在线精品国自产拍蜜月| 亚洲人成网站高清观看| 亚洲自拍偷在线| 男人和女人高潮做爰伦理| 国产精品人妻久久久久久| 男女下面进入的视频免费午夜| 久久久久九九精品影院| 成人精品一区二区免费| 国产不卡一卡二| 一本一本综合久久| 欧美激情国产日韩精品一区| 一边摸一边抽搐一进一小说| 久久国产精品人妻蜜桃| 亚洲中文字幕日韩| 精品久久国产蜜桃| 中文字幕av在线有码专区| 亚洲va日本ⅴa欧美va伊人久久| 日本 av在线| 一个人观看的视频www高清免费观看| 高清在线国产一区| 国产 一区 欧美 日韩| 国产一区二区三区av在线 | 韩国av在线不卡| 免费无遮挡裸体视频| 久久久成人免费电影| 精品无人区乱码1区二区| 少妇高潮的动态图| 欧美性猛交╳xxx乱大交人| 十八禁国产超污无遮挡网站| 成人国产一区最新在线观看| 日韩av在线大香蕉| 亚洲av成人精品一区久久| 婷婷色综合大香蕉| 18禁裸乳无遮挡免费网站照片| 午夜免费男女啪啪视频观看 | 色哟哟哟哟哟哟| 欧美性猛交╳xxx乱大交人| 精品午夜福利视频在线观看一区| .国产精品久久| 噜噜噜噜噜久久久久久91| 成年女人毛片免费观看观看9| 一级黄色大片毛片| 九色成人免费人妻av| 内地一区二区视频在线| 日本-黄色视频高清免费观看| 欧美日韩乱码在线| 国产又黄又爽又无遮挡在线| 色在线成人网| 欧美成人a在线观看| 成人综合一区亚洲| 少妇被粗大猛烈的视频| 无人区码免费观看不卡| 欧美日韩国产亚洲二区| 国产一区二区亚洲精品在线观看| 天堂网av新在线| 精品人妻视频免费看| 无人区码免费观看不卡| av在线亚洲专区| 国产精品国产三级国产av玫瑰| 中文字幕精品亚洲无线码一区| 九九在线视频观看精品| 国产一区二区激情短视频| 尤物成人国产欧美一区二区三区| 日韩欧美 国产精品| 亚洲五月天丁香| 亚洲在线自拍视频| 他把我摸到了高潮在线观看| 欧美极品一区二区三区四区| 1000部很黄的大片| 亚洲最大成人手机在线| 国产精品98久久久久久宅男小说| 嫩草影院新地址| 久久久久久久亚洲中文字幕| 精品国产三级普通话版| x7x7x7水蜜桃| 欧美不卡视频在线免费观看| or卡值多少钱| 高清在线国产一区| 日本与韩国留学比较| 久久午夜福利片| 亚洲无线观看免费| 国产私拍福利视频在线观看| 欧美中文日本在线观看视频| 在线观看免费视频日本深夜| 国产伦人伦偷精品视频| 国产大屁股一区二区在线视频| 成人av一区二区三区在线看| а√天堂www在线а√下载| 久99久视频精品免费| 伊人久久精品亚洲午夜| 国产一区二区在线av高清观看| 午夜福利成人在线免费观看| 男人舔奶头视频| 成人无遮挡网站| 国产精品亚洲美女久久久| 亚州av有码| 亚洲精品乱码久久久v下载方式| 91久久精品电影网| 美女cb高潮喷水在线观看| 直男gayav资源| 成人毛片a级毛片在线播放| 欧美性感艳星| 久久久久久久亚洲中文字幕| 午夜影院日韩av| 亚洲四区av| 精品一区二区三区av网在线观看| 久久亚洲精品不卡| 一级黄片播放器| h日本视频在线播放| 国产黄色小视频在线观看| 少妇的逼水好多| 国产伦精品一区二区三区视频9| 久久久国产成人精品二区| 人妻夜夜爽99麻豆av| 日本一二三区视频观看| h日本视频在线播放| 狂野欧美激情性xxxx在线观看| 淫妇啪啪啪对白视频| av国产免费在线观看| 女人被狂操c到高潮| 91午夜精品亚洲一区二区三区 | 亚洲第一电影网av| 91午夜精品亚洲一区二区三区 | 久久久久国产精品人妻aⅴ院| 天堂√8在线中文| 久久久精品大字幕| av国产免费在线观看| 一卡2卡三卡四卡精品乱码亚洲| 色在线成人网| 亚洲成人久久爱视频| ponron亚洲| 精品一区二区三区视频在线| 熟女电影av网| 成人无遮挡网站| av专区在线播放| 他把我摸到了高潮在线观看| 国产69精品久久久久777片| 波多野结衣高清无吗| 热99re8久久精品国产| 国产精品一区二区性色av| 亚洲精品久久国产高清桃花| 熟妇人妻久久中文字幕3abv| 俄罗斯特黄特色一大片| av在线亚洲专区| 噜噜噜噜噜久久久久久91| 美女被艹到高潮喷水动态| 97超视频在线观看视频| 国产在线精品亚洲第一网站| 久久中文看片网| 免费黄网站久久成人精品| 热99在线观看视频| 国产精品久久久久久久久免| 国产不卡一卡二| 精品久久久久久成人av| 色在线成人网| 性欧美人与动物交配| 深夜a级毛片| 成年女人毛片免费观看观看9| 无遮挡黄片免费观看| 狂野欧美白嫩少妇大欣赏| 亚洲乱码一区二区免费版| 两人在一起打扑克的视频| 午夜福利在线在线| 免费看av在线观看网站| 亚洲精品一卡2卡三卡4卡5卡| 国产精品久久视频播放| av天堂中文字幕网| 69人妻影院| 黄色日韩在线| 丝袜美腿在线中文| 久久香蕉精品热| 真人一进一出gif抽搐免费| 精品一区二区三区人妻视频| 国产人妻一区二区三区在| 亚洲中文字幕日韩|