[PICMI] Add ionization model configuration support

docs/source/usage/picmi/intro.rst

@@ -48,16 +48,26 @@ After you have installed the dependencies you must include the PIConGPU PICMI im
 .. note::
    Above, we used ``$PICSRC`` as a short hand for the path to picongpu's source code directory, provided from your shell environment if a pre-configured profile is used.
 
-After you have installed all PICMI dependencies, simply create a user script, see :ref:`here <example_PICMI_setup>`, and generate a picongpu setup, see :ref:`generating a PIConGPU setup with PICMI <generating_setups_with_PICMI>`.
+After you have installed all PICMI dependencies, simply create a user script, see the :ref:`warm plasma <example_PICMI_setup_warm_plasma>` and :ref:`laser wakefield <example_PICMI_setup_lwfa>` examples, and generate a picongpu setup, see :ref:`generating a PIConGPU setup with PICMI <generating_setups_with_PICMI>`.
 
 Example User Script for a warm plasma setup:
 --------------------------------------
-.. _example_PICMI_setup:
+.. _example_PICMI_setup_warm_plasma:
 
 .. literalinclude:: ../../../../share/picongpu/pypicongpu/examples/warm_plasma/main.py
    :language: python
 
-Creates a directory ``generated_input``, where you can run ``pic-build`` and subsequently ``tbg``.
+Creates a directory ``warm_plasma``, where you can run ``pic-build`` and subsequently ``tbg``.
+You can find this and more elaborate examples in `share/picongpu/pypicongpu/examples`.
+
+Example User Script for a laser wakefield setup:
+--------------------------------------
+.. _example_PICMI_setup_lwfa:
+
+.. literalinclude:: ../../../../share/picongpu/pypicongpu/examples/laser_wakefield/main.py
+   :language: python
+
+Creates a directory ``LWFA``, where you can run ``pic-build`` and subsequently ``tbg``.
 
 Generation of PIConGPU setups with PICMI
 ----------------------------------------
@@ -113,21 +123,49 @@ Parameters/Methods prefixed with ``picongpu_`` are PIConGPU-exclusive.
 
 - **Simulation**
 
-  - ``__init__(..., picongpu_template_dir)``:
-    Specify the template dir to use for code generation,
+  not supported methods:
+
+  - ``add_interaction(self, interaction)``:
+    The PIConGPU PICMI interface does not support the PICMI interaction specification, due to PICMI standard ambiguities.
+    Instead you must use the PIConGPU specific ``Interaction`` interface described below.
+
+  additional constructor/configuration options:
+
+  - ``picongpu_template_dir``:
+    Specify the template directory to use for code generation,
     please refer to :ref:`the documentation on the matter for details <picmi-custom-generation>`
-  - ``__init__(...,  picongpu_typical_ppc)`` typical ppc to be used for normalization in PIConGPU
+  - ``picongpu_typical_ppc``:
+    typical particle per cell(ppc) to be used for normalization in PIConGPU, if not set explicitly, PIConGPU will use the median ppc of all defined species
+  - ``picongpu_moving_window_move_point``:
+    portion of the simulation window a light ray reaches from the time of the start of the simulation until the simulation window begins to move.
+
+    .. warning::
+
+      If the moving window is active, one gpu row in y direction is reserved for initializing new spaces, thereby reducing the simulation window size accordingly
+
+  - ``picongpu_moving_window_stop_iteration``:
+    iteration at which to stop moving the simulation window
+  - ``picongpu_interaction``:
+    ``Interaction`` object specifying all interactions of the simulation, i.e. all ionization models and their configurations and so on.
+    This replaces the PICMI ``add_interaction`` method.
+
+  additional method arguments:
+
   - ``write_input_file(..., pypicongpu_simulation)``:
     use a :ref:`PyPIConGPU simulation<PyPIConGPU_Intro>` object instead of an PICMI- simulation object to generate a PIConGPU input.
+
+  additional methods:
+
   - ``get_as_pypicongpu()``:
     convert the PICMI simulation object to an equivalent :ref:`PyPIConGPU <PyPIConGPU_Intro>` simulation object.
   - ``picongpu_get_runner()``:
-    Retrieve a :ref:`PyPIConGPU Runner <pypicongpu-running>` for running a PIConGPU simulation from Python, **not recommended**
+    Retrieve a :ref:`PyPIConGPU Runner <pypicongpu-running>` for running a PIConGPU simulation from Python, **not recommended, see :ref:`PICMI setup generation <generating_setups_with_PICMI>`**.
   - ``picongpu_add_custom_user_input()``:
     pass custom user input to the code generation.
     This may be used in conjunction with custom templates to change the code generation.
     See :ref:`PICMI custom code generation<picmi-custom-generation>` for the documentation on using custom input.
 
+
 - **Grid**
 
   - ``picongpu_n_gpus``:
@@ -150,17 +188,28 @@ Parameters/Methods prefixed with ``picongpu_`` are PIConGPU-exclusive.
 
 - **Species**
 
-  - ``picongpu_ionization_electrons``:
-    Electron species to use for ionization.
-    Optional, will be guessed if possible.
-  - ``picongpu_fully_ionized``:
-    When defining an element (using ``particle_type``) it may or may not be ionizable
+  - ``picongpu_fixed_charge``:
+    When defining an ion species using ``particle_type`` it may or may not be ionizable
+
+    - to **enable** ionization add an ionization model to the Interaction object of the simulation and set the initial charge state using ``charge_state``.
+    - to **disable** ionization set ``picongpu_fixed_charge=True``, this will fix the charge of particles of this species for entire simulation.
 
-    - to **enable** ionization simulation set ``charge_state`` to an integer
-    - to **disable** ionization (ions are only core without electrons) set ``picongpu_fully_ionized=True``
+    ``picongpu_fixed_charge`` maybe combined with ``charge_state`` to control which charge state is to used for the ion species
 
     If neither is set a warning is printed prompting for either of the options above.
 
+Ionization:
+^^^^^^^^^^^
+The PIConGPU PICMI interface currently supports the configuration of ionization only through a picongpu specific PICMI extension, not the in the PICMI standard defined interface, due to the lack of standardization of ionization algorithm names in the PICMI standard.
+
+Use the **Interaction** interface
+
+
+- **Interaction**
+  picongpu specific configuration of PIC-algorithm extensions.
+
+    - ``__init__(ground_state_ionizaion_model_list= <list of ionization models>)``
+
 Output
 ^^^^^^
 Output is currently **not configurable** for picongpu using the PICMI interface.

lib/python/picongpu/picmi/__init__.py

@@ -1,7 +1,6 @@
 """
 PICMI for PIConGPU
 """
-
 from .simulation import Simulation
 from .grid import Cartesian3DGrid
 from .solver import ElectromagneticSolver
@@ -10,9 +9,10 @@
 from .layout import PseudoRandomLayout
 from . import constants
 
-from .distribution import FoilDistribution
-from .distribution import UniformDistribution
-from .distribution import GaussianDistribution
+from .distribution import FoilDistribution, UniformDistribution, GaussianDistribution
+from .interaction import Interaction
+from .interaction.ionization.fieldionization import ADK, ADKVariant, BSI, BSIExtension, Keldysh
+from .interaction.ionization.electroniccollisionalequilibrium import ThomasFermi
 
 import picmistandard
 
@@ -27,12 +27,20 @@
     "GaussianLaser",
     "Species",
     "PseudoRandomLayout",
+    "constants",
     "FoilDistribution",
     "UniformDistribution",
     "GaussianDistribution",
-    "constants",
+    "ADK",
+    "ADKVariant",
+    "BSI",
+    "BSIExtension",
+    "Keldysh",
+    "ThomasFermi",
+    "Interaction",
 ]
 
+
 codename = "picongpu"
 """
 name of this PICMI implementation

lib/python/picongpu/picmi/distribution/Distribution.py

@@ -9,6 +9,7 @@
 
 import typing
 import pydantic
+import typeguard
 
 """
 note on rms_velocity:
@@ -32,6 +33,7 @@
 """
 
 
+@typeguard.typechecked
 class Distribution(pydantic.BaseModel):
     rms_velocity: typing.Tuple[float, float, float] = (0, 0, 0)
     """thermal velocity spread [m/s]"""

lib/python/picongpu/picmi/distribution/FoilDistribution.py

@@ -42,8 +42,8 @@ def picongpu_get_rms_velocity_si(self) -> typing.Tuple[float, float, float]:
 
     def get_as_pypicongpu(self) -> species.operation.densityprofile.DensityProfile:
         util.unsupported("fill in", self.fill_in)
-        util.unsupported("lower bound", self.lower_bound, [None, None, None])
-        util.unsupported("upper bound", self.upper_bound, [None, None, None])
+        util.unsupported("lower bound", self.lower_bound, (None, None, None))
+        util.unsupported("upper bound", self.upper_bound, (None, None, None))
 
         foilProfile = species.operation.densityprofile.Foil()
         foilProfile.density_si = self.density

lib/python/picongpu/picmi/distribution/GaussianDistribution.py

@@ -61,8 +61,8 @@ def get_as_pypicongpu(self) -> species.operation.densityprofile.DensityProfile:
         util.unsupported("fill in not active", self.fill_in, True)
 
         # @todo support bounds, Brian Marre, 2024
-        util.unsupported("lower bound", self.lower_bound, [None, None, None])
-        util.unsupported("upper bound", self.upper_bound, [None, None, None])
+        util.unsupported("lower bound", self.lower_bound, (None, None, None))
+        util.unsupported("upper bound", self.upper_bound, (None, None, None))
 
         gaussian_profile = species.operation.densityprofile.Gaussian()
 

lib/python/picongpu/picmi/distribution/UniformDistribution.py

@@ -44,8 +44,8 @@ def picongpu_get_rms_velocity_si(self) -> typing.Tuple[float, float, float]:
 
     def get_as_pypicongpu(self) -> species.operation.densityprofile.DensityProfile:
         util.unsupported("fill in", self.fill_in)
-        util.unsupported("lower bound", self.lower_bound, [None, None, None])
-        util.unsupported("upper bound", self.upper_bound, [None, None, None])
+        util.unsupported("lower bound", self.lower_bound, (None, None, None))
+        util.unsupported("upper bound", self.upper_bound, (None, None, None))
 
         profile = species.operation.densityprofile.Uniform()
         profile.density_si = self.density

lib/python/picongpu/picmi/interaction/__init__.py

@@ -0,0 +1,4 @@
+from .interaction import Interaction
+from . import ionization
+
+__all__ = ["Interaction", "ionization"]

lib/python/picongpu/picmi/interaction/interaction.py

@@ -0,0 +1,155 @@
+"""
+This file is part of PIConGPU.
+Copyright 2024 PIConGPU contributors
+Authors: Brian Edward Marre
+License: GPLv3+
+"""
+
+from ... import pypicongpu
+
+from .ionization.groundstateionizationmodel import GroundStateIonizationModel, IonizationModel
+
+import picmistandard
+
+import typeguard
+import pydantic
+
+
+@typeguard.typechecked
+class Interaction(pydantic.BaseModel):
+    """
+    Common interface of Particle-In-Cell particle interaction extensions
+
+    e.g. collisions, ionization, nuclear reactions
+
+    This interface is only a semantic interface for typing interactions for storage in the simulation object.
+    It does not specify interface requirements for sub classes, since they differ too much.
+    """
+
+    ground_state_ionization_model_list: list[GroundStateIonizationModel]
+    """
+    list of all interaction models that change the charge state of ions
+
+    e.g. field ionization, collisional ionization, ...
+
+    """
+
+    # @todo add Collisions as elastic interaction model, Brian Marre, 2024
+
+    @staticmethod
+    def update_constant_list(
+        existing_list: list[pypicongpu.species.constant.Constant],
+        new_list: list[pypicongpu.species.constant.Constant],
+    ) -> None:
+        """check if dicts may be merged without overwriting previously set values"""
+
+        new_constant_list = []
+
+        for constant_new in new_list:
+            exists_already = False
+            for constant in existing_list:
+                if type(constant) == type(constant_new):
+                    # constant_new already exists in existing constants list
+                    exists_already = True
+
+                    if constant != constant_new:
+                        # same type of constant but conflicting values
+                        raise ValueError(f"Constants {constant} and {constant_new} conflict with each other.")
+
+            if not exists_already:
+                new_constant_list.append(constant_new)
+            # ignore already existing constants
+
+        # update constant_list
+        existing_list.extend(new_constant_list)
+
+    def get_interaction_constants(
+        self, picmi_species: picmistandard.PICMI_Species
+    ) -> tuple[
+        list[pypicongpu.species.constant.Constant],
+        dict[IonizationModel, pypicongpu.species.constant.ionizationmodel.IonizationModel],
+    ]:
+        """get list of all constants required by interactions for the given species"""
+
+        has_ionization = False
+        constant_list = []
+        ionization_model_conversion = {}
+        for model in self.ground_state_ionization_model_list:
+            if model.ion_species == picmi_species:
+                has_ionization = True
+                model_constants = model.get_constants()
+                Interaction.update_constant_list(constant_list, model_constants)
+                ionization_model_conversion[model] = model.get_as_pypicongpu()
+
+        if has_ionization:
+            # add GroundStateIonization constant for entire species
+            constant_list.append(
+                pypicongpu.species.constant.GroundStateIonization(
+                    ionization_model_list=ionization_model_conversion.values()
+                )
+            )
+
+        # add additional interaction sub groups needing constants here
+        return constant_list, ionization_model_conversion
+
+    def fill_in_ionization_electron_species(
+        self,
+        pypicongpu_by_picmi_species: dict[picmistandard.PICMI_Species, pypicongpu.species.Species],
+        ionization_model_conversion_by_type_and_species: dict[
+            picmistandard.PICMI_Species,
+            None | dict[IonizationModel, pypicongpu.species.constant.ionizationmodel.IonizationModel],
+        ],
+    ) -> None:
+        """
+        add ionization models to pypicongpu species
+
+        In PICMI ionization is defined as a list ionization models owned by an interaction object which in turn is a
+        member of the simulation, with each ionization model storing its PICMI ion and PICMI ionization electron
+        species.
+
+        In contrast in PyPIConGPU each ion PyPIConGPU species owns a list of ionization models, each storing its
+        PyPIConGPU ionization electron species.
+
+        This creates the problem that upon translation of the PICMI species to an PyPIConGPU species the PyPIConGPU
+        ionization electron species might not exist yet.
+
+        Therefore we leave the ionization electron unspecified upon species creation and fill it in from the PICMI
+        simulation ionization model list later.
+
+        (And because python uses pointers, this will be applied to the existing species objects passed in
+            pypicongpu_by_picmi_species)
+        """
+
+        # ground state ionization model
+        for species, ionization_model_conversion in ionization_model_conversion_by_type_and_species.items():
+            if ionization_model_conversion is not None:
+                for picmi_ionization_model, pypicongpu_ionization_model in ionization_model_conversion.items():
+                    try:
+                        pypicongpu_ionization_electron_species = pypicongpu_by_picmi_species[
+                            picmi_ionization_model.ionization_electron_species
+                        ]
+                    except KeyError:
+                        raise ValueError(
+                            f"Ionization electron species of {picmi_ionization_model} not known to simulation.\n"
+                            + f"Please add species {picmi_ionization_model.ionization_electron_species.name} to"
+                            + " the simulation."
+                        )
+                    pypicongpu_ionization_model.ionization_electron_species = pypicongpu_ionization_electron_species
+
+    def __has_ground_state_ionization(self, species) -> bool:
+        """does at least one ground state ionization model list species as ion species?"""
+
+        for ionization_model in self.ground_state_ionization_model_list:
+            if species == ionization_model.ion_species:
+                return True
+        return False
+
+    def has_ionization(self, species) -> bool:
+        """does at least one ionization model list species as ion species?"""
+        from ..species import Species
+
+        assert isinstance(species, Species)
+
+        # add additional groups of ionization models here
+        ionization_configured = self.__has_ground_state_ionization(species)
+        return ionization_configured

lib/python/picongpu/picmi/interaction/ionization/__init__.py

@@ -0,0 +1,11 @@
+from .ionizationmodel import IonizationModel
+from .groundstateionizationmodel import GroundStateIonizationModel
+from . import fieldionization
+from . import electroniccollisionalequilibrium
+
+__all__ = [
+    "IonizationModel",
+    "GroundStateIonizationModel",
+    "fieldionization",
+    "electroniccollisionalequilibrium",
+]

lib/python/picongpu/picmi/interaction/ionization/electroniccollisionalequilibrium/__init__.py

@@ -0,0 +1,3 @@
+from .thomasfermi import ThomasFermi
+
+__all__ = ["ThomasFermi"]