Fermi推定の推定 -

Fermiが書いた文書My Observations During the Explosion at Trinity on July 16, 1945では、1945年7月に行われたトリニティ実験に於ける原爆のエネルギー推定を行っているが、どういう計算をしたのか詳細は書かれていない。以下に、数値情報が含まれている部分を抜粋する。

On the morning of the 16th of July, I was stationed at the Base Camp at Trinity in a position about ten miles from the site of the explosion.

After a few seconds the rising flames lost their brightness and ...

About 40 seconds after the explosion the air blast reached me. I tried to estimate its strength by dropping from about six feet small pieces of paper before, during and after the passage of the blast wave. Since at the time, there was no wind I could observe very distinctly and actually measure the displacement of the pieces of paper that were in the process of falling while the blast was passing. The shift was about 2 1/2 meters, which, at the time, I estimated to correspond to the blast that would be produced by ten thousand tons of T.N.T.

Fermiがどういう計算をしたか疑問に思う人はいるようで、How did Enrico Fermi calculate the classical Fermi Problem?で、回答している人もいるけど、(私が判断する限り)正しいかどうかは疑わしい。

そもそも、Fermiが計算しようとした量は、原理的には以下のようなものだと思われる。

$E = \displaystyle \int_{0}^{R(t)} \left( \dfrac{1}{2}\rho u^2 + \dfrac{p-p_{0}}{\gamma-1} \right) 2\pi r^2 dr$

$r$ は、爆発地点からの距離。 $R=R(t)$ は、時刻 $t$ に於いて爆風が到達した爆発地点からの距離。 $\rho(r,t) , u(r,t) ,p(r,t) , p_{0} , \gamma$ は、それぞれ、大気の密度、流速、圧力、爆発前の大気圧、比熱比。

爆発は、半球状に広がり、地表からの反射があるので、正確には球対称性はないと思われるが、ここでは方向依存性は無視した。 $E$ は時刻 $t$ ごとに計算されるが、エネルギー保存則によって一定と考えて差し支えないだろう(実際は、粘性や地表との摩擦によって減衰するだろうが)。爆発がなければ、被積分関数は恒等的に0(流速は0で、圧力は至る所大気圧と等しい)で、どこまで積分しても0になる。

上の式には、輻射に関する項が含まれてない。今日、物理のいくつかの分野では、輻射を伴う流体力学を使うことがあって、輻射流体力学(radiation hydrodynamics)と呼ばれることもある。輻射流体シミュレーションの基礎という解説記事には『現在，世界的に見ると，輻射流体コードと呼ばれるものは宇宙物理と慣性核融合のコミュニティを中心として多数存在する．』と書いてる。

輻射流体力学の研究が本格的に進展したのは戦後のことだと思うが、Betheらによる報告書BLAST WAVE(PDF)でも、輻射に関する考察は、かなり行われている。しかし、私の知る限り、輻射の影響を簡単に見積もる方法は確立されてない。報告書の中で、Betheは、

In the Trinity test it was measured that about 15 to 20 percen t of the total energy was emitted in the form of radiation which could penetrate the air to a distance of several miles.

と書いている。このことの妥当性はともかく、今の問題はFermiの推定方法で、輻射に関する見積もりをできたとは考えにくいので、輻射を無視することにする。

G.I. Taylorは、論文The formation of a blast wave by a very intense explosion. - II. The atomic explosion of 1945で、輻射を無視した爆発エネルギーの大きさを見積もり、それは上式で与えられるものだった(正確には、上式は半球上のみの伝播を考えているのに対して、G.I.Taylorは全方位に伝播する球対称解を考えているので、上式の2倍)。

Taylorは、今日、Sedov-Taylor自己相似解と呼ばれる理論解で仮定される相似則が(少なくとも最初の数十ミリ秒の間は)実際に成立していることを、写真の解析から示し、爆発エネルギーを、 $\rho(r,t) , u(r,t) ,p(r,t)$ の理論解に基づいて計算した結果は、上エネルギーを16800TNTトン相当だとした。

Sedov-Taylor自己相似解では、爆発エネルギー $E$ は、

$E \propto \dfrac{R^5 \rho_0}{t^2}$

となり、比例係数は、比熱比に依存するが、大体1くらいである(比例係数は無次元の量になる)。 $\rho_0$ は標準状態での大気の密度。

G.I. Taylorの論文TABLE Iのデータを使うと、、例えば、t=62(msec)の時、R=185.0(m)とあり、空気の密度を1.293(kg/m^3)とすると、

E ≒ (185**5)*1.293/(62/1000.0)**2/4.184e9 = 17421(tons of TNT)

などと計算でき、G.I.Taylorの見積もりに近い値が出る。

長崎や広島の原子爆弾は、高度600メートル付近で爆発させたようなので、トリニティ実験も、同程度の高さで爆発させたとすれば、圧力波が600メートルほど広がるまでは、球対称解は十分良い近似だったと考えられ、R=185.0(m)で、この理論を当てはめるのは問題ないように思う。

G.I.Taylorとvon Neumannは、1941年には、Sedov-Taylor自己相似解を独立に得ていたようだが、公開されたのは、いずれも終戦後らしい。Fermiは、マンハッタン計画に参加していたので、このような結果は知ってたと思われる。が、Fermiのいた地点で、この式を適用しても法外な値が出てくる。実際、t=40+few secondsで、R=16(km)を代入すると、E=1.3億TNTトンというような値になってしまう。それに、これだと、わざわざ紙を落とした意味がない。

Fermiが紙を落として計測しようとした量は、多分、衝撃波が通り過ぎた直後の風速じゃないかと思う。高さ6フィート(1.8m)の地点から、紙片を落として、落下地点が2.5mずれたとあるので、落下時の空気抵抗を無視すれば、落下に要する時間が0.6秒で、平均風速4(m/s)程度と計算される。正しい風速の計算は、紙片の重量や形状に依存して厄介だが、軽い物体なので、殆ど瞬時に風速まで加速され、一定速度に達したと見なしたのだろう。

爆風が通過した直後の風速は、Rankine–Hugoniotの理論で、よくparticle velocityと呼ばれてる量である。衝撃波面後方の圧力、密度、内部エネルギーを $p_1,\rho_1,e_1$ として、衝撃波の速度を $u_s$ 、particle velocityを $u_p$ とすると

$\rho_0 u_s = \rho_1 (u_s-u_p)$

$p_1-p_0 = \rho_0 u_s u_p$

$e_1-e_0 = \dfrac{1}{2}(p_1+p_0)(\dfrac{1}{\rho_0} - \dfrac{1}{\rho_1})$

が成り立つ。 $\rho_0,p_0,e_0$ は、標準状態での大気の密度、圧力、内部エネルギー。ちょっと頑張ると、以下のような式が得られる。

$\dfrac{\rho_1}{\rho_0} = \dfrac{(\gamma-1)p_0 + (\gamma+1)p_1}{(\gamma+1)p_0+(\gamma-1)p_1} = \dfrac{2\gamma p_0 + (\gamma+1)(p_1-p_0)}{2\gamma p_0+(\gamma-1)(p_1-p_0)}$

$u_p = \dfrac{c_0 p_1}{\gamma p_0} \dfrac{1}{\sqrt{1+\dfrac{\gamma+1}{2\gamma}\dfrac{p_1-p_0}{p_0}}}$

$c_0$ は音速。物理の教科書では、あんまり見ない気がする式だけど、 $(\Delta p_{S}) = p_1-p_0$ は最大静止過圧に相当する量になっている。で、最大動圧 $q = \dfrac{1}{2} \rho_1 u_p^2$ が

$q = \dfrac{ (\Delta p_{S})^2 }{2\gamma p_0 + (\gamma-1)(\Delta p_{S})}$

と書ける。Fermiが生身で観測して平気だったということは、 $(\Delta p_{S}) << p_0$ だったと思っていい。そうすると、 $\rho_1 \approx \rho_0$ で、 $u_p$ が4(m/s)だった時、 $q$ は約10(N/m^2)で、 $(\Delta p_{S})$ は1670(N/m^2)くらい。0.02気圧以下なので、この程度の爆風は、人体に損傷を与えない

Taylorのエネルギー式に戻ると、 $(\Delta p_{S})$ は $r$ が $R$ より少し小さい地点での $p-p_0$ に相当する量になっている。爆風が通り過ぎると、比較的素早く大気圧に戻ると予想されるので $r=R$ に近い領域で積分すれば十分で、また、Taylorのエネルギーは第二項の方が支配的になる(これについては、Taylorの論文The formation of a blast wave by a very intense explosion I. Theoretical discussionの式(31)と(32)を参照)。両者を勘案して、積分領域は狭いので、積分領域の $p(r,t)-p_0$ を一定値 $(\Delta p_{S})$ に置き換えてみると、以下のようになる。

$E \approx \dfrac{(\Delta p_{S})}{\gamma-1} \displaystyle \int_{R-\epsilon}^{R} 2\pi r^2 dr \approx \dfrac{2 \pi (\Delta p_{S}) R^2}{\gamma-1} \epsilon$

$\epsilon$ は適当なカットオフで、 $r \lt R-\epsilon$ では、圧力 $p$ は、十分、大気圧 $p_0$ に等しいと見なせるくらいの数値に取る。 $R$ は16(km)という大きな距離なので、 $\epsilon \lt\lt R$ としている。

$\epsilon$ は、どれくらいに取るべきだろうか。オーダー計算なので、1000(m)なのか100(m)なのか、10(m)なのかという話。そもそも、紙片が落下するのに0.6秒くらいかかってて、その間、風が吹いてたなら、その間に爆風が伝播した距離340(m/s)×0.6(s)= 204(m)くらいに取るべきという気もする(衝撃波の伝播速度 $u_s$ はほぼ音速に等しい)。

一方、実際の $p(r,t)-p_{0}$ は、 $r$ が $r=R$ から少し離れるだけで、急激に小さくなる(通常、過圧は、時間の指数関数オーダーで減衰するとされる)が、 $(\Delta p_{S})$ をそのまま掛けてるので、 $\epsilon$ は小さめの10(m)とかにするのが正しそうにも思える。 $\epsilon$ を10(m)とした場合、 $E$ は、およそ16000(tons of TNT)になる。

あるいは、以下のように考えてみることもできる。

$p(r) \sim p_0 + (\Delta p_{S}) \exp(-B(r-R)/\epsilon)$

という形を仮定してみる(衝撃波が通り過ぎた時、圧力は一旦上昇した後、大気圧以下まで下がってから戻るので、これは、あんまり正確ではなく積分を過大評価しそうであるが)と、主要な項は

$\dfrac{ 2 \pi (\Delta p_{S})}{\gamma-1} \displaystyle \int_{R-\epsilon}^{R} r^2 e^{-B(r-R)/\epsilon}dr \sim \dfrac{2 \pi (\Delta p_{S}) R^2}{\gamma-1} \dfrac{\epsilon}{B}$

と書ける。 $B$ は減衰の早さを表すが、 $r\lt R-\epsilon$ では、 $p(r)$ は $p_0$ に十分近くなっているとしていた。 $(p(R-\epsilon)-p_0)/(p(R)-p_0) = e^{-B}$ であるから、 $B=1$ では、不十分だろう。 $B=5$ なら、比率は1%以下になる。

所詮、オーダー計算なので、あまり深く悩んでも仕方ない。 $\epsilon$ が100(m)で、 $B$ を10にすれば、 $E$ は、およそ16000(tons of TNT)になる。

G.I.Taylorの見積もりを知ってるから、それらしい数字を選べるけれど、正直なところ、一点に於ける最大静止過圧が分かっても、積分を正しく見積もるのは、結構難しいように思える。爆発は、半球状に広がり、方向依存性がないとしたのも、妥当な近似かどうかは疑わしい。

マンハッタン計画では、当然、原子爆弾の威力についても、爆風の伝播についても、事前に、かなりの検討が行われたはずだから、Fermiは、それらの知識を動員した上で、見積もりを行ったのかもしれない。何らかの実験データに基づいた数値をFermiが使ったということがあれば、見積もりの詳細を書かなかったことも理解できる。

Betheらによる報告書BLAST WAVEのChapter 5は、爆発地点から遠く離れた場所での爆風の伝播を考察していて、式(5.144)及び(5.144a)に、数値計算(IBM calculationなどと書かれてるが、マンハッタン計画では、電子計算機以前のIBMの"計算機"が使用されてたらしいので、それのことだろう)の結果に基づいて

$\dfrac{(\Delta p_{S})}{p_{0}} = \dfrac{260.7}{R\sqrt{\ln(R/260.1)}}$

$\dfrac{(\Delta p_{S})}{p_{0}} = \dfrac{278.6}{R\sqrt{\ln(R/129.1)}}$

などの(経験)式が書かれている。Rの単位はメートルで、10000~13500TNTトン相当のエネルギーが爆発によって解放されたと仮定した計算(エネルギーは保存せず増加していったらしい)だったようだ。この数値計算も、圧力波伝播の方向依存性はないものと仮定した計算になっている。当時の計算機の性能では、空間次元を減らすことは、とても重要だったのだろう。

Fermiのいた地点が、R=16000(m)であるから、どちらの式を使っても、最大静止過圧/大気圧比 $\frac{(\Delta p_{S})}{p_{0}}$ は、0.008程度。一方、既に書いたように、最大静止過圧は、風速から見積もることができ、 $u_p$ が4(m/s)とした時、 $(\Delta p_{S})$ は1670(N/m^2)くらいで、 $(\Delta p_{S})/p_{0}$ は、0.0167程度。数値計算とは2倍くらいずれているが、爆発のエネルギーが数値計算と同程度のオーダーだったと考えても差し支えないだろう。

Fermiは、単に、この結果を知っていて、10000TNTトン相当と結論付けたのかもしれない。