Conference Room SAD
[thread display] [new arrival display] [word search] [past log] [管理用]

Subject emittance計算
Date: 2012/11/27(Tue) 14:44:53
ContributorKentaro Harada

  PF 制御計算機が今年夏に更新されたのですが、最近、各種制御の為に SAD を連続で呼び出すと
だんだん計算機が重くなり、仕舞に落とす(ユーザープロセスで制御計算機本体が落ちていいのか
どうかはまあ、置いといて……)という現象が生じ、制御グループの方が SAD を更新しようとして
いるのですが、古い SAD と新しい SAD で計算結果が違います。
  具体的には、design_example.sad で、emittance の計算結果が違います。関連部分のみ
抜粋すると、

MOMENTUM = 1GEV;
DRIFT L1 =(L =1 )
;
BEND B =(L =2 ANGLE =.2617993877991494 )
;
QUAD QRD =(L =.5 K1 =-.1078621179953663 )
QRF =(L =1 K1 =.19457865291109586 )
QSD =(L =1 K1 =-.21812755206494103 )
QSF =(L =1 K1 =.15519425429347913 )
QD =(L =1 K1 =-.1560624058480439 )
QF =(L =1 K1 =.1269625766602043 )
;
SEXT SF =(L =1 K2 =.12116257080688167 )
SD =(L =1 K2 =-.20401978160297035 )
;
CAVI CA1 =(L =1 VOLT =1e+06 HARM =100 )
;
LINE

sample = (IP1 -QRD -L1 -L1 -L1 -CA1 -L1
-L1 -L1 -QRF -L1 -L1 -L1 -B
-L1 -L1 -L1 -QSD -L1 -L1 -L1
-L1 -L1 -L1 -L1 -L1 -QSF L1
SF L1 B L1 L1 L1 QD
L1 SD L1 B L1 L1 L1
QF L1 SF L1 B L1 L1
L1 QD L1 SD L1 B L1
L1 L1 QF L1 SF L1 B
L1 L1 L1 QD L1 SD L1
B L1 L1 L1 QF L1 SF
L1 B L1 L1 L1 QD L1
SD L1 B L1 L1 L1 QF
L1 SF L1 B L1 L1 L1
QD L1 SD L1 B L1 L1
L1 QSF L1 L1 L1 L1 L1
L1 L1 L1 QSD L1 L1 L1
B L1 L1 L1 QRF L1 L1
L1 CA1 L1 L1 L1 QRD -QRD
-L1 -L1 -L1 -CA1 -L1 -L1 -L1
-QRF -L1 -L1 -L1 -B -L1 -L1
-L1 -QSD -L1 -L1 -L1 -L1 -L1
-L1 -L1 -L1 -QSF -L1 -L1 -L1
-B -L1 -SD -L1 -QD -L1 -L1
-L1 -B -L1 -SF -L1 -QF -L1
-L1 -L1 -B -L1 -SD -L1 -QD
-L1 -L1 -L1 -B -L1 -SF -L1
-QF -L1 -L1 -L1 -B -L1 -SD
-L1 -QD -L1 -L1 -L1 -B -L1
-SF -L1 -QF -L1 -L1 -L1 -B
-L1 -SD -L1 -QD -L1 -L1 -L1
-B -L1 -SF -L1 -QF -L1 -L1
-L1 -B -L1 -SD -L1 -QD -L1
-L1 -L1 -B -L1 -SF -L1 QSF
L1 L1 L1 L1 L1 L1 L1
L1 QSD L1 L1 L1 B L1
L1 L1 QRF L1 L1 L1 CA1
L1 L1 L1 QRD -IP1 )
;

ON LOG;
ON RAD FLUC COD RFSW RADCOD;
FFS;
$FORM = "15.10";
PageWidth = 1999;

USE sample;
cell calc;
emit;
disp $$$;

bpos = LINE["POSITION","B*"];

GetIntegralH[num_]:=(
IH = 0;
! Print["Middle Parameters of ",LINE["NAME",num]];
bex = Twiss["BX",num+0.5];
alx = Twiss["AX",num+0.5];
gax = (1+alx*alx)/bex;
etx = Twiss["EX",num+0.5];
epx = Twiss["EPX",num+0.5];
L = LINE["L",num]/2;
tha = LINE["ANGLE",num]/2;
rho = L/tha;
! Print["BX = ", bex," AX = ",alx," GX = ", gax];
! Print["EX = ", etx," EPX= ",epx];
! Print["rho= ", rho," L = ",L];
aaa = 2*etx*epx*tha+2*epx*rho*(Sin[tha]-tha)+2*(etx-rho)*(1-Cos[tha])+rho/2*(1-Cos[2*tha]);
bbb = (tha+2*epx*epx*tha-Sin[2*tha]/2+4*epx*(1-Cos[tha]))/2;
ccc = (2*etx*etx*tha+4*etx*rho*(Sin[tha]-tha)+rho*rho*(3*tha-4*Sin[tha]+Sin[2*tha]/2))/2;
IH = IH + aaa * alx + bbb * bex + ccc * gax;
bex = Twiss["BX",num+0.5];
alx = -1*Twiss["AX",num+0.5];
gax = (1+alx*alx)/bex;
etx = Twiss["EX",num+0.5];
epx = -1*Twiss["EPX",num+0.5];
aaa = 2*etx*epx*tha+2*epx*rho*(Sin[tha]-tha)+2*(etx-rho)*(1-Cos[tha])+rho/2*(1-Cos[2*tha]);
bbb = (tha+2*epx*epx*tha-Sin[2*tha]/2+4*epx*(1-Cos[tha]))/2;
ccc = (2*etx*etx*tha+4*etx*rho*(Sin[tha]-tha)+rho*rho*(3*tha-4*Sin[tha]+Sin[2*tha]/2))/2;
IH = IH + aaa * alx + bbb * bex + ccc * gax;
! Print[IH];
IH
);


IntH = 0;
Do[
IntH = IntH + GetIntegralH[bpos[i]];
,{i,Length[bpos]}];

Print[IntH];
cq = 3.832/10000;
gamma = 1000/0.51099906;
emittance = cq * gamma * gamma * IntH /2 /Pi /rho;
Print[emittance];

stop;
stop;

です。下側はシンクロトロン積分の解析式でエミッタンスを計算する積分式です。
(BMの中央から、上流側、下流側の両方向へ積分式を計算し、和を求めます。)

  古い SAD では、emit の結果は、
*** Welcome to SAD Ver.1.0.10.2b1 built at 2007-09-10 11:46:34 JST ***
Emittance X = 2.45855E-7 m
ですが、afsad、amsad、PF でコンパイルし直した新しい sad では、
*** Welcome to SAD Revision 3995 built at 2012-11-22 20:42:37 +0900 ***
Emittance X = 4.58933E-7 m
です。上記の計算式では、244.9820119666 nmrad (古い方とほぼ同じ)と出ます。
afsad、amsadでは、新しい方とほぼ同じ結果になります。
 
  ちなみに、PF のラティスで計算すると、emittance の計算結果にこんなに
大きな差は出ません。
旧: Emittance X = 3.69803E-8 m
新: Emittance X = 3.71388E-8 m
解析式:37.1172052922 nmrad
 
  design example のラティスは、強いて言えば若干和共鳴に近いですが、
他に何か不運なケースなのでしょうか? 何かノウハウがありましたらお願い
致します。


- 関連一覧ツリー (Click ▼ to display all articles in a thread.)