基于鞍點法與互易性的遠探測波場模擬?

許家旗 胡恒山

(哈爾濱工業大學力學系 哈爾濱 150001)

0 引言

聲波遠探測技術是近年來發展較快的技術，可以對井外范圍數十米的裂縫、斷層、地質界面進行探測，具有巨大的應用前景。多年來，國內外對聲波遠探測方法開展了一系列研究，Hornby[1]在1989年從陣列波數據中提取反射波，并利用偏移的方法得到了井旁地層結構變化的圖像。國內薛梅[2]、喬文孝等[3]以及何峰江[4]針對單極源條件下的反射波成像問題展開了一系列研究。但是，單極子聲波遠探測無法識別反射體的方位，且單極子源較高的聲源頻率造成其探測范圍有限。針對單極子聲源存在的不足與問題，Tang[5]首次提出了用偶極子聲源來進行聲波遠探測。隨后，Tang等[6]利用四分量偶極數據對井旁裂縫和鹽丘內部的精細構造進行了成像。Wei等[7]采用有限差分模擬了偶極源遠探測中彈性波的輻射、反射以及接收系數，發現SH波有較大的輻射范圍以及較強的接收敏感度，適合反射波成像；在慢速地層中，P波也有較強的接收敏感度，也可用來反射波成像。在聲波遠探測中，快速準確地進行波場的正演模擬對于成像具有重要的意義。而當井外存在反射體時，結構是非軸對稱的，很難獲得解析解，通常采用有限差分進行數值模擬[7?9]。在有限差分中，需要設置較小的網格尺寸來模擬井孔，當反射界面距離井孔較遠時要耗費大量的時間。為了滿足現場測井數據實時處理的需求，快速準確地進行波場計算具有重要的意義。

聲波遠探測中井孔的存在對彈性波的輻射與接收的影響不可忽略。對于彈性波的輻射，Lee等[10]計算了低頻時井內點源激發的遠場波場。Gibson[11]計算了聲波測井頻率范圍內井內聲源激發的遠場波場。對于彈性波的接收，Schoenberg[12]計算了低頻時井外入射波在井內激發的波場，Lee[13]研究了井外平面P波、SV波和SH波斜入射時井壁上的位移分布。Lovell等[14]研究了井外平面波以不同角度入射時引起的井內聲壓。Peng[15]在Schoenberg理論的基礎上計算了聲波測井頻率范圍內井外入射波在井內激發的位移與位移場。但是采用Schoenberg理論計算井外入射波激發的井內響應是復雜的，Tang等[16]利用互易性將井外入射波激發的井內響應問題轉為井內聲源的遠場輻射問題，給出了偶極源聲波遠探測中反射波位移的漸近解可以快速模擬聲波遠探測問題。Wang等[9]給出了流固分層介質中多極源之間的互易關系，本文將Wang等的互易關系應用在聲波遠探測中，旨在獲得聲波遠探測中井內波場的漸近解。

本文首先通過鞍點法(也稱最速下降積分法)獲得井內單極源和偶極源遠場輻射場波數域積分的漸近表達式，并將漸近解與實軸積分獲得的精確解進行比較，驗證漸近解的準確性。井內聲源的輻射波遇到界面發生反射，界面反射波可以等效為集中力的輻射波，將反射波激發井內波場問題轉換為井外集中力激發井內波場問題。最后，利用井外集中力源與井內聲源之間的互易關系獲得了反射波場激發的井內響應的漸近解。本文的計算流程如圖1所示。

圖1 計算流程圖Fig.1 Computational flow chart

1 井中聲源的輻射場

圖2 井中聲源輻射示意圖Fig.2 The model of radiation of the source in the borehole

在本節中首先考慮充液井孔中的井軸處放置單極源或者偶極源，如圖2所示，井孔半徑為a，R是輻射波傳播的距離，φ是方位角，θ是輻射波相對于井軸的傾角。聲源用V0G(t)來描述，其中V0是體積變化量，G(t)是聲源的時間函數。在本研究里，傅里葉變換表達如下：

流體中位移可表示為

其中，?d是直達波的位移勢，?r是井內反射波的位移勢，?是柱坐標系下的拉普拉斯算子。

地層中的位移可表示為

其中，?s是地層P波位移勢函數，χ是地層SH波位移勢函數，Γ是地層SV波位移勢函數。

井內的位移勢函數可表示為[17]

井外的位移勢函數可表示為

其中，f2=k2?ω2/α2f，p2=k2?ω2/α2，s2=k2? ω2/β2。上標n=0代表單極源，n=1 代表偶極源；In和Kn分別是第一類和第二類貝塞爾函數；εn是紐曼因子，當n=0時εn=1，其余情況εn=2；r0是偶極源的半徑；k、f、p和s分別是軸向波數、流體徑向波數、縱波徑向波數、橫波徑向波數；αf、α和β分別是流體聲速、地層的縱波與橫波波速。系數A(n)、B(n)、C(1)和D(n)可以通過井壁上的邊界條件來獲得[18]。

1.1 最速積分法獲得遠場輻射場的漸近解

地層中的勢函數涉及到復雜的波數域積分，傳統的計算方法是采用實軸積分方法進行波數域積分的計算[19]，該方法計算較為繁瑣。當反射界面離井孔位置較遠時，可以采用最速下降積分法獲得井孔中聲源遠場輻射場的漸近表達式[20]：

1.2 最速積分法的漸近解與實軸積分法的比較

為了驗證遠場輻射波的漸近解的準確性，將采用式(6)得到的漸近結果與實軸積分計算的結果進行比較，計算采用的參數由表1給出。地層分別為慢速地層和快速地層2，井孔半徑為0.1 m，井外觀測點的位置為r=5 m,z=1 m,φ=0?，聲源為中心頻率3 kHz、半帶寬為2 kHz的余弦包絡脈沖函數。圖3給出了井外地層為快速地層時漸近解與實軸積分精確解的對比結果，圖4給出了井外地層為慢速地層時漸近解與實軸積分精確解的對比結果。從圖3和圖4中可以看出，無論是慢速地層還是快速地層，采用最速下降積分法得到的漸近結果與實軸積分的精確結果符合較好，這也說明了在遠場位置，漸近解是準確可信的，且鞍點的計算與頻率無關[6]。

表1 井中流體以及地層參數Table 1 The parameters of the borehole fluid and formations

將式(6)代入式(3)，并忽略高階小項，得到井內單極源的遠場輻射位移表達式：

井內偶極源的遠場輻射位移表達式為

圖3 快速地層條件下采用實軸積分與最速下降積分計算的井外輻射波Fig.3 The radiated wave calculated using real-axis integration and steepest descent integration in the fast formation

圖4 慢速地層條件下采用實軸積分與最速下降積分計算的井外輻射波Fig.4 The radiated wave calculated using real-axis integration and steepest descent integration in the slow formation

其中，

2 反射場的接收

下面研究聲源為偶極源時，井內接收到來自井外界面的反射波場。由式(8)可以看出，井內偶極源的遠場輻射波具有球面的形式，當遇到界面發生反射。當反射界面的尺寸遠大于波長時，可以將界面當成一個無限大平面。當反射界面距離聲源較遠時，可以將球面波入射可以視為平面波入射，球面波反射系數可以用平面波反射系數代替[20]，反射波可以寫作

圖5 井外存在反射體時的反射波示意圖Fig.5 The reflected wave from the reflector outside the borehole

其中，GSH、GSV和GP分別是SH波、SV波和P波的反射系數。反射波uφ、uθ和uR的偏振方向如圖5所示。

從界面反射回來的波由井中的接收器接收，在遠探測的聲波頻率范圍內，井孔的存在對波場的影響不可忽略。采用Schoenberg理論[12]計算井內的波場比較復雜，在本文中，在已知井內輻射場的基礎上，利用互易定理求解反射波在井內激發的位移以及聲壓。

2.1 反射波場的等效虛源

無限大均勻彈性介質中，集中力在該力的作用線上的位置激發且只激發出球面縱波[21]，

其中，F(ω)為集中力的大小，e為集中力的方向，是觀測位置與集中力之間的距離。在垂直于該力方向并經過其作用點的直線上集中力激發且只激發出球面橫波us，

將反射波式(9)代入式(10)和式(11)，可以得到與反射波等效的集中力源：

式(12)中，P0R0表示在距離聲源R0的壓強為P0，T(ω)為聲壓源的頻譜函數[22]。如圖6所示，等效集中力FR、Fφ和Fθ作用的位置與聲源關于反射界面對稱，其中FR沿反射P波的偏振方向，Fφ沿反射SH波的偏振方向，Fθ沿反射SV波的偏振方向。這樣就將求解反射波激發的井內響應問題轉化為求解井外虛源激發井內響應問題，虛源由三個相互正交的集中力構成。注意到在圖6中，在反射波等效為虛源的輻射波后，反射界面也不再存在了。

圖6 反射波的等效集中力Fig.6 The concentrated force equivalent to the reflected wave

2.2 利用互易性獲得井中反射波場的響應

為了獲得井外集中力激發的井內響應，利用互易定理求解。Wang等[9]建立了多極聲源以及集中力源在流固分層介質中的互易關系。首先求井外集中力激發的井內聲壓，如圖7所示。

圖7 井中接收聲壓信號示意圖Fig.7 The reception of pressure in the borehole

為了獲得井中O2處的聲壓，在O2處放置一個單極源膨脹點源，其在井外遠場O1處的輻射位移ui(i=R,θ,φ)由式(7)給出。單極源與集中力之間的互易關系為[9]

將式(7)和式(12)代入式(13)，得到反射波激發的井內聲壓p，

其中，

L是波傳播的總距離，θ1是輻射波傳播方向相對于井軸的傾角，θ2是反射波傳播方向相對于井軸的傾角。

類似地，如圖8所示，為獲得井外集中力激發的沿x方向井內位移，在接收位置O2處放置一個沿x方向的偶極源D，其在井外O1處的輻射位移ui(i=R,θ,φ)由式(8)給出。單極源與集中力之間的互易關系為[9]

其中，a=?ω2u′是O2處沿x方向的加速度。將式(8)和式(12)代入式(15)，獲得反射波激發的井內沿x方向井內位移，

其中，

圖8 井中接收位移信號示意圖Fig.8 The reception of displacement in the borehole

為了驗證漸近解的準確性，采用三維有限差分模擬聲波遠探測中反射波激發的聲壓與位移，如圖9所示，反射體距離聲源3 m，傾角為15?，計算參數在表1給出，聲源的中心頻率為3 kHz，聲源為偶極源并可以沿井軸旋轉。

當偶極源方向平行于反射界面，接收沿x方向的位移，采用漸近解與數值解的對比結果如圖10所示。從圖10中可以看出，漸近解模擬的反射SH波與有限差分吻合較好，當偶極源方向平行于反射體走向時，只有SH波的位移且SH波位移達到最大，除此之外，全波中最先達到的是彎曲波，幅度在全波中最大。

圖9 偶極聲波遠探測示意圖Fig.9 The model of dipole single-well imaging

當偶極源方向指向反射界面，接收聲壓，采用漸近解與數值解的對比結果如圖11所示。從圖11中可以看出，漸近解與數值解的結果吻合較好。漸近解與數值解之間的差異來自于P-SV轉換的貢獻，隨著收發距離的增大，轉換波的貢獻隨之增大，這是因為收發距離增大時入射角也變大，導致反射波的P-SV轉換系數增大。顯然，當反射界面距離聲源較遠時，入射角變小，轉換波減小，漸近解與數值解之間的誤差也會變小。值得注意的是，當采用偶極聲源發射，單極接收的收發方式，接收的全波中沒有井孔彎曲波，這也是與Kurkjian等[17]的結論是一致的。這種混合模式的收發方式，理論上可以達到消除井孔波而凸顯出反射波信號的目的，這將有利于后續的信號提取以及反演。另外，從式(14)可以看出，當反射體方位角改變180?時，接收的聲壓信號的極性發生偏轉。目前遠探測中存在的反射體方位角180?不確定性問題[5?6]，而本文的研究發現：采用偶極發射-單極接收的收發模式，當方位角改變180?時，接收的信號極性發生偏轉，這可以用來消除反射體方位角不確定性來唯一確定反射體方位角。

圖10 井內位移漸近解與數值解之間的對比Fig.10 Comparison between asymptotic solutions to horizontal displacement at the borehole axis and exact numerical results to displacement

圖11 井內聲壓漸近解與數值解之間的對比Fig.11 Comparison between asymptotic solutions to pressure at the borehole axis and exact numerical results to pressure

通過與有限差分結果的比較，驗證了漸近解的正確性。而采用有限差分要耗費大量的時間，同時對計算機的性能要求很高，而采用漸近解可以快速獲得波場的正確結果，對于實際聲波遠探測中數據實時處理具有重要的意義。

3 結論

本文采用解析法給出了求解聲波遠探測中非軸對稱波場的具體方法，獲得了以下結論和成果：

(1)采用鞍點法計算了充液井孔中偶極源和單極源激發的遠場輻射波場的漸近解，并與采用實軸積分法計算的精確結果相一致，證明了采用鞍點法計算井外輻射場的準確性。

(2)當井外存在傾向的反射界面時，來自界面的反射波可以等效為三個集中力源的輻射波，利用井內單極源與井外集中力之間的互易關系獲得了井內聲壓信號，利用井內偶極源與井外集中力之間的互易關系獲得了井內位移信號。此解析法計算結果與有限差分的結果一致，且解析法顯著提高了計算效率，對于實際測井中數據的實時處理具有重要意義。

(3)傳統的聲波遠探測儀器采用同一類型的聲源與接收器，例如聲源和接收器都為單極或者偶極，這種收發模式接收到的波場信號中含有較強的井孔模式波且幅度較大，不利于反射波的提取。采用偶極源發射，居中單極接收的收發模式，理論上可以消除井孔模式波來凸顯反射波，可提高實際測井數據實時處理的效率。

(4)在偶極發射-單極接收的收發模式下，當反射體方位角改變180?時，接收的反射波波形極性會發生偏轉，可以用來消除目前聲波遠探測中反射體方位角多解性問題來唯一確定反射體方位角。