5.6. [ Cell ] セクション

このセクションでは、 [surface] セクションで定義した面を使ってセル(小部屋)を定義します。 基本的にはセルは閉じた空間として定義し、様々な形状のセルを組み合わせることにより粒子を輸送する仮想空間を構築します。 PHITSでは、粒子が飛ぶことが可能な内部空間とそれ以外の空間を分けて定義する必要があり、後者を外部ボイド(真空)として明示的に指定しなければなりません。

コメント文字として $ の他に、行頭の c (直後に半角スペースがある場合)も使用できます。 このセクションではセル定義に # を用いるため、PHITS入力ファイルでのコメント文字 # は利用できません。 ファイルのインクルード文や定数定義などはセクションの途中でも用いることができ、継続行には先頭に5桁の空白が必要です。

書式は以下の通りです。

[ Cell ]
  セル番号    物質番号    物質密度    セル定義文    セルパラメータ

また、LIKE n BUT セルパラメータ の表式や lattice 構造も利用できます。 似たような構造のセルを複数並べる場合に効率良く設定できます。利用方法の詳細は 5.6.4 章 をご覧ください。 セルパラメータの種類については 表 5.6.2 にまとめています。

表 5.6.1 セル記述書式

項目

説明

セル番号

1 から 999999 まで使用できます。

物質番号

[material] セクションで指定した物質番号を入力します。内部ボイド(真空)のときは 0、外部ボイドのときは -1 を指定します。

物質密度

そのセルに含まれる物質の密度を与えます。正の場合は粒子密度(\(10^{24}\) atoms/cm3 )、負の場合は質量密度(g/cm3 )です。ここで指定した粒子密度が計算に使用されますが、物質の組成比は [material] セクションで指定した密度を基に与えられます。従って、同じ組成で密度の異なる物質を定義可能となります。その際に、異なる物質番号を付加するパラメータ matadd を加えました。セルがボイドのとき、すなわち物質番号が 0 か -1 の場合はこの項目は入力しません。

セル定義文

[surface] セクションで定義された面番号や集合代数の演算子を用いてセルの幾何形状を記述します。演算子は空白(積 and)、: (和 or)、# (否定 not)で、() も使用できます。記述方法は 5.6.1 章 を参照してください。

LIKE n BUT

セル番号 n とほとんど同じで、BUT 以下のセルパラメータだけが違うセルを定義します。

セルパラメータ

keyword=value 形式で指定します。VOLTMPTRCLULATFILL が使用できます。LIKE n BUT の場合はこれらに加えて MATRHO も使えます。

表 5.6.2 セルパラメータ

項目

説明

VOL

そのセルの体積(cm3 )を与えます。

TMP

そのセルの物質の温度(MeV)を与えます。

TRCL

そのセルに座標変換を適用します。書式は [transform] セクションと同じです。

U

Universe番号。1 から 999999 まで使用できます。Universe構造を利用する場合に定義します。利用方法は 5.6.2 章 を参照してください。

LAT

Lattice番号。LAT=1 の場合は四角柱、LAT=2 の場合は六角柱を定義します。これらを単位とした繰り返し構造を利用する場合に指定します。利用方法は 5.6.3 章 を参照してください。

FILL

U で定義した universe 番号を入力し、そのセルをその universe で満たします。利用方法は 5.6.2 章5.6.3 章 を参照してください。

MAT

LIKE n BUT MAT=m の表式で用います。物質番号だけ違うセルを定義できます。

RHO

LIKE n BUT RHO=x の表式で用います。密度だけ違うセルを定義できます。

同じ物質番号を指定して別の密度をもつセルを複数定義した場合は注意してください。 2番目以降のセル番号の物質番号として別の数字が割り振られます。 また、セル定義文において、積を意味する空白の方が、和を意味する : より先に演算されます。

5.6.1. セルの記述方法

任意の形状の構造物を仮想空間に構築する場合、最初にその構造物の表面を座標空間において定義し、次に対応する領域をその面で閉じることによってセルを作ります。 この方法は面の表側と裏側を区別して表現し、さらに集合代数の演算子を利用することで各セルを定義します。

考えている領域が定義した面の表側と裏側のどちらに属するかを区別するために、面の方程式 \(f(x,y,z)=0\) を利用します。 ある領域中の座標 \((x_0,y_0,z_0)\) を代入した \(f(x_0,y_0,z_0)\) が正であればその領域はプラス側、負であればマイナス側です。

リスト 5.6.1 例題(1)
   1:   [ Cell ]
   2:     1   0  -10
   3:     2  -1   10
   4:   [ Surface ]
   5:     10  SZ  3  5

この例では、面番号10で中心 \((0,0,3)\)、半径5 cm の球面を定義しています。 5行目の SZ を用いて定義される面番号10の面は中心を \((0,0,3)\)、半径を 5 cm とする球面であるため、方程式は \(f(x,y,z)=x^2+y^2+(z-3)^2-5^2=0\) となります。 球の内はマイナス側、外はプラス側と定義されますが、これらは球の中心の座標 \((0,0,3)\) や球の外の適当な座標 \((0,0,10)\) を代入して、 \(f(0,0,3)=-25<0\)\(f(0,0,10)=24>0\) となることから確認できます。 [cell] セクションではこのような領域の区切り方を行った上で、必要に応じて各領域にセル番号を付与します。 例題(1)の2行目では、セル番号1のセルがボイド(真空)であることと、その領域が球面のマイナス側であることが定義されており、 図 5.6.1 で示す仮想空間を構築します。 図 5.6.1 は、例題(1)でつくられる領域を \(xz\) 平面で切った断面図として見たものです。本来は球ですが、2次元平面で見た結果円となっています。 例題(1)の3行目では外部ボイドが明示的に定義されており、そのセル番号を2と指定しています。

cell example 1

図 5.6.1 例題(1)の空間を \(xz\) 平面で切り取った断面図。セル番号1の領域はボイド(真空)。

リスト 5.6.2 例題(2)
   1:   [ Cell ]
   2:     1   0  11  -12  13  -14  15  -16
   3:     2  -1  #1
   4:   [ Surface ]
   5:     11  PX  -6
   6:     12  PX   6
   7:     13  PY  -6
   8:     14  PY   6
   9:     15  PZ  -6
  10:     16  PZ   6

例題(1)の球面のように1つの面だけで閉じた空間が定義できる場合は簡単ですが、[surface] セクションで定義できる面の多くは半無限領域をつくるため、幾つかの領域を組み合わせる必要があります。 これを行うにあたって、PHITSでは集合代数の考え方を導入しており、演算子は空白、 :# の3つで、それぞれ積(and)、和(or)、否定(not)を表します。 また、 ( ) で囲んだ範囲は1つの領域とみなします。ただし、空白と : は面番号同士を演算するのに対して、 # はセル番号に掛かります。 また、 #all とした場合は、outer void および他の universe に定義されたものを除く全てのセルを除外するので、シミュレーション空間全体を定義するのに便利です。

この例題では半無限領域をつくる PXPYPZ の面を使用して、一辺が12 cmの立方体をつくっています。 例えば面番号11,12の面は、それぞれ座標 \((-6,0,0)\), \((6,0,0)\) を通り \(yz\) 平面に平行な面となります。 これらの方程式はそれぞれ \(x+6=0\), \(x-6=0\) ですから、原点 \((0,0,0)\) を含む領域は、面番号11に対してはプラス側、面番号12に対してはマイナス側となります。 よって、面番号11と12で囲まれた領域は「11のプラス側かつ12のマイナス側」となり、空白を用いて 11 -12 と記述します。 面番号13から16の面も同様ですから、 \(x,y,z\) の3方向に関して閉じた空間を表現するには、2行目のように6つの面を空白でつなげた書き方をします。 また、2行目ではセルに掛かる演算子 # を用いて外部ボイドを定義しています。これはセル番号1の領域以外の範囲を対象にするという意味になります。

cell example 2

図 5.6.2 例題(2)の空間を \(xz\) 平面で切り取った断面図。セル番号1の領域はボイド(真空)。

リスト 5.6.3 例題(3)
   1:   [ Cell ]
   2:     1   0  -10 : (11  -12  13  -14  15  -16)
   3:     2  -1  #1
   4:   [ Surface ]
   5:     10  SZ   3  5
   6:     11  PX  -6
   7:     12  PX   6
   8:     13  PY  -6
   9:     14  PY   6
  10:     15  PZ  -6
  11:     16  PZ   6

2行目の ( ) で囲んだ部分が例題(2)のセル番号1の領域に相当しており、これと例題(1)の球面の内側を合わせた領域がこの例題ではセル番号1として定義されています。 図 5.6.3 に示したのが例題(3)の領域を \(xz\) 平面で切った断面図です。本来は立方体から球の一部が飛び出た形状をもつ立体ですが、2次元平面では図のようになっています。 本例題のように、領域を足し合わせる場合に使用する演算子が : です。

cell example 3

図 5.6.3 例題(3)の空間を \(xz\) 平面で切り取った断面図。

リスト 5.6.4 例題(4)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:   [ Cell ]
   4:     1   0       -10
   5:     2   1  -1.0   10 (11  -12  13  -14  15  -16)
   6:     3  -1        #1 #2
   7:   [ Surface ]
   8:     10  SZ   3  5
   9:     11  PX  -6
  10:     12  PX   6
  11:     13  PY  -6
  12:     14  PY   6
  13:     15  PZ  -6
  14:     16  PZ   6

[surface] セクションは例題(3)と同じです。5行目で「立方体の内側でかつ球の外側」の領域を空白を用いて記述しており、これをセル番号2と定義しています。 また、本例題では [material] セクションで water を物質番号1として定義しており、セル番号2の領域には粒子密度 \(1.0\times10^{24}\) atoms/cm3 の water が満たされています。 セル番号1は4行目で定義しており、この球の中はボイドとしています。外部ボイドの定義は6行目で行っており、セル番号1と2の領域を除いた他の全てが外部となっています。

幾何学形状を定義する際、何もない領域(真空や空気)も他の領域と二重定義にならないように定義する必要があります。 そのような場合、 #all と定義することにより、定義された領域を全て除外するので便利です。 ただし、外部ボイドや他の universe に定義された領域は除外されません。また、1つのインプットファイルに #all は1つしか定義できません。

cell example 4

図 5.6.4 例題(4)の空間を \(xz\) 平面で切り取った断面図。

5.6.2. Universe構造

PHITSでは、セルパラメータ U を用いて複数の universe を定義できます。 各 universe に個別の仮想空間を構築することで、ある universe の一部として別の universe の内部構造を引用できます。 このように、多数の universe をそれぞれ定義し相互参照する仕組みを、ここでは universe 構造と呼びます。 ただし、実際に粒子輸送の舞台となる空間そのものは1つだけであり、この空間の一部を別の universe の内部構造で満たす(置き換える)場合に universe 構造を利用します。 単純に一部分だけを置き換える場合はあまり意味がありませんが、特に 5.6.4 章 の繰り返し幾何形状を利用するときに有効となる機能です。

universe structure 1

図 5.6.5 Universe構造の模式図 (1)。立方体を \(yz\) 平面で分割した2つの直方体。

universe structure 2

図 5.6.6 Universe構造の模式図 (2)。円柱内部に水を満たし、それ以外はボイドとした universe 1。

universe structure 3

図 5.6.7 Universe構造の模式図 (3)。鉄の円柱の周りを水で満たした universe 2。

リスト 5.6.5 例題(5)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:     mat[2]  Fe 1
   4:   [ Cell ]
   5:       1  0        11  -12  13  -14  15  -17  FILL=1
   6:       2  0        11  -12  13  -14  17  -16  FILL=2
   7:     101  1  -1.0  -10  13  -14  U=1
   8:     102  0       #101  U=1
   9:     201  2 -10.0  -10  13  -14  U=2
  10:     202  1  -1.0  #201  U=2
  11:       9 -1       #1 #2
  12:   [ Surface ]
  13:     10  CY   5
  14:     11  PX  -6
  15:     12  PX   6
  16:     13  PY  -6
  17:     14  PY   6
  18:     15  PZ  -6
  19:     16  PZ   6
  20:     17  PZ   0

7,8行目で universe 1 を、9,10行目で universe 2 を定義しています。これらは共に面番号10,13,14を使って定義した円柱を中心に置いた宇宙となっています。 Universe 1 では、セル番号101の円柱内部に物質番号1の水を入れ、その外部をボイドにしてセル番号102の領域と定義しています。 Universe 2 では、円柱内部をセル番号201として鉄柱にし、セル番号202の外部には水を満たしています。 7から10行目の最後にある U=1, 2 がそのセルがどの宇宙に属するかを指定するパラメータです。 図 5.6.5 のセル番号1,2の領域に、 図 5.6.6図 5.6.7 の universe をそれぞれ入れる設定になっています。5,6行目ではセル番号1,2を定義しており、最後の FILL=1, 2 によってそのセルをどの universe によって満たすかを決定します。この例題の結果を示したのが 図 5.6.8 です。 \(xz\) 平面の断面図ですが、セル番号1の領域が universe 1 に、セル番号2が universe 2 に置き換えられているのがわかります。

cell example 5 result

図 5.6.8 例題(5)の空間を \(xz\) 平面で切り取った断面図。セル番号1の領域に universe 1 のセル番号101と102が、セル番号2の領域に universe 2 のセル番号201と202が入っている。

shifted universe example

図 5.6.9 例題(5)において、セル番号1と2の範囲を \(x\) 方向に関してずらした場合。

Universe構造を利用する際の注意点として、未定義の領域を引用しないことと、座標系が全ての宇宙で一致していることが挙げられます。例えば前者については、例題(5)の8行目でセル番号102の領域をボイドとして定義する必要があり、これがないと5行目で引用した際にセル番号102の部分を適切に設定できません。各 universe の全ての領域を定義する必要はありませんが、引用する領域については何らかの物質、あるいはボイドを指定しておかなければなりません。またその際、置き換える元の領域よりも引用する領域を少し広く定義しておかないと、エラーが出る場合がありますので注意してください。

次に後者については、座標原点の位置、 \(x,y,z\) 軸の方向、空間のスケールがどの universe でも同じです。これは同じ [surface] セクションで定義した面を利用しているためです。別の宇宙を引用する場合は、対応する座標を確認してください。例えば、例題(5)の14,15行目で PX の値を変えると、 図 5.6.9 のように円柱の一部が立方体の中に入らなくなります。

5.6.3. Lattice構造

本節では、lattice 構造の仕組みと定義の仕方、および基本的な利用方法について説明します。この機能は似た形状を多数繰り返す場合に非常に有効です。より実際的な設定方法については 5.6.4 章 をご覧ください。

lattice 構造は、四角形と六角形の2種類の形状を基本単位として、格子状に繰り返し並べる機能です。 U パラメータで指定して、基本格子を無限に並べた宇宙をつくります。 ただし、その内部は別の U パラメータを用いて定義した宇宙で満たす必要があります。 LAT=1 が四角柱、 LAT=2 が六角柱です。 図 5.6.10 に示した各辺の番号は、各図形を定義する際の面の順番を意味しており、 [surface] セクションで定義した面番号を図の数字の順に並べます。1と2、3と4、5と6の面はそれぞれが平行な面として定義されている必要があります。もちろん、プラス側・マイナス側のどちらであるかも指定します。ただし、1と2の長さと3と4の長さを変えることは可能で、正方形や正六角形でなくても構いません。

lattice numbering

図 5.6.10 Lattice構造の基本単位として使用できる四角形(左)と六角形(右)。各辺の番号は、 LAT パラメータを定義する際の面の順番。

リスト 5.6.6 例題(6)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:   [ Cell ]
   4:       1   0        11  -12   13  -14   15  -16  FILL=1
   5:     101   0       -26   25  -24   23  -22   21  LAT=1  U=1  FILL=2
   6:     201   1  -1.0  -90  U=2
   7:       2  -1       #1
   8:   [ Surface ]
   9:     11  PX  -6
  10:     12  PX   6
  11:     13  PY  -6
  12:     14  PY   6
  13:     15  PZ  -6
  14:     16  PZ   6
  15:     21  PX  -2
  16:     22  PX   2
  17:     23  PY  -2
  18:     24  PY   2
  19:     25  PZ  -2
  20:     26  PZ   2
  21:     90  BOX  -10 -10 -10  20 0 0  0 20 0  0 0 20

5行目で lattice の基本格子を定義しています。 LAT=1 で扱うことができる形状は四角形ですが、3次元の形状を考えると実際には四角柱となります。 本例題で定義しているのは断面が一辺4 cmの正方形となる無限の高さをもつ四角柱であり、これが隙間なく並んでいる宇宙を universe 1 としています。 6行目では、面記号 BOX を用いて一辺20 cmの立方体が原点にある宇宙を定義しており、これは universe 2 としています。このセル番号201の立方体には水が入っており、 FILL=2 を用いることにより、セル番号101中の格子一つ一つが水で満たされることになります。 4行目でセル番号1を universe 1 で満たした領域として定義しており、これが一辺12 cmの立方体であることから、本例題の仮想空間は 図 5.6.12 のようになります。 これらは universe 1 においては高さが無限大ですが、セル番号1の領域が \(y\) 軸方向にも制限されているため、結果的に有限の高さをもつ四角柱となっています。もし、有限の高さをもった四角柱を lattice の基本単位としたい場合は、例えば面番号23,24による領域の制限を5行目で加える必要があります。

cell example 6 3d

図 5.6.11 例題(6)の図 (1)。立体図。

cell example 6 section

図 5.6.12 例題(6)の図 (2)。 \(xz\) 平面で切り取った断面図。各格子は、中に示された座標によって区別される。

図 5.6.12 の右図で示したように、各円柱には lattice 座標 \((i,j,k)\) がふられており、各々の領域がセル番号101の中のどれであるかを区別できるようになっています。ただし、この座標は必ずしも通常の座標 \((x,y,z)\) とは対応しておらず、単位格子を定義する際の面の順番で決まります。すなわち、 図 5.6.10 の左図において、辺2 \(\to\) 1の方向が \(i\) 、4 \(\to\) 3 が \(j\) の方向となります。残りの \(k\) の向きは、 \(i,j\) と同様にセル定義文の順番で決定されます。例題(6)の5行目では、 \(z\) 軸方向、 \(x\) 軸方向の順番で定義しているため、 \(i\)\(z\) 軸、 \(j\)\(x\) 軸に対応付けがなされています。実際に 図 5.6.12 の右図でも、 \(z\) が増える方向に \(i\) が増加し、 \(x\) が増える方向に \(j\) の値が増えていることがわかります。各タリーにおいて mesh=reg により領域指定する場合は、 (201 < 101[-1 0 0]) の書式を用います。これは「セル番号101の lattice 座標 \((-1,0,0)\) に含まれているセル番号201の領域」という意味です。この書式に関しては 06章「タリー形式」の領域メッシュの節を参照してください。なお lattice 座標 \((i,j,k)\) は、 [t-gshow] タリーにおいて output=7 or 8 を指定することで表示させることができます。

リスト 5.6.7 例題(7)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:   [ Cell ]
   4:       1   0        11 -12  13 -14  15 -16  FILL=1
   5:     101   0       -31  32 -33  34 -35  36 -24  23  LAT=2  U=1  FILL=2
   6:     201   1  -1.0  -90  U=2
   7:       2  -1       #1
   8:   [ Surface ]
   9:     11  PX  -6
  10:     12  PX   6
  11:     13  PY  -6
  12:     14  PY   6
  13:     15  PZ  -6
  14:     16  PZ   6
  15:     23  PY  -2
  16:     24  PY   2
  17:     set: c1[2]
  18:     31  PZ  [ c1*cos(pi/6)]
  19:     32  PZ  [-c1*cos(pi/6)]
  20:     33  P   1  0  [ 1/tan(pi/3)]  [ c1]
  21:     34  P   1  0  [ 1/tan(pi/3)]  [-c1]
  22:     35  P   1  0  [-1/tan(pi/3)]  [ c1]
  23:     36  P   1  0  [-1/tan(pi/3)]  [-c1]
  24:     90  BOX  -10 -10 -10  20 0 0  0 20 0  0 0 20

17から23行目で正六角形をつくるための面の定義を行っています。これらを組み合わせると、 \(xz\) 平面で見た場合に正六角形となり、 \(y\) 軸に垂直な六角柱となります。 ここで、正六角形の一辺の長さを2 cmとし、これを c1 として定数定義しています。また、組み込み定数である pi および面記号 PZ, P を用いて、60度ずつ回転した平面を定義しています。Lattice構造の基本格子は5行目で定義しており、セル定義文中の面番号の順番は 図 5.6.10 の右図にある数字の通りとなっています。ただし、最後に面番号23,24による \(y\) 軸方向の制限を加えています。 したがって本例題では、高さが4 cm、一辺が2 cmの六角柱を無限に並べた宇宙が universe 1 となります。各六角柱には6行目で定義した universe 2 が入っており、実質的には水で満たされたセル番号201が入っていることになります。4行目で、先の例題と同様に、セル番号1の領域を universe 1 で満たしており、本例題の結果は 図 5.6.14 となります。左図が立体図、右図がこれを \(xz\) 平面で切った断面図です。左図からわかるように、六角柱が詰められた universe 1 から立方体の領域だけを切り取った形状になっています。また、 \(y\) 軸方向にも3列が並んでおりますが、その高さは揃えられており、いわゆる3次元のハニカム構造はつくりません。右図では各領域に lattice 座標 \((i,j,k)\) が付けられています。図を見てわかるように、 \(i\) は左から右、 \(j\) は左下から右上に向かって増加しており、これらはやはり5行目のセル定義文のところに並べる順番に依存します。 図 5.6.10 の右図において、辺2 \(\to\) 1の方向が \(i\) 、4 \(\to\) 3 が \(j\) の方向となります。 \(k\) についても同様で、面番号24 \(\to\) 23の順番であれば \(y\) 軸方向に、23 \(\to\) 24であれば \(-y\) 軸方向に \(k\) の値は増加します。各タリーにおいて mesh=reg により領域指定する場合は、 (201 < 101[-2 0 0]) の書式を用います。これは「セル番号101の lattice 座標 \((-2,0,0)\) に含まれているセル番号201の領域」という意味です。この書式の詳細は 06章「タリー形式」の領域メッシュの節を参照してください。なお lattice 座標 \((i,j,k)\) は、 [t-gshow] タリーにおいて output=7 or 8 を指定することで表示させることができます。

cell example 7 3d

図 5.6.13 例題(7)の図 (1)。立体図。

cell example 7 section

図 5.6.14 例題(7)の図 (2)。 \(xz\) 平面で切り取った断面図。

5.6.4. 繰り返し幾何形状

多数の同じ、あるいは少し違うだけの構造のセルを並べる場合に、当然1つ1つのセルを定義することは可能ですが、基本となるセルを1つ用意しこれを繰り返すことで効率良く設定できます。 5.6.3 章 のlattice構造の利用もその1つであり、他には LIKE n BUT セルパラメータ の表式も利用できます。

5.6.4.1. LIKE n BUT セルパラメータ

セル番号 n のセルと同様であるが、BUT 以下の内容だけが違うセルを定義する方法です。 BUT の後に続くセルパラメータには 表 5.6.2 に挙げたものがあります。以下の例題では、 TRCLMAT を使用した場合を考えています。

リスト 5.6.8 例題(8)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:     mat[2]  Fe 1
   4:   [ Cell ]
   5:       1   0       -10  13 -14  #2 #3 #4
   6:       2   1  -1.0   11 -12  13 -14  15 -16
   7:       3   LIKE 2 BUT  TRCL=1
   8:       4   LIKE 2 BUT  TRCL=2  MAT=2
   9:       5  -1        #(-10  13  -14)
  10:   [ Surface ]
  11:     10  CY   10
  12:     11  PX  -2
  13:     12  PX   2
  14:     13  PY  -2
  15:     14  PY   2
  16:     15  PZ  -2
  17:     16  PZ   2
  18:   [ Transform ]
  19:   *tr1   3 0 -5
  20:   *tr2  0 0 6  0 30 0  0 0 0  0 0 0  2

7,8行目で LIKE n BUT セルパラメータ の表式を用いてセルを定義しています。この元となるセル番号2は6行目で定義しており、中に水が入った一辺4 cmの立方体です。 セル番号2は原点を中心として配置されていますが、これを19行目で定義した座標変換番号1で座標変換したのがセル番号3です。 次のセル番号4は TRCLMAT のセルパラメータを使用しており、座標変換番号2による平行移動および回転を行い、その上で物質番号2に中の物質を変更しています。 5行目では、セル番号1を半径10 cm、高さ4 cmの円柱と定義しています。ただし、セル番号2,3,4の領域が含まれないように、 # を用いてこれらを除いています。 また、これに伴い外部ボイドの定義も9行目に示すような書き方が必要です。単純に #1 のみではセル番号2,3,4の領域も外部ボイドとなってしまいます。

cell example 8

図 5.6.15 例題(8)の空間を \(xz\) 平面で切り取った断面図。

5.6.4.2. Latticeを用いた階層構造

5.6.3 章 で見た四角柱や六角柱の lattice をより複雑なものにするために、universe を用いて階層構造を構築することができます。考え方は universe 構造に基づいており、universe 1 の一部を universe 2 で満たし、更に universe 2 の一部を universe 3 で満たす、といった手続きで各階層を定義していくことが可能です。PHITSでは基本的に最大10( param.inc 中の mxlv で変更できます)の階層構造をつくることができます。

リスト 5.6.9 例題(9)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:     mat[2]  Fe 1
   4:   [ Cell ]
   5:       1   0        11  -12   13  -14  15  -16  FILL=1
   6:     101   0       -26   25  -22   21  LAT=1  U=1
   7:                   FILL=-1:1  -1:1  0:0
   8:                   2 2 3  2 3 2  3 2 2
   9:     201   1  -1.0  -90  U=2
  10:     301   2 -10.0  -10  U=3
  11:     302   0        10  U=3
  12:       2  -1       #1
  13:   [ Surface ]
  14:     10  CY   1.5
  15:     11  PX  -6
  16:     12  PX   6
  17:     13  PY  -6
  18:     14  PY   6
  19:     15  PZ  -6
  20:     16  PZ   6
  21:     21  PX  -2
  22:     22  PX   2
  23:     25  PZ  -2
  24:     26  PZ   2
  25:     90  BOX  -10 -10 -10  20 0 0  0 20 0  0 0 20

まず5,6行目にある四角柱の lattice 構造とそれを引用するセル番号1の定義は例題(6)と同じです。ただし、7,8行目におけるセルパラメータ FILL の使い方が違っており、ここでは各格子の状態を個別に指定しています。7行目では取り扱う範囲を決めており、各格子に割り振られた lattice 座標 \((i,j,k)\) を参照して FILL= の後ろで決定しています。本例題では、 \(i\) に関して -1 から1まで、 \(j\) に関して -1 から1まで、 \(k\) に関して0から0までをセル番号101の lattice 構造と定義しています。合計 \(3\times3\times1=9\) の格子が対象となり、それぞれを満たす universe は8行目で指定しています。この順番は lattice 座標 \((i,j,k)\) で表現して、 \((-1,-1,0), (0,-1,0), (1,-1,0), (-1,0,0), \cdots\) のようになっており、2であれば universe 2 で、3であれば universe 3 でその格子を満たすことを意味します。Universe 2 を水が入った一辺20 cmの立方体を原点付近に配置した宇宙として、universe 3 を原点に半径1.5 cmの鉄柱を置いた宇宙として定義しており、9つの四角柱のうち3つが特別な構造をもつことになります。結果を示したのが 図 5.6.16 で、 \(xz\) 平面に関する断面図です。この図から \((i,j,k)=(1,-1,0), (0,0,0), (-1,1,0)\) の格子が内部構造をもっていることがわかります。各タリーにおいて mesh=reg により領域指定する場合は、 (302 < 101[0 0 0] < 1) の書式を用います。これは「セル番号1を満たしたセル番号101の lattice 座標(0,0,0)に含まれているセル番号302の領域」という意味で、 図 5.6.16 では中心の四角柱にある鉄柱周辺のボイドの領域を指します。このような書式については 06章「タリー形式」の領域メッシュの節を参照してください。

cell example 9

図 5.6.16 例題(9)の空間を \(xz\) 平面で切り取った断面図。

リスト 5.6.10 例題(10)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:     mat[2]  Fe 1
   4:   [ Cell ]
   5:       1   0        11  -12   13  -14  15  -16  FILL=1
   6:     101   0       -26   25  -22   21  LAT=1  U=1
   7:                   FILL=-1:1  -1:1  0:0
   8:                   2 2 3(1 0 1)  2 3(1 0 1) 2  3(1 0 1) 2 2
   9:     201   1  -1.0  -90  U=2
  10:     301   0       -36   35  -32   31  LAT=1  U=3
  11:                   FILL=-1:0  -1:0  0:0
  12:                   4 2  2 4
  13:     401   2 -10.0  -10  U=4
  14:     402   0        10  U=4
  15:       2  -1       #1
  16:   [ Surface ]
  17:     10  CY   0.5
  18:     11  PX  -6
  19:     12  PX   6
  20:     13  PY  -6
  21:     14  PY   6
  22:     15  PZ  -6
  23:     16  PZ   6
  24:     21  PX  -2
  25:     22  PX   2
  26:     25  PZ  -2
  27:     26  PZ   2
  28:     31  PX  -1
  29:     32  PX   1
  30:     35  PZ  -1
  31:     36  PZ   1
  32:     90  BOX  -10 -10 -10  20 0 0  0 20 0  0 0 20

この例題でつくられる仮想空間は、 図 5.6.17 のようになります。ただし、 \(xz\) 平面で切った断面を示しています。まず、一辺12 cmの立方体を9つの四角柱に分け、そのうちの3つに内部構造をもたせた点は例題(9)と同じです。違うのは universe 3 として定義した宇宙を更に内部構造をもった格子で埋めたことで、9から11行目でその状態を指定しています。ここでは \(2\times2\times1=4\) の四角柱を定義しており、そのうちの2つが鉄柱をもつ格子となっています。また、8行目の (1 0 1) は座標変換を表しており、引用する際に universe 3 を \(x\) 軸方向に1 cm、 \(z\) 軸方向に1 cm平行移動しなさい、という意味となります。セル番号301は、原点 \((0,0,0)\) を中心にもつ四角柱を基本単位としているため、図のような配置にする場合は少しずらす必要があります。また、lattice を用いて階層構造を作成した場合の mesh=reg による領域指定は、階層構造を増やした書式により行います。例えば 図 5.6.17 において \(x=-1, z=-1\) を中心とする微小四角柱のボイド部分は、 (402 < 301[-1 -1 0] < 101[0 0 0] < 1) と記述します。もし対象とする階層が更に深くなった場合は、 < と lattice 座標を用いて階層を増やすことで各領域を指定できます。この書式の詳細については 06章「タリー形式」の領域メッシュの節を参照してください。

cell example 10

図 5.6.17 例題(10)の空間を \(xz\) 平面で切り取った断面図。

5.6.4.3. ボクセルファントム(voxel phantom)の利用

PHITSでは lattice 構造を利用してボクセルファントムを仮想空間として設定できます。ボクセルファントムとは微小立方体を積み重ねて生物などの複雑な構造物を表現したものです。CT等の画像データを元にして作成するボクセルデータを使用します。本節では、非常に簡単なボクセルデータを用意し、これを用いた設定方法について説明します。

基本的な考え方として、まず粒子輸送を行うためのある一定の大きさをもつ立方体、あるいは直方体を用意し、その中に微小立方体、すなわちボクセルを並べます。その際、lattice 構造の LAT=1 を利用し、各ボクセルの座標が指定できるようにセルパラメータ FILL を使います。 FILL を使用した場合、各ボクセルがどの状態か、すなわちどの universe で満たされているかが指定できるため、水や鉄、あるいはボイドで満たされた universe を引用することで、様々な構造物をつくることができます。したがって、ボクセルデータは各ボクセルの位置座標と構成物質の情報をもっている必要があり、ここでは \(L_{000}~L_{100}~L_{200}~\cdots~L_{stu}~\cdots\) の順番で並べたデータ群を想定しています。ここで、 \(L_{stu}\) は lattice 座標 \((i,j,k)\) にあるボクセルの構成物質、すなわち universe 番号です。

リスト 5.6.11 例題(11)
   1:   [ Material ]
   2:     mat[1]  1H 2  16O 1
   3:     mat[2]  Fe 1
   4:   [ Cell ]
   5:       1   0   11  -12   13  -14  15  -16  FILL=1
   6:     101   0  -20  LAT=1  U=1
   7:              FILL=-2:2  -2:2 -2:2
   8:              2 2 2 2 2  2 2 2 2 2  2 2 3 2 2  2 2 2 2 2  2 2 2 2 2
   9:              2 2 2 2 2  2 3 3 2 2  2 3 4 3 2  2 3 3 2 2  2 2 2 2 2
  10:              2 2 2 2 2  2 3 3 3 2  3 4 4 4 3  2 3 3 3 2  2 2 2 2 2
  11:              2 2 2 2 2  2 2 3 3 2  2 3 4 3 2  2 2 3 3 2  2 2 2 2 2
  12:              2 2 2 2 2  2 2 2 2 2  2 2 3 2 2  2 2 2 2 2  2 2 2 2 2
  13:     201   0       -90  U=2
  14:     301   2 -10.0  -90  U=3
  15:     401   1  -1.0  -90  U=4
  16:       2  -1       #1
  17:   [ Surface ]
  18:     11  PX  -5
  19:     12  PX   5
  20:     13  PY  -5
  21:     14  PY   5
  22:     15  PZ  -5
  23:     16  PZ   5
  24:     20  BOX  -1 -1 -1  2 0 0  0 2 0  0 0 2
  25:     90  BOX  -10 -10 -10  20 0 0  0 20 0  0 0 20

5行目でセル番号1の領域を一辺が10 cmの立方体として定義しています。これが次の行のセル番号101の lattice で満たされており、基本単位は面番号20で定義した一辺2 cmの立方体です。ここでは面記号 BOX を用いて微小立方体、すなわちボクセルを定義しています。7行目でボクセルを lattice 座標 \((i,j,k)\) のそれぞれに5つずつ並べるとしており、8から12行目で各ボクセルを満たす universe 番号を指定しています。ボクセルは \((-2,-2,-2), (-1,-2,-2), \cdots, (2,2,2)\) の順番で並んでおり、ボイドの場合は2、鉄は3、水は4の数字が対応しています。ボイドや鉄、水が入った空間の定義は13から15行目で行っており、それぞれが一辺20 cmの立方体に入った宇宙を universe 2,3,4 としています。本例題のボクセルで表現しているのは多少いびつな形状の鉄の箱に水が入っている構造物です。結果を示したのが 図 5.6.19 です。左図が立体的にそのまま見た結果、右図が表面の一部を透明にして見た結果です。少しずれていますが、図の左から右が \(i\) 軸、下から上が \(j\) 軸、奥から手前が \(k\) 軸の方向に対応します。見てわかるように、でこぼこした鉄の箱の中には水が入っており、この構造物を対象に粒子輸送計算を行うことができます。各タリーにおいて mesh=reg により領域指定する場合は、 (401 < 101[0 0 0] < 1) の書式を用います。あるいは、セル番号1中の水の部分全てを指定したい場合は、 (401 < 1) とすることができます。これらの書式については 06章「タリー形式」の領域メッシュの節を参照してください。ただし、 (301 < 101[-2:2 -2:2 -2:2] < 1) のような指定方法はできません。セル番号 101[-2:2 -2:2 -2:2] の全てにセル番号301の領域が含まれているわけではないためです。

voxel phantom 3d

図 5.6.18 例題(11)の図 (1)。立体図。

voxel phantom inside

図 5.6.19 例題(11)の図 (2)。内部構造を見た図。

例題(11)の8から12行目の universe の羅列について、バージョン3.09以降では短縮形式での記述が可能になりました。ボクセルファントムのようなデータでは、同じ universe 番号が連続することが多いため、この短縮形式を利用することで、大幅に記述を省略することができます。これにより、膨大なボクセルデータに対するファイルサイズの削減や読み書きに掛かる計算時間の短縮につながります。短縮形式では、連続する universe 番号について、一つ目の universe 番号はこれまでどおり記述し、二つ目には universe 番号ではなく、これ以降に同じ universe 番号が連続する数を負の符号を付けて記述します。

例題11d

 8:             2 -11 3 2 -17 3 -1 2 -2 3
 9:             4 3 2 -1 3 -1 2 -12 3 -2
10:             2 3 4 -2 3 2 3 -2 2 -12 3
11:             -1 2 -1 3 4 3 2 -2 3 -1 2
12:             -17 3 2 -11

例えば例題(11)では、universe番号の2が12個並んだ後に、universe番号の3が記述されているのに対し、例題(11')では、universe番号の2の後に同じ universe 番号がさらに11個並んでいるとの意味で、2 -11 のように記述できます。

PHITSでは、最初に読み取ったインプットデータを一時的にバイナリ形式で書き出し、それを再度読み込んで粒子輸送計算を行います。例題(11)で扱ったのは簡単なボクセルですが、一般的にボクセルデータは大規模なものとなり、それをPHITSの実行の度に読み書きするのは非常に時間がかかります。これを緩和するために、 [parameters] セクションの ivoxel パラメータを使う方法があります。 ivoxel は、ボクセルデータ、実質的には例題(11)の8から12行目をバイナリとして出力させて保存し、2度目からはそのバイナリを読んで計算するという機能をもちます。PHITS実行時に ivoxel=2 とすると file(18) にデータを書き出してそのまま終了します。ただし、バージョン2.30以前では終了せずに引き続き輸送計算を行います。以降 ivoxel=1 と設定することで該当ファイルからバイナリデータを読み込んで計算できます。また、大量のボクセルデータをPHITSのインプットファイルとは別のファイルに保存して計算するときは、 infl を使用します。この場合の ivoxel の使い方も同じです。

5.6.4.4. 連続四面体形状の利用

本節では、連続四面体形状の仕組みと定義の仕方、および基本的な利用方法について説明します。連続四面体形状はポリゴンの一種で、様々な大きさの四面体を組み合わせることで、複雑な形状をもつ物体を表現できるようになります。なお、四面体メッシュ体系内では反射面が動作しません。四面体メッシュ体系内で反射面を使いたい場合は、物理的に無視できる薄い四面体メッシュ体系外領域を内部に定義し、その中に反射面を定義すれば、反射を実現できます。

PHITSでは二種類の連続四面体形状の形式に対応しています。一つは四面体メッシュ分割ソフトウェアTetGenによって生成されるNodeファイルとElementファイルで四面体形状を表現する形式で、もう一つは有限要素法CAEソフトウェアNASTRANのバルクデータ形式です。

TetGenでの形式は、

  • Nodeファイル: 各四面体の各々の頂点の xyz 座標

  • Elementファイル: 各四面体要素を構成する4つの頂点

の2つのファイルで表現されます。NodeファイルおよびElementファイルは、共通のファイル名に .node および .ele という拡張子が付いた形を想定しています。PHITSのインプットファイルでは、この共通ファイル名を指定することになります。この機能の具体的な利用方法については、/phits/utility/TetraGEOM フォルダにある資料やサンプルインプットを参照してください。

リスト 5.6.12 Nodeファイルの例
   1:   # Simple two elements
   2:   5   3  0  0
   3:   #  pointID x y z
   4:   1    0.0  -2.5  -4.0
   5:   2   -4.0  -2.5   0.0
   6:   3    4.0  -2.5   0.0
   7:   4    0.0  -2.5   4.0
   8:   5    0.0   2.5   0.0

Nodeファイルでは、まず始めにこのファイルに含まれる全ての頂点の数と扱っている次元数を設定します(2行目)。PHITSで扱える次元は3次元のみです。2行目の最後の2つの0はPHITSでは使われずに無視されます。

その下に(4行目以下)、

[頂点番号] [x] [y] [z]

の順で頂点番号とその頂点の xyz 座標を指定します。

Elementファイルでは、まず始めにこのファイルに含まれる全要素数と1つの要素を構成する頂点の数を設定します。PHITSでは1次要素の四面体しか扱えないので、4となります。また、要素に付加する情報の数はPHITSでは1です(2行目)。

その下に(4行目以下)、

[要素番号] [頂点番号1] [頂点番号2] [頂点番号3] [頂点番号4] [要素universe番号]

の順で各要素の情報を指定します。要素番号は各四面体のID番号、頂点番号1から4はその四面体を構成する頂点の番号でNodeファイル内で定義した番号です。最後の要素 universe 番号で、各四面体をどの universe 番号の空間で満たすかを指定します。各 universe の定義はPHITSのインプットファイルで行います(例題(12)の13,14行目参照)。なお、NodeファイルおよびElementファイル内では、任意の行で # から始まるコメント行を挿入することができます。

リスト 5.6.13 Elementファイルの例
   1:   # Simple two elements
   2:   2  4  1
   3:   # elementID point(1:4) universe
   4:   1     1  2  3  5    5001
   5:   2     5  2  3  4    5002

NodeファイルおよびElementファイルは四面体メッシュ分割ソフトウェアTetGenを使って作成することができます。多少技術が必要ですが、一般的なポリゴンデータもこのソフトウェアを使って連続四面体に変換することが可能です。TetGenは次のウェブサイトから無料でダウンロード可能です。

http://wias-berlin.de/software/tetgen/

TetGenの使い方については、TetGenのマニュアルを参照してください。

NASTRANのバルクデータ形式は、構造解析計算や流体解析計算など様々なソフトウェアで一般的に出力可能な形式です。PHITSではバルクデータ形式の中で、頂点情報 GRID と四面体ソリッド要素情報 CTETRA で構成されるBDFファイルを読み込むことができます。通常、何らかのソフトウェアにより体系を作成し、それをバルクデータ形式で出力することが想定されることから、バルクデータ形式に関する詳細は省略します。詳しくはNASTRANバルクデータで検索してください。構造計算や流体計算と同一の形式を採用することで、これらの計算と放射線輸送計算をシームレスに繋いだ解析が可能となりました。具体的な利用方法については、/phits/utility/FLUENT フォルダにある資料を参照してください。

次に、PHITSのインプットファイルにおける指定の仕方について説明します。

リスト 5.6.14 例題(12)
   1:   [Material]
   2:   mat[1] 14N 78.1  16O 20.9  40Ar 0.93
   3:   mat[2] 1H 2  16O 1
   4:   mat[3] 56Fe 1
   5:   [Surface]
   6:   10  rpp   -5.0  5.0  -3.0  3.0  -5.0  5.0
   7:   20  rpp   -7.0  7.0  -5.0  5.0  -7.0  7.0
   8:   90  so     500.0
   9:   [Cell]
  10:   101   1 -0.001205 -20  U=1  LAT=3  tfile=Tetra  TSFAC=1.0
  11:     1   0           -10  FILL=1
  12:     2  -1            10
  13:   201   2 -1.0      -90  U=5000
  14:   202   3 -7.874    -90  U=5001

PHITSで連続四面体形状を取り扱うためには、直方体、すなわち面記号 RPP で指定する領域を定義し、その中に四面体連続形状を展開します。全ての四面体の頂点はこの直方体の中に含まれる必要があります。ただし、直方体の大きさを連続四面体形状の大きさからかけ離れて大きく取り過ぎると計算の都合で問題が生じるので、連続四面体形状が全て収まる適当な大きさの直方体としてください。 10行目では、セル番号101の領域として、面番号20で定義される14 cm \(\times\) 10 cm \(\times\) 14 cmの直方体の中にファイル名 Tetra で指定される連続四面体形状を展開しています。 LAT=3 は連続四面体形状を使用するオプションの指定で、 tfile に与えられるファイル名に対応するNodeファイルおよびElementファイルが読み込まれることになります。NASTRANバルクデータ形式のファイルを使用する場合は、例題(13)のように tfile の代わりに nfile でファイル名を指定してください。 なお、 nfile を指定する際、拡張子をもつ場合は、それまで含めたファイル名を指定してください。 [1] PHITS用HDF5ファイルフォーマットph5のファイルを使用する場合は、 リスト 5.6.16 の様にTFILEの代わりにHFILEでファイル名を拡張子も含めて指定してください。 Windows以外のOSではファイル名に大文字小文字の区別があるので、ファイル名の指定はどちらの形式でも大文字小文字を正しく入力してください。 TSFAC の指定により、Nodeファイルで与えられる四面体頂点の座標に対する規格化定数を与えることができます。ここでは TSFAC=1.0 なので、Nodeファイルで与えられる座標をそのまま使用します。 セル番号101の領域に対するマテリアルの指定は、領域101内で連続四面体を除く領域で使用されることになります。11行目でセル番号1の領域を10 cm \(\times\) 6 cm \(\times\) 10 cmの直方体として定義しています。 FILL の指定により、この領域が10行目のセル番号101で定義される連続四面体を含む領域によって満たされることになります。このように連続四面体形状、すなわち LAT=3 を使用する場合も、Lattice構造の指定と同じように、UniverseとFillを利用した入れ子構造を定義してください。Elementファイルで指定される各四面体の要素 universe に対応する領域の指定はLattice構造での指定と同様です。この例題では、四面体要素1には物質番号2の水が、四面体要素2には物質番号3の鉄が含まれることになります。

リスト 5.6.15 例題(13)
  10:   101   1 -0.001205 -20  U=1  LAT=3  nfile=Tetra.bdf  TSFAC=1.0

こちらはNASTRANバルクデータ形式、すなわちBDFを nfile で読み込む例です。

tetra mesh example

図 5.6.20 連続四面体形状の例。

リスト 5.6.16 例題(14)
  10: 101 1 -0.001205 -20 U=1 LAT=3 hfile=Tetra.ph5 TSFAC=1.0

連続四面体形状の構成物質が増えた場合は、各四面体の要素 universe に対応する領域の数が増えます。この手作業の負担を減らすため、自動的にこの領域を作成するオプションを追加しました。この自動領域追加オプションを使用するには、パラメータセクションで itetauto=1 と指定してください。連続四面体形状に含まれる全ての universe に対応する領域を自動で作成します。この際に、面番号5000番の非常に大きな球面とセル番号5001番から5000+n番を新たに追加します。ここで、nは連続四面体形状に含まれる universe の数です。ですので、このオプションを使う際にはこれらの面番号およびセル番号を使用しないようにしてください。NASTRANバルクデータ形式ファイルを使用する場合は、バルクデータに含まれる PSOLID タブおよび MAT タブの情報を読み取り、自動的に追加領域の密度が設定されます。TetGen形式ファイルを使用する場合は、これらの情報を外部ファイルで与える必要があります。同じファイル名で拡張子 .txt のファイルを用意し、例題(16)のように各 universe 番号と密度を1行ごとに記述してください。追加されたセル番号5001番から5000+n番の領域には自動的にセル番号と同じマテリアル番号が割り振られます。対応するマテリアルが存在しない状態で一度PHITSを実行すると、例題(17)のようなエラーメッセージが表示されますので、このメッセージを参考に連続四面体形状の universe に対応するマテリアルを5001番から5000+n番に定義してください。phits.out に含まれるインプットエコーに追加された面および領域が記述されるので、正しく指定が行えたかどうかを確認することができます。

リスト 5.6.17 txtファイルの例
   1:   5001 5001 -1.0
   2:   5002 5002 -7.874
リスト 5.6.18 エラーメッセージの例
*** ERROR : undefined material
 TETRA material (MID):  5001
   should be defined as material ID number:  5001
*** ERROR : undefined material
 TETRA material (MID):  5002
   should be defined as material ID number:  5002

エラーメッセージを参考に、連続四面体形状の universe に対応するマテリアルを追加定義してください。

LAT=3 を指定した cell には TRCL を指定することができません。連続四面体形状で読み込んだ体系に座標変換を加えるには、fill / universe を利用した階層構造を用意し、 LAT=3 を指定して連続四面体形状を読み込む直方体の定義よりも上の階層に TRCL を指定するようにしてください。具体的な利用方法については、/phits/utility/TetraGEOM フォルダにある資料やサンプルインプットを参照してください。 連続四面体形状を用いたより詳細な設定や変換の実例については、/phits/utility/TetraGEOM および /phits/utility/FLUENT の資料も参照してください。