In a similar way, the acceleration update profits from the structure of the electromagnetic forces acting on the particles. In the nonrelativistic Vlasov equation, the action of the Lorentz force transforms the phase space distribution inside one simulation time step like a solid rotator in velocity space: the magnetic field causes a gyration of the distribution function by the Larmor motion
while the electric field causes (1) acceleration and (2) drift motion perpendicular to the magnetic field direction.
Hence, the overall velocity space transformation is a rotation around the local magnetic field direction, with the gyration center given by the drift velocity. This transformation can be decomposed into three successive shear motions along the coordinate axes, thus repeating the same fundamental algorithmic structure as the spatial translation update. The acceleration update affects the velocity space completely locally within each spatial simulation cell.
In Vlasiator 5, the default solver uses a 5th order scheme for the acceleration of f w.r.t. , a 3rd order scheme for the translation of f w.r.t. . The Strang splitting scheme in time is 2nd order.
Actually this part is not hard to implement. Maybe there are some tricks for multi-dimensions.
Interestingly, the same idea appears in theoretical plasma physics for studying the phase space evolution. An extremely hard to solve problem can become quite easy and straightforward by moving in the phase space and shift your calculation completely to another spatial-temporal location.
The magnetic field is propagated using the algorithm by Londrillo and del Zanna (2004), which is a 2nd-order accurate (both in time and space) upwind constrained transport method. All reconstructions between volume-, face- and edge-averaged values follow Balsara+ (2009). The face-averaged values of are propagated forward in time using the integral form of Faraday's law
where the integral on the LHS is taken over the cell face and the line integral is evaluated along the contour of that face. Once has been computed based on Ohm's law, it is easy to propagate using the discretized form of the above. When computing each component of on an edge, the solver computes the candidate values on the four neighboring cells of the edge. In the supermagnetosonic case, when the plasma velocity exceeds the speed of the fast magnetosonic wave mode, the upwinded value from one of the cells is used. In the submagnetosonic case the value is computed as a weighted average of the electric field on the four cells and a diffusive flux is added to stabilize the scheme.
For the time integration a order Runge–Kutta method is used. To propagate the field from to we need the and values at both and . With the leap-frog algorithm there are no real values for the distribution function at these times. A order accurate interpolation is used to compute the required values:
The field solver also contributes to the dynamic computation of the timestep. With the simplest form of Ohm's law the fastest characteristic speed is the speed of the fast magnetosonic wave mode and we use that speed to compute the maximum timestep allowed by the field solver. For the field solver Courant numbers is used, as higher values cause numerical instability of the scheme.
In Vlasiator the magnetic field has been split into a perturbed field updated during the simulations, and a static background field. The electric field is computed based on the total magnetic field and all changes to the magnetic field are only added to the perturbed part of the magnetic field. The background field must be curl-free and thus the Hall term can be computed based on the perturbed part only. This avoids numerical integration errors arising from strong background field gradients. In magnetospheric simulations the background field consists of the Earth's dipole, as well as a constant IMF in all cells.
When adapting to AMR, Vlasiator developers made the decision of first keeping the same field solver on the finest level only. This makes the efficient implementation easier, at the cost of memory usage.
A typical ion population is localized in velocity space, for example a Maxwellian distribution, and a large portion of velocity space is effectively empty. A key technique in saving memory is the sparse velocity grid storage. The velocity grid is divided into velocity blocks comprising 4x4x4 velocity cells. The sparse representation is done at this level: either a block exists with all of its 64 velocity cells, or it does not exist at all. In the sparse representation we define that a block has content if any of its 64 velocity cells has a density above a user-specified threshold value. A velocity block exists if it has content, or if it is a neighbor to a block with content in any of the six dimensions. In the velocity space all 26 nearest neighbors are included, while in ordinary space blocks within a distance of 2 in the 6 normal directions and a distance of 1 in nodal directions are included.
There are two sources of loss when propagating the distribution function:
fluxes that flow out of the velocity grid;
distributions that go below the storage threshold.
Usually the term dominates.
|||Eulerian and Lagrangian descriptions of field appear most notably in fluid mechanics.|
|||This can be easily observed with a small velocity space, where the moments calculated from the distribution functions deviates from the analytical values.|