opalcycl.tex 28.3 KB
 snuverink_j committed Sep 07, 2017 1 \chapter{\textit{OPAL-cycl}}  snuverink_j committed Sep 06, 2017 2 3 4 5 6 7 8 \label{chp:opalcycl} \index{OPAL-cycl} % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Introduction}  snuverink_j committed Sep 07, 2017 9 \textit{OPAL-cycl}, as one of the flavors of the \textit{OPAL} framework, is a fully three-dimensional parallel beam dynamics simulation program dedicated to future high intensity cyclotrons and FFAG. it tracks multiple particles which takes into account the space charge effects. For the first time in the cyclotron community, \textit{OPAL-cycl} has the capability of tracking multiple bunches simultaneously  snuverink_j committed Sep 06, 2017 10 11 12 and take into account the beam-beam effects of the radially neighboring bunches (we call it neighboring bunch effects for short) by using a self-consistent numerical simulation model.  snuverink_j committed Sep 07, 2017 13 Apart from the multi-particle simulation mode, \textit{OPAL-cycl} also has two other serial tracking modes for conventional cyclotron machine design. One mode is the single particle tracking mode, which is a useful tool for the preliminary design of a new cyclotron. It allows one to compute basic parameters, such as reference orbit, phase shift history, stable region, and matching phase ellipse. The other one is the tune calculation mode, which can be used to compute the betatron oscillation frequency.  snuverink_j committed Sep 06, 2017 14 15 This is useful for evaluating the focusing characteristics of a given magnetic field map.  snuverink_j committed Sep 07, 2017 16 In addition, the widely used plugin elements, including collimator, radial profile probe, septum, trim-coil field and charge stripper, are currently implemented in \textit{OPAL-cycl}. These functionalities are very useful for designing, commissioning and upgrading of cyclotrons and FFAGs.  snuverink_j committed Sep 06, 2017 17 18  \section{Tracking modes}  snuverink_j committed Sep 06, 2017 19 According to the number of particles defined by the argument \texttt{npart} in \texttt{beam} see~Chapter~\ref{beam} ,  snuverink_j committed Sep 07, 2017 20 \textit{OPAL-cycl} works in one of the following three operation modes automatically.  snuverink_j committed Sep 06, 2017 21 22 23 24  \subsection{Single Particle Tracking mode}  snuverink_j committed Sep 07, 2017 25  In this mode, only one particle is tracked, either with acceleration or not. Working in this mode, \textit{OPAL-cycl}  snuverink_j committed Sep 06, 2017 26 27  can be used as a tool during the preliminary design phase of a cyclotron.  snuverink_j committed Sep 07, 2017 28  The 6D parameters of a single particle in the initial local frame must be read from a file. To do this, in the \textit{OPAL} input file,  snuverink_j committed Sep 06, 2017 29  the command line \texttt{DISTRIBUTION} see~Chapter~\ref{distribution} should be defined like this:  snuverink_j committed Sep 06, 2017 30 \begin{verbatim}  snuverink_j committed Sep 06, 2017 31  Dist1: DISTRIBUTION, TYPE=fromfile, FNAME="PartDatabase.dat";  snuverink_j committed Sep 06, 2017 32 \end{verbatim}  snuverink_j committed Sep 08, 2017 33  where the file \textit{PartDatabase.dat} should have two lines:  snuverink_j committed Sep 06, 2017 34 \begin{verbatim}  snuverink_j committed Sep 06, 2017 35 36  1 0.001 0.001 0.001 0.001 0.001 0.001  snuverink_j committed Sep 06, 2017 37 \end{verbatim}  snuverink_j committed Sep 06, 2017 38 The number in the first line is the total number of particles.  snuverink_j committed Sep 06, 2017 39 In the second line the data represents $x, p_x, y,$$p_y, z, p_z$ in the local reference frame. Their units are described in Section~\ref{variablesopalcycl}.  snuverink_j committed Sep 06, 2017 40   snuverink_j committed Sep 11, 2017 41 Please don't try to run this mode in parallel environment. You should believe that a single processor of the 21^{st} century is capable of doing  snuverink_j committed Sep 06, 2017 42 43 44 45 46 the single particle tracking. \subsection{Tune Calculation mode} In this mode, two particles are tracked, one with all data is set to zero is the reference particle and another one is an off-centering particle  snuverink_j committed Sep 07, 2017 47  which is off-centered in both $r$ and $z$ directions. Working in this mode, \textit{OPAL-cycl} can calculate the betatron oscillation frequency $\nu_r$ and $\nu_z$ for different energies to evaluate the focusing characteristics  snuverink_j committed Sep 06, 2017 48 49 50 51 52  for a given magnetic field. Like the single particle tracking mode, the initial 6D parameters of the two particles in the local reference frame must be read from a file, too. In the file should has three lines:  snuverink_j committed Sep 06, 2017 53 \begin{verbatim}  snuverink_j committed Sep 06, 2017 54 55 56  2 0.0 0.0 0.0 0.0 0.0 0.0 0.001 0.0 0.0 0.0 0.001 0.0  snuverink_j committed Sep 06, 2017 57 \end{verbatim}  snuverink_j committed Sep 06, 2017 58 59  When the total particle number equals 2, this mode is triggered automatically.  snuverink_j committed Sep 06, 2017 60 Only the element \texttt{CYCLOTRON} in the beam line is used and other elements are omitted if any exists.  snuverink_j committed Sep 06, 2017 61 62 63 64 65 66 67 68 69 70 71  Please don't try to run this mode in parallel environment, either. \subsection{Multi-particle tracking mode} In this mode, large scale particles can be tracked simultaneously, either with space charge or not, either single bunch or multi-bunch, either serial or parallel environment, either reading the initial distribution from a file or generating a typical distribution, either running from the beginning or restarting from the last step of a former simulation.  snuverink_j committed Sep 07, 2017 72  Because this is the main mode as well as the key part of \textit{OPAL-cycl},  snuverink_j committed Sep 06, 2017 73  we will describe this in detail in Section~\ref{opalcycl:MultiBunch}.  snuverink_j committed Sep 06, 2017 74 75 76 77 78  % -- -- -- -- -- -- Section -- -- -- -- -- --  snuverink_j committed Sep 07, 2017 79 \section{Variables in \textit{OPAL-cycl}}  snuverink_j committed Sep 06, 2017 80 81 82 83 84 85 86 \label{sec:variablesopalcycl} \index{OPAL-cycl:Variables} \label{sec:opalcycl:canon} \index{Canonical Variables} \index{Variables!Canonical}  snuverink_j committed Sep 07, 2017 87 \textit{OPAL-cycl} uses the following canonical variables to describe the motion of particles:  snuverink_j committed Sep 06, 2017 88 89 90  \begin{description} \item[X]  snuverink_j committed Sep 08, 2017 91  Horizontal position $x$ of a particle in given global Cartesian coordinates [{m}].  snuverink_j committed Sep 06, 2017 92 93  \item[PX]  snuverink_j committed Sep 08, 2017 94  Horizontal canonical momentum [{eV/c}].  snuverink_j committed Sep 06, 2017 95 96  \item[Y]  snuverink_j committed Sep 08, 2017 97  Longitudinal position $y$ of a particle in global Cartesian coordinates [{m}].  snuverink_j committed Sep 06, 2017 98 99  \item[PY]  snuverink_j committed Sep 08, 2017 100  Longitudinal canonical momentum [{eV/c}].  snuverink_j committed Sep 06, 2017 101 102  \item[Z]  snuverink_j committed Sep 08, 2017 103  Vertical position $z$ of a particle in global Cartesian coordinates [{m}].  snuverink_j committed Sep 06, 2017 104 105  \item[PZ]  snuverink_j committed Sep 08, 2017 106  Vertical canonical momentum [{eV/c}].  snuverink_j committed Sep 06, 2017 107 108 109  \end{description}  snuverink_j committed Sep 08, 2017 110 The independent variable is: \textbf{t} [{s}].  snuverink_j committed Sep 06, 2017 111 112   snuverink_j committed Sep 07, 2017 113 \subsubsection{NOTE: unit conversion of momentum in \textit{OPAL-t} and \textit{OPAL-cycl}}  snuverink_j committed Sep 11, 2017 114 Convert $\beta_x \gamma$ [dimensionless] to [{mrad}],  snuverink_j committed Sep 06, 2017 115 116 117 118 119 120 121  \label{eq:betagamma1} (\beta\gamma)_{\text{ref}}=\frac{P}{m_0c}=\frac{Pc}{m_0c^2}, \label{eq:betagamma2}  snuverink_j committed Sep 11, 2017 122 P_x[{mrad}]=1000\times\frac{(\beta_x\gamma)}{(\beta\gamma)_{\text{ref}}}.  snuverink_j committed Sep 06, 2017 123 124   snuverink_j committed Sep 08, 2017 125 Convert from [{eV/c}] to $\beta_x\gamma$ [dimensionless],  snuverink_j committed Sep 06, 2017 126 127  \label{eq:eVtobetagamma}  snuverink_j committed Sep 08, 2017 128 (\beta_x\gamma)=\sqrt{(\frac{P_x[{eV/c}]}{m_0c}+1)^2-1}.  snuverink_j committed Sep 06, 2017 129 130 131 132 133 134 135  This may be deduced by analogy for the $y$ and $z$ directions. \subsection{The initial distribution in the local reference frame } \label{sec:opalcycl:localframe} The initial distribution of the bunch,  snuverink_j committed Sep 06, 2017 136 either read from file or generated by a distribution generator see~Chapter~\ref{distribution},  snuverink_j committed Sep 07, 2017 137 is specified in the local reference frame of the \textit{OPAL-cycl} Cartesian coordinate system see~Section~\ref{opalcycl:canon}.  snuverink_j committed Sep 06, 2017 138 At the beginning of the run, the 6 phase space variables $$(x, y, z, p_x, p_y, p_z)$$  snuverink_j committed Sep 06, 2017 139 140 are transformed to the global Cartesian coordinates using the starting coordinates $r_0$ (\texttt{RINIT}), $\phi_0$ (\texttt{PHIINIT}), and $z_0$ (\texttt{ZINIT}), and the starting momenta $p_{r0}$ (\texttt{PRINIT}), and $p_{z0}$ (\texttt{PZINIT}) of the reference particle, defined  snuverink_j committed Sep 06, 2017 141 in the \texttt{CYCLOTRON} element see~Section~\ref{cyclotron}. Note that $p_{\phi 0}$ is calculated  snuverink_j committed Sep 06, 2017 142 143 144 automatically from $p_{total}$, $p_{r0}$, and $p_{z0}$. \begin{align*}  snuverink_j committed Sep 11, 2017 145 146 X &= x\cos(\phi_0) - y\sin(\phi_0) + r_0\cos(\phi_0) \\ Y &= x\sin(\phi_0) + y\cos(\phi_0) + r_0\sin(\phi_0) \\  snuverink_j committed Sep 06, 2017 147 148 149 Z &= z + z_0 \end{align*} \begin{align*}  snuverink_j committed Sep 11, 2017 150 151 PX &= (p_x+p_{r0})\cos(\phi_0) - (p_y+p_{\phi 0})\sin(\phi_0) \\ PY &= (p_x+p_{r0})\sin(\phi_0) + (p_y+p_{\phi 0})\cos(\phi_0) \\  snuverink_j committed Sep 06, 2017 152 153 154 155 156 157 158 PZ &= p_z + p_{z0} \end{align*} % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Field Maps} \label{sec:opalcycl:fieldmap}  snuverink_j committed Sep 07, 2017 159 In \textit{OPAL-cycl}, the magnetic field on the median plane is read from an ASCII type file. The field data should be stored in the cylinder coordinates frame (because the field map on the median plane of the cyclotron is usually measured in this frame).  snuverink_j committed Sep 06, 2017 160 161 162 163 164  There are two possible situations. One is the real field map on median plane of the exist cyclotron machine using measurement equipment. Limited by the narrow gaps of magnets, in most cases with cyclotrons, only vertical field $B_z$ on the median plane ($z=0$) is measured. Since the magnetic field data off the median plane field components is necessary for those particles with $z \neq 0$, the field need to be expanded in $Z$ direction. According to the approach given by Gordon and Taivassalo, by using a magnetic potential and measured $B_z$ on the median plane,  snuverink_j committed Sep 11, 2017 165 at the point $(r,\theta, z)$ in cylindrical polar coordinates, the third order field can be written as  snuverink_j committed Sep 06, 2017 166 \begin{eqnarray}\label{eq:Bfield}  snuverink_j committed Sep 11, 2017 167 168  B_r(r,\theta, z) & = & z\frac{\partial B_z}{\partial r}-\frac{1}{6}z^3 C_r, \\ B_\theta(r,\theta, z) & = & \frac{z}{r}\frac{\partial B_z}{\partial \theta}-\frac{1}{6}\frac{z^3}{r} C_{\theta}, \\  snuverink_j committed Sep 11, 2017 169  B_z(r,\theta, z) & = & B_z-\frac{1}{2}z^2 C_z,  snuverink_j committed Sep 06, 2017 170 171 172 \end{eqnarray} where $B_z\equiv B_z(r, \theta, 0)$ and \begin{eqnarray}\label{eq:Bcoeff}  snuverink_j committed Sep 11, 2017 173 174 175 176 177  C_r & = & \frac{\partial^{3}{B_z}}{\partial r^{3}} + \frac{1}{r}\frac{\partial^{2}{B_z}}{\partial r^{2}} - \frac{1}{r^2}\frac{\partial{B_z}}{\partial r} + \frac{1}{r^2}\frac{\partial^{3}{B_z}}{{\partial r}{\partial \theta^2}} - 2\frac{1}{r^3}\frac{\partial^{2}{B_z}}{\partial \theta^{2}}, \\ C_{\theta} & = & \frac{1}{r}\frac{\partial^{2}{B_z}}{{\partial r}{\partial \theta}} + \frac{\partial^{3}{B_z}}{{\partial r^2}{\partial \theta}} + \frac{1}{r^2}\frac{\partial^{3}{B_z}}{\partial \theta^{3}}, \\ C_z & = & \frac{1}{r}\frac{\partial{B_z}}{\partial r} + \frac{\partial^{2}{B_z}}{\partial r^{2}} + \frac{1}{r^2}\frac{\partial^{2}{B_z}}{\partial \theta^{2}}.  snuverink_j committed Sep 06, 2017 178 179 \end{eqnarray}  snuverink_j committed Sep 07, 2017 180 All the partial differential coefficients are on the median plane and can be calculated by interpolation. \textit{OPAL-cycl} uses Lagrange's 5-point formula.  snuverink_j committed Sep 06, 2017 181 182 183 184 185 186 187  The other situation is to calculate the field on the median plane or the 3D fields of the working gap for interesting region numerically by creating 3D model using commercial software, such as TOSCA, ANSOFT and ANSYS during the design phase of a new machine. If the field on the median plane is calculated, the field off the median plane can be obtained using the same expansion approach as the measured field map as described above. However, the 3D fields of the entire working gap should be more accurate than the expansion approach especially at the region not so close to the median plane in $Z$ direction. In the current version, we implemented the three specific type field-read functions \emph{Cyclotron::getFieldFromFile()} of the median plane fields.  snuverink_j committed Sep 06, 2017 188 That which function is used is controlled by the parameters \texttt{TYPE} of \texttt{CYCLOTRON} see~Section~\ref{cyclotron} in the input file.  snuverink_j committed Sep 06, 2017 189 190  \subsection{CARBONCYCL type}  snuverink_j committed Sep 06, 2017 191 If \texttt{TYPE=CARBONCYCL}, the program requires the $B_z$ data which is stored in a sequence shown in Figure~\ref{CYCLField}.  snuverink_j committed Sep 06, 2017 192 193 194 195 196 197 198 \begin{figure}[ht] \begin{center} \includegraphics[origin=bl,height=40mm]{./figures/cyclotron/CarbonFieldFormat.pdf} \caption{2D field map on the median plane with primary direction corresponding to the azimuthal direction, secondary direction to the radial direction} \label{fig:CYCLField} \end{center} \end{figure}  snuverink_j committed Sep 11, 2017 199 We need to add 6 parameters at the header of a plain $B_z$ [{kG}] data file, namely,  snuverink_j committed Sep 11, 2017 200 $r_{min}$ [{mm}], $\Delta r$ [{mm}], $\theta_{min}$ [{^{\circ}}], $\Delta \theta$ [{^{\circ}}],  snuverink_j committed Sep 06, 2017 201 $N_\theta$ (total data number in each arc path of azimuthal direction) and $N_r$ (total path number along radial direction).  snuverink_j committed Sep 08, 2017 202 If $\Delta r$ or $\Delta \theta$ is decimal, one can set its negative opposite number. For instance, if $\Delta \theta = \frac{1}{3}{^{\circ}}$, the fourth line of the header should be set to -3.0.  snuverink_j committed Sep 06, 2017 203 Example showing the above explained format:  snuverink_j committed Sep 06, 2017 204 \begin{verbatim}  snuverink_j committed Sep 06, 2017 205 206 207 208 209 210 211 212 213 214 215 3.0e+03 10.0 0.0 -3.0 300 161 1.414e-03 3.743e-03 8.517e-03 1.221e-02 2.296e-02 3.884e-02 5.999e-02 8.580e-02 1.150e-01 1.461e-01 1.779e-01 2.090e-01 2.392e-01 2.682e-01 2.964e-01 3.245e-01 3.534e-01 3.843e-01 4.184e-01 4.573e-01 ......  snuverink_j committed Sep 06, 2017 216 \end{verbatim}  snuverink_j committed Sep 06, 2017 217 218 219  \subsection{CYCIAE type}  snuverink_j committed Sep 08, 2017 220 If \texttt{TYPE=CYCIAE}, the program requires data format given by ANSYS\,10.0. This function is originally for the {100}{MeV} cyclotron of CIAE,  snuverink_j committed Sep 06, 2017 221 222 223  whose isochronous fields is numerically computed by by ANSYS. The median plane fields data is output by reading the APDL (ANSYS Parametric Design Language) script during the post-processing phase (you may need to do minor changes to adapt your own cyclotron model):  snuverink_j committed Sep 06, 2017 224 \begin{verbatim}  snuverink_j committed Sep 06, 2017 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 /post1 resume,solu,db csys,1 nsel,s,loc,x,0 nsel,r,loc,y,0 nsel,r,loc,z,0 PRNSOL,B,COMP CSYS,1 rsys,1 *do,count,0,200 path,cyc100_Ansys,2,5,45 ppath,1,,0.01*count,0,,1 ppath,2,,0.01*count/sqrt(2),0.01*count/sqrt(2),,1 pdef,bz,b,z paget,data,table *if,count,eq,0,then /output,cyc100_ANSYS,dat *STATUS,data,,,5,5 /output *else /output,cyc100_ANSYS,dat,,append *STATUS,data,,,5,5 /output *endif *enddo finish  snuverink_j committed Sep 06, 2017 255 \end{verbatim}  snuverink_j committed Sep 06, 2017 256   snuverink_j committed Sep 08, 2017 257 By running this in ANSYS, you can get a fields file with the name \textit{cyc100\_ANSYS.data}.  snuverink_j committed Sep 06, 2017 258 You need to put 6 parameters at the header of the file, namely,  snuverink_j committed Sep 11, 2017 259 $r_{min}$ [{mm}], $\Delta r$ [{mm}], $\theta_{min}$[{^{\circ}}], $\Delta \theta$[{^{\circ}}],  snuverink_j committed Sep 06, 2017 260 261 $N_\theta$(total data number in each arc path of azimuthal direction) and $N_r$(total path number along radial direction). If $\Delta r$ or $\Delta \theta$ is decimal,one can set its negative opposite number. This is useful is the decimal is unlimited.  snuverink_j committed Sep 08, 2017 262 For instance,if $\Delta \theta = \frac{1}{3} {^{\circ}}$, the fourth line of the header should be -3.0.  snuverink_j committed Sep 06, 2017 263 264 In a word, the fields are stored in the following format:  snuverink_j committed Sep 06, 2017 265 \begin{verbatim}  snuverink_j committed Sep 06, 2017 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 0.0 10.0 0.0e+00 1.0e+00 90 201 PARAMETER STATUS- DATA ( 336 PARAMETERS DEFINED) (INCLUDING 17 INTERNAL PARAMETERS) LOCATION VALUE 1 5 1 0.537657876 2 5 1 0.538079473 3 5 1 0.539086731 ...... 44 5 1 0.760777286 45 5 1 0.760918663 46 5 1 0.760969074 PARAMETER STATUS- DATA ( 336 PARAMETERS DEFINED) (INCLUDING 17 INTERNAL PARAMETERS) LOCATION VALUE 1 5 1 0.704927299 2 5 1 0.705050993 3 5 1 0.705341275 ......  snuverink_j committed Sep 06, 2017 292 \end{verbatim}  snuverink_j committed Sep 06, 2017 293 294  \subsection{BANDRF type}  snuverink_j committed Sep 06, 2017 295 296 If \texttt{TYPE=BANDRF}, the program requires the $B_z$ data format which is same as \texttt{CARBONCYCL}. But with \texttt{BANDRF} type, the program can also read in the 3D electric field(s).  snuverink_j committed Sep 06, 2017 297 For the detail about its usage, please see Section~\ref{cyclotron}.  snuverink_j committed Sep 06, 2017 298 299  \subsection{Default PSI format}  snuverink_j committed Sep 08, 2017 300 If the value of \texttt{TYPE} is other string rather than above mentioned, the program requires the data format like PSI format field file \textit{ZYKL9Z.NAR} and \textit{SO3AV.NAR}, which was given by the measurement.  snuverink_j committed Sep 06, 2017 301 We add 4 parameters at the header of the file, namely,  snuverink_j committed Sep 11, 2017 302 $r_{min}$ [{mm}], $\Delta r$ [{mm}], $\theta_{min}$[{^{\circ}}], $\Delta \theta$[{^{\circ}}],  snuverink_j committed Sep 06, 2017 303 If $\Delta r$ or $\Delta \theta$ is decimal,one can set its negative opposite number. This is useful is the decimal is unlimited.  snuverink_j committed Sep 08, 2017 304 For instance,if $\Delta \theta = \frac{1}{3}{^{\circ}}$, the fourth line of the header should be -3.0.  snuverink_j committed Sep 06, 2017 305   snuverink_j committed Sep 06, 2017 306 \begin{verbatim}  snuverink_j committed Sep 06, 2017 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 1900.0 20.0 0.0 -3.0 LABEL=S03AV CFELD=FIELD NREC= 141 NPAR= 3 LPAR= 7 IENT= 1 IPAR= 1 3 141 135 30 8 8 70 LPAR= 1089 IENT= 2 IPAR= 2 0.100000000E+01 0.190000000E+04 0.200000000E+02 0.000000000E+00 0.333333343E+00 0.506500015E+02 0.600000000E+01 0.938255981E+03 0.100000000E+01 0.240956593E+01 0.282477260E+01 0.340503168E+01 0.419502926E+01 0.505867147E+01 0.550443363E+01 0.570645094E+01 0.579413509E+01 0.583940887E+01 0.586580372E+01 0.588523722E+01 ......  snuverink_j committed Sep 06, 2017 322 \end{verbatim}  snuverink_j committed Sep 06, 2017 323 324 325  \subsection{3D field-map}  snuverink_j committed Sep 06, 2017 326 % NOTE: \texttt{SBEND3D}, \texttt{RINGDEFINITION} in elements.tex and \ubsection {3D fieldmap} in  snuverink_j committed Sep 06, 2017 327 328 329 % opalcycl.tex all refer to each other - if updating one check for update on % others to keep them consistent.  snuverink_j committed Sep 07, 2017 330 It is additionally possible to load 3D field-maps for tracking through \textit{OPAL-cycl}.  snuverink_j committed Sep 06, 2017 331 3D field-maps are loaded by sequentially adding new field elements to a line,  snuverink_j committed Sep 07, 2017 332 as is done in \textit{OPAL-t}. It is not possible to add  snuverink_j committed Sep 06, 2017 333 RF cavities while operating in this mode. In order to define ring parameters  snuverink_j committed Sep 06, 2017 334 335 such as initial ring radius a \texttt{RINGDEFINITION} type is loaded into the line, followed by one or more \texttt{SBEND3D} elements.  snuverink_j committed Sep 06, 2017 336   snuverink_j committed Sep 06, 2017 337 \begin{verbatim}  snuverink_j committed Sep 06, 2017 338 339 340 341 342 343 344 ringdef: RINGDEFINITION, HARMONIC_NUMBER=6, LAT_RINIT=2350.0, LAT_PHIINIT=0.0, LAT_THETAINIT=0.0, BEAM_PHIINIT=0.0, BEAM_PRINIT=0.0, BEAM_RINIT=2266.0, SYMMETRY=4.0, RFFREQ=0.2; triplet: SBEND3D, FMAPFN="fdf-tosca-field-map.table", LENGTH_UNITS=10., FIELD_UNITS=-1e-4; l1: Line = (ringdef, triplet, triplet);  snuverink_j committed Sep 06, 2017 345 \end{verbatim}  snuverink_j committed Sep 06, 2017 346   snuverink_j committed Sep 08, 2017 347 The field-map with file name \textit{fdf-tosca-field-map.table} is loaded, which is a  snuverink_j committed Sep 06, 2017 348 file like  snuverink_j committed Sep 06, 2017 349 \begin{verbatim}  snuverink_j committed Sep 06, 2017 350 351 352 353 354 355 356 357 358 359 360 361  422280 422280 422280 1 1 X [LENGU] 2 Y [LENGU] 3 Z [LENGU] 4 BX [FLUXU] 5 BY [FLUXU] 6 BZ [FLUXU] 0 194.01470 0.0000000 80.363520 0.68275932346E-07 -5.3752492577 0.28280706805E-07 194.36351 0.0000000 79.516210 0.42525693524E-07 -5.3827955117 0.17681348191E-07 194.70861 0.0000000 78.667380 0.19766168358E-07 -5.4350026348 0.82540823165E-08  snuverink_j committed Sep 06, 2017 362 363 364 \end{verbatim} The header parameters are ignored - user supplied parameters \texttt{LENGTH\_UNITS} and \texttt{FIELD\_UNITS} are used. Fields are supplied on points in a grid in  snuverink_j committed Sep 06, 2017 365 366 367 368 369 370 371 $(r, y, \phi)$. Positions and field elements are specified by Cartesian coordinates $(x, y, z)$. \subsection{User's own field-map} { You should revise the function or write your own function according to the instructions in the code to match your own field format if it is different to above types.}  snuverink_j committed Sep 06, 2017 372 For more detail about the parameters of \texttt{CYCLOTRON}, please refer to Section~\ref{cyclotron}.  snuverink_j committed Sep 06, 2017 373 374 375 376 377 378 379 380  % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{RF field} \label{sec:opalcycl:rffieldmap} \subsection{Read RF voltage profile} The RF cavities are treated as straight lines with infinitely narrow gaps and the electric field is a $\delta$ function plus a transit time correction. the two-gap cavity is treated as two independent single-gap cavities. the spiral gap cavity is not implemented yet.  snuverink_j committed Sep 06, 2017 381 For more detail about the parameters of cyclotron cavities, see Section~\ref{cavity-cycl}.  snuverink_j committed Sep 06, 2017 382 383  The voltage profile of a cavity gap is read from ASCII file. Here is an example:  snuverink_j committed Sep 06, 2017 384 \begin{verbatim}  snuverink_j committed Sep 06, 2017 385 386 387 388 389 390 391 6 0.00 0.15 2.40 0.20 0.65 2.41 0.40 0.98 0.66 0.60 0.88 -1.59 0.80 0.43 -2.65 1.00 -0.05 -1.71  snuverink_j committed Sep 06, 2017 392 \end{verbatim}  snuverink_j committed Sep 06, 2017 393 394 395 396 The number in the first line means 6 sample points and in the following lines the three values represent the normalized distance to the inner edge of the cavity, the normalized voltage and its derivative respectively. \subsection{Read 3D RF field-map}  snuverink_j committed Sep 06, 2017 397 The 3D RF field-map can be read from H5hut type file. This is useful for modeling the central region electric fields which usually has complicate shapes. For the detail about its usage, please see Section~\ref{cyclotron}.  snuverink_j committed Sep 06, 2017 398   snuverink_j committed Sep 06, 2017 399 Please note that in this case, the E field is treated as a part of \texttt{CYCLOTRON} element, rather than a independent \texttt{RFCAVITY} element.  snuverink_j committed Sep 06, 2017 400 401 402 403 % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Particle Tracking and Acceleration} The precision of the tracking methods is vital for the entire simulation process, especially for long distance tracking jobs.  snuverink_j committed Sep 11, 2017 404 \textit{OPAL-cycl} uses a fourth order Runge-Kutta algorithm and the second order Leap-Frog scheme. The fourth order Runge-Kutta algorithm needs four external magnetic field evaluation in each time step $\tau$ .  snuverink_j committed Sep 06, 2017 405 During the field interpolation process, for an arbitrary given point the code first interpolates Formula $B_z$  snuverink_j committed Sep 06, 2017 406 for its counterpart on the median plane and then expands to this give point using Equation~\ref{Bfield}.  snuverink_j committed Sep 06, 2017 407 408 409 410 411 412 413 414 415 416 417  After each time step $i$, the code detects whether the particle crosses any one of the RF cavities during this step. If it does, the time point $t_c$ of crossing is calculated and the particle return back to the start point of step $i$. Then this step is divided into three sub-steps: first, the code tracks this particle for $t_1 = \tau - (t_c-t_{i-1})$; then it calculates the voltage and adds momentum kick to the particle and refreshes its relativistic parameters $\beta$ and $\gamma$; and then tracks it for $t_2 = \tau - t_1$. % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Space Charge}  snuverink_j committed Sep 07, 2017 418 \textit{OPAL-cycl} uses the same solvers as \textit{OPAL-t} to calculate the space charge effects see~Chapter~\ref{fieldsolver}.  snuverink_j committed Sep 06, 2017 419   snuverink_j committed Sep 11, 2017 420 Typically, the space charge field is calculated once per time step. This is no surprise for the second order Boris-Buneman time integrator (leapfrog-like scheme) which has per default only one force evaluation per step. The fourth order Runge-Kutta integrator keeps the space charge field constant for one step, although there are four external field evaluations. There is an experimental multiple-time-stepping (MTS) variant of the Boris-Buneman/leapfrog-method, which evaluates space charge only every N^{th} step, thus greatly reducing computation time while usually being still accurate enough.  snuverink_j committed Sep 06, 2017 421 422 423 424 425 426 427 428 429 430  % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Multi-bunch Mode} \label{sec:opalcycl:MultiBunch} The neighboring bunches problem is motivated by the fact that for high intensity cyclotrons with small turn separation, single bunch space charge effects are not the only contribution. Along with the increment of beam current, the mutual interaction of neighboring bunches in radial direction becomes more and more important, especially at large radius where the distances between neighboring bunches get increasingly small and even they  snuverink_j committed Sep 11, 2017 431 can overlap each other. One good example is PSI {590}{MeV} Ring cyclotron with a current of about {2}{mA} in  snuverink_j committed Sep 08, 2017 432 433 434 CW operation and the beam power amounts to {1.2}{MW}. An upgrade project for Ring is in process with the goal of {1.8}{MW} CW on target by replacing four old aluminum resonators by four new copper cavities with peak voltage increasing from about {0.7}{MW} to above {0.9}{MW}. After upgrade, the total turn  snuverink_j committed Sep 06, 2017 435 436 number is reduced from 200 turns to less than 170 turns. Turn separation is increased a little bit, but still are at the same order  snuverink_j committed Sep 11, 2017 437 of magnitude as the radial size of the bunches. Hence once the beam current increases from {2}{mA} to {3}{mA}, the mutual space  snuverink_j committed Sep 06, 2017 438 439 440 charge effects between radially neighboring bunches can have significant impact on beam dynamics. In consequence, it is important to cover neighboring bunch effects in the simulation to quantitatively study its impact on the beam dynamics.  snuverink_j committed Sep 07, 2017 441 In \textit{OPAL-cycl}, we developed a new fully consistent algorithm of multi-bunch simulation. We implemented two working modes, namely ,  snuverink_j committed Sep 06, 2017 442 \texttt{AUTO} and \texttt{FORCE}. In the first mode, only a single bunch is tracked and accelerated at the beginning,  snuverink_j committed Sep 06, 2017 443 444 445 446 447 448 449 450 451 452 453 454 455 until the radial neighboring bunch effects become an important factor to the bunches' behavior. Then the new bunches will be injected automatically to take these effects into account. In this way, we can save time and memory, and more important, we can get higher precision for the simulation in the region where neighboring bunch effects are important. In the other mode, multi-bunch simulation starts from the injection points. In the space charge calculation for multi-bunch, the computation region covers all bunches. Because the energy of the bunches is quite different, it is inappropriate to use only one particle rest frame and a single Lorentz transformation any more. So the particles are grouped into different energy bins and in each bin the energy spread is relatively small. We apply Lorentz transforming, calculate the space charge fields and apply the back Lorentz transforming for each bin separately. Then add the field data together. Each particle has a ID number to identify which energy bin it belongs to. % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Input}  snuverink_j committed Sep 07, 2017 456 All the three working modes of \textit{OPAL-cycl} use an input file written in \texttt{mad} language which will be described in detail in the following chapters.  snuverink_j committed Sep 06, 2017 457 458  For the {\bfseries Tune Calculation mode}, one additional auxiliary file with the following format is needed.  snuverink_j committed Sep 06, 2017 459 \begin{verbatim}  snuverink_j committed Sep 06, 2017 460 461 462 463 464 465 466  72.000 2131.4 -0.240 74.000 2155.1 -0.296 76.000 2179.7 -0.319 78.000 2204.7 -0.309 80.000 2229.6 -0.264 82.000 2254.0 -0.166 84.000 2278.0 0.025  snuverink_j committed Sep 06, 2017 467 \end{verbatim}  snuverink_j committed Sep 06, 2017 468 In each line the three values represent energy $E$, radius $r$ and $P_r$ for the SEO (Static Equilibrium Orbit)  snuverink_j committed Sep 11, 2017 469 at starting point respectively and their units are {MeV}, {mm} and {mrad}.  snuverink_j committed Sep 06, 2017 470   snuverink_j committed Sep 08, 2017 471 A bash script \textit{tuning.sh} is shown on the next page, to execute \textit{OPAL-cycl} for tune calculations.  snuverink_j committed Sep 06, 2017 472 \examplefromfile{examples/tuning.sh}  snuverink_j committed Sep 08, 2017 473 474 To start execution, just run \textit{tuning.sh} which uses the input file \textit{testcycl.in} and the auxiliary file \textit{FIXPO\_ SEO}. The output file is \textit{plotdata} from which one can plot the tune diagram.  snuverink_j committed Sep 06, 2017 475 476 477 478 479 480 481  % -- -- -- -- -- -- Section -- -- -- -- -- -- \section{Output} \subsubsection{Single Particle Tracking mode} The intermediate phase space data is stored in an ASCII file which can be used to  snuverink_j committed Sep 08, 2017 482 the plot the orbit. The file's name is combined by input file name (without extension) and \textit{-trackOrbit.dat}.  snuverink_j committed Sep 06, 2017 483 The data are stored in the global Cartesian coordinates.  snuverink_j committed Sep 06, 2017 484 The frequency of the data output can be set using the option \texttt{SPTDUMPFREQ} of \texttt{OPTION} statement see~Section~\ref{option}  snuverink_j committed Sep 06, 2017 485   snuverink_j committed Sep 06, 2017 486 The phase space data per \texttt{STEPSPERTURN} (a parameter in the \texttt{TRACK} command) steps is stored in an ASCII file.  snuverink_j committed Sep 08, 2017 487 The file's name is a combination of input file name (without extension) and \textit{-afterEachTurn.dat}.  snuverink_j committed Sep 06, 2017 488 The data is stored in the global cylindrical coordinate system.  snuverink_j committed Sep 06, 2017 489 Please note that if the field map is ideally isochronous, the reference particle of a given energy take exactly one revolution in \texttt{STEPPERTURN} steps;  snuverink_j committed Sep 08, 2017 490 Otherwise, the particle may not go through a full {360}{^{\circ}} in \texttt{STEPPERTURN} steps.  snuverink_j committed Sep 06, 2017 491 492  There are 3 ASCII files which store the phase space data around $0$, $\pi/8$ and $\pi/4$ azimuths.  snuverink_j committed Sep 08, 2017 493 Their names are the combinations of input file name (without extension) and \textit{-Angle0.dat}, \textit{-Angle1.dat} and \textit{-Angle2.dat} respectively.  snuverink_j committed Sep 06, 2017 494 495 496 497 The data is stored in the global cylindrical coordinate system, which can be used to check the property of the closed orbit. \subsubsection{Tune calculation mode}  snuverink_j committed Sep 08, 2017 498  The tunes $\nu_r$ and $\nu_z$ of each energy are stored in a ASCII file with name \textit{tuningresult}.  snuverink_j committed Sep 06, 2017 499 500 501 502 503  \subsubsection{Multi-particle tracking mode} The intermediate phase space data of all particles and some interesting parameters, including RMS envelop size, RMS emittance, external field, time, energy, length of path, number of bunches and  snuverink_j committed Sep 08, 2017 504 505 tracking step, are stored in the H5hut file-format \ref{bib:howison2010} and can be analyzed using the H5root \ref{bib:schietinger}.  snuverink_j committed Sep 06, 2017 506 The frequency of the data output can be set using the \texttt{PSDUMPFREQ} option of \texttt{OPTION} statement see~Section~\ref{option}.  snuverink_j committed Sep 08, 2017 507 The file is named like the input file but the extension is \textit{.h5}.  snuverink_j committed Sep 06, 2017 508 509  The intermediate phase space data of central particle (with ID of 0) and an off-centering particle (with ID of 1)  snuverink_j committed Sep 08, 2017 510 are stored in an ASCII file. The file's name is combined by the input file name (without extension) and \textit{-trackOrbit.dat}.  snuverink_j committed Sep 06, 2017 511 The frequency of the data output can be set using the \texttt{SPTDUMPFREQ} option of \texttt{OPTION} statement see~Section~\ref{option}.  snuverink_j committed Sep 06, 2017 512 513 514 515 516  \section{Matched Distribution} In order to run matched distribution simulation one has to specify a periodic accelerator. The function call also needs  snuverink_j committed Sep 11, 2017 517 the symmetry of the machine as well as a field map. The user then specifies the emittance $\pi{mmmrad}$.  snuverink_j committed Sep 06, 2017 518 519 520 521 522 523 524 525 526 527 528 529 530 \begin{verbatim} /* * specify periodic accelerator */ l1 = ... /* * try finding a matched distribution */ Dist1:DISTRIBUTION, TYPE=GAUSSMATCHED, LINE=l1, FMAPFN=..., MAGSYM=..., EX = ..., EY = ..., ET = ...; \end{verbatim} \subsection{Example}  snuverink_j committed Sep 11, 2017 531 Simulation of the PSI Ring Cyclotron at {580}{MeV} and current {2.2}{mA}. The program finds a matched  snuverink_j committed Sep 06, 2017 532 533 534 535 536 537 538 539 540 541 542 distribution followed by a bunch initialization according to the matched covariance matrix.\\ The matched distribution algorithm works with normalized emittances, i.e. normalized by the lowest energy of the machine. The printed emittances, however, are the geometric emittances. In addition, it has to be paid attention that the computation is based on $(x,x',y,y',z,\delta)$ instead of $(x,p_{x},y,p_{y},z,p_{z})$. Since the particles are represented in the latter coordinate system, the corresponding transformation has to be applied to obtain the rms emittances that are given in the output. \subsubsection{Input file} \examplefromfile{examples/matchedDistribution.in} \subsubsection{Output} \examplefromfile{examples/matchedDistribution.output}