Internal

Public APIs

TestParticle.BiMaxwellianMethod
 BiMaxwellian(B::Vector{U}, u0::Vector{T}, 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.

source
TestParticle.MaxwellianMethod
 Maxwellian(u0::AbstractVector{T}, 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.

source
TestParticle.get_gcMethod
 get_gc(param::Union{TPTuple, FullTPTuple})

Get three functions for plotting the orbit of guiding center.

For example:

param = prepare(E, B; species=Proton)
gc = get_gc(param)
# 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_plot(x,y,z,vx,vy,vz) = (gc([x,y,z,vx,vy,vz])...,)
lines!(ax, sol, idxs=(gc_plot, 1, 2, 3, 4, 5, 6))
source
TestParticle.prepareMethod
 prepare(grid::CartesianGrid, E, B; kwargs...) -> (q2m, E, B)

Return a tuple consists of particle charge-mass ratio for a prescribed species and interpolated EM field functions.

keywords

  • order::Int=1: order of interpolation in [1,2,3].
  • bc::Int=1: type of boundary conditions, 1 -> NaN, 2 -> periodic.
  • species::Species=Proton: particle species.
  • q::AbstractFloat=1.0: particle charge. Only works when Species=User.
  • m::AbstractFloat=1.0: particle mass. Only works when Species=User.
 prepare(grid::CartesianGrid, E, B, F; species=Proton, q=1.0, m=1.0) -> (q, m, E, B, F)

Return a tuple consists of particle charge, mass for a prescribed species of charge q and mass m, interpolated EM field functions, and external force F.

 prepare(x::AbstractRange, y::AbstractRange, z::AbstractRange, E, B; kwargs...) -> (q2m, E, B)
 prepare(x, y, E, B; kwargs...) -> (q2m, E, B)

 prepare(x::AbstractRange, E, B; kwargs...) -> (q2m, E, B)

1D grid. An additional keyword dir is used for specifying the spatial direction, 1 -> x, 2 -> y, 3 -> z.

Direct range input for uniform grid in 2/3D is also accepted.

 prepare(E, B; kwargs...) -> (q2m, E, B)

Return a tuple consists of particle charge-mass ratio for a prescribed species of charge q and mass m and analytic EM field functions. Prescribed species are Electron and Proton; other species can be manually specified with species=Ion/User, q and m.

 prepare(E, B, F; kwargs...) -> (q, m, E, B, F)

Return a tuple consists of particle charge, mass for a prescribed species of charge q and mass m, analytic EM field functions, and external force F.

source
TestParticle.sampleMethod
 sample(vdf::Maxwellian)

Sample a 3D velocity from a Maxwellian distribution vdf using the Box-Muller method.

 sample(vdf::BiMaxwellian)

Sample a 3D velocity from a BiMaxwellian distribution vdf using the Box-Muller method.

source
TestParticle.trace!Method
 trace!(dy, y, p::TPTuple, t)
 trace!(dy, y, p::FullTPTuple, t)

ODE equations for charged particle moving in static EM field with in-place form.

ODE equations for charged particle moving in static EM field and external force field with in-place form.

source
TestParticle.traceMethod
 trace(y, p::TPTuple, t) -> SVector{6, Float64}
 trace(y, p::FullTPTuple, t) -> SVector{6, Float64}

ODE equations for charged particle moving in static EM field with out-of-place form.

ODE equations for charged particle moving in static EM field and external force field with out-of-place form.

source
TestParticle.trace_gc!Method
 trace_gc!(dy, y, p::TPTuple, t)

Guiding center equations for nonrelativistic charged particle moving in static 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).

source
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.

source
TestParticle.trace_normalized!Method
 trace_normalized!(dy, y, p::TPNormalizedTuple, t)

Normalized ODE equations for charged particle moving in static 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.

source
TestParticle.trace_relativistic!Method
 trace_relativistic!(dy, y, p::TPTuple, t)

ODE equations for relativistic charged particle (x, γv) moving in static EM field with in-place form.

source
TestParticle.trace_relativisticMethod
 trace_relativistic(y, p::TPTuple, t) -> SVector{6}

ODE equations for relativistic charged particle (x, γv) moving in static EM field with out-of-place form.

source

Private types and methods

TestParticle.FieldType
 Field{itd, F} <: AbstractField{itd}

A representation of a field function f, defined by:

time-independent field

\[\mathbf{F} = F(\mathbf{x}),\]

time-dependent field

\[\mathbf{F} = F(\mathbf{x}, t).\]

Arguments

  • field_function::Function: the function of field.
  • itd::Bool: whether the field function is time dependent.
  • F: the type of field_function.
source
TestParticle.dipole_fieldlineFunction
 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.

source
TestParticle.getB_CS_harrisFunction
 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.

source
TestParticle.getB_bottleMethod
 getB_bottle(x, y, z, distance, a, b, I1, I2) -> StaticVector{Float64, 3}

Get magnetic field from a magnetic bottle. Reference: https://en.wikipedia.org/wiki/Magneticmirror#Magneticbottles

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].

source
TestParticle.getB_current_loopMethod
 getB_current_loop(x, y, z, cl::Currentloop) -> StaticVector{Float64, 3}

Get magnetic field at [x, y, z] from a magnetic mirror generated from two coils.

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].
  • I1::Float: current in the solenoid times number of windings in side coils.
source
TestParticle.getB_mirrorMethod
 getB_mirror(x, y, z, distance, a, I1) -> StaticVector{Float64, 3}

Get magnetic field at [x, y, z] from a magnetic mirror generated from two coils.

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].
  • I1::Float: current in the solenoid times number of windings in side coils.
source
TestParticle.getB_tokamak_coilMethod
 getB_tokamak_coil(x, y, z, a, b, ICoils, IPlasma) -> StaticVector{Float64, 3}

Get the magnetic field from a Tokamak topology consists of 16 coils. Original: https://github.com/BoschSamuel/Simulation-of-a-Tokamak-Fusion-Reactor/blob/master/Simulation2.m

Arguments

  • x,y,z::Float: location in [m].
  • a::Float: radius of each coil in [m].
  • b::Float: radius of central region in [m].
  • ICoil::Float: current in the coil times number of windings in [A].
  • IPlasma::Float: current of the plasma in [A].
source
TestParticle.getB_tokamak_profileMethod
 getB_tokamak_profile(x, y, z, q_profile, a, R₀, Bζ0) -> StaticVector{Float64, 3}

Reconstruct the magnetic field distribution from a safe factor(q) profile. Reference: Tokamak, 4th Edition, John Wesson.

Arguments

  • x,y,z::Float: location in [m].
  • q_profile::Function: profile of q. The variable of this function must be the normalized radius.
  • a::Float: minor radius [m].
  • R₀::Float: major radius [m].
  • Bζ0::Float: toroidal magnetic field on axis [T].
source
TestParticle.get_rotation_matrixMethod
 get_rotation_matrix(axis::AbstractVector, angle::Real) --> 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̂, θ)
source
TestParticle.getchargemassMethod
 getchargemass(species::Species, q::AbstractFloat, m::AbstractFloat)

Return charge and mass for species. if species = User, input q and m are returned.

source
TestParticle.getinterpFunction
 getinterp(A, gridx, gridy, gridz, order::Int=1, bc::Int=1)
 getinterp(A, gridx, gridy, order::Int=1, bc::Int=1)
 getinterp(A, gridx, order::Int=1, bc::Int=1, dir::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.
source
TestParticle.getinterp_scalarFunction
 getinterp_scalar(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.
source
TestParticle.guiding_centerMethod
 guiding_center(xu, param::Union{TPTuple, FullTPTuple})

Calculate the coordinates of the guiding center according to the phase space coordinates of a particle. Reference: https://en.wikipedia.org/wiki/Guiding_center

A simple definition:

\[\mathbf{X}=\mathbf{x}-m\frac{\mathbf{v}\times\mathbf{B}}{qB}\]

source
TestParticle.set_axes_equalMethod
 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.

source
TestParticle.solveFunction
 solve(prob::TraceProblem; trajectories::Int=1,
	 savestepinterval::Int=1, isoutofdomain::Function=ODE_DEFAULT_ISOUTOFDOMAIN)

Trace particles using the Boris method with specified prob.

keywords

  • trajectories::Int: number of trajectories to trace.
  • savestepinterval::Int: saving output interval.
  • isoutofdomain::Function: a function with input of position and velocity vector xv that determines whether to stop tracing.
source
TestParticle.update_velocity!Method
 update_velocity!(xv, paramBoris, param, dt, t)

Update velocity using the Boris method, Birdsall, Plasma Physics via Computer Simulation. Reference: https://apps.dtic.mil/sti/citations/ADA023511

source