### Luminosity leveling in python/MAD-X G. Sterbini Thanks to Axel, Ilias, Gianni, Sofia, Stefania and Yannis --- ### Resources - References - [M. A. Furman, 2003](https://www.osti.gov/servlets/purl/836235) - [W. Herr and B. Muratori, 2006](https://cds.cern.ch/record/941318) - [M. Sypher, 2007](https://home.fnal.gov/~syphers/Accelerators/tevPapers/LumiCalc.pdf) - [R. Calaga, 2012](https://espace.cern.ch/acc-tec-sector/Chamonix/Chamx2012/papers/RC_9_04.pdf) - [S. Fartoukh, 2014](https://journals.aps.org/prab/pdf/10.1103/PhysRevSTAB.17.111001) - [S. Papadopoulou et al., 2020](https://arxiv.org/pdf/1806.07317.pdf) - Documentation - [Luminosity and its python implementation](http://lhcmaskdoc.web.cern.ch/ipynbs/luminosity_formula/luminosity_formula/) - Example - [MAD-X example](http://lhcmaskdoc.web.cern.ch/ipynbs/python_example_simple/Pythonic_approach_simple/) --- ### Introduction * **AIM**: he aim is to have a quite general formula of the luminosity to be used in MAD-X for the luminosity leveling module. * **NEEDS**: maintain some intrinsic aspect: crossing angle, hourglass effect, crab cavities, offsets, flat-optics, B1/2 asymmetry... --- ### Luminosity \begin{equation} \mathcal{L}= n_b f_r N_1 N_2 \underbrace{\sqrt{(\vec{v_1}-\vec{v_2})^2-\frac{(\vec{v_1} \times \vec{v_2})^2}{c^2}}}_{Moeller~factor~\approx~2 c} \int_{-\infty}^{\infty} dt \int_{-\infty}^{\infty} \rho^{B1} \rho^{B2} dx dy dz \end{equation} ::: info $\mathcal{L}$ is a time-integral of the overlapping space-integral. ::: ---- #### Normalization \begin{equation} \int_{-\infty}^{\infty} \rho^{B} dx dy dz =1,\quad {\rm with }~B\in{B1,B2} \end{equation} ::: info $\rho_i^B$ is in [$\frac{1}{\rm m^3}$]$\rightarrow \mathcal{L}$ is in [$\frac{\rm Hz}{\rm m^2}$]. ::: ---- #### HP1: factorization \begin{equation} \rho^{B}=\prod_{i\in{x,y,z}} \rho_i^B,\quad {\rm with }~B\in{B1,B2} \end{equation} ::: danger We are not considering $xy$-tilted distribution (i.e., strong coupling). ::: ---- #### HP2: transverse dependence wrt z \begin{eqnarray} \rho_x^{B}&=&\rho_x^{B}(x,\color{red}{z},t),\quad &{\rm with }~B\in{B1,B2} \\ \rho_y^{B}&=&\rho_y^{B}(y,\color{red}{z},t),\quad &{\rm with }~B\in{B1,B2} \\ \rho_z^{B}&=&\rho_z^{B}(z,t),\quad &{\rm with }~B\in{B1,B2} \\ \end{eqnarray} ::: info Paraxial approximation. ::: ---- #### HP3: normal distribution \begin{equation} \rho_i^B= N_i^B(i, \mu_{i,B}, \sigma_{i,B})=\frac{1}{\sqrt{2\pi}\sigma_{i,B}} e^{-\frac{(i-\mu_{i,B})^2}{2\sigma_{i,B}^2}},\\ \quad {\rm with }~B\in{B1,B2}~{\rm and}~i\in{x,y,z} \end{equation} ::: info Normal distributions: it can be generalized (see later). ::: --- ### How? The problem to compute $\mathcal{L}$ therefore reduces - in integrating (analytically) the general Gaussian expression of the luminosity in $x$ and $y$. - in expressing the parameters $\mu_{x,y}$ and $\sigma_{x,y}$ as function of $z$ and $t$. - in integrating (numerically) in $z$ and $t$. --- ### That is \begin{equation} \mathcal{L} \approx \frac{n_b f_r N_1 N_2}{2 \pi}~ ~ 2 c \int \prod_i^{1,2}\frac{e^{- \frac{\left(z- \beta_{ri} c t -\mu_{zi}\right)^{2}}{2 \sigma_{zi}^{2}}}}{\sqrt{2 \pi}\sigma_{zi}} \prod_j^{x,y}\frac{e^{-\frac{\left(\color{red}{\mu_{j1}}- \color{red}{\mu_{j2}}\right)^2}{2\left(\color{red}{\sigma_{j1}}^{2} + \color{red}{\sigma_{j2}}^{2}\right)}}}{\sqrt{\color{red}{\sigma_{j1}}^{2} + \color{red}{\sigma_{j2}}^{2}}} dzdt \end{equation} ::: info Longitudinal distribution are not integrated. ::: :::info Let us focus only on $\color{red}{\mu}$'s and $\color{red}{\sigma}$'s of the $\color{red}{xy}$-plane. ::: --- ### $\color{red}{\mu}$'s and $\color{red}{\sigma}$'s of the $\color{red}{xy}$-plane - Hourglass effect $\rightarrow$ $\color{red}{\sigma}(z)$ - Dispersion + momentum spread $\rightarrow$ $\color{red}{\sigma}(z)$ - Crossing angle and/or separation $\rightarrow$ $\color{red}{\mu}(z)$ - Crab crossing $\rightarrow$ $\color{red}{\mu}(z,\color{green}{t})$ ---- #### Example I \begin{equation} \sigma_{x1}(z)=\sqrt{\frac{\epsilon_{n, x1}\beta_{x1}(\color{red}{z})}{\beta_{r1}\gamma_{r1}}+D_{x1}^2(\color{red}{z}) \left(\frac{\Delta p}{p_{0}}\right)_1^2} \end{equation} ::: info $\rightarrow$ usual quadratic sum of betatronic and dispersive contribution ::: ---- #### Example II \begin{equation} \mu_{x1}(z, t)= \mu_{x1}(0) + R_{12} \theta_{CC}(\color{red}{z}-c\color{red}{t})+ (\theta_{x1} + R_{22} \theta_{CC}(\color{red}{z}-c\color{red}{t})) \color{red}{z} \end{equation} ::: info $\alpha_{IP}\approx0$ and $\Delta \mu_{CC\rightarrow IP}\approx \frac{\pi}{2}$ $\rightarrow R_{22}\approx0$ ::: ::: info Full crabbing compensation for $t=0$ $\rightarrow$ $\boxed{R_{12} \theta_{CC}(z)\approx - \theta_{x1} z, ~ {\rm assuming} ~R_{22}\approx0 }$ ::: --- ### (Cheap) Generalization to $\mathcal{M}$-Gaussian distribution \begin{equation} \rho_{x}^{B1}=\sum_{i=1}^\mathcal{M} c_{x,i}^{B1} N(x,\mu_{x,i}^{B1}, \sigma_{x,i}^{B1})=\sum_{i=1}^\mathcal{M} c_{x,i}^{B1} N_{x,i}^{B1} \end{equation} \begin{equation} {\rm with}~\sum_{i=1}^\mathcal{M} c_{x,i}^{B1}=1. \end{equation} :::info Let's consider the $\rho$ of the $xy$-plane as a sum of Normal distributions. In the z-plane they are anyhow numerically integrated. ::: ---- #### Polynominals multiplication of $\mathcal{L}$ \begin{equation} \mathcal{L_M}=\sum_{i=1}^{\mathcal{M}} \sum_{j=1}^{\mathcal{M}} \sum_{k=1}^{\mathcal{M}} \sum_{l=1}^{\mathcal{M}} c_{x,i}^{B1} c_{y,j}^{B1} c_{x,k}^{B2} c_{y,l}^{B2}\ \color{red}{\mathcal{L}}(N_{x,i}^{B1}, N_{y,j}^{B1}, N_{x,k}^{B2}, N_{y,l}^{B2}). \end{equation} :::success Just a linear combination (4 nested ```for``` loops) of something that we know, $\color{red}{\mathcal{L}}$. ::: --- ### Questions?  ---  --- ### Python implementation - Detailed explanation on [Luminosity and its python implementation](http://lhcmaskdoc.web.cern.ch/ipynbs/luminosity_formula/luminosity_formula/) - Code on github: [luminosity.py](https://github.com/sterbini/madxp/blob/master/madxp/luminosity.py) ```python= def L(f, nb, N1, N2, x_1, x_2, y_1, y_2, px_1, px_2, py_1, py_2, energy_tot1, energy_tot2, deltap_p0_1, deltap_p0_2, epsilon_x1, epsilon_x2, epsilon_y1, epsilon_y2, sigma_z1, sigma_z2, beta_x1, beta_x2,beta_y1, beta_y2, alpha_x1, alpha_x2,alpha_y1, alpha_y2, dx_1, dx_2, dy_1, dy_2, dpx_1, dpx_2, dpy_1, dpy_2, CC_V_x_1=0, CC_f_x_1=0, CC_phase_x_1=0, CC_V_x_2=0, CC_f_x_2=0, CC_phase_x_2=0, CC_V_y_1=0, CC_f_y_1=0, CC_phase_y_1=0, CC_V_y_2=0, CC_f_y_2=0, CC_phase_y_2=0, R12_1=0, R22_1=0, R34_1=0, R44_1=0, R12_2=0, R22_2=0, R34_2=0, R44_2=0, verbose=False, sigma_integration=3) ``` --- ### What about the leveling? :::info You can build yourself a target function and zero it with the ```least_squares``` method. ::: ```python= L_target=2e+32 starting_guess=1e-5 def function_to_minimize(delta): aux=L(f=B1.freq0*1e6, nb=parameter_dict['par_nco_IP8'], N1=B1.npart, N2=B2.npart, energy_tot1=B1.energy, energy_tot2=B2.energy, deltap_p0_1=B1.sige, deltap_p0_2=B2.sige, epsilon_x1=B1.exn, epsilon_x2=B2.exn, epsilon_y1=B1.eyn, epsilon_y2=B2.eyn, sigma_z1=B1.sigt, sigma_z2=B2.sigt, beta_x1=B1_IP.betx, beta_x2=B2_IP.betx, beta_y1=B1_IP.bety, beta_y2=B2_IP.bety, alpha_x1=B1_IP.alfx, alpha_x2=B2_IP.alfx, alpha_y1=B1_IP.alfy, alpha_y2=B2_IP.alfy, dx_1=B1_IP.dx, dx_2=B2_IP.dx, dpx_1=B1_IP.dpx, dpx_2=B2_IP.dpx, dy_1=B1_IP.dy, dy_2=B2_IP.dy, dpy_1=B1_IP.dpy, dpy_2=B2_IP.dpy, x_1=B1_IP.x, x_2=B2_IP.x, px_1=B1_IP.px, px_2=B2_IP.px, y_1=delta, y_2=-delta, py_1=B1_IP.py, py_2=B2_IP.py, verbose=False) return aux-L_target aux=least_squares(function_to_minimize, starting_guess, bounds=(0, 5e-5)) ``` ---- #### A real-world example in MAD-X  Click [here](http://lhcmaskdoc.web.cern.ch/ipynbs/python_example_simple/Pythonic_approach_simple/). --- ### Todo - Checks, checks, checks,... - Systematic checks of CC - Taking into account $\mu_{z1}$ and $\mu_{z2}$ (RF beam loading neglected) - Generalization with q-Gaussian Do not hesitate to leave your issues/suggestions on the [repository](https://github.com/sterbini/madxp). --- ### Thank you.