Parameters: Python

This documents on how to use ImpactX as a Python script (python3 run_script.py).

General

class impactx.ImpactX

This is the central simulation class.

property particle_shape

Control the particle B-spline order.

The order of the shape factors (splines) for the macro-particles along all spatial directions: 1 for linear, 2 for quadratic, 3 for cubic. Low-order shape factors result in faster simulations, but may lead to more noisy results. High-order shape factors are computationally more expensive, but may increase the overall accuracy of the results. For production runs it is generally safer to use high-order shape factors, such as cubic order.

Parameters

order (int) – B-spline order 1, 2, or 3

property n_cell

The number of grid points along each direction on the coarsest level.

property domain

The physical extent of the full simulation domain, relative to the reference particle position, in meters. When set, turns dynamic_size to False.

Note: particles that move outside the simulation domain are removed.

property prob_relative

The field mesh is expanded beyond the physical extent of particles by this factor. When set, turns dynamic_size to True.

property dynamic_size

Use dynamic (True) resizing of the field mesh or static sizing (False).

property space_charge

Enable (True) or disable (False) space charge calculations (default: True).

Whether to calculate space charge effects. This is in-development. At the moment, this flag only activates coordinate transformations and charge deposition.

property diagnostics

Enable (True) or disable (False) diagnostics generally (default: True). Disabling this is mostly used for benchmarking.

property slice_step_diagnostics

Enable (True) or disable (False) diagnostics every slice step in elements (default: True).

By default, diagnostics is performed at the beginning and end of the simulation. Enabling this flag will write diagnostics every step and slice step.

property diag_file_min_digits

The minimum number of digits (default: 6) used for the step number appended to the diagnostic file names.

init_grids()

Initialize AMReX blocks/grids for domain decomposition & space charge mesh.

This must come first, before particle beams and lattice elements are initialized.

add_particles(charge_C, distr, npart)

Generate and add n particles to the particle container. Note: Set the reference particle properties (charge, mass, energy) first.

Will also resize the geometry based on the updated particle distribution’s extent and then redistribute particles in according AMReX grid boxes.

Parameters
  • charge_C (float) – bunch charge (C)

  • distr – distribution function to draw from (object from impactx.distribution)

  • npart (int) – number of particles to draw

particle_container()

Access the beam particle container (impactx.ParticleContainer).

property lattice

Access the elements in the accelerator lattice. See impactx.elements for lattice elements.

property abort_on_warning_threshold

(optional) Set to “low”, “medium” or “high”. Cause the code to abort if a warning is raised that exceeds the warning threshold.

property abort_on_unused_inputs

Set to 1 to cause the simulation to fail after its completion if there were unused parameters. (default: 0 for false) It is mainly intended for continuous integration and automated testing to check that all tests and inputs are adapted to API changes.

property always_warn_immediately

If set to 1, ImpactX immediately prints every warning message as soon as it is generated. (default: 0 for false) It is mainly intended for debug purposes, in case a simulation crashes before a global warning report can be printed.

evolve()

Run the main simulation loop for a number of steps.

class impactx.Config

Configuration information on ImpactX that were set at compile-time.

property have_mpi

Indicates multi-process/multi-node support via the message-passing interface (MPI). Possible values: True/False

Note

Particle beam particles are not yet dynamically load balanced. Please see the progress in issue 198.

property have_gpu

Indicates GPU support. Possible values: True/False

property gpu_backend

Indicates the available GPU support. Possible values: None, "CUDA" (for Nvidia GPUs), "HIP" (for AMD GPUs) or "SYCL" (for Intel GPUs).

property have_omp

Indicates multi-threaded CPU support via OpenMP. Possible values: True/False`

Set the environment variable OMP_NUM_THREADS to control the number of threads.

Note

Not yet implemented. Please see the progress in issue 195.

Warning

By default, OpenMP spawns as many threads as there are available virtual cores on a host. When MPI and OpenMP support are used at the same time, it can easily happen that one over-subscribes the available physical CPU cores. This will lead to a severe slow-down of the simulation.

By setting appropriate environment variables for OpenMP, ensure that the number of MPI processes (ranks) per node multiplied with the number of OpenMP threads is equal to the number of physical (or virtual) CPU cores. Please see our examples in the high-performance computing (HPC) on how to run efficiently in parallel environments such as supercomputers.

Particles

class impactx.ParticleContainer

Beam Particles in ImpactX.

This class stores particles, distributed over MPI ranks.

add_n_particles(lev, x, y, z, px, py, pyz, qm, bchchg)

Add new particles to the container

Note: This can only be used after the initialization (grids) have

been created, meaning after the call to ImpactX.init_grids() has been made in the ImpactX class.

Parameters
  • lev – mesh-refinement level

  • x – positions in x

  • y – positions in y

  • z – positions in z

  • px – momentum in x

  • py – momentum in y

  • pz – momentum in z

  • qm – charge over mass in 1/eV

  • bchchg – total charge within a bunch in C

ref_particle()

Access the reference particle (impactx.RefPart).

Returns

return a data reference to the reference particle

Return type

impactx.RefPart

set_ref_particle(refpart)

Set reference particle attributes.

Parameters

refpart (impactx.RefPart) – a reference particle to copy all attributes from

class impactx.RefPart

This struct stores the reference particle attributes stored in impactx.ParticleContainer.

property s

integrated orbit path length, in meters

property x

horizontal position x, in meters

property y

vertical position y, in meters

property z

longitudinal position y, in meters

property t

clock time * c in meters

property px

momentum in x, normalized to proper velocity

property py

momentum in y, normalized to proper velocity

property pz

momentum in z, normalized to proper velocity

property pt

energy deviation, normalized by rest energy

property gamma

Read-only: Get reference particle relativistic gamma.

property beta

Read-only: Get reference particle relativistic beta.

property beta_gamma

Read-only: Get reference particle beta*gamma

property qm_qeeV

Read-only: Get reference particle charge to mass ratio (elementary charge/eV)

set_charge_qe(charge_qe)

Write-only: Set reference particle charge in (positive) elementary charges.

set_mass_MeV(massE)

Write-only: Set reference particle rest mass (MeV/c^2).

set_energy_MeV(energy_MeV)

Write-only: Set reference particle energy.

load_file(madx_file)

Load reference particle information from a MAD-X file.

Parameters

madx_file – file name to MAD-X file with a BEAM entry

Initial Beam Distributions

This module provides particle beam distributions that can be used to initialize particle beams in an impactx.ParticleContainer.

class impactx.distribution.Gaussian(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A 6D Gaussian distribution.

Parameters
  • sigx – for zero correlation, these are the related RMS sizes (in meters)

  • sigy – see sigx

  • sigt – see sigx

  • sigpx – RMS momentum

  • sigpy – see sigpx

  • sigpt – see sigpx

  • muxpx – correlation length-momentum

  • muypy – see muxpx

  • mutpt – see muxpx

class impactx.distribution.Kurth4D(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A 4D Kurth distribution transversely + a uniform distribution it t + a Gaussian distribution in pt.

class impactx.distribution.Kurth6D(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A 6D Kurth distribution.

R. Kurth, Quarterly of Applied Mathematics vol. 32, pp. 325-329 (1978) C. Mitchell, K. Hwang and R. D. Ryne, IPAC2021, WEPAB248 (2021)

class impactx.distribution.KVdist(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A K-V distribution transversely + a uniform distribution it t + a Gaussian distribution in pt.

class impactx.distribution.None

This distribution does nothing.

class impactx.distribution.Semigaussian(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A 6D Semi-Gaussian distribution (uniform in position, Gaussian in momentum).

class impactx.distribution.Waterbag(sigx, sigy, sigt, sigpx, sigpy, sigpt, muxpx=0.0, muypy=0.0, mutpt=0.0)

A 6D Waterbag distribution.

Lattice Elements

This module provides elements for the accelerator lattice.

class impactx.elements.KnownElementsList

An iterable, list-like type of elements.

clear()

Clear the list to become empty.

extend(list)

Add a list of elements to the list.

append(element)

Add a single element to the list.

load_file(madx_file, nslice=1)

Load and append an accelerator lattice description from a MAD-X file.

Parameters
  • madx_file – file name to MAD-X file with beamline elements

  • nslice – number of slices used for the application of space charge

class impactx.elements.ConstF(ds, kx, ky, kt, nslice=1)

A linear Constant Focusing element.

Parameters
  • ds – Segment length in m.

  • kx – Focusing strength for x in 1/m.

  • ky – Focusing strength for y in 1/m.

  • kt – Focusing strength for t in 1/m.

  • nslice – number of slices used for the application of space charge

class impactx.elements.DipEdge(psi, rc, g, K2)

Edge focusing associated with bend entry or exit

This model assumes a first-order effect of nonzero gap. Here we use the linear fringe field map, given to first order in g/rc (gap / radius of curvature).

References: * K. L. Brown, SLAC Report No. 75 (1982). * K. Hwang and S. Y. Lee, PRAB 18, 122401 (2015).

Parameters
  • psi – Pole face angle in rad

  • rc – Radius of curvature in m

  • g – Gap parameter in m

  • K2 – Fringe field integral (unitless)

class impactx.elements.Drift(ds, nslice=1)

A drift.

Parameters
  • ds – Segment length in m

  • nslice – number of slices used for the application of space charge

class impactx.elements.Multipole(multipole, K_normal, K_skew)

A general thin multipole element.

Parameters
  • multipole – index m (m=1 dipole, m=2 quadrupole, m=3 sextupole etc.)

  • K_normal – Integrated normal multipole coefficient (1/meter^m)

  • K_skew – Integrated skew multipole coefficient (1/meter^m)

class impactx.elements.NonlinearLens(knll, cnll)

Single short segment of the nonlinear magnetic insert element

A thin lens associated with a single short segment of the nonlinear magnetic insert described by V. Danilov and S. Nagaitsev, PRSTAB 13, 084002 (2010), Sect. V.A. This element appears in MAD-X as type NLLENS.

Parameters
  • knll – integrated strength of the nonlinear lens (m)

  • cnll – distance of singularities from the origin (m)

class impactx.elements.Quad(ds, k, nslice=1)

A Quadrupole magnet.

Parameters
  • ds – Segment length in m.

  • k – Quadrupole strength in m^(-2) (MADX convention) = (gradient in T/m) / (rigidity in T-m) k > 0 horizontal focusing k < 0 horizontal defocusing

  • nslice – number of slices used for the application of space charge

class impactx.elements.Sbend(ds, rc, nslice=1)

An ideal sector bend.

Parameters
  • ds – Segment length in m.

  • rc – Radius of curvature in m.

  • nslice – number of slices used for the application of space charge

class impactx.elements.ShortRF(V, k)

A short RF cavity element at zero crossing for bunching.

Parameters
  • V – Normalized RF voltage drop V = Emax*L/(c*Brho)

  • k – Wavenumber of RF in 1/m