Implementation of Bremsstrahlung collisions

Docs/source/usage/parameters.rst

@@ -2061,17 +2061,21 @@ Details about the collision models can be found in the :ref:`theory section <mul
     - ``background_stopping`` for slowing of ions due to collisions with electrons or ions.
       This implements the approximate formulae as derived in Introduction to Plasma Physics,
       from Goldston and Rutherford, section 14.2.
+    - ``bremsstrahlung`` for slowing of electrons due to Bremsstrahlung collisions with ions.
+      This uses the cross sections as given by `Seltzer and Berger <https://doi.org/10.1016/0092-640X(86)90014-8>`__.
 
 * ``<collision_name>.species`` (`strings`)
-    If using ``dsmc``, ``pairwisecoulomb`` or ``nuclearfusion``, this should be the name(s) of the species,
+    If using ``dsmc``, ``pairwisecoulomb``, ``nuclearfusion``, or ``bremsstrahlung``, this should be the name(s) of the species,
     between which the collision will be considered. (Provide only one name for intra-species collisions.)
+    With ``bremsstrahlung``, the electron species must be given first, followed by the target species.
     If using ``background_mcc`` or ``background_stopping`` type this should be the name of the
     species for which collisions with a background will be included.
     In this case, only one species name should be given.
 
 * ``<collision_name>.product_species`` (`strings`)
-    Only for ``dsmc`` and ``nuclearfusion``. The name(s) of the species in which to add
+    Only for ``dsmc``, ``nuclearfusion``, and ``bremsstrahlung``. The name(s) of the species in which to add
     the new macroparticles created by the reaction.
+    For ``bremsstrahlung``, the product species must be of type photon.
 
 * ``<collision_name>.ndt`` (`int`) optional
     Execute collision every # time steps. The default value is 1.
@@ -2197,6 +2201,22 @@ Details about the collision models can be found in the :ref:`theory section <mul
     Only for ``dsmc`` with impact ionization. This specifies which one of the
     colliding particles is ionized.
 
+* ``<collision_name>.Z`` (`integer`)
+    Only for ``bremsstrahlung``. The atomic number of the target ion species.
+    Currently, only the values 1, 2, 5, 6 are supported.
+
+* ``<collision_name>.multiplier`` (`float`)
+    Only for ``bremsstrahlung``. Multiplier for the collision probability.
+    Any resulting photons will have the electron weight divided the multiplier.
+    The default is 1. This must be greater than or equal to 1.
+
+* ``<collision_name>.create_photons`` (`integer`)
+    Only for ``bremsstrahlung``. Whether photons will be created, defaults to 1 (true).
+
+* ``<collision_name>.koT1_cut`` (`float`)
+    Only for ``bremsstrahlung``. Minimum energy of the photons created.
+    This is relative to the electron energy, defaulting to 1.e-4.
+
 .. _running-cpp-parameters-numerics:
 
 Numerics and algorithms

Examples/Tests/collision/CMakeLists.txt

@@ -60,3 +60,13 @@ add_warpx_test(
     "analysis_default_regression.py --path diags/diag1000150 --skip-particles"  # checksum
     OFF  # dependency
 )
+
+add_warpx_test(
+    test_1d_collision_z_Bremsstrahlung  # name
+    1  # dims
+    1  # nprocs
+    inputs_test_1d_collision_z_Bremsstrahlung  # inputs
+    "analysis_collision_1d_Bremsstrahlung.py diags/diag1000100"  # analysis
+    "analysis_default_regression.py --path diags/diag1000100"  # checksum
+    OFF  # dependency
+)

Examples/Tests/collision/analysis_collision_1d_Bremsstrahlung.py

@@ -0,0 +1,109 @@
+#!/usr/bin/env python3
+
+# Copyright 2025 David Grote
+#
+#
+# This file is part of WarpX.
+#
+# License: BSD-3-Clause-LBNL
+#
+# This is a script that analyses the simulation results from the script `inputs_test_1d_collision_z_Bremsstrahlung`.
+# run locally: python analysis_collision_1d_Bremsstrahlung.py diags/diag1000600/
+#
+# This is a 1D Bremsstrahlung collision of electrons on Borons ions
+# producing photons. This checks that the appropriate number of
+# photons are created with the correct distribution.
+#
+import os
+import sys
+
+import numpy as np
+from scipy import constants
+
+# this will be the name of the plot file
+last_fn = sys.argv[1]
+
+particle_energy = np.loadtxt(
+    os.path.join("diags", "reducedfiles", "particle_energy.txt"), skiprows=1
+)
+particle_momentum = np.loadtxt(
+    os.path.join("diags", "reducedfiles", "particle_momentum.txt"), skiprows=1
+)
+particle_number = np.loadtxt(
+    os.path.join("diags", "reducedfiles", "particle_number.txt"), skiprows=1
+)
+
+total_energy = particle_energy[:, 2]
+electron_energy = particle_energy[:, 3]
+ion_energy = particle_energy[:, 4]
+photon_energy = particle_energy[:, 5]
+
+total_momentum = particle_momentum[:, 2]
+electron_momentum = particle_momentum[:, 3]
+ion_momentum = particle_momentum[:, 4]
+photon_momentum = particle_momentum[:, 5]
+
+energy_tolerance = 1.0e-11
+
+print(f"initial total energy = {total_energy[0]}")
+print(f"final total energy = {total_energy[-1]}")
+
+dE_total = np.abs(total_energy[-1] - total_energy[0]) / total_energy[0]
+print(f"change in total energy = {dE_total}")
+assert dE_total < energy_tolerance
+
+momentum_tolerance = 1.0e-12
+
+print(f"initial total momentum = {total_momentum[0]}")
+print(f"final total momentum = {total_momentum[-1]}")
+
+dP_total = np.abs(total_momentum[-1] - total_momentum[0]) / total_momentum[0]
+print(f"change in total momentum = {dP_total}")
+assert dP_total < momentum_tolerance
+
+print()
+
+dt = 1.0e-2 * 1.0e-15
+Z = 5
+n_i = 5.47e31
+n_e = 5.47e30
+L = 1.0e-6  # 1 micron
+T1 = 1.0e6  # 1 MeV
+N_e = n_e * L  # number of electrons
+m_e_eV = constants.m_e * constants.c**2 / constants.e
+gamma = T1 / m_e_eV + 1.0
+gamma_beta = np.sqrt(gamma**2 - 1.0)
+beta = gamma_beta / gamma
+
+phirad = 6.761  # from Seltzer and Berger for 1 MeV electron and Boron
+
+
+dEdx_simulation = (
+    (particle_energy[0, 3] - particle_energy[1:, 3])
+    / particle_energy[1:, 0]
+    / (beta * constants.c * dt)
+    / constants.e
+    / N_e
+)
+
+Boron_weight = 20065.0 * constants.m_e
+r_e = (
+    1.0
+    / (4.0 * constants.pi * constants.epsilon_0)
+    * (constants.e**2 / (constants.m_e * constants.c**2))
+)
+dEdx = n_i * constants.alpha * r_e**2 * Z**2 * (T1 + m_e_eV) * phirad
+
+print(f"dE/dx analytic = {dEdx}")
+print(f"dE/dx simulated = {dEdx_simulation[-1]}")
+print(f"dE/dx simulated/analytic = {dEdx_simulation[-1] / dEdx}")
+assert np.abs(dEdx_simulation[-1] / dEdx - 1.0) < 0.03
+
+sigma_total = 1.818e-28  # Calculated from table with k cutoff=1.e-4
+N_photon = n_e * n_i * L * beta * constants.c * sigma_total * dt
+print(
+    f"New photons per step simulated/analytic = {particle_number[-1, -1] / (particle_energy[-1, 0] * N_photon)}"
+)
+assert (
+    np.abs(particle_number[-1, -1] / (particle_energy[-1, 0] * N_photon) - 1.0) < 0.02
+)

Examples/Tests/collision/inputs_test_1d_collision_z_Bremsstrahlung

@@ -0,0 +1,106 @@
+#################################
+########## CONSTANTS ############
+#################################
+
+my_constants.n_e = 5.47e30     # electron density, m^-3
+my_constants.Np_e = 1024       # electrons per cell
+my_constants.T1 = 1.e6         # 1 MeV, electron energy
+my_constants.m_e_eV = m_e*clight**2/q_e # electron rest mass in eV
+my_constants.gamma = T1/m_e_eV + 1. # electron gamma factor
+my_constants.gamma_beta = sqrt(gamma**2 - 1.)
+
+my_constants.n_i = 5.47e31     # ion density, m^-3, for rho = 1000 g/cm^3
+my_constants.Np_i = 1024       # ions per cell
+my_constants.T_i = 2.0         # ion temperature, eV
+
+my_constants.m_B11 = 20065.0*m_e  # 11.01*amu/me - 5, Boron
+my_constants.q_B11 = 5.*q_e
+
+#################################
+####### GENERAL PARAMETERS ######
+#################################
+max_step = 100
+amr.n_cell = 100
+amr.max_level = 0
+amr.blocking_factor = 4
+geometry.dims = 1
+geometry.prob_lo = 0.
+geometry.prob_hi = 1.e-6
+
+#################################
+###### Boundary Condition #######
+#################################
+boundary.field_lo = periodic
+boundary.field_hi = periodic
+
+#################################
+############ NUMERICS ###########
+#################################
+warpx.serialize_initial_conditions = 1
+warpx.verbose = 1
+warpx.const_dt = 1.e-2*1.e-15
+warpx.use_filter = 0
+
+# Do not evolve the E and B fields
+algo.maxwell_solver = none
+
+# Order of particle shape factors
+algo.particle_shape = 1
+
+#################################
+############ PLASMA #############
+#################################
+particles.species_names = electrons ions photons
+particles.photon_species = photons
+
+electrons.species_type = electron
+electrons.do_not_deposit = 1
+
+electrons.injection_style = "NUniformPerCell"
+electrons.num_particles_per_cell_each_dim = Np_e
+electrons.profile = "constant"
+electrons.density = n_e
+electrons.momentum_distribution_type = "constant"
+electrons.uz = gamma_beta
+
+ions.charge = q_B11
+ions.mass = m_B11
+ions.do_not_deposit = 1
+
+ions.injection_style = "NUniformPerCell"
+ions.num_particles_per_cell_each_dim = Np_i
+ions.profile = "constant"
+ions.density = n_i
+ions.momentum_distribution_type = "gaussian"
+ions.ux_th = sqrt(T_i*q_e/m_B11)/clight
+ions.uy_th = sqrt(T_i*q_e/m_B11)/clight
+ions.uz_th = sqrt(T_i*q_e/m_B11)/clight
+
+photons.species_type = photon
+photons.injection_style = none
+
+#################################
+############ COLLISION ##########
+#################################
+collisions.collision_names = bremsstrahlung
+bremsstrahlung.type = bremsstrahlung
+bremsstrahlung.species = electrons ions
+bremsstrahlung.Z = 5
+bremsstrahlung.product_species = photons
+bremsstrahlung.multiplier = 100.
+
+# Diagnostics
+diagnostics.diags_names = diag1
+diag1.intervals = 100
+diag1.diag_type = Full
+
+warpx.reduced_diags_names = particle_energy particle_momentum particle_number
+particle_energy.type = ParticleEnergy
+particle_energy.intervals = 10
+particle_energy.precision = 18
+particle_momentum.type = ParticleMomentum
+particle_momentum.intervals = 10
+particle_momentum.precision = 18
+particle_number.type = ParticleNumber
+particle_number.intervals = 10
+particle_number.precision = 18

Regression/Checksum/benchmarks_json/test_1d_collision_z_Bremsstrahlung.json

@@ -0,0 +1,34 @@
+{
+  "lev=0": {
+    "Bx": 0.0,
+    "By": 0.0,
+    "Bz": 0.0,
+    "Ex": 0.0,
+    "Ey": 0.0,
+    "Ez": 0.0,
+    "jx": 0.0,
+    "jy": 0.0,
+    "jz": 0.0
+  },
+  "ions": {
+    "particle_momentum_x": 6.244058618052987e-18,
+    "particle_momentum_y": 6.266839026611908e-18,
+    "particle_momentum_z": 6.237937377390938e-18,
+    "particle_position_x": 0.05119999887292794,
+    "particle_weight": 5.470000000000001e+25
+  },
+  "electrons": {
+    "particle_momentum_x": 1.7754689380349235e-26,
+    "particle_momentum_y": 1.8019938095335607e-26,
+    "particle_momentum_z": 7.780447460289813e-17,
+    "particle_position_x": 0.05119966127601074,
+    "particle_weight": 5.470000000000001e+24
+  },
+  "photons": {
+    "particle_momentum_x": 8.902000149233514e-25,
+    "particle_momentum_y": 9.022633314982659e-25,
+    "particle_momentum_z": 1.271533261332079e-18,
+    "particle_position_x": 0.014344002048764832,
+    "particle_weight": 1.5380635742187505e+22
+  }
+}