{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "5ddab7ea",
   "metadata": {},
   "source": [
    "# Ionization states in an interplanetary coronal mass ejection"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "31f467b1",
   "metadata": {},
   "source": [
    "[plasmapy.particles]: https://docs.plasmapy.org/en/stable/particles/index.html\n",
    "\n",
    "The ionization state distribution for an element refers to the fractions of that element at each ionic level. For example, the ionization state of helium in the solar wind might be 10% He$^{0+}$, 70% He$^{1+}$, and 20% He$^{2+}$. This notebook introduces the data structures in [plasmapy.particles] for representing the ionization state of a plasma."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6e9deb4a",
   "metadata": {},
   "outputs": [],
   "source": [
    "import astropy.units as u\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "from plasmapy.particles import IonizationState, IonizationStateCollection"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "279b1b3a",
   "metadata": {},
   "source": [
    "## The ionization state of a single element"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "39966b76",
   "metadata": {},
   "source": [
    "[IonizationState]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state.IonizationState.html#plasmapy.particles.ionization_state.IonizationState\n",
    "\n",
    "Let's create an [IonizationState] object for helium using the ionic fractions described above. We'll specify the total number density of helium via the `n_elem` keyword argument."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c4419185",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states = IonizationState(\"He-4\", [0.1, 0.7, 0.2], n_elem=1e13 * u.m**-3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b28f004a",
   "metadata": {},
   "source": [
    "[ionic_fractions]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state.IonizationState.html#plasmapy.particles.ionization_state.IonizationState.ionic_fractions\n",
    "[IonizationState]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state.IonizationState.html#plasmapy.particles.ionization_state.IonizationState\n",
    "\n",
    "The ionization state distribution is stored in the [ionic_fractions] attribute of [IonizationState]."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "06f5eaa6",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.ionic_fractions"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9d12a5a7",
   "metadata": {},
   "source": [
    "We can get the symbols for each ionic level in `He_states` too."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9ff206c4",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.ionic_symbols"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "938e60c9",
   "metadata": {},
   "source": [
    "Because we provided the number density of the element as a whole, we can get back the number density of each ionic level."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d7edcf41",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.number_densities"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7cec7b16",
   "metadata": {},
   "source": [
    "We can also get the electron number density required to balance the positive charges for ions of this element."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "22744f13",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.n_e"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3c7ff668",
   "metadata": {},
   "source": [
    "[IonicLevel]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state.IonicLevel.html\n",
    "[IonizationState]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state.IonizationState.html#plasmapy.particles.ionization_state.IonizationState\n",
    "\n",
    "We can provide an [IonizationState] with a charge number as an index to get an [IonicLevel] object that contains most of these attributes, but for a single ionic level (like He$^{1+}$).  This capability is useful if we wish to iterate over the ions of an element."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2c67356d",
   "metadata": {},
   "outputs": [],
   "source": [
    "for Z in range(3):\n",
    "    print(He_states[Z])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3c7995f0",
   "metadata": {},
   "source": [
    "We can get information about the average charge state via `Z_mean`, `Z_most_abundant`, and `Z_rms`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "efc11927",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.Z_mean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "edb4b146",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.Z_most_abundant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7fbd84bf",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.Z_rms"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ccda6535",
   "metadata": {},
   "source": [
    "We can calculate the properties of the average ionic level."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "847e998f",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.average_ion()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a70042f7",
   "metadata": {},
   "source": [
    "We can use the `summarize()` method to get information about the ionization state."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19a9ac38",
   "metadata": {},
   "outputs": [],
   "source": [
    "He_states.summarize()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6ac6a936",
   "metadata": {},
   "source": [
    "## Ionization states of multiple elements"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3f7c7d6b",
   "metadata": {},
   "source": [
    "[*Advanced Composition Explorer*]: https://en.wikipedia.org/wiki/Advanced_Composition_Explorer\n",
    "[Gilbert et al. (2012)]: https://doi.org/10.1088/0004-637X/751/1/20\n",
    "\n",
    "Now let's look at some actual average hourly densities for ions of C, O, and Fe during an interplanetary coronal mass ejection (ICME) observed by the [*Advanced Composition Explorer*] (*ACE*) near 1 AU.  The data were estimated from Figure 4 in [Gilbert et al. (2012)].  This data set is noteworthy because there is information from very low charge states to very high charge states for several elements."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1dea5963",
   "metadata": {},
   "source": [
    "[electron ionization]: https://en.wikipedia.org/wiki/Electron_ionization\n",
    "[radiative recombination]: https://en.wikipedia.org/wiki/Plasma_recombination\n",
    "\n",
    "The [electron ionization] and [radiative recombination] rates for a given temperature are proportional to $n_i n_e$, where $n_i$ is the ion number density and $n_e$ is the electron number density.  Ionization and recombination are fast at high densities and slow at low densities. High density plasma can reach *ionization equilibrium* (when the recombination and ionization rates balance each other out) more quickly than low density plasma.  Plasma undergoing rapid heating or cooling will be in *non-equilibrium ionization* because ionization and recombination can't keep up with the temperature changes.\n",
    "\n",
    "Quiescent plasma near the sun is typically close to ionization equilibrium because the density is high. As plasma moves away from the sun, the ionization and recombination rates drop rapidly because of the decreasing number density. The ionization states freeze out at several solar radii as the ionization and recombination time scales begin to exceed the time it takes for plasma to move from the sun to 1 AU. \n",
    "\n",
    "The ionization states of the solar wind are powerful diagnostics of the thermodynamic history of solar wind and ICME plasma.  The ionization states observed by *ACE* at 1 AU are essentially the same as at ∼$5R_☉$.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6f57b224",
   "metadata": {},
   "outputs": [],
   "source": [
    "number_densities = {\n",
    "    \"C\": [0, 5.7e-7, 4.3e-5, 3.6e-6, 2.35e-6, 1e-6, 1.29e-6] * u.cm**-3,\n",
    "    \"O\": [0, 1.2e-7, 2.2e-4, 7.8e-6, 8.8e-7, 1e-6, 4e-6, 1.3e-6, 1.2e-7] * u.cm**-3,\n",
    "    \"Fe\": [\n",
    "        0,\n",
    "        0,\n",
    "        1.4e-8,\n",
    "        1.1e-7,\n",
    "        2.5e-7,\n",
    "        2.2e-7,\n",
    "        1.4e-7,\n",
    "        1.2e-7,\n",
    "        2.1e-7,\n",
    "        2.1e-7,\n",
    "        1.6e-7,\n",
    "        8e-8,\n",
    "        6.3e-8,\n",
    "        4.2e-8,\n",
    "        2.5e-8,\n",
    "        2.3e-8,\n",
    "        1.5e-8,\n",
    "        3.1e-8,\n",
    "        6.1e-9,\n",
    "        2.3e-9,\n",
    "        5.3e-10,\n",
    "        2.3e-10,\n",
    "        0,\n",
    "        0,\n",
    "        0,\n",
    "        0,\n",
    "        0,\n",
    "    ]\n",
    "    * u.cm**-3,\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7b2d878e",
   "metadata": {},
   "source": [
    "[IonizationStateCollection]: https://docs.plasmapy.org/en/stable/api/plasmapy.particles.ionization_state_collection.IonizationStateCollection.html#plasmapy.particles.ionization_state_collection.IonizationStateCollection\n",
    "\n",
    "Let's use this information as an input for [IonizationStateCollection]: a data structure for the ionization states of multiple elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ed4da71a",
   "metadata": {},
   "outputs": [],
   "source": [
    "states = IonizationStateCollection(number_densities)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6e769fc0",
   "metadata": {},
   "source": [
    "We can index this to get an `IonizationState` for one of the elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ebce50d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "states[\"C\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ac7be92",
   "metadata": {},
   "source": [
    "We can get the relative abundances of each of the elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4b302f1e",
   "metadata": {},
   "outputs": [],
   "source": [
    "states.abundances"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5a7baf6c",
   "metadata": {},
   "outputs": [],
   "source": [
    "states.log_abundances"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1350ddb8",
   "metadata": {},
   "source": [
    "We can get the number densities as a `dict` (like what we provided) and the electron number density (assuming quasineutrality, but only for the elements contained in the data structure)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "566fab0c",
   "metadata": {},
   "outputs": [],
   "source": [
    "states.number_densities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "402686e8",
   "metadata": {},
   "outputs": [],
   "source": [
    "states.n_e"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cabc5516",
   "metadata": {},
   "source": [
    "We can summarize this information too, but let's specify the minimum ionic fraction to print."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "48e6f4ec",
   "metadata": {},
   "outputs": [],
   "source": [
    "states.summarize(minimum_ionic_fraction=0.02)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3d5a1cff",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, figsize=(10, 3), tight_layout=True)\n",
    "\n",
    "for state, ax in zip(states, axes, strict=False):\n",
    "    ax.bar(state.charge_numbers, state.ionic_fractions)\n",
    "    ax.set_title(f\"Ionization state for {state.base_particle}\")\n",
    "    ax.set_xlabel(\"ionic level\")\n",
    "    ax.set_ylabel(\"log$_{10}$ of ionic fraction\")\n",
    "    ax.set_yscale(\"log\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a282993c",
   "metadata": {},
   "source": [
    "[filament]: https://en.wikipedia.org/wiki/Solar_prominence\n",
    "\n",
    "The wide range of average ionic levels for each element is strong evidence that the plasma observed by *ACE* originated from a wide range of temperatures.  The lowest charge states are evidence of cool [filament] plasma while the high charge states are evidence of rapidly heated plasma."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}