Skip to content

Commit

Permalink
create new page for handling of kinetic particles, move push details …
Browse files Browse the repository at this point in the history 
…to that page
  • Loading branch information
@roelof-groenewald
roelof-groenewald committed Sep 13, 2024
1 parent 2e1b867 commit 367f97e
Show file tree
Hide file tree
Showing 3 changed files with 129 additions and 114 deletions.
17 changes: 16 additions & 1 deletion Docs/source/theory/intro.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,20 @@
Models & Algorithms
===================

In the *particle-in-cell method* :cite:p:`i-Birdsalllangdon`,
the electromagnetic fields are solved on a discretized grid while particles are
evolved in continuous space. A high-level schematic of the method is shown in the
figure below with details of the different field solvers, particle handling algorithms,
additional physics modules, etc. described in the sections linked further down.

.. _fig-pic:

.. figure:: PIC.png
:alt: [fig:PIC] The Particle-In-Cell (PIC) method follows the evolution of a collection of charged macro-particles (positively charged in blue on the left plot, negatively charged in red) that evolve self-consistently with their electromagnetic (or electrostatic) fields. The core PIC algorithm involves four operations at each time step: 1) evolve the velocity and position of the particles using the Newton-Lorentz equations, 2) deposit the charge and/or current densities through interpolation from the particles distributions onto the grid, 3) evolve Maxwell’s wave equations (for electromagnetic) or solve Poisson’s equation (for electrostatic) on the grid, 4) interpolate the fields from the grid onto the particles for the next particle push. Additional “add-ons” operations are inserted between these core operations to account for additional physics (e.g. absorption/emission of particles, addition of external forces to account for accelerator focusing or accelerating component) or numerical effects (e.g. smoothing/filtering of the charge/current densities and/or fields on the grid).

The Particle-In-Cell (PIC) method follows the evolution of a collection of charged macro-particles (positively charged in blue on the left plot, negatively charged in red) that evolve self-consistently with their electromagnetic (or electrostatic) fields. The core PIC algorithm involves four operations at each time step: 1) evolve the velocity and position of the particles using the Newton-Lorentz equations, 2) deposit the charge and/or current densities through interpolation from the particles distributions onto the grid, 3) evolve Maxwell’s wave equations (for electromagnetic) or solve Poisson’s equation (for electrostatic) on the grid, 4) interpolate the fields from the grid onto the particles for the next particle push. Additional “add-ons” operations are inserted between these core operations to account for additional physics (e.g. absorption/emission of particles, addition of external forces to account for accelerator focusing or accelerating component) or numerical effects (e.g. smoothing/filtering of the charge/current densities and/or fields on the grid).


.. figure:: Plasma_acceleration_sim.png
:alt: Plasma acceleration diagram

Expand All @@ -25,7 +39,7 @@ Field Solvers
.. toctree::
:maxdepth: 1

pic
maxwell_solvers
kinetic_fluid_hybrid_model

Grid & Geometries
Expand All @@ -45,6 +59,7 @@ Species Representations
.. toctree::
:maxdepth: 1

kinetic_particles
cold_fluid_model

Multiphysics Processes
Expand Down
110 changes: 110 additions & 0 deletions Docs/source/theory/kinetic_particles.rst
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,110 @@
.. _theory-kinetic-particles:

Kinetic Particles
=================

.. _theory-kinetic-particles-push:

Particle push
-------------

A centered finite-difference discretization of the Newton-Lorentz
equations of motion is given by

.. math::
\frac{\mathbf{x}^{i+1}-\mathbf{x}^{i}}{\Delta t} = \mathbf{v}^{i+1/2},
:label: leapfrog_x
.. math::
\frac{\gamma^{i+1/2}\mathbf{v}^{i+1/2}-\gamma^{i-1/2}\mathbf{v}^{i-1/2}}{\Delta t} = \frac{q}{m}\left(\mathbf{E}^{i}+\mathbf{\bar{v}}^{i}\times\mathbf{B}^{i}\right).
:label: leapfrog_v
In order to close the system, :math:`\bar{\mathbf{v}}^{i}` must be
expressed as a function of the other quantities. The two implementations that have become the most popular are presented below.

.. _theory-kinetic-particles-push-boris:

Boris relativistic velocity rotation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The solution proposed by Boris :cite:p:`kp-BorisICNSP70` is given by

.. math::
\mathbf{\bar{v}}^{i} = \frac{\gamma^{i+1/2}\mathbf{v}^{i+1/2}+\gamma^{i-1/2}\mathbf{v}^{i-1/2}}{2\bar{\gamma}^{i}}
:label: boris_v
where :math:`\bar{\gamma}^{i}` is defined by :math:`\bar{\gamma}^{i} \equiv (\gamma^{i+1/2}+\gamma^{i-1/2} )/2`.

The system (:eq:`leapfrog_v`, :eq:`boris_v`) is solved very
efficiently following Boris’ method, where the electric field push
is decoupled from the magnetic push. Setting :math:`\mathbf{u}=\gamma\mathbf{v}`, the
velocity is updated using the following sequence:

.. math::
\begin{aligned}
\mathbf{u^{-}} & = \mathbf{u}^{i-1/2}+\left(q\Delta t/2m\right)\mathbf{E}^{i}
\\
\mathbf{u'} & = \mathbf{u}^{-}+\mathbf{u}^{-}\times\mathbf{t}
\\
\mathbf{u}^{+} & = \mathbf{u}^{-}+\mathbf{u'}\times2\mathbf{t}/(1+\mathbf{t}^{2})
\\
\mathbf{u}^{i+1/2} & = \mathbf{u}^{+}+\left(q\Delta t/2m\right)\mathbf{E}^{i}
\end{aligned}
where :math:`\mathbf{t}=\left(q\Delta t/2m\right)\mathbf{B}^{i}/\bar{\gamma}^{i}` and where
:math:`\bar{\gamma}^{i}` can be calculated as :math:`\bar{\gamma}^{i}=\sqrt{1+(\mathbf{u}^-/c)^2}`.

The Boris implementation is second-order accurate, time-reversible and fast. Its implementation is very widespread and used in the vast majority of PIC codes.

.. _theory-kinetic-particles-push-vay:

Vay Lorentz-invariant formulation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

It was shown in :cite:t:`kp-Vaypop2008` that the Boris formulation is
not Lorentz invariant and can lead to significant errors in the treatment
of relativistic dynamics. A Lorentz invariant formulation is obtained
by considering the following velocity average

.. math::
\mathbf{\bar{v}}^{i} = \frac{\mathbf{v}^{i+1/2}+\mathbf{v}^{i-1/2}}{2}.
:label: new_v
This gives a system that is solvable analytically (see :cite:t:`kp-Vaypop2008`
for a detailed derivation), giving the following velocity update:

.. math::
\mathbf{u^{*}} = \mathbf{u}^{i-1/2}+\frac{q\Delta t}{m}\left(\mathbf{E}^{i}+\frac{\mathbf{v}^{i-1/2}}{2}\times\mathbf{B}^{i}\right),
:label: pusher_gamma
.. math::
\mathbf{u}^{i+1/2} = \frac{\mathbf{u^{*}}+\left(\mathbf{u^{*}}\cdot\mathbf{t}\right)\mathbf{t}+\mathbf{u^{*}}\times\mathbf{t}}{1+\mathbf{t}^{2}},
:label: pusher_upr
where

.. math::
\begin{align}
\mathbf{t} & = \boldsymbol{\tau}/\gamma^{i+1/2},
\\
\boldsymbol{\tau} & = \left(q\Delta t/2m\right)\mathbf{B}^{i},
\\
\gamma^{i+1/2} & = \sqrt{\sigma+\sqrt{\sigma^{2}+\left(\boldsymbol{\tau}^{2}+w^{2}\right)}},
\\
w & = \mathbf{u^{*}}\cdot\boldsymbol{\tau},
\\
\sigma & = \left(\gamma'^{2}-\boldsymbol{\tau}^{2}\right)/2,
\\
\gamma' & = \sqrt{1+(\mathbf{u}^{*}/c)^{2}}.
\end{align}
This Lorentz invariant formulation
is particularly well suited for the modeling of ultra-relativistic
charged particle beams, where the accurate account of the cancellation
of the self-generated electric and magnetic fields is essential, as
shown in :cite:t:`kp-Vaypop2008`.

.. bibliography::
:keyprefix: kp-
116 changes: 3 additions & 113 deletions Docs/source/theory/pic.rst → Docs/source/theory/maxwell_solvers.rst
View file
Original file line number Diff line number Diff line change
@@ -1,14 +1,7 @@
.. _theory-pic:
.. _theory-maxwell-solvers:

Particle-in-Cell Method
=======================

.. _fig-pic:

.. figure:: PIC.png
:alt: [fig:PIC] The Particle-In-Cell (PIC) method follows the evolution of a collection of charged macro-particles (positively charged in blue on the left plot, negatively charged in red) that evolve self-consistently with their electromagnetic (or electrostatic) fields. The core PIC algorithm involves four operations at each time step: 1) evolve the velocity and position of the particles using the Newton-Lorentz equations, 2) deposit the charge and/or current densities through interpolation from the particles distributions onto the grid, 3) evolve Maxwell’s wave equations (for electromagnetic) or solve Poisson’s equation (for electrostatic) on the grid, 4) interpolate the fields from the grid onto the particles for the next particle push. Additional “add-ons” operations are inserted between these core operations to account for additional physics (e.g. absorption/emission of particles, addition of external forces to account for accelerator focusing or accelerating component) or numerical effects (e.g. smoothing/filtering of the charge/current densities and/or fields on the grid).

The Particle-In-Cell (PIC) method follows the evolution of a collection of charged macro-particles (positively charged in blue on the left plot, negatively charged in red) that evolve self-consistently with their electromagnetic (or electrostatic) fields. The core PIC algorithm involves four operations at each time step: 1) evolve the velocity and position of the particles using the Newton-Lorentz equations, 2) deposit the charge and/or current densities through interpolation from the particles distributions onto the grid, 3) evolve Maxwell’s wave equations (for electromagnetic) or solve Poisson’s equation (for electrostatic) on the grid, 4) interpolate the fields from the grid onto the particles for the next particle push. Additional “add-ons” operations are inserted between these core operations to account for additional physics (e.g. absorption/emission of particles, addition of external forces to account for accelerator focusing or accelerating component) or numerical effects (e.g. smoothing/filtering of the charge/current densities and/or fields on the grid).
Maxwell solvers
===============

In the *electromagnetic particle-in-cell method* :cite:p:`pt-Birdsalllangdon,pt-HockneyEastwoodBook`,
the electromagnetic fields are solved on a grid, usually using Maxwell’s
Expand Down Expand Up @@ -51,109 +44,6 @@ on the grid from the particles’ positions and velocities, while the
electric and magnetic field components are interpolated from the grid
to the particles’ positions for the velocity update.

.. _theory-pic-push:

Particle push
-------------

A centered finite-difference discretization of the Newton-Lorentz
equations of motion is given by

.. math::
\frac{\mathbf{x}^{i+1}-\mathbf{x}^{i}}{\Delta t} = \mathbf{v}^{i+1/2},
:label: leapfrog_x
.. math::
\frac{\gamma^{i+1/2}\mathbf{v}^{i+1/2}-\gamma^{i-1/2}\mathbf{v}^{i-1/2}}{\Delta t} = \frac{q}{m}\left(\mathbf{E}^{i}+\mathbf{\bar{v}}^{i}\times\mathbf{B}^{i}\right).
:label: leapfrog_v
In order to close the system, :math:`\bar{\mathbf{v}}^{i}` must be
expressed as a function of the other quantities. The two implementations that have become the most popular are presented below.

.. _theory-pic-push-boris:

Boris relativistic velocity rotation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The solution proposed by Boris :cite:p:`pt-BorisICNSP70` is given by

.. math::
\mathbf{\bar{v}}^{i} = \frac{\gamma^{i+1/2}\mathbf{v}^{i+1/2}+\gamma^{i-1/2}\mathbf{v}^{i-1/2}}{2\bar{\gamma}^{i}}
:label: boris_v
where :math:`\bar{\gamma}^{i}` is defined by :math:`\bar{\gamma}^{i} \equiv (\gamma^{i+1/2}+\gamma^{i-1/2} )/2`.

The system (:eq:`leapfrog_v`, :eq:`boris_v`) is solved very
efficiently following Boris’ method, where the electric field push
is decoupled from the magnetic push. Setting :math:`\mathbf{u}=\gamma\mathbf{v}`, the
velocity is updated using the following sequence:

.. math::
\begin{aligned}
\mathbf{u^{-}} & = \mathbf{u}^{i-1/2}+\left(q\Delta t/2m\right)\mathbf{E}^{i}
\\
\mathbf{u'} & = \mathbf{u}^{-}+\mathbf{u}^{-}\times\mathbf{t}
\\
\mathbf{u}^{+} & = \mathbf{u}^{-}+\mathbf{u'}\times2\mathbf{t}/(1+\mathbf{t}^{2})
\\
\mathbf{u}^{i+1/2} & = \mathbf{u}^{+}+\left(q\Delta t/2m\right)\mathbf{E}^{i}
\end{aligned}
where :math:`\mathbf{t}=\left(q\Delta t/2m\right)\mathbf{B}^{i}/\bar{\gamma}^{i}` and where
:math:`\bar{\gamma}^{i}` can be calculated as :math:`\bar{\gamma}^{i}=\sqrt{1+(\mathbf{u}^-/c)^2}`.

The Boris implementation is second-order accurate, time-reversible and fast. Its implementation is very widespread and used in the vast majority of PIC codes.

.. _theory-pic-push-vay:

Vay Lorentz-invariant formulation
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

It was shown in :cite:t:`pt-Vaypop2008` that the Boris formulation is
not Lorentz invariant and can lead to significant errors in the treatment
of relativistic dynamics. A Lorentz invariant formulation is obtained
by considering the following velocity average

.. math::
\mathbf{\bar{v}}^{i} = \frac{\mathbf{v}^{i+1/2}+\mathbf{v}^{i-1/2}}{2}.
:label: new_v
This gives a system that is solvable analytically (see :cite:t:`pt-Vaypop2008`
for a detailed derivation), giving the following velocity update:

.. math::
\mathbf{u^{*}} = \mathbf{u}^{i-1/2}+\frac{q\Delta t}{m}\left(\mathbf{E}^{i}+\frac{\mathbf{v}^{i-1/2}}{2}\times\mathbf{B}^{i}\right),
:label: pusher_gamma
.. math::
\mathbf{u}^{i+1/2} = \frac{\mathbf{u^{*}}+\left(\mathbf{u^{*}}\cdot\mathbf{t}\right)\mathbf{t}+\mathbf{u^{*}}\times\mathbf{t}}{1+\mathbf{t}^{2}},
:label: pusher_upr
where

.. math::
\begin{align}
\mathbf{t} & = \boldsymbol{\tau}/\gamma^{i+1/2},
\\
\boldsymbol{\tau} & = \left(q\Delta t/2m\right)\mathbf{B}^{i},
\\
\gamma^{i+1/2} & = \sqrt{\sigma+\sqrt{\sigma^{2}+\left(\boldsymbol{\tau}^{2}+w^{2}\right)}},
\\
w & = \mathbf{u^{*}}\cdot\boldsymbol{\tau},
\\
\sigma & = \left(\gamma'^{2}-\boldsymbol{\tau}^{2}\right)/2,
\\
\gamma' & = \sqrt{1+(\mathbf{u}^{*}/c)^{2}}.
\end{align}
This Lorentz invariant formulation
is particularly well suited for the modeling of ultra-relativistic
charged particle beams, where the accurate account of the cancellation
of the self-generated electric and magnetic fields is essential, as
shown in :cite:t:`pt-Vaypop2008`.

.. _theory-pic-mwsolve:

Field solve
Expand Down

0 comments on commit 367f97e

Please sign in to comment.