十二、中国农历—定朔望（下）

福至运成

月球绕地球自西向东在星空中移动，与地球的自转方向相同，其视轨道称为月球的“白道”。白道面与黄道面之间平均有大约5º8'42.7”的夹角，朔望月的平均周期是29天12小时44分2.8016秒。月球受其它天体的影响，白道面与黄道面之间夹角在4º57'到5º19'之间变化，朔望月周期在29.272天至29.833天之间，这些变化是不规则的。要准确预测“朔”的时间几无可能，但是根据现代科学模型，可以进行非常精确地计算和预测。

• 黄道：以地球为基准，太阳在天空中穿行的视路径大圆。

• 白道：以地球为基准，月球在天空中穿行的视路径大圆。

• 恒星月（sidereal ）：以恒星为参照，月球在星空视运行一圈的周期，平均27.321 661天。

• 近点月（anomalistic ）：月球连续两次经过近地点的时间间隔，平均27.554 551天。

• 交点月（draconic ）：月球连续两次经过升交点的时间间隔，平均27.212 220天。每个恒星月，白道面和黄道面有两个交点。月球从南向北穿过黄道，叫“升交点”；月球从北向南穿过黄道，叫“降交点”。

• 朔望月（synodic ）：月相变化周期，连续两次朔的时间间隔，平均29.530 587 981天。

经验公式

月球绕地球公转是一个椭圆轨道。受角动量守恒的约束，月球的运行速度不是一个常数，是变化的。可以用经验公式进行近似计算：

M = 1.6 + 29.530588 m + 0.4 sin (1 - 0.45059 m)

M是从1900年开始的第m个朔日的天数，sin是正弦函数，m周期系数0.45059 = 2π(29.530588-27.554551)/27.554551，0.4为正弦函数的振幅；1.6和1分别是以1900年为起点的常数，分别代表m = 0时的“积日”和月球相对近地点的轨道位置。经对1900年到2200年期间的3723朔日验算发现，若对其常数稍加修改则与实际朔日误差会更小。如下：

M = 1.5925 + 29.530588 m + 0.41 sin (0.94 - 0.45059 m)

用修改后经验公式计算结果与实际朔日时间比较，最大偏差在（-0.2073, 0.2028）天之间，平均偏差在（-0.1097, 0.1103）天范围。在1900年到2200年期间的3723朔日中，有412个朔日的日期需要修正。也就是说，若用修改后经验公式计算，大约每9个朔日则会遇到需要修改朔日日期的情况。

显然，需要更加精确的方法，参见下面的半解析算法。

半解析算法

这里引用的是Jean Meeus编著《Astronomical Algorithms》一书（Willman-Bell出版社，1999，第二版）第49章第349页的算法。这个算法的方便之处是，只使用一个参数K（朔望月周期数），就可以达到分钟之内的精度。K是整数时，月球黄经在朔（0°新月），半整数在望（平满月），1/4上弦（平半满上弦月），3/4下弦（平半满下弦月），其它值没有意义。2000 年第一个朔日（ 1 月 7 日）相应的K值是0，计算 2000 年以前的月相，K取负数。

一儒略年（365.25天）中，平均有12.368531个朔望月。以百儒略年为单位的时间T = K/1236.8531，T2是其平方，T3表示三次方，T4为四次方。平月相的时间TD等于:

5.59766+29.530588853*K + 0.0001337*T2 -0.000000150*T3 + 0.00000000073*T4

结果已修正太阳光行时差及月球光行时差。TD是以日为单位的力学时间，是从2000年1月1日0时0分0秒起算的实数。

周期修正项

真正的月相时间需要对平月相时间TD加一些修正项，就是对一系列的正弦表达式A sin(θ)求和。A是振幅，θ是辐角（°）。θ由以下4个辐角参数求得：

太阳平近点角：
M = 2.5534 + 29.10535670k - 0.0000014T2 - 0.00000011T3

月亮的平近点角：
M′ = 201.5643 + 385.81693528K + 0.0107582T2 + 0.00001238T3 - 0.000000058T4

月亮的纬度参数：
F = 160.7108 + 390.67050284K - 0.0016118T2 - 0.00000227T3 + 0.000000011T4

月亮轨道升交点经度：
Ω = 124.7746 - 1.56375588K + 0.0020672T2 + 0.00000215T3

因为M的值与地球公转轨道的离心率有关，导致振幅A是个变量，需要根据时间进行修正。当辐角中包含了M或-M时，需要乘以系数E，若包含了2M或-2M，则需要乘以系数E的平方E2。E = 1 - 0.002516 T - 0.0000074 T2。这里仅介绍朔和望的计算参数, 如下表：

朔振幅A 望振幅A θ（°）
-0.40720 -0.40614 M′
0.17241E 0.17302E M
0.01608 0.01614 2M′
0.01039 0.01043 2F
0.00739E 0.00734E M′-M
-0.00514E -0.00514E M′+M
0.00208E2 0.00209E2 2M
-0.00111 -0.00111 M′-2F
-0.00057 -0.00057 M′+2F
0.00056E 0.00056E 2M′+M
-0.00042 -0.00042 3M′
0.00042E 0.00042E M+2F
0.00038E 0.00038E M-2F
-0.00024E -0.00024E 2M′-M
-0.00017 -0.00017 Ω
-0.00007 -0.00007 M′+2M
0.00004 0.00004 2M′-2F
0.00004 0.00004 3M
0.00003 0.00003 M′+M-2F
0.00003 0.00003 2M′+2F
-0.00003 -0.00003 M′+M+2F
0.00003 0.00003 M′-M+2F
-0.00002 -0.00002 M′-M-2F
-0.00002 -0.00002 3M′+M
0.00002 0.00002 4M′
行星影响

一天有86400秒，1分钟约等于0.0007天。如果仅希望达到分钟量级的精度，可以忽略一些影响较小的项，比如上表中Ω以下的项。如果希望达到比较好（分钟以内）的精度，则需要再对一些影响较大的行星影响做附加矫正，共有14个。也是对其正弦表达式A sin(θ)求和：

振幅A θ（°）
0.000325 299.77 + 0.107408K - 0.009173T2
0.000165 251.88 + 0.016321K
0.000164 251.83 + 26.651886K
0.000126 349.42 + 36.412478K
0.000110 84.66 + 18.206239K
0.000062 141.74 + 53.303771K
0.000060 207.14 + 2.453732K
0.000056 154.84 + 7.306860K
0.000047 34.52 + 27.261239K
0.000042 207.19 + 0.121824K
0.000040 291.34 + 1.844379K
0.000037 161.72 + 24.198154K
0.000035 239.56 + 25.513099K
0.000023 331.55 + 3.592518K
协调世界时转换

协调世界时（Coordinated Universal Time，UTC），又称世我们日常生活中，通常是以基于地球自转的世界时（Universal Time，UT或UT1）表示时间的。北京标准时间处于世界时区的东八区（UT+8），比世界时早8小时。世界时，目前使用的是协调世界时（Coordinated Universal Time，UTC），与格林威治时间（Greenwich Mean Time， GMT，UT0）一样，是以本初子午线（地球上的零度经线）的平子夜起算的平太阳时。

受地球自转速度长期变慢、季节性变化，以及不规则因素的影响，世界时不是一种均匀的时间系统。在进行天文计算时，使用的是地球时间（Terrestrial Time，TT，1991。之前使用的是力学时（Terrestrial Dynamical Time，TDT或TD，1984）。力学时间取代了再之前非相对论的星历（即历书）时（Ephemeris Time，ET，1952）。由于历书时和力学时的测定精度较低，从1967年起使用的是原子时（International Atomic Time，TAI）。

原子时钟以铯-133原子基态能级精细结构之间的光子能量吸收频率为标准，即1秒钟时间等于9 192 631 770周期。原子时是一个国际单位制（International System of Units， SI），它是一种均匀的时间计量系统，其误差只有百万亿分之二秒，即每一百四十万年才有一秒的误差。

因此，需要对基于力学时间的计算结果TD进行时间修正，变换到我们日常使用的世界时（UT）：ΔT=TD-UT。然而，由于世界时数值是不均匀、不规则的，导致ΔT的数值无法准确预测。只能根据天文观测，经过数据处理后，再给出具体数值。通常情况下，可以通过基于ΔT真实数据拟合的经验公式计算获得。

以下就是从美国国家航空航天局（NASA）网站获得的部分经验公式。这些公式中的Y是年份，是一个实数，月份M可以转换成年份的小数部分：(M - 0.5)/12。T是以相关年份为参照的时间，T2是其平方，T3是三次方，T4是四次方，T5是五次方。

1900-1920: T=Y-1900

ΔT(秒)=-2.79+1.494119T–0.0598939T2+0.0061966T3-0.000197T4

ΔT(秒)举例：

1900：-3

1910：11
1919： 21
1920-1941: T=Y-1920

ΔT(秒)=21.20+0.84493T-0.076100T2+0.0020936T3
ΔT(秒)举例：

1921：22
1931：24

1941-1961: T=Y-1950
ΔT(秒)=29.07+0.407T-T2/233+T3/2547
ΔT(秒)举例：

1941：25
1950：29

1961-1986: T=Y-1975
ΔT(秒)=45.45+1.067T-T2/260–T3/718
ΔT(秒)举例：

1961：34
1980：51

1986-2005: T=Y-2000
ΔT(秒)=63.86+0.3345T-0.060374T2+0.0017275T3+0.000651814T4+0.00002373599T5
ΔT(秒)举例：

1985：54
1990：57
1995：61
2000：64
2005：65

2005-2050: T=Y-2000
ΔT(秒)=62.92+0.32217T+0.005589T2
ΔT(秒)举例：

2010：67
2025：75

2050-2150: T=(Y-1820) /100
ΔT(秒)=-20+32T2-0.5628(2150-y)
ΔT(秒)举例：

2050：93
2100：203

2150以后: T=(Y-1820) /100
ΔT(秒)=-20+32T2
ΔT(秒)举例：

2150：330
2200：442

朔计算举例

这里比较使用的朔精确时间取之于《寿星天文历（V5.05Plus）》，有多处网站可以获得，比如：http://www.nongli.net/sxwnl/。

根据年份Y估算所需要计算的朔日K，可以通过如下公式进行近似计算：

K = 12.3685（Y-2000）- 0.2

Y是一个实数，可从月份M获得其小数部分：(M - 0.5)/12。0.2是因为2000年的第一个朔日在1月7日2:13:38。

例1，2022年第一个朔日，月份为1，Y小数部分为0.02。相应的朔日K约等于272.16 = (12.3685x22.02)–0.2，所以朔时刻K应为整数272，相应的T = 0.219912939。其平朔的时间TD等于= 8037.917838。进一步得到：

E = 0.999446341
M = 7919.2104 => 359.2104
M′ = 105143.7712 => 23.7712
F = 106423. 0875 => 223. 0875
Ω = -300.5669 => 59.4331

第一组周期参数累加结果是-0.143142，14个额外修正值是-0.000614。因此朔真实时间是：

TD = 8037.917838 - 0.143142 - 0.000614 = 8037.774081

这相当于2022年1月2日18:34:40。加8小时获得北京时间，2022年1月3日2:34:40。这时的ΔT=TD-UT是72.7秒，因此计算所得的最终结果时间是2022年1月3日2:33:27。比2022年的第一个朔日实际时间2:33:20偏大7秒钟。

例2，2022年最后一个朔日，月份为12，Y小数部分为0.96。相应的朔日K约等于283.78 = (12.3685x22.96)–0.2，所以朔时刻K应为整数284，相应的T = 0.229614980。其平朔的时间TD等于= 8392.284905。进一步得到：

E = 0.999421899
M = 8268.4747 => 348. 4747
M′ = 109773.5745 => 333. 5745
F = 111111.1335 => 231. 1335
Ω = -319.3320 => 40.6680

第一组周期参数累加结果是0.144953，14个额外修正值是-0.000687。因此朔真实时间是：

TD = 8392.284905 + 0. 144953 - 0. 000687= 8392.429170

这相当于2022年12月23日10:18:00。加8小时获得北京时间，2022年12月23日18:18:00。这时的ΔT=TD-UT是73.3秒，因此计算所得的最终结果时间是2022年12月23日18:16:47。比2022年的最后一个朔日实际时间18:16:41偏大6秒钟。

下面给出更多几个实例的简单计算结果比较：

1900年的第一个朔日在1月01日21:51:58，计算结果为21:51:55偏小2.8秒；
1950年的第一个朔日在1月18日15:59:30，计算结果为15:59:32偏大2.8秒；
1975年的第一个朔日在1月12日18:19:36，计算结果为18:19:43偏大7.5秒；
2000年的第一个朔日在1月07日02:13:38，计算结果为02:13:42 偏大4.4秒；
2020年的第一个朔日在1月25日05:41:52，计算结果为05:42:00偏大8.4秒；
2040年的第一个朔日在1月14日11:24:40，计算结果为11:25:00偏大20.8秒；
2060年的第一个朔日在1月04日00:39:28，计算结果为00:39:59偏大31.7秒；
2080年的第一个朔日在1月22日09:54:59，计算结果为09:55:17偏大18.5秒；
2100年的第一个朔日在1月10日20:55:17，计算结果为20:55:30偏大13.8秒；
2125年的第一个朔日在1月05日00:56:58，计算结果为00:57:06偏大8.4秒；
2150年的第一个朔日在1月29日03:15:42，计算结果为03:15:31偏小10.5秒；
2200年的第一个朔日在1月16日16:20:03，计算结果为16:19:49偏小13.8秒；
2200年最后一个朔日在12月7日18:13:01，计算结果为18:12:46偏小14.8秒。

1900年到2200年之间的所有朔日计算，没有发现朔日日期需要修正的例外情况。这些计算表明，20世纪的最大误差是14.5秒，平均误差是3.5秒；21世纪的最大误差44.9秒，平均19.6秒；22世纪的最大误差32.3秒，平均10.7秒。21世纪和22世纪误差偏大的原因，很可能是因为力学时转换为世界时引入了较大误差。原因是目前还没有实际的实测数据，仅仅利用经验公式外推，会带来较大的不确定性。 参阅 《二十四节气》
参阅 《定气法（上）》
参阅 《定气法（下）》
参阅 《定朔望（上）》
参阅 《历法编算》 返回 《中国农历简介》