From 5a6f5c85cdc8cfb5d09705744b7df5e80ff4fb35 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Mon, 16 Sep 2024 12:00:02 -0400 Subject: [PATCH 1/6] refactor result return type --- docs/reference.rst | 10 ++ ebtelplusplus/__init__.py | 4 +- ebtelplusplus/high_level.py | 161 ++++++++++++++++++------ ebtelplusplus/models.py | 27 ++-- ebtelplusplus/tests/test_compare_idl.py | 18 +-- ebtelplusplus/tests/test_solver.py | 1 + 6 files changed, 153 insertions(+), 68 deletions(-) diff --git a/docs/reference.rst b/docs/reference.rst index dc46966..49123c0 100644 --- a/docs/reference.rst +++ b/docs/reference.rst @@ -3,6 +3,16 @@ API Reference ============= +ebtelplusplus +------------- + .. automodapi:: ebtelplusplus + :no-heading: + :no-main-docstr: + :no-inheritance-diagram: + +ebtelplusplus.models +-------------------- .. automodapi:: ebtelplusplus.models + :no-heading: diff --git a/ebtelplusplus/__init__.py b/ebtelplusplus/__init__.py index 3b2b232..7437c69 100644 --- a/ebtelplusplus/__init__.py +++ b/ebtelplusplus/__init__.py @@ -1,11 +1,11 @@ """ ebtelplusplus: Zero-dimensional hydrodynamics of coronal loops """ -from ebtelplusplus.high_level import build_configuration, run +from ebtelplusplus.high_level import EbtelResult, run try: from ebtelplusplus._version import __version__ except ImportError: __version__ = "unknown" -__all__ = ["run", "build_configuration", "__version__"] +__all__ = ["run", "EbtelResult", "__version__"] diff --git a/ebtelplusplus/high_level.py b/ebtelplusplus/high_level.py index eb5ca23..5612e8e 100644 --- a/ebtelplusplus/high_level.py +++ b/ebtelplusplus/high_level.py @@ -2,45 +2,29 @@ High level interface to running an ebtel++ simulation """ import astropy.units as u +import dataclasses +import typing from astropy.utils.data import get_pkg_data_path -from collections import namedtuple import ebtelplusplus._low_level -from ebtelplusplus.models import DemModel, PhysicsModel, SolverModel +from ebtelplusplus.models import DemModel, HeatingModel, PhysicsModel, SolverModel -__all__ = ["run", "build_configuration"] - - -_UNITS_MAPPING = { - 'time': 's', - 'electron_temperature': 'K', - 'ion_temperature': 'K', - 'density': 'cm-3', - 'electron_pressure': 'dyne cm-2', - 'ion_pressure': 'dyne cm-2', - 'velocity': 'cm s-1', - 'heat': 'erg cm-3 s-1', - 'dem_temperature': 'K', - 'dem_tr': 'cm-5 K-1', - 'dem_corona': 'cm-5 K-1', - 'electron_thermal_conduction': 'erg cm-2 s-1', - 'ion_thermal_conduction': 'erg cm-2 s-1', - 'radiative_loss': 'erg cm-3 s-1', - 'tr_corona_radiative_loss_ratio': '', -} - - -EbtelResult = namedtuple('EbtelResult', - _UNITS_MAPPING.keys(), - defaults=[None,]*len(_UNITS_MAPPING)) +__all__ = ["run", "EbtelResult"] @u.quantity_input -def run(total_time: u.s, loop_length: u.cm, heating, physics=None, solver=None, dem=None): +def run(total_time: u.s, loop_length: u.cm, heating=None, physics=None, solver=None, dem=None): """ - Run an ebtel++ simulation + Run an ebtelplusplus simulation + + To run an ebtelplusplus simulation, at a minimum, you must specify a total simulation time + and a loop half length. Additional input parameters are configured through high-level objects + as listed below. Most users will likely want to modify the ``heating`` input by constructing + a `~ebtelplusplus.models.HeatingModel` that specifies the amount of energy injected into the + simulation as a function of time. Additional inputs are documented in the docstring of each + respective input model. Parameters ---------- @@ -48,35 +32,42 @@ def run(total_time: u.s, loop_length: u.cm, heating, physics=None, solver=None, Total duration of the simulation loop_length: `~astropy.units.Quantity` Loop half length - heating_model: `~ebtelplusplus.models.HeatingModel` + heating: `~ebtelplusplus.models.HeatingModel`, optional Configuration of heating model - physics_model: `~ebtelplusplus.models.PhysicsModel`, optional + physics: `~ebtelplusplus.models.PhysicsModel`, optional Configuration parameters related to the physics of the simulation - solver_model: `~ebtelplusplus.models.SolverModel`, optional + solver: `~ebtelplusplus.models.SolverModel`, optional Configuration parameters related to the numerical solver - dem_model: `~ebtelplusplus.models.DemModel`, optional + dem: `~ebtelplusplus.models.DemModel`, optional Configuration parameters related to the DEM calculation Returns ------- - results: `EbtelResult` - Data structure holding the results of the ebtel++ calculation + results: `~ebtelplusplus.EbtelResult` + Data structure holding the results of the ebtelplusplus calculation """ + if heating is None: + heating = HeatingModel() if solver is None: solver = SolverModel() if physics is None: physics = PhysicsModel() if dem is None: dem = DemModel() - config = build_configuration(total_time, loop_length, heating, physics, solver, dem) + config = _build_configuration(total_time, loop_length, heating, physics, solver, dem) results = ebtelplusplus._low_level.run(config) - return EbtelResult(**{k: u.Quantity(v, _UNITS_MAPPING[k]) for k,v in results.items()}) + return EbtelResult(results, config) @u.quantity_input -def build_configuration(total_time:u.s, loop_length:u.cm, heating, physics, solver, dem): +def _build_configuration(total_time:u.s, loop_length:u.cm, heating, physics, solver, dem): """ - Helper function for building a dictionary of ebtel++ configuration options + Helper function for building a dictionary of ebtelplusplus configuration options + + .. note:: + It is not necessary to use this function if you are running a simulation using + `~ebtelplusplus.run`. However, it may be useful to use this function directly + if you want to store the resulting inputs a dictionary. Parameters ---------- @@ -102,3 +93,95 @@ def build_configuration(total_time:u.s, loop_length:u.cm, heating, physics, solv **dem.to_dict(), 'heating': heating.to_dict(), } + + +@dataclasses.dataclass +class EbtelResult: + """ + Result of an ebtelplusplus simulation + + .. note:: + + This class is not meant to be instantiated directly. + Rather, it is meant to be returned by `~ebtelplusplus.run`. + + Parameters + ---------- + time: `~astropy.units.Quantity` + Simulation time at each time step. Has shape ``(N,)`` + electron_temperature: `~astropy.units.Quantity` + The temperature of the electrons at each time step + ion_temperature: `~astropy.units.Quantity` + The temperature of the ions at each time step + density: `~astropy.units.Quantity` + The density of the electrons and ions at each time step + electron_pressure: `~astropy.units.Quantity` + The pressure of the electrons at each time step + ion_pressure: `~astropy.units.Quantity` + The pressure of the ions at each time step + total_pressure: `~astropy.units.Quantity` + The total pressure at each time step. This is the sum of + the electron and ion pressures and is equivalent to the + pressure returned for single-fluid models. + velocity: `~astropy.units.Quantity` + The velocity at each timestep. This parameter should be + used cautiously as the EBTEL does not model the actual + transport of material through the spatial extent of the + loop. + heat: `~astropy.units.Quantity` + The energy released into the loop at each time step. + This is the amount of heating as specified in the + input `~ebtelplusplus.models.HeatingModel`. + electron_thermal_conduction: `~astropy.units.Quantity` + Electron thermal conductive flux at the TR-corona interface + at each time step. + ion_thermal_conduction: `~astropy.units.Quantity` + Ion thermal conductive flux at the TR-corona interface + at each time step. + radiative_loss: `~astropy.units.Quantity` + Energy lost due to radiation at each time step + tr_corona_radiative_loss_ratio: `~astropy.units.Quantity` + The ratio of the average radiative losses in the TR and the + corona. In the EBTEL model, this is often referred to as :math:`c_1`. + dem_temperature: `~astropy.units.Quantity` + The temperature grid on which the differential emission measure + distribution (DEM) is computed. Has shape ``(M,)``. These can be + interpreted as the centers of the temperature bins. + dem_tr: `~astropy.units.Quantity` + The TR DEM. Has shape ``(N,M)``. For more details of how this is + calculated, see Section 3 and the Appendix of :cite:t:`klimchuk_highly_2008`. + dem_corona: `~astropy.units.Quantity` + The coronal DEM. Has shape ``(N,M)``. + inputs: `dict` + All model inputs used to run the simulation. + """ + inputs: dict + time: u.Quantity[u.s] + electron_temperature: u.Quantity[u.K] + ion_temperature: u.Quantity[u.K] + density: u.Quantity[u.cm**(-3)] + electron_pressure: u.Quantity[u.dyne*u.cm**(-2)] + ion_pressure: u.Quantity[u.dyne*u.cm**(-2)] + velocity: u.Quantity[u.cm/u.s] + heat: u.Quantity[u.erg*u.cm**(-3)*u.s**(-1)] + electron_thermal_conduction: u.Quantity[u.erg*u.cm**(-2)*u.s**(-1)] + ion_thermal_conduction: u.Quantity[u.erg*u.cm**(-2)*u.s**(-1)] + radiative_loss: u.Quantity[u.erg*u.cm**(-3)*u.s**(-1)] + tr_corona_radiative_loss_ratio: u.Quantity[u.dimensionless_unscaled] + dem_temperature: u.Quantity[u.K] = None + dem_tr: u.Quantity[u.K**(-1)*u.cm**(-5)] = None + dem_corona: u.Quantity[u.K**(-1)*u.cm**(-5)] = None + # derived fields + total_pressure: u.Quantity[u.dyne*u.cm**(-2)] = dataclasses.field(init=None) + + def __init__(self, results, inputs): + # Use type hinting above to assign units to outputs; this is so + # we only have to specify the units of the outputs in one place + type_hints = type(self).__annotations__ + for k, v in results.items(): + _, unit = typing.get_args(type_hints[k]) + setattr(self, k, u.Quantity(v, unit)) + self.inputs = inputs + + def __post_init__(self): + self.total_pressure = self.electron_pressure + self.ion_pressure diff --git a/ebtelplusplus/models.py b/ebtelplusplus/models.py index c510534..0df793a 100644 --- a/ebtelplusplus/models.py +++ b/ebtelplusplus/models.py @@ -1,5 +1,8 @@ """ -Models for configuring ebtel++ runs +Classes for configuring inputs to `ebtelplusplus.run`. +All model inputs have default values such that each model only need be instantiated +using the values of the parameters a user wants to change. +The default values of each input are listed below. """ import astropy.units as u import dataclasses @@ -18,7 +21,7 @@ @dataclasses.dataclass class PhysicsModel: """ - ebtel++ input parameters related to the physics of the simulation + ebtelplusplus input parameters related to the physics of the simulation Parameters ---------- @@ -87,7 +90,7 @@ def to_dict(self): @dataclasses.dataclass class SolverModel: """ - ebtel++ input parameters related to the numerical solver + ebtelplusplus input parameters related to the numerical solver Parameters ---------- @@ -121,11 +124,11 @@ def to_dict(self): @dataclasses.dataclass class DemModel: """ - ebtel++ input parameters related to differential emission measure (DEM) calculation. + ebtelplusplus input parameters related to differential emission measure (DEM) calculation. - Optionally, ebtel++ can can also calculate the differential emission - measure (DEM) in both the transition region and the corona. See sections 2. - 2 and 3 of :cite:t:`klimchuk_highly_2008` for the details of this + Optionally, ebtelplusplus can can also calculate the differential emission + measure (DEM) in both the transition region and the corona. See sections 2.2 + and 3 of :cite:t:`klimchuk_highly_2008` for the details of this calculation. Note that this will result in much longer computation times. Parameters @@ -164,7 +167,7 @@ class HeatingEvent: Each event has a linear rise phase, a constant phase, and a linear decay phase. Using this format, it is easy to specify either - symmetric or asymmetric events of many different shapes + symmetric or asymmetric events of many different shapes. Parameters ---------- @@ -255,9 +258,9 @@ def __init__(self, rise_start, duration, rate): @dataclasses.dataclass class HeatingModel: """ - ebtel++ input parameters for the time-dependent heating + ebtelplusplus input parameters for the time-dependent heating - The ebtel++ time-dependent heating model is parameterized by a + The ebtelplusplus time-dependent heating model is parameterized by a series of discrete events (`HeatingEvent`) combined with a constant background heating rate. Furthermore, this energy can be injected into either the electrons or the ions or some admixture of the two. @@ -274,8 +277,8 @@ class HeatingModel: List of `HeatingEvent` objects that parameterize the energy injected into the loop by a series of discrete heating events. """ - background: u.Quantity[u.erg/(u.cm**3*u.s)] - partition: float = 1.0 + background: u.Quantity[u.erg/(u.cm**3*u.s)] = 1e-6 * u.Unit('erg cm-3 s-1') + partition: float = 0.5 events: list[HeatingEvent] = None def __post_init__(self): diff --git a/ebtelplusplus/tests/test_compare_idl.py b/ebtelplusplus/tests/test_compare_idl.py index 1c8701f..0dce7d3 100644 --- a/ebtelplusplus/tests/test_compare_idl.py +++ b/ebtelplusplus/tests/test_compare_idl.py @@ -7,7 +7,6 @@ import ebtelplusplus from ebtelplusplus.models import ( - DemModel, HeatingModel, PhysicsModel, SolverModel, @@ -57,16 +56,11 @@ def test_compare_idl_single_event(physics_model, heating_model, physics=physics_model, solver=solver_model) - config_dict = ebtelplusplus.build_configuration(total_time, - loop_length, - heating_model, - physics_model, - solver_model, - DemModel()) + r_idl = read_idl_test_data( 'idl_single_event.txt', ebtel_idl_path, - config_dict + r_cpp.inputs, ) if plot_idl_comparisons: plot_comparison(r_cpp, r_idl) @@ -98,13 +92,7 @@ def test_compare_idl_area_expansion( heating_model, physics=physics_model, solver=solver_model) - config_dict = ebtelplusplus.build_configuration(total_time, - loop_length, - heating_model, - physics_model, - solver_model, - DemModel()) - r_idl = read_idl_test_data(f'idl_area_expansion_{A_c=}_{A_0=}_{A_tr=}.txt', ebtel_idl_path, config_dict) + r_idl = read_idl_test_data(f'idl_area_expansion_{A_c=}_{A_0=}_{A_tr=}.txt', ebtel_idl_path, r_cpp.inputs) if plot_idl_comparisons: plot_comparison(r_cpp, r_idl) assert u.allclose(r_cpp.electron_temperature, r_idl['temperature'], rtol=RTOL) diff --git a/ebtelplusplus/tests/test_solver.py b/ebtelplusplus/tests/test_solver.py index 227fdd1..853b35a 100644 --- a/ebtelplusplus/tests/test_solver.py +++ b/ebtelplusplus/tests/test_solver.py @@ -99,6 +99,7 @@ def test_insufficient_heating(bad_heating, static_solver,): def test_NaNs_in_solver(use_adaptive_solver): solver = SolverModel(use_adaptive_solver=use_adaptive_solver) heating = HeatingModel( + partition=1.0, background=1e-6*u.Unit('erg cm-3 s-1'), events=[TriangularHeatingEvent(0.0*u.s, 200*u.s, -10*u.Unit('erg cm-3 s-1'))] ) From 66f0475970d7232a29a61a21fd10120f95c351f8 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Mon, 16 Sep 2024 12:00:50 -0400 Subject: [PATCH 2/6] add linux boost install instructions --- docs/development.rst | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/docs/development.rst b/docs/development.rst index 0417046..d9cca04 100644 --- a/docs/development.rst +++ b/docs/development.rst @@ -29,6 +29,12 @@ There are various ways to install this package, including: * - `conda-forge `__ - any - ``conda install -c conda-forge libboost-devel`` + * - APT + - Linux (Debian-based) + - ``apt install libboost-all-dev`` + * - YUM + - Linux (RHEL) + - ``yum install boost-devel`` * - `Homebrew `__ - macOS - ``brew install boost`` From 04e0ad40a043075f87e13f7e5b5b84a5f8344bf5 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Tue, 17 Sep 2024 02:04:15 -0400 Subject: [PATCH 3/6] minor touchups --- docs/development.rst | 2 +- ebtelplusplus/models.py | 2 +- pyproject.toml | 2 ++ 3 files changed, 4 insertions(+), 2 deletions(-) diff --git a/docs/development.rst b/docs/development.rst index d9cca04..bf0ce70 100644 --- a/docs/development.rst +++ b/docs/development.rst @@ -16,7 +16,7 @@ Next, clone your fork, git clone https://github.com//ebtelPlusPlus.git -Because ``ebtelplusplus`` implements the actual simulation code in C++ and thus is not a pure Python package, installation of this package requires first compiling the C++ and building the needed binaries. +Because ``ebtelplusplus`` implements the actual simulation code in C++ and thus is not a pure Python package, installation of this package requires first compiling the C++ code and building the needed binaries. To do this, you will first need to have the `Boost `__ package installed (at least v1.53). There are various ways to install this package, including: diff --git a/ebtelplusplus/models.py b/ebtelplusplus/models.py index 0df793a..ad2f313 100644 --- a/ebtelplusplus/models.py +++ b/ebtelplusplus/models.py @@ -258,7 +258,7 @@ def __init__(self, rise_start, duration, rate): @dataclasses.dataclass class HeatingModel: """ - ebtelplusplus input parameters for the time-dependent heating + ebtelplusplus input parameters for time-dependent heating The ebtelplusplus time-dependent heating model is parameterized by a series of discrete events (`HeatingEvent`) combined with a constant diff --git a/pyproject.toml b/pyproject.toml index 205d172..f7f14bc 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -135,6 +135,8 @@ lines_between_types = 1 [tool.codespell] skip = "*.fts,*.fits,venv,*.pro,*.bib,*.asdf,*.json" +ignore-words-list = "ebtel,hydrad" + [tool.ruff] target-version = "py310" From 0ffa07835090c90a52440ac88186b6eced9aa6e9 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Tue, 17 Sep 2024 02:05:12 -0400 Subject: [PATCH 4/6] add derivation topic guide;revamp landing page --- docs/conf.py | 2 +- docs/index.rst | 75 ++++--- docs/references.bib | 14 ++ docs/{ => topic_guides}/comparison.rst | 2 +- docs/topic_guides/derivation.rst | 271 +++++++++++++++++++++++++ docs/topic_guides/index.rst | 10 + 6 files changed, 344 insertions(+), 30 deletions(-) rename docs/{ => topic_guides}/comparison.rst (99%) create mode 100644 docs/topic_guides/derivation.rst create mode 100644 docs/topic_guides/index.rst diff --git a/docs/conf.py b/docs/conf.py index b6aa721..d1a37b7 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -92,6 +92,7 @@ html_theme = "pydata_sphinx_theme" html_theme_options = { + "show_nav_level": 2, "logo": { "text": f"ebtelplusplus {version}", }, @@ -119,7 +120,6 @@ html_sidebars = { "bibliography*": [], "development*": [], - "comparison*": [], } # Render inheritance diagrams in SVG graphviz_output_format = "svg" diff --git a/docs/index.rst b/docs/index.rst index b9e0326..3c657a2 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -6,23 +6,62 @@ ebtelplusplus :hidden: generated/gallery/index - comparison + topic_guides/index development reference bibliography ``ebtelplusplus`` (or ``ebtel++``) is a two-fluid version of the original enthalpy-based thermal evolution of loops (EBTEL) model implemented in C++ and wrapped in Python. -The EBTEL model, originally developed by :cite:t:`klimchuk_highly_2008`, efficiently computes spatially-averaged, time-dependent plasma parameters ( e.g. temperature, pressure, density) of dynamically-heated coronal loops. -It is often desirable to compute solutions for a large number of coronal loops, but the spatial and temporal scales needed to solve the full 1D-hydrodynamic equations lead to long compute times for even 1D hydrodynamic codes. -EBTEL computes quick and accurate solutions for spatially-averaged quantities, allowing efficient insight into how these monolithic structures evolve. -EBTEL also calculates the differential emission measure (DEM) for both the transition region and the corona. -Details regarding this formulation can be found in :cite:t:`klimchuk_highly_2008`. + +.. grid:: 1 2 2 2 + :gutter: 2 + + .. grid-item-card:: Examples + :link: generated/gallery + :text-align: center + + :material-outlined:`palette;8em;sd-text-secondary` + + Examples of how to set up a simulation + + .. grid-item-card:: Topic Guides + :link: ebtelplusplus-topic-guide + :link-type: ref + :text-align: center + + :material-outlined:`school;8em;sd-text-secondary` + + Detailed explanations about the EBTEL model + + .. grid-item-card:: Contributing + :link: ebtelplusplus-development + :link-type: ref + :text-align: center + + :material-outlined:`code;8em;sd-text-secondary` + + Instructions for how to contribute to `ebtelplusplus` + + .. grid-item-card:: Reference + :link: ebtelplusplus-reference + :link-type: ref + :text-align: center + + :material-outlined:`menu_book;8em;sd-text-secondary` + + Technical description of the inputs, outputs, and behavior of each component of ebteplusplus .. note:: If you are looking for the original IDL implementation, the repository for that code can be found `here `__. +The EBTEL model, originally developed by :cite:t:`klimchuk_highly_2008`, efficiently computes spatially-averaged, time-dependent plasma parameters ( e.g. temperature, pressure, density) of dynamically-heated coronal loops. +It is often desirable to compute solutions for a large number of coronal loops, but the spatial and temporal scales needed to solve the full 1D-hydrodynamic equations lead to long compute times for even 1D hydrodynamic codes. +EBTEL computes quick and accurate solutions for spatially-averaged quantities, allowing efficient insight into how these monolithic structures evolve. +EBTEL also calculates the differential emission measure (DEM) for both the transition region and the corona. +Details regarding this formulation can be found in :cite:t:`klimchuk_highly_2008`. + Installation ------------ @@ -40,38 +79,18 @@ Once you've successfully installed ``ebtelplusplus``, see the Example Gallery fo Citation -------- -If you use ``ebtelplusplus`` in any published work, please include citations to the following papers: +If you use ``ebtelplusplus`` in any published work, we ask that you please include citations to the following papers: * :cite:t:`klimchuk_highly_2008` * :cite:t:`cargill_enthalpy-based_2012` * :cite:t:`cargill_enthalpy-based_2012-1` * :cite:t:`barnes_inference_2016` +If you make use of the variable abundance feature in the radiative losses, please also cite :cite:t:`reep_modeling_2024`. Additionally, please add the following line within your methods, conclusion or acknowledgements sections: *This research used version X.Y.Z (software citation) of the ebtelplusplus open source software package.* The software citation should be the specific `Zenodo DOI`_ for the version used within your work. -The citation for most current version on Zenodo is, - -.. code:: bibtex - - @software{Barnes2024, - author = {Barnes, Will and - Klimchuk, James and - Cargill, Peter and - Bradshaw, Stephen and - Reep, Jeffrey and - Schonfeld, Samuel and - Collazzo, Rowan}, - title = {rice-solar-physics/ebtelPlusPlus: v0.2}, - month = jul, - year = 2024, - publisher = {Zenodo}, - version = {v0.2}, - doi = {10.5281/zenodo.12675386}, - url = {https://doi.org/10.5281/zenodo.12675386} - } - .. _Zenodo DOI: https://zenodo.org/records/12675374 diff --git a/docs/references.bib b/docs/references.bib index 53c6aa8..e7fd4f8 100644 --- a/docs/references.bib +++ b/docs/references.bib @@ -106,3 +106,17 @@ @article{reep_modeling_2024 annotation = {ADS Bibcode: 2024ApJ...970L..41R} } +@article{reep_geometric_2022, + title = {Geometric {{Assumptions}} in {{Hydrodynamic Modeling}} of {{Coronal}} and {{Flaring Loops}}}, + author = {Reep, Jeffrey W. and {Ugarte-Urra}, Ignacio and Warren, Harry P. and Barnes, Will T.}, + year = {2022}, + month = jul, + journal = {The Astrophysical Journal}, + volume = {933}, + pages = {106}, + issn = {0004-637X}, + doi = {10.3847/1538-4357/ac7398}, + urldate = {2023-08-29}, + annotation = {ADS Bibcode: 2022ApJ...933..106R} +} + diff --git a/docs/comparison.rst b/docs/topic_guides/comparison.rst similarity index 99% rename from docs/comparison.rst rename to docs/topic_guides/comparison.rst index 1409c2d..a5230b2 100644 --- a/docs/comparison.rst +++ b/docs/topic_guides/comparison.rst @@ -1,4 +1,4 @@ -.. _ebtelplusplus-comparison: +.. _ebtelplusplus-topic-guide-comparison: Why ebtelplusplus? ====================== diff --git a/docs/topic_guides/derivation.rst b/docs/topic_guides/derivation.rst new file mode 100644 index 0000000..00b9bcd --- /dev/null +++ b/docs/topic_guides/derivation.rst @@ -0,0 +1,271 @@ +.. _ebtelplusplus-topic-guide-derivation: + +Deriving the EBTEL Equations +============================ + +This page will briefly summarize the derivation of the version of the EBTEL +model used in `ebtelplusplus` which includes cross-sectional area expansion +and two-fluid effects. +For a more detailed explanation of the EBTEL model, including the physical +motivation behind the assumptions, see :cite:t:`klimchuk_highly_2008`. +Additionally, :cite:t:`barnes_inference_2016` provides a detailed discussion +of the two-fluid implementation while :cite:t:`cargill_static_2022` discusses +the inclusion of cross-sectional area expansion. + +The 1D Hydrodynamic Equations +----------------------------- + +Deriving the EBTEL model involves taking spatial integrals, +over both the transition region (TR) and corona, of the one-dimensional (1D) +hydrodynamic equations that govern the evolution of the dynamically-heated +solar atmosphere confined along a magnetic field line. +This yields a set of so-called "zero-dimensional" (0D) equations which +describe the evolution in time of the coronally-averaged temperature, density, +and pressure. +We begin by listing the relevant 1D hydrodynamic equations. +The 1D density, electron energy, and ion energy equations are given by, + +.. math:: + :label: density-1d + + \frac{\partial n}{\partial t} + \frac{1}{A}\frac{\partial}{\partial s}(Anv) = 0, + +.. math:: + :label: electron-energy-1d + + \frac{\partial E_e}{\partial t} + \frac{1}{A}\frac{\partial}{\partial s}(A(E_e + p_e)v) = &-\frac{1}{A}\frac{\partial}{\partial s}(AF_e) + v\frac{\partial p_e}{\partial s} - n^2\Lambda \\ + &+ \frac{k_Bn}{\gamma-1}\nu_{ie}(T_i - T_e) + Q_e, + +.. math:: + :label: ion-energy-1d + + \frac{\partial E_i}{\partial t} + \frac{1}{A}\frac{\partial}{\partial s}(A(E_i + p_i)v) = &-\frac{1}{A}\frac{\partial}{\partial s}(AF_i) - v\frac{\partial p_e}{\partial s} + m_i n v g_{\parallel} + \frac{1}{A}\left(\frac{4A\eta v}{3}\frac{\partial v}{\partial s}\right) \\ + &+ \frac{k_Bn}{\gamma-1}\nu_{ie}(T_e - T_i) + Q_i, + +where :math:`s` is the field-aligned spatial coordinate, :math:`n` is the density, +:math:`v` is the velocity, :math:`A` is the cross-sectional area, :math:`F_e` and +:math:`F_i` are the electron and ion thermal conductive fluxes, respectively, +:math:`T_e` and :math:`T_i` are the electron and ion temperatures respectively, +:math:`\Lambda` describes the energy lost to radiation, and :math:`Q_e,Q_i` describe +the *ad hoc* heating applied to the electrons and ions, respectively.These equations +are closed by set of equations of state for the electron and ion energies, :math:`E_e,E_i`, +and pressures, :math:`p_e,p_i`, + +.. math:: + :label: eqs-of-state + + E_e &= \frac{p_e}{\gamma - 1}, \\ + E_i &= \frac{p_i}{\gamma - 1} + \frac{1}{2}m_inv^2, \\ + p_e &= nk_BT_e, \\ + p_i &= nk_BT_i. + +For additional details about these 1D hydrodynamic equations, see :cite:t:`reep_geometric_2022`. + +Using the above equations of state coupled with the assumptions that (1) all flows are +subsonic and (2) gravity is negligible for loops of length :math:`<150` Mm +:cite:p:`{see Section 2 of}klimchuk_highly_2008`, we can simplify :eq:`electron-energy-1d` +and :eq:`ion-energy-1d` to, + +.. math:: + :label: electron-energy-1d-simple + + \frac{A}{\gamma-1}\frac{\partial p_e}{\partial t} + \frac{\gamma}{\gamma-1}\frac{\partial}{\partial s}(Ap_ev) = &-\frac{\partial}{\partial s}(AF_e) + Av\frac{\partial p_e}{\partial s} - An^2\Lambda \\ + &+ Ak_Bn\nu_{ie}(T_i - T_e) + AQ_e, + +.. math:: + :label: ion-energy-1d-simple + + \frac{A}{\gamma-1}\frac{\partial p_i}{\partial t} + \frac{\gamma}{\gamma-1}\frac{\partial}{\partial s}(Ap_iv) = -\frac{\partial}{\partial s}(AF_i) - Av\frac{\partial p_e}{\partial s} + Ak_Bn\nu_{ie}(T_e - T_i) + AQ_i, + +We can now apply the combined methodology of both :cite:t:`barnes_inference_2016` and +:cite:t:`cargill_static_2022` to Eqs. :eq:`density-1d`, :eq:`electron-energy-1d-simple`, +and :eq:`ion-energy-1d-simple` to derive the EBTEL equations including both two-fluid +effects and cross-sectional expansion. + +The EBTEL Electron Pressure Equation +------------------------------------ + +To derive the EBTEL electron pressure equation, we begin by taking a spatial integral +over Eq. :eq:`electron-energy-1d-simple` from the base of the corona to the apex of a +semi-circular loop that is symmetric about the apex, + +.. math:: + :label: e-energy-coronal-integral + + \frac{A_cL_c}{\gamma-1}\frac{d p_{e,c}}{dt} - \frac{\gamma}{\gamma-1}(Ap_ev)_0 = (AF_e)_0 + A_cL_cQ_{e,c} + A_c\psi_c - A_cR_c, + +where :math:`c` denotes an average taken over the coronal portion of the loop, :math:`0` +denotes evaluation at the TR-corona interface, :math:`L_c` is the coronal portion of the +loop half-length, and, + +.. math:: + :label: psi-corona + + \psi_c = \frac{1}{A_c}\int_c\mathrm{d}s\,Av\frac{\partial p_e}{\partial s} + \frac{1}{A_c}\int_c\mathrm{d}s\,\frac{Ak_Bn\nu_{ie}}{\gamma-1}(T_i - T_e), + +.. math:: + :label: losses-corona + + R_c = \frac{1}{A_c}\int_c\mathrm{d}s\,An^2\Lambda. + +Note that the coronal integral :math:`\int_c\mathrm{d}s=\int_{s=L_{TR}}^{s=L}\mathrm{d}s`, +where :math:`L_{TR}` is the length of the TR and :math:`L=L_{TR}+L_c` is the total loop +half-length from the bottom of the TR to the apex of the loop. Additionally, because the +loop is assumed symmetric about the apex and isolated from the lower atmosphere, the +velocity and heat flux terms vanish at those locations. + +Similarly, we can integrate Eq. :eq:`electron-energy-1d-simple` over the TR, + +.. math:: + :label: e-energy-tr-integral + + \frac{A_{TR}L_{TR}}{\gamma-1}\frac{d p_{e,TR}}{dt} + \frac{\gamma}{\gamma-1}(Ap_ev)_0 = &-(AF_e)_0 + A_{TR}L_{TR}Q_{e,TR} \\ + &+ A_{TR}\psi_{TR} - A_{TR}R_{TR}, + +where :math:`TR` denotes an average taken over the TR portion of the loop, :math:`\psi_{TR}` +and :math:`R_{TR}` have the same form as Eqs. :eq:`psi-corona` and :eq:`losses-corona`, respectively, +and the TR integral :math:`\int_{TR}\mathrm{d}s=\int_{s=0}^{s=L_{TR}}`. + +Following the approach of :cite:t:`cargill_static_2022`, we add Eqs. +:eq:`e-energy-coronal-integral` and :eq:`e-energy-tr-integral` and let :math:`p_e=p_{e,c}=p_{e,TR}` +and :math:`Q_e=Q_{e,c}=Q_{e,TR}` to get the EBTEL electron pressure equation, + +.. math:: + :label: ebtel-electron-pressure + + \boxed{\frac{1}{\gamma-1}\frac{dp_e}{dt} = Q_e + \frac{\psi_c}{L_*}\left(1+\frac{A_{TR}\psi_{TR}}{A_c\psi_c}\right) - \frac{R_c}{L_*}\left(1+c_1\frac{A_{TR}}{A_c}\right)}, + +where :math:`L_* = L_c + (A_{TR}/A_c)L_{TR}` and :math:`c_1=R_{TR}/R_c`. Eq. :eq:`ebtel-electron-pressure` describes the time-evolution of the spatially-averaged electron pressure. + +The EBTEL Ion Pressure Equation +------------------------------------ + +To derive the EBTEL ion pressure equation, we apply the same procedure as above to Eq. +:eq:`ion-energy-1d-simple`. The spatial integral of Eq. :eq:`ion-energy-1d-simple` over +the coronal portion of the loop is, + +.. math:: + :label: i-energy-c-integral + + \frac{A_cL_c}{\gamma-1}\frac{d p_{i,c}}{dt} - \frac{\gamma}{\gamma-1}(Ap_iv)_0 = (AF_i)_0 + A_cL_cQ_{i,c} - A_c\psi_c, + +and for the TR portion of the loop is, + +.. math:: + :label: i-energy-tr-integral + + \frac{A_{TR}L_{TR}}{\gamma-1}\frac{d p_{i,TR}}{dt} + \frac{\gamma}{\gamma-1}(Ap_iv)_0 = -(AF_i)_0 + A_{TR}L_{TR}Q_{i,TR} - A_{TR}\psi_{TR}. + +As above, we can add Eqs. :eq:`i-energy-c-integral` and :eq:`i-energy-tr-integral` together +and let :math:`p_i=p_{i,c}=p_{i,TR}` and :math:`Q_i=Q_{i,c}=Q_{i,TR}` to obtain the EBTEL +ion pressure equation, + +.. math:: + :label: ebtel-ion-pressure + + \boxed{\frac{1}{\gamma-1}\frac{dp_i}{dt} = Q_i - \frac{\psi_c}{L_*}\left(1 + \frac{A_{TR}\psi_{TR}}{A_c\psi_c}\right)}. + +Eq. :eq:`ebtel-ion-pressure` describes the time-evolution of the spatially-averaged ion pressure. + +The EBTEL Density Equation +------------------------------------ + +Lastly, we derive the EBTEL density equation. +We begin by taking a spatial integral of Eq. :eq:`density-1d` over the coronal portion of the loop, + +.. math:: + + A_cL_c\frac{dn}{dt} - (Anv)_0 = 0 + +and using the equation of state for :math:`p_e` from Eq. :eq:`eqs-of-state`, + +.. math:: + :label: density-coronal-integral + + A_cL_c\frac{dn}{dt} = \frac{(Ap_ev)_0}{k_BT_{e,0}}. + +The quantity :math:`(Ap_ev)_0` is the area-weighted electron enthalpy flux at the TR-corona interface. +We can derive an expression for this term by substituting Eq. :eq:`ebtel-electron-pressure` into +Eq. :eq:`e-energy-tr-integral` and doing a lot of algebra, + +.. math:: + :label: electron-enthalpy-flux + + \frac{\gamma}{\gamma-1}(Ap_ev)_0 = -\frac{A_{TR}L_cR_c}{L_*}\left(c_1 - \frac{L_{TR}}{L_c}\right) + \frac{A_{TR}L_c\psi_c}{L_*}\left(\frac{\psi_{TR}}{\psi_c} - \frac{L_{TR}}{L_c}\right) - (AF_e)_0. + +Similarly, we can derive an expression for the area-weighted ion enthalpy flux by substituting +Eq. :eq:`ebtel-ion-pressure` into Eq. :eq:`i-energy-tr-integral`, + +.. math:: + :label: ion-enthalpy-flux + + \frac{\gamma}{\gamma-1}(Ap_iv)_0 = -\frac{A_{TR}L_c\psi_c}{L_*}\left(\frac{\psi_{TR}}{\psi_c} - \frac{L_{TR}}{L_c}\right) - (AF_i)_0. + +Adding Eqs. :eq:`electron-enthalpy-flux` and :eq:`ion-enthalpy-flux`, + +.. math:: + :label: total-enthalpy-flux + + \frac{\gamma}{\gamma-1}(Ap_ev)_0 + \frac{\gamma}{\gamma-1}(Ap_iv)_0 = -\frac{A_{TR}L_cR_c}{L_*}\left(c_1 - \frac{L_{TR}}{L_c}\right) - (AF_e)_0 - (AF_i)_0. + +Additionally, we define :math:`\xi\equiv T_e/T_i` and again use the electron and ion equations +of state from Eq. :eq:`eqs-of-state` to find, + +.. math:: + :label: temperature-ratio + + \xi = \frac{T_e}{T_i} = \frac{T_{e,0}}{T_{i,0}} = \frac{A_0n_0k_BT_{e,0}v_0}{A_0n_0k_BT_{i,0}v_0} = \frac{(Ap_ev)_0}{(Ap_iv)_0}. + +Substituting this into Eq. :eq:`total-enthalpy-flux` gives, + +.. math:: + :label: electron-enthalpy-flux-simple + + (Ap_ev)_0 = -\frac{\xi(\gamma - 1)}{\gamma(\xi + 1)}\left(\frac{A_{TR}L_cR_c}{L_*}\left(c_1 - \frac{L_{TR}}{L_c}\right) + (AF_e)_0 - (AF_i)_0\right). + +Finally, substituting Eq. :eq:`electron-enthalpy-flux-simple` into Eq. :eq:`density-coronal-integral` +yields the EBTEL density equation, + +.. math:: + :label: ebtel-density + + \boxed{\frac{dn}{dt} = -\frac{(\gamma - 1)\xi c_2}{(\xi + 1)\gamma c_3 k_B L_c T_e}\left(\frac{A_{TR}L_c}{A_cL_*}R_c\left(c_1 - \frac{L_{TR}}{L_c}\right) + \frac{A_0}{A_c}(F_{e,0} + F_{i,0})\right)}, + +where :math:`c_2,c_3` are constants relating the base temperature to the spatially-averaged coronal +temperature. Eq. :eq:`ebtel-density` describes the time-evolution of the spatial averaged density. +In summary, Eqs. :eq:`ebtel-electron-pressure`, :eq:`ebtel-ion-pressure`, and :eq:`ebtel-density` +comprise the two-fluid equations including cross-sectional area expansion. + +Limiting Behavior +----------------- + +Below, we briefly describe how Eqs. :eq:`ebtel-electron-pressure`, :eq:`ebtel-ion-pressure`, and +:eq:`ebtel-density` reduce to the other versions of the EBTEL model. + +Constant Cross-section +++++++++++++++++++++++ + +Under the assumption of constant cross-section, :math:`A_c=A_{TR}=A_0` and :math:`L_{TR}\ll L_c` +such that :math:`L_*\approx L_c = L`. As such, the EBTEL equations simplify to, + +.. math:: + + \frac{dp_e}{dt} &= (\gamma-1)Q_e + \frac{\gamma-1}{L}(\psi_c+\psi_{TR} - R_c(1+c_1)), \\ + \frac{dp_i}{dt} &= (\gamma-1)Q_i - \frac{\gamma-1}{L}(\psi_c + \psi_{TR}), \\ + \frac{dn}{dt} &= -\frac{(\gamma - 1)\xi c_2}{(\xi + 1)\gamma c_3 k_B L_c T_e}(c_1R_c + F_{e,0} + F_{i,0}). + +Note that these expressions are equivalent to the two-fluid EBTEL equations as given in :cite:t:`barnes_inference_2016`. + +Single-fluid +++++++++++++ + +Under the single-fluid assumption, :math:`T_e=T_i` at all times. +Using this and adding Eqs. :eq:`ebtel-electron-pressure` and :eq:`ebtel-ion-pressure`, + +.. math:: + + \frac{dp}{dt} &= (\gamma - 1)\left(Q - \frac{R_c}{L_*}\left(1 + \frac{A_{TR}}{A_c}c_1\right)\right), \\ + \frac{dn}{dt} &= -\frac{(\gamma - 1)c_2}{2\gamma c_3 k_B L_c T}\left(\frac{A_{TR}L_c}{A_cL_*}R_c\left(c_1 - \frac{L_{TR}}{L_c}\right) + \frac{A_0}{A_c}F_0\right). + +Note that these expressions are equivalent to the expanding cross-section EBTEL equations given +in :cite:t:`cargill_static_2022`. diff --git a/docs/topic_guides/index.rst b/docs/topic_guides/index.rst new file mode 100644 index 0000000..348dafd --- /dev/null +++ b/docs/topic_guides/index.rst @@ -0,0 +1,10 @@ +.. _ebtelplusplus-topic-guide: + +Topic Guides +============ + +.. toctree:: + :maxdepth: 1 + + comparison + derivation From c94f60b2f94a814e04c4e75a39a86893b9caba9b Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Tue, 17 Sep 2024 02:25:36 -0400 Subject: [PATCH 5/6] add reference to derivation page --- docs/index.rst | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/docs/index.rst b/docs/index.rst index 3c657a2..fcc70a7 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -11,7 +11,7 @@ ebtelplusplus reference bibliography -``ebtelplusplus`` (or ``ebtel++``) is a two-fluid version of the original enthalpy-based thermal evolution of loops (EBTEL) model implemented in C++ and wrapped in Python. +``ebtelplusplus`` (or ``ebtel++``) is an implementation of the enthalpy-based thermal evolution of loops (EBTEL) model implemented in C++ and wrapped in Python. .. grid:: 1 2 2 2 :gutter: 2 @@ -61,6 +61,7 @@ It is often desirable to compute solutions for a large number of coronal loops, EBTEL computes quick and accurate solutions for spatially-averaged quantities, allowing efficient insight into how these monolithic structures evolve. EBTEL also calculates the differential emission measure (DEM) for both the transition region and the corona. Details regarding this formulation can be found in :cite:t:`klimchuk_highly_2008`. +For a more detailed overview of the model, see :ref:`ebtelplusplus-topic-guide-derivation`. Installation ------------ From d763fc37ef01e8e4e57e9de6073d6ce265dda591 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Tue, 17 Sep 2024 02:37:00 -0400 Subject: [PATCH 6/6] touchup topic guide --- docs/topic_guides/derivation.rst | 1 + 1 file changed, 1 insertion(+) diff --git a/docs/topic_guides/derivation.rst b/docs/topic_guides/derivation.rst index 00b9bcd..597f364 100644 --- a/docs/topic_guides/derivation.rst +++ b/docs/topic_guides/derivation.rst @@ -46,6 +46,7 @@ where :math:`s` is the field-aligned spatial coordinate, :math:`n` is the densit :math:`v` is the velocity, :math:`A` is the cross-sectional area, :math:`F_e` and :math:`F_i` are the electron and ion thermal conductive fluxes, respectively, :math:`T_e` and :math:`T_i` are the electron and ion temperatures respectively, +:math:`\nu_{ie}` is the electron-ion collision frequency, :math:`\Lambda` describes the energy lost to radiation, and :math:`Q_e,Q_i` describe the *ad hoc* heating applied to the electrons and ions, respectively.These equations are closed by set of equations of state for the electron and ion energies, :math:`E_e,E_i`,