Merzbild.jl public API reference

Particle{D}

A structure to store information about a single particle with a D-dimensional position vector.

Fields

• w: the computational weight of the particle
• v: the 3-dimensional velocity vector of the particle
• x: the D-dimensional position of the particle

ParticleVector{D}

The structure used to store particles, sort and keep track of particle indices, and keep track of unused particles. The lengths of the particles, index, cell, and buffer vectors are all the same (and stay the same during resizing of a ParticleVector instance). Only the first nbuffer elements of the buffer vector store indices of the actually unused particles. D determines the dimension of the particles' position vectors.

Accessing ParticleVector{D}[i] will return a Particle{D} instance, with the actual particle returned being ParticleVector.particles[ParticleVector.index[i]].

Fields

• particles: the vector of particles of a single species
• index: the vector of indices of the particles (these are sorted in grid sorting, not the particles themselves)
• cell: the vector storing information in which cell a particle is located (used in grid sorting routines)
• buffer: a last-in-first-out (LIFO) queue keeping track of pre-allocated but unused particles
• nbuffer: the number of elements in the buffer

ParticleVector{D}(np) where D

Create an empty ParticleVector{D} instance of length np (all vectors will have length np), filled with particles with weight 0, velocity 0, position 0.

Positional arguments

• np: the length of the ParticleVector{D} instance to create

Base.getindex(pv::ParticleVector{D}, i) where D

Returns the underlying particle in a ParticleVector instance with index i.

Is usually called as ParticleVector[i].

Positional arguments

• pv: ParticleVector instance
• i: the index of the particle to be selected

Base.setindex!(pv::ParticleVector{D}, p::Particle{D}, i::Integer) where D

Set the underlying particle in a ParticleVector instance with index i to a new particle.

Is usually called as ParticleVector[i] = p.

Positional arguments

• pv: ParticleVector instance
• p: the Particle instance to write
• i: the index of the particle to be written to

Base.length(pv::ParticleVector{D}) where D

Returns the length of a ParticleVector instance.

Is usually called as length(ParticleVector).

Positional arguments

• pv: ParticleVector instance

Base.resize!(pv::ParticleVector{D}, n::Integer) where D

Resize a ParticleVector instance, taking care of the indices, buffer, and creating placeholder new particles with weight 0, velocity 0, and position 0.

Positional arguments

• pv: ParticleVector instance
• n: the new length of the ParticleVector instance (i.e. the length of all the vector fields of the instance)

ParticleIndexer

The structure used to index particles of a given species in a given cell. It is assumed that the particle indices are contiguous; they may be split across two groups of contiguous indices.

Fields

• n_local: the number of particles in the cell
• start1: the first index in the first group of particle indices
• end1: the last index in the first group of particle indices
• n_group1: the number of particles in the first group (n_group1 = end1 - start1 + 1)
• start2: the first index in the second group of particle indices, if no particles are present in the group, it should be <= 0
• end2: the last index in the second group of particle indices
• n_group2: the number of particles in the second group (n_group2 = end2 - start2 + 1, unless no particles are present in the second group)

ParticleIndexer()

Create an empty ParticleIndexer. end1 and end2 are set to -1, the other fields are set to 0.

ParticleIndexer(n_particles)

Create a ParticleIndexer given a number of particles. All particles are in group 1, with indices starting from 1.

Positional arguments

• n_particles: the number of particles

ParticleIndexerArray

A structure to store an array of ParticleIndexer instances for each species in each cell. It also stores information about the whether the ParticleIndexer instances for each species are "contiguous". A set of ParticleIndexer instances (for a specific species) are "contiguous" if the following conditions are fulfilled (assuming all cells contain particles and in all cells have both n_group1>0 and n_group2>0, for a more detailed explanation, one is referred to the documentation on contiguous indexing):

• pia.indexer[cell,species].end1 + 1 == pia.indexer[cell+1,species].start1, cell = 1,...,n_cells - 1
• pia.indexer[n_cells,species].end1 + 1 == pia.indexer[1,species].start2
• pia.indexer[cell,species].end2 + 1 == pia.indexer[cell+1,species].start2, cell = 1,...,n_cells - 1

Fields

• n_cells: number of cells being indexed
• n_species: number of species being indexed
• indexer: the array of size (n_cells, n_species) (number of grid cells * number of species in the simulation) storing the ParticleIndexer instances
• n_total: vector of length n_species storing the total number of particles of each species
• contiguous: vector of length n_species storing a boolean flag whether the ParticleIndexer instances for a species are "contiguous"
• index_last: vector of length n_species storing the last valid index in the particle array for each species (tracks the highest index used, allowing for non-contiguous indexing when particles are deleted from the middle)

ParticleIndexerArray(indexer_arr::Array{ParticleIndexer,2}, n_total)

Creates a ParticleIndexerArray from a 2-D array of ParticleIndexer instances.

Positional arguments

• indexer_arr: the 2-D array of ParticleIndexer instances
• n_total: the vector of the total number of particles of each species

ParticleIndexerArray(n_cells::Integer, n_species::Integer)

Create an empty ParticleIndexerArray given the number of cells and species

Positional arguments

• n_cells: the number of grid cells
• n_species: the number of species

ParticleIndexerArray(n_particles::Integer)

Create a single-species/single-cell ParticleIndexerArray and set up the indexing for n_particles in group 1 in cell 1.

Positional arguments

• n_particles: the number of particles (integer number)

ParticleIndexerArray(n_particles::T) where T<:AbstractVector

Create a multi-species/single-cell ParticleIndexerArray and set up the indexing for n_particles[species] in group 1 in cell 1 for all species.

Positional arguments

• n_particles: the number of particles of each species (vector-like)

ParticleIndexerArray(grid, species_data::Array{Species})

Create an empty multi-species/multi-cell ParticleIndexerArray.

Positional arguments

• grid: the simulation grid
• species_data: array of Species data for all of the species in the simulation

squash_pia!(pv::ParticleVector{D}, pia, species) where D

Restore the continuity of indices in a ParticleVector and associated ParticleIndexerArray instance for a specific species. If for this species the instance has contiguous == true, nothing will be done.

Positional arguments

• pv: the ParticleVector
• pia: the ParticleIndexerArray instance
• species: the index of the species for which to restore continuity of indices

squash_pia!(particles, pia)

Restore the continuity of indices in a list of ParticleVectors and the associated ParticleIndexerArray instance for all species. If for a specific species the instance has contiguous == true, nothing will be done.

Positional arguments

• particles: the list of ParticleVectors for all species in the flow
• pia: the ParticleIndexerArray instance

restore_particle_ordering!(pv::ParticleVector{D}, inv_map::Vector{Int64}) where D

Restore the ordering of particles in a ParticleVector so that pv.index[i] == i. This function ensures that:

1. The index array is restored to identity mapping (pv.index[i] == i)
2. The buffer contains only unused particle indices in descending order
3. Used particles are placed in the first (length(pv.index) - pv.nbuffer) positions

Positional arguments

• pv: ParticleVector instance to restore ordering for
• inv_map: Vector of integers to store inverse map (inv_map[pv.index[i]] == i)

restore_particle_ordering!(pv::ParticleVector{D}, pia, species, inv_map::Vector{Int64}) where D

Restore the ordering of particles in a ParticleVector and update the associated ParticleIndexerArray. This function calls restore_particle_ordering! to restore the particle ordering and then updates pia.index_last[species] to reflect the new contiguous indexing.

Positional arguments

• pv: ParticleVector instance to restore ordering for
• pia: ParticleIndexerArray instance to update
• species: the index of the species for which to update index_last
• inv_map: Vector of integers to store inverse map (inv_map[pv.index[i]] == i)

count_disordered_particles(pv::ParticleVector{D}, pia, species) where D

Count number of particles of species for which pv.index[i] != i + offset(cell). The larger this count, the less orderly the layout of particles in memory, which can potentially lead to decreased performance. This only considers particles present in the simulation, i.e. i<=pia.n_total[species], assumes contiguous indexing, and that particles are only pointed to by the first group of a ParticleIndexer, which is the case after particles have been sorted.

In case use_offset is false, the value of offset is always 0, so the function simply counts all particles where pv.index[i] != i. If use_offset is set to true, for each cell cell, the offset is whilst iterating over the particles in a cell. So in case the indexing in a cell is simply offset by a constant value, or a subset of indices are offset by a constant value, and particles are still laid out continuously in memory, this provides a more accurate measurement of the degree of fragmentation of the particle array.

Positional arguments

• pv: the ParticleVector instance for which to compute the metric
• pia: the ParticleIndexerArray instance
• species: the index of the species for which the metric is computed

Keyword arguments

• use_offset: whether to account for cell-specific offsets in the indexing

Returns

The number of particles where the index to the index array and the value of the index array do not coincide.

check_unique_index(pv::ParticleVector{D}, pia, species) where D

Test that for all particles in a simulation, no two indices are the same, i.e. no two particles i and j, i!=j point to the same underlying particle. This function allocates a temporary array and is thus intended for debugging/verifying code and not for efficient simulations. It also checks that any particles present in the buffer are not pointed to by the indexing.

Positional arguments

• pv: the ParticleVector instance for which to check indexing
• pia: the ParticleIndexerArray instance
• species: the species for which to check indexing

Returns

If a particle exists to which more than 1 index is pointing, and particles in the buffer ARE NOT pointed to by the indicies, returns (false, n_index), where n_index is the number of indices pointing to the same particle.

If a particle exists to which more than 1 index is pointing, and particles in the buffer ARE pointed to by the indicies, returns (false, -n_index), where n_index is the number of indices pointing to the same particle.

If no particles exist to which more than 1 index is pointing, but a particle in a buffer is pointed to by the indexing, returns (false, -1).

If indexing is correct, returns (true, 0).

check_pia_is_correct(pia, species)

Check that a ParticleIndexerArray instance entries are correct. This means that for each cell for each species, the following should hold: * n_local == n_group1 + n_group2 * if n_group1 > 0, then n_group1 == end1 - start1 + 1 * if n_group1 == 0, then start1 == 0, end1 == -1 * if n_group2 > 0, then n_group2 == end2 - start2 + 1 * if n_group2 == 0, then start2 == 0, end2 == -1 * pia.n_total[species] == sum([pia.indexer[cell, species].n_local for cell in 1:n_cells]) and index_last[species] points to the last index used.

Positional arguments

• pia: the ParticleIndexerArray instance for which to check consistency
• species: the species for which to check indexing

Returns

If indexing is incorrect in a cell i, returns (false, i).

If the total number of particles as given by cell-wise particle indices is not equal to pia.n_total[species], returns (false, 0).

If index_last[species] < pia.n_total[species], returns (false, -1).

If index_last[species] is not the largest value of the indices pointed to by the groups, returns (false, -2).

If indexing is correct, returns (true, 0).

pretty_print_pia(pia, species)

Display a ParticleIndexerArray instance by showing the starting/ending indices of the groups over all cells for a specific species.

Positional arguments

• pia: the ParticleIndexerArray instance
• species: the index of the species for which the indices are displayed

check_unique_buffer(pv::ParticleVector{D}) where D

Test that the buffer elements in the active part of a buffer in a ParticleVector instance are unique (active part is pv.buffer[1:pv.nbuffer]). This function allocates a temporary dictionary and is thus intended for debugging/verifying code and not for efficient simulations.

Positional arguments

• pv: the ParticleVector instance for which to check buffer

Returns

If a non-unique index is in the active part of the buffer, the function returns (false, i), where i is the first location at which the non-unique index is encountered for the second time (so if it is present at locations i1, i2, i3, etc., i2 will be returned).

If buffer is correct, returns (true, 0).

MERZBILD_DATA_PATH

Absolute path to the data directory bundled with Merzbild.jl, holding the built-in species (particles.toml) and interaction (vhs.toml, pseudo_maxwell.toml, vss.toml) data files. Use it to load bundled data independently of the current working directory, e.g. load_species_data(joinpath(MERZBILD_DATA_PATH, "particles.toml"), "Ar").

Species

A structure to store information about a chemical species.

Fields

• name: the name of the species
• mass: the molecular mass of the species
• charge: the charge of the species in terms of elementary charge (i.e. 1, -1, etc.)
• charge_div_mass: the charge of the species divided by its mass, C/kg

Interaction

Structure to store interaction parameters for a 2-species interaction. The VHS model uses the following power law: $\sigma_{VHS} = C g^(1 - 2 \omega_{VHS})$, where $omega$ is the exponent of the VHS potential, and $C$ is the pre-computed factor: $C = \pi D_{VHS}^2 (2 T_{ref,VHS}/m_r)^{(\omega_{VHS} - 0.5)} \frac{1}{\Gamma(2.5 - \omega_{VHS})}$. The VSS model uses the same power law for the total cross-section, but a different scattering law.

The elastic scattering model is fixed for a species pair and is stored in the model field as a Merzbild.ScatteringModel enum value, so that the array of Interaction instances stays concretely typed.

Fields

• m_r: collision-reduced mass
• μ1: relative mass of the first species
• μ2: relative mass of the second species
• vhs_d: diameter for the VHS potential
• vhs_o: exponent for the VHS potential (omega)
• vhs_exp: pre-computed 1 - 2 vhs_o
• vhs_Tref: reference temperature for the VHS potential
• vhs_muref: reference viscosity for the VHS potential
• vhs_factor: pre-computed factor for calculation of the VHS cross-section
• vss_alpha: exponent of the VSS scattering law (equal to 1.0 for the models with isotropic scattering)
• vss_inv_alpha: pre-computed 1 / vss_alpha
• model: the ScatteringModel enum value of the elastic scattering model of the species pair

Interaction(m1::Float64, m2::Float64, vhs_d::Float64, vhs_o::Float64, vhs_Tref::Float64)

Construct an Interaction instance using the VHS model from the masses of the two species and the VHS parameters.

Positional arguments

• m1: molecular mass of the first species
• m2: molecular mass of the second species
• vhs_d: VHS diameter
• vhs_o: VHS exponent (omega)
• vhs_Tref: reference VHS temperature

Interaction(::VHS, m1::Float64, m2::Float64, vhs_d::Float64, vhs_o::Float64, vhs_Tref::Float64)

Construct an Interaction instance using the VHS model from the masses of the two species and the VHS parameters.

Positional arguments

• m1: molecular mass of the first species
• m2: molecular mass of the second species
• vhs_d: VHS diameter
• vhs_o: VHS exponent (omega)
• vhs_Tref: reference VHS temperature

Interaction(::VSS, m1::Float64, m2::Float64, vhs_d::Float64, vhs_o::Float64, vhs_Tref::Float64, vss_alpha::Float64)

Construct an Interaction instance using the VSS model from the masses of the two species, the VHS parameters (which define the total cross-section), and the VSS exponent (which defines the scattering law). The reference viscosity is corrected by the factor $(\alpha + 1)(\alpha + 2) / (6\alpha)$ accounting for the anisotropic scattering, see Merzbild.compute_vss_mu_ref_factor.

Positional arguments

• m1: molecular mass of the first species
• m2: molecular mass of the second species
• vhs_d: VHS diameter
• vhs_o: VHS exponent (omega)
• vhs_Tref: reference VHS temperature
• vss_alpha: VSS exponent (alpha)

load_species_data(species_filename, species_names)

Load a vector of species data (mass, charge, etc.) from a TOML file.

Positional arguments

• species_filename: the path to the TOML file containing the data
• species_names: a list of the names of the species for which to load the data

Returns

Vector of Species filled with data loaded from the file.

load_species_data(species_filename, species_name::String)

Load a vector of species data (mass, charge, etc.) from a TOML file for a single species.

Positional arguments

• species_filename: the path to the TOML file containing the data
• species_name: the name of the species for which to load the data

Returns

Vector of Species filled with data loaded from the file.

load_interaction_data(interactions_filename, species_data)

Load interaction data from a TOML file given a list of species' data (list of Species instances). It will load interaction data for all possible pair-wise interactions of the species in the list.

The resulting 2-D array has the interaction data for Species[i] with Species[k] in position [i,k]. It is not symmetric, as the relative collision masses μ1 and μ2 are swapped when comparing the Interaction instances in positions [i,k] and [k,i]. If no data is found, the function throws an error.

The elastic scattering model of a species pair is set by the optional model key of the pair's entry in the file, defaulting to the VHS model, see Merzbild.interaction_from_toml.

Positional arguments

• interactions_filename: the path to the TOML file containing the data
• species_data: list of Species instances for which to search for the interaction data

Returns

• 2-dimensional array of Interaction instances of size (n_species, n_species)

Throws

KeyError if interaction data not found in the file.

load_interaction_data(interactions_filename, species_data)

Load interaction data from a TOML file given a list of species' data (list of Species instances), filling in dummy VHS data in case no entry is found in the TOML file. Useful for interactions where the VHS model doesn't make sense, for example electron-neutral interactions.

It will load interaction data for all possible pair-wise interactions of the species in the list. The resulting 2-D array has the interaction data for Species[i] with Species[k] in position [i,k]. It is not symmetric, as the relative collision masses μ1 and μ2 are swapped when comparing the Interaction instances in positions [i,k] and [k,i].

The elastic scattering model of a species pair is set by the optional model key of the pair's entry in the file, defaulting to the VHS model, see Merzbild.interaction_from_toml; the dummy data is always filled in using the VHS model.

Positional arguments

• interactions_filename: the path to the TOML file containing the data
• species_data: list of Species instances for which to search for the interaction data
• dummy_vhs_d: value to use for the VHS diameter if no interaction data found in the file
• dummy_vhs_o: value to use for the VHS exponent if no interaction data found in the file
• dummy_vhs_Tref: value to use for the VHS reference temperature if no interaction data found in the file

Returns

• 2-dimensional array of Interaction instances of size (n_species, n_species)

load_interaction_data_with_dummy(interactions_filename, species_data)

Load interaction data from a TOML file given a list of species' data (list of Species instances), filling in dummy VHS data in case no entry is found in the TOML file. Useful for interactions where the VHS model doesn't make sense, for example electron-neutral interactions. Uses a value of 1e-10 for the dummy VHS diameter, 1.0 for the dummy VHS exponent, and 273.0 for the dummy VHS reference temperature.

It will load interaction data for all possible pair-wise interactions of the species in the list. The resulting 2-D array has the interaction data for Species[i] with Species[k] in position [i,k]. It is not symmetric, as the relative collision masses μ1 and μ2 are swapped when comparing the Interaction instances in positions [i,k] and [k,i].

Positional arguments

• interactions_filename: the path to the TOML file containing the data
• species_data: list of Species instances for which to search for the interaction data

Returns

• 2-dimensional array of Interaction instances of size (n_species, n_species)

load_species_and_interaction_data(species_filename, interactions_filename, species_names; fill_dummy=true)

Given a list of species' names, load the species and interaction data (filling with dummy data if needed).

Positional arguments

• species_filename: the path to the TOML file containing the species' data
• interactions_filename: the path to the TOML file containing the interaction data
• species_name: the name of the species for which to load the data

Keyword arguments

• fill_dummy: if true, fill interaction data with computed dummy values if no entry found for a species pair in the interaction file

Returns

• Vector of Species instances
• 2-dimensional array of Interaction instances of size (n_species, n_species)

AbstractScatteringModel

Abstract type for the elastic scattering models used by the DSMC/SWPM collision routines.

Each concrete model is a zero-field singleton type (VHS, VSS) that is used as a compile-time tag: the total cross-section Merzbild.sigma and the computation of the post-collision velocities Merzbild.scatter! are dispatched on it. As the tag is a zero-size value of a concrete type, it is not stored anywhere and the inner collision loops stay free of branching, dynamic dispatch, and allocations.

The model of a species pair is stored in the corresponding Interaction instance as a Merzbild.ScatteringModel enum value (so that the array of Interaction instances stays concretely typed); the collision routines convert it to the singleton tag exactly once per call, see Merzbild.@scattering_barrier.

VHS

Singleton type tag for the Variable Hard Sphere model: the total cross-section follows the power law $\sigma = C g^{1 - 2\omega}$ and the scattering is isotropic.

VSS

Singleton type tag for the Variable Soft Sphere model: the total cross-section is the same as for the VHS model, but the scattering is anisotropic, with the cosine of the deflection angle sampled as $\cos\chi = 2 R^{1/\alpha} - 1$, where $R$ is a uniformly distributed random number and $\alpha$ is the VSS exponent ($\alpha = 1$ recovers isotropic scattering).

References

• K. Koura, H. Matsumoto, Variable soft sphere molecular model for inverse-power-law or Lennard-Jones potential. Phys. Fluids A, 1991.

ScatteringModel ScatteringVHS=1 ScatteringVSS=2

Enum of the elastic scattering models, stored in an Interaction instance to define the model used for the species pair in question. ScatteringVHS corresponds to VHS, ScatteringVSS to VSS.

maxwellian(vx, vy, vz, m, T)

Evaluate the Maxwell distribution with temperature T for a species with mass m at a velocity (vx, vy, vz)

Positional arguments

• vx: x velocity
• vy: y velocity
• vz: z velocity
• m: species' mass
• T: temperature

bkw(vx, vy, vz, m, T, scaled_time)

Evaluate the Bobylev-Krook-Wu (BKW) distribution with temperature T for a species with mass m at a velocity (vx, vy, vz) and scaled time scaled_time

Positional arguments

• vx: x velocity
• vy: y velocity
• vz: z velocity
• m: species' mass
• T: temperature
• scaled_time: the scaled_time

sample_on_grid!(rng, vdf_func, particles::ParticleVector{D}, nv, m, T, n_total,
                xlo, xhi, ylo, yhi, zlo, zhi; v_mult=3.5, cutoff_mult=3.5, noise=0.0,
                v_offset=[0.0, 0.0, 0.0])

Sample particles by evaluating a distribution on a discrete velocity grid, considering only points inside a sphere of a given radius (the value of the VDF at points outside of the sphere will be 0.0). The values of the VDF at the grid points will then be the computational weights of the particles, and the particles velocities are taken to be the velocities of the corresponding grid nodes (with additional uniformly distributed noise). Note: this can produce a large amount of particles for fine grids (as the number of grid nodes scales as nv^3.) The grid is assumed to have the same number of nodes nv in each direction, and the extent is computed as v_mult * v_thermal, where v_thermal is the thermal velocity $\sqrt(2kT/m)$, and v_mult is a user-defined parameter. The positions of the particles are assumed to be randomly distributed in a cuboid.

Positional arguments

• rng: The random number generator
• vdf_func: the distribution function to be evaluated which takes the x, y, and z velocities as parameters
• particles: the Vector-like structure holding the particles
• nv: the number of grid nodes in each direction
• m: the molecular mass of the species
• T: the temperature used to compute the thermal velocity
• n_total: the total computational weight (number of physical particles) to be sampled
• xlo: the lower bound of the x coordinates of the particles
• xhi: the upper bound of the x coordinates of the particles
• ylo: the lower bound of the y coordinates of the particles
• yhi: the upper bound of the y coordinates of the particles
• zlo: the lower bound of the z coordinates of the particles
• zhi: the upper bound of the z coordinates of the particles

Keyword arguments

• v_mult: the value by which the thermal velocity is multiplied to compute the extent of the velocity grid
• cutoff_mult: the value by which the thermal velocity is multiplied to compute the radius for the sphere used to cut-off the higher velocities
• noise: controls the amount of noise added to the particle velocities (the noise is uniformly distributed on the interval [-noise*dv, noise*dv], where dv is the grid spacing)
• v_offset: the streaming velocity vector to be added to the particle velocities

Returns

• The number of particles created

sample_maxwellian_on_grid!(rng, particles::ParticleVector{D}, nv, m, T, n_total,
                           xlo, xhi, ylo, yhi, zlo, zhi; v_mult=3.5, cutoff_mult=3.5, noise=0.0,
                           v_offset=[0.0, 0.0, 0.0])

Sample particles by evaluating a Maxwellian on a discrete velocity grid, considering only points inside a sphere of a given radius (the value of the VDF at points outside of the sphere will be 0.0). The values of the VDF at the grid points will then be the computational weights of the particles, and the particles velocities are taken to be the velocities of the corresponding grid nodes (with additional uniformly distributed noise). Note: this can produce a large amount of particles for fine grids (as the number of grid nodes scales as nv^3.) The grid is assumed to have the same number of nodes nv in each direction, and the extent is computed as v_mult * v_thermal, where v_thermal is the thermal velocity $\sqrt(2kT/m)$, and v_mult is a user-defined parameter. The positions of the particles are assumed to be randomly distributed in a cuboid.

Positional arguments

• rng: The random number generator
• particles: the Vector-like structure holding the particles
• nv: the number of grid nodes in each direction
• m: the molecular mass of the species
• T: the temperature used to compute the thermal velocity
• n_total: the total computational weight (number of physical particles) to be sampled
• xlo: the lower bound of the x coordinates of the particles
• xhi: the upper bound of the x coordinates of the particles
• ylo: the lower bound of the y coordinates of the particles
• yhi: the upper bound of the y coordinates of the particles
• zlo: the lower bound of the z coordinates of the particles
• zhi: the upper bound of the z coordinates of the particles

Keyword arguments

• v_mult: the value by which the thermal velocity is multiplied to compute the extent of the velocity grid
• cutoff_mult: the value by which the thermal velocity is multiplied to compute the radius for the sphere used to cut-off the higher velocities
• noise: controls the amount of noise added to the particle velocities (the noise is uniformly distributed on the interval [-noise*dv, noise*dv], where dv is the grid spacing)
• v_offset: the streaming velocity vector to be added to the particle velocities

Returns

The function returns the number of particles created

sample_particles_equal_weight!(rng, particles, pia, cell, species,
                                    nparticles, m, T, Fnum, xlo, xhi, ylo, yhi, zlo, zhi;
                                    distribution=:Maxwellian, vx0=0.0, vy0=0.0, vz0=0.0)

Sample equal-weight particles of a specific species in a specific cell from a distribution. The positions of the particles are assumed to be randomly distributed in a cuboid. Note: this does not work if applied twice in a row to the same cell.

Positional arguments

• rng: the random number generator
• particles: the Vector-like structure holding the particles
• pia: the ParticleIndexerArray
• cell: the cell index
• species: the species index
• nparticles: the number of particles to sample
• m: species' mass
• T: temperature
• Fnum: the computational weight of the particles
• xlo: the lower bound of the x coordinates of the particles
• xhi: the upper bound of the x coordinates of the particles
• ylo: the lower bound of the y coordinates of the particles
• yhi: the upper bound of the y coordinates of the particles
• zlo: the lower bound of the z coordinates of the particles
• zhi: the upper bound of the z coordinates of the particles

Keyword arguments

• distribution: the distribution to sample from (either :Maxwellian or :BKW)
• vx0: the x-velocity offset to add to the particle velocities
• vy0: the y-velocity offset to add to the particle velocities
• vz0: the z-velocity offset to add to the particle velocities

sample_particles_equal_weight!(rng, grid1duniform, particles, pia, species, species_data, ppc::Integer, T, Fnum)

Sample particles from a Maxwellian distribution in each cell of a 1-D uniform grid given the number of particles per cell.

Positional arguments

• rng: the random number generator
• grid1duniform: the 1-D uniform grid
• particles: the ParticleVector of particles
• pia: the ParticleIndexerArray
• species: the index of the species to be sampled for
• species_data: Vector of Species data
• ppc: number of particles per cell to be sampled
• T: the temperature
• Fnum: the computational weight of the particles

sample_particles_equal_weight!(rng, grid1duniform, particles, pia, species, species_data, ppc::Integer, T, Fnum, cell_chunk)

Sample particles from a Maxwellian distribution in a sequentially ordered subset of cells of a 1-D uniform grid given the number of particles per cell.

Positional arguments

• rng: the random number generator
• grid1duniform: the 1-D uniform grid
• particles: the ParticleVector of particles
• pia: the ParticleIndexerArray
• species: the index of the species to be sampled for
• species_data: Vector of Species data
• ppc: number of particles per cell to be sampled
• T: the temperature
• Fnum: the computational weight of the particles
• cell_chunk: the list of cell indices or range in which to sample particles, should be ordered in increasing order

sample_particles_equal_weight!(rng, grid1duniform, particles, pia, species, species_data, ndens::Float64, T, Fnum)

Sample particles from a Maxwellian distribution in each cell in a sequentially ordered subset of cells of a 1-D uniform grid given the target number density. If the computed number of particles is not an integer value, the fractional remainder is used to probabilistically sample an extra particle, so that on average, the expected number density is achieved.

Positional arguments

• rng: the random number generator
• grid1duniform: the 1-D uniform grid
• particles: the ParticleVector of particles
• pia: the ParticleIndexerArray
• species: the index of the species to be sampled for
• species_data: Vector of Species data
• ndens: target number density
• T: the temperature
• Fnum: the computational weight of the particles

sample_particles_equal_weight!(rng, grid1duniform, particles, pia, species, species_data, ndens::Float64, T, Fnum, cell_chunk)

Sample particles from a Maxwellian distribution in each cell of 1-D uniform grid given the target number density. If the computed number of particles is not an integer value, the fractional remainder is used to probabilistically sample an extra particle, so that on average, the expected number density is achieved.

Positional arguments

• rng: the random number generator
• grid1duniform: the 1-D uniform grid
• particles: the ParticleVector of particles
• pia: the ParticleIndexerArray
• species: the index of the species to be sampled for
• species_data: Vector of Species data
• ndens: target number density
• T: the temperature
• Fnum: the computational weight of the particles
• cell_chunk: the list of cell indices or range in which to sample particles, should be ordered in increasing order

sample_particles_phase_box_weighted!(rng, particles::ParticleVector{D}, pia, cell, species,
                                     nparticles, m, T, n_total, xlo, xhi, ylo, yhi, zlo, zhi;
                                     v_mult=3.5, vx0=0.0, vy0=0.0, vz0=0.0)

Sample variable-weight particles of a specific species in a specific cell from a Maxwellian. The velocities are sampled in a uniform box in phase space with extent ±v_mult * v_thermal in each direction, where v_thermal is the thermal velocity $\sqrt(2kT/m)$, and v_mult is a user-defined parameter, and their weights are computed proportional to the value of a Maxwell distribution at temperature T at those velocities. A streaming velocity component can be added to all of the sampled particles' velocities. The positions of the particles are assumed to be randomly distributed in a cuboid.

Positional arguments

• rng: the random number generator
• particles: the Vector-like structure holding the particles
• pia: the ParticleIndexerArray
• cell: the cell index
• species: the species index
• nparticles: the number of particles to sample
• m: species' mass
• T: temperature
• n_total: the total computational weight (number of physical particles) to be sampled
• xlo: the lower bound of the x coordinates of the particles
• xhi: the upper bound of the x coordinates of the particles
• ylo: the lower bound of the y coordinates of the particles
• yhi: the upper bound of the y coordinates of the particles
• zlo: the lower bound of the z coordinates of the particles
• zhi: the upper bound of the z coordinates of the particles

Keyword arguments

• v_mult: the value by which the thermal velocity is multiplied to compute the extent of the velocity box
• vx0: the x-velocity offset to add to the particle velocities
• vy0: the y-velocity offset to add to the particle velocities
• vz0: the z-velocity offset to add to the particle velocities

PhysProps

Structure to store computed physical properties in a physical cell.

Fields

• ndens_not_Np: whether the n field stores number density (if true) and not the number of physical particles in a cell (if false)
• n_cells: number of physical cells
• n_species: number of species
• lpa: length of the particle array (vector of length n_species)
• np: number of particles (array of shape (n_cells, n_species))
• n: number density or number of physical particles in a cell (array of shape (n_cells, n_species))
• v: per-species flow velocity in a cell (array of shape (3, n_cells, n_species))
• T: per-species temperature in a cell (array of shape (n_cells, n_species))

PhysProps(n_cells, n_species; ndens_not_Np=false)

Construct physical properties given the number of cells and species.

Positional arguments

• n_cells: number of cells
• n_species: number of species

Keyword arguments

• ndens_not_Np: whether the n field stores number density (if true) and not the number of physical particles in a cell (if false)

PhysProps(pia::ParticleIndexerArray; ndens_not_Np=false)

Construct physical properties given a ParticleIndexerArray instance,.

Positional arguments

• pia: the ParticleIndexerArray instance

Keyword arguments

• ndens_not_Np: whether the n field stores number density (if true) and not the number of physical particles in a cell (if false)

SurfProps

Structure to store computed surface properties.

Fields

• n_elements: number of surface elements
• n_species: number of species
• areas: the vector of surface element areas
• inv_areas: the vector of the inverse surface element areas
• normals: the array of surface element normals with shape (3, n_elements)
• np: number of particles impacting the surface elements, array of shape (n_elements, n_species)
• flux_incident: incident mass flux per surface element, array of shape (n_elements, n_species)
• flux_reflected: reflected mass flux per surface element, array of shape (n_elements, n_species)
• force: force acting per surface element, array of shape (3, n_elements, n_species)
• normal_pressure: normal pressure per surface element, array of shape (n_elements, n_species)
• shear_pressure: shear pressure per surface element, array of shape (3, n_elements, n_species)
• kinetic_energy_flux: kinetic energy flux per surface element, array of shape (n_elements, n_species)

SurfProps(n_elements, n_species, areas, normals)

Construct a struct for holding computed surface properties for n_elements surface elements, n_species species. The surface elements' areas and normals are given by areas and normals, respectively.

Positional arguments

• n_elements: number of surface elements
• n_species: number of species
• areas: vector of areas of surface elements
• normals: vector of normals to surface elements

SurfProps(pia, grid::Grid1DUniform)

Create a SurfProps struct for a 1-D grid, with element 1 corresponding to the left wall and element 2 corresponding to the right wall. The areas of the wall are assumed to be equal to 1, the normals are parallel to the x axis.

Positional arguments

• pia: the ParticleIndexerArray instance
• grid: the Grid1DUniform grid

FluxProps

Structure to store computed flux densities in a physical cell. The kinetic energy flux density of a species is computed as $\frac{m}{2V} \sum_i w_i \mathbf{C}_i C_i^2$, where $\mathbf{C}_i$ is the peculiar velocity of a particle, $m$ is the molecular mass of the species, $V$ is the cell volume. The $k,l$ component of the the momentum flux density tensor is computed as $\frac{m}{2V} \sum_i w_i \C_{i,k} \C_{i,l}$. The $xx$, $yy$, $zz$ components thereof are stored in diagonal_momentum_flux. The $xy$, $xz$, $yz$ components thereof are stored in off_diagonal_momentum_flux.

Fields

• n_cells: number of physical cells
• n_species: number of species
• kinetic_energy_flux: per-species kinetic energy flux in a cell (array of shape (3, n_cells, n_species))
• diagonal_momentum_flux: diagonal components of the per-species momentum flux tensor in a cell (array of shape (3, n_cells, n_species))
• off_diagonal_momentum_flux: off-diagonal components of the per-species momentum flux tensor in a cell (array of shape (3, n_cells, n_species))

FluxProps(n_cells, n_species)

Construct a FluxProps instance given the number of cells and species.

Positional arguments

• n_cells: number of cells
• n_species: number of species

FluxProps(pia)

Construct a FluxProps instance given a ParticleIndexerArray instance.

Positional arguments

• pia: the ParticleIndexerArray instance

ElectrostaticFieldProps

Structure to store the electrostatic field quantities on the nodes of a grid. The quantities are node-centered: for a 1D uniform grid with n_cells cells, the nodes are located at $x_j = (j-1) \Delta x$, $j = 1 \dots n_{\textrm{cells}}+1$, so all the vectors are of length n_cells + 1. For a periodic simulation, the last node is the periodic image of the first one, and the values stored in it are mirrored from the first node.

The structure is grid-independent: it stores no reference to a grid or to a Poisson solver.

Fields

• n_nodes: number of nodes
• charge_density: the charge density $\rho$ in the nodes, C/m³ (before a call to normalize_charge_density! this holds the deposited charge in C and not a density)
• potential: the electrostatic potential $\phi$ in the nodes, V
• electric_field: the x-component of the electric field $E_x$ in the nodes, V/m
• net_charge_density: the mean charge density subtracted in the periodic case, 0.0 otherwise, C/m³

ElectrostaticFieldProps(n_nodes::Integer)

Construct an ElectrostaticFieldProps instance given the number of nodes (and not the number of cells, in contrast to PhysProps(n_cells, n_species)). For a 1-D uniform grid with n_cells cells, the number of nodes is n_cells + 1.

Positional arguments

• n_nodes: number of nodes

ElectrostaticFieldProps(grid::Grid1DUniform)

Construct an ElectrostaticFieldProps instance for a 1-D uniform grid, using n_nodes = grid.n_cells + 1 nodes.

Positional arguments

• grid: the Grid1DUniform grid

clear_charge_density!(field_props::ElectrostaticFieldProps)

Clear the charge density stored in an ElectrostaticFieldProps instance, leaving the potential and the electric field computed at the previous timestep untouched. This is to be called before the charge of the particles is deposited on the nodes at a new timestep.

Positional arguments

• field_props: the ElectrostaticFieldProps instance for which the charge density is cleared

compute_props!(particles, pia, species_data, phys_props)

Compute the physical properties of all species in all cells and store the result in a PhysProps instance. This function does not compute the total moments, even if phys_props.n_moments > 0.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored

compute_props_sorted!(particles, pia, species_data, phys_props, cell_chunk)

Compute the physical properties of all species in a subset of cells and store the result in a PhysProps instance, assuming the particles are sorted. This function does not compute the total moments, even if phys_props.n_moments > 0. Currently this does not compute the length of the particle array.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored
• cell_chunk: the list of cell indices or range of cell indices in which to compute the properties

compute_props_sorted!(particles, pia, species_data, phys_props)

Compute the physical properties of all species in all cells and store the result in a PhysProps instance, assuming the particles are sorted. This function does not compute the total moments, even if phys_props.n_moments > 0. Currently this does not compute the length of the particle array.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored

compute_props_sorted!(particles::Vector{ParticleVector{D}}, pia, species_data, phys_props, grid::G, cell_chunk) where {G<:AbstractGrid,D}

Compute the physical properties of all species in a subset of cells and store the result in a PhysProps instance, assuming the particles are sorted. This function does not compute the total moments, even if phys_props.n_moments > 0. If ndens_not_Np is true, the number density will be computed based on the volumes of the grid cells; otherwise, the number of physical particles in each cell will be computed. Currently this does not compute the length of the particle array.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored
• grid: the physical grid
• cell_chunk: the list of cell indices or range of cell indices in which to compute the properties

compute_props_sorted!(particles::Vector{ParticleVector{D}}, pia, species_data, phys_props, grid::G) where {G<:AbstractGrid,D}

Compute the physical properties of all species in all cells and store the result in a PhysProps instance, assuming the particles are sorted. This function does not compute the total moments, even if phys_props.n_moments > 0. If ndens_not_Np is true, the number density will be computed based on the volumes of the grid cells; otherwise, the number of physical particles in each cell will be computed. Currently this does not compute the length of the particle array.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored
• grid: the physical grid

compute_flux_props!(particles::Vector{ParticleVector{D}}, pia, species_data, phys_props::PhysProps, flux_props::FluxProps, grid::G) where {G<:AbstractGrid,D}

Compute the fluxes of all species in all cells and store the result in a FluxProps instance. This uses the pre-computed species-wise mean velocities from a PhysProps instance, which needs to be computed at the same timestep before calling this function. Particles of a cell can be distributed across group1 and group2 of indices pointed to by the ParticleIndexerArray instance, i.e. this function can be used to compute fluxes immediately after collisions (but before convection). In case particles are sorted, compute_flux_props_sorted! will be more efficient.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance computed at the same timestep
• flux_props: the FluxProps instance in which the computed fluxes are stored
• grid: the physical grid

compute_flux_props_sorted!(particles::Vector{ParticleVector{D}}, pia, species_data, phys_props, flux_props, grid::G, cell_chunk) where {G<:AbstractGrid, D}

Compute the flux densities of all species in a subset of cells and store the result in a FluxProps instance, assuming the particles are sorted. This uses the pre-computed species-wise mean velocities from a PhysProps instance, which needs to be computed at the same timestep before calling this function for the same subset of cells.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance computed at the same timestep for the same subset of cells
• flux_props: the FluxProps instance in which the computed fluxes are stored
• grid: the physical grid
• cell_chunk: the list of cell indices or range of cell indices in which to compute the properties

compute_flux_props_sorted!(particles::Vector{ParticleVector{D}}, pia, species_data, phys_props, flux_props, grid::G) where {G<:AbstractGrid,D}

Compute the flux densities of all species in all cells and store the result in a FluxProps instance, assuming the particles are sorted. This uses the pre-computed species-wise mean velocities from a PhysProps instance, which needs to be computed at the same timestep before calling this function.

Positional arguments

• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance computed at the same timestep
• flux_props: the FluxProps instance in which the computed fluxes are stored
• grid: the physical grid

avg_props!(phys_props_avg::PhysProps, phys_props::PhysProps, n_avg_timesteps)

Used to time-average computed physical properties, not including the total moments. For each instantaneous value of a property computed and stored in phys_props, it is divided by n_avg_timesteps and added to phys_props_avg.

Positional arguments

• phys_props_avg: the PhysProps instance used to store the time-averaged properties
• phys_props: the PhysProps instance holding the current values of the properties to be used for the averaging at the current timestep
• n_avg_timesteps: the number of timesteps over which the averaging is performed

Throws

ErrorException if phys_props_avg computes number density and phys_props computes number of physical particles, or vice versa.

avg_props!(flux_props_avg::FluxProps, flux_props::FluxProps, n_avg_timesteps)

Used to time-average computed flux densities. For each instantaneous value of a property computed and stored in flux_props, it is divided by n_avg_timesteps and added to flux_props_avg.

Positional arguments

• flux_props_avg: the FluxProps instance used to store the time-averaged properties
• flux_props: the FluxProps instance holding the current values of the properties to be used for the averaging at the current timestep
• n_avg_timesteps: the number of timesteps over which the averaging is performed

avg_props!(surf_props_avg::SurfProps, surf_props::SurfProps, n_avg_timesteps)

Used to time-average computed surface properties. For each instantaneous value of a property computed and stored in surf_props, it is divided by n_avg_timesteps and added to surf_props_avg.

Positional arguments

• surf_props_avg: the SurfProps instance used to store the time-averaged properties
• surf_props: the SurfProps instance holding the current values of the properties to be used for the averaging at the current timestep
• n_avg_timesteps: the number of timesteps over which the averaging is performed

clear_props!(phys_props::PhysProps)

Clear all data from PhysProps, for use when physical properties are averaged over timesteps and averaging over a new set of timesteps needs to be started.

Positional arguments

• phys_props: the PhysProps instance to be cleared

clear_props!(flux_props::FluxProps)

Clear all data from FluxProps, for use when flux densities are averaged over timesteps and averaging over a new set of timesteps needs to be started.

Positional arguments

• flux_props: the FluxProps instance to be cleared

clear_props!(field_props::ElectrostaticFieldProps)

Clear all data from an ElectrostaticFieldProps instance (the charge density, the potential, and the electric field). For clearing only the charge density in a time loop, use clear_charge_density!.

Positional arguments

• field_props: the ElectrostaticFieldProps instance to be cleared

clear_props!(surf_props::SurfProps)

Clear all data from a SurfProps instance, either at the start of a new convection step, or when physical properties are averaged over timesteps and averaging over a new set of timesteps needs to be started.

Positional arguments

• surf_props: the SurfProps instance to be cleared

compute_moment_scaling!(moment_scaling, moment_powers, species, species_data, Tref)

Fill the vector moment_scaling with the appropriate scaling based on the Maxwell-Boltzmann distribution for a specific species.

Positional arguments

• moment_scaling: vector to be filled with the moment scaling factors (length n_moments)
• moment_powers: vector of moment powers (length n_moments)
• species: the species for which the mixed moment is being computed
• species_data: the Vector of SpeciesData
• Tref: reference temperature used to compute the total moments of a Maxwellian distribution which are used to scale the total moments

compute_moments!(moment_values, moment_scaling, moment_powers, particles, pia, cell, species, species_data, phys_props)

Fill the vector moment_values with the computed scaled moment values for a single species in a single cell. The values of n (number density or number of physical particles) and v (mean velocity) are taken from phys_props. It's important that function compute_props! must be called before this function, so that the correct nandv`` values are available.

Positional arguments

• moment_values: vector to be filled with the computed moment values (length n_moments)
• moment_scaling: vector of precomputed moment scaling factors (length n_moments)
• moment_powers: vector of moment powers (length n_moments)
• particles: the Vector of ParticleVectors containing all the particles in a simulation
• pia: the ParticleIndexerArray instance
• cell: the cell in which the moment is being computed
• species: the species for which the mixed moment is being computed
• species_data: the Vector of SpeciesData
• phys_props: the PhysProps instance in which the computed physical properties are stored

mean_free_path(interaction, species, n, T)

Compute elastic VHS mean free path for single-species collisions. This also holds for a VSS interaction, as the VSS model uses the same total cross-section as the VHS model.

Positional arguments

• interaction: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species: the species for which to compute the mean free path
• n: the number density
• T: the temperature

Returns

Mean free path.

References

• Eqn. (4.65) in "Molecular Gas Dynamics and the Direct Simulation of Gas Flows" or (D.21) in "Nonequilibrium Gas Dynamics and Molecular Simulation".

mean_collision_frequency(interaction, species, species_data, n, T)

Compute elastic VHS mean collision frequency for single-species collisions. This also holds for a VSS interaction, as the VSS model uses the same total cross-section as the VHS model.

Positional arguments

• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species: the species for which to compute the mean free path
• species_data: the vector of Species instances of the species in the flow
• n: the number density
• T: the temperature

Returns

Mean collision frequency.

References

• Eqn. (4.64) in "Molecular Gas Dynamics and the Direct Simulation of Gas Flows"

debye_length(n, T)

Compute the Debye length $\lambda_D = \sqrt{\varepsilon_0 k_B T / (n q_e^2)}$ of a plasma with singly charged particles with density $n$ at a temperature $T$.

Positional arguments

• n: the number density of the charged particles
• T: the temperature of the charged particles

Returns

Debye length.

plasma_frequency(species, species_data, n)

Compute the plasma frequency $\omega_p = \sqrt{n q^2 / (\varepsilon_0 m)}$ of a species with charge $q$ and mass $m$ with a number density $n$.

Positional arguments

• species: the species for which to compute the plasma frequency
• species_data: the vector of Species instances of the species in the flow
• n: the number density of the species

Returns

Plasma frequency in rad/s.

CollisionData

Structure to store temporary collision data for a specific collision (relative velocity, collision energy, post-collision velocities, etc.)

Fields

• v_com: vector of center-of-mass velocity
• g: magnitude of relative velocity
• E_coll: relative translational energy of the colliding particles
• E_coll_eV: relative translational energy of the colliding particles in electron-volt
• g_vec: vector of pre-collisional relative velocity
• g_vec_new: vector of post-collisional relative velocity
• g_new_1: magnitude of post-collisional relative velocity of the first particle in the collision pair
• g_new_2: magnitude of post-collisional relative velocity of the second particle in the collision pair

CollisionData()

Create an empty CollisionData instance.

CollisionFactors

Structure to store NTC-related collision factors for collisions between particles of two species in a given cell.

Fields

• n1: the number of particles of the first species in the cell
• n2: the number of particles of the second species in the cell
• sigma_g_w_max: estimate of the $(\sigma g w)_{max}$ ($\sigma$ is the total collision cross-section, $g$ is the relative collision velocity, $w$ is the computational weight of the particles)
• n_coll: number of collisions to be tested
• n_coll_performed: number of collisions actually performed
• n_eq_w_coll_performed: number of collisions between particles with equal weights actually performed

CollisionFactors()

Create an empty CollisionFactors instance (all values set to 0).

CollisionFactorsSWPM

Structure to store SWPM-related collision factors for collisions between particles of two species in a given cell.

Fields

• n1: the number of particles of the first species in the cell
• n2: the number of particles of the second species in the cell
• sigma_g_max: estimate of the $(\sigma g)_{max}$ ($\sigma$ is the total collision cross-section, $g$ is the relative collision velocity
• n_coll: number of collisions to be tested
• n_coll_performed: number of collisions actually performed

CollisionFactorsSWPM()

Create an empty CollisionFactorsSWPM instance (all values set to 0).

CollisionDataFP

Structure to store temporary collision data for the particle Fokker-Planck approach.

Fields

• vel_ave: average velocity of the particles in the cell
• mean: mean value of the sampled velocities
• stddev: standard deviation of the sampled velocities
• xvel_rand: pre-allocated storage for sampled x-velocity components
• yvel_rand: pre-allocated storage for sampled y-velocity components
• zvel_rand: pre-allocated storage for sampled z-velocity components

CollisionDataFP()

Create an empty CollisionDataFP instance.

CollisionDataFP(n_particles_in_cell)

Create an empty CollisionDataFP instance, pre-allocating the arrays for sampled normal variables for n_particles_in_cell.

Positional arguments

• n_particles_in_cell: estimate of expected maximum number of particles in cell

create_collision_factors_array(n_species)

Create a 3-dimensional array of collision factors for all interaction pairs for a 0-D case (1 spatial cell), with shape (n_species,n_species,1).

Positional arguments

• n_species: number of species in the flow

Returns

3-dimensional array of CollisionFactors instances with shape (n_species,n_species,1).

create_collision_factors_array(n_species, n_cells)

Create a 3-dimensional array of collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells).

Positional arguments

• n_species: number of species in the flow
• n_cells: number of cells in the simulation

Returns

3-dimensional array of CollisionFactors instances with shape (n_species,n_species,n_cells).

create_collision_factors_array(pia)

Create a 3-dimensional array of collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells).

Positional arguments

• pia: the ParticleIndexerArray instance

Returns

3-dimensional array of CollisionFactors instances with shape (n_species,n_species,n_cells).

create_collision_factors_array(pia, interactions, species_data, T_list, Fnum::Real; mult_factor=1.0)

Create a 3-dimensional array of collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells). This will fill the array with the estimates $(\sigma g w)_{max}$ for all species in all cells, assuming a constant particle computational weight Fnum, a VHS cross-section, and that the temperature of each species is constant across all cells.

Positional arguments

• pia: the ParticleIndexerArray instance
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T_list: the list of temperatures of the species
• Fnum: the constant computational weight of the particles

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

3-dimensional array of CollisionFactors instances with shape (n_species,n_species,n_cells) filled with estimated values of $(\sigma g w)_{max}$.

create_collision_factors_array(pia, interactions, species_data, T::Real, Fnum::Real; mult_factor=1.0)

Create a 3-dimensional array of collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells). This will fill the array with the estimates $(\sigma g w)_{max}$ for all species in all cells, assuming a constant particle computational weight Fnum, a VHS cross-section, and that all species have a single temperature that is constant across all cells.

Positional arguments

• pia: the ParticleIndexerArray instance
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T: the temperatures of the flow
• Fnum: the constant computational weight of the particles

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

3-dimensional array of CollisionFactors instances with shape (n_species,n_species,n_cells) filled with estimated values of $(\sigma g w)_{max}$.

create_collision_factors_swpm_array(n_species)

Create a 3-dimensional array of SWPM collision factors for all interaction pairs for a 0-D case (1 spatial cell), with shape (n_species,n_species,1).

Positional arguments

• n_species: number of species in the flow

Returns

3-dimensional array of CollisionFactorsSWPM instances with shape (n_species,n_species,1).

create_collision_factors_swpm_array(n_species, n_cells)

Create a 3-dimensional array of SWPM collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells).

Positional arguments

• n_species: number of species in the flow
• n_cells: number of cells in the simulation

Returns

3-dimensional array of CollisionFactorsSWPM instances with shape (n_species,n_species,n_cells).

create_collision_factors_swpm_array(pia)

Create a 3-dimensional array of SWPM collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells).

Positional arguments

• pia: the ParticleIndexerArray instance

Returns

3-dimensional array of CollisionFactorsSWPM instances with shape (n_species,n_species,n_cells).

create_collision_factors_swpm_array(pia, interactions, species_data, T_list; mult_factor=1.0)

Create a 3-dimensional array of SWPM collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells). This will fill the array with the estimates $(\sigma g)_{max}$ for all species in all cells, assuming a VHS cross-section, and that the temperature of each species is constant across all cells.

Positional arguments

• pia: the ParticleIndexerArray instance
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T_list: the list of temperatures of the species

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

3-dimensional array of CollisionFactorsSWPM instances with shape (n_species,n_species,n_cells) filled with estimated values of $(\sigma g)_{max}$.

create_collision_factors_swpm_array(pia, interactions, species_data, T::Real; mult_factor=1.0)

Create a 3-dimensional array of SWPM collision factors for all interaction pairs for all cells in the simulation, with shape (n_species,n_species,n_cells). This will fill the array with the estimates $(\sigma g)_{max}$ for all species in all cells, assuming a VHS cross-section, and that all species have a single temperature that is constant across all cells.

Positional arguments

• pia: the ParticleIndexerArray instance
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T: the temperatures of the flow

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

3-dimensional array of CollisionFactorsSWPM instances with shape (n_species,n_species,n_cells) filled with estimated values of $(\sigma g)_{max}$.

create_computed_crosssections(electron_neutral_interactions)

Create a vector of ComputedCrossSections instances for the electron-neutral interactions.

Positional arguments

• electron_neutral_interactions: the ElectronNeutralInteractions instance for which the cross-sections will be computed

Returns

Vector of ComputedCrossSections of length electron_neutral_interactions.n_neutrals.

estimate_sigma_g_w_max(interaction, species1, species2, T1, T2, Fnum; mult_factor=1.0)

Estimate $(\sigma g w)_{max}$ for a two-species interaction, assuming a constant particle computational weight Fnum and a VHS cross-section. The relative velocity $g$ is estimated as $g = 0.5 (\sqrt{2T_1 k_B / m_1} + \sqrt{2T_2 k_B / m_2})$, where $T_1$ and $m_1$ are the temperature and mass of the first species (species1), and $T_2$ and $m_2$ are the temperature and mass of the second species (species2). This relative velocity estimate is then plugged into the VHS cross-section model to compute $\sigma$. The result is then multiplied by Fnum and an (optional) factor mult_factor. This is also the correct estimate for the VSS model, as it uses the same total cross-section as the VHS model.

Positional arguments

• interaction: the Interaction instance for the interacting species
• species1: the Species instance of the first interacting species
• species2: the Species instance of the second interacting species
• T1: the temperature of the first interacting species
• T2: the temperature of the second interacting species
• Fnum: the constant computational weight of the particles

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

The estimate of $(\sigma g w)_{max}$.

estimate_sigma_g_w_max(interaction, species, T, Fnum; mult_factor=1.0)

Estimate $(\sigma g w)_{max}$ for a single-species interaction, assuming a constant particle computational weight Fnum and a VHS cross-section. Uses the same methodology as the estimate for a two-species interaction.

Positional arguments

• interaction: the Interaction instance for the interacting species
• species: the Species instance of the interacting species
• T: the temperature of the interacting species
• Fnum: the constant computational weight of the particles

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

Returns

• the estimate of $(\sigma g w)_{max}$

estimate_sigma_g_w_max!(collision_factors, interactions, species_data, T_list, Fnum; mult_factor=1.0)

Estimate $(\sigma g w)_{max}$ for all species in all cells, assuming a constant particle computational weight Fnum, a VHS cross-section, and that each species' temperature is constant across all cells. Uses the same methodology as the estimate for a two-species interaction.

Positional arguments

• collision_factors: 3-dimensional array of CollisionFactors of shape (n_species, n_species, n_cells)
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T_list: the list of temperatures of the species
• Fnum: the constant computational weight of the particles

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

estimate_sigma_g_w_max!(collision_factors, interactions, species_data, T_list, Fnum; mult_factor=1.0)

Estimate $(\sigma g w)_{max}$ for all species in all cells, assuming species-specific computational weights Fnum, a VHS cross-section, and that each species' temperature is constant across all cells. Uses the same methodology as the estimate for a two-species interaction.

Positional arguments

• collision_factors: 3-dimensional array of CollisionFactors of shape (n_species, n_species, n_cells)
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T_list: the list of temperatures of the species
• Fnum: a list of the computational weights of the species

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

estimate_sigma_g_w_max_ntc_n_e!(rng, collision_factors, collision_data, interaction,
                                n_e_interactions, n_e_cs, particles_n, particles_e,
                                pia, cell, species_n, species_e, Δt, V; min_coll=5, n_loops=3)

Estimate $(\sigma g w)_{max}$ for an electron-neutral interaction by stochastically choosing particle pairs multiple times and computing $(\sigma g w)$ for each pair. The number of collisions is computed using the standard variable-weight NTC formula, the value of min_coll is added to this number, and particles are randomly sampled. The whole procedure is repeated n_loops times, so that an increased value $(\sigma g w)_{max}$ can have an impact on the computed number of pairs to select during the next loop iteration.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• n_e_interactions: the ElectronNeutralInteractions instance
• n_e_cs: a vector of ComputedCrossSections instances
• particles_n: ParticleVector of the particles of neutral species
• particles_e: ParticleVector of the particles of the electron species
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species_n: the index of the neutral species
• species_e: the index of the electron species
• Δt: timestep
• V: cell volume

Keyword arguments:

• min_coll: the minimum number of pairs to test
• n_loops: the number of loops to perform (in each loop the number of collisions is computed using the estimated value of $(\sigma g w)_{max}$)
• extend: enum of CSExtend type that sets how out-of-range energy values are treated when computing cross-sections

estimate_sigma_g_max!(collision_factors::Array{CollisionFactorsSWPM}, interactions, species_data, T_list; mult_factor=1.0)

Estimate $(\sigma g)_{max}$ for all species in all cells, assuming a constant particle computational weight Fnum, a VHS cross-section, and that each species' temperature is constant across all cells. Uses the same methodology as the estimate for a two-species interaction.

Positional arguments

• collision_factors: 3-dimensional array of CollisionFactorsSWPM of shape (n_species, n_species, n_cells)
• interactions: the 2-dimensional array of Interaction instances (of shape (n_species, n_species)) of all the pair-wise interactions
• species_data: the vector of Species instances of the species in the flow
• T_list: the list of temperatures of the species

Keyword arguments

• mult_factor: a factor by which to multiply the result (default value is 1.0)

ntc!(rng, collision_factors, collision_data, interaction, particles::ParticleVector{D}, pia,
     cell, species, Δt, V; dw_tol=1e-16) where D

Perform elastic collisions between particles of same species using the NTC algorithm. The elastic scattering model is taken from the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• Δt: timestep
• V: cell volume

Keyword arguments

• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

ntc!(rng, model::AbstractScatteringModel, collision_factors, collision_data, interaction,
     particles::ParticleVector{D}, pia, cell, species, Δt, V; dw_tol=1e-16) where D

Perform elastic collisions between particles of same species using the NTC algorithm and the elastic scattering model model, overriding the model stored in the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• model: the AbstractScatteringModel singleton for the scattering model
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• Δt: timestep
• V: cell volume

Keyword arguments

• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

ntc!(rng, collision_factors, collision_data, interaction,
     particles_1::ParticleVector{D}, particles_2::ParticleVector{D}, pia,
     cell, species1, species2, Δt, V; dw_tol=1e-16) where D

Perform elastic collisions between particles of different species using the NTC algorithm. The elastic scattering model is taken from the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles_1: ParticleVector of the particles of the first species being collided
• particles_2: ParticleVector of the particles of the second species being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species1: the index of the first species for which collisions are performed
• species2: the index of the second species for which collisions are performed
• Δt: timestep
• V: cell volume

Keyword arguments

• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

ntc!(rng, model::AbstractScatteringModel, collision_factors, collision_data, interaction,
     particles_1::ParticleVector{D}, particles_2::ParticleVector{D}, pia,
     cell, species1, species2, Δt, V; dw_tol=1e-16) where D

Perform elastic collisions between particles of different species using the NTC algorithm and the elastic scattering model model, overriding the model stored in the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• model: the AbstractScatteringModel singleton for the scattering model
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles_1: ParticleVector of the particles of the first species being collided
• particles_2: ParticleVector of the particles of the second species being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species1: the index of the first species for which collisions are performed
• species2: the index of the second species for which collisions are performed
• Δt: timestep
• V: cell volume

Keyword arguments

• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

ntc_equal_weight!(rng, collision_factors, collision_data, interaction, particles::ParticleVector{D}, pia,
                  cell, species, Δt, V) where D

Perform elastic collisions between particles of same species using the NTC algorithm. Particle weights are assumed to be equal, and no weight checks/splitting is performed. The elastic scattering model is taken from the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• Δt: timestep
• V: cell volume

References

• G.A. Bird, Molecular gas dynamics and the direct simulation of gas flows, Clarendon Press, Oxford, 1994.

ntc_equal_weight!(rng, model::AbstractScatteringModel, collision_factors, collision_data, interaction,
                  particles::ParticleVector{D}, pia, cell, species, Δt, V) where D

Perform elastic collisions between particles of same species using the NTC algorithm and the elastic scattering model model, overriding the model stored in the Interaction instance of the species pair. Particle weights are assumed to be equal, and no weight checks/splitting is performed.

Positional arguments

• rng: the random number generator
• model: the AbstractScatteringModel singleton for the scattering model
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• Δt: timestep
• V: cell volume

References

• G.A. Bird, Molecular gas dynamics and the direct simulation of gas flows, Clarendon Press, Oxford, 1994.

ntc_equal_weight!(rng, collision_factors, collision_data, interaction,
                  particles_1, particles_2, pia,
                  cell, species1, species2, Δt, V)

Perform elastic collisions between particles of different species using the NTC algorithm. Particle weights are assumed to be equal, and no weight checks/splitting is performed. The elastic scattering model is taken from the Interaction instance of the species pair.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles_1: ParticleVector of the particles of the first species being collided
• particles_2: ParticleVector of the particles of the second species being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species1: the index of the first species for which collisions are performed
• species2: the index of the second species for which collisions are performed
• Δt: timestep
• V: cell volume

References

• G.A. Bird, Molecular gas dynamics and the direct simulation of gas flows, Clarendon Press, Oxford, 1994.

ntc_equal_weight!(rng, model::AbstractScatteringModel, collision_factors, collision_data, interaction,
                  particles_1, particles_2, pia,
                  cell, species1, species2, Δt, V)

Perform elastic collisions between particles of different species using the NTC algorithm and the elastic scattering model model, overriding the model stored in the Interaction instance of the species pair. Particle weights are assumed to be equal, and no weight checks/splitting is performed.

Positional arguments

• rng: the random number generator
• model: the AbstractScatteringModel singleton for the scattering model
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles_1: ParticleVector of the particles of the first species being collided
• particles_2: ParticleVector of the particles of the second species being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species1: the index of the first species for which collisions are performed
• species2: the index of the second species for which collisions are performed
• Δt: timestep
• V: cell volume

References

• G.A. Bird, Molecular gas dynamics and the direct simulation of gas flows, Clarendon Press, Oxford, 1994.

ntc_n_e!(rng, collision_factors, collision_data, interaction,
         n_e_interactions, n_e_cs, particles_n, particles_e, particles_ion,
         pia, cell, species_n, species_e, species_ion, Δt, V; extend::CSExtend=CSExtendConstant, dw_tol=1e-16)

Perform electron-neutral elastic scattering and electron-impact ionization collisions.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• n_e_interactions: the ElectronNeutralInteractions instance
• n_e_cs: a vector of ComputedCrossSections instances
• particles_n: ParticleVector of the particles of neutral species
• particles_e: ParticleVector of the particles of the electron species
• particles_ion: ParticleVector of the particles of the ion species
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species_n: the index of the neutral species
• species_e: the index of the electron species
• species_ion: the index of the ion species
• Δt: timestep
• V: cell volume

Keyword arguments

• extend: enum of CSExtend type that sets how out-of-range energy values are treated when computing cross-sections
• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

ntc_n_e_es!(rng, collision_factors, collision_data, interaction,
         n_e_interactions, n_e_cs, particles_n, particles_e, particles_ion,
         pia, cell, species_n, species_e, species_ion, Δt, V; extend::CSExtend=CSExtendConstant, dw_tol=1e-16)

Perform electron-neutral elastic scattering and electron-impact ionization collisions using the event splitting method.

Positional arguments

• rng: the random number generator
• collision_factors: the CollisionFactors for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• n_e_interactions: the ElectronNeutralInteractions instance
• n_e_cs: a vector of ComputedCrossSections instances
• particles_n: ParticleVector of the particles of neutral species
• particles_e: ParticleVector of the particles of the electron species
• particles_ion: ParticleVector of the particles of the ion species
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species_n: the index of the neutral species
• species_e: the index of the electron species
• species_ion: the index of the ion species
• Δt: timestep
• V: cell volume

Keyword arguments

• extend: enum of CSExtend type that sets how out-of-range energy values are treated when computing cross-sections
• dw_tol: if weights of particles differ by less than this amount, an equal-weight collision is assumed

and no particle splitting is performed

References

• G. Oblapenko, D. Goldstein, P. Varghese, C. Moore, Hedging direct simulation Monte Carlo bets via event splitting. J. Comput. Phys, 2022.
• D.P. Schmidt, C.J. Rutland, A New Droplet Collision Algorithm. J. Comput. Phys, 2000.

swpm!(rng, collision_factors_swpm, collision_data, interaction, particles::ParticleVector{D}, pia,
      cell, species, G, Δt, V)

Perform elastic collisions between variable-weight particles of same species using the SWPM algorithm. The elastic scattering model is taken from the Interaction instance of the species pair. During a collision of particles with weights $w_i$, $w_j$, the weights are depleted by $\min(w_i, w_j) / (1+G)$, where $G \geq 0$ is a user-defined parameter.

Positional arguments

• rng: the random number generator
• collision_factors_swpm: the CollisionFactorsSWPM for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• G: non-negative value defining the weight transfer function
• Δt: timestep
• V: cell volume

References

• S. Rjasanow, W.Wagner, Stochastic numerics for the Boltzmann equation. Springer Berlin, Heidelberg, 2005.

swpm!(rng, model::AbstractScatteringModel, collision_factors_swpm, collision_data, interaction,
      particles::ParticleVector{D}, pia, cell, species, G, Δt, V)

Perform elastic collisions between variable-weight particles of same species using the SWPM algorithm and the elastic scattering model model, overriding the model stored in the Interaction instance of the species pair. During a collision of particles with weights $w_i$, $w_j$, the weights are depleted by $\min(w_i, w_j) / (1+G)$, where $G \geq 0$ is a user-defined parameter.

Positional arguments

• rng: the random number generator
• model: the AbstractScatteringModel singleton tag of the scattering model
• collision_factors_swpm: the CollisionFactorsSWPM for the species in question in the cell
• collision_data: CollisionData instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• G: non-negative value defining the weight transfer function
• Δt: timestep
• V: cell volume

References

• S. Rjasanow, W.Wagner, Stochastic numerics for the Boltzmann equation. Springer Berlin, Heidelberg, 2005.

fp_linear!(rng, collision_data_fp, interaction, particles::ParticleVector{D}, pia, cell, species, species_data, Δt, V) where D

Model single-species elastic collisions using a linear Fokker-Planck approximation.

Positional arguments

• rng: the random number generator
• collision_data_fp: CollisionDataFP instance used for storing collisional quantities
• interaction: 2-dimensional array of Interaction instances for all possible species pairs
• particles: ParticleVector of the particles being collided
• pia: the ParticleIndexerArray
• cell: the index of the cell in which collisions are performed
• species: the index of the species for which collisions are performed
• species_data: the vector of SpeciesData
• Δt: timestep
• V: cell volume

References

• M.H. Gorji, M. Torrilhon, P. Jenny, Fokker-Planck model for computational studies of monatomic rarefied gas flows. J. Fluid Mech., 2011.

ElectronNeutralInteractions

Structure to hold data on electron-neutral interactions. The neutral_indexer field is used to obtain the index of the species inside the structure, given an index of the neutral species in the full list of species in the simulation. For example, if we have a following list of species in the simulation: [e-, He, Ar+, He+, Ar], n_neutrals=2, the ElectronNeutralInteractions instance stores data for interactions of electrons with [He, Ar], and neutral_indexer[2] = 1, neutral_indexer[5] = 2.

Fields

• n_neutrals: number of neutral species for which the data has been loaded
• mass_ratios: ratio of electron mass to mass of neutrals; this uses the electron mass as given in the Species instance, so heavy electrons are possible
• neutral_indexer: array that maps indices of neutral species in the full list of species to local indices in the ElectronNeutralInteractions structure
• elastic: an array of ElasticScattering instances (of length n_neutrals) holding the cross-section data on elastic scattering for each neutral species
• ionization: an array of Ionization instances (of length n_neutrals) holding the cross-section data on electron-impact ionization for each neutral species
• excitation_sink: an array of ExcitationSink instances (of length n_neutrals) holding the cross-section data on electron-impact electronic excitation for each neutral species

ElectronNeutralInteractions(species_data, filename, databases, scattering_laws, energy_splits)

Load electron-neutral interaction data from an LXCAT format XML file for a set of given neutral species and instantiate a ElectronNeutralInteractions instance.

Positional arguments

• species_data: vector of Species data for the neutral species
• filename: path to XML file
• databases: dictionary of (species.name => database) pairs, specifying the name

of the cross-section database in the XML file to use for the species

• scattering_laws: vector of ScatteringLaw instances to use for each species
• energy_splits: vector of ElectronEnergySplit instances to use for each species

Returns

ElectronNeutralInteractions structure containing the electron-neutral interaction data.

Throws

DataMissingException if data not found or not all required data present.

ComputedCrossSections

Structure to hold data on computed cross-sections of electron-neutral interactions for a specific neutral species.

Fields

• n_excitations: number of electron-impact excitation reactions
• cs_total: the computed total cross-section (sum of cross-sections of all processes)
• cs_elastic: the computed elastic scattering cross-section
• cs_ionization: the computed electron-impact ionization cross-section
• cs_excitation: the computed electron-impact electronic excitation cross-section
• prob_vec: a vector of probabilities of the processes (of length 2+n_excitations). prob_vec[1] is the probability of elastic scattering, prob_vec[2] is the probability of electron-impact ionization, prob_vec[3:2+n_excitations] are the probabilities of the electron-impact electronic excitation processes
• cdf_prob_vec: a vector of cumulative probabilities of the processes (of length 2+n_excitations), used for sampling a specific process: cfd_prob_vec[1] = 0.0, cfd_prob_vec[n] = cfd_prob_vec[n-1] + prob_vec[n-1], n>1

ElectronEnergySplit ElectronEnergySplitEqual=1 ElectronEnergySplitZeroE=2

Enum for various splittings of electron energy in electron-impact ionization reactions. ElectronEnergySplitEqual corresponds to energy being shared equally amongst the electrons, ElectronEnergySplitZeroE corresponds to the one-takes-all sharing model.

ScatteringLaw ScatteringIsotropic=1 ScatteringOkhrimovskyy=2

Enum for various scattering laws in electron-neutral interactions. ScatteringIsotropic corresponds to isotropic scattering, ScatteringOkhrimovskyy to the scattering model of A. Okhrimovskyy et al., 2002.

CSExtend CSExtendConstant=1 CSExtendZero=2

Enum defining how to extend cross-sections in case energy is outside of tabulated range.

Possible values:

• CSExtendConstant: continue with closest value in array
• CSExtendZero: continue with zero

GridN2Merge{D}

Struct for merging particles with D-dimensional position vectors using a grid in velocity space. Particles in each cell are merged down. Particles outside of the grid are merged based on the octant they are in. So for an N:2 merge one would expect at most 2*(Nx*Ny*Nz+8) post-merge particles, where Nx, Ny, Nz are the number of grid cells in each velocity direction. The grid bounds in each direction can be computed using the mean thermal velocity of the particles and the mean streaming velocity: [v0-extent_multiplier*sqrt(2*k_B*T/m),v0+extent_multiplier*sqrt(2*k_B*T/m)]. Here T is the temperature of the species in question, m is the molecular mass, v0 is the mean velocity and extent_multiplier is a user-defined parameter (3.5 is a reasonable choice) defining the extent of the grid.

Fields

• Nx: number of grid cells in x velocity direction
• Ny: number of grid cells in y velocity direction
• Nz: number of grid cells in z velocity direction
• NyNz: product of Ny and Nz
• Ntotal: total number of grid cells (equal to Nx*Ny*Nz+8, as we account for the external octants)
• extent_multiplier: the vector of factors by which to multiply the thermal velocity to determine the grid bounds in each velocity direction
• extent_v_lower: the lower bounds of the velocity grid in each velocity direction
• extent_v_upper: the upper bounds of the velocity grid in each velocity direction
• extent_v_mid: (extent_v_lower + extent_v_upper)/2
• Δv: the grid cell size in each velocity direction
• Δv_inv: the inverse grid cell size in each velocity direction
• direction_vec: used to store randomly sampled direction signs
• direction_vecD: used to store randomly sampled direction signs of length D
• cells: vector of GridCell instances for each grid cell, as well as the external octants

GridN2Merge{D}(Nx::Int, Ny::Int, Nz::Int, extent_multiplier::T) where {D, T <: AbstractArray}

Create velocity grid-based merging for particles with D-dimensional position vectors.

Positional arguments

• Nx: number of cells in vx direction
• Ny: number of cells in vy direction
• Nz: number of cells in vz direction
• extent_multiplier: the vector of factors by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge{D}(N::Int, extent_multiplier::T) where {D, T <: AbstractArray}

Create velocity grid-based merging for particles with D-dimensional position vectors with an equal number of cells in each direction.

Positional arguments

• N: number of cells in each velocity direction
• extent_multiplier: the vector of factors by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge{D}(Nx::Int, Ny::Int, Nz::Int, extent_multiplier::Float64) where D

Create velocity grid-based merging for particles with D-dimensional position vectors with equal extent multipliers in each direction.

Positional arguments

• Nx: number of cells in vx direction
• Ny: number of cells in vy direction
• Nz: number of cells in vz direction
• extent_multiplier: the factor by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge{D}(Nx::Int, Ny::Int, Nz::Int,
            extent_multiplier_x::Float64,
            extent_multiplier_y::Float64,
            extent_multiplier_z::Float64) where D

Create velocity grid-based for merging particles with D-dimensional position vectors

Positional arguments

• Nx: number of cells in vx direction
• Ny: number of cells in vy direction
• Nz: number of cells in vz direction
• extent_multiplier_x: the factor by which to multiply the thermal velocity to determine the grid bounds

in the x-velocity direction

• extent_multiplier_y: the factor by which to multiply the thermal velocity to determine the grid bounds

in the y-velocity direction

• extent_multiplier_z: the factor by which to multiply the thermal velocity to determine the grid bounds

in the z-velocity direction

GridN2Merge{D}(N::Int, extent_multiplier::Float64) where D

Create velocity grid-based merging for particles with D-dimensional position vectors with an equal number of cells in each direction and equal extent multipliers in each direction.

Positional arguments

• N: number of cells in each velocity direction
• extent_multiplier: the factor by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge(N::Int, extent_multiplier::T) where T <: AbstractArray

Create velocity grid-based merging for particles with 3-dimensional position vectors with an equal number of cells in each direction.

Positional arguments

• N: number of cells in each velocity direction
• extent_multiplier: the vector of factors by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge(Nx::Int, Ny::Int, Nz::Int, extent_multiplier::Float64)

Create velocity grid-based merging for particles with 3-dimensional position vectors with equal extent multipliers in each direction.

Positional arguments

• Nx: number of cells in vx direction
• Ny: number of cells in vy direction
• Nz: number of cells in vz direction
• extent_multiplier: the factor by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

GridN2Merge(Nx::Int, Ny::Int, Nz::Int,
            extent_multiplier_x::Float64,
            extent_multiplier_y::Float64,
            extent_multiplier_z::Float64)

Create velocity grid-based merging for particles with 3-dimensional position vectors.

Positional arguments

• Nx: number of cells in vx direction
• Ny: number of cells in vy direction
• Nz: number of cells in vz direction
• extent_multiplier_x: the factor by which to multiply the thermal velocity to determine the grid bounds

in the x-velocity direction

• extent_multiplier_y: the factor by which to multiply the thermal velocity to determine the grid bounds

in the y-velocity direction

• extent_multiplier_z: the factor by which to multiply the thermal velocity to determine the grid bounds

in the z-velocity direction

merge_grid_based!(rng, merging_grid::GridN2Merge{D}, particles::ParticleVector{D}, pia, cell, species, species_data, phys_props::PhysProps) where D

Merge particles using a velocity grid-based merging approach. A Cartesian grid in velocity space is used to group particles together (particles outside of the grid are group by velocity octant), and in each cell/octant, particles are merged down to 2 particles. The extent of the grid is based on the temperature for the species in question in the physical grid cell being considered, as stored in the phys_props parameter.

Positional arguments:

• rng: the random number generator instance
• merging_grid: the grid merging (GridN2Merge) instance defining the velocity space grid
• particles: the ParticleVector instance of the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the cell index
• species: the species index
• species_data: the array of Species data
• phys_props: the PhysProps instance containing the computed temperature

References

• M. Vranic, T. Grismayer, J.L. Martins, R.A. Fonseca, L.O. Silva, Particle merging algorithm for PIC codes. Comput. Phys. Comm., 2015.
• G. Oblapenko, D. Goldstein, P. Varghese, C. Moore, A velocity space hybridization-based Boltzmann equation solver. J. Comput. Phys, 2020.

merge_grid_based!(rng, merging_grid::GridN2Merge{D}, particles::ParticleVector{D}, pia, cell, species, species_data, vx_extent, vy_extent, vz_extent)

Merge particles using a velocity grid-based merging approach. A Cartesian grid in velocity space is used to group particles together (particles outside of the grid are group by velocity octant), and in each cell/octant, particles are merged down to 2 particles. The extent of the grid is specified explicitly.

Positional arguments:

• rng: the random number generator instance
• merging_grid: the grid merging (GridN2Merge) instance defining the velocity space grid
• particles: the ParticleVector instance of the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the cell index
• species: the species index
• species_data: the array of Species data
• vx_extent: lower and upper bounds of the grid extent in the x velocity direction
• vy_extent: lower and upper bounds of the grid extent in the y velocity direction
• vz_extent: lower and upper bounds of the grid extent in the z velocity direction

References

• M. Vranic, T. Grismayer, J.L. Martins, R.A. Fonseca, L.O. Silva, Particle merging algorithm for PIC codes. Comput. Phys. Comm., 2015.
• G. Oblapenko, D. Goldstein, P. Varghese, C. Moore, A velocity space hybridization-based Boltzmann equation solver. J. Comput. Phys, 2020.

merge_grid_based!(rng, merging_grid::GridN2Merge{D}, particles::ParticleVector{D}, pia, cell, species, species_data, phys_props::PhysProps, grid::Grid1DUniform) where D

Merge particles using a velocity grid-based merging approach. A Cartesian grid in velocity space is used to group particles together (particles outside of the grid are group by velocity octant), and in each cell/octant, particles are merged down to 2 particles. The extent of the grid is based on the temperature for the species in question in the physical grid cell being considered, as stored in the phys_props parameter. If particle positions end up outside of the simulation domain, the particles are placed back into the domain.

Positional arguments:

• rng: the random number generator instance
• merging_grid: the grid merging (GridN2Merge) instance defining the velocity space grid
• particles: the ParticleVector instance of the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the cell index
• species: the species index
• species_data: the array of Species data
• phys_props: the PhysProps instance containing the computed temperature
• grid: the Grid1DUniform grid

References

• M. Vranic, T. Grismayer, J.L. Martins, R.A. Fonseca, L.O. Silva, Particle merging algorithm for PIC codes. Comput. Phys. Comm., 2015.
• G. Oblapenko, D. Goldstein, P. Varghese, C. Moore, A velocity space hybridization-based Boltzmann equation solver. J. Comput. Phys, 2020.

merge_grid_based!(rng, merging_grid::GridN2Merge{D}, particles::ParticleVector{D}, pia, cell, species, species_data, vx_extent, vy_extent, vz_extent, grid::Grid1DUniform) where D

Merge particles using a velocity grid-based merging approach. A Cartesian grid in velocity space is used to group particles together (particles outside of the grid are group by velocity octant), and in each cell/octant, particles are merged down to 2 particles. The extent of the grid is specified explicitly. If particle positions end up outside of the simulation domain, the particles are placed back into the domain.

Positional arguments:

• rng: the random number generator instance
• merging_grid: the grid merging (GridN2Merge) instance defining the velocity space grid
• particles: the ParticleVector instance of the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the cell index
• species: the species index
• species_data: the array of Species data
• vx_extent: lower and upper bounds of the grid extent in the x velocity direction
• vy_extent: lower and upper bounds of the grid extent in the y velocity direction
• vz_extent: lower and upper bounds of the grid extent in the z velocity direction
• grid: the Grid1DUniform grid

References

• M. Vranic, T. Grismayer, J.L. Martins, R.A. Fonseca, L.O. Silva, Particle merging algorithm for PIC codes. Comput. Phys. Comm., 2015.
• G. Oblapenko, D. Goldstein, P. Varghese, C. Moore, A velocity space hybridization-based Boltzmann equation solver. J. Comput. Phys, 2020.

GridN2Merge(N::Int, extent_multiplier::Float64)

Create velocity grid-based merging for particles with 3-dimensional position vectors with an equal number of cells in each direction and equal extenr multipliers in each direction.

Positional arguments

• N: number of cells in each velocity direction
• extent_multiplier: the factor by which to multiply the thermal velocity to determine the grid bounds

in each velocity direction

NNLSMerge{D}

Struct for keeping track of merging-related quantities for NNLS-based merging of particles with D-dimensional position vectors.

Fields

• v0: vector of the mean velocity of the particles
• x0: D-dimensional vector of the mean position of the particles
• vref: reference velocity magnitude for scaling
• inv_vref: inverse of reference velocity magnitude for scaling
• Ev: vector of standard deviation of velocities of particles
• Ex: vector of standard deviation of positions of particles
• w_total: total computational of the particles
• scalev: vector of values used to scale the velocity of particles
• scalex: vector of values used to scale the positions of particles
• n_total_conserved: total number of moments conserved
• n_moments_vel: number of velocity moments to preserve
• rhs_vector: vector of computed moments
• row_scale: scratch vector of the row-wise scaling factors applied to the LHS matrix
• mim: vector of 3-tuples of multi-indices for the velocity moments to preserve
• tot_order: vector of total orders of the velocity moments to preserve
• vel_powers: scratch table of powers of the centered velocity components of a particle
• n_moments_pos: number of spatial moments to preserve
• mim_pos: vector of 3-tuples of multi-indices for the spatial moments to preserve
• tot_order_pos: vector of total orders of the spatial moments to preserve
• pos_powers: scratch table of powers of the centered position components of a particle
• pos_i_x: index of the spatial moment corresponding to preservation of the center of mass in the x direction
• pos_i_y: index of the spatial moment corresponding to preservation of the center of mass in the y direction
• pos_i_z: index of the spatial moment corresponding to preservation of the center of mass in the z direction
• lhs_matrix_ncols_start: the number of columns in the first pre-allocated matrix
• lhs_matrix_ncols_end: the number of columns in the last pre-allocated matrix
• column_norms: vector of vectors of the column-wise inverse norms of the LHS matrices
• vel_pos_matrices: vector of matrices of size (3+D)xNp that store the velocities and positions of the particles
• column_norms_scratch: column-wise inverse norms used when the number of columns is outside the pre-allocated range
• vel_pos_matrix_scratch: velocities and positions used when the number of columns is outside the pre-allocated range
• work: Vector of NNLSWorkspace instances

NNLSMerge{D}(multi_index_moments, init_np; rate_preserving=false, multi_index_moments_pos=[], matrix_ncol_nprealloc=0) where D

Create NNLS-based merging. Mass, momentum, directional energy are always conserved: if not in the list multi_index_moments of moments to preserve, the corresponding moment multi-indices will be added automatically. These indices are (0,0,0) for mass, (1,0,0), (0,1,0), (0,0,1) for momentum, and (2,0,0), (0,2,0), (0,0,2) for directional energies. Spatial moments can be preserved by setting multi_index_moments_pos. If multi_index_moments_pos is non-empty, it is currently left to the user to include any relevant 1st order moments (corresponding center of mass conservation) i.e. (1, 0, 0), (0, 1, 0), (0, 0, 1); otherwise the corresponding spatial moments including higher-order ones will not be conserved. By default, a single NNLSWorkspace is pre-allocated for the system. The number of columns in the LHS matrix is the number of particles to be merged (+ any fictitious particles); so it cannot be fixed in advance. By setting matrix_ncol_nprealloc to a value larger than 0, one can pre-allocate a range of workspaces with a fixed number of columns spanning [init_np, init_np+matrix_ncol_nprealloc]. In this case, matrix_ncol_nprealloc+2 NNLSWorkspace instances are pre-allocated, with the last one being used in case the number of columns is not in the range (its LHS matrix is then re-allocated on the fly, but only when the number of columns differs from the previous call). Vectors column_norms and vel_pos_matrices are also then pre-allocated, holding respectively the inverses of the column-wise norms of the LHS matrices used for their scaling, and the particle velocities and positions. Outside the pre-allocated range the column_norms_scratch / vel_pos_matrix_scratch buffers are used instead; these grow as needed and are re-used across calls.

Positional arguments

• multi_index_moments: vector of mixed moments to preserve of the form [(i1, j1, k1), (i2, j2, k2), ...]`
• init_np: assumption on pre-merge number of particles to pre-allocate memory for

Keyword arguments:

• rate_preserving: used for rate-preserving merging of electrons, preserves approximate elastic collision and ionization rates
• multi_index_moments_pos: list of spatial moments to preserve
• matrix_ncol_nprealloc: number of NNLS workspaces with a fixed number of columns to pre-allocate

Throws

• ArgumentError: if any of moment powers is negative

NNLSMerge(multi_index_moments, init_np; rate_preserving=false, multi_index_moments_pos=[], matrix_ncol_nprealloc=0)

Create NNLS-based merging for particles with a 3-dimensional position vector. Mass, momentum, directional energy are always conserved: if not in the list multi_index_moments of moments to preserve, the corresponding moment multi-indices will be added automatically. These indices are (0,0,0) for mass, (1,0,0), (0,1,0), (0,0,1) for momentum, and (2,0,0), (0,2,0), (0,0,2) for directional energies. Spatial moments can be preserved by setting multi_index_moments_pos. If multi_index_moments_pos is non-empty, it is currently left to the user to include any relevant 1st order moments (corresponding center of mass conservation) i.e. (1, 0, 0), (0, 1, 0), (0, 0, 1); otherwise the corresponding spatial moments including higher-order ones will not be conserved. By default, a single NNLSWorkspace is pre-allocated for the system. The number of columns in the LHS matrix is the number of particles to be merged (+ any fictitious particles); so it cannot be fixed in advance. By setting matrix_ncol_nprealloc to a value larger than 0, one can pre-allocate a range of workspaces with a fixed number of columns spanning [init_np, init_np+matrix_ncol_nprealloc]. In this case, matrix_ncol_nprealloc+2 NNLSWorkspace instances are pre-allocated, with the last one being used in case the number of columns is not in the range (its LHS matrix is then re-allocated on the fly, but only when the number of columns differs from the previous call). Vectors column_norms and vel_pos_matrices are also then pre-allocated, holding respectively the inverses of the column-wise norms of the LHS matrices used for their scaling, and the particle velocities and positions. Outside the pre-allocated range the column_norms_scratch / vel_pos_matrix_scratch buffers are used instead; these grow as needed and are re-used across calls.

Positional arguments

• multi_index_moments: vector of mixed moments to preserve of the form [(i1, j1, k1), (i2, j2, k2), ...]`
• init_np: assumption on pre-merge number of particles to pre-allocate memory for

Keyword arguments:

• rate_preserving: used for rate-preserving merging of electrons, preserves approximate elastic collision and ionization rates
• multi_index_moments_pos: list of spatial moments to preserve
• matrix_ncol_nprealloc: number of NNLS workspaces with a fixed number of columns to pre-allocate

Throws

• ArgumentError: if any of moment powers is negative

compute_multi_index_moments(n)

Compute all mixed moment multi-indices of total order up to n, i.e. all 3-tuples (i,j,k)such thati+j+k <= n`.

Positional arguments

• n: maximum total order

Returns

Vector of 3-tuples of moment multi-indices.

merge_nnls_based!(rng, nnls_merging::NNLSMerge{D}, particles::ParticleVector{D}, pia, cell, species;
                  vref=1.0, scaling=:variance,
                  max_err=1e-11, iteration_mult=2, w_threshold=0.0) where D

Perform NNLS-based merging. The NNLS system is scaled to improve numerical stability, the scaling algorithm is set by the scaling parameter. Even if scaling is done using the computed variances, vref might be used in case those variances are small.

Positional arguments

• rng: the random number generator instance
• nnls_merging: the NNLSMerge instance
• particles: the ParticleVector instance containing the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the index of the grid cell in which particles are being merged
• species: the index of the species being merged

Keyword arguments

• vref: the reference velocity used to scale the velocities
• scaling: how to scale entries in the LHS and RHS of the NNLS system - either based on the reference velocity vref (scaling=:vref) or on the computed variances in each direction (scaling=:variance)
• max_err: maximum allowed value of the residual of the NNLS system
• iteration_mult: the number by which the number of columns of the NNLS system matrix is multiplied, this gives the maximum number of iterations of the NNLS algorithm
• w_threshold: any particles with a relative weight smaller than this value will be discarded (and the weight of the remaining particles re-scaled)

Returns

If the residual exceeds max_err or the number of non-zero (or smaller than w_threshold) elements in the solution vector is equal to the original number of particles, -1 is returned to signify a failure of the merging algorithm.

References

• G. Oblapenko, M. Torrilhon, Moment-preserving particle merging via non-negative least squares. arXiv preprint, 2026.

merge_nnls_based_rate_preserving!(rng, nnls_merging::NNLSMerge{D},
                                  interaction, electron_neutral_interactions, computed_cs,
                                  particles::ParticleVector{D}, pia, cell, species, neutral_species_index,
                                  ref_cs_elastic, ref_cs_ion; scaling=:variance,
                                  vref=1.0,  max_err=1e-11,
                                  iteration_mult=2,
                                  extend::CSExtend=CSExtendConstant) where D

Perform NNLS-based merging of electrons that conserves approximate elastic scattering and electron-impact ionization rates. The NNLS system is scaled to improve numerical stability, the scaling algorithm is set by the scaling parameter. Even if scaling is done using the computed variances, vref might be used in case those variances are small.

The reference velocity is also used in conjunction with the reference cross-sections to scale the parts of the NNLS matrix and RHS corresponding to conservation of electron-neutral collision rates. The reference rate is computed as $\sigma_{r,ref} v_{ref}$, where $\sigma_{r,ref}$ is the reference process cross-section. In case scaling==:variance, the reference velocity for the computation of reference rates is computed as $v_{ref} = \sqrt{E_x^2 + E_y^2 + E_z^2}$, where $E_x$, $E_y$, $E_z$ are the variances of the velocity in the corresponding directions.

Positional arguments

• rng: the random number generator instance
• nnls_merging: the NNLSMerge instance
• interaction: the Interaction instance describing the electron-neutral interaction being considered
• electron_neutral_interactions: the ElectronNeutralInteractions instance storing the tabulated cross-section data used to compute the rates
• computed_cs: the vector of ComputedCrossSection instances in which the computed values will be stored
• particles: the ParticleVector instance containing the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the index of the grid cell in which particles are being merged
• species: the index of the species being merged
• neutral_species_index: the index of the neutral species which is the collision partner in the electron-neutral collisions for which approximate rates are being preserved.
• ref_cs_elastic: the reference elastic scattering cross-section used to scale the rates
• ref_cs_ion: the reference electron-impact ionization cross-section used to scale the rates

Keyword arguments

• vref: the reference velocity used to scale the velocities in the case of scaling=:vref
• scaling: how to scale entries in the LHS and RHS of the NNLS system - either based on the reference velocity vref (scaling=:vref) or on the computed variances in each direction (scaling=:variance)
• max_err: maximum allowed value of the residual of the NNLS system
• iteration_mult: the number by which the number of columns of the NNLS system matrix is multiplied, this gives the maximum number of iterations of the NNLS algorithm
• w_threshold: any particles with a relative weight smaller than this value will be discarded (and the weight of the remaining particles re-scaled)
• extend: enum of CSExtend type that sets how out-of-range energy values are treated when computing cross-sections

Returns

If the residual exceeds max_err or the number of non-zero (or smaller than w_threshold) elements in the solution vector is equal to the original number of particles, -1 is returned to signify a failure of the merging algorithm.

References

• G. Oblapenko, M. Torrilhon, Moment-preserving particle merging via non-negative least squares. arXiv preprint, 2026.

merge_nnls_based_rate_preserving!(rng, nnls_merging::NNLSMerge{D},
                                  interaction, electron_neutral_interactions, computed_cs,
                                  particles::ParticleVector{D}, particles_neutral::ParticleVector{D}, pia, cell, species, neutral_species_index,
                                  ref_cs_elastic, ref_cs_ion; vref=1.0, scaling=:variance,
                                  max_err=1e-11,
                                  iteration_mult=2, w_threshold=0.0,
                                  extend::CSExtend=CSExtendConstant) where D

Perform NNLS-based merging of electrons that conserves exact elastic scattering and electron-impact ionization rates for one specific neutral species. The NNLS system is scaled to improve numerical stability, the scaling algorithm is set by the scaling parameter. Even if scaling is done using the computed variances, vref might be used in case those variances are small.

The reference velocity is also used in conjunction with the reference cross-sections to scale the parts of the NNLS matrix and RHS corresponding to conservation of electron-neutral collision rates. The reference rate is computed as $\sigma_{r,ref} v_{ref}$, where $\sigma_{r,ref}$ is the reference process cross-section. In case scaling==:variance, the reference velocity for the computation of reference rates is computed as $v_{ref} = \sqrt{E_x^2 + E_y^2 + E_z^2}$, where $E_x$, $E_y$, $E_z$ are the variances of the velocity in the corresponding directions.

Positional arguments

• rng: the random number generator instance
• nnls_merging: the NNLSMerge instance
• interaction: the Interaction instance describing the electron-neutral interaction being considered
• electron_neutral_interactions: the ElectronNeutralInteractions instance storing the tabulated cross-section data used to compute the rates
• computed_cs: the vector of ComputedCrossSection instances in which the computed values will be stored
• particles: the ParticleVector instance containing the particles to be merged
• particles_neutral: the ParticleVector instance containing the neutral collision partner particles (they are not affected by the merge)
• pia: the ParticleIndexerArray instance
• cell: the index of the grid cell in which particles are being merged
• species: the index of the species being merged
• neutral_species_index: the index of the neutral species which is the collision partner in the electron-neutral collisions for which approximate rates are being preserved.
• ref_cs_elastic: the reference elastic scattering cross-section used to scale the rates
• ref_cs_ion: the reference electron-impact ionization cross-section used to scale the rates

Keyword arguments

• vref: the reference velocity used to scale the velocities in the case of scaling=:vref
• scaling: how to scale entries in the LHS and RHS of the NNLS system - either based on the reference velocity vref (scaling=:vref) or on the computed variances in each direction (scaling=:variance)
• max_err: maximum allowed value of the residual of the NNLS system
• iteration_mult: the number by which the number of columns of the NNLS system matrix is multiplied, this gives the maximum number of iterations of the NNLS algorithm
• w_threshold: any particles with a relative weight smaller than this value will be discarded (and the weight of the remaining particles re-scaled)
• extend: enum of CSExtend type that sets how out-of-range energy values are treated when computing cross-sections

Returns

If the residual exceeds max_err or the number of non-zero (or smaller than w_threshold) elements in the solution vector is equal to the original number of particles, -1 is returned to signify a failure of the merging algorithm.

References

• G. Oblapenko, M. Torrilhon, Moment-preserving particle merging via non-negative least squares. arXiv preprint, 2026.

OctreeBinSplit OctreeBinMidSplit=1 OctreeBinMeanSplit=2 OctreeBinMedianSplit=3

Enum defining how the velocity along which the bin is split is chosen.

Possible values:

• OctreeBinMidSplit: the bin is split along the middle velocity
• OctreeBinMeanSplit: the bin is split along the mean velocity of the particles in the bin
• OctreeBinMedianSplit: the bin is split along the median velocity of the particles in the bin

OctreeInitBin OctreeInitBinMinMaxVel=1 OctreeInitBinMinMaxVelSym=2 OctreeInitBinC=3

Enum defining how the bounds of the initial bin are computed.

Possible values:

• OctreeInitBinMinMaxVel: the minimum and maximum velocities of the particles being merged are used to compute the bounds
• OctreeInitBinMinMaxVelSym: the minimum and maximum velocities of the particles being merged are used to compute the bounds, but the bounds are then symmetrized in each velocity direction: [-max(abs(min_v), abs(max_v)), max(abs(min_v), abs(max_v))]
• OctreeInitBinC: the initial bounds are set to [-c, c]in each direction, wherec`` speed of light

OctreeBinBounds OctreeBinBoundsInherit=1 OctreeBinBoundsRecompute=2

Enum defining how the bounds of a split sub-octant bin are computed.

Possible values:

• OctreeBinBoundsInherit: the splitting velocity and the appropriate bounds of the parent bin are inherited
• OctreeBinBoundsRecompute: the bounds are recomputed based on the particles in the bin

OctreeMerge{D,M}

Struct for N:M Octree merging for particles with D-dimensional position vectors.

Fields

• max_Nbins: maximum possible number of bins
• Nbins: number of bins currently used
• bins: Vector of OctreeCell instances used to compute the properties required for bin refinement
• full_bins: Vector of OctreeFullCell{D,M} instances used to compute the post-merge particles in each bin
• n_particles: total number of particles being merged
• bin_start: denotes start of indices of particles in bin i in the particle_indexes_sorted array
• bin_end: denotes end of indices of particles in bin i in the particle_indexes_sorted array
• particle_indexes_sorted: Vector of particle indices of the particles being merged
• particle_octants: Vector of particle octants for each particle used during radix sort
• particles_sort_output: Vector of integer indices used to store particle indices during radix sort
• particle_in_bin_counter: MVector of size 8, stores the number of particles in each bin
• nonempty_counter: MVector of size 8, stores the number of particles in each non-empty bin (the octants to which these bins correspond to are in nonempty_bins)
• nonempty_bins: MVector of size 8, a sequential list of non-empty octants
• ndens_counter: MVector of size 8, used in bin splitting, stores number density in each (non-empty) bin
• bin_bounds_compute: enum of OctreeBinBounds type defining whether bin bounds are fully defined by the parent bin and splitting velocity (vel_middle), or whether they are recomputed for each new sub-octant bin
• split: enum of OctreeBinSplit type defining how bins are split
• vel_middle: used to store the velocity along which a bin is split into octants
• v_min_parent: used in bin splitting to store the vector of the per-component lower bounds of the velocities in the cell
• v_max_parent: used in bin splitting to store the vector of the per-component upper bounds of the velocities in the cell
• direction_vec3: used to store randomly sampled direction signs, of length 3
• direction_vecD: used to store randomly sampled direction signs, of length D
• init_bin_bounds: enum of OctreeInitBin type defining how the bounds of the top-level bin are set
• max_depth: maximum allowed depth of a bin
• total_post_merge_np: used to keep track of number of post-merge particles
• v_mean: used to store mean velocity for conservative N:1 merging
• x_mean: used to store mean position for conservative N:1 merging
• v_var_before: used to store pre-merge variance of velocity for conservative N:1 merging
• x_var_before: used to store pre-merge variance of position for conservative N:1 merging
• v_mean_post: used to store post-merge mean velocity for conservative N:1 merging
• x_mean_post: used to store post-merge mean position for conservative N:1 merging
• v_var_post: used to store post-merge variance of velocity for conservative N:1 merging
• x_var_post: used to store post-merge variance of position for conservative N:1 merging

OctreeMerge{D,M}(split::OctreeBinSplit; init_bin_bounds=OctreeInitBinMinMaxVel, bin_bounds_compute=OctreeBinBoundsInherit,
          max_Nbins=4096, max_depth=10)

Create an Octree N:M merging instance for particles with D-dimensional position vectors.

Positional arguments:

• split: a enum of OctreeBinSplit type which tells how to split a bin into sub-bins

Keyword arguments:

• init_bin_bounds: a enum of OctreeInitBin type which defines how the bounds of the top-level bin are set
• bin_bounds_compute: a enum of OctreeBinBounds type which defines whether the bounds of sub-bins are recomputed based on the minimum/maximum velocities of the particles in those sub-bins, or the bounds are inherited from the bin that was split
• max_Nbins: maximum number of bins allowed (this only counts leaf-level bins)
• max_depth: maximum depth of a sub-bin starting from the top-level bin containing all particles (which has a depth of 0)

Returns: OctreeMerge instance with everything set to 0.

OctreeMerge(split::OctreeBinSplit; init_bin_bounds=OctreeInitBinMinMaxVel, bin_bounds_compute=OctreeBinBoundsInherit,
          max_Nbins=4096, max_depth=10)

Create an Octree N:2 merging instance for particles with 3-dimensional position vectors.

Positional arguments:

• split: a enum of OctreeBinSplit type which tells how to split a bin into sub-bins

Keyword arguments:

• init_bin_bounds: a enum of OctreeInitBin type which defines how the bounds of the top-level bin are set
• bin_bounds_compute: a enum of OctreeBinBounds type which defines whether the bounds of sub-bins are recomputed based on the minimum/maximum velocities of the particles in those sub-bins, or the bounds are inherited from the bin that was split
• max_Nbins: maximum number of bins allowed (this only counts leaf-level bins)
• max_depth: maximum depth of a sub-bin starting from the top-level bin containing all particles (which has a depth of 0)

Returns: OctreeMerge instance with everything set to 0.

merge_octree!(rng, octree::OctreeMerge{D,M}, particles::ParticleVector{D}, pia, cell, species, target_np) where {D,M}

Perform octree N:M merging without checking whether particle positions end up outside of the simulation domain.

Positional arguments

• rng: the random number generator instance
• octree: the OctreeMerge instance
• particles: the ParticleVector instance containing the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the index of the grid cell in which particles are being merged
• species: the index of the species being merged
• target_np: the target post-merge number of particles; the post-merge number of particles will not exceed this value but may be not exactly equal to it

References

• R.S. Martin, J.-L. Cambier, Octree particle management for DSMC and PIC simulations. J. Comput. Phys., 2016.

merge_octree!(rng, octree::OctreeMerge{D,M}, particles::ParticleVector{D}, pia, cell, species, target_np, grid::Grid1DUniform) where {D,M}

Perform octree N:M merging, checking whether particle positions end up outside of the simulation domain, and placing them back into the domain if needed.

Positional arguments

• rng: the random number generator instance

• octree: the OctreeMerge instance

• particles: the ParticleVector instance containing the particles to be merged

• pia: the ParticleIndexerArray instance

• cell: the index of the grid cell in which particles are being merged

• species: the index of the species being merged

• target_np: the target post-merge number of particles; the post-merge number of particles will not exceed this value but may be not exactly equal to it

• grid: the Grid1DUniform grid

• R.S. Martin, J.-L. Cambier, Octree particle management for DSMC and PIC simulations. J. Comput. Phys., 2016.

merge_roulette!(rng, particles::ParticleVector{D}, pia, cell, species, target_np) where D

Perform roulette merging - delete random particles until target number of particles is reached, and re-weight remaining particles to conserve number density.

Positional arguments

• rng: the random number generator instance
• particles: the ParticleVector instance containing the particles to be merged
• pia: the ParticleIndexerArray instance
• cell: the index of the grid cell in which particles are being merged
• species: the index of the species being merged
• target_np: the target post-merge number of particles

Keyword arguments

• conservative: if true, the post-merge particles' velocities will be corrected to ensure conservation of momentum and energy

References

• J. Watrous, D.B. Seidel, C.H. Moore, W. McDoniel, Improvements to Particle Merge Algorithms for Sandia National Laboratories Plasma Physics Modeling Code, EMPIRE. Presentation, 2023.

AbstractGrid

Abstract grid type

Grid1DUniform

1-D Uniform grid for a domain $[0,L]$

Fields

• L: length of the domain
• n_cells: number of cells
• Δx: cell size
• inv_Δx: inverse of cell size
• cells: Vector of Cell1D elements
• min_x: minimum allowed x coordinate for particles (slightly larger than $0$)
• max_x: maximum allowed x coordinate for particles (slightly smaller than $L$)
• surface_normals: vector of surface normals (1st element corresponds to the left wall, 2nd to the right wall)

Grid1DUniform(L, nx; wall_offset=1e-12)

Create 1-D uniform grid for a domain $[0, L]$ with nx cells

Positional arguments

• L: length of the domain
• nx: number of cells

Keyword arguments

• wall_offset: specifies a small relative offset from the walls so that particles never end up with a coordinate of exactly $0$ or $L$, otherwise some sorting routines may produce cell indices outside of the 1:nx range. The offset is computed as Δx * wall_offset, where Δx is the cell size.

GridSortInPlace

Struct for in-place sorting of particles.

The occ_lo/occ_hi fields record the first and last cell holding particles, as recorded by the most recent sorting call performed with this instance, to be used for multi-threaded simulations with particle exchange between threads. They describe that call only, and are meant to be consumed right after it, by update_occupancy_bounds!; a GridSortInPlace shared between species or ParticleIndexerArray instances therefore holds the bounds of whichever was sorted last. They are 1/n_cells before the first sort.

Fields

• cell_counts: scratch vector used to count the number of particles in each cell and to compute the resulting offsets
• sorted_indices: vector to store sorted particle indices
• occ_lo: first cell holding particles during the last sort
• occ_hi: last cell holding particles during the last sort

GridSortInPlace(n_cells::Integer, n_particles::Integer)

Create a GridSortInPlace instance given a number of grid cells and number of particles.

Positional arguments

• n_cells: the number of grid cells
• n_particles: the (expected) number of particles in the simulation (to pre-allocate the sorted_indices vector) - it is recommended to set this to the maximum expected number of particles in the simulation to avoid resizing of arrays during a simulation

GridSortInPlace(grid::G, n_particles::Integer) where {G<:AbstractGrid}

Create a GridSortInPlace instance given a grid and number of particles.

Positional arguments

• grid: the grid on which to sort the particles
• n_particles: the (expected) number of particles in the simulation (to pre-allocate the sorted_indices vector) - it is recommended to set this to the maximum expected number of particles in the simulation to avoid resizing of arrays during a simulation

sort_particles!(gridsort::GridSortInPlace, grid, particles::ParticleVector{D}, pia, species) where D

Sort particles on a grid using an in-place sorting algorithm. The pia instance is allowed to have non-contiguous indices (arising for example from merging). This function assumes that at the start of the sorting, it is not known in which cell each particle is located, and therefore the cell for each particle has to be determined (by calling get_cell).

The first and last cell holding particles are recorded in the GridSortInPlace instance, for use in update_occupancy_bounds!.

Positional arguments

• gridsort: the GridSortInPlace structure
• grid: the grid (should have an n_cells field, and a get_cell function has to be defined for the grid type)
• particles: the ParticleVector of particles to be sorted
• pia: the ParticleIndexerArray instance
• species: the index of the species being sorted

sort_particles!(gridsort::GridSortInPlace, particles, pia, species)

Sort particles on a grid using an in-place sorting algorithm. The pia instance is allowed to have non-contiguous indices (arising for example from merging). This function assumes that at the start of the sorting, it is known in which cell each particle is located.

The first and last cell holding particles are recorded in the GridSortInPlace instance, for use in update_occupancy_bounds!.

Positional arguments

• gridsort: the GridSortInPlace structure
• particles: the ParticleVector of particles to be sorted
• pia: the ParticleIndexerArray instance
• species: the index of the species being sorted

convect_particles!(rng, grid::Grid1DUniform, bc_list, particles::ParticleVector{D}, pia, species, species_data, Δt) where D

Convect particles on a 1-D uniform grid.

Positional arguments

• rng: the random number generator
• grid: the grid on which the convection is performed
• bc_list: the Tuple of boundary conditions (left and right wall)
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• species_data: the vector of Species data
• Δt: the convection timestep

convect_particles!(rng, grid::Grid1DUniform, bc_list, surf_props::SurfProps, particles::ParticleVector{D}, pia, species, species_data, Δt) where D

Convect particles on a 1-D uniform grid, computing surface properties if particles hit a surface.

Positional arguments

• rng: the random number generator
• grid: the grid on which the convection is performed
• bc_list: the Tuple of boundary conditions (left and right wall)
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• species_data: the vector of Species data
• surf_props: the SurfProps struct where the computed surface properties will be stored
• Δt: the convection timestep

convect_particles_periodic!(grid::Grid1DUniform, particles::ParticleVector{D}, pia, species, Δt) where D

Convect particles on a 1-D uniform grid assuming a periodic grid.

Positional arguments

• grid: the grid on which the convection is performed
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• Δt: the convection timestep

convect_particles_and_compute_cell!(rng, grid::Grid1DUniform, bc_list, particles::ParticleVector{D}, pia, species, species_data, Δt) where D

Convect particles on a 1-D uniform grid and write post-convection cell index to particles.cell.

Positional arguments

• rng: the random number generator
• grid: the grid on which the convection is performed
• bc_list: the Tuple of boundary conditions (left and right wall)
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• species_data: the vector of Species data
• Δt: the convection timestep

convect_particles_and_compute_cell!(rng, grid::Grid1DUniform, bc_list, surf_props::SurfProps, particles::ParticleVector{D}, pia, species, species_data, Δt) where D

Convect particles on a 1-D uniform grid and write post-convection cell index to particles.cell, computing surface properties if particles hit a surface.

Positional arguments

• rng: the random number generator
• grid: the grid on which the convection is performed
• bc_list: the Tuple of boundary conditions (left and right wall)
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• species_data: the vector of Species data
• surf_props: the SurfProps struct where the computed surface properties will be stored
• Δt: the convection timestep

convect_particles_and_compute_cell_periodic!(grid::Grid1DUniform, particles::ParticleVector{D}, pia, species, Δt) where D

Convect particles on a 1-D uniform grid and write post-convection cell index to particles.cell assuming a periodic grid.

Positional arguments

• grid: the grid on which the convection is performed
• particles: the ParticleVector of particles to be convected
• pia: the ParticleIndexerArray instance
• species: the index of the species being convected
• Δt: the convection timestep

MaxwellWallBC1D <: AbstractBC

A struct to hold information about a Maxwell reflecting wall (mixture of specular and diffuse scattering) orthogonal to the x-axis. This is species-specific, as reflection_velocity_sq is dependent on the species' mass.

Fields

• T: temperature
• v: wall velocity vector
• accommodation: accommodation coefficient (a value of 0 corresponds to specular reflection, a value of 1 corresponds to purely diffuse reflection) - these are implemented in separate more efficient functions as well
• reflection_velocity_sq: pre-computed squared thermal reflection velocity

MaxwellWallBC1D(species, species_data, T::Float64, v, accommodation::Float64)

Construct a MaxwellWallBC1D instance for a given species.

Positional arguments

• species: index of the species for which to create the BC
• species_data: the list of Species data
• T: wall temperature
• v: wall velocity vector
• accommodation: accommodation coefficient

FullyDiffuseBC1D <: AbstractBC

A struct to hold information about a fully diffuse reflecting wall orthogonal to the x-axis. This is species-specific, as reflection_velocity_sq is dependent on the species' mass.

Fields

• T: temperature
• v: wall velocity vector
• reflection_velocity_sq: pre-computed squared thermal reflection velocities

FullyDiffuseBC1D(species, species_data, T::Float64, v)

Construct a FullyDiffuseBC1D instance for a given species.

Positional arguments

• species: index of the species for which to create the BC
• species_data: the list of Species data
• T: wall temperature
• v: wall velocity vector

FullyDiffuseBC1D <: AbstractBC

A struct to hold information about a fully specularly reflecting wall orthogonal to the x-axis. Since such reflection simply flips the sign of the particle's x-velocity, no data is stored in the struct.

write_grid(nc_filename, grid::Grid1DUniform; global_attributes=Dict{Any,Any}())

Write grid info to a NetCDF file

Positional arguments

• nc_filename: path to the NetCDF file
• grid: the grid to be written out to a file

Keyword arguments

• global_attributes: dictionary of additional global attributes to be written to the netCDF file

IOSkipList

Struct that holds track of which variables are not to be written to NetCDF file for physical properties computed on a grid. If the field value is true, the corresponding physical grid property will not be output to the file.

Fields

• skip_length_particle_array: whether the length of the particle array should be skipped
• skip_number_of_particles: whether the output of the number of particles should be skipped
• skip_number_density: whether the output of the number density/number of physical particles should be skipped
• skip_velocity: whether the output of the velocity should be skipped
• skip_temperature: whether the output of the temperature should be skipped

IOSkipList(list_of_variables_to_skip)

Construct an IOSkipList from a list of variable names. The possible names are: length_particle_array, np or nparticles, ndens, v, T.

Positional arguments

• list_of_variables_to_skip: list of variable names to skip

IOSkipList()

Construct an empty IOSkipList.

IOSkipListSurf

Struct that holds track of which variables are not to be written to NetCDF file for computed surface properties. If the field value is true, the corresponding surface property will not be output to the file.

Fields

• skip_number_of_particles: whether the output of the number of particles should be skipped
• skip_fluxes: whether the output of the incident/reflected fluxes should be skipped
• skip_force: whether the output of the force should be skipped
• skip_normal_pressure: whether the output of the normal pressure should be skipped
• skip_shear_pressure: whether the output of the shear pressure should be skipped
• skip_kinetic_energy_flux: whether the output of the kinetic energy flux should be skipped

IOSkipListSurf(list_of_variables_to_skip)

Construct an IOSkipListSurf from a list of variable names. The possible names are: np or nparticles, fluxes, force, normal_pressure, shear_pressure, kinetic_energy_flux.

Positional arguments

• list_of_variables_to_skip: list of variable names to skip

IOSkipListSurf()

Construct an empty IOSkipListSurf

IOSkipListFlux

Struct that holds track of which variables are not to be written to NetCDF file for computed fluxes. If the field value is true, the corresponding flux will not be output to the file.

Fields

• skip_kinetic_energy_flux: whether the output of the kinetic energy flux should be skipped
• skip_diagonal_momentum_flux: whether the output of the diagonal components of the momentum flux tensor should be skipped
• skip_off_diagonal_momentum_flux: whether the output of the off-diagonal components of the momentum flux tensor should be skipped

IOSkipListFlux(list_of_variables_to_skip)

Construct an IOSkipListFlux from a list of variable names. The possible names are: kinetic_energy_flux, diagonal_momentum_flux, off_diagonal_momentum_flux.

Positional arguments

• list_of_variables_to_skip: list of variable names to skip

IOSkipListFlux()

Construct an empty IOSkipListFlux

IOSkipListField

Struct that holds track of which variables are not to be written to NetCDF file for the electrostatic field quantities computed on the nodes of a grid. If the field value is true, the corresponding field quantity will not be output to the file.

Fields

• skip_charge_density: whether the output of the charge density should be skipped
• skip_potential: whether the output of the electrostatic potential should be skipped
• skip_electric_field: whether the output of the electric field should be skipped

IOSkipListField(list_of_variables_to_skip)

Construct an IOSkipListField from a list of variable names. The possible names are: rho or charge_density, phi or potential, E or electric_field.

Positional arguments

• list_of_variables_to_skip: list of variable names to skip

IOSkipListField()

Construct an empty IOSkipListField.

NCDataHolder <: AbstractNCDataHolder

Struct that holds NetCDF-output related data for physical properties (grid properties) I/O.

Fields

• filehandle: handle to the open NetCDF file
• ndens_not_Np: whether the number density or the number of physical particles is being output
• timestep_dim: timestep dimension that used to keep track of the number of output steps
• v_timestep: variable to hold the simulation timestep number (dimension time)
• v_lpa: variable to hold lengths of particle arrays (dimension n_species x time)
• v_np: variable to hold number of particles (dimension n_cells x n_species x time)
• v_ndens: variable to hold number density or the number of physical particles (dimension n_cells x n_species x time)
• v_v: variable to hold velocity (dimension 3 x n_cells x n_species x time)
• v_T: variable to hold temperature (dimension n_cells x n_species x time)
• n_species_1: constant vector [n_species, 1] (used for offsets during I/O)
• n_cells_n_species_1: constant vector [n_cells, n_species, 1] (used for offsets during I/O)
• n_v_n_cells_n_species_1: constant vector [3, n_cells, n_species, 1] (used for offsets during I/O)
• currtimesteps: vector [n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1: vector [1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1: vector [1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1_1: vector [1, 1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• timestep: vector storing the current simulation timestep
• skip_list: IOSkipList instance of variables to skip during output

NCDataHolder(nc_filename, names_skip_list, species_data, phys_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolder instance with a list of variables to skip.

Positional arguments

• nc_filename: filename to write output to
• names_skip_list: list of variable names to skip, see IOSkipList for more details
• species_data: the vector of Species data for the species in the simulation
• phys_props: the PhysProps instance which will be used for the computation and output of physical properties

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolder(nc_filename, species_data, phys_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolder instance with an empty list of variable to skip,

Positional arguments

• nc_filename: filename to write output to
• species_data: the vector of Species data for the species in the simulation
• phys_props: the PhysProps instance which will be used for the computation and output of physical properties

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderSurf

Struct that holds NetCDF-output related data for surface properties I/O.

Fields

• filehandle: handle to the open NetCDF file
• timestep_dim: timestep dimension that used to keep track of the number of output steps
• v_timestep: variable to hold the simulation timestep number (dimension time)
• v_np: variable to hold number of particles that hit the surface (dimension n_elements x n_species x time)
• v_flux_incident: variable to hold incident mass flux (dimension n_elements x n_species x time)
• v_flux_reflected: variable to hold reflected mass flux (dimension n_elements x n_species x time)
• v_force: variable to hold force (dimension 3 x n_elements x n_species x time)
• v_normal_pressure: variable to hold normal pressure (dimension n_elements x n_species x time)
• v_shear_pressure: variable to hold shear pressure (dimension 3 x n_elements x n_species x time)
• v_kinetic_energy_flux: variable to hold kinetic energy flux (dimension n_elements x n_species x time)
• n_species_1: constant vector [n_species, 1] (used for offsets during I/O)
• n_elements_n_species_1: constant vector [n_elements, n_species, 1] (used for offsets during I/O)
• n_v_n_elements_n_species_1: constant vector [3, n_elements, n_species, 1] (used for offsets during I/O)
• currtimesteps: vector [n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1: vector [1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1: vector [1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1_1: vector [1, 1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• timestep: vector storing the current simulation timestep
• skip_list: IOSkipListSurf instance of variables to skip during output

NCDataHolderSurf(nc_filename, names_skip_list, species_data, surf_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderSurf instance with a list of variables to skip.

Positional arguments

• nc_filename: filename to write output to
• names_skip_list: list of variable names to skip, see IOSkipListSurf for more details
• species_data: the vector of Species data for the species in the simulation
• surf_props: the SurfProps instance which will be used for the computation and output of surface properties

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderSurf(nc_filename, species_data, surf_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderSurf instance with an empty list of variable to skip.

Positional arguments

Positional arguments

• nc_filename: filename to write output to
• species_data: the vector of Species data for the species in the simulation
• surf_props: the SurfProps instance which will be used for the computation and output of surface properties

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderFlux

Struct that holds NetCDF-output related data for fluxes I/O.

Fields

• filehandle: handle to the open NetCDF file
• timestep_dim: timestep dimension that used to keep track of the number of output steps
• v_timestep: variable to hold the simulation timestep number (dimension time)
• v_kinetic_energy_flux: variable to hold kinetic energy flux (dimension 3 x n_elements x n_species x time)
• v_diagonal_momentum_flux: variable to the diagonal components of the momentum flux tensor (dimension 3 x n_elements x n_species x time)
• v_off_diagonal_momentum_flux: variable to the off-diagonal components of the momentum flux tensor (dimension 3 x n_elements x n_species x time)
• n_v_n_elements_n_species_1: constant vector [3, n_elements, n_species, 1] (used for offsets during I/O)
• currtimesteps: vector [n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1_1: vector [1, 1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• timestep: vector storing the current simulation timestep
• skip_list: IOSkipListFlux instance of variables to skip during output

NCDataHolderFlux(nc_filename, names_skip_list, species_data, flux_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderFlux instance with a list of variables to skip.

Positional arguments

• nc_filename: filename to write output to
• names_skip_list: list of variable names to skip, see IOSkipListSurf for more details
• species_data: the vector of Species data for the species in the simulation
• flux_props: the FluxProps instance which will be used for the computation and output of fluxes

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderFlux(nc_filename, species_data, flux_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderFlux instance with an empty list of variable to skip.

Positional arguments

Positional arguments

• nc_filename: filename to write output to
• species_data: the vector of Species data for the species in the simulation
• flux_props: the FluxProps instance which will be used for the computation and output of fluxes

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderMoments <: AbstractNCDataHolder

Struct that holds NetCDF-output related data for moments I/O.

Fields

• filehandle: handle to the open NetCDF file
• timestep_dim: timestep dimension that used to keep track of the number of output steps
• v_timestep: variable to hold the simulation timestep number (dimension time)
• v_mompows: variable to hold list of total moment powers (dimension n_moments)
• v_moments: variable to hold total moments (dimension n_moments x n_cells x n_species x time)
• n_moments: number of moments
• n_cells: number of cells
• n_species: number of species
• n_moments_n_cells_n_species_1: constant vector [n_moments, n_cells, n_species, 1] (used for offsets during I/O)
• currtimesteps: vector [n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1_1_1: vector [1, 1, 1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• timestep: vector storing the current simulation timestep

NCDataHolderMoments(nc_filename, species_data, n_cells, n_species, moment_powers; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderMoments instance.

Positional arguments

• nc_filename: filename to write output to
• species_data: the vector of Species data for the species in the simulation
• n_cells: number of cells
• n_species: number of species
• moment_powers: vector of moment powers (e.g., Int8[2, 4, 6])

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, can use NC_NETCDF4 for NetCDF4 format)

NCDataHolderField <: AbstractNCDataHolder

Struct that holds NetCDF-output related data for the I/O of the electrostatic field quantities stored in the nodes of a grid.

Fields

• filehandle: handle to the open NetCDF file
• timestep_dim: timestep dimension that used to keep track of the number of output steps
• v_timestep: variable to hold the simulation timestep number (dimension time)
• v_charge_density: variable to hold the charge density (dimension n_nodes x time)
• v_potential: variable to hold the electrostatic potential (dimension n_nodes x time)
• v_electric_field: variable to hold the x-component of the electric field (dimension n_nodes x time)
• n_nodes: number of nodes
• n_nodes_1: constant vector [n_nodes, 1] (used for offsets during I/O)
• currtimesteps: vector [n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• currtimesteps_1: vector [1, n_t_output], where n_t_output is the current output timestep (i.e. how many times the properties have already been output, not the simulation timestep) (used for offsets during I/O)
• timestep: vector storing the current simulation timestep
• skip_list: IOSkipListField instance of variables to skip during output

NCDataHolderField(nc_filename, names_skip_list, field_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderField instance with a list of variables to skip.

Positional arguments

• nc_filename: filename to write output to
• names_skip_list: list of variable names to skip, see IOSkipListField for more details
• field_props: the ElectrostaticFieldProps instance which will be used for the output of the field quantities

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, NC_NETCDF4 for NetCDF4 format)

NCDataHolderField(nc_filename, field_props; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Construct a NCDataHolderField instance with an empty list of variables to skip.

Positional arguments

• nc_filename: filename to write output to
• field_props: the ElectrostaticFieldProps instance which will be used for the output of the field quantities

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute
• mode: NetCDF file format mode (default: NC_64BIT_OFFSET for older and faster format, NC_NETCDF4 for NetCDF4 format)

 write_netcdf(ds, phys_props::PhysProps, timestep; sync_freq=0)

Write computed PhysProps to NetCDF file and synchronize file to disk if necessary.

Positional arguments

• ds: the NCDataHolder for the file to which the output will be written
• phys_props: the PhysProps instance containing the computed properties
• timestep: the simulation timestep

Keyword arguments

• sync_freq: if larger than 0 and if the number of timesteps output is proportional to sync_freq, the data will be synchronized to disk. If set to 1, will sync data to disk at every timestep at which data is written to the file.

Throws

ErrorException if the NCDataHolder expects number density and the phys_props holds the number of physical particles, or vice versa.

write_netcdf(ds, surf_props::SurfProps, timestep; sync_freq=0)

Write SurfProps to a NetCDF file and synchronize file to disk if necessary.

Positional arguments

• ds: the NCDataHolderSurf for the file to which the output will be written
• surf_props: the SurfProps instance containing the computed properties
• timestep: the simulation timestep

Keyword arguments

• sync_freq: if larger than 0 and if the number of timesteps output is proportional to sync_freq, the data will be synchronized to disk. If set to 1, will sync data to disk at every timestep at which data is written to the file.

write_netcdf(ds, flux_props::FluxProps, timestep; sync_freq=0)

Write FluxProps to a NetCDF file and synchronize file to disk if necessary.

Positional arguments

• ds: the NCDataHolderFlux for the file to which the output will be written
• flux_props: the FluxProps instance containing the computed properties
• timestep: the simulation timestep

Keyword arguments

• sync_freq: if larger than 0 and if the number of timesteps output is proportional to sync_freq, the data will be synchronized to disk. If set to 1, will sync data to disk at every timestep at which data is written to the file.

write_netcdf(nc_filename, particles::Vector{ParticleVector{D}}, pia, species_data, species_ids; global_attributes=Dict{Any,Any}(), mode=NC_64BIT_OFFSET)

Write particles of species listed in species_ids to a NetCDF file. The particles are written cell-wise, so the ordering is not preserved in case particles are present in the set of indices pointed to by group2 indices. If D==0 (particle position is not tracked), neither the particle position nor the cell index is written to the file.

Positional arguments

• nc_filename: filename to write output to
• particles: the Vector of ParticleVector instances of particles to be written
• pia: the ParticleIndexerArray instance
• species_data: the vector of Species data for the species in the simulation
• species_ids: list of species ids (in range 1:n_species) of species for which to write particle data

Keyword arguments

• global_attributes: dictionary of any additional attributes to write to the netCDF file as a global attribute

write_netcdf(ds::NCDataHolderMoments, moments, timestep; sync_freq=0)

Write 'moments' to a NetCDF file and synchronize file to disk if necessary.

Positional arguments

• ds: the NCDataHolderMoments for the file to which the output will be written
• moments: the moments array (dimension n_moments x n_cells x n_species)
• timestep: the simulation timestep

Keyword arguments

• sync_freq: if larger than 0 and if the number of timesteps output is proportional to sync_freq, the data will be synchronized to disk. If set to 1, will sync data to disk at every timestep at which data is written to the file.

write_netcdf(ds, field_props::ElectrostaticFieldProps, timestep; sync_freq=0)

Write the electrostatic field quantities stored in an ElectrostaticFieldProps instance to a NetCDF file and synchronize file to disk if necessary.

Positional arguments

• ds: the NCDataHolderField for the file to which the output will be written
• field_props: the ElectrostaticFieldProps instance containing the field quantities
• timestep: the simulation timestep

Keyword arguments

• sync_freq: if larger than 0 and if the number of timesteps output is proportional to sync_freq, the data will be synchronized to disk. If set to 1, will sync data to disk at every timestep at which data is written to the file.

close_netcdf(ds::T) where {T<:AbstractNCDataHolder}

Close NetCDF file.

Positional arguments

• ds: an AbstractNCDataHolder instance to close.

ChunkExchanger

Struct used to organize exchange of particles between independent ParticleVector instances for chunked multi-threaded simulations. It is assumed that the cell indices within each chunk are contiguous, i.e. chunk[i+1] = chunk[i]+1.

The indexing for the exchanged groups is stored as arrays of shape (n_chunks, n_cells), and [chunk_id, cell] describes the particles that belong to cell cell and came from particle chunk chunk_id. Group 1 holds the particles that arrived via swapping, group 2 the particles that arrived via pushing, i.e. those added to the end of the particle array. The end of a group is not stored, as it is given by start + n_group - 1 (which is -1 for an empty group, matching the convention used by ParticleIndexer).

The table is self-clearing: every entry written during an exchange is read exactly once, in the same timestep, by the chunk owning the cell, and sort_particles_after_exchange! zeroes it as it reads it. It is therefore already empty at the start of each timestep and does not need to be reset! between timesteps.

The occ_lo/occ_hi fields hold, per chunk, the first and last cell in which the chunk holds any particles. They let exchange_particles! reject a pair of chunks that cannot have anything to exchange without scanning the cells of either chunk. They are updated by update_occupancy_bounds! and start out as the whole grid, i.e. a simulation that never calls it is still correct, just slower.

Fields

• n_chunks: number of chunks used in the simulation
• n_cells: number of grid cells in the simulation
• start1: index of the first particle of group 1
• n_group1: number of particles in group 1
• start2: index of the first particle of group 2
• n_group2: number of particles in group 2
• occ_lo: first cell in which a chunk holds particles
• occ_hi: last cell in which a chunk holds particles

ChunkExchanger(chunks, n_cells)

Create a ChunkExchanger for length(chunks) chunks and n_cells cell.

Positional arguments

• chunks: list of cell chunks
• n_cells: total number of cells in the simulation

update_occupancy_bounds!(chunk_exchanger, gridsort, pia, chunk_id, species)

Update the first and last cell in which chunk chunk_id holds particles of the given species, so that exchange_particles! can reject pairs of chunks with nothing to exchange without scanning any cells.

Must be called after sort_particles! and before exchange_particles!, as it simply copies over the bounds recorded by the sort. If not called, one should set occ_lo to 1 and occ_hi to n_cells in each chunk.

Positional arguments

• chunk_exchanger: the ChunkExchanger instance in which to store the bounds
• gridsort: the GridSortInPlace used to sort the chunk's particles
• pia: the ParticleIndexerArray instance associated with the chunk
• chunk_id: the chunk for which to update the bounds
• species: the particle species for which the bounds are computed

exchange_particles!(chunk_exchanger, particles_chunks::Vector{Vector{ParticleVector{D}}}, pia_chunks, cell_chunks, species, i, j) where D

Redistribute particles between chunks i and j based on their spatial cell ownership.

This function ensures each particle resides in the chunk responsible for its current cell. It performs symmetric swaps when possible, and pushes remaining particles if needed. The indexing metadata (chunk_exchanger and pia_chunks) is updated accordingly, and particles pushed to another chunk (not swapped) are added to the buffer for future re-use. The particles before the start of the re-distribution need to be sorted, so that no particles are indexed by the start2:end2 part of a ParticleIndexer. After the operation, the n_total[species] value of pia_chunks[chunk_id] will not include particles that were pushed to another chunk (the appropriate n_group1, start1, end1 values will be set to 0, 0, -1). However indexing should not be relied on until particles are re-sorted, see (sort_particles_after_exchange!)[@ref].

Positional arguments

• chunk_exchanger: the ChunkExchanger instance to track post-swap and post-push indices
• particles_chunks: Vector of Vector of ParticleVector (per chunk and per species, i.e. particles_chunks[chunk_id][species] is the correct order of access)
• pia_chunks: Vector of ParticleIndexerArray instances for each chunk
• cell_chunks: cell ownership list for each chunk, i.e. cell_chunks[chunk_id] is a list of cells belonging to chunk chunk_id; the cells within cell_chunks[chunk_id] should be ordered in increasing order and be continuous: i.e. cell_chunks[chunk_id][i] == cell_chunks[chunk_id][i-1] + 1
• species: the particle species being redistributed
• i: index of first chunk
• j: index of second chunk

exchange_particles!(chunk_exchanger, particles_chunks::Vector{Vector{ParticleVector{D}}}, pia_chunks, cell_chunks, species) where D

Redistribute particles between chunks based on their spatial cell ownership.

This function ensures each particle resides in the chunk responsible for its current cell. It performs symmetric swaps when possible, and pushes remaining particles if needed. The indexing metadata (chunk_exchanger and pia_chunks) is updated accordingly, and particles pushed to another chunk (not swapped) are added to the buffer for future re-use. The particles before the start of the re-distribution need to be sorted, so that no particles are indexed by the start2:end2 part of a ParticleIndexer. After the operation, the n_total[species] value of pia_chunks[chunk_id] will not include particles that were pushed to another chunk (the appropriate n_group1, start1, end1 values will be set to 0, 0, -1). However indexing should not be relied on until particles are re-sorted, see (sort_particles_after_exchange!)[@ref].

Positional arguments

• chunk_exchanger: the ChunkExchanger instance to track post-swap and post-push indices
• particles_chunks: Vector of Vector of ParticleVector (per chunk and per species, i.e. particles_chunks[chunk_id][species] is the correct order of access)
• pia_chunks: Vector of ParticleIndexerArray instances for each chunk
• cell_chunks: cell ownership list for each chunk, i.e. cell_chunks[chunk_id] is a list of cells belonging to chunk chunk_id; the cells within cell_chunks[chunk_id] should be ordered in increasing order and be continuous: i.e. cell_chunks[chunk_id][i] == cell_chunks[chunk_id][i-1] + 1
• species: the particle species being redistributed