Transcript 応用実習1：便利な機能

PHITS
Multi-Purpose Particle and Heavy Ion Transport code System
PHITS 応用実習①：
便利な機能
2014年8月改訂
Title
1
実習内容
1. [transform]
2. [magnetic field]
3. [counter]
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
入力形式
ただし、毎回方向余弦を用いて回転行列を設定するのは
面倒なので、下記のサンプルを利用する。
[Transf
set: c10[0]
set: c20[0]
set: c30[0]
orm]
$ angle of around Z (degree)
$ angle of around Y (degree)
$ angle of around X (degree)
Z軸, Y軸, X軸の周りに
回転させる角度
tr1 0 0 0
平行移動のX,Y,Z成分
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
[transform]
5
例題
PHITSを実行
transform.inp （最初は体系を確認）
[Material]
mat[1] 1H 2 16O 1
[Surface]
10 so
500.
11 cz
5.
12 pz
0.
13 pz
10.
[Cell]
100 -1
101
1 -1.
110
0
track_xy.eps
Y軸の周りに30度回
転させてみよう！
10
-11 12 -13
-10 #101
半径5cm,
高さ10cm
の円柱
track_xz.eps
[transform]
6
例題（回転）
transform.inp
[Transform]
set: c10[0] $ angle of around Z (degree)
set: c20[0]
c20[30] $$angle
angleofofaround
aroundYY(degree)
(degree)
set: c30[0] $ angle of around X (degree)
tr1 0 0 0
・・・・・・・・・・・・・・・・・・
[Surface]
10 so
500.
11 cz
5.
Z軸,
Y軸,
12 pz
0.X軸周りに回転
13 pz
10. 101番のセルに座標変換
させる角度
番号1の変換を適用する。
[Cell]
平行移動のX,Y,Z成分
100 -1
10
101
1 -1.
-11 12 -13 trcl=1
110
0
-10 #101
Y軸の周りに
30度回転させ
てみよう！
track_xz.eps
[transform]
7
例題（平行移動）
transform.inp
[Transform]
set: c10[0] $ angle of around Z (degree)
set: c20[30] $ angle of around Y (degree)
set: c30[0] $ angle of around X (degree)
Z軸, Y軸, X軸周りに回転
させる角度
tr1 0 0 0
50
平行移動のX,Y,Z成分
・・・・・・・・・・・・・・・・・・
M=1の場合、回転した後に平行移動を行う。
Z軸の正の方向に50cm
平行移動させてみよう！
track_xz.eps
[transform]
8
課題1
transform.inp
Y軸の周りに45度回転させ、X軸
の正の方向に10cm、Z軸の正の
方向に45cm平行移動させよう！
[Transform]
set: c10[0] $ angle of around Z (degree)
set: c20[45]
c20[30] $ angle of around Y (degree)
set: c30[0] $ angle of around X (degree)
tr1 10
0 00 50
45
・・・・・・・・・・・・・・・・・・
track_xz.eps
[transform]
9
課題2
更に3つのセルを追加して、右下の図のような体系を作ってみよう！
• [surface]セクションは変更しない。
• 変更するのは[cell]セクションと[transform]セクションのみ。
（ [transform]セクションでは、set: c10[*]の行からコピーして、座
標変換番号2,3,4の変換を設定する。）
[transform]
10
実習内容
1. [transform]
2. [magnetic field]
3. [counter]
Table of Contents
11
[magnetic field]セクション
磁場(2重極場及び4重極場)を任意の領域に設定して
荷電粒子の軌道に対する影響を調べる。
2重極磁場：
1方向に一様な磁場
F=ev×B
[magnetic field]
12
[magnetic field]セクション
磁場(2重極場及び4重極場)を任意の領域に設定して
荷電粒子の軌道に対する影響を調べる。
4重極磁場：
広がった荷電粒子のビームを収束させる
Y軸
N極
S極
X軸
S極
N極
4重極磁石（参考）
[magnetic field]
13
入力形式
reg: 磁場を置く領域
typ: =2; 2重極磁場, =4; 4重極磁場
gap: 磁石の間の距離の半分[cm]
（2重極の時は使わないが値を入れる）
mgf: 磁場の強さ[kG]
trcl: 座標変換番号
[parameters]セクションにおいてimagnf=1
2重極磁場はy軸の正方向
→正の荷電粒子がz軸方向に進む時、
x軸の負の方向に動く
4重極磁場は
→正の荷電粒子がz軸方向に進む時、
x軸方向に収束、y軸方向に分散する
（z軸を中心軸とし無限に4重極磁石を配
置した時の磁場として定義されるため、
体系にあわせて座標変換が必要。）
[magnetic field]
14
例題
magfield.inp
体系を確認した後、icntl=0として輸送計算を実行。
[Source]
s-type = 1
proj = proton
dir = 1
r0 = 2.5
z0 = -10.
z1 = -10.
e0 = 290
[Surface]
10 so
500.
11 cz
5.
12 pz
0.
13 pz
10.
[Cell]
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
[magnetic field]
track_xz.eps
15
例題（2重極磁場）
magfield.inp
[Parameters]
icntl = 0
maxcas = 100
maxbch = 10
c imagnf
imagnf==11
file(6) = phits.out
[magnetic field]
reg typ gap mgf
104
2 10 100
track_xz.eps
[magnetic field]
16
例題（4重極磁場）
magfield.inp
[Parameters]
icntl = 0
maxcas = 100
maxbch = 10
imagnf = 1
file(6) = phits.out
[magnetic field]
reg typ gap mgf
104
2
4 10 100
track_xz.eps
[magnetic field]
17
課題1
4重極磁場の場合、x軸方向に
は収束するが、y軸方向には分
散する。これを確認してみよう！
magfield.inp
[T-TRACK]
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
1
30.0
40.0 50.0 60.0
-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
実習内容
1. [transform]
2. [magnetic field]
3. [counter]
Table of Contents
20
Counterとは？
PHITSでは、各粒子がエネルギーや位置情報など以外にカウンター
という値を持っています。カウンター値を特定のイベント（領域に入る、
領域から出る、その領域で散乱（反応）する、境界で反射する）が起
きる度に変化させることにより、タリーした粒子がどのような経緯で
そこに到達したかを調べることができます。
注目している領域からの寄与のみを
タリーしたい
→ 注目している領域で発生した粒子
のカウンター値を増やす
この領域で散乱した
場合にcounter値を+1
[counter]
タリー
21
入力形式
カウンター番号(3つまで)
カウンターをつける粒
子を指定
(*partはカウンターを
つけない粒子を指定)
reg: カウンターを回す領域
• in: 領域に入った時
• out: 領域から出た時
• coll: 領域で散乱（反応）した時
• ref: 領域の境界で反射した時
• 0: 変化なし
• 10000: ゼロセット
• それ以外: 指定した数だけ増減
[counter]
22
例題
counter.inp
体系を確認した後、icntl=0として輸送計算を実行。
[Material]
MAT[ 1 ]
1H 2 16O 1
MAT[ 2 ]
14N 8 16O 2
MAT[ 3 ]
63Cu 0.6915 65Cu 0.3085
[Cell]
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
・・・・・・
[Counter]
counter = 1
part = proton
reg
in
out
102
0
0
coll
1
ref
0
track_xz.eps
[counter]
23
例題
counter.inp
[T-Cross]
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
領域110から
unit = 1
101に入る粒
axis = eng
子をタリーする
file = cross.dat
output = current
part = proton
angel = ymin[1.E-06] ymax[1.E-02]
epsout = 1
cross.eps
[counter]
24
例題
タリーする粒子にカウンター値に
関する条件を加えてみよう！
•ctmin(1)=1, ctmax(1)=1
counter.inp
[Counter]
counter = 1
part = proton
reg in out coll ref
102 0
0
1
0
[T-Cross]
・・・・・・
・・・・・・
・・・・・・
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]
counter = 1
part = proton
reg in out coll ref
102 0
0
1
0
タリーする粒子に中性子(neutron)を加え、
中性子の振る舞いを調べてみよう！
[T-Cross]
・・・・・・
・・・・・・
・・・・・・
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
領域102の場合と同様に、103と104のそ
れぞれで陽子が反応したイベントをタリー
してみよう！
•領域103と104で起こるイベントに対応す
る様に、カウンター2と3を追加する
•[t-cross]を2回コピー&ペーストし、それぞ
れカウンター2と3がタリー条件となるよう
に設定する（出力file名も変更）
タリー
[counter]
27
課題2の答え
counter.inp
[Counter]
・・・・・・
counter = 2
part = proton
reg
in
out coll
ref
103
0
0
1
0
counter = 3
part = proton
reg
in
out coll
ref
104
0
0
1
0
[T-Cross]
・・・・・・
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