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

考慮速度趨膚效應與接觸電阻影響的樞軌界面電流分布特性

2022-10-14 06:32:50王增基陳立學尤彭昊蘭昕宇葛一凡
電工技術學報 2022年19期
關鍵詞:界面模型

王增基 陳立學 尤彭昊 蘭昕宇 葛一凡

考慮速度趨膚效應與接觸電阻影響的樞軌界面電流分布特性

王增基1,2陳立學1,2尤彭昊1,2蘭昕宇1,2葛一凡3

(1. 強電磁工程與新技術國家重點實驗室(華中科技大學電氣與電子工程學院) 武漢 430074 2. 脈沖功率技術教育部重點實驗室(華中科技大學) 武漢 430074 3. 武漢工程大學 武漢 430205)

該文采用有限元法,仿真研究了軌道炮的電流分布特性。根據有限元模型中影響因素的不同,建立了三種模型,即速度趨膚效應模型、接觸電阻模型以及兩者都包含的綜合模型。首先在運動電樞上建立坐標系,在麥克斯韋方程組中增加速度項,實現考慮速度趨膚效應的運動電磁場計算;其次接觸電阻模型考慮了電樞與軌道的粗糙表面特征以及樞軌間接觸壓力;最后將兩種模型結合建立了綜合模型。該文通過對比分析各個模型中的電流在樞軌界面上的分布特性,表明綜合模型能更準確地揭示電流分布特性?;诖?,以電感梯度為例,該文進一步分析了不同模型中電感梯度變化規律,研究表明,速度趨膚效應與接觸電阻都是影響電磁軌道炮性能的關鍵因素。

軌道炮 速度趨膚效應 接觸電阻 電流分布

0 引言

軌道炮中的電流分布是進行溫度分布、界面損傷、結構優化等研究的基礎。必須實現電流分布的準確定量計算,才能實現進一步的可靠研究。電流分布的鄰近效應、趨膚效應都可以通過對麥克斯韋方程組的瞬態求解來實現。所以,對于軌道炮電流分布的研究有兩個重點:如何考慮接觸電阻以及速度趨膚效應對電流分布的影響。

目前,軌道炮的電流分布研究主要分為以下三種模型:僅考慮接觸電阻的電磁場模型[1-6](接觸電阻模型);僅考慮速度趨膚效應的電磁場模型[6-14](速度趨膚效應模型);理想接觸且不考慮速度趨膚效應[14-23](簡單模型)。

目前最為常用的是簡單模型,采用商業軟件進行電磁場瞬態求解即可獲得軌道炮內部的電流分布[14-23]。其次是接觸電阻模型,常用的建模方法有電樞分段建模法[1-3]、Holm接觸電阻公式[4-5]及CMY(cooper-mikic-yoranovich)接觸電阻公式[6]。接觸電阻公式提出了額外的力學仿真要求,而電樞分段建模法則主要取決于研究人員的經驗。對速度趨膚效應的研究也有多種方法。許多學者多采用磁擴散方程進行二維模型的速度趨膚效應的研究[6-7,14]。金龍文采用LS-DYNA實現了速度趨膚效應的三維模型計算,意味著采用有限元法是可行的[8-9]。李白、李湘平等采用有限體積法實現了三維模型的速度趨膚效應的計算[11-12]。林慶華則通過自編程序對速度趨膚效應進行研究[13]。

一般認為,在軌道炮的起始階段時速度較低,可以忽略速度趨膚效應。而當速度較高時,速度趨膚效應占據主導,可以忽略接觸電阻。事實上兩者計算得到的電流分布截然不同,接觸電阻會將電流限制在高接觸壓力區域附近,而速度趨膚效應將導致電流集中在電樞尾翼端部。兩種影響電流分布的機制是如何協同作用的,目前尚未開展相關研究。

本文首先給出了考慮速度的麥克斯韋方程組,用于三維仿真模型中實現速度趨膚效應的計算;其次在接觸電阻模型中,認為接觸斑點總是塑性變形的,采用CMY接觸電阻計算方法進行電磁場-結構場模型的雙向耦合計算;最后,結合以上的控制方程和邊界條件,建立綜合模型,并對比分析計算得到的電流分布。

1 速度趨膚效應模型

模型采用如圖1所示電流波形。上升沿階段從0開始線性增加至0.5ms,保持200kA至3ms,隨后下降,至8ms為0。模型中采用20mm口徑軌道炮,其幾何模型如圖2所示,其中電樞與軌道為完全接觸。

圖1 電流波形

圖2 軌道炮的幾何模型

模型中軌道高度為30mm(方向),軌道厚度為20mm(方向),電樞高度為(方向)19.4mm,樞軌接觸面長度為(方向)20mm。電樞位于軌道的中間部位。為實現速度趨膚效應而不引入幾何模型的運動,本文將坐標系建立在電樞上,則軌道在發射過程中不斷向后方移動。將速度對電磁場的影響添加到麥克斯韋方程組中,有

式中,為電樞的發射速度;為磁場強度;為電位移矢量;為磁通密度;為電場強度;為電流密度。相比于一般的麥克斯韋方程組,式(1)中添加了動生電動勢產生的電流。由于電樞不動,電樞中控制方程的速度項恒為0,而軌道中的速度項與實際發射過程中電樞的速度大小相等方向相反。

本文根據圖2所示幾何模型中的對稱關系,建立了軌道炮的1/4模型。基于COMSOL軟件進行仿真計算可得圖2中研究區域對應的樞軌界面電流密度分布,如圖3所示。速度趨膚效應模型參數見表1。

表1 速度趨膚效應模型參數

Tab.1 The parameters of velocity skin effect model

圖3 速度趨膚效應模型中樞軌界面的電流密度分布

當速度較大時,電流幾乎全部集中于電樞尾翼的邊緣角落,此處的電流密度遠超附近區域。本文進一步提取出樞軌界面邊緣線上的電流分布,如圖4所示??梢钥闯?,在各個時刻尾翼末端區域都是電流集中的區域。

圖4 速度趨膚效應模型中樞軌界面邊緣處的電流分布

在0.3ms、0.7ms、5ms時刻的速度分別為20m/s、150m/s、1 370m/s。在0.3ms,由于速度較低,瞬態趨膚效應占據主導地位,速度趨膚效應對電流分布的作用不明顯。此時電流集中于電樞尾翼頭部與末端。而在0.7ms時刻,電流波形進入平頂階段,速度趨膚效應的影響開始突顯。尾翼末端的電流密度幅值遠超頭部區域。而在5ms時刻,速度的大幅提高沒有導致電流分布形狀在對數坐標系下的顯著變化。但在靠近電樞尾翼末端的位置(=0),電流的數值梯度大幅增加。

由此可以看出,在電樞速度為150m/s時,速度趨膚效應已經非常顯著。在電流的下降階段由于電樞進入高速階段,電流變化產生的感生電動勢已經無法觀察到。

2 接觸電阻模型

本文建立的接觸電阻模型采用尾翼向兩側張開的電樞,其單邊過盈量為0.9mm。首先對裝配后的狀態進行穩態計算,得到接觸壓力分布。在發射的過程中,在電磁力的作用下,樞軌界面的接觸壓力將發生變化。

在樞軌界面上采用CMY接觸電阻計算方法[24],對于任意一點的接觸壓力,其對應的接觸電導率c為

表2 接觸電阻模型參數

Tab.2 The parameters of contact resistance model

在穩態計算后進行接觸電阻模型的瞬態計算。瞬態計算時采用電磁場與結構場的雙向耦合模型。模型以不含速度項的麥克斯韋方程組作為控制方程。其多物理場耦合方法如圖5所示。每一次迭代中,電磁場向結構場輸出洛侖茲力作為載荷,結構場輸出根據接觸壓力計算到的接觸電阻。

圖5 多物理場耦合方法

考慮接觸電阻作為樞軌界面的邊界條件,可以得到樞軌界面各時刻的電流分布如圖6所示,接觸壓力分布如圖7所示。樞軌界面邊緣處的電流分布如圖8所示。

0時刻接觸壓力即為樞軌界面的預接觸壓力。此時高接觸壓力區域主要集中在邊緣區域。以后在電磁力的作用下,接觸壓力逐漸變得均勻。0.5~3ms是電流的平頂階段,電流分布呈現出向中央擴散的趨勢。而電樞邊緣處的電流密度也略微下降,呈現出均勻化的趨勢。相較于接觸壓力的均勻化過程,電流密度的均勻化過程較為滯后。此時,由于電流幅值沒有變化,所以接觸壓力的分布幾乎沒有變化。在5ms時刻,電流處于下降階段,受到感生渦流的影響,電流更趨于向中央匯集。邊緣的電流密度在電流幅值降低、感生電動勢的雙重作用下,下降幅度明顯。在此階段,樞軌界面中央區域的接觸壓力降低,接觸壓力向邊緣處集中。根據接觸電阻模型的計算結果,在整個發射過程中,電樞尾翼末端的接觸壓力總為0,意味著沒有與軌道接觸。

圖6 接觸電阻模型中樞軌界面的電流分布

圖8 接觸電阻模型中樞軌界面邊緣處的電流分布

3 綜合模型

對比上述兩個模型的電流分布云圖與電樞邊緣電流密度的曲線圖可以發現,兩者具有極為顯著的區別。因此,為建立能夠更加準確描述實際軌道炮發射過程的模型,本文結合速度趨膚效應模型與接觸電阻模型,建立綜合模型進行電流分布研究。

綜合模型首先進行樞軌裝配過程的穩態計算,獲得樞軌界面的初始接觸電阻。然后進行電磁-結構雙向耦合的瞬態計算。綜合模型的結構場計算方法與接觸電阻模型相同。在電磁場計算上采用與接觸電阻模型相同的邊界條件以及材料屬性,而使用速度趨膚效應模型的控制方程,即式(1)。

綜合模型中樞軌界面的電流分布如圖9所示,接觸壓力分布如圖10所示。

圖9 綜合模型中樞軌界面的電流分布

圖10 綜合模型中樞軌界面的接觸壓力分布

綜合模型合并考慮了接觸電阻模型和速度趨膚效應模型的特點。在0.3ms時刻,電流集中于樞軌界面邊緣靠近頭部的某點,此處為高接觸壓力區域的上邊緣。這與接觸電阻模型中得到的電流分布是一致的。而隨著速度的增加,速度趨膚效應逐漸顯現,使得電流集中區域不斷向電樞尾翼末端移動。

盡管速度趨膚效應模型中電流密度的峰值位于電樞尾翼末端,但在綜合模型中距離尾翼末端的區域內電流密度卻為0,如圖9所示。原因在于,此區域接觸壓力為0,如圖10所示。速度趨膚效應能夠使得電流不斷向后方移動,但產生的電磁力仍然不足以推動尾翼末端與軌道接觸。電樞尾翼末端與軌道失接觸,電流無法通過。樞軌界面邊緣處不同時刻的電流分布如圖11所示,直觀地展示了綜合模型中電流分布特性。

圖11 綜合模型中樞軌界面邊緣處的電流分布

簡單模型、速度趨膚效應模型、接觸電阻模型以及綜合模型中得到的電感梯度在發射過程中的變化如圖12所示。由于7ms以后電流趨于0造成了較大的數值波動,且在實際的發射過程中為了避免轉捩,往往在下降階段開始不久即出膛。所以不考慮7~8ms之間的數據。

四種模型在起始時刻的電感梯度區別不大。接觸電阻模型與簡單模型相差1%左右。這說明軌道結構與電樞結構才是決定電感梯度的根本因素。

在發射的主要過程中,電感梯度最大的為接觸電阻模型,其次為簡單模型、綜合模型,最小的為速度趨膚效應模型。在速度趨膚效應模型中,軌道電流被強烈地限制在電樞尾翼末端。軌道電流距離電樞電流較遠,提供的電磁推力也較小,如圖13a所示。且速度越大,電樞尾翼末端區域電流占比越大,因此電感梯度呈現下降的趨勢。

圖12 四種模型中的電感梯度

圖13 兩種模型中的電流分布與樞軌相互作用示意圖

對于接觸電阻模型而言,接觸壓力主要集中于樞軌界面的中央區域,導致此處的接觸電阻小,所以電流主要經過中央區域。這使得電樞與軌道的作用距離大為降低,所以電感梯度較大,如圖13b所示。而在磁擴散的作用下,樞軌間的作用距離會進一步減小,所以電感梯度呈現增加的趨勢。

而綜合模型的初期階段在接觸電阻的作用下,電流被限制在樞軌界面邊緣的中間區域,從而電感梯度在初期與接觸電阻模型近似。但隨著速度的提高,速度趨膚效應逐漸顯著,使得電感梯度也呈現出下降趨勢,但接觸電阻的作用仍不可忽略,所以下降速度比速度趨膚效應模型小一些。

最后,簡單模型的磁擴散過程在平頂、下降階段中發展較為緩慢,所以整體的電感梯度幾乎是均勻的。0~1ms時刻的電感梯度上升是由于電流線性上升,磁擴散較為迅速造成的。

基于上述分析可知,綜合模型能夠準確地反映兩種機制對電流分布的影響,也是目前最貼近實際的電流分布計算方法。

4 結論

電磁軌道發射中,速度趨膚效應和接觸電阻對電流分布的影響機制不同,導致電流分布特性也差異顯著。本文提出了包含接觸電阻和速度趨膚效應的綜合模型。通過與其他模型的對比分析,證明了兩種效應在綜合模型中均有顯著作用,可真實反映電磁特性。綜合模型克服了傳統模型僅適用于某一個階段的局限性,為軌道炮的全程發射仿真提供了精確的電流分布計算方法。

[1] 劉勇, 國偉, 張濤, 等. 雙層電樞結構設計及其電磁力和電流密度分布[J]. 彈箭與制導學報, 2020, 40(1): 160-164.

Liu Yong, Guo Wei, Zhang Tao, et al. Geometry design of double-layer armature and its electromagnetic force and current density distribution[J]. Journal of Projectiles, Rockets, Missiles and Guidance, 2020, 40(1): 160-164.

[2] 申澤軍, 左鵬, 袁建生. 電磁軌道炮電樞與軌道接觸面大小對電流密度的影響分析[J]. 高電壓技術, 2014, 40(4): 1084-1090.

Shen Zejun, Zuo Peng, Yuan Jiansheng. Influence of contact area size between armature and rails in railguns on current density[J]. High Voltage Engineering, 2014, 40 (4): 1084-1090.

[3] Li Chengxian, Chen Lixue, Xia Shengguo, et al. Simulations on saddle armature with concave arc surface in small caliber railgun[J]. IEEE Transactions on Plasma Science, 2019, 47(5): 2347-2353.

[4] 高翔, 張暉輝, 劉暢, 等. 考慮接觸電阻的電磁發射裝置溫度場分析[J]. 兵器材料科學與工程, 2020, 43(3): 59-64.

Gao Xiang, Zhang Huihui, Liu Chang, et al. Temperature field analysis of electromagnetic launcher under consideration of contact resistance[J]. Ordnance Materials Science and Engineering, 2020, 43(3): 59-64.

[5] 王雪軍, 顧金良, 羅紅娥, 等. 電磁軌道炮內彈道接觸電阻分析[J]. 南京理工大學學報, 2019, 43(3): 306-311.

Wang Xuejun, Gu Jinliang, Luo Honge, et al. Analysis of internal ballistic contact resistance of electromagnetic railgun[J]. Journal of Nanjing University of Science and Technology, 2019, 43(3): 306-311.

[6] 殷強, 張合, 李豪杰, 等. 考慮電樞與導軌實際接觸狀態的電磁軌道炮膛內磁場分析[J]. 兵工學報, 2019, 40(3): 464-472.

Yin Qiang, Zhang He, Li Haojie, et al. Analysis of in-bore magnetic field in electromagnetic railgun considering the realistic armature-rail contact status[J]. Acta Armamentarii, 2019, 40(3): 464-472.

[7] 殷強, 張合, 李豪杰. 動態條件下電磁軌道炮膛內磁場和電場分析[J]. 兵工學報, 2017, 38(6): 1059-1066.

Yin Qiang, Zhang He, Li Haojie. Analysis of in-bore magnetic and electric fields in electromagnetic railgun under dynamic condition[J]. Acta Armamentarii, 2017, 38(6): 1059-1066.

[8] Jin Longwen, Li Jun, Lei Bin. Approximate field scaling of railgun launcher under the condition of matching projectile dynamic parameters[J]. IEEE Transactions on Plasma Science, 2015, 43(9): 3286-3292.

[9] Jin Longwen, Lei Bin, Zhang Qian, et al. Electromechanical performance of rails with different cross-sectional shapes in railgun[J]. IEEE Transactions on Plasma Science, 2015, 43(5): 1220-1224.

[10] Dwight L, Sikhanda S. Eddy current effects in the laminated containment structure of railguns[J]. IEEE Transactions on Magnetics, 2007, 43(1): 150-156.

[11] 李湘平, 魯軍勇, 譚賽, 等. 基于Fluent二次開發的電磁軌道發射運動磁場仿真[J]. 中國電機工程學報, 2020, 40(19): 6364-6371.

Li Xiangping, Lu Junyong, Tan Sai, et al. Simulation on moving magnetic field of electromagnetic rail launch based on fluent secondary development[J]. Proceedings of the CSEE, 2020, 40(19): 6364-6371.

[12] 李白, 魯軍勇, 譚賽, 等. 有限體積法在滑動電接觸問題中的應用[J]. 海軍工程大學學報, 2019, 31(6): 23-28.

Li Bai, Lu Junyong, Tan Sai, et al. Application of finite volume method in analyzing sliding electrical contact problem[J]. Journal of Naval University of Engineering, 2019, 31 (6): 23-28.

[13] 林慶華, 栗保明. 基于瞬態多物理場求解器的電磁軌道炮發射過程建模與仿真[J]. 兵工學報, 2020, 41(9): 1697-1707.

Lin Qinghua, Li Baoming. Modeling and simulation of electromagnetic railgun launching process based on a transient multi-physical field solver[J]. Acta Armamentarii, 2020, 41(9): 1697-1707.

[14] 湯鈴鈴, 李豪杰. 電磁軌道炮膛內磁場環境仿真分析[J]. 計算機仿真, 2014, 31(11): 1-5, 46.

Tang Lingling, Li Haojie. Simulation analysis of railgun in-bore high magnetic field[J]. Computer Simulation, 2014, 31(11): 1-5, 46.

[15] 任先進, 張春. 靜止條件下電磁軌道炮膛內磁場環境仿真分析[J]. 火控雷達技術, 2018, 47(2): 82-84, 90.

Ren Xianjin, Zhang Chun. Simulation analysis of in-bore magnetic field environment of electromagnetic rail-gun at static condition[J]. Fire Control Radar Technology, 2018, 47(2): 82-84, 90.

[16] 呂慶敖, 王維剛, 邢彥昌, 等. 電磁軌道炮鐵磁材料對銅帶內電流分布的影響[J]. 強激光與粒子束, 2015, 27(10): 268-271.

Lü Qingao, Wang Weigang, Xing Yanchang, et al. Effect of ferromagnetism material on current distribution in copper strips for electromagnetic railguns[J]. High Power Laser and Particle Beams, 2015, 27(10): 268-271.

[17] 邢彥昌, 呂慶敖, 陳建偉, 等. 軌道炮不同激勵電流下的發射特性對比分析[J]. 火力與指揮控制, 2019, 44(11): 58-60, 66.

Xing Yanchang, Lü Qingao, Chen Jianwei, et al. Comparative analysis of launching characteristic for railgun with different excitation current[J]. Fire Control & Command Control, 2019, 44(11): 58-60, 66.

[18] Zuo Peng, Geng Yiqing, Li Jun, et al. An approach for eddy-current calculation in railguns based on the finite-element method[J]. IEEE Transactions on Plasma Science, 2015, 43(5): 1592-1596.

[19] Mohammad S, Asghar K. Study of the current distribution, magnetic field, and inductance gradient of rectangular and circular railguns[J]. IEEE Transactions on Plasma Science, 2013, 41(5): 1376-1381.

[20] 李鶴, 李治源, 雷彬, 等. 電磁軌道炮不同截面軌道的特性分析[J]. 火力與指揮控制, 2014, 39(4): 45-48, 53.

Li He, Li Zhiyuan, Lei Bin, et al. Analysis on rail performance of EM railgun with different cross sections[J]. Fire Control & Command Control, 2014, 39 (4): 45-48, 53

[21] 杜佩佩, 魯軍勇, 馮軍紅, 等. 電磁軌道發射器電磁結構耦合動態發射過程數值模擬[J]. 電工技術學報, 2020, 35(18): 3802-3810.

Du Peipei, Lu Junyong, Feng Junhong, et al. Numerical simulation of dynamic launching process of electromagnetic rail launcher with electromagnetic and structural coupling[J]. Transactions of China Electrotechnical Society, 2020, 35(18): 3802-3810.

[22] 古剛, 吳立周, 耿昊, 等. 基于電磁-流場耦合的軌道冷卻仿真分析[J]. 電工技術學報, 2020, 35(17): 3601-3608.

Gu Gang, Wu Lizhou, Geng Hao, et al. Simulation and analysis of rail cooling based on electromagnetic and fluid field coupling[J]. Transactions of China Electrotechnical Society, 2020, 35(17): 3601-3608.

[23] 張威, 李海元, 栗保明. 基于有限元仿真的不同截面軌道炮特性分析[J]. 計算機仿真, 2021, 38(8): 18-22, 171.

Zhang Wei, Li Haiyuan, Li Baoming. Characteristics analysis of railgun with different cross section based on finite element simulation[J]. Computer Simulation, 2021, 38 (8): 18-22, 171.

[24] Yovanovich M. Four decades of research on thermal contact, gap, and joint resistance in microelectronics[J]. IEEE Transactions on Components and Packaging Technologies, 2005, 28(2): 182-206.

Current Distribution Characteristics of Armature-Rail Interface under Velocity Skin Effect and Contact Resistance

Wang Zengji1,2Chen Lixue1,2You Penghao1,2Lan Xinyu1,2Ge Yifan3

(1. State Key Laboratory of Advanced Electromagnetic Engineering and Technology School of Electrical and Electronic Engineering Huazhong University of Science and Technology Wuhan 430074 China 2. Key Laboratory of Pulsed Power Technology Ministry of Education Huazhong University of Science and Technology Wuhan 430074 China 3.Wuhan Institute of Technology Wuhan 430205 China)

In this paper, the current distribution of railgun is calculated by finite element method. In this paper, according to the different factors in the model, three models are established: the velocity skin effect model, the contact resistance model and the composite model. Firstly, by setting up the coordinate system on the armature and adding the velocity term to the Maxwell equations, the velocity skin effect caused by the armature motion is calculated. Secondly, the rough surface characteristics of armature and rail and the contact pressure between armature and rail are considered in the contact resistance model. Finally, a composite model is established by combining above models. By comparing and analyzing the current distribution characteristics of each model on the armature-rail interface, this paper shows that the composite model is the most accurate. It is shown that the velocity skin effect and contact resistance are the key factors affecting the performance of electromagnetic railgun.

Railgun, velocity skin effect, contact resistance, current distribution

10.19595/j.cnki.1000-6753.tces.211767

TM89

國家自然科學基金(51877096)和重點實驗室基金(6142217200507)資助項目。

2021-11-03

2022-01-22

王增基 男,1997年生,博士研究生,研究方向為電磁發射。E-mail:1379130580@qq.com

陳立學 男,1984年生,教授,博士生導師,研究方向為電磁發射、開關電器等。E-mail:chenlixue@hust.edu.cn(通信作者)

(編輯 赫蕾)

猜你喜歡
界面模型
一半模型
重要模型『一線三等角』
國企黨委前置研究的“四個界面”
當代陜西(2020年13期)2020-08-24 08:22:02
重尾非線性自回歸模型自加權M-估計的漸近分布
基于FANUC PICTURE的虛擬軸坐標顯示界面開發方法研究
空間界面
金秋(2017年4期)2017-06-07 08:22:16
電子顯微打開材料界面世界之門
人機交互界面發展趨勢研究
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
主站蜘蛛池模板: 四虎影院国产| 日韩在线成年视频人网站观看| 国产精品成人一区二区| 欧美中日韩在线| 东京热一区二区三区无码视频| 亚洲熟妇AV日韩熟妇在线| 日韩精品一区二区三区大桥未久| 国产流白浆视频| 亚洲国产成人自拍| 在线日韩一区二区| 久草青青在线视频| 国内自拍久第一页| 成人午夜福利视频| 亚洲婷婷丁香| 亚洲精品无码高潮喷水A| 亚洲视频无码| 中国毛片网| 国产高清毛片| 高清欧美性猛交XXXX黑人猛交| 欧美一区二区三区国产精品| 久久精品无码专区免费| 91青青草视频在线观看的| 国产黄网永久免费| 国产小视频a在线观看| 亚洲最新在线| 国产丰满大乳无码免费播放| 国产国模一区二区三区四区| 毛片手机在线看| 区国产精品搜索视频| 午夜在线不卡| 亚洲天堂久久新| 亚洲熟妇AV日韩熟妇在线| 国产无码性爱一区二区三区| 日韩高清一区 | 久久精品视频亚洲| 国产中文在线亚洲精品官网| 国产欧美自拍视频| 亚洲一区二区日韩欧美gif| 国产欧美日本在线观看| 亚洲精品麻豆| 久久熟女AV| 99久久精品免费看国产电影| 成人国产一区二区三区| 国产精品v欧美| 99精品一区二区免费视频| 青青青视频免费一区二区| 免费毛片视频| 久久久国产精品免费视频| 美女国产在线| 狠狠色丁香婷婷综合| 国产区人妖精品人妖精品视频| 久久精品国产91久久综合麻豆自制| 日韩天堂在线观看| 中文字幕 91| 国产欧美在线视频免费| 亚洲精品动漫| 国产手机在线小视频免费观看| 全部免费特黄特色大片视频| 欧美日韩专区| 黄色免费在线网址| 国产成人亚洲综合A∨在线播放| 亚洲人成网站在线播放2019| 国产精品男人的天堂| 少妇露出福利视频| 精品伊人久久大香线蕉网站| 国产视频久久久久| 日韩av高清无码一区二区三区| 手机在线免费不卡一区二| 亚洲欧洲AV一区二区三区| 国产成人在线无码免费视频| 久热re国产手机在线观看| 国产又爽又黄无遮挡免费观看| 伊人天堂网| 亚洲三级成人| 欧美在线免费| 亚洲美女高潮久久久久久久| 免费大黄网站在线观看| 国产欧美日韩在线一区| 2021国产精品自产拍在线观看| 日韩天堂网| 99资源在线| 欧美精品亚洲日韩a|