物理大地测量学笔记（二）

4 边值问题

边值问题，即在边界面上已知某些函数值，且函数值满足一定条件，据此来求的外部或内部的调和、并且在无穷远处正则的函数。

根据给定的边值条件不同，分为三种边值问题：

• 第一边值问题：Dirichlet问题——已知边界面上引力位 $V$
• 第二边值问题：Neumann问题——已知边界面上引力位的法向导数 $\frac{\partial V}{\partial n}$
• 第三边值问题：Robin问题——已知边界面上引力位法向导数和引力位的线性组合 $hV+h\frac{\partial V}{\partial n}$

4.1 格林公式求解边值问题

4.1.1 格林公式

内部第一格林公式

记$D(U,V)$为：
$$D(U,V)=\frac{\partial U}{\partial x}\frac{\partial V}{\partial x}+\frac{\partial U}{\partial y}\frac{\partial V}{\partial y}+\frac{\partial U}{\partial z}\frac{\partial V}{\partial z}$$
则内部第一格林公式表示为：
$$\underset{v}{\iiint }U\Delta Vdv+\underset{v}{\iiint }D(U,V)dv=\underset{s}{\iint }U\frac{\partial V}{\partial n}dS$$

内部第二格林公式

将内部第一格林公式中的 $U,V$ 互换：
$$\underset{v}{\iiint }V\Delta Udv+\underset{v}{\iiint }\left(\frac{\partial U}{\partial x}\frac{\partial V}{\partial x}+\frac{\partial U}{\partial y}\frac{\partial V}{\partial y}+\frac{\partial U}{\partial z}\frac{\partial V}{\partial z}\right)dv=\underset{s}{\iint }V\frac{\partial U}{\partial n}dS$$
和第一恒等式相减：
$$\underset{v}{\iiint }\left(U\Delta V-V\Delta U\right)dv=\underset{s}{\iint }\left(U\frac{\partial V}{\partial n}-V\frac{\partial U}{\partial n}\right)dS$$
即得到内部第二格林公式

外部第一格林公式

在外部时，为了使得 $n$ 仍然为向外的法线，需要将 $\frac{\partial }{\partial n}$ 的符号反号，即得到：
$$\underset{v}{\iiint }U\Delta Vdv+\underset{v}{\iiint }D(U,V)dv=-\underset{s}{\iint }U\frac{\partial V}{\partial n}dS$$

外部第二格林公式

同理得：
$$\underset{v}{\iiint }\left(U\Delta V-V\Delta U\right)dv=-\underset{s}{\iint }\left(U\frac{\partial V}{\partial n}-V\frac{\partial U}{\partial n}\right)dS$$

4.1.2 格林公式的应用

4.1.2.1 计算地球质量

令 $U=1$ ，代入内部第二格林公式：
$$\underset{v}{\iiint }\Delta Vdv=\underset{s}{\iint }\frac{\partial V}{\partial n}dS\text{\hspace{1em}(1)}$$
令 $V=W(重力位)$，则有
$$\Delta W=-4\pi G\rho +2{\omega }^2$$
设 $v$ 为地球体积，$S$ 为地球实际表面，则$\frac{\partial W}{\partial n}=-g_n$，它是地球表面 $S$ 的法线方向重力分量。

将以上条件代入（1）式：
$$\begin{gathered} \underset{v}{\iiint }\left(-4\pi G\rho +2{\omega }^2\right)dv=-\underset{s}{\iint }{g}_{n}dS \\ \Rightarrow -4\pi G\underset{v}{\iiint }\rho dv+\underset{v}{\iiint }2{\omega }^{2}dv=-\underset{s}{\iint }{g}_{n}dS \\ \Rightarrow -4\pi GM+2{\omega }^{2}v=-\underset{s}{\iint }{g}_{n}dS \\ \Rightarrow M=\frac{1}{4\pi G}\left(\underset{s}{\iint }{g}_{n}dS+2{\omega }^{2}v\right) \end{gathered}$$
通过这一式子，即可根据地球表面重力来计算地球质量。

4.1.3 第一边值问题（Dirichlet问题）

令
$$U=\frac{1}{l}，V为质体内部的位函数$$
代入内部第二格林公式：
$$\underset{v}{\iiint }\left(\frac{1}{l}\Delta V-V\Delta \frac{1}{l}\right)dv=\underset{s}{\iint }\left(\frac{1}{l}\frac{\partial V}{\partial n}-V\frac{\partial }{\partial n}\left(\frac{1}{l}\right)\right)dS$$
之前证明过， $\frac{1}{l}$是一个谐函数，即$\Delta \frac{1}{l}=0$。又根据泊松公式：$\Delta V=-4\pi G\rho$ 。因此化简得：
$$-4\pi G\underset{v}{\iiint }\frac{\rho }{l}dv=\underset{s}{\iint }\left(\frac{1}{l}\frac{\partial V}{\partial n}-V\frac{\partial }{\partial n}\left(\frac{1}{l}\right)\right)dS$$
根据位函数基本公式：
$$V=G\underset{v}{\iiint }\frac{\rho }{l}dv$$
代入得：
$$\underset{s}{\iint }\left(\frac{1}{l}\frac{\partial V}{\partial n}-V\frac{\partial }{\partial n}\left(\frac{1}{l}\right)\right)dS=-4\pi V\text{\hspace{1em}(1)}$$
又令 $u$ 为任意一个在外部调和、无穷远处正则的函数，代入外部第二格林公式得：
$$\underset{v}{\iiint }\left(u\Delta V-V\Delta u\right)dv=-\underset{s}{\iint }\left(u\frac{\partial V}{\partial n}-V\frac{\partial u}{\partial n}\right)dS$$
在外部 $V$ 同样是调和函数，因此式子左边为0，即：
$$\underset{s}{\iint }\left(u\frac{\partial V}{\partial n}-V\frac{\partial u}{\partial n}\right)dS=0\text{\hspace{1em}(2)}$$
（1）（2）两式相减得：
$$\underset{s}{\iint }\left(\left(u-\frac{1}{l}\right)\frac{\partial V}{\partial n}-V\frac{\partial }{\partial n}\left(u-\frac{1}{l}\right)\right)dS=4\pi V\text{\hspace{1em}(3)}$$
引入格林函数：
$$f_G=u-\frac{1}{l}$$
假定 $u$ 在 $S$ 面上等于 $\frac{1}{l}$，则代入（3）式可得：
$$4\pi V=-\underset{s}{\iint }\left(V\frac{\partial f_G}{\partial n}\right)dS\text{\hspace{1em}(3)}$$
可以看到，如果已知物体表面上 $v$ 的值、边界面以及 $f_G$ 的值，可以求出外部引力位。

4.1.2 第二边值问题（Neumann问题)

类似地，假定 $u$ 使得格林函数在边界面上的法向导数为0，即：
$$\frac{\partial f_G}{\partial n}{|}_s=0$$
依然是代入（3）式：
$$4\pi V=\underset{s}{\iint }\left(\frac{\partial v}{\partial n}{f}_G\right)dS$$
如果已知物体表面上 $\frac{\partial v}{\partial n}$ 的值、边界面以及 $f_G$ 的值，可以求出外部引力位。

4.2.3 第三边值问题（Robin问题）

一样的方法，根据边值条件来假定u：
$$\left[\alpha \left(u-\frac{1}{l}\right)+\beta \frac{\partial }{\partial n}\left(u-\frac{1}{l}\right)\right]{|}_s=0$$
这里推导较为繁琐，懒得打了，直接从PPT里截了个图

4.2 球函数求解边值问题

4.2.1 第一边值问题（Dirichlet问题）

假设边界面是一个单位球面，位函数表示为$V(1,\theta ,\lambda )$，此时 $r=1$，则可以将其展开为面球谐函数：
$$V(1,\theta ,\lambda )=\sum _{n=0}^{\infty }{Y}_n(\theta ,\lambda )$$
在球面之外，容易得到函数解：
$$V_e(r,\theta ,\lambda )=\sum _{n=0}^{\infty }\frac{1}{r^{n+1}}{Y}_n(\theta ,\lambda )$$
球面之内则有函数解：
$$V_i(r,\theta ,\lambda )=\sum _{n=0}^{\infty }{r}^n{Y}_n(\theta ,\lambda )$$
这几个函数都收敛，并且都是谐函数，因此满足条件。

那么对于半径为 $R$ 的球体，同样先展开为：
$$V(R,\theta ,\lambda )=\sum _{n=0}^{\infty }{Y}_n(\theta ,\lambda )$$
书上直接给出了面谐函数的公式（也可以通过分解公式代入求得，但我不会）：
$$Y_n(\theta ,\lambda )=\frac{2n+1}{4\pi }{\int }_{{\lambda }^{\prime }=0}^{2\pi }{\int }_{{\theta }^{\prime }=0}^{\pi }V(R,{\theta }^{\prime },{\lambda }^{\prime })P_n(\cos \psi )\sin {\theta }^{\prime }d{\theta }^{\prime }d{\lambda }^{\prime }$$
这里可以简单理解成将单位球面等比例缩放，将原来的$1$替换为 $\frac{r}{R}$ ，解得面外和面内的函数为：
$$\begin{gathered} V_e(r,\theta ,\lambda )=\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{Y}_n(\theta ,\lambda ) \\ {V}_i(r,\theta ,\lambda )=\sum _{n=0}^{\infty }{\left(\frac{r}{R}\right)}^n{Y}_n(\theta ,\lambda ) \end{gathered}$$
仅考虑面外部分，代入面谐函数公式：
$$\begin{gathered} V_e(r,\theta ,\lambda )=\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}\frac{2n+1}{4\pi }{\int }_{{\lambda }^{\prime }=0}^{2\pi }{\int }_{{\theta }^{\prime }=0}^{\pi }V(R,{\theta }^{\prime },{\lambda }^{\prime })P_n(\cos \psi )\sin {\theta }^{\prime }d{\theta }^{\prime }d{\lambda }^{\prime } \\ =\frac{1}{4\pi }{\int }_{{\lambda }^{\prime }=0}^{2\pi }{\int }_{{\theta }^{\prime }=0}^{\pi }V(R,{\theta }^{\prime },{\lambda }^{\prime })\left[\sum _{n=0}^{\infty }(2n+1){\left(\frac{R}{r}\right)}^{n+1}{P}_n(\cos \psi )\right]\sin {\theta }^{\prime }d{\theta }^{\prime }d{\lambda }^{\prime } \end{gathered}$$
接下来对方括号内的值进行化简：
$$\sum _{n=0}^{\infty }(2n+1){\left(\frac{R}{r}\right)}^{n+1}{P}_n(\cos \psi )$$
根据距离倒数的展开式：
$$\frac{1}{l}=\frac{1}{\sqrt{r^2+R^2-2rR\cos \psi }}=\sum _{n=0}^{\infty }\frac{R^n}{r^{n+1}}{P}_n\left(\cos \psi \right)\text{\hspace{1em}(1)}$$
左右同时对 $r$ 微分得：
$$\begin{gathered} -\frac{1}{2}\frac{2r-2R\cos \psi }{{\left(\sqrt{r^2+R^2-2rR\cos \psi }\right)}^3}=-(n+1)\sum _{n=0}^{\infty }\frac{R^n}{r^{n+2}}{P}_n\left(\cos \psi \right) \\ \Rightarrow \frac{r-R\cos \psi }{l^3}=(n+1)\sum _{n=0}^{\infty }\frac{R^n}{r^{n+2}}{P}_n\left(\cos \psi \right)\text{\hspace{1em}(2)} \end{gathered}$$
为了构造出$(2n+1)$的形式，并且让 $r,R$ 的次数都为 $n+1$，

将$(1)$式乘以$-R$ 得：
$$-\frac{R}{l}=-\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{P}_n\left(\cos \psi \right)$$
将$(2)$式乘以 $2rR$ 得：
$$\frac{2rR(r-R\cos \psi )}{l^3}=(2n+2)\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{P}_n\left(\cos \psi \right)$$
两式相加得：
$$\begin{gathered} \frac{-R{l}^2+2rR(r-R\cos \psi )}{l^3}=(2n+1)\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{P}_n\left(\cos \psi \right) \\ \Rightarrow \frac{-R(r^2+R^2-2rR\cos \psi )+2rR(r-R\cos \psi )}{l^3}=(2n+1)\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{P}_n\left(\cos \psi \right) \\ \Rightarrow \frac{R(r^2-R^2)}{l^3}=(2n+1)\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}{P}_n\left(\cos \psi \right) \end{gathered}$$
至此就已经完成了化简。将式子代入：
$$V_e(r,\theta ,\lambda )=\frac{R(r^2-R^2)}{4\pi }{\int }_{{\lambda }^{\prime }=0}^{2\pi }{\int }_{{\theta }^{\prime }=0}^{\pi }\frac{V(R,{\theta }^{\prime },{\lambda }^{\prime })}{l^3}\sin {\theta }^{\prime }d{\theta }^{\prime }d{\lambda }^{\prime }$$
这就是泊松积分式。也是Dirichlet问题的直接解。

4.2.2 第二边值问题（Neumann问题)

对于第二边值问题，基本思路和第一边值问题没有区别。简单来说，给定了什么函数值，就展开什么函数值。

第二边值问题中给定 $S$ 面上法向导数$\frac{\partial V}{\partial n}$，所以先将法向导数展开为面球谐函数：
$${\left(\frac{\partial V}{\partial n}\right)}_{r=R}=\sum _{n=0}^{\infty }{Y}_n(\theta ,\lambda )$$
同时注意到，在球面上法向导数和径向梯度是相等的，即：
$${\left(\frac{\partial V}{\partial n}\right)}_{r=R}={\left(\frac{\partial V}{\partial r}\right)}_{r=R}$$
因此Neumann问题变成了寻找一个表达式，使得它对 $r$ 微分后等于面球谐函数展开式。

直接给出结果：
$$V_e(r,\theta ,\lambda )=-R\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}\frac{Y_n(\theta ,\lambda )}{n+1}$$
其正确性很容易通过求微分来验证。

4.2.3 第三边值问题（Robin问题）

第三边值问题可以用于通过重力异常确定大地水准面起伏。

仍然是将给定条件展开：
$$hV+k\frac{\partial V}{\partial n}=\sum _{n=0}^{\infty }{Y}_n(\theta ,\lambda )$$
直接给出结果：
$$V_e(r,\theta ,\lambda )=\sum _{n=0}^{\infty }{\left(\frac{R}{r}\right)}^{n+1}\frac{Y_n(\theta ,\lambda )}{h-(k/R)(n+1)}$$

4.3 Stokes定理和Stokes问题

4.3.1 Stokes定理

若已知⼀个水准面形状S、S面上的位W0（或它内部所包含的物质的总质量M）及该物体绕某⼀固定轴的旋转⻆速度ω，则该水准面上及其外部空间任意点的重力位都可唯一确定，并且不需要知道物体内的质量分布情况。

4.3.2 Stokes问题

已知水准面上的重力 $g$ 和重⼒位 $W$（或地球的总质量 $M$），以及地球的自转角速度 $\omega$，需求定水准面的形状 $S$ 及其外部的重力位。

4.4 边值问题解的唯一性

外部第一格林公式：
$$\underset{v_e}{\int }\left[u_e\Delta v_e+D\left(u_e,v_e\right)\right]dv=-\underset{\sigma }{\int }{u}_e\frac{\partial v_e}{\partial n}d\sigma$$
假设解不是唯一的，则在 $\sigma$ 外有两个调和并在无穷远处正则的函数 $V_1$ 和 $V_2$ 同时满足 $\sigma$ 上的边界条件，记 $T=V_1-V_2$ ，则 $T$ 也是在 $\sigma$ 外调和且在无穷远处正则的函数（$\Delta T\equiv 0$ ）。

将外部第一格林公式应用于 $T$，并在该式中令 $u=v=T$ ，则有
$$\underset{v_e}{\int }\left[T\Delta T+D\left(T,T\right)\right]dv=-\underset{\sigma }{\int }T\frac{\partial T}{\partial n}d\sigma$$
根据条件，在 $v_e$ 内 $\Delta T\equiv 0$，在 $\sigma$ 面上 $T_{\sigma }^{}=0$，则上式为：
$$\underset{v_e}{\int }D\left(T,T\right)dv=0\Rightarrow \underset{v_e}{\int }\left[{\left(\frac{\partial T}{\partial x}\right)}^2+{\left(\frac{\partial T}{\partial y}\right)}^2+{\left(\frac{\partial T}{\partial z}\right)}^2\right]dv=0$$
要使上式成立，则必须在 $\sigma$ 外任意点上都满足
$$\frac{\partial T}{\partial x}=\frac{\partial T}{\partial y}=\frac{\partial T}{\partial z}=0$$
即 $T$ 为常数。又因为 $T$ 是在无穷远处的正则函数，则 $T$ 只能等于零，亦即 $V_1=V_2$ ，证明解是唯一的。

5 地球重力场

5.1 重力和重力位

5.1.1 重力位的推导

第二章讨论了引力和引力位，而重力是引力和离心力的合力。因此先推导离心力和离心力位：

易得 向量 $p=(x,y,0)$，因此有离心力 $f$ 等于：
$$f={\omega }^{2}p=({\omega }^{2}x,{\omega }^{2}y,0)$$
离心力位应满足梯度等于离心力，直接给出离心力位的表达式：
$$\Phi =\frac{1}{2}{\omega }^2(x^2+y^2)$$
则重力位表达为：
$$W=V+\Phi =G\underset{v}{\iiint }\frac{\rho }{l}dv+\frac{1}{2}{\omega }^2(x^2+y^2)$$
还是老套路，求一下重力位的 $Laplace$ 方程值：
$$\Delta \Phi =\frac{\partial \Phi }{\partial x^2}+\frac{\partial \Phi }{\partial y^2}+\frac{\partial \Phi }{\partial z^2}=2{\omega }^2$$
由 $Possion$ 方程（这里考虑质体内部）
$$\Delta V=-4\pi G\rho$$
得到重力位的拉普拉斯方程值：
$$\Rightarrow \Delta W=\Delta V+\Delta \Phi =-4\pi G\rho +2{\omega }^2$$
该式称为广义 Poisson 公式。
比较明显的结论是：重力位在外部连续但不调和（因为离心力的存在）

5.1.2 重力的大小和方向

重力单位为伽 $(gal)$，在赤道上大小约为 $978gal$，在两极大小约为$983gal$。之前更常用的单位是牛顿$(N)$，其中$$

重力的方向称为铅垂方向，和该点所在的重力等位面垂直，因此铅垂线呈现为不规则的曲线

5.1.3 水准面

水准面也就是重力等位面，在水准面上，重力大小处处相等，即 $W(x,y,z)=W_0$。

大地水准面是特殊的水准面，它和全球无潮平均海水面密合，作为海拔的起算面。

虽然重力和等位面垂直似乎是个显而易见的结论，还是简单证明一下：

重力位的梯度，等于重力在该方向的分量，用微分式来表示：
$$\begin{array}{c}dW=\frac{\partial W}{\partial x}dx+\frac{\partial W}{\partial y}dy+\frac{\partial W}{\partial z}dz\\ \\ =\left(\frac{\partial W}{\partial x},\frac{\partial W}{\partial y},\frac{\partial W}{\partial z}\right)\cdot {\left(dx,dy,dz\right)}^T\\ \\ l=\mathbf{g}\cdot \mathrm{d}\mathbf{x}\end{array}$$
令 $dx$ 和水准面相切，$W$为常数，则 $dW=0\Rightarrow \mathbf{g}=0$，说明重力在等位面的切线方向没有分量，因此重力方向垂直于等位面。

5.1.4 正高

从大地水准面起，沿铅垂线方向至某点的距离成为该点的正高。

从前面的经验来看，无论什么物理量，最终都归结到位函数上。正高也不例外。

沿铅垂线增高的方向取矢量 $d\mathbf{x}$，则有$\left|dx\right|=dH$，$dH$ 和 重力 $\mathbf{g}$ 方向相反。对重力位取微分，则有：
$$dW=\mathbf{g}\cdot d\mathbf{x}=\mathbf{g}\cdot dH\cdot \cos (\mathbf{g},d\mathbf{x})=-\mathbf{g}dH$$
或表示为：

$$\{\begin{array}{c}g=-\frac{\partial W}{\partial H}\\ \\ dH=-\frac{dW}{g}\end{array}$$
这就建立了重力位和正高之间的数学联系。

从这个式子可以得出水准面的一些重要性质：

• 由于重力随纬度有变大的趋势，因此两个水准面之间的距离随纬度有缩小的趋势
• 同一水准面上重力位相等，但重力大小不相等，因此 两个水准面之间的$dH$ 不为常数，也就是水准面之间不平行
• 因为重力 $\mathbf{g}$ 大小不为0，因此 $dH$ 不为0，即水准面之间不相交

5.1.5 水准面曲率

上一部分建立了位差和高差的数学关系，在几何概念和物理概念之间进行了转换。这一部分则介绍了Bruns公式的推导，它建立了垂直重力梯度和水准面平均曲率的关系。

首先对于曲线 $y=f(x)$，有曲率公式：
$$\kappa =\frac{1}{\rho }=\frac{y^{\prime \prime }}{(1+y^{\prime 2}{)}^{3/2}}$$
当切线平行于$x$轴时，$y^{\prime }=0$，则有：
$$\kappa =\frac{1}{\rho }=y^{\prime \prime }=\frac{d^{2}y}{d{x}^2}$$
为了求水准面的平均曲率，首先建立空间直角坐标系：

其中，原点为 $P$ ，在大地水准面上，图中画出了 $x-z$ 平面和水准面的交线和铅垂线，$x-y$ 平面实际上就是水准面在 $P$ 点的切面。

考虑此时 $x-z$ 平面内的曲率，根据上面结论，就等于二阶导。将 $W(x,y,z)=W_0$ 对 $x$ 微分得：
$$W_x+W_z\frac{dz}{dx}=0$$
二次微分得（涉及复合函数求微分）：
$$W_{xx}+(W_{xz}+W_{zz}\frac{dz}{dx})\frac{dz}{dx}+W_z\frac{d^{2}z}{d{x}^2}=0$$
因为在P点有 $\frac{dz}{dx}=0$，代入二次微分式得：
$$\frac{d^{2}z}{d{x}^2}=\frac{W_{xx}}{W_z}=K_1$$
由于 $z$ 轴在 $P$ 点垂直于水准面，$W_z=\frac{\partial W}{\partial z}=\frac{\partial W}{\partial H}=-g$，得到交线的曲率为：
$$K_1=-\frac{W_{xx}}{g}$$
同理，水准面和 $y-z$平面的交线的曲率为
$$K_2=-\frac{W_{yy}}{g}$$
定义水准面上 $P$ 点的曲率为两条交线的曲率平均数，则有
$$J=\frac{1}{2}(K_1+K_2)=-\frac{W_{xx}+W_{yy}}{2g}\text{\hspace{1em}(1)}$$
根据前面推导的广义 $Poission$ 公式
$$\Delta W=W_{xx}+W_{yy}+W_{zz}=-4\pi G\rho +2{\omega }^2\text{\hspace{1em}(2)}$$
其中考虑到：
$$W_{zz}=\frac{\partial W_z}{\partial z}=-\frac{\partial g}{\partial z}=-\frac{\partial g}{\partial H}\text{\hspace{1em}(3)}$$
联立(1)(2)(3)三式，可以表示出重力梯度公式：
$$\frac{\partial g}{\partial H}=-2gJ+4\pi G\rho -2{\omega }^2$$
上式建立了垂直重力梯度和水准面平均曲率之间的联系。

5.2 引力位的球谐展开

5.2.1 引力位的球谐展开式

回忆一下母函数：
$$\frac{1}{l}=\sum _{n=0}^{\infty }\frac{r^{\prime n}}{r^{n+1}}{P}_n(\cos \psi )$$
根据分解公式：
$$P_n\left(\cos \psi \right)=P_n\left(\cos \theta \right)P_n\left(\cos {\theta }^{\prime }\right)+2\sum _{m=1}^n\frac{\left(n-m\right)!}{\left(n+m\right)!}\left[R_{nm}\left(\theta ,\lambda \right)R_{nm}\left({\theta }^{\prime },{\lambda }^{\prime }\right)+S_{nm}\left(\theta ,\lambda \right)S_{nm}\left({\theta }^{\prime },{\lambda }^{\prime }\right)\right]$$
联立得：
$$\frac{1}{l}=\sum _{n=0}^{\infty }\frac{r^{\prime n}}{r^{n+1}}\left\{P_n\left(\cos \theta \right)P_n\left(\cos {\theta }^{\prime }\right)+2\sum _{m=1}^n\frac{\left(n-m\right)!}{\left(n+m\right)!}\left[R_{nm}\left(\theta ,\lambda \right)R_{nm}\left({\theta }^{\prime },{\lambda }^{\prime }\right)+S_{nm}\left(\theta ,\lambda \right)S_{nm}\left({\theta }^{\prime },{\lambda }^{\prime }\right)\right]\right\}$$
将第n阶引力位用勒让德函数来表示：
$$\begin{array}{c}{V}_n=\frac{1}{r^{n+1}}\left[A_n{P}_n(\cos \theta )+\sum _{m=1}^n(A_{nm}\cos m\lambda +B_{nm}\sin m\lambda )P_{nm}(\cos \theta )\right]\\ \\ \text{或者：}{V}_n=\frac{1}{r^{n+1}}\sum _{m=0}^n\left[a_{nm}{P}_{nm}\left(\cos \theta \right)\cos m\lambda +b_{nm}{P}_{nm}\left(\cos \theta \right)\sin m\lambda \right]\end{array}$$
两种表达方式没有区别，因为 $m=0$ 时 $\sin m\lambda$ 为0，也就不存在对应的系数 $B_n$

对于式中的几个参数，有几个基本的概念：

• $P_n(\cos \theta )$：主球谐函数，或带球谐函数，即勒让德多项式
• $P_{nm}(\cos \theta )$：缔和勒让德函数
• $P_{nm}(\cos \theta )\cos m\lambda$ 和 $P_{nm}(\cos \theta )\sin m\lambda$：缔和球谐函数，或面球谐函数
• $A_n,A_{nm},B_{nm}$：Stokes系数

这里给出Stokes系数的计算公式（可利用面谐函数正交性求解）：
$$\begin{array}{c}{A}_{n0}=G{\iiint }_{earth}{r}^{\prime n}{P}_n(\cos \theta )dM=G{\iiint }_{earth}{z}^{\prime }dM\\ \\ A_{nm}=2\frac{(n-m)!}{(n+m)!}G{\iiint }_{earth}{r}^{\prime n}{R}_{nm}({\theta }^{\prime },{\lambda }^{\prime })dM=G{\iiint }_{earth}{z}^{\prime }dM\\ \\ B_{nm}=2\frac{(n-m)!}{(n+m)!}G{\iiint }_{earth}{r}^{\prime n}{S}_{nm}({\theta }^{\prime },{\lambda }^{\prime })dM=G{\iiint }_{earth}{z}^{\prime }dM\end{array}$$

5.2.2 矩

矩是物理中的概念，这里引入是为了表达引力位展开式低阶项的物理意义。

定义：质体的质量与（到某点、某轴或某平面）距离d的k次方的乘积的物理量统称为矩。

将质体的k阶矩表示为：
$$\underset{v}{\iiint }{d}^{k}dm$$
值得注意的是，这里的距离d可以是到点或线或面

而低阶矩本身就具有一些性质：

• 零阶矩
  $$\underset{v}{\iiint }{d}^{0}dm=\underset{v}{\iiint }dm=M$$
• 到坐标平面的一阶矩
  $$\begin{array}{c}\underset{v}{\iiint }{x}^{1}dm=\underset{v}{\iiint }x\cdot dm=x_{c}M\\ \\ \text{其中}{x}_c表示质体质心的横坐标，对于y,z同理\end{array}$$
  二阶矩则和转动惯量等有关，这里不再赘述。重点在于引力位低阶项的物理含义。

5.2.3 引力位低阶项物理含义

• 零阶项：取 $n=0$，得
  $$V_0=\frac{A_0}{r}=\frac{1}{r}G{\iiint }_{earth}{r}^{\prime 0}{P}_0\left(\cos {\theta }^{\prime }\right)dM=\frac{GM}{r}$$
  这表示了将地球质量全部集中在地球质心上所产生的引力位，因此零阶项与地球质量有关

• 一阶项：取 $n=1$，得
  $$V_1=\frac{1}{r^2}\left[A_1{P}_1\left(\cos \theta \right)+A_{11}{R}_{11}+B_{11}{S}_{11}\right]=\frac{1}{r^2}\left[A_1\cos \theta +A_{11}\sin \theta \cos \lambda +B_{11}\sin \theta \sin \lambda \right]$$
  其中
  $$\begin{array}{c}{A}_1=G{\iiint }_{earth}{r}^{\prime 1}{P}_1\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{z}^{\prime }dM\\ \\ A_{11}=G{\iiint }_{earth}{r}^{\prime 1}{R}_{11}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{x}^{\prime }dM\\ \\ B_{11}=G{\iiint }_{earth}{r}^{\prime 1}{S}_{11}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{y}^{\prime }dM\end{array}$$
  容易看到，一阶项和地球质心坐标有关。因此，若将坐标原点选在重心，则三个系数全为零，即一阶项为0

• 二阶项：取 $n=2$，得
  $$V_2=\frac{1}{r^3}\left[A_2{P}_2\left(\cos \theta \right)+A_{21}{R}_{21}\left(\cos \theta \right)+B_{21}{S}_{21}\left(\cos \theta \right)+A_{22}{R}_{22}\left(\cos \theta \right)+B_{22}{S}_{22}\left(\cos \theta \right)\right]$$
  二阶项系数为：
  $$\begin{array}{c}{A}_2=G{\iiint }_{earth}{r}^{\prime 2}{P}_1\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}(-x^{\prime 2}-y^{\prime 2}+2z^{\prime 2})dM\\ \\ A_{21}=\frac{1}{3}G{\iiint }_{earth}{r}^{\prime 2}{R}_{21}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{x}^{\prime }{z}^{\prime }dM\\ \\ B_{21}=\frac{1}{3}G{\iiint }_{earth}{r}^{\prime 2}{S}_{21}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{y}^{\prime }{z}^{\prime }dM\\ \\ A_{22}=\frac{1}{12}G{\iiint }_{earth}{r}^{\prime 2}{R}_{22}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}(x^{\prime 2}-y^{\prime 2})dM\\ \\ B_{22}=\frac{1}{12}G{\iiint }_{earth}{r}^{\prime 2}{S}_{21}\left(\cos {\theta }^{\prime }\right)dM=G{\iiint }_{earth}{x}^{\prime }{y}^{\prime }dM\end{array}$$
  其中xy，xz，yz形式的积分，是一种二阶矩惯性积。当坐标轴和主惯性轴重合，它们就为零。

但对于地球而言，需要比较多的条件才能使得它们全为零。一般来说，在地球坐标系下，所有的一阶谐函数和二阶一次谐函数在展开式中都默认为零，也被称为禁止的或不容许的谐函数。

5.3 正常重力

这一块内容相对轻松，主要以概念为主。

5.3.1 引入正常重力场的意义

引⼊⼀个近似的地球重⼒位，它函数关系简单，⾮常接近真实的地球重⼒位，称为正常重⼒位，记为 $U$。这样就将地球重力场的求解归结为扰动场或异常重力场（微小量）的求解，保证了其解的存在性，并方便求解。

5.3.2 确定正常重力位的方法

• Laplace方法
  • 将地球重力位 $W$ 展开为球谐级数
  • 保留前面最大的几项作为正常重力
  • 选取一个正常重力位水准面，假设它是产生正常重力位的质体的表面，则正常重力场就理解为该质体产生的重力场
• Stokes方法
  • 选择一个形状大小已知的质体（一般为旋转椭球体）
  • 根据质体的总质量和旋转角速度，解其形成的外部重力位
  • 以此重力位为正常重力位

5.3.3 正常椭球

• 基本要求
  • 质量和角速度等于地球总质量和旋转角速度，旋转轴和地球自转轴重合
  • 椭球表面为水准面 $W=W_0$，且外部没有物质存在
  • 椭球中心与地球质心重合
• 基本参数（任选四个可以确定正常重力场）
  • $U_0$：椭球表面上的正常重力位
  • $A_0=GM$
  • $A_2=G(A-C)$
  • $\omega$：自转角速度
  • $f=\frac{a-b}{a}$：椭球扁率
  • $f^{\ast }=\frac{{\gamma }_b-{\gamma }_a}{{\gamma }_a}$：正常重力扁率
  • ${\gamma }_a$：赤道上正常重力

5.4 扰动重力

5.4.1 定义

真实重力位和正常重力位之差，表示为 $T$
$$W(x,y,z)=U(x,y,z)+T(x,y,z)$$

5.4.2 大地水准面高

将大地水准面上一点 $P$ 沿参考椭球面的法线 $n^{\prime }$ 投影到椭球面（正常重力）上的 $Q$ 点，距离 $PQ$ 就称为大地水准面高，用 $N$ 表示，如图所示：

5.4.3 重力异常

$P$ 点的重力矢量 ${\mathbf{g}}_P$和 $Q$ 点的正常重力矢量 ${\gamma }_Q$ 之差称为重力异常矢量

两者的大小（模）的差称为重力异常

5.4.4 垂线偏差

观测点处垂线与椭球面法线的夹角称为垂线偏差

计算公式：
$$\{\begin{array}{c}\xi =\Phi -\phi \\ \eta =\left(\Lambda -\lambda \right)\cos \phi \end{array}$$

5.4.5 重力扰动

重力扰动和重力异常不同之处在于：重力扰动在同一点上比较。定义为：

在大地水准面上同一点 $P$ 比较重力矢量 ${\mathbf{g}}_P$ 和正常重力矢量 ${\gamma }_P$，得到重力扰动矢量

两者大小之差为重力扰动

重力扰动矢量的方向和垂线偏差一样

5.5 Bruns公式

Bruns公式表示为：
$$N=\frac{T}{\gamma }$$

5.5.1 推导

还是这张图。由正常重力的定义，对于 $P$ 点有：
$$W_P=U_P+T_P$$
而 $PQ$ 的正常重力位之差可以表示为正常重力在 $N$ 上做的功。因此有关系式：
$$U_P=U_Q+{\left(\frac{\partial U}{\partial n}\right)}_{Q}N=U_Q-\gamma N$$
又根据水准面定义：
$$W_P=U_Q=W_0$$
将三式联立得到Bruns公式：
$$N=\frac{U_Q-U_P}{\gamma }=\frac{W_P-U_P}{\gamma }=\frac{T}{\gamma }$$

5.6 物理大地测量学基本方程式

对于 $P$ 点，其扰动重力定义为：
$$\delta g=g_P-{\gamma }_P\text{\hspace{1em}(1)}$$
和Bruns公式推导一样，其正常重力可以表示为：
$${\gamma }_P={\gamma }_Q+\frac{\partial \gamma }{\partial h}N\text{\hspace{1em}(2)}$$
联立(1)(2)得：
$$\delta g=g_P-{\gamma }_Q-\frac{\partial \gamma }{\partial h}N\text{\hspace{1em}(3)}$$
又有重力异常的定义为：
$$\Delta g=g_P-{\gamma }_Q\text{\hspace{1em}(4)}$$
代入(3)式得：
$$\delta g=\Delta g-\frac{\partial \gamma }{\partial h}N\text{\hspace{1em}(5)}$$
根据Bruns公式 $N=\frac{T}{\gamma }$，将式子中的 $N$ 替换，得到
$$\delta g=\Delta g-\frac{\partial \gamma }{\partial h}\frac{T}{\gamma }\text{\hspace{1em}(6)}$$
将重力扰动表示为扰动位的梯度，并近似处理：
$$\delta g=-\frac{\partial T}{\partial n}\approx -\frac{\partial T}{\partial h}\text{\hspace{1em}(7)}$$
代入(6)式得到物理大地测量学基本方程式：
$$-\frac{\partial T}{\partial h}=\Delta g-\frac{\partial \gamma }{\partial h}\frac{T}{\gamma }\text{\hspace{1em}(6)}$$
该式建立了扰动位、正常重力和重力异常之间的关联。

5.7 Possion积分（内涵）

Possion积分公式说明了空间点上谐函数的值可以由边界面上的值来确定，因此是一个将谐函数由边界向上延拓的方法

5.8 改进的Poisson积分（内涵）

改进的Poisson积分公式去掉了零阶项和一阶项，使得计算更加简便。

对于重力异常来说，重力异常并不是一个谐函数，需要乘以r构造成谐函数 $r\Delta g$
将此谐函数展开后，应用改进的Poisson积分，即可由地面重力异常求出地球外部重力异常

5.9 Stokes公式

$$N=\frac{R}{4\pi {\gamma }_0}\underset{\sigma }{\iint }\Delta gS(\psi )d\sigma$$
Stokes公式是目前为止物理大地测量最重要的公式，它可以根据重力异常来确定大地水准面。

6 重力归算

尽管Stokes公式似乎把问题变得很简单，但想要确定重力异常，需要的是大地水准面上的重力值。现实情况下，我们只能得到地表的重力值。另外，Stokes公式中大地水准面外必须无质量，而事实上并非如此。

重力归算的任务就是将大地水准面外部的质量去掉（或移到大地水准面以内），再将重力点从地面归算到大地水准面。

6.1 空间改正

6.1.1 定义

空间改正，就是将海拔高程为 $H$ 的重力点的重力观测值 $g$，归算到大地水准面的重力观测值 $g_0$，归算时不考虑地面和大地水准面间的质量影响，只考虑高度对重力的改正

6.1.2 推导

空间改正中的不考虑地面和大地水准面间的质量影响，指的是将所有质量压入到大地水准面内。

对于地面上 $P$ 点的引力：
$$g=G\frac{M}{(R+H{)}^2}\text{\hspace{1em}(1)}$$
对于大地水准面上对应的 $P_0$点：
$$g_0=G\frac{M}{R^2}\text{\hspace{1em}(2)}$$
注意一下符号问题：地面上观测值＋改正值=大地水准面上值。因此改正值为：
$${\Delta }_{1}g=g_0-g=G\frac{M}{R^2}-G\frac{M}{(R+H{)}^2}=GM\left[\frac{1}{R^2}-\frac{1}{(R+H{)}^2}\right]=\frac{GM}{R^2}\left[1-\frac{1}{(1+\frac{H}{R}{)}^2}\right]\text{\hspace{1em}(3)}$$
可以看到式子中能够提出$\frac{GM}{R^2}$，记作 $\gamma$，并将 $\frac{1}{(1+\frac{H}{R}{)}^2}$ 级数展开得：
$${\Delta }_{1}g=\gamma \left[1-\left(1-\frac{2H}{R}+\frac{3H^2}{R^2}\right)\right]=2\frac{\gamma }{R}H-3\frac{\gamma }{R^2}{H}^2=0.3086H-0.72\times 10^{-7}{H}^2(mGal)\text{\hspace{1em}(4)}$$
另一种方法推导更为简单：
根据泰勒级数展开，并保留线性部分得：
$$g_0=g-\frac{\partial g}{\partial H}H$$
根据近似关系：
$$-\frac{\partial g}{\partial H}\approx -\frac{\partial \gamma }{\partial h}=+0.3086H(mGal)$$
即完成了推导。

6.1.3 空间重力异常

将空间重力异常记作：
$$\Delta g_{\text{空}}=g+{\Delta }_{1}g-\gamma$$

6.2 层间改正

6.2.1 定义

层间改正又称中间层改正。设地球表面和大地水准面均为平面。那么两平面之间的质量对两个面上点的重力存在着影响、去掉这部分质量引起的重力改正称为层间改正

6.2.2 推导

这里的过程不过多展开，实际上就是求圆柱体对轴线上点的引力。需要注意的是，圆柱体的半径a在这里趋近于无穷大，从而接近地球表面。

布格片引力表示为：
$$A_B=2\pi G\rho H$$
代入标准密度 $\rho =2.67g/c{m}^3$ 时，得：
$$A_B=0.1119H(mGal)$$

6.2.3 符号

在层间改正中，去掉了布格片质量影响，使得重力值减小，因此改正值为负。得到层间改正表达式为：
$${\Delta }_{2}g=-0.1119H(mGal)$$

6.3 地形改正

6.3.1 定义

在层间改正的基础上，考虑地形不是平面的情况，需要将 $I$ 区的质量补充进 $II$ 区，以补偿被多余去掉的质量。

将局部地形改正定义为：计算点周围地形起伏的质量对计算点重力值的影响，记作 ${\Delta }_{3}g$

6.3.2 符号

在局部地形改正中，仅考虑局部地形的影响。在挖空的 $II$ 区补进质量，会使得重力值变大；在多余的 $I$ 区挖去质量，同样会使得重力值变大。因此符号为正。

6.4 各种改正的定义

• 简单布格改正：空间改正+层间改正
• 布格改正：空间改正+层间改正+局部地形改正
• 法耶改正：空间改正+局部地形改正

对应的有布格重力异常、简单布格重力异常、法耶重力异常，不再赘述，通用的公式是重力+对应改正-正常重力。

6.5 地形均衡理论

在前面的改正中，都假设了地球的密度均匀，导致布格异常在地形起伏的地方存在系统的误差。因此针对地形对密度的影响，提出了普拉特均衡模型和艾里均衡模型。

6.5.1 普拉特均衡模型

6.5.1.1 主要思想

认为在地下某深度有一个补偿面（海水面和补偿面距离几乎处处相等），在补偿面之上，将地壳分割成截面相同的柱体，超出海平面的山区柱体密度小，低于海平面的海底柱体密度大，但柱体质量相等。

6.5.1.2 密度计算

设海平面到补偿面距离为 $D$，山体柱体高出海平面的高度记为 $H$，海底柱体低于海平面的高度记为 $H^{\prime }$。当 $H=0$时柱体密度记为 ${\rho }_0=2.67g/c{m}^3$（也就是标准密度的值）

对于陆地：
$$(D+H)\rho =D{\rho }_0$$
对于海洋（额外考虑海水密度 ${\rho }_w$）：
$$(D-H^{\prime })\rho +H^{\prime }{\rho }_w=D{\rho }_0$$

PPT上还介绍了一种将海平面以上和以下的密度区分开的模型（如下图），但书上没有，这里就不写了。

6.5.2 艾里均衡模型

6.5.2.1 主要思想

艾里均衡模型认为地壳由厚度不同的轻的岩石组成，各个柱体漂浮在密度较大的岩浆之上。
每个柱体的密度相同，并且露出岩浆和陷入岩浆的部分是相对应的。凸起部分越高，陷入部分也越深

6.5.2.2 密度计算

假定山体的密度为 ${\rho }_0=2.67g/c{m}^3$，下层岩浆密度为 ${\rho }_1=3.27g/c{m}^3$

对于山体，其陷入岩浆部分长度为 $t$，则根据浮力平衡原理：
$$t{\rho }_1=H{\rho }_0+t{\rho }_0\Rightarrow t\Delta \rho =H{\rho }_0$$
因此得到 $t$ 的值为：
$$t=\frac{H{\rho }_0}{\Delta \rho }=4.45H$$
对于海洋，则有
$$t^{\prime }\Delta \rho =H^{\prime }({\rho }_0-{\rho }_w)\Rightarrow t^{\prime }=\frac{{\rho }_0-{\rho }_w}{{\rho }_1-{\rho }_0}{H}^{\prime }=2.73H^{\prime }$$

6.5.3 均衡改正

6.5.3.1 定义

依据某种均衡模型（可以是普拉特也可以是艾里-海斯卡涅）调整地壳，最后结果是一个想象的密度为 ${\rho }_0$ 的均匀地壳，它并不像布格改正一样完全去掉地形质量，而是将这些质量移到大地水准面内部，从而弥补山下的质量亏损，使得密度均衡。

均衡改正包括三个步骤：

1. 移去地形部分
2. 移去补偿部分
3. 加上空间改正，归化到大地水准面上

6.5.3.2 均衡重力异常

均衡重力异常=重力值+空间改正+层间改正+局部地形改正+均衡改正-正常重力

6.6 重力归算总结

6.6.1 重力归算的目的

• 求定大地水准面
• 内插和外推重力值
• 研究地壳

6.6.2 重力归算的基本要求

• 大地水准面的外部没有质量
• 不改变地球质心的位置，即满足椭球体和大地水准面质心重合的条件
• 地球的总质量不变
• 不改变大地水准面的形状

6.6.3 重力归算的步骤

• 移去地形质量部分（层间改正、地形改正）
• 移去补偿质量部分（均衡改正）
• 加空间改正归化到大地水准面上（空间改正）