Adding EPAC2004 space charge benchmarks

docs/source/usage/examples.rst

@@ -30,6 +30,7 @@ This section allows you to **download input files** that correspond to different
    examples/kicker/README.rst
    examples/thin_dipole/README.rst
    examples/aperture/README.rst
+   examples/epac2004_benchmarks/README.rst
    examples/pytorch_surrogate_model/README.rst
    examples/apochromatic/README.rst
 

docs/source/usage/parameters.rst

@@ -236,6 +236,20 @@ Initial Beam Distributions
         * ``beam.muypy`` (``float``, dimensionless, default: ``0``) correlation Y-Py
         * ``beam.mutpt`` (``float``, dimensionless, default: ``0``) correlation T-Pt
 
+    * ``thermal`` for a 6D stationary thermal or bithermal distribution.
+      This distribution type is described, for example in:
+      R. D. Ryne et al, "A Test Suite of Space-Charge Problems for Code Benchmarking", in Proc. EPAC2004, Lucerne, Switzerland.
+      C. E. Mitchell et al, "ImpactX Modeling of Benchmark Tests for Space Charge Validation", in Proc. HB2023, Geneva, Switzerland.
+      With additional parameters:
+
+        * ``beam.k`` (``float``, in inverse meters) external focusing strength
+        * ``beam.kT`` (``float``, dimensionless) temperature of core population
+        * ``beam.kT_halo`` (``float``, dimensionless, default ``kT``) temperature of halo population
+        * ``beam.normalize`` (``float``, dimensionless) normalizing constant for core population
+        * ``beam.normalize_halo`` (``float``, dimensionless) normalizing constant for halo population
+        * ``beam.halo`` (``float``, dimensionless) fraction of charge in halo
+
+
 .. _running-cpp-parameters-lattice:
 
 Lattice Elements

docs/source/usage/python.rst

@@ -427,6 +427,10 @@ This module provides particle beam distributions that can be used to initialize
 
    A 6D Waterbag distribution.
 
+.. py:class:: impactx.distribution.Thermal(k, kT, kT_halo, normalize, normalize_halo, halo)
+
+   A 6D stationary thermal or bithermal distribution.
+
 
 Lattice Elements
 ----------------

examples/CMakeLists.txt

@@ -548,7 +548,7 @@ add_impactx_test(positron_channel.py
     OFF  # no plot script yet
 )
 
-# Cyclotron ############################################################
+# Cyclotron ###################################################################
 #
 # w/o space charge
 add_impactx_test(cyclotron
@@ -566,7 +566,7 @@ add_impactx_test(cyclotron.py
     OFF  # no plot script yet
 )
 
-# Combined-function bend ########################################################
+# Combined-function bend ######################################################
 #
 # w/o space charge
 add_impactx_test(cfbend
@@ -584,7 +584,7 @@ add_impactx_test(cfbend.py
     OFF  # no plot script yet
 )
 
-# Ballistic compression ########################################################
+# Ballistic compression #######################################################
 #
 # w/o space charge
 add_impactx_test(compression
@@ -602,7 +602,7 @@ add_impactx_test(compression.py
     OFF  # no plot script yet
 )
 
-# Kicker test ##########################################################
+# Kicker test #################################################################
 #
 # w/o space charge
 add_impactx_test(kicker
@@ -640,7 +640,47 @@ add_impactx_test(hvkicker_madx.py
         OFF  # no plot script yet
 )
 
-# IOTA s-dependent nonlinear lens test ##########################################################
+# FODO + RF EPAC2004 ##########################################################
+#
+# with space charge
+add_impactx_test(fodo_rf_sc
+    examples/epac2004_benchmarks/input_fodo_rf_SC.in
+      ON  # ImpactX MPI-parallel
+      OFF  # ImpactX Python interface
+    examples/epac2004_benchmarks/analysis_fodo_rf_SC.py
+    OFF  # no plot script yet
+)
+add_impactx_test(fodo_rf_sc.py
+    examples/epac2004_benchmarks/run_fodo_rf_SC.py
+      OFF  # ImpactX MPI-parallel
+      ON   # ImpactX Python interface
+    examples/epac2004_benchmarks/analysis_fodo_rf_SC.py
+    OFF  # no plot script yet
+)
+
+# Thermal Beam EPAC2004 #######################################################
+#
+# with space charge
+add_impactx_test(thermal
+    examples/epac2004_benchmarks/input_thermal.in
+      ON  # ImpactX MPI-parallel
+      OFF  # ImpactX Python interface
+    examples/epac2004_benchmarks/analysis_thermal.py
+    OFF  # no plot script yet
+)
+
+# Bithermal Beam EPAC2004 #####################################################
+#
+# with space charge
+add_impactx_test(bithermal
+    examples/epac2004_benchmarks/input_bithermal.in
+      ON  # ImpactX MPI-parallel
+      OFF  # ImpactX Python interface
+    examples/epac2004_benchmarks/analysis_bithermal.py
+    examples/epac2004_benchmarks/plot_bithermal.py
+)
+
+# IOTA s-dependent nonlinear lens test ########################################
 #
 # w/o space charge
 add_impactx_test(IOTA_nll
@@ -658,7 +698,7 @@ add_impactx_test(IOTA_nll.py
     OFF  # no plot script yet
 )
 
-# IOTA nonlinear lattice test ##########################################################
+# IOTA nonlinear lattice test #################################################
 #
 # w/o space charge
 add_impactx_test(IOTA_lattice
@@ -676,7 +716,7 @@ add_impactx_test(IOTA_lattice.py
     OFF  # no plot script yet
 )
 
-# Thin dipole ########################################################
+# Thin dipole #################################################################
 #
 # w/o space charge
 add_impactx_test(thin_dipole
@@ -694,7 +734,7 @@ add_impactx_test(thin_dipole.py
     OFF  # no plot script yet
 )
 
-# Aperture collimation ############################################################
+# Aperture collimation ########################################################
 #
 # w/o space charge
 add_impactx_test(aperture
@@ -712,7 +752,7 @@ add_impactx_test(aperture.py
     OFF  # no plot script yet
 )
 
-# Apochromat example ########################################################
+# Apochromat example ##########################################################
 #
 # w/o space charge
 add_impactx_test(apochromat

examples/epac2004_benchmarks/README.rst

@@ -0,0 +1,206 @@
+.. _examples-fodo-rf-sc:
+
+Cold Beam in a FODO Channel with RF Cavities (and Space Charge)
+===============================================================
+
+This example is based on the subsection of the same name in:
+R. D. Ryne et al, "A Test Suite of Space-Charge Problems for Code Benchmarking", in Proc. EPAC2004, Lucerne, Switzerland.
+
+See additional documentation in:
+C. E. Mitchell et al, "ImpactX Modeling of Benchmark Tests for Space Charge Validation", in Proc. HB2023, Geneva, Switzerland.
+
+A cold (zero momentum spread), uniform density, 250 MeV, 143 pC proton bunch propagates in a FODO lattice with 700 MHz RF
+cavities added for longitudinal confinement.  The on-axis profile of the RF electric field is given by:
+
+.. math::
+
+   E(z)=\exp(-(4z)^4)\cos(\frac{5\pi}{2}\tanh(5z)).
+
+The beam is matched to the 3D focusing, with space charge, using the rms envelope equations.
+
+The particle distribution should remain unchanged, to within the level expected due to numerical particle noise.
+This is tested using the second moments of the distribution.
+
+In this test, the initial and final values of :math:`\sigma_x`, :math:`\sigma_y`, :math:`\sigma_t`, :math:`\epsilon_x`, :math:`\epsilon_y`, and :math:`\epsilon_t` must agree with nominal values.
+
+
+Run
+---
+
+This example can be run as a Python script (``python3 run_fodo_rf_SC.py``) or with an app with an input file (``impactx input_fodo_rf_SC.in``).
+Each can also be prefixed with an `MPI executor <https://www.mpi-forum.org>`__, such as ``mpiexec -n 4 ...`` or ``srun -n 4 ...``, depending on the system.
+
+.. tab-set::
+
+   .. tab-item:: Python Script
+
+       .. literalinclude:: run_fodo_rf_SC.py
+          :language: python3
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/run_fodo_rf_SC.py``.
+
+   .. tab-item:: App Input File
+
+       .. literalinclude:: input_fodo_rf_SC.in
+          :language: ini
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/input_fodo_rf_SC.in``.
+
+
+Analyze
+-------
+
+We run the following script to analyze correctness:
+
+.. dropdown:: Script ``analysis_fodo_rf_SC.py``
+
+   .. literalinclude:: analysis_fodo_rf_SC.py
+      :language: python3
+      :caption: You can copy this file from ``examples/epac2004_benchmarks/analysis_fodo_rf_SC.py``.
+
+
+
+.. _examples-thermal-beam:
+
+Thermal Beam in a Constant Focusing Channel (with Space Charge)
+===============================================================
+
+This example is based on the subsection of the same name in:
+R. D. Ryne et al, "A Test Suite of Space-Charge Problems for Code Benchmarking", in Proc. EPAC2004, Lucerne, Switzerland.
+
+See additional documentation in:
+C. E. Mitchell et al, "ImpactX Modeling of Benchmark Tests for Space Charge Validation", in Proc. HB2023, Geneva, Switzerland.
+
+This example illustrates a stationary solution of the Vlasov-Poisson equations with spherical symmetry (in the beam
+rest frame).  The distribution represents a thermal equilibrium of the form:
+
+.. math::
+
+   f=C\exp(-H/kT),
+
+where :math:`C` and :math:`kT` are constants, and :math:`H` denotes the self-consistent Hamiltonian with space charge.
+
+In this example, a 0.1 MeV, 143 pC proton bunch with :math:`kT=36\times 10^{-6}` propagates in a constant focusing lattice
+with 3D isotropic focusing.  (The isotropy is exact in the beam rest frame.)
+
+The particle distribution should remain unchanged, to within the level expected due to numerical particle noise.
+This is tested using the second moments of the distribution.
+
+In this test, the initial and final values of :math:`\sigma_x`, :math:`\sigma_y`, :math:`\sigma_t`, :math:`\epsilon_x`, :math:`\epsilon_y`, and :math:`\epsilon_t` must agree with nominal values.
+
+
+Run
+---
+
+This example can be run as a Python script (``python3 run_thermal.py``) or as an app with an input file (``impactx input_thermal.in``).
+Each can also be prefixed with an `MPI executor <https://www.mpi-forum.org>`__, such as ``mpiexec -n 4 ...`` or ``srun -n 4 ...``, depending on the system.
+
+.. tab-set::
+
+   .. tab-item:: Python Script
+
+       .. literalinclude:: run_thermal.py
+          :language: python3
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/run_thermal.py``.
+
+
+   .. tab-item:: App Input File
+
+       .. literalinclude:: input_thermal.in
+          :language: ini
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/input_thermal.in``.
+
+
+Analyze
+-------
+
+We run the following script to analyze correctness:
+
+.. dropdown:: Script ``analysis_thermal.py``
+
+   .. literalinclude:: analysis_thermal.py
+      :language: python3
+      :caption: You can copy this file from ``examples/epac2004_benchmarks/analysis_thermal.py``.
+
+
+
+.. _examples-bithermal-beam:
+
+Bithermal Beam in a Constant Focusing Channel (with Space Charge)
+=================================================================
+
+This example is based on the subsection of the same name in:
+R. D. Ryne et al, "A Test Suite of Space-Charge Problems for Code Benchmarking", in Proc. EPAC2004, Lucerne, Switzerland.
+
+See additional documentation in:
+C. E. Mitchell et al, "ImpactX Modeling of Benchmark Tests for Space Charge Validation", in Proc. HB2023, Geneva, Switzerland.
+
+This example illustrates a stationary solution of the Vlasov-Poisson equations with spherical symmetry (in the beam rest frame).
+It provides a self-consistent model of a 3D bunch with a nontrivial core-halo distribution.
+
+The distribution represents a bithermal stationary distribution of the form:
+
+.. math::
+
+   f=c_1\exp(-H/kT_1)+c_2\exp(-H/kT_2),
+
+where :math:`c_j`, :math:`kT_j` :math:`(j=1,2)` are constants, and :math:`H` denotes the self-consistent Hamiltonian with space charge.
+
+In this example, a 0.1 MeV, 143 pC proton bunch with :math:`kT_1=36\times 10^{-6}` and :math:`kT_1=900\times 10^{-6}` propagates in a constant focusing lattice
+with 3D isotropic focusing.
+(The isotropy is exact in the beam rest frame.)
+5% of the total charge lies in the second (halo) population.
+
+The particle distribution should remain unchanged, to within the level expected due to numerical particle noise.
+This is tested using the second moments of the distribution.
+
+In this test, the initial and final values of :math:`\sigma_x`, :math:`\sigma_y`, :math:`\sigma_t`, :math:`\epsilon_x`, :math:`\epsilon_y`, and :math:`\epsilon_t` must agree with nominal values.
+
+
+Run
+---
+
+This example can be run as a Python script (``python3 run_bithermal.py``) or as an app with an input file (``impactx input_bithermal.in``).
+Each can also be prefixed with an `MPI executor <https://www.mpi-forum.org>`__, such as ``mpiexec -n 4 ...`` or ``srun -n 4 ...``, depending on the system.
+
+.. tab-set::
+
+   .. tab-item:: Python Script
+
+       .. literalinclude:: run_bithermal.py
+          :language: python3
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/run_bithermal.py``.
+
+
+   .. tab-item:: App Input File
+
+       .. literalinclude:: input_bithermal.in
+          :language: ini
+          :caption: You can copy this file from ``examples/epac2004_benchmarks/input_bithermal.in``.
+
+
+Analyze
+-------
+
+We run the following script to analyze correctness:
+
+.. dropdown:: Script ``analysis_bithermal.py``
+
+   .. literalinclude:: analysis_bithermal.py
+      :language: python3
+      :caption: You can copy this file from ``examples/epac2004_benchmarks/analysis_bithermal.py``.
+
+
+Visualize
+---------
+
+You can run the following script to visualize the initial and final beam distribution:
+
+.. dropdown:: Script ``plot_bithermal.py``
+
+   .. literalinclude:: plot_bithermal.py
+      :language: python3
+      :caption: You can copy this file from ``examples/fodo/plot_bithermal.py``.
+
+.. figure:: https://user-images.githubusercontent.com/1353258/293794130-9aaee337-d810-4221-8f6d-1d7d2134c1b7.png
+   :alt: Initial and final beam distribution.
+
+   Initial and final beam distribution.