Hide code cell source
from map2D import map2D
import numpy as np
import matplotlib.pyplot as plt
#from matplotlib import animation
#from IPython.display import HTML
%matplotlib widget
plt.rcParams.update({'font.size': 12})

4. Perturbation of transverse motion#

The Hamiltonian of an unperturbed motion is

(4.1)#\[\begin{align} H=\frac{p_x^2}{2}+\frac{p_y^2}{2}+\left(\frac{1}{2\rho^2}+\frac{k(s)}{2} \right) x^2-\frac{k(s)}{2}y^2 \end{align}\]

Now we assume that the motion is perturbed by \(\Delta \mathbf{B}(s) = \Delta B_x(s) \hat{x}+\Delta B_y(s) \hat{y}\), where

(4.2)#\[\begin{align} \Delta B_y+i\Delta B_x &= B_0 \sum_{n=0}^\infty \left(\Delta b_n + i \Delta a_n\right) \left(x+iy\right)^{n} \end{align}\]

The vector potential becomes:

(4.3)#\[\begin{align} \Delta A_s=-B_0 \Re \sum_{n=0}^\infty \frac {\Delta b_n + i \Delta a_n}{n+1} \left(x+iy\right)^{n+1} \end{align}\]

The Hamiltonian becomes:

(4.4)#\[\begin{align} H=\frac{p_x^2}{2}+\frac{p_y^2}{2}+\left(\frac{1}{2\rho^2}+\frac{k(s)}{2} \right) x^2-\frac{k(s)}{2}y^2 - \frac{1}{B\rho} \Delta A_s \end{align}\]

Then the Hill’s equation becomes:

(4.5)#\[\begin{align} x''+\left(\frac{1}{\rho^2}+k(s)\right)x&= -\frac{\Delta B_y}{B\rho} \\ y''-k(s)y&=\frac{\Delta B_x}{B\rho} \end{align}\]

If \(\Delta B_y\) is not a function of \(y\), and \(\Delta B_x\) is not a function of \(x\), these two equation remain decoupled. The linear effect of uncoupled case will be introduce in the next two sub-sections. The linear coupled case will be introduced in the last sub-section.

4.1. Perturbed orbit#

The perturbed orbit rises from the dipole errors, which can be contributed from a dipole strength error or a quadrupole misalignment.

4.1.1. Green function of Hill’s equation#

At the present of the dipole error, the general Hill’s equation becomes:

(4.6)#\[\begin{equation} x''\left(s\right)+\kappa\left(s\right)x\left(s\right)=-\frac{\Delta B_y(s)}{B_0\rho} \end{equation}\]

The solution can be retrieved from the Green’s function method. The green function is the solution of:

(4.7)#\[\begin{equation} x''\left(s\right)+\kappa\left(s\right)x\left(s\right)=\theta\delta\left(s-s_0\right) \end{equation}\]

This shows that a delta function dipole changes the orbit by \(\theta\) at location \(s_0\)

Integrate the function at the vicinity of the \(s_0\) gives:

(4.8)#\[\begin{equation} x'\left(s_{0+}\right)-x'\left(s_{0-}\right)=\theta %%\label{eq:angle_kick} \end{equation}\]

4.1.1.1. Linac accelerator#

From the linear accelerator, the beam only pass \(s_0\) once. From the transfer matrix from \(s_0\) to \(s_1\)

(4.9)#\[\begin{split}\begin{equation} M\left(s_{1}\mid s_{0}\right) =\left(\begin{array}{cc} \sqrt{{\frac{\beta_{1}}{\beta_{0}}}}\left(\cos\psi+\alpha_{0}\sin\psi\right) & \sqrt{\beta_{0}\beta_{1}}\sin\psi\\ -\frac{1+\alpha_{0}\alpha_{1}}{\sqrt{\beta_{0}\beta_{1}}}\sin\psi+\frac{\alpha_{0}-\alpha_{1}}{\sqrt{\beta_{0}\beta_{1}}}\cos\psi & \sqrt{{\frac{\beta_{0}}{\beta_{1}}}}\left(\cos\psi-\alpha_{1}\sin\psi\right) \end{array}\right) %%\label{eq:tmatrix_s0s1} \end{equation}\end{split}\]

where \(\psi(s)=\psi(s_1)-\psi(s_0)\). We get the response of the delta dipole error, or the Green’s function (normalized by \(\theta\)) as:

(4.10)#\[\begin{equation} G(s, s_0)=\frac{x(s)}{\theta}=\sqrt{\beta_{s}\beta_{s_0}}\sin\left(\psi\left(s\right)-\psi\left(s_0\right)\right) H\left(s-s_0\right) \end{equation}\]

Here, \(H(x)\) is the Heaviside step function.

Then for an arbitrary dipole error distribution \(\frac{\Delta B(s)}{B_0\rho}\), the centroid orbit becomes:

(4.11)#\[\begin{align} x(s) &=\int_{-\infty}^{\infty}\frac{\Delta B(s_0)}{B_0\rho}G(s,s_0)ds_0 \nonumber\\ &=\sqrt{\beta(s)}\int_{-\infty}^{s}\frac{\sqrt{\beta(s_0)}\Delta B(s_0)}{B_0\rho}\sin\left(\psi\left(s\right)-\psi\left(s_0\right)\right) ds_0 \end{align}\]

4.1.1.2. Ring accelerator#

In ring accelerator, the beam will pass \(s_0\) many times. There exist a stable orbit, deviate from the reference orbit, due to the dipole error. This is called closed orbit.

In every turn, the beam will have a kick at \(s_0\) described in (4.8). Combined with the one-turn matrix \(M(s_0)\),

(4.12)#\[\begin{equation} M\left(s_{0}\right)= \left(\begin{array}{cc} \cos\Phi+\alpha_{0}\sin\Phi & \beta_{0}\sin\Phi\\ -\frac{1+\alpha_{0}^{2}}{\beta_{0}}\sin\Phi & \cos\Phi-\alpha_{0}\sin\Phi \end{array}\right) \end{equation}\]

Here \(\Phi=2\pi\nu\). we have a closed orbit condition:

(4.13)#\[\begin{equation} \left(\begin{array}{c} x\\ x'_{-}=x'_{+}-\theta \end{array}\right)= \left(\begin{array}{cc} \cos\Phi+\alpha_{0}\sin\Phi & \beta_{0}\sin\Phi\\ -\frac{1+\alpha_{0}^{2}}{\beta_{0}}\sin\Phi & \cos\Phi-\alpha_{0}\sin\Phi \end{array}\right) \left(\begin{array}{c} x\\ x'_{+} \end{array}\right) \end{equation}\]

Solve the equation, we have:

(4.14)#\[\begin{equation} \left(I-M\right) \left(\begin{array}{c} x\\ x'_{+} \end{array}\right) = \left(\begin{array}{c} 0\\ \theta \end{array}\right) \end{equation}\]

We know that \(I-M\) is invertabile when \(\det(I-M)=2-2\cos\Phi \ne 0\). Then the closed orbit at \(s_0\) can be easily found as:

(4.15)#\[\begin{align} x(s_0)&=\frac{\beta\left(s_0\right)\sin\Phi}{2-2\cos\Phi}\theta=\frac{\beta\left(s_0\right)\cos\pi\nu}{2\sin\pi\nu}\theta \\ x'(s_{0+})&=\frac{1-\cos\Phi-\alpha\sin\Phi}{2-2\cos\Phi}\theta=\frac{\sin\pi\nu-\alpha\left(s_0\right)\cos\pi\nu}{2\sin\pi\nu}\theta \end{align}\]

Here, the betatron tune is \(\nu\). Then using the transfer matrix (4.9), we get the Green’s function:

(4.16)#\[\begin{align} G\left(s,s_{0}\right) & =\frac{x(s)}{\theta}\nonumber\\ & =\frac{\sqrt{\beta\left(s\right)\beta\left(s_{0}\right)}}{{2\sin\pi\nu}}\left(\cos\pi\nu\left(\cos\psi+\alpha\left(s_{0}\right)\sin\psi\right)+\sin\psi\left(\sin\pi\nu-\alpha\left(s_{0}\right)\cos\pi\nu\right)\right)\nonumber\\ & =\frac{\sqrt{\beta\left(s\right)\beta\left(s_{0}\right)}}{{2\sin\pi\nu}}\cos\left(\psi\left(s\right)-\psi\left(s_{0}\right)-\pi\nu\right) \end{align}\]

Therefore, for an arbitrary dipole error distribution \(\frac{\Delta B(s)}{B_0\rho}\), the closed orbit is simply:

(4.17)#\[\begin{align} x_{co}(s)&=\int_{s}^{s+C}\frac{\Delta B(s_0)}{B_0\rho}G\left(s,s_{0}\right) ds_0 \\ &=\frac{\sqrt{\beta\left(s\right)}}{2\sin\pi\nu}\int_{s}^{s+C}\sqrt{\beta\left(s_0\right)}\frac{\Delta B(s_0)}{B_0\rho} \cos\left(\psi\left(s\right)-\psi\left(s_{0}\right)-\pi\nu\right) ds_0 \end{align}\]

Let’s define the ‘rotating angle’ to \(\phi(s)=\psi(s)/\nu\), then

(4.18)#\[\begin{align} d\phi=\frac{1}{\nu}\frac{d\psi}{ds}ds=\frac{1}{\nu \beta}ds \end{align} \]

Therefore:

(4.19)#\[\begin{align} x_{co}(s)&=\frac{\sqrt{\beta\left(s\right)}}{2\sin(\pi\nu)}\int_{\phi(s)}^{\phi(s)+2\pi}\nu\sqrt{\beta\left(\phi\right)^3}\frac{\Delta B(\phi)}{B_0\rho} \cos\nu\left(\phi\left(s\right)-\phi-\pi\right) d\phi \end{align}\]

The advantage of writing the integral using the ‘rotating angle’ is to perform Fourier analysis easily. We can expand:

(4.20)#\[\begin{align} \nu\beta^{3/2}(\phi)\frac{\Delta B(\phi)}{B_0\rho}=\sum_{m=-\infty}^{\infty}f_m e^{im\phi} \end{align}\]

Then the closed orbit is written as:

(4.21)#\[\begin{align} x_{co}&=\sqrt{\beta\left(s\right)}\sum_{m=-\infty}^{\infty} \frac{\nu f_m}{\nu^2-m^2} e^{im\phi(s)} \end{align}\]

Note

We can see clearly that only the harmonic of the error close to the betatron tune will contribute the most on the closed orbit. Therefore, integer tune is not stable under the dipole errors. The above expression can be used to evaluate the stop-band of the integer tune.

We can see clearly that only the harmonic of the error close to the betatron tune will contribute the most on the closed orbit. The following picture shows the comparision of integer and half-integer tune under dipole error. integer_tune

Example: Dipole errors with one-turn map.

A transverse kick is applied after the one turn map. We may observe the effects of:

  • Dependence on tune

  • Conservation of beam size and emittance

  • Effect of energy deviation, will be detaled later.

dipole_error=map2D(npart=10000, twiss=[2,1], twiss_beam=[2,1],tune=0.02, chrom=3, espr=0e-4)
avex,avep,sizex,sizep,emit=dipole_error.statistics()
emitlist=[]
sizelist=[]
avelist=[]
N_turn=10000


def evolve_func(turns, kick_turn_start=2000, kick_angle=5e-4,
               ):
    for i in range(turns):
        if i>=kick_turn_start:
            dipole_error.coor2D[1,:]+=kick_angle
            
        
        dipole_error.propagate()
        
        avex,avep,sizex,sizep,emit=dipole_error.statistics()
        avelist.append(avex)
        sizelist.append(sizex)
        emitlist.append(emit)
        yield dipole_error.coor2D
        
evolve=evolve_func(N_turn+2)
    
#fig,ax=plt.subplots()
for i in range(N_turn):
    arr=next(evolve)
    
    
    
fig1,(ax1,ax2,ax3)=plt.subplots(3,1, sharex= True)
ax3.set_xlabel("Turns")
ax1.set_ylabel("Centroid [mm]")
ax2.set_ylabel("Beam size [mm]")
ax3.set_ylabel("Emit. [mm mrad]")
ax1.plot(np.array(avelist)*1e3)  
ax1.set_ylim([-10,10])
ax2.plot(np.array(sizelist)*1e3) 
ax2.set_ylim([0,3])
ax3.set_ylim([0,4])
ax3.plot(np.array(emitlist)*1e6) 
fig1.tight_layout()

4.2. Perturbed Tune#

The gradient error of the quadrupole contribute to the perturbation of the betatron tune. Similar to the dipole errors, we start from the delta function response:

(4.22)#\[\begin{equation} x''\left(s\right)+\kappa\left(s\right)x\left(s\right)=-(kds)\delta(s-s_0)x(s) \end{equation}\]

We can easily find the tune change from the perturbed matrix as

(4.23)#\[\begin{align} M_{perturbed}&= M\left(s_{0}\right) \left(\begin{array}{cc} 1 & 0\\ -kds & 1 \end{array}\right)\nonumber\\ &= \left(\begin{array}{cc} \cos\Phi+\alpha_{0}\sin\Phi & \beta_{0}\sin\Phi\\ -\frac{1+\alpha_{0}^{2}}{\beta_{0}}\sin\Phi & \cos\Phi-\alpha_{0}\sin\Phi \end{array}\right) \left(\begin{array}{cc} 1 & 0\\ -kds & 1 \end{array}\right) \nonumber\\ &=\left(\begin{array}{cc} \cos\Phi+\alpha_{0}\sin\Phi-kds\beta_{0}\sin\Phi & \beta_{0}\sin\Phi\\ -\frac{1+\alpha_{0}^{2}}{\beta_{0}}\sin\Phi-kds(\cos\Phi-\alpha_{0}\sin\Phi)& \cos\Phi-\alpha_{0}\sin\Phi \end{array}\right) \end{align}\]

Treat the tune change to be small, therefore the phase advance change is given by

(4.24)#\[\begin{align} \Delta \Phi &\approx-\frac{\cos\Phi_n-\cos\Phi}{\sin\Phi} \nonumber\\ &=\frac{1}{2}\beta_0 k ds \end{align}\]

Note

Note that the phase advance change is actually nonliear, especially when the perturned tune is approaching 0 or 0.5 as shown below.

Hide code cell source
tunelist=[0.05,0.1,0.2,0.3, 0.4, 0.45]
kdsbeta=np.linspace(0,3,100)
fig, ax=plt.subplots(figsize=(8,6))
for tune in tunelist:
    trace = 2.0 * np.cos(2.0 * np.pi * tune) - kdsbeta * np.sin(2.0 * np.pi * tune)
    mask= np.abs(trace/2.0)<=1
    tuneshift = np.arccos(trace[mask]/2.0) / (2.0 * np.pi) - tune
    line=ax.plot(kdsbeta[mask], tuneshift, label=f"Unperturned Tune={tune}")
    trace = 2.0 * np.cos(2.0 * np.pi * tune) + kdsbeta * np.sin(2.0 * np.pi * tune)
    mask= np.abs(trace/2.0)<=1
    tuneshift = np.arccos(trace[mask]/2.0) / (2.0 * np.pi) - tune
    ax.plot(kdsbeta[mask], tuneshift, color=line[0].get_color(), linestyle="--", )

ax.set_xlabel(r"Kick strength $\beta k ds$")
ax.set_ylabel(r"Tune shift $\Delta \nu $")
ax.set_xlim([0,3])
ax.set_ylim([-0.2,0.2])
l=ax.legend()

Therefore the tune shift of entire ring is just:

(4.25)#\[\begin{align} \Delta \nu= \frac{1}{4\pi}\oint \beta(s)k(s)ds \end{align}\]

Consider the location \(s_1\), the one turn matrix can be found as:

(4.26)#\[\begin{align} M(s_1)=M(s_1,s_0)M(s|_0)M(s_0,s_1) \end{align}\]

The change of \(m_{12}\) element of \(M(s_1)\) is:

(4.27)#\[\begin{align} \Delta m_{12}=-kds\beta_0\beta_1\sin(\psi_0-\psi_1) \sin(2\pi\nu_0+\psi_1-\psi_0) \end{align}\]

The change of the beta function at \(s_1\) is simply:

(4.28)#\[\begin{align} \Delta \beta_1&=\frac{1}{\sin\Phi}\left(\Delta m_{12}-\beta_1\cos\Phi\Delta\Phi \right) \end{align}\]

Therefore the beta function error gives:

(4.29)#\[\begin{align} \frac{\Delta \beta_1}{\beta_1}&=\frac{-1}{2\sin\Phi}\int_{s}^{s+C} k(s_0)\beta(s_0) \cos\left(\psi\left(s\right)-\psi\left(s_{0}\right)+2\pi\nu\right) ds_0 \end{align}\]

Follow the same procedure as in the orbit error with transformation (4.18), we can define the Fourier sequence of the quad errors:

(4.30)#\[\begin{align} \nu\beta^{2}(\phi)k(\phi)=\sum_{m=-\infty}^{\infty}f_m e^{im\phi} \end{align}\]

we have beta beat as:

(4.31)#\[\begin{align} \frac{\Delta \beta}{\beta}=-\frac{1}{2}\sum_{m=-\infty}^{\infty} \frac{\nu f_m}{\nu^2-(m/2)^2} e^{im\phi(s)} \end{align}\]

Note

Half integer tune is not stable under the quadrupole errors. The above expression can be used to evaluate the stop-band of the half integer tune.

Example: Quad errors with one-turn map.

A quadrupole error is applied after the one turn map.

Hide code cell source
quad_error=map2D(npart=10000, twiss=[1,1], twiss_beam=[1,1],tune=0.525, chrom=3.0, espr=0.0e-4)
avex,avep,sizex,sizep,emit=quad_error.statistics()
emitlist=[]
sizelist=[]
avelist=[]
N_turn=6000


def evolve_func(turns, kick_turn_start=1000, df=4*np.pi*1.0e-2):
    for i in range(turns):
        if i>=kick_turn_start:
            quad_error.coor2D[1,:]+=df*quad_error.coor2D[0,:]
            
            
        
        quad_error.propagate()
        
        avex,avep,sizex,sizep,emit=quad_error.statistics()
        avelist.append(avex)
        sizelist.append(sizex)
        emitlist.append(emit)
        yield quad_error.coor2D
        
evolve=evolve_func(N_turn+2)
    
for i in range(N_turn):
    arr=next(evolve)
    
fig1,(ax1,ax2,ax3)=plt.subplots(3,1, sharex= True, figsize=(8,10))
ax3.set_xlabel("Turns")
ax1.set_ylabel("Centroid [mm]")
ax2.set_ylabel("Beam size [mm]")
ax3.set_ylabel("Emit. [mm mrad]")
ax1.plot(np.array(avelist)*1e3)  
ax1.set_ylim([-2,2])
ax2.plot(np.array(sizelist)*1e3) 

ax2.set_ylim([0,10])
ax3.set_ylim([0,2])
ax3.plot(np.array(emitlist)*1e6) 
fig1.tight_layout()