• 
    

    
    

      99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

      負(fù)荷勒夫數(shù)漸近值解析表達(dá)式

      2021-10-20 06:33:50周江存潘爾年
      地球物理學(xué)報(bào) 2021年10期
      關(guān)鍵詞:勒夫格林引力

      周江存,潘爾年

      1 中國(guó)科學(xué)院精密測(cè)量科學(xué)與技術(shù)創(chuàng)新研究院,大地測(cè)量與地球動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,武漢 430077 2 College of Engineering, Cleveland State University, Cleveland, OH 44115, USA

      0 引言

      全球水循環(huán)是地球表層質(zhì)量遷移的重要形式,它引起全球水質(zhì)量在時(shí)間和空間的重新分布,會(huì)進(jìn)一步導(dǎo)致地球的變形.該變形主要有三種形式:(1)大氣負(fù)荷,主要是大氣質(zhì)量的變化對(duì)地球的引力作用與壓強(qiáng)變化對(duì)地表的負(fù)載作用(Merriam,1992;Van Dam et al.,1994;Sun and Luo,1998;Guo et al.,2004);(2)潮汐和非潮汐海洋負(fù)荷,主要是海水質(zhì)量的時(shí)空變化對(duì)地球的引力作用和對(duì)地表的負(fù)載作用(Longman,1962,1963;Farrell,1972;Goad,1980;許厚澤等,1982);(3)陸地水負(fù)荷,主要是陸地水的質(zhì)量變化對(duì)地球的引力作用和對(duì)地表的負(fù)載作用(Zhou et al.,2009).這三種效應(yīng)可以引起地表幾厘米(cm)的位移和數(shù)十微伽(10-8m·s-2)的重力變化.在目前高精度的大地測(cè)量學(xué)研究中,這三種效應(yīng)都是必須要加以考慮的.對(duì)這些效應(yīng)的模擬,目前已有比較成熟的理論,通常采用負(fù)荷質(zhì)量與負(fù)荷格林函數(shù)的褶積積分計(jì)算獲得,其中格林函數(shù)是負(fù)荷勒夫數(shù)或其組合與Legendre函數(shù)(或稱為球函數(shù))或其導(dǎo)數(shù)的乘積的無(wú)窮級(jí)數(shù)求和(Longman,1963;Farrell,1972;Guo et al.,2004).鑒于球函數(shù)有解析的遞推公式,因此計(jì)算格林函數(shù)的基礎(chǔ)就是計(jì)算出負(fù)荷勒夫數(shù).

      由于格林函數(shù)是一個(gè)無(wú)窮級(jí)數(shù)的求和,對(duì)于一個(gè)無(wú)窮級(jí)數(shù)求和的問(wèn)題來(lái)說(shuō),高階勒夫數(shù)的高效計(jì)算以及格林函數(shù)的收斂性如何是數(shù)值計(jì)算的關(guān)鍵.為了計(jì)算高階的負(fù)荷勒夫數(shù),Longman(1963)提出了對(duì)變量進(jìn)行無(wú)量綱化的處理;汪漢勝等(1996)通過(guò)觀察各變量的變化特征提出了r冪因子法;特別地,Pan等(2015)基于合理的假設(shè)提出了計(jì)算負(fù)荷勒夫數(shù)的解析解方法;Chen等(2018)基于解析解方法提出了一種稱為剛度矩陣的傳播矩陣方法來(lái)計(jì)算超高階負(fù)荷勒夫數(shù),克服傳統(tǒng)傳播矩陣方法的溢出問(wèn)題;Zhou等(2019a)運(yùn)用DVP(dual variable and position)傳播矩陣法獲得了任意球諧階數(shù)的勒夫數(shù),并提供了相應(yīng)的MATLAB程序計(jì)算不同類型的勒夫數(shù)(包括潮汐、負(fù)荷、剪切和位錯(cuò)).實(shí)際上,這些表征地球在不同力源(如引潮力、地表正應(yīng)力或剪切應(yīng)力、以及內(nèi)部斷層的等效體力等)作用下產(chǎn)生變形的勒夫數(shù)之間并不完全獨(dú)立,它們之間的關(guān)系可參考潘爾年和丁中一(1986)和Okubo(1993).

      為了保證格林函數(shù)的收斂,F(xiàn)arrell(1972)給出了負(fù)荷勒夫數(shù)的漸近值,即當(dāng)球諧階數(shù)n趨于無(wú)窮大時(shí)勒夫數(shù)的值,然后通過(guò)Kummer變換可有效降低截?cái)嚯A數(shù)(約10000階)并加速格林函數(shù)的收斂,從而為能夠準(zhǔn)確地模擬負(fù)荷引起的地球形變奠定了基礎(chǔ).可見(jiàn),負(fù)荷勒夫數(shù)漸近值對(duì)負(fù)荷格林函數(shù)來(lái)說(shuō)是至關(guān)重要的.除了Farrell(1972)之外,Okubo和Endo(1985)、Guo等(2004)也對(duì)如何計(jì)算負(fù)荷勒夫數(shù)漸近值給出了不一樣的理論計(jì)算方法.我們最近提出了計(jì)算位錯(cuò)勒夫數(shù)的漸近值的方法,在忽略自引力的情形下,漸近值可根據(jù)傳播矩陣的特性很容易地求出(Zhou et al.,2021a).位錯(cuò)問(wèn)題和負(fù)荷問(wèn)題的差別在于邊界條件的不同,因此該方法可直接應(yīng)用于求取負(fù)荷勒夫數(shù)的漸近值.

      本文首先對(duì)目前已有的上述三種計(jì)算漸近值的方法進(jìn)行簡(jiǎn)單的回顧,然后著重介紹一種新的能夠獲取負(fù)荷勒夫數(shù)漸近值解析表達(dá)式的方法.最后將我們的方法與已有的計(jì)算方法進(jìn)行了比較,并給出了定量的數(shù)值結(jié)果展示新的解析表達(dá)式在表達(dá)負(fù)荷勒夫數(shù)漸近特征方面的優(yōu)勢(shì).

      1 已有方法介紹

      對(duì)于地球受力變形,我們采用如下的描述地球變形的微分方程組邊值問(wèn)題.其中,微分方程組為

      (1)

      其中λ,μ是拉梅常數(shù),ρ是密度,g是重力加速度(習(xí)慣上稱之為重力),它們都是球坐標(biāo)系的徑向分量r的函數(shù),G是萬(wàn)有引力常數(shù),UL、UM、Φ、TL、TM分別是垂向位移、水平位移、附加引力位、垂向正應(yīng)力、垂向剪應(yīng)力的球諧展開(kāi)的第n階的系數(shù)(Pan et al.,2015),Q滿足

      (2)

      因此,這6個(gè)變量要分別對(duì)不同的球諧階數(shù)n求解.該方程是研究地球在不同力源作用下發(fā)生靜態(tài)變形的基礎(chǔ)微分方程組,對(duì)于不同的力源需采用相應(yīng)的邊界條件(Love,1911;Farrell,1972;Okubo and Endo,1985;Sun and Okubo,1993).

      對(duì)于地表質(zhì)量負(fù)荷問(wèn)題,通常將負(fù)荷看成點(diǎn)源,并且將該點(diǎn)源放置于北極點(diǎn).由于對(duì)稱性,在負(fù)荷作用下的變形與經(jīng)度無(wú)關(guān),即球諧展開(kāi)的次數(shù)m=0,此時(shí)球函數(shù)退化為L(zhǎng)egendre多項(xiàng)式.

      下面介紹邊界條件.為方便起見(jiàn),我們定義兩個(gè)矢量

      (3)

      在地心處正則,即

      U(r=0)=0.

      (4)

      在內(nèi)部邊界連續(xù),以及在地表

      (5)

      其中,ga是地表重力,a是地球半徑.在推導(dǎo)該式中的邊界值時(shí),認(rèn)為負(fù)荷點(diǎn)源的質(zhì)量M滿足單位引力位假定(即GM/a=1,可參考Sun and Sj?berg,1999),也就是說(shuō)負(fù)荷質(zhì)量M=a/G=me/aga,其中me為地球質(zhì)量.也有一些研究采用地球質(zhì)量假定,即M=me(如Guo et al.,2004).不管用哪種定義,在計(jì)算格林函數(shù)時(shí)都要除以負(fù)荷質(zhì)量以得到單位質(zhì)量負(fù)荷的結(jié)果,因此只要推導(dǎo)過(guò)程自洽都不影響最終的結(jié)果.

      負(fù)荷勒夫數(shù)定義為

      (6)

      其中,帶下標(biāo)n的負(fù)荷勒夫數(shù)hn、ln、kn表示對(duì)應(yīng)于球諧展開(kāi)的n階的地表結(jié)果.因?yàn)?1)式是關(guān)于各個(gè)物理量的球諧展開(kāi)系數(shù)的,所以不同的球諧階數(shù)對(duì)應(yīng)不同的勒夫數(shù).相應(yīng)的負(fù)荷格林函數(shù)為

      (7)

      其中,u,v和φ分別為垂直位移、切向位移和附加引力位,Pn是Legendre多項(xiàng)式,θ是計(jì)算點(diǎn)余緯.可以看到格林函數(shù)是一個(gè)無(wú)窮級(jí)數(shù)求和,需要截?cái)嗵幚?但是對(duì)于負(fù)荷問(wèn)題需要截?cái)嗟胶芨叩碾A數(shù),這給計(jì)算帶來(lái)很大的時(shí)間成本.為此Farrell(1972)根據(jù)Boussinesq問(wèn)題的解(見(jiàn)下文)給出了負(fù)荷勒夫數(shù)的漸近值,從而可有效降低截?cái)嗟碾A數(shù),并加速格林函數(shù)計(jì)算的收斂.

      應(yīng)用Kummer變換,(7)式改寫(xiě)為

      (8)

      其中,帶下標(biāo)∞的負(fù)荷勒夫數(shù)h∞、l∞、k∞表示相應(yīng)的漸近值.由于負(fù)荷勒夫數(shù)漸近值具有嚴(yán)格解析表達(dá)式,右端的第二項(xiàng)通常具有用三角函數(shù)表示的嚴(yán)格解析形式(Farrell,1972;Zhou et al.,2019b),而第一項(xiàng)求和的收斂要比(7)式中的求和快得多,因此截?cái)嚯A數(shù)N可以取得足夠小(對(duì)于位移與附加引力位,N=3000;對(duì)于傾斜或應(yīng)變N=10000足以保證計(jì)算的收斂).可見(jiàn)負(fù)荷勒夫數(shù)漸近值在計(jì)算負(fù)荷格林函數(shù)時(shí)發(fā)揮著至關(guān)重要的作用.

      對(duì)于當(dāng)n很大時(shí)的情形,我們可以忽略液核的影響,因?yàn)橐汉说拇嬖谥挥绊戄^低階的負(fù)荷勒夫數(shù),所以可以采用純固體的地球模型.實(shí)際上,正如下文將要討論的那樣,只有地球最外層的參數(shù)才會(huì)影響負(fù)荷勒夫數(shù)的漸近值,這也是為什么可以直接用均勻球(采用地球最外層結(jié)構(gòu)參數(shù)和地表重力值,如Okubo,1988)或僅僅采用最外層參數(shù)(如郭俊義,2000;Guo et al.,2004)來(lái)解算負(fù)荷勒夫數(shù)的漸近值.

      下面分別介紹目前已有的三種求解負(fù)荷勒夫數(shù)漸近值的方法.

      1.1 Farrell方法

      在球坐標(biāo)系下,當(dāng)角距離(θ)非常小的時(shí)候,問(wèn)題就等效于平面的Boussinesq問(wèn)題,因?yàn)榇藭r(shí)地球的曲率效應(yīng)可忽略(Farrell,1972).對(duì)于Boussinesq問(wèn)題,采用柱坐標(biāo)及柱函數(shù)(Bessel函數(shù))系(用z,r,φ表示,將負(fù)荷點(diǎn)源放置于坐標(biāo)原點(diǎn),由于對(duì)稱性,變形與φ無(wú)關(guān).注意此時(shí)的r表示距離z軸的距離),所有變量都利用Hankel變換(潘爾年等,1986),即

      (9)

      其中,ξ表示波數(shù).為了簡(jiǎn)便,仍用與球諧展開(kāi)系數(shù)相同的符號(hào)表示變換后的系數(shù).在此情形下,球函數(shù)趨近于柱函數(shù),即

      (10)

      并且

      (11)

      因此,(7)式的求和可以表示成如下的積分:

      (12)

      根據(jù)(11)式,(12)式變成

      (13)

      對(duì)于Boussinesq問(wèn)題,一個(gè)點(diǎn)負(fù)荷作用下的解為(在Hankel變換域)

      (14)

      其中ρa(bǔ)表示地球模型最上層的密度,即地表密度,〈ρ〉表示地球平均密度.將(14)式代入(9)式,并與(13)式進(jìn)行比較,可得

      (15)

      在此過(guò)程中,要運(yùn)用到(11)式中的n與ξ的近似關(guān)系n=aξ,因此其結(jié)果精度為O(1/n).我們稱(15)式的結(jié)果為一級(jí)漸近值,以區(qū)別于Guo等(2004)給出的更高一級(jí)的漸近值.

      1.2 Okubo和Endo方法

      Okubo和Endo(1985)、Okubo(1988)推導(dǎo)負(fù)荷勒夫數(shù)漸近值的方法是基于Love(1911)給出的均勻球變形的通解表達(dá)式,用X=[ULUMΦTLTMQ]T表示由6個(gè)變量表示的這個(gè)通解矢量,其中上標(biāo)T表示轉(zhuǎn)置,其解可表示為

      (16)

      其中,c是由6個(gè)常數(shù)標(biāo)量組成的待定系數(shù)矢量,且

      [D(r)]=

      (17)

      (18)

      (19)

      式中

      +[2f±+n(n+1)]zn(k±r)},

      (20)

      (21)

      (22)

      (23)

      α和β分別為P波和S波速度,jn(x)是球Bessel函數(shù),n為球諧階數(shù).利用邊界條件(5)式可確定c,從而求得地表的位移與附加引力位,即

      U(a)=D(a)[E(a)]-1T(a).

      (24)

      由于高階負(fù)荷勒夫數(shù)只與最上層的結(jié)構(gòu)有關(guān),而如果采用密度為地表密度的均勻球,則該均勻球表面的重力將不等于真實(shí)地球的地表重力.Okubo(1988)指出,可用真實(shí)地球的地表重力值代替均勻球表面的重力,這樣得到的負(fù)荷勒夫數(shù)漸近值就是真實(shí)地球的結(jié)果.因此地球模型的參數(shù)代入(24)式,并將各量展開(kāi)成n的級(jí)數(shù)后,就可得到如(15)式所示的負(fù)荷勒夫數(shù)漸近值.

      1.3 Guo方法

      郭俊義等(郭俊義,2000;Guo et al.,2004)根據(jù)地表邊界條件考察了每個(gè)變量與球諧階數(shù)n的量級(jí)關(guān)系,得到

      (25)

      例如,TL隨著n的增大而增大,與n呈線性關(guān)系,而UL不隨n的增大而增大.然后通過(guò)變量代換

      (26)

      使得所有的新變量zi具有關(guān)于n相同的量級(jí),并在微分方程中忽略高階小量,從而得到

      (27)

      該微分方程組存在解析解,可以求出負(fù)荷勒夫數(shù)h和l的漸近值.但是k的漸近值無(wú)法求出,原因是舍去的項(xiàng)太多了.因此他們繼續(xù)保留高一階小量,并求解一個(gè)非齊次的微分方程組,從而求出了k的漸近值.

      在此過(guò)程中,Guo等(2004)未介紹如何求得(27)式的解析解,而是直接給出了結(jié)果.因?yàn)樵摲匠探M實(shí)際上是變系數(shù)的,盡管將n/r移到了方程組的左邊.所以其是否存在解析解,并不是那么直觀地能夠看出來(lái).實(shí)際上,因?yàn)?27)式是歐拉型,所以通過(guò)變量代換,

      r=eξ,

      (28)

      可將該方程組化成常系數(shù)的微分方程組,從而有解析解.Guo等(2004)以一級(jí)漸近值為基礎(chǔ),求解關(guān)于高一級(jí)小量的非齊次微分方程組,給出了更高一級(jí)的漸近值.但是對(duì)于負(fù)荷問(wèn)題,更高一級(jí)的漸近值在計(jì)算負(fù)荷格林函數(shù)時(shí)是沒(méi)有必要的,因?yàn)槔?15)式的漸近值已經(jīng)能夠很高效地獲得收斂的格林函數(shù)(參考下文的圖2).

      圖2 負(fù)荷格林函數(shù)隨角距離的變化(規(guī)格化系數(shù)為1012aθ)“asym.with g”表示采用考慮自引力的50000階數(shù)值結(jié)果作為漸近值,“asym.without g”表示采用不考慮自引力的50000階數(shù)值結(jié)果作為漸近值.

      2 解析解方法

      本文介紹的推導(dǎo)負(fù)荷勒夫數(shù)漸近值的新方法是以Pan等(2015)提出的解析法和DVP傳播矩陣方法為基礎(chǔ).其中解析法的思想是使得地球變形的通解具有解析的形式,而DVP傳播矩陣法是建立高階地表負(fù)荷勒夫數(shù)與地表邊界條件之間的直接聯(lián)系.下面首先簡(jiǎn)單介紹解析法和DVP傳播矩陣法,具體內(nèi)容可參考(Pan et al.,2015),然后給出推導(dǎo)負(fù)荷勒夫數(shù)漸近值的過(guò)程.

      2.1 解析解構(gòu)造

      Pan等(2015)提出了求解(1)式的解析法.對(duì)于一個(gè)分層地球模型,(1)式是沒(méi)有解析解的,但是通過(guò)合理的假設(shè)可以構(gòu)造出近似的解析解.首先,在地核中假設(shè)每層密度是常數(shù),重力是r的線性函數(shù); 在地幔和地殼中,重力是常數(shù),密度是r的反函數(shù),如在第j層中,重力和密度滿足

      (29)

      (30)

      其中,b和c是常數(shù),rj-1表示第j層的下邊界,rj表示第j層的上邊界,rc表示核幔邊界的半徑,ρ0j是第j層的平均密度,r0j是使該層質(zhì)量守恒的等效半徑,它滿足

      (31)

      以上的假設(shè)是否合理,一方面可以考察假設(shè)后的模型與原模型的近似程度,另一方面取決于計(jì)算結(jié)果是否準(zhǔn)確.Pan等(2015)的數(shù)值結(jié)果已經(jīng)證明了結(jié)果的可靠性以及假設(shè)的合理性,實(shí)際上上述假設(shè)正是基于PREM模型的參數(shù)分布特征作出的,假設(shè)后的模型與原模型的符合程度是非常高的,差異僅為0.6%.并且這種方法也成功運(yùn)用到求解位錯(cuò)以及地表熱負(fù)荷引起的變形問(wèn)題中,也被證明是非常高效的(Zhou et al.,2019a,b,2020,2021a,b).

      此外,作一個(gè)變量代換

      (32)

      顯然,對(duì)于j層的下邊界有ξ=0.再以rTL,rTM和rQ為新變量,則(1)式變?yōu)?/p>

      (33)

      如果再假設(shè)第j層內(nèi)的彈性模量為常數(shù)(目前的研究基本都是這么做的),那么可以看出上述微分方程組的系數(shù)都是常數(shù),這樣就存在解析解.為簡(jiǎn)便計(jì),把(33)式表示成矢量和矩陣的形式,即

      (34)

      其中,

      (35)

      為了方便,我們用球諧階數(shù)n對(duì)變量作了規(guī)格化處理使各個(gè)變量量級(jí)相當(dāng).相應(yīng)的矩陣A可從(33)式導(dǎo)出(系數(shù)作相應(yīng)的調(diào)整,乘上與n有關(guān)的系數(shù)).此時(shí),通解可以表示成(如在地球模型的第j層)

      (36)

      (36)式右邊第一個(gè)矩陣是系數(shù)矩陣A的特征向量組成的矩陣,下標(biāo)j表示它由第j層的參數(shù)計(jì)算獲得,第二個(gè)矩陣中的s1和s2分別是由具有正實(shí)部(上標(biāo)+)和負(fù)實(shí)部(上標(biāo)-)的三個(gè)特征值組成,即

      (37)

      從而,由傳播矩陣法建立第j層上下邊界處(ξ=ξj和ξ=0)的兩組解之間的關(guān)系:

      (38)

      其中上標(biāo)u和l分別表示第j層的上下邊界,上標(biāo)-1表示矩陣求逆.然而,Pan等(2015)的數(shù)值結(jié)果表明,由于數(shù)值溢出,普通的傳播矩陣方法,即(38)式,無(wú)法正確算出高階(如n>6000)結(jié)果.下面要介紹的DVP傳播矩陣法可以很好地解決這個(gè)問(wèn)題.

      2.2 DVP傳播矩陣法

      該方法的核心思想是避免導(dǎo)致溢出的大數(shù)產(chǎn)生,其技巧僅在于對(duì)普通的傳播矩陣法進(jìn)行小的改動(dòng),即

      1)將(36)式改寫(xiě)為

      由于ξj是常數(shù),所以最終只是改變了待定系數(shù),從(36)式中的c1變成(39)式的c3,并不會(huì)改變最終的解.

      2)在第j層中,由上下邊界的解我們可以得到與(38)式不一樣的關(guān)系:

      (40)

      其中,

      (41)

      需要指出的是,這種傳播矩陣方法的特別之處在于變量的位置關(guān)系與傳統(tǒng)的傳播矩陣方法是不一樣的(對(duì)比(38)與(40)兩式).從(41)式可以看出,在傳播矩陣的計(jì)算中不涉及大的正數(shù)的指數(shù),從而避免大數(shù)的產(chǎn)生,解決了溢出的問(wèn)題.

      2.3 漸近值推導(dǎo)

      基于上面的方法原理,我們就可以推導(dǎo)出負(fù)荷勒夫數(shù)的漸近值.特別地,當(dāng)球諧階數(shù)n很大時(shí),由于特征值接近于±n,所以(41)式中的〈-s1ξj〉和〈s2ξj〉都接近于0,當(dāng)n為無(wú)窮大時(shí),它們就等于0.因此(41)式變成

      (42)

      由(40)式進(jìn)一步得到

      (43)

      這說(shuō)明,當(dāng)球諧階數(shù)很大時(shí),U和T是相關(guān)的,這種相關(guān)性會(huì)導(dǎo)致傳統(tǒng)的龍格-庫(kù)塔數(shù)值積分法無(wú)法準(zhǔn)確計(jì)算出非常高階的負(fù)荷勒夫數(shù); 然而這種相關(guān)性可以使我們能夠直接求得表面的負(fù)荷勒夫數(shù)漸近值.對(duì)于地球的最外層(設(shè)為第p層),我們可由(43)式直接得出

      (44)

      其中,方括號(hào)的下標(biāo)p表示方括號(hào)內(nèi)的矩陣由第p層的參數(shù)求得.這說(shuō)明,只要我們知道最外層的系數(shù)矩陣A的特征向量,就可以得到當(dāng)n趨于無(wú)窮大時(shí)的U,從而得到負(fù)荷勒夫數(shù)的漸近值.

      通常情況下,矩陣A的特征向量是無(wú)法解析地表達(dá)出來(lái)的.但是我們知道,對(duì)于地表質(zhì)量負(fù)荷問(wèn)題,在計(jì)算高階的負(fù)荷勒夫數(shù)時(shí),地球的自引力效應(yīng)是可以忽略的,球諧階數(shù)越高,自引力效應(yīng)的影響越不明顯(見(jiàn)下文圖1).實(shí)際上Farrell(1972)從Boussinesq近似得到的漸近值也是忽略了自引力影響的.

      圖1 負(fù)荷勒夫數(shù)隨球諧階數(shù)的變化

      對(duì)于忽略自引力的情形,在(33)式關(guān)于應(yīng)力的微分方程中就可以舍去與密度或重力有關(guān)的項(xiàng),從而位移和應(yīng)力是與附加引力位無(wú)關(guān)的,因此可以首先解出(如Wason and Singh,1972;Fang et al.,2014),然而附加引力位是與位移相關(guān)的,所以繼續(xù)求附加引力位時(shí)需要求解一個(gè)非齊次的微分方程組,這也是1.3節(jié)中Guo等(2004)求負(fù)荷勒夫數(shù)k的漸近值采用的方法,但這反而增加了難度.實(shí)際上,我們可以直接運(yùn)用2.1節(jié)的方法,此時(shí)矩陣A的特征向量可以很容易地求出(如通過(guò)MATLAB的符號(hào)運(yùn)算,具體表達(dá)式可參考Zhou et al.,2021a).

      因此根據(jù)矩陣A的特征向量,以及(5)式和(44)式可以很容易求得

      (45)

      (45)式是高階負(fù)荷勒夫數(shù)的嚴(yán)格解析表達(dá)式,給出了負(fù)荷勒夫數(shù)與球諧階數(shù)n的關(guān)系.特別地,當(dāng)n趨于無(wú)窮大時(shí), 取極限就得到

      (46)

      由于我們是求負(fù)荷勒夫數(shù)在n趨于無(wú)窮大時(shí)的漸近值,我們可以將第p層設(shè)得足夠薄,因此有

      ρ0pr0p=ρa(bǔ)a,

      (47)

      即地表密度與地球半徑的乘積.這樣就可以得到與(15)式完全一致的結(jié)果.實(shí)際上對(duì)于最外層來(lái)說(shuō),地球模型通常是均勻的,因此ρ0p就是這個(gè)平均值,而由(31)式計(jì)算得到的r0p實(shí)際已經(jīng)非常接近地球半徑a了(對(duì)于PREM模型來(lái)說(shuō),誤差小于萬(wàn)分之三),就算直接使用r0p也可以確保k∞的準(zhǔn)確性.

      2.4 數(shù)值結(jié)果

      圖1給出了根據(jù)(45)式確定的負(fù)荷勒夫數(shù)隨球諧階數(shù)的變化, 并與Guo等(2004)考慮高一階小量的漸近值進(jìn)行了比較.其中地球模型采用PREM(Dziewonski and Anderson, 1981),其最外層的海洋層替換為固體層,其參數(shù)設(shè)置如下:λ=3.43×1010Pa,μ=2.66×1010Pa,ρ=2600 kg·m-3.從圖中可以看出,隨著球諧階數(shù)n的增大,負(fù)荷勒夫數(shù)趨近于一個(gè)常數(shù).由于我們推導(dǎo)高階負(fù)荷勒夫數(shù)解析公式時(shí)未考慮自引力的影響,因此所得結(jié)果在高階部分與不考慮自引力的數(shù)值結(jié)果比較接近,在n=5000左右時(shí)二者就已經(jīng)比較一致;但是在低階部分與考慮自引力的數(shù)值結(jié)果仍然存在很大的差異.通常計(jì)算負(fù)荷格林函數(shù)時(shí)只要考慮10000階以下的負(fù)荷勒夫數(shù),因此在研究負(fù)荷形變時(shí)自引力的影響一般是需要加以考慮的.但是從圖中可以看出,隨著n的增大,自引力的影響越來(lái)越小,因此在研究負(fù)荷勒夫數(shù)漸近特征時(shí),自引力是可以忽略的,這可從圖1中不同結(jié)果在n很大時(shí)具有的一致性得到驗(yàn)證.

      此外,我們也將本文給出的漸近值與Guo等(2004)的結(jié)果進(jìn)行了比較.結(jié)果說(shuō)明,在n<20000左右時(shí),Guo等(2004)的漸近值與考慮自引力的數(shù)值結(jié)果更加接近,這是因?yàn)樗麄兺茖?dǎo)漸近值時(shí)部分地考慮了自引力的影響;但是對(duì)于n>20000的部分,我們的結(jié)果不管是和考慮自引力還是和不考慮自引力的數(shù)值結(jié)果都更加接近,對(duì)于位移的勒夫數(shù),這種情況更加明顯.

      為了顯示漸近值在計(jì)算格林函數(shù)中的作用,作為例子,我們根據(jù)(8)式計(jì)算了垂直位移、水平位移和附加引力位的格林函數(shù)隨角距離θ的變化關(guān)系,結(jié)果如圖2所示.為了更好的顯示結(jié)果,格林函數(shù)值都乘上了一個(gè)規(guī)格化系數(shù)1012aθ.首先,通過(guò)對(duì)比采用不同漸近值獲得的格林函數(shù),發(fā)現(xiàn)漸近值之間的微小差異不會(huì)影響格林函數(shù)的計(jì)算.并且我們發(fā)現(xiàn)負(fù)荷勒夫數(shù)漸近值在附加引力位的計(jì)算中并不是必需的,不采用Kummer變換同樣可以獲得收斂的格林函數(shù).圖2也說(shuō)明了如果不采用Kummer變換,即使截?cái)嚯A數(shù)N=50000,位移格林函數(shù)仍然是不收斂的,說(shuō)明漸近值對(duì)于位移格林函數(shù)的計(jì)算來(lái)說(shuō)是必不可少的,這從另一個(gè)側(cè)面證明了漸近值的重要性.

      3 討論

      從2.3節(jié)的推導(dǎo)我們可以看出,負(fù)荷勒夫數(shù)的漸近值只和地球模型地表的參數(shù)有關(guān),而不論地球內(nèi)部的結(jié)構(gòu)如何.其實(shí),這一點(diǎn)早就為我們所熟知,即高階的負(fù)荷勒夫數(shù)與地殼和上地幔介質(zhì)密切相關(guān),因此采用龍格-庫(kù)塔方法計(jì)算高階負(fù)荷勒夫數(shù)時(shí)可以直接從地幔開(kāi)始積分.然而,通過(guò)我們的推導(dǎo)可以從本質(zhì)上了解其原因:地球每個(gè)分層內(nèi)的U和T之間的相關(guān)性隨著球諧階數(shù)的增大而增加,從而地表的負(fù)荷勒夫數(shù)可由地表的邊界條件直接確定.

      相比于Okubo和Endo方法,我們給出的新方法不需要作均勻球的假設(shè)(或基于均勻球模型開(kāi)展相關(guān)推導(dǎo)),不需要用真實(shí)地球的地表重力代替均勻球表面的重力.目前基于均勻球的地球受力變形的數(shù)值模擬必須要替換地表重力才能得到正確的結(jié)果,如Sun(2004)對(duì)位錯(cuò)引起的地表位移場(chǎng)變化的模擬.與Guo等(2004)的方法一樣,我們的方法實(shí)際上也只需要考慮地球模型的最外層即可.

      從Guo等(2004)給出的負(fù)荷勒夫數(shù)漸近值的更高階的結(jié)果可以看出,在位移(h和l)結(jié)果中耦合進(jìn)了重力的影響.但是將(44)式展開(kāi)成n的羅朗級(jí)數(shù)后高階解中沒(méi)有重力項(xiàng),原因在于我們忽略了自引力的影響.對(duì)于地表負(fù)荷問(wèn)題,忽略自引力不會(huì)影響負(fù)荷勒夫數(shù)漸近值的一級(jí)結(jié)果,即O(1)量級(jí)的常數(shù)項(xiàng),但是會(huì)影響高一級(jí)的結(jié)果,但是這種影響隨著球諧階數(shù)的增大越來(lái)越弱,并且高一級(jí)的結(jié)果在格林函數(shù)的計(jì)算中是不必要的.在2.3節(jié)中我們忽略自引力的目的是為了獲得矩陣A的特征向量的解析形式,從而獲得漸近值的解析形式,這樣就可以運(yùn)用Kummer變換.Kummer變換的關(guān)鍵有兩點(diǎn):一方面,漸近值必須與實(shí)際的勒夫數(shù)比較接近,從而它們的差是一個(gè)高階小量,并且這個(gè)差異隨著球諧階數(shù)的增大而減小,這樣(8)式右端的第一個(gè)級(jí)數(shù)求和能夠很快地收斂;另一方面,漸近值需要具有嚴(yán)格解析式,從而(8)式右端的第二個(gè)級(jí)數(shù)求和有簡(jiǎn)單的解析式,能夠用簡(jiǎn)單的函數(shù)表示.Chen等(2018)采用的剛度矩陣法或者2.2節(jié)介紹的DVP都能夠計(jì)算任意球諧階數(shù)的負(fù)荷勒夫數(shù),但是要獲得收斂的負(fù)荷格林函數(shù),具有解析表達(dá)式的漸近值對(duì)于簡(jiǎn)化負(fù)荷格林函數(shù)的計(jì)算仍是非常必要的,因?yàn)榻財(cái)嗟胶芨叩那蛑C階數(shù)是很費(fèi)時(shí)低效的.采用擬合的方法確定當(dāng)n非常大時(shí)的勒夫數(shù)隨n變化的解析式(即漸近值)也是一個(gè)非常有效的方法(如Pan et al.,2015;Zhou et al.,2019b).

      本文的方法除了采用Pan等(2015)的解析解方法獲得常微分方程組的系數(shù)矩陣A的特征向量解析表達(dá)式外,另一個(gè)關(guān)鍵之處是采用了DVP傳播矩陣方法,(42)式顯示了其獨(dú)特的性質(zhì),階數(shù)n越大,其計(jì)算方面的優(yōu)越性越明顯.DVP方法只是在傳統(tǒng)的傳播矩陣方法基礎(chǔ)上作了小的改動(dòng),變換了變量的位置,就解決了數(shù)值溢出的問(wèn)題,并且計(jì)算也非常的高效.

      最后需要強(qiáng)調(diào)的是,從表面上看,本文的方法涉及構(gòu)造解析解的理論方法和DVP傳播矩陣方法.實(shí)際上推導(dǎo)負(fù)荷勒夫數(shù)的過(guò)程非常簡(jiǎn)單,可以分成以下兩個(gè)步驟:首先計(jì)算常微分方程組的系數(shù)矩陣A的特征向量從而得到(44)式右端方括號(hào)中的矩陣;然后根據(jù)地表邊界條件由(44)式通過(guò)矩陣運(yùn)算直接得到需要的漸近值的解析表達(dá)式.該過(guò)程只涉及矩陣特征向量的計(jì)算及相關(guān)的矩陣運(yùn)算,原理簡(jiǎn)單,也非常容易編程實(shí)現(xiàn)(Pan,2019;Zhou et al.,2019a).

      4 結(jié)論

      本文基于求解地球變形的解析解方法和DVP傳播矩陣方法,導(dǎo)出了高階負(fù)荷勒夫數(shù)的解析表達(dá)式,進(jìn)而確定了負(fù)荷勒夫數(shù)的漸近值.漸近值的推導(dǎo)只與地球模型的最外層有關(guān),推導(dǎo)過(guò)程只涉及矩陣特征向量計(jì)算和簡(jiǎn)單的矩陣運(yùn)算.本文所得結(jié)果與已有方法給出的結(jié)果是一致的,并且更加準(zhǔn)確地描述了負(fù)荷勒夫數(shù)的漸近特征.

      利用DVP傳播矩陣法對(duì)于n非常大時(shí)所具備的特點(diǎn),我們不但可以看出變量U和T的關(guān)聯(lián)性,這也是導(dǎo)致傳統(tǒng)的龍格-庫(kù)塔數(shù)值積分法失效的本質(zhì)原因,還可以直接根據(jù)地表邊界條件求出負(fù)荷勒夫數(shù)的漸近值.這說(shuō)明DVP傳播矩陣法在數(shù)值計(jì)算方面有其獨(dú)特的優(yōu)勢(shì),可以廣泛應(yīng)用于分層地球模型的變形研究中,如負(fù)荷形變、地震形變等靜態(tài)變形問(wèn)題,以及地球簡(jiǎn)正?;蚝铣傻卣饒D等動(dòng)態(tài)變形問(wèn)題中.

      致謝本文第二作者于1982年至1984年在北京大學(xué)地質(zhì)系師承王仁教授攻讀碩士學(xué)位,學(xué)位論文題目為《成層地球模型對(duì)體勢(shì)力及表載的響應(yīng)》.該論文利用球諧函數(shù)與傳播矩陣方法研究分層球形地球在日月引潮力和地表負(fù)載作用下的位移場(chǎng)和應(yīng)力場(chǎng)變化,并采用數(shù)值積分方法計(jì)算勒夫數(shù).在王仁教授的教誨下,科研能力得到了很大的提升,為后來(lái)從事力學(xué)和數(shù)學(xué)相關(guān)的研究奠定了堅(jiān)實(shí)的基礎(chǔ),對(duì)現(xiàn)在的研究工作仍有很大的幫助.

      猜你喜歡
      勒夫格林引力
      麻辣老師
      我喜歡小狼格林
      小讀者(2020年4期)2020-06-16 03:34:04
      綠毛怪格林奇
      電影(2018年12期)2018-12-23 02:19:00
      凱文·勒夫感動(dòng)
      NBA特刊(2018年15期)2018-11-21 08:47:32
      KEVIN LOVE'S OUTLET PASSES 凱文·勒夫 快攻發(fā)動(dòng)機(jī)
      NBA特刊(2018年7期)2018-06-08 05:48:45
      引力
      初中生(2017年3期)2017-02-21 09:17:40
      勒夫:美學(xué)家,戰(zhàn)略家,世界冠軍
      出版廣角(2016年15期)2016-10-18 00:47:17
      感受引力
      格林的遺憾
      山東青年(2016年1期)2016-02-28 14:25:24
      A dew drop
      邯郸县| 东城区| 神木县| 长丰县| 乐昌市| 双辽市| 垫江县| 昌都县| 法库县| 景德镇市| 津市市| 滨海县| 淄博市| 封开县| 莲花县| 台山市| 浦北县| 合肥市| 无极县| 哈尔滨市| 定南县| 肥城市| 昭通市| 建水县| 正安县| 江门市| 牡丹江市| 凤庆县| 咸丰县| 水富县| 射阳县| 甘谷县| 宜川县| 自贡市| 安福县| 静安区| 平泉县| 吉首市| 雷州市| 长阳| 曲麻莱县|