确定球近似大地水准面的积分Stokes方法

何风勇+杨吉明

（济南市勘察测绘研究院，山东 济南 250013）

摘要：文章介绍了计算求解大地水准面的Stokes公式的一种新方法—积分Stokes方法。给出了积分Stokes函数，通过该函数计算大地水准面时，在计算点处不存在奇异性。作为示例文章采用了连续分布的重力异常函数进行了试算，结果表明，利用该方法计算所得的结果与普通的数值积分方法比较具有较高的精度。

关键词：大地水准面；Stokes公式积分；Stokes函数；计算点处

中图分类号：P223文献标识码：A文章编号：1009-2374（2014）22-0061-02确定地球形状是大地测量工作者的一项基本任务，通常采用大地水准面来表征地球形状，它是覆盖全球的一个重力等位面，其上的位与正常椭球表面的位相等。1849年，Stokes导出了计算大地水准面的积分公式—Stokes公式，为确定地球形状提供了有力的工具。由于该公式是对全球进行积分，围绕计算的精度与速度，不同学者进行了广泛的研究，提出了许多计算方法，如远近区相结合的方法和谱分析方法，以及Stokes函数的改化。1980年Goad在计算海潮负荷时提出了积分格林函数方法，我国学者吴庆鹏也对此做了更进一步的研究，他们的研究表明该方法具有显著的优点，特别是能够使得计算点处的奇异性降低1阶。本文提出的方法即是将这一思想运用于由重力异常分布确定大地水准面的计

算中。

1计算理论

Stokes公式具有非常简洁的形式，

（1）

其中是地球半径，是平均椭球体表面的正常重力的平均值，是重力异常，被称为Stokes函数，它是测站与计算点的角距的函数（如图1），有如下的公式：

（2）

其中：

（3）

分别为计算点的余纬与经度，分别为测站的余纬与经度，表示半径为的单位球，且有：

（4）

上式中，A为方位角。

图1 Stokes函数

从（1）式可以看出，计算大地水准面差距也就是计算重力异常与Stokes函数的乘积在整个球面上的积分，因此在某种意义上，Stokes函数可以认为是一种格林函数，即任一测站的重力异常对大地水准面差距的贡献的响应函数。以下我们所称的格林函数即为Stokes

函数。

我们定义积分格林函数为：

（5）

δ为积分Stokes函数所处的区间间隔大小，即以该区间的中点代替积分Stokes函数的自变量。通过推导可得积分Stokes函数的解析表达式为：

（6）

其中：

（7）

显然，在计算点处，y=0，即t=0，由（2）知，Stokes函数表达式中的分母为0，因此产生计算时产生奇异，有许多文献讨论了避免奇异性的方法，如将该点处的重力异常展开为泰勒级数。

而（6）中包含对数函数的部分的值为

（8）

上式中运用了求极限的罗必塔法则。因此，积分Stokes函数在计算点是没有奇异性的。图2给出了其随角距离的变化，可以看出，积分Stokes函数的变化是比较平缓的，而不像Stokes函数那样在近区的变化比较剧烈（注意图1中横轴取了对数），图中红色的线表示由（6）式计算的结果，而黑色的线表示由（5）式计算的结果（积分步长为0.001°），很明显二者具有非常好的一致性，说明（6）式是正确的。

图2 积分Stokes函数

令DA为沿方位角方向的网格划分，那么数值积分时，认为重力异常在任一网格内为常数，则将（5）式代入（1）式可得：

（9）

2试验算例

由于我们仅仅考虑我们提出的方法的优劣，因此采用有解析结果的算例进行验证。本试验中，令重力异常在全球的分布满足

（10）

那么对于北极点（y=q）的大地水准面差距来说，其具有严格的解析解：

（11）

这样我们就可以用来检核不同方法计算结果的正确性与精度。对于北极此时我们就可以利用重力异常分布的对称性，并且在这两点的大地水准面差距具有严格的解析解。

从图1中的Stokes函数可以看出，在角距离较小时（近区），函数变化比较剧烈，因此在积分时，近区的步长取得比较小，而远区的步长可以大一些。在计算中，具体采用的步长见表1。

由于我们假设的重力异常分布与经度无关，显然与经度有关的积分值为2p，因此在计算结果的比较中我们忽略这一常数，也就是只对纬度进行积分，另外我们也忽略（9）式中求和符号前的常数，这并不影响结果的比较，计算结果在表2中列出。

由表中的结果可以看出，积分Stokes方法计算结果的精度比普通的数值积分方法的精度提高了1个数量级，这在确定高精度的大地水准面的计算中就显得尤为重要了。另外我们减小了近区的积分步长，但是对结果的影响甚微，表明所采用的步长已经可以达到最佳的

效果。

3结语

我们提出了一种计算大地水准面差距的新方法—积分Stokes方法，该方法显著的优点是去除了计算点处的奇异性，并且在计算效率上也优于普通的数值积分方法，该方法使得计算结果的精度提高了一个数量级，这在实际的工作中将发挥重要的作用。

参考文献

[1] 管泽霖，管铮，黄谟涛，翟国君．局部重力场逼近理

论和方法[M]．北京：测绘出版社，1997．

[2] Goad，C，C．Gravimetric tidal loading computed from

integrated Greens function，J．G．Res．，85，

2679-2683，1980．

[3] Wu Qingpeng， Zhu Wenlu．， The load integrated

Greens function （LIGF） for spherical elastic earth

models and their applications． Proceedings of the

twelfth international symposium on earth tides．Science

press， Beijing， New York， 1995．

[4] 吴庆鹏．重力学与固体潮[M]．北京：地震出版社，

1997．

endprint

（济南市勘察测绘研究院，山东 济南 250013）

摘要：文章介绍了计算求解大地水准面的Stokes公式的一种新方法—积分Stokes方法。给出了积分Stokes函数，通过该函数计算大地水准面时，在计算点处不存在奇异性。作为示例文章采用了连续分布的重力异常函数进行了试算，结果表明，利用该方法计算所得的结果与普通的数值积分方法比较具有较高的精度。

关键词：大地水准面；Stokes公式积分；Stokes函数；计算点处

中图分类号：P223文献标识码：A文章编号：1009-2374（2014）22-0061-02确定地球形状是大地测量工作者的一项基本任务，通常采用大地水准面来表征地球形状，它是覆盖全球的一个重力等位面，其上的位与正常椭球表面的位相等。1849年，Stokes导出了计算大地水准面的积分公式—Stokes公式，为确定地球形状提供了有力的工具。由于该公式是对全球进行积分，围绕计算的精度与速度，不同学者进行了广泛的研究，提出了许多计算方法，如远近区相结合的方法和谱分析方法，以及Stokes函数的改化。1980年Goad在计算海潮负荷时提出了积分格林函数方法，我国学者吴庆鹏也对此做了更进一步的研究，他们的研究表明该方法具有显著的优点，特别是能够使得计算点处的奇异性降低1阶。本文提出的方法即是将这一思想运用于由重力异常分布确定大地水准面的计

算中。

1计算理论

Stokes公式具有非常简洁的形式，

（1）

其中是地球半径，是平均椭球体表面的正常重力的平均值，是重力异常，被称为Stokes函数，它是测站与计算点的角距的函数（如图1），有如下的公式：

（2）

其中：

（3）

分别为计算点的余纬与经度，分别为测站的余纬与经度，表示半径为的单位球，且有：

（4）

上式中，A为方位角。

图1 Stokes函数

从（1）式可以看出，计算大地水准面差距也就是计算重力异常与Stokes函数的乘积在整个球面上的积分，因此在某种意义上，Stokes函数可以认为是一种格林函数，即任一测站的重力异常对大地水准面差距的贡献的响应函数。以下我们所称的格林函数即为Stokes

函数。

我们定义积分格林函数为：

（5）

δ为积分Stokes函数所处的区间间隔大小，即以该区间的中点代替积分Stokes函数的自变量。通过推导可得积分Stokes函数的解析表达式为：

（6）

其中：

（7）

显然，在计算点处，y=0，即t=0，由（2）知，Stokes函数表达式中的分母为0，因此产生计算时产生奇异，有许多文献讨论了避免奇异性的方法，如将该点处的重力异常展开为泰勒级数。

而（6）中包含对数函数的部分的值为

（8）

上式中运用了求极限的罗必塔法则。因此，积分Stokes函数在计算点是没有奇异性的。图2给出了其随角距离的变化，可以看出，积分Stokes函数的变化是比较平缓的，而不像Stokes函数那样在近区的变化比较剧烈（注意图1中横轴取了对数），图中红色的线表示由（6）式计算的结果，而黑色的线表示由（5）式计算的结果（积分步长为0.001°），很明显二者具有非常好的一致性，说明（6）式是正确的。

图2 积分Stokes函数

令DA为沿方位角方向的网格划分，那么数值积分时，认为重力异常在任一网格内为常数，则将（5）式代入（1）式可得：

（9）

2试验算例

由于我们仅仅考虑我们提出的方法的优劣，因此采用有解析结果的算例进行验证。本试验中，令重力异常在全球的分布满足

（10）

那么对于北极点（y=q）的大地水准面差距来说，其具有严格的解析解：

（11）

这样我们就可以用来检核不同方法计算结果的正确性与精度。对于北极此时我们就可以利用重力异常分布的对称性，并且在这两点的大地水准面差距具有严格的解析解。

从图1中的Stokes函数可以看出，在角距离较小时（近区），函数变化比较剧烈，因此在积分时，近区的步长取得比较小，而远区的步长可以大一些。在计算中，具体采用的步长见表1。

由于我们假设的重力异常分布与经度无关，显然与经度有关的积分值为2p，因此在计算结果的比较中我们忽略这一常数，也就是只对纬度进行积分，另外我们也忽略（9）式中求和符号前的常数，这并不影响结果的比较，计算结果在表2中列出。

由表中的结果可以看出，积分Stokes方法计算结果的精度比普通的数值积分方法的精度提高了1个数量级，这在确定高精度的大地水准面的计算中就显得尤为重要了。另外我们减小了近区的积分步长，但是对结果的影响甚微，表明所采用的步长已经可以达到最佳的

效果。

3结语

我们提出了一种计算大地水准面差距的新方法—积分Stokes方法，该方法显著的优点是去除了计算点处的奇异性，并且在计算效率上也优于普通的数值积分方法，该方法使得计算结果的精度提高了一个数量级，这在实际的工作中将发挥重要的作用。

参考文献

[1] 管泽霖，管铮，黄谟涛，翟国君．局部重力场逼近理

论和方法[M]．北京：测绘出版社，1997．

[2] Goad，C，C．Gravimetric tidal loading computed from

integrated Greens function，J．G．Res．，85，

2679-2683，1980．

[3] Wu Qingpeng， Zhu Wenlu．， The load integrated

Greens function （LIGF） for spherical elastic earth

models and their applications． Proceedings of the

twelfth international symposium on earth tides．Science

press， Beijing， New York， 1995．

[4] 吴庆鹏．重力学与固体潮[M]．北京：地震出版社，

1997．

endprint

（济南市勘察测绘研究院，山东 济南 250013）

摘要：文章介绍了计算求解大地水准面的Stokes公式的一种新方法—积分Stokes方法。给出了积分Stokes函数，通过该函数计算大地水准面时，在计算点处不存在奇异性。作为示例文章采用了连续分布的重力异常函数进行了试算，结果表明，利用该方法计算所得的结果与普通的数值积分方法比较具有较高的精度。

关键词：大地水准面；Stokes公式积分；Stokes函数；计算点处

中图分类号：P223文献标识码：A文章编号：1009-2374（2014）22-0061-02确定地球形状是大地测量工作者的一项基本任务，通常采用大地水准面来表征地球形状，它是覆盖全球的一个重力等位面，其上的位与正常椭球表面的位相等。1849年，Stokes导出了计算大地水准面的积分公式—Stokes公式，为确定地球形状提供了有力的工具。由于该公式是对全球进行积分，围绕计算的精度与速度，不同学者进行了广泛的研究，提出了许多计算方法，如远近区相结合的方法和谱分析方法，以及Stokes函数的改化。1980年Goad在计算海潮负荷时提出了积分格林函数方法，我国学者吴庆鹏也对此做了更进一步的研究，他们的研究表明该方法具有显著的优点，特别是能够使得计算点处的奇异性降低1阶。本文提出的方法即是将这一思想运用于由重力异常分布确定大地水准面的计

算中。

1计算理论

Stokes公式具有非常简洁的形式，

（1）

其中是地球半径，是平均椭球体表面的正常重力的平均值，是重力异常，被称为Stokes函数，它是测站与计算点的角距的函数（如图1），有如下的公式：

（2）

其中：

（3）

分别为计算点的余纬与经度，分别为测站的余纬与经度，表示半径为的单位球，且有：

（4）

上式中，A为方位角。

图1 Stokes函数

从（1）式可以看出，计算大地水准面差距也就是计算重力异常与Stokes函数的乘积在整个球面上的积分，因此在某种意义上，Stokes函数可以认为是一种格林函数，即任一测站的重力异常对大地水准面差距的贡献的响应函数。以下我们所称的格林函数即为Stokes

函数。

我们定义积分格林函数为：

（5）

δ为积分Stokes函数所处的区间间隔大小，即以该区间的中点代替积分Stokes函数的自变量。通过推导可得积分Stokes函数的解析表达式为：

（6）

其中：

（7）

显然，在计算点处，y=0，即t=0，由（2）知，Stokes函数表达式中的分母为0，因此产生计算时产生奇异，有许多文献讨论了避免奇异性的方法，如将该点处的重力异常展开为泰勒级数。

而（6）中包含对数函数的部分的值为

（8）

上式中运用了求极限的罗必塔法则。因此，积分Stokes函数在计算点是没有奇异性的。图2给出了其随角距离的变化，可以看出，积分Stokes函数的变化是比较平缓的，而不像Stokes函数那样在近区的变化比较剧烈（注意图1中横轴取了对数），图中红色的线表示由（6）式计算的结果，而黑色的线表示由（5）式计算的结果（积分步长为0.001°），很明显二者具有非常好的一致性，说明（6）式是正确的。

图2 积分Stokes函数

令DA为沿方位角方向的网格划分，那么数值积分时，认为重力异常在任一网格内为常数，则将（5）式代入（1）式可得：

（9）

2试验算例

由于我们仅仅考虑我们提出的方法的优劣，因此采用有解析结果的算例进行验证。本试验中，令重力异常在全球的分布满足

（10）

那么对于北极点（y=q）的大地水准面差距来说，其具有严格的解析解：

（11）

这样我们就可以用来检核不同方法计算结果的正确性与精度。对于北极此时我们就可以利用重力异常分布的对称性，并且在这两点的大地水准面差距具有严格的解析解。

从图1中的Stokes函数可以看出，在角距离较小时（近区），函数变化比较剧烈，因此在积分时，近区的步长取得比较小，而远区的步长可以大一些。在计算中，具体采用的步长见表1。

由于我们假设的重力异常分布与经度无关，显然与经度有关的积分值为2p，因此在计算结果的比较中我们忽略这一常数，也就是只对纬度进行积分，另外我们也忽略（9）式中求和符号前的常数，这并不影响结果的比较，计算结果在表2中列出。

由表中的结果可以看出，积分Stokes方法计算结果的精度比普通的数值积分方法的精度提高了1个数量级，这在确定高精度的大地水准面的计算中就显得尤为重要了。另外我们减小了近区的积分步长，但是对结果的影响甚微，表明所采用的步长已经可以达到最佳的

效果。

3结语

我们提出了一种计算大地水准面差距的新方法—积分Stokes方法，该方法显著的优点是去除了计算点处的奇异性，并且在计算效率上也优于普通的数值积分方法，该方法使得计算结果的精度提高了一个数量级，这在实际的工作中将发挥重要的作用。

参考文献

[1] 管泽霖，管铮，黄谟涛，翟国君．局部重力场逼近理

论和方法[M]．北京：测绘出版社，1997．

[2] Goad，C，C．Gravimetric tidal loading computed from

integrated Greens function，J．G．Res．，85，

2679-2683，1980．

[3] Wu Qingpeng， Zhu Wenlu．， The load integrated

Greens function （LIGF） for spherical elastic earth

models and their applications． Proceedings of the

twelfth international symposium on earth tides．Science

press， Beijing， New York， 1995．

[4] 吴庆鹏．重力学与固体潮[M]．北京：地震出版社，

1997．

endprint