不同震源對地震波場的影響

不同震源對地震波場的影響

在彈性波場模擬中，源模擬是一個非常重要的方面。源加載方式和源模擬的質(zhì)量直接影響數(shù)值模擬結(jié)果。文獻(xiàn)利用傅氏變換法對脹縮源、剪切源、方向震源在均勻各向同性介質(zhì)中進(jìn)行了數(shù)值模擬;文獻(xiàn)對地震勘探中的炸藥震源藥量進(jìn)行了理論和實驗分析;文獻(xiàn)通過計算爆炸震源作用下的位移公式,分析了影響球腔震源子波的多個因素;本文從基于交錯網(wǎng)格高階有限差分的波動方程出發(fā),模擬了集中力源、脹縮源及剪切力源在各向異性介質(zhì)中的波場,通過模擬結(jié)果充分地說明震源的加載方式對數(shù)值模擬的重要性,同時該模擬結(jié)果對野外觀測系統(tǒng)中的震源研究、波場分析及橫波分裂的研究有很好的借鑒意義。1方法原理1.1階交錯網(wǎng)格速度—一階交錯網(wǎng)格速度——應(yīng)力波動方程有限差分法由于其算法簡便、對解決近、遠(yuǎn)場、復(fù)雜介質(zhì)結(jié)構(gòu)和邊界的波動問題具有較大的適應(yīng)性,故在地震波場模擬中得到廣泛的應(yīng)用。原理是根據(jù)波動方程的差分形式,計算離散網(wǎng)格點上的地震波場,這種算法有一個缺點是會產(chǎn)生數(shù)值頻散,而且傳統(tǒng)差分法的精度也較低。近10余年來發(fā)展了一種時間和空間差分精度可達(dá)任意階的高階交錯網(wǎng)格有限差分法,以其計算精度高、提高計算效率和節(jié)省存儲空間等優(yōu)點越來越受到人們的重視。通過Hooke定律、運動平衡方程和Cauchy公式等基本方程推導(dǎo)可得到ρ?Vx?t=?σxx?x+?τxz?z;(1a)ρ?Vz?t=?τxz?x+?σzz?z;(1b)?σxx?t=c11?Vx?x+c13?Vz?z;(1c)?σzz?t=c13?Vx?x+c33?Vz?z;(1d)?σxz?t=c44(?Vz?x+?Vx?z);(1e)ρ?Vy?t=?τxy?x+?τyz?z;(1f)?σyz?t=c44?Vy?z;(1g)?σxy?t=c66?Vy?x。(1h)式(1)為橫向各向同性介質(zhì)中二維三分量一階交錯網(wǎng)格速度——彈性波方程。Vx,Vy,Vz是位移速度分量,σxx,σzz,σxz,σyz,σyx是應(yīng)力分量,c11,c13,c33,c44,c66是彈性系數(shù)。用U,V,R,T,H,W,S,Q分別表示Vx,Vz,σxx,σzz,σxz,Vy,σyz,σyx對應(yīng)的離散值,k表示時間節(jié)點,i,j表示x,z方向網(wǎng)格節(jié)點,將位移速度分量Vx,Vy,Vz定義在離散時間k+12和k-12上,應(yīng)力分量σxx,σzz,σxz,σyz,σyx定義在離散時間k和k-1上,可得到各項精度為O(Δt4+Δx8)離散形式的差分格式Uk+12i,j=Uk-12i,j+Δtρi,j{1Δx4∑n=1C(4)n[Rki+(2n-1)2,j-Rki-(2n-1)2,j]}+Δtρi,j{1Δz4∑n=1C(4)n[Ηki,j+(2n-1)2-Ηki,j-(2n-1)2]}+Δt324ρ2i,jΡU；(2a)Vk+12i+12,j+12=Vk-12i+12,j+12+Δtρi+12,j+12{1Δx4∑n=1C(4)n[Ηki+n,j+12-Ηki-(n-1),j+12]}+Δtρi+12,j+12{1Δz4∑n=1C(4)n[Τki+12,j+n-Τki+12,j-(n-1)]}+Δt324ρ2i+12,j+12ΡV；(2b)Rki+12,j=Rk-1i+12,j+Δt*{C11Δx4∑n=1C(4)n[Uk-12i+n,j-Uk-12i-(n-1),j]}+Δt*{C12Δz4∑n=1C(4)n[Vk-12i+12,j+(2n-1)2-Vk-12i+12,j-(2n-1)2]}+Δt324ρi+12,jΡR；(2c)Τki+12,j=Τk-1i+12,j+Δt{C12Δx4∑n=1C(4)n[Uk-12i+n,j-Uk-12i-(n-1),j]}+Δt{C22Δz4∑n=1C(4)n[Vk-12i+12,j+(2n-1)2-Vk-12i+12,j-(2n-1)2]}+Δt324ρi+12,jΡΤ；(2d)Ηki,j+12=Ηk-1i,j+12+Δt{C55Δx4∑n=1C4n[Vk-12i+(2n-1)2,j+12-Vk-12i-(2n-1)2,j+12]}+Δt{C55Δz4∑n=1C(4)n[Uk-12i,j+n-Uk-12i,j-(n-1)]}+Δt324ρi,j+12ΡΗ；(2e)Wi,jk+12=Wi,jk-12+Δtρi,j{1Δx∑n=14Cn(4)[Qi+(2n-1)2,jk-Qi-(2n-1)2,jk]}+Δtρi,j{1Δz∑n=14Cn(4)[Si,j+(2n-1)2k-Si,j-(2n-1)2k]}+Δt324ρi,j2ΡW；(2f)Si,j+12k=Si,j+12k-1+Δt*{C44Δz∑n=14Cn(4)[Wi,j+nk-12-Wi,j-(n-1)k-12]}+Δt324ρi,j+12ΡS；(2g)Qi+12,jk=Qi+12,jk-1+Δt{C55Δx∑n=14Cn(4)[Wi+n,jk-12-Wi-(n-1),jk-12]}+Δt324ρi+12,jΡQ(2h)其中PU,PV,PR,PT,PH,PW,PS,PQ均是空間導(dǎo)數(shù)的高階項,考慮到文章的篇幅,就不詳細(xì)地列出。采用速度——應(yīng)力交錯網(wǎng)格的離散剖分后,方程(2)中應(yīng)力和位移速度的空間節(jié)點位置如圖1所示;應(yīng)力在時間的整數(shù)節(jié)點上計算,位移速度在時間的半數(shù)節(jié)點上計算。采用交錯網(wǎng)格的優(yōu)點是計算應(yīng)力(位移速度)時,用到應(yīng)力節(jié)點兩邊半個網(wǎng)格節(jié)點上的位移速度(應(yīng)力)進(jìn)行計算,計算精度要比規(guī)則網(wǎng)格高。1.2吸收邊界條件用一階彈性波方程的交錯網(wǎng)格高階差分解法模擬地震波場,精度雖然提高了,但同傳統(tǒng)差分解法一樣,應(yīng)考慮因數(shù)值參數(shù)選擇不合理而引起的頻散問題,文中穩(wěn)定性參考文獻(xiàn),這里不作詳述。進(jìn)行數(shù)值模擬時,往往引進(jìn)人工邊界以限制計算范圍。對于彈性波方程來講,就是在邊界上人工反射應(yīng)盡可能地小,使之對區(qū)域內(nèi)部解的影響在允許誤差的范圍之內(nèi)。文獻(xiàn)提出一種非常有效的吸收邊界條件—完全匹配吸收邊界,即PML(Perfectlymatchedlayer)吸收邊界,這種方法目前已被廣泛地應(yīng)用于有限差分和有限元求電磁波方程和彈性動力學(xué)方程的數(shù)值模擬,但在有些情況下,PML吸收邊界條件卻不盡人意,模擬時會出現(xiàn)數(shù)值不穩(wěn)定的現(xiàn)象。本文引入了MPML(Multi-axialPML)邊界條件,使得在兩個正交方向用互成比例的衰減系數(shù)來衰減匹配層內(nèi)的波,獲得了很好的效果。傳統(tǒng)的分裂PML的衰減系數(shù)(以右邊界為例)dx=dx(x),dy=0,dz=0。衰減只是x的函數(shù),只能對垂直于邊界方向的地震波進(jìn)行衰減。MPML則可對各個方向的地震波進(jìn)行衰減,對應(yīng)的衰減系數(shù)為dx=dx(x)(x);dy=Ρyxdx(x)(x);dz=Ρzxdx(x)(x)。(3)式(3)中衰減系數(shù)dx(x)(x)=(n+1)V2δ?ln(1R),其中:n為控制衰減快慢的參數(shù),本文取n=3;V為波速;δ為吸收層厚;R為反射系數(shù);Ρzx是MPML介質(zhì)穩(wěn)定性的微調(diào)因子。2種介質(zhì)中的波場震源模擬就是在差分網(wǎng)格上施加力源,為時間層和遞推提供初始能量。從震源的加載方式講,可以是模擬錘擊的集中力源和炸藥源的爆炸源,也可以是純剪切力源。下面詳細(xì)討論這3種震源在各向異性介質(zhì)中的波場。文中震源模擬的模型大小為500m×500m,主頻為40Hz,震源子波為ricker子波,水平方向和垂直方向的網(wǎng)格大小均2m,震源的有效半徑為10m,數(shù)據(jù)的采樣間隔為0.2ms。橫向各向異性介質(zhì)參數(shù)為:ρ=2150m/s,C11=25.24GPa,C13=9.85GPa,C33=12.9GPa,C44=3.62GPa,C44=9.94GPa,震源坐標(biāo)為(250m,250m)。2.1向有角度下的qp波和準(zhǔn)s波集中力源模擬是差分網(wǎng)格某個結(jié)點上施加一個隨時間變化的作用力,力的作用可以是水平力、垂直力,也可以是傾斜力,所以在彈性介質(zhì)中激發(fā)的彈性波場也具有方向效應(yīng)。式(4)中只有ρfx體力分量時為水平集中力,只有ρfz體力分量時為垂向集中力,若沿水平或垂向有夾角時則為傾斜集中力ρ?Vx?t=?σxx?x+?τxz?z+ρfx;(4a)ρ?Vy?t=?τxy?x+?τyz?z+ρfy;(4b)ρ?Vz?t=?τxz?x+?σzz?z+ρfz。(4c)從圖2可看出集中力源在各向異性介質(zhì)中不僅激發(fā)P波(縱波),也激發(fā)S波(橫波)。由于各向異性介質(zhì)中波的傳播方向不一定和極化方向一致或垂直,所以稱為準(zhǔn)P波和準(zhǔn)S波。波前圖不再如各向同性介質(zhì)中為兩個同心圓,而像個橢圓或似橢圓。x分量圖沿著力源方向和垂直于力源方向qP波和qS波能量均為零,沿著與力源方向斜交的方向上qP波和qS波均有較強(qiáng)的能量出現(xiàn);z分量圖沿著力源方向qS波能量在z方向為零,x方向能量最強(qiáng),qP波能量在x方向為零,z方向能量最強(qiáng)。與各向同性介質(zhì)中激發(fā)的能量分配方式基本相似,但強(qiáng)弱不同。2.2水平和垂直方向應(yīng)力模擬亦稱為爆炸震源,其作用效果相當(dāng)球形空腔源,該震源產(chǎn)生的是無旋場,這種震源在各向同性介質(zhì)中只產(chǎn)生縱波,常用來模擬縱波源;本文通過在網(wǎng)格點處采用同時加水平方向應(yīng)力和垂直方向應(yīng)力來模擬脹縮源,而不是賦在位移處,這樣就很好地避開了震源處的無窮速度問題。由圖3可明顯地看出,在各向異性介質(zhì)中加載脹縮震源,既產(chǎn)生縱波也產(chǎn)生橫波,橫波的能量相對較弱些。x分量的qP波能量在x方向最強(qiáng),qS波能量主要集中在與x方向斜交的方向上;z分量的qP波能量在z方向最強(qiáng),橫波能量也集中在與x方向斜交的方向上,但斜交的角度與x分量的橫波顯然不同。2.3分量方程簡化計算剪切力源在各向同性介質(zhì)中只產(chǎn)生橫波,常用來模擬橫波震源。自20世紀(jì)80年代以來,地球物理學(xué)家就已經(jīng)通過檢測地震資料中的橫波分裂現(xiàn)象提取裂縫參數(shù)。由于橫波激發(fā)地震比縱波激發(fā)地震復(fù)雜而且耗資較高,至今沒有得到廣泛的應(yīng)用。隨著油氣田勘探開發(fā)的需要,提高采樣率迫切需要橫波技術(shù)的推進(jìn)。本文二維三分量方程在加載了剪切震源簡化后的公式為ρ?Vx?t=?σxx?x+??z(τxz-ρSy);(5a)ρ?Vy?t=??x(τxy-ρSz)+??z(τyz+ρSx);(5b)ρ?Vz?t=??x(τxz+ρSy)+?σzz?z。(5c)其中S為震源函數(shù)。圖4中波場顯然表明剪切力源在各向異性介質(zhì)中既激發(fā)qP波也激發(fā)qS波,但qP波能量很弱。x分量的qP波能量沿著x方向最強(qiáng),沿著z方向最弱,qS波能量在z方向最強(qiáng),而在x方向為零。Z分量的qP波和qS波能量分配與x分量的相反。3震源在各向異性介質(zhì)中的釋放特性通過用交錯網(wǎng)格高階有限差法對地震波震源的數(shù)值模擬,得出具有方向性的垂直集中力源在各向異性介質(zhì)中均可激發(fā)出qP波和qS波,說明能量分配和偏振方向會隨集中力的方向有明顯地變化。脹縮震源在各向異性介質(zhì)中除了產(chǎn)生qP波也產(chǎn)生qS波,但