HallThruster.Air
HallThruster.Argon
HallThruster.Bismuth
HallThruster.Electron
HallThruster.Krypton
HallThruster.LOOKUP_ZS
HallThruster.Mercury
HallThruster.NA
HallThruster.R0
HallThruster.Xenon
HallThruster.e
HallThruster.kB
HallThruster.me
HallThruster.AnomalousTransportModel
HallThruster.Bohm
HallThruster.Braginskii
HallThruster.Config
HallThruster.ExcitationLookup
HallThruster.GKElectronNeutral
HallThruster.Gas
HallThruster.Gas
HallThruster.Grid1D
HallThruster.IonizationLookup
HallThruster.LANDMARK_conductivity
HallThruster.LandmarkExcitationLookup
HallThruster.LandmarkIonizationLookup
HallThruster.Mitchner
HallThruster.MultiLogBohm
HallThruster.NoAnom
HallThruster.NoExcitation
HallThruster.ShiftedGaussianBohm
HallThruster.ShiftedMultiBohm
HallThruster.ShiftedTwoZoneBohm
HallThruster.Species
HallThruster.ThermalConductivityModel
HallThruster.Thruster
HallThruster.TwoZoneBohm
HallThruster.WallLossModel
HallThruster.EvenGrid
HallThruster.UnevenGrid
HallThruster.allocate_anom_variables
HallThruster.backward_diff_coeffs
HallThruster.backward_difference
HallThruster.central_diff_coeffs
HallThruster.central_difference
HallThruster.channel_area
HallThruster.channel_perimeter
HallThruster.channel_width
HallThruster.compute_current
HallThruster.coulomb_logarithm
HallThruster.downwind_diff_coeffs
HallThruster.electron_mobility
HallThruster.forward_diff_coeffs
HallThruster.forward_difference
HallThruster.freq_electron_electron
HallThruster.freq_electron_ion
HallThruster.freq_electron_neutral
HallThruster.generate_grid
HallThruster.interpolation_coeffs
HallThruster.lerp
HallThruster.load_reactions
HallThruster.load_reactions
HallThruster.maximum_charge_state
HallThruster.nodes_from_density
HallThruster.num_anom_variables
HallThruster.rate_coeff
HallThruster.read_restart
HallThruster.run_from_setup
HallThruster.run_simulation
HallThruster.second_deriv_central_diff
HallThruster.second_deriv_coeffs
HallThruster.sheath_potential
HallThruster.supported_gases
HallThruster.time_average
HallThruster.upwind_diff_coeffs
HallThruster.write_restart
HallThruster.Air
— ConstantAir::Gas
Earth air at standard temperature and pressure
HallThruster.Argon
— ConstantArgon::Gas
Argon gas
HallThruster.Bismuth
— ConstantBismuth::Gas
Bismuth vapor
HallThruster.Electron
— ConstantElectron::Species
Electron
HallThruster.Krypton
— ConstantKrypton::Gas
Krypton gas
HallThruster.LOOKUP_ZS
— ConstantLookup for thermal conductivity coefficients, from Table 1 of S. I. Braginskii, in Reviews of Plasma Physics (1965), Vol. 1, p. 205.
HallThruster.Mercury
— ConstantMercury::Gas
Mercury vapor
HallThruster.NA
— ConstantNA
Number of atoms in a kg-mol (6.02214076e26 / kmol)
HallThruster.R0
— ConstantR0
Universal gas constant (8314.46261815324 J / kmol K)
HallThruster.Xenon
— ConstantXenon::Gas
Xenon gas
HallThruster.e
— Constante
Electron charge (1.602176634e-19 Coulomb)
HallThruster.kB
— ConstantkB
Boltzmann constant (1.380649e-23 J/K)
HallThruster.me
— Constantme
Electron mass (9.10938356e-31 kilograms)
HallThruster.AnomalousTransportModel
— TypeAnomalousTransportModel
The abstract supertype of all types of anomalous transport models. Subtype this to define your own model.
HallThruster.Bohm
— TypeBohm(c) <: AnomalousTransportModel
Model where the anomalous collision frequency scales with the electron cyclotron frequency ωce times some scaling factor c
HallThruster.Braginskii
— TypeBraginskii conductivity model for fully-ionized plasmasa
S. I. Braginskii, in Reviews of Plasma Physics (1965), Vol. 1, p. 205.
HallThruster.Config
— Typestruct 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}
HallThruster.ExcitationLookup
— TypeExcitationLookup(;[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
HallThruster.GKElectronNeutral
— TypeElectron neutral collision model as defined by Eq. 3.6-13, from Fundamentals of Electric Propulsion, Goebel and Katz, 2008.
HallThruster.Gas
— TypeGas
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
HallThruster.Gas
— MethodGas(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
HallThruster.Grid1D
— MethodGrid1D(geometry, z_edge)
Given 1-D edge coordinates and thruster geometry, compute cell centers and cell volumes for grid
HallThruster.IonizationLookup
— TypeIonizationLookup(;[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
HallThruster.LANDMARK_conductivity
— TypeLANDMARK, uses 10/9 μnϵ
HallThruster.LandmarkExcitationLookup
— TypeLandmarkExcitationLookup()
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.
HallThruster.LandmarkIonizationLookup
— TypeLandmarkIonizationLookup()
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.
HallThruster.Mitchner
— TypeMitchner-Kruger conductivity model for partially-ionized plasmas
M. Mitchner and C. H. Kruger, Jr., Partially Ionized Gases (1973). Pg. 94
HallThruster.MultiLogBohm
— TypeMultiLogBohm(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.
HallThruster.NoAnom
— TypeNoAnom <: AnomalousTransportModel
No anomalous collision frequency included in simulation
HallThruster.NoExcitation
— TypeNoExcitation <: ExcitationModel
Model for neglecting excitation energy losses
HallThruster.ShiftedGaussianBohm
— TypeShiftedGaussianBohm(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 parametertrough_max
: the maximum Hall parametertrough_location
: the axial position (in meters) of the mean of the Gaussian troughtrough_width
: the standard deviation (in meters) of the Gaussian troughz0
: the furthest upstream displacement permitted at high back-pressures, relative totrough_location
dz
: the maximum allowable amount of axial displacementpstar
: the background pressure at which the shift upstream halts/plateausalpha
: the slope of the pressure response curve, with a higher value corresponding to a steeper pressure response
HallThruster.ShiftedMultiBohm
— TypeShiftedMultiBohm(zs, cs, z0, dz, alpha, pstar) <: AnomalousTransportModel
A version of MultiLogBohm where the coefficients are shifted depending on the background pressure.
HallThruster.ShiftedTwoZoneBohm
— TypeShiftedTwoZoneBohm(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.
HallThruster.Species
— TypeSpecies
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+
HallThruster.ThermalConductivityModel
— TypeThermal_Conductivity
The abstract supertype of all types of thermal conductivity models. Subtype this to define your own model.
HallThruster.Thruster
— TypeThruster
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
HallThruster.TwoZoneBohm
— TypeTwoZoneBohm(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.
HallThruster.WallLossModel
— Typeabstract 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. SeeWallMaterial
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.
HallThruster.EvenGrid
— MethodEvenGrid(n)
Construct an evenly-spaced grid with n cells.
HallThruster.UnevenGrid
— FunctionUnevenGrid(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
HallThruster.allocate_anom_variables
— Methodallocate_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
HallThruster.backward_diff_coeffs
— Methodbackward_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)
HallThruster.backward_difference
— Methodbackward_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
HallThruster.central_diff_coeffs
— Methodcentral_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)
HallThruster.central_difference
— Methodcentral_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
HallThruster.channel_area
— Methodchannel_area(outer_radius, inner_radius)
Compute the cross-sectional area of a Hall thruster channel from its dimensions
HallThruster.channel_perimeter
— Methodchannel_perimeter(outer_radius, inner_radius)
Compute the perimeteter of the thruster channel, equal to the sum of the inner and outer circumferences
HallThruster.channel_width
— Methodchannel_width(outer_radius, inner_radius)
Compute the thruster channel width
HallThruster.compute_current
— Functioncompute_current(sol, location)
compute current at anode or cathode = outflow in 1D code.
HallThruster.coulomb_logarithm
— Functioncoulomb_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.
HallThruster.downwind_diff_coeffs
— Methoddownwind_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)
HallThruster.electron_mobility
— Methodelectron_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.
HallThruster.forward_diff_coeffs
— Methodforward_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)
HallThruster.forward_difference
— Methodforward_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
HallThruster.freq_electron_electron
— Methodfreq_electron_electron(ne, Tev)
Effective frequency at which electrons are scattered due to collisions with other electrons
HallThruster.freq_electron_ion
— Methodfreq_electron_ion(ne, Tev, Z)
Effective frequency at which electrons are scattered due to collisions with ions
HallThruster.freq_electron_neutral
— Methodfreq_electron_neutral(model::ElectronNeutralModel, nn, Tev)
Effective frequency of electron scattering caused by collisions with neutrals
HallThruster.generate_grid
— Methodgenerate_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.
HallThruster.interpolation_coeffs
— Methodinterpolation_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
HallThruster.lerp
— Methodlerp(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
HallThruster.load_reactions
— MethodCharge exchange fits from Hause, Prince and Bemish, 2013
HallThruster.load_reactions
— Methodload_reactions(model::ReactionModel, species)::Vector{IonizationReaction}
Load ionization reactions for the provided species and ionization model
HallThruster.maximum_charge_state
— Methodmaximum_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.
HallThruster.nodes_from_density
— Methodnodes_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.
HallThruster.num_anom_variables
— Methodnum_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.
HallThruster.rate_coeff
— MethodBy default, rate_coeff looks for a lookup table stored in the reaction struct
HallThruster.read_restart
— Methodread_restart(path::AbstractString)
Load a JLD2 restart file from path
.
HallThruster.run_from_setup
— Methodrun_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.
HallThruster.run_simulation
— Methodrun_simulation(config; duration, nsave, verbose, kwargs...)
Run a Hall thruster simulation using the provided Config object.
Arguments
config
: aConfig
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.
HallThruster.second_deriv_central_diff
— Methodsecond_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
HallThruster.second_deriv_coeffs
— Methodsecond_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)
HallThruster.sheath_potential
— Methodsheath_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.
HallThruster.supported_gases
— Methodsupported_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.
HallThruster.time_average
— Functiontime_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.
HallThruster.upwind_diff_coeffs
— Methodupwind_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)
HallThruster.write_restart
— Methodwrite_restart(path::AbstractString, sol)
Write a JLD2 restart file to path
`.
This can be reloaded to resume a simulation.