-
Notifications
You must be signed in to change notification settings - Fork 39
PTC_TWISS equaltwiss
Calculation of optical parameters from A_ is done by equaltwiss subroutine.
- Linear parameters beta, alpha and gamma are calculated as follows:
beta11 = A_(1,1)*A_(1,1) + A_(1,2)*A_(1,2)
alfa11 = -(A_(1,1)*A_(2,1) + A_(1,2)*A_(2,2))
gama11 = A_(2,1)*A_(2,1) + A_(2,2)*A_(2,2)
For simplicity the matrix notation is used above. The first index corresponds to polynomial number, the second to coefficient next to variable. In the code beta11 would be
A_(1).sub.'100000'*A_(1).sub.'100000' +
A_(1).sub.'010000'*A_(1).sub.'010000'
It is implemented using these 2 loops
do i=1,c_%nd2-2*ndel
do jj=i,c_%nd2-2*ndel
do k=1,c_%nd-ndel
J(2*k-1)=1
lat(i,jj,k) = (A_script(i).sub.J)*(A_script(jj).sub.J)
J(2*k-1)=0
J(2*k)=1
lat(i,jj,k) = lat(i,jj,k) + (A_script(i).sub.J)*(A_script(jj).sub.J)
lat(jj,i,k) = lat(i,jj,k)
J(2*k)=
enddo
enddo
enddo
do i=1,c_%nd
do k=1,c_%nd
s1%beta(k,i)= lat(2*k-1,2*k-1,i)
s1%alfa(k,i)=-lat(2*k-1,2*k ,i)
s1%gama(k,i)= lat(2*k ,2*k ,i)
enddo
enddo
-
Linear dispersion: there different algorithms for 5D and 6D.
For 4D of course there is no dispersion.- 5D
The 5th coordinate is a parameter, so coefficient of the 5th variable is dispersion:
DX == disp1 = A_(1).sub.'000010' DPX == disp2 = A_(2).sub.'000010'
- 6D
Here it is more complicated because the 5th coordinate is changing due to the synchrotron motion and one needs to average over all the momenta. Analytical algorithm does not exist. It is approximation using the so called stroboscopic average. The algorithm is due to Chao and Sands See the following docs: http://iopscience.iop.org/1748-0221
It is implemented in subroutine dispersion6D It uses 3 Iaa matrices that are precomputed in ptc_twiss initialisation in routine initIaaMatrix
- A_script linear part is copied to a regular matrix amatrix
- It is inverted and stored in amatrix_inv
- H matrices are computed as
H(n) = (amatrix*Iaa(n))*amatrix_inv
where n = 1,2,3
Dispersion than is
disp(n) = Ha(n,5,3)/Ha(5,5,3)
-
Phase advances are calculated as,
mu1 = atan2(A_(1,2),A_(1,1))/twopi,
mu2 = atan2(A_(3,2),A_(3,1))/twopi,
mu2 = atan2(A_(6,2),A_(6,1))/twopi.
The way the A_ is initialised does not assure any initial value of phase advance. Therefore, delta of phase needs to be evaluated at each step and added to the phase variable. The calculation goes as follows:
- if result of atan2 is negative then it is made positive by adding 1
- difference w.r.t. phase at the previous element (kept in array testold) is computed
- the difference is added to phase
- If user requested deltap dependence of the Twiss parameters, they are evaluated in computeDeltapDependency. It is possible only for icase=5. The programmed equations correspond to derivatives of given variable. For example,
beta11 = A_(1,1)*A_(1,1) + A_(1,2)*A_(1,2)
d(beta11)/dp = d(A_(1,1))/dp * A_(1,1) + A_(1,1)*d(A_(1,1))/dp + ...
d(A_(1,1))/dp == A_(1).sub.'100010'
Result is stored in array beta_p, alfa_p, gama_p.
For dispersion up to order 4 is available, which in the icase=5 it trivial.
For example, 3rd order dispersion, which is d(X)^4/d^4p, is
s1%disp_p3(1) = 4.0*3.0*2.0*(A_script(1).sub.'000040')
4 3 2 of course is to take out factorial to obtain value of the partial derivative.
There is an implementation of the 2nd order dispersion for 6D, that worked for old PTC, but I am afraid that with the new one it is not correct anymore.
In fact this routine is never called for 56D or 6D, so the present code is dummy.