Intro

PICMI is a Python interface for configuring and running Particle-In-Cell simulations, defined via the Python reference implementation available from github and pypip.

PICMI allows a user to configure a simulation in a Python script, called user script, by building a simulation Python object from different standardised building blocks.

../../_images/PICMI_structure_DetailMiddle.png

Fig. 2 general overview of PICMI interface, see the reference implementation for details

From this simulation object a user may then generate input files for all PIC-simulation codes supporting PICMI, but some features may be code specific or not supported by all codes.

Usage Quick-start

To use the PICMI interface you need a working PIConGPU environment, see install instructions and the Setup part of PIConGPU in 5 Minutes on Hemera for instructions.

In addition you need to install the Python dependencies of the PIConGPU PICMI implementation to your Python environment.

To install the Python dependencies you may either run the command below

pip install -r $PICSRC/lib/python/picongpu/picmi/requirements.txt

or install all the requirements listed in

  • $PICSRC/lib/python/picongpu/picmi/requirements.txt

  • and $PICSRC/lib/python/picongpu/pypicongpu/requirements.txt

After you have installed the dependencies you must include the PIConGPU PICMI implementation in your PYHTONPATH environment variable, for example by

export PYTHONPATH=$PICSRC/lib/python:$PYTHONPATH

Note

If you are using one of our pre-configured profiles, this is done automatically when you source a profile

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 the warm plasma and laser wakefield examples, and generate a picongpu setup, see generating a PIConGPU setup with PICMI.

Example User Script for a warm plasma setup:

"""
This file is part of PIConGPU.
Copyright 2021-2024 PIConGPU contributors
Authors: Hannes Troepgen
License: GPLv3+
"""

from picongpu import picmi

OUTPUT_DIRECTORY_PATH = "warm_plasma"

boundary_conditions = ["periodic", "periodic", "periodic"]
grid = picmi.Cartesian3DGrid(
    # note: [x] * 3 == [x, x, x]
    number_of_cells=[192] * 3,
    lower_bound=[0, 0, 0],
    upper_bound=[0.0111152256] * 3,
    # delta {x, y, z} is computed implicitly
    # lower & upper boundary conditions must be equal
    lower_boundary_conditions=boundary_conditions,
    upper_boundary_conditions=boundary_conditions,
)
solver = picmi.ElectromagneticSolver(method="Yee", grid=grid)

profile = picmi.UniformDistribution(
    density=1e20,
    # most probable E_kin = 5 mc^2
    # approx. 9000 keV for electrons
    # must be equal for all three components
    rms_velocity=[4.18 * picmi.constants.c] * 3,
)
electron = picmi.Species(
    name="e",
    # openPMD particle type
    particle_type="electron",
    initial_distribution=profile,
)

sim = picmi.Simulation(
    time_step_size=9.65531e-14,
    max_steps=1024,
    solver=solver,
)

layout = picmi.PseudoRandomLayout(n_macroparticles_per_cell=25)
sim.add_species(electron, layout)

if __name__ == "__main__":
    sim.write_input_file(OUTPUT_DIRECTORY_PATH)

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:

"""
This file is part of PIConGPU.
Copyright 2024-2024 PIConGPU contributors
Authors: Masoud Afshari, Brian Edward Marre
License: GPLv3+
"""

from picongpu import picmi
from picongpu import pypicongpu
import numpy as np

"""
@file PICMI user script reproducing the PIConGPU LWFA example

This Python script is example PICMI user script reproducing the LaserWakefield example setup, based on 8.cfg.
"""

# generation modifiers
ENABLE_IONS = True
ENABLE_IONIZATION = True
ADD_CUSTOM_INPUT = True
OUTPUT_DIRECTORY_PATH = "LWFA"

numberCells = np.array([192, 2048, 192])
cellSize = np.array([0.1772e-6, 0.4430e-7, 0.1772e-6])  # unit: meter)

# Define the simulation grid based on grid.param
grid = picmi.Cartesian3DGrid(
    picongpu_n_gpus=[2, 4, 1],
    number_of_cells=numberCells.tolist(),
    lower_bound=[0, 0, 0],
    upper_bound=(numberCells * cellSize).tolist(),
    lower_boundary_conditions=["open", "open", "open"],
    upper_boundary_conditions=["open", "open", "open"],
)

gaussianProfile = picmi.distribution.GaussianDistribution(
    density=1.0e25,
    center_front=8.0e-5,
    sigma_front=8.0e-5,
    center_rear=10.0e-5,
    sigma_rear=8.0e-5,
    factor=-1.0,
    power=4.0,
    vacuum_cells_front=50,
)

solver = picmi.ElectromagneticSolver(
    grid=grid,
    method="Yee",
)

laser = picmi.GaussianLaser(
    wavelength=0.8e-6,
    waist=5.0e-6 / 1.17741,
    duration=5.0e-15,
    propagation_direction=[0.0, 1.0, 0.0],
    polarization_direction=[1.0, 0.0, 0.0],
    focal_position=[float(numberCells[0] * cellSize[0] / 2.0), 4.62e-5, float(numberCells[2] * cellSize[2] / 2.0)],
    centroid_position=[float(numberCells[0] * cellSize[0] / 2.0), 0.0, float(numberCells[2] * cellSize[2] / 2.0)],
    picongpu_polarization_type=pypicongpu.laser.GaussianLaser.PolarizationType.CIRCULAR,
    a0=8.0,
    picongpu_phase=0.0,
)

random_layout = picmi.PseudoRandomLayout(n_macroparticles_per_cell=2)

# Initialize particles  based on speciesInitialization.param
# simulation schema : https://github.com/BrianMarre/picongpu/blob/2ddcdab4c1aca70e1fc0ba02dbda8bd5e29d98eb/share/picongpu/pypicongpu/schema/simulation.Simulation.json

# for particle type see https://github.com/openPMD/openPMD-standard/blob/upcoming-2.0.0/EXT_SpeciesType.md
species_list = []
if not ENABLE_IONIZATION:
    interaction = None

    electron_placed = picmi.Species(particle_type="electron", name="electron", initial_distribution=gaussianProfile)
    species_list.append((electron_placed, random_layout))

    if ENABLE_IONS:
        hydrogen_fully_ionized = picmi.Species(
            particle_type="H", name="hydrogen", picongpu_fixed_charge=True, initial_distribution=gaussianProfile
        )
        species_list.append((hydrogen_fully_ionized, random_layout))
else:
    if not ENABLE_IONS:
        raise ValueError("Ions species required for ionization")

    hydrogen_with_ionization = picmi.Species(
        particle_type="H", name="hydrogen", charge_state=0, initial_distribution=gaussianProfile
    )
    species_list.append((hydrogen_with_ionization, random_layout))

    electron_not_placed = picmi.Species(particle_type="electron", name="electron", initial_distribution=None)
    species_list.append((electron_not_placed, None))

    adk_ionization_model = picmi.ADK(
        ADK_variant=picmi.ADKVariant.CircularPolarization,
        ion_species=hydrogen_with_ionization,
        ionization_electron_species=electron_not_placed,
        ionization_current=None,
    )

    bsi_effectiveZ_ionization_model = picmi.BSI(
        BSI_extensions=[picmi.BSIExtension.EffectiveZ],
        ion_species=hydrogen_with_ionization,
        ionization_electron_species=electron_not_placed,
        ionization_current=None,
    )

    interaction = picmi.Interaction(
        ground_state_ionization_model_list=[adk_ionization_model, bsi_effectiveZ_ionization_model]
    )

sim = picmi.Simulation(
    solver=solver,
    max_steps=4000,
    time_step_size=1.39e-16,
    picongpu_moving_window_move_point=0.9,
    picongpu_interaction=interaction,
)
for species, layout in species_list:
    sim.add_species(species, layout=layout)

sim.add_laser(laser, None)

# additional non standardized custom user input
# only active if custom templates are used

# for generating setup with custom input see standard implementation,
#  see https://picongpu.readthedocs.io/en/latest/usage/picmi/custom_template.html
if ADD_CUSTOM_INPUT:
    min_weight_input = pypicongpu.customuserinput.CustomUserInput()
    min_weight_input.addToCustomInput({"minimum_weight": 10.0}, "minimum_weight")
    sim.picongpu_add_custom_user_input(min_weight_input)

    output_configuration = pypicongpu.customuserinput.CustomUserInput()
    output_configuration.addToCustomInput(
        {
            "png_plugin_data_list": "['Ex', 'Ey', 'Ez', 'Bx', 'By', 'Bz', 'Jx', 'Jy', 'Jz']",
            "png_plugin_SCALE_IMAGE": 1.0,
            "png_plugin_SCALE_TO_CELLSIZE": True,
            "png_plugin_WHITE_BOX_PER_GPU": False,
            "png_plugin_EM_FIELD_SCALE_CHANNEL1": 7,
            "png_plugin_EM_FIELD_SCALE_CHANNEL2": -1,
            "png_plugin_EM_FIELD_SCALE_CHANNEL3": -1,
            "png_plugin_CUSTOM_NORMALIZATION_SI": "5.0e12 / constants.c, 5.0e12, 15.0",
            "png_plugin_PRE_PARTICLE_DENS_OPACITY": 0.25,
            "png_plugin_PRE_CHANNEL1_OPACITY": 1.0,
            "png_plugin_PRE_CHANNEL2_OPACITY": 1.0,
            "png_plugin_PRE_CHANNEL3_OPACITY": 1.0,
            "png_plugin_preParticleDensCol": "colorScales::grayInv",
            "png_plugin_preChannel1Col": "colorScales::green",
            "png_plugin_preChannel2Col": "colorScales::none",
            "png_plugin_preChannel3Col": "colorScales::none",
            "png_plugin_preChannel1": "field_E.x() * field_E.x();",
            "png_plugin_preChannel2": "field_E.y()",
            "png_plugin_preChannel3": "-1.0_X * field_E.y()",
            "png_plugin_period": 100,
            "png_plugin_axis": "yx",
            "png_plugin_slicePoint": 0.5,
            "png_plugin_species_name": "electron",
            "png_plugin_folder_name": "pngElectronsYX",
        },
        "png plugin configuration",
    )

    output_configuration.addToCustomInput(
        {
            "energy_histogram_species_name": "electron",
            "energy_histogram_period": 100,
            "energy_histogram_bin_count": 1024,
            "energy_histogram_min_energy": 0.0,
            "energy_histogram_maxEnergy": 1000.0,
            "energy_histogram_filter": "all",
        },
        "energy histogram plugin configuration",
    )

    output_configuration.addToCustomInput(
        {
            "phase_space_species_name": "electron",
            "phase_space_period": 100,
            "phase_space_space": "y",
            "phase_space_momentum": "py",
            "phase_space_min": -1.0,
            "phase_space_max": 1.0,
            "phase_space_filter": "all",
        },
        "phase space plugin configuration",
    )

    output_configuration.addToCustomInput(
        {"openPMD_period": 100, "openPMD_file": "simData", "openPMD_extension": "bp"}, "openPMD plugin configuration"
    )

    output_configuration.addToCustomInput(
        {"checkpoint_period": 100, "checkpoint_backend": "openPMD", "checkpoint_restart_backend": "openPMD"},
        "checkpoint configuration",
    )

    output_configuration.addToCustomInput(
        {"macro_particle_count_period": 100, "macro_particle_count_species_name": "electron"},
        "macro particle count plugin configuration",
    )
    sim.picongpu_add_custom_user_input(output_configuration)

if __name__ == "__main__":
    sim.write_input_file(OUTPUT_DIRECTORY_PATH)

Creates a directory LWFA, where you can run pic-build and subsequently tbg.

Generation of PIConGPU setups with PICMI

The recommended way to use the generated simulations is to

  1. create the simulation in the PICMI

  2. call simulation.write_input_file(DIR)

  3. use the normal PIConGPU toolchain (pic-build, tbg) on the generated PIConGPU setup

Note

Rationale: PICMI does not (yet) support enough parameters to meaningfully control the execution process.

Additionally, the following methods work (but are not recommended):

  • call Simulation.step(NUM)

    • directly builds and runs the simulation

    • NUM must be the maximum number of steps

    • has no diagnostic output (i.e. console hangs without output)

  • call Simulation.picongpu_run()

    • equivalent to Simulation.step() with the maximum number of steps

  • use the PyPIConGPU runner

PICMI Reference

The full PICMI standard interface reference is available upstream.

PIConGPU specifics

PIConGPU has it’s own implementation of the PICMI standard with picongpu specific extensions. Therefore all PICMI implementations for use with PIConGPU are located in the picongpu.pimci namespace instead of the usual picmistandard namespace.

In addition names of classes differ a little from the standard names, specifically we usually strip the PICMI_ prefix, and sometimes a class may provide additional options, see Extensions below.

Extensions

Parameters/Methods prefixed with picongpu_ are PIConGPU-exclusive.

Warning

We strive to quickly contribute these parameters to PICMI upstream, so this list is to be considered volatile.

  • Simulation

    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 the documentation on the matter for details

    • 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 PyPIConGPU simulation 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 PyPIConGPU simulation object.

    • picongpu_get_runner(): Retrieve a PyPIConGPU Runner 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 PICMI custom code generation for the documentation on using custom input.

  • Grid

    • picongpu_n_gpus: list of a 1 or 3 integers, greater than zero, describing GPU distribution in space 3-integer list: [N_gpu_x, N_gpu_y, N_gpu_z] 1-integer list: [1, N_gpu_y, 1] Default is None equal to [1, 1, 1]

  • Gaussian Laser

    • Laguerre Modes (picongpu_laguerre_modes and picongpu_laguerre_phases): Two lists of float, passed to PIConGPU laser definition to use laguerre modes for laser description

    • picongpu_polarization_type: configuration of polarization of the laser, either linear or circular, default is linear.

    • picongpu_phase: phase offset of the laser, default is 0

    • picongpu_huygens_surface_positions configuration of the position of the huygens surfaces used by PIConGPU for laser feed in, in cells

  • Species

    • 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.

      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.

Warning

This is subject to change.

If you are using the the default templates some output is automatically enabled, including PNGs. For this the period is chosen that the output is generated (approx.) 100 times over the entire simulation duration.

To configure output you must change the generated files by using use custom user input and custom templates.

Unsupported Features

The PIConGPU PICMI interface currently does not support the entire PICMI interface due to standard ambiguities and/or incomplete implementation. If you try to use an unsupported feature, you will be alerted by either a warning printed to stdout or an error thrown (including because a class does not exist).

In this case read the error message to fix this.

For reference you can see how the tests in $PICSRC/test/python/picongpu/quick/picmi use the interface.

Note

If a feature is not out of the box supported by the PIConGPU PICMI interface but supported by PIConGPU it may still be configured in PICMI through custom code generation.

In general everything that may be configured using the extensive PIConGPU .param files , may by configured via PICMI using custom input and custom templates.

Note

Please consider contributing all custom template features and custom user input back to PIConGPU to allow further improvements of the standard, making your life and everybody else’ easier.

PyPIConGPU

In addition to the PICMI interface PIConGPU has a second Python interface called PyPIConGPU. This is the native Python interface of PIConGPU to which PICMI inputs are converted to.

This interface offers additional configuration options above and beyond the PICMI interface and may be used instead of PICMI to configure a PIConGPU simulation.

To generate a PIConGPU setup from a PyPIConGPU simulation object use the following.

Note

The PyPIConGPU interface may be used stand-alone or together with the PICMI interface.

In the latter case configure the PICMI simulation first, generate a PyPIConGPU simulation from the PICMI simulation and then continue configuring the PyPIConGPU object.

Changes of the PICMI object after creation of the PyPIConGPU simulation will not be reflected in the PyPIConGPU simulation object. Both simulations representations are independent of each other, after generation.