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

    軸對稱指向性球面波的界面反射波*

    2022-04-15 07:33:20段韻達胡恒山
    物理學報 2022年7期
    關(guān)鍵詞:界面

    段韻達 胡恒山

    (哈爾濱工業(yè)大學航天科學與力學系,哈爾濱 150001)

    無限大平面剛性障板中圓形活塞源的聲輻射場可近似為軸對稱指向性球面波,前人只給出了活塞面與界面平行時軸對稱指向性球面波的界面響應(yīng)表達式,本文針對活塞面與界面不平行的情況,推導了軸對稱指向性球面波的錐面波展開式,并進一步導出了其界面反射波的表達式.在源距遠大于聲波波長的情況下通過鞍點法將界面反射波的表達式化簡為了簡化表達式.簡化式不僅計算上簡潔,而且物理含義清楚:軸對稱指向性球面波的界面反射波可視為鏡像活塞源激發(fā)的軸對稱指向性球面波與反射系數(shù)的乘積.計算表明,當活塞半徑小于聲波波長時,反射波對活塞與界面的夾角和接收點的環(huán)向方位角不太敏感,反射波的指向性較弱;當活塞半徑大于聲波波長時,反射波對活塞與界面的夾角和接收點的環(huán)向方位角很敏感,反射波的指向性很強.增加活塞與界面的夾角,反射波先增加后減小,反射波的指向性先增強后減弱;當活塞與界面的夾角等于活塞中心鏡像點與接收點的連線與界面法線的夾角時,反射波最大,反射波的指向性最強.

    1 引言

    聲波是海洋中探測目標和傳遞信息的有效工具,而聲納是利用聲波獲取水下信息的儀器,被廣泛用于海洋軍事和地質(zhì)勘探[1,2].聲納的重要組成之一是換能器,換能器類型多種多樣,其中復合棒換能器是目前應(yīng)用最為普遍的換能器之一,由于它產(chǎn)生的是一種類似于活塞式的振動,因而也稱為活塞式換能器或縱向換能器[1,3].聲源通常嵌入到一堅硬面板中,使聲源的輻射場只存在于面板的單側(cè)空間,并使該側(cè)空間內(nèi)的聲輻射場比無面板時更大,此時該面板稱為障板[4].障板可以是球形、圓柱或平面,其中平面障板是最基本的障板模型,當聲波波長遠小于障板尺度時,有限平面障板模型可由無限大平面障板模型描述[5,6].多個活塞安裝在障板上形成陣列,可提高聲源的發(fā)射功率,使聲輻射集中于同一方向,更有利于定位水下目標[7].

    無限大平面剛性障板中圓形活塞源的聲輻射是經(jīng)典的聲輻射模型[8-10],其中活塞面內(nèi)所有的面元從零時刻開始以相同的速度沿活塞軸線方向振動,障板與活塞面共面,障板在任意時刻都不會發(fā)生變形和運動.當活塞中心到接收點的距離遠大于活塞半徑時,無限大平面剛性障板中圓形活塞(簡稱為障板活塞)的聲輻射場可近似為一種軸對稱指向性球面波(軸對稱是指波場關(guān)于活塞的軸線對稱),即[4]

    式中第2 個等號的左邊為障板活塞源的輻射場,右邊為軸對稱指向性球面波;pf是流體中的聲壓(流體可以是水或者空氣),a為活塞半徑,kf是流體中的聲波波數(shù);RMP為接收點M到活塞面上任意點P的距離,SP為圓形活塞面,RME為接收點M到活塞中心E的距離;D(θME) 為指向性函數(shù),有

    其中J1為一階第一類貝塞爾函數(shù),θME為接收點M與活塞中心E的連線與活塞軸線的夾角,D(θME)的量綱為1;A(ω) 為聲壓源的頻譜(單位為Pa),有

    其中v0(ω)為活塞振動速度,ρf為流體密度,ω為角頻率,i 為單位虛數(shù).這種軸對稱指向性球面波的輻射阻抗和指向特征對指導水下?lián)Q能器設(shè)計[6,11]、水雷探測器設(shè)計[12]、隧道微壓波計算[13]和醫(yī)學超聲渦旋聲場設(shè)計[14]有著重要意義.

    前面介紹了軸對稱指向性球面波在聲輻射方面的諸多研究,但軸對稱指向性球面波的界面響應(yīng)卻少見報道.Amédin 等[15]針對孔隙介質(zhì)夾層介于障板活塞與空氣之間的情況,推導了障板活塞源的聲輻射場對孔隙介質(zhì)-空氣界面的響應(yīng)表達式,當孔隙介質(zhì)夾層厚度遠大于活塞半徑時,該障板活塞源的界面響應(yīng)式等于軸對稱指向性球面波的界面響應(yīng)式.Wang和Cho[16]考慮了各向異性孔隙介質(zhì)夾層介于障板活塞與空氣之間的情況,當孔隙介質(zhì)夾層厚度遠大于活塞半徑時,Wang和Cho 給出的障板活塞源的界面響應(yīng)式等于軸對稱指向性球面波的界面響應(yīng)式.Schakel 等[17-19]為了檢驗孔隙介質(zhì)動電效應(yīng)理論[20,21]的有效性,推導了指向性球面聲波入射流體-孔隙介質(zhì)界面的響應(yīng)式.文獻[15-19]都是基于人為設(shè)計的實驗?zāi)P徒o出軸對稱指向性球面波的界面響應(yīng)式,由于實驗?zāi)P椭谢钊媾c界面平行,所以他們只考慮活塞面與界面平行這一特殊情況,但從理論模型的角度看,活塞面與界面不平行的情況比兩者平行的情況更為一般,然而目前未見活塞與界面不平行時軸對稱指向性球面波的界面響應(yīng)式的文獻報道.從應(yīng)用的角度看,海洋地震勘探和水下目標探測時活塞面與界面通常不平行[5,22],在沒有理論公式的情況下,雖然可通過數(shù)值仿真模擬軸對稱指向性球面波的界面反射波,然而當目標界面遠離活塞面和接收點時,數(shù)值仿真將面臨很大的計算量和很長的計算時間,這不利于海底儲層和水下目標的探測.因此,本文以軸對稱指向性球面波從流體入射固體界面為例,推導活塞面不平行于界面時軸對稱指向性球面波的界面反射波的表達式,完善障板活塞源的聲輻射場的界面響應(yīng)理論;推導軸對稱指向性球面波的界面反射波的簡化表達式,為障板活塞型聲源快速探測水下界面提供思路.

    2 軸對稱指向性球面波的錐面波展開

    欲得均勻球面波的界面反射波需先導出均勻球面波的錐面波展開式[23],類似地,欲得軸對稱指向性球面波的界面反射波,需先導出軸對稱指向性球面波的錐面波展開式.由(1)式知RME ?a時障板活塞源輻射場近似為軸對稱指向性球面波,因而本節(jié)從障板活塞源輻射場出發(fā),先導出障板活塞源輻射場的錐面波展開式,再令RME ?a,則該錐面波展開式近似為軸對稱指向性球面波的錐面波展開式.圖1 是無限大平面剛性障板中的圓形活塞與界面的示意圖,坐標系建立的方式為:從圓形活塞中心E向界面作垂線,交點為O,以EO的連線為z軸,方向向下;使活塞軸線與界面相交,該交點與O點的連線為x軸,方向向右,Oxyz構(gòu)成右手直角坐標系;使Oxy面與界面重合.活塞中心E到界面的距離為dE,活塞軸線與界面法線的夾角為α,亦即活塞面與界面的夾角為α.本節(jié)只考慮入射波的情況,下節(jié)再考慮界面反射波的情況.

    圖1 無限大平面剛性障板中的圓形活塞與界面的示意圖Fig.1.Schematic diagram of a piston in an infinite plane rigid baffle and an interface.

    2.1 障板活塞源輻射場的錐面波展開式

    無限大平面剛性障板中圓形活塞的模型如圖1所示,圖中P是活塞面內(nèi)的點,M為接收點,由(1)式可知障板活塞源的聲壓輻射場為[4]

    (4)式表明障板活塞源的輻射場可視為活塞面SP內(nèi)頻譜為的點源輻射場在M點處疊加(點源指均勻球面波的點源,下文均為該含義),所以可利用均勻球面波的錐面波展開式處理(4)式.均勻球面波的錐面波展開式為[23]

    由于徑向波數(shù)kr及kr積分的上下限均與空間坐標無關(guān),所以可以交換(6)式中的積分次序,得

    由圖1 知,圓形活塞的面元 dSP與其在Oxy面的投影面元 dSF滿足

    這表明圓形活塞面SP在Oxy面的投影SF是一個橢圓,并且該橢圓的周線方程在直角坐標系和極坐標系下可分別寫為

    記P點在Oxy面上的投影點為F,由(8)式可知

    其中dE為活塞中心E到Oxy面的距離,-acosα≤xF≤acosα(因xF位于橢圓SF內(nèi)).由于P在Oxy面上的投影為F,所以rMPrMF,將rMPrMF,(8)式和(10)式代入(7)式得

    其中

    由貝塞爾函數(shù)加法定理可得[24]

    其中εn為紐曼因子,滿足ε01,而n≥1時εn2 ;rM和rF分別為M點和F點的徑向坐標,φM和φF分別為M點和F點的環(huán)向坐標(見圖1).將(13)式代入(12)式,并交換積分和求和的次序得

    其中

    利用三角關(guān)系,可將(15)式改寫為

    觀察(16)式第二行的積分式可知,n≥1 時Jn(krrF)和 exp(ikzxFtanα)均關(guān)于φF0對稱,而sin(nφF)關(guān)于φF0反對稱,又橢圓面SF關(guān)于φF0 對稱(見(9)式),所以n≥1 時(16)式第二行的積分式為零.此外,n0 時 sin(nφF)0,所以綜上知n≥0時(16)式第二行的積分項為零,得

    其中

    將(17)式和(14)式代入(11)式得障板活塞源輻射場的錐面波展開式

    其中Gn由(18)式表達.(19)式表明障板活塞源的輻射場可視為頂點位于z軸zM-dE處的多階錐面波在M點處的疊加.因(4)式推倒到(19)式的過程中沒有用到任何近似條件,而(4)式適用于任意場點[4],所以(19)式也適用于任意場點.

    2.2 軸對稱指向性球面波的錐面波展開式

    由(1)式可知,當活塞中心E到接收點M的距離遠大于活塞半徑時(RME ?a),障板活塞源輻射場(4)式近似為軸對稱指向性球面波,所以RME ?a時障板活塞源輻射場的錐面波展開式(19)近似為軸對稱指向性球面波的錐面波展開式.RME ?a時(19)式不僅能用來推導軸對稱指向性球面波的界面反射波,也能用來推導界面透射波,本文關(guān)注的是界面反射波.

    2.3 算例檢驗

    前文看到,RME ?a時(19)式為軸對稱指向性球面波的錐面波展開式,而RME ?a時(1)式為軸對稱指向性球面波(即(1)式中的exp(ikfRME)),所以RME ?a時(19)式應(yīng)等于(1)式.因此,這里為了檢驗(19)式,在RME ?a的條件下比較(19)式和(1)式.

    由軸對稱指向性球面波的表達式(1)式和(2)式知,參數(shù)kfa控制了障板活塞源的指向能力(kf為流體中聲波波長,a為活塞半徑),因而本文考慮兩類活塞,一類是活塞半徑小于聲波波長的小活塞,即a<λ,而另一類是活塞半徑大于聲波波長的大活塞,即a >λ.(19)式中的A為聲源頻譜,設(shè)聲源波形為余弦包絡(luò)脈沖函數(shù),脈沖長度為0.5 ms,中心頻率為6 kHz,頻率范圍為2—10 kHz(高頻海洋地震勘探的聲源頻段[25]),則有λ0.15—0.75 m,所以設(shè)小活塞的半徑為a0.1 m,大活塞的半徑為a1 m .設(shè)距離小活塞源1 m 處的聲壓為 105Pa,距大活塞源1 m 處的聲壓為 107Pa.對于流體中的入射波而言,只需考慮流體的參數(shù)(見表1),設(shè)流體的品質(zhì)因子為100.(19)式和(1)式乘以聲源頻譜函數(shù)可得頻域波形,再經(jīng)傅里葉變換可得時域波形.

    表1 流體和固體參數(shù)Table 1.Parameters of fluid and solid.

    設(shè)活塞中心E位于z軸上z-dE-25 m處,活塞面與Oxy面的夾角α30°,接收點M位于 (rM,φM,zM)(21.65,0,-12.5),其中rM,φM和zM的單位分別為米、度和米,則有RME25 m?1 ≥a,因而RME ?a條件滿足.由M點與E點的位置可知(參考圖1),M點位于Oxz面,線段ME與z軸的夾角為 60°,ME與活塞軸線的夾角為 30°,因而M點與O點關(guān)于活塞的軸線對稱,則M與O點的聲壓應(yīng)相等.記錐面波展開式(19)在M點處的計算結(jié)果為pf(M),(19)式在O點處的計算結(jié)果為pf(O),(1)式在RME25 m且θME30°下的計算結(jié)果為pf0.圖2 給出了pf(M) ,pf(O)和pf0的計算結(jié)果,其中圖2(a)和圖2(b)分別對應(yīng)a0.1 m和a1 m,即分別對應(yīng)小活塞和大活塞.圖2(a)和圖2(b)表明pf(M),pf(O)和pf0幾乎沒有區(qū)別,因此RME ?a時軸對稱指向性球面波的錐面波展開式(19)與軸對稱指向性球面波(1)式一致.

    圖2 pf(M) ,pf(O) 和 pf0的時域波形 (a) a=0.1 m ;(b)a=1 mFig.2.Time-domain waveforms of pf(M),pf(O)and pf0 :(a) a=0.1 m;(b) a=1 m .

    3 軸對稱指向性球面波的界面反射波的表達式

    前文表明RME ?a時障板活塞源輻射場的錐面波展開式(19)近似為軸對稱指向性球面波的錐面波展開式,由此可知RME ?a時障板活塞源輻射場的界面反射波表達式近似為軸對稱指向性球面波的界面反射波表達式,因此本節(jié)先推導障板活塞源輻射場的界面反射波表達式,再令RME ?a即得軸對稱指向性球面波的界面反射波表達式.圖3 給出了活塞面關(guān)于界面的鏡像,圖中活塞SP1為活塞SP的鏡像,E1為E的鏡像,P1為P的鏡像,M為反射波的接收點,M2為M在Oxy面上的投影點.

    圖3 活塞面關(guān)于界面的鏡像Fig.3.Mirror image of the piston surface with respect to the interface.

    參考錐面波展開式(19)的形式,設(shè)障板活塞源輻射場的界面反射波為

    其中M點為反射波的接收點,BnBn(kr) 為待定系數(shù),需通過邊界條件確定,具體過程為:參考(19)式的形式先列出透射縱波和橫波的表達式,再將它們和(20)式代入流體-固體界面邊界連續(xù)條件可得Bn滿足的線性方程組,最后求解該線性方程組即得Bn.然而,我們并不打算通過上述方式求解Bn以獲得界面反射波的表達,因為由(1)式注意到障板活塞源的輻射場等于活塞面SP上點源輻射場的疊加,則障板活塞源輻射場的界面反射波等于活塞面SP上點源輻射場的界面反射波的疊加.對于活塞面SP上的單個點源P而言,其界面反射波為[23]

    由于徑向波數(shù)kr及kr的積分上下限均與空間坐標無關(guān),所以可以交換(22)式中的積分次序,得

    由于P1為P的鏡像,所以由(10)式得

    由 于SP在Oxy面的投影為SF,而SP和關(guān) 于Oxy面對稱,所以SP1面在Oxy面的投影也是SF(參照圖3),則有·cosαdSF和rMF,再將這兩式代入(23)式得

    將(24)式代入(25)式得

    其中

    其中Gn由(18)式表達;εn為紐曼因子,滿足ε01,而n≥1時εn2 .至此,我們導出了活塞面與界面不平行時障板活塞源的界面反射波表達式(28)式,令?a則(28)式近似為活塞面與界面不平行時軸對稱指向性球面波的界面反射波表達式,相比前人只考慮活塞面與界面平行的情況而言[15-19],(28)式完善了軸對稱指向性球面波的界面響應(yīng)理論.比較(28)式和(20)式可知BnεnB(kr),而B(kr)已由前人給出[23],所以無需求解Bn滿足的線性方程組,通過點源疊加的方式就可以獲得軸對稱指向性球面波的界面反射波.此外,前文的算例表明RME ?a時軸對稱指向性球面波的錐面波展開式(19)正確,而這里從(19)式到(28)式的過程中用的都是恒等式,因而RME ?a時軸對稱指向性球面波的界面反射波(28)式也正確,從而無需再通過算例檢驗(28)式.

    4 軸對稱指向性球面波的界面反射波的簡化表達式

    當活塞半徑遠小于聲波波長 (a ?λ) 時有kfa ?1,代入(2)式得D(θME)→1,即指向性函數(shù)幾乎不隨M點的方位改變,則由(1)式知軸對稱指向性球面波退化為均勻球面波,從而其界面反射波退化為均勻面波的界面反射波,因此a ?λ的情況無需研究,而對于a不遠小于λ的一般情況,本文研究軸對稱指向性球面波的界面反射波的簡化表達式.海洋地震勘探[25]和水下目標超聲探測[26,27]時,活塞中心E關(guān)于界面的鏡像E1到接收點M的距離通常遠大于聲波波長(?λ),則鞍點法的使用條件?1 得到滿足,從而可用鞍點法化簡軸對稱指向性球面波的界面反射波的表達式(28).

    4.1 基于鞍點法的簡化表達式

    由圖3 可知 (zM -dE)-|zM -dE|,將其代入(28)式,再將n寫為偶數(shù)和奇數(shù)階的形式得

    由貝塞爾函數(shù)教材可知[28]J2m(x)J2m(-x)和J2m+1(x)-J2m+1(-x),其中m0,1,2,3,···,再代入(18)式得

    由貝塞爾函數(shù)教材可知[28]

    利用(30)式和(31)式,并考慮到kz和B(kr) 是關(guān)于kr的偶函數(shù)[23],則(29)式可改寫為

    再將(33)式代入(32)式得

    其中

    其中Gn(kr) 由(18)式表達.由附錄A 可知,最速下降路徑和鞍點是由e 指數(shù)函數(shù)中的宗量決定,因而需確定(34)式中所有的e 指數(shù)函數(shù).(34)式中B(kr)是平面波的反射系數(shù),不是e 指數(shù)函數(shù)[23],而G*(kr) 也不是e 指數(shù)函數(shù)(見附錄B),因此(34)式中只有 exp[i(krrM+kz|zM -dE|)] 是e 指數(shù)函數(shù),則(34)式可改寫為

    由圖3 知M2是M在Oxy面的投影,E1是E在Oxy面的投影,則MM2與OE1平行,OM2垂直于OE1,又rM,所以

    其中θME為M與E的連線與活塞SP軸線的夾角.由(42)式可得

    其中 (rM,φM,zM) 為M點的柱坐標(見圖3).比較反射波簡化式(44)和入射波(1)式可知,軸對稱指向性球面波的反射波可視為鏡像源的波乘以反射系數(shù)B,其中鏡像源的波指鏡像活塞面激發(fā)的軸對稱指向性球面波.因此簡化式(44)的物理含義比(28)式更為清楚.此外,與(28)式相比,簡化式(44)無需計算貝塞爾函數(shù)在橢圓面上積分、不同階錐面波的求和及波數(shù)域的積分,所以簡化式(44)省去了巨大的計算量,減少了大量計算時間.另外,需要說明的是,(1)式表明當活塞中心到接收點的距離遠大于活塞半徑時(RME ?a),障板活塞源輻射場近似為軸對稱指向性球面波,所以我們考慮軸對稱指向性球面波時應(yīng)默認RME ?a成立,從而考慮其界面反射波時有?a成立.因此簡化式(44)的成立條件是?a和?1同時成立,其中?1 是鞍點法的使用條件,而?a是默認成立的條件.

    值得一提的是,簡化式(44)也可用另一種更為簡單的方式導出.因考慮軸對稱指向性球面波的界面反射波時默認?a成立,則由?1條件可導出?1,即對于活塞面上任意一點P鞍點法的使用條件成立,所以可通過鞍點法將活塞面SP上單個點源P的輻射場的界面反射波(21)式化簡為[23]

    再利用(1)式可得

    最后將(49)式代入(48)式即得簡化式(44).

    4.2 算例檢驗

    本節(jié)通過算例檢驗簡化式(44)的正確性,即比較(44)式和(28)式的一致性.聲源設(shè)定與前文算例相同,流體和固體參數(shù)見表1,地層縱波和橫波的品質(zhì)因子分別為80和60,流體中聲波的品質(zhì)因子為 100 .

    參照圖3,設(shè)活塞面SP與界面Oxy的夾角為α30°,活塞中心E與界面的距離為dE5 m,則E點位于z軸上z-dE-5 m 處;設(shè)接收點M的位置為,其中rM,φM和zM的單位分別為米(m)、度(°)和米(m),此后默認為該單位.由M點和E1點的坐標可知20 m?1 m,而小活塞和大活塞的半徑分別為a0.1 m和a1 m,所以有?a;由于流體中聲波波長為λ0.15—0.75 m,則/λ≥26.7 m?1 m,所以鞍點法的使用條件?1成立.綜上可得?a和?1 成立,則簡化式(44)可以用來計算軸對稱指向性球面波的界面反射波.

    記(28)式和簡化式(44)計算幅度譜分別為Pref和Pref_jian.圖4 為Pref與Pref_jian的比較,其中圖4(a)對應(yīng)a0.1 m 時快速地層,圖4(b)對應(yīng)a0.1 m時慢速地層,圖4(c)對應(yīng)a1 m 時快速地層,圖4(d)對應(yīng)a1 m 時慢速地層.圖4(a)—圖4(d)表明,無論快速地層還是慢速地層,無論小活塞還是大活塞,均有Pref_jian與Pref一致,因此簡化式(44)與(28)式一致,即簡化式(44)正確.另外,因簡化式(44)是通過鞍點法導出的,所以鞍點法的使用條件?1 決定了簡化式(44)的適用范圍.上述算例中聲源頻率為高頻海洋地震勘探頻率(2—10 kHz),并且有20 m,所以簡化式(44)適用于≥20 m 的高頻海洋地震勘探問題,由于活塞中心E距地層界面 10 m 以上還是容易滿足的,所以簡化式(44)有著較為廣泛的適用范圍.另外,超聲頻率下(f > 20 kHz),≥1 m 時簡化式(44)也是適用的,但這里不再給出具體算例.

    圖4 Pref與 Pref_jian的比較 (a) a=0.1 m且快速地層;(b) a=0.1 m且慢速地層;(c) a=1 m且快速地層;(d) a=1 m 且慢速地層Fig.4.Comparison of Prefand Pref_jian:(a) a=0.1 mand fast formation;(b) a=0.1 mand slow formation;(c) a=1 m and fast formation;(d) a=1 m and slow formation.

    5 位置參數(shù)對界面反射波的影響

    本節(jié)基于簡化式(44),研究活塞面與界面的夾角α、接收點M的環(huán)向方位角φM和線段ME1與界面法線的夾角對軸對稱指向性球面波的界面反射波的影響,其中接收點M的位置由柱坐標(rM,φM,zM)描述(見圖3),活塞中心E位于(0,0,-dE),則活塞中心鏡像E1位于 (0,0,dE) .由于活塞源嵌入無限大障板中,所以障板背面不存在活塞源的輻射場,例如當活塞與界面的夾角為α90°時,由圖3 可知源的輻射場(即入射波)只存在于x≥0的半空間,而當α0°時源的輻射場只存在于z≥-dE的半空間;因此對于α0°—90°的活塞源而言,設(shè)界面反射波的接收點M位于x≥0 且-dE≤z≤0的空間.聲源設(shè)定與第4 節(jié)一致,流體和固體地層參數(shù)見表1,流體中聲波品質(zhì)因子為100,固體中縱波和橫波的品質(zhì)因子分別為80和60.

    5.1 活塞與界面的夾角的影響

    記過活塞軸線且與界面垂直的面為η面,則φM0°時M位于η面,φM90°時M位于與η面垂直的面.設(shè)活塞中心E位于z軸上z-dE-10 m處,反射波接收點M位于(rM,φM,zM)0,-10),其中rM,φM和zM的單位分別為米(m)、度(°)和米(m),此后默認為該單位.參考圖3 可知,M位于η面,并且βME1∠ME1O60°.地層取為快速地層,見表1.圖5 給出了不同夾角α下在M點處的反射波,其中圖5(a)和圖5(b)分別對應(yīng)a0.1 m和a1 m .由圖5(a)可知,小活塞源(a<λ)的輻射場的界面反射波隨α變化不大,波形隨α增加出現(xiàn)先增后減的現(xiàn)象,并且與活塞面和界面平行時(α0°)的結(jié)果相差不大.由圖5(b)可知,大活塞源 (a >λ) 的輻射場的界面反射波隨α顯著改變,波形隨著α的增加出現(xiàn)先增后減的現(xiàn)象,并且α30°和α90°下的波形相等.慢速地層下的結(jié)果類似,結(jié)論同上.接下來對圖5 的結(jié)果做出解釋,由于a<λ時入射波的指向性函數(shù)D隨θ的改變不太敏感,所以小活塞源的輻射場的指向性較弱,因而其界面反射波隨α的改變不大(見圖5(a)),而a >λ時情況相反(見圖5(b)).(44)式表明軸對稱指向性球面波的界面反射波可以視為鏡像活塞源激發(fā)的波(下文簡稱鏡像源的波)與反射系數(shù)B的乘積,由于圖5 中變化的位置參數(shù)只有α,則由(45)式和(2)式可知鏡像源的波受影響,而反射系數(shù)B不受影響,所以鏡像源的波隨α的變化規(guī)律即界面反射波隨α的變化規(guī)律.對于圖5 而言,α60°時鏡像活塞軸線經(jīng)過M點,得(見圖3),則D取得最大值1,此時M點處鏡像源的波取得最大值;由于圖5中不變,則從 0°開始增加α直至 90°的過程中,鏡像活塞軸線先接近M點再遠離M點,所以M點處鏡像源的波先增大后減小.因此,由鏡像源的波隨α的變化規(guī)律可知,圖5 中界面反射波隨α的增加而先增大后減小,并且在α60°處取得最大值.此外,圖5 中α30°和α90°下的鏡像活塞軸線關(guān)于ME1對稱,則α30°和α90°下的界面反射波也關(guān)于ME1對稱,又因M點接收的是聲壓,所以圖5 中α30°和α90°下的界面反射波相等.

    圖5 不同夾角α下在M 點(200,-10)處的反射波(a)a=0.1m;(b)a=1 mFig.5.Reflected wave under different included angleαfor pointMat (200,-10):(a)a=0.1 m;(b) a=1 m .

    前面算例表明,在給定M點 (rM,φM,zM)和E1點 (0,0,dE) 的情況下,改變α則鏡像活塞的軸線與M點的距離d改變,又因不隨α改變,所以最小的d值對應(yīng)鏡像源的波最大(見(44)式、(45)式和(2)式).因此,需要找到使d最小的α值,從而判斷鏡像源的波隨α的變化規(guī)律,繼而判斷界面反射波隨α的變化規(guī)律.由圖3 知鏡像活塞的軸線與M點的距離d等于,則由(45)式知

    在給定M點 (rM,φM,zM)和E1點 (0,0,dE) 的情況下,求d的導數(shù)零點(d對α的導數(shù))可得

    由于d的二階導數(shù)在αα*處大于零,所以αα*時鏡像活塞軸線與M點的距離d最小,此時M點處鏡像源的波最大.由于反射系數(shù)不隨α改變,則鏡像源的波隨α的變化規(guī)律即界面反射波隨α的變化規(guī)律,所以αα*時M點處的界面反射波最大.因此 0°≤α<α*時界面反射隨α增加而增大,α*≤α≤90°時界面反射隨α增加而減小.

    圖6 給出了不同夾角α下在M點45,-10)處(此時M點離開η面)的反射波,其中圖6(a)和圖6(b)分別對應(yīng)a0.1 m和a1 m .由圖6(a)可知,小活塞源 (a<λ) 的輻射場的界面反射波隨α變化較小,與活塞面和界面平行時(α0°)的結(jié)果相差較小.由圖6(b)可知大活塞源 (a >λ) 的輻射場的界面反射波隨α改變明顯,波形隨著α的增加出現(xiàn)先增后減的現(xiàn)象.慢速地層下的結(jié)果類似,結(jié)論同上.圖7 給出了不同夾角α下在M點90,-10)處(此時M點位于與η面垂直的平面)的反射波,其中圖7(a)和圖7(b)分別對應(yīng)a0.1 m和a1 m .由圖7(a)可知,小活塞源 (a<λ) 的輻射場的界面反射波幾乎不隨α變化,與活塞面和界面平行時(α0°)的結(jié)果幾乎沒有差別.由圖7(b)可知大活塞源 (a >λ) 的輻射場的界面反射波隨α改變較為明顯,波形隨著α的增加而減小,不同于圖6 中波形隨著α的增加而先增后減.慢速地層下的結(jié)果類似,結(jié)論同上.接下來對上述結(jié)果做出解釋,由(51)式知圖6 對應(yīng)α*50.8°,則界面反射波在αα*50.8°處最大,所以波形隨α的增加而先增后減;由(51)式知圖7對應(yīng)α*0°,則界面反射波在αα*0°處最大,所以波形隨α增加而減小.

    圖6 不同夾角α下在M點45,-10)處的反射波(a)a=0.1m;(b)a=1 mFig.6.Reflected wave under different included angleαfor point M at 45,-10):(a) a=0.1 m;(b) a=1 m .

    圖7 不同夾角α 下在M 點 90,-10)處的反射波 (a) a=0.1 m;(b)a=1 mFig.7.Reflected wave under different included angle α for point M at 90,-10):(a) a=0.1 m;(b) a=1 m.

    5.2 接收點的環(huán)向方位角的影響

    由于障板背面無輻射場存在,則有意義的方位角范圍為φM-90°—90°.由圖3 知活塞和界面均關(guān)于Oxz面對稱,又M點接收的是流體中的聲壓,所以φM0°—90°范圍的界面反射波與φM-90°—0°范圍的界面反射波對稱且相等.設(shè)活塞中心E位于 (0,0,-10),則E的鏡像E1位于 (0,0,10),設(shè)接收點M的位置為φM,-10),得40 m和60°,并且和不隨φM改變.地層取為快速地層,見表1.

    圖8 給出了α45°時不同φM下的界面反射波和其指向圖,其中指向圖是向φM0°處的波形的最大幅度歸一,圖8(a)對應(yīng)a0.1 m 時的反射波,圖8(b)對應(yīng)a0.1 m 時的指向圖,圖8(c)對應(yīng)a1 m時的反射波,圖8(d)對應(yīng)a1 m 時的指向圖.由圖8(a)和圖8(c)可知,反射波在φM0°處最大,并且隨著φM的增加而減小.因不隨φM改變,則反射系數(shù)也不隨φM改變(見(45)式),所以鏡像源的波隨φM的變化規(guī)律決定了反射波隨φM的變化規(guī)律;由于鏡像活塞的軸線位于φ0°面(即Oxz面,見圖3),且不變,則增加φM導致M點與鏡像活塞軸線的距離增加,所以鏡像源的波減小,因而反射波在φM0°處最大,并且隨φM的增加而減小.另外,比較圖8(a)和圖8(c)可知,小活塞情況下反射波隨φM的增加緩慢減小,而大活塞情況下反射波隨φM的增加迅速減小,所以大活塞情況下反射波對φM更為敏感.此外,比較圖8(b)和圖8(d)可知,大活塞情況下的反射波的指向性強于小活塞,這是因為在反射系數(shù)不變的情況下,鏡像活塞源的輻射場的指向性決定了反射波的指向性,而由(2)式可知大活塞源的輻射場的指向性強于小活塞源,那么其反射波的指向性也強于小活塞.

    圖8 α=45°時不同 φM下的反射波和其指向圖 (a) a=0.1 m時的反射波;(b) a=0.1 m時的指向圖;(c) a=1 m 時的反射波;(d) a=1 m 時的指向圖Fig.8.Reflected waves and directional diagrams at different azimuths φMfor α=45°:(a) Reflected waves for a=0.1 m ;(b) directional diagram for a=0.1 m;(c) reflected waves for a=1 m;(d) directional diagram for a=1 m .

    圖9 給出了不同α下的反射波的指向圖,其中圖9(a)和圖9(b)分別對應(yīng)a0.1 m和a1 m .由圖9(a)和圖9(b)可知有兩個特征存在:小活塞情況下反射波的指向性對α的改變不太敏感,而大活塞情況下反射波的指向性對α的改變很敏感;反射波的指向性隨α的增加出現(xiàn)先增后減的現(xiàn)象.第一個特征的出現(xiàn)是因為α的改變不影響反射系數(shù),則鏡像活塞源的輻射場的指向性決定了反射波的指向性,由于鏡像小活塞源的指向性較弱,所以反射波的指向性對α的改變不太敏感,由于鏡像大活塞源的指向性強,所以反射波的指向性對α的改變很敏感.第二個特征的出現(xiàn)是因為且φM0°時鏡像活塞的軸線通過M點,對應(yīng)鏡像源的波最大,又60°(由M點E1點的位置確定,但不隨φM改變),則α60°且φM0°時的鏡像源的波最大,所以在α60°下從φM0°開始增加或減小φM將導致鏡像源的波劇烈變化,因而α60°下反射波的指向性最強,則反射波的指向性隨α的增加出現(xiàn)先增后減的現(xiàn)象.

    圖9 不同夾角α 下的反射波的指向圖 (a) a=0.1 m;(b)a=1 mFig.9.Directional diagrams under different angle α:(a) a=0.1 m;(b) a=1 m .

    5.3 線段ME1 與界面法線的夾角的影響

    線段ME1與界面法線的夾角即,設(shè)活塞中心E位于z軸上z-dE-10 m 處,當M點位于 (7.279,φM,-10)時有20°,而當M點位于 (16.782,φM,-10)時有40°.圖10 給出了a1 m 時不同α下的反射波的指向圖,其中圖10(a)和圖10(b)分別對應(yīng)20°和40°.圖10(a)表明反射波的指向性隨著α的增加而先增后減,在α20°處達到最大值;圖10(b)表明反射波的指向性隨著α的增加而先增后減,在α40°處達到最大值;因此,界面反射波的指向性隨著α的增加而先增后減,在α處達到最大值.

    圖10 a=1 m時不同夾角α 下的指向圖 (a) =20°;(b)=40°Fig.10.Directional diagrams under different angle α for a=1 m:(a)=20°;(b) =40° .

    由(44)式和(45)式可知,βME1直接影響反射系數(shù)B,所以在給定鏡像源的波和地層參數(shù)的情況下,決定了界面反射波的大小.圖11 給出了反射系數(shù)的絕對值|B|隨的變化曲線(B的量綱為1),可以看出,快速地層下的反射系數(shù)大于慢速地層下的反射系數(shù),因而快速地層下的界面反射波總是大于慢速地層,而圖4 恰好說明了這一點.

    圖11 |B| 隨的變化曲線Fig.11.Variation curve of|B|with .

    另外,由(44)式可知,軸對稱指向性球面波的反射波中含有項,所以反射波的到時td等于鏡像活塞中心到接收點的距離除以流體中聲速vf,即td其中zM和rM分別為接收點的軸向和徑向坐標.因此,利用反射波的到時td可推算活塞源與界面的距離dE,而圖5—圖7 中反射波的到時都符合td的上述表達式.

    6 結(jié)論

    本文導出了無限大平面剛性障板中圓形活塞源(簡稱障板活塞源)輻射場的錐面波展開式(19),該式由多階錐面波構(gòu)成,適用于任意場點.前人在活塞中心到接收點的距離遠大于活塞半徑的條件下,將障板活塞源輻射場近似為軸對稱指向性球面波(1)式,而本文算例表明該條件下障板活塞源輻射場的錐面波展開式(19)近似為軸對稱指向性球面波的錐面波展開式.利用(19)式導出了活塞面與界面呈任意夾角時障板活塞源輻射場的界面反射波(28)式,并且當活塞中心的鏡像點到接收點的距離遠大活塞半徑時,(28)式近似為軸對稱指向性球面波的界面反射波表達式.(28)式相比前人只考慮活塞面與界面平行的情況而言,完善了無限大平面剛性障板中圓形活塞源的界面響應(yīng)理論.

    海洋地震勘探和水下目標超聲探測中,活塞中心的鏡像點到接收點的距離通常遠大于聲波波長,從而鞍點法的使用條件得到滿足.本文利用鞍點法和(1)式將軸對稱指向性球面波的反射波(28)式化簡為了(44)式.簡化式(44)適用于活塞中心的鏡像點到接收點的距離同時遠大于聲波波長和活塞半徑的情況.簡化式(44)數(shù)學形式簡潔,物理含義清楚,即軸對稱指向性球面波的反射波等于鏡像活塞源激發(fā)的軸對稱指向性球面波與反射系數(shù)的乘積.

    本文基于簡化式(44),通過數(shù)值算例研究了活塞面與界面的夾角α、接收點M的環(huán)向方位角φM和線段ME1與界面法線的夾角對軸對稱指向性球面波的界面反射波的影響.當活塞半徑小于流體中聲波波長時(a<λ),反射波對α不太敏感,當活塞半徑大于流體中聲波波長時(a >λ),反射波對α很敏感.增加活塞與界面的夾角α,反射波先增加后減小,而且α時取得最大值.因障板背面無波場存在,則接收點M有意義的方位角范圍為φM-90°—90°,又反射波場關(guān)于φM0°平面對稱,且M點接收的是聲壓,所以φM0°—90°與φM-90°—0°的反射波對稱且相等.在φM-90°—90°上,反射波隨方位角φM的增加而先增后減,在φM0°處取得最大值.a<λ時反射波對φM不太敏感,反射波的指向性較弱,a >λ時反射波對φM很敏感,反射波的指向性很強(相比a<λ而言).a<λ時反射波的指向性對α不太敏感,a >λ時反射波的指向性對α很敏感;增加活塞與界面的夾角α,反射波的指向性先增強后減小,當α時反射波的指向性最強.直接影響反射系數(shù),0°—90°時快速度地層的反射系數(shù)大于慢速地層的反射系數(shù),因而快速地層下的反射波大于慢速地層下的反射波.

    附錄A 鞍點法

    這里給出鞍點法的簡要介紹[23].對于積分式

    其中F(p) 和f(p) 是關(guān)于復數(shù)p的解析函數(shù),x是正實數(shù),C是積分路徑(通常為實軸積分路徑).當x是足夠大的正實數(shù)時,I(p)的值由 exp(xRef(p)) 控制,若將積分路徑C變換為最速下降路徑(即 Imf(p)=const 且經(jīng)過鞍點,而鞍點滿足f′(pan)=0),則I(p) 主要由鞍點附近的最速下降路徑積分決定,并可以寫為

    其中χ是鞍點附近的最速下降路徑與正 Re(p) 軸的夾角,設(shè)φ是f′′(pan)的相位角(φ∈(-π,π)),當積分路徑是從左往右走時則χ=反之則

    由(35)式可知,G*(kr)中只有Gn(kr) 可能是e 指數(shù)函數(shù).Gn(kr)是貝塞爾函數(shù)在橢圓面上的積分,與活塞面與界面的夾角α有關(guān),而與接收點的位置無關(guān).當活塞面與界面平行時α=0,由(18)式得Gn(kr)=0 (n≥1) 和G0(kr)=又因所以α=0 時Gn(kr) 不是e 指數(shù)函數(shù).因此我們推測當α=0時,Gn(kr) 也不是e 指數(shù)函數(shù),而數(shù)值計算表明確實如此,但這里不再給出詳細結(jié)果.綜上可知Gn(kr)不是e 指數(shù)函數(shù),所以G*(kr) 中不含e 指數(shù)函數(shù).

    猜你喜歡
    界面
    聲波在海底界面反射系數(shù)仿真計算分析
    微重力下兩相控溫型儲液器內(nèi)氣液界面仿真分析
    國企黨委前置研究的“四個界面”
    當代陜西(2020年13期)2020-08-24 08:22:02
    基于FANUC PICTURE的虛擬軸坐標顯示界面開發(fā)方法研究
    西門子Easy Screen對倒棱機床界面二次開發(fā)
    空間界面
    金秋(2017年4期)2017-06-07 08:22:16
    鐵電隧道結(jié)界面效應(yīng)與界面調(diào)控
    電子顯微打開材料界面世界之門
    人機交互界面發(fā)展趨勢研究
    手機界面中圖形符號的發(fā)展趨向
    新聞傳播(2015年11期)2015-07-18 11:15:04
    不卡视频在线观看欧美| 男的添女的下面高潮视频| 免费在线观看黄色视频的| 啦啦啦在线观看免费高清www| 成人18禁高潮啪啪吃奶动态图| 国产在视频线精品| 丝瓜视频免费看黄片| 亚洲,欧美精品.| 亚洲第一区二区三区不卡| 观看av在线不卡| 久久久精品94久久精品| 久久国内精品自在自线图片| 精品亚洲乱码少妇综合久久| 看十八女毛片水多多多| 免费人妻精品一区二区三区视频| 午夜激情久久久久久久| 国产精品久久久久久精品电影小说| 18禁观看日本| 免费大片黄手机在线观看| av网站免费在线观看视频| 亚洲一级一片aⅴ在线观看| xxxhd国产人妻xxx| 久久久久久久国产电影| 国产欧美亚洲国产| 日韩伦理黄色片| 黄色怎么调成土黄色| 亚洲五月色婷婷综合| 亚洲三级黄色毛片| 亚洲av日韩在线播放| 久久久久久伊人网av| 日本爱情动作片www.在线观看| 侵犯人妻中文字幕一二三四区| 丰满少妇做爰视频| 大片免费播放器 马上看| 亚洲精品,欧美精品| 91aial.com中文字幕在线观看| 成人黄色视频免费在线看| 亚洲精品色激情综合| 中文字幕最新亚洲高清| 最近中文字幕高清免费大全6| 午夜免费鲁丝| 好男人视频免费观看在线| 欧美 亚洲 国产 日韩一| 国产免费一区二区三区四区乱码| 一级,二级,三级黄色视频| 蜜臀久久99精品久久宅男| 欧美日韩成人在线一区二区| 日本色播在线视频| 五月伊人婷婷丁香| 桃花免费在线播放| 国产精品免费大片| 久久人人爽av亚洲精品天堂| 精品少妇久久久久久888优播| av国产精品久久久久影院| 最后的刺客免费高清国语| 2021少妇久久久久久久久久久| 久久久久精品性色| 97人妻天天添夜夜摸| 美女视频免费永久观看网站| 丝袜喷水一区| 久热这里只有精品99| 大香蕉久久网| 人人妻人人澡人人看| 亚洲成人av在线免费| 国产熟女午夜一区二区三区| 毛片一级片免费看久久久久| 18禁动态无遮挡网站| 国产毛片在线视频| 亚洲色图 男人天堂 中文字幕 | 亚洲精品国产色婷婷电影| 下体分泌物呈黄色| 午夜免费观看性视频| 美女大奶头黄色视频| 午夜激情久久久久久久| 97人妻天天添夜夜摸| 人妻 亚洲 视频| 1024视频免费在线观看| 天堂8中文在线网| 免费看av在线观看网站| 深夜精品福利| 如日韩欧美国产精品一区二区三区| 97在线视频观看| 女性被躁到高潮视频| 高清在线视频一区二区三区| 国产av码专区亚洲av| 国产麻豆69| 丝袜美足系列| 精品一区二区三卡| 免费人妻精品一区二区三区视频| 久久热在线av| 999精品在线视频| 欧美+日韩+精品| 边亲边吃奶的免费视频| 久久久久人妻精品一区果冻| 飞空精品影院首页| 一级毛片黄色毛片免费观看视频| 人成视频在线观看免费观看| 青春草亚洲视频在线观看| 黄网站色视频无遮挡免费观看| 欧美变态另类bdsm刘玥| 欧美精品一区二区免费开放| 精品少妇内射三级| av在线app专区| 亚洲欧美精品自产自拍| 国产精品偷伦视频观看了| 国产成人av激情在线播放| 久久久精品区二区三区| 免费av不卡在线播放| 99热网站在线观看| a级毛片黄视频| 日韩一区二区视频免费看| 男女下面插进去视频免费观看 | 少妇的逼水好多| 最新的欧美精品一区二区| 欧美性感艳星| 国产69精品久久久久777片| 国产成人一区二区在线| av福利片在线| 亚洲四区av| 又黄又爽又刺激的免费视频.| 插逼视频在线观看| 高清av免费在线| 国产精品99久久99久久久不卡 | 青春草国产在线视频| 99re6热这里在线精品视频| 欧美人与性动交α欧美精品济南到 | 亚洲av欧美aⅴ国产| www.熟女人妻精品国产 | 亚洲美女黄色视频免费看| 啦啦啦中文免费视频观看日本| 婷婷成人精品国产| 亚洲中文av在线| 亚洲综合色惰| av不卡在线播放| 日韩精品免费视频一区二区三区 | 99九九在线精品视频| 十八禁高潮呻吟视频| av.在线天堂| 又粗又硬又长又爽又黄的视频| 久久韩国三级中文字幕| 国产69精品久久久久777片| 久久久久精品人妻al黑| 777米奇影视久久| 97精品久久久久久久久久精品| 成人综合一区亚洲| 亚洲国产日韩一区二区| 精品午夜福利在线看| av卡一久久| 一个人免费看片子| 美女大奶头黄色视频| 国产综合精华液| 自线自在国产av| 日本免费在线观看一区| 欧美精品av麻豆av| 久久人人97超碰香蕉20202| 一级毛片黄色毛片免费观看视频| 亚洲综合色网址| 99热6这里只有精品| 免费大片黄手机在线观看| 日日撸夜夜添| 九色亚洲精品在线播放| 日韩,欧美,国产一区二区三区| freevideosex欧美| 亚洲综合精品二区| 亚洲在久久综合| 日本vs欧美在线观看视频| 欧美国产精品va在线观看不卡| 男女无遮挡免费网站观看| 黄片无遮挡物在线观看| 国产老妇伦熟女老妇高清| 蜜桃在线观看..| 中文字幕人妻熟女乱码| 欧美xxxx性猛交bbbb| 国产毛片在线视频| 午夜激情久久久久久久| 久久久久久人妻| 亚洲欧美一区二区三区黑人 | 下体分泌物呈黄色| 亚洲经典国产精华液单| 免费久久久久久久精品成人欧美视频 | 精品一区二区三卡| 视频区图区小说| 日本欧美国产在线视频| 久久久欧美国产精品| 亚洲天堂av无毛| 最近手机中文字幕大全| 日韩成人伦理影院| 桃花免费在线播放| 免费看av在线观看网站| 亚洲欧美精品自产自拍| 啦啦啦中文免费视频观看日本| 亚洲少妇的诱惑av| 最新的欧美精品一区二区| 久久韩国三级中文字幕| freevideosex欧美| 午夜91福利影院| 亚洲人成网站在线观看播放| 亚洲av成人精品一二三区| 下体分泌物呈黄色| 九九在线视频观看精品| 男女高潮啪啪啪动态图| 久久精品国产鲁丝片午夜精品| 亚洲国产精品999| 国产精品 国内视频| 精品久久蜜臀av无| a级毛片黄视频| 少妇精品久久久久久久| 亚洲欧美精品自产自拍| 国产精品人妻久久久影院| 精品少妇久久久久久888优播| 成年人午夜在线观看视频| 老熟女久久久| 看十八女毛片水多多多| 久久精品国产自在天天线| √禁漫天堂资源中文www| 婷婷色综合www| 亚洲伊人色综图| 久久久久国产网址| 2022亚洲国产成人精品| 天堂8中文在线网| 人妻系列 视频| h视频一区二区三区| 狠狠婷婷综合久久久久久88av| 91aial.com中文字幕在线观看| 亚洲综合色网址| 99久久精品国产国产毛片| 国产精品蜜桃在线观看| 欧美最新免费一区二区三区| 制服诱惑二区| 亚洲精品av麻豆狂野| 亚洲av综合色区一区| kizo精华| 国产 一区精品| 午夜久久久在线观看| 亚洲国产欧美日韩在线播放| 久久99一区二区三区| 日本爱情动作片www.在线观看| 精品亚洲成国产av| 99香蕉大伊视频| 美女国产高潮福利片在线看| videosex国产| 看免费av毛片| 亚洲精品一二三| 亚洲国产精品专区欧美| 国产乱来视频区| 免费女性裸体啪啪无遮挡网站| 欧美少妇被猛烈插入视频| 丝瓜视频免费看黄片| 午夜日本视频在线| 曰老女人黄片| 欧美日韩综合久久久久久| 九色成人免费人妻av| 九九爱精品视频在线观看| 国产成人精品久久久久久| 丰满迷人的少妇在线观看| 日本91视频免费播放| 亚洲色图 男人天堂 中文字幕 | 看十八女毛片水多多多| 在线观看国产h片| 又黄又爽又刺激的免费视频.| 一区二区三区精品91| 我要看黄色一级片免费的| 亚洲在久久综合| 少妇的逼好多水| 丝袜在线中文字幕| 国产黄频视频在线观看| 中文乱码字字幕精品一区二区三区| 精品少妇黑人巨大在线播放| 激情五月婷婷亚洲| 夜夜爽夜夜爽视频| 国产精品国产三级国产av玫瑰| 韩国高清视频一区二区三区| 22中文网久久字幕| 欧美日韩国产mv在线观看视频| 午夜视频国产福利| 中国三级夫妇交换| av线在线观看网站| 国产一区二区在线观看日韩| 亚洲精品成人av观看孕妇| 在线天堂最新版资源| 成人18禁高潮啪啪吃奶动态图| 精品人妻一区二区三区麻豆| 国产极品粉嫩免费观看在线| 另类精品久久| 22中文网久久字幕| 春色校园在线视频观看| 国产一区二区激情短视频 | 国产精品99久久99久久久不卡 | 18禁国产床啪视频网站| 国产国拍精品亚洲av在线观看| 亚洲国产欧美日韩在线播放| 黄片无遮挡物在线观看| 麻豆精品久久久久久蜜桃| 中文字幕人妻丝袜制服| 黄色毛片三级朝国网站| 国产探花极品一区二区| 久久久精品94久久精品| 日本欧美视频一区| 亚洲精品色激情综合| 九九爱精品视频在线观看| 亚洲av中文av极速乱| 欧美日韩国产mv在线观看视频| 久久av网站| 日韩av免费高清视频| 亚洲av在线观看美女高潮| 日本av免费视频播放| 国产成人aa在线观看| 狠狠婷婷综合久久久久久88av| xxxhd国产人妻xxx| 亚洲精品美女久久久久99蜜臀 | 另类亚洲欧美激情| 日韩av免费高清视频| 日本91视频免费播放| 国产国语露脸激情在线看| 丰满饥渴人妻一区二区三| a级毛片在线看网站| 99久久精品国产国产毛片| 看免费成人av毛片| 成人漫画全彩无遮挡| 久久国产精品男人的天堂亚洲 | 精品视频人人做人人爽| 2018国产大陆天天弄谢| 91成人精品电影| 成人国产麻豆网| 国产精品一二三区在线看| 亚洲av日韩在线播放| 最后的刺客免费高清国语| 捣出白浆h1v1| 两个人看的免费小视频| 亚洲精品一二三| 色哟哟·www| 男人舔女人的私密视频| 国产极品天堂在线| 美女内射精品一级片tv| 久久久a久久爽久久v久久| 男女免费视频国产| 校园人妻丝袜中文字幕| 精品亚洲成a人片在线观看| 男的添女的下面高潮视频| 在线观看美女被高潮喷水网站| 另类精品久久| 永久网站在线| 美女xxoo啪啪120秒动态图| 麻豆精品久久久久久蜜桃| 两个人免费观看高清视频| 亚洲熟女精品中文字幕| 男人爽女人下面视频在线观看| 日韩一区二区视频免费看| 亚洲精品,欧美精品| av国产久精品久网站免费入址| 国产一区有黄有色的免费视频| 欧美精品人与动牲交sv欧美| 久久国产亚洲av麻豆专区| 日韩电影二区| 日韩三级伦理在线观看| 免费黄频网站在线观看国产| 久久久久视频综合| 在线观看一区二区三区激情| 日韩三级伦理在线观看| 欧美精品高潮呻吟av久久| 最近最新中文字幕免费大全7| 三上悠亚av全集在线观看| 日本av免费视频播放| 超碰97精品在线观看| 黄片无遮挡物在线观看| 女的被弄到高潮叫床怎么办| 国产一区二区三区av在线| 建设人人有责人人尽责人人享有的| 国产成人精品无人区| 日本欧美国产在线视频| 国产av码专区亚洲av| 一区二区三区乱码不卡18| 精品国产一区二区三区久久久樱花| 91精品国产国语对白视频| 女人久久www免费人成看片| 国产av国产精品国产| 国产在线视频一区二区| 国产成人免费无遮挡视频| 伦理电影大哥的女人| 免费观看性生交大片5| 捣出白浆h1v1| 精品少妇内射三级| 中文字幕另类日韩欧美亚洲嫩草| 夜夜骑夜夜射夜夜干| 香蕉丝袜av| 婷婷色综合www| 欧美xxⅹ黑人| 国产黄色视频一区二区在线观看| 欧美精品人与动牲交sv欧美| 99热全是精品| 欧美日韩成人在线一区二区| 国产精品久久久久久久久免| 哪个播放器可以免费观看大片| 啦啦啦啦在线视频资源| 如日韩欧美国产精品一区二区三区| 黄网站色视频无遮挡免费观看| 视频在线观看一区二区三区| 久久精品久久久久久久性| 欧美精品一区二区免费开放| 狠狠精品人妻久久久久久综合| 日本av免费视频播放| 国产精品久久久久久久久免| av福利片在线| 涩涩av久久男人的天堂| 在线观看三级黄色| 久久久欧美国产精品| av一本久久久久| 少妇被粗大的猛进出69影院 | 午夜福利视频在线观看免费| 一级毛片 在线播放| 一级爰片在线观看| 男女高潮啪啪啪动态图| 在线精品无人区一区二区三| 国产成人免费无遮挡视频| 国产高清不卡午夜福利| 国语对白做爰xxxⅹ性视频网站| 色94色欧美一区二区| 亚洲婷婷狠狠爱综合网| 插逼视频在线观看| 99久国产av精品国产电影| 中文字幕另类日韩欧美亚洲嫩草| 免费日韩欧美在线观看| 天天躁夜夜躁狠狠久久av| 狂野欧美激情性bbbbbb| 国产精品无大码| 婷婷色av中文字幕| 精品国产一区二区久久| 免费少妇av软件| 国产精品免费大片| 高清视频免费观看一区二区| 91精品国产国语对白视频| 成人亚洲欧美一区二区av| 王馨瑶露胸无遮挡在线观看| 亚洲第一区二区三区不卡| 国产精品熟女久久久久浪| 精品人妻偷拍中文字幕| 在线免费观看不下载黄p国产| 日韩一本色道免费dvd| 一二三四在线观看免费中文在 | 国产老妇伦熟女老妇高清| 亚洲精品乱码久久久久久按摩| 亚洲五月色婷婷综合| 男人操女人黄网站| 亚洲国产精品国产精品| 成人国语在线视频| 成人午夜精彩视频在线观看| 两个人看的免费小视频| 人人妻人人澡人人爽人人夜夜| 精品少妇内射三级| 午夜福利网站1000一区二区三区| 水蜜桃什么品种好| 欧美人与性动交α欧美精品济南到 | 水蜜桃什么品种好| 男人操女人黄网站| 丁香六月天网| 欧美少妇被猛烈插入视频| 99re6热这里在线精品视频| 国产精品久久久久久av不卡| 香蕉丝袜av| 最近2019中文字幕mv第一页| 免费人妻精品一区二区三区视频| 国产成人欧美| 一级毛片 在线播放| h视频一区二区三区| 黑人欧美特级aaaaaa片| 日韩av免费高清视频| 十八禁高潮呻吟视频| 高清av免费在线| 国产伦理片在线播放av一区| 亚洲精品av麻豆狂野| 中文字幕另类日韩欧美亚洲嫩草| 中文乱码字字幕精品一区二区三区| 久久久久久久久久久久大奶| 另类精品久久| 桃花免费在线播放| 久久久欧美国产精品| 人人妻人人添人人爽欧美一区卜| 女性生殖器流出的白浆| 精品一区二区三区视频在线| 日韩人妻精品一区2区三区| 久久精品久久久久久久性| 草草在线视频免费看| 日本91视频免费播放| 伦理电影免费视频| 成人亚洲精品一区在线观看| 亚洲欧美中文字幕日韩二区| 免费人妻精品一区二区三区视频| 三上悠亚av全集在线观看| 日韩三级伦理在线观看| 另类亚洲欧美激情| av卡一久久| 日本-黄色视频高清免费观看| 久久久久精品性色| 波多野结衣一区麻豆| av黄色大香蕉| 国产亚洲最大av| 黄色配什么色好看| 婷婷色综合www| 2022亚洲国产成人精品| 亚洲国产精品国产精品| 国产男女内射视频| 久久久久国产网址| 欧美人与性动交α欧美软件 | 中文字幕制服av| 黑人高潮一二区| 成年女人在线观看亚洲视频| 色视频在线一区二区三区| 国产熟女欧美一区二区| 狂野欧美激情性bbbbbb| 国产国语露脸激情在线看| 日韩电影二区| 亚洲欧美日韩卡通动漫| 在线看a的网站| 999精品在线视频| 各种免费的搞黄视频| 国产精品人妻久久久久久| 亚洲情色 制服丝袜| 亚洲人成77777在线视频| 狂野欧美激情性bbbbbb| 一区二区av电影网| 亚洲精品久久成人aⅴ小说| 久久97久久精品| 免费观看在线日韩| 下体分泌物呈黄色| 国产欧美日韩综合在线一区二区| 涩涩av久久男人的天堂| 亚洲国产精品一区三区| 国产精品一区二区在线观看99| 亚洲,一卡二卡三卡| 精品酒店卫生间| 国产片特级美女逼逼视频| 精品一品国产午夜福利视频| 国产免费视频播放在线视频| 亚洲精品aⅴ在线观看| 亚洲国产欧美日韩在线播放| 欧美精品国产亚洲| 在线精品无人区一区二区三| 亚洲av免费高清在线观看| 18在线观看网站| 亚洲精品一区蜜桃| 中国国产av一级| 人人妻人人澡人人爽人人夜夜| av女优亚洲男人天堂| 久久99精品国语久久久| 久久久a久久爽久久v久久| 狂野欧美激情性xxxx在线观看| av有码第一页| 国产精品久久久久久精品电影小说| 日本av免费视频播放| 日本欧美视频一区| 免费观看无遮挡的男女| 日韩一本色道免费dvd| 天天躁夜夜躁狠狠躁躁| av一本久久久久| 99国产综合亚洲精品| 青春草国产在线视频| 各种免费的搞黄视频| 91精品三级在线观看| 在线 av 中文字幕| 日韩精品有码人妻一区| 在线观看免费视频网站a站| 久久综合国产亚洲精品| 精品一品国产午夜福利视频| 精品第一国产精品| 建设人人有责人人尽责人人享有的| 国产探花极品一区二区| av国产精品久久久久影院| 久久婷婷青草| 18禁裸乳无遮挡动漫免费视频| 国产色婷婷99| 亚洲国产欧美在线一区| 免费黄频网站在线观看国产| 中文精品一卡2卡3卡4更新| 国产不卡av网站在线观看| 亚洲精品av麻豆狂野| 肉色欧美久久久久久久蜜桃| 亚洲国产精品成人久久小说| 国产一区二区三区av在线| 久久精品国产亚洲av天美| 精品久久蜜臀av无| 国产日韩欧美亚洲二区| 啦啦啦在线观看免费高清www| 在线观看免费视频网站a站| 久久综合国产亚洲精品| 捣出白浆h1v1| 欧美日韩国产mv在线观看视频| 欧美日韩综合久久久久久| 人妻人人澡人人爽人人| 制服丝袜香蕉在线| 欧美成人午夜精品| 日韩中字成人| 丝瓜视频免费看黄片| 自拍欧美九色日韩亚洲蝌蚪91| 在线观看免费视频网站a站| 黑人欧美特级aaaaaa片| 90打野战视频偷拍视频| 我要看黄色一级片免费的| 1024视频免费在线观看| 欧美日韩亚洲高清精品| 亚洲成国产人片在线观看| 两个人免费观看高清视频| 久久精品久久精品一区二区三区| 久久精品人人爽人人爽视色| 成年人午夜在线观看视频| 街头女战士在线观看网站| 伊人久久国产一区二区| 久久久国产一区二区| 国产成人一区二区在线| 日日摸夜夜添夜夜爱| 波野结衣二区三区在线| 国产亚洲精品久久久com| 视频在线观看一区二区三区| 亚洲精品国产av蜜桃| 777米奇影视久久| 欧美人与善性xxx|