impedances Package¶
impedance
Module¶
Module to compute intensity effects
- Authors
Juan F. Esteban Mueller, Danilo Quartullo, Alexandre Lasheen, Markus Schwarz
-
class
blond.impedances.impedance.
InducedVoltageFreq
(Beam, Profile, impedance_source_list, frequency_resolution=None, multi_turn_wake=False, front_wake_length=0, RFParams=None, mtw_mode=None)¶ Bases:
blond.impedances.impedance._InducedVoltage
Induced voltage derived from the sum of several impedances
- Parameters
Beam (object) – Beam object
Profile (object) – Profile object
impedance_source_list (list) – Impedance sources list (e.g. list of Resonator objects)
frequency_resolution (float, optional) – Frequency resolution of the impedance [Hz]
multi_turn_wake (boolean, optional) – Multi-turn wake enable flag
front_wake_length (float, optional) – Lenght [s] of the front wake (if any) for multi-turn wake mode
RFParams (object, optional) – RFStation object for turn counter and revolution period
mtw_mode (boolean, optional) – Multi-turn wake mode can be ‘freq’ or ‘time’ (default)
-
impedance_source_list
¶ Impedance sources list (e.g. list of Resonator objects)
- Type
list
-
total_impedance
¶ Total impedance array of all sources in* \(\Omega\)
- Type
float array
-
front_wake_length
¶ Lenght [s] of the front wake (if any) for multi-turn wake mode
- Type
float
-
process
()¶ Reprocess the impedance contributions. To be run when profile change
-
sum_impedances
(freq)¶ Summing all the wake contributions in one total impedance.
-
class
blond.impedances.impedance.
InducedVoltageResonator
(Beam, Profile, Resonators, timeArray=None)¶ Bases:
blond.impedances.impedance._InducedVoltage
Calculates the induced voltage of several resonators for arbitrary line density. It does so by linearily interpolating the line density and solving the convolution integral with the resonator impedance analytically. The line density need NOT be sampled at equidistant points. The times where the induced voltage is calculated need to be the same where the line density is sampled. If no timeArray is passed, the induced voltage is evaluated at the points of the line density. This is nececassry of compatability with other functions that calculate the induced voltage. Currently, it requires the all quality factors :math:`Q>0.5` Currently, only works for single turn.
- Parameters
Beam (object) – Beam object
Profile (object) – Profile object
Resonators (object) – Resonators object
timeArray (float array, optional) – Array of time values where the induced voltage is calculated. If left out, the induced voltage is calculated at the times of the line density.
-
beam
¶ Copy of the Beam object in order to access the beam info.
- Type
object
-
profile
¶ Copy of the Profile object in order to access the line density.
- Type
object
-
tArray
¶ array of time values where the induced voltage is calculated. If left out, the induced voltage is calculated at the times of the line density
- Type
float array
-
atLineDensityTimes
¶ flag indicating if the induced voltage has to be computed for timeArray or for the line density
- Type
boolean
-
n_time
¶ length of tArray
- Type
int
-
R, omega_r, Q
Resonators parameters
- Type
lists of float
-
n_resonators
¶ Number of resonators
- Type
int
-
induced_voltage
¶ Computed induced voltage [V]
- Type
float array
-
Heaviside
(x)¶ Heaviside function, which returns 1 if x>1, 0 if x<0, and 1/2 if x=0
-
induced_voltage_1turn
(beam_spectrum_dict={})¶ Method to calculate the induced voltage through linearily interpolating the line density and applying the analytic equation to the result.
-
process
()¶ Reprocess the impedance contributions. To be run when slicing changes
-
class
blond.impedances.impedance.
InducedVoltageTime
(Beam, Profile, wake_source_list, wake_length=None, multi_turn_wake=False, RFParams=None, mtw_mode=None)¶ Bases:
blond.impedances.impedance._InducedVoltage
Induced voltage derived from the sum of several wake fields (time domain)
- Parameters
Beam (object) – Beam object
Profile (object) – Profile object
wake_source_list (list) – Wake sources list (e.g. list of Resonator objects)
wake_length (float, optional) – Wake length [s]
multi_turn_wake (boolean, optional) – Multi-turn wake enable flag
RFParams (object, optional) – RFStation object for turn counter and revolution period
mtw_mode (boolean, optional) – Multi-turn wake mode can be ‘freq’ or ‘time’ (default)
-
wake_source_list
¶ Wake sources list (e.g. list of Resonator objects)
- Type
list
-
total_wake
¶ Total wake array of all sources in \(\Omega / s\)
- Type
float array
-
process
()¶ Reprocess the impedance contributions. To be run when profile changes
-
sum_wakes
(time_array)¶ Summing all the wake contributions in one total wake.
-
class
blond.impedances.impedance.
InductiveImpedance
(Beam, Profile, Z_over_n, RFParams, deriv_mode='gradient')¶ Bases:
blond.impedances.impedance._InducedVoltage
Constant imaginary Z/n impedance
- Parameters
Beam (object) – Beam object
Profile (object) – Profile object
Z_over_n (float array) – Constant imaginary Z/n program in* \(\Omega\).
RFParams (object) – RFStation object for turn counter and revolution period
deriv_mode (string, optional) – Derivation method to compute induced voltage
-
Z_over_n
¶ Constant imaginary Z/n program in* \(\Omega\).
- Type
float array
-
deriv_mode
¶ Derivation method to compute induced voltage
- Type
string, optional
-
induced_voltage_1turn
(beam_spectrum_dict={})¶ Method to calculate the induced voltage through the derivative of the profile. The impedance must be a constant Z/n.
-
class
blond.impedances.impedance.
TotalInducedVoltage
(Beam, Profile, induced_voltage_list)¶ Bases:
object
Object gathering all the induced voltage contributions. The input is a list of objects able to compute induced voltages (InducedVoltageTime, InducedVoltageFreq, InductiveImpedance). All the induced voltages will be summed in order to reduce the computing time. All the induced voltages should have the same slicing resolution.
- Parameters
Beam (object) – Beam object
Profile (object) – Profile object
induced_voltage_list (object list) – List of objects for which induced voltages have to be calculated
-
beam
¶ Copy of the Beam object in order to access the beam info
- Type
object
-
profile
¶ Copy of the Profile object in order to access the profile info
- Type
object
-
induced_voltage_list
¶ List of objects for which induced voltages have to be calculated
- Type
object list
-
induced_voltage
¶ Array to store the computed induced voltage [V]
- Type
float array
-
time_array
¶ Time array corresponding to induced_voltage [s]
- Type
float array
-
induced_voltage_sum
()¶ Method to sum all the induced voltages in one single array.
-
induced_voltage_sum_packed
()¶ Method to sum all the induced voltages in one single array.
-
reprocess
()¶ Reprocess the impedance contributions. To be run when profile changes
-
track
()¶ Track method to apply the induced voltage kick on the beam.
-
track_ghosts_particles
(ghostBeam)¶
impedance_sources
Module¶
Module to describe classes for the calculation of wakes and impedances to be used either alone or by InducedVoltage objects described in impedance.py. The module consists of a parent class called _ImpedanceObject and several child classes, as for example InputTable, Resonators and TravelingWaveCavity.
- Authors
Danilo Quartullo, Alexandre Lasheen,
Juan F. Esteban Mueller
-
class
blond.impedances.impedance_sources.
InputTable
(input_1, input_2, input_3=None)¶ Bases:
blond.impedances.impedance_sources._ImpedanceObject
Intensity effects from impedance and wake tables. If the constructor takes just two arguments, then a wake table is passed; if it takes three arguments, then an impedance table is passed. Be careful that if you input a wake, the input wake for W(t=0) should be already divided by two (beam loading theorem) ; and that if you input impedance, only the positive frequencies of the impedance is needed (the impedance will be assumed to be Hermitian (Real part symmetric and Imaginary part antisymmetric).Note that we add the point (f, Z(f)) = (0, 0) to the frequency and impedance arrays derived from the table.
- Parameters
input_1 (float array) – Time array of the wake in s or frequency array of the impedance in Hz
input_2 (float array) – Wake array in \(\Omega / s\) or real part of impedance in \(\Omega\)
input_3 (float array) – Imaginary part of impedance in \(\Omega\)
-
time_array
¶ Input time array of the wake in s
- Type
float array
-
wake_array
¶ Input wake array in \(\Omega / s\)
- Type
float array
-
frequency_array_loaded
¶ Input frequency array of the impedance in Hz
- Type
float array
-
Re_Z_array_loaded
¶ Input real part of impedance in \(\Omega\)
- Type
float array
-
Im_Z_array_loaded
¶ Input imaginary part of impedance in \(\Omega\)
- Type
float array
-
impedance_loaded
¶ Input impedance array in \(\Omega + j \Omega\)
- Type
complex array
Examples
>>> ##### WAKEFIELD TABLE >>> time = np.array([1.1,2.2,3.3]) >>> wake = np.array([4.4,5.5,6.6]) >>> table = InputTable(time, wake) >>> new_time_array = np.array([7.7,8.8,9.9]) >>> table.wake_calc(new_time_array) >>> >>> ##### IMPEDANCE TABLE >>> frequency = np.array([1.1,2.2,3.3]) >>> real_part = np.array([4.4,5.5,6.6]) >>> imaginary_part = np.array([7.7,8.8,9.9]) >>> table2 = InputTable(frequency, real_part, imaginary_part) >>> new_freq_array = np.array([7.7,8.8,9.9]) >>> table2.imped_calc(new_freq_array)
-
imped_calc
(new_frequency_array)¶ The impedance from the table is interpolated using the new frequency array.
- Parameters
new_frequency_array (float array) – frequency array in \(\Omega\)
-
frequency_array
¶ frequency array in \(\Omega\)
- Type
float array
-
Re_Z_array
¶ Real part of the interpolated impedance in \(\Omega\)
- Type
float array
-
Im_Z_array
¶ Imaginary part of the interpolated impedance in \(\Omega\)
- Type
float array
-
impedance
¶ Output interpolated impedance array in \(\Omega + j \Omega\)
- Type
complex array
-
class
blond.impedances.impedance_sources.
ResistiveWall
(pipe_radius, pipe_length, resistivity=None, conductivity=None)¶ Bases:
blond.impedances.impedance_sources._ImpedanceObject
Impedance contribution from resistive wall for a cilindrical beam pipe
The model is the following:
\[Z(f) = \frac{Z_0 c L}{ \pi } \frac{ 1 }{ \left[1 - i \sign{f}\right] 2 b c \sqrt{ \frac{\sigma_c Z_0 c }{ 4 \pi |f| } + i 2 \pi b^2 f }\]- Parameters
pipe_radius (float) – Beam pipe radius in m
pipe_length (float) – Beam pipe length in m
resistivity (float) – Beam pipe resistivity in \(m / s\)
conductivity (float) – Beam pipe conductivity in \(s / m\)
-
pipe_radius
¶ Beam pipe radius in m
- Type
float
-
pipe_length
¶ Beam pipe length in m
- Type
float
-
conductivity
¶ Beam pipe conductivity in \(s / m\)
- Type
float
-
Z0
¶ Characteristic impedance of vacuum in* \(\Omega\)
- Type
float
Examples
>>> pipe_radius = 1 >>> pipe_length = 2 >>> resistivity = 3 >>> rw = ResistiveWall(pipe_radius, pipe_length, resistivity) >>> frequency = np.array(1,2,3) >>> rw.imped_calc(frequency)
-
property
conductivity
¶
-
imped_calc
(frequency_array)¶ Impedance calculation method as a function of frequency.
- Parameters
frequency_array (float array) – Input frequency array in Hz
-
frequency_array
¶ Input frequency array in Hz
- Type
float array
-
impedance
¶ Output impedance in \(\Omega + j \Omega\)
- Type
complex array
-
property
resistivity
¶
-
class
blond.impedances.impedance_sources.
Resonators
(R_S, frequency_R, Q, method='c++')¶ Bases:
blond.impedances.impedance_sources._ImpedanceObject
Impedance contribution from resonators, analytic formulas for both wake and impedance. The resonant modes (and the corresponding R and Q) can be inputed as a list in case of several modes. The model is the following:
\[Z(f) = \frac{R}{1 + j Q \left(\frac{f}{f_r}-\frac{f_r}{f}\right)}\]\[ \begin{align}\begin{aligned}W(t>0) = 2\alpha R e^{-\alpha t}\left(\cos{\bar{\omega}t} - \frac{\alpha}{\bar{\omega}}\sin{\bar{\omega}t}\right)\\W(0) = \alpha R\end{aligned}\end{align} \]\[ \begin{align}\begin{aligned}\omega_r = 2 \pi f_r\\\alpha = \frac{\omega_r}{2Q}\\\bar{\omega} = \sqrt{\omega_r^2 - \alpha^2}\end{aligned}\end{align} \]- Parameters
R_S (float list) – Shunt impepdance in \(\Omega\)
frequency_R (float list) – Resonant frequency in Hz
Q (float list) – Quality factor
method (string) – It defines which algorithm to use to calculate the impedance (C++ or Python)
-
R_S
¶ Shunt impepdance in \(\Omega\)
- Type
float array
-
frequency_R
¶ Resonant frequency in Hz
- Type
float array
-
Q
¶ Quality factor
- Type
float array
-
n_resonators
¶ number of resonators
- Type
int
Examples
>>> R_S = [1, 2, 3] >>> frequency_R = [1, 2, 3] >>> Q = [1, 2, 3] >>> resonators = Resonators(R_S, frequency_R, Q) >>> time = np.array(1,2,3) >>> resonators.wake_calc(time) >>> frequency = np.array(1,2,3) >>> resonators.imped_calc(frequency)
-
property
frequency_R
¶
-
property
omega_R
¶
-
class
blond.impedances.impedance_sources.
TravelingWaveCavity
(R_S, frequency_R, a_factor)¶ Bases:
blond.impedances.impedance_sources._ImpedanceObject
Impedance contribution from travelling wave cavities, analytic formulas for both wake and impedance. The resonance modes (and the corresponding R and a) can be inputed as a list in case of several modes.
The model is the following:
\[\begin{split}Z &= Z_+ + Z_- \\ Z_-(f) &= R \left[\left(\frac{\sin{\frac{a(f-f_r)}{2}}}{\frac{a(f-f_r)}{2}}\right)^2 - 2i \frac{a(f-f_r) - \sin{a(f-f_r)}}{\left(a(f-f_r)\right)^2}\right] \\ Z_+(f) &= R \left[\left(\frac{\sin{\frac{a(f+f_r)}{2}}}{\frac{a(f+f_r)}{2}}\right)^2 - 2i \frac{a(f+f_r) - \sin{a(f+f_r)}}{\left(a(f+f_r)\right)^2}\right]\end{split}\]\[\begin{split}W(0<t<\tilde{a}) &= \frac{4R}{\tilde{a}}\left(1-\frac{t}{\tilde{a}}\right)\cos{\omega_r t} \\ W(0) &= \frac{2R}{\tilde{a}}\end{split}\]\[a = 2 \pi \tilde{a}\]- Parameters
R_S (float list) – Shunt impepdance in \(\Omega\)
frequency_R (float list) – Resonant frequency in Hz
a_factor (float list) – Damping time a in s
-
R_S
¶ Shunt impepdance in \(\Omega\)
- Type
float array
-
frequency_R
¶ Resonant frequency in Hz
- Type
float array
-
a_factor
¶ Damping time a in s
- Type
float array
-
n_twc
¶ number of resonant modes
- Type
int
Examples
>>> R_S = [1, 2, 3] >>> frequency_R = [1, 2, 3] >>> a_factor = [1, 2, 3] >>> twc = TravelingWaveCavity(R_S, frequency_R, a_factor) >>> time = np.array(1,2,3) >>> twc.wake_calc(time) >>> frequency = np.array(1,2,3) >>> twc.imped_calc(frequency)
induced_voltage_analytical
Module¶
- Authors
Danilo Quartullo, Joel Repond
-
blond.impedances.induced_voltage_analytical.
analytical_gaussian_resonator
(sigma_t, Q, R_s, omega_r, tau_array, n_particles)¶ It gives the analytical induced voltage for a gaussian bunch and a resonator.
- Parameters
sigma_t (float) – Rms of the beam longitudinal coordinates [s].
Q (float) – Quality factor of the resonator [1].
R_s (float) – Shunt impedance of the resonator [\(\Omega\)].
omega_r (float) – Angular resonant frequency [rad/s].
tau_array (float array) – Time array [s] for which the induced voltage should be calculated.
n_particles (int) – Beam intensity [1].
- Returns
induced_voltage – Induced voltage [V] (multiplied by -1 according to BLonD conventions).
- Return type
float array
Notes
The formula can be derived using the following lines of codes in Mathematica:
lambda = Exp[-(tau - t)^2 / (2 * sigma^2)] / ((2*Pi)^0.5 * sigma) wake = 2 * alpha * Rs * Exp[-alpha * t] * (Cos[ombar * t] - alpha / ombar
Sin[ombar * t])
- output = Integrate[wake * lambda, {t, 0, Infinity},
Assumptions -> {sigma > 0, alpha > 0, Rs > 0, ombar > 0, tau in Reals}]
- outputreal = Simplify[ComplexExpand[Re[output]],
{sigma > 0, alpha > 0, Rs > 0, ombar > 0, tau in Reals}]
The imaginary part of output is identically equal to zero even if Mathematica cannot see that. To verify that the expression is correct launch:
output2 = Integrate[wake * lambda, t] outputreal2 = Simplify[ComplexExpand[Re[output2]],
{sigma > 0, alpha > 0, Rs > 0, ombar > 0, tau in Reals, t > 0}]
d/dt output2; Simplify[ComplexExpand[Re[%]],
{sigma > 0, alpha > 0, Rs > 0, ombar > 0, tau in Reals, t > 0}]
and check that applying the fundamental calculus theorem to outputreal2, one obtains back outputreal; the formula [5] in
H. E. Salzer, “Formulas for Calculating the Error Function of a Complex Variable”, Mathematical Tables and Other Aids to Computation, Vol. 5, No. 34 (Apr., 1951),pp. 67-70
can be useful.
music
Module¶
- Authors
Danilo Quartullo, Konstantinos Iliakis
-
class
blond.impedances.music.
Music
(Beam, resonator, n_macroparticles, n_particles, t_rev)¶ Bases:
object
Implementation of the MuSiC algorithm in C++ to calculate the exact induced voltage generated by resonant modes in time domain without using slices, cost = O(n). The corresponding methods in Python are kept for reference. The method track_classic, which calculates in time domain the exact voltage with the O(n^2) algorithm used in the usual voltage definition, is kept just for reference.
- Parameters
Beam (object) – Beam object.
resonator (float list) – List of the resonator parameters: [shunt impedance [\(\Omega\)], angular resonant frequency [rad/s], quality factor [1]].
n_macroparticles (int) – Number of macro-particles [1].
n_particles (float) – Beam intensity [1].
t_rev (float) – Revolution period [s]
-
beam
¶ Beam object.
- Type
object
-
R_S
¶ shunt impedance [\(\Omega\)]
- Type
float
-
omega_R
¶ angular resonant frequency [rad/s]
- Type
float
-
Q
¶ quality factor [1]
- Type
float
-
n_macroparticles
¶ Number of macro-particles [1].
- Type
int
-
n_particles
¶ Beam intensity [1].
- Type
float
-
alpha
¶ Definition dependent on previously defined attributes.
- Type
float
-
omega_bar
¶ Definition dependent on previously defined attributes.
- Type
float
-
const
¶ Definition dependent on previously defined attributes.
- Type
float
-
induced_voltage
¶ Output induced voltage [V] (multiplied by -1 for BLonD conventions)
- Type
float array
-
coeff1
¶ Definition dependent on previously defined attributes.
- Type
float
-
coeff2
¶ Definition dependent on previously defined attributes.
- Type
float
-
coeff3
¶ Definition dependent on previously defined attributes.
- Type
float
-
coeff4
¶ Definition dependent on previously defined attributes.
- Type
float
-
input_first_component
¶ First component of vertical array in MuSiC algorithm
- Type
float
-
input_second_component
¶ Second component of vertical array in MuSiC algorithm
- Type
float
-
t_rev
¶ Revolution period [s]
- Type
float
-
last_dt
¶ Last longitudinal coordinate of the beam [s]
- Type
float
-
array_parameters
¶ Array gathering four attributes already defined to be used in the C++ algorithm.
- Type
float array
Notes
The energies dE of the particles in the beam object are updated after the induced voltage calculation.
See also
The MuSiC algorithm is described in: M. Migliorati, L. Palumbo, ‘Multibunch and multiparticle simulation code with an alternative approach to wakefield effects’, Phys. Rev. ST Accel. Beams 18, 2015.
-
track_classic
()¶ Voltage in time domain using the basic definition (Python code)
-
track_cpp
()¶ Voltage in time domain (single-turn) using MuSiC (C++ code). Note: this method should also be called at turn number 1 when multi-turn voltage computations are needed.
Examples
>>> import impedances.music as musClass >>> from setup_cpp import libblond >>> >>> music_cpp = musClass.Music(my_beam, [R_S, 2*np.pi*frequency_R, Q], >>> n_macroparticles, n_particles, t_rev) >>> music_cpp.track_cpp()
-
track_cpp_multi_turn
()¶ Voltage in time domain (multi-turn) using MuSiC (C++ code). Note: this method should be called from turn number 2 onwards when multi-turn voltage computations are needed..
Examples
>>> import impedances.music as musClass >>> from setup_cpp import libblond >>> >>> music_cpp = musClass.Music(my_beam, [R_S, 2*np.pi*frequency_R, Q], >>> n_macroparticles, n_particles, t_rev) >>> music_cpp.track_cpp() >>> for i in range(2, n_turns): >>> music_cpp.track_cpp_multi_turn()
-
track_py
()¶ Voltage in time domain (single-turn) using MuSiC (Python code). Note: this method should also be called at turn number 1 when multi-turn voltage computations are needed.
Examples
>>> import impedances.music as musClass >>> >>> music_cpp = musClass.Music(my_beam, [R_S, 2*np.pi*frequency_R, Q], >>> n_macroparticles, n_particles, t_rev) >>> music_cpp.track_py()
-
track_py_multi_turn
()¶ Voltage in time domain (multi-turn) using MuSiC (Python code). Note: this method should be called from turn number 2 onwards when multi-turn voltage computations are needed..
Examples
>>> import impedances.music as musClass >>> >>> music_cpp = musClass.Music(my_beam, [R_S, 2*np.pi*frequency_R, Q], >>> n_macroparticles, n_particles, t_rev) >>> music_cpp.track_py() >>> for i in range(2, n_turns): >>> music_cpp.track_py_multi_turn()