FINDSYNCORBIT finds closed orbit, synchronous with the RF cavity
and momentum deviation dP (first 5 components of the phase space vector)
by numerically solving for a fixed point
of the one turn map M calculated with LINEPASS
(X, PX, Y, PY, dP2, CT2 ) = M (X, PX, Y, PY, dP1, CT1)
under constraints CT2 - CT1 = dCT = C(1/Frev - 1/Frev0) and dP2 = dP1 , where
Frev0 = Frf0/HarmNumber is the design revolution frequency
Frev = (Frf0 + dFrf)/HarmNumber is the imposed revolution frequency
IMPORTANT!!! FINDSYNCORBIT imposes a constraint (CT2 - CT1) and
dP2 = dP1 but no constraint on the value of dP1, dP2
The algorithm assumes time-independent fixed-momentum ring
to reduce the dimensionality of the problem.
To impose this artificial constraint in FINDSYNCORBIT
PassMethod used for any element SHOULD NOT
1. change the longitudinal momentum dP (cavities , magnets with radiation)
2. have any time dependence (localized impedance, fast kickers etc).
FINDSYNCORBIT(RING,dCT) is 5x1 vector - fixed point at the
entrance of the 1-st element of the RING (x,px,y,py,dP)
FINDSYNCORBIT(RING,dCT,REFPTS) is 5-by-Length(REFPTS)
array of column vectors - fixed points (x,px,y,py,dP)
at the entrance of each element indexed in REFPTS array.
REFPTS is an array of increasing indexes that select elements
from the range 1 to length(RING)+1.
See further explanation of REFPTS in the 'help' for FINDSPOS
FINDSYNCORBIT(RING,dCT,REFPTS,GUESS) - same as above but the search
for the fixed point starts at the initial condition GUESS
Otherwise the search starts from [0 0 0 0 0 0]'.
GUESS must be a 6-by-1 vector;
[ORBIT, FIXEDPOINT] = FINDSYNCORBIT( ... )
The optional second return parameter is
a 6-by-1 vector of initial conditions
on the close orbit at the entrance to the RING.
See also FINDORBIT4, FINDORBIT6.0001 function [orbit, fixedpoint] = findsyncorbit(RING,dCT,varargin) 0002 %FINDSYNCORBIT finds closed orbit, synchronous with the RF cavity 0003 % and momentum deviation dP (first 5 components of the phase space vector) 0004 % by numerically solving for a fixed point 0005 % of the one turn map M calculated with LINEPASS 0006 % 0007 % (X, PX, Y, PY, dP2, CT2 ) = M (X, PX, Y, PY, dP1, CT1) 0008 % 0009 % under constraints CT2 - CT1 = dCT = C(1/Frev - 1/Frev0) and dP2 = dP1 , where 0010 % Frev0 = Frf0/HarmNumber is the design revolution frequency 0011 % Frev = (Frf0 + dFrf)/HarmNumber is the imposed revolution frequency 0012 % 0013 % IMPORTANT!!! FINDSYNCORBIT imposes a constraint (CT2 - CT1) and 0014 % dP2 = dP1 but no constraint on the value of dP1, dP2 0015 % The algorithm assumes time-independent fixed-momentum ring 0016 % to reduce the dimensionality of the problem. 0017 % To impose this artificial constraint in FINDSYNCORBIT 0018 % PassMethod used for any element SHOULD NOT 0019 % 1. change the longitudinal momentum dP (cavities , magnets with radiation) 0020 % 2. have any time dependence (localized impedance, fast kickers etc). 0021 % 0022 % 0023 % FINDSYNCORBIT(RING,dCT) is 5x1 vector - fixed point at the 0024 % entrance of the 1-st element of the RING (x,px,y,py,dP) 0025 % 0026 % FINDSYNCORBIT(RING,dCT,REFPTS) is 5-by-Length(REFPTS) 0027 % array of column vectors - fixed points (x,px,y,py,dP) 0028 % at the entrance of each element indexed in REFPTS array. 0029 % REFPTS is an array of increasing indexes that select elements 0030 % from the range 1 to length(RING)+1. 0031 % See further explanation of REFPTS in the 'help' for FINDSPOS 0032 % 0033 % FINDSYNCORBIT(RING,dCT,REFPTS,GUESS) - same as above but the search 0034 % for the fixed point starts at the initial condition GUESS 0035 % Otherwise the search starts from [0 0 0 0 0 0]'. 0036 % GUESS must be a 6-by-1 vector; 0037 % 0038 % [ORBIT, FIXEDPOINT] = FINDSYNCORBIT( ... ) 0039 % The optional second return parameter is 0040 % a 6-by-1 vector of initial conditions 0041 % on the close orbit at the entrance to the RING. 0042 % 0043 % See also FINDORBIT4, FINDORBIT6. 0044 % 0045 if ~iscell(RING) 0046 error('First argument must be a cell array'); 0047 end 0048 0049 d = 1e-6; % step size for numerical differentiation 0050 dps = 1e-12; % convergence threshold 0051 %dps=eps; % convergence threshold 0052 max_iterations = 20; 0053 0054 if nargin >= 4 % Check if guess argument was supplied 0055 if isnumeric(varargin{2}) && all(eq(size(varargin{2}),[6,1])) 0056 Ri=varargin{2}; 0057 else 0058 error('The last argument GUESS must be a 6-by-1 vector'); 0059 end 0060 else 0061 Ri = zeros(6,1); 0062 end 0063 0064 D5 = [d*eye(5) zeros(5,1); zeros(1,6)]; 0065 %D5 = [0.5*d*eye(5) -0.5*d*eye(5) zeros(5,1); zeros(1,11)]; 0066 theta5 = [0 0 0 0 dCT]'; 0067 0068 args={}; 0069 change=Inf; 0070 itercount = 0; 0071 while (change > dps) && (itercount < max_iterations) 0072 RMATi= Ri(:,ones(1,6)) + D5; 0073 %RMATi= Ri(:,ones(1,11)) + D5; 0074 RMATf = linepass(RING,RMATi,args{:}); 0075 Rf = RMATf(:,end); 0076 % compute the transverse part of the Jacobian 0077 J5 = (RMATf([1:4,6],1:5)-RMATf([1:4,6],6*ones(1,5)))/d; 0078 %J5 = (RMATf([1:4,6],1:5)-RMATf([1:4,6],6:10))/d; 0079 Ri_next = Ri + [(diag([1 1 1 1 0]) - J5)\(Rf([1:4,6])-[Ri(1:4);0]-theta5);0]; 0080 change = norm(Ri_next - Ri); 0081 Ri = Ri_next; 0082 itercount = itercount+1; 0083 args={'KeepLattice'}; 0084 end 0085 0086 if itercount == max_iterations 0087 warning('Maximum number of iterations reached. Possible non-convergence') 0088 end 0089 0090 if (nargin<3) || (isscalar(varargin{1}) && (varargin{1}==(length(RING)+1))) 0091 % return only the fixed point at the entrance of RING{1} 0092 orbit=Ri(1:5,1); 0093 else % 3-rd input argument - vector of reference points along the Ring 0094 % is supplied - return orbit 0095 orb6 = linepass(RING,Ri,varargin{1},'KeepLattice'); 0096 orbit = orb6(1:5,:); 0097 end 0098 0099 if nargout==2 0100 fixedpoint=Ri; 0101 end