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

    彈性楔形體入水砰擊載荷及結(jié)構(gòu)響應的理論計算與數(shù)值模擬研究*

    2021-12-03 09:06:30王一雯吳衛(wèi)國
    爆炸與沖擊 2021年11期
    關(guān)鍵詞:楔形邊界條件極值

    王一雯,鄭 成,吳衛(wèi)國

    (武漢理工大學綠色智能江海直達船舶與郵輪游艇研究中心,湖北 武漢 430063)

    結(jié)構(gòu)物入水砰擊問題在船艏艉入水、多體船濕甲板砰擊、空投魚雷入水、水上飛機降落、救生艇拋落等情況下廣泛存在,針對入水砰擊問題的研究對船舶以及航空航天等領(lǐng)域具有重要的工程意義。在船舶、海洋平臺以及相關(guān)航空結(jié)構(gòu)物的設計中,對砰擊載荷及其結(jié)構(gòu)響應進行合理可靠的評估是極為重要的。由于砰擊載荷的瞬態(tài)特性以及較高的幅值特性,對于較小斜升角的薄板結(jié)構(gòu)以及復合材料夾層結(jié)構(gòu)需考慮砰擊載荷對彈性結(jié)構(gòu)的水彈性效應,結(jié)構(gòu)變形使相對砰擊速度降低并且影響流場特性。而在準靜態(tài)結(jié)構(gòu)入水砰擊載荷及其結(jié)構(gòu)響應的研究中,采用解耦的方法進行計算,在對剛體結(jié)構(gòu)砰擊載荷計算的基礎(chǔ)上,將砰擊載荷加載到彈性結(jié)構(gòu)體以求解結(jié)構(gòu)響應。與準靜態(tài)結(jié)構(gòu)入水砰擊載荷及結(jié)構(gòu)響應計算方法相比,水彈性方法可考慮結(jié)構(gòu)彈性效應與砰擊載荷的耦合作用,并可直接預報彈性結(jié)構(gòu)入水砰擊過程的結(jié)構(gòu)響應。

    彈性結(jié)構(gòu)物入水砰擊問題的水彈性分析方法可分為基于砰擊載荷模型的解析方法和基于計算流體力學的流固耦合數(shù)值模擬方法。Kvalsvold 等[1]和Faltinsen 等[2]在雙體船的濕甲板砰擊問題中率先針對彈性加筋板的砰擊載荷及結(jié)構(gòu)響應進行理論和試驗研究,將Wagner 入水理論推廣到楔形體左右為正交異性板的入水砰擊情況,分別分析斜升角與入水速度對砰擊載荷及結(jié)構(gòu)響應的影響。Khabakhpasheva等[3]將歐拉梁模型與Wagner 模型進行流固耦合,并提出彈性楔形體的理論分析方法,分別針對簡支與線彈簧相連的邊界條件進行討論。Shams 等[4]將Wagner 模型推廣至不同邊界條件情況,采用混合邊界值法計及結(jié)構(gòu)變形的影響。Lu 等[5]基于線性化離散的Bernoulli 方程,通過對有限元法和邊界元法進行流固耦合,求解彈性楔形體結(jié)構(gòu)響應的水彈性效應。Stenius 等[6]采用任意拉格朗日-歐拉(arbitrary Lagrangian-Eulerian, ALE)流固耦合算法開展入水砰擊流固耦合計算,對彈性楔形體入水砰擊載荷及結(jié)構(gòu)響應進行討論,分析不同斜升角、邊界條件以及入水速度的影響,并與準靜態(tài)計算結(jié)果進行對比。Panciroli 等[7-8]針對彈性楔形體結(jié)構(gòu)進行了水彈性落體砰擊模型實驗,并通過高速攝像技術(shù)捕捉入水過程中彈性結(jié)構(gòu)及自由液面的變形,采用粒子圖像測試技術(shù)捕捉自由液面及流場的運動;除此之外,對0°~50°斜升角的鋁板及玻璃纖維彈性楔形板在不同入水速度下的結(jié)構(gòu)應變進行測試,并與有限元-無網(wǎng)格光滑粒子流固耦合方法進行對比;發(fā)現(xiàn)在結(jié)構(gòu)完全浸入液面時間小于自身固有周期時長的情況下水彈性效應較顯著。由于復合材料具有較好的抗沖擊性,被廣泛運用到高性能船舶、多體船以及艦船中,金屬夾層板作為非傳統(tǒng)結(jié)構(gòu)形式,其砰擊載荷特性及結(jié)構(gòu)響應近年來也得到研究者們的關(guān)注。Hassoon等[9-11]針對不同厚度的乙烯基編織玻璃纖維層合板以及玻璃纖維/乙烯基面板與聚氯乙烯泡沫夾層板的10°楔形體模型,對其勻速入水過程中的砰擊載荷及結(jié)構(gòu)響應開展了實驗研究,并且采用無網(wǎng)格光滑粒子數(shù)值模擬方法和Abaqus 中的歐拉-拉格朗日耦合算法進行了數(shù)值計算和對比分析,討論了復合材料材料特性、結(jié)構(gòu)剛度、入水速度以及不同邊界條件對砰擊載荷及其結(jié)構(gòu)動態(tài)變形響應的影響以及不同的損傷模式,并建立了復合材料連續(xù)損傷模型,包括層合板Hashine 失效準則以及夾層板的Christensen 失效準則。綜合以上分析可知,基于計算流體力學的流固耦合數(shù)值模擬方法需對模型進行精細化處理,合理選取參數(shù),并且求解耗時較長;而基于砰擊載荷模型解析計算方法可高效求解,在結(jié)構(gòu)方案優(yōu)化階段具有明顯優(yōu)勢。

    本文中,分別采用ALE 流固耦合數(shù)值模擬方法和理論模型解析方法對彈性楔形體結(jié)構(gòu)的砰擊載荷及其結(jié)構(gòu)響應開展計算并對比分析,實現(xiàn)彈性楔形體結(jié)構(gòu)入水砰擊過程結(jié)構(gòu)響應的快速預報;與剛體結(jié)構(gòu)進行對比,分析水彈性對砰擊載荷作用下結(jié)構(gòu)響應特性的影響;與Lu 等[5]所開展的邊界元計算結(jié)果進行對比,驗證本文中所提出的彈性結(jié)構(gòu)砰擊載荷及結(jié)構(gòu)響應預報方法的高精度和可靠性。

    1 二維彈性楔形體砰擊載荷和結(jié)構(gòu)響應的理論計算方法

    船底結(jié)構(gòu)由船底板、龍骨、縱骨、縱桁以及實肋板等組成,針對等厚度以及有限寬度梁結(jié)構(gòu)而言,其結(jié)構(gòu)沿橫向方向可簡化為楔形體結(jié)構(gòu),該結(jié)構(gòu)入水砰擊過程可簡化為楔形體砰擊靜水面的水動力沖擊問題。在二維不可壓縮理想流體域內(nèi),對稱彈性楔形體以勻速垂直進入靜水面,楔形體和流場均關(guān)于Oz軸對稱,如圖1 所示,圖中為楔形體斜升角,L為楔形體物面長度,b為楔形體板厚,c(t)為浸濕半寬,t為時間。

    圖1 彈性楔形體砰擊載荷及結(jié)構(gòu)響應示意圖Fig. 1 Schematics of slamming load and structural response of a flexible wedge

    由于二維彈性楔形體結(jié)構(gòu)形式及入水運動過程均關(guān)于Oz軸對稱,可針對其一半結(jié)構(gòu)模型及流體域進行分析,楔形體結(jié)構(gòu)可簡化為歐拉梁,在不計及剪切和轉(zhuǎn)動慣量的情況下,其自由振動方程為:

    式中:E為彈性模量,I為剖面慣性矩,m為長度質(zhì)量,w為梁的撓度,t為時間,x′為梁局部縱向坐標,p(x′,w,t)為梁局部縱向坐標x′處的局部水動力壓力。

    基于模態(tài)疊加法,可通過前n階的主坐標an求得梁擾度w(x′,t)為[3]:

    在彈性楔形體入水砰擊過程中,根據(jù)MLM (modified Logvinovich model)、GWM (generalized Wagner model)和OLM (original Logvinovich model)砰擊載荷模型解析解[12],可知沿楔形體物面各處的砰擊載荷可表示為:

    針對上述不同砰擊載荷模型進行對比分析,分別基于MLM、GWM 和OLM 砰擊載荷模型對彈性楔形體結(jié)構(gòu)響應開展研究,探討基于不同模型的理論解析解對結(jié)構(gòu)變形的影響。長L=0.4 m、板厚b=8 mm 的彈性楔形體采用船體鋼材,其密度為7 850 kg/m3,彈性模量為210 GPa,泊松比為0.3;同時水的密度為1 000 kg/m3。入水過程中忽略流體表面張力及重力的影響,同時流體為無旋無黏的不可壓縮理想流體,并且流場域為無限水深,起始時刻流體為靜止自由狀態(tài)。在1 m/s 入水速度下,對不同斜升角的彈性楔形體采用不同理論模型,所得x′=0.5L處的結(jié)構(gòu)變形w時間歷程結(jié)果如圖2 所示。

    從圖2 可看出:不同斜升角下,這3 種模型的計算結(jié)果具有相同的變化趨勢;在結(jié)構(gòu)入水階段,3 種模型的計算結(jié)果較接近,而隨著結(jié)構(gòu)入水運動過程中射流的分離,其偏差增大。在x′=0.5L處結(jié)構(gòu)變形時間歷程結(jié)果中,基于MLM 模型的理論計算結(jié)果均最大,而基于GWM 模型的理論計算結(jié)果則最小。在10°的斜升角下,三者的偏差較小;而隨著斜升角的增大,三者間的偏差顯著增大。在45°斜升角楔形體中,基于GWM 模型計算的結(jié)構(gòu)變形極值均遠小于其他兩種模型的計算結(jié)果,主要是由于GWM 模型中將飛濺高度線性化處理所致的,并且偏差也隨著斜升角的增大而增大。可見,GWM 模型并不適用于較大斜升角的彈性結(jié)構(gòu)砰擊載荷分析。而OLM 模型中忽略了結(jié)構(gòu)物面斜升角高階項的影響,在物面射流根部附近的砰擊壓力峰值顯著低于MLM 模型結(jié)果,因此在大斜升角情況下低估了結(jié)構(gòu)變形。綜合分析可知,在砰擊載荷作用下的彈性結(jié)構(gòu)變形分析中,MLM 模型可適用于分析斜升角范圍較大的彈性楔形體入水砰擊過程。

    圖2 采用不同模型得到的板厚為8mm、以 1 m/s 的速度垂直入水的彈性楔形體 x′=處結(jié)構(gòu)變形的時間歷程Fig. 2 Time series of the structural response at x′=of the flexible wedges with different deadrise angles and the plate thickness of 8 mm at the vertical water-entry velocity of 1 m/s calculated by different models

    10°斜升角彈性楔形體在t=20, 40 ms 時沿物面的結(jié)構(gòu)變形分布情況如圖3 所示。由圖3 可見,基于MLM、GWM 以及OLM 模型的理論解析解中結(jié)構(gòu)變形分布總體變化趨勢一致,但是GWM 整體量值偏低,與圖2 中所得結(jié)果一致。而隨著結(jié)構(gòu)入水繼續(xù)深入自t=20 ms 至t=40 ms 時刻,其結(jié)構(gòu)變形最大值所處位置則從x′/L=0.46 處移至x′/L=0.50 處,即結(jié)構(gòu)變形最大值所處位置隨結(jié)構(gòu)入水深入而上移。結(jié)合圖2 可知,當t=40 ms 時結(jié)構(gòu)變形升至其極值并且位于結(jié)構(gòu)x′/L=0.50處,因此需對結(jié)構(gòu)中點x′/L=0.50 處的結(jié)構(gòu)變形極值進行對比。

    圖3 t=20, 40 ms 時沿彈性楔形體物面的結(jié)構(gòu)變形分布(10°, b=8 mm)Fig. 3 Response distribution along the wedge structure at t=20, 40 ms (=10°, b=8 mm)

    綜合分析可知,在砰擊載荷作用下的彈性楔形體結(jié)構(gòu)響應分析中,MLM 砰擊載荷模型適用斜升角范圍較大。本文后續(xù)將采用基于MLM砰擊載荷模型的彈性結(jié)構(gòu)理論計算方法進行分析,并與ALE 流固耦合模擬方法進行對比。

    2 彈性楔形體的砰擊載荷

    通過LS-DYNA 軟件采用任意拉格朗日-歐拉流固耦合模擬方法對二維彈性楔形體結(jié)構(gòu)入水砰擊過程進行研究,討論模型精細化程度對砰擊載荷的影響。并且分別對不同斜升角、板厚以及邊界條件的彈性楔形體結(jié)構(gòu)在不同入水速度砰擊靜水面的情況開展數(shù)值模擬,以探究結(jié)構(gòu)水彈性效應對砰擊載荷及結(jié)構(gòu)響應特性的影響。

    2.1 任意拉格朗日-歐拉耦合算法及計算模型

    由于結(jié)構(gòu)物砰擊入水過程中伴隨著流體大變形以及大位移等特點,任意拉格朗日-歐拉耦合(ALE)算法兼具了Lagrange 算法和Euler 算法的優(yōu)勢,可在計算網(wǎng)格不發(fā)生畸變的情況下,通過網(wǎng)格重映射數(shù)值處理方法實現(xiàn)Lagrange 算法和Euler 算法相互結(jié)合,能夠有效追蹤物質(zhì)結(jié)構(gòu)的邊界并可對自由邊界和運動邊界進行求解,合理準確模擬結(jié)構(gòu)物入水砰擊問題。

    任意拉格朗日-歐拉流固耦合模擬方法中:彈性楔形體結(jié)構(gòu)單元為拉格朗日彈性單元,通過關(guān)鍵字*MAT_ELASTIC 描述;空氣域和水域則采用多物質(zhì)歐拉單元,采用關(guān)鍵字*MAT_NULL 描述;空氣域采用線性多項式狀態(tài)方程描述,水域則適用于Grüneisen 狀態(tài)方程描述。線性多項式狀態(tài)方程中氣體滿足γ 定律狀態(tài)方程,其壓力值p與體積的關(guān)系為:

    式中:C0、C1、C2、C3、C4、C5和C6為常數(shù);C4=C5=γ?1;γ 為比熱比;μ=1/V?1,V為相對體積;E0為初始體積內(nèi)能。

    壓縮狀態(tài)下的Grüneisen 狀態(tài)方程可通過沖擊速度定義為:

    針對本文中所述的空氣域和水域的狀態(tài)方程,參數(shù)見表1。

    表1 空氣域及水域狀態(tài)方程參數(shù)Table 1 Parameters for equations of state of air and water

    通過任意拉格朗日-歐拉算法和罰函數(shù)約束算法進行流固耦合,流體邊界處理為無反射邊界條件以實現(xiàn)無界流域來消除邊界的影響。在入水砰擊主要作用區(qū)域內(nèi)采用密集均勻網(wǎng)格,在較遠區(qū)域內(nèi)采用漸變型網(wǎng)格,如圖4 所示,可在保障計算精度的基礎(chǔ)上縮短計算時間及降低計算成本。水域和空氣域的密集區(qū)域分別為L4×L6=0.4 m×0.6 m 和L3×L6=0.2 m×0.6 m,整個水域和空氣域模型尺寸分別為L5×L2=1.35 m×0.9 m 和L5×L1=1.35 m×0.4 m。約束楔形體結(jié)構(gòu)的z軸方向位移,使其僅于xOy平面內(nèi)二維運動。在對稱面yOz平面處建立對稱邊界條件,則僅可建立一半結(jié)構(gòu)模型及流體域,可顯著縮短計算時間及降低計算成本。

    圖4 砰擊入水計算模型Fig. 4 The water-entry impact computation model with mesh generation

    2.2 網(wǎng)格驗證

    由于砰擊入水這類強非線性流固耦合問題的計算受模型精細化程度的影響極大,精細網(wǎng)格模型的計算精度較高但計算量過大從而影響計算效率,而網(wǎng)格尺寸過大的模型不僅不能模擬流體濺射以及自由液面變化的非線性過程,而且其計算精度顯著降低。因此,在結(jié)構(gòu)物入水砰擊模擬研究中,確定合理的網(wǎng)格尺寸是保障數(shù)值模擬結(jié)果可靠、有效及高精度的前提。

    針對長L=0.4 m 的彈性楔形體進行砰擊入水流固耦合數(shù)值模擬分析,分別對網(wǎng)格尺寸為4、2、1 mm 的10°彈性楔形體剖面進行數(shù)值計算并開展參數(shù)分析。不同網(wǎng)格尺寸模型的詳細信息如表2 所示,對比了其單元數(shù)量以及計算時長,3個模型中的時間步長均采用50 μs。同時其彈性結(jié)構(gòu)所受的無量綱砰擊力CF=F/(0.5ρv2Lsinβ)如圖5 所示,其中F為結(jié)構(gòu)所受的長度砰擊力,橫坐標vt/(Lsinβ)為無量綱時間。由圖5 可見:當模型網(wǎng)格尺寸為4 mm 時,入水初始階段存在高頻振蕩;而網(wǎng)格尺寸分別為1 mm 及2 mm 時,其砰擊壓力時間歷程較接近,呈較好的一致性,射流與物面分離后其剖面所受的砰擊力則迅速下降;當模型尺寸為1 mm時,相較于2 mm 網(wǎng)格尺寸,其模擬耗時顯著增長至11.69 倍??梢?,當網(wǎng)格尺寸為2 mm 時,可在保障計算精度的前提下顯著提高計算效率。因此,后續(xù)研究中均選用網(wǎng)格尺寸為2 mm。

    表2 不同網(wǎng)格尺寸模型信息Table 2 Three models with different mesh sizes

    圖5 不同網(wǎng)格尺寸的楔形體模型所受結(jié)構(gòu)砰擊力Fig. 5 Comparison of slamming forces of wedge models with different mesh sizes

    2.3 結(jié)構(gòu)砰擊力

    彈性結(jié)構(gòu)物入水砰擊過程相較于剛體入水情況更復雜,還需要考慮結(jié)構(gòu)物不同邊界條件的影響。分別對兩端簡支和兩端固支邊界條件進行分析,探討邊界條件對結(jié)構(gòu)砰擊力的影響。10°斜升角的彈性楔形體以2、4、6、8 m/s 速度的砰擊入水情況中,剛性楔形體和不同邊界條件的彈性楔形體所受砰擊力時間歷程如圖6 所示。圖6 中c 與s 分別表示兩端固支以及兩端簡支邊界條件。

    圖6 剛性和彈性楔形體所受到的無量綱砰擊力 (β=10°)Fig. 6 Dimensionless slamming forces on rigid and elastic wedges with β=10°

    可見,在入水初期,剛性楔形體和彈性楔形體所受砰擊力均成線性增大且較接近;而在入水砰擊中期,彈性體結(jié)構(gòu)所受砰擊力相較剛性體而言顯著偏低,并且隨著入水速度的升高,兩者間的差距隨之增大。而在射流分離后砰擊力均達峰值后則迅速降低,彈性體的砰擊壓力峰值以及增長速率均高于剛性體的并且時間滯后,隨著入水速度的升高,彈性體所受結(jié)構(gòu)砰擊力峰值明顯高于剛性體結(jié)構(gòu)。可見,此時結(jié)構(gòu)彈性變形效應對流場已產(chǎn)生顯著影響。兩種不同邊界條件下的彈性楔形體結(jié)構(gòu)砰擊力變化趨勢總體相同,在較低入水速度下兩者間差異較小并不顯著??梢姡噍^于結(jié)構(gòu)邊界條件而言,結(jié)構(gòu)砰擊力對入水速度的敏感程度更高。

    針對30°斜升角的剛性體和彈性體在不同入水速度下結(jié)構(gòu)所受砰擊力進行分析,并與10°斜升角楔形體的無量綱砰擊力結(jié)果進行對比,探究斜升角對彈性結(jié)構(gòu)所受砰擊力的影響。不同入水速度下30°斜升角剛性以及彈性楔形體結(jié)構(gòu)所受到的砰擊力如圖7 所示。與10°斜升角的楔形體結(jié)構(gòu)砰擊力變化趨勢相同的是,在較高入水速度作用下,彈性體結(jié)構(gòu)所受的砰擊力峰值均大于剛性體結(jié)構(gòu)所受的,隨著入水速度自2 m/s 升高至8 m/s,兩種結(jié)構(gòu)所受砰擊力峰值的差值自2.4%增大至10.9%。但30°斜升角楔形體結(jié)構(gòu)所受砰擊力均大幅降低,并且剛性楔形體結(jié)構(gòu)的砰擊載荷曲線與彈性楔形體結(jié)構(gòu)的砰擊載荷曲線更接近,且在低速入水砰擊過程中剛性體與彈性體間的砰擊力差異較小。直至入水速度升高6 m/s時,彈性體所受到的砰擊壓力峰值明顯高于剛性體所受到的砰擊力峰值。但相較于10°斜升角楔形體結(jié)構(gòu),剛體與彈性體結(jié)構(gòu)所受到的砰擊力峰值差值從134%縮小至110.9%。兩端簡支彈性楔形體在6 m/s入水速度下的無量綱砰擊力峰值自0.85 降至0.059??梢姡噍^于入水速度以及邊界條件,彈性體結(jié)構(gòu)所受砰擊力受物面斜升角的影響更顯著。

    圖7 剛性和彈性楔形體所受到的無量綱砰擊力 (β=30°)Fig. 7 Dimensionless slamming forces on rigid and elastic wedges with β=30°

    2.4 砰擊壓力及自由液面

    10°及30°斜升角的彈性楔形體入水過程中自由液面變化以及流場壓力分布如圖8~9 所示??梢?,在入水砰擊初期階段,射流沿物面濺射,彈性楔形體所受砰擊壓力峰值集中于射流根部,物面底部區(qū)域附近所受砰擊壓力則顯著降低。相較于30°斜升角,10°斜升角楔形體的砰擊壓力分布則更為集中于射流根部且幅值更高。隨著結(jié)構(gòu)物完全入水后,濺射射流與物面分離而后濺起水花回落至自由液面,物面所受砰擊壓力則迅速降低。相較于30°斜升角楔形體,10°斜升角楔形體底部區(qū)域所受砰擊壓力則顯著高于頂部區(qū)域。

    10°彈性楔形體入水過程中,t=4 ms 時刻浸濕液面位于結(jié)構(gòu)物面中部并未升至物面邊界,而當t=8 ms 時刻楔形體結(jié)構(gòu)完全入水。由圖8 可見,入水砰擊過程中結(jié)構(gòu)的彈性變形影響整個流域內(nèi)的壓力分布及自由液面形態(tài),從而增大了結(jié)構(gòu)物所受的砰擊力峰值及延長了各處的砰擊壓力作用時間,結(jié)構(gòu)彈性變形的影響使得楔形體底部局部斜升角增大,而楔形體頂端局部斜升角則減小。而在30°斜升角的彈性楔形體結(jié)構(gòu)入水過程中,由于斜升角的增大其結(jié)構(gòu)所受砰擊力顯著降低,無量綱砰擊力峰值自0.850 降至0.059,從而使得結(jié)構(gòu)變形響應明顯降低。由于斜升角的增大對其結(jié)構(gòu)彈性變形、流場壓力分布以及自由液面均有較大影響,針對物面不同位置處的砰擊壓力時間歷程開展對比分析,對物面L/4、L/2、3L/4 以及L處(分別記為點1、點2、點3 及點4)的砰擊壓力時間歷程對比如圖10~11 所示,圖中縱坐標Cp=p/(0.5ρv2)為無量綱砰擊壓力,橫坐標為無量綱時間。

    圖8 10°斜升角彈性楔形體入水過程中自由液面變化及流場壓力分布(v=8 m/s)Fig. 8 Fluid evolution and pressure distribution during the process of the elastic wedge with β=10° entering the water (v=8 m/s)

    圖9 30°斜升角彈性楔形體入水過程中自由液面變化及流場壓力分布(v=8 m/s)Fig. 9 Fluid evolution and pressure distribution during the process of the elastic wedge with β=30° entering the water (v=8 m/s)

    圖10 剛性和彈性楔形體點1~4 處的無量綱砰擊壓力-時間歷程(β=10°,v=6 m/s)Fig. 10 Dimensionless slamming pressure-time history curves at points 1?4 of the rigid and elastic wedges (β=10°, v=6 m/s)

    與剛性楔形體結(jié)構(gòu)相比,砰擊壓力起始作用時間、砰擊壓力峰值以及變化趨勢受結(jié)構(gòu)彈性效應影響較大。彈性結(jié)構(gòu)由于結(jié)構(gòu)變形的影響自由液面抬高速度較慢,使砰擊壓力起始時間滯后于剛性體。相較于兩端固支邊界條件的彈性楔形體,兩端簡支邊界條件的彈性體楔形體遭遇砰擊壓力的起始時間稍許滯后。相對于其他位置而言,P1 點處兩者時滯相差則較小,這是由于在入水砰擊初期,在較短時間內(nèi)砰擊壓力作用下結(jié)構(gòu)變形對流場的影響有限。而隨著結(jié)構(gòu)入水深入,結(jié)構(gòu)水彈性效應對結(jié)構(gòu)物上部區(qū)域處的砰擊壓力影響則加劇。

    而對于剛性體結(jié)構(gòu),砰擊壓力峰值均處于遭遇砰擊載荷的起始階段,物面各點砰擊壓力均在遭遇流體沖擊瞬時達到峰值,靠近物面頂端砰擊壓力衰減較快。然而對于彈性體的點1 及點2 處,砰擊壓力雖在遭遇砰擊載荷瞬時到達峰值,但其砰擊載荷并未隨著結(jié)構(gòu)物入水過程而衰減,并在入水的中后期出現(xiàn)次峰值,該兩點處固支邊界條件結(jié)構(gòu)的Cp分別為34.5 和35.4,分別在vt/(Lsinβ)為0.79 和0.75 下達到次峰值30.5 和30.8,分別占主峰值的89.2%和87.0%。而簡支邊界條件結(jié)構(gòu)的無量綱砰擊壓力峰值分別為38.8 和43.0,分別在vt/(Lsinβ)為0.78 和0.70 下達到次峰值30.3 和30.6,分別占主峰值的78.1%和71.2%。除此之外,簡支邊界條件的砰擊載荷峰值均稍大于固支邊界條件,并且峰值時刻均滯后于固支邊界條件。相對于剛體結(jié)構(gòu)的砰擊載荷峰值而言,點1 處兩彈性體模型均偏小且分別為剛性體的73.3%和82.4%。而點4 處兩彈性體模型的砰擊壓力峰值則均大于剛性體且分別為剛性體的112.2%和114.1%??梢?,水彈性效應和結(jié)構(gòu)邊界條件對砰擊壓力峰值及其沿物面分布情況均產(chǎn)生顯著影響,且結(jié)構(gòu)邊界條件的影響更大。在靠近結(jié)構(gòu)底端的點1 處,由于彈性結(jié)構(gòu)變形的影響使局部斜升角變大,因此彈性結(jié)構(gòu)所受砰擊壓力峰值較剛體結(jié)構(gòu)有所降低且存在明顯的次峰值現(xiàn)象。而在靠近結(jié)構(gòu)頂端的點2 處,由于彈性結(jié)構(gòu)變形的影響使得局部斜升角變小,因此彈性結(jié)構(gòu)所受砰擊壓力峰值較剛體結(jié)構(gòu)顯著增大。

    圖11 剛性和彈性楔形體點1~4 處的無量綱砰擊壓力時間歷程(β=30°,v=6 m/s)Fig. 11 Dimensionless slamming pressure-time history curves at points 1?4 of the rigid and elastic wedges (β=30°, v=6 m/s)

    與10°斜升角結(jié)果對比可知,30°斜升角的剛性體與彈性體的砰擊壓力時域曲線的變化趨勢較接近。相較于剛形體而言,彈性體物面點3 及點4 處的砰擊壓力峰值均偏大。相較于兩端固支邊界條件而言,兩端簡支邊界條件下點3 及點4 處的砰擊壓力峰值均偏大。值得注意的是,楔形體物面中點(點2)處,由于彈性體結(jié)構(gòu)變形的 影響,遭遇的砰擊壓力峰值持續(xù)時間則遠長于剛性體,直至楔形體入水砰擊后期迅速衰減。

    3 砰擊載荷作用下彈性楔形體的結(jié)構(gòu)響應

    3.1 板厚的影響

    不同板厚的30°斜升角彈性楔形體物面中心點處的結(jié)構(gòu)變形如圖12 所示,不同斜升角和板厚的彈性楔形體結(jié)構(gòu)響應峰值見表3。不同板厚的結(jié)構(gòu)變形時間歷程變化趨勢近似,在砰擊入水初期理論解析解稍高于ALE 水彈性模擬結(jié)果以及BEM 計算結(jié)果,隨著板厚的增大,不同模擬結(jié)果間的差異隨之減小,且BEM 模擬計算時間歷程曲線的局部振動頻率也隨之升高。射流分離后由于結(jié)構(gòu)所受砰擊力削弱,結(jié)構(gòu)變形在砰擊力升至峰值后則迅速回落。板厚由11 mm 分別減至8 mm 和5 mm 后,ALE 水彈性模擬計算中結(jié)構(gòu)變形極值則分別由39 μm 激增至100 μm 和420 μm,板厚分別減至73%和45%后導致其結(jié)構(gòu)響應極值增大至2.57 倍和12.13 倍。板厚b=5 mm 的彈性楔形體結(jié)構(gòu)變形時間歷程結(jié)果中,理論解析解與ALE 水彈性數(shù)值計算結(jié)果初始階段吻合較好。砰擊入水初期理論解析解稍偏大于ALE 水彈性模擬解,但在入水砰擊中后期理論解析解極值則偏小,主要是由于理論解析解中對基于MLM 砰擊載荷速度勢高階項的簡化導致的。而BEM 模擬計算結(jié)果則均偏小,其偏差主要源于結(jié)構(gòu)水彈性效應及砰擊載荷的差異。隨著板厚的增大,結(jié)構(gòu)變形的理論解析解在射流分離之前均偏大于其他模擬方法的計算結(jié)果,但其結(jié)構(gòu)變形極值與ALE 水彈性模擬結(jié)果則吻合較好。盡管不同板厚的彈性結(jié)構(gòu)所受砰擊力較接近,但板厚的削弱對局部結(jié)構(gòu)的安全較不利,在抗砰擊結(jié)構(gòu)設計中需結(jié)合輕量化以及安全性的考慮,綜合分析砰擊載荷幅值以及結(jié)構(gòu)局部響應水平。

    表 3 兩端簡支的不同板厚的彈性楔形體處的結(jié)構(gòu)變形峰值Table 3 The maximum structural responses atof the supported wedges with different plate thickesses

    表 3 兩端簡支的不同板厚的彈性楔形體處的結(jié)構(gòu)變形峰值Table 3 The maximum structural responses atof the supported wedges with different plate thickesses

    b/mm wmax/mm理論解析解ALEBEM[5]β=10°β=30°β=45°β=10°β=30°β=45°β=30°β=45°5 6.9700.4200.2004.3200.4100.2100.3800.150 8 1.7100.1000.0491.5300.1000.0510.0950.035 110.6600.0390.0190.5900.4000.0200.0360.015

    圖12 不同板厚的楔形體x′ =0.5L 處的結(jié)構(gòu)變形 (β=30°)Fig. 12 Structural responses at x′=0.5Lof the wedges with different thicknesses (β=30°)

    3.2 斜升角的影響

    采用ALE 流固耦合數(shù)值模擬和理論計算方法分別對10°和45°斜升角的彈性楔形體結(jié)構(gòu)進行研究,并與上節(jié)中30°斜升角的結(jié)構(gòu)變形分析結(jié)果進行對比,探究不同斜升角對砰擊載荷作用下的結(jié)構(gòu)動態(tài)響應的影響。10°和45°斜升角彈性楔形體以2 m/s 勻速入水過程中,在砰擊載荷作用下的結(jié)構(gòu)動態(tài)響應如圖13~14 所示。

    圖13 不同板厚的楔形體x′=0.5L 處的結(jié)構(gòu)變形(β=10°)Fig. 13 Structural responses at x′=0.5Lof the wedges with different plate thickness (β=10°)

    圖14 不同板厚的楔形體x′=0.5L 處的結(jié)構(gòu)變形(β=45°)Fig. 14 Structural responses at x′=0.5Lof the wedges with different plate thicknesses (β=45°)

    在10°及30°斜升角結(jié)果中各種研究方法得到結(jié)果較一致,隨著斜升角增大至45°后理論解析結(jié)果與其他兩種數(shù)值方法間差異較顯著,結(jié)構(gòu)變形曲線增長速率增大,但結(jié)構(gòu)變形極值與ALE 流固耦合模擬結(jié)果吻合較好,其偏差主要是由于MLM 模型中并未考慮斜升角的高階項,非線性程度有限。而Lu 等[5]的BEM 方法模擬結(jié)果僅針對入水砰擊初期,并未對結(jié)構(gòu)完全入水的中后期階段結(jié)構(gòu)變形進行分析,其分析結(jié)果與ALE 流固耦合數(shù)值模擬結(jié)果較接近。8 mm 板厚的彈性楔形體斜升角自10°分別增大至30°和45°后,其結(jié)構(gòu)變形極值顯著減小至6.5%和3.3%,而10°斜升角楔形體的板厚自8 mm 增大至11 mm后,結(jié)構(gòu)變形極值僅減小了61.7%。對比分析可知,結(jié)構(gòu)動態(tài)響應受斜升角變化的影響程度更高,通過增大斜升角可有效減小局部結(jié)構(gòu)變形極值。

    3.3 邊界條件的影響

    對兩端固支邊界條件的彈性楔形體進行分析,與上述兩端簡支邊界條件情況進行對比,探究不同邊界條件對砰擊載荷作用下的彈性楔形體結(jié)構(gòu)動態(tài)響應的影響。板厚b分別為5、8、11 mm 的兩端固支彈性楔形體物面中心點處的結(jié)構(gòu)變形如圖15 所示,不同邊界條件的彈性楔形體結(jié)構(gòu)變形峰值見表4。

    表4 的不同板厚的彈性楔形體 x′=處的結(jié)構(gòu)變形峰值Table 4 The maximum structural responses atof the wedges with and different plate thicknesses

    表4 的不同板厚的彈性楔形體 x′=處的結(jié)構(gòu)變形峰值Table 4 The maximum structural responses atof the wedges with and different plate thicknesses

    b/mm wmax/μm理論解析ALE兩端簡支兩端固支兩端簡支兩端固支 542084.041096.0 810021.010024.0 11 397.94009.1

    圖15 兩端固支的不同板厚的彈性楔形體x′=0.5L 處的結(jié)構(gòu)變形(β=30°)Fig. 15 Structural responses at x′=0.5Lof the clamped wedges with different plate thicknesses (β=30°)

    可見,兩種研究方法在入水砰擊初期可得到較好的吻合,但至射流分離后ALE 流固耦合模擬結(jié)果中的結(jié)構(gòu)變形極值均偏高。相較于理論解析結(jié)果,不同板厚下的ALE 計算結(jié)果分別偏高14.7%、15.01%和15.46%,并且時間分別滯后19、16、16 ms,其偏差主要是來源于理論解析解中對于射流分離后速度勢高階項簡化的影響。相較于圖12 中兩端簡支邊界條件結(jié)果,盡管結(jié)構(gòu)變形變化趨勢一致,但由于邊界條件的影響,相同板厚下的結(jié)構(gòu)變形極值大幅減小至20%左右??梢?,結(jié)構(gòu)邊界條件也是影響砰擊載荷作用下結(jié)構(gòu)動態(tài)響應的重要因素之一。

    3.4 入水速度的影響

    分別采用ALE 流固耦合數(shù)值模擬和理論計算方法對不同入水速度砰擊下的彈性楔形體結(jié)構(gòu)變形極值進行研究,探究入水速度對結(jié)構(gòu)響應極值的影響。不同板厚以及不同邊界條件的30°斜升角彈性楔形體,在不同入水速度砰擊載荷作用下的結(jié)構(gòu)響應極值如圖16~17 所示。

    圖16 不同入水速度下兩端簡支楔形體的結(jié)構(gòu)變形極值 (β=30°)Fig. 16 The maximum responses of the supported wedges at different impact velocities (β=30°)

    在1 m/s 和2 m/s 入水速度下,兩種計算結(jié)果間的偏差較小,而隨著入水速度的升高其偏差更顯著,尤其是在板厚較小的情況下,理論解析結(jié)果顯著偏高于ALE 流固耦合計算結(jié)果。該誤差主要是由于小板厚以及高速砰擊情況下,理論計算方法對結(jié)構(gòu)變形對流場的影響考慮有限。而在低速入水及板厚增大至11 mm 后,結(jié)構(gòu)變形極值顯著減小并且兩者間誤差較小。相同板厚條件下,兩端固支邊界條件下與兩端固支邊界條件具有類似的變化趨勢,但是其結(jié)構(gòu)變形極值大幅降低至20%。可見,結(jié)構(gòu)邊界條件對砰擊載荷作用下結(jié)構(gòu)的水彈性效應具有較顯著的影響。除板厚b=5 mm 的彈性楔形體在入水速度v=8 m/s 情況下兩者偏差明顯,其他板厚及入水速度情況下兩種方法預報的結(jié)構(gòu)變形極值則吻合較好,相對誤差低至9.6%。可見,對于兩端固支邊界條件的彈性結(jié)構(gòu),采用理論計算方法可得到精度較高的結(jié)構(gòu)變形極值。

    圖17 不同入水速度下兩端固支楔形體的結(jié)構(gòu)變形極值(β=30°)Fig. 17 The maximum responses of the clamped wedges at different impact velocities (β=30°)

    表5 不同板厚的楔形體在不同入水速度下的水彈性因數(shù)Table 5 The hydroelastic factors of the wedges with different plate thicknesses at different water-enter velocities

    4 結(jié) 論

    針對彈性楔形體入水的砰擊載荷及結(jié)構(gòu)響應開展研究,提出了彈性楔形體在砰擊載荷作用下的理論解析計算模型,評估其結(jié)構(gòu)動態(tài)響應特性,基于模態(tài)疊加法實現(xiàn)了彈性楔形體入水過程中結(jié)構(gòu)響應的高精度快速預報。并且采用ALE 流固耦合數(shù)值模擬方法評估砰擊載荷作用下的結(jié)構(gòu)動態(tài)響應,并與BEM 邊界元數(shù)值計算結(jié)果進行對比,驗證其理論模型解析計算方法的適用性及有效性。分別針對不同的斜升角及板厚的彈性結(jié)構(gòu)進行分析,討論不同入水速度、斜升角、板厚和不同邊界約束條件對砰擊載荷特性和結(jié)構(gòu)動態(tài)響應特性的影響。

    相較于板厚而言,彈性結(jié)構(gòu)所受砰擊力和結(jié)構(gòu)動態(tài)響應對物面斜升角和入水速度更敏感,通過增大斜升角可有效降低結(jié)構(gòu)砰擊載荷和局部結(jié)構(gòu)變形極值,加強結(jié)構(gòu)邊界約束也可顯著降低結(jié)構(gòu)動態(tài)響應。在結(jié)構(gòu)斜升角較小和入水速度較高的情況下,其彈性效應引起的砰擊載荷變化更顯著。在考慮結(jié)構(gòu)水彈性效應時,彈性體所受到的砰擊壓力峰值顯著高于剛性體結(jié)構(gòu)所受到的砰擊壓力峰值并且到達峰值的時間滯后,隨著入水速度的升高,彈性體與剛性體的砰擊壓力峰值比值相較亦顯著增大。在滿足結(jié)構(gòu)輕量化要求下,提高斜升角可顯著降低結(jié)構(gòu)所受的砰擊載荷及結(jié)構(gòu)響應的幅值,有效保障抗砰擊結(jié)構(gòu)的安全性。

    猜你喜歡
    楔形邊界條件極值
    極值點帶你去“漂移”
    極值點偏移攔路,三法可取
    一類帶有Stieltjes積分邊界條件的分數(shù)階微分方程邊值問題正解
    帶有積分邊界條件的奇異攝動邊值問題的漸近解
    History of the Alphabet
    鋼絲繩楔形接頭連接失效分析與預防
    Eight Surprising Foods You’er Never Tried to Grill Before
    一類“極值點偏移”問題的解法與反思
    腹腔鏡下胃楔形切除術(shù)治療胃間質(zhì)瘤30例
    匹配數(shù)為1的極值2-均衡4-部4-圖的結(jié)構(gòu)
    久久久久久久大尺度免费视频| 成年av动漫网址| 1024视频免费在线观看| 天天躁狠狠躁夜夜躁狠狠躁| 三上悠亚av全集在线观看| 水蜜桃什么品种好| 日日摸夜夜添夜夜爱| 国产人伦9x9x在线观看 | 国精品久久久久久国模美| 校园人妻丝袜中文字幕| 国产精品女同一区二区软件| 777米奇影视久久| 性少妇av在线| 一级,二级,三级黄色视频| 在线天堂最新版资源| 国产午夜精品一二区理论片| 在线观看人妻少妇| 天天躁夜夜躁狠狠躁躁| 一级黄片播放器| 亚洲精品第二区| 日韩一区二区三区影片| 不卡av一区二区三区| 在线观看美女被高潮喷水网站| 久久综合国产亚洲精品| 啦啦啦在线观看免费高清www| 国产一区二区在线观看av| 最近最新中文字幕免费大全7| 日产精品乱码卡一卡2卡三| 精品酒店卫生间| 国产在线视频一区二区| 成人午夜精彩视频在线观看| 国产乱人偷精品视频| 丰满饥渴人妻一区二区三| 日韩免费高清中文字幕av| 亚洲欧美清纯卡通| 日韩欧美精品免费久久| 天美传媒精品一区二区| 国产视频首页在线观看| 国产片内射在线| 久久ye,这里只有精品| 精品久久蜜臀av无| 亚洲综合色惰| 中文字幕人妻熟女乱码| 一区在线观看完整版| 午夜av观看不卡| 91成人精品电影| 99久国产av精品国产电影| 王馨瑶露胸无遮挡在线观看| 大片免费播放器 马上看| 亚洲色图综合在线观看| 中国国产av一级| 精品少妇久久久久久888优播| 青春草亚洲视频在线观看| 久久热在线av| 人人妻人人添人人爽欧美一区卜| 91久久精品国产一区二区三区| 搡老乐熟女国产| 国产日韩一区二区三区精品不卡| 亚洲国产精品一区三区| 国产日韩一区二区三区精品不卡| xxxhd国产人妻xxx| 精品少妇久久久久久888优播| xxxhd国产人妻xxx| 99九九在线精品视频| 日韩免费高清中文字幕av| 少妇猛男粗大的猛烈进出视频| 精品卡一卡二卡四卡免费| 丝袜美足系列| 好男人视频免费观看在线| 丝袜美足系列| 嫩草影院入口| 五月天丁香电影| 97精品久久久久久久久久精品| 麻豆乱淫一区二区| 麻豆精品久久久久久蜜桃| 亚洲五月色婷婷综合| 婷婷色av中文字幕| av在线观看视频网站免费| 黄网站色视频无遮挡免费观看| 美女福利国产在线| 国产国语露脸激情在线看| 亚洲国产看品久久| 中文精品一卡2卡3卡4更新| 丝袜喷水一区| 成人18禁高潮啪啪吃奶动态图| 亚洲欧美精品自产自拍| 大话2 男鬼变身卡| 午夜福利在线免费观看网站| 国产精品二区激情视频| 亚洲中文av在线| 天堂中文最新版在线下载| 久久久久久久久久久免费av| 黑丝袜美女国产一区| 久久久久久久久久人人人人人人| 国产伦理片在线播放av一区| 尾随美女入室| 日韩精品免费视频一区二区三区| 日韩欧美一区视频在线观看| 亚洲精品一二三| 国产av一区二区精品久久| 国产成人一区二区在线| 亚洲人成77777在线视频| 成年美女黄网站色视频大全免费| 亚洲成人av在线免费| 香蕉丝袜av| 黑丝袜美女国产一区| 99精国产麻豆久久婷婷| 精品国产露脸久久av麻豆| 男女免费视频国产| 国产精品 国内视频| 亚洲激情五月婷婷啪啪| 街头女战士在线观看网站| 男人操女人黄网站| 一区二区av电影网| 国产成人免费观看mmmm| 国产精品秋霞免费鲁丝片| 免费在线观看黄色视频的| 亚洲国产av影院在线观看| 777米奇影视久久| 国产精品三级大全| www.熟女人妻精品国产| 久久久精品94久久精品| 成人18禁高潮啪啪吃奶动态图| 嫩草影院入口| 国语对白做爰xxxⅹ性视频网站| 欧美人与性动交α欧美软件| 夜夜骑夜夜射夜夜干| 国产 精品1| 天天躁夜夜躁狠狠躁躁| 亚洲伊人色综图| 日日啪夜夜爽| 国产日韩欧美在线精品| 人人妻人人爽人人添夜夜欢视频| 国产男女超爽视频在线观看| 最近最新中文字幕大全免费视频 | 不卡av一区二区三区| 亚洲欧美精品综合一区二区三区 | 久久精品熟女亚洲av麻豆精品| 丝袜喷水一区| 欧美老熟妇乱子伦牲交| 日本av手机在线免费观看| 欧美日韩亚洲国产一区二区在线观看 | 叶爱在线成人免费视频播放| 久久人妻熟女aⅴ| 永久免费av网站大全| 日产精品乱码卡一卡2卡三| 亚洲av.av天堂| 一区二区av电影网| 国产男女内射视频| 青春草亚洲视频在线观看| 啦啦啦视频在线资源免费观看| 日韩欧美精品免费久久| 亚洲欧美精品综合一区二区三区 | 久久国内精品自在自线图片| 男女边吃奶边做爰视频| 精品国产乱码久久久久久小说| 一本—道久久a久久精品蜜桃钙片| 日韩伦理黄色片| 成人国语在线视频| 最近最新中文字幕免费大全7| 丰满迷人的少妇在线观看| 婷婷色综合www| 亚洲五月色婷婷综合| 波多野结衣av一区二区av| 狠狠精品人妻久久久久久综合| 一区二区三区激情视频| videossex国产| 777久久人妻少妇嫩草av网站| 999精品在线视频| 精品99又大又爽又粗少妇毛片| 天天躁夜夜躁狠狠躁躁| 80岁老熟妇乱子伦牲交| 欧美激情极品国产一区二区三区| 天堂中文最新版在线下载| 观看av在线不卡| 自拍欧美九色日韩亚洲蝌蚪91| 国产片内射在线| 欧美人与性动交α欧美软件| 国产精品国产三级国产专区5o| 亚洲av在线观看美女高潮| 久久精品国产a三级三级三级| h视频一区二区三区| 18+在线观看网站| 叶爱在线成人免费视频播放| 男女边摸边吃奶| 久久久久久久久久久久大奶| av免费观看日本| 亚洲欧美成人综合另类久久久| 免费av中文字幕在线| 丰满迷人的少妇在线观看| 午夜91福利影院| 又黄又粗又硬又大视频| 在线免费观看不下载黄p国产| 黑人巨大精品欧美一区二区蜜桃| 青春草视频在线免费观看| 亚洲精品久久成人aⅴ小说| 黄色配什么色好看| 亚洲av日韩在线播放| 男人操女人黄网站| 国产男人的电影天堂91| 人人澡人人妻人| 国产精品嫩草影院av在线观看| 精品亚洲乱码少妇综合久久| 人妻一区二区av| www.av在线官网国产| 欧美亚洲 丝袜 人妻 在线| 极品人妻少妇av视频| 午夜激情久久久久久久| 男女高潮啪啪啪动态图| 国产日韩一区二区三区精品不卡| 91精品国产国语对白视频| 久久久精品区二区三区| 天天躁夜夜躁狠狠躁躁| 1024香蕉在线观看| 女性生殖器流出的白浆| 亚洲精品,欧美精品| 午夜av观看不卡| 亚洲欧美成人综合另类久久久| 97在线视频观看| 色婷婷久久久亚洲欧美| 高清在线视频一区二区三区| 亚洲国产毛片av蜜桃av| 波多野结衣av一区二区av| 亚洲中文av在线| 一级,二级,三级黄色视频| 久久久久精品久久久久真实原创| 国产男女内射视频| 久久久久人妻精品一区果冻| 日韩免费高清中文字幕av| 国产一级毛片在线| 久久毛片免费看一区二区三区| 人成视频在线观看免费观看| 纵有疾风起免费观看全集完整版| 性少妇av在线| 欧美日韩视频高清一区二区三区二| videossex国产| 自线自在国产av| 欧美激情极品国产一区二区三区| 国产精品成人在线| 中国三级夫妇交换| 国产av国产精品国产| 久久久久久免费高清国产稀缺| 久久久精品免费免费高清| 日本91视频免费播放| 亚洲成av片中文字幕在线观看 | 波多野结衣一区麻豆| 色哟哟·www| av不卡在线播放| 欧美成人午夜免费资源| 亚洲精品美女久久av网站| 日本av免费视频播放| 国产精品三级大全| 91在线精品国自产拍蜜月| 婷婷成人精品国产| 尾随美女入室| 成人毛片60女人毛片免费| 成人国语在线视频| 黄色毛片三级朝国网站| 亚洲国产成人一精品久久久| 国产伦理片在线播放av一区| 亚洲精品视频女| 成年女人毛片免费观看观看9 | 人妻人人澡人人爽人人| 精品99又大又爽又粗少妇毛片| 日本黄色日本黄色录像| 街头女战士在线观看网站| 久久久久久久精品精品| 亚洲国产欧美日韩在线播放| 久久久久国产网址| 亚洲天堂av无毛| 久久综合国产亚洲精品| 日韩欧美精品免费久久| 久久精品国产a三级三级三级| 成人亚洲欧美一区二区av| 精品一区二区三区四区五区乱码 | 91成人精品电影| 中文字幕人妻丝袜制服| 99久久中文字幕三级久久日本| 久久精品久久精品一区二区三区| 亚洲三区欧美一区| 中文字幕人妻丝袜一区二区 | av.在线天堂| 亚洲精品视频女| 2018国产大陆天天弄谢| 欧美国产精品一级二级三级| 高清欧美精品videossex| 亚洲精品一区蜜桃| 亚洲欧美精品自产自拍| 在线观看免费日韩欧美大片| 人人妻人人添人人爽欧美一区卜| 三级国产精品片| 欧美黄色片欧美黄色片| 两性夫妻黄色片| 美女福利国产在线| 国产亚洲av片在线观看秒播厂| 免费女性裸体啪啪无遮挡网站| 日韩三级伦理在线观看| 五月天丁香电影| 精品人妻在线不人妻| 制服诱惑二区| 亚洲 欧美一区二区三区| 晚上一个人看的免费电影| 国产白丝娇喘喷水9色精品| 美女主播在线视频| 国产成人免费无遮挡视频| 成年av动漫网址| 国产白丝娇喘喷水9色精品| 精品国产乱码久久久久久小说| 日韩欧美精品免费久久| 久久国产精品男人的天堂亚洲| 午夜福利网站1000一区二区三区| 亚洲内射少妇av| 日韩 亚洲 欧美在线| 人妻 亚洲 视频| av在线观看视频网站免费| 黄色毛片三级朝国网站| 晚上一个人看的免费电影| 天堂中文最新版在线下载| 国产亚洲午夜精品一区二区久久| 激情视频va一区二区三区| 久久久久国产网址| 汤姆久久久久久久影院中文字幕| 久久久国产一区二区| 日韩 亚洲 欧美在线| 国产欧美亚洲国产| 十分钟在线观看高清视频www| 亚洲精品日本国产第一区| a级片在线免费高清观看视频| 在线观看免费日韩欧美大片| 国产成人免费无遮挡视频| 国产免费又黄又爽又色| 亚洲激情五月婷婷啪啪| 亚洲 欧美一区二区三区| 国产乱人偷精品视频| 久久久国产精品麻豆| 搡女人真爽免费视频火全软件| 老司机亚洲免费影院| 男人舔女人的私密视频| 欧美日韩综合久久久久久| 精品国产一区二区久久| 国产免费一区二区三区四区乱码| 日日啪夜夜爽| 叶爱在线成人免费视频播放| 2018国产大陆天天弄谢| 久久久久国产精品人妻一区二区| 精品人妻熟女毛片av久久网站| 欧美精品一区二区免费开放| 久久午夜综合久久蜜桃| 纯流量卡能插随身wifi吗| 九草在线视频观看| 十分钟在线观看高清视频www| 曰老女人黄片| 欧美老熟妇乱子伦牲交| 五月伊人婷婷丁香| 欧美 亚洲 国产 日韩一| 久久婷婷青草| av在线观看视频网站免费| 日日摸夜夜添夜夜爱| 国产又爽黄色视频| 两性夫妻黄色片| 国产免费一区二区三区四区乱码| 欧美av亚洲av综合av国产av | 亚洲,欧美,日韩| 亚洲伊人久久精品综合| 欧美国产精品va在线观看不卡| 国产午夜精品一二区理论片| 免费在线观看黄色视频的| 色播在线永久视频| 国产成人一区二区在线| 黄片无遮挡物在线观看| 日日爽夜夜爽网站| 男人添女人高潮全过程视频| 2022亚洲国产成人精品| 91精品伊人久久大香线蕉| 国产xxxxx性猛交| 国产精品无大码| 高清欧美精品videossex| 国产又色又爽无遮挡免| 日韩伦理黄色片| 777久久人妻少妇嫩草av网站| 欧美人与善性xxx| 日本vs欧美在线观看视频| 男人爽女人下面视频在线观看| av视频免费观看在线观看| 欧美日韩精品成人综合77777| 黑人欧美特级aaaaaa片| 国产爽快片一区二区三区| 交换朋友夫妻互换小说| xxxhd国产人妻xxx| 久久久a久久爽久久v久久| 只有这里有精品99| 精品国产国语对白av| 欧美变态另类bdsm刘玥| 国产精品秋霞免费鲁丝片| 狂野欧美激情性bbbbbb| 成年人免费黄色播放视频| 成人国语在线视频| 亚洲欧美日韩另类电影网站| 国产白丝娇喘喷水9色精品| 久久久精品94久久精品| 精品亚洲成a人片在线观看| 丝袜人妻中文字幕| 美国免费a级毛片| 亚洲精品一二三| 成人国语在线视频| 999精品在线视频| 久久99一区二区三区| 国产综合精华液| 免费少妇av软件| 国产精品成人在线| 99久国产av精品国产电影| 免费大片黄手机在线观看| 久久久久国产一级毛片高清牌| 欧美成人精品欧美一级黄| 纵有疾风起免费观看全集完整版| 亚洲精品久久午夜乱码| 亚洲国产看品久久| 少妇人妻 视频| 日韩中文字幕欧美一区二区 | 天堂俺去俺来也www色官网| 最近中文字幕2019免费版| 黄色怎么调成土黄色| 久久国产精品男人的天堂亚洲| 成人国产麻豆网| 黄网站色视频无遮挡免费观看| 亚洲人成77777在线视频| 我的亚洲天堂| h视频一区二区三区| av片东京热男人的天堂| 国产无遮挡羞羞视频在线观看| 欧美精品一区二区大全| 免费观看无遮挡的男女| 久热这里只有精品99| 国产精品亚洲av一区麻豆 | 亚洲图色成人| 91aial.com中文字幕在线观看| 中文乱码字字幕精品一区二区三区| 夫妻性生交免费视频一级片| 综合色丁香网| 一区二区日韩欧美中文字幕| 亚洲av成人精品一二三区| 久久精品久久久久久久性| 欧美精品一区二区大全| 久久精品国产综合久久久| 久久av网站| 欧美精品国产亚洲| 新久久久久国产一级毛片| 日本av免费视频播放| 最近的中文字幕免费完整| 一级毛片黄色毛片免费观看视频| 黄色怎么调成土黄色| 1024香蕉在线观看| 一边亲一边摸免费视频| 日韩中字成人| 国产黄频视频在线观看| 免费日韩欧美在线观看| 性高湖久久久久久久久免费观看| 在线观看国产h片| 久久这里有精品视频免费| 各种免费的搞黄视频| 18禁动态无遮挡网站| 91国产中文字幕| 最近的中文字幕免费完整| 久久精品国产亚洲av高清一级| 日日摸夜夜添夜夜爱| 黄片小视频在线播放| 高清av免费在线| 国产精品麻豆人妻色哟哟久久| 国产麻豆69| 中文字幕制服av| 国产精品成人在线| 午夜激情av网站| 成人毛片60女人毛片免费| 国语对白做爰xxxⅹ性视频网站| 性高湖久久久久久久久免费观看| 精品第一国产精品| 精品一品国产午夜福利视频| 十分钟在线观看高清视频www| 乱人伦中国视频| 久久精品国产自在天天线| 久久久精品区二区三区| 男女边摸边吃奶| 丝袜脚勾引网站| 亚洲一码二码三码区别大吗| 精品国产国语对白av| 黄色配什么色好看| 精品第一国产精品| 久热这里只有精品99| 亚洲国产欧美在线一区| 丁香六月天网| 精品一区二区三区四区五区乱码 | 国产女主播在线喷水免费视频网站| 最近手机中文字幕大全| 国产成人精品在线电影| 国产综合精华液| 久久亚洲国产成人精品v| √禁漫天堂资源中文www| 极品少妇高潮喷水抽搐| 国产一区二区在线观看av| 在线亚洲精品国产二区图片欧美| 欧美精品一区二区免费开放| 人妻系列 视频| 久久99一区二区三区| 精品少妇内射三级| 性高湖久久久久久久久免费观看| 亚洲男人天堂网一区| √禁漫天堂资源中文www| 巨乳人妻的诱惑在线观看| 国产精品亚洲av一区麻豆 | 伊人久久国产一区二区| 国产极品粉嫩免费观看在线| 蜜桃国产av成人99| 日韩av在线免费看完整版不卡| 亚洲国产精品成人久久小说| 欧美日韩综合久久久久久| 大片免费播放器 马上看| 伦理电影大哥的女人| 免费观看无遮挡的男女| 中文字幕人妻丝袜一区二区 | www.精华液| 国产xxxxx性猛交| 久久久久国产一级毛片高清牌| 丁香六月天网| 色吧在线观看| 午夜影院在线不卡| 一本—道久久a久久精品蜜桃钙片| 高清在线视频一区二区三区| 日韩精品有码人妻一区| 999久久久国产精品视频| 少妇被粗大猛烈的视频| 韩国高清视频一区二区三区| 婷婷色麻豆天堂久久| 国产日韩欧美视频二区| 电影成人av| 亚洲欧美日韩另类电影网站| 一级片免费观看大全| 亚洲欧美一区二区三区国产| 亚洲av欧美aⅴ国产| 男人爽女人下面视频在线观看| 午夜福利在线观看免费完整高清在| 亚洲国产最新在线播放| 亚洲人成77777在线视频| 国产免费一区二区三区四区乱码| 成人亚洲欧美一区二区av| 精品少妇内射三级| 亚洲精品美女久久av网站| 国产1区2区3区精品| 国产欧美亚洲国产| av在线播放精品| 亚洲经典国产精华液单| 91精品伊人久久大香线蕉| 成人免费观看视频高清| 美女福利国产在线| 国产高清不卡午夜福利| 国产在视频线精品| 丝袜美腿诱惑在线| 婷婷色综合www| 欧美xxⅹ黑人| 9热在线视频观看99| 男女边摸边吃奶| 韩国高清视频一区二区三区| 制服诱惑二区| 久热久热在线精品观看| 美女主播在线视频| 老熟女久久久| 中文字幕色久视频| 精品视频人人做人人爽| 亚洲精品国产一区二区精华液| 午夜免费观看性视频| 中文字幕最新亚洲高清| 日本黄色日本黄色录像| 国产成人免费无遮挡视频| 精品国产乱码久久久久久男人| 免费在线观看视频国产中文字幕亚洲 | 天堂中文最新版在线下载| 午夜福利影视在线免费观看| 免费女性裸体啪啪无遮挡网站| 亚洲av免费高清在线观看| 大话2 男鬼变身卡| 国产精品免费视频内射| 国产黄频视频在线观看| 在线亚洲精品国产二区图片欧美| 久久久精品国产亚洲av高清涩受| 在线观看人妻少妇| 中文精品一卡2卡3卡4更新| 久久国产精品男人的天堂亚洲| 多毛熟女@视频| 尾随美女入室| 国产日韩一区二区三区精品不卡| 1024香蕉在线观看| 亚洲精品久久午夜乱码| 成人亚洲精品一区在线观看| 亚洲综合精品二区| 久久久亚洲精品成人影院| 亚洲三区欧美一区| 精品国产乱码久久久久久小说| 国产日韩欧美亚洲二区| 亚洲精品国产一区二区精华液| 国产伦理片在线播放av一区| 免费观看无遮挡的男女| 你懂的网址亚洲精品在线观看| 黄色视频在线播放观看不卡| 国产精品成人在线| 韩国精品一区二区三区| 各种免费的搞黄视频| 国产av一区二区精品久久| 校园人妻丝袜中文字幕| 亚洲人成77777在线视频| 亚洲一区中文字幕在线| 丝袜人妻中文字幕| 久久久久久久国产电影| 国产97色在线日韩免费| 另类亚洲欧美激情| 久久韩国三级中文字幕| 精品人妻偷拍中文字幕|