---
title: "用已发表的能流表复算一条 EMRI 频率演化"
description: "运行一个 Python 脚本，把 Fujita–Tagoshi 能流表换成波周数和频率随时间的变化，再故意改三处，看结果怎么变。"
canonical_url: "https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/recompute-a-frequency-evolution/"
language: "zh-CN"
page_type: "story"
alternate_languages:
  en: "https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/recompute-a-frequency-evolution/"
last_modified: "2026-11-12"
---

[极端质量比旋近，从轨道读到探测预测](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/)

# 用已发表的能流表复算一条 EMRI 频率演化

Sandbox Physics · 2026-11-12

这份教程复算圆轨道 EMRI 的波频率怎样随时间升高。能流不用四极公式，而用一份已发表 Teukolsky 能流表的冻结副本。所需的东西都在一个文件夹里：输入、只用标准库的 Python 脚本、预期输出和校验和。运行不到一秒。

## 开始前需要什么

Python 3，以及本课程讲能量平衡的那一章。除标准库外不需要任何包。

## 要复算的结果

三条圆赤道旋近，$q=10^{-5}$，$M=10^6\,M_\odot$：零自转从 $r=10$ 到 ISCO 外侧，自转0.9从 $r=10$ 到 $r=6$，自转 $-0.9$ 从 $r=20$ 到 $r=10$。每一条，脚本都给出用表中能流和用四极能流的时间与波周数，以及九个点的频率轨迹。

| 情形 | 周数（表中能流） | 年数（表中能流） | 周数（四极能流） | 比值 |
| --- | --- | --- | --- | --- |
| $a=0$，$10\to6.003$ | $112\,981$ | $1.417$ | $112\,045$ | $0.9917$ |
| $a=0.9$，$10\to6$ | $260\,673$ | $3.157$ | $222\,719$ | $0.8544$ |
| $a=-0.9$，$20\to10$ | $971\,669$ | $31.47$ | $970\,777$ | $0.9991$ |

## 第一步：冻结的能流表

Fujita 和 Tagoshi 的表VII给出圆赤道轨道到 $\ell=6$ 为止的引力波光度，自转从 $-0.9$ 到 $0.9$ 共七种。输入文件收录 $r=6$ 到 $r=20$ 各行，是从 arXiv 第1版手工抄录的，表注原文一并保存。

表里没有可以直接引用的归一化说明，所以“按 $(\mu/M)^2$ 归一、到达无穷远的能流”这个解读是我们的推断。有两项检验支持它：$r=6$ 处逐个模式的表III加起来，与表VII的值相差不到 $5\times10^{-5}$；在 $r\ge50$ 的远处，表中数值与 Tagoshi–Sasaki 后牛顿无穷远能流相差不到 $2\times10^{-4}$。

对应阅读：Fujita & Tagoshi (2004) · Tagoshi & Sasaki (1994)

## 第二步：插值一个平缓的修正因子，而不是能流本身

整张表里能流跨了三个数量级，所以脚本插值的是表中能流与四极能流之比。这个比值始终在0.8到2之间，随半径平缓变化。脚本用最近的四个表列半径做三次插值，6到20以外的半径一律拒绝。

能量斜率 $E'(r)$ 由 Kerr 圆轨道精确能量做四阶有限差分得到。两个积分都在 $\ln r$ 上用辛普森公式，分2000段。

$$
t=\int_{\ln r_1}^{\ln r_0}\frac{r\,E'(r)}{\mathcal F(r)}\,d\ln r,\qquad N=\int_{\ln r_1}^{\ln r_0}\frac{\Omega_\phi(r)}{\pi}\,\frac{r\,E'(r)}{\mathcal F(r)}\,d\ln r,\qquad \mathcal F=\tfrac{32}{5}\,q\,\Omega_\phi^{10/3}\times\frac{\mathcal F_{\rm table}}{\mathcal F_{\rm quad}}
$$

在 $\ln r$ 上积分，能把取样点均匀地放在被积函数变化最快的地方。

## 第三步：运行并比较

下载文件夹，先用校验和确认文件完整，用输入文件和一个新的输出文件名运行脚本，再对照你的输出和预期输出。在生成预期文件的机器上（macOS，Python 3.9.6），重新运行逐字节一致。换一台机器，预计相对差在万亿分之一左右。

时间用 $GM_\odot/c^3=4.925490947\times10^{-6}\,\mathrm s$ 换成秒，用儒略年31,557,600秒换成年。轨迹里的波频率是 $f=2\Omega_\phi/(2\pi)$，单位赫兹，对应 $M=10^6\,M_\odot$。

## 第四步：故意改坏

输出里记了三处故意的改动，以及每处让零自转周数变了多少。修正因子改用线性插值：−0.54%。辛普森公式从2000段减到4段：+0.088%。把 Kerr 能量换成牛顿束缚能：+99%。

最后一个数最值得记住：靠近 ISCO 时，束缚能和能流一样要紧。要公平地检验一个能流模型，轨道本身必须是精确的。

## 结果说明不了的事

拒绝外推。自转0.9的 ISCO 在 $2.32$，远低于表的第一行，所以这一情形停在 $r=6$。输出里把 $r=4$ 的请求列为已拒绝。

缺少的多极。$r=6$ 处零自转能流里，$\ell=6$ 模式占0.29%，略去的 $\ell\ge7$ 模式那里大约占0.1%。这是顺着趋势的外推，不是计算。

缺少的物理。视界能流、偏心和倾斜轨道、共振，以及所有后绝热项都没有算进去。后绝热相位修正约一弧度，所以这些周数只在质量比的领头阶上是对的。

对应阅读：Fujita & Tagoshi (2004) · Barack & Pound (2019)

## 下载并运行

把六个文件放进同一个文件夹，运行下面两条命令。第二条会在预期输出旁边写出你自己的结果。

`shasum -a 256 -c SHA256SUMS python3 recompute_frequency_evolution.py inputs.json my-output.json`

- [recompute\_frequency\_evolution.py](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/recompute_frequency_evolution.py) · 脚本，只用 Python 标准库
- [inputs.json](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/inputs.json) · 情形、单位与冻结的能流表
- [expected-output.json](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/expected-output.json) · 预期输出，Python 3.9.6
- [README.txt](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/README.txt) · 单位、模型与失败情形
- [SHA256SUMS](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/SHA256SUMS) · 每个文件的 SHA-256 校验和
- [LICENSE.txt](https://sandboxphysics.com/stories/extreme-mass-ratio-inspirals/materials/frequency-evolution/LICENSE.txt) · 使用许可与署名

## 在实验室里试一试

### [与实验室的四极周数对照](https://sandboxphysics.com/zh/labs/emri-spacetime-mapping/)

1. 在实验室里设自转为零、偏心率为零、$p=10$。
2. 读出“剩余波周数”。实验室保留三位有效数字，与上表第一行四极列的112,045比较。

实验室从你设定的轨道大小开始计数，只用四极能流。

## 练习与答案

哪个情形里表中能流与四极能流差得最多？为什么是它？

自转0.9、从 $r=10$ 到 $r=6$，比值0.854。在这段范围里四极能流处处偏大15%到19%，误差不变号，也就没法抵消。零自转时误差在 $r=8$ 附近变号。

把输入里零自转的终止半径改成5.9再运行，会发生什么？

脚本报错停下，指出这个半径并拒绝外推，因为 $r=5.9$ 低于表的第一行。即使表更大，它也落在零自转 ISCO $r=6$ 以内，那里没有稳定的圆轨道。

## 一手来源，以及每篇读了哪些部分

1. [Barack & Pound (2019), Self-force and radiation reaction in general relativity, Rep. Prog. Phys., doi:10.1088/1361-6633/aae552](https://arxiv.org/abs/1805.10385v2)

   arXiv 第2版，2018年11月19日。第5.2节与式(72)至(74)：双时间尺度相位展开，以及绝热阶与后绝热阶各需要什么；第5.3节瞬时共振；第6.1节与式(77)至(78)，按轨道平均的平衡律，包括为负的视界能流。

2. [Fujita & Tagoshi (2004), New numerical methods to evaluate homogeneous solutions of the Teukolsky equation, Prog. Theor. Phys. 112, 415](https://arxiv.org/abs/gr-qc/0410018v1)

   arXiv 第1版。表VII，七种自转下圆赤道轨道到 l = 6 的引力波光度；表III，r = 6M 处逐个模式的光度。单位归一化以及“到达无穷远的能流”这一解读是我们的推断，已用四极与后牛顿极限检验。

3. [Tagoshi & Sasaki (1994), Post-Newtonian expansion of gravitational waves from a particle in circular orbit around a Schwarzschild black hole, Prog. Theor. Phys. 92, 745](https://arxiv.org/abs/gr-qc/9405062v1)

   arXiv 第1版。式(43)：到四阶后牛顿的无穷远能流与四极值之比；只用作能流表弱场一端的独立检验。

## 继续阅读

极端质量比旋近 · 阅读路径

### 小黑洞为大黑洞测绘

1. [文章 **绕上几十万圈的轨道，能读出大黑洞的什么**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/what-hundreds-of-thousands-of-orbits-reveal/)
2. [文章 **自转黑洞周围的轨道怎样进动，又在哪里失去稳定**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/geodesics-precession-and-the-last-stable-orbit/)
3. [文章 **能量流失怎样决定最后稳定轨道前还剩多少圈**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/radiation-reaction-and-the-adiabatic-inspiral/)
4. [文章 **偏心轨道的引力波为什么是一排谱线**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/fundamental-frequencies-and-harmonics/)
5. [当前页面 **用已发表的能流表复算一条 EMRI 频率演化**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/recompute-a-frequency-evolution/)
6. [文章 **怎样读 EMRI 的探测预测、自力、共振与多极矩检验文献**](https://sandboxphysics.com/zh/stories/extreme-mass-ratio-inspirals/reading-emri-forecasts-self-force-resonances-multipoles/)
