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

    潛艇圍殼線型優(yōu)化抑制脈動壓力與流激噪聲的數(shù)值模擬研究

    2014-06-22 05:08:10呂世金沈泓萃
    船舶力學(xué) 2014年4期
    關(guān)鍵詞:大渦聲級馬蹄

    張 楠,呂世金,沈泓萃,謝 華

    (中國船舶科學(xué)研究中心,江蘇 無錫 214082)

    潛艇圍殼線型優(yōu)化抑制脈動壓力與流激噪聲的數(shù)值模擬研究

    張 楠,呂世金,沈泓萃,謝 華

    (中國船舶科學(xué)研究中心,江蘇 無錫 214082)

    潛艇脈動壓力與流激噪聲的數(shù)值預(yù)報與分析已經(jīng)成為流聲耦合研究領(lǐng)域的重要課題。文章對于潛艇圍殼進行了前緣加裝填角以及三維座艙型圍殼的線型優(yōu)化,并采用大渦模擬與聲學(xué)類比方法對于原圍殼、帶填角圍殼與三維圍殼的渦旋流場(馬蹄渦)和聲學(xué)特征進行了數(shù)值計算評估,分析了圍殼線型優(yōu)化對于壁面脈動壓力和流激噪聲的定量影響與抑制效果,研究表明直立型圍殼加裝填角或者進行三維流線型設(shè)計,可以明顯改善圍殼區(qū)域特別是圍殼與艇體交接部的流動品質(zhì),使得脈動壓力與流激噪聲顯著下降。文中工作有益于潛艇流聲耦合學(xué)術(shù)研究以及未來新型潛艇設(shè)計。

    圍殼;線型優(yōu)化;馬蹄渦;脈動壓力;流激噪聲;大渦模擬;聲學(xué)類比

    1 引 言

    當(dāng)火箭、飛機、導(dǎo)彈等飛行器在空氣介質(zhì)中運動,或艦船、潛艇、魚雷等航行體在水介質(zhì)中運動時,它們的外表面都會形成邊界層。邊界層內(nèi)的脈動速度誘導(dǎo)產(chǎn)生物面上的脈動壓力(pressure fluctuations或者wall pressure fluctuations),進而產(chǎn)生遠(yuǎn)場流激輻射噪聲。邊界層中的脈動壓力特性是湍流研究中的基礎(chǔ)性問題。由于脈動壓力是湍流非定常特性的重要表征,而且是流激噪聲的重要來源,所以在流體誘發(fā)振動與噪聲的許多工程應(yīng)用問題中脈動壓力都備受關(guān)注[1-5]。一般認(rèn)為,在水下噪聲問題中,邊界層脈動壓力的重要性主要表現(xiàn)在兩方面[4-5]:第一、它是聲納設(shè)備自噪聲的主要來源。第二、它激勵載體的壁面和其他構(gòu)件振動,一方面向周圍介質(zhì)輻射噪聲,另一方面也在載體內(nèi)部產(chǎn)生艙內(nèi)噪聲。

    由于以上原因,研究湍流邊界層脈動壓力的特性及其抑制方法具有十分重要的實際意義。在早期,大量有關(guān)湍流脈動壓力的研究都是在平板邊界層或管路邊界層中開展的,試驗測量是一個主要的研究手段,當(dāng)時要做到精確的測量也比較困難,主要在于脈動壓力對于傳感器的形式和尺寸較敏感,而且脈動壓力具有非常廣泛的時間與空間尺度,目前,由于壓力傳感器的更新?lián)Q代,試驗研究又有了廣泛且深入的進展。

    利用數(shù)值模擬手段計算壁面脈動壓力與流激噪聲是當(dāng)前該領(lǐng)域的研究熱點,對于簡單幾何形式的物體,國外已有應(yīng)用直接數(shù)值模擬的研究,但限于計算條件,多年來進展并不明顯,對于復(fù)雜幾何形式的物體,還是趨向于大渦模擬(LES)方法[6-14]。LES方法本身的特點就在于可以求解非定常的精細(xì)湍流特征,這種非定常的內(nèi)涵是LES與RANS最大的區(qū)別。2007年劍橋大學(xué)出版社出版了一部由三十位學(xué)者合作完成的專著—《聲學(xué)中的大渦模擬》[15],在這本書中,作者們明確指出,用大渦模擬方法計算非定常流場,在此基礎(chǔ)上,進而計算壁面脈動壓力與流激輻射噪聲,是大渦模擬方法真正的用武之地。目前,大渦模擬在聲學(xué)計算中已經(jīng)嶄露頭角,將來更要發(fā)揮舉足輕重的作用。而要建立有工程實用價值的大渦模擬方法,就必須需要可靠的試驗數(shù)據(jù)來做驗證,近年來,中國船舶科學(xué)研究中心在拖曳水池、低速風(fēng)洞、循環(huán)水槽、低噪聲直流小風(fēng)洞和小型多功能高速空泡水筒中進行了大量湍流脈動壓力的模型試驗,測試模型涵蓋平板、翼型、回轉(zhuǎn)體、小突體等基礎(chǔ)模型以及潛艇、艦船、流水孔、聲吶罩等實際模型,同時在循環(huán)水槽中也進行了大量流激噪聲測量,為深入分析湍流脈動壓力與流激噪聲特性與規(guī)律,以及驗證CFD數(shù)值模擬方法在計算流動非定常特征方面的能力提供了堅實的技術(shù)支撐。

    潛艇指揮臺圍殼無疑是潛艇脈動壓力和流激噪聲的最主要來源之一。潛艇指揮臺圍殼是一流線型的導(dǎo)流罩,其對那些需要經(jīng)常接觸水面的作戰(zhàn)系統(tǒng)起到保護作用,這些被圍殼包圍在內(nèi)的作戰(zhàn)系統(tǒng)包括:潛望鏡、通信天線、通氣管設(shè)備、聲吶陣、艦橋出入艙口、導(dǎo)航設(shè)備以及其他裝置等。自從上世紀(jì)四十年代開始,潛艇圍殼形狀與尺寸就已基本固定,一直沿用至今,并無根本性改變,其設(shè)計思路為:先選擇一個低阻力的母型,將所需作戰(zhàn)系統(tǒng)全部包圍在內(nèi),然后不斷改型,盡量達(dá)到最小體積最小阻力的要求,這種設(shè)計約束直接導(dǎo)致翼型圍殼的誕生,并首先在美國潛艇上應(yīng)用,被其他各國沿襲至今。

    圍殼外型與尺寸對潛艇阻力、流場、脈動壓力與流激噪聲有重要影響。多年來,為了進一步提高潛艇的水動力性能與噪聲性能,人們一直致力于改變圍殼外型與尺寸,甚至嘗試過移除圍殼的探索研究。在早期,圍殼外型優(yōu)化焦點集中在如何降低圍殼尾流對潛艇槳盤面入流的影響,改善槳盤面入流品質(zhì),進而降低螺旋槳低頻線譜噪聲。近年來的研究熱點集中在如何抑制或消除圍殼與主艇體接合部馬蹄渦及其所誘發(fā)的流激噪聲。美國海浪級與弗吉尼亞級潛艇都在圍殼前緣加裝了填角,流場品質(zhì)得到很大改善,研究表明在圍殼前緣與主艇體結(jié)合部加裝不同形式的填角對抑制圍殼結(jié)合部馬蹄渦都有較好的作用。近年來,在復(fù)合材料和復(fù)雜構(gòu)型技術(shù)的支撐下,國外還興起了研究三維座艙型圍殼的熱潮,旨在降低潛艇阻力和流激噪聲并提高潛艇槳盤面入流均勻性,美國水面武器研究中心的Carderock分部(NSWCCD)致力于潛艇水動力、結(jié)構(gòu)與噪聲的綜合研究,專門設(shè)立了先進圍殼研發(fā)計劃(Advanced Sail Project)[16],從流體動力、結(jié)構(gòu)強度、水動力噪聲、設(shè)計技術(shù)、復(fù)合材料等全方位對三維圍殼進行深入研究。

    為了加快研究進度并減小設(shè)計成本,1997年美國先進潛艇研發(fā)部(Advanced Submarine R&D Office)專門成立了CFD計算團隊來用數(shù)值模擬手段承擔(dān)三維圍殼設(shè)計研制中的大量評估問題,支撐項目代號NAVSEA(SEA93R)。迄今為止,這一團隊發(fā)揮了巨大的作用,節(jié)省了大量物理模型試驗,并給出了可靠的精細(xì)流場與流噪聲信息[16]。數(shù)值模擬研究主要分兩步走,第一步,對于弗吉尼亞級潛艇模型進行系統(tǒng)的阻力、流場和噪聲試驗,為驗證CFD方法提供詳細(xì)的模型試驗數(shù)據(jù),利用CFD工具對于上述物理量進行全面模擬,并與試驗結(jié)果對比,不斷改進計算方法,提高計算精度,直到達(dá)到設(shè)計要求的可接受的計算精度為止;第二步,利用已經(jīng)驗證的CFD方法來指導(dǎo)三維圍殼設(shè)計,從1998年至今已經(jīng)進行了大量計算,在此基礎(chǔ)上,產(chǎn)生了大尺度模型(LSV1),在實艇制造之前進行相關(guān)各項性能的試驗驗證。LSV1是自主式1/4尺度的海浪級(SSN21)潛艇模型,采用三維圍殼,NSWCCD聲學(xué)研究分部(Acoustic Research Detachment)對此模型進行了大量試驗,并利用大渦模擬方法(LES)進行大量計算評估,試驗與計算研究表明三維圍殼可以顯著降低脈動壓力和流激噪聲。可以講,目前潛艇線型設(shè)計優(yōu)化重點已經(jīng)由降低阻力和改善槳盤面流場逐步走向降低脈動壓力和流激噪聲,數(shù)值優(yōu)化方法也由RANS發(fā)展為LES。

    脈動壓力與流激噪聲屬于流體力學(xué)與聲學(xué)的學(xué)科交叉研究領(lǐng)域,是流聲耦合研究的關(guān)鍵所在,近年來,課題組進行了大量系列數(shù)值計算[10-14],目前已基本建立了脈動壓力與流激噪聲的數(shù)值預(yù)報方法。在潛艇線型優(yōu)化與流動控制需求的推動下,本文對于某假想潛艇圍殼進行了前緣加裝填角以及三維座艙型圍殼的線型優(yōu)化,并采用大渦模擬方法進行了數(shù)值計算評估,分析了圍殼線型優(yōu)化對于馬蹄渦、脈動壓力與流激噪聲的定量影響與抑制效果,為新型潛艇的聲學(xué)設(shè)計與水動力設(shè)計提供技術(shù)支撐。

    2 數(shù)值模擬方法

    2.1 大渦模擬方程

    濾波的連續(xù)性方程和NS方程可以表示為:

    其中:σij為分子粘性引起的應(yīng)力張量,τij為亞格子應(yīng)力,需用亞格子渦模型進行模擬:

    其中:第一項稱作Leonard項,也稱作外散射項(outscatter term),代表兩個大渦間的相互作用,以產(chǎn)生小尺度湍流;第二項稱作交叉項(cross term),代表大、小渦間的相互作用,其間能量可以從大渦向小渦傳遞,也可以反向傳遞,從小渦傳向大渦,但總體平均起來,能量還是以從大渦向小渦傳遞為主;第三項稱作反散射項或逆散射項(backscatter term),代表小渦間的相互作用以產(chǎn)生大渦,并帶來能量從小渦到大渦的傳遞。過去曾認(rèn)為,既然各項的物理意義不同,應(yīng)該分別用不同模型去近似它們,但是由于?;夹g(shù)尚未完善,分別?;幢販?zhǔn)確,所以意義不大,還是傾向于合在一起作總體上的?;痆17]。

    本文采用動態(tài)Smagorinsky模型對亞格子應(yīng)力進行模擬。這種模型是Germano在1991年提出的[18],通過動態(tài)計算渦粘性系數(shù)來盡可能地反映實際流動情況。通過對最小可解尺度的信息進行采樣,然后利用這些信息來模擬亞格子尺度應(yīng)力。此模型在接近壁面邊界時給出了正確的漸近特性,而且并不需要阻尼函數(shù)或者間歇函數(shù)。此模型還能夠考慮逆散射的影響。Lilly(1992)[19]利用最小二乘法對此模型進行了改進。詳細(xì)的公式推導(dǎo)與含義請見文獻[11]。

    2.2 聲學(xué)類比方程

    Ffowcs Williams與Hawkings聲學(xué)類比方程,即FW-H方程表達(dá)如下:

    其中:ui為xi方向的流體速度分量;un為垂直于物體表面(f=0)的流體速度分量;vi為xi方向的物體表面速度分量;vn為物體表面法向速度分量;δ(f)為Dirac delta函數(shù);H(f)為Heaviside階躍函數(shù)。p′為遠(yuǎn)場聲壓。f=0表示物體表面,f>0表示外部無界的自由空間。ni為物面外法線,指向流體內(nèi)部。c0為遠(yuǎn)場聲速,ρ0為遠(yuǎn)場密度,Tij為 Lighthill應(yīng)力張量。

    利用自由空間格林函數(shù)δ(g)/4πr,可得遠(yuǎn)場解,本文實際使用的是結(jié)合Kirchhoff積分與可滲透表面的滲流FW-H方法。遠(yuǎn)場聲中包含單極子噪聲pT′(x,t),偶極子噪聲pL′(x,t)以及四極子噪聲pQ′(x,t)。

    2.3 流激噪聲計算特點

    Crighton(1988)[22]曾經(jīng)指出,相對于一般的計算流體力學(xué)問題而言,計算流動產(chǎn)生的聲要面臨更多的挑戰(zhàn)。近年來,這一問題又被Colonius與Lele(2004)[23]加以詳細(xì)地討論。

    首先,產(chǎn)生噪聲的流動一定是非定常的,因此不能使用定常計算方法,而非定常的RANS方法又難以達(dá)到足夠的計算要求,因而對于非定常流動發(fā)聲問題的模擬而言,應(yīng)采用諸如DNS、LES、DES這樣的現(xiàn)代湍流模擬技術(shù),但它們的計算花費也是巨大的,一般都要在高性能并行計算機上進行。LES和DES都包含不同尺度流動的模擬與近似,這種處理對于流動噪聲預(yù)報的影響,是一個具有相當(dāng)理論深度的難題,迄今為止,沒人能給出透徹的檢驗與解釋。

    其次,流體力學(xué)與聲學(xué)的擾動幅值和特征尺度存在巨大的差異,除了包含激波的高速流動等少數(shù)情況之外,只有很小部分的流體能可以在遠(yuǎn)場轉(zhuǎn)化為聲能。在直接方法(CAA)中,聲與流動是同時被求解的,這就對數(shù)值差分格式有了非常嚴(yán)格的要求。在混合方法(CFD+聲學(xué)類比)中,流動與聲是分別計算的,數(shù)值精度也要達(dá)到合理的水平,而且所用聲源的表達(dá)式必須忠實于其實際的輻射特性。

    流與聲的尺度分離是流動聲學(xué)中的突出問題,它既帶來了計算上的挑戰(zhàn),同時也帶來了計算上的簡化,關(guān)鍵取決于研究者所選擇的方法。在無限區(qū)域,聲波與流體動力源的時間尺度是相匹配的,其波長λ與源(渦)尺度l之間通過馬赫數(shù)M而近似表達(dá)為λ=l/M。在水中低馬赫數(shù)流動中,M<<1,這就產(chǎn)生了l與λ之間的巨大差異,使得直接計算聲變得十分困難。另一方面,尺度上的巨大差異卻使混合方法成為合理而簡便的實用方法,即假設(shè)聲能很小,難以對流體運動產(chǎn)生有效的影響,所以可以先解算出流動信息,然后再代入聲學(xué)方程解算出聲。

    2.4 邊界條件

    速度入口:潛艇首部向前2L,設(shè)定來流速度的大小與方向。

    壓力出口:潛艇尾部向后4L,設(shè)定相對于參考壓力點的流體靜壓值。

    壁 面:潛艇外表面,設(shè)定無滑移粘附條件。

    外 場:距潛艇表面約2L,速度為未受擾動的主流區(qū)速度。

    其中:L指潛艇長度。聲學(xué)邊界處理方法取自文獻[20]。

    2.5 數(shù)值求解

    時間項采用二階隱式格式離散,動量方程采用限界中心差分格式離散,壓力速度耦合采用SIMPLE算法。利用代數(shù)多重網(wǎng)格方法加速收斂。計算中時間步長Δt=1×10-5s。壁面y+≈1。采用FFT結(jié)合Hanning窗處理非定常信號時間序列。

    Hanning窗表達(dá)式為:

    本文所用FFT算法簡述如下:

    設(shè)x(n)為N點序列,N為2的整數(shù)次方(4,8,16,……,1 024,……),將序列按偶數(shù)點與奇數(shù)點分解為兩個N/2點的序列:

    這樣分組后,計算過程變?yōu)橄扔嬎銉蓚€N/2長度的序列DFT,然后求其線性組合。后兩個公式表示了一種蝶形運算的流程。分組后的計算量大約是N2/2次復(fù)乘,比原來減少一半。分組后的兩部分,還可再各自分為兩組進行類似的簡化,最后分解到每組只有兩個點為止。計算量大為減少,計算效率大為提高。上面介紹的FFT算法是對序列分解到只有兩個點,逐級采用蝶形運算,故又稱為以2為基的FFT算法。為了提高計算速度和減少所需存儲容量,人們還研究了其他基以及混合基的各種算法,其中被公認(rèn)為有效的快速算法是由Winograd提出來的,該算法是將點數(shù)N分解為互為質(zhì)數(shù)的因子的乘積,該算法效率高、效果好,本文就采用這種算法。

    3 計算方法驗證

    驗證算法所用計算模型為用于基礎(chǔ)性研究的圍殼,長約2 m。對于水速V=6 m/s的脈動壓力與流激噪聲進行了計算,基于圍殼長度雷諾數(shù)為Re=1.2×107,并將計算結(jié)果與試驗結(jié)果對比,驗證了計算方法有效性。計算網(wǎng)格與脈動壓力測點位置見圖1。流激噪聲已折算為等效聲源級,全計算域都采用結(jié)構(gòu)化網(wǎng)格離散,網(wǎng)格總數(shù)為9.6×106。水中參考聲壓pref_water=1×10-6Pa。本文所有算例都是在無錫超級計算中心采用并行計算完成的,總的計算時間為2個月。

    圖1 模型表面網(wǎng)格與脈動壓力傳感器位置Fig.1 Mesh on model surface and position of sensor

    圖2 w1點脈動壓力計算與試驗對比Fig.2 Comparison of computed pressure fluctuations with measurement(w1)

    圖3 w2點脈動壓力計算與試驗對比Fig.3 Comparison of computed pressure fluctuations with measurement(w2)

    圖4 w3點脈動壓力計算與試驗對比Fig.4 Comparison of computed pressure fluctuations with measurement(w3)

    圖5 w4點脈動壓力計算與試驗對比Fig.5 Comparison of computed pressure fluctuations with measurement(w4)

    圖6 流激噪聲計算與試驗對比Fig.6 Comparison of computed flow induced noise with measurement

    圍殼表面四個測點位置的湍流脈動壓力計算結(jié)果與試驗結(jié)果的對比見圖2~5。首先,從譜型來看,計算得到的脈動壓力譜型與試驗吻合,隨頻率增加,脈動壓力衰減明顯;其次,從幅值來看,在測試頻率范圍內(nèi),計算誤差在2~6 dB之內(nèi),高頻計算誤差比低頻計算誤差大;最后,計算分析總聲級可知,圍殼脈動壓力總聲級計算誤差在1.5~3.2 dB之內(nèi)。綜上所述,圍殼脈動壓力計算結(jié)果是可靠的。

    流經(jīng)圍殼的三維流動會產(chǎn)生寬帶噪聲,渦旋的形成與演化以及其與壁面的相互作用是流動輻射噪聲產(chǎn)生的根源,請見文獻[14]。圖6給出了圍殼流激噪聲計算與試驗對比。首先,從譜型來看,計算得到的流激噪聲譜型與試驗吻合,隨頻率增加,噪聲降低,8 kHz以上,譜型趨于穩(wěn)定無明顯變化;其次,從幅值來看,在測試頻率范圍內(nèi),計算誤差在3~6 dB之內(nèi),高頻計算誤差比低頻計算誤差稍大;最后,計算分析總聲級可知,圍殼流激噪聲總聲級計算誤差為3.1 dB。綜上所述,圍殼流激噪聲計算結(jié)果與試驗結(jié)果也是吻合較好的,這也證明數(shù)值計算方法是可靠的。

    4 圍殼線型優(yōu)化計算結(jié)果分析

    本文針對某假想潛艇的直立型圍殼進行線型優(yōu)化,首先設(shè)計了前緣帶填角圍殼與三維座艙型圍殼,以期用這些措施來消弱結(jié)合部馬蹄渦并抑制脈動壓力,然后采用已經(jīng)過驗證的大渦模擬方法計算三種類型圍殼的渦流場,主要研究馬蹄渦的變化情況,同時分析了直立型圍殼、前緣帶填角圍殼與三維座艙型圍殼對于壁面脈動壓力和流激噪聲的定量影響,評估了線型優(yōu)化效果。各圍殼模型主尺度與網(wǎng)格數(shù)見表1。計算模型表面網(wǎng)格劃分與8個脈動壓力采樣點位置見圖7,三個圍殼上采樣點相對于主體的位置相同,其上的脈動壓力變化情況具有可比性。模型主體長L≈5 m,取水速V=6 m/s,模型長度雷諾數(shù)為 Re≈3.0×107。

    圖7 模型表面網(wǎng)格與脈動壓力傳感器位置Fig.7 Mesh on model surface and position of sensor

    表1 圍殼模型主尺度與網(wǎng)格數(shù)量Tab.1 Characteristic of fairwater model and number of cell

    4.1 馬蹄渦計算結(jié)果分析

    三種圍殼與主體結(jié)合部馬蹄渦計算結(jié)果見圖8,圖中渦面是利用Q渦判據(jù)表達(dá)的等渦面(Q=3 000),Q取其他數(shù)值不改變?nèi)N圍殼馬蹄渦之間的相互關(guān)系。

    圖8 三種圍殼與主體結(jié)合部馬蹄渦計算結(jié)果Fig.8 Computed horse-shoe vortex around three fairwater models

    從圖8可知,原方案直立型圍殼馬蹄渦呈現(xiàn)明顯的“U”型結(jié)構(gòu),在向圍殼下游發(fā)展的過程中渦面保持穩(wěn)定且耗散較慢,能夠維持到圍殼尾流之中,馬蹄渦的渦腿是“粗壯”的,而且通過計算發(fā)現(xiàn)其表面上的脈動壓力幅值較大;圍殼加裝填角之后,填角起到了很好的導(dǎo)流作用,與原方案相比,馬蹄渦明顯得到弱化,渦腿比原方案縮減很多,且馬蹄渦在向圍殼下游發(fā)展的過程中渦面耗散較快,文獻[21]曾指出填角能夠?qū)⒔诿娓邷u量流動向外排擠和抬升,從而起到抑制渦旋作用,其表面上的脈動壓力幅值比原圍殼有所減小;三維座艙型圍殼對馬蹄渦的抑制作用最為明顯,與原方案相比,馬蹄渦顯著得到弱化,渦腿比原方案縮減很多,馬蹄渦的渦腿是最為“纖細(xì)”的,且馬蹄渦在向圍殼下游發(fā)展的過程中渦面耗散最快,其表面上的脈動壓力幅值最小。綜上所述,直立型圍殼加裝填角或者進行三維流線型設(shè)計,可以明顯改善圍殼區(qū)域流動品質(zhì),使得馬蹄渦得到有效抑制。

    4.2 脈動壓力計算結(jié)果分析

    三種圍殼表面8個采樣點位置上的脈動壓力三分之一倍頻程計算結(jié)果見圖9。由圖可知,與原方案圍殼相比,帶填角圍殼與三維座艙型圍殼的脈動壓力頻譜幅值都有不同程度的下降,而且首部P1、P5兩點下降最多,頻譜低頻段比高頻段減幅更明顯。圍殼加裝填角之后,P1點脈動壓力降低1~26.7 dB,P2點降低1~11.5 dB,P3點脈動壓力降低2~10.7 dB,P4點脈動壓力降低1~8.8 dB,P5點脈動壓力降低1~31.5 dB,P6點脈動壓力降低1~16.6 dB,P7點脈動壓力降低2~24.4 dB,P8點脈動壓力降低1~3.8 dB;三維座艙型圍殼與原圍殼相比,P1點脈動壓力降低1~25.6 dB,P2點降低1~17.5 dB,P3點脈動壓力降低1~21.6 dB,P4點脈動壓力降低2~22.5 dB,P5點脈動壓力降低2~31.8 dB,P6點脈動壓力降低 1~14.8 dB,P7點脈動壓力降低 5~25.9 dB,P8點脈動壓力降低 1~27.8 dB。

    經(jīng)計算12.5 Hz~10 kHz頻段內(nèi)三種圍殼脈動壓力總聲級可知,與原方案圍殼相比,圍殼加裝填角之后,脈動壓力總聲級下降2~14 dB,三維圍殼使得脈動壓力總聲級下降2~21 dB。這說明,圍殼的優(yōu)化設(shè)計改善了流動特征,弱化了馬蹄渦,從而收到了抑制脈動壓力的良好效果。

    圖9 不同圍殼表面脈動壓力計算結(jié)果對比Fig.9 Comparison of computed pressure fluctuations of three fairwater models

    4.3 流激噪聲計算結(jié)果分析

    在流激噪聲計算中,聲壓接收器取在遠(yuǎn)場,然后將計算結(jié)果折算為距離等效聲中心1 m處的聲源級,再進行不同圍殼之間的相互比較。計算得到的三分之一倍頻程聲源級見圖10。由圖可知,與原方案圍殼相比,帶填角圍殼與三維座艙型圍殼的流激噪聲頻譜幅值都有不同程度的下降,在100 Hz~1 kHz頻段,帶填角圍殼下降2~7 dB,三維座艙型圍殼下降 3~9 dB;在 1 kHz~10 kHz頻段,帶填角圍殼下降1~8 dB,三維座艙型圍殼下降2~5 dB。經(jīng)計算100 Hz~10 kHz頻段內(nèi)三種圍殼流激噪聲總聲級可知,與原方案圍殼相比,圍殼加裝填角之后,流激噪聲總聲級下降4 dB,三維圍殼使得流激噪聲總聲級下降6 dB。計算方法在前面章節(jié)得到驗證,三個算例的計算精度基本相同,降噪效果為三個圍殼之間相互比較的結(jié)果,應(yīng)是可信的。

    綜上所述,直立型圍殼加裝填角或者進行三維流線型設(shè)計,可以明顯改善圍殼區(qū)域流動品質(zhì),使得脈動壓力總聲級與頻譜幅值以及流激噪聲總聲級與頻譜幅值都顯著下降。

    圖10 不同圍殼流激噪聲計算結(jié)果對比Fig.10 Comparison of computed flow induced noise of three fairwater models

    5 結(jié) 語

    利用數(shù)值模擬方法來對潛艇進行優(yōu)化設(shè)計已是近年來的研究熱點。在潛艇線型優(yōu)化與流動控制需求的推動下,作者從抑制馬蹄渦并降低脈動壓力與流激噪聲的優(yōu)化目的出發(fā),對于潛艇圍殼進行了前緣加裝填角以及三維圍殼的線型優(yōu)化,并采用大渦模擬方法對于原圍殼、帶填角圍殼與三維座艙型圍殼的渦旋流場與壁面脈動壓力以及流激噪聲進行了數(shù)值計算評估,分析了圍殼線型優(yōu)化對于脈動壓力與流激噪聲的定量影響與抑制效果,研究表明直立型圍殼加裝填角或者進行三維流線型設(shè)計,可以明顯改善圍殼區(qū)域流動品質(zhì),使得馬蹄渦得到有效抑制,并使得脈動壓力總聲級與頻譜幅值顯著下降,圍殼加裝填角之后,脈動壓力總聲級下降2~14 dB,三維圍殼使得脈動壓力總聲級下降2~21 dB,圍殼加裝填角之后,流激噪聲總聲級下降4 dB,三維圍殼使得流激噪聲總聲級下降6 dB。本文可為未來新型潛艇的聲學(xué)設(shè)計與水動力設(shè)計提供技術(shù)支撐。另外,脈動壓力與流激噪聲的研究屬于流體力學(xué)與聲學(xué)的交叉領(lǐng)域,開展相應(yīng)的理論分析、模型試驗和數(shù)值計算研究是十分必要的,研究脈動壓力與流動發(fā)聲之間的定量關(guān)系是今后重要發(fā)展方向。

    [1]Blake W K.Mechanics of flow-induced sound and vibration[M].Academic Press.INC,New York,1986.

    [2]Ross D.Mechanics of underwater noise[M].Peninsula Publishing,Los Altos,California,U.S.A.,1987.

    [3]Skelton E A,James J H.Theoretical acoustics of underwater structures[M].Imperial College Press,1997.

    [4]湯渭霖.水下噪聲學(xué)原理[G].上海交通大學(xué)講義,2004.

    [5]朱錫清.船舶流體動力噪聲原理[G].無錫:中國船舶科學(xué)研究中心講義,2006.

    [6]Chisachi Kato,Yoshinobu Yamade et al.Numerical prediction of sound generated from flows with a low Mach number[J].Computers&Fluids,2007,36:53-68.

    [7]Khalighi Y,Mani A,Ham F,Moin P.Prediction of sound generated by complex flows at low Mach number regimes[G].Center for Turbulence Research,Stanford University,Annual Research Briefs 2008,2008:313-324.

    [8]Ask J,Davidson L.A numerical investigation of the flow past a generic side mirror and its impact on sound generation[J].Journal of Fluids Engineering,2009,131:061102-1~061102-12.

    [9]Gand F,Deck S.Flow dynamics past a simplified wing body junction[J].Physics of Fluids.2010,22:115111.

    [10]Zhang Nan,Shen Hongcui,Yao Huizhi.Numerical investigation of the flow induced noise of submarine[C]//8th International Conference on Hydrodynamics(ICHD2008).Nantes,France,2008:591-599.

    [11]張 楠.孔腔流動和流激噪聲機理及耦合計算方法研究[D].無錫:中國船舶科學(xué)研究中心,2010.

    [12]Zhang Nan,Shen Hongcui,Yao Huizhi.Numerical simulation of cavity flow induced noise By LES And FW-H acoustic analogy[C]//9th International Conference on Hydrodynamics(ICHD2010).Shanghai,China,2010:242-247.

    [13]張 楠,沈泓萃,朱錫清,姚惠之,謝 華.三維孔腔流激噪聲的大渦模擬與聲學(xué)類比預(yù)報與驗證研究[J].船舶力學(xué),2010,14(1-2):181-190.

    Zhang Nan,Shen Hongcui,et al.Validation and prediction of flow induced naise of 3-dimensional cavity with large eddy simulation and acoustic analogy[J].Journal of Ship Mechanics,2010,14(1-2):181-190.

    [14]張 楠,張勝利,沈泓萃,謝 華.翼/板結(jié)合部渦旋流動結(jié)構(gòu)與壁面脈動壓力的大渦模擬研究[J].船舶力學(xué),2013,17(7):729-740.Zhang Nan,Zhang Shengli,et al.Large eddy simulation of vortical flow structure and wall pressure fluctuations around wing-plate junction[J].Journal of Ship Mechanics,2013,17(7):729-740.

    [15]Wagner C,Huttl T,Sagaut P.Large eddy simulation for acoustics[M].Cambridge,UK:Cambridge University Press,2007.

    [16]Stout M C,Dozier D F.Advanced Submarine Sail[G].Carderock Division,NSWC-Technical digest,2011.

    [17]劉應(yīng)中,繆國平.高等流體力學(xué)[M].上海:上海交通大學(xué)出版社,2000.

    [18]Germano M,Piomelli U,Moin P,Cabot W H.A dynamic subgrid-scale eddy viscosity model[J].Phys.Fluids A.1991,3(7):1760-1765.

    [19]Lilly D K.A proposed modification of the Germano Subgrid-scale closure method[J].Phys.Fluids A.1992,4(3):633-635.

    [20]Brun C,et al.(Eds.)Numerical simulation of turbulent flows&noise generation[M].Springer-Verlag Berlin Heidelberg.2009.

    [21]張 楠,沈泓萃,姚惠之.潛艇阻力與流場的數(shù)值模擬與驗證及艇型的數(shù)值優(yōu)化研究[J].船舶力學(xué),2005,9(1):1-13.

    Zhang Nan,Shen Hongcui,et al.Validation of numerical simulation on resistance and flow field of submarine and numerical optimization of submarine hull form[J].Journal of Ship Mechanics,2005,9(1):1-13.

    [22]Crighton D G.Goals for computational aeroacoustics[M].In Computational Acoustics:Algorithms and Applications,Amsterdam:Elsevier-North Holland.1988:3-20.

    [23]Colonius T,Lele S K.Computational aeroacoustics:Progress on nonlinear problems of sound generation[J].Prog.Aerospace Sci.,2004,40(6):345-416.

    Numerical simulation on the effect of fairwater optimization to suppress the wall pressure fluctuations and flow induced noise

    ZHANG Nan,Lü Shi-jin,SHEN Hong-cui,XIE Hua

    (China Ship Scientific Research Center,Wuxi 214082,China)

    Numerical simulation of wall pressure fluctuations and flow induced noise of submarine has become vigorous in the field of flow-acoustic coupling.The fairwater of a submarine is optimized by two shapes,one is the fairwater with fillet,and another is three-dimensional canopy-like fairwater.Numerical investigation by LES and acoustic analogy is carried out to analyze the vortical flow fields(horse-shoe vortex)and acoustic characteristics of the three fairwaters,and the effect of optimization is studied to describe the quantitative change of pressure fluctuations and flow induced noise.This study shows that the optimization of fairwater is able to improve the flow quality around fairwater especially in the junction of fairwater and hull,and the wall pressure fluctuations and flow induced noise can be suppressed obviously.The work is beneficial to the research in the field of flow-acoustic coupling and the design of the submarine with a new type.

    fairwater;shape optimization;horse-shoe vortex;pressure fluctuations;flow induced noise;large eddy simulation(LES);acoustic analogy

    U661.3

    A

    10.3969/j.issn.1007-7294.2014.04.013

    1007-7294(2014)04-0448-11

    2014-03-04

    自然科學(xué)基金資助項目(51079133);江蘇省自然科學(xué)基金資助項目(BK2010162)

    張 楠(1977-),男,博士,中國船舶科學(xué)研究中心高級工程師,E-mail:zn_nan@sina.com;

    呂世金(1973-),男,中國船舶科學(xué)研究中心研究員。

    猜你喜歡
    大渦聲級馬蹄
    一種等效40 h聲級的計算方法探討
    馬蹄子與北海道男人的選擇
    吃“馬蹄”
    基于壁面射流的下?lián)舯┝鞣欠€(wěn)態(tài)風(fēng)場大渦模擬
    軸流風(fēng)機葉尖泄漏流動的大渦模擬
    馬蹄島(男中音獨唱)
    草原歌聲(2016年1期)2016-11-15 05:46:54
    火車站列車運行噪聲測試與分析
    河口動物 馬蹄蟹
    基于大渦模擬的旋風(fēng)分離器錐體結(jié)構(gòu)影響研究
    AWA6228型聲級計監(jiān)測室內(nèi)噪聲的問題探討
    男女无遮挡免费网站观看| 午夜免费鲁丝| 边亲边吃奶的免费视频| 日韩电影二区| 热99国产精品久久久久久7| 一级毛片电影观看| 国内少妇人妻偷人精品xxx网站| 一级爰片在线观看| 久久久久精品性色| 日韩中文字幕视频在线看片| av网站免费在线观看视频| 麻豆成人av视频| 欧美成人精品欧美一级黄| 女性被躁到高潮视频| 精品一区在线观看国产| 日本午夜av视频| 成人亚洲欧美一区二区av| 18+在线观看网站| 精品99又大又爽又粗少妇毛片| 插阴视频在线观看视频| 自拍欧美九色日韩亚洲蝌蚪91 | 一本—道久久a久久精品蜜桃钙片| 国产精品久久久久久久电影| 在线观看免费视频网站a站| 蜜桃久久精品国产亚洲av| 妹子高潮喷水视频| 春色校园在线视频观看| 视频中文字幕在线观看| 午夜激情久久久久久久| 久久久久久伊人网av| 婷婷色综合大香蕉| 国内少妇人妻偷人精品xxx网站| 国产男女超爽视频在线观看| 久久久久久久久久成人| 又黄又爽又刺激的免费视频.| 日日啪夜夜爽| 少妇人妻精品综合一区二区| 国产精品一区www在线观看| 国产精品熟女久久久久浪| 天堂8中文在线网| 天堂8中文在线网| 国产欧美另类精品又又久久亚洲欧美| 一区在线观看完整版| 在线观看av片永久免费下载| 国产熟女欧美一区二区| 久久精品国产亚洲av涩爱| 日韩一区二区视频免费看| 在线亚洲精品国产二区图片欧美 | 成人漫画全彩无遮挡| 丰满人妻一区二区三区视频av| 免费高清在线观看视频在线观看| 肉色欧美久久久久久久蜜桃| 人妻制服诱惑在线中文字幕| 99re6热这里在线精品视频| 80岁老熟妇乱子伦牲交| 男人添女人高潮全过程视频| 国产深夜福利视频在线观看| 国产成人aa在线观看| 欧美日韩国产mv在线观看视频| 王馨瑶露胸无遮挡在线观看| av在线播放精品| 日韩成人av中文字幕在线观看| 欧美性感艳星| 久久久久视频综合| 国产男女内射视频| 91午夜精品亚洲一区二区三区| 久久久久久久久久人人人人人人| 午夜免费男女啪啪视频观看| 国产乱人偷精品视频| 欧美3d第一页| 久久久久精品性色| 一本大道久久a久久精品| 成人特级av手机在线观看| www.色视频.com| 亚洲精品日韩av片在线观看| 天天躁夜夜躁狠狠久久av| 日韩熟女老妇一区二区性免费视频| 秋霞伦理黄片| 亚洲精品乱码久久久v下载方式| 男人狂女人下面高潮的视频| av不卡在线播放| 久久韩国三级中文字幕| 一级二级三级毛片免费看| 男女边吃奶边做爰视频| 久久久久久久久大av| 色吧在线观看| 看非洲黑人一级黄片| 精品国产一区二区久久| 搡女人真爽免费视频火全软件| 亚洲欧美成人综合另类久久久| 女人久久www免费人成看片| 久久亚洲国产成人精品v| 人妻 亚洲 视频| 免费黄网站久久成人精品| 国产在视频线精品| av在线播放精品| 日本av手机在线免费观看| 午夜激情久久久久久久| 少妇的逼好多水| 国产精品三级大全| 日本猛色少妇xxxxx猛交久久| 色婷婷av一区二区三区视频| 精华霜和精华液先用哪个| 亚洲精品视频女| 制服丝袜香蕉在线| 天堂中文最新版在线下载| 岛国毛片在线播放| 蜜桃在线观看..| 国产免费又黄又爽又色| 久久久国产欧美日韩av| 国产精品麻豆人妻色哟哟久久| 久久久a久久爽久久v久久| 久久久久久久精品精品| 欧美成人精品欧美一级黄| 极品教师在线视频| 一级av片app| 亚洲国产日韩一区二区| 妹子高潮喷水视频| 一本一本综合久久| 午夜91福利影院| av国产久精品久网站免费入址| 91久久精品国产一区二区三区| av播播在线观看一区| 久久久久久久久久久久大奶| 国产精品一区二区在线观看99| 国产精品国产av在线观看| 晚上一个人看的免费电影| 97精品久久久久久久久久精品| 下体分泌物呈黄色| 中文字幕精品免费在线观看视频 | 99热这里只有精品一区| 偷拍熟女少妇极品色| 国产亚洲91精品色在线| 日韩精品有码人妻一区| 亚洲av国产av综合av卡| 亚洲激情五月婷婷啪啪| 黑丝袜美女国产一区| 亚洲在久久综合| 国产色婷婷99| av在线app专区| 欧美日韩亚洲高清精品| 嘟嘟电影网在线观看| 午夜激情福利司机影院| 少妇 在线观看| 秋霞在线观看毛片| 国产亚洲一区二区精品| 国内精品宾馆在线| 毛片一级片免费看久久久久| 久久久久人妻精品一区果冻| 久久人人爽人人爽人人片va| 99久久中文字幕三级久久日本| 少妇被粗大猛烈的视频| 精华霜和精华液先用哪个| 久久国产乱子免费精品| 精品国产一区二区久久| 亚洲精品视频女| av在线观看视频网站免费| 国产永久视频网站| 天天操日日干夜夜撸| 一区二区三区精品91| 久久99蜜桃精品久久| 国产熟女午夜一区二区三区 | 免费观看在线日韩| 边亲边吃奶的免费视频| 久久国产乱子免费精品| 十分钟在线观看高清视频www | 在线观看免费高清a一片| 亚洲国产最新在线播放| 99久国产av精品国产电影| 少妇 在线观看| 少妇猛男粗大的猛烈进出视频| 嫩草影院入口| 高清av免费在线| 国产精品女同一区二区软件| 综合色丁香网| 99九九在线精品视频 | 夫妻性生交免费视频一级片| 久久国产精品大桥未久av | 亚洲精品456在线播放app| 国产69精品久久久久777片| 久久精品国产鲁丝片午夜精品| 国产精品麻豆人妻色哟哟久久| 高清视频免费观看一区二区| 亚洲国产日韩一区二区| 黄色日韩在线| av网站免费在线观看视频| 国产精品偷伦视频观看了| 99热国产这里只有精品6| 我的女老师完整版在线观看| 夫妻午夜视频| 人人妻人人添人人爽欧美一区卜| 99久久精品国产国产毛片| 欧美日韩av久久| 免费看av在线观看网站| 在线亚洲精品国产二区图片欧美 | 性色av一级| 中文字幕制服av| 黄色欧美视频在线观看| 亚洲久久久国产精品| 亚洲精华国产精华液的使用体验| 日日摸夜夜添夜夜添av毛片| 久久久久国产网址| 伊人久久精品亚洲午夜| 国产精品.久久久| 亚洲欧美一区二区三区黑人 | 国产午夜精品一二区理论片| 丝袜喷水一区| 在线观看免费视频网站a站| 在线精品无人区一区二区三| 国产精品一区二区三区四区免费观看| 国产精品一区二区在线观看99| 成人毛片60女人毛片免费| 自拍偷自拍亚洲精品老妇| 七月丁香在线播放| 男女边吃奶边做爰视频| 热re99久久国产66热| 中文乱码字字幕精品一区二区三区| 久久99热6这里只有精品| 精品99又大又爽又粗少妇毛片| 亚洲av成人精品一区久久| 色5月婷婷丁香| 欧美人与善性xxx| 高清午夜精品一区二区三区| 久久99精品国语久久久| 又大又黄又爽视频免费| 美女主播在线视频| 五月玫瑰六月丁香| 国产高清国产精品国产三级| 国产黄片视频在线免费观看| 欧美日韩av久久| 欧美成人午夜免费资源| 精品亚洲乱码少妇综合久久| 一级毛片aaaaaa免费看小| 欧美+日韩+精品| 成人免费观看视频高清| 欧美xxⅹ黑人| 国产 一区精品| 精品午夜福利在线看| 国产欧美另类精品又又久久亚洲欧美| 一级二级三级毛片免费看| 99九九线精品视频在线观看视频| 免费高清在线观看视频在线观看| 欧美亚洲 丝袜 人妻 在线| 久久久久精品久久久久真实原创| 一区二区三区精品91| av播播在线观看一区| 国产一区二区三区av在线| 免费高清在线观看视频在线观看| 午夜福利网站1000一区二区三区| 91午夜精品亚洲一区二区三区| 国产极品粉嫩免费观看在线 | 99热国产这里只有精品6| 亚洲三级黄色毛片| 久久毛片免费看一区二区三区| 国产伦理片在线播放av一区| 人体艺术视频欧美日本| 汤姆久久久久久久影院中文字幕| 一级,二级,三级黄色视频| 少妇人妻久久综合中文| 黄色日韩在线| 日本av免费视频播放| 日韩精品有码人妻一区| 中文精品一卡2卡3卡4更新| 亚洲一级一片aⅴ在线观看| 国产极品天堂在线| 哪个播放器可以免费观看大片| 久久久国产一区二区| 五月玫瑰六月丁香| 精品国产一区二区久久| 久久久国产一区二区| 成人18禁高潮啪啪吃奶动态图 | 97在线视频观看| 妹子高潮喷水视频| h日本视频在线播放| av福利片在线| 青青草视频在线视频观看| 国产一区二区三区综合在线观看 | 国产成人精品婷婷| 亚洲精品色激情综合| 久久ye,这里只有精品| 一区二区三区四区激情视频| 丝袜脚勾引网站| 国产日韩欧美亚洲二区| 亚洲,欧美,日韩| 国产黄色免费在线视频| 色吧在线观看| 亚洲,一卡二卡三卡| 51国产日韩欧美| 三级国产精品片| 久久久久久人妻| 日韩电影二区| av.在线天堂| 日韩av不卡免费在线播放| 婷婷色麻豆天堂久久| 国产极品粉嫩免费观看在线 | 水蜜桃什么品种好| 熟妇人妻不卡中文字幕| 美女中出高潮动态图| 久久精品夜色国产| 人妻少妇偷人精品九色| 欧美日韩亚洲高清精品| 亚洲av成人精品一区久久| 国产白丝娇喘喷水9色精品| 欧美精品一区二区大全| 免费观看a级毛片全部| 久久久久人妻精品一区果冻| 水蜜桃什么品种好| 久久国产精品男人的天堂亚洲 | 成人综合一区亚洲| 一个人看视频在线观看www免费| 欧美三级亚洲精品| 一区二区av电影网| 国产精品免费大片| 国精品久久久久久国模美| 亚洲av免费高清在线观看| 三级国产精品片| 美女视频免费永久观看网站| 亚洲性久久影院| 国产精品偷伦视频观看了| 亚洲av福利一区| 国模一区二区三区四区视频| 婷婷色av中文字幕| 国产av精品麻豆| 精品一区二区免费观看| 我要看黄色一级片免费的| tube8黄色片| 欧美最新免费一区二区三区| 成人国产麻豆网| 视频中文字幕在线观看| 另类精品久久| 久久鲁丝午夜福利片| 国产一区二区三区av在线| 一区二区三区四区激情视频| 成年av动漫网址| 最近中文字幕2019免费版| 日韩中字成人| 久久影院123| 久久久久久久大尺度免费视频| 成人免费观看视频高清| 久久久久久久久大av| 99久久精品国产国产毛片| 少妇被粗大的猛进出69影院 | 亚洲色图综合在线观看| 色视频www国产| 国产成人精品福利久久| 视频区图区小说| 春色校园在线视频观看| a级一级毛片免费在线观看| 亚洲一区二区三区欧美精品| 菩萨蛮人人尽说江南好唐韦庄| 99视频精品全部免费 在线| 久久精品久久久久久久性| 亚洲av男天堂| 热re99久久精品国产66热6| 中文字幕久久专区| 日本色播在线视频| 免费观看的影片在线观看| 久久国产精品男人的天堂亚洲 | 亚洲高清免费不卡视频| 日本91视频免费播放| 在线看a的网站| 国产亚洲最大av| 秋霞在线观看毛片| 大香蕉久久网| 观看av在线不卡| 婷婷色综合www| 99久久人妻综合| 成人亚洲欧美一区二区av| 国产男人的电影天堂91| 少妇人妻精品综合一区二区| 免费黄色在线免费观看| 在线观看人妻少妇| 两个人的视频大全免费| 偷拍熟女少妇极品色| 2021少妇久久久久久久久久久| 97精品久久久久久久久久精品| 亚洲无线观看免费| 日产精品乱码卡一卡2卡三| xxx大片免费视频| 精品久久久精品久久久| 高清毛片免费看| 97超碰精品成人国产| 中文欧美无线码| 亚洲欧洲日产国产| 97精品久久久久久久久久精品| 国产精品嫩草影院av在线观看| 18禁动态无遮挡网站| 国产精品人妻久久久久久| 五月伊人婷婷丁香| 久久99热这里只频精品6学生| 建设人人有责人人尽责人人享有的| 亚洲无线观看免费| 天天躁夜夜躁狠狠久久av| 午夜激情久久久久久久| 亚洲欧美中文字幕日韩二区| 精品一区二区三卡| 国产探花极品一区二区| 男女免费视频国产| 草草在线视频免费看| 国产一区二区在线观看日韩| 亚洲情色 制服丝袜| 日韩成人av中文字幕在线观看| 极品教师在线视频| 精品一品国产午夜福利视频| av网站免费在线观看视频| 99热网站在线观看| 夜夜骑夜夜射夜夜干| 成年女人在线观看亚洲视频| 亚洲熟女精品中文字幕| 免费久久久久久久精品成人欧美视频 | 国产成人a∨麻豆精品| 精品亚洲成a人片在线观看| 高清午夜精品一区二区三区| 九色成人免费人妻av| 九九爱精品视频在线观看| 在线观看免费视频网站a站| 男女边摸边吃奶| 国产伦精品一区二区三区视频9| 一级毛片我不卡| 天堂俺去俺来也www色官网| 亚洲精品一二三| 国产黄频视频在线观看| 黄色怎么调成土黄色| 啦啦啦啦在线视频资源| 老司机亚洲免费影院| 日韩精品免费视频一区二区三区 | 超碰97精品在线观看| 国产一区有黄有色的免费视频| 午夜福利在线观看免费完整高清在| 最近最新中文字幕免费大全7| 国产 精品1| 黑丝袜美女国产一区| 色吧在线观看| 女的被弄到高潮叫床怎么办| 国产男女超爽视频在线观看| 91aial.com中文字幕在线观看| 国模一区二区三区四区视频| 日韩精品免费视频一区二区三区 | 午夜福利,免费看| 久久99一区二区三区| 欧美少妇被猛烈插入视频| 久热久热在线精品观看| 亚洲精品日韩av片在线观看| 在线观看免费日韩欧美大片 | 边亲边吃奶的免费视频| 好男人视频免费观看在线| 中国美白少妇内射xxxbb| 2018国产大陆天天弄谢| 中文字幕亚洲精品专区| 最新的欧美精品一区二区| 欧美 亚洲 国产 日韩一| 成人免费观看视频高清| 国产极品天堂在线| 日日爽夜夜爽网站| 尾随美女入室| 国产亚洲欧美精品永久| 亚洲不卡免费看| 欧美精品人与动牲交sv欧美| 中文欧美无线码| 男人爽女人下面视频在线观看| 在线亚洲精品国产二区图片欧美 | 久久国产精品男人的天堂亚洲 | 国产午夜精品久久久久久一区二区三区| 日韩免费高清中文字幕av| 免费看不卡的av| 少妇人妻 视频| 日韩欧美精品免费久久| 亚洲国产精品999| 亚洲美女视频黄频| 日韩亚洲欧美综合| 国产欧美亚洲国产| 狂野欧美激情性xxxx在线观看| 乱人伦中国视频| 街头女战士在线观看网站| 在线观看免费高清a一片| 亚洲国产毛片av蜜桃av| 欧美精品一区二区大全| 六月丁香七月| 少妇的逼水好多| 亚洲综合色惰| 晚上一个人看的免费电影| 国产男女超爽视频在线观看| 国产精品一区二区在线不卡| 亚洲精品一区蜜桃| 汤姆久久久久久久影院中文字幕| 99九九在线精品视频 | 观看美女的网站| 亚洲人成网站在线播| 亚洲国产日韩一区二区| 一级二级三级毛片免费看| 亚洲av在线观看美女高潮| 欧美亚洲 丝袜 人妻 在线| 精品一区二区三区视频在线| 欧美激情极品国产一区二区三区 | 黄色视频在线播放观看不卡| 欧美高清成人免费视频www| 精品久久久久久久久亚洲| 国产精品久久久久久精品电影小说| 亚洲性久久影院| 五月玫瑰六月丁香| 老司机亚洲免费影院| av福利片在线观看| 免费播放大片免费观看视频在线观看| 久久99热6这里只有精品| 国产在线免费精品| 久久午夜综合久久蜜桃| 中文字幕亚洲精品专区| 精品一区二区三卡| 国产极品天堂在线| av国产久精品久网站免费入址| 另类亚洲欧美激情| 国产国拍精品亚洲av在线观看| 久热久热在线精品观看| 三级国产精品欧美在线观看| 精品久久久久久久久av| 内地一区二区视频在线| 边亲边吃奶的免费视频| 午夜免费鲁丝| 晚上一个人看的免费电影| 80岁老熟妇乱子伦牲交| 亚洲精品一二三| 在线免费观看不下载黄p国产| 一本一本综合久久| 卡戴珊不雅视频在线播放| 黄片无遮挡物在线观看| 国产亚洲午夜精品一区二区久久| 国产精品国产三级国产专区5o| 黄色一级大片看看| 亚洲精品久久午夜乱码| 久久精品熟女亚洲av麻豆精品| 亚洲综合色惰| 纯流量卡能插随身wifi吗| 寂寞人妻少妇视频99o| 国产精品女同一区二区软件| 国产日韩欧美在线精品| 97超视频在线观看视频| 国产亚洲一区二区精品| 99视频精品全部免费 在线| 亚洲av福利一区| 国产免费又黄又爽又色| 亚洲欧美精品专区久久| 日韩欧美精品免费久久| 91精品一卡2卡3卡4卡| 国国产精品蜜臀av免费| 熟女人妻精品中文字幕| 亚洲成色77777| 嫩草影院新地址| 人人妻人人爽人人添夜夜欢视频 | 少妇高潮的动态图| 国产欧美另类精品又又久久亚洲欧美| 国产亚洲av片在线观看秒播厂| 成人黄色视频免费在线看| 黄色视频在线播放观看不卡| 高清黄色对白视频在线免费看 | 欧美成人午夜免费资源| 九九久久精品国产亚洲av麻豆| 丝袜脚勾引网站| 免费观看在线日韩| 一区二区三区乱码不卡18| 91精品一卡2卡3卡4卡| 精品午夜福利在线看| 老熟女久久久| 精品国产乱码久久久久久小说| 欧美 亚洲 国产 日韩一| av福利片在线观看| 国产成人a∨麻豆精品| 一级爰片在线观看| 日本wwww免费看| 99热6这里只有精品| 天天躁夜夜躁狠狠久久av| 久久99热这里只频精品6学生| 精品人妻一区二区三区麻豆| 一区二区三区精品91| 最近中文字幕高清免费大全6| 日韩一区二区视频免费看| 亚洲国产欧美在线一区| 97在线人人人人妻| 久久久久久伊人网av| 亚洲av成人精品一二三区| 精品人妻一区二区三区麻豆| 色婷婷av一区二区三区视频| 中文字幕制服av| 亚洲欧美成人综合另类久久久| 亚洲国产精品国产精品| 亚洲一级一片aⅴ在线观看| 欧美一级a爱片免费观看看| 九色成人免费人妻av| 国产精品人妻久久久久久| 精品少妇久久久久久888优播| 大香蕉久久网| 精品久久久精品久久久| 男女免费视频国产| 水蜜桃什么品种好| 国产精品偷伦视频观看了| 久久久久久久大尺度免费视频| 久久久久精品久久久久真实原创| 极品少妇高潮喷水抽搐| xxx大片免费视频| 成年女人在线观看亚洲视频| 搡老乐熟女国产| 国产在线免费精品| 一本久久精品| av网站免费在线观看视频| 人妻少妇偷人精品九色| 亚洲av中文av极速乱| 国产精品一二三区在线看| 免费观看a级毛片全部| 亚洲av在线观看美女高潮| 亚洲精品乱久久久久久| 80岁老熟妇乱子伦牲交| 国产精品女同一区二区软件| 黄色配什么色好看| 能在线免费看毛片的网站| 日本-黄色视频高清免费观看| 黄色视频在线播放观看不卡|