HallThruster.LOOKUP_ZSConstant

Lookup for thermal conductivity coefficients, from Table 1 of S. I. Braginskii, in Reviews of Plasma Physics (1965), Vol. 1, p. 205.

source
HallThruster.BohmType
Bohm(c) <: AnomalousTransportModel

Model where the anomalous collision frequency scales with the electron cyclotron frequency ωce times some scaling factor c

source
HallThruster.BraginskiiType
Braginskii conductivity model for fully-ionized plasmasa
S. I. Braginskii, in Reviews of Plasma Physics (1965), Vol. 1, p. 205.
source
HallThruster.ConfigType
struct Config{A<:HallThruster.AnomalousTransportModel, TC<:HallThruster.ThermalConductivityModel, W<:HallThruster.WallLossModel, IZ<:HallThruster.IonizationModel, EX<:HallThruster.ExcitationModel, EN<:HallThruster.ElectronNeutralModel, HET<:HallThruster.Thruster, S_N, S_IC, S_IM, S_ϕ, S_E, IC<:HallThruster.InitialCondition, HS<:HallThruster.HyperbolicScheme}

Hall thruster configuration struct. Only four mandatory fields: discharge_voltage, thruster, anode_mass_flow_rate, and domain.

Fields

  • discharge_voltage::Float64

  • cathode_potential::Float64

  • anode_Te::Float64

  • cathode_Te::Float64

  • wall_loss_model::HallThruster.WallLossModel

  • neutral_velocity::Float64

  • neutral_temperature::Float64

  • implicit_energy::Float64

  • propellant::HallThruster.Gas

  • ncharge::Int64

  • ion_temperature::Float64

  • anom_model::HallThruster.AnomalousTransportModel

  • conductivity_model::HallThruster.ThermalConductivityModel

  • ionization_model::HallThruster.IonizationModel

  • excitation_model::HallThruster.ExcitationModel

  • electron_neutral_model::HallThruster.ElectronNeutralModel

  • electron_ion_collisions::Bool

  • min_number_density::Float64

  • min_electron_temperature::Float64

  • transition_length::Float64

  • initial_condition::HallThruster.InitialCondition

  • magnetic_field_scale::Float64

  • source_neutrals::Any

  • source_ion_continuity::Any

  • source_ion_momentum::Any

  • source_potential::Any

  • source_energy::Any

  • scheme::HallThruster.HyperbolicScheme

  • thruster::HallThruster.Thruster

  • domain::Tuple{Float64, Float64}

  • LANDMARK::Bool

  • anode_mass_flow_rate::Float64

  • ion_wall_losses::Bool

  • background_pressure::Float64

  • background_neutral_temperature::Float64

  • neutral_ingestion_multiplier::Float64

  • anode_boundary_condition::Symbol

  • anom_smoothing_iters::Int64

  • solve_plume::Bool

  • apply_thrust_divergence_correction::Bool

  • electron_plume_loss_scale::Float64

  • reaction_rate_directories::Vector{String}

source
HallThruster.ExcitationLookupType
ExcitationLookup(;[directories::Vector{String} = String[]])

Default excitation model for HallThruster.jl. Reads excitation rate coefficients from file. Looks (preferentially) in provided directories and in the reactions subfolder for rate coefficient files

source
HallThruster.GasType
Gas

A chemical element in the gaseous state. Container for element properties used in fluid computations.

Fields

name::String Full name of gas (i.e. Xenon)

short_name::String Short name/symbol (i.e. Xe for Xenon)

γ::Float64 Specific heat ratio / adiabatic index

M::Float64 Molar mass (grams/mol) or atomic mass units

m::Float64 Mass of atom in kg

cp::Float64 Specific heat at constant pressure in J / kg / K

cv::Float64 Specific heat at constant volume in J / kg / K

R::Float64 Gas constant in J / kg / K

source
HallThruster.GasMethod
Gas(name::String, short_name::String; γ::Float64, M::Float64)

Instantiate a new Gas, providing a name, short name, the adiabatic index, and the molar mass. Other gas properties, including gas constant, specific heats at constant pressure/volume, and mass of atom/molecule in kg will are then computed.

julia> Gas("Xenon", "Xe", γ = 5/3, M = 83.798)
Xenon
source
HallThruster.Grid1DMethod
Grid1D(geometry, z_edge)

Given 1-D edge coordinates and thruster geometry, compute cell centers and cell volumes for grid

source
HallThruster.IonizationLookupType
IonizationLookup(;[directories::Vector{String} = String[]])

Default ionization model for HallThruster.jl. Reads ionization rate coefficients from file. Looks (preferentially) in provided directories and in the reactions subfolder for rate coefficient files

source
HallThruster.LandmarkExcitationLookupType
LandmarkExcitationLookup()

Excitation model for the LANDMARK benchmark.

Reads excitation rate coefficients from the landmark/landmark_rates.csv file in the HallThruster.jl main directory. Supports only singly-charged Xenon.

source
HallThruster.LandmarkIonizationLookupType
LandmarkIonizationLookup()

Ionization model for the LANDMARK benchmark.

Reads ionization rate coefficients from the landmark/landmark_rates.csv file in the HallThruster.jl main directory. Supports only singly-charged Xenon.

source
HallThruster.MitchnerType
Mitchner-Kruger conductivity model for partially-ionized plasmas
M. Mitchner and C. H. Kruger, Jr., Partially Ionized Gases (1973). Pg. 94
source
HallThruster.MultiLogBohmType
MultiLogBohm(zs, cs) <: AnomalousTransportModel

Model similar to that employed in Hall2De, where the mobility is Bohm-like (i.e. νan(z) = c(z) * ωce(z)) and z is in meters.

The function c(z) is defined by a sequence of nodes (z, c) provided by the user. At z = z[1], c(z) = c[1], and so forth.

At z[i] < z < z[i+1], log(c) is defined by linearly interpolating between log(c[i]) and log(c[i+1]).

For z < z[1], c = c[1] and for z > z[end], c(z) = c[end].

The user may also provide a single array of [z[1], z[2], ..., z[end], c[1], c[2], ..., c[end]]. The number of z values must be equal to the number of c values.

source
HallThruster.ShiftedGaussianBohmType
ShiftedGaussianBohm(trough_location, trough_width, trough_depth, z0, dz, alpha, pstar) <: AnomalousTransportModel

Model in which the anomalous collision frequency is Bohm-like (νan ~ ωce), except in a Gaussian-shaped region defined by the parameters troughlocation, troughwidth, troughmax, and troughmin where the collision frequency is lower. The location of the trough is based on the background pressure and the user-provided coefficients.

Arguments

  • trough_min: the minimum Hall parameter
  • trough_max: the maximum Hall parameter
  • trough_location: the axial position (in meters) of the mean of the Gaussian trough
  • trough_width: the standard deviation (in meters) of the Gaussian trough
  • z0: the furthest upstream displacement permitted at high back-pressures, relative to trough_location
  • dz: the maximum allowable amount of axial displacement
  • pstar: the background pressure at which the shift upstream halts/plateaus
  • alpha: the slope of the pressure response curve, with a higher value corresponding to a steeper pressure response
source
HallThruster.ShiftedMultiBohmType
ShiftedMultiBohm(zs, cs, z0, dz, alpha, pstar) <: AnomalousTransportModel

A version of MultiLogBohm where the coefficients are shifted depending on the background pressure.

source
HallThruster.ShiftedTwoZoneBohmType
ShiftedTwoZoneBohm(coeffs, z0, dz, alpha, pstar) <: AnomalousTransportModel

Model where the anomalous collision frequency has two values: c1 * ωce before some transition location and c2 * ωce after. Takes two arguments: c1 and c2. The transition between these values can be smoothed by the user-provided transition function. The location of the transition is based on the background pressure and the user-provided coefficients.

source
HallThruster.SpeciesType
Species

Represents a gas with a specific charge state. In a plasma, different ionization states of the same gas may coexist, so we need to be able to differentiate between these.

julia> Species(Xenon, 0)
Xe

julia> Species(Xenon, 1)
Xe+

julia> Species(Xenon, 3)
Xe3+
source
HallThruster.ThrusterType
Thruster

Struct containing information about a Hall thruster. This includes a name, geometry (a Geometry1D object), magnetic_field (radial magnetic field along centerline, a function which takes z in meters and outputs B in Tesla), and a shielded (a flag indicating whether the thruster is magnetically-shielded).

Fields

  • name::String

  • geometry::HallThruster.Geometry1D

  • magnetic_field::Any

  • shielded::Bool

source
HallThruster.TwoZoneBohmType
TwoZoneBohm(c1, c2) <: AnomalousTransportModel

Model where the anomalous collision frequency has two values: c1 * ωce inside the channel and c2 * ωce outside of the channel. Takes two arguments: c1 and c2. The transition between these values can be smoothed by the user-provided transition function.

source
HallThruster.WallLossModelType
abstract type WallLossModel

Abstract type for wall loss models in the electron energy equation. Types included with HallThruster.jl are:

  • NoWallLosses No electron energy losses to the wall.
  • ConstantSheathPotential(sheath_potential, inner_loss_coeff, outer_loss_coeff) Power to the walls scales with α * exp(-sheath_potential)), where α = inner_loss_coeff inside the channel and α = outer_loss_coeff outside.
  • WallSheath(material::WallMaterial) Power to the walls is computed self-consistently using a sheath model, dependent on the secondary electron emission yield of the provided material. See WallMaterial for material options.

Users implementing their own WallLossModel will need to implement at least three methods 1) freq_electron_wall(model, params, i): Compute the electron-wall momentum transfer collision frequency in cell i 2) wall_power_loss!(Q, model, params): Compute the electron power lost to the walls in array Q

A third method, wall_electron_current(model, params, i), will compute the electron current to the walls in cell i. If left unimplemented, it defaults to Ie,w = e ne νew Vcell where Vcell is the cell volume.

A fourth method, wall_ion_current(model, params, i, Z), for computing the current of ions of charge Z to the walls in cell i, may also be implemented. If left unimplemented, it will default to computing the current assuming Ie,w = Ii,w.

source
HallThruster.UnevenGridFunction
UnevenGrid(n, density = HallThruster.default_density)

Construct an unevenly-spaced grid according to provided density function. Defaults to twice as many grids inside of channel than outside. Provided density functions must have a signature of (z, z0, z1, Lch) where z is the location, (z0, z1) are the extents of the domain and Lch is the channel length

source
HallThruster.allocate_anom_variablesMethod
allocate_anom_variables(::AnomalousTransportModel, ncells)

Allocate arrays for anomalous transport state variables. ncells is the length of the arrays to be allocated. These anomalous transport variables are then stored in params.cache.anom_variables

source
HallThruster.backward_diff_coeffsMethod
backward_diff_coeffs(x0, x1, x2)

Generate finite difference coefficients for a backward first derivative approximation at the point x2 on a three-point stencil at points x0, x1, and x2

julia> backward_diff_coeffs(-2//1, -1//1, 0//1)
(1//2, -2//1, 3//2)
julia> backward_diff_coeffs(-3//2, -1//1, 0//1)
(4//3, -3//1, 5//3)
source
HallThruster.backward_differenceMethod
backward_difference(f0, f1, f2, x0, x1, x2)

Given three points x0, x1, and x2, and the function values at those points, f0, f1, f2, compute the second-order approximation of the derivative at x2

f(x) = x^4
x0, x1, x2 = 1.9999998, 1.9999999, 2
bd = backward_difference(f(x0), f(x1), f(x2), x0, x1, x2)
bd ≈ 32

# output

true
source
HallThruster.central_diff_coeffsMethod
central_diff_coeffs(x0, x1, x2)

Generate finite difference coefficients for a central first derivative approximation at the point x1 on a three-point stencil at points x0, x1, and x2

julia> central_diff_coeffs(-1//1, 0//1, 1//1)
(-1//2, 0//1, 1//2)
julia> central_diff_coeffs(-1//2, 0//1, 1//1)
(-4//3, 1//1, 1//3)
source
HallThruster.central_differenceMethod
central_difference(f0, f1, f2, x0, x1, x2)

Given three points x0, x1, and x2, and the function values at those points, f0, f1, f2, compute the second-order approximation of the derivative at x1

f(x) = x^4
x0, x1, x2 = 1.9999999, 2, 2.0000001
cd = central_difference(f(x0), f(x1), f(x2), x0, x1, x2)
cd ≈ 32

# output

true
source
HallThruster.channel_areaMethod
channel_area(outer_radius, inner_radius)

Compute the cross-sectional area of a Hall thruster channel from its dimensions

source
HallThruster.channel_perimeterMethod
channel_perimeter(outer_radius, inner_radius)

Compute the perimeteter of the thruster channel, equal to the sum of the inner and outer circumferences

source
HallThruster.coulomb_logarithmFunction
coulomb_logarithm(ne, Tev, Z = 1)

calculate coulomb logarithm for electron-ion collisions as a function of ion charge state Z, electron number density in m^-3, and electron temperature in eV.

source
HallThruster.downwind_diff_coeffsMethod
downwind_diff_coeffs(x0, x1, x2)

Generate finite difference coefficients for a downwind first derivative approximation at the point x2 on a three-point stencil at points x0, x1, and x2 (uses only points x1 and x2)

julia> downwind_diff_coeffs(-1//1, 0//1, 2//1)
(0//1, -1//2, 1//2)
source
HallThruster.electron_mobilityMethod
electron_mobility(νe, B)

calculates electron transport according to the generalized Ohm's law as a function of sum of the classical and anomalous collision frequencies and the magnetic field.

source
HallThruster.forward_diff_coeffsMethod
forward_diff_coeffs(x0, x1, x2)

Generate finite difference coefficients for a forward first derivative approximation at the point x0 on a three-point stencil at points x0, x1, and x2

julia> forward_diff_coeffs(1.0, 2.0, 3.0)
(-1.5, 2.0, -0.5)
julia> forward_diff_coeffs(0//1, 1//2, 3//2)
(-8//3, 3//1, -1//3)
source
HallThruster.forward_differenceMethod
forward_difference(f0, f1, f2, x0, x1, x2)

Given three points x0, x1, and x2, and the function values at those points, f0, f1, f2, compute the second-order approximation of the derivative at x0

f(x) = x^4
x0, x1, x2 = 2.0, 2.000001, 2.000002
fd = forward_difference(f(x0), f(x1), f(x2), x0, x1, x2)
fd ≈ 32

# output

true
source
HallThruster.generate_gridMethod
generate_grid(geometry, ncells)

Generate a one-dimensional uniform grid on the domain specified in the geomety. Returns number of cells, coordinates of cell centers (plus ghost cells face coordinates), interface/edges and volume of a cell for number density calculations.

source
HallThruster.interpolation_coeffsMethod
interpolation_coeffs(x, x0, x1, y0, y1)

Compute the coefficients for interpolation between two points (x0, y0) and (x1, y1) such that y = c0 * y0 + c1 * y1 ```jldoctest;setup = :(using HallThruster: itpcoeffs, lerp) julia> c0, c1 = interpolationcoeffs(0.5, 0.0, 1.0, 0.0, 2.0) (0.5, 0.5) julia> c0 * 0.0 + c1 * 2.0 == lerp(0.5, 0.0, 1.0, 0.0, 2.0) true

source
HallThruster.lerpMethod
lerp(x, x0, x1, y0, y1)

Interpolate between two points (x0, y0) and (x1, y1) ```jldoctest;setup = :(using HallThruster: lerp) julia> lerp(0.5, 0.0, 1.0, 0.0, 2.0) 1.0

source
HallThruster.load_reactionsMethod
load_reactions(model::ReactionModel, species)::Vector{IonizationReaction}

Load ionization reactions for the provided species and ionization model

source
HallThruster.maximum_charge_stateMethod
maximum_charge_state(model::ReactionModel)::Int

Return the maximum supported charge state for a given reaction model. If 0 is returned, then no charge state restriction is applied.

source
HallThruster.nodes_from_densityMethod
nodes_from_density(density, x0, x1, N)

Given bounds x0, x1, a number of points N, and a density function density(x), generate N nodes betweeen x0 and x1 spaced according to the provided desity function using inverse CDF transformation.

source
HallThruster.num_anom_variablesMethod
num_anom_variables(::AnomalousTransportModel)::Int

The number of variable arrays that should be allocated for the provided anomalous transport model. These arrays are used to save state beyond the anomalous collision frequency, and are useful for defining more complex anomalous transport models. If not defined by the user, this defaults to zero.

source
HallThruster.run_from_setupMethod
run_from_setup(U, params; duration, nsave, verbose)

Given a state vector U and params struct generated by setup_simulation, run for provided duration, saving nsave snapshots.

source
HallThruster.run_simulationMethod
run_simulation(config; duration, nsave, verbose, kwargs...)

Run a Hall thruster simulation using the provided Config object.

Arguments

  • config: a Config containing simulation parameters.
  • dt: The timestep, in seconds. Typical values are O(10 ns) (1e-8 seconds).
  • duration: How long to run the simulation, in seconds (simulation time, not wall time). Typical runtimes are O(1 ms) (1e-3 seconds).
  • ncells: How many cells to use. Typical values are 100 - 1000 cells.
  • nsave: How many frames to save.
source
HallThruster.second_deriv_central_diffMethod
second_deriv_central_diff(f0, f1, f2, x0, x1, x2)

Given three points x0, x1, and x2, and the function values at those points, f0, f1, f2, compute the second-order approximation of the second derivative at x1

f(x) = x^4
x0, x1, x2 = 1.9999, 2.0, 2.0001
sd = second_deriv_central_diff(f(x0), f(x1), f(x2), x0, x1, x2)
sd ≈ 48

# output

true
source
HallThruster.second_deriv_coeffsMethod
second_deriv_coeffs(x0, x1, x2)

Generate finite difference coefficients for a central second derivative approximation at the point x1 on a three-point stencil at points x0, x1, and x2

julia> second_deriv_coeffs(-2//1, 0//1, 2//1)
(1//4, -1//2, 1//4)
julia> second_deriv_coeffs(-1//2, 0//1, 1//1)
(8//3, -4//1, 4//3)
source
HallThruster.sheath_potentialMethod
sheath_potential(Tev, γ, mi))

compute wall sheath to be used for radiative losses and loss to wall. Goebel Katz equ. 7.3-29, 7.3-44. Assumed nₑuₑ/nᵢuᵢ ≈ 0.5 Sheath potentials are positive by convention in HallThruster.jl.

source
HallThruster.supported_gasesMethod
supported_gases(model::ReactionModel)::Vector{HallThruster.Gas}

Check which gases are supported by a given reaction model. If an empty vector is provided, then there are no restrictions on what gases can be used.

source
HallThruster.time_averageFunction
time_average(sol, tstampstart)

compute time-averaged solution, input Solution type and the frame at which averaging starts. Returns a Solution object with a single frame.

source
HallThruster.upwind_diff_coeffsMethod
upwind_diff_coeffs(x0, x1, x2)

Generate finite difference coefficients for a upwind first derivative approximation at the point x2 on a three-point stencil at points x0, x1, and x2 (uses only points x0 and x1)

julia> upwind_diff_coeffs(-3//1, 0//1, 2//1)
(-1//3, 1//3, 0//1)
source
HallThruster.write_restartMethod
write_restart(path::AbstractString, sol)

Write a JLD2 restart file to path`.

This can be reloaded to resume a simulation.

source