trident.add_ion_fields

trident.add_ion_fields(ds, ions, ftype='gas', ionization_table=None, field_suffix=False, line_database=None, sampling_type='auto', particle_type=None)[source]

Preferred method for adding ion fields to a yt dataset.

Select ions based on the selection indexing set up in parse_subset_to_ions function, that is, by specifying a list of strings where each string represents an ion or line. Strings are of one of three forms:

  • <element>
  • <element> <ion state>
  • <element> <ion state> <line_wavelength>

If a line_database is selected, then the ions chosen will be a subset of the ions present in the equivalent LineDatabase, nominally located in trident.__path__/data/line_lists.

For each ion species selected, four fields will be added (example for Mg II):

  • Ion fraction field. e.g. (ftype, ‘Mg_p1_ion_fraction’)
  • Number density field. e.g. (ftype, ‘Mg_p1_number_density’)
  • Density field. e.g. (ftype, ‘Mg_p1_density’)
  • Mass field. e.g. (ftype, ‘Mg_p1_mass’)

This function is the preferred method for adding ion fields to one’s dataset, but for more fine-grained control, one can also employ the add_ion_fraction_field, add_ion_number_density_field, add_ion_density_field, add_ion_mass_field functions individually.

Fields are added assuming collisional ionization equilibrium and photoionization in the optically thin limit from a redshift-dependent metagalactic ionizing background using the ionization_table specified.

WARNING: The “ftype” must match the field type that you’re using for the field interpolation. So for particle-based codes, this must be the ftype of the gas particles (e.g., PartType0, Gas). Using the default of gas in this instance will interpolate on the grid-based fields, which will give the wrong answers for particle-based codes, since the ion field interpolation will take place on the already deposited grid-based fields.

Parameters

Ds:

yt dataset object

This is the dataset to which the ion fraction field will be added.

Ions:

list of strings

List of strings matching possible lines. Strings can be of the form: * Atom - Examples: “H”, “C”, “Mg” * Ion - Examples: “H I”, “H II”, “C IV”, “Mg II” * Line - Examples: “H I 1216”, “C II 1336”, “Mg II 1240”

If set to ‘all’, creates all ions for the first 30 elements: (ie hydrogen to zinc). If set to ‘all’ with line_database keyword set, then creates all ions associated with the lines specified in the equivalent LineDatabase.

Ftype:

string, optional

The field type of the field to add. it is the first string in the field tuple e.g. “gas” in (“gas”, “O_p5_ion_fraction”) ftype must correspond to the ftype of the ‘density’, and ‘temperature’ fields in your dataset you wish to use to generate the ion field. Default: “gas”

Ionization_table:
 

string, optional

Path to an appropriately formatted HDF5 table that can be used to compute the ion fraction as a function of density, temperature, metallicity, and redshift. When set to None, it uses the table specified in ~/.trident/config Default: None

Field_suffix:

boolean, optional

Determines whether or not to append a suffix to the field name that indicates what ionization table was used. Useful when using generating ion_fields that already exist in a dataset.

Line_database:

string, optional

Ions are selected out of the set of ions present in the line_database constructed from the line list filename specified here. See LineDatabase for more information.

Sampling_type:

string, optional

Set to ‘particle’ if the field should be for particles. Set to ‘cell’ if the field should be for grids/cells. Set to ‘auto’ for this to be determined automatically. Default: ‘auto’

Particle_type:

boolean, optional

This is deprecated in favor of ‘sampling_type’. Set to True if you are adding ion fields to particles, as specified by the ‘ftype’. Set to False if you are not. Set to ‘auto’, if you want the code to autodetermine if the field specified by the ‘ftype’ is particle or not. Default: ‘auto’

Example

To add ionized hydrogen, doubly-ionized Carbon, and all of the Magnesium species fields to a dataset, you would run:

>>> import yt
>>> import trident
>>> ds = yt.load('path/to/file')
>>> trident.add_ion_fields(ds, ions=['H II', 'C III', 'Mg'])