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

淺海中圓柱殼的聲輻射特性分析

2017-06-19 19:00:47繆宇躍李天勻朱翔王鵬張冠軍
哈爾濱工程大學學報 2017年5期

繆宇躍, 李天勻, 朱翔, 王鵬, 張冠軍

(1.中國艦船研究設計中心,湖北 武漢 430064; 2.華中科技大學 船舶與海洋工程學院,湖北 武漢 430074; 3.船舶與海洋水動力湖北省重點實驗室,湖北 武漢 430074; 4.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心, 上海 200240)

?

淺海中圓柱殼的聲輻射特性分析

繆宇躍1, 李天勻2,3,4, 朱翔2,3,4, 王鵬2,3, 張冠軍2,3

(1.中國艦船研究設計中心,湖北 武漢 430064; 2.華中科技大學 船舶與海洋工程學院,湖北 武漢 430074; 3.船舶與海洋水動力湖北省重點實驗室,湖北 武漢 430074; 4.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心, 上海 200240)

根據聲學邊界元理論,本文在考慮海面和吸聲海底的聲反射作用基礎上,對格林函數進行修正,建立波導域聲輻射模型,并分析了淺海中圓柱殼的聲輻射特性。研究發(fā)現吸聲和剛性海底邊界條件下圓柱殼的輻射聲場有明顯差別;可采用波導中點源格林函數預測圓柱殼的聲輻射特性。結果表明邊界反射作用隨著圓柱殼離開海面或海底而減弱,隨著場點離開圓柱殼而增強;且水深方向上聲壓呈周期性波動,其周期規(guī)律與海底邊界類型無關。本文用聲學點源的輻射疊加原理解釋研究現象的產生機理,可為工程應用提供理論參考。

邊界元;淺海;圓柱殼;聲輻射;吸聲海底;波導;格林函數;點源

有限空間中結構的聲輻射特性往往不可避免地受到邊界影響,其中半空間聲輻射問題的研究較為成熟。黎勝等利用半空間聲學邊界元法研究了自由液面和剛性壁面對脈動球聲輻射的影響,結果表明剛性壁面對脈動球輻射聲功率的影響遠小于自由液面對脈動球輻射聲功率的影響[1-2]。鄒元杰等建立了半無限流體域中結構流固耦合振動方程, 探討了自由液面和剛性壁面對結構的固有頻率、振動響應和有關聲學物理量的影響,結果表明邊界影響明顯存在并隨結構離開邊界而減弱[3]。另外還有一些與半空間聲學邊界元法相關的文獻,其重點在于流固耦合算法和半空間快速多極聲學邊界元法的研究[4-9]。陳爐云等對垂直和水平壁面同時存在的四分之一和八分之一空間聲輻射問題展開研究,建立了多虛源的邊界元模型,拓展了半空間聲學邊界元法的應用范圍[10-12]。對類似平行波導的有限空間中結構聲輻射問題的研究并不多見。鄒元杰等建立淺水域聲學邊界元方程和相應的有限元/邊界元流固耦合振動方程,研究平行邊界對結構振動和聲輻射的影響,其水底邊界近似為剛性水底[13]。Chen等利用波疊加法研究了半空間和平行波導中結構的聲輻射問題[14];Zou等采用三維聲彈性理論分析了淺海環(huán)境對結構聲輻射的影響,其中水底邊界的反射系數為掠射角的函數[15-16]。白振國等采用解析法建立了淺水環(huán)境中二維圓柱殼的振動聲輻射數學物理模型,研究了水深和浸深對圓柱殼聲振特性的影響,結果表明圓柱殼離開水面一段小距離后其振動響應和表面聲壓與在無限水域中一致,邊界對輻射聲場的影響十分顯著[17]。在淺海環(huán)境中(水深一般為幾十到幾百米),海底吸聲作用不可忽略,本文首先對波導中格林函數進行修正,提出可考慮邊界吸聲作用的平行波導空間聲學邊界元法,然后研究吸聲和剛性海底邊界條件下淺海波導域中圓柱殼的聲輻射特性及其規(guī)律,最后運用聲學點源的輻射疊加原理解釋輻射規(guī)律的產生機理。

1 海底吸聲邊界條件下的波導聲學邊界元法

(1)

圖1 淺海平行波導空間中結構聲輻射示意圖Fig.1 The sketch of the acoustical radiation of structures in the shallow sea

平行波導空間中三維聲學Helmholtz方程的基本解(格林函數)由實源項和虛源項疊加而成[13],可表示為

(2)

式中:第一項代表實源,第二項代表虛源,R0為P和實源Q的距離,RUn為P和虛源QUn的距離,RLn為P和虛源QLn的距離。FUn和FLn分別為QUn和QLn在海面和海底的反射系數,對于自由液面FUn=-1,對于剛性海底FLn=1。k0表示海水中的波數,m是與反射次數n相關的系數,m=j/2+[1+(-1)j+1]/4。

由于實源和虛源在一條直線上,QUn和QLn與Q的x坐標關系如下

(3)

式中:hU和hL分別為海面和海底到坐標原點o的距離,x0為Q的x坐標,h=hU+hL表示淺海深度。

實際上海底并非剛性,其對聲波的吸收是很明顯的。由于基本解中的三項e-ik0R0/(4πR0)、e-ik0RUn/(4πRUn)和e-ik0RLn/4πRLn分別表示位于Q、QUn和QLn的點源作用于場點P的聲壓,故式(1)表示分布于結構表面不同強度點源及其鏡像在場點P的聲壓疊加[14-16]。在吸聲邊界條件下,點源格林函數中海底反射系數可表示為[19]

FLn=F+B

(4)

其中,

海水密度和聲速為ρ0=1 025kg/m3和c0=1 500m/s,海底沉積層密度ρ1=1 400kg/m3,考慮沉積層對聲波吸收的聲速[20]c1=1 530(1-iη),其中η=0.021。

在結構表面劃分邊界單元,通過插值求得x0并代入式(3)得到xLn,則掠射角θn為

(5)

式中 xp是P的x坐標。

由式(4)和(5)求得反射系數FLn,再聯合式(1)~(3)便可求得吸聲海底邊界條件下的淺海平行波導空間聲場。在FLn中,當k0RLn足夠大時可將B忽略,則FLn≈F。因為RLn隨n的增大而增大,當hL=50m時,RLn≥50m。令k0=1rad/m,RLn分別取50、500、1 000m,其他參數同上。反射系數隨掠射角的變化曲線如圖2所示,其中反射系數由FLn和F的絕對值表示。

圖2 反射系數隨掠射角的變化曲線Fig.2 The variation of reflectances with the glancing angle

圖2中FLn曲線隨著RLn的增大而靠近F曲線,當RLn=500m時,兩條曲線已經十分接近,RLn=1 000m時,兩條曲線基本重合。同時可看出當掠射角超過30°后幾條曲線幾乎重合,所以在大角度掠射情況下k0RLn對反射系數的影響較小。由于在40°~90°區(qū)間內,反射系數下降非常緩慢,為方便收斂性分析及聲場計算,此區(qū)間內反射系數可近似看作常數。

2 收斂性分析和算例驗證

在Matlab中編程分析格林函數G的收斂性,設頻率f=200 Hz和800 Hz,hU=10 m,hL=50 m,Q點坐標為(0, 0, 0),P點坐標為(0, 50 m, 0),格林函數幅值|G|隨反射次數n的變化曲線如圖3所示。

圖3 格林函數收斂性曲線Fig.3 The convergence curves of green functions

由圖3可看出吸聲海底條件下的格林函數很快收斂,而剛性海底條件下的格林函數波動明顯,收斂緩慢,這是因為兩種海底邊界的反射能力有很大差別,導致疊加效果不同。為保證收斂,在下文的計算中兩種海底條件下的n分別取5和20。

在海底吸聲邊界條件下,以淺海平行波導空間中脈動球輻射解析解驗證提出的聲學邊界元法的正確性。在脈動球表面劃分若干個四邊形單元,球心位于坐標原點o,f=200 Hz,場點分布為y軸上10~50 m的一條直線。此時FLn近似為常數,根據疊加原理[1]可知場點聲壓解析解為

(6)

式中:vn=1×10-6m/s,脈動球半徑a=1m,r0為球心到場點的距離,rUn和rLn分別是球心在海面以上和海底以下的虛源到場點的距離。

圖4展示了吸聲海底條件下聲學邊界元法與解析法計算的聲壓級(SPL)曲線,可看出兩種方法的結果十分吻合,證明所提方法是正確的。

3 淺海中圓柱殼的聲輻射特性

簡支在半無限障板上的圓柱殼模型如圖5所示,圓柱殼長Lc=1.284 m,半徑Rc=0.18 m,厚度為0.003 m,材料的楊氏模量為2.1×1011Pa,泊松比為0.3,密度為7 850 kg/m3。圓柱殼軸線與z軸重合,圓柱殼中心位于坐標系原點o,場點P在平面xy上,P到o的距離為r,γ為周向角。在圓柱殼下表面沿x軸施加頻率為f,幅值為1 N的法向簡諧激勵力,激勵力指向o。利用有限元法計算圓柱殼流固耦合模型的表面振速,當h=50 m時其均方振速在圓柱殼離開邊界五倍半徑距離后不再變化,這種情況與文獻[17]研究結果一致。當h增大,邊界影響將更小,故在h≥50 m情況下,圓柱殼離開邊界五倍半徑距離后可采用無限流域中圓柱殼的振速輸入聲學邊界元Matlab程序中計算不同水深情況下圓柱殼的聲輻射特性。

圖4 脈動球場點聲壓級對比Fig.4 The comparison of the sound pressure level

圖5 圓柱殼模型Fig.5 The model of cylindrical shell

設f=200 Hz,hU=10 m,hL=50 m,不同半徑下的聲壓級指向性曲線如圖6所示,其中SPLu、SPLa和SPLr分別表示忽略海底只考慮海面邊界、海面+吸聲海底和海面+剛性海底邊界條件下的聲壓級。從圖6可看出兩種海底邊界條件下的聲壓級指向性曲線總體上存在明顯差別,SPLa曲線比SPLr曲線更接近SPLu曲線,這是因為吸聲海底對聲波的反射能力較弱。當hU和hL不變,r越大,SPLa和SPLr曲線與SPLu曲線的差別越大,這是因為不同場點處的疊加效果不同。

圖6 不同半徑下的聲壓級指向性曲線Fig.6 The directivity of the sound pressure level with different radius

令f=200 Hz,hL=100 m,r=15 m,浸深hU分別為50 m和200 m,且考慮靜水壓力作用,不同浸深下的聲壓級指向性曲線如圖7所示,其中SPLi為無限流域的聲壓級。從圖7可看出隨著圓柱殼離開海面,SPLa曲線越來越接近SPLi曲線,這是因為圓柱殼離海底足夠遠,其反射作用可忽略,而海面的反射作用隨浸深增加而減弱。再令hU=hL=50 m,聲壓級沿x軸的變化如圖8所示。從圖8看出,在圓柱殼下方三條曲線差別顯著,SPLu曲線單調下降,而SPLa和SPLr曲線波動下降,并且SPLr曲線波動的劇烈程度遠高于SPLa曲線;在圓柱殼上方三條曲線差別較小,波動情況相同。進一步發(fā)現,聲壓在任何情況下都表現出一種周期波動性,表明邊界對聲波的反射造成深度方向上出現駐波。波峰位置和波峰間距L并不因邊界條件改變而改變,當f=200 Hz時波長λ=7.5 m,L=3.79 m,λ/L=1.98;當f=800 Hz時波長λ=1.875 m,L=0.94 m,λ/L=1.99,波長與波峰間距之比近似為一常數。

當hU=10 m,hL=50 m時,兩種海底邊界條件下的聲壓級分布云圖如圖9所示。從圖9可見兩種海底邊界條件下的聲壓級分布有明顯不同,吸聲海底邊界條件下的聲壓變化較小。在整個空間都有駐波分布,頻率越高,駐波越密集,疊加聲場中波峰波谷交錯出現。

4 圓柱殼聲輻射特性的產生機理

圓柱殼輻射聲場中P點聲壓可以看作是無數個不同強度點源的輻射聲壓和反射聲壓在P點的疊加,式(2)中格林函數第一項幅值為G1,代表輻射聲,第二項幅值為G2,代表反射聲。當f=200 Hz時,設hU=10 m,hL=50 m,Q點坐標為(0, 0, 0),P點坐標為(0,y, 0);再設hL=100 m,P點坐標為(0, 15 m, 0),G1和G2隨P點位置和浸深變化如圖10所示。從圖10看到G2/G1隨著P點離開圓柱殼而增大,隨著圓柱殼離開海面而減小,其根本原因是反射聲與輻射聲的成分比例隨場點到虛源距離的變化而變化。

圖7 不同浸深下的聲壓級指向性曲線Fig.7 The directivity of the sound pressure level at different submerged depths

圖8 聲壓級沿x軸的變化Fig.8 The variation of the sound pressure level in the direction of x-axis

圖9 聲壓級分布云圖Fig.9 The nephogram of the sound pressure level

圖10 G1和G2變化曲線Fig.10 The variation of G1 and G2

在結構聲振問題中,半無限域和無限域指在某些條件下可以忽略某些邊界的聲波反射作用的情況。在淺海波導域中,當P點離圓柱殼較近和圓柱殼距離海面或海底足夠遠時,有可能忽略邊界反射作用,在此情況下布置測點得到圓柱殼的輻射聲壓等同于半無限域或無限域測試結果。為尋找滿足聲學半無限域或無限域條件的P點,可采用點源格林函數進行預測,不同參數下的格林函數幅值如圖11所示,其中Ga、Gu和Gi表示海面+吸聲海底、忽略海底只考慮海面和無限流域中的格林函數幅值。從圖11可見格林函數幅值的匹配情況與圖6和圖7中聲壓級的匹配情況完全一致,這點進一步解釋了圖6和圖7所展示現象的產生原因,說明采用點源格林函數預測淺海中圓柱殼聲輻射特性是可行的。

圖11 格林函數幅值Fig.11 The amplitude of the green function

聲壓在深度方向上的周期性波動現象同樣源自點源的聲壓疊加性。圖1中實源和無限個虛源組成一條點源鏈,其中點源相位周期性變化。以構成一個周期的四個點源QU1、Q、QL1、QL2為例,設點源強度為A,則該周期的點源在P點的疊加聲壓為

(7)

當P點下降ΔR=λ/2時,P點的疊加聲壓變?yōu)?/p>

(8)

由式(7)、(8)和指數函數的周期性質可知,p0和p1具有相同的周期性,當| p0|為極大值時| p1|也為極大值,距離和邊界條件的變化只改變幅值大小而不改變周期性,這便是λ/L=2的原因。

當只考慮海面反射作用時,實源和虛源組成相位相反的偶極子,P點的疊加聲壓為

(9)

當P點在實源和海面之間下降ΔR=λ/2時,P點的疊加聲壓變?yōu)?/p>

(10)

同樣,p2和p3具有相同的周期性,當時| p2|為極大值時| p3|也為極大值,在實源和海面之間λ/L=2。

當P點在實源下方下降時,λ/L=2規(guī)律的產生條件便不再滿足。此時RU1=R0+2hU恒成立,P點的疊加聲壓為:

(11)

可得到:

(12)

由式(12)可知| p4|隨著P點下降而減小,故圖8中SPLu在圓柱殼下方單調衰減。對應圖8的格林函數幅值如圖12所示(Gr表示海面+剛性海底邊界條件下的格林函數幅值),從圖12可看出格林函數幅值表現出與聲壓級相同的變化規(guī)律。以上現象是由流域邊界對聲波的反射作用產生的,而這點是早期圓柱殼聲輻射問題研究文獻并未關注的[21]。

圖12 格林函數幅值沿x軸的變化Fig.12 The variation of the green function in the direction of x axis

5 結論

1) 剛性和吸聲海底邊界條件下圓柱殼的輻射聲壓幅值有明顯差異,吸聲海底的聲場疊加效果總體弱于剛性海底的聲場疊加效果。

2) 可采用點源格林函數預測淺海波導域中圓柱殼的聲輻射特性和尋找滿足聲學半無限域或無限域條件的場點,邊界反射作用隨著圓柱殼離開海面或海底而減弱,隨著場點離開圓柱殼而增強。當圓柱殼與海面和海底距離滿足一定條件時,其近場聲輻射特性與半無限域或無限域中情況相同。

3)海面和海底對聲波的多次反射在淺海波導域中產生復雜的疊加聲場,波峰波谷交錯出現,頻率越高,駐波越密集,其中深度方向上的駐波遵循波峰間距為波長一半的規(guī)律,此規(guī)律不受邊界類型影響。

[1]黎勝, 趙德有. 半空間內結構聲輻射研究[J]. 船舶力學, 2004, 8(1): 106-112.

LI Sheng, ZHAO Deyou. Research on acoustic radiation in a three-dimensional half space[J]. Journal of ship mechanics, 2004, 8(1): 106-112.

[2]SEYBERT A F, SOENARKO B. Radiation and scattering of acoustic waves from bodies of arbitrary shape in a three-dimensional half space[J]. Journal of vibration and acoustics, 1988, 110(1): 112-117.

[3]鄒元杰, 趙德有, 黎勝. 自由液面和剛性壁面對結構振動聲輻射的影響[J]. 聲學學報, 2005,30(1): 89-96.

ZOU Yuanjie, ZHAO Deyou, LI Sheng. Impact of soft surface and hard plane on structural vibration and acoustic radiation[J]. Acta acustica, 2005, 30(1): 89-96.

[4]ZHOU Q, JOSEPH P F. A numerical method for the calculation of dynamic response and acoustic radiation from an underwater structure[J]. Journal of sound and vibration, 2005, 283:853-873.

[5]BRUNNER D, JUNGE M, CABOS C, et al. Vibroacoustic simulation of partly immersed bodies by a coupled fast BE-FE approach[J]. Journal of the acoustical society of america, 2008, 123(5): 3418.

[6]JUNGE M, BECKER J, BRUNNER D, et al. FE-model reduction for FE-BE Coupling with large fluid structure interfaces[J]. Journal of the acoustical society of America, 2008, 123(5): 3726.

[7]CHEN P T, LIN C S, YANG T. Responses of partially immersed elastic structures using a symmetric formulation for coupled boundary element and finite element methods[J]. Journal of the acoustical society of America, 2002, 112(3): 866-875.

[8]JUNGE M, BRUNNER D, BECKER J, et al. Interface-reduction for the craig-bampton and rubin method applied to FE-BE coupling with a large fluid-structure interface[J]. International journal for numerical methods in engineering, 2009, 77(12): 1731-1752.

[9]ZHENG C J, CHEN H B, CHEN L L. A wideband fast multipole boundary element method for half-space/plane-symmetric acoustic wave problems[J]. Acta mechanica sinica, 2013, 29(2): 219-232.

[10]陳爐云, 王德禹, 張立軍. 直角域內的結構聲輻射特性研究[J]. 船舶力學, 2011, 15(1): 175-181.

CHEN Luyun, WANG Deyu, ZHANG Lijun. Research on acoustic radiation characteristic in a three-dimensional right-angle space[J]. Journal of ship mechanics, 2011, 15(1): 175-181.

[11]CHEN L, ZHANG Y. Acoustic radiation analysis based on essential solution of Green’s function[J]. Journal of Shanghai Jiaotong University: Science, 2013, 18: 409-417.

[12]CHEN L, LIANG X, YI H. Acoustic radiation analysis for a control domain based on Green′s function[J]. Applied mathematical modelling, 2016, 40(4): 2514-2528.

[13]鄒元杰, 趙德有. 結構在淺水中的振動和聲輻射特性研究[J]. 振動工程學報, 2005, 17(3): 269-274.

ZOU Yuanjie, ZHAO Deyou. A vibro-acoustic study on structures in shallow water[J]. Journal of vibration engineering, 2005, 17(3): 269-274.

[14]CHEN H, LI Q, SHANG D. Fast prediction of acoustic radiation from a hemi-capped cylindrical shell in waveguide[J]. Journal of marine science and application, 2014, 13(4): 437-448.

[15]ZOU M, WU Y, LIU Y, et al. A three-dimensional hydroelasticity theory for ship structures in acoustic field of shallow sea[J]. Journal of hydrodynamics, Ser. B, 2013, 25(6): 929-937.

[16]ZOU M S, WU Y S, LIU Y M. The application of three-dimensional hydroelastic analysis of ship structures in Pekeris hydro-acoustic waveguide environment[J]. Acta mechanica sinica, 2014, 30(1): 59-66.

[17]白振國, 吳文偉, 左成魁,等. 有限水深環(huán)境圓柱 殼聲輻射及傳播特性[J]. 船舶力學, 2014, (1): 178-190.

BAI Zhenguo, WU Wenwei, ZUO Chengkui, et al. Sound radiation and spread characteristics of cylindrical shell in finite depth water[J]. Journal of ship mechanics, 2014, (1): 178-190.

[18]CHEN L H, SCHWEIKERT D G. Sound radiation from an arbitrary body[J]. The journal of the acoustical society of America, 1963, 35(10): 1626-1632.

[19]ZAMPOLLI M, TESEI A, CANEPA G, et al. Computing the far field scattered or radiated by objects inside layered fluid media using approximate Green’s functions[J]. The journal of the acoustical society of America, 2008, 123(6): 4051-4058.

[20]汪德昭, 尚爾昌. 水聲學[M]. 北京:科學出版社, 1981: 175-178.

WANG Dezhao, SANG Erchang. Hydroacoustics[M]. Beijing: science press, 1981.175-178.

[21]周鋒, 駱東平, 蔡敏波, 等. 有限長環(huán)肋圓柱殼低階模態(tài)聲輻射性能分析[J]. 應用科技, 2004, 31(9): 38-41.

ZHOU Feng, LUO Dongping, CAI Minbo, et al. Sound radiation analysis of low order modes from the ring-stiffened cylindrical shells in fluid medium[J]. Applied science and technology, 2004, 31(9): 38-41.

本文引用格式:

繆宇躍, 李天勻, 朱翔,等. 淺海中圓柱殼的聲輻射特性分析[J]. 哈爾濱工程大學學報, 2017, 38(5): 719-726.

MIAO Yuyue, LI Tianyun, ZHU Xiang, et al. Research on the acoustical radiation characteristics of cylindrical shells in ashallow sea [J]. Journal of Harbin Engineering University, 2017, 38(5): 719-726.

Research on the acoustical radiation characteristics of cylindrical shells in ashallow sea

MIAO Yuyue1, LI Tianyun2,3,4, ZHU Xiang2,3,4, WANG Peng2,3, ZHANG Guanjun2,3

(1.Ship Development and Design Center, Wuhan 430064, China; 2.School of Naval Architecture and Ocean Engineering,Huazhong University of Science & Technology, Wuhan 430074, China; 3.Hubei Key Laboratory of Naval Architecture & Ocean Engineering Hydrodynamics, Wuhan 430074, China; 4.Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration, Shanghai 200240, China)

According to the acoustical boundary element theory on acoustics and considering the sound reflection effect of the sea surface and sound-absorbing seabed, the Green Function was modified, and an acoustical radiation model in the waveguide domain was established. In addition, the acoustical radiation properties of a cylindrical shell in the shallow sea was analyzed. It is found that the sound absorption is apparently different from the radiative sound field of the cylindrical shell under the conditions of a rigid seafloor. The Green Function on the point source in the parallel waveguide space can be used to forecast the acoustical radiation of the cylindrical shell. When a cylindrical shell leaves the sea surface or seafloor, the boundary reflection effect will weaken; when a field point leaves the cylindrical shell, the effect will strengthen. In the direction of the water’s depth, the sound pressure fluctuates periodically, its periodic law is irrelevant to the type of seabed boundary. The radiation superposition principle of the acoustical point source is used to explain the formation theory of such a phenomenon, which can provide a theoretical reference for engineering applications.

boundary element method; shallow sea; cylindrical shell; acoustical radiation; sound-absorbing seabed; waveguide; Green function; point source

2016-01-15.

日期:2017-04-26.

國家自然科學基金項目(51379083,51479079,51579109);高等學校博士學科點專項科研基金項目(20120142110051).

繆宇躍(1988-), 男, 博士研究生; 李天勻(1969-), 男, 教授, 博士生導師.

李天勻,E-mail:ltyz801@ hust.edu.cn

10.11990/jheu.201601057

TB532

A

1006-7043(2017)05-0719-08

網絡出版地址:http://www.cnki.net/kcms/detail/23.1390.u.20170426.1040.014.html

主站蜘蛛池模板: 日韩第一页在线| 久久精品免费看一| 国产在线无码av完整版在线观看| 波多野衣结在线精品二区| 国产美女叼嘿视频免费看| 国产成人欧美| 久久99久久无码毛片一区二区| 免费国产一级 片内射老| 自拍偷拍一区| 亚洲人成影视在线观看| 久久精品国产电影| 波多野结衣AV无码久久一区| 精品人妻一区二区三区蜜桃AⅤ| 精品久久久久成人码免费动漫| 大乳丰满人妻中文字幕日本| 午夜无码一区二区三区| 91亚洲精品第一| 丁香五月亚洲综合在线| 国产杨幂丝袜av在线播放| 99热国产这里只有精品9九| 片在线无码观看| 五月婷婷伊人网| 5388国产亚洲欧美在线观看| a毛片免费在线观看| 日韩AV手机在线观看蜜芽| 国产SUV精品一区二区| 精品免费在线视频| 中文字幕久久亚洲一区| aa级毛片毛片免费观看久| 欧美成人午夜视频| av在线人妻熟妇| 国禁国产you女视频网站| 久青草网站| 中国美女**毛片录像在线 | 中文一区二区视频| av免费在线观看美女叉开腿| 自拍中文字幕| 91日本在线观看亚洲精品| 五月天香蕉视频国产亚| 1769国产精品视频免费观看| 亚洲五月激情网| 欧美成人午夜视频免看| 精品人妻系列无码专区久久| 手机在线免费不卡一区二| 亚洲第一中文字幕| 国产成人精品男人的天堂下载| 国产成人乱无码视频| 一本久道热中字伊人| 日韩乱码免费一区二区三区| 日韩av无码DVD| 亚洲自拍另类| 国产精品任我爽爆在线播放6080 | 国产91高清视频| 精品视频在线一区| 日韩第九页| 中文字幕在线观| 一区二区偷拍美女撒尿视频| 久久情精品国产品免费| 欧美高清三区| 亚洲欧美极品| 国产在线观看99| 久久无码免费束人妻| 在线视频精品一区| 亚洲欧美精品一中文字幕| 99人体免费视频| 99这里只有精品6| 无码内射中文字幕岛国片 | 亚洲日本中文字幕乱码中文| 色欲国产一区二区日韩欧美| www.av男人.com| 日韩a级毛片| 无码免费的亚洲视频| 欧美日韩国产精品va| 日韩欧美中文字幕在线韩免费| 亚洲欧美成人网| 特级欧美视频aaaaaa| 99久久性生片| 无码免费试看| 动漫精品啪啪一区二区三区| 亚洲天堂视频在线观看免费| 1级黄色毛片| 国产精品福利社|