Multi-Purpose Particle and Heavy Ion Transport code System PHITS Multi-Purpose Particle and Heavy Ion Transport code System 便利な機能 PHITS講習会　入門実習 2015年10月改訂 Title 1

実習内容 [transform] [magnetic field] [counter] Table of Contents 2 PHITS講習会　入門実習 Table of Contents 2

[transform]セクション ソース、surfaceやcellの定義、タリーのr-zやxyzメッシュ、磁場の定義等の際に、平行移動や回転を行うことが可能。 Z Y 回転 X XYZ座標系 平行移動 [transform] 3

入力形式 nは座標変換番号 O1,O2,O3は平行移動を表すX,Y,Z成分 B1～B9は回転行列の各成分 Mで変換式を選択 （ここではM=1のみ扱う） [transform] 4

入力形式 ただし、毎回方向余弦を用いて回転行列を設定するのは 面倒なので、下記のサンプルを利用する。 [transform] 5 set: c10[0] $ angle of around Z (degree) set: c20[0] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 0 cos(c10/180*pi)*cos(c20/180*pi) sin(c10/180*pi)*cos(c30/180*pi)+cos(c10/180*pi)*sin(c20/180*pi)*sin(c30/180*pi) sin(c10/180*pi)*sin(c30/180*pi)-cos(c10/180*pi)*sin(c20/180*pi)*cos(c30/180*pi) -sin(c10/180*pi)*cos(c20/180*pi) cos(c10/180*pi)*cos(c30/180*pi)-sin(c10/180*pi)*sin(c20/180*pi)*sin(c30/180*pi) cos(c10/180*pi)*sin(c30/180*pi)+sin(c10/180*pi)*sin(c20/180*pi)*cos(c30/180*pi) sin(c20/180*pi) -cos(c20/180*pi)*sin(c30/180*pi) cos(c20/180*pi)*cos(c30/180*pi) 1 Z軸, Y軸, X軸の周りに回転させる角度 平行移動のX,Y,Z成分 [transform] 5

例題 transform.inp [transform] 6 PHITSを実行 （最初は体系を確認） track_xy.eps [ M a t e r i a l ] mat[1] 1H 2 16O 1 [ S u r f a c e ] 10 so 500. 11 cz 5. 12 pz 0. 13 pz 10. [ C e l l ] 100 -1 10 101 1 -1. -11 12 -13 110 0 -10 #101 Y軸の周りに30度回転させてみよう！ track_xy.eps 半径5cm, 高さ10cm の円柱 track_xz.eps [transform] 6

例題（回転） transform.inp [transform] 7 Z軸, Y軸, X軸周りに回転させる角度 [ S u r f a c e ] 10 so 500. 11 cz 5. 12 pz 0. 13 pz 10. [ C e l l ] 100 -1 10 101 1 -1. -11 12 -13 trcl=1 110 0 -10 #101 [ T r a n s f o r m ] set: c10[0] $ angle of around Z (degree) set: c20[30] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 0 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ [ T r a n s f o r m ] set: c10[0] $ angle of around Z (degree) set: c20[0] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 0 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ Z軸, Y軸, X軸周りに回転させる角度 101番のセルに座標変換番号1の変換を適用する。 平行移動のX,Y,Z成分 Y軸の周りに30度回転させてみよう！ track_xz.eps [transform] 7

例題（平行移動） transform.inp [transform] 8 Z軸, Y軸, X軸周りに回転させる角度 平行移動のX,Y,Z成分 set: c10[0] $ angle of around Z (degree) set: c20[30] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 50 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ [ T r a n s f o r m ] set: c10[0] $ angle of around Z (degree) set: c20[30] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 0 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ Z軸, Y軸, X軸周りに回転させる角度 平行移動のX,Y,Z成分 M=1の場合、回転した後に平行移動を行う。 Z軸の正の方向に50cm平行移動させてみよう！ track_xz.eps [transform] 8

課題1 transform.inp [transform] 9 Y軸の周りに45度回転させ、X軸の正の方向に10cm、Z軸の正の方向に45cm平行移動させよう！ [ T r a n s f o r m ] set: c10[0] $ angle of around Z (degree) set: c20[45] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 10 0 45 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ [ T r a n s f o r m ] set: c10[0] $ angle of around Z (degree) set: c20[30] $ angle of around Y (degree) set: c30[0] $ angle of around X (degree) tr1 0 0 50 ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ track_xz.eps [transform] 9

課題2 [transform] 10 更に3つのセルを追加して、右下の図のような体系を作ってみよう！ [surface]セクションは変更しない。 変更するのは[cell]セクションと[transform]セクションのみ。 （ [transform]セクションでは、set: c10[*]の行からコピーして、座標変換番号2,3,4の変換を設定する。） [transform] 10

実習内容 [transform] [magnetic field] [counter] Table of Contents 11 PHITS講習会　入門実習 Table of Contents 11

[magnetic field]セクション 磁場(2重極場及び4重極場)を任意の領域に設定して荷電粒子の軌道に対する影響を調べる。 2重極磁場： 1方向に一様な磁場 F=ev×B [magnetic field] 12

[magnetic field]セクション 磁場(2重極場及び4重極場)を任意の領域に設定して荷電粒子の軌道に対する影響を調べる。 4重極磁場： 広がった荷電粒子のビームを収束させる N極 S極 X軸 Y軸 4重極磁石（参考） [magnetic field] 13

入力形式 [magnetic field] 14 2重極磁場はy軸の正方向 reg: 磁場を置く領域 →正の荷電粒子がz軸方向に進む時、 typ: =2; 2重極磁場, =4; 4重極磁場 gap: 磁石の間の距離の半分[cm] （2重極の時は使わないが値を入れる） mgf: 磁場の強さ[kG] trcl: 座標変換番号 2重極磁場はy軸の正方向 →正の荷電粒子がz軸方向に進む時、 　x軸の負の方向に動く 4重極磁場は 　x軸方向に収束、y軸方向に分散する （z軸を中心軸とし無限に4重極磁石を配置した時の磁場として定義されるため、体系にあわせて座標変換が必要。） [parameters]セクションにおいてimagnf=1 [magnetic field] 14

例題 magfield.inp [magnetic field] 15 体系を確認した後、icntl=0として輸送計算を実行。 [ S o u r c e ] s-type = 1 proj = proton dir = 1 r0 = 2.5 z0 = -10. z1 = -10. e0 = 290 [ S u r f a c e ] 10 so 500. 11 cz 5. 12 pz 0. 13 pz 10. [ C e l l ] 100 -1 10 101 1 -1. -11 12 -13 trcl=1 102 1 -1. -11 12 -13 trcl=2 103 1 -1. -11 12 -13 trcl=3 104 1 -1. -11 12 -13 trcl=4 110 0 -10 #101 #102 #103 #104 track_xz.eps [magnetic field] 15

例題（2重極磁場） magfield.inp [magnetic field] 16 track_xz.eps [ P a r a m e t e r s ] icntl = 0 maxcas = 100 maxbch = 10 imagnf = 1 file(6) = phits.out [magnetic field] reg typ gap mgf 104 2 10 100 [ P a r a m e t e r s ] icntl = 0 maxcas = 100 maxbch = 10 c imagnf = 1 file(6) = phits.out [magnetic field] reg typ gap mgf 104 2 10 100 track_xz.eps [magnetic field] 16

例題（4重極磁場） magfield.inp [magnetic field] 17 track_xz.eps [ P a r a m e t e r s ] icntl = 0 maxcas = 100 maxbch = 10 imagnf = 1 file(6) = phits.out [magnetic field] reg typ gap mgf 104 4 10 100 [ P a r a m e t e r s ] icntl = 0 maxcas = 100 maxbch = 10 imagnf = 1 file(6) = phits.out [magnetic field] reg typ gap mgf 104 2 10 100 track_xz.eps [magnetic field] 17

課題1 magfield.inp [magnetic field] 18 4重極磁場の場合、x軸方向には収束するが、y軸方向には分散する。これを確認してみよう！ magfield.inp [ T - T R A C K ] mesh = xyz x-type = 2 nx = 100 xmin = -25. xmax = 25. y-type = 2 ny = 100 ymin = -25. ymax = 25. z-type = 1 nz = 3 30.0 40.0 50.0 60.0 part = proton ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ file = xy_track_proton.dat [ T - T R A C K ] mesh = xyz x-type = 2 nx = 100 xmin = -25. xmax = 25. y-type = 2 ny = 100 ymin = -25. ymax = 25. z-type = 1 nz = 1 -5.0 5.0 part = proton ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ ・ file = xy_track_proton.dat track_xz.eps track_xy.eps [magnetic field] 18

課題2 2重極と4重極の磁場を使って、領域101を図の様に収束したビームが通るように設定してみよう！（ただし、y軸方向の分散は考えない） 磁場を配置する領域を2つ用意する 磁場の方向や強さを変えてビームをコントロールする [magnetic field] 19

実習内容 [transform] [magnetic field] [counter] Table of Contents 20 PHITS講習会　入門実習 Table of Contents 20

Counterとは？ PHITSでは、各粒子がエネルギーや位置情報など以外にカウンターという値を持っています。カウンター値を特定のイベント（領域に入る、領域から出る、その領域で散乱（反応）する、境界で反射する）が起きる度に変化させることにより、タリーした粒子がどのような経緯でそこに到達したかを調べることができます。 タリー 注目している領域からの寄与のみをタリーしたい →　注目している領域で発生した粒子のカウンター値を増やす この領域で散乱した場合にcounter値を+1 [counter] 21

入力形式 [counter] 22 カウンター番号(3つまで) カウンターをつける粒子を指定 (*partはカウンターをつけない粒子を指定) reg: カウンターを回す領域 in: 領域に入った時 out: 領域から出た時 coll: 領域で散乱（反応）した時 ref: 領域の境界で反射した時 0: 変化なし 10000: ゼロセット それ以外: 指定した数だけ増減 [counter] 22

例題 counter.inp [counter] 23 体系を確認した後、icntl=0として輸送計算を実行。 track_xz.eps [ M a t e r i a l ] MAT[ 1 ] 1H 2 16O 1 MAT[ 2 ] 14N 8 16O 2 MAT[ 3 ] 63Cu 0.6915 65Cu 0.3085 [ C e l l ] 100 -1 10 101 1 -1. -11 12 -13 trcl=1 102 1 -1. -11 12 -13 trcl=2 103 2 -1.20e-03 -11 12 -13 trcl=3 104 3 -8.93 -11 12 -13 trcl=4 110 0 -10 #101 #102 #103 #104 ・ ・ ・ ・ ・ ・ [ C o u n t e r ] counter = 1 part = proton reg in out coll ref 102 0 0 1 0 track_xz.eps [counter] 23

例題 counter.inp [counter] 24 領域110から101に入る粒子をタリーする cross.eps [ T - C r o s s ] title = Particle current using [T-cross] tally mesh = reg reg = 1 non r-in r-out area 1 110 101 1.0 e-type = 3 emin = 1.0 emax = 500. ne = 30 unit = 1 axis = eng file = cross.dat output = current part = proton angel = ymin[1.E-06] ymax[1.E-02] epsout = 1 領域110から101に入る粒子をタリーする cross.eps [counter] 24

例題 counter.inp [counter] 25 タリーする粒子にカウンター値に関する条件を加えてみよう！ ctmin(1)=1, ctmax(1)=1 counter.inp [ C o u n t e r ] counter = 1 part = proton reg in out coll ref 102 0 0 1 0 [ T - C r o s s ] ・ ・ ・ ・ ・ ・ file = cross-c1.dat output = current part = proton angel = ymin[1.E-06] ymax[1.E-02] epsout = 1 ctmin(1) = 1 ctmax(1) = 1 cross-c1.eps タリーする粒子のカウンター値の下限と上限 領域102で散乱される陽子のみを区別して評価した [counter] 25

課題1 counter.inp [counter] 26 タリーする粒子に中性子(neutron)を加え、中性子の振る舞いを調べてみよう！ part = proton reg in out coll ref 102 0 0 1 0 [ T - C r o s s ] ・ ・ ・ ・ ・ ・ file = cross-c1.dat output = current part = proton angel = ymin[1.E-06] ymax[1.E-02] epsout = 1 ctmin(1) = 1 ctmax(1) = 1 [ C o u n t e r ] counter = 1 part = proton reg in out coll ref 102 0 0 1 0 [ T - C r o s s ] ・ ・ ・ ・ ・ ・ file = cross-c1.dat output = current part = proton neutron angel = ymin[1.E-06] ymax[1.E-02] epsout = 1 ctmin(1) = 1 ctmax(1) = 1 cross-c1.eps 中性子も同様に区別して評価することが可能 [counter] 26

課題2 [counter] 27 領域102の場合と同様に、103と104のそれぞれで陽子が反応したイベントをタリーしてみよう！ 領域103と104で起こるイベントに対応する様に、カウンター2と3を追加する [t-cross]を2回コピー&ペーストし、それぞれカウンター2と3がタリー条件となるように設定する（出力file名も変更） タリー [counter] 27

課題2の答え counter.inp [counter] 28 cross-c2.eps cross-c3.eps ・ ・ ・ ・ ・ ・ counter = 2 part = proton reg in out coll ref 103 0 0 1 0 counter = 3 104 0 0 1 0 [ T - C r o s s ] file = cross-c2.dat output = current part = proton neutron angel = ymin[1.E-06] ymax[1.E-02] epsout = 1 ctmin(2) = 1 ctmax(2) = 1 cross-c2.eps cross-c3.eps 空気の領域ではほとんど反応しない 銅の領域から散乱される陽子は少ない [counter] 28

まとめ [transform]セクションを利用することにより、複雑な ジオメトリーやソース、タリーの設定が容易になる。 [magnetic field]セクションを用いて磁場を設定する ことにより、加速器施設等でのビームに対する磁場 の影響をシミュレーションできる。 [counter]セクションで適切な設定を行うことにより、 シミュレーションの重要な利点である結果の詳細な 分析が可能となる。 《休憩はさむ》 まとめ Summary 29