diff --git a/.gitignore b/.gitignore
index 125cd925d3..c8dcbd0fff 100644
--- a/.gitignore
+++ b/.gitignore
@@ -37,3 +37,10 @@ docs/build/
# Visual Studio configuration and output files
out
/.vs
+
+# clangd files
+.cache/clangd
+compile_commands.json
+
+# Zed configuration files
+.zed
diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml
index 8dd51e3d75..9da022926d 100644
--- a/.pre-commit-config.yaml
+++ b/.pre-commit-config.yaml
@@ -39,7 +39,7 @@ repos:
hooks:
- id: enforce-ascii
types_or: [file]
- exclude_types: [rst, tex, tsv, jupyter, markdown]
+ exclude_types: [rst, tex, tsv, jupyter, markdown, svg]
exclude: .zenodo.json
- repo: local
hooks:
diff --git a/docs/source/models/media/picongpu-main-loop.svg b/docs/source/models/media/picongpu-main-loop.svg
new file mode 100644
index 0000000000..0189625337
--- /dev/null
+++ b/docs/source/models/media/picongpu-main-loop.svg
@@ -0,0 +1,4 @@
+
+
+
+
diff --git a/docs/source/models/pic.rst b/docs/source/models/pic.rst
index 7792d6fdbe..da6c245c3a 100644
--- a/docs/source/models/pic.rst
+++ b/docs/source/models/pic.rst
@@ -3,14 +3,15 @@
The Particle-in-Cell Algorithm
==============================
-.. sectionauthor:: Axel Huebl, Klaus Steiniger
+.. sectionauthor:: Axel Huebl, Klaus Steiniger, Fabia Dietrich
Please also refer to the textbooks [BirdsallLangdon]_, [HockneyEastwood]_, our :ref:`latest paper on PIConGPU ` and the works in [Huebl2014]_ and [Huebl2019]_ .
-System of Equations
+Maxwell's Equations
-------------------
.. math::
+ :label: Maxwell
\nabla \cdot \mathbf{E} &= \frac{1}{\varepsilon_0}\sum_s \rho_s
@@ -25,65 +26,104 @@ for multiple particle species :math:`s`.
Except for normalization of constants, PIConGPU implements the governing equations in SI units.
-Relativistic Plasma Physics
+Vlasov--Maxwell equation
---------------------------
-The 3D3V particle-in-cell method is used to describe many-body systems such as a plasmas.
-It approximates the Vlasov--Maxwell--Equation
+The 3D3V particle-in-cell algorithm is used to describe many-body systems such as plasmas.
+It simulates its behaviour by calculating the motion of its electrons and ions in the present electromagnetic fields with the Vlasov--Maxwell equation
.. math::
:label: VlasovMaxwell
- \partial_t f_s(\mathbf{x},\mathbf{v},t) + \mathbf{v} \cdot \nabla_x f_s(\mathbf{x},\mathbf{v},t) + \frac{q_s}{m_s} \left[ \mathbf{E}(\mathbf{x},t) + \mathbf{v} \times \mathbf{B}(\mathbf{x},t) \right] \cdot \nabla_v f_s(\mathbf{x},\mathbf{v},t) = 0
+ \partial_t f_s(\mathbf{x},\mathbf{p},t) + \mathbf{v} \cdot \nabla_x f_s(\mathbf{x},\mathbf{p},t) + q_s \left[ \mathbf{E}(\mathbf{x},t) + \frac{\mathbf{v}}{c} \times \mathbf{B}(\mathbf{x},t) \right] \cdot \nabla_p f_s(\mathbf{x},\mathbf{p},t) = 0
-with :math:`f_s` as the distribution function of a particle species :math:`s`, :math:`\mathbf{x},\mathbf{v},t` as position, velocity and time and :math:`\frac{q_s}{m_s}` the charge to mass-ratio of a species.
-The momentum is related to the velocity by :math:`\mathbf{p} = \gamma m_s \mathbf{v}`.
-
-The equations of motion are given by the Lorentz force as
+with :math:`f_s` as the distribution function of a particle species :math:`s`, :math:`\mathbf{x},\mathbf{p},t` as position, momentum and time and :math:`{q_s},{m_s}` the charge and mass of a species.
+The velocity is related to the momentum by
.. math::
+ \mathbf{v} = \frac{c \mathbf{p}}{\sqrt{(m_s c)^2 + \abs{p}^2}}.
- \frac{\mathrm{d}}{\mathrm{d}t} \mathbf{V_s}(t) &= \frac{q_s}{m_s} \left[ \mathbf{E}(\mathbf{X_s}(t),t) + \mathbf{V_s}(t) \times \mathbf{B}(\mathbf{X_s}(t),t) \right]\\
- \frac{\mathrm{d}}{\mathrm{d}t} \mathbf{X_s}(t) &= \mathbf{V_s}(t) .
-
-.. attention::
- TODO: write proper relativistic form
+Electro-Magnetic PIC Method
+---------------------------
-:math:`\mathbf{X}_s = (\mathbf x_1, \mathbf x_2, ...)_s` and :math:`\mathbf{V}_s = (\mathbf v_1, \mathbf v_2, ...)_s` are vectors of *marker* positions and velocities, respectively, which describe the ensemble of particles belonging to species :math:`s`.
+The distribution of **particles** is described by the distribution function :math:`f_s(\mathbf{x},\mathbf{p},t)` and tracked in a continuous 6D phase space (or 5D, for 2D3V simulations) :math:`(\mathbf{r},\mathbf{p})`.
+As the number of particles involved in typical plasma simulations is usually extremely large, they are combined to so-called *macroparticles*, containing up to millions of particles of a species :math:`s`.
+Those macroparticles sample the distribution function :math:`f_s` by representing blobs of incompressible phase fluid with finite size and velocity moving in phase space.
+Since the Lorentz force only depends on the charge-to-mass ratio :math:`q_s / m_s`, a macroparticle follows the same trajectory as an individual particle of the same species would.
.. note::
Particles in a particle species can have different charge states in PIConGPU.
In the general case, :math:`\frac{q_s}{m_s}` is not required to be constant per particle species.
-Electro-Magnetic PIC Method
----------------------------
+The temporal evolution of :math:`f_s` is simulated by advancing the macroparticles over time with :eq:`VlasovMaxwell`.
+
+The charge of the macroparticles is assigned to neighboring mesh points with help of an *assignment function* :math:`W(\mathbf{x})`.
+These assignment functions are called *shapes* in PIConGPU.
+See :ref:`section Hierarchy of Charge Assignment Schemes ` for a list of available macroparticle shapes in PIConGPU.
+
+The **fields** :math:`\mathbf{E}(t)` and :math:`\mathbf{B}(t)` are discretized on a regular mesh in Eulerian frame of reference (= static, see [EulerLagrangeFrameOfReference]_).
+Their temporal evolution is computed with the :ref:`Finite-Difference Time-Domain (FDTD) method `, meaning that the partial space
+and time derivatives occurring in Maxwell’s equations are approximated by centered finite differences.
+In this approach, the derivatives of field values are computed at positions intermediate to those where the field values are explicitly known.
+This results in a staggered grid arrangement, where the electric and magnetic field are offset by half a grid cell in space and half a time step.
+
+The basic PIC cycle, representing a single simulation time step :math:`(n \to n + 1)`, comprises four consecutive steps, which are computed for all macroparticles involved in the simulation:
+
+#. *Field to particle:* The fields :math:`\mathbf{E}^n`, :math:`\mathbf{B}^n` are interpolated onto the particle position math:`\mathbf{x}^n` with help of the assignment function :math:`W`.
+#. *Particle push:* By integrating the equation of motion obtained from the Lorentz force acting on the macroparticle, its change in position :math:`\mathbf{x}` and velocity :math:`\mathbf{u} = \gamma \mathbf{v}` is obtained:
+
+ .. math::
+ \mathbf{u}^{n+3/2} = \mathbf{u}^{n+1/2} + \Delta t \frac{q}{m} (\mathbf{E}^{n} + \bar{\mathbf{v}}^{n} \times \mathbf{B}^{n})
+ .. math::
+ \mathbf{x}^{n+1} = \mathbf{x}^n + \Delta t ~\mathbf{u}^{n+1/2}
+ with :math:`\bar{\mathbf{v}}^{n}` as an approximated velocity at :math:`t = n`.
+ There exist various particle pushers, differing in accuracy and performance and offering different approximation schemes for :math:`\bar{\mathbf{v}}^{n}`. Please refer to [Ripperda2018]_ for further information.
+ The particle pusher can be chosen and configured in :ref:`pusher.param `.
+#. *Current deposition:* The current :math:`\mathbf{J}` caused by the macroparticle movement is computed on the grid by solving the continuity equation of electrodynamics:
+
+ .. math::
+ \rho^{n+1} = \rho^n + \Delta t \cdot \nabla \mathbf{J}^{n + 1/2}
+ The charge densities :math:`\rho^{n}, \rho^{n+1}` are obtained from the macroparticle position before and after movement with help of the assignment function :math:`W`.
+ In PIConGPU, the user can choose between Esirkepov's current deposition method [Esirkepov2001]_ and the performance-increased EZ method [Steiniger2023]_.
+ See :ref:`here ` for further details.
+#. *Field update*: By inserting the current :math:`\mathbf{J}`, the electromagnetic fields can be updated by applying the third and fourth Maxwell equation :eq:`Maxwell`
+ using the :ref:`FDTD method `.
+ Thereby, the magnetic field update is split in two parts of each half a timestep, so that electric and magnetic fields are known at the same time.
+ This is necessary to caclulate the correct Lorentz force in the next run of the PIC cycle.
-**Fields** such as :math:`\mathbf{E}(t), \mathbf{B}(t)` and :math:`\mathbf{J}(t)` are discretized on a regular mesh in Eulerian frame of reference (see [EulerLagrangeFrameOfReference]_).
-See :ref:`section Finite-Difference Time-Domain Method ` describing how Maxwell's equations are discretized on a mesh in PIConGPU.
+By running the cycle multiple times, a longer simulation duration is achieved.
-The distribution of **Particles** is described by the distribution function :math:`f_s(\mathbf{x},\mathbf{v},t)`.
-This distribution function is sampled by *markers*, which are commonly referred to as *macroparticles*.
-These markers represent blobs of incompressible phase fluid moving in phase space.
-The temporal evolution of the distribution function is simulated by advancing the markers over time according to the Vlasov--Maxwell--Equation in Lagrangian frame (see eq. :eq:`VlasovMaxwell` and [EulerLagrangeFrameOfReference]_).
-A marker has a finite-size and a velocity, such that it can be regarded as a cloud of particles, whose center of mass is the marker's position and whose mean velocity is the marker's velocity.
-The cloud shape :math:`S^n(x)` of order :math:`n` of a marker describes its charge density distribution.
-See :ref:`section Hierarchy of Charge Assignment Schemes ` for a list of available marker shapes in PIConGPU.
+The field update requires only information from neighbouring or next-to neighbouring cells, depending on the FDTD scheme involved.
+This locality allows a parallelized computation of the algorithm by sub-dividing the simulation volume and mapping the underlying mesh onto local memory subsets, which can then be assigned to a single processing unit.
+Data transfer between processing units is required at every time step to transfer field data between neighbouring cells and when particles cross cell boundaries.
+
+The following flowchart provides a summary and visualization of the implementation of the PIC cycle within PIConGPU.
+
+.. image:: media/picongpu-main-loop.svg
+ :width: 700
+ :alt: PIC cycle in PIConGPU
References
----------
-.. [EulerLagrangeFrameOfReference]
- *Eulerian and Lagrangian specification of the flow field.*
- https://en.wikipedia.org/wiki/Lagrangian_and_Eulerian_specification_of_the_flow_field
-
.. [BirdsallLangdon]
C.K. Birdsall, A.B. Langdon.
*Plasma Physics via Computer Simulation*,
McGraw-Hill (1985),
ISBN 0-07-005371-5
+.. [EulerLagrangeFrameOfReference]
+ *Eulerian and Lagrangian specification of the flow field.*
+ https://en.wikipedia.org/wiki/Lagrangian_and_Eulerian_specification_of_the_flow_field
+
+.. [Esirkepov2001]
+ T.Zh. Esirkepov,
+ *Exact charge conservation scheme for particle-in-cell simulation with an arbitrary form-factor*,
+ Computer Physics Communications 135.2 (2001): 144-153,
+ `DOI:10.1016/S0010-4655(00)00228-9 `_
+
.. [HockneyEastwood]
R.W. Hockney, J.W. Eastwood.
*Computer Simulation Using Particles*,
@@ -101,3 +141,15 @@ References
*PIConGPU: Predictive Simulations of Laser-Particle Accelerators with Manycore Hardware*,
PhD Thesis at TU Dresden & Helmholtz-Zentrum Dresden - Rossendorf (2019),
`DOI:10.5281/zenodo.3266820 `_
+
+.. [Ripperda2018]
+ B. Ripperda et. al.
+ *A Comprehensive Comparison of Relativistic Particle Integrators*,
+ The American Astronomical Society (2018),
+ `DOI:10.3847/1538-4365/aab114 `_
+
+.. [Steiniger2023]
+ K. Steiniger et. al.
+ *EZ: An efficient, charge conserving current deposition algorithm for electromagnetic Particle-in-Cell simulations*,
+ Computer Physics Communications (2023)
+ `DOI:10.1016/j.cpc.2023.108849 `_