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

    不同狀態(tài)方程對(duì)固體介質(zhì)中斜激波反射的影響1)

    2017-11-11 01:55:24
    力學(xué)學(xué)報(bào) 2017年5期
    關(guān)鍵詞:馬赫狀態(tài)方程激波

    黃 蕭 于 鑫

    (北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所,北京100094)

    不同狀態(tài)方程對(duì)固體介質(zhì)中斜激波反射的影響1)

    黃 蕭 于 鑫2)

    (北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所,北京100094)

    相比氣體,固體介質(zhì)在高壓下的狀態(tài)方程更為復(fù)雜,形式也多種多樣.現(xiàn)有關(guān)于固體介質(zhì)中激波反射的理論研究,一般直接采用某種狀態(tài)方程,缺乏對(duì)采用不同狀態(tài)方程得到的結(jié)果的對(duì)比.本項(xiàng)工作采用激波極曲線的理論分析方法,選擇4種不同組合形式的狀態(tài)方程(一次沖擊激波采用線性的沖擊波速度與粒子速度關(guān)系式,二次沖擊激波采用Gr¨uneisen狀態(tài)方程;一次沖擊和二次沖擊激波均采用沖擊波速度與粒子速度關(guān)系式;一次沖擊激波采用線性沖擊波速度與粒子速度關(guān)系式,二次沖擊激波采用剛性氣體狀態(tài)方程;以及一次沖擊激波和二次沖擊激波均采用剛性氣體狀態(tài)方程),研究固體介質(zhì)中的斜激波反射,比較了采用不同組合形式的狀態(tài)方程對(duì)反射激波波后壓力的影響.利用量綱分析方法討論了簡(jiǎn)化狀態(tài)方程達(dá)到較高精度的條件.此外,用ANSYS/LS-DYNA軟件,對(duì)激波極曲線理論給出的結(jié)果進(jìn)行了驗(yàn)證.本項(xiàng)工作可為固體介質(zhì)中激波反射問題狀態(tài)方程的選取提供一定的指導(dǎo).

    斜激波反射,激波極曲線,狀態(tài)方程,固體介質(zhì)

    引言

    工程研究中經(jīng)常涉及固體介質(zhì)里的激波反射問題,如高能炸藥的爆轟波與固體結(jié)構(gòu)的相互作用、彈丸侵徹固體靶片、兩物體高速碰撞產(chǎn)生射流等[16].爆轟波在固體介質(zhì)中誘導(dǎo)產(chǎn)生的斜激波,或者物體高速斜碰撞產(chǎn)生的斜激波,在遇到固壁、交界面時(shí)會(huì)發(fā)生反射,產(chǎn)生較高的壓力.急劇升高的壓力可能在固體介質(zhì)中引發(fā)特殊的力學(xué)現(xiàn)象[78],因此開展相關(guān)研究對(duì)國(guó)防應(yīng)用及民用經(jīng)濟(jì)都具有重要的意義.

    氣體里的激波反射現(xiàn)象很早就受到關(guān)注[9].馬赫最早在實(shí)驗(yàn)中發(fā)現(xiàn)并記錄了激波反射現(xiàn)象,他記錄了兩種不同的激波反射結(jié)構(gòu),一種是正規(guī)反射,另一種是馬赫反射.之后馮·諾依曼針對(duì)正規(guī)反射和馬赫反射分別提出了“兩激波理論”和“三激波理論”,求解得到反射點(diǎn)附近各區(qū)域的流動(dòng)參數(shù).由于激波反射求解中邊界條件的約束是轉(zhuǎn)角和壓強(qiáng),在實(shí)際應(yīng)用中更為有效的是壓力–轉(zhuǎn)角極曲線[1014].例如反射激波極曲線與縱軸的交點(diǎn)是正規(guī)反射的解,而馬赫桿極曲線與反射激波極曲線的交點(diǎn)即為馬赫反射的解.通常根據(jù)激波的運(yùn)動(dòng)狀態(tài)又可以分為定常的反射結(jié)構(gòu)和非定常的反射結(jié)構(gòu).氣體中超音速來流繞楔角的流動(dòng)會(huì)產(chǎn)生駐定的激波反射結(jié)構(gòu),而運(yùn)動(dòng)激波沿楔面的運(yùn)動(dòng)會(huì)產(chǎn)生準(zhǔn)定常的反射結(jié)構(gòu),準(zhǔn)定常的反射形式比定常的更為豐富,例如有可能產(chǎn)生雙馬赫反射[15].

    相比氣體,固體介質(zhì)中斜激波反射的系統(tǒng)研究工作在公開發(fā)表的文獻(xiàn)中鮮有出現(xiàn).原因之一是固體介質(zhì)在高壓下的狀態(tài)方程較為復(fù)雜,可采用的狀態(tài)方程形式也多種多樣,例如沖擊波速度與波后粒子速度關(guān)系、Gr¨uneisen狀態(tài)方程、剛性氣體狀態(tài)方程等[1619].近年關(guān)于固體介質(zhì)中激波反射的研究,一般直接采用某種狀態(tài)方程,例如陳大偉等[20]采用剛性氣體狀態(tài)方程研究了凝聚介質(zhì)里的斜激波反射問題,Brown和Ravichandran[21]采用Gr¨uneisen狀態(tài)方程研究了固體介質(zhì)界面間的激波折射.可以看到這類研究缺乏對(duì)不同形式狀態(tài)方程得到的結(jié)果的對(duì)比,若能給出采用不同狀態(tài)方程得到結(jié)果之間的差異,對(duì)以后的理論或數(shù)值模擬具有重要的意義,例如能對(duì)采用某種狀態(tài)方程得到結(jié)果的精度有定量的認(rèn)識(shí).因此本工作基于激波極曲線的理論分析方法,比較在不同狀態(tài)方程下固體介質(zhì)中的斜激波反射,探討采用不同狀態(tài)方程對(duì)反射激波波后壓力精度的影響.

    1 理論模型

    1.1 狀態(tài)方程

    圖1描述了一個(gè)典型的馬赫反射結(jié)構(gòu).入射激波(OI)和馬赫桿(OT)是一次沖擊激波,而反射激波(OR)是二次沖擊激波.采用不同的狀態(tài)方程,得到的一次沖擊激波關(guān)系式和二次沖擊激波關(guān)系式會(huì)不同.對(duì)于一次沖擊,大量的實(shí)驗(yàn)研究表明,許多固體材料的沖擊波速度與粒子速度呈現(xiàn)線性關(guān)系[2223]

    圖1 馬赫反射結(jié)構(gòu)示意圖Fig.1 Schematic diagram of a Mach re fl ection

    其中,常壓下的聲速c0和λ均為材料參數(shù).對(duì)于二次沖擊,并沒有直接的實(shí)驗(yàn)性的經(jīng)驗(yàn)關(guān)系式.本文研究4種不同組合形式狀態(tài)方程下的斜激波反射問題.一次沖擊激波采用線性的沖擊波速度與粒子速度關(guān)系式,二次沖擊激波采用Gr¨uneisen狀態(tài)方程,用1Du2Gr¨uneisen表示;一次沖擊激波采用線性沖擊波速度與粒子速度關(guān)系式,二次沖擊激波采用二次的沖擊波速度與粒子速度關(guān)系式,用1Du2Du表示;一次沖擊激波采用線性沖擊波速度與粒子速度關(guān)系式、二次沖擊激波采用剛性氣體狀態(tài)方程,用1Du2Gama表示;一次沖擊激波和二次沖擊激波均采用剛性氣體狀態(tài)方程,用 1Gama2Gama表示.由于固體介質(zhì)參數(shù)的測(cè)定利用Gr¨uneisen狀態(tài)方程,一次沖擊利用沖擊波速度與粒子速度關(guān)系式.因此本文假設(shè)理論上1Du2Gr¨uneisen狀態(tài)方程得到的結(jié)果與實(shí)驗(yàn)結(jié)果是一致的,在后續(xù)的理論分析中,以1Du2Gr¨uneisen得到的結(jié)果作為參考值,比較另外3種狀態(tài)方程與其差異.

    1.2 激波極曲線關(guān)系式

    1.2.1 采用1Du2Gr¨uneisen狀態(tài)方程

    對(duì)于一次沖擊采用線性的沖擊波速度與粒子速度關(guān)系式,激波極曲線關(guān)系式為

    其中 p0和 p1分別表示入射激波波前和波后的壓力,對(duì)于固體介質(zhì),通常令p0=0.ρ0和q0分別表示入射激波波前的密度和來流速度(參考系固定在激波上),D1和u1分別表示入射激波速度和波后的粒子速度,θ1為入射激波波前、波后速度的轉(zhuǎn)角.對(duì)于二次沖擊激波,Brown和Ravichandran[21]給出采用Gr¨uneisen狀態(tài)方程時(shí)的壓力與比容(p-v)關(guān)系如下

    其中v0為初始狀態(tài)下的比容,v2和p2分別表示反射激波波后的比容和壓力,Γ為材料參數(shù),為簡(jiǎn)便起見,取 Γ/v= Γ0/v0及 Γ0=2λ?2/3[24].反射激波波前、波后速度的轉(zhuǎn)角θ2有如下關(guān)系式

    其中v1和q1分別表示反射激波波前的比容和來流速度(參考系固結(jié)在激波上).式(3)和式(4)構(gòu)成反射激波采用Gr¨uneisen狀態(tài)方程時(shí)的激波極曲線關(guān)系式.

    1.2.2 采用1Du2Du狀態(tài)方程

    一次沖擊的激波極曲線關(guān)系同式(2),二次沖擊的沖擊波速度和粒子速度的關(guān)系可設(shè)為

    其中 a1,b1和 c1為待定參數(shù),D2nr和 u2nr分別表示相對(duì)于u1n的反射激波速度和波后粒子速度,即D2nr=D2n?u1n和u2nr=u2n?u1n,其中u1n為入射激波波后粒子速度u1在垂直于反射激波方向上的分量.相對(duì)一次沖擊時(shí)的沖擊波速度與粒子速度關(guān)系式(1),式(5)的參數(shù)較難通過實(shí)驗(yàn)進(jìn)行標(biāo)定,Neal[25]和經(jīng)福謙[26]給出了一種把壓力--比容關(guān)系轉(zhuǎn)換為沖擊波速度--粒子速度關(guān)系的方法,基于此,待定參數(shù)a1,b1和c1可表示為

    若采用Gr¨uneisen狀態(tài)方程來表示二次沖擊時(shí)的p-v關(guān)系,結(jié)合式(3)和(6)可以得到a1,b1和c1的具體表達(dá)式,從而得到反射激波采用沖擊波速度與粒子速度關(guān)系式時(shí)的激波極曲線關(guān)系式

    實(shí)際上,此時(shí)反射激波采用沖擊波速度與粒子速度關(guān)系式等效于采用Gr¨uneisen狀態(tài)方程的二階近似.

    1.2.3 采用1Du2Gama狀態(tài)方程

    一次沖擊的激波極曲線關(guān)系同式(2),二次沖擊采用剛性狀態(tài)方程時(shí)的激波極曲線關(guān)系可表示為

    其中 γ為絕熱指數(shù).為方便不同狀態(tài)方程間的比較,絕熱指數(shù)γ取為等熵的常數(shù)近似γ=2λ?1+

    1.2.4 采用1Gama2Gama狀態(tài)方程

    此時(shí)所有激波關(guān)系式均采用剛性氣體狀態(tài)方程,對(duì)于一次沖擊的入射激波和馬赫桿,有如下激波極曲線關(guān)系

    對(duì)于二次沖擊的反射激波,激波極曲線關(guān)系同式(8).

    2 理論分析結(jié)果和討論

    2.1 不同形式狀態(tài)方程對(duì)激波極曲線的影響

    以金屬鈹(Be)為例,給定材料參數(shù):常壓下的聲速 c0=0.799cm/μs、初始密度 ρ0=1.85g/cm3以及系數(shù)λ=1.13,加載條件給定入射激波強(qiáng)度,即入射激波波后壓力p1=100GPa.圖2給出了不同入射激波傾角α0下的激波極曲線圖,當(dāng)α0=35°時(shí),從圖2(a)可以看出,此時(shí)發(fā)生正規(guī)反射,采用1Du2Du狀態(tài)方程的解和1Du2Gr¨uneisen最為接近,而1Du2Gama偏離1Du2Gr¨uneisen較多.

    當(dāng)入射激波傾角α0=45°時(shí),激波極曲線圖2(b)顯示了4種不同狀態(tài)方程下馬赫桿極曲線和反射激波極曲線的交點(diǎn),可以看出此時(shí)發(fā)生的是非正規(guī)反射,而且采用不同狀態(tài)方程得到的結(jié)果差異明顯.這種差異不僅反映在反射激波波后壓力上,也表現(xiàn)在反射模式上.當(dāng)采用1Du2Gama和1Gama2Gama狀態(tài)方程時(shí),交點(diǎn)在反射激波極曲線頂點(diǎn)的左側(cè),發(fā)生的是馬赫反射(MR),而采用1Du2Gr¨uneisen和1Du2Du狀態(tài)方程時(shí),交點(diǎn)在反射激波極曲線頂點(diǎn)的右側(cè),此時(shí)發(fā)生的是馮·諾依曼反射(vNR)[11,2931].

    圖2 不同入射激波傾角下的激波極曲線Fig.2 Shock polar under di ff erent incident angles

    2.2 不同形式狀態(tài)方程對(duì)反射激波波后壓力的影響

    改變?nèi)肷浼げ▋A角 α0,可以得到在給定激波強(qiáng)度條件下,反射波波后壓力 p2隨 α0的變化.圖3描繪了不同狀態(tài)方程對(duì)反射激波波后壓力的影響.從圖中可以看出,采用1Du2Du狀態(tài)方程的解和1Du2Gr¨uneisen最為接近,而1Du2Gama的解和1Du2Gr¨uneisen偏離最遠(yuǎn).而且,總體上正規(guī)反射階段4種狀態(tài)方程得出的結(jié)果差異較小,非正規(guī)反射階段得出的結(jié)果差異較為明顯.

    圖3 不同狀態(tài)方程對(duì)反射激波波后壓力的影響Fig.3 E ff ect of di ff erent equations of state on the pressure behind the re fl ected shock

    此外,圖3表明隨著α0增大到某一角度時(shí),p2的值會(huì)發(fā)生突變,此時(shí)對(duì)應(yīng)的角度即為正規(guī)反射過渡到非正規(guī)反射的臨界角.當(dāng)入射激波強(qiáng)度較大(如p1=100GPa)時(shí),采用不同狀態(tài)方程得出的臨界角在37°到42°之間,而強(qiáng)度較小(如 p1=15GPa)的入射激波給出的臨界角在47°到57°之間,此時(shí)采用不同狀態(tài)方程得到的結(jié)果之間差距更大.

    2.3 不同形式狀態(tài)方程對(duì)反射類型的影響

    根據(jù)反射激波脫體的判據(jù),由激波極曲線理論可以得到不同入射激波強(qiáng)度下從正規(guī)反射過渡到非正規(guī)反射的臨界角,由此可將激波反射類型劃分為兩個(gè)區(qū)域:正規(guī)反射區(qū)和非正規(guī)反射區(qū),如圖4所示.橫軸為代表入射激波強(qiáng)度的入射激波波后壓力,縱軸為入射激波與來流的角度.從圖中可以看出,入射激波強(qiáng)度越小,從正規(guī)反射過渡到非正規(guī)反射的臨界角越大.此外,圖中表明采用 4種不同狀態(tài)方程得到的區(qū)域分界線差異較為明顯,由1Du2Du和1Du2Gr¨uneisen得到的結(jié)果基本一致,由1Gama2Gama得到的結(jié)果與1Du2Gr¨uneisen有所差異,1Du2Gama得到的結(jié)果與1Du2Gr¨uneisen差異最大.由此可見,采用不同狀態(tài)方程對(duì)反射激波波后壓力以及反射類型的判定都有較大的影響,在實(shí)際應(yīng)用中需要加以注意.

    圖4 激波反射區(qū)域類型的劃分Fig.4 Domains of di ff erent types of shock wave re fl ection

    2.4 量綱分析

    前面的小節(jié)給出了特定材料參數(shù)下的結(jié)果,接下來通過量綱分析研究更一般的結(jié)果,并基于此給出簡(jiǎn)化形式的剛性氣體狀態(tài)方程的適用條件.對(duì)于本文研究的問題,反射激波波后壓力是以下獨(dú)立參數(shù)的函數(shù)

    即p2只與初始密度ρ0、聲速c0、一次沖擊D?u關(guān)系式的系數(shù)λ、代表入射激波強(qiáng)度的p1以及入射激波的傾角α0有關(guān).這5個(gè)參數(shù)只存在兩個(gè)具有獨(dú)立量綱的參數(shù),不失一般性,選為ρ0和c0,其他參數(shù)的量綱可以由這兩個(gè)基本參數(shù)的量綱導(dǎo)出,利用量綱分析中的π定理,式(10)可以重新寫為如下的無量綱函數(shù)形式

    即無量綱化的反射激波波后壓力只與無量綱化的入射激波波后壓力、λ以及α0相關(guān),無量綱化后可將原本需要討論的5個(gè)參數(shù)降為3個(gè),很大程度上簡(jiǎn)化了問題.

    圖5 正規(guī)反射階段不同狀態(tài)方程解的誤差Fig.5 Solution error of regular re fl ection under di ff erent equations of state

    當(dāng)入射激波的傾角 α0較小時(shí),系統(tǒng)產(chǎn)生正規(guī)反射,在正規(guī)反射階段反射激波波后壓力隨傾角增大變化不大.給定 α0=20°,討論不同的 λ值下,隨的變化,并比較4種不同形式狀態(tài)方程的影響,此時(shí)結(jié)果可近似表示正規(guī)反射階段不同狀態(tài)方程對(duì)無量綱化的反射激波波后壓力的影響.以1Du2Gr¨uneisen得到的值為標(biāo)準(zhǔn),比較另外3種狀態(tài)方程得到的值相對(duì)其的誤差η,η隨的變化趨勢(shì)如圖5所示.從圖中可以看出,對(duì)于1Du2Du,大多數(shù)情況下得到的結(jié)果都很精確,因?yàn)槠涫?Du2Gr¨uneisen解的二階近似.對(duì)于1Du2Gama,在λ為1.5~1.8時(shí)得到的結(jié)果與1Du2Gr¨uneisen符合得較好,λ為其他值時(shí)精度較差.而對(duì)于1Gama2Gama,在λ為1.2~1.3時(shí)得到的結(jié)果與1Du2Gr¨uneisen符合得較好;λ為其他值時(shí)只有

    圖5 正規(guī)反射階段不同狀態(tài)方程解的誤差(續(xù))Fig.5 Solution error of regular re fl ection under di ff erent equations of state(continued)

    圖6 非正規(guī)反射階段不同狀態(tài)方程解的誤差Fig.6 Solution error of irregular re fl ection under di ff erent equations of state

    當(dāng)入射激波傾角較大時(shí),系統(tǒng)有可能出現(xiàn)非正規(guī)反射,在正規(guī)反射和非正規(guī)反射過渡的階段,不同狀態(tài)方程對(duì)結(jié)果的影響較為劇烈.圖6(a)和圖6(b)分別描繪了在入射激波傾角α0=35°及α0=40°條件下的誤差隨的變化趨勢(shì),此時(shí)給定λ=1.2.從圖中可以看出,隨著無量綱化入射激波強(qiáng)度增大,不同狀態(tài)方程得到解的誤差都會(huì)出現(xiàn)突然增大,而在一定區(qū)間后又會(huì)突然減小的現(xiàn)象.此區(qū)間對(duì)應(yīng)著系統(tǒng)的解從正規(guī)反射到非正規(guī)反射的過渡區(qū),該過渡區(qū)也可以從圖4加以解釋,即為以α0等于定值的直線交于不同狀態(tài)方程得到的區(qū)域分界線的區(qū)間段,從圖中可以看出1Du2Gama的過渡區(qū)最大.在過渡區(qū)內(nèi)另外3種狀態(tài)方程的相對(duì)誤差都較大(η接近50%).而當(dāng)?shù)闹敌∮诖藚^(qū)間時(shí),系統(tǒng)的解為正規(guī)反射,從圖中可以看出1Du2Du得到的結(jié)果很精確,而1Du2Gama的結(jié)果誤差最大,約為10%.當(dāng)?shù)闹荡笥诖藚^(qū)間時(shí),系統(tǒng)的解為非正規(guī)反射,從圖中可以看出1Du2Du得到的結(jié)果很精確,而1Gama2Gama的結(jié)果誤差最大,約為10%.

    綜上所述,對(duì)于大多數(shù)固體(λ=1.2左右),無論是求解正規(guī)反射還是非正規(guī)反射問題,采用簡(jiǎn)化形式的剛性狀態(tài)方程(1Gama2Gama或1Du2Gama)均會(huì)帶來一定的誤差.采用1Gama2Gama得到的解誤差在正規(guī)反射問題時(shí)比1Du2Gama小,誤差在5%以內(nèi);而采用1Du2Gama得到的解誤差在非正規(guī)反射問題時(shí)比1Gama2Gama小,誤差在5%以內(nèi).在正規(guī)反射到非正規(guī)反射的轉(zhuǎn)變階段,1Du2Gama得到的解誤差在過渡區(qū)最大,建議相關(guān)理論研究工作慎重采用.

    3 數(shù)值模擬

    為驗(yàn)證激波極曲線理論得出的結(jié)果,利用有限元程序ANSYS/LS-DYNA對(duì)不同入射激波傾角下的激波反射進(jìn)行模擬.由于固體介質(zhì)中的斜激波反射一般非定常,因此建立運(yùn)動(dòng)激波過楔面發(fā)生反射的數(shù)值模型.計(jì)算模型如圖7所示,左側(cè)給定入射激波波后壓力的邊界條件,上邊、下邊和斜邊施加滑移固壁邊界條件,模型采用固定的歐拉網(wǎng)格.固體介質(zhì)采用Gr¨uneisen狀態(tài)方程,材料參數(shù)同3.1節(jié).圖7中α為入射激波和楔面的夾角,由于非正規(guī)反射中來流方向平行于三波點(diǎn)軌跡線,而非平行于固壁,因此當(dāng)運(yùn)動(dòng)激波發(fā)生非正規(guī)反射時(shí),激波極曲線理論需要相應(yīng)地考慮三波點(diǎn)軌跡角的影響,此時(shí)有α=α0+δ.

    圖7 激波反射數(shù)值模擬的有限元模型Fig.7 Finite element model for numerical simulations of shock wave re fl ection

    隨著入射激波和楔面夾角 α的增大,數(shù)值模擬得到的斜激波反射波后狀態(tài)及反射模式將發(fā)生改變,圖 8(a)~ 圖 8(c)描繪了 α =35°,α =45°和α=66°時(shí)數(shù)值模擬給出的斜激波反射模式以及相應(yīng)的激波極曲線理論解.

    圖8 不同入射激波傾角下的密度等值線圖及相應(yīng)的激波極曲線Fig.8 Density contours and the corresponding shock polar under di ff erent incident angles

    當(dāng)α=35°時(shí),理論分析中反射激波極曲線和縱軸有交點(diǎn),說明此時(shí)發(fā)生正規(guī)反射,數(shù)值模擬的圖8(a)給出了由入射激波和反射激波組成的二波結(jié)構(gòu),即正規(guī)反射結(jié)構(gòu).當(dāng)α=45°時(shí),理論分析中馬赫桿極曲線和反射激波極曲線的交點(diǎn)位于反射激波極曲線頂點(diǎn)的左側(cè),說明此時(shí)發(fā)生馬赫反射.從圖8(b)可以清晰地看到由入射激波、弧形的反射激波和馬赫桿組成的三波結(jié)構(gòu),因此結(jié)果為馬赫反射.繼續(xù)增大α到66°,激波極曲線理論表明狀態(tài)方程采用1Du2Gr¨uneisen和1Du2Du得到的交點(diǎn)位于反射激波極曲線頂點(diǎn)的右側(cè),即 vNR 反射,而采用1Du2Gama和1Gama2Gama狀態(tài)方程得到的交點(diǎn)位于反射激波極曲線頂點(diǎn)的左側(cè),因此仍然為馬赫反射.從圖8(c)數(shù)值模擬的結(jié)果可以看出,此時(shí)發(fā)生的是vNR反射,相比馬赫反射,vNR反射結(jié)構(gòu)中馬赫桿和入射激波之間光滑連接,而且三波點(diǎn)附近的反射波是連續(xù)的壓縮波.?dāng)?shù)值模擬的結(jié)果表明與采用1Du2Gr¨uneisen和1Du2Du狀態(tài)方程得到的結(jié)果更為一致.

    圖8 不同入射激波傾角下的密度等值線圖及相應(yīng)的激波極曲線(續(xù))Fig.8 Density contours and the corresponding shock polar under di ff erent incident angles(continued)

    4 結(jié)論

    本文利用激波極曲線的理論分析方法,研究不同形式狀態(tài)方程對(duì)固體介質(zhì)中斜激波反射的影響.結(jié)果表明:

    (1)不同的狀態(tài)方程對(duì)反射激波波后壓力以及對(duì)反射類型的判別都有較大的影響,在實(shí)際應(yīng)用中需要加以注意.

    (2)采用沖擊波速度和波后粒子速度關(guān)系或Gr¨uneisen狀態(tài)方程得到的結(jié)果比采用剛性氣體狀態(tài)方程更為精確.

    (3)通過量綱分析討論了采用剛性氣體狀態(tài)方程能達(dá)到較高精度的條件:對(duì)于大多數(shù)固體(λ=1.2左右),在求解正規(guī)反射問題時(shí),一次和二次激波均采用剛性氣體狀態(tài)方程得到的解誤差比一次采用線性沖擊波速度與粒子速度關(guān)系式、二次采用剛性氣體狀態(tài)方程小,誤差在5%以內(nèi);而在求解非正規(guī)反射問題時(shí),一次激波采用線性沖擊波速度與粒子速度關(guān)系式、二次激波采用剛性氣體狀態(tài)方程得到的解誤差比一次和二次激波均采用剛性氣體狀態(tài)方程小,誤差在5%以內(nèi).

    1 Loomis E,Swift D.Oblique shock waves incident on an interface between two materials for general equations of state.Journal of Applied Physics,2008,103(2):023518

    2馬天寶,任會(huì)蘭,李健等.爆炸與沖擊問題的大規(guī)模高精度計(jì)算.力學(xué)學(xué)報(bào),2016,48(3):599-608(Ma Tianbao,Ren Huilan,Li Jian,et al.Large scale high precision computation for explosion and impactproblems.ChineseJournalofTheoreticalandAppliedMechanics,2016,48(3):599-608(in Chinese))

    3王成,王萬軍,寧建國(guó).聚能裝藥對(duì)混凝土靶板的侵徹研究.力學(xué)學(xué)報(bào),2015,47(4):672-686(Wang Cheng,Wang Wanjun,Ning Jianguo.Investigation on shaped charge penetration into concretetargets.Chinese Journal of Theoretical and Applied Mechanics,2015,47(4):672-686(in Chinese))

    4于明,孫宇濤,劉全.爆轟波在炸藥--金屬界面上的折射分析.物理學(xué)報(bào),2015,64(11):282-291(Yu Ming,Sun Yutao,Liu Quan.Analysis on refraction of detonation wave at the explosive-metal interface.Acta Physica Sinica,2015,47(4):672-686(in Chinese))

    5于明,劉全.凝聚炸藥爆轟波在高聲速材料界面上的折射現(xiàn)象分析.物理學(xué)報(bào),2016,65(2):235-244(Yu Ming,Liu Quan.Refraction of detonation wave at interface between condensed explosives and high sound-speed material.Acta Physica Sinica,2016,65(2):235-244(in Chinese))

    6魏晨慧,朱萬成,白羽等.不同節(jié)理角度和地應(yīng)力條件下巖石雙孔爆破的數(shù)值模擬.力學(xué)學(xué)報(bào),2016,48(4):926-935(Wei Chenhui,Zhu Wancheng,Bai Yu,et al.Numerical simulation on two-hole blasting of rock under di ff erent joint angles and in-situ stress conditions.Chinese Journal of Theoretical and Applied Mechanics,2016,48(4):926-935(in Chinese))

    7 Singh M,Suneja HR,Bola MS,et al.Dynamic tensile deformation and fracture of metal cylinders at high strain rates.International Journal of Impact Engineering,2002,27(9):939-954

    8張崇玉,胡海波,李慶忠等.爆轟波對(duì)碰驅(qū)動(dòng)下平面鉛飛層對(duì)碰區(qū)動(dòng)載行為實(shí)驗(yàn)研究.高壓物理學(xué)報(bào),2009,23(4):283-287(Zhang Chongyu,Hu Haibo,Li Qingzhong,et al.Experimental study on dynamic behavior of lead plate driven by two head-on colliding detonation waves.Chinese Journal of High Pressure Physics,2009,23(4):283-287(in Chinese))

    9 Hornung H.Regular and Mach re fl ection of shock waves.Annual Review of Fluid Mechanics,1986,18:33-58

    10 Henderson LF,Colella P,Puckett EG.On the refraction of shock waves at a slow–fast gas interface.Journal of Fluid Mechanics,1991,224:1-27

    11 Ben-Dor G.Shock Wave Re fl ection Phenomena. New York:Springer,2007

    12 Gvozdeva L,Gavrenkov S,Nesterov A.A study of slipstreams in triple shock wave con fi gurations.Shock Waves,2015,25(3):283-291

    13 Matheis J,Hickel S.On the transition between regular and irregular shock patterns of shock-wave/boundary-layer interactions.Journal of Fluid Mechanics,2015,776:200-234

    14 Kobayashi S,Adachi T.Shock-tube experiments on the stability of regular re fl ection in the dual-solution domain.Shock Waves,2016,26(3):263-267

    15 Li H,Ben-Dor G.Reconsideration of pseudo-steady shock wave re fl ections and the transition criteria between them.Shock Waves,1995,5(1):59-73

    16 Meniko ffR,Plohr BJ.The Riemann problem for fl uid fl ow of real materials.Reviews of Modern Physics,1989,61(1):75-130

    17王繼海.二維非定常流和激波.北京:科學(xué)出版社,1994(Wang Jihai.Two-Dimensional Unsteady Flow and Shock Waves.Beijing:Science Press,1994(in Chinese))

    18 Davison L.Fundamentals of Shock Wave Propagation in Solids.Berlin:Springer,2008

    19 McCoy CA,Gregor MC,Polsin DN,et al.Shock-wave equationof-state measurements in fused silica up to 1600GPa.Journal of Applied Physics,2016,119(21):215901

    20陳大偉,秦承森,王裴等.凝聚介質(zhì)中斜激波的反射.計(jì)算物理,2011,28(6):791-796(Chen Dawei,Qin Chengsen,Wang Pei,et al.Oblique shock wave re fl ection in condensed matter.Chinese Journal of Computational Physics,2011,28(6):791-796(in Chinese))

    21 Brown JL,Ravichandran G.Analysis of oblique shock waves in solids using shock polars.Shock Waves,2014,24(4):403-413

    22 Ruo ff AL.Linear shock-velocity-particle-velocity relationship.Journal of Applied Physics,1967,38(13):4976-4980

    23 Neal T.Mach waves and re fl ected rarefactions in aluminum.Journal of Applied Physics,1975,46(6):2521-2527

    24譚華.實(shí)驗(yàn)沖擊波物理導(dǎo)引.北京:國(guó)防工業(yè)出版社,2007(Tan Hua.Introduction to Experimental Shock-Wave Physics.Beijing:National Defense Industry Press,2007(in Chinese))

    25 Neal T.Second Hugoniot relationship for solids.Journal of Physics and Chemistry of Solids,1977,38(3):225-231

    26經(jīng)福謙.沖擊波速度–粒子速度關(guān)系式的一個(gè)簡(jiǎn)單推導(dǎo)及其直線表達(dá)式適用范圍的討論.爆炸與沖擊,1982(3):17-24(Jing Fuqian. A simpli fi ed derivation of the shock-velocity-particlevelocity relationship and discussion on the applicable range of its linearized expression.Explosion and Shock Waves,1982(3):17-24(in Chinese))

    27李維新.凝聚介質(zhì)的簡(jiǎn)化狀態(tài)方程.爆炸與沖擊,1983(2):30-38(Li Weixin.A simpli fi ed equation of state in condensed medium.Explosion and Shock Waves,1983(2):30-38(in Chinese))

    28李維新.一維不定常流與沖擊波.北京:國(guó)防工業(yè)出版社,2003(Li Weixin.One-dimensional Nonsteady Flow and Shock Waves.Beijing:National Defence Industry Press,2003(in Chinese))

    29 Ben-Dor G,Takayama K.The phenomena of shock wave re fl ection—A review of unsolved problems and future research needs.Shock Waves,1992,2(4):211-223

    30 Kobayashi S,Adachi T,Suzuki T.Non-self-similar behavior of the von Neumann re fl ection.Physics of Fluids,2000,12(7):1869-1877

    31 Wu F,Dai H,Kong D.Mechanism for the transition from a regular re fl ection to a mach re fl ection or a von neumann re fl ection.Acta Mathematica Scientia,2016,36(3):931-944

    EFFECTS OF DIFFERENT EQUATIONS OF STATE ON THE OBLIQUE SHOCK WAVE REFLECTION IN SOLIDS1)

    Huang Xiao Yu Xin2)
    (Institute of Applied Physics and Computational Mathematics,Beijing 100094,China)

    The equations of state of solids under high pressure are more complicated than that of gases in a variety of forms.While the existing investigations on the oblique shock wave re fl ection usually take one of the equations of state,lacking of the comparisons among them.Therefore,this paper aims at the oblique shock wave re fl ection in solids through shock polar methodology under four di ff erent forms of equations of state(principal shock taking with linear shock–particle velocity relationship and second shock taking with Gr¨uneisen equation of state,principal and second shock both taking with shock–particle velocity relationship,principal shock taking with linear shock–particle velocity relationship and second shock taking with sti ff ened gas equation of state,and principal and second shock both taking with sti ff ened gas equation of state).The e ff ects of di ff erent equations of state on the pressure behind the re fl ected shock wave are discussed.By conducting the dimensional analysis,we provide an applicable condition for employing a simpli fi ed equation of state to achieve high accuracy.Moreover,numerical simulations performed by ANSYS/LS-DYNA software are conducted to validate the results through shock polar methodology.The results of this paper could be helpful for the decision of the equation of state on the oblique shock wave re fl ection in solids.

    oblique shock wave re fl ection,shock polar,equation of state,solid

    O357.4+2

    A

    10.6052/0459-1879-17-015

    2017–01–09收稿,2017–08–23 錄用,2017–08–23 網(wǎng)絡(luò)版發(fā)表.

    1)國(guó)家自然科學(xué)基金資助項(xiàng)目(11602027).

    2)于鑫,研究員,主要研究方向:爆炸與沖擊動(dòng)力學(xué).E-mail:yummuy@iapcm.ac.cn

    黃蕭,于鑫.不同狀態(tài)方程對(duì)固體介質(zhì)中斜激波反射的影響.力學(xué)學(xué)報(bào),2017,49(5):1145-1153

    Huang Xiao,Yu Xin.E ff ects of di ff erent equations of state on the oblique shock wave re fl ection in solids.Chinese Journal of Theoretical and Applied Mechanics,2017,49(5):1145-1153

    猜你喜歡
    馬赫狀態(tài)方程激波
    東風(fēng)風(fēng)行T5馬赫版
    汽車觀察(2022年12期)2023-01-17 02:19:58
    LKP狀態(tài)方程在天然氣熱物性參數(shù)計(jì)算的應(yīng)用
    煤氣與熱力(2021年6期)2021-07-28 07:21:30
    一種基于聚類分析的二維激波模式識(shí)別算法
    基于HIFiRE-2超燃發(fā)動(dòng)機(jī)內(nèi)流道的激波邊界層干擾分析
    穿越“馬赫谷”
    27馬赫,刺破蒼穹
    斜激波入射V形鈍前緣溢流口激波干擾研究
    適于可壓縮多尺度流動(dòng)的緊致型激波捕捉格式
    基于隨機(jī)與區(qū)間分析的狀態(tài)方程不確定性比較
    用狀態(tài)方程模擬氨基酸水溶液的熱力學(xué)性質(zhì)
    每晚都被弄得嗷嗷叫到高潮| 18禁黄网站禁片午夜丰满| 亚洲不卡免费看| 9191精品国产免费久久| 悠悠久久av| 黄片小视频在线播放| 69人妻影院| 日本撒尿小便嘘嘘汇集6| 观看免费一级毛片| 亚洲 欧美 日韩 在线 免费| 五月玫瑰六月丁香| 精品一区二区三区视频在线 | 波多野结衣高清无吗| 中文字幕人妻丝袜一区二区| 毛片女人毛片| 国产色爽女视频免费观看| 波野结衣二区三区在线 | 亚洲成人久久爱视频| 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 日韩国内少妇激情av| 人人妻人人澡欧美一区二区| 免费看日本二区| 18+在线观看网站| 12—13女人毛片做爰片一| 欧美黄色片欧美黄色片| 午夜精品久久久久久毛片777| 99热这里只有是精品50| 亚洲电影在线观看av| 九九久久精品国产亚洲av麻豆| 日日干狠狠操夜夜爽| 首页视频小说图片口味搜索| 亚洲最大成人中文| 男人舔女人下体高潮全视频| 欧美xxxx黑人xx丫x性爽| 亚洲人成网站在线播| 在线十欧美十亚洲十日本专区| 国产真实伦视频高清在线观看 | 欧美激情在线99| 久久久国产成人精品二区| 午夜日韩欧美国产| 精品久久久久久,| 深爱激情五月婷婷| 欧美成人免费av一区二区三区| 无限看片的www在线观看| 无遮挡黄片免费观看| 偷拍熟女少妇极品色| 国产视频内射| 757午夜福利合集在线观看| 亚洲国产精品999在线| 少妇裸体淫交视频免费看高清| 欧美色欧美亚洲另类二区| 午夜免费男女啪啪视频观看 | 亚洲人成网站在线播| 激情在线观看视频在线高清| 18禁黄网站禁片免费观看直播| 香蕉丝袜av| 美女 人体艺术 gogo| 久久久久久久久中文| 此物有八面人人有两片| 久久精品综合一区二区三区| 欧美黄色淫秽网站| 亚洲色图av天堂| 亚洲片人在线观看| 亚洲精品影视一区二区三区av| 国产精品久久电影中文字幕| 国产精品一区二区三区四区免费观看 | 最近在线观看免费完整版| 最新美女视频免费是黄的| 中文亚洲av片在线观看爽| 一进一出好大好爽视频| 高清毛片免费观看视频网站| 中文字幕人妻熟人妻熟丝袜美 | 热99re8久久精品国产| 麻豆一二三区av精品| 两个人的视频大全免费| 国产精品自产拍在线观看55亚洲| 最新在线观看一区二区三区| 国产精品1区2区在线观看.| 亚洲美女视频黄频| a在线观看视频网站| 观看免费一级毛片| 男人舔奶头视频| 极品教师在线免费播放| svipshipincom国产片| 在线观看66精品国产| 国产欧美日韩一区二区三| 18禁裸乳无遮挡免费网站照片| 最近最新中文字幕大全电影3| 欧美三级亚洲精品| 乱人视频在线观看| 国产精品 欧美亚洲| 精品不卡国产一区二区三区| 老汉色av国产亚洲站长工具| 女人被狂操c到高潮| www日本在线高清视频| 国产不卡一卡二| 亚洲精华国产精华精| 最新美女视频免费是黄的| 人妻夜夜爽99麻豆av| 麻豆国产97在线/欧美| 国产麻豆成人av免费视频| 精品久久久久久久人妻蜜臀av| 午夜免费成人在线视频| 亚洲aⅴ乱码一区二区在线播放| 亚洲avbb在线观看| 麻豆久久精品国产亚洲av| 少妇的逼好多水| 精品无人区乱码1区二区| 老司机深夜福利视频在线观看| 久久精品国产自在天天线| 国产精品 国内视频| 国产aⅴ精品一区二区三区波| 午夜福利在线观看免费完整高清在 | 久久久国产精品麻豆| av片东京热男人的天堂| 一级黄色大片毛片| 90打野战视频偷拍视频| 亚洲久久久久久中文字幕| 三级男女做爰猛烈吃奶摸视频| 免费电影在线观看免费观看| 国产日本99.免费观看| 99国产精品一区二区三区| 婷婷六月久久综合丁香| 国产亚洲精品久久久com| 国产成人av激情在线播放| 久久精品人妻少妇| 99久久无色码亚洲精品果冻| 国产精品电影一区二区三区| av中文乱码字幕在线| 精品免费久久久久久久清纯| 久久人妻av系列| 久久久久久九九精品二区国产| 亚洲av日韩精品久久久久久密| 最新在线观看一区二区三区| 亚洲精品在线美女| 亚洲av一区综合| 蜜桃久久精品国产亚洲av| 精品一区二区三区视频在线观看免费| 国产精品女同一区二区软件 | 国产69精品久久久久777片| 免费看光身美女| 国产私拍福利视频在线观看| 色视频www国产| 老熟妇乱子伦视频在线观看| 国产三级黄色录像| 亚洲人成网站在线播| 国产99白浆流出| 国产一区二区三区在线臀色熟女| 国产av一区在线观看免费| 三级国产精品欧美在线观看| 久久久久久久精品吃奶| 尤物成人国产欧美一区二区三区| 欧美一区二区精品小视频在线| 免费观看人在逋| 成人18禁在线播放| 精品人妻1区二区| 91久久精品国产一区二区成人 | 久久99热这里只有精品18| 国产亚洲精品久久久com| 中文字幕av在线有码专区| 美女大奶头视频| 欧美成人一区二区免费高清观看| 中文资源天堂在线| 97超视频在线观看视频| 少妇人妻精品综合一区二区 | 精品一区二区三区视频在线观看免费| 露出奶头的视频| 欧美av亚洲av综合av国产av| 最近在线观看免费完整版| 好看av亚洲va欧美ⅴa在| 男女做爰动态图高潮gif福利片| 免费观看的影片在线观看| 午夜免费男女啪啪视频观看 | 丁香六月欧美| 97人妻精品一区二区三区麻豆| 亚洲五月婷婷丁香| 国产免费男女视频| 午夜免费观看网址| 精品无人区乱码1区二区| 成年版毛片免费区| a级毛片a级免费在线| 欧美极品一区二区三区四区| 欧美+亚洲+日韩+国产| 亚洲七黄色美女视频| 日本一本二区三区精品| 国产不卡一卡二| 亚洲国产色片| av福利片在线观看| 亚洲av电影在线进入| 熟女少妇亚洲综合色aaa.| 麻豆久久精品国产亚洲av| 午夜亚洲福利在线播放| 免费人成在线观看视频色| 在线视频色国产色| 亚洲成人久久爱视频| 少妇高潮的动态图| www.999成人在线观看| 美女高潮的动态| 中文亚洲av片在线观看爽| 少妇的逼好多水| 成人亚洲精品av一区二区| 12—13女人毛片做爰片一| 91在线精品国自产拍蜜月 | 久久久久久久午夜电影| 国产在视频线在精品| 此物有八面人人有两片| tocl精华| 久久国产精品影院| 国产高清有码在线观看视频| 少妇的逼水好多| 久久亚洲精品不卡| 老司机午夜十八禁免费视频| 老汉色∧v一级毛片| 性色av乱码一区二区三区2| 国产成人a区在线观看| 国产三级在线视频| 给我免费播放毛片高清在线观看| 成人国产综合亚洲| 亚洲性夜色夜夜综合| 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 午夜a级毛片| av女优亚洲男人天堂| 色老头精品视频在线观看| 久久精品影院6| 成熟少妇高潮喷水视频| 男女午夜视频在线观看| 熟女人妻精品中文字幕| 精品一区二区三区视频在线观看免费| av视频在线观看入口| 午夜精品一区二区三区免费看| 搡老岳熟女国产| 久久久久九九精品影院| 婷婷亚洲欧美| 制服丝袜大香蕉在线| 看免费av毛片| 97超视频在线观看视频| 欧美国产日韩亚洲一区| 亚洲国产精品sss在线观看| 真人一进一出gif抽搐免费| 国产毛片a区久久久久| 亚洲精品国产精品久久久不卡| 精品日产1卡2卡| 国产极品精品免费视频能看的| 男女之事视频高清在线观看| 9191精品国产免费久久| 成年女人永久免费观看视频| 99热这里只有是精品50| 欧美午夜高清在线| 成人鲁丝片一二三区免费| 三级国产精品欧美在线观看| 嫁个100分男人电影在线观看| 免费在线观看日本一区| 久久中文看片网| 狂野欧美激情性xxxx| 亚洲电影在线观看av| 女人十人毛片免费观看3o分钟| 国内精品久久久久精免费| 热99在线观看视频| 久久草成人影院| 日韩欧美 国产精品| 欧美成人a在线观看| 欧美乱码精品一区二区三区| 女生性感内裤真人,穿戴方法视频| 久久久久精品国产欧美久久久| 精品久久久久久,| 亚洲人成电影免费在线| 999久久久精品免费观看国产| 午夜亚洲福利在线播放| 国产午夜精品论理片| 九色成人免费人妻av| svipshipincom国产片| 无遮挡黄片免费观看| 日韩有码中文字幕| 精品国产亚洲在线| 99riav亚洲国产免费| 成人高潮视频无遮挡免费网站| 淫秽高清视频在线观看| 可以在线观看毛片的网站| 久久婷婷人人爽人人干人人爱| 女人被狂操c到高潮| 国内少妇人妻偷人精品xxx网站| 久久久久国内视频| 亚洲色图av天堂| 超碰av人人做人人爽久久 | 精品国产美女av久久久久小说| 又爽又黄无遮挡网站| 亚洲 国产 在线| 久久午夜亚洲精品久久| 午夜福利在线观看免费完整高清在 | 老鸭窝网址在线观看| 国产伦一二天堂av在线观看| or卡值多少钱| 午夜视频国产福利| 国产三级黄色录像| 免费一级毛片在线播放高清视频| 免费高清视频大片| 亚洲av成人精品一区久久| 黄色女人牲交| 亚洲精品乱码久久久v下载方式 | xxx96com| 午夜免费观看网址| 中文在线观看免费www的网站| 国产精品香港三级国产av潘金莲| 99久国产av精品| 2021天堂中文幕一二区在线观| 最近在线观看免费完整版| 国产亚洲精品久久久久久毛片| 日韩av在线大香蕉| 欧美日韩国产亚洲二区| 黄色视频,在线免费观看| 国产激情偷乱视频一区二区| 88av欧美| 高潮久久久久久久久久久不卡| 国产高清激情床上av| 亚洲第一欧美日韩一区二区三区| 亚洲精品一区av在线观看| 好男人在线观看高清免费视频| 国产亚洲精品久久久久久毛片| 老鸭窝网址在线观看| 亚洲最大成人中文| 欧美激情在线99| 国产精品亚洲美女久久久| 19禁男女啪啪无遮挡网站| 两个人看的免费小视频| 美女高潮喷水抽搐中文字幕| 天堂网av新在线| 18禁国产床啪视频网站| 久久久久精品国产欧美久久久| 一区福利在线观看| 99热这里只有精品一区| 一级a爱片免费观看的视频| 嫩草影院精品99| 免费人成在线观看视频色| 美女黄网站色视频| 观看免费一级毛片| 亚洲18禁久久av| 国产 一区 欧美 日韩| 黄色女人牲交| 国产 一区 欧美 日韩| 国产精品久久久久久精品电影| 内地一区二区视频在线| 久久亚洲精品不卡| 人妻久久中文字幕网| 少妇人妻一区二区三区视频| 99久国产av精品| 成人三级黄色视频| 国产v大片淫在线免费观看| 国产真实伦视频高清在线观看 | 午夜精品久久久久久毛片777| 国内精品美女久久久久久| 给我免费播放毛片高清在线观看| 老汉色∧v一级毛片| 欧美成狂野欧美在线观看| 全区人妻精品视频| 国产亚洲av嫩草精品影院| 国产视频一区二区在线看| 国模一区二区三区四区视频| 在线天堂最新版资源| 18美女黄网站色大片免费观看| 国语自产精品视频在线第100页| 亚洲欧美日韩卡通动漫| 亚洲色图av天堂| 老鸭窝网址在线观看| 国产乱人伦免费视频| 啪啪无遮挡十八禁网站| 麻豆成人av在线观看| 18美女黄网站色大片免费观看| 在线观看66精品国产| 18美女黄网站色大片免费观看| 三级男女做爰猛烈吃奶摸视频| 噜噜噜噜噜久久久久久91| 一卡2卡三卡四卡精品乱码亚洲| 中文亚洲av片在线观看爽| 90打野战视频偷拍视频| 特大巨黑吊av在线直播| 亚洲一区二区三区色噜噜| 国产三级黄色录像| 国产成人av教育| 老司机在亚洲福利影院| 久久久久久久久久黄片| 99精品久久久久人妻精品| av天堂中文字幕网| 嫁个100分男人电影在线观看| 一区二区三区激情视频| av女优亚洲男人天堂| 91久久精品电影网| 很黄的视频免费| 国产av在哪里看| 最新在线观看一区二区三区| 久久精品国产自在天天线| 法律面前人人平等表现在哪些方面| 在线观看一区二区三区| 黄色丝袜av网址大全| 国产精品,欧美在线| 日本一本二区三区精品| 国产精品三级大全| 国产亚洲精品一区二区www| 观看美女的网站| 成人高潮视频无遮挡免费网站| 久久久久久久精品吃奶| 欧美区成人在线视频| 国产亚洲av嫩草精品影院| 精品一区二区三区人妻视频| 国产成人aa在线观看| 大型黄色视频在线免费观看| 成人鲁丝片一二三区免费| 五月伊人婷婷丁香| 日韩欧美国产在线观看| 亚洲熟妇中文字幕五十中出| 日韩欧美在线乱码| 免费看美女性在线毛片视频| 亚洲最大成人中文| 日韩欧美精品免费久久 | 国产精品久久久人人做人人爽| 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 九九热线精品视视频播放| 韩国av一区二区三区四区| 18禁美女被吸乳视频| 国产精品一区二区三区四区免费观看 | 99久久99久久久精品蜜桃| 偷拍熟女少妇极品色| 成人高潮视频无遮挡免费网站| 日本 欧美在线| 国产精品亚洲一级av第二区| 久9热在线精品视频| 精品人妻偷拍中文字幕| 欧美中文日本在线观看视频| 婷婷六月久久综合丁香| 日韩大尺度精品在线看网址| av在线天堂中文字幕| 两性午夜刺激爽爽歪歪视频在线观看| 国产精品 国内视频| 亚洲国产高清在线一区二区三| 麻豆成人午夜福利视频| 午夜影院日韩av| 成人亚洲精品av一区二区| 哪里可以看免费的av片| 麻豆一二三区av精品| 国内毛片毛片毛片毛片毛片| 国语自产精品视频在线第100页| 男插女下体视频免费在线播放| 一区二区三区高清视频在线| 国产精品亚洲av一区麻豆| 国内揄拍国产精品人妻在线| 最近最新中文字幕大全电影3| 天堂√8在线中文| 中文字幕av在线有码专区| 观看美女的网站| 国产精品98久久久久久宅男小说| 麻豆国产av国片精品| 欧美+亚洲+日韩+国产| 又爽又黄无遮挡网站| 草草在线视频免费看| xxx96com| 欧美黄色片欧美黄色片| 最新美女视频免费是黄的| 岛国视频午夜一区免费看| 麻豆一二三区av精品| 哪里可以看免费的av片| 免费观看精品视频网站| 久久久久久国产a免费观看| 熟妇人妻久久中文字幕3abv| 亚洲精华国产精华精| 色哟哟哟哟哟哟| 美女 人体艺术 gogo| 18禁黄网站禁片免费观看直播| 精品人妻一区二区三区麻豆 | 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 国产精品女同一区二区软件 | 欧美绝顶高潮抽搐喷水| a在线观看视频网站| 在线观看av片永久免费下载| 特大巨黑吊av在线直播| 99精品欧美一区二区三区四区| 午夜视频国产福利| 免费人成视频x8x8入口观看| 亚洲一区二区三区色噜噜| 日本撒尿小便嘘嘘汇集6| 久久亚洲精品不卡| 国产激情欧美一区二区| 日日干狠狠操夜夜爽| 成人特级黄色片久久久久久久| 日韩有码中文字幕| 午夜福利欧美成人| 亚洲欧美日韩高清在线视频| 国内毛片毛片毛片毛片毛片| 99视频精品全部免费 在线| 亚洲成人久久爱视频| 五月伊人婷婷丁香| 人妻丰满熟妇av一区二区三区| 久久久精品欧美日韩精品| 网址你懂的国产日韩在线| 久9热在线精品视频| 欧美av亚洲av综合av国产av| 黑人欧美特级aaaaaa片| 偷拍熟女少妇极品色| 国产老妇女一区| 热99re8久久精品国产| 男人舔女人下体高潮全视频| 看免费av毛片| 最近最新中文字幕大全免费视频| 精品一区二区三区人妻视频| 久久草成人影院| 亚洲av免费在线观看| 天堂av国产一区二区熟女人妻| av视频在线观看入口| 亚洲国产高清在线一区二区三| 久久久久性生活片| 高清日韩中文字幕在线| 午夜日韩欧美国产| 免费看十八禁软件| 婷婷精品国产亚洲av在线| svipshipincom国产片| 久久久久久国产a免费观看| 色老头精品视频在线观看| 精品欧美国产一区二区三| 在线视频色国产色| 亚洲在线观看片| 天美传媒精品一区二区| 不卡一级毛片| 老司机午夜福利在线观看视频| 深夜精品福利| 日本熟妇午夜| 久久久国产成人免费| 欧美在线黄色| 国产欧美日韩精品一区二区| 12—13女人毛片做爰片一| 又爽又黄无遮挡网站| 岛国视频午夜一区免费看| 国产午夜精品论理片| 男插女下体视频免费在线播放| 日本五十路高清| 中文在线观看免费www的网站| 国产精品自产拍在线观看55亚洲| 午夜两性在线视频| 女人被狂操c到高潮| 久久中文看片网| 国产成+人综合+亚洲专区| 俺也久久电影网| 久久人人精品亚洲av| 最新美女视频免费是黄的| 国产精品女同一区二区软件 | 国产毛片a区久久久久| 国产真实乱freesex| 欧美bdsm另类| av专区在线播放| 亚洲aⅴ乱码一区二区在线播放| 老熟妇仑乱视频hdxx| 亚洲久久久久久中文字幕| 国产亚洲欧美在线一区二区| 激情在线观看视频在线高清| 国产麻豆成人av免费视频| 九九久久精品国产亚洲av麻豆| 国产伦精品一区二区三区四那| 欧美zozozo另类| 日韩人妻高清精品专区| 欧美中文综合在线视频| 欧美日韩中文字幕国产精品一区二区三区| 黄色日韩在线| 免费在线观看亚洲国产| 日韩欧美三级三区| 亚洲欧美日韩东京热| 国产淫片久久久久久久久 | 97超视频在线观看视频| 久久精品国产亚洲av涩爱 | 午夜视频国产福利| www.www免费av| 亚洲美女黄片视频| 中文字幕久久专区| 99国产精品一区二区蜜桃av| 成人国产一区最新在线观看| 高清在线国产一区| 免费大片18禁| 国产一区二区三区视频了| 色综合欧美亚洲国产小说| 19禁男女啪啪无遮挡网站| 欧美最新免费一区二区三区 | 蜜桃久久精品国产亚洲av| a级毛片a级免费在线| 国产精品久久久久久久电影 | 五月玫瑰六月丁香| 色老头精品视频在线观看| 美女 人体艺术 gogo| 丰满人妻一区二区三区视频av | 欧美性猛交黑人性爽| 国产精品亚洲一级av第二区| 少妇的丰满在线观看| 国产欧美日韩精品亚洲av| 亚洲av五月六月丁香网| 午夜激情欧美在线| 美女高潮喷水抽搐中文字幕| 少妇高潮的动态图| 日本撒尿小便嘘嘘汇集6| 国产精品一区二区三区四区久久| 亚洲av五月六月丁香网| 97超级碰碰碰精品色视频在线观看| 黄片小视频在线播放| 日韩精品中文字幕看吧| 久久久国产成人免费| 日韩欧美一区二区三区在线观看| 亚洲一区二区三区不卡视频| 别揉我奶头~嗯~啊~动态视频| 免费av毛片视频| 特级一级黄色大片| 国产真实乱freesex| 亚洲激情在线av| 久久精品国产亚洲av涩爱 | 色综合婷婷激情| 俺也久久电影网| 国产不卡一卡二| 亚洲国产日韩欧美精品在线观看 | 女人高潮潮喷娇喘18禁视频| 啦啦啦观看免费观看视频高清| 99久久无色码亚洲精品果冻| 午夜影院日韩av|