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

艦船集防區對外通道結構變形及氣密性仿真分析

2016-11-03 05:25:28楊光趙尚輝李天勻朱翔
中國艦船研究 2016年5期
關鍵詞:有限元變形結構

楊光,趙尚輝,李天勻,朱翔

1華中科技大學船舶與海洋工程學院,湖北武漢430074

2船舶與海洋水動力湖北省重點實驗室,湖北武漢430074

3中國船舶及海洋工程設計研究院,上海200011

艦船集防區對外通道結構變形及氣密性仿真分析

楊光1,2,趙尚輝3,李天勻1,2,朱翔1,2

1華中科技大學船舶與海洋工程學院,湖北武漢430074

2船舶與海洋水動力湖北省重點實驗室,湖北武漢430074

3中國船舶及海洋工程設計研究院,上海200011

艦船集防區對外通道中,氣密門是氣密周界上主要的泄漏點之一。目前氣密性評估一般采用實驗檢測的方法。結合有限元仿真方法,在集體防護系統建造前對氣密門的氣密性進行初步評估。采用設計波法結合國內外軍民船規范,利用有限元軟件將長期預報得到的設計波波浪動水壓力值、重力、慣性力、靜水壓力等施加到全艦有限元模型上,全面反映真實受載下艦船在設計波下的響應,進而預報全艦變形情況。選擇變形較大的密封結構,分析變形條件下密封結構的密封性能。研究結果表明:借鑒全艦結構強度有限元直接計算方法,能有效分析對外通道的局部變形特征;氣密門處的局部變形將導致橡膠密封圈接觸壓力分布不均。結合局部變形的特點,可對密封結構的優化設計提供參考。

集體防護系統;變形;設計波;密封結構;氣密性

網絡出版地址:http://www.cnki.net/kcms/detail/42.1755.TJ.20160921.1338.020.html期刊網址:www.ship-research.com

引用格式:楊光,趙尚輝,李天勻,等.艦船集防區對外通道結構變形及氣密性仿真分析[J].中國艦船研究,2016,11(5):71-77.

YANG Guang,ZHAO Shanghui,LI Tianyun,et al.Numerical simulation of the deformation and airtightness of the external passagewaystructureofnavalships'collectiveprotectionsystem[J].ChineseJournalof ShipResearch,2016,11(5):71-77.

0 引言

艦船集體防護(簡稱“集防”)系統是艦船核生化防護能力最核心的組成部分。如果不具備有效的防護手段,在面臨核生化武器突襲的情況下,艦船會在短時間內即喪失生存力和戰斗力[1]。目前,國外艦船裝備的集體防護系統主要分為兩大類:一類是德國的全船式集體防護系統,即DSK系統;另一類是美國的分區式集體防護系統,即CPS[2]。最早用于對抗核生化威脅的集體防護系統DSK由德國海軍于上世紀60年代中期研發設計成功,該防護系統將空氣凈化和空氣調節整合成一個單元系統[3]。在“漢堡”級驅逐艦上首次裝備的DSK防護系統是結合在全船加熱、通風、空調系統上的濾毒增壓通風系統,具備數小時的核生化防護能力[4]。美國海軍于上世紀70年代中后期至80年代初期,在水面艦船上嘗試建立了有限集體防護系統[5]。此后,美國海軍第一艘具備全艦全時集體防護能力的主戰艦船DDG-51型導彈驅逐艦誕生,表明艦船集防系統得到了迅速發展[6]。

我國軍事裝備集體防護的研究工作始于上世紀50年代,直到上世紀60~70年代集體防護裝備的國產化批量生產才得以實現,而這些裝備當時也僅配備在主戰坦克和步兵戰車上。針對水面艦船集體防護系統的研究則相對起步較晚,直至21世紀初才在水面艦船成功配備了集體防護系統[7]。無論是小規模的個人防護,還是大規模的集體防護,其防護原理是一樣的,即與外部隔斷,依靠核生化防護過濾系統與外部通風,因此,集體防護區中對外通道部分的氣密性至關重要[8-9]。

艦船在水面航行時,其結構的變形是不可避免的。當船體變形的量級與船用密封結構的指標精度相當或超出時,可能會使得密封結構無法正常使用。因此,抑制、降低船體變形對密封結構的密封性能具有重要意義。本文的研究對象為有密性要求的艦船集防區對外通道,將主要分析其變形特征以及變形條件下的密封性能。

在對外通道中,氣密門是氣密周界上的主要泄漏點之一。目前氣密性評估一般采用實驗檢測的手段,可采用的常規氣密性檢測方法有涂抹煤油檢驗法、粉檢驗法、沖氣檢驗法等。本文將結合有限元仿真方法對氣密門密封性能進行初步評估,以便在其建造前對氣密性能進行初步評估,為密封結構的設計提供參考。

由于艙段模型邊界條件對變形的影響較大[10],且對外通道的密封結構大量存在于上層建筑中,故本文將根據軍民船規范[11-12],結合設計波法對全船有限元模型進行變形計算,預報全船的總體變形和局部變形特征,彌補艙段模型的不足,最終提取出變形后的密封結構以分析其密封性能。

1 全船有限元直接計算

1.1設計波法簡介

設計波法是通過三維勢流理論求解頻率響應函數,然后通過譜分析方法對艦船剖面的波浪載荷進行長期預報,根據艦船實際運行地區的海況和艦船設計要求,利用長期預報中得到的一定超越概率水平的波浪載荷(通常取概率水平為10-8),將其認為是該艦船使用壽期內所遇到的波浪載荷極值,然后選取一規則波,使其產生的載荷與長期預報所得的相同。也即設計一個波浪,用來模擬艦船所遇到的最危險海況,然后將艦船放置在該波浪下,觀測其響應,以此來觀察艦船是否安全,是否能滿足規范要求。

上述頻率響應函數,是指艦船在單位波幅的規則波中的響應(如6個自由度的運動、剪力、彎矩等),其相應的波浪頻率曲線如圖1所示。從圖1中可以很清楚地得到艦船在何處載荷達到極值,因此便可以得知處于該位置時相應的艦船狀態,這樣便可以知道一個規則波中的響應極值。而通過對多個不同的規則波使用相同方法便最終能夠確定出艦船在規則波中運動響應的極值[13-14]。

圖1 頻率響應函數圖Fig.1Frequency response function curve

設計波法的具體步驟可以歸納如下:

1)利用三維勢流理論對艦船進行水動力參數的分析,并計算艦船在各浪向角、各遭遇頻率的規則波中的響應值,包括艦船6個自由度的位移值,船體所受的彎矩、剪力和扭矩等,還有波浪作用于艦船濕表面上的動壓力分布等。

2)利用艦船在規則波中的響應值,通過海浪譜資料,對波浪誘導船體載荷進行長期預報,從而確定出艦船在給定時間內航行于實際海況中的波浪載荷的統計特性。在確定了艦船的運行海域和使用壽期內可能遭遇到的最大波浪載荷的概率水平(通常取10-8)后,便可得到波浪載荷的極值,也就是船舶在遭遇波浪數為N的使用壽期內,平均可能出現一次的最危險的海況所遭遇到的最大波浪載荷。找出與波浪載荷極值相對應的規則波和波高,此即為設計波。

1.2船體載荷及邊界條件

靜水載荷:靜水載荷將作為全艦直接強度分析的一部分施加到船體梁上進行分析,靜水載荷包括靜水壓力和重力兩部分:在MSC/PATRAN里根據工況通過定義一個與垂向坐標z線性相關的域(field)函數,即可施加船體濕表面各點處的靜水壓力;重力則是根據空船重量、壓載水以及設備和甲板載荷的具體分配位置,改變相應位置板的密度來模擬該處重力的變化。

波浪動壓力計算與施加[15]:利用設計波法可得到船體濕表面動水壓力,然而計算船體濕表面的動水壓力值的網格和船體有限元外殼模型網格不匹配,計算水動壓力的濕表面網格一般較粗,而有限元網格尺寸則相對小得多,因此需編制PATRAN內部的PCL程序實現加載,即PCL加載程序。該程序的編制流程如下:

1)讀取各濕表面網格中心點的三坐標值和水動力壓力值。

2)找出與有限元模型表面任意點最近的濕表面網格中心點,則作用在該船體外表面有限元板網格的法向動壓力為pk。

3)賦值,實現加載。

利用該PCL程序定義域(field),可實現板殼網格表面壓力的加載。圖2所示為船體濕表面在極限波中的水動力壓力分布示意圖。

圖2 船體濕表面的水動力壓力分布(軸側圖)Fig.2Hydrodynamic pressure distribution on wet surface of ship hull(isometric view)

本文中船體邊界條件采用NASTRAN中的慣性釋放功能。

1.3結果分析

此艦船主船體結構均采用AH36鋼,屈服應力為355 MPa;上層建筑采用普通鋼,屈服應力為235 MPa。鋼材的楊氏模量為2.06×105MPa,泊松比為0.3。有限元模型采用板(Shell)單元模擬構件中的板殼結構,如甲板、外板、縱桁等以及各種構件之間的連接肘板。采用梁(Beam)單元模擬各種構件上承受側向壓力的縱骨、扶強材等,并按照實際情況考慮梁的截面和偏心。主要構件的面板和加強筋用桿單元模擬,如非水密的縱桁、肋板、平面艙壁析材上的面板和加強筋等。全艦有限元模型網格尺寸最大為700 mm,節點數為239 894,單元數為442 457。

對要求的不同工況分別進行分析,應力和變形的計算結果可在PATRAN后處理模塊實現顯示。根據有限元計算結果,分析集防區對外通道的應力、變形特征,包括通道結構的變形和密封構件(門、窗等)處開口的變形特征。整船的應力云圖和變形示意圖如圖3所示。

圖3 全艦應力變形云圖Fig.3The contours of stress and deformation of a whole ship

2 氣密門密封性能評估

氣密門的性能評估存在非線性材料和接觸的問題,對網格的要求也與整船模型不同,故不宜在整船模型中建立氣密門模型進行計算。由于氣密門對門框附近的結構剛度影響不大,故在整船有限元模型中,氣密門處用開口處理并建立門框,變形相對實際情況偏危險。

在密封條與門框的接觸分析中,由于密封條在接觸過程中變形很大,可用變形接觸體描述,而門框是由鋼板焊接而成的,在接觸過程中變形很小,可以忽略,可用剛體描述。在密封條與門框的接觸中,密封條變形的力與位移通過與之相接觸的門框的運動產生[16-18]。變形后的門框從全艦結構有限元計算結果中提取。

2.1門框的變形

由于要考慮變形對密封性能的影響,將從全艦結構有限元計算結果中提取出變形后的門框,相當于利用變形結果預制變形門框,進行氣密性評估。選取全艦應力變形最大的艙壁為研究對象,目標艦中該艙壁為主甲板上的上層建筑后端壁,由于該處應力集中、變形較大,同時為了避免接觸分析中相鄰面元間出現大的轉角,需要在全艦模型中對該處的艙壁網格進行細化。細化后的有限元模型如圖4所示。

圖4 門框處有限元模型(局部)Fig.4Finite element model of the doorframe(partial)

在PATRAN中利用變形的結果建立位移場,通過位移場將所有節點移動到變形后的位置,然后使用“Create surface from shell elements”功能將變形后的門框作有限元模型生成幾何模型導出備用。

由于在與橡膠密封圈的接觸分析中,門框僅接觸的一側參與計算,故僅采用該側門框進行后續接觸分析。門框幾何模型如圖5所示。

圖5 變形門框幾何模型Fig.53D model of the deformed doorframe

由于目前沒有相應的變形衡準,且有限元軟件中得到的變形為總體變形。故采取以門框左下角節點3為參考點,選擇門框邊緣節點相對參考點的位移來反映門框處的相對變形情況。部分節點的相對位移如表1所示。x方向為船體縱向,y方向為船體橫向,z方向為船體垂向。

表1 門框節點的相對位移Tab.1Relative displacement of doorframe at different nodes

可以看出門框邊界的相對變形主要集中在x方向,z方向次之,y方向最小,這與整船變形主要由總縱彎曲引起的原因相符。由于x方向是氣密門開閉的方向,也是氣密橡膠密封圈主要的壓縮方向,因而對氣密性影響較大。

2.2橡膠Mooney-Rivlin模型

由于橡膠材料的體積彈性模量非常高,因此可假定橡膠材料具有不可壓縮性,而且在未應變狀態下是各向同性的。本文采用應用廣泛的2參數橡膠Mooney-Rivlin本構方程模擬橡膠材料,其應變函數如下式所示:

式中:W為修正的應變勢能;C01,C10為材料常數;I1,I2為應力張量的第1、第2不變量。

在不進行復雜的材料實驗測試情況下,可通過經驗公式求得Mooney-Rivlin模型材料常數。由于橡膠的靜剪切模量是橡膠元件設計中最基本的參數之一,其與橡膠硬度及成分有關(其中最主要的決定因素是橡膠硬度)。對于硬度相同成分不同的橡膠材料,其值之差不超過10%。因此,可以根據橡膠硬度,并利用經驗公式確定Mooney-Rivlin模型的材料常數[19]。

測得橡膠材料的邵氏硬度HA,代入下式:

在小變形的情況下,根據橡膠材料的不可壓縮性,泊松比μ=0.5。橡膠材料的剪切模量和楊氏模量與Mooney-Rivlin模型常數有如下關系:

根據法國PAULSTRA公司給出的不同橡膠硬度下支座的載荷—變形曲線進行有限元建模,并與實測值對比,以確定不同硬度下材料常數的最佳取值。3種橡膠硬度C01/C10值如表2所示。對3個C01/C10值分段進行線性擬合,可以得到各硬度下C01/C10值。

表2 不同硬度下橡膠材料C01/C10的值Tab.2TheC01/C10of rubber materials with different hardness

由式(2)、式(3)和表2計算得到不同橡膠硬度下Mooney-Rivlin模型的材料常數C01和C10如表3所示。

表3 不同硬度下橡膠材料Mooney-Rivlin常數Tab.3The Mooney-Rivlin constant of rubber materials with different hardness

2.3接觸分析

影響橡膠密封性能的因素很多,判別密封性能的參數也有很多選擇,如:壓縮量、接觸面積、接觸壓力等。本文以密封圈與門框的接觸壓力為判別準則,當接觸壓力小于安全值時認為密封失效。由于本艦氣密門位于鎖氣室,而為了防止外側空氣進入鎖氣室,鎖氣室內設有氣閘保證內外存在壓差,一般壓差取P=250 Pa。若σn為接觸應力,則認為當σn<P時,密封失效。

本算例中,門框模型為變形的Shell單元生成的幾何模型,密封圈使用簡單的O形密封圈,直徑D=30 mm,采用實體單元劃分網格,共得到約96 000個單元。由于橡膠接觸問題,變形較大,對單元質量要求較高,模型中所有單元均為六面體單元。

橡膠材料為丁腈橡膠(NBR),采用Mooney-Rivlin兩參數模型模擬,由提供的橡膠硬度HA=43插值得到模型材料常數,由表3可知:C10=0.289,C01=0.028。忽略門框與密封圈之間的摩擦。接觸分析有限元模型如圖6所示。

圖6 三維接觸分析有限元模型Fig.63D finite element model for contact analysis

在接觸分析中,O形密封圈作為變形體,密封圈與門相連處暫以固支處理,門框作為剛體向密封圈移動,壓縮量取10 mm進行接觸分析。計算結果如圖7所示。

圖7 密封圈接觸壓應力云圖Fig.7The contours of contact stress distribution on seal ring

沿順時針選取橡膠密封圈與門框接觸一側的最外層節點。由于該處節點最先與門框接觸,接觸應力最大,能體現該橫截面處的最大接觸應力。由PATRAN后處理導出沿該路徑的接觸應力分布曲線,如圖8(a)所示。圖8(b)為使用未變形門框的計算結果。

在圖8(a)所有節點中選取接觸應力最小的節點node 79370,接觸應力σmin=23.01 MPa,因此任意截面都滿足σn≥P,認為密封性能良好。

由圖8(a)與圖8(b)對比可以看出,由于門框的變形,接觸壓力波動較大,這可能會引起2種危險結果:1)當壓縮量較小時,有部分區域接觸壓力小于安全值,引起泄露;2)當壓縮量過大時,部分橡膠材料應力超過屈服極限,容易出現裂紋,造成壽命降低,甚至有斷裂的危險。

圖8 接觸壓應力分布曲線Fig.8Contact stress distribution curves

3 優化措施

對于上述情況可采取2種思路防止危險的發生:一是從優化結構設計上考慮,減小門框變形;二是從橡膠密封圈上考慮。優化措施如下:

1)增加艙壁板厚,或在該處圍壁與上、下甲板間增設肘板等結構加強方案,改善該處由于總縱彎曲產生的沿船長方向的相對變形。

2)在該處縱橫艙壁交界處設計一些彈性連接,減小門框結構參與總縱彎曲的程度,從而減小門框處的相對變形。

3)改變橡膠密封圈與門框的接觸方式,避免沿船長方向的相對變形對密封性能的影響,如采用橡膠密封圈與門框內側壁接觸。

對于其他工況、不同位置的氣密門,可以采用類似的方法進行分析,選取最不利工況,結合各自變形的特點進行結構設計。

4 結語

本文開展了評估艦船集防區對外通道結構變形及其密封結構氣密性的方法研究。借鑒整船結構強度的有限元直接計算,結合設計波法,預報全船總體變形和局部變形,提取變形后的密封結構進行有限元仿真,分析其密封性。在建造前對氣密門的氣密性進行仿真分析,可以為密封結構的優化設計提供參考。

針對集防區對外通道中門框變形的特點,提出了合理的密封改善措施。在結構優化方面,認為應在保證密封性能的前提下對結構進行優化,減小結構變形,從而提高對外通道密封設計的可靠性。在改變結構、調整局部剛度的同時,還需要兼顧結構強度的變化,以及結構重量的因素,進行多方面綜合考慮,最終得到合理的優化方案。

[1]ELDRIDGE J.Maritime CBRN countermeasures-a clink in the armour[J].Jane's Defense Weekly,2005(2):5-11.

[2]DAVIS R.Chemical and biological warfare:use of collective protection on vehicles,aircraft,and ship:GAO/NSIAD-91-273 FS[R].Washington:US General Accounting Office,1991:11.

[3]HAUPT R E.Nuclear biological chemical(NBC)warfare protection in the German navy[J].Naval Engineers Journal,1977,89(5):26-32.

[4]劉虹,劉飛,王斌.海上核化生安全威脅與水面艦艇集體防護[J].艦船科學技術,2011,33(7):150-155. LIU Hong,LIU Fei,WANG Bin.Nautical nuclear chemical and biological security threat and collective protection of surface ship[J].Ship Science and Technology,2011,33(7):150-155.

[5]FEENEY J J.Chemical contamination survivability of equipment[J].Naval Engineers Journal,1989,101(2):71-75.

[6]HORNING H K.Performance assessmentoftheDDG-51 collective protection system[J].Naval Engineers Journal,1994,106(2):118-125.

[7]程代云,李戰國,劉江歌.化學防護裝備研制50年回顧[J].防化研究,2003(1):6-12.

[8]李旭霞,張鄆城.外軍艦船對核生化武器的防護[J].中華航海醫學與高氣壓醫學雜志,2004,11(2):124-126.

[9]中國船舶工業總公司第七研究院第七一八研究所.艦船對核武器、生物武器和化學武器的防護要求:GJB 2327-95[S].北京:國防科學工業技術委員會,1995.

[10]耿保陽,謝永和,劉瑩瑩,等.3 600 DWT超規范油船艙段結構有限元分析[J].船舶工程,2014,36(4):13-16,20. GENG Baoyang,XIE Yonghe,LIU Yingying,et al. FEA of structural strength for cabin of 3 600 DWT supergauge tanker[J].Ship Engineering,2014,36(4):13-16,20.

[11]海軍規范所,船總701所,船總719所,等.艦船通用規范:GJB 4000-2000[S].北京:總裝備部軍標出版發行部,2000.

[12]中國船級社.鋼質海船入級規范[S].北京:人民交通出版社,2012.

[13]戴仰山,沈進威,宋競正.船舶波浪載荷[M].北京:國防工業出版社,2007.

[14]顧永寧,滕曉青,戴立廣,等.大開口船波浪載荷長期預報和彎扭強度整船有限元分析[J].中國造船,1998(2):63-70. GU Yongning,TENG Xiaoqing,DAI Liguang,et al. Long term prediction of wave loading and bending-torsion strength analysis of open hatch ship on entire ship FE model[J].Shipbuilding of China,1998(2):63-70.

[15]馮國慶,劉相春,任慧龍.基于PCL語言的波浪壓力自動加載方法[J].船舶力學,2006,10(5):107-112. FENG Guoqing,LIU Xiangchun,REN Huilong.A PCL based automatic loading approach for wave pressures[J].Journal of Ship Mechanics,2006,10(5):107-112.

[16]劉偉,劉志群,趙曉博,等.飛機艙門密封結構的氣密可靠性研究[J].材料工程,2012(4):47-51. LIUWei,LIUZhiqun,ZHAOXiaobo,etal. Air-tight reliability analysis of aero-port seal structure[J].Journal of Materials Engineering,2012(4):47-51.

[17]王勇.轎車門體彈性仿真及密封性研究[D].武漢:武漢理工大學,2008.

[18]王偉,趙樹高.橡膠O形密封圈的接觸變形及應力分析[J].彈性體,2006,15(4):28-31. WANG Wei,ZHAO Shugao.Analysis on contact deformation and stress for rubber O-sealing ring[J]. China Elastomerics,2006,15(4):28-31.

[19]周明帥,李天勻,朱翔,等.水潤滑橡膠艉軸承橡膠軸瓦硬度分區取值后的接觸性能研究[J].中國艦船研究,2014,9(2):62-67. ZHOU Mingshuai,LI Tianyun,ZHU Xiang,et al. The contact properties of water lubricated rubber stern bearings after the differentiation of rubber bush hardness[J].Chinese Journal of Ship Research,2014,9(2):62-67.

Numerical simulation of the deformation and airtightness of the external passageway structure of naval ships'collective protection system

YANG Guang1,2,ZHAO Shanghui3,LI Tianyun1,2,ZHU Xiang1,2
1 School of Naval Architecture and Ocean Engineering,Huazhong University of Science and Technology,Wuhan 430074,China
2 Hubei Key Laboratory of Naval Architecture and Ocean Engineering Hydrodynamics,Wuhan 430074,China
3 Marine Design and Research Institute of China,Shanghai 200011,China

In this paper,an evaluation method for the deformation and airtightness of the external passageway of a ship collective protection system is presented,where the seal structure like doors and windows is selected as the research object.Firstly,long-term prediction of ship's cross-sectional wave-included loads is conducted based on spectral analysis.Next,the parameters of designed waves are determined,and with FEM software MSC/PATRAN,a FEM model of the ship is established.The designed wave pressure is automatically applied to the ship's wet surface elements through PCL language.After obtaining the result of the deformation of the external passageway,a FEM model of sealed structures is built,where the 3D model of deformed doorframe are calculated based on the result of the general ship FEM model.Finally,contact analysis is conducted between sealing strip and the doorframe,and the result under normal contact stress is compared with the given value to evaluate the airtightness of the sealed structure.It can be concluded that the proposed method serves as a good reference to the ship design that has similar requirements.

collective protection system;deformation;design wave approach;sealed structure;airtightness

U667.3

A

10.3969/j.issn.1673-3185.2016.05.011

2015-10-29網絡出版時間:2016-9-21 13:38

楊光,男,1990年生,碩士生。研究方向:船舶結構力學。E-mail:yang1990824@163.com

李天勻(通信作者),男,1969年生,博士,教授,博士生導師。研究方向:結構振動與噪聲控制。E-mail:ltyz801@mail.hust.edu.cn

猜你喜歡
有限元變形結構
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
談詩的變形
中華詩詞(2020年1期)2020-09-21 09:24:52
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
“我”的變形計
例談拼圖與整式變形
會變形的餅
論《日出》的結構
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
磨削淬硬殘余應力的有限元分析
基于SolidWorks的吸嘴支撐臂有限元分析
主站蜘蛛池模板: 国产无码高清视频不卡| 国产成人三级| 亚洲日韩久久综合中文字幕| 中文字幕亚洲专区第19页| 777国产精品永久免费观看| 成人中文在线| 无遮挡国产高潮视频免费观看| 欧美性精品不卡在线观看| 夜夜操狠狠操| 欧美精品色视频| 国产综合亚洲欧洲区精品无码| 一级在线毛片| 1级黄色毛片| 亚洲综合激情另类专区| 亚洲精品成人片在线观看| 在线精品亚洲国产| 搞黄网站免费观看| 自拍偷拍欧美日韩| 青青草国产在线视频| 亚洲高清无码久久久| 四虎永久在线精品国产免费 | 国产在线精品人成导航| 成人年鲁鲁在线观看视频| 中文字幕色在线| 网久久综合| 福利小视频在线播放| 老司机精品99在线播放| 亚洲一区网站| 国产精品久久久久久久久kt| 茄子视频毛片免费观看| 日韩午夜伦| 美女潮喷出白浆在线观看视频| 日韩一区二区三免费高清| 在线观看热码亚洲av每日更新| 国产91线观看| www.99精品视频在线播放| 97青青青国产在线播放| 日本伊人色综合网| 久久久成年黄色视频| 91精品专区| 99热这里只有免费国产精品 | 国产精品美女网站| 波多野结衣一区二区三区四区| 国产偷倩视频| 波多野结衣无码AV在线| 三区在线视频| 国产农村妇女精品一二区| 成人在线亚洲| 色偷偷男人的天堂亚洲av| 国产福利在线免费| 国产精品永久不卡免费视频| 国产丰满大乳无码免费播放| 中文字幕伦视频| 日本精品影院| 久久精品国产精品青草app| 亚洲精品中文字幕无乱码| 亚洲精品777| 国产美女精品人人做人人爽| 白丝美女办公室高潮喷水视频| 国产高清在线精品一区二区三区| 欧美黄色a| 国产一区二区三区免费观看| 69av在线| 综合色亚洲| 影音先锋丝袜制服| 人与鲁专区| 国产精品太粉嫩高中在线观看 | 国产成人精品18| 97亚洲色综久久精品| 国产日韩欧美在线视频免费观看 | 日韩天堂网| 久久频这里精品99香蕉久网址| 亚洲大学生视频在线播放| 国产精品一区在线麻豆| 国产精品部在线观看| 全部免费毛片免费播放| 精品亚洲欧美中文字幕在线看| 国产在线无码av完整版在线观看| 91成人试看福利体验区| 久久久久久久97| 99999久久久久久亚洲| 日韩高清中文字幕|