入力
電磁場取り込み
進行波型加速管
User elemetns
Python
GPT進行波型加速管の計算方法SUPERFISH の電磁場を取り込み計算
General Particle Tracer(GPT)を使い進行波型加速管のビームトラッキング計算を行う方法を示します．GPT には，進行波型加速管を取り扱うためのエレメントは標準で用意されています．しかし，その電磁場のモデルには問題があります．特に，バンチャー加速管でのトラッキング計算ではエラーが大きいと思われます．そこで，ここでは進行波型加速管内の電磁場を SUPERFISH で計算し，それをGPTに適用する方法を示します．この方法を使うと，GPT の標準のエレメントを使う場合に比べ，計算精度が良いはずです．
目次
はじめに
進行波型の加速管(空胴)を計算するために，GPTには3つのビルトインエレメント(trwcell, trwlinac, trwlinbm)が用意されています．これらのエレメントでは，進行波型加速管の電磁場は次のようにして取り扱っています．
\begin{align} E_z&=E_0\sin(\omega t+\phi-k_zz)I_0(k_rr) \nonumber \\ E_r&=E_0\cos(\omega t+\phi-k_zz)I_1(k_rr) \\ B_\phi&= E_r/c \nonumber \end{align}ここで，$(E_z,E_r,B_\phi)$は電磁場．
基本的な考え方
GPTで定在波型の空胴の計算は容易です．コマンド sf7 と fish2gdf を使い，SUPERFISH の計算結果から GPT で読み込み可能な 電磁場を記述した gdf ファイルを作成します．これを直接，GPTの入力ファイルで指定します．それに対して，進行波型加速管の場合はかなりやっかいです．SUPERFISHでは進行波の電磁場が計算できないことと，GPTでは進行波の電磁場を直接取り込めないことが，取り扱いを難しくします．しかし，工夫をすれば，進行波の取り扱いは可能になります．ここで，筆者がいつも使っているテクニックを紹介します．特に電子銃から出射される低エネルギー(数十 — 数百 [keV])のビームを進行波型のバンチャーでバンチするシミュレーションに適していると考えています．
GPTで進行波型加速管内でのビームトラッキングを行う場合，その電磁場は図1のようにします．両端の入出力カップラーの半分は定在波(SW)，残りの部分は進行波(TW)にします．このようにすることで，カップラー空胴を除いた部分では，ほぼ実際の加速管内の電磁場となります．カップラー部の電場もほぼ実際と同じようになります．後述しますが，磁場の方は少し問題があります．バンチャー管を計算する場合では磁場よりも電場の方が重要なので，磁場の精度を犠牲にします．
このモデルでは，進行波型加速管を完全な軸対象として取り扱います．従って，SUPERFISH を使い電磁場を計算することができます．
|
| 進行波型加速管を計算する場合の領域の分け方．SWは定在波，TWは進行波の領域を表す． |
進行波型加速管の電磁場の計算方法
ここでは，先に示した進行波型加速管のモデルの電磁場をSUPERFISHで計算する具体的な方法を述べます．進行波と定在波のそれぞれの領域について，説明します．
進行波領域
電磁場の計算方法
よく知られているように，境界条件が異なる二つの電磁場を重ね合わせることにより，進行波を作ることができます．具体的な加速管では，図2に示す形状で計算します．1.5セル(1.5周期)の領域で，両端の境界条件が異なっています．一方は両端とも電気的短絡面，もう一方は磁気的短絡面です．電磁場(固有関数)は異なりますが，周波数(固有値)は同じです．これらの蓄積エネルギーを同一にし(規格化)，両者を π/2 [rad]位相を変えて重ね合わせると，進行波になります．
|
|
| 両端の境界条件がショート(SS) | 両端の境界条件がオープン(OO) |
もう少し，具体的な式で考えます．進行波の電場と磁場は， \begin{align} \vm{E}_{TW}(\vm{r},t)&=\left[\vm{E}_{SS}(\vm{r},\omega)+i\vm{E}_{OO}(\vm{r},\omega)\right]e^{i\omega t}\\ \vm{H}_{TW}(\vm{r},t)&=\left[\vm{H}_{OO}(\vm{r},\omega)-i\vm{H}_{SS}(\vm{r},\omega)\right]e^{i\omega t} \nonumber \\ \end{align} となります．ETWとHTWは，進行波の電磁場を表します．ESSとHSSは，SUPERFISHで計算された両端の境界条件がショート(図2)の場合の電磁場です．EOOとHOOは，SUPERFISHで計算された両端の境界条件がオープン(図3)の場合の電磁場です．式(2)のようにすると，SUPERFISHで計算できる定在波から進行波の電磁場が計算できます．
エネルギーの流れ
次に電磁場のエネルギーの流れを考えます．加速管の左端のセル(入力カップラー側)の進行波領域の左側の電磁場の電力は，加速管のインプットパワーになります．この加速管内の電磁場のエネルギーは，進行波の領域でインプットカップラーから出力カップラーに移動します．空洞の有限なQの影響により，その移動に伴い電磁場のエネルギーは減衰します．その減衰は，群速度とQ値により決まります．加速管のu進行方向をzとすると，RF電力(P)の変化は \begin{align} \ddiffA{P}{z}=-\cfrac{\omega P}{v_gQ_0} \end{align} となります．従って，加速管の1セル移動する毎に，RF電力が \(\exp(-\omega D/v_gQ)\) の割合で減衰します．また，蓄積エネルギーと郡速度には，次の関係 \begin{align} \cfrac{v_gU}{D}=P \end{align} があります．図1の加速管のモデルでは，隣の進行波の領域(TW)に移るときに，その領域の蓄積エネルギーをこの減衰を考慮します．
定在波領域
カップラーセルの半分は定在波にします．実際の加速管でも，両端付近は定在波になっています．カップラー空洞の真中で，定在波から進行波，あるいはその逆の変化が生じます．もちろん，これは実際の加速管とは異なります．これらの境界では，軸上電場を連続になるように，定在波の領域の電磁場を調整します．すると，境界領域にわたり，ほぼ電場は連続とすることができます．一方，磁場は不連続が生じます．図1のようなモデルを使う限り，これは致し方ないことです．
ビームダイナミックスに，TMモードの磁場が影響を与えるのは，粒子の速度(v)が高速より小さい場合です．特に，電子リニアックのバンチャー加速管のインプット空洞付近が影響を受けます．しかし，実際の加速管ではここでのビームの半径はそれほど大きくなく，その影響はそれほど大きくありません．この磁場に比べ，電場の方がはるかにビームダイナミックスに影響を与えます．したがって，先に示したように進行波と定在波の境界で，電場を連続にすることは，磁場を連続にするよりも良い近似になります．
|
|
| 左側のカップラーの半分(定在波) | 左側のカップラーの半分(定在波) |
自動計算のための Perl プログラム
図に示したSUPERFISHのインプットデータをひとつずつ手入力していたら，データの作成に大変時間がかかります．それを避けるためには，自動的に計算するプログラムが必要になります．ここでは，筆者が作成したプログラムを示し，そのソースコードを公開します．
著作権は放棄しますので，勝手に使ってください．修正も可能です．このプログラムを使いその結果を発表するときに，このwebページを引用してもらえるとうれしいです．引用は強制しませんので，自由にしてください．
プログラムの内容
効率よく計算するためには，進行波型加速管の電磁場の計算とGPTのインプットファイルの自動作成が必要です．そのために，筆者は Perl で記述した二つのプログラムを作成しました．GPTの電磁場の入力ファイルとなる「*.gdf」を作成する mk_TW_gdf_02B.pl と，それに応じた GPT の入力ファイル「*.in」の一部を作成する mk_GPTin_file.pl です．
mk_TW_gdf_02B.pl は入力ファイルに従い，SUPERFISHを起動し，加速管内の電磁場を計算します．このとき，自動的に空胴の直径を調整し，周波数を2856[MHz]に合わせます．進行波領域では，図2や図3の 1.5セルを計算し，必要な1セル分の電磁場を GPT用のファイル(*.gdf)に変換して保管します．同じように定在波領域(図4, 図5)も計算し，電磁場をファイルに保管します．また，郡速度や蓄積エネルギー，Q値などこの後必要な加速セルのパラメーターをファイル(*.tw, *.sw)に保管します．
mk_GPTin_file.pl は入力ファイルと mk_TW_gdf_02B.pl の出力ファイル(*.tw, *.sw)に従い，GPTの入力ファイル(*.in)の一部，進行波型加速管を記述する行を作成します．このプログラムを実行すると，二つのテキストファイル(*.para, *.forGPT)が作成されます．前者(*.para)は加速管に関するパラメーターが記述されます．後者(*.forGPT)は GPT のインプットファイル(*.in)のうち進行波型加速管の部分です．これは，GPTのインプットファイルにコピーし，貼り付けます．
以上の説明の通り，進行波型加速管のインプットデータを作るためには，二つのプログラムをmk_TW_gdf_02B.pl → mk_GPTin_file.plと順番に実行する必要があります．これらの詳細については，以降に述べます．
使い方
プログラムと実行方法
進行波型加速管の電磁場を計算し GPT のインプットファイルを作るためには，次のプログラム/インプットファイルを使います．
- GPT の gdf ファイルを作成する Perl のプログラム(mk_TW_gdf_03B.pl)．
- GPT のインプットファイル(*.in)を作成する Perl のプログラム(mk_GPTin_file.pl)．
- GPT で必要なファイルを一度に作成するバッチファイル(mk_file_TW.bat)．
- 進行波型加速管のインプットファイルの例(TW_structure.dat)
これらのプログラムは，Windows のコマンドプロンプトから，次のようにして実行します．これらの二つのプログラムのインプットファイル，ここでは「TW_structure.dat」は同一でなくてはなりません．
> perl mk_TW_gdf_02B.pl TW_structure.dat
> perl mk_GPTin_file.pl TW_structure.dat
バッチファイルを次のように実行することも可能です．
> mk_file_TW.bat
入力ファイル
この二つのプログラムには，進行波型加速管を定義する入力ファイルが必要です．以下の例(TW_structure.dat) を用いて，入力ファイルの書き方を説明します．
入力ファイルの例(TW_structure.dat)
001 #------------------------------------------------------ 002 # GPT の *.in ファイルから見た gdf のディレクトリー 003 #------------------------------------------------------ 004 $directory 005 ../TW/ 006 $end 007 #------------------------------------------------------ 008 $input_RF 009 5.0 #input RF power [MW] 010 $end 011 #------------------------------------------------------ 012 $sw 013 #pipeL 2a 014 30.000 20.00 015 30.000 20.00 016 $end 017 #------------------------------------------------------ 018 $cell 019 #beta 2a 2b t 020 0.700 22.000 83.158 5.000 021 0.800 21.900 83.515 5.000 022 0.900 21.800 83.493 5.000 023 1.000 21.700 83.470 5.000 024 1.000 21.600 83.470 5.000 025 1.000 21.500 83.470 5.000 026 1.000 21.400 83.470 5.000 027 $end
入力ファイル中のナンバー記号「#」はコメントを表します．この記号があると，それ以降はコメントとして取り扱われ，プログラム中では無視されます．行頭にあると，その行は全て無視されます．
入力ファイルは，4つのセクションから構成されます．これらのセクションは入力ファイルには必須です．以下にその説明を行いますので，上の入力ファイルを参照しながら，入力データの内容を理解してください．
- ディレクトリーセクション
- 「$directory」から「$end」までが，ディレクトリーセクションになります．
- GPTのインプットファイル(*.in)から見て，電磁場を記述したgdfファイルのディレクトリーを示します．ここでのディレクトリーの指定は，mk_GPTin_file.plのアウトプットファイル(*.para, *.forGPT)に記述される gdf ファイルのディレクトリーにのみ影響を与えます．
- 上の入力ファイルの例では，「../TW/」と指定していますので，親ディレクトリー「..」のサブディレクトリー「TW」に，gdfファイルが格納されていることを示します．もし，GPTのインプットファイルとgdfファイルが同じディレクトリーならば，「./」と指定します．
- 入力RF電力セクション
- 「$input_RF」から「$end」までが，入力RF電力セクションになります．
- 加速管の入力電力[MW]を指定します．ここでのRF電力の指定は，mk_GPTin_file.plのアウトプットファイル(*.para, *.forGPT)に記述される gdf ファイルの増倍計数(Multication factor)，すなわちGPTのエレメント「map25D_TM」の ffac に影響を与えます．
- 定在波セクション
- 「$sw」から「$end」までが，定在波セクションになります．
- 入力カップラーの定在波領域と出力カップラーの定在波領域の形状を二行に分けて記述します．最初の行が入力カップラーで次の行が出力カップラーです．
- 行内はタブ区切りで，パイプ部の長さ，引き続きパイプ部の直径を記述します(図6)．
- パイプ部と空洞部の接続部のコーナーRは，隣接するディスクのコーナーRと同一になります．
- 空洞部の直径は，隣り合う空胴の直径と同一になります．
- 空洞部の長さは，隣接する空胴(隣の進行波の1セル)の半分になります．
- 進行波セクション
- 「$cell」から「$end」までが，定在波セクションになります．
- 各行が進行波の空胴の1セルになります．
- 行内はタブ区切りで順に，\(\beta=v_p/c\)，ディスク孔径(2a[mm])，空胴直径(2b[mm])，ディスク厚さ(t[mm])を記述します(図6)．
- ディスクの先端Rは，ディスク厚さの半分になります．
- 空洞直径(2b)は\(2\pi/3\)モードが 2856[MHz]になる近似値とします．プログラム「mk_TW_gdf_02B.pl」では，周波数が\(10^{-6}\)の精度で2856[MHz]になるように空胴直径を再計算しています．ここで与えられた空胴直径は計算の出発点に過ぎません．
|
| インプットファイルに記載する寸法． |
出力ファイル
しかるべき入力ファイルを作成した後，プログラムと実行方法で示した方法で，二つのプログラム(mk_TW_gdf_02B.plとmk_GPTin_file.pl)を実行すると，さまざまなのファイルが作成されます．作成されたファイルとその使い方について説明します．
mk_TW_gdf_02B.plはSUPERFISHを自動起動し，図1の領域の電磁場のファイル(*.gdf)を作成します．それとともに，SUPERFISHの計算結果をまとめたファイル(*.tw と *.sw)が作成されます．具体例でアウトプットファイルについて説明します．例えば，インプットファイル名が「TW_structure.dat」とします．そして，コマンド
> perl mk_TW_gdf_03B.pl TW_structure.dat
を実行すると，以下のファイルが作成されます．
- 進行波領域の電磁場ファイル (TW_structure_XX_SS.gdfとTW_structure_XX_OO.gdf)
- これらは，図2や図3の左側から1セル分の電磁場のファイルです．径方向は，ディスクの半径+1[mm]です．
- XXは，セルの番号を表します．
- これらの gdf ファイルは後の GPT のエレメント「map25D_TM」で指定するマップファイルになります．
- GPTがインストールされていればファイルをダブルクリックするとGPTが起動し，電磁場の様子が分かります．
- 定在波領域の電磁場ファイル (TW_structure_inSW.gdf と TW_structure_outSW.gdf)
- これらは，図4や図5で示したカップラーの半分，すなわち定在波領域の電磁場のファイルです．径方向は，パイプ部の半径+1[mm]です．
- これらの gdf ファイルは後の GPT のエレメント「map25D_TM」で指定するマップファイルになります．
- GPTがインストールされていればファイルをダブルクリックするとGPTが起動し，電磁場の様子が分かります．
- SUPERFISHの計算結果のサマリーファイル (TW_structure.tw と TW_structure.sw)
- これらは，SUPERFISH計算された各種パラメーターをまとめたテキストファイルです．
- 拡張子が「tw」は進行波領域，「sw」は定在波領域の計算結果をまとめたものです．
- これらのファイルは，次のGPTのインプットファイルを作成するプログラム「mk_GPTin_file.pl」で使われます．
mk_GPTin_file.pl は，このインプットファイル (TW_structure.dat)と，mk_TW_gdf_02B.plのアウトプットファイル (TW_structure.tw と TW_structure.sw)から，GPTのインプットファイルとそのサマリーファイルを作成します．具体的には，次のコマンド
> perl mk_GPTin_file.pl TW_structure.dat
を実行すると，以下のファイルが作成されます．
- GPTのインプットファイルの一部になる (TW_structure.forGPT )
- これは，GPTのインプットファイル (*.in)の一部で，これをコピー/ペーストして使います．
- テキストファイルなので，適当なエディターで開き，内容の確認ができます．
- 計算結果のサマリーファイル (TW_structure.para )
- 進行波部分の計算結果のサマリーが書かれています．
- テキストファイルなので，適当なエディターで開き，内容の確認ができます．
GPTへの反映
できあがった進行波型加速管のファイルを使い図7のようなモデルを計算します．この計算過程で，GPTのインプットファイル中にこの進行波型加速管の反映の仕方を説明します．
|
| GPTでのシミュレーションモデル． |
このモデルのトラッキング計算を行う手順は，以下の通りです．前節で示した進行波型加速管のGPTの計算に必要なファイル( *.gdf と *.fotGPT ) は完成しているとします．
- GPTの計算を行うディレクトリーを作成します．仮に，このディレクトリー名を「GPT_dir」としましょう．
- この GPT の計算ディレクトリの親ディレクトリーに移動し，そこにサブディレクトリー「TW」を作成します．これは，リスト1:入力ファイルの 005 行目で，gdf ファイルのディレクトリーを「../TW/」としたためです．
- 前節で作成した全ての gdf フィアイルをこのディレクトリー「TW」に移動させます．
- 進行波型の加速管部を除いたの GPT インプットファイル ( *.in ) をディレクトリー「GPT_dir」に作成します．その例をリスト2に示します．このリストの 001 行 — 017 行と 067 行 — 070 行が，それに対応します．
- mk_GPTin_file.pl により作成されたファイル ( *.forGPT )，ここでは，「TW_structure.forGPT」を開き，GPTのインプットファイルにコピーします．リスト2 中の 024 行 — 064 が，コピーした部分です．
- その他，必要な情報を GPT のインプットファイル(リスト2)に書き込みます．
- 長さの単位を表す変数「mm」を001 行のように記述します．
- 加速管の角振動数を変数「omega」に定義します．ここでは，002 行と019 行 — 020のようにしました．
- 加速管の入力カップラーの中心座標を変数「StrZ」に定義します(021 行)．単位は，[mm]です．
- 加速管の初期位相の変数「StrP」に定義します(022 行)．単位は，[deg]です．
- そして，適当な名前，ここでは「TW_Buncher.in」を付けて保存します．
GPTのインプットファイル(TW_Buncher.in)
001 mm = 1e-3; 002 MHz = 1e6; 003 004 #======================================================== 005 # difinition: initial condition of electorn beam 006 #======================================================== 007 Energy=100e+3; # beam energy [eV] 008 009 G=1-qe*Energy/(me*c^2); 010 011 setparticles("beam",50000,me,qe,-1e-9); 012 setrxydist("beam","u", 1.25*mm, 2.5*mm); 013 setphidist("beam","u", 0, 2*pi); 014 setzdist("beam","u", -287*mm/2, 287*mm); 015 setGBrxydist("beam","u",0, 0); 016 setGBphidist("beam","u",0, 0); 017 setGdist("beam","u",G, 0); 018 019 facc = 2856*MHz; # Operatin freq. [Hz] 020 omega = 2*pi*facc; 021 StrZ = 100; # start point [mm] 022 StrP = 30; # initial phase [deg] 023 024 #========================================================= 025 # TW Acelerating Structure 026 # 2014.10.27 19:29:58 027 # data dir : ./ 028 # input data : TW_structure.dat 029 #========================================================== 030 031 #----- SW region ------- 032 map25D_TM("wcs", "z", (StrZ-39.746)*mm, "../TW/TW_structure_inSW.gdf", "R", "Z", "Er", "Ez", "H", -3.730944e+000, 0, (StrP+0.00)/deg, omega); 033 rmax("wcs", "z", (StrZ-24.746)*mm, 10.000*mm, 30.000*mm); 034 #----- cell x0701 ------- 035 map25D_TM("wcs", "z", (StrZ+0.000)*mm, "../TW/TW_structure_01_SS.gdf", "R", "Z", "Er", "Ez", "H", 6.132265e+000, 0, (StrP+0.00)/deg, omega); 036 map25D_TM("wcs", "z", (StrZ+0.000)*mm, "../TW/TW_structure_01_OO.gdf", "R", "Z", "Er", "Ez", "H", -5.796206e+000, 0, (StrP-90.00)/deg, omega); 037 rmax("wcs", "z", (StrZ+12.246)*mm, 11.000*mm, 5.000*mm); 038 #----- cell x0702 ------- 039 map25D_TM("wcs", "z", (StrZ+24.493)*mm, "../TW/TW_structure_02_SS.gdf", "R", "Z", "Er", "Ez", "H", 6.577072e+000, 0, (StrP-120.00)/deg, omega); 040 map25D_TM("wcs", "z", (StrZ+24.493)*mm, "../TW/TW_structure_02_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.156720e+000, 0, (StrP-210.00)/deg, omega); 041 rmax("wcs", "z", (StrZ+38.489)*mm, 10.950*mm, 5.000*mm); 042 #----- cell x0703 ------- 043 map25D_TM("wcs", "z", (StrZ+52.485)*mm, "../TW/TW_structure_03_SS.gdf", "R", "Z", "Er", "Ez", "H", 6.953617e+000, 0, (StrP-240.00)/deg, omega); 044 map25D_TM("wcs", "z", (StrZ+52.485)*mm, "../TW/TW_structure_03_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.456409e+000, 0, (StrP-330.00)/deg, omega); 045 rmax("wcs", "z", (StrZ+68.231)*mm, 10.900*mm, 5.000*mm); 046 #----- cell x0704 ------- 047 map25D_TM("wcs", "z", (StrZ+83.976)*mm, "../TW/TW_structure_04_SS.gdf", "R", "Z", "Er", "Ez", "H", 7.275100e+000, 0, (StrP-360.00)/deg, omega); 048 map25D_TM("wcs", "z", (StrZ+83.976)*mm, "../TW/TW_structure_04_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.706165e+000, 0, (StrP-450.00)/deg, omega); 049 rmax("wcs", "z", (StrZ+101.471)*mm, 10.850*mm, 5.000*mm); 050 #----- cell x0705 ------- 051 map25D_TM("wcs", "z", (StrZ+118.966)*mm, "../TW/TW_structure_05_SS.gdf", "R", "Z", "Er", "Ez", "H", 7.297176e+000, 0, (StrP-480.00)/deg, omega); 052 map25D_TM("wcs", "z", (StrZ+118.966)*mm, "../TW/TW_structure_05_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.726994e+000, 0, (StrP-570.00)/deg, omega); 053 rmax("wcs", "z", (StrZ+136.461)*mm, 10.800*mm, 5.000*mm); 054 #----- cell x0706 ------- 055 map25D_TM("wcs", "z", (StrZ+153.956)*mm, "../TW/TW_structure_06_SS.gdf", "R", "Z", "Er", "Ez", "H", 7.318745e+000, 0, (StrP-600.00)/deg, omega); 056 map25D_TM("wcs", "z", (StrZ+153.956)*mm, "../TW/TW_structure_06_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.744383e+000, 0, (StrP-690.00)/deg, omega); 057 rmax("wcs", "z", (StrZ+171.451)*mm, 10.750*mm, 5.000*mm); 058 #----- cell x0707 ------- 059 map25D_TM("wcs", "z", (StrZ+188.946)*mm, "../TW/TW_structure_07_SS.gdf", "R", "Z", "Er", "Ez", "H", 7.339642e+000, 0, (StrP-720.00)/deg, omega); 060 map25D_TM("wcs", "z", (StrZ+188.946)*mm, "../TW/TW_structure_07_OO.gdf", "R", "Z", "Er", "Ez", "H", -6.761152e+000, 0, (StrP-810.00)/deg, omega); 061 rmax("wcs", "z", (StrZ+206.441)*mm, 10.700*mm, 5.000*mm); 062 #----- SW region ------- 063 map25D_TM("wcs", "z", (StrZ+223.936)*mm, "../TW/TW_structure_outSW.gdf", "R", "Z", "Er", "Ez", "H", -5.199715e+000, 0, (StrP-840.00)/deg, omega); 064 rmax("wcs", "z", (StrZ+253.931)*mm, 10.000*mm, 30.000*mm); 065 066 067 #======================================================== 068 # beam control 069 #======================================================== 070 tout(0, 11/facc, 1/(12*facc));
以上で，すべてのGPTのインプットファイルができあがります．後は，適当なバッチファイル(例: run.bat )を作り，GPTを実行させます．
計算結果
図8や図9が GPT の計算結果をプロットした例です．バンチャーでのバンチングの様子が分かります．何も考えずに加速管のパラメーターを決めたので，バンチャーとしてはあまり良くないですね．GPTではアニメーションにすることもできますが，ファイルが大きいので，ここでは載せていません．
|
|
| 計算結果(r vs z)． | 電子密度プロット(z方向)． |
プログラムの変更
ここで公開しているプログラムの著作権は放棄しますので，必要に応じて変更は可能です．図6に示すように，公開しているプログラム「mk_TW_gdf_02B.pl」(リスト3)は，クラッシックな進行波型加速管を計算します．少しばかり，プログラムを変更するともう少しモダンな加速管の計算に対応させることもできます．実際，あいちシンクロトロンの加速管の設計では，このプログラムを少し変更したものを使いました．
俺はこんなクラッシックな加速管は嫌だ，もう少し精度良く計算したい—とのたまう人のために，変更箇所を簡単に説明します．
- SUPERFISHの計算メッシュサイズは変数 $SFPara->dx です．0157 行で定義しています．
- 加速管の周波数を表す変数は $freq です．0166 行で定義しています．Sバンド(2856[MHz])以外を計算する場合は，ここを変更します．
- 空胴にコーナーRを付けるような加速管に形状変更する場合には，プログラムの変更が必要です．主な変更点は，以下の通りです．
- 形状を決めるパラメーターに応じて，入力ファイルのフォーマットを変更する必要があるでしょう．
- また，空胴形状を表す構造体 Cell_Shape や $SW_Region の変更が必要です．
- データを読み込むサブルーチン read_cell_data や init_SW_regin の変更が必要です．
- SUPERFISHの計算用のファイルを作成するサブルーチン write_sf_data や mk_sw_region の変更が必要です．
- SUPERFISHでQ値を計算するときに必要な両端の境界のセグメント番号の再設定が必要です．0172 行と 0173 行です．
SUPERFISHを使いgdfファイルを作成するプログラム (mk_TW_gdf_03B.pl)
0001 use strict; 0002 use Class::Struct; 0003 use File::Find; 0004 use File::Basename; 0005 use Math::Trig; 0006 0007 our($file); 0008 our($SFPara); 0009 our($mm, $sec, $MHz); 0010 our($freq, $light, $lambda); 0011 our($it_max, $error); 0012 our($z_min, $z_max); 0013 our($cal); 0014 our(@input, @tw_result); 0015 our($seg_Wall_L, $seg_Wall_R); 0016 our($U_SS, @rL_SS, @rR_SS, @HL_SS, @HR_SS, $N_divR); 0017 our($U_OO, @rL_OO, @rR_OO, @EL_OO, @ER_OO); 0018 our($inSW, $outSW, @sw_result); 0019 0020 struct File => { 0021 input => '$', 0022 base => '$', 0023 af => '$', 0024 seg => '$', 0025 T35 => '$', 0026 in7 => '$', 0027 SFO => '$', 0028 tw => '$', 0029 GPT => '$', 0030 sf7 => '$', 0031 gdf => '$', 0032 sw => '$' 0033 }; 0034 0035 struct forSF => { 0036 title => '$', 0037 conv => '$', 0038 dx => '$' 0039 }; 0040 0041 struct INPUT => { 0042 shape => 'CELL_Shape', 0043 }; 0044 0045 struct ResutSFO =>{ 0046 freq => '$', 0047 U => '$', 0048 Q => '$', 0049 ro_v_Q => '$', 0050 Rs => '$', 0051 Ez_S => '$', 0052 Ez_E => '$' 0053 }; 0054 0055 struct Result => { 0056 i => '$', 0057 N => '$', 0058 shape => 'CELL_Shape', 0059 vg => '$', 0060 SS => 'ResutSFO', 0061 gdf_SS => '$', 0062 OO => 'ResutSFO', 0063 gdf_OO => '$' 0064 }; 0065 0066 struct CELL_Shape =>{ 0067 beta => '$', 0068 a2 => '$', 0069 b2 => '$', 0070 t => '$', 0071 D => '$' 0072 }; 0073 0074 struct CalMode => { 0075 i => '$', 0076 mode => '$', 0077 N => '$' 0078 }; 0079 0080 struct SW_Region => { 0081 pipe_L => '$', 0082 a2 => '$', 0083 ra => '$', 0084 b2 => '$', 0085 D_half => '$' 0086 }; 0087 0088 struct SW_Result => { 0089 i => '$', 0090 shape => 'SW_Region', 0091 SS => 'ResutSFO', 0092 gdf => '$', 0093 }; 0094 0095 print "SUPERFISH has been started.\n"; 0096 0097 0098 #------- cell shape parameter ----------- 0099 0100 &init($ARGV[0]); 0101 &read_cell_data; # read cell data form input cell geometory file 0102 0103 for($cal->i(1); $cal->i <= $cal->N; $cal->i($cal->i+1)){ 0104 printf "\n------------ Cell %d ---------------", $cal->i; 0105 printf "\n\2a: %.3f\tD: %.3f\tt: %.3f\n", $input[$cal->i]->shape->a2, $input[$cal->i]->shape->b2, $input[$cal->i]->shape->t; 0106 &tune_b_acc; # Tune b which is acc structure 0107 &cal_SF_two; # Calculate two modes (SO, OS) by SF 0108 &cal_vg; 0109 &write_result; 0110 } 0111 &init_SW_region; 0112 &mk_sw_region; 0113 &rm_files; 0114 0115 #============================================================================= 0116 # initialize 0117 #============================================================================= 0118 sub init{ 0119 my $input = $_[0]; 0120 my($path, $fname, $ext); 0121 0122 fileparse_set_fstype("Unix"); 0123 ($fname, $path, $ext) = fileparse($input, qr/\.[^\.]+$/); 0124 0125 printf "\n"; 0126 printf "input file: %s\n", $input; 0127 printf "file name : %s\n", $fname; 0128 printf "path : %s\n", $path; 0129 printf "ext : %s\n", $ext; 0130 printf "\n"; 0131 0132 if($ext ne ".dat"){ 0133 printf "\nfile extension error!!! file extension should be \".dat\". \n"; 0134 exit(0); 0135 } 0136 0137 chdir $path; 0138 0139 $file = File->new(); 0140 0141 $file -> input($fname.$ext); 0142 $file -> base($fname); 0143 $file -> af($fname.".af"); 0144 $file -> seg($fname.".seg"); 0145 $file -> T35($fname.".T35"); 0146 $file -> in7($fname.".in7"); 0147 $file -> SFO($fname.".SFO"); 0148 $file -> tw( $fname.".tw"); 0149 $file -> GPT($fname.".gpt"); 0150 $file -> sf7("OUTSF7.TXT"); 0151 $file -> gdf($fname.".gdf"); 0152 $file -> sw( $fname.".sw"); 0153 0154 $SFPara = forSF->new(); 0155 $SFPara -> title("S-Band Accelerator Structure"); 0156 $SFPara -> conv(0.1); 0157 $SFPara -> dx(0.25); 0158 0159 $cal = CalMode->new(); 0160 $cal->i(1); 0161 0162 $mm=1; 0163 $sec=1; 0164 $MHz=1e+6; 0165 0166 $freq=2856*$MHz; 0167 $light=2.99792458e+8; 0168 $lambda=$light/($freq); 0169 0170 $it_max=10; # 2b調整時の最大イタレーション回数 0171 $error=1e-6; # 2b調整時の周波数エラーの許容値 0172 $seg_Wall_L = 3; # 左端の壁のセグメント番号 0173 $seg_Wall_R = 12; # 右 0174 $N_divR = 1000; # ポインティングベクトル計算の壁の分割数 0175 0176 0177 unlink($file->af, $file->T35, $file -> in7, $file->SFO); 0178 unlink($file->tw, $file->sw); 0179 0180 opendir(DIR,"./") or die "could not open directory: $!"; 0181 while ($fname = readdir(DIR)){ 0182 if($fname =~ /\.gdf$/){unlink($fname);next;} 0183 if($fname =~ /\.TXT$/){unlink($fname);next;} 0184 if($fname =~ /\.TBL$/){unlink($fname);next;} 0185 if($fname =~ /\.log$/){unlink($fname);next;} 0186 if($fname =~ /\.INF$/){unlink($fname);next;} 0187 if($fname =~ /\.PMI$/){unlink($fname);next;} 0188 } 0189 0190 open(CELL_DATA,">".$file->tw) or die "open: $!"; 0191 printf CELL_DATA "#i\t"; 0192 printf CELL_DATA "beta\t"; 0193 printf CELL_DATA "2a\t"; 0194 printf CELL_DATA "2b\t"; 0195 printf CELL_DATA "D\t"; 0196 printf CELL_DATA "t\t"; 0197 printf CELL_DATA "vg\t"; 0198 printf CELL_DATA "f(SS)\t"; 0199 printf CELL_DATA "U\t"; 0200 printf CELL_DATA "Q\t"; 0201 printf CELL_DATA "r/Q\t"; 0202 printf CELL_DATA "Rs\t"; 0203 printf CELL_DATA "Ez_S\t"; 0204 printf CELL_DATA "Ez_E\t"; 0205 printf CELL_DATA "gdf\t"; 0206 printf CELL_DATA "f(OO)\t"; 0207 printf CELL_DATA "U\t"; 0208 printf CELL_DATA "Q\t"; 0209 printf CELL_DATA "r/Q\t"; 0210 printf CELL_DATA "Rs\t"; 0211 printf CELL_DATA "Ez_S\t"; 0212 printf CELL_DATA "Ez_E\t"; 0213 printf CELL_DATA "gdf\n"; 0214 0215 close CELL_DATA; 0216 0217 open(SW_DATA,">".$file->sw) or die "open: $!"; 0218 printf SW_DATA "#i\t"; 0219 printf SW_DATA "pipe_L\t"; 0220 printf SW_DATA "2a\t"; 0221 printf SW_DATA "ra\t"; 0222 printf SW_DATA "2b\t"; 0223 printf SW_DATA "Dh\t"; 0224 printf SW_DATA "f(SS)\t"; 0225 printf SW_DATA "U\t"; 0226 printf SW_DATA "Q\t"; 0227 printf SW_DATA "r/Q\t"; 0228 printf SW_DATA "Rs\t"; 0229 printf SW_DATA "Ez_S\t"; 0230 printf SW_DATA "Ez_E\t"; 0231 printf SW_DATA "gdf\n"; 0232 0233 close SW_DATA; 0234 } 0235 0236 0237 #============================================================================= 0238 # make seg_file for calculating Q factor 0239 #============================================================================= 0240 sub mk_seg{ 0241 my($i); 0242 my($seg_start, $seg_end); 0243 0244 $seg_start = $_[0]; 0245 $seg_end = $_[1]; 0246 0247 open(SEG,">".$file->seg) or die "open: $!"; 0248 0249 printf SEG "FieldSegments\n"; 0250 for($i=$seg_start; $i<$seg_end; $i++){ 0251 printf SEG "%d ",$i; 0252 } 0253 printf SEG "%d\n",$i; 0254 printf SEG "EndData\n"; 0255 printf SEG "End\n"; 0256 0257 close(SEG); 0258 0259 } 0260 0261 #============================================================================= 0262 # make seg_file for calculating vg 0263 #============================================================================= 0264 sub mk_seg2{ 0265 my($segL, $segR); 0266 0267 $segL = $_[0]; 0268 $segR = $_[1]; 0269 0270 open(SEG,">".$file->seg) or die "open: $!"; 0271 0272 printf SEG "FieldSegments\n"; 0273 printf SEG "%d %d\n", $segL, $segR; 0274 printf SEG "EndData\n"; 0275 printf SEG "End\n"; 0276 0277 close(SEG); 0278 0279 } 0280 0281 #============================================================================= 0282 # read data 0283 #============================================================================= 0284 sub read_cell_data{ 0285 my($i); 0286 my(@in_dat); 0287 my $find = 0; 0288 0289 $i=1; 0290 0291 open(CELL_DATA,"<".$file->input) or die "open: $!"; 0292 0293 while(<CELL_DATA>){ 0294 chomp; 0295 if(/^#/){next;} x07 先頭に'x07'があるとコメント文扱い 0296 if(/#/){s/(^.+)x07.*/$1/} x07 'x07'以降はコメント文扱い 0297 if($find or /^\$cell/){ # 先頭に$cellがあるまでスキップ 0298 if($find == 0){$find = 1;next;} 0299 }else{ 0300 next; 0301 } 0302 if(/^\$end/){last;} # データの終わり 0303 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 0304 0305 printf "%d\t%s\n",$i, $_; 0306 0307 @in_dat = split(/\s+/,$_); 0308 0309 @input[$i] = INPUT->new(); 0310 @input[$i]->shape(CELL_Shape->new()); 0311 0312 $input[$i] -> shape -> beta($in_dat[0]); 0313 $input[$i] -> shape -> a2( $in_dat[1]); 0314 $input[$i] -> shape -> b2( $in_dat[2]); 0315 $input[$i] -> shape -> t( $in_dat[3]); 0316 $input[$i] -> shape -> D( $in_dat[0]*$lambda/3.0*1e3); 0317 0318 $i++ 0319 } 0320 0321 close(CELL_DATA); 0322 0323 $cal->N($i-1); 0324 0325 &check_data; 0326 0327 } 0328 0329 #============================================================================= 0330 # tunning the accelerating cavity radius 0331 #============================================================================= 0332 sub tune_b_acc{ 0333 my(@tune_b_acc, @sf_frequency, $i, $j); 0334 0335 $i = $cal->i; 0336 0337 printf "\n"; 0338 $tune_b_acc[1]=$input[$i]->shape->b2; # 初期値 0339 for($j=1; $j<=$it_max; $j++){ 0340 write_sf_data(0, 0); 0341 system "autofish ".$file->af; 0342 $sf_frequency[$j]=&read_fr; 0343 printf "\t%d\t%f\t%f\n",$j,$tune_b_acc[$j],$sf_frequency[$j]/$MHz; 0344 if(abs(($sf_frequency[$j]-$freq)/$freq)<$error){last}; 0345 0346 0347 if($j==1){ 0348 $tune_b_acc[2]=$sf_frequency[1]*$tune_b_acc[1]/$freq; 0349 }else{ 0350 $tune_b_acc[$j+1]= 0351 ($tune_b_acc[$j]-$tune_b_acc[$j-1]) 0352 /($sf_frequency[$j]-$sf_frequency[$j-1]) 0353 *($freq-$sf_frequency[$j]) 0354 +$tune_b_acc[$j]; 0355 } 0356 $input[$i]->shape->b2($tune_b_acc[$j+1]); 0357 } 0358 } 0359 0360 #============================================================================= 0361 # calculate OS and SO modes 0362 #============================================================================= 0363 sub cal_SF_two{ 0364 my($i); 0365 my($freq, $U, $Q, $r, $r_ov_Q, $Rs); 0366 my($f_error); 0367 my($Ez_S, $Ez_E); 0368 my($gdf_file); 0369 0370 &mk_seg($seg_Wall_L+1, $seg_Wall_R-1); 0371 0372 $i = $cal->i; 0373 $tw_result[$i]=Result->new(); 0374 0375 $tw_result[$i]->i($i); 0376 0377 $tw_result[$i]->shape(CELL_Shape->new()); 0378 $tw_result[$i]->shape($input[$i]->shape); 0379 0380 #--------------- calculate SS mode ------------------ 0381 write_sf_data(1, 1); 0382 system "autofish ".$file->af; 0383 printf "\t ---- Short - Short mode ----\n"; 0384 ($freq, $U, $Q, $r_ov_Q, $Rs)=&read_SFO; 0385 ($Ez_S, $Ez_E) = &side_Ez(0.0, $tw_result[$i]->shape->D); 0386 0387 printf "\t Frequency \t%.5f\n", $freq/$MHz; 0388 printf "\t Stored Energy\t%e\n", $U; 0389 printf "\t Q value \t%.2f\n", $Q; 0390 printf "\t r/Q \t%e\n", $r_ov_Q; 0391 printf "\t Rs \t%e\n", $Rs; 0392 printf "\t Ez(start) \t%e\n", $Ez_S; 0393 printf "\t Ez(end) \t%e\n", $Ez_E; 0394 0395 $tw_result[$i]->SS(ResutSFO->new()); 0396 $tw_result[$i]->SS->freq($freq); 0397 $tw_result[$i]->SS->U($U); 0398 $tw_result[$i]->SS->Q($Q); 0399 $tw_result[$i]->SS->ro_v_Q($r_ov_Q); 0400 $tw_result[$i]->SS->Rs($Rs); 0401 $tw_result[$i]->SS->Ez_S($Ez_S); 0402 $tw_result[$i]->SS->Ez_E($Ez_E); 0403 $gdf_file = &mk_gdf_file_name("SS"); 0404 $tw_result[$i]->gdf_SS($gdf_file); 0405 &make_gdf($gdf_file); # make a gdf file 0406 0407 0408 #--------------- calculate OO mode ------------------ 0409 write_sf_data(0, 0); 0410 system "autofish ".$file->af; 0411 printf "\t ---- Open - Open mode ----\n"; 0412 ($freq, $U, $Q, $r_ov_Q, $Rs)=&read_SFO; 0413 ($Ez_S, $Ez_E)=&side_Ez(0.0, $tw_result[$i]->shape->D); 0414 0415 printf "\t Frequency\t%.5f\n", $freq/$MHz; 0416 printf "\t Stored Energy\t%e\n", $U; 0417 printf "\t Q value\t%.2f\n", $Q; 0418 printf "\t r/Q\t%e\n", $r_ov_Q; 0419 printf "\t Rs\t%e\n", $Rs; 0420 printf "\t Ez(start) \t%e\n", $Ez_S; 0421 printf "\t Ez(end) \t%e\n", $Ez_E; 0422 0423 $tw_result[$i]->OO(ResutSFO->new()); 0424 $tw_result[$i]->OO->freq($freq); 0425 $tw_result[$i]->OO->U($U); 0426 $tw_result[$i]->OO->Q($Q); 0427 $tw_result[$i]->OO->ro_v_Q($r_ov_Q); 0428 $tw_result[$i]->OO->Rs($Rs); 0429 $tw_result[$i]->OO->Ez_S($Ez_S); 0430 $tw_result[$i]->OO->Ez_E($Ez_E); 0431 $gdf_file = &mk_gdf_file_name("OO"); 0432 $tw_result[$i]->gdf_OO($gdf_file); 0433 &make_gdf($gdf_file); # make a gdf file 0434 0435 #------------- frequency check ------------ 0436 $f_error = abs($tw_result[$i]->OO->freq - $tw_result[$i]->SS->freq); 0437 if(100*$error< $f_error/$tw_result[$i]->OO->freq){ 0438 printf "Frequency error !!\n"; 0439 printf "\tSS mode : %f\n",$tw_result[$i]->SS->freq/$MHz; 0440 printf "\tOO mode : %f\n",$tw_result[$i]->OO->freq/$MHz; 0441 exit; 0442 } 0443 0444 } 0445 0446 #============================================================================= 0447 # calculate group velocity 0448 #============================================================================= 0449 sub cal_vg{ 0450 my($i); 0451 my($freq_SS, $Q_SS, $r_SS, $r_ov_Q_SS, $Rs_SS); 0452 my($freq_OO, $Q_OO, $r_OO, $r_ov_Q_OO, $Rs_OO); 0453 my($f_err, $U_err); # error 0454 my($a, $b, $Lz); 0455 my($U_tw, $RF_pwr_L, $RF_pwr_R); # Stored Energy and RF power 0456 my($vg_L, $vg_R); 0457 0458 $i = $cal->i; 0459 0460 $a = $tw_result[$i]->shape->a2/2.0; 0461 $b = $tw_result[$i]->shape->b2/2.0; 0462 $Lz =$tw_result[$i]->shape->D*1.5; 0463 0464 #--------------- calculate SS mode ------------------ 0465 write_sf_data(1, 1); 0466 system "autofish ".$file->af; 0467 printf "\t ---- Short - Short mode for vg cal ----\n"; 0468 ($freq_SS, $U_SS, $Q_SS, $r_ov_Q_SS, $Rs_SS)=&read_SFO; 0469 printf "\t Frequency\t%f.5\n", $freq_SS/$MHz; 0470 # --- mode check --------- 0471 $f_err = abs($freq_SS - $tw_result[$i]->SS->freq); 0472 $U_err = abs($U_SS - $tw_result[$i]->SS->U); 0473 if(1e-6 < $f_err/$freq_SS || 1e-6 < $U_err/$U_SS){ 0474 printf "mode check error at cal_vg SS mode !!\n"; 0475 exit(0); 0476 } 0477 0478 &set_field("Ht", 0.0, 0.0, $b, \@rL_SS, \@HL_SS); 0479 &set_field("Ht", $Lz, 0.0, $a, \@rR_SS, \@HR_SS); 0480 0481 #--------------- calculate OO mode ------------------ 0482 write_sf_data(0, 0); 0483 system "autofish ".$file->af; 0484 printf "\t ---- Open - Open mode for vg cal ----\n"; 0485 ($freq_OO, $U_OO, $Q_OO, $r_ov_Q_OO, $Rs_OO)=&read_SFO; 0486 printf "\t Frequency\t%f.5\n", $freq_OO/$MHz; 0487 # --- mode check --------- 0488 $f_err = abs($freq_OO - $tw_result[$i]->OO->freq); 0489 $U_err = abs($U_OO - $tw_result[$i]->OO->U); 0490 if(1e-6 < $f_err/$freq_OO || 1e-6 < $U_err/$U_OO){ 0491 printf "mode check error at cal_vg OO mode !!\n"; 0492 exit(0); 0493 } 0494 0495 &set_field("Er", 0.0, 0.0, $b, \@rL_OO, \@EL_OO); 0496 &set_field("Er", $Lz, 0.0, $a, \@rR_OO, \@ER_OO); 0497 0498 #------------- frequency_check ------------ 0499 $f_err = abs($freq_SS - $freq_OO); 0500 if(100*$error< $f_err/$freq_SS){ 0501 printf "Frequency error at cal_vg !!\n"; 0502 printf "\tSS mode : %f\n",$freq_SS/$MHz; 0503 printf "\tOO mode : %f\n",$freq_OO/$MHz; 0504 exit; 0505 } 0506 0507 ($U_tw, $RF_pwr_L, $RF_pwr_R) = &cal_RF_pwr; 0508 printf("\tStoroed Energy in 1.5 cell : %e [Joules]\n", $U_tw); 0509 0510 $vg_L = $Lz*1e-3*$RF_pwr_L/$U_tw; 0511 printf("\tRF power flow(left wall) : %e [Watt]\n", $RF_pwr_L); 0512 printf("\tGrout velocity(left wall) : %e [m/sec]\n", $vg_L); 0513 printf("\t\tvg/c(left wall) : %e\n", $vg_L/$light); 0514 0515 $vg_R = $Lz*1e-3*$RF_pwr_R/$U_tw; 0516 printf("\tRF power flow(right wall) : %e [Watt]\n", $RF_pwr_R); 0517 printf("\tGrout velocity(right wall) : %e [m/sec]\n", $vg_R); 0518 printf("\t\tvg/c(rith wall) : %e\n", $vg_R/$light); 0519 0520 $tw_result[$i]->vg($vg_L/$light); 0521 0522 } 0523 0524 0525 0526 #============================================================================= 0527 # check data 0528 #============================================================================= 0529 sub check_data{ 0530 my($i); 0531 0532 printf "\n=============== input data check ==============\n"; 0533 0534 for($i=1; $i<=$cal->N; $i++){ 0535 printf "Cell Number %d\n", $i; 0536 printf "beta \t%f\n",$input[$i]->shape->beta; 0537 printf "2a \t%f\n",$input[$i]->shape->a2; 0538 printf "2b \t%f\n",$input[$i]->shape->b2; 0539 printf "D \t%f\n",$input[$i]->shape->D; 0540 printf "t \t%f\n",$input[$i]->shape->t; 0541 printf "\n"; 0542 } 0543 0544 } 0545 0546 0547 #============================================================================= 0548 # write result 0549 #============================================================================= 0550 sub write_result{ 0551 my($i); 0552 $i = $cal->i; 0553 0554 open(TW_DATA,">>".$file->tw) or die "open: $!"; 0555 0556 printf TW_DATA "%d\t", $tw_result[$i]->i; 0557 printf TW_DATA "%.2f\t", $tw_result[$i]->shape->beta; 0558 printf TW_DATA "%.3f\t", $tw_result[$i]->shape->a2; 0559 printf TW_DATA "%.3f\t", $tw_result[$i]->shape->b2; 0560 printf TW_DATA "%.3f\t", $tw_result[$i]->shape->D; 0561 printf TW_DATA "%.3f\t", $tw_result[$i]->shape->t; 0562 printf TW_DATA "%e\t", $tw_result[$i]->vg; 0563 printf TW_DATA "%.5f\t", $tw_result[$i]->SS->freq/$MHz; 0564 printf TW_DATA "%e\t", $tw_result[$i]->SS->U; 0565 printf TW_DATA "%.2f\t", $tw_result[$i]->SS->Q; 0566 printf TW_DATA "%e\t", $tw_result[$i]->SS->ro_v_Q; 0567 printf TW_DATA "%e\t", $tw_result[$i]->SS->Rs; 0568 printf TW_DATA "%e\t", $tw_result[$i]->SS->Ez_S; 0569 printf TW_DATA "%e\t", $tw_result[$i]->SS->Ez_E; 0570 printf TW_DATA "%s\t", $tw_result[$i]->gdf_SS; 0571 printf TW_DATA "%.5f\t", $tw_result[$i]->OO->freq/$MHz; 0572 printf TW_DATA "%e\t", $tw_result[$i]->OO->U; 0573 printf TW_DATA "%.2f\t", $tw_result[$i]->OO->Q; 0574 printf TW_DATA "%e\t", $tw_result[$i]->OO->ro_v_Q; 0575 printf TW_DATA "%e\t", $tw_result[$i]->OO->Rs; 0576 printf TW_DATA "%e\t", $tw_result[$i]->OO->Ez_S; 0577 printf TW_DATA "%e\t", $tw_result[$i]->OO->Ez_E; 0578 printf TW_DATA "%s\n", $tw_result[$i]->gdf_OO; 0579 0580 close TW_DATA; 0581 0582 } 0583 0584 #============================================================================= 0585 # write superfish data 0586 #============================================================================= 0587 sub write_sf_data { 0588 my($nbslf, $nbsrt); 0589 my($beta, $a, $b, $D, $t); 0590 my($dx); 0591 my($xdri, $ydri); 0592 my($i); 0593 my(@xreg, @kreg, $kmax); 0594 my(@yreg, @lreg, $lmax); 0595 0596 $i = $cal->i; 0597 0598 $nbslf = $_[0]; 0599 $nbsrt = $_[1]; 0600 $beta = $input[$i]->shape->beta; 0601 $a = ($input[$i]->shape->a2)/2.0; 0602 $b = ($input[$i]->shape->b2)/2.0; 0603 $D = $input[$i]->shape->D; 0604 $t = $input[$i]->shape->t; 0605 $dx = $SFPara->dx; 0606 $xdri = 1.0*$D; 0607 $ydri = 0.7*$b; 0608 0609 $xreg[0] = 0.5*$D-0.5*$t; 0610 $xreg[1] = 0.5*$D+0.5*$t; 0611 $xreg[2] = 1.5*$D-0.5*$t; 0612 $kreg[0] = 0; 0613 $kreg[1] = int((0.5*$D-0.5*$t)/$dx); 0614 $kreg[2] = int((0.5*$D+0.5*$t)/$dx); 0615 $kreg[3] = int((1.5*$D-0.5*$t)/$dx); 0616 $kmax = int(1.5*$D/$dx); 0617 0618 $yreg[0] = $a; 0619 $yreg[1] = $a+0.5*$t; 0620 $lreg[0] = 0; 0621 $lreg[1] = int($a/$dx); 0622 $lreg[2] = int(($a+0.5*$t)/$dx); 0623 $lmax = int($b/$dx); 0624 0625 open(SF_DATA,">".$file->af) or die "open: $!"; 0626 0627 printf SF_DATA "%s\n", $SFPara -> title; 0628 printf SF_DATA " \$reg kprob=1, "; 0629 printf SF_DATA "dx=%f ,", $SFPara->dx; 0630 printf SF_DATA "conv=%f ,", $SFPara->conv; 0631 printf SF_DATA "kprob=1,\n"; 0632 printf SF_DATA " xreg=%f, %f, %f,\n", $xreg[0], $xreg[1], $xreg[2]; 0633 printf SF_DATA " kreg=%d, %d, %d, %d,\n", $kreg[0], $kreg[1], $kreg[2], $kreg[3]; 0634 printf SF_DATA " kmax=%d,\n", $kmax; 0635 printf SF_DATA " yreg=%f, %f,\n", $yreg[0], $yreg[1]; 0636 printf SF_DATA " lreg=%d, %d, %d,\n", $lreg[0], $lreg[1], $lreg[2]; 0637 printf SF_DATA " lmax=%d,\n", $lmax; 0638 printf SF_DATA " xdri=$xdri, ydri=$ydri,\n"; 0639 printf SF_DATA " nbsup=1, nbslo=0, nbslf=$nbslf, nbsrt=$nbsrt\n"; 0640 printf SF_DATA " freq=%f,\n",$freq/$MHz; 0641 printf SF_DATA " kmethod=1,\n"; 0642 printf SF_DATA " beta=%f\$\n",$beta; 0643 0644 0645 # ------ accelerating cell (half half) ----------- 0646 printf SF_DATA "!---------- half half cell -----------\n"; 0647 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 0648 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, $b; 0649 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.5*$D-0.5*$t, $b; 0650 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.5*$D-0.5*$t, $a+0.5*$t; 0651 printf SF_DATA 0652 " \$po nt=2, r=%f, theta=270, x0=%f, y0=%f \$\n", 0653 0.5*$t, 0.5*$D, $a+0.5*$t; 0654 0655 # ------ accelerating cell (full cell) ----------- 0656 printf SF_DATA "!---------- full cell -----------\n"; 0657 printf SF_DATA 0658 " \$po nt=2, r=%f, theta=0, x0=%f, y0=%f \$\n", 0659 0.5*$t, 0.5*$D, $a+0.5*$t; 0660 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.5*$D+0.5*$t, $b; 0661 printf SF_DATA " \$po x=%f, y=%f \$\n", 1.5*$D-0.5*$t, $b; 0662 printf SF_DATA " \$po x=%f, y=%f \$\n", 1.5*$D-0.5*$t, $a+0.5*$t; 0663 printf SF_DATA 0664 " \$po nt=2, r=%f, theta=270, x0=%f, y0=%f \$\n", 0665 0.5*$t, 1.5*$D, $a+0.5*$t; 0666 0667 printf SF_DATA " \$po x=%f, y=%f \$\n", 1.5*$D, 0.0; 0668 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 0669 0670 close SF_DATA; 0671 0672 } 0673 0674 #============================================================================= 0675 # READ frequency (***.SFO) 0676 #============================================================================= 0677 sub read_fr{ 0678 0679 my($resonance_fr); 0680 my(@oneline); 0681 0682 open(SFOUT,"<".$file->SFO) or die "open:$file $!";; 0683 while(<SFOUT>){ 0684 chomp; 0685 0686 if(/^Frequency\s+=/){ 0687 @oneline=split(/\s+/,$_); 0688 $resonance_fr=@oneline[2]; 0689 next; 0690 } 0691 } 0692 0693 close SFOUT; 0694 return $resonance_fr*$MHz; 0695 } 0696 0697 #============================================================================= 0698 # READ frequency and Rs etc (***.SFO) 0699 #============================================================================= 0700 sub read_SFO{ 0701 0702 my($f, $U, $Q, $r, $r_ov_Q); 0703 my(@oneline); 0704 0705 $f = 0; 0706 $U = 0; 0707 $Q = 0; 0708 $r = 0; 0709 $r_ov_Q = 0; 0710 0711 open(sfout,"<".$file->SFO) or die "open:$file $!"; 0712 0713 while(<sfout>){ 0714 chomp; 0715 0716 if(/^All calculated values below refer to the mesh geometry only./){ 0717 last; 0718 } 0719 } 0720 0721 while(<sfout>){ 0722 chomp; 0723 0724 if(/^Frequency\s+=/){ 0725 @oneline=split(/\s+/,$_); 0726 $f=$oneline[2]*$MHz; 0727 } 0728 0729 if(/^Stored\s+energy/){ 0730 @oneline=split(/\s+/,$_); 0731 $U=$oneline[3]; 0732 } 0733 0734 0735 if(/^Q\s+=/){ 0736 @oneline=split(/\s+/,$_); 0737 $Q=$oneline[2]; 0738 } 0739 0740 if(/Z\*T\*T/){ 0741 @oneline=split(/\s+/,$_); 0742 $r=$oneline[6]; 0743 } 0744 0745 if(/r\/Q/){ 0746 @oneline=split(/\s+/,$_); 0747 $r_ov_Q=$oneline[2]; 0748 } 0749 0750 if(/^Wall segments:/){last;} 0751 } 0752 0753 close sfout; 0754 return ($f, $U, $Q, $r_ov_Q, $r); 0755 } 0756 0757 0758 #============================================================================= 0759 # read Ez at side 0760 #============================================================================= 0761 sub side_Ez{ 0762 my($nz); 0763 my($d, $min_z, $max_z, $r); 0764 my($i); 0765 my($mode); 0766 my($Ez_S, $Ez_E); 0767 my($line, @data); 0768 0769 $min_z = $_[0]; 0770 $max_z = $_[1]; 0771 $r = 0.0; 0772 $nz=1; 0773 0774 open(IN7,">".$file->in7) or die "open: $!"; 0775 0776 printf IN7 "line noscreen\n"; 0777 printf IN7 "%.3f\t%.3f\t%.3f\t%.3f\n", $min_z, $r ,$max_z, $r; 0778 printf IN7 "%d\t%d\n", $nz; 0779 printf IN7 "end"; 0780 0781 close(IN7); 0782 0783 system "sf7"; 0784 0785 open(SFOUT,"<".$file->sf7); 0786 0787 while(<SFOUT>){ 0788 chomp; 0789 s/^\s*(.*?)\s*$/$1/; 0790 if(/^\(mm\)\s+\(mm\)/){last;} 0791 } 0792 0793 $line = <SFOUT>; 0794 chomp($line); 0795 $line =~ s/^\s*(.*?)\s*$/$1/; # 前後の空白の削除 0796 @data = split(/\s+/,$line); 0797 $Ez_S = $data[2]; 0798 0799 $line = <SFOUT>; 0800 chomp($line); 0801 $line =~ s/^\s*(.*?)\s*$/$1/; 0802 @data = split(/\s+/,$line); # 前後の空白の削除 0803 $Ez_E = $data[2]; 0804 0805 close(SFOUT); 0806 0807 return($Ez_S, $Ez_E); 0808 0809 } 0810 0811 0812 #============================================================================= 0813 # set Er or Ht field 0814 #============================================================================= 0815 sub set_field{ 0816 0817 my($fld, $z, $rmin, $rmax, $r, $EH); 0818 my($Nr, $dr); 0819 my($idx, $fac); 0820 my($i, $line, @data); 0821 0822 $fld = $_[0]; # Er or Ht 0823 $z = $_[1]; # H_theta を格納する z座標 0824 $rmin = $_[2]; # r座標(最小) 0825 $rmax = $_[3]; # r座標(最大) 0826 $r = $_[4]; # r座標の配列の参照渡し 0827 $EH = $_[5]; # 磁場を格納する配列の参照渡し 0828 0829 if($fld eq "Er"){ 0830 $idx = 3; 0831 $fac = 1e+6; 0832 }elsif($fld eq "Ht"){ 0833 $idx = 5; 0834 $fac = 1.0; 0835 }else{ 0836 printf("Error occured at set_field !!\n"); 0837 printf("1st arg of set_field should be \"Er\" or\" Ht\".\n"); 0838 printf("\t arg:%s\n", $fld); 0839 exit(0); 0840 } 0841 0842 $dr = ($rmax - $rmin)/$N_divR; 0843 0844 open(IN7,">".$file->in7) or die "open: $!"; 0845 printf IN7 "line noscreen\n"; 0846 printf IN7 "%.3f\t%.3f\t%.3f\t%.3f\n", $z, $rmin ,$z, $rmax; 0847 printf IN7 "%d\n", $N_divR; 0848 printf IN7 "end"; 0849 close(IN7); 0850 0851 system "sf7"; 0852 0853 open(OUTSF7,"<".$file->sf7) or die "at set_field:$file $!"; 0854 0855 while(<OUTSF7>){ 0856 chomp; 0857 s/^\s*(.*?)\s*$/$1/; 0858 if(/^\(mm\)\s+\(mm\)\s+\(MV\/m\)\s+\(MV\/m\)\s+\(MV\/m\)\s+\(A\/m\)/){last;} 0859 } 0860 0861 for($i=0; $i<=$N_divR; $i++){ 0862 $line = <OUTSF7>; 0863 chomp($line); 0864 $line =~ s/^\s*(.*?)\s*$/$1/; # 前後の空白の削除 0865 if($line=~/^\s*$/){last;} # 空行ならば読み込み終了 0866 @data = split(/\s+/,$line); 0867 if(0.001 < abs($data[0]-$z) || 0.001 < abs($data[1]-$i*$dr)){ 0868 printf("Error occured at set_field !!\n"); 0869 printf("\t\t%d\t\tZ = %f\t\t R = %f\n", $i, $data[0], $data[1]); 0870 exit(0); 0871 } 0872 $$r[$i] = $data[1]; 0873 $$EH[$i] = $fac*$data[$idx]; 0874 } 0875 0876 if(0.001<abs($$r[0]-$rmin) || 0.001<abs($$r[$N_divR]-$rmax)){ 0877 printf("Error occured at set_field !!\n"); 0878 printf("r[0]:%f\t\trmin:%f\n", $$r[0], $rmin); 0879 printf("r[%d]:%f\t\trmax:%f\n", $N_divR, $$r[$N_divR], $rmax); 0880 exit(0); 0881 } 0882 0883 close OUTSF7; 0884 0885 } 0886 0887 #============================================================================= 0888 # calculate RF power case of stored energy 2 [Joules] in 1.5 cells 0889 #============================================================================= 0890 sub cal_RF_pwr{ 0891 my($U_unit); 0892 my($fac_H, $fac_E); 0893 my($E, $H, $r, $dr); 0894 my($E, $H); 0895 my($pwrL, $pwrR); 0896 my($i); 0897 0898 $U_unit = 1; # " [Joules]" 0899 0900 $fac_H = sqrt($U_unit/$U_SS); 0901 $fac_E = sqrt($U_unit/$U_OO); 0902 0903 #========= left side calculation ======== 0904 $pwrL = 0.0; 0905 for($i=0; $i<$N_divR; $i++){ 0906 $E = $fac_E*($EL_OO[$i]+$EL_OO[$i+1])/2.0; 0907 $H = $fac_H*($HL_SS[$i]+$HL_SS[$i+1])/2.0; 0908 if(0.001<abs($rL_SS[$i]-$rL_OO[$i])){ 0909 printf("error occured at cal_RF_pwr !!\n"); 0910 printf("\t left side r coodinate.\n"); 0911 printf("\t%d\trL_SS:%f\trL_OO:%f\n", $rL_SS[$i], $rL_OO[$i]); 0912 exit(0); 0913 }else{ 0914 $r = ($rL_SS[$i]+$rL_SS[$i+1])/2.0*1e-3; # r [m] 0915 $dr = ($rL_SS[$i+1]-$rL_SS[$i])*1e-3; # dr [m] 0916 } 0917 $pwrL += 2*pi*$r*$dr*$E*$H; 0918 } 0919 $pwrL /= 2.0; 0920 0921 0922 #========= Right side calculation ======== 0923 $pwrR = 0.0; 0924 for($i=0; $i<$N_divR; $i++){ 0925 $E = $fac_E*($ER_OO[$i]+$ER_OO[$i+1])/2.0; 0926 $H = $fac_H*($HR_SS[$i]+$HR_SS[$i+1])/2.0; 0927 if(0.001<abs($rR_SS[$i]-$rR_OO[$i])){ 0928 printf("error occured at cal_RF_pwr !!\n"); 0929 printf("\t right side r coodinate.\n"); 0930 printf("\t%d\trR_SS:%f\trR_OO:%f\n", $rR_SS[$i], $rR_OO[$i]); 0931 exit(0); 0932 }else{ 0933 $r = ($rR_SS[$i]+$rR_SS[$i+1])/2.0*1e-3; # r [m] 0934 $dr = ($rR_SS[$i+1]-$rR_SS[$i])*1e-3; # dr [m] 0935 } 0936 0937 $pwrR += 2.0*pi*$r*$dr*$E*$H; 0938 0939 } 0940 $pwrR /= 2.0; 0941 0942 0943 return($U_unit*2.0, $pwrL, $pwrR); 0944 0945 } 0946 0947 0948 #============================================================================= 0949 # make a gdf file 0950 #============================================================================= 0951 sub make_gdf{ 0952 my($nz, $nr); 0953 my($d, $min_z, $max_z, $min_r, $max_r); 0954 my($i); 0955 my($gdf_file); 0956 0957 $gdf_file = $_[0]; 0958 0959 $i = $cal->i; 0960 0961 $d=0.1; 0962 $min_z = 0.0; 0963 $max_z = $tw_result[$i]->shape->D; 0964 $min_r = 0.0; 0965 $max_r = $tw_result[$i]->shape->a2/2.0+1.0; 0966 0967 $nz=int(($max_z-$min_z)/$d); 0968 $nr=int(($max_r-$min_r)/$d); 0969 0970 open(IN7,">".$file->in7) or die "open: $!"; 0971 0972 printf IN7 "rect noscreen\n"; 0973 printf IN7 "%.3f\t%.3f\t%.3f\t%.3f\n", $min_z, $min_r ,$max_z, $max_r; 0974 printf IN7 "%d\t%d\n", $nz, $nr; 0975 printf IN7 "end"; 0976 0977 close(IN7); 0978 0979 system "sf7"; 0980 system "fish2gdf -o ".$gdf_file." ". $file->sf7; 0981 0982 } 0983 0984 #============================================================================= 0985 # initialize standing wave region 0986 #============================================================================= 0987 sub mk_gdf_file_name{ 0988 my($mode, $file_name, $add, $ext, $gdf_file); 0989 my($i); 0990 0991 $mode = $_[0]; 0992 $i = $cal->i; 0993 0994 $file_name = $file->gdf; 0995 0996 $file_name =~ s/^(.*)[.].*$/$1/; # 拡張子削除 0997 $ext = $file->gdf; 0998 $ext =~ s/^(.*)[.](.*)$/$2/; # 拡張子を取り出し 0999 $add = sprintf("%02d", $i); 1000 $gdf_file = $file_name.'_'.$add.'_'.$mode.'.'.$ext; 1001 1002 return($gdf_file); 1003 1004 } 1005 1006 #============================================================================= 1007 # initialize standing wave region 1008 #============================================================================= 1009 sub init_SW_region{ 1010 my($N); 1011 my(@in_dat); 1012 my $find = 0; 1013 my $i = 0; 1014 1015 $inSW = SW_Region -> new(); 1016 $outSW = SW_Region -> new(); 1017 1018 $inSW -> ra(($input[1]->shape->t)/2.0); 1019 $inSW -> b2($input[1]->shape->b2); 1020 $inSW -> D_half(($input[1]->shape->D - $input[1]->shape->t)/2.0); 1021 1022 $N = $cal->N; 1023 1024 $outSW -> ra(($input[$N]->shape->t)/2.0); 1025 $outSW -> b2($input[$N]->shape->b2); 1026 $outSW -> D_half(($input[$N]->shape->D - $input[$N]->shape->t)/2.0); 1027 1028 open(CELL_DATA,"<".$file->input) or die "open: $!"; 1029 1030 while(<CELL_DATA>){ 1031 chomp; 1032 if(/^#/){next;} x07 先頭に'x07'があるとコメント文扱い 1033 if(/#/){s/(^.+)x07.*/$1/} x07 'x07'以降はコメント文扱い 1034 if($find or /^\$sw/){ # 先頭に$cellがあるまでスキップ 1035 if($find == 0){$find = 1;next;} 1036 }else{ 1037 next; 1038 } 1039 if(/^\$end/){last;} # データの終わり 1040 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 1041 1042 @in_dat = split(/\s+/,$_); 1043 1044 if($i==0){ 1045 $inSW -> pipe_L($in_dat[0]); 1046 $inSW -> a2($in_dat[1]); 1047 }else{ 1048 $outSW -> pipe_L($in_dat[0]); 1049 $outSW -> a2($in_dat[1]); 1050 } 1051 1052 $i++; 1053 } 1054 1055 close(CELL_DATA); 1056 1057 } 1058 1059 1060 #============================================================================= 1061 # write superfish data for input SW region 1062 #============================================================================= 1063 sub mk_sw_region{ 1064 my($freq, $U, $Q, $r, $r_ov_Q, $Rs); 1065 my($Ez_S, $Ez_E); 1066 my($gdf_file); 1067 1068 printf "\n\n------------ RF input coupler SW region ---------------\n"; 1069 &mk_seg(2, 5); 1070 &write_sf_inSW; 1071 system "autofish ".$file->af; 1072 ($freq, $U, $Q, $r_ov_Q, $Rs)=&read_SFO; 1073 ($Ez_S, $Ez_E) = &side_Ez(0, $inSW->pipe_L + $inSW->D_half); 1074 1075 printf "\t Frequency \t%f.5\n", $freq/$MHz; 1076 printf "\t Stored Energy\t%e\n", $U; 1077 printf "\t Q value \t%.2f\n", $Q; 1078 printf "\t r/Q \t%e\n", $r_ov_Q; 1079 printf "\t Rs \t%e\n", $Rs; 1080 printf "\t Ez(start) \t%e\n", $Ez_S; 1081 printf "\t Ez(end) \t%e\n", $Ez_E; 1082 1083 $sw_result[1] = SW_Result->new(); 1084 $sw_result[1] -> i(1); 1085 1086 $sw_result[1] -> shape(SW_Region->new()); 1087 $sw_result[1] -> shape($inSW); 1088 1089 $sw_result[1]->SS(ResutSFO->new()); 1090 $sw_result[1]->SS->freq($freq); 1091 $sw_result[1]->SS->U($U); 1092 $sw_result[1]->SS->Q($Q); 1093 $sw_result[1]->SS->ro_v_Q($r_ov_Q); 1094 $sw_result[1]->SS->Rs($Rs); 1095 $sw_result[1]->SS->Ez_S($Ez_S); 1096 $sw_result[1]->SS->Ez_E($Ez_E); 1097 $gdf_file = &mk_gdf_file_name_sw("inSW"); 1098 $sw_result[1]->gdf($gdf_file); 1099 &make_gdf_sw($sw_result[1]); 1100 1101 printf "\n\n------------ RF output coupler SW region ---------------\n"; 1102 &mk_seg(4, 7); 1103 &write_sf_outSW; 1104 system "autofish ".$file->af; 1105 ($freq, $U, $Q, $r_ov_Q, $Rs)=&read_SFO; 1106 ($Ez_S, $Ez_E) = &side_Ez(0, $outSW->pipe_L + $outSW->D_half); 1107 1108 printf "\t Frequency \t%f.5\n", $freq/$MHz; 1109 printf "\t Stored Energy\t%e\n", $U; 1110 printf "\t Q value \t%.2f\n", $Q; 1111 printf "\t r/Q \t%e\n", $r_ov_Q; 1112 printf "\t Rs \t%e\n", $Rs; 1113 printf "\t Ez(start) \t%e\n", $Ez_S; 1114 printf "\t Ez(end) \t%e\n", $Ez_E; 1115 1116 $sw_result[2] = SW_Result->new(); 1117 $sw_result[2] -> i(2); 1118 1119 $sw_result[2] -> shape(SW_Region->new()); 1120 $sw_result[2] -> shape($outSW); 1121 1122 $sw_result[2]->SS(ResutSFO->new()); 1123 $sw_result[2]->SS->freq($freq); 1124 $sw_result[2]->SS->U($U); 1125 $sw_result[2]->SS->Q($Q); 1126 $sw_result[2]->SS->ro_v_Q($r_ov_Q); 1127 $sw_result[2]->SS->Rs($Rs); 1128 $sw_result[2]->SS->Ez_S($Ez_S); 1129 $sw_result[2]->SS->Ez_E($Ez_E); 1130 $gdf_file = &mk_gdf_file_name_sw("outSW"); 1131 $sw_result[2]->gdf($gdf_file); 1132 &make_gdf_sw($sw_result[2]); 1133 1134 &write_sw_result; 1135 } 1136 1137 1138 #============================================================================= 1139 # write superfish data for input SW region 1140 #============================================================================= 1141 sub write_sf_inSW { 1142 my($nbslf, $nbsrt); 1143 my($pl, $a, $ra, $b, $D); 1144 my($dx); 1145 my($xdri, $ydri); 1146 my(@xreg, @kreg, $kmax); 1147 my(@yreg, @lreg, $lmax); 1148 1149 $pl = $inSW-> pipe_L; 1150 $a = ($inSW -> a2)/2; 1151 $ra = $inSW -> ra; 1152 $b = ($inSW -> b2)/2; 1153 $D = $inSW -> D_half; 1154 $dx = $SFPara->dx; 1155 1156 $nbslf = 1; 1157 $nbsrt = 1; 1158 $xdri = 0.99*($pl+$D); 1159 $ydri = 0.99*$b; 1160 1161 $xreg[0] = $pl-$ra; 1162 $xreg[1] = $pl; 1163 $kreg[0] = 0; 1164 $kreg[1] = int($xreg[0]/$dx); 1165 $kreg[2] = int($xreg[1]/$dx); 1166 $kmax = int(($pl+$D)/$dx); 1167 1168 $yreg[0] = $a; 1169 $yreg[1] = $a+$ra; 1170 $lreg[0] = 0; 1171 $lreg[1] = int($a/$dx); 1172 $lreg[2] = int(($a+$ra)/$dx); 1173 $lmax = int($b/$dx); 1174 1175 open(SF_DATA,">".$file->af) or die "open: $!"; 1176 1177 printf SF_DATA "%s\n", $SFPara -> title; 1178 printf SF_DATA " \$reg kprob=1, "; 1179 printf SF_DATA "dx=%f ,", $SFPara->dx; 1180 printf SF_DATA "conv=%f ,", $SFPara->conv; 1181 printf SF_DATA "kprob=1,\n"; 1182 printf SF_DATA " xreg=%f, %f,\n", $xreg[0], $xreg[1]; 1183 printf SF_DATA " kreg=%d, %d, %d,\n", $kreg[0], $kreg[1], $kreg[2]; 1184 printf SF_DATA " kmax=%d,\n", $kmax; 1185 printf SF_DATA " yreg=%f, %f,\n", $yreg[0], $yreg[1]; 1186 printf SF_DATA " lreg=%d, %d, %d,\n", $lreg[0], $lreg[1], $lreg[2]; 1187 printf SF_DATA " lmax=%d,\n", $lmax; 1188 1189 printf SF_DATA " xdri=$xdri, ydri=$ydri,\n"; 1190 printf SF_DATA " nbsup=1, nbslo=0, nbslf=$nbslf, nbsrt=$nbsrt\n"; 1191 printf SF_DATA " freq=%f,\n",$freq/$MHz; 1192 printf SF_DATA " kmethod=1,\n"; 1193 printf SF_DATA " beta=%f\$\n",1; 1194 1195 printf SF_DATA "!---------- input SW region -----------\n"; 1196 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 1197 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, $a; 1198 printf SF_DATA " \$po x=%f, y=%f \$\n", $pl-$ra, $a; 1199 printf SF_DATA 1200 " \$po nt=2, r=%f, theta=0, x0=%f, y0=%f \$\n", 1201 $ra, $pl-$ra, $a+$ra; 1202 printf SF_DATA " \$po x=%f, y=%f \$\n", $pl, $b; 1203 printf SF_DATA " \$po x=%f, y=%f \$\n", $pl+$D, $b; 1204 printf SF_DATA " \$po x=%f, y=%f \$\n", $pl+$D, 0.0; 1205 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 1206 1207 close SF_DATA; 1208 1209 } 1210 1211 #============================================================================= 1212 # write superfish data for output SW region 1213 #============================================================================= 1214 sub write_sf_outSW { 1215 my($nbslf, $nbsrt); 1216 my($pl, $a, $ra, $b, $D); 1217 my($dx); 1218 my($xdri, $ydri); 1219 my(@xreg, @kreg, $kmax); 1220 my(@yreg, @lreg, $lmax); 1221 1222 $pl = $outSW-> pipe_L; 1223 $a = ($outSW -> a2)/2; 1224 $ra = $outSW -> ra; 1225 $b = ($outSW -> b2)/2; 1226 $D = $outSW -> D_half; 1227 $dx = $SFPara->dx; 1228 1229 $nbslf = 1; 1230 $nbsrt = 1; 1231 $xdri = 0.0; 1232 $ydri = 0.99*$b; 1233 1234 $xreg[0] = $D; 1235 $xreg[1] = $D+$ra; 1236 $kreg[0] = 0; 1237 $kreg[1] = int($xreg[0]/$dx); 1238 $kreg[2] = int($xreg[1]/$dx); 1239 $kmax = int(($pl+$D)/$dx); 1240 1241 $yreg[0] = $a; 1242 $yreg[1] = $a+$ra; 1243 $lreg[0] = 0; 1244 $lreg[1] = int($a/$dx); 1245 $lreg[2] = int(($a+$ra)/$dx); 1246 $lmax = int($b/$dx); 1247 1248 open(SF_DATA,">".$file->af) or die "open: $!"; 1249 1250 printf SF_DATA "%s\n", $SFPara -> title; 1251 printf SF_DATA " \$reg kprob=1, "; 1252 printf SF_DATA "dx=%f ,", $SFPara->dx; 1253 printf SF_DATA "conv=%f ,", $SFPara->conv; 1254 printf SF_DATA "kprob=1,\n"; 1255 printf SF_DATA " xreg=%f, %f,\n", $xreg[0], $xreg[1]; 1256 printf SF_DATA " kreg=%d, %d, %d,\n", $kreg[0], $kreg[1], $kreg[2]; 1257 printf SF_DATA " kmax=%d,\n", $kmax; 1258 printf SF_DATA " yreg=%f, %f,\n", $yreg[0], $yreg[1]; 1259 printf SF_DATA " lreg=%d, %d, %d,\n", $lreg[0], $lreg[1], $lreg[2]; 1260 printf SF_DATA " lmax=%d,\n", $lmax; 1261 printf SF_DATA " xdri=$xdri, ydri=$ydri,\n"; 1262 printf SF_DATA " nbsup=1, nbslo=0, nbslf=$nbslf, nbsrt=$nbsrt\n"; 1263 printf SF_DATA " freq=%f,\n",$freq/$MHz; 1264 printf SF_DATA " kmethod=1,\n"; 1265 printf SF_DATA " beta=%f\$\n",1; 1266 1267 printf SF_DATA "!---------- output SW region -----------\n"; 1268 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 1269 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, $b; 1270 printf SF_DATA " \$po x=%f, y=%f \$\n", $D, $b; 1271 printf SF_DATA " \$po x=%f, y=%f \$\n", $D, $a+$ra; 1272 printf SF_DATA 1273 " \$po nt=2, r=%f, theta=270, x0=%f, y0=%f \$\n", 1274 $ra, $D+$ra, $a+$ra; 1275 printf SF_DATA " \$po x=%f, y=%f \$\n", $D+$pl, $a; 1276 printf SF_DATA " \$po x=%f, y=%f \$\n", $pl+$D, 0.0; 1277 printf SF_DATA " \$po x=%f, y=%f \$\n", 0.0, 0.0; 1278 1279 close SF_DATA; 1280 1281 } 1282 1283 1284 #============================================================================= 1285 # write result 1286 #============================================================================= 1287 sub write_sw_result{ 1288 my($i); 1289 1290 $i=1; 1291 1292 open(SW_DATA,">>".$file->sw) or die "open: $!"; 1293 1294 for($i=1; $i<=2; $i++){ 1295 printf SW_DATA "%d\t", $sw_result[$i]->i; 1296 printf SW_DATA "%.3f\t", $sw_result[$i]->shape->pipe_L; 1297 printf SW_DATA "%.3f\t", $sw_result[$i]->shape->a2; 1298 printf SW_DATA "%.3f\t", $sw_result[$i]->shape->ra; 1299 printf SW_DATA "%.3f\t", $sw_result[$i]->shape->b2; 1300 printf SW_DATA "%.3f\t", $sw_result[$i]->shape->D_half; 1301 printf SW_DATA "%.5f\t", $sw_result[$i]->SS->freq/$MHz; 1302 printf SW_DATA "%e\t", $sw_result[$i]->SS->U; 1303 printf SW_DATA "%.2f\t", $sw_result[$i]->SS->Q; 1304 printf SW_DATA "%e\t", $sw_result[$i]->SS->ro_v_Q; 1305 printf SW_DATA "%e\t", $sw_result[$i]->SS->Rs; 1306 printf SW_DATA "%e\t", $sw_result[$i]->SS->Ez_S; 1307 printf SW_DATA "%e\t", $sw_result[$i]->SS->Ez_E; 1308 printf SW_DATA "%s\n", $sw_result[$i]->gdf; 1309 } 1310 1311 close SW_DATA; 1312 1313 } 1314 1315 #============================================================================= 1316 # make a gdf file 1317 #============================================================================= 1318 sub make_gdf_sw{ 1319 my($nz, $nr); 1320 my($d, $min_z, $max_z, $min_r, $max_r); 1321 my($gdf_file); 1322 my($sw); 1323 my($command); 1324 1325 $sw = SW_Result->new(); 1326 $sw = $_[0]; 1327 1328 $d=0.1; 1329 $min_z = 0.0; 1330 $max_z = $sw->shape->D_half + $sw->shape->pipe_L; 1331 $min_r = 0.0; 1332 $max_r = ($sw->shape->a2)/2.0+1; 1333 1334 $nz=int(($max_z-$min_z)/$d); 1335 $nr=int(($max_r-$min_r)/$d); 1336 1337 # --- make .in7 file ----------- 1338 open(IN7,">".$file->in7) or die "open: $!"; 1339 printf IN7 "rect noscreen\n"; 1340 printf IN7 "%.3f\t%.3f\t%.3f\t%.3f\n", $min_z, $min_r ,$max_z, $max_r; 1341 printf IN7 "%d\t%d\n", $nz, $nr; 1342 printf IN7 "end"; 1343 close(IN7); 1344 1345 # --- make gdf file --------- 1346 $gdf_file = $sw->gdf; 1347 system "sf7"; 1348 $command = sprintf("fish2gdf -o %s %s",$gdf_file, $file->sf7); 1349 system $command; 1350 1351 } 1352 1353 #============================================================================= 1354 # make a gdf file name for SW region 1355 #============================================================================= 1356 sub mk_gdf_file_name_sw{ 1357 my($file_name, $ext, $gdf_file); 1358 my($mode); 1359 1360 $mode = $_[0]; 1361 $file_name = $file->gdf; 1362 $file_name =~ s/^(.*)[.].*$/$1/; # 拡張子削除 1363 $ext = $file->gdf; 1364 $ext =~ s/^(.*)[.](.*)$/$2/; # 拡張子を取り出し 1365 $gdf_file = $file_name.'_'.$mode.'.'.$ext; 1366 1367 return($gdf_file); 1368 } 1369 1370 1371 #============================================================================= 1372 # delete files 1373 #============================================================================= 1374 sub rm_files{ 1375 my($fname); 1376 1377 opendir(DIR,"./") or die "could not open directory: $!"; 1378 1379 unlink($file->af); 1380 unlink($file->seg); 1381 unlink($file->T35); 1382 unlink($file->in7); 1383 unlink($file->SFO); 1384 unlink($file->sf7); 1385 unlink("OUTAUT.TXT"); 1386 unlink("OUTFIS.TXT"); 1387 unlink("TAPE35.INF"); 1388 unlink($file->base.".PMI"); 1389 1390 }
自動計算のための Python プログラム
今時は Python だろう，ということで Python の自動計算プログラムを作成しました．GPT で進行波型加速管を計算するインプットファイルを作成します．ビームローディングも考慮した定常状態のRFを計算します．ビームローディングの計算方法については，時間のあるときに公表します．これは，定電界型や定インピーダンス型も含んだどんな電場分布でも計算可能です．
計算プログラムは「auto_superfish」です．右クリックでダウロード可能です．解凍すると，三個のフォルダー (programFiles, input_data, manual) が表示されます．programFiles に含まれる run.bat を実行すると，プログラムが起動し，実行されます．加速管他のパラメーターは，フォルダー: input_data の中にテキストファイルで記載します．
時間ができたら，詳しい説明を書きますので，悪しからず．
|
| Python バージョンの計算パラメーター． |
位相プロット
バンチャー加速管の設計では，粒子の位置とRF位相の関係が重要です．しかし，これをプロットするための機能は，GPT にはありません．そこで，GPT の出力ファイル(*.gdf)と進行波型加速管をつくるデータファイル(*.dat：mk_TW_gdf_03B.pl や mk_GPTin_file.plの入力)から，位相プロットを作成します．ここでは，そのプロットの作成方法を示します．
粒子情報ファイルの作成
位相プロットを作成するためには，粒子の情報をz座標軸毎に出力します．具体的には，リスト 2 の 071 行を削除し，そこに以下を追加します．
位相プロットのためのGPTのインプットファイル(TW_Buncher.in)
068 #=========================================== 069 # beam control 070 #=========================================== 071 zminmax("wcs", "I", -300*mm, 410*mm); 072 rmax("wcs","I", 30*mm); 073 screen("wcs", "I", 0.5*mm, 400*mm, 1*mm);
071と072行でスクリーンに到達しない，あるいは到達するまでに多大に時間がかかる粒子を削除します．ここで削除される粒子は通常であれば，壁に衝突して失われるため，計算精度に影響を与えません．このエレメントを追加しない場合，GPTは全ての粒子が screen を通過するまで計算を続けます．場合によっては膨大な CPU 時間が必要になるので，避けるべきです．073行で粒子情報をファイルに出力する z 座標を指定します．
このインプットファイル(TW_Buncher.in)で GPT を実行すると，指定した z 座標毎の粒子情報が書かれたバイナリーファイル(ex: TW-Buncher.gdf)が作成されます．このバイナリファイルには，粒子の位置と時刻が書かれています．RFの情報があれば，それを位相に変換することが可能です．
位相情報ファイルの作成
先に作成した GPT のアウトプットファイル(ex: TW-Buncher.gdf)には位相情報が含まれていません．そこで，最初に述べたように，進行波型加速管をつくるデータファイル(ex: TW_structure.dat)を用いて，位相情報を作成します．これら二つのファイルから，位相プロットのためのz座標毎の全ての粒子の位相情報が書かれたファイルを作成します．そのための Perl のプログラムをリスト5に示します．
GPTのインプットファイル(mk_phase_plot.pl)
001 use strict; 002 use Class::Struct; 003 use Math::Trig; 004 use File::Find; 005 use File::Basename; 006 use FileHandle; 007 008 struct cell_struct=>{ 009 beta => '$', 010 a2 => '$', 011 b2 => '$', 012 t => '$', 013 zs => '$', 014 zc => '$', 015 ze => '$', 016 D => '$', 017 p0 => '$' 018 }; 019 020 struct structure_struct=>{ 021 start_str => '$', # z-position of 1st cell center [mm] *.in 022 BaPhase => '$', # initial phase of buncher [deg] *.bat 023 N => '$' # number of cells *.dat 024 }; 025 026 027 our($structure, @cell); 028 our($structure_file, $GDF_in, $GDF_out, $temp_time, $temp_phase); 029 our($frequency, $omega, $omega_deg, $lambda, $k_deg, $c); 030 031 $structure = structure_struct->new(); 032 033 034 &initialize; 035 &read_cell_data_file; 036 &set_cell_pzD; 037 &print_cell_zp; 038 &print_structure_data; 039 040 &gdf2a; 041 &time2phase; 042 &asci2gdf; 043 &delete_temp; 044 045 #===================================================== 046 # initialize parameters 047 #===================================================== 048 sub initialize{ 049 050 $structure_file = "../../TW/TW_structure.dat"; 051 $GDF_in = "../TW_Buncher.gdf"; 052 $GDF_out = "phase.gdf"; 053 $temp_time = "temp_time.txt"; 054 $temp_phase = "temp_phase.txt"; 055 056 $frequency = 2.856e9; 057 $omega = 2*pi*$frequency; 058 $omega_deg = $frequency*360; 059 060 $c = 2.99792458e+8; # speed of light [m/sec] 061 $lambda = $c/$frequency; # wavelength [m] 062 $k_deg = 360.0/$lambda; 063 064 $structure->start_str(100); # start point: center of 1st cell(input cpl) 065 $structure->BaPhase(30); # initial phase [deg] 066 067 } 068 069 #===================================================== 070 # read input file 071 #===================================================== 072 sub read_cell_data_file{ 073 074 my($fh); 075 my(@in_dat); 076 077 $fh = new FileHandle; 078 079 printf "read input data form %s\n", $structure_file; 080 081 $fh->open("<".$structure_file) or die "open: $!"; 082 083 while(<$fh>){ 084 chomp; 085 if(/^#/){next;} x07 先頭に'x07'があるとコメント文 086 if(/#/){s/(^.+)x07.*/$1/} x07 'x07' 以降はコメント 087 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 088 089 if(/^\$directory/){&read_dir($fh);} 090 if(/^\$input_RF/){&read_RF($fh);} 091 if(/^\$cell/){&read_cell($fh);} 092 } 093 094 $fh->close; 095 096 } 097 098 #===================================================== 099 # read directory 100 #===================================================== 101 sub read_dir{ 102 my $fh; 103 104 $fh = $_[0]; 105 106 while(<$fh>){ 107 chomp; 108 if(/^#/){next;} x07 先頭に'x07'があるとコメント文 109 if(/#/){s/(^.+)x07.*/$1/} x07 'x07' 以降はコメント 110 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 111 if(/^\$end/){last;}; 112 113 printf "\tdirectrory:%s\n", $_; 114 115 last; 116 } 117 } 118 119 #===================================================== 120 # read RF 121 #===================================================== 122 sub read_RF{ 123 my $fh; 124 125 $fh = $_[0]; 126 127 while(<$fh>){ 128 chomp; 129 if(/^#/){next;} x07 先頭に'x07'があるとコメント文 130 if(/#/){s/(^.+)x07.*/$1/} x07 'x07' 以降はコメント 131 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 132 if(/^\$end/){last;}; 133 134 printf "\tRF Power:%.3f [MW]\n", $_; 135 136 last; 137 } 138 } 139 140 #===================================================== 141 # read cell 142 #===================================================== 143 sub read_cell{ 144 my $fh; 145 my @input; 146 my $i; 147 148 $fh = $_[0]; 149 150 printf "\tbeta\t2a\t2b\tt\trb\tn\n"; 151 152 $i=0; 153 while(<$fh>){ 154 chomp; 155 if(/^#/){next;} x07 先頭に'x07'があるとコメント文 156 if(/#/){s/(^.+)x07.*/$1/} x07 'x07' 以降はコメント 157 s/^\s*(.*?)\s*$/$1/; # 前後の空白を削除 158 if(/^\$end/){last;}; 159 @input = split(/\s+/,$_); 160 161 our $one_cell = cell_struct->new(); 162 163 $one_cell->beta($input[0]); 164 $one_cell->a2($input[1]); 165 $one_cell->b2($input[2]); 166 $one_cell->t($input[3]); 167 168 $i++; 169 $cell[$i]=$one_cell; 170 171 printf "\t%.3f\t%.3f\t%.3f\t%.3f\n", 172 $cell[$i]->beta, $cell[$i]->a2, $cell[$i]->b2, $cell[$i]->t; 173 } 174 175 $structure->N($i); 176 } 177 178 #===================================================== 179 # check every cell data 180 #===================================================== 181 sub set_cell_pzD{ 182 my($D, $i, $N); 183 184 #------------ sw region --------- 185 our $one_cell = cell_struct->new(); 186 187 $D=$cell[1]->beta * $lambda /3.0; 188 $one_cell->D($D); 189 $one_cell->beta($cell[1]->beta); 190 $one_cell->p0($structure->BaPhase); 191 $one_cell->ze($structure->start_str*1e-3); 192 $one_cell->zs($one_cell->ze-$D/2.0); 193 $one_cell->zc(($one_cell->zs + $one_cell->ze)/2.0); 194 195 $cell[0] = $one_cell; 196 197 #------------ in tw region --------- 198 199 for($i=1; $i<=$structure->N; $i++){ 200 $D=$cell[$i]->beta * $lambda /3.0; 201 202 $cell[$i]->D($D); 203 $cell[$i]->p0($structure->BaPhase-($i-1)*120); 204 $cell[$i]->zs($cell[$i-1]->ze); 205 $cell[$i]->ze($cell[$i]->zs+$D); 206 $cell[$i]->zc(($cell[$i]->zs + $cell[$i]->ze)/2.0); 207 } 208 209 #------------- sw region --------- 210 $N=$structure->N; 211 our $one_cell = cell_struct->new(); 212 213 $D=$cell[$N]->beta * $lambda /3.0; 214 $one_cell->D($D); 215 $one_cell->beta($cell[$N]->beta); 216 $one_cell->p0($cell[$N]->p0-120); 217 $one_cell->zs($cell[$N]->ze); 218 $one_cell->ze($cell[$N]->ze+$D/2.0); 219 $one_cell->zc(($one_cell->zs + $one_cell->ze)/2.0); 220 221 $cell[$N+1] = $one_cell; 222 223 } 224 225 #===================================================== 226 # print cell posiotion initial phase 227 #===================================================== 228 sub print_cell_zp{ 229 my($i); 230 231 printf "\n-----------------------------------------------------\n"; 232 printf "\t#\tzs\t\tzc\t\tze\t\tphase\n"; 233 printf "=====================================================\n"; 234 235 236 for($i=0; $i<=$structure->N+1; $i++){ 237 printf "\t%d\t", $i; 238 printf "%10.3f\t", $cell[$i]->zs * 1e3; 239 printf "%10.3f\t", $cell[$i]->zc * 1e3; 240 printf "%10.3f\t", $cell[$i]->ze * 1e3; 241 printf "%7.1f\n", $cell[$i]->p0; 242 } 243 244 printf "\n-----------------------------------------------------\n"; 245 246 } 247 248 #===================================================== 249 # print structure data 250 #===================================================== 251 sub print_structure_data{ 252 253 printf "\n-------------------------------------------------------\n"; 254 printf "\t 1st cell center position : %.3f[mm]\n", $structure->start_str; 255 printf "\t Buncher initial phase : %.3f[deg]\n", $structure->BaPhase; 256 printf "\t Numbe of cells : %d\n", $structure->N; 257 printf "-------------------------------------------------------\n\n"; 258 259 } 260 261 262 263 264 #===================================================== 265 # convert gdf file to ascii file 266 #===================================================== 267 sub gdf2a{ 268 269 system("gdf2a -w 15 -o $temp_time $GDF_in z t"); 270 271 } 272 273 274 #===================================================== 275 # convert to phase from time 276 #===================================================== 277 sub time2phase{ 278 279 my($z, $time); 280 my(@input); 281 my($n, $dp); 282 my $print_zp = 0; 283 284 open(TIME, "<$temp_time") or die "file:$temp_time\n $!"; 285 open(PHAS, ">$temp_phase") or die "file:$temp_phase\n $!"; 286 287 foreach(<TIME>){ 288 if(/^\s*\d/){ 289 chomp; 290 s/^\s+//; # 行頭の空白を削除 291 @input = split(/\s+/,$_); 292 $z = $input[0]; 293 $time = $input[1]; 294 if( ($z< $cell[0]->zs) || ($cell[$structure->N+1]->ze < $z)){next;} 295 ($n, $dp) = &z2phase($z, $time); 296 printf PHAS " %f %f\n", $z , $dp; 297 }elsif(/^\s*z\s+t/){ 298 if($print_zp !=1){ 299 print PHAS " z phase\n"; 300 } 301 $print_zp =1; 302 }else{ 303 print PHAS "#$_"; 304 } 305 } 306 307 close TIME; 308 close PHAS; 309 310 } 311 #===================================================== 312 # find cell 313 #===================================================== 314 sub z2phase{ 315 316 my $z = $_[0]; 317 my $t = $_[1]; 318 my($dz); 319 my($n, $dp, $p); 320 321 $n = &find_cell($z); 322 $dz = $z-$cell[$n]->zs; 323 324 if($n==0 || $n==$structure->N+1){ 325 $p = $cell[$n]->p0 + 360*$frequency*$t; 326 }else{ 327 $p = $cell[$n]->p0 + 360*($frequency*$t - $dz/($cell[$n]->beta * $lambda)); 328 } 329 330 $dp = $p-int($p/360.0)*360; 331 332 if($dp<-180){$dp += 360.0;} 333 if(180<$dp) {$dp -= 360.0;} 334 335 return ($n, $dp); 336 337 } 338 339 #===================================================== 340 # find cell 341 #===================================================== 342 sub find_cell{ 343 my $z = $_[0]; 344 my($min, $max, $tst); 345 346 $min = 0; 347 $max = $structure->N+1; 348 349 do{ 350 $tst = int(($max + $min)/2); 351 352 if($cell[$tst]->ze < $z){ 353 $min = $tst+1; 354 }else{ 355 $max = $tst; 356 } 357 358 }while($max != $min); 359 360 return $max; 361 } 362 363 364 #===================================================== 365 # convert gdf file to ascii file 366 #===================================================== 367 sub asci2gdf{ 368 369 system("asci2gdf -o $GDF_out $temp_phase"); 370 371 } 372 373 #===================================================== 374 # delete temp file 375 #===================================================== 376 sub delete_temp{ 377 378 system("rm -f $temp_time"); 379 system("rm -f $temp_phase"); 380 381 }
このプログラムを実行するためには，計算モデルに応じて以下の修正が必要です．
- 進行波型加速管の情報が書かれている「mk_TW_gdf_03B.pl」 や 「mk_GPTin_file.pl」の入力ファイルの指定．050 行で相対パスで指定します．
- GPTの計算結果である z座標毎に粒子が書かれた gdf ファイルの指定．051 行で相対パスで指定します．
- この Perl プログラムの計算結果である位相情報が書かれた gdf ファイル名の指定．052 行でファイル名を指定します．
- 加速管の運転周波数の指定．056 行で指定します．
- 第１セル(入力カップラーセル)の中心の z座標の指定．これは定在波と進行波の境の座標で，リスト2 の029 行の$StrZです．位相プロットを作成するプログラムの064 行にこの値を記述します．
- RFの初期位相の指定．リスト2 の030 行の$StrPの値です．位相プロットを作成するプログラムの065 行にこの値を記述します．
以上の修正を行いこの Perl プログラムを実行すると，z 座標毎に全ての粒子の位相情報が書かれた gdf ファイル(ex:phase.gdf)が作成されます．
位相プロットの作成
位相情報が書かれた gdf ファイルをダブルクリックすると，GPTが起動し，通常のプロット同様に位相プロットが作成できます．図10にその結果を示します．
|
| 位相プロット |
位相情報は加速管のRFがあるところのみ色づけされています．両端の定在波の部分では，波が止まっているので粒子の移動に従い位相が変化しています．
公開プログラム
バッチファイル (mk_GPTin_file.pl) を実行すると，(1) インプットファイル (TW_structure.dat) に記述された進行波型加速管の電磁場を計算し，(2) GPTのインプットファイルを作成します．作成されるファイルは，電磁場のマップファイル (*_OO.gdf, *_SS.gdf) と GPT のインプットファイル (*.forGPT)，SUPERFISH の計算結果 (*.para, *.sw, *.tw) です．
最新バージョン (ver 03B)
メッシュ作成時にエラーにより計算が停止する問題を修正しました．旧バージョンでは，セルの寸法やメッシュサイズに依存して，メッシュエラー作成時にエラーが発生することがありました．それが起きないように，メッシュの切り方を変えました．
| 進行波型加速管のインプットファイルの例 | TW_structure.dat |
| gdf ファイルを作成するプログラム | mk_TW_gdf_03B.pl |
| インプットファイル (*.in) を作成するプログラム | mk_GPTin_file.pl |
| GPT 入力データを作成するするバッチファイル | make_file_TW.bat |
旧バージョン
ver 02B
郡速度も計算できるようにしました．
| 進行波型加速管のインプットファイルの例 | TW_structure.dat |
| gdf ファイルを作成するプログラム | mk_TW_gdf_02B.pl |
| インプットファイル (*.in) を作成するプログラム | mk_GPTin_file.pl |
| GPT 入力データを作成するするバッチファイル | make_file_TW.bat |
ver 01A
初期バージョンです．
| 進行波型加速管のインプットファイルの例 | TW_structure.dat |
| gdf ファイルを作成するプログラム | mk_TW_gdf.pl |
| インプットファイル (*.in) を作成するプログラム | mk_GPTin_file.pl |
| GPT 入力データを作成するするバッチファイル | make_file_TW.bat |
ページ作成情報
参考資料
- Gneral Particle Tracer — Programmer's Reference Version 2.80 —
更新履歴
| 2013年7月 | ページの新規作成 |
| 2013年08月15日 | 群速度の計算を可能にした． |
| 2013年11月10日 | 位相プロットのプログラム追加． |
| 2018年07月22日 | メッシュ計算時のエラーの修正 (ver_03B)． |
| 2022年08月06日 | Python バージョンのプログラムの追加． |