forked from lhcopt/lhctoolkit
-
Notifications
You must be signed in to change notification settings - Fork 0
/
FineCouplingCorrection.madx
119 lines (95 loc) · 3 KB
/
FineCouplingCorrection.madx
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
! Fine tuning of coupling correction
! S.Fartoukh March 2009
! 2011/11/21 From SLHCV3.0 R. De Maria
! 2016/04/27 R. De Maria & F.F Van der Veken
! work around mad-x bug on integer part of the tune in presence of coupling for ctap calculation
if(mylhcbeam==1)
{kqtf=kqtf.b1;kqtf0=kqtf.b1;
kqtd=kqtd.b1;kqtd0=kqtd.b1;
kqtf.b1:=kqtf;
kqtd.b1:=kqtd;};
if(mylhcbeam > 1)
{kqtf=kqtf.b2;kqtf0=kqtf.b2;
kqtd=kqtd.b2;kqtd0=kqtd.b2;
kqtf.b2:=kqtf;
kqtd.b2:=kqtd;};
! closest tune
qmid=(qx0-qx00+qy0-qy00)*0.5;
match;
global, q1=qx00+qmid,q2=qy00+qmid;
vary, name=kqtf, step=1.E-9;
vary, name=kqtd, step=1.E-9;
lmdif, calls=50, tolerance=1.E-10;
endmatch;
! Quick minimization based on linear machine
cmrskew0 = cmrskew;
cmiskew0 = cmiskew;
twiss; qx=table(summ,q1); qy=table(summ,q2);
cta0 = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
closest0 = cta0;
!value,closest0;
cmrskew = cmrskew0+cta0/2.;
twiss; qx=table(summ,q1); qy=table(summ,q2);
ctap = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,ctap;
cmrskew = cmrskew0-cta0/2.;
twiss; qx=table(summ,q1); qy=table(summ,q2);
ctam = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,ctam;
cmrskew = cmrskew0+(ctam^2-ctap^2)/2./cta0;
twiss; qx=table(summ,q1); qy=table(summ,q2);
cta0 = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,cta0;
cmiskew = cmiskew0+cta0/2.;
twiss; qx=table(summ,q1); qy=table(summ,q2);
ctap = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,ctap;
cmiskew = cmiskew0-cta0/2.;
twiss; qx=table(summ,q1); qy=table(summ,q2);
ctam = abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,ctam;
cmiskew = cmiskew0+(ctam^2-ctap^2)/2./cta0;
twiss; qx=table(summ,q1); qy=table(summ,q2);
closest1 =abs(2*(qx-qy)-round(2*(qx-qy)))/2;
!value,closest1;
!Empirical minimisation
match;
global, q1=qx00+qmid, q2=qy00+qmid;
vary, name=kqtf, step=1.E-9;
vary, name=kqtd, step=1.E-9;
lmdif, calls=100,tolerance=1.E-6;
endmatch;
match;
global, q1=qx00+qmid, q2=qy00+qmid;
vary, name=cmrskew, step=1.E-9;
vary, name=cmiskew, step=1.E-9;
lmdif, calls=150, tolerance=2.E-6;
endmatch;
match;
global, q1=qx00+qmid, q2=qy00+qmid;
vary, name=kqtf, step=1.E-9;
vary, name=kqtd, step=1.E-9;
lmdif, calls=100, tolerance=1.E-7;
endmatch;
match;
global, q1=qx00+qmid, q2=qy00+qmid;
vary, name=cmrskew, step=1.E-9;
vary, name=cmiskew, step=1.E-9;
lmdif, calls=150, tolerance=2.E-7;
endmatch;
match;
global, q1=qx00+qmid, q2=qy00+qmid;
vary, name=kqtf, step=1.E-9;
vary, name=kqtd, step=1.E-9;
lmdif, calls=100, tolerance=1.E-7;
endmatch;
twiss; qx=table(summ,q1); qy=table(summ,q2);
closest2=abs(2*(qx-qy)-round(2*(qx-qy)))/2;
value, closest0, closest1, closest2, cmrskew, cmiskew;
if(mylhcbeam==1)
{kqtf.b1=kqtf0;
kqtd.b1=kqtd0;};
if(mylhcbeam > 1)
{kqtf.b2=kqtf0;
kqtd.b2=kqtd0;};
return;