DSMC: Add impact ionization (#5654)

This PR extends the work already done by @RemiLehe and @oshapoval in #5524. Some changes were made to how the ionization process is initialized, for example, the user must now specify the `target_species` (i.e., the species which undergoes ionization) as well as the `product_species`. The product species can include the colliding species (for example electron + neutral -> 2 x electron + ion), but does not have to (for example H$^+$ + D -> H$^+$ + D$^+$ + electron). The test created by @archermarx is now passing (at least early on): ![image](https://github.com/user-attachments/assets/344476b8-bc63-4395-92c1-d795183048b9) Todo: - [x] test implementation - [x] fix scattering process to conserve momentum - [x] clean up code --------- Signed-off-by: roelof-groenewald <[email protected]> Co-authored-by: Remi Lehe <[email protected]> Co-authored-by: Olga Shapoval <[email protected]> Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
ECP-WarpX · Mar 4, 2025 · 645a7a3 · 645a7a3
1 parent b858e36
commit 645a7a3
Show file tree

Hide file tree

Showing 18 changed files with 905 additions and 149 deletions.
diff --git a/Docs/source/usage/parameters.rst b/Docs/source/usage/parameters.rst
@@ -2070,7 +2070,7 @@ Details about the collision models can be found in the :ref:`theory section <mul
     In this case, only one species name should be given.
 
 * ``<collision_name>.product_species`` (`strings`)
-    Only for ``nuclearfusion``. The name(s) of the species in which to add
+    Only for ``dsmc`` and ``nuclearfusion``. The name(s) of the species in which to add
     the new macroparticles created by the reaction.
 
 * ``<collision_name>.ndt`` (`int`) optional
@@ -2185,14 +2185,18 @@ Details about the collision models can be found in the :ref:`theory section <mul
     represent the kinetic energy of the colliding particles in the center-of-mass frame.
 
 * ``<collision_name>.<scattering_process>_energy`` (`float`)
-    Only for ``background_mcc``. If the scattering process is either
+    Only for ``dsmc`` and ``background_mcc``. If the scattering process is either
     ``excitationX`` or ``ionization`` the energy cost of that process must be given in eV.
 
 * ``<collision_name>.ionization_species`` (`float`)
     Only for ``background_mcc``. If the scattering process is ``ionization`` the
     produced species must also be given. For example if argon properties is used
     for the background gas, a species of argon ions should be specified here.
 
+* ``<collision_name>.ionization_target_species`` (`string`)
+    Only for ``dsmc`` with impact ionization. This specifies which one of the
+    colliding particles is ionized.
+
 .. _running-cpp-parameters-numerics:
 
 Numerics and algorithms

diff --git a/Examples/Physics_applications/capacitive_discharge/analysis_dsmc.py b/Examples/Physics_applications/capacitive_discharge/analysis_dsmc.py
@@ -6,39 +6,39 @@
 
 # fmt: off
 ref_density = np.array([
-    1.27942709e+14, 2.23579371e+14, 2.55384387e+14, 2.55660663e+14,
-    2.55830911e+14, 2.55814337e+14, 2.55798906e+14, 2.55744891e+14,
-    2.55915585e+14, 2.56083194e+14, 2.55942354e+14, 2.55833026e+14,
-    2.56036175e+14, 2.56234141e+14, 2.56196179e+14, 2.56146141e+14,
-    2.56168022e+14, 2.56216909e+14, 2.56119961e+14, 2.56065167e+14,
-    2.56194764e+14, 2.56416398e+14, 2.56465239e+14, 2.56234337e+14,
-    2.56234503e+14, 2.56316003e+14, 2.56175023e+14, 2.56030269e+14,
-    2.56189301e+14, 2.56286379e+14, 2.56130396e+14, 2.56295225e+14,
-    2.56474082e+14, 2.56340375e+14, 2.56350864e+14, 2.56462330e+14,
-    2.56469391e+14, 2.56412726e+14, 2.56241788e+14, 2.56355650e+14,
-    2.56650599e+14, 2.56674748e+14, 2.56642480e+14, 2.56823508e+14,
-    2.57025029e+14, 2.57110614e+14, 2.57042364e+14, 2.56950884e+14,
-    2.57051822e+14, 2.56952148e+14, 2.56684016e+14, 2.56481130e+14,
-    2.56277073e+14, 2.56065774e+14, 2.56190033e+14, 2.56411074e+14,
-    2.56202418e+14, 2.56128368e+14, 2.56227002e+14, 2.56083004e+14,
-    2.56056768e+14, 2.56343831e+14, 2.56443659e+14, 2.56280541e+14,
-    2.56191572e+14, 2.56147304e+14, 2.56342794e+14, 2.56735473e+14,
-    2.56994680e+14, 2.56901500e+14, 2.56527131e+14, 2.56490824e+14,
-    2.56614730e+14, 2.56382744e+14, 2.56588214e+14, 2.57160270e+14,
-    2.57230435e+14, 2.57116530e+14, 2.57065771e+14, 2.57236507e+14,
-    2.57112865e+14, 2.56540177e+14, 2.56416828e+14, 2.56648954e+14,
-    2.56625594e+14, 2.56411003e+14, 2.56523754e+14, 2.56841108e+14,
-    2.56856368e+14, 2.56757912e+14, 2.56895134e+14, 2.57144419e+14,
-    2.57001944e+14, 2.56371759e+14, 2.56179404e+14, 2.56541905e+14,
-    2.56715727e+14, 2.56851681e+14, 2.57114458e+14, 2.57001739e+14,
-    2.56825690e+14, 2.56879682e+14, 2.56699673e+14, 2.56532841e+14,
-    2.56479582e+14, 2.56630989e+14, 2.56885996e+14, 2.56694637e+14,
-    2.56250819e+14, 2.56045278e+14, 2.56366075e+14, 2.56693733e+14,
-    2.56618530e+14, 2.56580918e+14, 2.56812781e+14, 2.56754216e+14,
-    2.56444736e+14, 2.56473391e+14, 2.56538398e+14, 2.56626551e+14,
-    2.56471950e+14, 2.56274969e+14, 2.56489423e+14, 2.56645266e+14,
-    2.56611124e+14, 2.56344324e+14, 2.56244156e+14, 2.24183727e+14,
-    1.27909856e+14
+    1.27939695e+14, 2.23589080e+14, 2.55400046e+14, 2.55652603e+14,
+    2.55810704e+14, 2.55816145e+14, 2.55810457e+14, 2.55743643e+14,
+    2.55908052e+14, 2.56076623e+14, 2.55948081e+14, 2.55841574e+14,
+    2.56029524e+14, 2.56320511e+14, 2.56608595e+14, 2.56755504e+14,
+    2.56699377e+14, 2.56700767e+14, 2.56497253e+14, 2.56481560e+14,
+    2.56832303e+14, 2.57064841e+14, 2.57023000e+14, 2.56614315e+14,
+    2.56368670e+14, 2.56370666e+14, 2.56227710e+14, 2.56240281e+14,
+    2.56673842e+14, 2.56837209e+14, 2.56625623e+14, 2.56729845e+14,
+    2.56975973e+14, 2.56801701e+14, 2.56491181e+14, 2.56516559e+14,
+    2.56468688e+14, 2.56251727e+14, 2.56243466e+14, 2.56484137e+14,
+    2.56637978e+14, 2.56448971e+14, 2.56140684e+14, 2.56117358e+14,
+    2.56274706e+14, 2.56233588e+14, 2.56047578e+14, 2.56087060e+14,
+    2.56365128e+14, 2.56357745e+14, 2.56269776e+14, 2.56419914e+14,
+    2.56392856e+14, 2.56202826e+14, 2.56363244e+14, 2.56572545e+14,
+    2.56351695e+14, 2.56393353e+14, 2.56759784e+14, 2.56767115e+14,
+    2.56700246e+14, 2.56618056e+14, 2.56234915e+14, 2.56237788e+14,
+    2.56606031e+14, 2.56520133e+14, 2.56316818e+14, 2.56184858e+14,
+    2.56246807e+14, 2.56626394e+14, 2.56747253e+14, 2.56630112e+14,
+    2.56518940e+14, 2.56358089e+14, 2.56249884e+14, 2.56271535e+14,
+    2.56420396e+14, 2.56704340e+14, 2.56912250e+14, 2.56823163e+14,
+    2.56694985e+14, 2.56822690e+14, 2.56736406e+14, 2.56438911e+14,
+    2.56359312e+14, 2.56356028e+14, 2.56415261e+14, 2.56408702e+14,
+    2.56267048e+14, 2.56274807e+14, 2.56494202e+14, 2.56789842e+14,
+    2.56939719e+14, 2.56875327e+14, 2.56831776e+14, 2.56827482e+14,
+    2.56698383e+14, 2.56712727e+14, 2.56879409e+14, 2.56629297e+14,
+    2.56322165e+14, 2.56377317e+14, 2.56277894e+14, 2.56112364e+14,
+    2.56171697e+14, 2.56370929e+14, 2.56855124e+14, 2.57621107e+14,
+    2.57656000e+14, 2.56760729e+14, 2.56449741e+14, 2.56716250e+14,
+    2.56721224e+14, 2.56506121e+14, 2.56236691e+14, 2.56270200e+14,
+    2.56745053e+14, 2.56940581e+14, 2.56539958e+14, 2.56403313e+14,
+    2.56600509e+14, 2.56776206e+14, 2.56884434e+14, 2.56755321e+14,
+    2.56558818e+14, 2.56400159e+14, 2.56223931e+14, 2.23879043e+14,
+    1.27601051e+14
 ])
 # fmt: on
 

diff --git a/Examples/Physics_applications/capacitive_discharge/inputs_base_1d_picmi.py b/Examples/Physics_applications/capacitive_discharge/inputs_base_1d_picmi.py
@@ -268,32 +268,55 @@ def setup_run(self):
         #######################################################################
 
         cross_sec_direc = "../../../../warpx-data/MCC_cross_sections/He/"
-        electron_colls = picmi.MCCCollisions(
-            name="coll_elec",
-            species=self.electrons,
-            background_density=self.gas_density,
-            background_temperature=self.gas_temp,
-            background_mass=self.ions.mass,
-            ndt=self.mcc_subcycling_steps,
-            scattering_processes={
-                "elastic": {
-                    "cross_section": cross_sec_direc + "electron_scattering.dat"
-                },
-                "excitation1": {
-                    "cross_section": cross_sec_direc + "excitation_1.dat",
-                    "energy": 19.82,
-                },
-                "excitation2": {
-                    "cross_section": cross_sec_direc + "excitation_2.dat",
-                    "energy": 20.61,
-                },
-                "ionization": {
-                    "cross_section": cross_sec_direc + "ionization.dat",
-                    "energy": 24.55,
-                    "species": self.ions,
-                },
+
+        electron_scattering_processes = {
+            "elastic": {"cross_section": cross_sec_direc + "electron_scattering.dat"},
+            "excitation1": {
+                "cross_section": cross_sec_direc + "excitation_1.dat",
+                "energy": 19.82,
             },
-        )
+            "excitation2": {
+                "cross_section": cross_sec_direc + "excitation_2.dat",
+                "energy": 20.61,
+            },
+            "ionization": {
+                "cross_section": cross_sec_direc + "ionization.dat",
+                "energy": 24.55,
+                "species": self.ions,
+            },
+        }
+        if self.dsmc:
+            ionization = {"ionization": electron_scattering_processes.pop("ionization")}
+            ionization["ionization"]["target_species"] = self.neutrals
+            ionization["ionization"].pop("species")
+            electron_colls_dsmc = picmi.DSMCCollisions(
+                name="coll_elec_dsmc",
+                species=[self.electrons, self.neutrals],
+                product_species=[self.ions, self.electrons],
+                ndt=4,
+                scattering_processes=ionization,
+            )
+            electron_colls_mcc = picmi.MCCCollisions(
+                name="coll_elec",
+                species=self.electrons,
+                background_density=self.gas_density,
+                background_temperature=self.gas_temp,
+                background_mass=self.ions.mass,
+                ndt=self.mcc_subcycling_steps,
+                scattering_processes=electron_scattering_processes,
+            )
+            electron_colls = [electron_colls_mcc, electron_colls_dsmc]
+        else:
+            electron_colls_mcc = picmi.MCCCollisions(
+                name="coll_elec",
+                species=self.electrons,
+                background_density=self.gas_density,
+                background_temperature=self.gas_temp,
+                background_mass=self.ions.mass,
+                ndt=self.mcc_subcycling_steps,
+                scattering_processes=electron_scattering_processes,
+            )
+            electron_colls = [electron_colls_mcc]
 
         ion_scattering_processes = {
             "elastic": {"cross_section": cross_sec_direc + "ion_scattering.dat"},
@@ -316,6 +339,7 @@ def setup_run(self):
                 ndt=self.mcc_subcycling_steps,
                 scattering_processes=ion_scattering_processes,
             )
+        ion_colls = [ion_colls]
 
         #######################################################################
         # Initialize simulation                                               #
@@ -325,7 +349,7 @@ def setup_run(self):
             solver=self.solver,
             time_step_size=self.dt,
             max_steps=self.max_steps,
-            warpx_collisions=[electron_colls, ion_colls],
+            warpx_collisions=electron_colls + ion_colls,
             verbose=self.test,
         )
         self.solver.sim = self.sim

diff --git a/Examples/Tests/CMakeLists.txt b/Examples/Tests/CMakeLists.txt
@@ -26,6 +26,7 @@ add_subdirectory(gaussian_beam)
 add_subdirectory(implicit)
 add_subdirectory(initial_distribution)
 add_subdirectory(initial_plasma_profile)
+add_subdirectory(ionization_dsmc)
 add_subdirectory(field_ionization)
 add_subdirectory(ion_stopping)
 add_subdirectory(langmuir)

diff --git a/Examples/Tests/ionization_dsmc/CMakeLists.txt b/Examples/Tests/ionization_dsmc/CMakeLists.txt
@@ -0,0 +1,12 @@
+# Add tests (alphabetical order) ##############################################
+#
+
+add_warpx_test(
+    test_3d_ionization_dsmc  # name
+    3  # dims
+    2  # nprocs
+    inputs_test_3d_ionization_dsmc  # inputs
+    "analysis_ionization_dsmc_3d.py"  # analysis
+    "analysis_default_regression.py --path diags/diag1000250"  # checksum
+    OFF  # dependency
+)
diff --git a/Examples/Tests/ionization_dsmc/analysis_default_regression.py b/Examples/Tests/ionization_dsmc/analysis_default_regression.py
@@ -0,0 +1 @@
+../../analysis_default_regression.py