天不生牛顿，万古如长夜。

牛顿在研究万有引力的时候，很是为一件事苦恼：均匀球体的引力是否精确地等于所有质量集中于球心时产生的引力？当他解决了这个问题以后，才确认自己走在正确的道路上。

这个问题可以分为三步。

首先，均匀球壳对内部质点的引力为零。几何方法最容易了。对于球壳内的任意一点P，过该点的任意平面与球壳相交为圆。过P点的两条直线与圆相交于A、B、C、D四点。三角形PAB与PCD相似（对顶角相等，同一段圆弧的张角相等），所以，AB段的引力与CD段大小相等、方向相反，二者的合力为零。因为P点、平面和两条直线都是任意选择的，所以，均匀球壳对内部质点的引力为零。注意，这个结果只对完美的球壳成立，对其他形状不成立，特别的，对于圆圈也是不成立的，因为引力与距离成平方反比关系，而球壳上AB段对应的质量是其平方值，而圆圈上AB段正比于长度。

其次，均匀球体对外部质点的引力精确地等于所有质量集中于球心时产生的引力。这个用积分算很容易。$OP=R > OA=OB=r$，根据对称性可知，引力方向位于OP，张角$\theta$的圆环对P的引力正比于

$\int ^{\pi}_0 2\pi r\sin \theta r d \theta r \frac{R-r\cos \theta}{(R^2+r^2-2Rr\cos \theta)^{3/2}}$

$=s^2\int ^1_{-1} dx\frac{1-sx}{(1+s^2-2sx)^{3/2}}$

其中，$s=r/R$，$x=\cos \theta$。把积分项的分子$1-sx$拆解为$1+s^2-2sx$和$1-s^2$之差，二者的原函数就是$\sqrt{1+s^2-2sx}/s$和$(1-s^2)\sqrt{1+s^2-2sx}/s$。定积分的结果就正比于$s^2$，由此得证。注意，这个结论仍然只对完美的球壳成立，对于其他形状乃至圆圈都不成立，因为积分里包含了$\sin \theta$。

第三步就简单了。均匀球体对球外质点的引力等于全部质量集中于球心时的引力。实际上，如果密度只依赖于半径，这个结论就仍然成立。

现在考虑岁差。岁差的计算是牛顿《自然哲学之数学原理》中最重要的结果之一，甚至可以把“之一”去掉。

首先，什么是岁差？简单地说，是这样的（力学教学笔记之曰若稽古）：

地球绕着太阳公转，地球轨道所在的平面就是黄道面（古人以为太阳绕着地球转，黄道面是太阳在天空中的运行面）。地球的自转轴与黄道面有个夹角（大约66.5度）。地球的公转轨道决定了一个平面，太阳和地球自转轴决定了另一个平面，随着地球的公转，这两个平面的夹角也在变化。当这两个平面垂直的时候，就是冬至和夏至（这两个位置是对径点）；当这两个面的夹角最小的时候，就是春分和秋分（这两个位置也叫作两分点，它们是对径点）。

因为地球不是完美的球形，而是扁率0.3%的旋转椭球，太阳和月亮对地球的引力会改变地球自转轴的方向，虽然自转轴与黄道面夹角不变。换句话说，地球自转轴在锥角为23.5度的圆锥面上进动，每年变化50”（这就是岁差，古希腊的喜帕恰斯和我国晋朝的虞喜都注意到了这一点），进动周期大约是26000年。

我们来计算这个岁差。

假设地球的密度是均匀的，其扁率是0.3%。日地连线与地球自转轴的夹角$\theta _0$为66.5度。简单地说，地球就是个陀螺，太阳对超出球体部分的引力不在球心上，这样就产生了一个力矩，使得地球进动，就像倾斜的陀螺在重力作用下的进动一样。

地球角动量的变化率是

$\frac{2}{5}{m}{r^2}\cdot \omega \cdot  \cos \theta \cdot \frac{2\pi}{T}$

第一项是地球的转动惯量，第二项是地球的自转角速度（每天转动一周），第三项表示到黄道面的投影，第四项表示变化的快慢，$T$就是进动周期。

（计算力矩的示意图：完全不合比例。蓝色为地球偏离球体的部分；红色为太阳。）

接下来估计力矩的大小。作用在非球形部分的力类似于潮汐力，正比于$\frac{GM \delta m}{R^2}\frac{r}{R}$，其中，$\delta $是地球的扁率，$r$是地球半径，$R$是日地距离，$M$和$m$分别是太阳和地球的质量。利用$GMm/R^3=m\Omega ^2$，其中，$\Omega $是地球公转角速度，就可以得到力矩

$ \alpha \cdot \delta \cdot  mr^2 \Omega ^2$

其中，$\alpha$是个待定系数。让力矩等于角动量的变化率，就可以得到进动的周期

$T=\frac{2/5}{\alpha \delta} \cdot \frac{2\pi}{\Omega} \cdot \frac{\omega}{\Omega}$

其中，$\frac{2\pi}{\Omega}$就是地球的公转周期（1年），$\frac{\omega}{\Omega}$是地球自转速度与公转速度的比值（365），$2/5$来自于转动惯量的系数，$\delta =0.003$，把这些数值带入以后，可以得到

$T=\frac{4.8}{\delta}\times 10^4 /$（年）

正好是几万年的量级。$\delta$的数值没办法估计了，我只能硬算了。对于太阳，可以得到

$\delta \approx \frac{9\pi}{32}\sin(2\times 23.5) \approx 0.63$

由此可知

$T=7.7\times 10^4$（年）

这个数值与地球进动的周期26000年相差甚远，因为我们还没有考虑月亮的影响。月亮的潮汐力是太阳的2.2倍，所以上面这个数值应该再除以3.2，也就是说

$T=2.4\times 10^4$（年）

很好的估计，误差不到10%。

考虑到地球的密度并不是均匀的，而且地球的扁率也没有采用非常准确的数值，这个估计值应该算是很好的了。但是还比不上牛顿的估计，他在《原理》中得到的岁差是每年50角秒，周期非常接近于26000年。牛顿采用的月亮与太阳的潮汐力比值是4.4815:1，所以，我觉得他可能也是硬往上凑的——当然，更有可能是他用岁差数据（以及潮汐数据）来反推月球质量等参数。

因为我看不懂《原理》，所以并不知道牛顿到底怎么做的。但是无论如何，《原理》开创了科学的新时代，其中的具体结果可能会随着科学的进步而显得过时，但是牛顿的研究方法和独创精神必然是“历千万祀，与天壤而同久，共三光而永光。”