Internal
Public APIs
TestParticle.AdaptiveBoris Method
AdaptiveBoris(; dtmin, dtmax, safety=0.1)Adaptive Boris method with adaptive time stepping based on local gyrofrequency. The time step is determined by dt = safety / |qB/m|, clamped by dtmin and dtmax.
TestParticle.CurrentLoop Type
CurrentLoop{T}A circular current loop.
Fields
radius::T: Radius of the loop [m].current::T: Current in the loop [A].center::SVector{3, T}: Position of the loop center [m].normal::SVector{3, T}: Unit normal vector of the loop.
TestParticle.BiKappa Method
BiKappa(args...; kw...)Construct a BiKappa distribution. Forwards to VelocityDistributionFunctions.BiKappa.
BiKappa(B, u0, ppar, pperp, n, kappa; m=mᵢ)Construct a BiKappa distribution with magnetic field B, bulk velocity u0, parallel thermal pressure ppar, perpendicular thermal pressure pperp, number density n, and spectral index kappa in SI units. The default particle is proton.
TestParticle.BiMaxwellian Method
BiMaxwellian(args...; kw...)Construct a BiMaxwellian distribution. Forwards to VelocityDistributionFunctions.BiMaxwellian.
BiMaxwellian(B, u0, ppar, pperp, n; m=mᵢ)Construct a BiMaxwellian distribution with magnetic field B, bulk velocity u0, parallel thermal pressure ppar, perpendicular thermal pressure pperp, and number density n in SI units. The default particle is proton.
TestParticle.Kappa Method
Kappa(args...; kw...)Construct a Kappa distribution. Forwards to VelocityDistributionFunctions.Kappa.
Kappa(u0, p, n, kappa; m=mᵢ)Construct a Kappa distribution with bulk velocity u0, thermal pressure p, number density n, and spectral index kappa in SI units. The default particle is proton.
TestParticle.Maxwellian Method
Maxwellian(args...; kw...)Construct a Maxwellian distribution. Forwards to VelocityDistributionFunctions.Maxwellian.
Maxwellian(u0, p, n; m=mᵢ)Construct a Maxwellian distribution with bulk velocity u0, thermal pressure p, and number density n in SI units. The default particle is proton.
TestParticle.getB_bottle Method
getB_bottle(x, y, z, distance, a, b, I1, I2) :: StaticVector{3}Get magnetic field from a magnetic bottle. Reference: wiki
Arguments
x,y,z::Float: particle coordinates in [m].distance::Float: distance between solenoids in [m].a::Float: radius of each side coil in [m].b::Float: radius of central coil in [m].I1::Float: current in the solenoid times number of windings in side coils in [A].I2::Float: current in the central solenoid times number of windings in the central loop in [A].
TestParticle.getB_loop Method
getB_loop(r, R, a, I, n)Calculate the magnetic field B [T] at point r from a current loop with current I [A], radius a [m], centered at R, and normal vector n.
TestParticle.getB_loop Method
getB_loop(r, loop::CurrentLoop)Calculate the magnetic field B [T] at point r from a CurrentLoop.
TestParticle.getB_mirror Method
getB_mirror(x, y, z, distance, a, I1) :: StaticVector{3}Get magnetic field at [x, y, z] from a magnetic mirror generated from two coils.
Arguments
x,y,z: particle coordinates in [m].distance: distance between solenoids in [m].a: radius of each side coil in [m].I1: current in the solenoid times number of windings in side coils.
TestParticle.getB_tokamak_coil Method
getB_tokamak_coil(x, y, z, a, b, ICoils, IPlasma) -> StaticVector{3}Get the magnetic field from a Tokamak topology consists of 16 coils. Original: Tokamak-Fusion-Reactor
Arguments
x,y,z: location in [m].a: radius of each coil in [m].b: radius of central region in [m].ICoil: current in the coil times number of windings in [A].IPlasma: current of the plasma in [A].
TestParticle.getB_zpinch Method
getB_zpinch(x, y, z, I, a) -> StaticVector{Float64, 3}Get magnetic field from a Z-pinch configuration. Reference: Z-pinch
Arguments
x,y,z::Float: particle coordinates in [m].I::Float: current in the wire [A].a::Float: radius of the wire [m].
TestParticle.get_gc Method
get_gc(xu, param)
get_gc(x, y, z, vx, vy, vz, bx, by, bz, q2m)Calculate the coordinates of the guiding center according to the phase space coordinates of a particle. Reference: wiki
Nonrelativistic definition:
TestParticle.get_gc_func Method
get_gc_func(param)Return the function for plotting the orbit of guiding center.
Example
param = prepare(E, B; species = Proton)
# The definitions of stateinit, tspan, E and B are ignored.
prob = ODEProblem(trace!, stateinit, tspan, param)
sol = solve(prob, Vern7(); dt = 2e-11)
f = Figure(fontsize = 18)
ax = Axis3(f[1, 1], aspect = :data)
gc = param |> get_gc_func
gc_plot(x, y, z, vx, vy, vz) = (gc(SA[x, y, z, vx, vy, vz])...,)
lines!(ax, sol, idxs = (gc_plot, 1, 2, 3, 4, 5, 6))TestParticle.get_gyrofrequency Function
get_gyrofrequency(B=5e-9; q=qᵢ, m=mᵢ)Return the gyrofrequency [rad/s].
Arguments
B: Magnetic field magnitude [T]. Default is 5 nT.q: Charge [C]. Default is proton charge.m: Mass [kg]. Default is proton mass.
TestParticle.get_gyroperiod Function
get_gyroperiod(B=5e-9; q=qᵢ, m=mᵢ)Return the gyroperiod [s].
Arguments
B: Magnetic field magnitude [T]. Default is 5 nT.q: Charge [C]. Default is proton charge.m: Mass [kg]. Default is proton mass.
TestParticle.get_gyroradius Method
get_gyroradius(V, B; q=qᵢ, m=mᵢ)Return the gyroradius [m].
Arguments
V: Velocity magnitude [m/s] (usually perpendicular to the magnetic field).B: Magnetic field magnitude [T].q: Charge [C]. Default is proton charge.m: Mass [kg]. Default is proton mass.
TestParticle.get_number_density_flux Method
get_number_density_flux(grid::CartesianGrid, sols, dt)Calculate the steady state particle number density flux on a uniform Cartesian grid. The flux is estimated by accumulating the number of particles in each cell at time steps dt, divided by the surface area of each cell.
Arguments
grid: ACartesianGridfromMeshes.jl.sols: Particle trajectory solutions (e.g.EnsembleSolution).dt: Time step for sampling particle positions.
TestParticle.prepare Function
prepare(args...; kwargs...) -> (q2m, m, E, B, F)
prepare(E, B, F = ZeroField(); kwargs...)
prepare(grid::CartesianGrid, E, B, F = ZeroField(); kwargs...)
prepare(x, E, B, F = ZeroField(); dir = 1, kwargs...)
prepare(x, y, E, B, F = ZeroField(); kwargs...)
prepare(x, y, z, E, B, F = ZeroField(); kwargs...)
prepare(B; E = ZeroField(), F = ZeroField(), kwargs...)Return a tuple consists of particle charge-mass ratio for a prescribed species of charge q and mass m, mass m for a prescribed species, analytic/interpolated EM field functions, and external force F.
Prescribed species are Electron and Proton; other species can be manually specified with m and q keywords or species = Ion(m̄, q̄), where m̄ and q̄ are the mass and charge numbers respectively.
Direct range input for uniform grid in 1/2/3D is supported. For 1D grid, an additional keyword dir is used for specifying the spatial direction, 1 -> x, 2 -> y, 3 -> z. For 3D grid, the default grid type is CartesianGrid. To use StructuredGrid (spherical) grid, an additional keyword gridtype is needed. For StructuredGrid (spherical) grid, dimensions of field arrays should be (Br, Bθ, Bϕ).
Keywords
order::Int=1: order of interpolation in [1,2,3].bc::Int=1: type of boundary conditions, 1 -> NaN, 2 -> periodic.species=Proton: particle species.q=nothing: particle charge.m=nothing: particle mass.gridtype:CartesianGrid,RectilinearGrid,StructuredGrid.
TestParticle.sample_unit_sphere Method
sample_unit_sphere()Sample a unit vector on a sphere uniformly.
sourceTestParticle.trace! Method
trace!(dy, y, p, t)ODE equations for charged particle moving in EM field and external force field with in-place form.
sourceTestParticle.trace Method
trace(y, p, t)::SVector{6}ODE equations for charged particle moving in EM field and external force field with out-of-place form.
sourceTestParticle.trace_fieldline! Method
trace_fieldline!(dx, x, p, s)Equation for tracing magnetic field lines with in-place form. The parameter p is the magnetic field function. Note that the independent variable s represents the arc length.
TestParticle.trace_fieldline Method
trace_fieldline(x, p, s)Equation for tracing magnetic field lines with out-of-place form.
sourceTestParticle.trace_fieldline Method
trace_fieldline(u0, p, tspan::Tuple{<:Real, <:Real}; mode::Symbol=:both, kw...)Helper function to create ODEProblems for tracing magnetic field lines.
Arguments
u0: Initial position.p: Parameter tuple containing the magnetic field. It can also be the magnetic field array (numerical) or function (analytical/interpolator).tspan: Time span (arc length) for tracing. typically(0.0, L).mode: Tracing mode, one of:forward,:backward,:both.kw: Keyword arguments passed topreparewhenpis an array.
Returns
If
modeis:forwardor:backward, returns a singleODEProblem.If
modeis:both, returns a vector of twoODEProblems (forward and backward).
TestParticle.trace_gc! Method
trace_gc!(dy, y, p, t)Guiding center equations for nonrelativistic charged particle moving in EM field with in-place form. Variable y = (x, y, z, u), where u is the velocity along the magnetic field at (x,y,z).
TestParticle.trace_gc_drifts! Method
trace_gc_drifts!(dx, x, p, t)Equations for tracing the guiding center using analytical drifts, including the grad-B drift, curvature drift, and ExB drift. Parallel velocity is also added. This expression requires the full particle trajectory p.sol.
TestParticle.trace_gc_exb! Method
trace_gc_exb!(dx, x, p, t)Equations for tracing the guiding center using the ExB drift and parallel velocity from a reference trajectory.
sourceTestParticle.trace_gc_flr! Method
trace_gc_flr!(dx, x, p, t)Equations for tracing the guiding center using the ExB drift with FLR corrections and parallel velocity.
sourceTestParticle.trace_normalized! Method
trace_normalized!(dy, y, p, t)Normalized ODE equations for charged particle moving in EM field with in-place form. If the field is in 2D X-Y plane, periodic boundary should be applied for the field in z via the extrapolation function provided by Interpolations.jl.
sourceTestParticle.trace_relativistic! Method
trace_relativistic!(dy, y, p, t)ODE equations for relativistic charged particle (x, γv) moving in EM field with in-place form.
sourceTestParticle.trace_relativistic Method
trace_relativistic(y, p, t) -> SVector{6}ODE equations for relativistic charged particle (x, γv) moving in static EM field with out-of-place form.
sourceTestParticle.trace_relativistic_normalized! Method
trace_relativistic_normalized!(dy, y, p, t)Normalized ODE equations for relativistic charged particle (x, γv) moving in EM field with in-place form.
sourceTestParticle.trace_relativistic_normalized Method
trace_relativistic_normalized(y, p, t)Normalized ODE equations for relativistic charged particle (x, γv) moving in EM field with out-of-place form.
sourcePrivate types and methods
TestParticle.Field Type
Field{itd, F} <: AbstractField{itd}A representation of a field function f, defined by:
time-independent field
time-dependent field
Arguments
field_function::Function: the function of field.itd::Bool: whether the field function is time dependent.F: the type offield_function.
TestParticle.dipole Method
Calculates the magnetic field from a dipole with magnetic moment M at r.
TestParticle.dipole_fieldline Function
dipole_fieldline(ϕ, L=2.5, nP=100)Creates nP points on one field line of the magnetic field from a dipole. In a centered dipole magnetic field model, the path along a given L shell can be described as r = L*cos²λ, where r is the radial distance (in planetary radii) to a point on the line, λ is its co-latitude, and L is the L-shell of interest.
TestParticle.getB_CS_harris Function
getB_CS_harris(B₀, L)Return the magnetic field at location r near a current sheet with magnetic strength B₀ and sheet length L. The current sheet is assumed to lie in the z = 0 plane.
TestParticle.getB_dipole Method
Analytic magnetic field function for testing. Return in SI unit.
sourceTestParticle.getB_tokamak_profile Method
getB_tokamak_profile(x, y, z, q_profile, a, R₀, Bζ0) :: StaticVector{3}Reconstruct the magnetic field distribution from a safe factor(q) profile. Reference: Tokamak, 4th Edition, John Wesson.
Arguments
x,y,z: location in [m].q_profile::Function: profile of q. The variable of this function must be the normalized radius.a: minor radius [m].R₀: major radius [m].Bζ0: toroidal magnetic field on axis [T].
TestParticle.get_interpolator Method
get_interpolator(A, gridx, gridy, gridz, order::Int=1, bc::Int=1)
get_interpolator(gridtype, A, grid1, grid2, grid3, order::Int=1, bc::Int=1)Return a function for interpolating field array A on the grid given by gridx, gridy, and gridz.
Arguments
gridtype:CartesianGrid,RectilinearGridorStructuredGrid.A: field array. For vector field, the first dimension should be 3.order::Int=1: order of interpolation in [1,2,3].bc::Int=1: type of boundary conditions, 1 -> NaN, 2 -> periodic, 3 -> Flat.
Notes
The input array A may be modified in-place for memory optimization.
TestParticle.get_number_density Method
get_number_density(sols, grid, t_start, t_end, dt)Calculate time-averaged particle number density from t_start to t_end with step dt.
TestParticle.get_number_density Method
get_number_density(sols, grid, t)Calculate particle number density at time t in a given grid.
TestParticle.get_rotation_matrix Method
get_rotation_matrix(axis::AbstractVector, angle) :: SMatrix{3,3}Create a rotation matrix for rotating a 3D vector around a unit axis by an angle in radians. Reference: Rotation matrix from axis and angle
Example
using LinearAlgebra
v = [-0.5, 1.0, 1.0]
v̂ = normalize(v)
θ = deg2rad(-74)
R = get_rotation_matrix(v̂, θ)TestParticle.getinterp Function
getinterp(::Type{<:CartesianGrid}, A, gridx, gridy, gridz, order::Int=1, bc::Int=1)Return a function for interpolating field array A on the grid given by gridx, gridy, and gridz.
Arguments
order::Int=1: order of interpolation in [1,2,3].bc::Int=1: type of boundary conditions, 1 -> NaN, 2 -> periodic, 3 -> Flat.dir::Int: 1/2/3, representing x/y/z direction.
Notes
The input array A may be modified in-place for memory optimization.
TestParticle.getinterp_scalar Function
getinterp_scalar(::Type{<:CartesianGrid}, A, gridx, gridy, gridz, order::Int=1, bc::Int=1)Return a function for interpolating scalar array A on the grid given by gridx, gridy, and gridz. Currently only 3D arrays are supported.
Arguments
order::Int=1: order of interpolation in [1,2,3].bc::Int=1: type of boundary conditions, 1 -> NaN, 2 -> periodic, 3 -> Flat.dir::Int: 1/2/3, representing x/y/z direction.
TestParticle.set_axes_equal Method
set_axes_equal(ax)Set 3D plot axes to equal scale for Matplotlib. Make axes of 3D plot have equal scale so that spheres appear as spheres and cubes as cubes. Required since ax.axis('equal') and ax.set_aspect('equal') don't work on 3D.
TestParticle.solve Function
solve(prob::TraceProblem; trajectories::Int=1, dt::AbstractFloat,
savestepinterval::Int=1, isoutofdomain::Function=ODE_DEFAULT_ISOUTOFDOMAIN,
n::Int=1)Trace particles using the Boris method with specified prob.
keywords
trajectories::Int: number of trajectories to trace.dt::AbstractFloat: time step.savestepinterval::Int: saving output interval.isoutofdomain::Function: a function with input of position and velocity vectorxvthat determines whether to stop tracing.n::Int: number of substeps for the Multistep Boris method. Default is 1 (standard Boris).save_start::Bool: save the initial condition. Default istrue.save_end::Bool: save the final condition. Default istrue.save_everystep::Bool: save the state at everysavestepinterval. Default istrue.
TestParticle.update_velocity! Method
update_velocity!(xv, paramBoris, param, dt, t)Update velocity using the Boris method, Birdsall, Plasma Physics via Computer Simulation. Reference: DTIC
sourceTestParticle.update_velocity_multistep! Method
update_velocity_multistep!(xv, paramBoris, param, dt, t, n)Update velocity using the Multistep Boris method. Reference: Zenitani & Kato 2025
source