Particle In Cell (physics basics) - openpixi/openpixi_pic GitHub Wiki

PIC codes are used to simulate charged, interacting particles (Plasmas). The main idea behind the PIC concept is that electric and magnetic fields are defined at discrete grid points, but the particle positions remain in the continuous domain. The most basic PIC code would have square grid cells where electric (E) and magnetic (B) fields are defined on the four surrounding points of each cell. A more advanced concept is the Finite-difference time-domain method (FDTD) introduced 1966 by K. Yee, where E and B fields are displaced in the spacial and time domains allowing for leapfrog integration. Modern PIC codes also employ adaptive and non orthogonal grids. A comprehensive introduction to plasma simulations can be found in "Plasma Physics via Computer Simulation", Birdsall (1991) a basic overview is given in Verboncoeur (2005).

A PIC code consists of four components

1) Particle mover 
2) Interpolator (particles to grid)
3) Field solver 
4) Interpolator (grid to particles)

###Particle mover Equations of motion need to be solved for each particle. Currently the folowing algorithms are implemented: Euler, Semi-implicit Euler, Euler Richardson, Leap frog and a modification of the latter which was introduced by Boris in 1970. There are also "damped" versions of the leap frog solvers which implement the velocity dependent drag force in a way so that the algorithm does not degrade to a lower accuracy level. This is only needed for the leap frog solvers because there , velocity and positions are defined at different times i.e. at execution time the algorithm knows x(t) and v(t+dt/2). In the non-damped versions the dragging is implemented in a more naive way. The relativistic solvers take care of the relativistic effects arising at high particle velocities.

An other important feature are the prepare and complete methods. These synchronize the velocity and the position of the particle in time. This is needed to perform collisions among particles, between particles and boundaries or in any other event where the velocity of the particle is needed at the same time as the particle position.

The Boris algorithm is second-order accurate in the presence of magnetic fields and is therefore the most widely used algorithm in PIC codes. Derivations and proofs can be found in "Plasma Physics via Computer Simulation", Birdsall (1991).

###Interpolator (particles to grid) Charged particles produce electrostatic fields when they are at rest. One can calculate these fields by solving the Poisson equation. Currently this is done with a Fast Fourier Transform.

Moving charged particles produce currents which in turn produce magnetic and electric fields. They can be calculated by solving Maxwell's equations, which are first order partial differential equations. Since the fields and currents (current desity on each grid point which corresponds to the current in the according cell) are defined in the descrete domain an interpolation from the continous domain of the particles in needed to couple the particle motion to the current. This can be done in several ways.

One can interpolate the particle movement in a way so that the equation of continuity, i.e. the current through cell boundaries is equal to the decrease of charge in that cell, is satisfied. A charge conserving interpolation method for square unit size charges is described in Villasenor (1991). A charge conserving method for general form factors is described in Esirkepov (2000). If these interpolation methods are used then there is no need for further corrections of the electric or magnetic fields, Maxwell's equations take care of that. Note that this still imposes certain restrictions on the grid and timestep size through the Courant-Friedrichs-Lewy condition.

An other method would be to calculate the currents in a non conserving way but subtract the "unphysical" part of the electric from the calculated one. To determine the "unphysical" part one has to solve the poisson equation at each timestep and therefore determine the "electrostatic" (the rotation free) part of the electric field. Further details on this can be found here. There one can also find information on the "Hyperbolic cleaning" method.

Except of the charge conserving algorithm there is an other method in the code called "cloud in cell", this method calculates the currents based on the current (either v(t) or v(t+dt/2)) velocity of the particle. This quantity is then weighted to the four surrounding grid points depending on the distance of the particle to each point.

###Field solver This component must be implemented according to the FDTD method. Currently this component is not completely finished because the E and B fields are calculated at the same time and not in a leapfrog fashion.

###Interpolator (grid to particles) Here it is important that the same form factor is used as in the other interpolator component otherwise the particles will exhibit a self force. It is also important to take into account that the E and B fields are defined at different grid points and different times (this in not implemented at the moment!). A description of an interpolation scheme can be found in Verboncoeur (2005).

###Further notes We use the central difference method throughout all calculations.

At the moment the boundary conditions are "hardcoded" in each component, to be able to switch between different boundary conditions drastic changes to some parts of the code will be needed.