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

降落傘充氣過程中傘衣外形及流場變化研究

2011-04-07 08:58:32張紅英劉衛華秦福德童明波
空氣動力學學報 2011年3期

張紅英,劉衛華,秦福德,童明波

(南京航空航天大學航空宇航學院,江蘇 南京 210016)

0 引言

在降落傘拉直、充氣及穩定下降等一系列動作中,充氣過程是影響其安全穩定性最為關鍵的環節,也是物理過程中最復雜的一個過程,因而成為了當前降落傘研究領域的熱點及難點問題之一。在對降落傘開傘過程的分析中,傳統上主要依靠試驗研究,但由于充氣時間一般很短,有時不到半秒就完成了充氣過程,導致在這樣的一個動態充氣過程中試驗狀態的模擬和各項參數的記錄非常困難。因此降落傘的風洞試驗研究一般偏重在開傘載荷與氣動力系數的研究上[1-4]。而事實上,由于降落傘是個柔性織物透氣體,它在開傘過程中經歷了急劇的結構大變形,因而傘衣內外的流場十分復雜,只有深刻地了解降落傘周圍的流場,才能深刻理解降落傘工作時的工作機理,也才能更好地從理論上解決我們關心的問題,解釋傘衣充氣時的諸多氣動力現象。這就必須依靠充氣過程的動態數值模擬來獲得傘衣形狀變化和流場變化之間的關系,從而提高降落傘的理論分析水平。

隨著計算機技術的發展和數值模擬技術精確度的提高,同時降落傘設計、研制所要求的周期縮短及日益高漲的試驗費用,近年來,已有不少研究人員發表了論文,探討如何采用計算流體力學(Computational Fluid Dynamics,CFD)方法來模擬降落傘開傘過程[5-10]。這些研究論文已見報道并給了我們很多啟示。但同時我們也發現,這些計算模型都不能模擬降落傘從任意初始形狀迅速展開直到完全充滿,即不能對降落傘初始充氣階段的開傘過程進行模擬。在文獻中,一些研究者使用了充分膨脹的形狀作為最初形狀[7,10];另一些研究者使用了軸對稱的假定來進行近似計算[5-6,8];有的甚至僅考慮穩定狀態的情況[9,11]。為此,本文以平面圓形傘為原型,對軸對稱降落傘整個充氣階段的流場進行了數值模擬和分析,以詳細了解降落傘整個開傘過程中傘衣外形變化及傘衣內外流場的變化情況。

1 充氣過程數學模型

1.1 結構動力學模型

本文以平面圓形傘為原型建立模型,對降落傘模型作如下假設:

1)忽略質點的重力影響,徑向加強帶的充氣過程是軸對稱的;

2)傘繩、加強帶和傘衣是彈性體,符合虎克定律。

降落傘傘衣由一定數量的相同傘衣幅彼此縫合連接而成,相鄰傘衣幅之間由加強帶相連。在對稱充氣情況下,選取兩對稱徑向加強帶來表示傘衣充氣過程的形狀變化,將研究對象(一根加強帶和兩邊的兩個半幅傘衣)離散為一系列用阻尼彈簧連接的質點,分布在加強帶上。質點與質點之間的受力示意圖如圖1所示,其中FP為傘衣質點受到的氣動壓力,α為氣動力與Y軸的夾角,Ti、Ti-1分別為i與i+1質點、i與i-1質點之間的彈性力和阻尼力之和,β為傘衣質點之間的張力和X軸的夾角,Fz為傘衣織物(收口繩、傘頂孔加強帶)緯向張力在傘坐標X軸方向上的分力(當兩加強帶之間的距離大于傘衣幅的寬度時才存在)。根據牛頓運動定律,傘衣各質點的結構動力學方程為:

圖1 質點間受力示意圖Fig.1 Force diagram of the free body of canopy

1.2 計算流體力學模型

降落傘傘衣的厚度很小,遠小于傘衣的幾何尺度。因此從流體力學的角度來考慮,可以把傘衣看成是薄膜材料,其厚度對流場本質不會產生影響,在對降落傘進行流體力學數值模擬時可以忽略傘衣厚度的影響。采用準定常假設,把充氣過程中的降落傘視為剛體,通過對不同充氣時刻的不同傘衣外形進行不可壓Navier-Stokes(N-S)方程求解來模擬降落傘的流場特性。選用二維雷諾平均N-S方程作為控制方程,采用有限體積結構網格數值格式進行流場模擬。

二維守恒型雷諾平均Navier-Stokes方程可寫為:式中:W為守恒通量,Re為雷諾數,E、F為對流通量,Ev,Fv為粘性通量。

在這基礎上對結構模型和流場模型進行耦合計算:先根據初始數據,在每一時間步長開始時,降落傘到達一個新位置,將傘面附近的流場網格點移動至傘面上;通過修改這些點相應的動量方程,更新源項,利用CFD程序計算得出交界處質點的壓差,將計算結果傳給結構動力學方程,進行表面力計算和進行傘面變形運動計算,得到下一形狀,依據下一形狀,流場計算程序進入下一個時間步長,開始新一輪的計算。

織物的透氣性是指在織物的兩側存在空氣壓力差時,空氣從織物的孔隙透過的性能。本文引用了多孔介質(Porous jump)邊界條件,并以602錦絲66綢為例進行計算模型驗證。Porous jump邊界條件是對多孔介質模型(porous media model)的一維簡化,用于模擬流場中已知“壓力損失-透氣速度”關系的“膜”。該邊界條件的“壓力損失-透氣速度”關系由以下公式給定:其中μ是流體粘性系數,ρ為流體密度,α是介質的透氣系數,C2是壓力躍變系數,v是垂直于介質表面的速度大小,△m是介質厚度。α、C2和△m是Porous jump邊界條件中需要設置的參數。

從公式(4)中可見,將壓力損失△p表示為透氣速度v的二次函數。在實際計算中,流體粘性系數μ、流體密度ρ和介質厚度△m為已知量,還有α和C2兩個參數需要確定。如果已知一組透氣速度和壓力損失的實驗數據,則可以通過對實驗數據進行二次多項式插值的方法確定α和C2兩個參數。下面以602錦絲66綢為例說明α和C2的確定過程。

表1為602錦絲66綢透氣速度與壓力損失的一組實驗數據[12],采用過原點的二次多項式對這組數據進行擬合(如圖2所示)可得:

比較式(5)和式(4),由對應常數項相等可得:

標準大氣條件下空氣密度ρ為1.225kg/m3,粘性系數μ為1.726×106kg/ms,設介質厚度△m為0.0005m,則由以上兩式計算可得:

表1 不同透氣速度下壓力損失的計算值與實驗值[12]對比Table1 The comparison of calculated and experimental values [12] of air speed and pressure loss

圖2 二次多項式曲線擬合結果Fig.2 Quadratic polynomial curve fitting results

為了驗證采用Porous jump邊界條件模擬降落傘傘衣透氣性的可靠性,制作了一個如圖3所示的計算模型。該模型的計算域為二維長方形“通道”,側邊指定為周期邊界,Porous jump邊界設置在“通道”中間位置,相當于在通道中間加了一層“透氣膜”,在速度入口指定空氣流速,壓力出口設置為標準大氣壓。對該模型進行計算直至收斂后,可以得到入口和出口之間的壓力差,此壓力差即為Porous jump邊界模擬的“透氣膜”在當前流速(即透氣速度)下的壓力損失。根據之前對602錦絲66綢透氣參數的計算,設定Porous jump邊界的參數為:△m=0.0005m,C2=435395.9,α =1.1270×10-10。將速度入口的空氣流速分別設定為表1中的速度,即可計算得到相應的壓力損失,結果如表1所示。從表1可看出計算值與實驗值的結果比較接近。本文的傘衣材料為411平紋綢,也可使用該方法模擬其透氣性。

圖3 透氣性驗證計算模型的邊界條件和網格示意圖Fig.3 Boundary condition and grids schematic for porosity verification calculation

2 試驗裝置及測試系統

充氣時間一般很短 ,有時不到半秒就完成了充氣過程,在如此短的時間內,傘衣發生了劇烈的變形,傘衣周圍流場變化十分復雜,在這樣的一個動態充氣過程中試驗狀態的模擬和各項參數的記錄非常困難。因此只對降落傘充滿時的繞流流場進行了定量測量,將該狀態的測量結果和計算結果進行對比分析,以檢測計算結果的可靠性。

試驗是在南京航空航天大學的非定常回流低速風洞進行的,降落傘在風洞中的照片如圖4所示。該風洞是國內首座非定常風洞,通過水平并列旁路加上非定常流動控制機構實現試驗段的非定常流場,在作為定常風洞使用時具有低湍流度,低噪聲的特點。實驗段開口為1.5m×1m的矩形,全長1.7m,紊流度小于0.05%,最大風速是40 m/s,最低穩定風速為0.5m/s。

圖4 降落傘在風洞中的照片Fig.4 Experimental apparatus

降落傘的內外流場十分復雜,氣流有較大的偏角,因此采用七孔探針來測量其周圍的流場,七孔探頭的結構如圖5所示。七孔探針可以得到流場中的三維速度及壓力信息,對氣流偏角為78°的大偏角流動,其測試精度為1%[9]。為提高測量效率,利用多根七孔探頭制成耙,可以同時測量空間多個點的氣動參數。氣孔探針在試驗前需要進行校準,由于氣孔探針的校準系數較多,實際試驗中七孔探針將流場中局部點的七個孔壓力值由多通道微壓變送測試儀轉換成模擬電壓量,經高精度采集模塊將電壓量轉換成數字量。經過計算機采集后按校準試驗得到的系數運算得到速度向量,局部總壓、靜壓值等測量結果。整個測試系統的配置如圖6所示。

圖5 七孔探頭的結構Fig.5 Seven-hole probe

圖6 測試系統配置圖Fig.6 Testing system

本試驗采用有傘頂空的平面圓形傘,其物理參數為:傘衣副數為8塊,傘頂孔直徑為100 mm,傘衣直徑為600mm。傘繩長465 mm,材料為3.2-120錦絲繩,其彈性模數為4000 N。傘衣名義面積為0.32m2,材料為411平紋綢,彈性模數為30000N/m。

3 結果分析

3.1 充氣過程傘衣外形變化及流場變化

計算模型與試驗用傘的物理模型一致,采用以第1節所介紹的數學模型對該模型進行了整個充氣階段的數值模擬,計算的初始狀態與實驗狀態一致,氣流來流速度為20m/s。通過計算,得到了傘衣外形變化和傘衣周圍的流場變化。

初始充氣階段是氣流從傘衣底部沖到傘頂的階段,圖7為初始充氣階段的傘衣外形變化及流場變化圖。圖7(a)為氣流剛進入傘衣口時的狀態,從圖中可看出,外擴狀“喇叭口”導致大的分離區,在“喇叭口”之后產生較大面積的局部低壓,該區域傘衣內外壓差較大,達到1.6左右,促進傘衣口迅速膨脹擴張。圖7(b)所示的狀態中,傘衣口外形擴張為“口杯”狀,流動分離區緊貼傘衣外表面,相比狀態7(a)明顯變窄,分離區傘衣的內外壓差除緊鄰傘衣底邊的小區域外,總體上小于狀態7(a)的情況,傘衣內存在渦流,但流速極低,基本不影響傘衣內滯止壓力的大小。圖7(c)所示的狀態中,氣流在傘衣內形成一對方向相反的漩渦。傘衣呈先擴后收的形狀,傘衣外基本為附壁流動。擴張部分傘衣向外排擠氣流,使氣流加速,傘衣外表面壓力降低,傘衣內外壓差系數在1.5左右。收縮部分氣流減速,傘衣外表面壓力逐漸恢復至大氣壓,傘衣內外壓差隨之逐漸減小。圖7(d)所示的狀態為傘頂孔完全沖開時的狀態,此時氣流從傘衣底邊流入,再從傘頂孔流出,沒有發生流動分離,可見傘衣的結構透氣性對流場會產生顯著影響。從壓力等值線圖來看,傘衣內的滯止壓力系數在頂端開口和傘衣透氣的條件下仍然達到最高0.85左右,傘衣外表面的壓力系數則因為平滑的傘衣外形和傘衣透氣性的影響降低到僅有-0.15左右。說明傘衣的透氣性對傘衣內滯止壓力的影響較小,對傘衣外流場結構的影響較大,從而對傘衣外表面的壓力影響較大。

圖7 初始充氣階段的傘衣外形變化及流場變化圖Fig.7 Canopy shape and flow-field vs.time in the initial inflation phase

圖8為主充氣階段的傘衣外形變化及流場變化圖。圖8(a)所示的狀態為主傘主充氣階段的初期,可以看出傘頂孔張開至更大,氣流在傘衣內的滯止有所減弱,且由于傘衣織物透氣的作用,氣流滯止過程從傘衣口至傘衣頂持續進行,傘衣內部壓力系數從傘衣口附近的0.76增加至傘衣頂附近的0.9以上。盡管傘衣外形在氣流方向前大后小,由于傘衣織物透氣的作用,傘衣外表面仍沒有發生明顯的氣流分離。傘衣尾部由于有傘頂孔“射流”的存在,也沒有發生明顯的流動分離。傘衣外表面壓力除傘衣口附近局部較低外,其余部分壓力接近大氣壓,傘衣尾部壓力也很快恢復至大氣壓。

圖8(b)所示的狀態中傘頂孔面積張至最大,空氣首先在傘衣頂部聚集,使傘衣頂部膨脹。沒有明顯的流動分離發生,傘衣內部流動滯止效應進一步減弱,壓力系數為0.74。與圖8(a)中的狀態不同的是,由于傘頂孔面積增大,結構透氣量增大,傘衣織物透氣性對氣流在傘衣內滯止的影響相對減弱,傘衣內壓力系數從上一個狀態的逐漸增大轉變為基本不變。另外,由于傘衣直徑在尾部有增大,形成局部低壓。

圖8 主充氣階段的傘衣外形變化及流場變化圖(前:壓力系數等值線圖,后:流線圖)Fig.8 Canopy shape and flow-field vs.time in the Inflation Phase

如圖8(c)所示,傘衣直徑繼續擴大,呈直筒狀。從圖中的流線圖可以看出,在傘衣外表面附近流線發生扭曲,但沒有形成漩渦,這是傘衣織物透氣的影響(如果沒有傘衣織物透氣性,將會有漩渦形成)。在傘衣尾部仍然沒有看到明顯的流動分離發生,流線以很大的曲率彎曲繞過傘衣尾部,這一方面是受傘衣透氣性的影響,另一方面是傘頂孔射流對周圍氣流有吸入作用。從壓力等值線圖來看,傘衣內部壓力系數從傘衣口的0.9左右向傘頂方向逐漸增加至0.96以上,逐漸滯止現象再次出現圖8(b)所示狀態中(該現象不明顯),這是因為傘頂孔面積張至最大保持不變以后,隨著傘衣的繼續擴張,傘衣展開面積增大,傘衣側面織物透氣的面積繼續增大,對氣流滯止的影響也相對減小。

如圖8(d)所示,相對8(c)中的狀態,傘衣形狀呈長寬比更小一些的直筒狀(從傘衣口到傘衣頂稍有擴張)。由于傘衣尾部面積的增大,氣流無法再順滑地繞過而發生分離漩渦。由于主流和傘頂孔射流分別在漩渦區的兩側都沿順流方向對旋渦區內空氣形成強剪切,因此在一個漩渦區內形成兩個相反旋轉方向的漩渦。從壓力等值線來看,漩渦區使得整個傘衣尾部的壓力有所降低,該區域壓力需要遠離傘衣尾部更遠的距離才完全恢復至大氣壓力。

在圖8(e)中,傘衣直徑繼續擴大,傘衣尾部已經呈現圓弧形。漩渦區隨著傘衣直徑的擴張而增大。主流對漩渦區氣流的剪切作用相對傘頂孔射流的剪切作用增大(傘衣直徑擴張分離區擴大,則主流對漩渦區的剪切面積增大,而傘頂孔面積不變,其射流的剪切面積也基本不變),因此漩渦區兩個漩渦中外側的一個漩渦變大。從整個漩渦區的尺度來看,已經與傘衣張開的尺度相當。從壓力等值線來看,傘衣尾部壓力隨著旋渦區的增大繼續降低,且低壓范圍也繼續擴大。傘衣內部壓力仍主要受氣流滯止作用影響。由于傘衣開口面積增大,傘衣透氣量相對進入傘衣的氣流量的比例有所下降,因而傘衣內滯止壓力有所上升,達到接近1.0。

在圖8(f)中,傘衣完全張滿。主流對漩渦區氣流的剪切面積更大,主要受主流剪切作用影響的漩渦(單側的兩個漩渦中靠外的一個)一面向外擴大尺度,一面向內擠壓傘頂孔射流和另一個漩渦。在傘衣尾流軸線上形成了較長區域的回流,傘頂孔射流的繞行現象更明顯,射流和內側漩渦被擠壓至離傘頂更近的區域。從壓力等值線來看,傘衣后出現了更大面積的低壓區,這是漩渦區擴大的結果。

3.2 計算與試驗的對比

圖9為試驗中三維測量結果在半幅中心面上的流場顯示,對比圖9(f)中計算所得到的流場可以看出,傘衣前、側流場結構均相似,從傘衣尾部看,計算場和試驗場的截面流線都出現了匯聚,從傘頂及繞傘衣底邊流出的流線都會包圍中間的紊流旋渦向尾部流去。但由于測量場的測量間距偏大,流場拓撲結構不明顯。再對壓力場進行比較分析,在圖9中顯示的是壓力,測量壓力單位為mmH2O,試驗的來流速度是20m/s。圖中顯示傘衣內最大壓力約為240Pa(24.37),換算壓力系數約為0.98。而傘衣外部漩渦中心區最低壓力約為-95Pa(-9.58 mmH2O),換算壓力系數約為-0.4。對比圖8(f)中的壓力系數,二者誤差不大。因此所采用的數學模型是可靠的,能通過其計算來定性分析繞傘衣流動的流場特性。

圖9 風洞流場測量結果Fig.9 The experimental flow-field around canopy

4 結論

本文以平面圓形傘為原型,對軸對稱降落傘整個充氣階段的流場進行了數值模擬和分析,以詳細了解降落傘整個開傘過程中傘衣外形變化及傘衣內外流場的變化情況。通過計算得到如下結論:

(1)初始充氣階段傘衣外形變化為:傘衣底邊張開后,氣流進入傘衣,在氣流作用下,傘衣折疊不分從下到上依次張開,直至氣流沖到傘衣頂部,傘頂孔打開。整個初部分從下到上依次張開,直至氣流沖到傘衣頂部,傘頂孔打開。整個初始充氣階段,傘衣展開部分外形基本保持較光滑的直筒形狀,而非喇叭形。在主充氣階段:空氣首先在傘衣頂部聚集,使傘衣頂部膨脹,然后膨脹部分向傘衣底邊擴展,直到傘衣完全張滿。

(2)對于頂部有傘頂孔的平面圓形傘,當傘頂孔被氣流沖開后,傘衣的結構透氣性對流場會產生顯著影響。但從內外壓力系數的變化來看,透氣性對傘衣內滯止壓力的影響較小,對傘衣外流場結構的影響較大,從而對傘衣外表面的壓力影響較大。

(3)對于此類有傘頂空的平面圓形傘,當傘衣充氣張開后,傘衣尾部出現氣流分離,由于主流和傘頂孔射流分別在漩渦區的兩側都沿順流方向對旋渦區內空氣形成強剪切,因此在一個漩渦區內形成兩個相反旋轉方向的漩渦。且隨著傘衣直徑擴張,分離區擴大,主流對漩渦區的剪切面積增大,因此漩渦區兩個漩渦中外側的一個漩渦增大,內測漩渦被擠壓至離傘頂更近的區域。

(4)對降落傘充滿時的繞流流場進行定量測量,對比計算結果發現,從流場拓撲結構來看,計算結果和實驗結果的拓撲結構相似。再從壓力場的比較可看出,計算與試驗得到的壓力系數誤差不大。因此所采用的數學模型是可靠的,能通過其計算來定性分析繞傘衣流動的流場特性。

[1]BARBER J,JOHARI H.Experimental investigation of personnel parachute designs using scale model wind tunnel testing[A].16th AIAA Aerodynamic Decelerator Systems Technology Conference[C].Boston:MA.AIAA 2001-2074,2001.

[2]SHANNON M P,GEORGIA A G.Experimentalanalysis ofthe pressure distribution on a 35-foot personnel parachute[C].A-merican Institute of Aeronautics&Ast-ronautics or Published with Permission of Author(s)and/or Author(s)'Sponsoring Organization[R].AIAA 2001-2008 ,2001:114-121.

[3]DESABRAIS K L,JOHARIH.The flow in the near wake of an inflating parachute canopy[R].AIAA 2001-2009,2001:122-130.

[4]DESABRAIS K J.Velocity field measurement s in the near wake of a parachute canopy[D].Doctor Thesis of Worcester Polytechnic Institute.2002.

[5]STEIN K R,BENNEY R J.Parachute inflation:a problem in aero-elasticity[R].US Army Tech.Report,NATICK/TR-94/015,1995.

[6]STEIN K R,BENNEY R,TEZDUYAR T,et al.3-D computation of parachute fluid-structure interactions:performance and control[R].AIAA 99-1714,1999.

[7]STEIN K R,BENNEY R,KALRO V,TEZDUYAR T E,LEONARD J,ACCORSI M.Parachute fluid-structure interactions:3-D computation[J].Compuer Methods in Applied Mechanics and Engineering,2000,190:373-386.

[8]HAUG E,LASRY D,Kermel P.Dynamic simulation of industrial membranes including their interaction with surrounding media[A].Third International Symposium of the SBF 230 Evolution of Natural Structures[C].University of Stuttgart,1994.

[9]SAHU J,COOPER G,BENNEY R.3-D parachute descent analysis using coupled parachute descent codes[A].13th AIAA Aerodynamic Decelerator Systems Technical Conference[C].AIAA paper 99-1580,1999.

[10]余莉,史獻林,明曉.降落傘充氣過程數值模擬[J].航空學報,2007,28(1):53-57.

[11]蔣崇文,曹義華,蘇文翰.軸對稱降落傘小迎角穩定下降時流場特性[J].中國空間科學技術,2007,28(2):59-65.

[12]馬衍富.降落傘織物的透氣性[J].產業用紡織品,1987,5(1):26-30.

主站蜘蛛池模板: 国产一二三区在线| 成人午夜网址| 色综合天天视频在线观看| 全部毛片免费看| 国产大片黄在线观看| 亚洲国产精品日韩av专区| 国产午夜小视频| 久久久亚洲色| 中国国语毛片免费观看视频| 欧美精品在线免费| 亚洲人精品亚洲人成在线| 色成人亚洲| 国产swag在线观看| 2020极品精品国产| 57pao国产成视频免费播放| 国产一国产一有一级毛片视频| 91精品福利自产拍在线观看| 国产精品永久在线| 亚洲成人一区在线| 国产三级视频网站| 亚洲激情99| 免费国产高清视频| 91在线高清视频| 色妞www精品视频一级下载| 精品乱码久久久久久久| 大陆精大陆国产国语精品1024| 亚洲成人高清无码| 国产99精品视频| 精品国产自在在线在线观看| 久久精品丝袜高跟鞋| 久久香蕉欧美精品| 啪啪国产视频| 国产精品亚洲专区一区| 一级黄色片网| 99尹人香蕉国产免费天天拍| 日本精品一在线观看视频| 精品久久香蕉国产线看观看gif| 午夜a级毛片| 国产探花在线视频| 影音先锋亚洲无码| 少妇精品网站| 亚洲黄网视频| 丰满的熟女一区二区三区l| 色播五月婷婷| 中文字幕在线日韩91| 91精品啪在线观看国产91| 九九久久精品免费观看| 亚洲伊人久久精品影院| 啪啪永久免费av| 美女无遮挡拍拍拍免费视频| 久久黄色视频影| 老司国产精品视频91| 5555国产在线观看| 午夜视频在线观看免费网站 | 综合色区亚洲熟妇在线| 日韩国产黄色网站| 亚洲无码一区在线观看| 亚洲视频三级| 国产成人在线无码免费视频| 高清亚洲欧美在线看| 五月六月伊人狠狠丁香网| 久久亚洲精少妇毛片午夜无码 | www.91在线播放| 又爽又大又光又色的午夜视频| 久久精品国产一区二区小说| 亚洲欧美不卡视频| 国产女人在线视频| 无码专区在线观看| 欧美亚洲一二三区| 综合人妻久久一区二区精品| 日本影院一区| 久久久久青草大香线综合精品| av性天堂网| 在线另类稀缺国产呦| 波多野结衣久久高清免费| 在线欧美国产| 精品无码人妻一区二区| 中文字幕在线播放不卡| 欧美亚洲国产精品第一页| 国产在线欧美| 91欧美在线| 亚洲天堂区|