引力不變量曲線坐標系及在定位導航上的應用

于錦海，徐 煥，萬曉云

1. 中國科學院計算地球動力學重點實驗室，北京 100049; 2. 中國科學院大學地球與行星科學學院，北京 100049; 3. 中國地質(zhì)大學(北京)土地科學與技術(shù)學院，北京 100083

地球重力場(或引力場)的研究工作目前已經(jīng)取得了較大進展。例如：綜合利用各類重力觀測數(shù)據(jù)建立了EGM2008模型，重力場的分辨率有了極大的提高。又如：衛(wèi)星重力技術(shù)的應用，特別是GRACE衛(wèi)星計劃的實施，應用時變重力場信息可以解釋地球質(zhì)量的遷移特征[1-2]。此外，重力測量技術(shù)也獲得了較大的發(fā)展，特別是衛(wèi)星測高技術(shù)與重力梯度測量技術(shù)的應用[3-6]，為解算重力場提供了豐富的觀測數(shù)據(jù)。所有這些關(guān)于重力場方面的進展，使得進一步地拓展重力場應用成為可能。

重力輔助導航是重力場理論的潛在應用方向之一[7-9]。自20世紀90年代開始，便有國外的學者進行了研究與模擬測試[10-14]，而我國也有學者開展了相應的研究工作[15-16]。若將上述研究工作進一步延伸，能否利用重力理論來實現(xiàn)自主定位與導航(即單純使用重力測量數(shù)據(jù)實現(xiàn)目標的定位與導航)。顯然，這是一項極有挑戰(zhàn)的工作，其中包含了利用重力進行定位與導航的理論問題、重力測量儀器的研制、高精度與高分辨率的重力場模型研發(fā)等諸多因素[17]。

如何建立重力場(引力場)數(shù)據(jù)與相應空間點位之間的直接關(guān)系，以及建立的關(guān)系有何應用價值，是能否利用重力數(shù)據(jù)進行自主定位與導航的基礎(chǔ)理論問題。本文的目標從引力場對應的不變量入手，引入一組曲線坐標系，并研究其性質(zhì)和探討潛在的應用。因為引入的曲線坐標系是由不變量構(gòu)成的，因此該坐標系直接描述引力場與空間點位之間直接的關(guān)系。

1 引力場的不變量

(1)

(2)

經(jīng)計算，這些不變量可寫為[19-20]

(3)

由于引力位v滿足Laplace方程，即I1=0，所以從引力梯度數(shù)據(jù)可構(gòu)造出函數(shù)I2和I3。因為I2和I3是不變量，這意味著I2和I3也僅是點位的函數(shù)，與引力梯度分量的方向無關(guān)。

至此，就從引力場自身引入了3個僅與點位相關(guān)的函數(shù)g、I2和I3。由于這3個函數(shù)都有各自的物理量綱，所以需進行去量綱化處理。

為了確保不變量坐標系(f1,f2,f3)構(gòu)成了曲線坐標系，需要在理論上證明相應的Jacobi行列式是非零的，即

(4)

事實上，因為Jacobi行列式J是不恒為零的，所以從三維空間分布看面，J=0僅可能構(gòu)成空間的曲面，圖1的結(jié)果表明，該曲面大致是赤道沿徑向向外延伸的曲面。

圖1 利用EGM2008前180階次模型計算的Jacobi行列式在地球平均球面上的分布Fig.1 The distribution of Jacobi’s determinant for EGM2008 model with the former 180 degree and order on the average sphere of the earth

2 不變量坐標系與空間點位之間的關(guān)系

本節(jié)將通過某些實例計算來驗證：假若事先給出了不變量坐標系的具體表達式，則利用該表達式就可以進行空間點的定位。正如地圖一樣，通過對若干參照物的觀測可以利用地圖作為工具來確定點位的坐標，因此在某種意義上看，不變量坐標系起著地圖的作用。

本節(jié)將以EGM2008引力場模型的前360階次作為實際地球引力場，由此便可以得到相應的引力不變量曲線坐標系(f1,f2,f3)。假設(shè)在地面上某點P進行引力和引力梯度觀測，并進行數(shù)據(jù)處理后得到了P點在不變量坐標系下的值(F1,F2,F3)，其目標是解算出P點的空間坐標(r,θ,λ)。事實上，從P點的觀測值(F1,F2,F3)，可以建立下列方程組

(5)

由于不變量坐標系(f1,f2,f3)是已知的，故通過解算上述方程組便能得到P點對應的球坐標。顯然方程組(5)關(guān)于待解算的變量(r,θ,λ)是非線性的，因此直接進行求解是難以實現(xiàn)的。理論上講，對于形式如方程組(5)的非線性方程而言，只要在某點處的Jacobi行列式非零，那么在該點附近是局部可解的，用數(shù)學語言描述，就是存在該點的某個鄰域，在該鄰域內(nèi)方程組是唯一可解，而且解還具有穩(wěn)定性[21]。至于求解方法通常是采用Newton迭代法，即線性化迭代解法。

下面就方程組(5)，給出具體的線性化迭代過程如下

(6)

這里j=0,1,2,…，而初始值(r0,θ0,λ0)的選擇與待計算的P點要盡可能地接近，其理由就是因為非線性方程組的解僅在局部具有存在性、唯一性、穩(wěn)定性。在方程組(6)中，待解的量是(rj+1,θj+1,λj+1)，因此方程組(6)關(guān)于解算量是線性的，所以方程組(6)稱為方程組(5)的線性化迭代形式。由于方程組(6)是線性的，所以其可解性完全取決于在(rj,θj,λj)處的Jocobi行列式J(rj,θj,λj)。如果能通過實際解算驗證方程組(6)關(guān)于(rj+1,θj+1,λj+1)是唯一可解的，那么便間接地證明了Jocobi行列式是非零的。此外，至于需要進行多少次迭代，這取決于事先給出的誤差。事實上，只要迭代前后解的差小于誤差要求，便可以停止迭代，從而得到所需的解。

為了驗證方程組(5)以及線性化迭代形式(6)的可解性，本文通過3個低、中、高緯度實際點位上的不變量坐標值對其進行了驗證，其解算結(jié)果匯總見表1—表3。

表1 低緯度點驗證結(jié)果

表2 中緯度點驗證結(jié)果

表3 高緯度點驗證結(jié)果

表1—表3中，第1列是點位坐標與實測不變量坐標值的記號，第2列是實際的點位與不變量坐標值，第3列是迭代形式的初始值(r0,θ0,λ0)以及對應的不變量坐標值，第4列是通過迭代解算式(6)后得到結(jié)果，第5與第6列則是初始結(jié)果和最終計算結(jié)果分別與實際值的差，而最下面的行(距離)則是誤差。表1表明，對于低緯度點，迭代初始值的點位距實際解算點位17.676 km，通過迭代計算后得到的最終點位坐標與實際點位的差僅有3.87 mm。表2表明，對于中緯度點，迭代初始值的點位距實際解算點位32.634 km，通過迭代計算后得到的最終點位坐標與實際點位的差僅有1.34 mm。表3表明，對于高緯度點，迭代初始值的點位距實際解算點位42.519 km，通過迭代計算后得到的最終點位坐標與實際點位的差僅有7.69 mm。為了清晰地反映迭代求解過程，本文將表1—表3中的計算迭代過程與次數(shù)匯制成圖2。從圖2能清晰可見關(guān)于低、中、高緯度點位的迭代解算過程，此外迭代次數(shù)分別是4、3、9次。

圖2 迭代求解過程Fig.2 Graphical processes of iteration

通過上述關(guān)于低、中、高緯度點的驗證，可得到如下結(jié)論，若事先能給出不變量坐標系，那么通過重力與重力梯度觀測是可以反解出空間點的位置坐標的，這為單純利用重力數(shù)據(jù)進行自主定位導航提供了新的途徑。

需要說明的是：①若出現(xiàn)Jacobi行列式為零的情況，那么本文給出的方法將無法反演出點位坐標，如圖1所示，在赤道附近有可能會出現(xiàn)這樣的情況；②諸如引力(重力)和梯度測量是可以連續(xù)實施的，因此在求解線性化方程組(6)時其初始值通常應該取為上一個點位的坐標，這樣就能保證初始值與定位點的距離差距不大。

3 結(jié)論與潛在的應用

本文利用地球引力場(或重力場)的性質(zhì)，描述了相應的引力不變量曲線坐標系。利用該坐標系，可以進行空間點位的反演計算。事實上，不變量坐標系在定位時就像地圖一樣，是需要事先確定的。不同的是，地圖是紙質(zhì)或數(shù)字形式表現(xiàn)出來的，而不變量坐標系是以函數(shù)形式表現(xiàn)出來的。若要將本文論述不變量曲線坐標系用于實際定位與導航應用，大致需要進行下列幾項工作：①通過各種引力(重力)測量方法來獲取相關(guān)數(shù)據(jù)，對其進行處理后給出不變量坐標系的表達形式；②在進行實際定位或?qū)Ш綍r，需要進行實時的引力與引力梯度測量，即需要精度滿足要求的加速度計和梯度測量儀器；③采用理論上可行的算法，利用事先建立的不變量坐標系與實時測量值來反演點位的坐標。本文工作主要是論證上述方法在理論上的可行性，并針對第3部分內(nèi)容進行分析與討論。

如何事先給出引力(重力)不變量曲線坐標系，這取決于地球引力場(重力場)模型的建立。如果建立了精度足夠高的引力場(重力場)模型，則利用該模型自然就可以生成不變量坐標系。目前EGM2008重力場模型已經(jīng)達到了2160階次，分辨率可達8 km左右，而國內(nèi)和國際上即將推出更高階次的重力場模型，這對于建立高精度的全球不變量曲線坐標系無疑有著直接的作用。對于局部問題(例如南海區(qū)域)可以通過在該區(qū)域進行高分辨率的重力測量，結(jié)合衛(wèi)星海洋測高數(shù)據(jù)解算出區(qū)域重力場模型，例如，可采用球冠諧級數(shù)模型或移去恢復法建立局部模型[22-26]，由此便可得到局部區(qū)域相應的不變量坐標系。總之，要實現(xiàn)事先能給出不變量曲線坐標系，需要進行大量的實際測量、數(shù)據(jù)處理等系列工作，正如要繪制出精確的地圖一樣，需要事先完成系列工程性的工作。顯然，這不是本文能夠解決的問題。

為何引入引力梯度不變量作為曲線坐標系，而不是采用梯度的分量。這是因為梯度分量不僅與點位有關(guān)，而且依賴于方向，這會導致更多參數(shù)的引入。單純從數(shù)學角度講，若使用梯度分量，則對應數(shù)組已不再是曲線坐標系的概念了。從實際應用看，更多參數(shù)的引入會增加理論與計算的復雜性，甚至會導致欠定問題，即給出的條件數(shù)少于待解未知數(shù)的個數(shù)。

在引力不變量曲線坐標系中，引力(重力)是可以直接測量的，而引力(重力)梯度不變量I2和I3則是梯度分量的組合。事實上，文獻[20]曾證明：計算I2和I3僅需梯度的3個對角分量vrr、vθθ和vλλ即可，理由是這3個對角分量是梯度分量的主項，其余分量的計算可利用重力場模型迭代計算來解決。對于該3個對角分量，由于徑向r與引力(重力)的反方向相差不大，所以可以利用垂線方向作為標定軸向，其余兩個軸向只要與標定軸垂直即可?？傊?，在計算不變量I2和I3時，梯度的3個對角分量是必不可少的。概括來說，若在某點給出了引力以及3個梯度對角分量(共4個觀測量)，則可以算出該點處的不變量坐標值，從而解出該點的空間坐標位置。如果不采用不變量坐標系，需要解算的量是坐標位置(3個分量)和梯度姿態(tài)(3個歐拉角)共6個未知數(shù)，因此將導致欠定問題，即4個觀測量，6個未知數(shù)。

至于第2項主要工作，即加速度計與梯度測量儀器的研制與應用，顯然是能否將本文提出的方法予以實施的關(guān)鍵。對于引力(重力)測量，微伽級精度的加速度計已經(jīng)用于實際測量，因此引力測量的條件是具備的。關(guān)于引力(重力)梯度測量，目前得到應用的最好成果當屬GOCE衛(wèi)星上搭載的梯度測量儀[27]，其精度能達到幾個mE(10-12s-2)量級[3]。此外，美國Maryland大學[28]在實驗室里研制出了精度達到0.14 mE的重力梯度儀，而我國也有研究機構(gòu)正致力于重力梯度儀的研制。期待不久的將來，能生產(chǎn)出可用于實際測量的重力梯度儀。

本文的目的是闡述單純利用引力場(或重力場)理論進行自主定位與導航的可行性。盡管學術(shù)界已經(jīng)提出過重力輔助匹配導航的思路，但是以曲線坐標系的形式進行引力場(重力場)獨立定位導航的思路尚屬首次提出，這對于拓展重力場的應用無疑是有益的。