PHITS2.52 の特徴
PHITS2.52 は，2012 年 12 月に公開された最新版 PHITS パッケージであり，従来，原子力コードセンターなどに登録されていた PHITS2.24 パッケージと比べて，計算機能や精度，ユーザーサポートツールなど，あらゆる面が改良されている。PHITS2.24と比較した PHITS2.52 の特徴をコード自体及びパッケージの特徴に分類して表１及び２に示す。以下，いくつかの項目に関して，その詳細を解説する。

1 最新の核内カスケードモデルの組込核反応モデルは，放射線挙動解析計算コードの「肝」とも言うべき部分であり，この精度が計算コード全体の優劣を決定すると言っても過言ではない。PHITS は，量子分子動力学模型に基づいて全ての核子間相互作用を解析する核反応モデル JQMD[1]や，核内カスケード模型に基づいて数 GeV 以上の共鳴状態を丁寧に再現する核反応モデル
JAM[2]を世界に先駆けて導入したため，重イオン入射や高エネルギー核子入射核反応に対する計算精度は良い。しかし，PHITS2.24 では，数 10MeV から数 GeV の核子入射に対して，古いタイプの核内カスケード模型 Bertini に改良を加えたモデル[3]を初期設定として採用していたため，その計算精度は最新のモデルと比較してあまり良くないことが報告されていた[4]。
そこで，PHITS2.52 には，中間領域エネルギーに対する最新の核内カスケード模型を２つ，それぞれ独立して組み込んだ。１つはフランス CEA が中心となって開発している INCL[5]で，もう１つは九州大学で開発している INC-ELF[6]である。どちらも，近年注目されているコアレッセンス模型を組み込んでおり，従来のカスケード模型では全く再現できなかった高エネルギーフラグメント生成を再現することができる。例として，JAM，INCL，INC-ELF モデルでそれぞれ計算した陽子入射に対する中性子生成及び重陽子生成反応断面積を図１に示す。図より，中性子生成に対しては，どのモデルもほぼ実験結果を再現可能だが，重陽子生成に関しては，コアレッセンス模型が組み込まれていない JAM モデルは大幅に過小評価してしまうことが分かる。また，INC-ELFの重陽子生成断面積が高エネルギー側で高いのは，複合粒子のノックアウト過程も考慮しているためである。
PHITS2.52 では，モデルの計算時間や汎用性を検討した結果，中間エネルギーの陽子・中性子・π粒子・重陽子・三重陽子・3He・α粒子入射反応に対する初期設定モデルとして INCL を採用することとした。これらの改良は，原子力機構と CEA 及び原子力機構と九州大学の共同研究による成果である。

2 統計マルチフラグメンテーションモデルの組込
PHITS では，JQMD，JAM，INCL など各核子の挙動を運動力学的に扱う動的過程

モデルと，動的過程を経て平衡状態に達した原子核からの核子放出や核分裂など統計学的に扱う静的過程モデルを組み合わせて核反応を模擬する。この静的過程モデルには，
通常，蒸発モデル GEM[7]が用いられる。しかし GEM は，核分裂は扱えるものの，高い励起エネルギーの場合に原子核が３つ以上に分裂するマルチフラグメンテーション反応を模擬できず，重核からの軽核生成反応断面積を過小評価してしまう問題点が指摘されていた[8]。
そこで，PHITS2.52 では，その静的過程シミュレーションに統計マルチフラグメン
テーションモデル SMM[9]を新たに組み込み，原子核が 3 つ以上に分裂する反応を模擬できるようにした。具体的には，動的過程直後の残留核励起エネルギーが 2MeV/u 以上の場合，SMM を起動して数多くある分裂パターンからエントロピー的に取り得る状態をモンテカルロ法により決定し，分裂後の各残留核の蒸発過程を GEM で再現するよう変更した。これにより，PHITS2.24 で見られた一部の残留核収率に対する過小評価が改善された（図２参照）。なお，SMM を導入することにより計算時間が従来の数倍程度になる場合があるため，PHITS2.52 の初期設定では SMM を起動しない設定とした。

3 巨大共鳴反応断面積の改訂
光子のエネルギーが中性子（もしくは陽子）放出のしきい値エネルギー（約 8MeV）を超えると，光核反応が起きる。光核反応には，光子が原子核全体と共鳴する巨大共鳴反応（8～25MeV），核内にある仮想的な重陽子と共鳴する準重陽子崩壊反応（25～150MeV），核内の核子単体と共鳴してπ粒子などを放出する核子共鳴反応（150MeV～）がある（括弧内は，支配的な光子エネルギー）。PHITS2.24 は，これらの反応機構の中で巨大共鳴のみ扱うことができたが，その反応断面積は，単純なローレンツカーブで表現されており，その精度に関する検証は十分に行われていなかった。
そこで PHITS2.52 では，光核反応断面積として，評価済核データライブラリJENDL/PD-2004[10]に格納された値を採用した。例として，PHITS2.24 及び 2.52 で採用した 6Li と 184W に対する光核反応断面積を図３に示す。図より，6Li に対しては新旧 PHITS の光核反応断面積が全エネルギー領域で大きく異なることが分かる。
これは，6Li など軽いターゲットに対する光核反応断面積が単純なローレンツカーブでは再現できないためである。一方，184W に対しては，低エネルギー側ではほぼ一致するものの25MeV 以上の高エネルギー領域で新旧 PHITS の光核反応断面積が異なる。これは，
PHITS2.52 で採用した核データライブラリ JENDL/PD-2004 が準重陽子崩壊反応による寄与も含んでいるためである。ただし，その反応機構そのものはまだ組み込まれていないため，PHITS2.52 を用いても高エネルギー光核反応を精度よく再現することはできない。今後，準重陽子崩壊，核子共鳴反応機構を順次，PHITS に組み込んでいく予定である。

4 新しい反応断面積計算モデルの組込
PHITS2.24 には，核子-原子核及び原子核-原子核反応断面積を計算するモデルとして，Pearlstein-Niita の式，Shen の式，NASA の式などが組み込まれている。これらに加え，PHITS2.52 では，Kurotama 模型[11]を組み込んだ。例として，NASA の式とKurotama 模型で計算した 12C-12C 反応断面積の入射エネルギー依存性を図４に示す。
図より，Kurotama 模型は，粒子線治療で用いられる数 100MeV/u エネルギー領域において，NASA モデルよりも実験値の再現性が良いことが分かる。ただし，初期設定で使うモデルは，PHITS2.52 でも PHITS2.24 と同じく Pearlstein-Niita の式（核子-原子核）と NASA の式（核子-核子）とした。本改良は，小濱洋央氏（理研），飯田圭氏（高知大学），親松和浩氏（愛知淑徳大学）らとの共同開発による成果である。

5 メモリ共有型並列計算機能の導入
並列計算には，MPI を使ったメモリ分散型並列と，OpenMP を使ったメモリ共有型並列がある。PHITS2.24 では，メモリ分散型並列のみ対応していた。しかし，メモリ分散型並列計算では，並列プロセス毎にシングル計算と同等のメモリを使用する（８並列で８倍のメモリ領域を使用する）ため，メモリ不足の問題から高分解能ボクセルファントムなど巨大なメモリを必要とする体系を並列計算で処理できない欠点があった。
そこで，PHITS2.52 では，メモリ共有型並列計算にも適応可能となるようソースコードを抜本的に改良した。また，メモリ分散型並列とメモリ共有型並列を組み合わせたハイブリッド並列計算も可能とした。その結果，ハイブリッド並列を基本とする「京」などのスーパーコンピュータのみならず，MPI プロトコルをインストールしていない汎用の PC でも並列計算が可能となった。例として，メモリ共有型並列計算でコア数を１～８まで増やしたときのスケーラビリティ（コア数が１と N のときの計算時間比：t1/tN）を図５に示す。図より，コア数を増やしてもスケーラビリティは理想値の 90%以上を保持し，コア数にほぼ反比例して計算時間が短縮できることが分かる。ただし，コア数が１の場合のメモリ共有型並列計算の計算時間は，シングル計算と比べて約２倍となるため，コア数が２以下の場合は，メモリ共有型並列計算を実行するメリットはほとんどない。
本成果は，次世代生命体統合シミュレーションソフトウェアの研究開発プロジェクト，理化学研究所戦略的研究展開事業，理科学研究所基礎科学特別研究員制度の支援によって得られたものであり，原子力機構と理化学研究所の共同研究による成果である。また，開発において「京」コンピュータ試験利用および理研情報基盤センターRICC システムを利用した。

6 統計誤差計算方法の改良と再開始計算機能の追加
統計誤差は，モンテカルロ計算精度を検証する上で重要な指標であるが，PHITS2.24

では，誤った計算方法で統計誤差を導出していた。そこで PHITS2.52 では，検出器の応答関数などヒストリー間の分散を計算するタリーでは相対標準誤差を，それ以外のタリーでは各ヒストリーもしくは各バッチ計算結果の相対標準偏差を統計誤差として出力するように全てのタリーを改良した。例として，PHITS の 2 次元プロットで出力した粒子フラックスとその統計誤差分布を図６に示す。図より，フラックスの小さい領域（左図寒色系）では，その相対誤差が大きい（右図暖色系）ことが分かる。このように統計誤差を可視化することにより，統計の足りない領域が直感的に分かるようになり，より効果的な計算が可能となった。
また，一度終了した計算からタリー出力や初期乱数を読み込んで計算を再開する機能を追加した。この機能を用いれば，あらかじめ設定していたヒストリー数やバッチ数では統計が十分でなかった場合に，これまでの計算結果を無駄にすることなく，より統計精度のよい結果を得ることができる。ただし，古いバージョンの PHITS で計算した結果からの再開始計算はできない。これらの改良は，原子力機構・システム計算科学セン
ター・原子力コード高速化作業の一環として，(株) 富士通システムズ・イーストの大日向大地氏らと実施した共同開発の成果である。

7 残留放射能計算機能の追加
残留放射能の時間変化の計算は，加速器の遮へい設計や粒子線治療の医療従事者の２次被ばくなどを評価する際，不可欠となる。そのためには，照射直後に生成された放射性残留核の収率のみならず，その時間減衰を評価する必要がある。しかし，PHITS2.24
では照射直後の残留核収率計算のみ可能で，その時間変化を追跡するためには，ユーザ
ー自身がその計算結果を基に DCHAIN-SP[12]など残留放射能計算コードの入力ファイルを作成する必要があった。
そこで，PHITS2.52 では新たなタリー[t-dchain]を導入し，照射時間や照射後の冷却時間などを指定すれば，PHITS から DCHAIN-SP 用の入力ファイルを直接出力できるようにした。また，PHITS との接続計算用に調整した DCHAIN-SP やそれを実行するためのバッチファイルを整備し，PHITS2.52 パッケージに組み込んだ。これらの成果により，PHITS と DCHAIN-SP を用いて簡単に残留放射能の時間変化を計算可能となった。例として，PHITS2.52 パッケージを用いて計算した 150MeV 陽子を水ファント
ムに 6 分間で 5Gy 照射したときの残留放射能の時間変化を図７に示す。図より，照射直後は 14O や 15O が多数生成されているが，数 10 分後にはそれらはほとんど崩壊し，11C や 13N のみが残ることが分かる。今後は，PHITS と DCHAIN-SP をソースレベルで統合し，一度計算した残留放射能分布からその周辺の被ばく線量を PHITS で再計算する機能を開発する予定である。

8 電子・光子輸送計算アルゴリズムの改良
PHITS2.24 では，電子のカットオフエネルギーを低くすると，高エネルギー電子の挙動が変化してしまう問題点があった。これは，PHITS2.24 で採用している電子の阻止能がカットオフエネルギーに依らず一定のため，大量の 2 次電子を放出する高エネルギー電子に対して阻止能と 2 次電子放出によるエネルギー損失の重複が起こり，その飛程が短くなっていたためである。そこで，PHITS2.52 では，電子の阻止能をカットオ
フエネルギーに依存して変化させ，上記重複が起きないように改良した。
また，PHITS2.24 では，光電効果やコンプトン散乱が起きた際，各イベントでエネルギーが保存しない問題点があった。そこで，PHITS2.52 では，中性子や陽子の場合と同じく，光子や電子の輸送計算でもイベントジェネレータとなるよう改良した。これ
により，付与エネルギーのヒストリー間での分散を正しく計算できるようになった。例
として，PHITS2.24 と PHITS2.52 で計算した 667keV 光子入射に対する CsI 検出器の応答関数を図８に示す。PHITS2.24 ではエネルギー保存が成立しないため光電ピークの幅が広くなっていたが，PHITS2.52 ではその問題が解決されている。

9 核データライブラリの改訂
PHITS で低エネルギー中性子，光子及び電子を輸送するためにはデータライブラリが必要となる。しかし PHITS2.24 パッケージには，JENDL-4.0[13]に基づいて整備した中性子に対する核データライブラリしか付録されておらず，光子や電子に対する原子データライブラリは含まれていなかった。そこで，JENDL-4.0 及び EEDL[14]を元に光子-原子及び電子-原子相互作用ライブラリをそれぞれ整備し，PHITS2.52 パッケージに付録した。また，中性子ライブラリに関しては，一部の核種で Kerma 係数が異常に大きい問題があったため，その評価方法を改訂した。さらに，水など 15 物質に対して熱中性子散乱則 S(α,β)データを収納した。
例として，新旧ライブラリ及び ENDF/V-II.1[15]に格納された 35Cl に対する中性子
Kerma 係数を図９に示す。図より，古いライブラリは，低エネルギー中性子に対するKerma 係数を新しいライブラリや ENDF/B-VII.1 と比較して約 1000 倍も過大評価していたことが分かる。古いライブラリを使って中性子による人体内の被ばく線量を計算した場合，人体に 0.1%程度しか含まれない 35Cl の有無により結果が大きく変わってしまう問題があったが，新しいライブラリを使えば，このような問題は生じない。ただし，このような大きな変化があった核種は希であり，ほとんどの核種に対しては，新旧ライブラリでほぼ同等の結果を与える。

参考文献
[1].

K. Niita et al., Phys. Rev. C52, 2620 (1995)[2].

Y. Nara et al., Phys. Rev. C61, 024901 (1990)

[3].

H. Takada et al., JAERI-Data/Code 98-005 (1998).

[4].

岩瀬広, 核データニュース, 96, 9-17 (2010)[5].

J. Cugnon et al., J. Korean Phys. Soc. 59, 955 (2011)[6].

Y. Sawada et al., Nucl. Instr. Meth. B 291, 38-44 (2012)[7].

S. Furihata et al., Nucl. Instr. Meth. B 171, 251-258 (2000)[8].

T. Ogawa et al., Prog. Nucl. Sci. Technol. in press.

[9].

J.P. Bondorf et al., Physics Reports, 257, 133 (1995)[10]. http://wwwndc.jaea.go.jp/ftpnd/jendl/jendl-pd-2004.html[11]. K. Iida et al., J. Phys. Soc. Japan 76, 044201 (2007)[12]. T. Kai et al., JAERI-Data/Code 2001-016 (2001)[13]. http://wwwndc.jaea.go.jp/jendl/j40/J40_J.html[14]. S. T. Perkins et al., Lawrence Livermore National Lab. UCRL-50400 Vol. 31(1991)[15]. http://www.nndc.bnl.gov/endf/b7.1/[16]. Y. Iwamoto et al., Nucl. Instr. Meth. B 274, 57-64 (2012)[17]. http://rcwww.kek.jp/research/egs/

表１

PHITS2.24 と比較した PHITS2.52 コードの特徴

項目

概要

最新の核内カスケードモデ 中間領域エネルギーの核反応に対する計算精度を高める
ルの組込

た め ， 最 新 の 核 内 カ ス ケ ー ド モ デ ル INCL[5] 及 び（詳細は 2.1 参照）

INC-ELF[6]を組み込んだ。

統計マルチフラグメンテー 核反応により生成する残留核収率の計算精度を高めるたションモデルの組込

め，静的過程に統計マルチフラグメンテーションモデル

（詳細は 2.2 参照）（SMM）[9]を組み込んだ。

巨大共鳴反応断面積の改訂

光核反応データライブラリ JENDL/PD-2004 [10]に基づ（詳細は 2.3 参照）

いて，光核反応断面積を改訂した。

新しい反応断面積計算モデ 理研などで開発された核子-原子核および原子核-原子核
ルの組込

反応断面積計算モデル Kurotama[11]を組み込んだ。

（詳細は 2.4 参照）
メモリ共有型並列計算機能 OpenMP を利用したメモリ共有型並列計算機能を導入の導入し，メモリ分散型並列と組み合わせたハイブリッド並列（詳細は 2.5 参照）計算や，汎用 PC での並列計算を可能とした。

統計誤差計算方法の改良と ヒストリー間の分散を計算するタリーでは相対標準誤差再開始計算機能の追加を，それ以外のタリーでは相対標準偏差を統計誤差とし（詳細は 2.6 参照）て出力するように全てのタリーを改良した。また，統計が足りなかった場合などに対応するため，一度終了した計算を再開始する機能を追加した。

残留放射能計算機能の追加

残留放射能計算コード DCHAIN-SP[12]を PHITS パッ（詳細は 2.7 参照）ケージに組み込むとともに，PHITS から DCHAIN-SPの入力ファイルを自動で作成する[t-dchain]タリーを作成した。

電子・光子輸送計算アルゴ 電子の阻止能をそのカットオフエネルギーに応じて変化リズムの改良させ，高エネルギー電子の挙動がカットオフエネルギー（詳細は 2.8 参照）に依存しないようにした。また，光子や電子の輸送計算でもイベントジェネレータとなるよう改良した。

ユーザー定義タリーの導入

多様なニーズに応じるため，ユーザーが各自の必要とする物理量を自由にスコアするための[t-UserDefined]タリーを作成した。

δ線生成機能の組込

荷電粒子の飛跡周辺にノックアウトされる電子（δ線）の生成を考慮できるように改良し，エネルギー付与の空間的な分布をより精度よく計算可能とした。

DPA 計算機能の改良

放射線による材料損傷の指標となる DPA を計算する際，（詳細は文献[16]）クーロン散乱の効果を適切に考慮できるように改良し，重イオン照射に対しても精度よく DPA を計算できるようにした。

バグ修正と高速化

等方照射線源などに関するバグを修正し，JQMD など一部のルーチンを高速化した。

表２

PHITS2.24 と比較した PHITS2.52 パッケージの特徴

項目

概要

インストーラの整備

Windows 及び Mac 用のインストーラを作成し，インストールフォルダを指定するのみでインストールが完了するようにした。

実行シェルの整備

Windows では右クリックの「送る」コマンドより，Mac
ではアイコンのドラッグ＆ドロップにより PHITS を実行可能とした。

核データライブラリの改訂

JENDL-4.0[13]に基づいて作成した中性子核データライ（詳細は 2.9 参照）ブラリを改訂するとともに，光子及び電子の原子データライブラリを整備した。

講習会用資料の改訂

基礎実習（ジオメトリ・タリー・パラメータ設定），応用（HP よりダウンロード可） 実習（便利な機能・分散低減・ボクセルファントム），総合実習用の資料をそれぞれ作成した。
奨励設定ファイルの整備

ユーザーが PHITS 入力ファイルを作る際の参考とする（HP よりダウンロード可） ため，様々な使用目的に合わせた PHITS パラメータの奨励設定を決定し，そのサンプル入力ファイルを作成した。
マニュアルの改訂

マニュアルの不備を大幅に改訂し，原子力分野以外のユ（HP よりダウンロード可） ーザーでも理解できるようにした。
便利なツールの整備

PHITS で作成した 3 次元体系や放射線挙動をアニメーシ（HP よりダウンロード可） ョン化するツールなど，PHITS を使った便利な機能の紹介とその使用方法の解説を作成した。
PHITS+EGS5 の暫定版を PHITS と EGS5[17]を統合したコードの暫定版実行ファ追加

イルを作成した。ただし，メモリ共有型並列計算や再開始計算機能など PHITS2.52 から組み込まれた様々な機能に未対応であり，光核反応も考慮しない。

56Fe(p,d)

208Pb(p,n)

JAM，INCL，INC-ELF モデルでそれぞれ計算した 208Pb-陽子(258 MeV)反応の

図１

中性子生成断面積（左図）及び 56Fe-陽子(558 MeV)反応の重陽子生成断面積（右図）

24Na

Cross‐section (mb)

1.E+02102

SMMあり
SMMなし

1.E+02102

Experiment(Ogawa et al.)Cross‐section (mb)

1.E+03103

1.E+00100

1.E+00100

1.E-0110‐100

図２Experiment(Ogawa et al.)

1.E+01101

1.E+0110110‐21.E-02

SMMあり
SMMなし

75Se

100100NatPb(C,x)

200200

Energy (MeV/u)

300300400400

1.E-0110‐100100100200200300300

Energy (MeV/u)

フラグメンテーション反応による 24Na と 75Se 生成断面積の入射エ

ネルギー依存性

図３400400

PHITS2.24 及び 2.52 で採用した 6Li と 184W に対する光核反応断面積

図４NASA の式と Kurotama0 で計算した 12C-12C 反応断面積の入射エネルギー依存性。Kurotama0 とは，高エネルギー側で有効な Kurotama 模型に基づく公式[18]を，低エネルギー側で NASA の式により補完したものである。

図５ メモリ共有型並列（OpenMP）計算でコア数を１～８まで増やしたときのスケー
ラビリティ（コア数が１と N のときの計算時間比：t1/tN）

図６PHITS の 2 次元プロットで出力した W ターゲットから発生する X 線フラック

ス（左図）とその統計誤差分布（右図）

Irradiation

図７Cooling

PHITS2.52 パッケージを用いて計算した 150MeV 陽子を水ファントムに 6 分間で 5Gy 照射したときの残留放射能の時間変化

Count / MeV / source

PHITS 2.24
PHITS 2.52
10100−20

図８0.20.4
Deposit energy (MeV)

0.6

PHITS2.24 と PHITS2.52 で計算した 667keV 光子入射に対する CsI 検出器の応答関数

図９ 新旧ライブラリ及び ENDF/V-II.1[22]に格納された 35Cl に対する中性子 Kerma係数