NEI
- class plasmapy_nei.nei.NEI(inputs, abundances: ~typing.Optional[~typing.Union[~typing.Dict, str]] = None, T_e: ~typing.Optional[~typing.Union[~typing.Callable, ~astropy.units.quantity.Quantity]] = None, n: ~typing.Optional[~typing.Union[~typing.Callable, ~astropy.units.quantity.Quantity]] = None, time_input: ~typing.Optional[~astropy.units.quantity.Quantity] = None, time_start: ~typing.Optional[~astropy.units.quantity.Quantity] = None, time_max: ~typing.Optional[~astropy.units.quantity.Quantity] = None, max_steps: ~typing.Union[int, ~numpy.integer] = 10000, tol: ~typing.Union[int, float] = 1e-15, dt: ~typing.Optional[~astropy.units.quantity.Quantity] = None, dt_max: ~astropy.units.quantity.Quantity = <Quantity inf s>, dt_min: ~astropy.units.quantity.Quantity = <Quantity 0. s>, adapt_dt: ~typing.Optional[bool] = None, safety_factor: ~typing.Union[int, float] = 1, verbose: bool = False)
Bases:
object
Perform and analyze a non-equilibrium ionization simulation.
- Parameters
inputs –
T_e (astropy.units.Quantity or callable) – The electron temperature, which may be a constant, an array of temperatures corresponding to the times in
time_input
, or a function that yields the temperature as a function of time.n (astropy.units.Quantity or callable) – The number density multiplicative factor. The number density of each element will be
n
times the abundance given inabundances
. For example, ifabundance['H'] = 1
, then this will correspond to the number density of hydrogen (including neutral hydrogen and protons). This factor may be a constant, an array of number densities over time, or a function that yields a number density as a function of time.time_input (astropy.units.Quantity, optional) – An array containing the times associated with
n
andT_e
in units of time.time_start (astropy.units.Quantity, optional) – The start time for the simulation. If density and/or temperature are given by arrays, then this argument must be greater than
time_input[0]
. If this argument is not supplied, thentime_start
defaults totime_input[0]
(if given) and zero seconds otherwise.time_max (astropy.units.Quantity) – The maximum time for the simulation. If density and/or temperature are given by arrays, then this argument must be less than
time_input[-1]
.max_steps (
int
) – The maximum number of time steps to be taken during a simulation.dt (astropy.units.Quantity) – The time step. If
adapt_dt
isFalse
, thendt
is the time step for the whole simulation.dt_max (astropy.units.Quantity) – The maximum time step to be used with an adaptive time step.
dt_min (astropy.units.Quantity) – The minimum time step to be used with an adaptive time step.
adapt_dt (
bool
) – IfTrue
, change the time step based on the characteristic ionization and recombination time scales and change in temperature. Not yet implemented.safety_factor (
float
orint
) – A multiplicative factor to multiply by the time step whenadapt_dt
isTrue
. Lower values improve accuracy, whereas higher values reduce computational time. Not yet implemented.tol (float) – The absolute tolerance to be used in comparing ionic fractions.
verbose (bool, optional) – A flag stating whether or not to print out information for every time step. Setting
verbose
toTrue
is useful for testing. Defaults toFalse
.abundances (dict) –
Examples
>>> import numpy as np >>> import astropy.units as u
>>> inputs = {'H': [0.9, 0.1], 'He': [0.9, 0.099, 0.001]} >>> abund = {'H': 1, 'He': 0.085} >>> n = u.Quantity([1e9, 1e8], u.cm**-3) >>> T_e = np.array([10000, 40000]) * u.K >>> time = np.array([0, 300]) * u.s >>> dt = 0.25 * u.s
The initial conditions can be accessed using the initial attribute.
>>> sim = NEI(inputs=inputs, abundances=abund, n=n, T_e=T_e, time_input=time, adapt_dt=False, dt=dt)
After having inputted all of the necessary information, we can run the simulation.
>>> results = sim.simulate()
The initial results are stored in the
initial
attribute.>>> sim.initial.ionic_fractions['H'] array([0.9, 0.1])
The final results can be access with the
final
attribute.>>> sim.final.ionic_fractions['H'] array([0.16665179, 0.83334821]) >>> sim.final.ionic_fractions['He'] array([0.88685261, 0.11218358, 0.00096381]) >>> sim.final.T_e <Quantity 40000. K>
Both
initial
andfinal
are instances of theIonizationStateCollection
class.Notes
The ionization and recombination rates are from Chianti version 8.7. These rates include radiative and dielectronic recombination. Photoionization is not included.
Attributes Summary
The temperature input.
Return the abundances.
Return
True
if the time step is set to be adaptive,False
if the time step is set to not be adapted, andNone
if this attribute was not set.Return the inputted time step.
The minimum time step.
Return a
dict
containingEigenData
instances for each element.A
list
of the elements.Return the ionization states of the plasma at the end of the simulation.
Return the ionization states of the plasma at the beginning of the simulation.
The maximum number of steps that a simulation will be allowed to take.
The number density factor input.
Return the
SimulationResults
class instance that corresponds to the simulation results.The multiplicative factor that the time step is to be multiplied by when using an adaptive time step.
The maximum time allowed for the simulation.
The start time of the simulation.
The tolerance for comparisons between different ionization states.
Return
True
if verbose output during a simulation is requested, andFalse
otherwise.Methods Summary
electron_temperature
(time)equil_ionic_fractions
([T_e, time])Return the equilibrium ionic fractions for a temperature or at a given time.
hydrogen_number_density
(time)in_time_interval
(time[, buffer])Return
True
if thetime
is betweentime_start - buffer
andtime_max + buffer
, andFalse
otherwise.index_to_time
(index)Returns the time value or array given the index/indices
save
([filename])Save the
NEI
instance to an HDF5 file.set_timestep
([dt])Set the time step for the next non-equilibrium ionization time advance.
simulate
()Perform a non-equilibrium ionization simulation.
Advance the simulation by one time step.
time_to_index
(time)Returns the closest index value or array for the given time(s)
Attributes Documentation
- T_e_input
The temperature input.
- abundances
Return the abundances.
- adapt_dt
Return
True
if the time step is set to be adaptive,False
if the time step is set to not be adapted, andNone
if this attribute was not set.
- dt_input
Return the inputted time step.
- dt_max
- dt_min
The minimum time step.
- final
Return the ionization states of the plasma at the end of the simulation.
- initial
Return the ionization states of the plasma at the beginning of the simulation.
- max_steps
The maximum number of steps that a simulation will be allowed to take.
- n_input
The number density factor input.
- results
Return the
SimulationResults
class instance that corresponds to the simulation results.
- safety_factor
The multiplicative factor that the time step is to be multiplied by when using an adaptive time step.
- time_input
- time_max
The maximum time allowed for the simulation.
- time_start
The start time of the simulation.
- tol
The tolerance for comparisons between different ionization states.
Methods Documentation
-
electron_temperature(time:
astropy.units.Quantity
)astropy.units.Quantity
-
equil_ionic_fractions(T_e: Optional[
astropy.units.Quantity
] = None, time: Optional[astropy.units.Quantity
] = None) Dict[str, ndarray] Return the equilibrium ionic fractions for a temperature or at a given time.
- Parameters
T_e (astropy.units.Quantity, optional) – The electron temperature in units that can be converted to kelvin.
time (astropy.units.Quantity, optional) – The time in units that can be converted to seconds.
- Returns
equil_ionfracs – The equilibrium ionic fractions for the elements contained within this class
- Return type
Notes
Only one of
T_e
andtime
may be included as an argument. If neitherT_e
ortime
is provided and the temperature for the simulation is given by a constant, the this method will assume thatT_e
is the temperature of the simulation.
-
hydrogen_number_density(time:
astropy.units.Quantity
)astropy.units.Quantity
- in_time_interval(time: Unit("s"), buffer: Unit("s") = <Quantity 1.e-09 s>)
Return
True
if thetime
is betweentime_start - buffer
andtime_max + buffer
, andFalse
otherwise.- Raises
TypeError – If
time
orbuffer
is not aastropy.units.Quantity
astropy.units.UnitsError – If
time
orbuffer
is not in units of time.
- index_to_time(index)
Returns the time value or array given the index/indices
- Parameters
index (array-like) – A value or array of values representing the index of the time array created by the simulation
- Returns
get_time – The time value associated with index input(s)
- Return type
-
set_timestep(dt: Optional[
astropy.units.Quantity
] = None) Set the time step for the next non-equilibrium ionization time advance.
- Parameters
dt (astropy.units.Quantity, optional) – The time step to be used for the next time advance.
Notes
If
dt
is notNone
, then the time step will be set todt
.If
dt
is not set and theadapt_dt
attribute of anNEI
instance isTrue
, then this method will calculate the time step corresponding to how long it will be until the temperature rises or drops into the next temperature bin. If this time step is betweendtmin
anddtmax
, thenIf
dt
is not set and theadapt_dt
attribute isFalse
, then this method will set the time step as what was inputted to theNEI
class upon instantiation in thedt
argument or through theNEI
class’sdt_input
attribute.- Raises
NEIError – If the time step cannot be set, for example if the
dt
argument is invalid or the time step cannot be adapted.
- simulate() SimulationResults
Perform a non-equilibrium ionization simulation.
- Returns
results – The results from the simulation (which are also stored in the
results
attribute of theNEI
instance this method was called from.- Return type
Simulation
- time_advance()
Advance the simulation by one time step.
- time_to_index(time)
Returns the closest index value or array for the given time(s)
- Parameters
time (array-like,) – A value or array of values representing the values of the time array created by the simulation
- Returns
index – The index value associated with the time input(s)
- Return type
int or array-like,