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

    高超聲速M(fèi)HD 球頭激波脫體距離理論求解1)

    2021-11-10 09:49:04李逸翔李進(jìn)平
    力學(xué)學(xué)報(bào) 2021年9期
    關(guān)鍵詞:駐點(diǎn)激波超聲速

    李逸翔 汪 球 羅 凱 李進(jìn)平 趙 偉

    * (中國科學(xué)院力學(xué)研究所高溫氣體動(dòng)力學(xué)國家重點(diǎn)實(shí)驗(yàn)室,北京 100190)

    ? (中國科學(xué)院大學(xué)工程科學(xué)學(xué)院,北京 100049)

    引言

    新一輪的高超聲速技術(shù)對(duì)飛行器在“極端”環(huán)境和“極端”動(dòng)力條件下的飛行提出了新的挑戰(zhàn),它也面臨降熱、減阻、控制等一系列的難題.當(dāng)飛行器以高超聲速飛行時(shí),頭部弓形會(huì)由于激波壓縮和黏性阻滯減速而產(chǎn)生高溫,其溫度高達(dá)幾千度甚至上萬度,高溫下氣體化學(xué)反應(yīng)產(chǎn)生的離子和電子在飛行器周圍形成一個(gè)等離子體鞘層.高超聲速M(fèi)HD流動(dòng)控制是利用飛行器內(nèi)部布置的磁場與電離氣體的相互作用來達(dá)到流動(dòng)控制的目的.它作為一種非接觸的主動(dòng)流動(dòng)控制方法,在飛行器氣動(dòng)力操控、熱環(huán)境管理、“電磁窗口”減弱“黑障”、“電磁舵面”飛行姿態(tài)控制、磁流體發(fā)電等領(lǐng)域均有廣闊的應(yīng)用前景,受到國內(nèi)外的廣泛關(guān)注[1-3].

    實(shí)驗(yàn)研究、數(shù)值模擬和理論分析是研究科學(xué)問題的常用手段,三者相輔相成.高超聲速M(fèi)HD 流動(dòng)控制的研究對(duì)實(shí)驗(yàn)條件要求較高,需要高超聲速設(shè)備能夠模擬高速飛行狀態(tài)下的來流以及產(chǎn)生高磁場強(qiáng)度的磁場發(fā)生裝置;同時(shí),電磁干擾以及脈沖設(shè)備較短的實(shí)驗(yàn)時(shí)間也會(huì)帶來測(cè)量方面的困難[4].由于當(dāng)前實(shí)驗(yàn)可進(jìn)行測(cè)量的物理量較為有限,所以研究主要是針對(duì)激波脫體距離[5-8]、模型氣動(dòng)力(阻力/磁阻力等)[8-12]以及模型表面熱流[13]等方向來開展,但總體來講,隨著高超聲速地面設(shè)備模擬能力及電磁技術(shù)的發(fā)展,越來越多的學(xué)者考慮通過地面實(shí)驗(yàn)來研究磁場流動(dòng)控制的效果及規(guī)律.

    數(shù)值模擬相較于實(shí)驗(yàn)研究而言,所需要的研究成本大為降低,通過數(shù)值模擬可以得到全流場信息,有助于直觀地了解流場規(guī)律.經(jīng)過近二三十年的發(fā)展,尤其是對(duì)計(jì)算格式的不斷修正[14-18],MHD 數(shù)值模擬發(fā)揮著越來越重要的作用.然而,MHD 涉及多種復(fù)雜物理現(xiàn)象,包括產(chǎn)生和維持帶電粒子的過程、帶電粒子與磁場的作用、非平衡能量輸運(yùn)機(jī)理以及電磁場中的霍爾效應(yīng)等,每一種效應(yīng)都包含復(fù)雜的規(guī)律,即使是全MHD 的數(shù)值模擬也是基于一些簡化假設(shè)(如電磁介質(zhì)各向同性、無磁極化、忽略粒子滑移等).因此,MHD 數(shù)值模擬也亟待可靠的MHD 流動(dòng)控制實(shí)驗(yàn)數(shù)據(jù)來對(duì)計(jì)算方法及結(jié)果進(jìn)行校驗(yàn).同時(shí),在相同的計(jì)算條件下,它往往要比CFD模擬耗費(fèi)更大的計(jì)算量;除了CFD 模擬存在的問題外,MHD 模擬還存在一些如磁場散度、全MHD 奇異性、導(dǎo)電壁面邊界等特有的問題[19].

    基于高超聲速M(fèi)HD 的理論研究雖然起步非常早,但相比于實(shí)驗(yàn)研究和數(shù)值模擬,其研究工作仍然非常有限,且成果主要集中于二十世紀(jì)五六十年代.較為經(jīng)典的是Bush[20]的研究工作,他在激波層內(nèi)常密度、常電導(dǎo)率、無黏、球形激波等假設(shè)下,采用分離變量的方法求解流動(dòng)方程和電磁方程,將復(fù)雜的偏微分方程組化簡為常微分方程組,以轉(zhuǎn)化為一個(gè)初值問題進(jìn)行離散求解.Smith 和Wu[21]、Porter和Cambel[22]沿用Bush 的方法,分別研究了黏性系數(shù)和Hall 效應(yīng)的影響.該領(lǐng)域其他方面的理論研究工作也在同時(shí)期展開,Ludford 和Murray[23]研究了不可壓縮無黏導(dǎo)電流體流經(jīng)磁化球的穩(wěn)態(tài)流動(dòng)問題.Meyer[24]采用Newton-Plus-Centrifugal 近似方法,發(fā)展了MHD回轉(zhuǎn)體穩(wěn)態(tài)流動(dòng)理論.Lykoudis[25]引入Newton-Busemann 壓力定律,通過對(duì)無黏動(dòng)量方程進(jìn)行積分得到激波脫體距離的表達(dá)式.Eriscon 和Maciulaitis[26]研究了高度?速度圖上的電離流與磁場相互作用的分布范圍,并進(jìn)一步評(píng)估了較強(qiáng)相互作用流態(tài)下的磁流體控制作用力.Hooks 和Lewis[27]基于激波坐標(biāo)系和物體坐標(biāo)系之間的關(guān)系,得到了駐點(diǎn)處的無黏速度梯度、物體周圍的靜壓分布和邊界層外沿的速度分布等.激波脫體距離作為高超聲速飛行器飛行問題中的一個(gè)重要參數(shù),是氣動(dòng)領(lǐng)域研究的關(guān)注重點(diǎn)之一.然而,上述理論工作對(duì)于磁場影響下激波脫體距離的研究相對(duì)較少,相關(guān)結(jié)果也較為有限.有學(xué)者采用的方法雖然可以得到駐點(diǎn)線附近的流場信息,但其求解過程和結(jié)果不能直觀快速的反映磁控效果,有一定局限性[20-22];另一方面,Lykoudis[25]、Hooks 和Lewis[27]對(duì)激波脫體距離解析形式的理論推導(dǎo)是基于激波層內(nèi)磁場方向處處垂直于流速方向的條件下進(jìn)行的,該類型磁場可與電離流充分作用,然而其分布形式較為理想,且積分形式的激波脫體距離表達(dá)式需要進(jìn)一步數(shù)值求解.因此,高超聲速M(fèi)HD 激波脫體距離的快速求解方法仍有待進(jìn)一步發(fā)展和完善,以達(dá)到從理論上快速評(píng)估磁控效果以及能夠反映MHD 流動(dòng)控制本質(zhì)的目標(biāo).

    本文通過理論分析的方法得到高超聲速M(fèi)HD無量綱激波脫體距離的解析表達(dá)式.求解過程主要分為兩部分,首先對(duì)連續(xù)方程進(jìn)行積分,將無量綱激波脫體距離作為待求量包含在內(nèi);進(jìn)一步,引入切向速度切向偏導(dǎo)數(shù)沿駐點(diǎn)線上的線性近似分布,利用徑向動(dòng)量方程建立磁場與激波脫體距離的聯(lián)系,尋找球頭駐點(diǎn)處切向速度沿切向的偏導(dǎo)數(shù)值.最后,通過求解關(guān)于激波脫體距離的一元二次方程,獲得包含磁場影響的解析表達(dá)式.

    1 無量綱激波脫體距離理論求解

    磁場添加后對(duì)流動(dòng)影響最直觀的體現(xiàn)是脫體激波距離或激波結(jié)構(gòu)的變化,學(xué)者在多種設(shè)備中均通過光學(xué)方法觀測(cè)到了磁場對(duì)脫體激波距離的影響,一般而言,磁場越強(qiáng),脫體激波距離越大.

    本文的理論研究是基于以下條件:

    (1)低磁雷諾數(shù)假設(shè),對(duì)于高超聲速飛行條件,激波層內(nèi)氣體電離度一般較低,磁雷諾數(shù)遠(yuǎn)小于1,此時(shí)誘導(dǎo)磁場對(duì)于外加磁場的影響可以忽略不記[1,28-29].Porter 和Cambel[29]詳細(xì)研究了不同磁雷諾數(shù)條件下的磁場分布,其結(jié)果也表明在低磁雷諾數(shù)條件下,電流導(dǎo)致的誘導(dǎo)磁場對(duì)外加磁場的影響很小,在激波層內(nèi)只考慮外加磁場即可;

    (2)激波層內(nèi)流體不可壓縮,盡管自由來流速度是高超聲速,但強(qiáng)激波后激波層內(nèi)駐點(diǎn)線附近馬赫數(shù)較小,該區(qū)域內(nèi)的流體可當(dāng)不可壓縮流體近似處理,Bush[20]、Smith 和Wu[21]在求解高超聲速M(fèi)HD問題時(shí)也采用了該假設(shè);

    (3)忽略黏性影響,高超聲速來流條件下,流體雷諾數(shù)較大,在關(guān)注脫體激波距離問題時(shí)可以忽略黏性項(xiàng),將問題簡化為無黏流動(dòng)問題[20,25];

    (4)激波層內(nèi)電導(dǎo)率為常值,電導(dǎo)率一般與溫度相關(guān),在無黏假設(shè)條件下,激波層內(nèi)不考慮邊界層,沿駐 點(diǎn)線溫度變化較小,可將其以常值近似處理[17,20,25].

    1.1 對(duì)連續(xù)方程的處理

    對(duì)連續(xù)方程的處理借鑒了Olivier[30]的方法,他通過對(duì)連續(xù)方程積分的方式獲得了包含激波脫體距離的方程,本文在其基礎(chǔ)上進(jìn)一步考慮添加磁場對(duì)脫體激波距離的影響,對(duì)該方法的具體應(yīng)用過程如下.球坐標(biāo)系下的連續(xù)方程為

    其中 ρ∞和u∞為來流密度和速度,ρs為激波層內(nèi)密度,rs和rb分別為球頭半徑和激波半徑,Δ 為激波脫體距離,各變量示意圖如圖1 所示.

    圖1 鈍頭體和激波結(jié)構(gòu)示意圖Fig.1 Schematic of the blunt body and the shock

    將式(2)沿著駐點(diǎn)線由駐點(diǎn)向激波處進(jìn)行積分,積分的目的是將無量綱激波脫體距離加入到方程中,將其作為待求量來進(jìn)一步求解,積分過程具體如下

    由于求解過程在駐點(diǎn)線附近,此時(shí)可以通過引入駐點(diǎn)線假設(shè)簡化問題的求解過程,此過程實(shí)質(zhì)上是去除方程中的高階小量,同時(shí)保留低階小量,從而保證了低階量的可求解性.根據(jù)駐點(diǎn)線假設(shè),引入如下近似條件[30]

    Olivier[30]近似地給出了/?θ在駐點(diǎn)線上的線性分布表達(dá)式,既可以使/?θ滿足激波處的邊界條件,同時(shí)提供了壁面處 (/?θ)b的待求量,通過 (/?θ)b可以建立激波脫體距離與磁場和來流的聯(lián)系.以下給出/?θ在駐點(diǎn)線附近的線性分布形式

    1.2 對(duì)動(dòng)量方程的處理

    動(dòng)量方程包含了待求問題的特征元素,即電磁作用項(xiàng),同時(shí)需要通過對(duì)動(dòng)量方程的處理,尋找激波脫體距離與磁場影響的聯(lián)系,因此這一部分較為關(guān)鍵.通過處理徑向動(dòng)量方程,得到駐點(diǎn)處/?θ的具體數(shù)值,即 (/?θ)b.此處選取徑向動(dòng)量方程求解,這是考慮到壓力項(xiàng)的量級(jí)在分析的過程中會(huì)被排除掉,有助于求解的簡化.同時(shí),需要額外地說明一下,對(duì)動(dòng)量方程取旋度的方式是不可行的,原因是分離變量后,在駐點(diǎn)處方程兩邊恒為0,對(duì)于問題的求解無作用.對(duì)徑向的動(dòng)量方程進(jìn)行處理具體步驟如下.包含電磁項(xiàng)的徑向動(dòng)量方程為

    其中σ為電導(dǎo)率.根據(jù)電磁感應(yīng)定律,變化的磁場會(huì)導(dǎo)致感生電場,但由于本文涉及的是穩(wěn)態(tài)磁場,所以電場需要滿足?×E=0,對(duì)于二維軸對(duì)稱問題,在駐點(diǎn)線上滿足E=0[20],于是式(9)可以進(jìn)一步整理

    對(duì)于磁場分布,選取較為常用的模型為均勻磁化球形磁體的磁場分布,該模型的磁場分布簡單,存在解析形式的表達(dá),便于實(shí)驗(yàn)和理論的研究[12,20,22],這種磁場的方向與激波層內(nèi)電離流速度方向趨近于垂直,可以較為充分地與電離流發(fā)生相互作用,其分布形式見圖2.由式(10)知,電流為環(huán)向,垂直于磁場和速度矢量所在的面,于是MHD 動(dòng)量方程中的電磁力項(xiàng)J×B會(huì)受到增益的影響.這里參考Porter給出的偶極子磁場分布形式[22]

    圖2 均勻磁化球體磁感線分布示意圖Fig.2 Distribution of magnetic induction lines of uniformly magnetized sphere

    根據(jù)量級(jí)分析,式(16) 的壓力項(xiàng)不包含小量O(θ2),所以當(dāng)θ →0 時(shí),等式兩端包含O(θ2)的項(xiàng)都趨近0,而此時(shí)壓力項(xiàng)為0,將該結(jié)果代入式(16),可以得到

    由于磁場會(huì)明顯影響流場,切向速度需要包含磁場以及來流所帶來的影響,所以這里選取第二個(gè)解.額外地,由式(16)還會(huì)得到于是可以發(fā)現(xiàn),當(dāng)出現(xiàn)磁場時(shí),駐點(diǎn)附近地切向速度有向駐點(diǎn)方向流動(dòng)的趨勢(shì),而該處的徑向速度可由微分形式的連續(xù)方程得到

    此結(jié)果說明,駐點(diǎn)附近徑向速度方向有與來流相反的趨勢(shì).

    根據(jù)激波前后質(zhì)量守恒關(guān)系,激波處駐點(diǎn)線上的速度有[19]

    2 結(jié)果與分析

    根據(jù)圖3,對(duì)于相同的相互作用系數(shù)Q,隨著ε的增加(來流速度降低)無量綱激波脫體距離變大.對(duì)于相同的ε,隨著Q增加激波無量綱脫體距離增加.從圖3 中也可以觀察到隨著ε的減小,無量綱激波脫體距離關(guān)于ε的靈敏度降低,在高超聲速M(fèi)HD實(shí)驗(yàn)或數(shù)值模擬開始之前,可以通過合理估計(jì)的ε值以及由初始條件構(gòu)成的Q來預(yù)測(cè)結(jié)果,這有助于快速評(píng)估磁控效果.

    圖3 無量綱激波脫體距離與磁相互作用系數(shù)Q 的關(guān)系曲線Fig.3 The curve of dimensionless shock stand-off distance vs.MHD interaction parameter Q

    圖4 對(duì)比了本文結(jié)果與Lefevre 等[31]在膨脹管中的實(shí)驗(yàn)結(jié)果.C2~ C6 代表Lefevre 選取的不同實(shí)驗(yàn)條件,每一個(gè)實(shí)驗(yàn)條件對(duì)應(yīng)著不同的相互作用系數(shù)Q以及激波前后氣流密度比ε,但由于實(shí)驗(yàn)時(shí)間很短,不足以形成穩(wěn)定流場,文章以20 μs 為間隔進(jìn)行數(shù)據(jù)取點(diǎn),Lefevre 同時(shí)對(duì)比了MHD 的CFD 模擬結(jié)果.從圖4 中可以看到,理論結(jié)果、實(shí)驗(yàn)和數(shù)值結(jié)果符合得較好.

    圖4 無量綱激波脫體距離比/NM 與磁相互作用系數(shù)Q 的關(guān)系曲線Fig.4 The curve of ratio of shock stand-off distance /NM vs.MHD interaction parameter Q

    Lykoudis[25]同樣在低磁雷諾數(shù)以及激波層內(nèi)常值密度、常值電導(dǎo)率、無黏的假設(shè)條件下推導(dǎo)了球頭激波脫體距離表達(dá)式,他采用了較為理想的磁場分布模型:J=σvBφ,即磁場處垂直于流線,這樣的磁場分布有助于問題簡化.在對(duì)動(dòng)量方程的處理過程中,他認(rèn)為壓力分布隨外加磁場的變化很小,于是通過引入Newton-Busemann 壓力分布近似公式進(jìn)一步補(bǔ)充了流場信息,結(jié)合流管上下游的質(zhì)量守恒關(guān)系,最終推導(dǎo)得到了激波脫體距離積分形式的解析表達(dá)式,分析結(jié)果反映了磁場對(duì)氣流流動(dòng)的影響效果.

    圖5 給出了本文與Lykoudis 理論結(jié)果的對(duì)比,兩者結(jié)果趨勢(shì)一致,都體現(xiàn)了磁場與導(dǎo)電流體相互作用產(chǎn)生的阻滯效應(yīng),但兩者也存在一定偏差,主要原因如下:Lykoudis 采用的磁場分布模型與本文存在一定的差別,他的推導(dǎo)是基于激波層內(nèi)磁場處處垂直于流速的條件,相較于偶極子磁場而言,Lykoudis使用的磁場可與流場進(jìn)行更加充分的相互作用,這也會(huì)產(chǎn)生更強(qiáng)的阻滯作用,使激波層厚度增益效果更加明顯,從而導(dǎo)致其結(jié)果在一定Q值范圍內(nèi)高于本文所得的結(jié)果.當(dāng)然,由于Lykoudis 采用的磁場過于理想,在現(xiàn)實(shí)中難以用磁場裝置產(chǎn)生該類磁場,在實(shí)驗(yàn)研究中難以采用.另外,需要注意的是,包括Lykoudis 文章在內(nèi)的一些文獻(xiàn)中的相互作用系數(shù)表達(dá)式為Qs=σB2rs/(ρ∞u∞),與本文使用的特征長度不一致,需要進(jìn)行轉(zhuǎn)換.

    圖5 本文結(jié)果與Lykoudis 結(jié)果對(duì)比圖[25]Fig.5 Comparison between results of this paper and the Lykoudis’s

    通常情況下,理論分析引入的假設(shè)條件會(huì)導(dǎo)致模型存在一定的局限性.根據(jù)線性近似分布假設(shè),/?θ 是僅關(guān)于r?的函數(shù),當(dāng)來流速度越高時(shí),對(duì)于相同磁相互作用系數(shù)Q條件下,慣性力的影響越強(qiáng),/?θ在駐點(diǎn)線上的分布越符合線性近似假設(shè),然而隨著激波脫體距離的增加,該假設(shè)會(huì)出現(xiàn)一定程度的失效,即/?θ 將與r?和兩個(gè)量相關(guān),/?θ的分布與激波層厚度有關(guān)聯(lián),于是式(5)的積分結(jié)果將會(huì)受到影響,導(dǎo)致無法將其化為關(guān)于激波脫體距離的一元二次方程形式.總之,在Q較大時(shí),/?θ分布偏離線性假設(shè)是限制本文模型應(yīng)用的關(guān)鍵點(diǎn).然而,在Q>1 時(shí)就會(huì)出現(xiàn)明顯的MHD 流動(dòng)控制效果,且當(dāng)前實(shí)驗(yàn)條件還難以達(dá)到較大的Q值,因此本文的模型仍然有較強(qiáng)的指導(dǎo)意義.當(dāng)然,后續(xù)會(huì)通過更多的理論和實(shí)驗(yàn)來對(duì)該模型在Q較大時(shí),/?θ分布偏離線性假設(shè)是限制本文模型應(yīng)用值下進(jìn)行改進(jìn)優(yōu)化.

    本文的結(jié)果中有一個(gè)值得關(guān)注的地方,當(dāng)磁場存在時(shí),駐點(diǎn)附近的切向速度有指向駐點(diǎn)的趨勢(shì),通過連續(xù)方程得到與來流相反的徑向速度,如圖6.這里認(rèn)為該結(jié)果是低磁條件假設(shè)造成的,原因是磁場的旋度方程中,如果考慮電流對(duì)磁場的影響,該影響將會(huì)是一個(gè)小量,這個(gè)小量會(huì)影響徑向動(dòng)量方程的量級(jí)分析過程,從而改變求解結(jié)果.然而這種流動(dòng)趨勢(shì)是實(shí)際存在的,由于電離流與磁場相互作用的影響,流體會(huì)受到與流動(dòng)相反方向的洛倫茲力,從而使激波層內(nèi)的流速降低,從而導(dǎo)致了激波脫體距離的增加.

    圖6 駐點(diǎn)處流動(dòng)示意圖Fig.6 Schematic of flow near stagnation point

    3 結(jié)論

    本文利用連續(xù)方程建立了關(guān)于激波脫體距離的待求方程,通過徑向動(dòng)量方程尋找 (?v/?θ)b的具體數(shù)值,從而使磁場影響與待求方程得以關(guān)聯(lián),將問題化為求解一個(gè)關(guān)于激波脫體距離的一元二次方程,得到了較為直觀簡潔的解析表達(dá)式.本文主要結(jié)論如下:

    (1)理論分析結(jié)果表明,當(dāng)磁相互作用系數(shù)增加時(shí),激波脫體距離增加,即磁控效果增強(qiáng).將本文結(jié)果與Lefevre 的實(shí)驗(yàn)數(shù)據(jù)、Lykoudis 的結(jié)果進(jìn)行了對(duì)比,三者吻合較好.本文理論模型對(duì)高超聲速M(fèi)HD實(shí)驗(yàn)設(shè)計(jì)和結(jié)果分析具有較好的指導(dǎo)意義;

    (2)本文模型的建立是基于特定的假設(shè)條件,存在一定的適用范圍,在Q較大時(shí),/?θ分布偏離線性假設(shè)會(huì)限制本文模型的應(yīng)用.通常情況下,當(dāng)Q>1 時(shí)即會(huì)出現(xiàn)明顯的MHD 流動(dòng)控制效果,本文結(jié)果的適用性可以達(dá)到并遠(yuǎn)大于該范圍.

    猜你喜歡
    駐點(diǎn)激波超聲速
    高超聲速出版工程
    高超聲速飛行器
    一種基于聚類分析的二維激波模式識(shí)別算法
    基于HIFiRE-2超燃發(fā)動(dòng)機(jī)內(nèi)流道的激波邊界層干擾分析
    斜激波入射V形鈍前緣溢流口激波干擾研究
    基于游人游賞行為的留園駐點(diǎn)分布規(guī)律研究
    中國園林(2018年7期)2018-08-07 07:07:48
    超聲速旅行
    適于可壓縮多尺度流動(dòng)的緊致型激波捕捉格式
    利用遠(yuǎn)教站點(diǎn),落實(shí)駐點(diǎn)干部帶學(xué)
    利用遠(yuǎn)教站點(diǎn),落實(shí)駐點(diǎn)干部帶學(xué)
    免费观看性生交大片5| 久久婷婷青草| 成人18禁高潮啪啪吃奶动态图 | 深夜a级毛片| 日韩,欧美,国产一区二区三区| 美女中出高潮动态图| 亚洲精品国产av蜜桃| 97精品久久久久久久久久精品| 国产精品免费大片| 又大又黄又爽视频免费| 久久人妻熟女aⅴ| 午夜视频国产福利| 欧美xxxx黑人xx丫x性爽| 一级毛片 在线播放| 少妇人妻久久综合中文| 欧美日韩视频精品一区| 又大又黄又爽视频免费| av又黄又爽大尺度在线免费看| 国国产精品蜜臀av免费| 精品少妇久久久久久888优播| 一级毛片久久久久久久久女| 女人十人毛片免费观看3o分钟| 日韩欧美精品免费久久| 国产在线视频一区二区| 国产69精品久久久久777片| 女的被弄到高潮叫床怎么办| 高清在线视频一区二区三区| 一个人免费看片子| 性高湖久久久久久久久免费观看| av福利片在线观看| 亚洲精品日韩在线中文字幕| 精华霜和精华液先用哪个| 午夜福利网站1000一区二区三区| 久久精品久久久久久噜噜老黄| 菩萨蛮人人尽说江南好唐韦庄| 一级毛片电影观看| 观看美女的网站| 嫩草影院入口| 最近的中文字幕免费完整| 一个人免费看片子| 中国美白少妇内射xxxbb| 久久女婷五月综合色啪小说| 免费观看在线日韩| 在线观看免费日韩欧美大片 | 久久精品夜色国产| 国产在线一区二区三区精| 久久av网站| 精品国产一区二区三区久久久樱花 | 晚上一个人看的免费电影| av黄色大香蕉| 欧美变态另类bdsm刘玥| 国产在线视频一区二区| 国产成人精品久久久久久| 亚洲成人中文字幕在线播放| 美女国产视频在线观看| 国内揄拍国产精品人妻在线| 午夜激情久久久久久久| 欧美日韩视频高清一区二区三区二| 国产黄片视频在线免费观看| 女人十人毛片免费观看3o分钟| 午夜福利视频精品| 美女cb高潮喷水在线观看| 日韩 亚洲 欧美在线| 如何舔出高潮| 大又大粗又爽又黄少妇毛片口| 99国产精品免费福利视频| 99视频精品全部免费 在线| 国产爱豆传媒在线观看| 国产 一区精品| 我的女老师完整版在线观看| 国产精品爽爽va在线观看网站| 亚洲成人手机| 一级av片app| 精品视频人人做人人爽| 久久久久久九九精品二区国产| 国产亚洲午夜精品一区二区久久| 99热全是精品| 亚洲性久久影院| 日本av手机在线免费观看| a级毛色黄片| 80岁老熟妇乱子伦牲交| 汤姆久久久久久久影院中文字幕| 成人18禁高潮啪啪吃奶动态图 | av线在线观看网站| 色综合色国产| www.色视频.com| 18禁在线无遮挡免费观看视频| 国产精品久久久久久精品古装| 麻豆成人av视频| 男人舔奶头视频| 又粗又硬又长又爽又黄的视频| 国模一区二区三区四区视频| 日日啪夜夜爽| 女性生殖器流出的白浆| 欧美日韩精品成人综合77777| 国产精品女同一区二区软件| 99热全是精品| 中国美白少妇内射xxxbb| 亚洲激情五月婷婷啪啪| 国产精品欧美亚洲77777| 日韩一区二区三区影片| 熟女av电影| 纵有疾风起免费观看全集完整版| 五月玫瑰六月丁香| 午夜老司机福利剧场| 成年av动漫网址| 免费看光身美女| 国产亚洲5aaaaa淫片| 国产亚洲欧美精品永久| 久久久色成人| 日韩视频在线欧美| 中文天堂在线官网| 丰满人妻一区二区三区视频av| 久久 成人 亚洲| 色婷婷久久久亚洲欧美| 亚洲国产精品一区三区| 免费高清在线观看视频在线观看| 久久99蜜桃精品久久| 菩萨蛮人人尽说江南好唐韦庄| 国产精品欧美亚洲77777| 91在线精品国自产拍蜜月| 我要看黄色一级片免费的| 午夜激情福利司机影院| 日韩免费高清中文字幕av| 日韩一区二区三区影片| 极品教师在线视频| 中文字幕免费在线视频6| 中文字幕久久专区| 夜夜看夜夜爽夜夜摸| 亚洲av成人精品一二三区| 又粗又硬又长又爽又黄的视频| 亚洲四区av| tube8黄色片| 熟女电影av网| 久久久久视频综合| 在线免费十八禁| 男人添女人高潮全过程视频| 国产免费一级a男人的天堂| 日韩三级伦理在线观看| 一本久久精品| 国产精品麻豆人妻色哟哟久久| 日本-黄色视频高清免费观看| 久久久久久伊人网av| 丰满少妇做爰视频| 伊人久久国产一区二区| 国产欧美亚洲国产| 成人国产麻豆网| 国产免费视频播放在线视频| 久久精品国产亚洲av天美| 免费av不卡在线播放| 亚洲av成人精品一二三区| 亚洲婷婷狠狠爱综合网| videos熟女内射| 欧美三级亚洲精品| 免费不卡的大黄色大毛片视频在线观看| 一区二区三区免费毛片| 亚洲美女搞黄在线观看| 成人亚洲精品一区在线观看 | 日韩精品有码人妻一区| 国产黄片视频在线免费观看| 国产欧美日韩一区二区三区在线 | 激情 狠狠 欧美| 99热6这里只有精品| 久久毛片免费看一区二区三区| av网站免费在线观看视频| 精品久久久久久久久av| 一级黄片播放器| 国产精品人妻久久久久久| 欧美精品一区二区免费开放| 91在线精品国自产拍蜜月| 国产91av在线免费观看| 夜夜看夜夜爽夜夜摸| 国产亚洲91精品色在线| 51国产日韩欧美| 久久久欧美国产精品| 男人舔奶头视频| 18禁在线播放成人免费| 最近的中文字幕免费完整| 80岁老熟妇乱子伦牲交| 国产精品国产三级专区第一集| 青春草亚洲视频在线观看| 天堂俺去俺来也www色官网| 国产毛片在线视频| 国产 一区精品| 91精品国产国语对白视频| 高清毛片免费看| 久久久久久久久久人人人人人人| xxx大片免费视频| 国产精品99久久99久久久不卡 | 国产精品欧美亚洲77777| 日韩一区二区三区影片| 国产探花极品一区二区| 男人舔奶头视频| 亚洲精品国产色婷婷电影| 国产亚洲一区二区精品| 网址你懂的国产日韩在线| 国产成人91sexporn| 大香蕉97超碰在线| 中国三级夫妇交换| 下体分泌物呈黄色| 在线观看一区二区三区激情| 大话2 男鬼变身卡| 黑人猛操日本美女一级片| 欧美丝袜亚洲另类| 人妻一区二区av| 国产一区二区三区av在线| 毛片一级片免费看久久久久| 中文字幕免费在线视频6| 最近中文字幕2019免费版| 日韩三级伦理在线观看| 国产精品99久久99久久久不卡 | 国产在视频线精品| 欧美变态另类bdsm刘玥| 国产成人免费观看mmmm| 亚洲成色77777| 中文天堂在线官网| videossex国产| 在线观看免费视频网站a站| 精华霜和精华液先用哪个| 亚洲精品日本国产第一区| 内射极品少妇av片p| 国产男人的电影天堂91| 女的被弄到高潮叫床怎么办| 中文天堂在线官网| 亚洲国产最新在线播放| 日日摸夜夜添夜夜添av毛片| 亚洲国产精品专区欧美| 久久精品久久久久久久性| 你懂的网址亚洲精品在线观看| 国产精品女同一区二区软件| 日韩电影二区| 国产伦精品一区二区三区视频9| 国产美女午夜福利| 亚洲国产高清在线一区二区三| 中国国产av一级| 精品人妻一区二区三区麻豆| 91aial.com中文字幕在线观看| 99久久精品国产国产毛片| 精品少妇黑人巨大在线播放| 免费久久久久久久精品成人欧美视频 | 亚洲三级黄色毛片| av国产免费在线观看| av女优亚洲男人天堂| 国产白丝娇喘喷水9色精品| 亚洲欧洲日产国产| av在线老鸭窝| 亚洲欧美成人综合另类久久久| 亚洲欧美一区二区三区国产| av.在线天堂| 免费观看性生交大片5| 秋霞伦理黄片| 国产成人免费观看mmmm| 18禁动态无遮挡网站| 男女啪啪激烈高潮av片| 最近手机中文字幕大全| 九草在线视频观看| 日韩强制内射视频| 国产大屁股一区二区在线视频| 国产精品偷伦视频观看了| 国产亚洲精品久久久com| 99re6热这里在线精品视频| 九色成人免费人妻av| 欧美高清性xxxxhd video| 日韩av不卡免费在线播放| 一区二区三区乱码不卡18| 色5月婷婷丁香| 国产一区亚洲一区在线观看| 亚洲不卡免费看| av在线观看视频网站免费| 免费观看的影片在线观看| 精品酒店卫生间| 国产高清国产精品国产三级 | 搡老乐熟女国产| 一个人免费看片子| 国产精品av视频在线免费观看| 99热6这里只有精品| 日本-黄色视频高清免费观看| 一二三四中文在线观看免费高清| 中国美白少妇内射xxxbb| 久久久国产一区二区| av黄色大香蕉| 亚洲美女视频黄频| 晚上一个人看的免费电影| 亚洲精品亚洲一区二区| 国产精品久久久久久精品电影小说 | 99九九线精品视频在线观看视频| 久热这里只有精品99| 日日撸夜夜添| 成年av动漫网址| 街头女战士在线观看网站| 成年人午夜在线观看视频| 黄色日韩在线| 啦啦啦中文免费视频观看日本| 99久久精品一区二区三区| 久久久久久久久大av| 18禁动态无遮挡网站| 国产亚洲5aaaaa淫片| 国产精品久久久久久精品电影小说 | 天堂8中文在线网| 亚洲av中文字字幕乱码综合| 国产亚洲91精品色在线| 99热这里只有是精品50| 99久久中文字幕三级久久日本| 午夜精品国产一区二区电影| 国产精品无大码| 99热6这里只有精品| 两个人的视频大全免费| 国产精品嫩草影院av在线观看| 国产黄片美女视频| 免费观看性生交大片5| 久久亚洲国产成人精品v| 国产精品三级大全| 在线观看免费视频网站a站| 国产精品国产三级国产av玫瑰| 18禁动态无遮挡网站| 街头女战士在线观看网站| 久久久久国产精品人妻一区二区| 欧美精品一区二区大全| 大香蕉97超碰在线| 欧美xxⅹ黑人| 春色校园在线视频观看| 一级毛片黄色毛片免费观看视频| 91精品伊人久久大香线蕉| 国产成人freesex在线| 国产精品一区二区在线不卡| 小蜜桃在线观看免费完整版高清| 麻豆国产97在线/欧美| 丝袜喷水一区| 日韩三级伦理在线观看| 亚洲人成网站高清观看| 亚洲一级一片aⅴ在线观看| 精品亚洲乱码少妇综合久久| 国产精品av视频在线免费观看| 欧美日韩亚洲高清精品| 亚洲天堂av无毛| 日本黄色片子视频| 国产一区二区三区av在线| a级毛色黄片| 男人舔奶头视频| 精品人妻熟女av久视频| 久久精品人妻少妇| 丝袜脚勾引网站| 午夜福利在线观看免费完整高清在| 一区二区三区免费毛片| 插逼视频在线观看| 国产精品久久久久成人av| 欧美最新免费一区二区三区| 亚洲成人手机| 国产av精品麻豆| 国产免费一区二区三区四区乱码| 欧美区成人在线视频| 久久午夜福利片| 日韩一本色道免费dvd| 最近手机中文字幕大全| 精品酒店卫生间| 国精品久久久久久国模美| 国产伦精品一区二区三区四那| 高清av免费在线| 精品一区二区三卡| 亚洲内射少妇av| 国产色婷婷99| 有码 亚洲区| 男人添女人高潮全过程视频| 亚洲国产欧美人成| 在线播放无遮挡| 日韩强制内射视频| 伦精品一区二区三区| 一级毛片久久久久久久久女| 看免费成人av毛片| 六月丁香七月| 成人亚洲欧美一区二区av| 国产精品久久久久久av不卡| 免费av不卡在线播放| 欧美xxxx黑人xx丫x性爽| 国产高清不卡午夜福利| 亚洲四区av| 三级国产精品欧美在线观看| 亚洲国产精品成人久久小说| 男人舔奶头视频| 韩国av在线不卡| 自拍偷自拍亚洲精品老妇| 热re99久久精品国产66热6| 亚洲精品久久久久久婷婷小说| 最近的中文字幕免费完整| 亚洲精品视频女| 韩国av在线不卡| 国产黄频视频在线观看| 中文字幕亚洲精品专区| 欧美bdsm另类| 久久精品国产a三级三级三级| 日韩不卡一区二区三区视频在线| 老师上课跳d突然被开到最大视频| 丰满人妻一区二区三区视频av| 国产色爽女视频免费观看| 成人国产av品久久久| 成人二区视频| 成人国产麻豆网| 亚洲国产精品国产精品| 十分钟在线观看高清视频www | 91精品一卡2卡3卡4卡| 丰满迷人的少妇在线观看| 在线精品无人区一区二区三 | 寂寞人妻少妇视频99o| 亚洲四区av| 边亲边吃奶的免费视频| 国产视频内射| 黄片wwwwww| 亚洲av在线观看美女高潮| 亚洲国产精品专区欧美| 青春草视频在线免费观看| 久久久国产一区二区| 联通29元200g的流量卡| 丰满乱子伦码专区| 亚洲性久久影院| 国精品久久久久久国模美| 黄色日韩在线| 久久久久网色| 国产精品福利在线免费观看| 熟妇人妻不卡中文字幕| 国产一区亚洲一区在线观看| 小蜜桃在线观看免费完整版高清| 校园人妻丝袜中文字幕| 黑丝袜美女国产一区| 国内精品宾馆在线| av网站免费在线观看视频| 麻豆乱淫一区二区| 国产 一区精品| 丝瓜视频免费看黄片| 十分钟在线观看高清视频www | 伊人久久精品亚洲午夜| av国产久精品久网站免费入址| 天堂中文最新版在线下载| 成年人午夜在线观看视频| 亚洲av国产av综合av卡| 嘟嘟电影网在线观看| 亚洲怡红院男人天堂| 99热这里只有是精品50| 国产精品爽爽va在线观看网站| 99国产精品免费福利视频| 91精品伊人久久大香线蕉| 日本免费在线观看一区| 精品人妻熟女av久视频| 久久国产乱子免费精品| 高清在线视频一区二区三区| 有码 亚洲区| 高清视频免费观看一区二区| 毛片一级片免费看久久久久| 亚洲欧美精品自产自拍| 麻豆精品久久久久久蜜桃| av免费观看日本| 极品少妇高潮喷水抽搐| 国产亚洲欧美精品永久| 欧美人与善性xxx| 三级国产精品片| 22中文网久久字幕| 国产一级毛片在线| 国产成人91sexporn| 久久国产精品男人的天堂亚洲 | 成人亚洲欧美一区二区av| 女性被躁到高潮视频| 欧美亚洲 丝袜 人妻 在线| 亚洲成人一二三区av| 草草在线视频免费看| 最近最新中文字幕大全电影3| 国产亚洲欧美精品永久| 国产精品99久久久久久久久| 18禁裸乳无遮挡免费网站照片| 成人毛片a级毛片在线播放| 日韩一区二区三区影片| 又大又黄又爽视频免费| 青春草视频在线免费观看| 欧美性感艳星| 大陆偷拍与自拍| 中文乱码字字幕精品一区二区三区| 亚洲av国产av综合av卡| 欧美高清性xxxxhd video| 大码成人一级视频| 欧美+日韩+精品| 美女内射精品一级片tv| 交换朋友夫妻互换小说| 欧美成人精品欧美一级黄| 成人国产av品久久久| 欧美3d第一页| 亚洲av.av天堂| 亚洲三级黄色毛片| 亚洲精品久久午夜乱码| 午夜福利在线观看免费完整高清在| 久久久久国产精品人妻一区二区| 国产一区二区在线观看日韩| 成人综合一区亚洲| 国产精品精品国产色婷婷| 欧美激情极品国产一区二区三区 | 99精国产麻豆久久婷婷| 最近的中文字幕免费完整| 国产一区二区在线观看日韩| 日本欧美视频一区| 黄色怎么调成土黄色| 99视频精品全部免费 在线| 国产成人aa在线观看| 免费观看a级毛片全部| 国产 一区 欧美 日韩| 免费av中文字幕在线| 久久鲁丝午夜福利片| 日本欧美国产在线视频| 亚洲美女黄色视频免费看| 极品少妇高潮喷水抽搐| 久久ye,这里只有精品| 九草在线视频观看| 久久久久精品性色| 3wmmmm亚洲av在线观看| 亚洲三级黄色毛片| 亚洲,一卡二卡三卡| 人妻夜夜爽99麻豆av| 亚洲精品456在线播放app| 精品一区二区免费观看| 亚洲av中文字字幕乱码综合| av卡一久久| av在线蜜桃| 嘟嘟电影网在线观看| 亚州av有码| 亚洲,欧美,日韩| 一级毛片我不卡| 日韩一本色道免费dvd| 久久久久久人妻| 精品久久久久久久久av| 成人免费观看视频高清| 简卡轻食公司| 又大又黄又爽视频免费| 国产伦精品一区二区三区四那| 高清不卡的av网站| 91久久精品电影网| 免费黄网站久久成人精品| 国产精品免费大片| 国产黄片美女视频| 一本色道久久久久久精品综合| 欧美高清成人免费视频www| av国产久精品久网站免费入址| 日韩av在线免费看完整版不卡| 久久精品国产亚洲av天美| 激情五月婷婷亚洲| 国产成人freesex在线| 五月天丁香电影| 久久精品国产自在天天线| 联通29元200g的流量卡| 女人十人毛片免费观看3o分钟| 欧美精品一区二区大全| 国产爽快片一区二区三区| 少妇的逼好多水| 看十八女毛片水多多多| 日韩免费高清中文字幕av| 91aial.com中文字幕在线观看| 最近最新中文字幕大全电影3| 超碰av人人做人人爽久久| 欧美日韩综合久久久久久| 久久99热6这里只有精品| 爱豆传媒免费全集在线观看| 久久精品国产亚洲av涩爱| 国产精品久久久久成人av| 这个男人来自地球电影免费观看 | 国产爽快片一区二区三区| 多毛熟女@视频| av天堂中文字幕网| 天堂俺去俺来也www色官网| 中国三级夫妇交换| 国产成人freesex在线| 简卡轻食公司| 日本wwww免费看| 亚洲第一av免费看| 欧美zozozo另类| 国产乱人视频| 老师上课跳d突然被开到最大视频| 一级毛片久久久久久久久女| 日本黄大片高清| 成年av动漫网址| 91久久精品国产一区二区三区| 成人亚洲精品一区在线观看 | av网站免费在线观看视频| 九九爱精品视频在线观看| 黑人猛操日本美女一级片| 亚洲成人av在线免费| 国产中年淑女户外野战色| 国产精品一及| 久久av网站| 精品一区二区三区视频在线| 日本wwww免费看| 亚洲性久久影院| 成年人午夜在线观看视频| 国产精品秋霞免费鲁丝片| 另类亚洲欧美激情| 午夜福利影视在线免费观看| 免费观看a级毛片全部| 国产av一区二区精品久久 | 亚洲久久久国产精品| 国产精品熟女久久久久浪| 日韩不卡一区二区三区视频在线| 亚洲精华国产精华液的使用体验| 亚洲va在线va天堂va国产| 国产精品无大码| 一级毛片aaaaaa免费看小| 久久 成人 亚洲| 日产精品乱码卡一卡2卡三| 欧美一级a爱片免费观看看| 亚洲成人手机| 中文欧美无线码| 日韩一区二区视频免费看| 欧美日韩国产mv在线观看视频 | 亚洲av成人精品一区久久| 欧美zozozo另类| 欧美亚洲 丝袜 人妻 在线| 精品久久久噜噜| 嫩草影院入口| 欧美日韩视频精品一区| 亚洲人成网站在线观看播放| 婷婷色综合www| 亚洲欧洲国产日韩|