From 2afb94e6365d828431dc184582c0d98c8d32f6c3 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Birk=20Karlsen-B=C3=A6ck?= <birk.karlsen.baeck@gmail.com>
Date: Thu, 31 Oct 2024 20:20:17 +0100
Subject: [PATCH 1/9] Added BLonD-Xsuite interface to BLonD along with
benchmarking examples
---
..._single_particle_with_oneturnmatrix_lhc.py | 226 ++++++++++
.../EX_24_single_particle_ramp_lhc.py | 236 +++++++++++
...X_25_single_particle_two_rfstations_lhc.py | 256 +++++++++++
...EX_26_single_particle_oneturnmatrix_psb.py | 232 ++++++++++
.../EX_27_single_particle_ramp_psb.py | 238 +++++++++++
__EXAMPLES/main_files/EX_28_phase_loop_lhc.py | 246 +++++++++++
blond/interfaces/__init__.py | 0
blond/interfaces/xsuite.py | 401 ++++++++++++++++++
8 files changed, 1835 insertions(+)
create mode 100644 __EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
create mode 100644 __EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
create mode 100644 __EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
create mode 100644 __EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
create mode 100644 __EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
create mode 100644 __EXAMPLES/main_files/EX_28_phase_loop_lhc.py
create mode 100644 blond/interfaces/__init__.py
create mode 100644 blond/interfaces/xsuite.py
diff --git a/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
new file mode 100644
index 00000000..988d40c9
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
@@ -0,0 +1,226 @@
+'''
+LHC simulation of a single particle.
+
+:Author: **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+#mpl.use('Agg')
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+
+# Interface
+from xsuite_blond_interface import (BlondElement, BlondObserver,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_01'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 26658.8832 # Machine circumference [m]
+p_s = 450e9 # Synchronous momentum [eV/c]
+p_f = 450.0e9 # Synchronous momentum, final
+h = 35640 # Harmonic number [-]
+alpha = 0.00034849575112251314 # First order mom. comp. factor [-]
+V = 5e6 # RF voltage [V]
+dphi = 0 # Phase modulation/offset [rad]
+
+# Bunch parameters
+N_bunches = 1 # Number of bunches [-]
+N_m = 1 # Number of macroparticles [-]
+N_p = 1.15e11 # Intensity
+blen = 1.25e-9 # Bunch length [s]
+energy_err = 100e6 # Beamenergy error [eV/c]
+
+# Simulation parameters
+N_t = 481 # Number of (tracked) turns [-]
+N_buckets = 2 # Number of buckets [-]
+dt_plt = 50 # Timestep between plots [-]
+input_dt = 2 * blen - 0.1e-9 # Input particles dt [s]
+input_dE = 0.0 # Input particles dE [eV]
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring(C, alpha, np.linspace(p_s, p_f, N_t + 1), Proton(), n_turns=N_t)
+rfstation = RFStation(ring, [h], [V], [dphi])
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particle ---
+beam.dt = np.array([input_dt])
+beam.dE = np.array([input_dE])
+
+# --- RF tracker object ---
+blond_tracker = RingAndRFTracker(rfstation, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation.omega_rf[0, 0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation.phi_s[0]) + (np.pi - 2 * rfstation.phi_s[0]) / 2 * np.sin(rfstation.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation.charge) * rfstation.voltage[0, 0] * rfstation.beta[0] * rfstation.beta[0]
+ * rfstation.energy[0] / (np.pi * rfstation.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation, beam,
+ results_directory + '/EX_01_output_data', Profile=None)
+
+format_options = {'dirname': results_directory + '/EX_01_fig'}
+plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=None, h5file=results_directory + '/EX_01_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity = BlondElement(blond_tracker, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Preparing Xsuite items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing Xsuite items...")
+# --- Setup matrix ---
+# Make First order matrix map (takes care of drift in Xsuite)
+matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C
+)
+
+# Create line
+line = xt.Line(elements=[matrix], element_names={'matrix'} )
+line['matrix'].length = C
+
+# Insert the BLonD elements
+line.insert_element(index=0, element=cavity, name='blond_cavity')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# blond_cavity 0 BlondElement False False None True
+# matrix 0 LineSegmentMap True False None False
+# blond_monitor 26658.9 BlondObserver False False None True
+# blond_plots 26658.9 BlondElement False False None True
+# _end_point 26658.9 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0],
+ omega_rf=rfstation.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n').format(
+ beam.dE[0], beam.dt[0], rfstation.momentum[0], rfstation.gamma[0], rfstation.beta[0],
+ rfstation.energy[0], mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0]))
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i],
+ rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(
+ dE[0], dt[0], rfstation.momentum[i], rfstation.gamma[i], rfstation.beta[i],
+ rfstation.energy[i], mon.x[0, i].T, mon.px[0, i].T, mon.y[0, i].T, mon.py[0, i].T)
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i]))
+
+with open(results_directory + '/EX_01_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_01_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T, mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
new file mode 100644
index 00000000..91e88d32
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
@@ -0,0 +1,236 @@
+'''
+LHC simulation of a single particle with a ramp.
+
+: Author: **Birk Emil Karlsen-Baeck**, **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+#mpl.use('Agg')
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+
+# Interface
+from xsuite_blond_interface import (BlondElement, BlondObserver, EnergyUpdate,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_02'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 26658.8832 # Machine circumference [m]
+p_s = 450e9 # Synchronous momentum [eV/c]
+p_f = 450.1e9 # Synchronous momentum, final
+h = 35640 # Harmonic number [-]
+alpha = 0.00034849575112251314 # First order mom. comp. factor [-]
+V = 5e6 # RF voltage [V]
+dphi = 0 # Phase modulation/offset [rad]
+
+# Bunch parameters
+N_bunches = 1 # Number of bunches [-]
+N_m = 1 # Number of macroparticles [-]
+N_p = 1.15e11 # Intensity
+blen = 1.25e-9 # Bunch length [s]
+energy_err = 100e6 # Beamenergy error [eV/c]
+
+# Simulation parameters
+N_t = 330 # Number of (tracked) turns [-]
+N_buckets = 2 # Number of buckets [-]
+dt_plt = 30 # Timestep between plots [-]
+input_dt = 2 * blen - 0.4e-9 # Input particles dt [s]
+input_dE = 0.0 # Input particles dE [eV]
+
+# Xsuite parameters
+cavity_name = 'acsca.d5l4.b1'
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring(C, alpha, np.linspace(p_s, p_f, N_t + 1), Proton(), n_turns=N_t)
+rfstation = RFStation(ring, [h], [V], [dphi])
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particle ---
+beam.dt = np.array([input_dt])
+beam.dE = np.array([input_dE])
+
+# --- RF tracker object ---
+blond_tracker = RingAndRFTracker(rfstation, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation.omega_rf[0, 0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation.phi_s[0]) + (np.pi - 2 * rfstation.phi_s[0]) / 2 * np.sin(rfstation.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation.charge) * rfstation.voltage[0, 0] * rfstation.beta[0] * rfstation.beta[0]
+ * rfstation.energy[0] / (np.pi * rfstation.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation, beam,
+ results_directory + '/EX_02_output_data', Profile=None)
+
+format_options = {'dirname': results_directory + '/EX_02_fig'}
+plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2*t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=None, h5file=results_directory + '/EX_02_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity = BlondElement(blond_tracker, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Preparing Xsuite items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing Xsuite items...")
+# --- Setup matrix ---
+# Make First order matrix map (takes care of drift in Xsuite)
+matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C)
+
+# Create line
+line = xt.Line(elements=[matrix], element_names={'matrix'})
+line['matrix'].length = C
+
+# Insert the BLonD elements
+line.insert_element(index=0, element=cavity, name='blond_cavity')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Insert energy ramp
+energy_update = EnergyUpdate(ring.momentum[0, :])
+line.insert_element(index='matrix', element=energy_update, name='energy_update')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# energy_update 0 EnergyUpdate False False None True
+# blond_cavity 0 BlondElement False False None True
+# matrix 0 LineSegmentMap True False None False
+# blond_monitor 26658.9 BlondObserver False False None True
+# blond_plots 26658.9 BlondElement False False None True
+# _end_point 26658.9 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0],
+ omega_rf=rfstation.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(beam.dE[0], beam.dt[0], rfstation.momentum[0],
+ rfstation.gamma[0], rfstation.beta[0], rfstation.energy[0],
+ mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0])
+)
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i],
+ rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(dE[0], dt[0], rfstation.momentum[i], rfstation.gamma[i],
+ rfstation.beta[i], rfstation.energy[i], mon.x[0, i].T,
+ mon.px[0, i].T, mon.y[0, i].T, mon.py[0, i].T)
+
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i])
+ )
+
+with open(results_directory + '/EX_02_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_02_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T,mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
new file mode 100644
index 00000000..45c325e3
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
@@ -0,0 +1,256 @@
+'''
+LHC simulation of a single particle with two RF stations.
+
+: Author: **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+import matplotlib.pyplot as plt
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+
+# Interface
+from xsuite_blond_interface import (BlondElement,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_03'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 26658.8832 # Machine circumference [m]
+p_s = 450e9 # Synchronous momentum [eV/c]
+p_f = 450.0e9 # Synchronous momentum, final
+h = 35640 # Harmonic number [-]
+alpha = 0.00034849575112251314 # First order mom. comp. factor [-]
+V = 5e6 # RF voltage [V]
+dphi = 0 # Phase modulation/offset [rad]
+
+# Bunch parameters
+N_bunches = 1 # Number of bunches [-]
+N_m = 1 # Number of macroparticles [-]
+N_p = 1.15e11 # Intensity
+blen = 1.25e-9 # Bunch length [s]
+energy_err = 100e6 # Beamenergy error [eV/c]
+
+# Simulation parameters
+N_t = 481 # Number of (tracked) turns [-]
+N_buckets = 2 # Number of buckets [-]
+dt_plt = 50 # Timestep between plots [-]
+input_dt = 2 * blen - 0.1e-9 # Input particles dt [s]
+input_dE = 0.0 # Input particles dE [eV]
+
+# Additional input
+C_part = 0.3
+V_part = 0.5
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring([C_part * C, (1 - C_part) * C],
+ [[alpha], [alpha]],
+ [p_s, p_s],
+ Proton(), n_turns=N_t, n_sections=2)
+rfstation_1 = RFStation(ring, [h], [V * V_part], [dphi], section_index=1)
+rfstation_2 = RFStation(ring, [h], [V * (1 - V_part)], [dphi], section_index=2)
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particle ---
+beam.dt = np.array([input_dt])
+beam.dE = np.array([input_dE])
+
+# --- RF tracker object ---
+rftracker_1 = RingAndRFTracker(rfstation_1, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+rftracker_2 = RingAndRFTracker(rfstation_2, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation_1.omega_rf[0,0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation_1.phi_s[0]) + (np.pi - 2*rfstation_1.phi_s[0])
+ / 2 * np.sin(rfstation_1.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation_1.charge) * V * rfstation_1.beta[0] * rfstation_1.beta[0] * rfstation_1.energy[0]
+ / (np.pi * rfstation_1.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation_1, beam,
+ results_directory + '/EX_03_output_data', Profile=None)
+
+format_options = {'dirname': results_directory + '/EX_03_fig'}
+plots = Plot(ring, rfstation_1, beam, dt_plt, N_t, 0, 2 * t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=None, h5file=results_directory + '/EX_03_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity_1 = BlondElement(rftracker_1, beam)
+cavity_2 = BlondElement(rftracker_2, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Xsuite initialization -----------------------------------------------------------------------------------------------
+print(f"\nXsuite initialization...")
+# --- Setup matrices ---
+# Make First order matrix maps (takes care of drift in Xsuite)
+matrix_1 = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C * C_part)
+matrix_2 = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C * (1 - C_part))
+
+# Create new line
+line = xt.Line(elements=[matrix_1], element_names={'matrix_1'} )
+line['matrix_1'].length = C_part * C
+line.append_element(element=matrix_2, name='matrix_2')
+line['matrix_2'].length = (1 - C_part) * C
+
+# Insert the BLonD elements
+line.insert_element(index='matrix_1', element=cavity_1, name='blond_cavity_1')
+line.insert_element(index='matrix_2', element=cavity_2, name='blond_cavity_2')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# blond_cavity_1 0 BlondElement False False None True
+# matrix_1 0 LineSegmentMap True False None False
+# blond_cavity_2 7997.66 BlondElement False False None True
+# matrix_2 7997.66 LineSegmentMap True False None False
+# blond_monitor 26658.9 BlondObserver False False None True
+# blond_plots 26658.9 BlondElement False False None True
+# _end_point 26658.9 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation_1.phi_s[0],
+ omega_rf=rfstation_1.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(beam.dE[0], beam.dt[0], rfstation_1.momentum[0],
+ rfstation_1.gamma[0], rfstation_1.beta[0], rfstation_1.energy[0],
+ mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0])
+)
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:,i].T, mon.ptau[:,i].T,
+ rfstation_1.beta[i],
+ rfstation_1.energy[i],
+ phi_s=rfstation_1.phi_s[i] - rfstation_1.phi_rf[0, 0],
+ omega_rf=rfstation_1.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(dE[0], dt[0], rfstation_1.momentum[i],
+ rfstation_1.gamma[i], rfstation_1.beta[i],
+ rfstation_1.energy[i], mon.x[0, i].T, mon.px[0, i].T,
+ mon.y[0, i].T, mon.py[0, i].T)
+
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i])
+ )
+
+with open(results_directory + '/EX_03_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_03_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T,mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
new file mode 100644
index 00000000..aac0e0e4
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
@@ -0,0 +1,232 @@
+'''
+PSB simulation of a single particle.
+
+:Author: **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+from scipy.constants import c, e, m_p
+import matplotlib.pyplot as plt
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+
+# Interface
+from xsuite_blond_interface import (BlondElement,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_04'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 2 * np.pi * 25.000057945260654 # Machine circumference, radius 25 [m]
+gamma_t = 4.11635447373496 # Transition gamma [-]
+alpha = 1./gamma_t/gamma_t # First order mom. comp. factor [-]
+h = 1 # Harmonic number [-]
+V = 8e3 # RF voltage [V]
+dphi = np.pi # Phase modulation/offset [rad]
+
+# Beam parameters
+N_m = 1 # Number of macroparticles [-]
+N_p = 1e11 # Intensity / number of particles
+sigma_dt = 180e-9 / 4 # [s]
+blen = sigma_dt * 4 # Initial bunch length, 4 sigma [s]
+kin_beam_energy = 1.4e9 # [eV]
+
+
+# Simulation parameters
+N_t = 6800 # Number of (tracked) turns [-]
+N_buckets = 1 # Number of buckets [-]
+dt_plt = 500 # Timestep between plots [-]
+input_dt = 5.720357923415153e-07 - sigma_dt # Input particles dt [s]
+input_dE = 0.0 # Input particles dE [eV]
+
+# Derived parameters
+E_0 = m_p * c**2 / e # [eV]
+tot_beam_energy = E_0 + kin_beam_energy # [eV]
+sync_momentum = np.sqrt(tot_beam_energy**2 - E_0**2) # [eV / c]
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring(C, alpha, np.linspace(sync_momentum, sync_momentum, N_t + 1), Proton(), n_turns=N_t)
+rfstation = RFStation(ring, [h], [V], [dphi])
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particle ---
+beam.dt = np.array([input_dt])
+beam.dE = np.array([input_dE])
+
+# --- RF tracker object ---
+blond_tracker = RingAndRFTracker(rfstation, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation.omega_rf[0, 0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation.phi_s[0]) + (np.pi - 2 * rfstation.phi_s[0]) / 2 * np.sin(rfstation.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation.charge) * rfstation.voltage[0, 0] * rfstation.beta[0] * rfstation.beta[0]
+ * rfstation.energy[0] / (np.pi * rfstation.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation, beam,
+ results_directory + '/EX_04_output_data', Profile=None)
+
+format_options = {'dirname': results_directory + '/EX_04_fig'}
+plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=None, h5file=results_directory + '/EX_04_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity = BlondElement(blond_tracker, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Preparing Xsuite items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing Xsuite items...")
+# --- Setup matrix ---
+# Make First order matrix map (takes care of drift in Xsuite)
+matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C)
+
+# Create line
+line = xt.Line(elements=[matrix], element_names={'matrix'})
+line['matrix'].length = C
+
+# Insert the BLonD elements
+line.insert_element(index=0, element=cavity, name='blond_cavity')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=sync_momentum, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# blond_cavity 0 BlondElement False False None True
+# matrix 0 LineSegmentMap True False None False
+# blond_monitor 157.08 BlondObserver False False None True
+# blond_plots 157.08 BlondElement False False None True
+# _end_point 157.08 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(beam.dE[0], beam.dt[0], rfstation.momentum[0],
+ rfstation.gamma[0], rfstation.beta[0], rfstation.energy[0],
+ mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0])
+)
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i],
+ rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(dE[0], dt[0], rfstation.momentum[i], rfstation.gamma[i],
+ rfstation.beta[i], rfstation.energy[i], mon.x[0, i].T,
+ mon.px[0, i].T, mon.y[0, i].T, mon.py[0, i].T)
+
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i])
+ )
+
+with open(results_directory + '/EX_04_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_04_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T,mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
new file mode 100644
index 00000000..18df366e
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
@@ -0,0 +1,238 @@
+'''
+PSB simulation of a single particle with a ramp.
+
+:Author: **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+from scipy.constants import c, e, m_p
+import matplotlib.pyplot as plt
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+
+# Interface
+from xsuite_blond_interface import (BlondElement, EnergyUpdate,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_05'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 2 * np.pi * 25.000057945260654 # Machine circumference, radius 25 [m]
+gamma_t = 4.11635447373496 # Transition gamma [-]
+alpha = 1./gamma_t/gamma_t # First order mom. comp. factor [-]
+h = 1 # Harmonic number [-]
+V = 8e3 # RF voltage [V]
+dphi = np.pi # Phase modulation/offset [rad]
+
+# Beam parameters
+N_m = 1 # Number of macroparticles [-]
+N_p = 1e11 # Intensity / number of particles
+sigma_dt = 180e-9 / 4 # [s]
+blen = sigma_dt * 4 # Initial bunch length, 4 sigma [s]
+kin_beam_energy = 1.4e9 # Kinetic energy [eV]
+
+
+# Simulation parameters
+N_t = 6000 # Number of (tracked) turns [-]
+N_buckets = 1 # Number of buckets [-]
+dt_plt = 500 # Timestep between plots [-]
+input_dt = 5.720357923415153e-07 - sigma_dt # Input particles dt [s]
+input_dE = 0.0 # Input particles dE [eV]
+
+# Derived parameters
+E_0 = m_p * c**2 / e # [eV]
+tot_beam_energy = E_0 + kin_beam_energy # [eV]
+sync_momentum = np.sqrt(tot_beam_energy**2 - E_0**2) # [eV / c]
+momentum_increase = 0.001e9 # [eV / c]
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring(C, alpha, np.linspace(sync_momentum, sync_momentum + momentum_increase, N_t + 1), Proton(), n_turns=N_t)
+rfstation = RFStation(ring, [h], [V], [dphi])
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particle ---
+beam.dt = np.array([input_dt])
+beam.dE = np.array([input_dE])
+
+# --- RF tracker object ---
+blond_tracker = RingAndRFTracker(rfstation, beam, Profile=None, interpolation=False,
+ CavityFeedback=None, BeamFeedback=None,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation.omega_rf[0, 0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation.phi_s[0]) + (np.pi - 2 * rfstation.phi_s[0]) / 2 * np.sin(rfstation.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation.charge) * rfstation.voltage[0, 0] * rfstation.beta[0] * rfstation.beta[0]
+ * rfstation.energy[0] / (np.pi * rfstation.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation, beam,
+ results_directory + '/EX_05_output_data', Profile=None)
+
+format_options = {'dirname': results_directory + '/EX_05_fig'}
+plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=None, h5file=results_directory + '/EX_05_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity = BlondElement(blond_tracker, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Preparing Xsuite items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing Xsuite items...")
+# --- Setup matrix ---
+# Make First order matrix map (takes care of drift in Xsuite)
+matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C)
+
+# Create line
+line = xt.Line(elements=[matrix], element_names={'matrix'})
+line['matrix'].length = C
+
+# Insert the BLonD elements
+line.insert_element(index=0, element=cavity, name='blond_cavity')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Insert energy ramp
+energy_update = EnergyUpdate(ring.momentum[0, :])
+line.insert_element(index=0, element=energy_update, name='energy_update')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=sync_momentum, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# energy_update 0 EnergyUpdate False False None True
+# blond_cavity 0 BlondElement False False None True
+# matrix 0 LineSegmentMap True False None False
+# blond_monitor 157.08 BlondObserver False False None True
+# blond_plots 157.08 BlondElement False False None True
+# _end_point 157.08 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(beam.dE[0], beam.dt[0], rfstation.momentum[0],
+ rfstation.gamma[0], rfstation.beta[0], rfstation.energy[0],
+ mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0])
+)
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i],
+ rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(dE[0], dt[0], rfstation.momentum[i], rfstation.gamma[i],
+ rfstation.beta[i], rfstation.energy[i], mon.x[0, i].T,
+ mon.px[0, i].T, mon.y[0, i].T, mon.py[0, i].T)
+
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i])
+ )
+
+with open(results_directory + '/EX_05_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_05_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T,mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
new file mode 100644
index 00000000..b605c4d7
--- /dev/null
+++ b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
@@ -0,0 +1,246 @@
+'''
+LHC simulation using the LHC phase loop.
+
+:Author: **Birk Emil Karlsen-Baeck**, **Thom Arnoldus van Rijswijk**
+'''
+
+# Imports
+import os
+import numpy as np
+import matplotlib.pyplot as plt
+
+# Xsuite imports
+import xpart as xp
+import xtrack as xt
+import xplt
+
+# BLonD objects
+from blond.beam.beam import Beam, Proton
+from blond.beam.profile import Profile, CutOptions, FitOptions
+from blond.beam.distributions import parabolic
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.rf_parameters import RFStation
+from blond.trackers.tracker import RingAndRFTracker
+from blond.llrf.beam_feedback import BeamFeedback
+
+# Interface
+from xsuite_blond_interface import (BlondElement, BlondObserver,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+# Monitor objects
+from blond.monitors.monitors import BunchMonitor
+from blond.plots.plot import Plot
+
+# Directory to save files
+this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+results_directory = this_directory + '../output_files/EX_06'
+os.makedirs(results_directory, exist_ok=True)
+
+# Parameters ----------------------------------------------------------------------------------------------------------
+# Accelerator parameters
+C = 26658.8832 # Machine circumference [m]
+p_s = 450e9 # Synchronous momentum [eV/c]
+p_f = 450.0e9 # Synchronous momentum, final
+h = 35640 # Harmonic number [-]
+alpha = 0.00034849575112251314 # First order mom. comp. factor [-]
+V = 5e6 # RF voltage [V]
+dphi = 0 # Phase modulation/offset [rad]
+
+# Bunch parameters
+N_bunches = 1 # Number of bunches [-]
+N_m = int(1e3) # Number of macroparticles [-]
+N_p = 1.15e11 # Intensity
+blen = 1.25e-9 # Bunch length [s]
+energy_err = 100e6 # Beamenergy error [eV/c]
+
+# Simulation parameters
+N_t = 2 * 481 # Number of (tracked) turns [-]
+N_buckets = 2 # Number of buckets [-]
+dt_plt = 50 # Timestep between plots [-]
+
+# Preparing BLonD items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing BLonD items...")
+# --- BLonD objects ---
+ring = Ring(C, alpha, np.linspace(p_s, p_f, N_t + 1), Proton(), n_turns=N_t)
+rfstation = RFStation(ring, [h], [V], [dphi])
+beam = Beam(ring, N_m, N_p)
+
+# --- Insert particles ---
+parabolic(ring, rfstation, beam, bunch_length=blen)
+beam.dE += energy_err
+
+# Create Profile
+profile = Profile(beam,
+ CutOptions(cut_left=-1.5 * rfstation.t_rf[0, 0],
+ cut_right=(N_buckets + 1.5) * rfstation.t_rf[0, 0],
+ n_slices=(N_buckets + 3) * 2**5),
+ FitOptions(fit_option='gaussian'))
+profile.track()
+
+# --- LHC beam-phase loop and synchronization loop ---
+PL_gain = 1 / (5 * ring.t_rev[0])
+SL_gain = PL_gain / 10
+
+bl_config = {'machine': 'LHC',
+ 'PL_gain': PL_gain,
+ 'SL_gain': SL_gain}
+BL = BeamFeedback(ring, rfstation, profile, bl_config)
+
+# --- RF tracker object ---
+blond_tracker = RingAndRFTracker(rfstation, beam, Profile=profile, interpolation=False,
+ CavityFeedback=None, BeamFeedback=BL,
+ TotalInducedVoltage=None,
+ with_xsuite=True)
+
+# --- Bunch monitor ---
+# Print parameters for plotting
+print(f"\n Calculating plotting parameters...")
+
+t_u = np.pi / rfstation.omega_rf[0, 0]
+Y_phi_s = np.sqrt(abs(-np.cos(rfstation.phi_s[0]) + (np.pi - 2 * rfstation.phi_s[0]) / 2 * np.sin(rfstation.phi_s[0])))
+dE_sep = np.sqrt(2 * abs(rfstation.charge) * rfstation.voltage[0, 0] * rfstation.beta[0] * rfstation.beta[0]
+ * rfstation.energy[0] / (np.pi * rfstation.harmonic[0, 0] * abs(ring.eta_0[0, 0]))) * Y_phi_s
+
+print('t_u: ' + str(2 * t_u))
+print('dE_sep,m: ' + str(dE_sep))
+
+# Make bunchmonitor
+bunchmonitor = BunchMonitor(ring, rfstation, beam,
+ results_directory + '/EX_06_output_data', Profile=profile, PhaseLoop=BL)
+
+format_options = {'dirname': results_directory + '/EX_06_fig'}
+plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2*t_u,
+ -1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
+ Profile=profile, PhaseLoop=BL, h5file=results_directory + '/EX_06_output_data',
+ format_options=format_options)
+
+# --- Creating interface elements ---
+# Blond tracker
+cavity = BlondElement(blond_tracker, beam)
+# BLonD monitor
+beam_monitor = BlondElement(bunchmonitor, beam)
+beam_profile = BlondElement(profile, beam)
+beam_plots = BlondElement(plots, beam)
+
+# Preparing Xsuite items -----------------------------------------------------------------------------------------------
+print(f"\nPreparing Xsuite items...")
+# --- Setup matrix ---
+# Make First order matrix map (takes care of drift in Xsuite)
+matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=alpha,
+ length=C)
+
+# Create line
+line = xt.Line(elements=[matrix], element_names={'matrix'})
+line['matrix'].length = C
+
+# Insert the BLonD elements
+line.insert_element(index=0, element=cavity, name='blond_cavity')
+line.append_element(element=beam_monitor, name='blond_monitor')
+line.append_element(element=beam_profile, name='blond_profile')
+line.append_element(element=beam_plots, name='blond_plots')
+
+# Add particles to line and build tracker
+line.particle_ref = xp.Particles(p0c=p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+line.build_tracker()
+
+# Show table
+line.get_table().show()
+# prints:
+
+# name s element_type isthick isreplica parent_name iscollective
+# blond_cavity 0 BlondElement False False None True
+# matrix 0 LineSegmentMap True False None False
+# blond_monitor 26658.9 BlondObserver False False None True
+# blond_plots 26658.9 BlondElement False False None True
+# _end_point 26658.9 False False None False
+
+# Simulating ----------------------------------------------------------------------------------------------------------
+print(f"\nSetting up simulation...")
+
+# --- Convert the initial BLonD distribution to xsuite coordinates ---
+zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0],
+ omega_rf=rfstation.omega_rf[0, 0])
+
+# --- Track matrix ---
+particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+mon = line.record_last_track
+
+# Saving turn-by-turn particle coordinates
+test_string = ''
+test_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'dE', 'dt', 'momentum', 'gamma', 'beta', 'energy', 'x', 'px', 'y', 'py')
+test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(beam.dE[0], beam.dt[0], rfstation.momentum[0],
+ rfstation.gamma[0], rfstation.beta[0], rfstation.energy[0],
+ mon.x[0, 0].T, mon.px[0, 0].T, mon.y[0, 0].T, mon.py[0, 0].T)
+
+# Saving turn-by-turn Xsuite particle coordinates
+xsuite_string = ''
+xsuite_string += '{:<17}\t{:<17}\t{:<17}\t{:<17}\t{:<17}\n'.format(
+ 'ptau', 'zeta', 'momentum', 'gamma', 'beta')
+xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, 0].T, mon.zeta[0, 0].T, np.float64(mon.p0c[0, 0]), np.float64(mon.gamma0[0, 0]),
+ np.float64(mon.beta0[0, 0])
+)
+
+# Convert the xsuite particle coordinates back to BLonD
+for i in range(N_t):
+ dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i],
+ rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i])
+
+ # Statistics
+ test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
+ '\t{:+10.10e}\t{:+10.10e}\n').format(dE[0], dt[0], rfstation.momentum[i], rfstation.gamma[i],
+ rfstation.beta[i], rfstation.energy[i], mon.x[0, i].T,
+ mon.px[0, i].T, mon.y[0, i].T, mon.py[0, i].T)
+
+ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.format(
+ mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
+ np.float64(mon.beta0[0, i])
+ )
+
+with open(results_directory + '/EX_06_output_data.txt', 'w') as f:
+ f.write(test_string)
+with open(results_directory + '/EX_06_output_data_xsuite.txt', 'w') as f:
+ f.write(xsuite_string)
+
+# Results --------------------------------------------------------------------------------------------------------------
+print(f"\nPlotting result...")
+
+# Plot Phasespace
+
+plt.figure()
+plt.plot(mon.zeta.T,mon.ptau.T)
+plt.scatter(np.mean(mon.zeta.T), np.mean(mon.ptau.T), label='mon')
+plt.grid()
+plt.title('Phasespace')
+plt.xlabel(r'$\zeta$')
+plt.ylabel(r'$p_{\tau}$')
+plt.savefig(results_directory + '/figure_Phasespace.png')
+
+# Use Xplt to plot Phasespace
+
+plot = xplt.PhaseSpacePlot(
+ mon
+)
+plot.fig.suptitle("Particle distribution for a single turn")
+plot.fig.savefig(results_directory + '/figure_ps_mon.png')
+
+# Show plots
+plt.show()
diff --git a/blond/interfaces/__init__.py b/blond/interfaces/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/blond/interfaces/xsuite.py b/blond/interfaces/xsuite.py
new file mode 100644
index 00000000..2a4260e1
--- /dev/null
+++ b/blond/interfaces/xsuite.py
@@ -0,0 +1,401 @@
+'''
+Functions and classes to interface BLonD with xsuite.
+
+:Authors: **Birk Emil Karlsen-Baeck**, **Thom Arnoldus van Rijswijk**, **Helga Timko**
+'''
+
+from __future__ import annotations
+from typing import TYPE_CHECKING
+import numpy as np
+from scipy.constants import c as clight
+
+from xtrack import ReferenceEnergyIncrease, ZetaShift
+from blond.trackers.tracker import RingAndRFTracker
+from blond.impedances.impedance import TotalInducedVoltage
+
+if TYPE_CHECKING:
+ from xtrack import Particles
+ from xtrack import Line
+ from blond.beam.beam import Beam
+ from blond.beam.profile import Profile
+ from typing import Any, Union
+
+
+def blond_beam_to_xsuite_coords(beam: Beam, beta0: float, energy0: float,
+ omega_rf: float, phi_s: float = 0):
+ r'''
+ Coordinate transformation from Xsuite to BLonD at a given turn or multiple turns if numpy arrays ar given.
+
+ :param beam: blond.beam.beam.Beam class.
+ BLonD beam class to extract the initial coordinates dE and dt.
+ :param beta0: float.
+ Synchronous beta.
+ :param energy0: float.
+ Synchronous energy.
+ :param omega_rf: float.
+ The rf angular frequency.
+ :param phi_s: float.
+ Synchronous phase in radians equivalent to Xsuite's phi_rf (below transition energy input should be phi_s-phi_rf).
+ The default value is 0.
+
+ :return zeta, ptau: as numpy-arrays (or single variable)
+ '''
+ ptau = beam.dE / (beta0 * energy0)
+ zeta = -(beam.dt - phi_s / omega_rf) * beta0 * clight
+ return zeta, ptau
+
+
+def xsuite_coords_to_blond_coords(zeta: Union[float, np.ndarray], ptau: Union[float, np.ndarray],
+ beta0: float, energy0: float, omega_rf: float, phi_s: float = 0):
+ r'''
+ Coordinate transformation from Xsuite to BLonD.
+
+ :param zeta: float or numpy-array.
+ The zeta coordinate as defined in Xsuite.
+ :param ptau: float or numpy-array.
+ The ptau coordinate as defined in Xsuite.
+ :param beta0: float.
+ The synchronous beta.
+ :param energy0: float.
+ The synchronous energy in eV.
+ :param omega_rf: float.
+ The rf angular frequency.
+ :param phi_s: float.
+ The synchronous phase in radians equivalent to Xsuite's phi_rf
+ (below transition energy input should be phi_s-phi_rf)
+
+ :return dt, dE: as numpy-arrays (or single variable)
+ '''
+
+ dE = ptau * beta0 * energy0
+ dt = -zeta / (beta0 * clight) + phi_s / omega_rf
+ return dt, dE
+
+
+class BlondElement:
+
+ def __init__(self, trackable: Any, beam: Beam, update_zeta: bool = False) -> None:
+ r"""
+ The BlondElement class contains a trackable object from the BLonD simulation suite and the Beam object from
+ BLonD. The class is intended to be an element to be added to the Line in XTrack instead of the default
+ RF cavities in XTrack.
+
+ The behavior of the tracking of the class depends on what trackable object is passed.
+ If the RingAndRFTracker class is passed then a coordinate transformation will be performed and the energy
+ deviation of the particles will be updated. If the TotalInducedVoltage class is passed then the
+ TotalInducedVoltage.induced_voltage_sum() method will be called. Lastly, if any other class is passed then
+ the tracking of BlondElement will track this class.
+
+ trackable : BLonD trackable object
+ BLonD object to be tracked, e.g. the RingAndRFTracker.
+ beam : blond.beam.beam.Beam class
+ BLonD Beam class used for tracking.
+ """
+
+ self.trackable = trackable
+ self.beam = beam
+ self.update_zeta = update_zeta
+
+ # Initialize time- and orbit-shift to BLonD coordinates
+ self._dt_shift = None
+ self.orbit_shift = ZetaShift(dzeta=0)
+
+ # Check what BLonD trackable has been passed and track the object accordingly
+ if isinstance(self.trackable, RingAndRFTracker):
+ self.track = getattr(self, "rf_track")
+ elif isinstance(self.trackable, TotalInducedVoltage):
+ self.track = getattr(self, "ind_track")
+ else:
+ self.track = getattr(self, "obs_track")
+
+ def rf_track(self, particles: Particles) -> None:
+ r"""
+ Tracking method which is called if the trackable BLonD class is the RingAndRFTracker.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ """
+
+ # Compute the shift to BLonD coordinates
+ self._get_time_shift()
+
+ # Convert the Xsuite coordinates to BLonD coordinates
+ self.xsuite_part_to_blond_beam(particles)
+
+ # Track with only the energy kick
+ self.trackable.track()
+
+ # Convert the BLonD energy coordinate to the equivalent Xsuite coordinate
+ self.blond_beam_to_xsuite_part(particles, self.update_zeta)
+
+ # Update the zeta shift due to potential frequency shifts in the RF
+ self._orbit_shift(particles)
+
+ def ind_track(self, particles: Particles) -> None:
+ r"""
+ Tracking method which is called if the trackable BLonD class is the TotalInducedVoltage.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ """
+ self.trackable.induced_voltage_sum()
+
+ def obs_track(self, particles: Particles) -> None:
+ r"""
+ Tracking method which is called if the trackable BLonD class is not the RingAndRFTracker or TotalInducedVoltage.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ """
+ self.trackable.track()
+
+ def xsuite_part_to_blond_beam(self, particles: Particles) -> None:
+ r"""
+ Coordinate transformation from Xsuite to BLonD.
+ It uses the initial particle coordinates and beam properties in stored in the Particles class from xtrack.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ """
+ # Convert Xsuite momentum to BLonD energy deviation
+ self.beam.dE[:] = particles.beta0 * particles.energy0 * particles.ptau
+
+ # Convert Xsuite zeta coordinate to BLonD time deviation
+ self.beam.dt[:] = -particles.zeta / particles.beta0 / clight + self._dt_shift
+
+ # Check what particles are still alive
+ self.beam.id[:] = np.int_(particles.state > 0)
+
+ def blond_beam_to_xsuite_part(self, particles: Particles, update_zeta: bool = False) -> None:
+ r"""
+ Coordinate transformation from BLonD to Xsuite.
+ It uses the particle coordinates stored in Beam class of BLonD
+ It uses the beam properties in stored in the Particles class from xtrack
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ """
+ # Subtract the given acceleration kick in BLonD, in Xsuite this is dealt with differently
+ if isinstance(self.trackable, RingAndRFTracker):
+ self.beam.dE = self.beam.dE - self.trackable.acceleration_kick[self.trackable.counter[0] - 1]
+
+ # Convert BLonD energy deviation to Xsuite momentum
+ particles.ptau = self.beam.dE / (particles.beta0 * particles.energy0)
+
+ # Convert BLonD time deviation to Xsuite zeta.
+ # This step is not needed usually because the BLonD simulation only does the kick, so dt is not changed.
+ if update_zeta:
+ particles.zeta = - (self.beam.dt - self._dt_shift) * particles.beta0 * clight
+
+ # Check what particles are still alive after the BLonD track
+ mask_lost = (self.beam.id <= 0) & particles.state > 0
+
+ # If the particle is lost its state is set to -500 by convention
+ particles.state[mask_lost] = -500
+
+ def _get_time_shift(self) -> None:
+ r'''
+ Computes the time-shift between the Xsuite and BLonD coordinate systems.
+ '''
+ # Get turn counter from the RingAndRFTracker
+ counter = self.trackable.rf_params.counter[0]
+
+ # Compute the time-shift based on the synchronous phase
+ self._dt_shift = ((self.trackable.rf_params.phi_s[counter] - self.trackable.rf_params.phi_rf[0, counter])
+ / self.trackable.rf_params.omega_rf[0, counter])
+
+ def _orbit_shift(self, particles) -> None:
+ r'''
+ Computes the radial steering due to rf periods which are not an integer multiple of the revolution period.
+ This is for example needed when tracking with global LLRF feedback loops.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ '''
+ # Get turn counter from the RingAndRFTracker
+ counter = self.trackable.counter[0]
+
+ # Compute the orbit shift due to the difference in rf frequency
+ dzeta = self.trackable.rf_params.ring_circumference
+ omega_rf = self.trackable.rf_params.omega_rf[:, counter]
+ omega_rf_design = (2 * np.pi * self.trackable.rf_params.harmonic[:, counter]
+ / self.trackable.rf_params.t_rev[counter])
+ domega = omega_rf - omega_rf_design
+
+ dzeta *= domega / omega_rf_design
+
+ # Apply the shift
+ self.orbit_shift.dzeta = dzeta
+
+ # Track the shift
+ self.orbit_shift.track(particles)
+
+
+class EnergyUpdate:
+
+ def __init__(self, momentum: Union[list, np.ndarray]) -> None:
+ r"""
+ Class to update the synchronous energy from the momentum program in BLonD.
+
+ :param momentum: list or numpy-array.
+ Momentum program from BLonD in units of eV
+ """
+
+ # Load momentum program
+ self.momentum = momentum
+
+ # Find initial momentum update
+ init_p0c = self.momentum[1] - self.momentum[0]
+
+ # Enter the initial momentum update in the ReferenceEnergyIncrease class in xsuite
+ self.xsuite_energy_update = ReferenceEnergyIncrease(Delta_p0c=init_p0c)
+
+ def track(self, particles: Particles) -> None:
+ r'''
+ Track method for the class to update the synchronous energy.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ '''
+ # Check for particles which are still alive
+ mask_alive = particles.state > 0
+
+ # Use the still alive particles to find the current turn momentum
+ p0c_before = particles.p0c[mask_alive]
+
+ # Find the momentum for the next turn
+ p0c_after = self.momentum[particles.at_turn[mask_alive][0]]
+
+ # Update the energy increment
+ self.xsuite_energy_update.Delta_p0c = p0c_after - p0c_before[0]
+
+ # Apply the energy increment to the particles
+ self.xsuite_energy_update.track(particles)
+
+
+class EnergyFrequencyUpdate:
+
+ def __init__(self, momentum: np.ndarray, f_rf: np.ndarray, line: Line, cavity_name: str) -> None:
+ r"""
+ Class to update energy of Particles class turn-by-turn with the ReferenceEnergyIncrease function
+ from xtrack. Additionally it updates the frequency of the xtrack cavity in the line.
+ Intended to be used without BLonD-Xsuite interface.
+
+ :param momentum: numpy-array
+ The momentum program from BLonD in eV
+ :param f_rf: numpy-array
+ The frequency program from BLonD in Hz.
+ particles : xtrack.Line
+ Line class from xtrack
+ :param cavity_name: string
+ Name of cavity to update frequency.
+ """
+
+ # Load the parameters
+ self.momentum = momentum
+ self.f_rf = f_rf
+ self.line = line
+ self.cavity_name = cavity_name
+
+ # Find initial momentum update
+ init_p0c = self.momentum[1] - self.momentum[0]
+
+ # Enter the initial momentum update in the ReferenceEnergyIncrease class in xsuite
+ self.xsuite_energy_update = ReferenceEnergyIncrease(Delta_p0c=init_p0c)
+
+ def track(self, particles: Particles) -> None:
+ r'''
+ Track-method from for the class. This method updates the synchronous momentum and the rf frequency.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ '''
+ # Check for particles which are still alive
+ mask_alive = particles.state > 0
+
+ # Use the still alive particles to find the current turn momentum
+ p0c_before = particles.p0c[mask_alive]
+
+ # Find the momentum for the next turn
+ p0c_after = self.momentum[particles.at_turn[mask_alive][0]]
+
+ # Update the energy increment
+ self.xsuite_energy_update.Delta_p0c = p0c_after - p0c_before[0]
+
+ # Apply the energy increment to the particles
+ self.xsuite_energy_update.track(particles)
+
+ # Update the rf frequency
+ self.line[self.cavity_name].frequency = self.f_rf[particles.at_turn[mask_alive][0]]
+
+
+class BlondObserver(BlondElement):
+
+ def __init__(self, trackable: Any, beam: Beam, blond_cavity: bool, update_zeta: bool = False,
+ profile: Profile = None) -> None:
+ r'''
+ Child-class of the BlondElement, except that it updates the coordinates
+ in BLonD when an observing element is used such as BunchMonitor.
+
+ :param trackable: BLonD trackable object
+ BLonD object to be tracked, e.g. the RingAndRFTracker.
+ :param beam: blond.beam.beam.Beam class
+ BLonD Beam class used for tracking.
+ :param blond_cavity: bool.
+ If there is no BlondCavity (bool = False), it updates its own turn-counter.
+ :param update_zeta: bool.
+ Boolean that decides whether zeta is converter back to dt after tracking object or not.
+ Usually not necessary so default is False.
+ :param profile: blond.beam.profile.Profile class
+ BLonD Profile class used for tracking.
+ '''
+ # Initialize the parent class
+ super().__init__(trackable, beam, update_zeta)
+
+ # Load the parameters
+ self.blond_cavity = blond_cavity
+ # For bunch monitoring we need the profile separately added
+ self.profile = profile
+
+ # Initializing arrays for storing some turn-by-turn properties of xtrack.Particles class
+ self.xsuite_ref_energy = np.zeros(self.trackable.rf_params.n_turns + 1)
+ self.xsuite_trev = np.zeros(self.trackable.rf_params.n_turns + 1)
+
+ def obs_track(self, particles: Particles) -> None:
+ r'''
+ observation tracker which performs the coordinate transformations.
+
+ particles : xtrack.Particles
+ Particles class from xtrack
+ '''
+ # Compute the shift to BLonD coordinates
+ self._get_time_shift()
+
+ # Convert the Xsuite coordinates to BLonD coordinates
+ self.xsuite_part_to_blond_beam(particles)
+
+ # Track profile if given
+ if self.profile is not None:
+ self.profile.track()
+
+ # Track
+ self.trackable.track()
+
+ if not self.blond_cavity:
+ # Updating the beam synchronous momentum etc.
+ turn = self.trackable.rf_params.counter[0]
+
+ self.beam.beta = self.trackable.rf_params.beta[turn + 1]
+ self.beam.gamma = self.trackable.rf_params.gamma[turn + 1]
+ self.beam.energy = self.trackable.rf_params.energy[turn + 1]
+ self.beam.momentum = self.trackable.rf_params.momentum[turn + 1]
+
+ # Update custom counter
+ self.trackable.rf_params.counter[0] += 1
+
+ # Convert the BLonD energy coordinate to the equivalent Xsuite coordinate
+ self.blond_beam_to_xsuite_part(particles, self.update_zeta)
+
+ # Track properties of xtrack.Particles
+ self.xsuite_ref_energy[self.trackable.rf_params.counter[0] - 1] = particles.energy0[0]
+ self.xsuite_trev[self.trackable.rf_params.counter[0] - 1] = particles.t_sim # Does not update per turn!
--
GitLab
From 6baa7a0e9e6f8b17dfa2f6e8390c103b1fb7a948 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Birk=20Karlsen-B=C3=A6ck?= <birk.karlsen.baeck@gmail.com>
Date: Fri, 1 Nov 2024 09:08:03 +0100
Subject: [PATCH 2/9] [examples][main files] Minor clean up of xsuite-blond
examples
---
..._single_particle_with_oneturnmatrix_lhc.py | 19 +++++++++----------
.../EX_24_single_particle_ramp_lhc.py | 18 +++++++++---------
...X_25_single_particle_two_rfstations_lhc.py | 18 +++++++++---------
...EX_26_single_particle_oneturnmatrix_psb.py | 18 +++++++++---------
.../EX_27_single_particle_ramp_psb.py | 18 +++++++++---------
__EXAMPLES/main_files/EX_28_phase_loop_lhc.py | 18 +++++++++---------
6 files changed, 54 insertions(+), 55 deletions(-)
diff --git a/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
index 988d40c9..9052aab0 100644
--- a/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
+++ b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
@@ -9,7 +9,6 @@ import os
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
-#mpl.use('Agg')
# Xsuite imports
import xpart as xp
@@ -23,9 +22,9 @@ from blond.input_parameters.rf_parameters import RFStation
from blond.trackers.tracker import RingAndRFTracker
# Interface
-from xsuite_blond_interface import (BlondElement, BlondObserver,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement, BlondObserver,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -33,7 +32,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_01'
+results_directory = this_directory + '../output_files/EX_23'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -91,12 +90,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation, beam,
- results_directory + '/EX_01_output_data', Profile=None)
+ results_directory + '/EX_23_output_data', Profile=None)
-format_options = {'dirname': results_directory + '/EX_01_fig'}
+format_options = {'dirname': results_directory + '/EX_23_fig'}
plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=None, h5file=results_directory + '/EX_01_output_data',
+ Profile=None, h5file=results_directory + '/EX_23_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -195,9 +194,9 @@ for i in range(N_t):
mon.ptau[0, i].T, mon.zeta[0, i].T, np.float64(mon.p0c[0, i]), np.float64(mon.gamma0[0, i]),
np.float64(mon.beta0[0, i]))
-with open(results_directory + '/EX_01_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_23_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_01_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_23_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
diff --git a/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
index 91e88d32..a36d3809 100644
--- a/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
+++ b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
@@ -23,9 +23,9 @@ from blond.input_parameters.rf_parameters import RFStation
from blond.trackers.tracker import RingAndRFTracker
# Interface
-from xsuite_blond_interface import (BlondElement, BlondObserver, EnergyUpdate,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement, BlondObserver, EnergyUpdate,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -33,7 +33,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_02'
+results_directory = this_directory + '../output_files/EX_24'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -94,12 +94,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation, beam,
- results_directory + '/EX_02_output_data', Profile=None)
+ results_directory + '/EX_24_output_data', Profile=None)
-format_options = {'dirname': results_directory + '/EX_02_fig'}
+format_options = {'dirname': results_directory + '/EX_24_fig'}
plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2*t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=None, h5file=results_directory + '/EX_02_output_data',
+ Profile=None, h5file=results_directory + '/EX_24_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -205,9 +205,9 @@ for i in range(N_t):
np.float64(mon.beta0[0, i])
)
-with open(results_directory + '/EX_02_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_24_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_02_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_24_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
diff --git a/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
index 45c325e3..41e6e88f 100644
--- a/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
+++ b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
@@ -21,9 +21,9 @@ from blond.input_parameters.rf_parameters import RFStation
from blond.trackers.tracker import RingAndRFTracker
# Interface
-from xsuite_blond_interface import (BlondElement,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -31,7 +31,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_03'
+results_directory = this_directory + '../output_files/EX_25'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -102,12 +102,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation_1, beam,
- results_directory + '/EX_03_output_data', Profile=None)
+ results_directory + '/EX_25_output_data', Profile=None)
-format_options = {'dirname': results_directory + '/EX_03_fig'}
+format_options = {'dirname': results_directory + '/EX_25_fig'}
plots = Plot(ring, rfstation_1, beam, dt_plt, N_t, 0, 2 * t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=None, h5file=results_directory + '/EX_03_output_data',
+ Profile=None, h5file=results_directory + '/EX_25_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -225,9 +225,9 @@ for i in range(N_t):
np.float64(mon.beta0[0, i])
)
-with open(results_directory + '/EX_03_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_25_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_03_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_25_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
diff --git a/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
index aac0e0e4..6ea893ad 100644
--- a/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
+++ b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
@@ -22,9 +22,9 @@ from blond.input_parameters.rf_parameters import RFStation
from blond.trackers.tracker import RingAndRFTracker
# Interface
-from xsuite_blond_interface import (BlondElement,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -32,7 +32,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_04'
+results_directory = this_directory + '../output_files/EX_26'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -95,12 +95,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation, beam,
- results_directory + '/EX_04_output_data', Profile=None)
+ results_directory + '/EX_26_output_data', Profile=None)
-format_options = {'dirname': results_directory + '/EX_04_fig'}
+format_options = {'dirname': results_directory + '/EX_26_fig'}
plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=None, h5file=results_directory + '/EX_04_output_data',
+ Profile=None, h5file=results_directory + '/EX_26_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -201,9 +201,9 @@ for i in range(N_t):
np.float64(mon.beta0[0, i])
)
-with open(results_directory + '/EX_04_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_26_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_04_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_26_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
diff --git a/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
index 18df366e..b4707f0f 100644
--- a/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
+++ b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
@@ -22,9 +22,9 @@ from blond.input_parameters.rf_parameters import RFStation
from blond.trackers.tracker import RingAndRFTracker
# Interface
-from xsuite_blond_interface import (BlondElement, EnergyUpdate,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement, EnergyUpdate,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -32,7 +32,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_05'
+results_directory = this_directory + '../output_files/EX_27'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -96,12 +96,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation, beam,
- results_directory + '/EX_05_output_data', Profile=None)
+ results_directory + '/EX_27_output_data', Profile=None)
-format_options = {'dirname': results_directory + '/EX_05_fig'}
+format_options = {'dirname': results_directory + '/EX_27_fig'}
plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2 * t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=None, h5file=results_directory + '/EX_05_output_data',
+ Profile=None, h5file=results_directory + '/EX_27_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -207,9 +207,9 @@ for i in range(N_t):
np.float64(mon.beta0[0, i])
)
-with open(results_directory + '/EX_05_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_27_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_05_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_27_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
diff --git a/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
index b605c4d7..cefbbabd 100644
--- a/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
+++ b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
@@ -24,9 +24,9 @@ from blond.trackers.tracker import RingAndRFTracker
from blond.llrf.beam_feedback import BeamFeedback
# Interface
-from xsuite_blond_interface import (BlondElement, BlondObserver,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+from blond.interfaces.xsuite import (BlondElement, BlondObserver,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -34,7 +34,7 @@ from blond.plots.plot import Plot
# Directory to save files
this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-results_directory = this_directory + '../output_files/EX_06'
+results_directory = this_directory + '../output_files/EX_28'
os.makedirs(results_directory, exist_ok=True)
# Parameters ----------------------------------------------------------------------------------------------------------
@@ -107,12 +107,12 @@ print('dE_sep,m: ' + str(dE_sep))
# Make bunchmonitor
bunchmonitor = BunchMonitor(ring, rfstation, beam,
- results_directory + '/EX_06_output_data', Profile=profile, PhaseLoop=BL)
+ results_directory + '/EX_28_output_data', Profile=profile, PhaseLoop=BL)
-format_options = {'dirname': results_directory + '/EX_06_fig'}
+format_options = {'dirname': results_directory + '/EX_28_fig'}
plots = Plot(ring, rfstation, beam, dt_plt, N_t, 0, 2*t_u,
-1.05 * dE_sep, 1.05 * dE_sep, xunit='s', separatrix_plot=True,
- Profile=profile, PhaseLoop=BL, h5file=results_directory + '/EX_06_output_data',
+ Profile=profile, PhaseLoop=BL, h5file=results_directory + '/EX_28_output_data',
format_options=format_options)
# --- Creating interface elements ---
@@ -215,9 +215,9 @@ for i in range(N_t):
np.float64(mon.beta0[0, i])
)
-with open(results_directory + '/EX_06_output_data.txt', 'w') as f:
+with open(results_directory + '/EX_28_output_data.txt', 'w') as f:
f.write(test_string)
-with open(results_directory + '/EX_06_output_data_xsuite.txt', 'w') as f:
+with open(results_directory + '/EX_28_output_data_xsuite.txt', 'w') as f:
f.write(xsuite_string)
# Results --------------------------------------------------------------------------------------------------------------
--
GitLab
From ba7846b90fd105f9f1941f4c953231b481506119 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Birk=20Karlsen-B=C3=A6ck?= <birk.karlsen.baeck@gmail.com>
Date: Fri, 1 Nov 2024 11:49:28 +0100
Subject: [PATCH 3/9] [unittests][interfaces] Made a first draft of the
unittests for the xsuite interface
---
unittests/interfaces/__init__.py | 0
unittests/interfaces/test_xsuite.py | 482 ++++++++++++++++++++++++++++
2 files changed, 482 insertions(+)
create mode 100644 unittests/interfaces/__init__.py
create mode 100644 unittests/interfaces/test_xsuite.py
diff --git a/unittests/interfaces/__init__.py b/unittests/interfaces/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/unittests/interfaces/test_xsuite.py b/unittests/interfaces/test_xsuite.py
new file mode 100644
index 00000000..be3770a7
--- /dev/null
+++ b/unittests/interfaces/test_xsuite.py
@@ -0,0 +1,482 @@
+# coding: utf8
+# Copyright 2014-2017 CERN. This software is distributed under the
+# terms of the GNU General Public Licence version 3 (GPL Version 3),
+# copied verbatim in the file LICENCE.md.
+# In applying this licence, CERN does not waive the privileges and immunities
+# granted to it by virtue of its status as an Intergovernmental Organization or
+# submit itself to any jurisdiction.
+# Project website: http://blond.web.cern.ch/
+
+"""
+Unittest for interfaces.xsuite
+
+:Authors: **Birk Emil Karlsen-Bæck**, **Thom Arnoldus van Rijswijk**
+"""
+
+import os
+import unittest
+
+import numpy as np
+from scipy.constants import c, e, m_p
+
+import xpart as xp
+import xtrack as xt
+
+from blond.beam.beam import Beam, Proton
+from blond.beam.profile import CutOptions, Profile
+from blond.input_parameters.rf_parameters import RFStation
+from blond.input_parameters.ring import Ring
+from blond.trackers.tracker import RingAndRFTracker, FullRingAndRF
+
+from blond.interfaces.xsuite import (BlondElement, BlondObserver,
+ EnergyUpdate,
+ blond_beam_to_xsuite_coords,
+ xsuite_coords_to_blond_coords)
+
+this_directory = os.path.dirname(os.path.realpath(__file__))
+
+# TODO: Make a test of the EnergyFrequencyUpdate and BlondObserver classes
+
+
+class TestXsuiteBLonDTransforms(unittest.TestCase):
+ # TODO: finish implementation of the test of the particle coordinates
+ def setUp(self):
+ pass
+
+
+class TestXsuiteLHC(unittest.TestCase):
+
+ def setUp(self):
+ # Accelerator parameters
+ self.C = 26658.8832 # Machine circumference [m]
+ self.p_s = 450e9 # Synchronous momentum [eV/c]
+ self.p_f = 450.1e9 # Synchronous momentum, final
+ self.h = 35640 # Harmonic number [-]
+ self.alpha = 0.00034849575112251314 # First order mom. comp. factor [-]
+ self.V = 5e6 # RF voltage [V]
+ self.dphi = 0 # Phase modulation/offset [rad]
+
+ def testSingleParticle(self):
+ r'''Test of a single particle simulation in the LHC at constant energy'''
+ rtol = 1e-2 # relative tolerance
+ atol = 0 # absolute tolerance
+
+ # ----- Interface Simulation -----
+ # Initialize the BLonD cavity
+ blond_tracker, beam = self.singleParticleBLonDSimulation()
+ cavity = BlondElement(blond_tracker, beam)
+
+ line = self.initializeXsuiteLine()
+
+ # Insert the BLonD elements
+ line.insert_element(index=0, element=cavity, name='blond_cavity')
+
+ dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+
+ # ----- Pure BLonD Simulation -----
+ blond_tracker, beam = self.singleParticleBLonDSimulation()
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+
+ # ----- Perform test -----
+ np.testing.assert_allclose(
+ dt_inter, dt_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticle the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dE_inter, dE_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticle the energy-coordinate differs'
+ )
+
+ def testSingleParticleRamp(self):
+ r'''Test of acceleration of a single particle simulation in the LHC'''
+ rtol = 1e-2 # relative tolerance
+ atol = 0 # absolute tolerance
+
+ # Initialize the BLonD cavity
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ cavity = BlondElement(blond_tracker, beam)
+
+ line = self.initializeXsuiteLine()
+
+ # Insert the BLonD elements
+ line.insert_element(index=0, element=cavity, name='blond_cavity')
+
+ # Insert energy ramp
+ energy_update = EnergyUpdate(blond_tracker.rf_params.momentum)
+ line.insert_element(index='matrix', element=energy_update, name='energy_update')
+
+ dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+
+ # ----- Pure BLonD Simulation -----
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+
+ # ----- Perform test -----
+ np.testing.assert_allclose(
+ dt_inter, dt_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleRamp the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dE_inter, dE_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleRamp the energy-coordinate differs'
+ )
+
+ def testSingleParticleTwoRFStations(self):
+ r'''Test of a single particle simulation in the LHC at constant energy with two RF stations'''
+ rtol = 1e-2 # relative tolerance
+ atol = 0 # absolute tolerance
+
+ # Initialize the BLonD cavities
+ blond_tracker_1, blond_tracker_2, beam = self.singleParticleBLonDSimulation(two_rfstations=True)
+ cavity_1 = BlondElement(blond_tracker_1, beam)
+ cavity_2 = BlondElement(blond_tracker_2, beam)
+
+ line = self.initializeXsuiteLine(two_rfstations=True)
+
+ # Insert the BLonD elements
+ line.insert_element(index='matrix', element=cavity_1, name='blond_cavity_1')
+ line.insert_element(index='matrix_2', element=cavity_2, name='blond_cavity_2')
+
+ dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker_1, beam)
+
+ # ----- Pure BLonD Simulation -----
+ blond_tracker, beam = self.singleParticleBLonDSimulation(two_rfstations=True)
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+
+ # ----- Perform test -----
+ np.testing.assert_allclose(
+ dt_inter, dt_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleTwoRFStations the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dE_inter, dE_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleTwoRFStations the energy-coordinate differs'
+ )
+
+ def testBunchWithPhaseLoop(self):
+ r'''Test of a single bunch simulation in the LHC at constant energy with the beam-phase loop'''
+ # TODO: finish implementation of the test with the beam-phase loop
+ line = self.initializeXsuiteLine()
+
+ def singleParticleBLonDSimulation(self, ramp: bool = False, two_rfstations: bool = False):
+ r'''
+ Method to generate a BLonD simulation in the LHC with a single particle with a time offset with respect to the
+ RF bucket.
+ '''
+ # Bunch parameters
+ N_m = 1 # Number of macroparticles [-]
+ N_p = 1.15e11 # Intensity
+
+ # Simulation parameters
+ N_t = 481 # Number of (tracked) turns [-]
+ dt_offset = 0.1 # Input particles dt [s]
+ v_part = 1
+
+ # --- BLonD objects ---
+ if ramp:
+ ring = Ring(self.C, self.alpha, np.linspace(self.p_s, self.p_f, N_t + 1), Proton(), n_turns=N_t)
+ else:
+ if two_rfstations:
+ # Additional input
+ c_part = 0.3
+ v_part = 0.5
+
+ ring = Ring(
+ [self.C * c_part, self.C * (1 - c_part)],
+ [[self.alpha], [self.alpha]],
+ [self.p_s, self.p_s],
+ Proton(), n_turns=N_t, n_sections=2
+ )
+ else:
+ ring = Ring(self.C, self.alpha, self.p_s, Proton(), n_turns=N_t)
+
+ beam = Beam(ring, N_m, N_p)
+ rfstation = RFStation(ring, [self.h], [self.V * v_part], [self.dphi])
+
+ # --- Insert particle ---
+ beam.dt = np.array([rfstation.t_rf[0, 0] * (1 - dt_offset)])
+ beam.dE = np.array([0])
+
+ if two_rfstations:
+ second_rfstation = RFStation(ring, [self.h], [self.V * (1 - v_part)],
+ [self.dphi], section_index=2)
+
+ return RingAndRFTracker(rfstation, beam), RingAndRFTracker(second_rfstation, beam), beam
+
+ return RingAndRFTracker(rfstation, beam), beam
+
+ def initializeXsuiteLine(self, two_rfstations: bool = False):
+ r'''
+ Method to generate an Xsuite Line using the one-turn matrix.
+ '''
+ c_part = 1
+
+ if two_rfstations:
+ c_part = 0.3
+ matrix_2 = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=self.alpha,
+ length=self.C * (1 - c_part)
+ )
+
+ matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=self.alpha,
+ length=self.C * c_part
+ )
+
+ # Create line
+ line = xt.Line(elements=[matrix], element_names={'matrix'})
+ line['matrix'].length = self.C * c_part
+
+ if two_rfstations:
+ line.append_element(element=matrix_2, name='matrix_2')
+ line['matrix_2'].length = (1 - c_part) * self.C
+
+ return line
+
+ def performSingleParticleInterfaceSimulation(self, line, blond_track, beam):
+ # Add particles to line and build tracker
+ line.particle_ref = xp.Particles(p0c=self.p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+ line.build_tracker()
+
+ # --- Convert the initial BLonD distribution to xsuite coordinates ---
+ zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=blond_track.rf_params.phi_s[0],
+ omega_rf=blond_track.rf_params.omega_rf[0, 0])
+
+ # --- Track matrix ---
+ N_t = len(blond_track.rf_params.phi_s) - 1
+ particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+ line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+ mon = line.record_last_track
+
+ dt_array = np.zeros(N_t)
+ dE_array = np.zeros(N_t)
+ # Convert the xsuite particle coordinates back to BLonD
+ for i in range(N_t):
+ dt_array[i], dE_array[i] = xsuite_coords_to_blond_coords(
+ mon.zeta[:, i].T, mon.ptau[:, i].T, blond_track.rf_params.beta[i], blond_track.rf_params.energy[i],
+ phi_s=blond_track.rf_params.phi_s[i] - blond_track.rf_params.phi_rf[0, 0],
+ omega_rf=blond_track.rf_params.omega_rf[0, i]
+ )
+
+ return dt_array, dE_array
+
+ @staticmethod
+ def performSingleParticleBLonDSimulation(blond_track, beam, second_blond_tracker: RingAndRFTracker = None):
+ # The number of turns to track
+ N_t = len(blond_track.rf_params.phi_s) - 1
+
+ # Save particle coordinates
+ dt_array = np.zeros(N_t)
+ dE_array = np.zeros(N_t)
+
+ if second_blond_tracker is not None:
+ full_tracker = FullRingAndRF([blond_track, second_blond_tracker])
+ else:
+ full_tracker = blond_track
+
+ for i in range(N_t):
+ full_tracker.track()
+ dt_array[i] = beam.dt[0]
+ dE_array[i] = beam.dE[0]
+
+ return dt_array, dE_array
+
+
+class TestXsuitePSB(unittest.TestCase):
+
+ def setUp(self):
+ # Accelerator parameters
+ self.C = 2 * np.pi * 25.00005794526065 # Machine circumference, radius 25 [m]
+ gamma_t = 4.11635447373496 # Transition gamma [-]
+ self.alpha = 1. / gamma_t / gamma_t # First order mom. comp. factor [-]
+ self.h = 1 # Harmonic number [-]
+ self.V = 8e3 # RF voltage [V]
+ self.dphi = np.pi # Phase modulation/offset [rad]
+
+ # Derived parameters
+ kin_energy = 1.4e9 # Kinetic energy [eV]
+ E_0 = m_p * c ** 2 / e # [eV]
+ tot_energy = E_0 + kin_energy # [eV]
+
+ self.p_s = np.sqrt(tot_energy ** 2 - E_0 ** 2)
+ self.p_f = self.p_s + 0.001e9
+
+ def testSingleParticle(self):
+ r'''Test of a single particle simulation in the PSB at constant energy'''
+ rtol = 1e-2 # relative tolerance
+ atol = 0 # absolute tolerance
+
+ # ----- Interface Simulation -----
+ # Initialize the BLonD cavity
+ blond_tracker, beam = self.singleParticleBLonDSimulation()
+ cavity = BlondElement(blond_tracker, beam)
+
+ line = self.initializeXsuiteLine()
+
+ # Insert the BLonD elements
+ line.insert_element(index=0, element=cavity, name='blond_cavity')
+
+ dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+
+ # ----- Pure BLonD Simulation -----
+ blond_tracker, beam = self.singleParticleBLonDSimulation()
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+
+ # ----- Perform test -----
+ np.testing.assert_allclose(
+ dt_inter, dt_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticle the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dE_inter, dE_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticle the energy-coordinate differs'
+ )
+
+ def testSingleParticleRamp(self):
+ r'''Test of acceleration of a single particle simulation in the PSB'''
+ rtol = 1e-2 # relative tolerance
+ atol = 0 # absolute tolerance
+
+ # Initialize the BLonD cavity
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ cavity = BlondElement(blond_tracker, beam)
+
+ line = self.initializeXsuiteLine()
+
+ # Insert the BLonD elements
+ line.insert_element(index=0, element=cavity, name='blond_cavity')
+
+ # Insert energy ramp
+ energy_update = EnergyUpdate(blond_tracker.rf_params.momentum)
+ line.insert_element(index='matrix', element=energy_update, name='energy_update')
+
+ dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+
+ # ----- Pure BLonD Simulation -----
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+
+ # ----- Perform test -----
+ np.testing.assert_allclose(
+ dt_inter, dt_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleRamp the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dE_inter, dE_blond, rtol=rtol, atol=atol,
+ err_msg='In testSingleParticleRamp the energy-coordinate differs'
+ )
+
+ def singleParticleBLonDSimulation(self, ramp: bool = False):
+ r'''
+ Method to generate a BLonD simulation in the LHC with a single particle with a time offset with respect to the
+ RF bucket.
+ '''
+ # Bunch parameters
+ N_m = 1 # Number of macroparticles [-]
+ N_p = 1e11 # Intensity
+
+ # Simulation parameters
+ N_t = 6000 # Number of (tracked) turns [-]
+ dt_offset = 0.1 # Input particles dt [s]
+
+ # --- BLonD objects ---
+ if ramp:
+ ring = Ring(self.C, self.alpha, np.linspace(self.p_s, self.p_f, N_t + 1), Proton(), n_turns=N_t)
+ else:
+ ring = Ring(self.C, self.alpha, self.p_s, Proton(), n_turns=N_t)
+
+ beam = Beam(ring, N_m, N_p)
+ rfstation = RFStation(ring, [self.h], [self.V], [self.dphi])
+
+ # --- Insert particle ---
+ beam.dt = np.array([rfstation.t_rf[0, 0] * (1 - dt_offset)])
+ beam.dE = np.array([0])
+
+ return RingAndRFTracker(rfstation, beam), beam
+
+ def initializeXsuiteLine(self):
+ r'''
+ Method to generate an Xsuite Line using the one-turn matrix.
+ '''
+ matrix = xt.LineSegmentMap(
+ longitudinal_mode='nonlinear',
+ qx=1.1, qy=1.2,
+ betx=1.,
+ bety=1.,
+ voltage_rf=0,
+ frequency_rf=0,
+ lag_rf=0,
+ momentum_compaction_factor=self.alpha,
+ length=self.C
+ )
+
+ # Create line
+ line = xt.Line(elements=[matrix], element_names={'matrix'})
+ line['matrix'].length = self.C
+
+ return line
+
+ def performSingleParticleInterfaceSimulation(self, line, blond_track, beam):
+ # Add particles to line and build tracker
+ line.particle_ref = xp.Particles(p0c=self.p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
+ line.build_tracker()
+
+ # --- Convert the initial BLonD distribution to xsuite coordinates ---
+ zeta, ptau = blond_beam_to_xsuite_coords(beam,
+ line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0],
+ phi_s=blond_track.rf_params.phi_s[0],
+ omega_rf=blond_track.rf_params.omega_rf[0, 0])
+
+ # --- Track matrix ---
+ N_t = len(blond_track.rf_params.phi_s) - 1
+ particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
+ line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+ mon = line.record_last_track
+
+ dt_array = np.zeros(N_t)
+ dE_array = np.zeros(N_t)
+ # Convert the xsuite particle coordinates back to BLonD
+ for i in range(N_t):
+ dt_array[i], dE_array[i] = xsuite_coords_to_blond_coords(
+ mon.zeta[:, i].T, mon.ptau[:, i].T, blond_track.rf_params.beta[i], blond_track.rf_params.energy[i],
+ phi_s=blond_track.rf_params.phi_s[i] - blond_track.rf_params.phi_rf[0, 0],
+ omega_rf=blond_track.rf_params.omega_rf[0, i]
+ )
+
+ return dt_array, dE_array
+
+ @staticmethod
+ def performSingleParticleBLonDSimulation(blond_track, beam):
+ # The number of turns to track
+ N_t = len(blond_track.rf_params.phi_s) - 1
+
+ # Save particle coordinates
+ dt_array = np.zeros(N_t)
+ dE_array = np.zeros(N_t)
+
+ full_tracker = blond_track
+
+ for i in range(N_t):
+ full_tracker.track()
+ dt_array[i] = beam.dt[0]
+ dE_array[i] = beam.dE[0]
+
+ return dt_array, dE_array
\ No newline at end of file
--
GitLab
From 104da833f958c11e4a3f2bcf24d660b5e617f87e Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Birk=20Karlsen-B=C3=A6ck?= <birk.karlsen.baeck@gmail.com>
Date: Fri, 1 Nov 2024 12:07:28 +0100
Subject: [PATCH 4/9] [unittests][interfaces][test xsuite] changed initial
energy to be float
---
unittests/interfaces/test_xsuite.py | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/unittests/interfaces/test_xsuite.py b/unittests/interfaces/test_xsuite.py
index be3770a7..dd8f26ff 100644
--- a/unittests/interfaces/test_xsuite.py
+++ b/unittests/interfaces/test_xsuite.py
@@ -195,7 +195,7 @@ class TestXsuiteLHC(unittest.TestCase):
# --- Insert particle ---
beam.dt = np.array([rfstation.t_rf[0, 0] * (1 - dt_offset)])
- beam.dE = np.array([0])
+ beam.dE = np.array([0.0])
if two_rfstations:
second_rfstation = RFStation(ring, [self.h], [self.V * (1 - v_part)],
@@ -407,7 +407,7 @@ class TestXsuitePSB(unittest.TestCase):
# --- Insert particle ---
beam.dt = np.array([rfstation.t_rf[0, 0] * (1 - dt_offset)])
- beam.dE = np.array([0])
+ beam.dE = np.array([0.0])
return RingAndRFTracker(rfstation, beam), beam
--
GitLab
From d47ae36d4b590e1f3084f76a1bff4afdd670dd1d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Birk=20Karlsen-B=C3=A6ck?= <birk.karlsen.baeck@gmail.com>
Date: Wed, 6 Nov 2024 23:08:29 +0100
Subject: [PATCH 5/9] [unittests][interfaces][test xsuite] Finished unit tests
of single particle in the LHC and PSB
---
unittests/interfaces/test_xsuite.py | 114 +++++++++++++++++++---------
1 file changed, 78 insertions(+), 36 deletions(-)
diff --git a/unittests/interfaces/test_xsuite.py b/unittests/interfaces/test_xsuite.py
index dd8f26ff..d9709927 100644
--- a/unittests/interfaces/test_xsuite.py
+++ b/unittests/interfaces/test_xsuite.py
@@ -27,6 +27,7 @@ from blond.beam.profile import CutOptions, Profile
from blond.input_parameters.rf_parameters import RFStation
from blond.input_parameters.ring import Ring
from blond.trackers.tracker import RingAndRFTracker, FullRingAndRF
+from blond.monitors.monitors import BunchMonitor
from blond.interfaces.xsuite import (BlondElement, BlondObserver,
EnergyUpdate,
@@ -58,8 +59,8 @@ class TestXsuiteLHC(unittest.TestCase):
def testSingleParticle(self):
r'''Test of a single particle simulation in the LHC at constant energy'''
- rtol = 1e-2 # relative tolerance
- atol = 0 # absolute tolerance
+ rtol = 4e-9 # relative tolerance
+ atol = 5e-9 # absolute tolerance
# ----- Interface Simulation -----
# Initialize the BLonD cavity
@@ -72,10 +73,14 @@ class TestXsuiteLHC(unittest.TestCase):
line.insert_element(index=0, element=cavity, name='blond_cavity')
dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+ dt_inter = dt_inter / np.max(dt_inter)
+ dE_inter = dE_inter / np.max(dE_inter)
# ----- Pure BLonD Simulation -----
blond_tracker, beam = self.singleParticleBLonDSimulation()
dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+ dt_blond = dt_blond / np.max(dt_blond)
+ dE_blond = dE_blond / np.max(dE_blond)
# ----- Perform test -----
np.testing.assert_allclose(
@@ -89,8 +94,8 @@ class TestXsuiteLHC(unittest.TestCase):
def testSingleParticleRamp(self):
r'''Test of acceleration of a single particle simulation in the LHC'''
- rtol = 1e-2 # relative tolerance
- atol = 0 # absolute tolerance
+ rtol = 4e-9 # relative tolerance
+ atol = 5e-9 # absolute tolerance
# Initialize the BLonD cavity
blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
@@ -103,13 +108,18 @@ class TestXsuiteLHC(unittest.TestCase):
# Insert energy ramp
energy_update = EnergyUpdate(blond_tracker.rf_params.momentum)
- line.insert_element(index='matrix', element=energy_update, name='energy_update')
+ line.insert_element(index='blond_cavity', element=energy_update, name='energy_update')
dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+ dt_inter = dt_inter / np.max(dt_inter)
+ dE_inter = dE_inter / np.max(dE_inter)
+ line.get_table().show()
# ----- Pure BLonD Simulation -----
- blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True, pure_blond=True)
dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+ dt_blond = dt_blond / np.max(dt_blond)
+ dE_blond = dE_blond / np.max(dE_blond)
# ----- Perform test -----
np.testing.assert_allclose(
@@ -123,8 +133,8 @@ class TestXsuiteLHC(unittest.TestCase):
def testSingleParticleTwoRFStations(self):
r'''Test of a single particle simulation in the LHC at constant energy with two RF stations'''
- rtol = 1e-2 # relative tolerance
- atol = 0 # absolute tolerance
+ rtol = 4e-9 # relative tolerance
+ atol = 5e-9 # absolute tolerance
# Initialize the BLonD cavities
blond_tracker_1, blond_tracker_2, beam = self.singleParticleBLonDSimulation(two_rfstations=True)
@@ -138,10 +148,14 @@ class TestXsuiteLHC(unittest.TestCase):
line.insert_element(index='matrix_2', element=cavity_2, name='blond_cavity_2')
dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker_1, beam)
+ dt_inter = dt_inter / np.max(dt_inter)
+ dE_inter = dE_inter / np.max(dE_inter)
# ----- Pure BLonD Simulation -----
- blond_tracker, beam = self.singleParticleBLonDSimulation(two_rfstations=True)
- dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+ blond_tracker_1, blond_tracker_2, beam = self.singleParticleBLonDSimulation(two_rfstations=True)
+ dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker_1, beam, blond_tracker_2)
+ dt_blond = dt_blond / np.max(dt_blond)
+ dE_blond = dE_blond / np.max(dE_blond)
# ----- Perform test -----
np.testing.assert_allclose(
@@ -153,12 +167,8 @@ class TestXsuiteLHC(unittest.TestCase):
err_msg='In testSingleParticleTwoRFStations the energy-coordinate differs'
)
- def testBunchWithPhaseLoop(self):
- r'''Test of a single bunch simulation in the LHC at constant energy with the beam-phase loop'''
- # TODO: finish implementation of the test with the beam-phase loop
- line = self.initializeXsuiteLine()
-
- def singleParticleBLonDSimulation(self, ramp: bool = False, two_rfstations: bool = False):
+ def singleParticleBLonDSimulation(self, ramp: bool = False, two_rfstations: bool = False,
+ pure_blond: bool = False):
r'''
Method to generate a BLonD simulation in the LHC with a single particle with a time offset with respect to the
RF bucket.
@@ -169,12 +179,17 @@ class TestXsuiteLHC(unittest.TestCase):
# Simulation parameters
N_t = 481 # Number of (tracked) turns [-]
- dt_offset = 0.1 # Input particles dt [s]
+ dt_offset = 0.15 # Input particles dt [s]
v_part = 1
# --- BLonD objects ---
if ramp:
- ring = Ring(self.C, self.alpha, np.linspace(self.p_s, self.p_f, N_t + 1), Proton(), n_turns=N_t)
+ cycle = np.linspace(self.p_s, self.p_f, N_t + 1)
+ if pure_blond:
+ cycle = np.concatenate(
+ (np.array([self.p_s]), np.linspace(self.p_s, self.p_f, N_t + 1)[:-1])
+ )
+ ring = Ring(self.C, self.alpha, cycle, Proton(), n_turns=N_t)
else:
if two_rfstations:
# Additional input
@@ -201,7 +216,9 @@ class TestXsuiteLHC(unittest.TestCase):
second_rfstation = RFStation(ring, [self.h], [self.V * (1 - v_part)],
[self.dphi], section_index=2)
- return RingAndRFTracker(rfstation, beam), RingAndRFTracker(second_rfstation, beam), beam
+ return (RingAndRFTracker(rfstation, beam),
+ RingAndRFTracker(second_rfstation, beam),
+ beam)
return RingAndRFTracker(rfstation, beam), beam
@@ -248,6 +265,9 @@ class TestXsuiteLHC(unittest.TestCase):
return line
def performSingleParticleInterfaceSimulation(self, line, blond_track, beam):
+ r'''
+ Perform simulation with the interface.
+ '''
# Add particles to line and build tracker
line.particle_ref = xp.Particles(p0c=self.p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
line.build_tracker()
@@ -262,7 +282,7 @@ class TestXsuiteLHC(unittest.TestCase):
# --- Track matrix ---
N_t = len(blond_track.rf_params.phi_s) - 1
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
- line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+ line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=False)
mon = line.record_last_track
dt_array = np.zeros(N_t)
@@ -279,6 +299,9 @@ class TestXsuiteLHC(unittest.TestCase):
@staticmethod
def performSingleParticleBLonDSimulation(blond_track, beam, second_blond_tracker: RingAndRFTracker = None):
+ r'''
+ Perform simulation in pure BLonD.
+ '''
# The number of turns to track
N_t = len(blond_track.rf_params.phi_s) - 1
@@ -292,9 +315,9 @@ class TestXsuiteLHC(unittest.TestCase):
full_tracker = blond_track
for i in range(N_t):
- full_tracker.track()
dt_array[i] = beam.dt[0]
dE_array[i] = beam.dE[0]
+ full_tracker.track()
return dt_array, dE_array
@@ -320,8 +343,8 @@ class TestXsuitePSB(unittest.TestCase):
def testSingleParticle(self):
r'''Test of a single particle simulation in the PSB at constant energy'''
- rtol = 1e-2 # relative tolerance
- atol = 0 # absolute tolerance
+ rtol = 4e-4 # relative tolerance
+ atol = 4.1e-4 # absolute tolerance
# ----- Interface Simulation -----
# Initialize the BLonD cavity
@@ -334,10 +357,14 @@ class TestXsuitePSB(unittest.TestCase):
line.insert_element(index=0, element=cavity, name='blond_cavity')
dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+ dt_inter = dt_inter / np.max(dt_inter)
+ dE_inter = dE_inter / np.max(dE_inter)
# ----- Pure BLonD Simulation -----
blond_tracker, beam = self.singleParticleBLonDSimulation()
dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+ dt_blond = dt_blond / np.max(dt_blond)
+ dE_blond = dE_blond / np.max(dE_blond)
# ----- Perform test -----
np.testing.assert_allclose(
@@ -351,8 +378,8 @@ class TestXsuitePSB(unittest.TestCase):
def testSingleParticleRamp(self):
r'''Test of acceleration of a single particle simulation in the PSB'''
- rtol = 1e-2 # relative tolerance
- atol = 0 # absolute tolerance
+ rtol = 4e-4 # relative tolerance
+ atol = 4.1e-4 # absolute tolerance
# Initialize the BLonD cavity
blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
@@ -368,10 +395,14 @@ class TestXsuitePSB(unittest.TestCase):
line.insert_element(index='matrix', element=energy_update, name='energy_update')
dt_inter, dE_inter = self.performSingleParticleInterfaceSimulation(line, blond_tracker, beam)
+ dt_inter = dt_inter / np.max(dt_inter)
+ dE_inter = dE_inter / np.max(dE_inter)
# ----- Pure BLonD Simulation -----
- blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True)
+ blond_tracker, beam = self.singleParticleBLonDSimulation(ramp=True, pure_blond=True)
dt_blond, dE_blond = self.performSingleParticleBLonDSimulation(blond_tracker, beam)
+ dt_blond = dt_blond / np.max(dt_blond)
+ dE_blond = dE_blond / np.max(dE_blond)
# ----- Perform test -----
np.testing.assert_allclose(
@@ -383,7 +414,7 @@ class TestXsuitePSB(unittest.TestCase):
err_msg='In testSingleParticleRamp the energy-coordinate differs'
)
- def singleParticleBLonDSimulation(self, ramp: bool = False):
+ def singleParticleBLonDSimulation(self, ramp: bool = False, pure_blond: bool = False):
r'''
Method to generate a BLonD simulation in the LHC with a single particle with a time offset with respect to the
RF bucket.
@@ -394,11 +425,16 @@ class TestXsuitePSB(unittest.TestCase):
# Simulation parameters
N_t = 6000 # Number of (tracked) turns [-]
- dt_offset = 0.1 # Input particles dt [s]
+ dt_offset = 0.15 # Input particles dt [s]
# --- BLonD objects ---
if ramp:
- ring = Ring(self.C, self.alpha, np.linspace(self.p_s, self.p_f, N_t + 1), Proton(), n_turns=N_t)
+ cycle = np.linspace(self.p_s, self.p_f, N_t + 1)
+ if pure_blond:
+ cycle = np.concatenate(
+ (np.array([self.p_s]), np.linspace(self.p_s, self.p_f, N_t + 1)[:-1])
+ )
+ ring = Ring(self.C, self.alpha, cycle, Proton(), n_turns=N_t)
else:
ring = Ring(self.C, self.alpha, self.p_s, Proton(), n_turns=N_t)
@@ -434,21 +470,24 @@ class TestXsuitePSB(unittest.TestCase):
return line
def performSingleParticleInterfaceSimulation(self, line, blond_track, beam):
+ r'''
+ Perform simulation with the interface.
+ '''
# Add particles to line and build tracker
line.particle_ref = xp.Particles(p0c=self.p_s, mass0=xp.PROTON_MASS_EV, q0=1.)
line.build_tracker()
# --- Convert the initial BLonD distribution to xsuite coordinates ---
- zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=blond_track.rf_params.phi_s[0],
- omega_rf=blond_track.rf_params.omega_rf[0, 0])
+ zeta, ptau = blond_beam_to_xsuite_coords(
+ beam, line.particle_ref.beta0[0], line.particle_ref.energy0[0],
+ phi_s=blond_track.rf_params.phi_s[0] - blond_track.rf_params.phi_rf[0, 0],
+ omega_rf=blond_track.rf_params.omega_rf[0, 0]
+ )
# --- Track matrix ---
N_t = len(blond_track.rf_params.phi_s) - 1
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
- line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=True)
+ line.track(particles, num_turns=N_t, turn_by_turn_monitor=True, with_progress=False)
mon = line.record_last_track
dt_array = np.zeros(N_t)
@@ -465,6 +504,9 @@ class TestXsuitePSB(unittest.TestCase):
@staticmethod
def performSingleParticleBLonDSimulation(blond_track, beam):
+ r'''
+ Perform simulation in pure BLonD.
+ '''
# The number of turns to track
N_t = len(blond_track.rf_params.phi_s) - 1
@@ -475,8 +517,8 @@ class TestXsuitePSB(unittest.TestCase):
full_tracker = blond_track
for i in range(N_t):
- full_tracker.track()
dt_array[i] = beam.dt[0]
dE_array[i] = beam.dE[0]
+ full_tracker.track()
return dt_array, dE_array
\ No newline at end of file
--
GitLab
From 838a615a1c78c7345f0f9afe0b2159586c1387c4 Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Thu, 7 Nov 2024 10:08:49 +0100
Subject: [PATCH 6/9] Added xsuite to the Dockerfile
---
.../docker/images/python310-cuda/Dockerfile | 33 +++++++++++++++++--
1 file changed, 31 insertions(+), 2 deletions(-)
diff --git a/developer_tools/docker/images/python310-cuda/Dockerfile b/developer_tools/docker/images/python310-cuda/Dockerfile
index 92547389..e7c492a5 100644
--- a/developer_tools/docker/images/python310-cuda/Dockerfile
+++ b/developer_tools/docker/images/python310-cuda/Dockerfile
@@ -1,9 +1,29 @@
-ENTRYPOINT["/bin/bash"]
FROM python:3.10
+# Create a non-root user
+RUN useradd -m non_root
+
+# Switch to root to install packages
+USER root
+
# Set the environment variable to non-interactive
ENV DEBIAN_FRONTEND=noninteractive
+#######################################################################################
+# Obtaining CERN rf-noise-cpp
+#######################################################################################
+# requires read api token !
+# 1. obtained token at https://gitlab.cern.ch/-/user_settings/personal_access_tokens
+# 2. export GITLAB_READ_TOKEN=<YOUR-TOKEN>
+# 3. Build your docker container
+RUN --mount=type=secret,id=gitlab_token \
+ GITLAB_READ_TOKEN=$(cat /run/secrets/gitlab_token) && \
+ git clone https://oauth2:$GITLAB_READ_TOKEN@gitlab.cern.ch/be-rf-cs/Tools-and-libs/rf-noise-cpp.git/ /usr/local/rf-noise-cpp
+RUN apt-get update && \
+ apt-get install -y libboost-all-dev
+ENV RF_NOISE_DIR=/usr/local/rf-noise-cpp
+#######################################################################################
+
RUN apt-get update && \
apt-get install -y gcc g++ build-essential mpich libmpich-dev
RUN wget https://developer.download.nvidia.com/compute/cuda/repos/debian12/x86_64/cuda-keyring_1.1-1_all.deb && \
@@ -13,6 +33,7 @@ RUN wget https://developer.download.nvidia.com/compute/cuda/repos/debian12/x86_6
apt-get install -y cuda-toolkit && \
rm -rf /var/lib/apt/lists/*
+
# Make sure NVCC is accessible
ENV CUDA_HOME=/usr/local/cuda
ENV PATH=${CUDA_HOME}/bin:${PATH}
@@ -20,5 +41,13 @@ ENV LD_LIBRARY_PATH=${CUDA_HOME}/lib64:$LD_LIBRARY_PATH
RUN bash -c 'python3 -m pip install --upgrade pip&& \
- python3 -m pip install cupy-cuda12x pytest_cov mpi4py'
+ python3 -m pip install cupy-cuda12x pytest_cov mpi4py xsuite'
+
+# Switch to the non-root user
+USER non_root
+
+# Set the working directory (optional, adjust as needed)
+WORKDIR /home/non_root
+# Set the entrypoint
+ENTRYPOINT ["/bin/bash"]
\ No newline at end of file
--
GitLab
From 6e2fc298498919dd78b04cdc38471927e4500c7b Mon Sep 17 00:00:00 2001
From: Birk Emil Karlsen-Baeck <birk.beck@cern.ch>
Date: Fri, 31 Jan 2025 13:00:41 +0100
Subject: [PATCH 7/9] Renamed the transformation functions and changed the
examples accordingly
---
..._single_particle_with_oneturnmatrix_lhc.py | 24 +++++++++---------
.../EX_24_single_particle_ramp_lhc.py | 25 ++++++++++---------
...X_25_single_particle_two_rfstations_lhc.py | 24 +++++++++---------
...EX_26_single_particle_oneturnmatrix_psb.py | 25 ++++++++++---------
.../EX_27_single_particle_ramp_psb.py | 25 ++++++++++---------
__EXAMPLES/main_files/EX_28_phase_loop_lhc.py | 25 ++++++++++---------
blond/interfaces/xsuite.py | 19 ++++++++------
7 files changed, 87 insertions(+), 80 deletions(-)
diff --git a/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
index 9052aab0..6edecffd 100644
--- a/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
+++ b/__EXAMPLES/main_files/EX_23_single_particle_with_oneturnmatrix_lhc.py
@@ -23,8 +23,8 @@ from blond.trackers.tracker import RingAndRFTracker
# Interface
from blond.interfaces.xsuite import (BlondElement, BlondObserver,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -149,11 +149,10 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation.phi_s[0],
- omega_rf=rfstation.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0], line.particle_ref.energy0[0],
+ phi_s=rfstation.phi_s[0], omega_rf=rfstation.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -179,11 +178,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
- rfstation.beta[i],
- rfstation.energy[i],
- phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i], rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
index a36d3809..cfaa4644 100644
--- a/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
+++ b/__EXAMPLES/main_files/EX_24_single_particle_ramp_lhc.py
@@ -24,8 +24,8 @@ from blond.trackers.tracker import RingAndRFTracker
# Interface
from blond.interfaces.xsuite import (BlondElement, BlondObserver, EnergyUpdate,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -157,11 +157,11 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation.phi_s[0],
- omega_rf=rfstation.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0], phi_s=rfstation.phi_s[0],
+ omega_rf=rfstation.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -188,11 +188,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
- rfstation.beta[i],
- rfstation.energy[i],
- phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i], rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
index 41e6e88f..9103a32f 100644
--- a/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
+++ b/__EXAMPLES/main_files/EX_25_single_particle_two_rfstations_lhc.py
@@ -22,8 +22,8 @@ from blond.trackers.tracker import RingAndRFTracker
# Interface
from blond.interfaces.xsuite import (BlondElement,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -176,11 +176,10 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation_1.phi_s[0],
- omega_rf=rfstation_1.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0], line.particle_ref.energy0[0],
+ phi_s=rfstation_1.phi_s[0], omega_rf=rfstation_1.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -207,11 +206,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:,i].T, mon.ptau[:,i].T,
- rfstation_1.beta[i],
- rfstation_1.energy[i],
- phi_s=rfstation_1.phi_s[i] - rfstation_1.phi_rf[0, 0],
- omega_rf=rfstation_1.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:,i].T, mon.ptau[:,i].T,
+ rfstation_1.beta[i], rfstation_1.energy[i],
+ phi_s=rfstation_1.phi_s[i] - rfstation_1.phi_rf[0, 0],
+ omega_rf=rfstation_1.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
index 6ea893ad..845ffa3c 100644
--- a/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
+++ b/__EXAMPLES/main_files/EX_26_single_particle_oneturnmatrix_psb.py
@@ -23,8 +23,8 @@ from blond.trackers.tracker import RingAndRFTracker
# Interface
from blond.interfaces.xsuite import (BlondElement,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -153,11 +153,11 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0], phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -184,11 +184,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
- rfstation.beta[i],
- rfstation.energy[i],
- phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i], rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
index b4707f0f..dd82280e 100644
--- a/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
+++ b/__EXAMPLES/main_files/EX_27_single_particle_ramp_psb.py
@@ -23,8 +23,8 @@ from blond.trackers.tracker import RingAndRFTracker
# Interface
from blond.interfaces.xsuite import (BlondElement, EnergyUpdate,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -159,11 +159,11 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0], phi_s=rfstation.phi_s[0] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -190,11 +190,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
- rfstation.beta[i],
- rfstation.energy[i],
- phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i], rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
index cefbbabd..bcb349ff 100644
--- a/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
+++ b/__EXAMPLES/main_files/EX_28_phase_loop_lhc.py
@@ -25,8 +25,8 @@ from blond.llrf.beam_feedback import BeamFeedback
# Interface
from blond.interfaces.xsuite import (BlondElement, BlondObserver,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# Monitor objects
from blond.monitors.monitors import BunchMonitor
@@ -167,11 +167,11 @@ line.get_table().show()
print(f"\nSetting up simulation...")
# --- Convert the initial BLonD distribution to xsuite coordinates ---
-zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=rfstation.phi_s[0],
- omega_rf=rfstation.omega_rf[0, 0])
+zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0], phi_s=rfstation.phi_s[0],
+ omega_rf=rfstation.omega_rf[0, 0]
+)
# --- Track matrix ---
particles = line.build_particles(x=0, y=0, px=0, py=0, zeta=np.copy(zeta), ptau=np.copy(ptau))
@@ -198,11 +198,12 @@ xsuite_string += '{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\n'.
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt, dE = xsuite_coords_to_blond_coords(mon.zeta[:, i].T, mon.ptau[:, i].T,
- rfstation.beta[i],
- rfstation.energy[i],
- phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
- omega_rf=rfstation.omega_rf[0, i])
+ dt, dE = xsuite_to_blond_transform(
+ mon.zeta[:, i].T, mon.ptau[:, i].T,
+ rfstation.beta[i], rfstation.energy[i],
+ phi_s=rfstation.phi_s[i] - rfstation.phi_rf[0, 0], # Correct for RF phase
+ omega_rf=rfstation.omega_rf[0, i]
+ )
# Statistics
test_string += ('{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}\t{:+10.10e}'
diff --git a/blond/interfaces/xsuite.py b/blond/interfaces/xsuite.py
index 2a4260e1..bd30be38 100644
--- a/blond/interfaces/xsuite.py
+++ b/blond/interfaces/xsuite.py
@@ -21,13 +21,15 @@ if TYPE_CHECKING:
from typing import Any, Union
-def blond_beam_to_xsuite_coords(beam: Beam, beta0: float, energy0: float,
- omega_rf: float, phi_s: float = 0):
+def blond_to_xsuite_transform(dt: Union[float, np.ndarray], de: Union[float, np.ndarray], beta0: float,
+ energy0: float, omega_rf: float, phi_s: float = 0):
r'''
Coordinate transformation from Xsuite to BLonD at a given turn or multiple turns if numpy arrays ar given.
- :param beam: blond.beam.beam.Beam class.
- BLonD beam class to extract the initial coordinates dE and dt.
+ :param dt: float or NDArray.
+ The deviation in time from the reference clock in BLonD.
+ :param de: float or NDArray.
+ The deviation in energy from the synchronous particle.
:param beta0: float.
Synchronous beta.
:param energy0: float.
@@ -40,13 +42,14 @@ def blond_beam_to_xsuite_coords(beam: Beam, beta0: float, energy0: float,
:return zeta, ptau: as numpy-arrays (or single variable)
'''
- ptau = beam.dE / (beta0 * energy0)
- zeta = -(beam.dt - phi_s / omega_rf) * beta0 * clight
+
+ ptau = de / (beta0 * energy0)
+ zeta = -(dt - phi_s / omega_rf) * beta0 * clight
return zeta, ptau
-def xsuite_coords_to_blond_coords(zeta: Union[float, np.ndarray], ptau: Union[float, np.ndarray],
- beta0: float, energy0: float, omega_rf: float, phi_s: float = 0):
+def xsuite_to_blond_transform(zeta: Union[float, np.ndarray], ptau: Union[float, np.ndarray],
+ beta0: float, energy0: float, omega_rf: float, phi_s: float = 0):
r'''
Coordinate transformation from Xsuite to BLonD.
--
GitLab
From 03a13d1b23f0ff7567295b9c5dbd05ed30b4419a Mon Sep 17 00:00:00 2001
From: Birk Emil Karlsen-Baeck <birk.beck@cern.ch>
Date: Fri, 31 Jan 2025 13:01:33 +0100
Subject: [PATCH 8/9] [unittests][interfaces][xsuite] Added unittests for the
coordinate transformations between BLonD and Xsuite
---
unittests/interfaces/test_xsuite.py | 135 ++++++++++++++++++++++++----
1 file changed, 120 insertions(+), 15 deletions(-)
diff --git a/unittests/interfaces/test_xsuite.py b/unittests/interfaces/test_xsuite.py
index d9709927..c5a0e57d 100644
--- a/unittests/interfaces/test_xsuite.py
+++ b/unittests/interfaces/test_xsuite.py
@@ -31,19 +31,124 @@ from blond.monitors.monitors import BunchMonitor
from blond.interfaces.xsuite import (BlondElement, BlondObserver,
EnergyUpdate,
- blond_beam_to_xsuite_coords,
- xsuite_coords_to_blond_coords)
-
-this_directory = os.path.dirname(os.path.realpath(__file__))
+ blond_to_xsuite_transform,
+ xsuite_to_blond_transform)
# TODO: Make a test of the EnergyFrequencyUpdate and BlondObserver classes
class TestXsuiteBLonDTransforms(unittest.TestCase):
- # TODO: finish implementation of the test of the particle coordinates
+
def setUp(self):
- pass
+ # Accelerator parameters
+ C = 26658.8832
+ p_s = 450e9
+ h = 35640
+ alpha = 0.00034849575112251314
+ V = 5e6
+ dphi = 0
+
+ # BLonD ring and RF station
+ self.ring = Ring(C, alpha, p_s, Proton(), n_turns=1)
+ self.rfstation = RFStation(self.ring, [h], [V], [dphi])
+
+ def testBlondToXsuiteTransform(self):
+ # Arrays of time and energy offsets
+ dts = np.linspace(
+ self.rfstation.t_rf[0, 0] * 1e-3,
+ self.rfstation.t_rf[0, 0], 10
+ )
+ dEs = np.linspace(-200e6, 200e6, 10)
+
+ # Transform BLonD coordinates to Xsuite
+ zeta, ptau = blond_to_xsuite_transform(
+ dts, dEs, self.ring.beta[0, 0], self.ring.energy[0, 0],
+ self.rfstation.omega_rf[0, 0], self.rfstation.phi_s[0]
+ )
+
+ # Correct values
+ act_zeta = np.array(
+ [0.37325428, 0.29022578, 0.20719727, 0.12416876, 0.04114025,
+ -0.04188826, -0.12491676, -0.20794527, -0.29097378, -0.37400229]
+ )
+ act_ptau = np.array(
+ [-4.44444444e-04, -3.45679012e-04, -2.46913580e-04, -1.48148148e-04, -4.93827160e-05,
+ 4.93827160e-05, 1.48148148e-04, 2.46913580e-04, 3.45679012e-04, 4.44444444e-04]
+ )
+
+ # Check if the coordinates are the same
+ np.testing.assert_allclose(
+ act_zeta, zeta,
+ err_msg='In testBlondToXsuiteTransform the zeta-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ act_ptau, ptau,
+ err_msg='In testBlondToXsuiteTransform the ptau-coordinate differs'
+ )
+
+ def testXsuiteToBlondTransform(self):
+ # Arrays of z and ptau from xsuite
+ zeta = np.linspace(
+ -0.37, 0.37, 10
+ )
+ ptau = np.linspace(
+ -4e-4, 4e-4, 10
+ )
+ # Transform BLonD coordinates to Xsuite
+ dt, de = xsuite_to_blond_transform(
+ zeta, ptau, self.ring.beta[0, 0], self.ring.energy[0, 0],
+ self.rfstation.omega_rf[0, 0], self.rfstation.phi_s[0]
+ )
+
+ # Correct values
+ act_dt = np.array(
+ [2.48172990e-09, 2.20746549e-09, 1.93320108e-09, 1.65893668e-09, 1.38467227e-09,
+ 1.11040786e-09, 8.36143453e-10, 5.61879046e-10, 2.87614638e-10, 1.33502300e-11]
+ )
+ act_de = np.array(
+ [-1.8e+08, -1.4e+08, -1.0e+08, -6.0e+07, -2.0e+07,
+ 2.0e+07, 6.0e+07, 1.0e+08, 1.4e+08, 1.8e+08]
+ )
+
+ # Check if the coordinates are the same
+ np.testing.assert_allclose(
+ act_dt, dt,
+ err_msg='In testXsuiteToBlondTransform the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ act_de, de,
+ err_msg='In testXsuiteToBlondTransform the energy-coordinate differs'
+ )
+
+ def testTransformInverse(self):
+ # Arrays of time and energy offsets
+ dts = np.linspace(
+ self.rfstation.t_rf[0, 0] * 1e-3,
+ self.rfstation.t_rf[0, 0], 20
+ )
+ dEs = np.linspace(-200e6, 200e6, 20)
+
+ # Perform transforms to xsuite coordinates and back
+ zeta, ptau = blond_to_xsuite_transform(
+ dts, dEs, self.ring.beta[0, 0], self.ring.energy[0, 0],
+ self.rfstation.omega_rf[0, 0], self.rfstation.phi_s[0]
+ )
+
+ _dts, _dEs = xsuite_to_blond_transform(
+ zeta, ptau, self.ring.beta[0, 0], self.ring.energy[0, 0],
+ self.rfstation.omega_rf[0, 0], self.rfstation.phi_s[0]
+ )
+
+ # Check if the coordinates are the same
+ np.testing.assert_allclose(
+ dts, _dts,
+ err_msg='In testTransformInverse the time-coordinate differs'
+ )
+ np.testing.assert_allclose(
+ dEs, _dEs,
+ err_msg='In testTransformInverse the energy-coordinate differs'
+ )
class TestXsuiteLHC(unittest.TestCase):
@@ -273,11 +378,11 @@ class TestXsuiteLHC(unittest.TestCase):
line.build_tracker()
# --- Convert the initial BLonD distribution to xsuite coordinates ---
- zeta, ptau = blond_beam_to_xsuite_coords(beam,
- line.particle_ref.beta0[0],
- line.particle_ref.energy0[0],
- phi_s=blond_track.rf_params.phi_s[0],
- omega_rf=blond_track.rf_params.omega_rf[0, 0])
+ zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0],
+ line.particle_ref.energy0[0], phi_s=blond_track.rf_params.phi_s[0],
+ omega_rf=blond_track.rf_params.omega_rf[0, 0]
+ )
# --- Track matrix ---
N_t = len(blond_track.rf_params.phi_s) - 1
@@ -289,7 +394,7 @@ class TestXsuiteLHC(unittest.TestCase):
dE_array = np.zeros(N_t)
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt_array[i], dE_array[i] = xsuite_coords_to_blond_coords(
+ dt_array[i], dE_array[i] = xsuite_to_blond_transform(
mon.zeta[:, i].T, mon.ptau[:, i].T, blond_track.rf_params.beta[i], blond_track.rf_params.energy[i],
phi_s=blond_track.rf_params.phi_s[i] - blond_track.rf_params.phi_rf[0, 0],
omega_rf=blond_track.rf_params.omega_rf[0, i]
@@ -478,8 +583,8 @@ class TestXsuitePSB(unittest.TestCase):
line.build_tracker()
# --- Convert the initial BLonD distribution to xsuite coordinates ---
- zeta, ptau = blond_beam_to_xsuite_coords(
- beam, line.particle_ref.beta0[0], line.particle_ref.energy0[0],
+ zeta, ptau = blond_to_xsuite_transform(
+ beam.dt, beam.dE, line.particle_ref.beta0[0], line.particle_ref.energy0[0],
phi_s=blond_track.rf_params.phi_s[0] - blond_track.rf_params.phi_rf[0, 0],
omega_rf=blond_track.rf_params.omega_rf[0, 0]
)
@@ -494,7 +599,7 @@ class TestXsuitePSB(unittest.TestCase):
dE_array = np.zeros(N_t)
# Convert the xsuite particle coordinates back to BLonD
for i in range(N_t):
- dt_array[i], dE_array[i] = xsuite_coords_to_blond_coords(
+ dt_array[i], dE_array[i] = xsuite_to_blond_transform(
mon.zeta[:, i].T, mon.ptau[:, i].T, blond_track.rf_params.beta[i], blond_track.rf_params.energy[i],
phi_s=blond_track.rf_params.phi_s[i] - blond_track.rf_params.phi_rf[0, 0],
omega_rf=blond_track.rf_params.omega_rf[0, i]
--
GitLab
From cf5d4ba55c7cff2e4233a4b3c678226aa1e7b073 Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Tue, 11 Feb 2025 15:01:55 +0100
Subject: [PATCH 9/9] Fixed CI to not require admin privileges
---
.gitlab-ci.yml | 1 -
1 file changed, 1 deletion(-)
diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index 364f6330..1d021d2a 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -12,7 +12,6 @@ variables:
.before:
before_script:
- - apt-get update
# Verify CUDA installation
- nvidia-smi
- gcc --version
--
GitLab