A Single Bend with ISR
This test illustrates the effects of incoherent synchrotron radiation (ISR) on a high-energy electron bunch in a bending dipole. The effects of CSR are turned off.
An initially cold (zero emittance), 100 GeV electron bunch with a 0.1 mm rms beam size propagates in an ideal sector bend, through a region of uniform magnetic field.
Due to the stochastic emission of synchrotron radiation, the bunch experiences a mean loss of energy as well as a growth in rms energy spread. Due to the presence of nonzero dispersion, this also results in a growth of the bunch rms emittance.
The resulting quantities are compared against the expressions given in:
N. Yampolsky and B. Carlsten, “Beam debunching due to ISR-induced energy diffusion,” Nucl. Instrum. and Methods in Phys. Res. A, 870, 156-162 (2017), equations (17-20).
In this test, the values of \(\langle{p_x\rangle}\), \(\langle{p_y\rangle}\), \(\langle{p_t\rangle}\), \(\sigma_{p_x}\), \(\sigma_{p_y}\), and \(\sigma_{p_t}\) must agree with nominal values.
In addition, the final values of \(\langle{p_t\rangle}\) and \(\sigma_{p_t}\) must agree with predicted values.
Run
This example can be run either as:
Python script:
python3 run_bend_isr.pyorImpactX executable using an input file:
impactx input_bend_isr.in
For MPI-parallel runs, prefix these lines with mpiexec -n 4 ... or srun -n 4 ..., depending on the system.
#!/usr/bin/env python3
#
# Copyright 2022-2023 ImpactX contributors
# Authors: Axel Huebl, Chad Mitchell
# License: BSD-3-Clause-LBNL
#
# -*- coding: utf-8 -*-
from impactx import ImpactX, distribution, elements
sim = ImpactX()
# set numerical parameters and IO control
sim.space_charge = False
sim.isr = True
sim.isr_order = 1
sim.slice_step_diagnostics = True
# domain decomposition & space charge mesh
sim.init_grids()
# load a 2 GeV electron beam with an initial
# unnormalized rms emittance of 2 nm
kin_energy_MeV = 100.0e3 # reference energy
bunch_charge_C = 1.0e-12 # used with space charge
npart = 100000 # number of macro particles
# reference particle
ref = sim.particle_container().ref_particle()
ref.set_charge_qe(-1.0).set_mass_MeV(0.510998950).set_kin_energy_MeV(kin_energy_MeV)
# particle bunch
distr = distribution.Gaussian(
lambdaX=1.0e-4,
lambdaY=1.0e-4,
lambdaT=1.0e-4,
lambdaPx=0.0,
lambdaPy=0.0,
lambdaPt=0.0,
muxpx=0.0,
muypy=0.0,
mutpt=0.0,
)
sim.add_particles(bunch_charge_C, distr, npart)
# add beam diagnostics
monitor = elements.BeamMonitor("monitor", backend="h5")
# design the accelerator lattice)
ns = 25 # number of slices per ds in the element
bend = elements.ExactSbend(
name="sbend_exact", ds=0.1, phi=0.5729577951308232, nslice=ns
)
lattice_line = [monitor, bend, monitor]
# define the lattice
sim.lattice.extend(lattice_line)
# run simulation
sim.track_particles()
# clean shutdown
sim.finalize()
###############################################################################
# Particle Beam(s)
###############################################################################
beam.npart = 100000
beam.units = static
beam.kin_energy = 100.0e3 # 100 GeV nominal energy
beam.charge = 1.0e-12 # 1 pC charge
beam.particle = electron
beam.distribution = gaussian
beam.lambdaX = 1.0e-4
beam.lambdaY = beam.lambdaX
beam.lambdaT = 1.0e-4
beam.lambdaPx = 0.0
beam.lambdaPy = 0.0
beam.lambdaPt = 0.0
beam.muxpx = 0.0
beam.muypy = 0.0
beam.mutpt = 0.0
###############################################################################
# Beamline: lattice elements and segments
###############################################################################
lattice.elements = monitor bend monitor
lattice.nslice = 25
monitor.type = beam_monitor
monitor.backend = h5
bend.type = sbend_exact
bend.ds = 0.1 #used only to parameterize reference arc length
bend.phi = 0.5729577951308232 #deg
#bend.type = sbend
#bend.ds = 0.1 #used only to parameterize reference arc length
#bend.rc = 10.0 #0.5729577951308232 deg
###############################################################################
# Algorithms
###############################################################################
algo.space_charge = false
algo.isr = true
algo.isr_order = 1
###############################################################################
# Diagnostics
###############################################################################
diag.slice_step_diagnostics = true
diag.eigenemittances = true
Analyze
We run the following script to analyze correctness:
Script analysis_bend_isr.py
#!/usr/bin/env python3
#
# Copyright 2022-2023 ImpactX contributors
# Authors: Axel Huebl, Chad Mitchell
# License: BSD-3-Clause-LBNL
#
import glob
import re
import numpy as np
import pandas as pd
def read_file(file_pattern):
for filename in glob.glob(file_pattern):
df = pd.read_csv(filename, delimiter=r"\s+")
if "step" not in df.columns:
step = int(re.findall(r"[0-9]+", filename)[0])
df["step"] = step
yield df
def read_time_series(file_pattern):
"""Read in all CSV files from each MPI rank (and potentially OpenMP
thread). Concatenate into one Pandas dataframe.
Returns
-------
pandas.DataFrame
"""
return pd.concat(
read_file(file_pattern),
axis=0,
ignore_index=True,
) # .set_index('id')
# read reduced diagnostics
rbc = read_time_series("diags/reduced_beam_characteristics.*")
s = rbc["s"]
px_mean = rbc["mean_px"]
py_mean = rbc["mean_py"]
pt_mean = rbc["mean_pt"]
sigma_px = rbc["sigma_px"]
sigma_py = rbc["sigma_py"]
sigma_pt = rbc["sigma_pt"]
px_meani = px_mean.iloc[0]
py_meani = py_mean.iloc[0]
pt_meani = pt_mean.iloc[0]
sig_pxi = sigma_px.iloc[0]
sig_pyi = sigma_py.iloc[0]
sig_pti = sigma_pt.iloc[0]
length = len(s) - 1
sf = s.iloc[length]
px_meanf = px_mean.iloc[length]
py_meanf = py_mean.iloc[length]
pt_meanf = pt_mean.iloc[length]
sig_pxf = sigma_px.iloc[length]
sig_pyf = sigma_py.iloc[length]
sig_ptf = sigma_pt.iloc[length]
print("Initial Beam:")
print(f" px_mean={px_meani:e} py_mean={py_meani:e} pt_mean={pt_meani:e}")
print(f" sig_px={sig_pxi:e} sig_py={sig_pyi:e} sig_pt={sig_pti:e}")
atol = 1.0e-6
print(f" atol={atol}")
assert np.allclose(
[px_meani, py_meani, pt_meani, sig_pxi, sig_pyi, sig_pti],
[0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
atol=atol,
)
# Physical constants:
re_classical = 2.8179403205e-15 # classical electron radius
lambda_compton_reduced = 3.8615926744e-13 # reduced Compton wavelength
# Problem parameters:
ds = 0.1
rc = 10.0
gamma = 195696.117901
num_particles = 10000
# Characteristic length of energy loss:
length_isr = 3.0 / 2.0 * rc**2 / (re_classical * gamma**3)
print("")
print("Length and characteristic length for energy loss [m]:")
print(f" Length={ds:e} Length_ISR={length_isr:e}")
# Predicted energy loss and energy spread:
dpt = 2.0 / 3.0 * re_classical * ds / (rc**2) * gamma**3
dsigpt2 = (
55.0
/ (24 * np.sqrt(3.0))
* lambda_compton_reduced
* re_classical
* ds
/ rc**3
* gamma**5
)
dsigpt = np.sqrt(dsigpt2)
print("")
print("Final Beam:")
print(f" pt_mean={pt_meanf:e} sig_pt={sig_ptf}")
print("Predicted (assuming that Length << Length_isr):")
print(f" pt_mean={dpt:e} sig_pt={dsigpt:e}")
print("")
rtol = 10.0 * num_particles**-0.5 # from random sampling of a smooth distribution
assert np.allclose(
[pt_meanf, sig_ptf],
[
dpt,
dsigpt,
],
rtol=rtol,
atol=atol,
)