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

wake_calc(new_time_array)

The wake from the table is interpolated using the new time array.

Parameters

new_time_array (float array) – Input time array in s

new_time_array

Input time array in s

Type

float array

wake

Output interpolated wake in \(\Omega / s\)

Type

float 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
wake_calc(time_array)

Wake calculation method as a function of time.

Parameters

time_array (float array) – Input time array in s

time_array

Input time array in s

Type

float array

wake

Output wake in \(\Omega / s\)

Type

float array

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)
imped_calc(frequency_array)

Impedance calculation method as a function of frequency.

Parameters

frequency_array (float array) – Input frequency array in Hz

frequency_array

nput frequency array in Hz

Type

float array

impedance

Output impedance in \(\Omega + j \Omega\)

Type

complex array

wake_calc(time_array)

Wake calculation method as a function of time.

Parameters

time_array (float array) – Input time array in s

time_array

Input time array in s

Type

float array

wake

Output wake in \(\Omega / s\)

Type

float array

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()