999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

爆炸荷載作用下覆土庫外部沖擊波傳播規(guī)律

2020-07-23 00:36:46董文學(xué)唐紅亮
工程爆破 2020年3期
關(guān)鍵詞:方向模型

榮 凱,楊 軍,董文學(xué),唐紅亮,崔 寧

(1.北京理工大學(xué)爆炸科學(xué)與技術(shù)國家重點(diǎn)實(shí)驗(yàn)室,北京 100081;2.中國兵器工業(yè)火炸藥工程與安全技術(shù)研究院,北京 100053)

覆土庫是庫體結(jié)構(gòu)表面覆蓋土的危險(xiǎn)品倉庫,一般用于存儲火炮、彈藥及軍用武器等。覆土庫由于其獨(dú)特的結(jié)構(gòu)特點(diǎn),在受內(nèi)部爆炸荷載作用下,對沖擊波的傳播具有一定的導(dǎo)向性。目前,許多學(xué)者對覆土庫及其他建筑結(jié)構(gòu)在內(nèi)部和外部爆炸沖擊波荷載作用下的動(dòng)態(tài)響應(yīng)以及沖擊波的傳播規(guī)律進(jìn)行了研究。

在受爆炸沖擊波荷載作用下結(jié)構(gòu)動(dòng)態(tài)響應(yīng)的研究方面,Yong Lu等[1]、洪武等[2]采用非線性動(dòng)力學(xué)有限元軟件ANSYS/LS-DYNA對地下混凝土結(jié)構(gòu)受外部爆炸荷載作用下的動(dòng)態(tài)響應(yīng)進(jìn)行模擬,比較了三維模型和二維模型模擬所得沖擊波峰值壓力和結(jié)構(gòu)加速度等參數(shù),得出了結(jié)構(gòu)的動(dòng)態(tài)響應(yīng)規(guī)律以及沖擊波的傳播過程??垫玫萚3]采用數(shù)值模擬的方法對地下拱形結(jié)構(gòu)受外部爆炸荷載作用下的動(dòng)態(tài)響應(yīng)進(jìn)行研究,將結(jié)構(gòu)的動(dòng)態(tài)響應(yīng)劃分為4個(gè)階段,得出了不同階段結(jié)構(gòu)動(dòng)態(tài)響應(yīng)能量占比的規(guī)律。Li Tian[4]對地下結(jié)構(gòu)受內(nèi)部爆炸時(shí),覆土對結(jié)構(gòu)的動(dòng)態(tài)響應(yīng)及沖擊波超壓傳播的影響進(jìn)行了研究。結(jié)果表明,覆土厚度對結(jié)構(gòu)的動(dòng)態(tài)響應(yīng)有一定的影響;結(jié)構(gòu)內(nèi)部爆炸時(shí)沖擊波的傳播十分復(fù)雜,沖擊波超壓隨與爆心比例距離的增大而減小。

在沖擊波傳播規(guī)律方面,仲倩等[5]、夏曼曼等[6]采用試驗(yàn)測試系統(tǒng)對炸藥在空中爆炸時(shí)近區(qū)爆炸沖擊波的傳播規(guī)律進(jìn)行了研究,提出了沖擊波峰值超壓與比例距離關(guān)系的修正經(jīng)驗(yàn)公式。劉偉等[7]采用試驗(yàn)和數(shù)值模擬的手段對近地面TNT爆炸時(shí)比例距離小于4 m/kg1/3區(qū)域內(nèi)沖擊波的傳播規(guī)律進(jìn)行研究,將模擬結(jié)果與試驗(yàn)測試結(jié)果進(jìn)行對比,驗(yàn)證了數(shù)值模擬算法及參數(shù)的準(zhǔn)確性。辛春亮等[8]采用不同算法對TNT空中爆炸進(jìn)行了沖擊波超壓計(jì)算,將數(shù)值模擬與經(jīng)驗(yàn)公式計(jì)算所得超壓值進(jìn)行比較,得出了2種不同算法的優(yōu)缺點(diǎn)。汪維等[9]采用數(shù)值模擬的方法對大當(dāng)量TNT在空氣中爆炸進(jìn)行數(shù)值模擬,將計(jì)算結(jié)果與經(jīng)驗(yàn)公式相比,得出了比例距離小于5 m/kg1/3時(shí)沖擊波超壓與沖擊波到達(dá)時(shí)間隨爆高的變化規(guī)律。

目前,學(xué)者們在結(jié)構(gòu)受內(nèi)部和外部爆炸時(shí)的動(dòng)態(tài)響應(yīng)方面研究較多,而在覆土庫結(jié)構(gòu)對沖擊波傳播影響方面研究較少。對于沖擊波傳播規(guī)律的研究,大部分學(xué)者側(cè)重于炸藥自由場爆炸情況下沖擊波的傳播,且研究范圍較小,并未對炸藥爆炸作用下覆土庫結(jié)構(gòu)外部較大范圍內(nèi)的沖擊波傳播規(guī)律進(jìn)行研究。筆者采用非線性動(dòng)力學(xué)有限元軟件ANSYS/LS- DYNA,結(jié)合將覆土庫結(jié)構(gòu)破壞與沖擊波傳播先后模擬的新手段,對覆土庫外部較大范圍內(nèi)爆炸沖擊波傳播規(guī)律進(jìn)行研究。分別建立覆土庫受炸藥內(nèi)部爆炸模型和沖擊波在空氣中傳播的數(shù)值模型,提取炸藥爆炸傳出覆土庫模型后不同方向單元內(nèi)的壓力時(shí)程,將其加載到?jīng)_擊波在空氣中傳播的數(shù)值模型中,計(jì)算并提取不同測點(diǎn)處沖擊波參數(shù),得到覆土庫結(jié)構(gòu)對沖擊波影響的傳播規(guī)律。

1 數(shù)值模型

1.1 模型的建立

本模型是用于研究覆土庫受1 t TNT內(nèi)部爆炸作用下比例距離小于15 m/kg1/3范圍內(nèi)空氣沖擊波的傳播規(guī)律。由于炸藥當(dāng)量較大,考慮空氣沖擊波傳播范圍較大,且涉及到覆土庫的破壞過程,如果將覆土庫的破壞與沖擊波的傳播在同一個(gè)有限元模型中進(jìn)行模擬,模型網(wǎng)格量巨大,難以對該問題進(jìn)行有效的數(shù)值模擬。因此提出了一種新的有限元模型建模方式,即將覆土庫的破壞和沖擊波的傳播分別建立有限元模型,并先后進(jìn)行模擬。覆土庫模型與不同方向沖擊波壓力測試線的位置關(guān)系如圖1所示,其中0°測線方向?yàn)楦餐翈烨皦Ψ较颉?/p>

圖1 數(shù)據(jù)測試線Fig.1 Data test line

覆土庫的破壞模擬完成后,提取覆土庫外0°、60°、90°、135°和180°方向單個(gè)空氣單元內(nèi)沖擊波壓力時(shí)程,將該壓力時(shí)程作為沖擊波傳播模型中的等效爆炸荷載,加載到5個(gè)不同測線方向的沖擊波傳播模型中,對各個(gè)方向沖擊波的傳播過程分別進(jìn)行計(jì)算。以0°測線方向沖擊波傳播計(jì)算為例,當(dāng)覆土庫的破壞模擬完成后,提取覆土庫外部0°方向沖擊波壓力時(shí)程(見圖2)。將0°方向沖擊波傳播模型起點(diǎn)處空氣單元(見圖5 BE處單元)設(shè)置壓力輸入屬性,并在K文件中將所提取的壓力時(shí)程數(shù)據(jù)點(diǎn)寫入關(guān)鍵字* DEFINE _CURVE中,從而完成沖擊波傳播模型等效爆炸荷載的施加。

圖2 0°方向等效爆炸荷載壓力時(shí)程Fig.2 Time history of equivalent explosion load pressure in the 0° direction

對于覆土庫受爆炸荷載破壞計(jì)算部分的數(shù)值模型,由于其具有對稱性,為減少計(jì)算時(shí)間,只建立1/2有限元模型(見圖3)。其中混凝土墻厚0.3 m。覆土庫的主體結(jié)構(gòu)是由波紋鋼搭建而成,并在波紋鋼拱頂部覆蓋0.50 m厚的覆土。

圖3 覆土庫模型Fig.3 ECM model

覆土庫受炸藥內(nèi)部爆炸荷載作用下破壞過程的數(shù)值模擬采用流固耦合算法進(jìn)行計(jì)算,需要建立流固耦合模型(見圖4)。其中,鋼筋與混凝土的耦合使用LS-DYNA中的關(guān)鍵字*CONSTRAINED_ LAGRANGE _IN_SOLID來實(shí)現(xiàn)。炸藥放置于庫體中心,其特征尺寸為0.85 m×0.85 m×0.425 m,空氣模型的尺寸略大于覆土庫結(jié)構(gòu)模型特征尺寸。

圖4 流固耦合模型Fig.4 Fluid-solid coupling model

對于空氣沖擊波傳播計(jì)算部分的數(shù)值模型,需要考慮沖擊波傳出覆土庫結(jié)構(gòu)后在空氣中的傳播規(guī)律進(jìn)行建模。文獻(xiàn)[10]提到,沖擊波在介質(zhì)中的衰減規(guī)律,隨距爆心比例距離的變化而變化。在比例距離小于5時(shí),沖擊波峰值超壓在介質(zhì)中呈二維衰減;在比例距離大于5 m/kg1/3時(shí),沖擊波峰值超壓在介質(zhì)中呈一維衰減。所以空氣沖擊波傳播計(jì)算模型在比例距離小于5 m/kg1/3時(shí)為扇形形狀;在比例距離大于5 m/kg1/3時(shí),呈條形形狀。

在炸藥爆炸近區(qū),測線方向覆土的有無會對沖擊波傳播計(jì)算模型近區(qū)的扇形形狀產(chǎn)生影響。土體對爆炸沖擊波具有很好的衰減作用,不同測線方向覆土的有無會對沖擊波傳出覆土庫的時(shí)間產(chǎn)生影響。由于覆土庫0°測線方向無覆土,此方向沖擊波最先傳出。在爆炸近區(qū),沖擊波不僅在豎直方向有壓力衰減,而且在水平方向也存在壓力衰減,故0°測線方向的沖擊波傳播計(jì)算模型在水平及豎直方向均采用扇形形狀。覆土庫結(jié)構(gòu)在其他測線方向均有覆土對沖擊波進(jìn)行衰減,沖擊波的傳出時(shí)間基本一致。在爆炸近區(qū),其他測線方向沖擊波傳播計(jì)算模型只在豎直方向采用扇形形狀。

以0°測線方向沖擊波傳播模型中參數(shù)的計(jì)算原理為例,說明各測線方向沖擊波傳播模型參數(shù)的計(jì)算方法(見圖5 ~ 圖6),其中O為爆心位置,A~L為沖擊波傳播模型各特征位置點(diǎn)編號。其中0°測線模型需考慮圖5、圖6中模型參數(shù)的計(jì)算,60°、90°、135°和180°測線方向只需要考慮圖5中模型參數(shù)的計(jì)算。

圖5 沖擊波傳播模型(側(cè)視)Fig.5 The shock wave propagation model (Side view)

圖6 沖擊波傳播模型(俯視)Fig.6 The shock wave propagation model (Top view)

如圖5所示,將多次數(shù)值模擬所得沖擊波峰值超壓值與經(jīng)驗(yàn)公式進(jìn)行比較,得出在比例距離小于4.4 m/kg1/3時(shí),沖擊波呈二維衰減;在比例距離大于4.4 m/kg1/3時(shí),沖擊波呈一維衰減。故AF離爆心O的水平比例距離為4.4 m/kg1/3(即44 m)。BE右側(cè)為覆土庫破壞計(jì)算部分。覆土庫破壞模型計(jì)算完成后,提取CD處一個(gè)單元內(nèi)的沖擊波壓力時(shí)程加載到?jīng)_擊波傳播計(jì)算模型的空氣域起點(diǎn)BE上。此計(jì)算方法需要保證提取沖擊波壓力時(shí)程的單元尺寸與加載沖擊波壓力時(shí)程的單元尺寸一致。在覆土庫破壞模型中,空氣的網(wǎng)格尺寸為10 cm,故沖擊波傳播計(jì)算模型BE處的網(wǎng)格尺寸也需要為10 cm。由于爆心高度為42.5 cm,取BE高50 cm,使BE處的高度正好為5個(gè)網(wǎng)格大小,以保證計(jì)算的精度。由BE高度和爆心的高度即可計(jì)算出扇形對應(yīng)圓心角∠AOF的度數(shù)。取HJ為一個(gè)網(wǎng)格尺寸10 cm,同理可以計(jì)算出垂直于測線方向扇形的角度∠GOK的度數(shù)。

根據(jù)圖5、圖6沖擊波傳播模型參數(shù)計(jì)算方法所得各方向測線扇形角度如表1所示。

1.2 網(wǎng)格劃分與邊界條件

在數(shù)值模擬中,網(wǎng)格尺寸的大小對計(jì)算結(jié)果影響很大,需要選取合適的網(wǎng)格尺寸進(jìn)行計(jì)算。覆土庫破壞模型中取網(wǎng)格邊長為10 cm。沖擊波傳播模型中取模型豎直方向網(wǎng)格尺寸為10 cm,網(wǎng)格尺寸隨模型高度的增大而增大;在模型水平方向網(wǎng)格尺寸為固定值10 cm。將數(shù)值模擬結(jié)果與經(jīng)驗(yàn)公式進(jìn)行比較,證明了該網(wǎng)格尺寸計(jì)算的準(zhǔn)確性。

數(shù)值模型的邊界條件是數(shù)值計(jì)算過程中的重要步驟。對于流固耦合模型,混凝土與波紋鋼底部采用全約束,取一半的對稱面設(shè)置對稱約束,其余表面不設(shè)置約束。空氣和炸藥模型在底部設(shè)置全約束,對稱面設(shè)置對稱約束,其余表面設(shè)置無反射邊界。對于沖擊波傳播模型,在模型底部設(shè)置全約束,壓力加載端不設(shè)置約束,模型遠(yuǎn)端設(shè)置無反射邊界,其余表面設(shè)置對稱約束。

1.3 材料參數(shù)

炸藥選取TNT,采用MAT_HIGH_EXPLOSIVE _BURN材料模型,并利用JWL狀態(tài)方程來實(shí)現(xiàn)爆炸荷載的施加[11],炸藥材料參數(shù)如表2[12]所示。

表2 炸藥材料參數(shù)

(1)

空氣采用MAT_NULL材料模型和LINEAR _POLYNOMIAL狀態(tài)方程來描述[11]:

p=C0+C1μ+C2μ2+C3μ3+(C4+C5μ+C6μ2)E

(2)

空氣材料參數(shù)見表3[13],其中ρ為材料密度,E0為初始內(nèi)能。

表3 空氣材料參數(shù)

土介質(zhì)采用MAT_SOIL_AND_FOAM_FAILURE材料模型,其屈服函數(shù)為[11]

(3)

式中:a0、a1、a2為動(dòng)力屈服常數(shù);φ為土介質(zhì)材料內(nèi)摩擦角;sij為應(yīng)力偏量;p為材料受到的壓應(yīng)力。

土介質(zhì)材料參數(shù)如表4所示[13]。

表4 土介質(zhì)材料參數(shù)

混凝土采用HJC材料模型,波紋鋼和鋼筋采用*MAT_PLASTIC_KINEMATIC材料模型,其材料參數(shù)如文獻(xiàn)[14]所示。

2 模型驗(yàn)證

土體對于爆炸沖擊波的衰減具有顯著的作用。由于覆土庫0°測線方向無覆土,炸藥在覆土庫內(nèi)部爆炸后從該方向傳出覆土庫的沖擊波峰值超壓與同等藥量炸藥自由場爆炸沖擊波峰值超壓基本一致。提取覆土庫0°測線方向不同比例距離處沖擊波峰值超壓,與Henrych、Mills和M.A.Sadovskyi經(jīng)驗(yàn)公式計(jì)算所得沖擊波峰值超壓進(jìn)行比較(見圖7)。

圖7 峰值超壓模擬結(jié)果與經(jīng)驗(yàn)公式計(jì)算結(jié)果Fig.7 The numerical and empirical formula results of peak overpressure

從圖7可知,0°測線方向模擬沖擊波峰值超壓與3種經(jīng)驗(yàn)公式在比例距離小于15 m/kg1/3時(shí),沖擊波峰值超壓大小基本一致,沖擊波峰值超壓均呈指數(shù)型衰減趨勢,驗(yàn)證了數(shù)值模型以及算法的正確性。

3 結(jié)果與討論

通過對覆土庫受內(nèi)部爆炸作用下結(jié)構(gòu)的破壞過程以及沖擊波傳播過程的模擬,對覆土庫結(jié)構(gòu)的破壞過程、沖擊波峰值超壓以及沖擊波到達(dá)時(shí)間等數(shù)據(jù)進(jìn)行分析,還原了覆土庫受內(nèi)部爆炸時(shí)的破壞過程,得出了沖擊波的傳播規(guī)律。

3.1 結(jié)構(gòu)破壞過程

覆土庫模型受炸藥內(nèi)部爆炸作用時(shí)的破壞過程如圖8所示。

圖8 覆土庫破壞過程Fig.8 The destroys process of ECM

從圖8可知,1.30 ms時(shí),爆炸沖擊波到達(dá)覆土庫門以及后墻的位置,庫門開始破壞;2.35 ms時(shí),爆炸沖擊波繼續(xù)對覆土庫結(jié)構(gòu)進(jìn)行沖擊,庫門及后墻破壞程度加劇,覆土庫頂部由中心開始破壞;7.60 ms時(shí),庫體側(cè)向及后部覆土庫破壞,覆土庫頂部完全破壞。

由于覆土庫0°測線方向無土體覆蓋,對爆炸沖擊波的抵抗力最弱,故最先破壞。由于庫體頂部覆土厚度相對于其他方向覆土厚度較小,對爆炸沖擊波的抵抗能力也較弱。當(dāng)爆炸沖擊波作用于土體時(shí),庫體頂部覆土最先開始破壞,且完全破壞的時(shí)間也早于其他方向土體。

3.2 沖擊波峰值超壓

覆土庫由于其結(jié)構(gòu)的特殊性,受內(nèi)部爆炸作用時(shí),對沖擊波的傳播具有一定的導(dǎo)向性,不同測線方向沖擊波的傳播規(guī)律不同。為研究沖擊波在不同測線方向上的傳播規(guī)律,提取不同測線方向沖擊波峰值超壓值(見圖9)進(jìn)行比較。

圖9 不同測線方向峰值超壓Fig.9 Peak overpressure in different directions

由圖9可知,不同測線方向,同一測點(diǎn)的沖擊波峰值超壓不同。0°測線方向沖擊波超壓最大,因?yàn)?°測線方向無覆土,對爆炸沖擊波的衰減能力最弱;135°和180°測線同一測點(diǎn)處峰值超壓比60°和90°測線方向峰值超壓略大,是因?yàn)?35°和180°測線方向覆土庫結(jié)構(gòu)的破壞略早于60°和90°測線方向。不同方向沖擊波峰值超壓均呈指數(shù)型衰減規(guī)律,隨著比例距離的增加,沖擊波峰值超壓趨于一致。測點(diǎn)距爆心比例距離越大,沖擊波峰值超壓受覆土庫結(jié)構(gòu)影響越小。

取不同比例距離處測點(diǎn)峰值超壓擬合值進(jìn)行比較(見表5)。當(dāng)覆土庫受內(nèi)部爆炸作用時(shí),含覆土測線方向的沖擊波峰值超壓明顯低于0°測線方向沖擊波峰值超壓;當(dāng)測點(diǎn)距爆心比例距離在1~15 m/kg1/3范圍時(shí),隨比例距離的增加,沖擊波峰值超壓衰減率逐漸降低。分別提取60°、90°和135°、180°相同比例距離處沖擊波峰值超壓求和取平均值,計(jì)算各測線方向?qū)_擊波峰值超壓的衰減率。在60°和90°測線方向,隨著測點(diǎn)比例距離的增加,沖擊波峰值超壓衰減率從87.63%降到26.39%;在135°和180°測線方向,隨著測點(diǎn)比例距離的增加,沖擊波峰值超壓衰減率從81.19%降到1.39%??梢钥闯?,60°和90°測線方向沖擊波峰值超壓最小,135°和180°測線方向沖擊波峰值超壓次之,0°測線方向沖擊波峰值超壓最大。在比例距離大于11 m/kg1/3時(shí),0°測線峰值超壓與135°和180°測線峰值超壓相差低于20%;在比例距離為15 m/kg1/3時(shí),0°測線峰值超壓與60°和90°測線峰值超壓相差25%。

表5 不同比例距離處峰值超壓擬合值

3.3 沖擊波到達(dá)時(shí)間

沖擊波到達(dá)不同比例距離處的時(shí)間是描述沖擊波傳播規(guī)律的重要參數(shù)。不同測線方向沖擊波到達(dá)時(shí)間的擬合曲線如圖10所示。

圖10 不同測線方向沖擊波到達(dá)時(shí)間Fig.10 Arrival time of shock wave in different directions

由圖10可以看出,不同測線方向沖擊波到達(dá)各測點(diǎn)的時(shí)間呈線性增長。0°測線方向沖擊波到達(dá)時(shí)間最早,135°和180°測線方向沖擊波到達(dá)時(shí)間次之,60°和90°測線方向沖擊波到達(dá)時(shí)間最晚。

4 結(jié)論

1)覆土庫受內(nèi)部炸藥爆炸作用時(shí),庫門最先破壞,庫頂次之,側(cè)向覆土最后破壞。135°和180°測線方向的覆土破壞時(shí)間要早于60°和90°方向的覆土破壞時(shí)間。

2)覆土對沖擊波峰值超壓的影響,隨測點(diǎn)比例距離的增大而減小。測點(diǎn)距爆心比例距離在1~15 m/kg1/3范圍內(nèi)時(shí),隨比例距離的增大,在60°和90°測線方向,沖擊波峰值超壓衰減率從87.63%降到26.39%;在135°和180°測線方向,沖擊波峰值超壓衰減率從81.19%降到1.39%。

3)沖擊波到達(dá)時(shí)間隨測點(diǎn)比例距離的增大,呈線性增大,0°測線方向沖擊波到達(dá)各測點(diǎn)時(shí)間最早。

猜你喜歡
方向模型
一半模型
2022年組稿方向
2022年組稿方向
2021年組稿方向
2021年組稿方向
2021年組稿方向
重要模型『一線三等角』
重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
主站蜘蛛池模板: 国产色爱av资源综合区| 亚洲一区精品视频在线| 亚洲欧美天堂网| 日本午夜影院| 亚洲午夜福利在线| 在线欧美国产| 四虎影视国产精品| 国产一国产一有一级毛片视频| 99色亚洲国产精品11p| 国产成人三级| 国产精品自在在线午夜| 免费无码网站| 午夜不卡视频| 黄色网站在线观看无码| 久久国产亚洲欧美日韩精品| 国产亚洲欧美日韩在线一区二区三区| 亚洲欧美成人影院| 免费视频在线2021入口| 国产精品 欧美激情 在线播放 | 欧美激情第一欧美在线| 久久网欧美| 毛片网站在线播放| 第九色区aⅴ天堂久久香| 亚洲成年人网| 精品久久高清| 四虎永久免费地址| 久久精品娱乐亚洲领先| 欧美精品另类| 在线观看国产一区二区三区99| 亚洲第一国产综合| 欧美一级高清免费a| 综合亚洲网| 国产女人综合久久精品视| 亚洲精品在线影院| 中国一级特黄视频| 欧美日韩一区二区三区四区在线观看| 欧美日韩中文字幕在线| 成人免费黄色小视频| 国产美女免费网站| 五月天久久综合| 自拍中文字幕| 日本午夜影院| 国产玖玖玖精品视频| 91区国产福利在线观看午夜| 9丨情侣偷在线精品国产| 中文字幕在线欧美| 国产毛片久久国产| 国产精品尤物铁牛tv| 成人精品在线观看| 免费看美女自慰的网站| 国产香蕉97碰碰视频VA碰碰看 | 久草热视频在线| 91免费精品国偷自产在线在线| 找国产毛片看| 九色在线观看视频| 久久影院一区二区h| 国产va欧美va在线观看| 免费毛片a| 好久久免费视频高清| 亚洲嫩模喷白浆| 久青草免费视频| 毛片久久网站小视频| 欧美日韩亚洲国产| 午夜人性色福利无码视频在线观看| 欧洲一区二区三区无码| 亚洲男人的天堂久久香蕉| 狠狠色噜噜狠狠狠狠奇米777| 中文字幕1区2区| 亚洲国产无码有码| 久久美女精品| 在线中文字幕日韩| 国产99久久亚洲综合精品西瓜tv| 91九色国产在线| 久草网视频在线| 久久午夜影院| 99热这里只有精品免费| P尤物久久99国产综合精品| 久久久久久久久久国产精品| 中文字幕永久在线看| 免费毛片网站在线观看| 国产香蕉在线视频| 国产福利一区视频|