diff --git a/CHANGELOG.md b/CHANGELOG.md index 895cfef6..f6ec0008 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -24,10 +24,15 @@ New: * Improved beam direction calculation to allow for natural broadening of the BES line shape due to beam divergence. (#414) * Add kwargs to invert_regularised_nnls to pass them to scipy.optimize.nnls. (#438) * StarkBroadenedLine now supports Doppler broadening and Zeeman splitting. (#393) +* Add thermal charge-exchange emission model. (#57) +* PECs for C VI spectral lines for n <= 5 are now included in populate(). Rerun populate() after upgrading to 1.5 to update the atomic data repository. +* All interpolated atomic rates now return 0 if plasma parameters <= 0, which matches the behaviour of emission models. (#450) Bug fixes: * Fix deprecated transforms being cached in LaserMaterial after laser.transform update (#420) * Fix IRVB calculate sensitivity method. +* Fix missing donor_metastable attribute in the core BeamCXPEC class (#411). +* **Fix the receiver ion density being passed to the BeamCXPEC instead of the total ion density in the BeamCXLine. Also fix incorrect BeamCXPEC dosctrings. Attention!!! The results of CX spectroscopy are affected by this change. (#441)** Release 1.4.0 (3 Feb 2023) ------------------- diff --git a/cherab/core/atomic/interface.pxd b/cherab/core/atomic/interface.pxd index 2b671176..1dacfb57 100644 --- a/cherab/core/atomic/interface.pxd +++ b/cherab/core/atomic/interface.pxd @@ -45,6 +45,8 @@ cdef class AtomicData: cpdef RecombinationPEC recombination_pec(self, Element ion, int charge, tuple transition) + cpdef ThermalCXPEC thermal_cx_pec(self, Element donor_ion, int donor_charge, Element receiver_ion, int receiver_charge, tuple transition) + cpdef TotalRadiatedPower total_radiated_power(self, Element element) cpdef LineRadiationPower line_radiated_power_rate(self, Element element, int charge) diff --git a/cherab/core/atomic/interface.pyx b/cherab/core/atomic/interface.pyx index 9e4e1113..52c396d6 100644 --- a/cherab/core/atomic/interface.pyx +++ b/cherab/core/atomic/interface.pyx @@ -1,6 +1,6 @@ -# Copyright 2016-2023 Euratom -# Copyright 2016-2023 United Kingdom Atomic Energy Authority -# Copyright 2016-2023 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -31,70 +31,124 @@ cdef class AtomicData: cpdef double wavelength(self, Element ion, int charge, tuple transition): """ - Returns the natural wavelength of the specified transition in nm. + The natural wavelength of the specified transition in nm. """ raise NotImplementedError("The wavelength() virtual method is not implemented for this atomic data source.") cpdef IonisationRate ionisation_rate(self, Element ion, int charge): + """ + Electron impact ionisation rate for a given species in m^3/s. + """ + raise NotImplementedError("The ionisation_rate() virtual method is not implemented for this atomic data source.") cpdef RecombinationRate recombination_rate(self, Element ion, int charge): + """ + Recombination rate for a given species in m^3/s. + """ + raise NotImplementedError("The recombination_rate() virtual method is not implemented for this atomic data source.") cpdef ThermalCXRate thermal_cx_rate(self, Element donor_ion, int donor_charge, Element receiver_ion, int receiver_charge): + """ + Thermal charge exchange effective rate coefficient for a given donor and receiver species in m^3/s. + """ + raise NotImplementedError("The thermal_cx_rate() virtual method is not implemented for this atomic data source.") cpdef list beam_cx_pec(self, Element donor_ion, Element receiver_ion, int receiver_charge, tuple transition): """ - Returns a list of applicable charge exchange emission rates in W.m^3. + A list of Effective charge exchange photon emission coefficient for a given donor (beam) in W.m^3. """ raise NotImplementedError("The cxs_rates() virtual method is not implemented for this atomic data source.") cpdef BeamStoppingRate beam_stopping_rate(self, Element beam_ion, Element plasma_ion, int charge): """ - Returns a list of applicable beam stopping coefficients in m^3/s. + Beam stopping coefficient for a given beam and target species in m^3/s. """ raise NotImplementedError("The beam_stopping() virtual method is not implemented for this atomic data source.") cpdef BeamPopulationRate beam_population_rate(self, Element beam_ion, int metastable, Element plasma_ion, int charge): """ - Returns a list of applicable dimensionless beam population coefficients. + Dimensionless Beam population coefficient for a given beam and target species. """ raise NotImplementedError("The beam_population() virtual method is not implemented for this atomic data source.") cpdef BeamEmissionPEC beam_emission_pec(self, Element beam_ion, Element plasma_ion, int charge, tuple transition): """ - Returns a list of applicable beam emission coefficients in W.m^3. + The beam photon emission coefficient for a given beam and target species + and a given transition in W.m^3. """ raise NotImplementedError("The beam_emission() virtual method is not implemented for this atomic data source.") cpdef ImpactExcitationPEC impact_excitation_pec(self, Element ion, int charge, tuple transition): + """ + Electron impact excitation photon emission coefficient for a given species in W.m^3. + """ + raise NotImplementedError("The impact_excitation() virtual method is not implemented for this atomic data source.") cpdef RecombinationPEC recombination_pec(self, Element ion, int charge, tuple transition): + """ + Recombination photon emission coefficient for a given species in W.m^3. + """ + raise NotImplementedError("The recombination() virtual method is not implemented for this atomic data source.") + cpdef ThermalCXPEC thermal_cx_pec(self, Element donor_ion, int donor_charge, Element receiver_ion, int receiver_charge, tuple transition): + """ + Thermal charge exchange photon emission coefficient for given donor and receiver species in W.m^3. + """ + + raise NotImplementedError("The thermal_cx_pec() virtual method is not implemented for this atomic data source.") + cpdef TotalRadiatedPower total_radiated_power(self, Element element): + """ + The total (summed over all charge states) radiated power + in equilibrium conditions for a given species in W.m^3. + """ + raise NotImplementedError("The total_radiated_power() virtual method is not implemented for this atomic data source.") cpdef LineRadiationPower line_radiated_power_rate(self, Element element, int charge): + """ + Line radiated power coefficient for a given species in W.m^3. + """ + raise NotImplementedError("The line_radiated_power_rate() virtual method is not implemented for this atomic data source.") cpdef ContinuumPower continuum_radiated_power_rate(self, Element element, int charge): + """ + Continuum radiated power coefficient for a given species in W.m^3. + """ + raise NotImplementedError("The continuum_radiated_power_rate() virtual method is not implemented for this atomic data source.") cpdef CXRadiationPower cx_radiated_power_rate(self, Element element, int charge): + """ + Charge exchange radiated power coefficient for a given species in W.m^3. + """ + raise NotImplementedError("The cx_radiated_power_rate() virtual method is not implemented for this atomic data source.") cpdef FractionalAbundance fractional_abundance(self, Element ion, int charge): + """ + Fractional abundance of a given species in thermodynamic equilibrium. + """ + raise NotImplementedError("The fractional_abundance() virtual method is not implemented for this atomic data source.") cpdef ZeemanStructure zeeman_structure(self, Line line, object b_field=None): + r""" + Wavelengths and ratios of :math:`\pi`-/:math:`\sigma`-polarised Zeeman components + for any given value of magnetic field strength. + """ + raise NotImplementedError("The zeeman_structure() virtual method is not implemented for this atomic data source.") cpdef tuple zeeman_triplet_parameters(self, Line line): diff --git a/cherab/core/atomic/rates.pxd b/cherab/core/atomic/rates.pxd index 1f844122..41871191 100644 --- a/cherab/core/atomic/rates.pxd +++ b/cherab/core/atomic/rates.pxd @@ -46,11 +46,13 @@ cdef class RecombinationPEC(_PECRate): pass -cdef class ThermalCXPEC(_PECRate): - pass +cdef class ThermalCXPEC: + cpdef double evaluate(self, double electron_density, double electron_temperature, double donor_temperature) except? -1e999 cdef class BeamCXPEC: + cdef readonly int donor_metastable + cpdef double evaluate(self, double energy, double temperature, double density, double z_effective, double b_field) except? -1e999 diff --git a/cherab/core/atomic/rates.pyx b/cherab/core/atomic/rates.pyx index f3a653db..37225e73 100644 --- a/cherab/core/atomic/rates.pyx +++ b/cherab/core/atomic/rates.pyx @@ -130,11 +130,27 @@ cdef class RecombinationPEC(_PECRate): pass -cdef class ThermalCXPEC(_PECRate): +cdef class ThermalCXPEC: """ Thermal charge exchange rate coefficient. """ - pass + + def __call__(self, double electron_density, double electron_temperature, donor_temperature): + """Returns a CX photon emissivity coefficient at the specified plasma conditions. + + This function just wraps the cython evaluate() method. + """ + return self.evaluate(electron_density, electron_temperature, donor_temperature) + + cpdef double evaluate(self, double electron_density, double electron_temperature, double donor_temperature) except? -1e999: + """Returns a CX photon emissivity coefficient at given conditions. + + :param electron_density: Electron density in m^-3. + :param electron_temperature: Electron temperature in eV. + :param donor_temperature: Donor temperature in eV. + :return: The effective CX PEC rate in W/m^3. + """ + raise NotImplementedError("The evaluate() virtual method must be implemented.") cdef class BeamCXPEC: @@ -144,8 +160,13 @@ cdef class BeamCXPEC: transition :math:`n\rightarrow n'` of ion :math:`Z^{(\alpha+1)+}` with electron donor :math:`H^0` in metastable state :math:`m_{i}`. Equivalent to :math:`q^{eff}_{n\rightarrow n'}` in `adf12 _`. + + :param donor_metastable: The metastable state of the donor species for which the rate data applies. """ + def __init__(self, int donor_metastable): + self.donor_metastable = donor_metastable + def __call__(self, double energy, double temperature, double density, double z_effective, double b_field): """Evaluates the Beam CX rate at the given plasma conditions. @@ -158,7 +179,7 @@ cdef class BeamCXPEC: :param float energy: Interaction energy in eV/amu. :param float temperature: Receiver ion temperature in eV. - :param float density: Receiver ion density in m^-3 + :param float density: Plasma total ion density in m^-3 :param float z_effective: Plasma Z-effective. :param float b_field: Magnetic field magnitude in Tesla. :return: The effective rate diff --git a/cherab/core/model/beam/charge_exchange.pxd b/cherab/core/model/beam/charge_exchange.pxd index cd0b48db..46e22891 100644 --- a/cherab/core/model/beam/charge_exchange.pxd +++ b/cherab/core/model/beam/charge_exchange.pxd @@ -36,7 +36,7 @@ cdef class BeamCXLine(BeamModel): object _lineshape_class, _lineshape_args, _lineshape_kwargs cdef double _composite_cx_rate(self, double x, double y, double z, double interaction_energy, - Vector3D donor_velocity, double receiver_temperature, double receiver_density) except? -1e999 + Vector3D donor_velocity, double receiver_temperature) except? -1e999 cdef double _beam_population(self, double x, double y, double z, Vector3D beam_velocity, list population_data) except? -1e999 diff --git a/cherab/core/model/beam/charge_exchange.pyx b/cherab/core/model/beam/charge_exchange.pyx index 50f34249..4b5c292e 100644 --- a/cherab/core/model/beam/charge_exchange.pyx +++ b/cherab/core/model/beam/charge_exchange.pyx @@ -61,6 +61,7 @@ cdef class BeamCXLine(BeamModel): :ivar Line line: The emission line object. .. code-block:: pycon + >>> from cherab.core.model import BeamCXLine >>> from cherab.core.atomic import carbon >>> from cherab.core.model import ParametrisedZeemanTriplet @@ -158,7 +159,7 @@ cdef class BeamCXLine(BeamModel): interaction_energy = ms_to_evamu(interaction_speed) # calculate the composite charge-exchange emission coefficient - emission_rate = self._composite_cx_rate(x, y, z, interaction_energy, donor_velocity, receiver_temperature, receiver_density) + emission_rate = self._composite_cx_rate(x, y, z, interaction_energy, donor_velocity, receiver_temperature) # spectral line emission in W/m^3/str radiance = RECIP_4_PI * donor_density * receiver_density * emission_rate @@ -169,7 +170,7 @@ cdef class BeamCXLine(BeamModel): @cython.wraparound(False) @cython.cdivision(True) cdef double _composite_cx_rate(self, double x, double y, double z, double interaction_energy, - Vector3D donor_velocity, double receiver_temperature, double receiver_density) except? -1e999: + Vector3D donor_velocity, double receiver_temperature) except? -1e999: """ Performs a beam population weighted average of the effective cx rates. @@ -188,23 +189,26 @@ cdef class BeamCXLine(BeamModel): :param interaction_energy: The donor-receiver interaction energy in eV/amu. :param donor_velocity: A Vector defining the donor particle velocity in m/s. :param receiver_temperature: The receiver species temperature in eV. - :param receiver_density: The receiver species density in m^-3 :return: The composite charge exchange rate in W.m^3. """ cdef: - double z_effective, b_field, rate, total_population, population, effective_rate + double ion_density, z_effective + double b_field + double rate, effective_rate + double population, total_population BeamCXPEC cx_rate list population_data - # calculate z_effective and the B-field magnitude + # calculate ion density, z_effective and the B-field magnitude + ion_density = self._plasma.ion_density(x, y, z) z_effective = self._plasma.z_effective(x, y, z) b_field = self._plasma.get_b_field().evaluate(x, y, z).get_length() # rate for the ground state (metastable = 1) rate = self._ground_beam_rate.evaluate(interaction_energy, receiver_temperature, - receiver_density, + ion_density, z_effective, b_field) @@ -219,7 +223,7 @@ cdef class BeamCXLine(BeamModel): effective_rate = cx_rate.evaluate(interaction_energy, receiver_temperature, - receiver_density, + ion_density, z_effective, b_field) diff --git a/cherab/core/model/plasma/__init__.pxd b/cherab/core/model/plasma/__init__.pxd index 6971f828..c398fd30 100644 --- a/cherab/core/model/plasma/__init__.pxd +++ b/cherab/core/model/plasma/__init__.pxd @@ -20,4 +20,5 @@ from cherab.core.model.plasma.bremsstrahlung cimport Bremsstrahlung from cherab.core.model.plasma.impact_excitation cimport ExcitationLine from cherab.core.model.plasma.recombination cimport RecombinationLine +from cherab.core.model.plasma.thermal_cx cimport ThermalCXLine from cherab.core.model.plasma.total_radiated_power cimport TotalRadiatedPower diff --git a/cherab/core/model/plasma/__init__.py b/cherab/core/model/plasma/__init__.py index 52766032..436060b5 100644 --- a/cherab/core/model/plasma/__init__.py +++ b/cherab/core/model/plasma/__init__.py @@ -20,4 +20,5 @@ from .bremsstrahlung import Bremsstrahlung from .impact_excitation import ExcitationLine from .recombination import RecombinationLine +from .thermal_cx import ThermalCXLine from .total_radiated_power import TotalRadiatedPower diff --git a/cherab/core/model/plasma/thermal_cx.pxd b/cherab/core/model/plasma/thermal_cx.pxd new file mode 100644 index 00000000..636fddaa --- /dev/null +++ b/cherab/core/model/plasma/thermal_cx.pxd @@ -0,0 +1,35 @@ +# Copyright 2016-2018 Euratom +# Copyright 2016-2018 United Kingdom Atomic Energy Authority +# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# +# Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the +# European Commission - subsequent versions of the EUPL (the "Licence"); +# You may not use this work except in compliance with the Licence. +# You may obtain a copy of the Licence at: +# +# https://joinup.ec.europa.eu/software/page/eupl5 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the Licence is distributed on an "AS IS" basis, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. +# +# See the Licence for the specific language governing permissions and limitations +# under the Licence. + +from cherab.core.atomic cimport Line +from cherab.core.plasma cimport PlasmaModel +from cherab.core.species cimport Species +from cherab.core.model.lineshape cimport LineShapeModel + + +cdef class ThermalCXLine(PlasmaModel): + + cdef: + Line _line + double _wavelength + Species _target_species + list _rates + LineShapeModel _lineshape + object _lineshape_class, _lineshape_args, _lineshape_kwargs + + cdef int _populate_cache(self) except -1 diff --git a/cherab/core/model/plasma/thermal_cx.pyx b/cherab/core/model/plasma/thermal_cx.pyx new file mode 100644 index 00000000..1325919c --- /dev/null +++ b/cherab/core/model/plasma/thermal_cx.pyx @@ -0,0 +1,163 @@ +# Copyright 2016-2018 Euratom +# Copyright 2016-2018 United Kingdom Atomic Energy Authority +# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# +# Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the +# European Commission - subsequent versions of the EUPL (the "Licence"); +# You may not use this work except in compliance with the Licence. +# You may obtain a copy of the Licence at: +# +# https://joinup.ec.europa.eu/software/page/eupl5 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the Licence is distributed on an "AS IS" basis, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. +# +# See the Licence for the specific language governing permissions and limitations +# under the Licence. + +from raysect.optical cimport Spectrum, Point3D, Vector3D +from cherab.core cimport Plasma, AtomicData +from cherab.core.atomic cimport ThermalCXPEC +from cherab.core.model.lineshape cimport GaussianLine, LineShapeModel +from cherab.core.utility.constants cimport RECIP_4_PI + + +cdef class ThermalCXLine(PlasmaModel): + r""" + Emitter that calculates spectral line emission from a plasma object + as a result of thermal charge exchange of the target species with the donor species. + + .. math:: + \epsilon_{\mathrm{CX}}(\lambda) = \frac{1}{4 \pi} n_{Z_\mathrm{i} + 1} + \sum_j{n_{Z_\mathrm{j}} \mathrm{PEC}_{\mathrm{cx}}(n_\mathrm{e}, T_\mathrm{e}, T_{Z_\mathrm{j}})} + f(\lambda), + + where :math:`n_{Z_\mathrm{i} + 1}` is the receiver species density, + :math:`n_{Z_\mathrm{j}}` is the donor species density, + :math:`\mathrm{PEC}_{\mathrm{cx}}` is the thermal CX photon emission coefficient + for the specified spectral line of the :math:`Z_\mathrm{i}` ion, + :math:`T_{Z_\mathrm{j}}` is the donor species temperature, + :math:`f(\lambda)` is the normalised spectral line shape, + + :param Line line: Spectroscopic emission line object. + :param Plasma plasma: The plasma to which this emission model is attached. Default is None. + :param AtomicData atomic_data: The atomic data provider for this model. Default is None. + :param object lineshape: Line shape model class. Default is None (GaussianLine). + :param object lineshape_args: A list of line shape model arguments. Default is None. + :param object lineshape_kwargs: A dictionary of line shape model keyword arguments. Default is None. + + :ivar Plasma plasma: The plasma to which this emission model is attached. + :ivar AtomicData atomic_data: The atomic data provider for this model. + """ + + def __init__(self, Line line, Plasma plasma=None, AtomicData atomic_data=None, object lineshape=None, + object lineshape_args=None, object lineshape_kwargs=None): + + super().__init__(plasma, atomic_data) + + self._line = line + + self._lineshape_class = lineshape or GaussianLine + if not issubclass(self._lineshape_class, LineShapeModel): + raise TypeError("The attribute lineshape must be a subclass of LineShapeModel.") + + if lineshape_args: + self._lineshape_args = lineshape_args + else: + self._lineshape_args = [] + if lineshape_kwargs: + self._lineshape_kwargs = lineshape_kwargs + else: + self._lineshape_kwargs = {} + + # ensure that cache is initialised + self._change() + + def __repr__(self): + return ''.format(self._line.element.name, self._line.charge, self._line.transition) + + cpdef Spectrum emission(self, Point3D point, Vector3D direction, Spectrum spectrum): + + cdef: + double ne, te, receiver_density, donor_density, donor_temperature, weighted_rate, radiance + Species species + ThermalCXPEC rate + + # cache data on first run + if self._target_species is None: + self._populate_cache() + + ne = self._plasma.get_electron_distribution().density(point.x, point.y, point.z) + if ne <= 0.0: + return spectrum + + te = self._plasma.get_electron_distribution().effective_temperature(point.x, point.y, point.z) + if te <= 0.0: + return spectrum + + receiver_density = self._target_species.distribution.density(point.x, point.y, point.z) + if receiver_density <= 0.0: + return spectrum + + # obtain composite CX PEC by iterating over all possible CX donors + weighted_rate = 0 + for species, rate in self._rates: + donor_density = species.distribution.density(point.x, point.y, point.z) + donor_temperature = species.distribution.effective_temperature(point.x, point.y, point.z) + weighted_rate += donor_density * rate.evaluate(ne, te, donor_temperature) + + # add emission line to spectrum + radiance = RECIP_4_PI * weighted_rate * receiver_density + return self._lineshape.add_line(radiance, point, direction, spectrum) + + cdef int _populate_cache(self) except -1: + + cdef: + int receiver_charge + Species species + ThermalCXPEC rate + + # sanity checks + if self._plasma is None: + raise RuntimeError("The emission model is not connected to a plasma object.") + if self._atomic_data is None: + raise RuntimeError("The emission model is not connected to an atomic data source.") + + if self._line is None: + raise RuntimeError("The emission line has not been set.") + + # locate target species + receiver_charge = self._line.charge + 1 + try: + self._target_species = self._plasma.composition.get(self._line.element, receiver_charge) + except ValueError: + raise RuntimeError("The plasma object does not contain the ion species for the specified CX line " + "(element={}, ionisation={}).".format(self._line.element.symbol, receiver_charge)) + + # obtain rate functions + self._rates = [] + # iterate over all posible electron donors in plasma composition + # and for each donor, cache the PEC rate function for the CX reaction with this receiver + for species in self._plasma.composition: + # exclude the receiver species from the list of donors and omit fully ionised species + if species != self._target_species and species.charge < species.element.atomic_number: + rate = self._atomic_data.thermal_cx_pec(species.element, species.charge, # donor + self._line.element, receiver_charge, # receiver + self._line.transition) + self._rates.append((species, rate)) + + # identify wavelength + self._wavelength = self._atomic_data.wavelength(self._line.element, self._line.charge, self._line.transition) + + # instance line shape renderer + self._lineshape = self._lineshape_class(self._line, self._wavelength, self._target_species, self._plasma, + *self._lineshape_args, **self._lineshape_kwargs) + + def _change(self): + + # clear cache to force regeneration on first use + self._target_species = None + self._wavelength = 0.0 + self._rates = None + self._lineshape = None diff --git a/cherab/core/tests/test_beam.py b/cherab/core/tests/test_beam.py index ef10126f..7d2ec203 100644 --- a/cherab/core/tests/test_beam.py +++ b/cherab/core/tests/test_beam.py @@ -55,32 +55,38 @@ def beam_stopping_rate(self, beam_ion, plasma_ion, charge): class TestBeam(unittest.TestCase): - atomic_data = MockAtomicData() - - world = World() - - plasma_density = 1.e19 - plasma_temperature = 1.e3 - plasma_species = [(deuterium, 1, plasma_density, plasma_temperature, Vector3D(0, 0, 0))] - plasma = build_constant_slab_plasma(length=1, width=1, height=1, electron_density=plasma_density, - electron_temperature=plasma_temperature, - plasma_species=plasma_species) - plasma.atomic_data = atomic_data - plasma.parent = world - - beam = Beam(transform=translate(0.5, 0, 0)) - beam.atomic_data = atomic_data - beam.plasma = plasma - beam.attenuator = SingleRayAttenuator(clamp_to_zero=True) - beam.energy = 50000 - beam.power = 1e6 - beam.temperature = 10 - beam.element = deuterium - beam.parent = world - beam.sigma = 0.2 - beam.divergence_x = 1. - beam.divergence_y = 2. - beam.length = 10. + def setUp(self): + + self.atomic_data = MockAtomicData() + + self.world = World() + + self.plasma_density = 1.e19 + self.plasma_temperature = 1.e3 + plasma_species = [(deuterium, 1, self.plasma_density, self.plasma_temperature, Vector3D(0, 0, 0))] + plasma = build_constant_slab_plasma(length=1, width=1, height=1, + electron_density=self.plasma_density, + electron_temperature=self.plasma_temperature, + plasma_species=plasma_species) + plasma.atomic_data = self.atomic_data + plasma.parent = self.world + + beam = Beam(transform=translate(0.5, 0, 0)) + beam.atomic_data = self.atomic_data + beam.plasma = plasma + beam.attenuator = SingleRayAttenuator(clamp_to_zero=True) + beam.energy = 50000 + beam.power = 1e6 + beam.temperature = 10 + beam.element = deuterium + beam.parent = self.world + beam.sigma = 0.2 + beam.divergence_x = 1. + beam.divergence_y = 2. + beam.length = 10. + + self.plasma = plasma + self.beam = beam def test_beam_density(self): diff --git a/cherab/core/tests/test_beamcxline.py b/cherab/core/tests/test_beamcxline.py index e8336710..a93cf946 100644 --- a/cherab/core/tests/test_beamcxline.py +++ b/cherab/core/tests/test_beamcxline.py @@ -36,7 +36,7 @@ class ConstantBeamCXPEC(BeamCXPEC): """ def __init__(self, donor_metastable, value): - self.donor_metastable = donor_metastable + super().__init__(donor_metastable) self.value = value def evaluate(self, energy, temperature, density, z_effective, b_field): @@ -76,25 +76,33 @@ def wavelength(self, ion, charge, transition): class TestBeamCXLine(unittest.TestCase): - world = World() - - atomic_data = MockAtomicData() - - plasma_species = [(deuterium, 1, 1.e19, 200., Vector3D(0, 0, 0))] - plasma = build_constant_slab_plasma(length=1, width=1, height=1, electron_density=1e19, electron_temperature=200., - plasma_species=plasma_species, b_field=Vector3D(0, 10., 0)) - plasma.atomic_data = atomic_data - plasma.parent = world - - beam = Beam(transform=translate(0.5, 0, 0)) - beam.atomic_data = atomic_data - beam.plasma = plasma - beam.attenuator = SingleRayAttenuator(clamp_to_zero=True) - beam.energy = 50000 - beam.power = 1e6 - beam.temperature = 10 - beam.element = deuterium - beam.parent = world + def setUp(self): + + self.world = World() + + self.atomic_data = MockAtomicData() + + plasma_species = [(deuterium, 1, 1.e19, 200., Vector3D(0, 0, 0))] + plasma = build_constant_slab_plasma(length=1, width=1, height=1, + electron_density=1e19, + electron_temperature=200., + plasma_species=plasma_species, + b_field=Vector3D(0, 10., 0)) + plasma.atomic_data = self.atomic_data + plasma.parent = self.world + + beam = Beam(transform=translate(0.5, 0, 0)) + beam.atomic_data = self.atomic_data + beam.plasma = plasma + beam.attenuator = SingleRayAttenuator(clamp_to_zero=True) + beam.energy = 50000 + beam.power = 1e6 + beam.temperature = 10 + beam.element = deuterium + beam.parent = self.world + + self.plasma = plasma + self.beam = beam def test_default_lineshape(self): # setting up the model diff --git a/cherab/core/tests/test_bremsstrahlung.py b/cherab/core/tests/test_bremsstrahlung.py index 5373776b..a77a1da5 100644 --- a/cherab/core/tests/test_bremsstrahlung.py +++ b/cherab/core/tests/test_bremsstrahlung.py @@ -34,13 +34,18 @@ class TestBremsstrahlung(unittest.TestCase): - world = World() - - plasma_species = [(deuterium, 1, 1.e19, 2000., Vector3D(0, 0, 0)), (nitrogen, 7, 1.e18, 2000., Vector3D(0, 0, 0))] - plasma = build_constant_slab_plasma(length=1, width=1, height=1, electron_density=1e19, electron_temperature=2000., - plasma_species=plasma_species) - plasma.parent = world - plasma.atomic_data = AtomicData() + def setUp(self): + + self.world = World() + + plasma_species = [(deuterium, 1, 1.e19, 2000., Vector3D(0, 0, 0)), + (nitrogen, 7, 1.e18, 2000., Vector3D(0, 0, 0))] + self.plasma = build_constant_slab_plasma(length=1, width=1, height=1, + electron_density=1e19, + electron_temperature=2000., + plasma_species=plasma_species) + self.plasma.parent = self.world + self.plasma.atomic_data = AtomicData() def test_bremsstrahlung_model(self): # setting up the model diff --git a/cherab/core/tests/test_line_emission.py b/cherab/core/tests/test_line_emission.py new file mode 100644 index 00000000..b55af871 --- /dev/null +++ b/cherab/core/tests/test_line_emission.py @@ -0,0 +1,324 @@ +# Copyright 2016-2023 Euratom +# Copyright 2016-2023 United Kingdom Atomic Energy Authority +# Copyright 2016-2023 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# +# Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the +# European Commission - subsequent versions of the EUPL (the "Licence"); +# You may not use this work except in compliance with the Licence. +# You may obtain a copy of the Licence at: +# +# https://joinup.ec.europa.eu/software/page/eupl5 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the Licence is distributed on an "AS IS" basis, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. +# +# See the Licence for the specific language governing permissions and limitations +# under the Licence. + +import unittest + +import numpy as np + +from raysect.core import Point3D, Vector3D, translate +from raysect.optical import World, Spectrum, Ray + +from cherab.core.atomic import Line, AtomicData, ImpactExcitationPEC, RecombinationPEC, ThermalCXPEC +from cherab.core.atomic import deuterium, carbon +from cherab.tools.plasmas.slab import build_constant_slab_plasma +from cherab.core.model import ExcitationLine, RecombinationLine, ThermalCXLine, GaussianLine, ZeemanTriplet + + +class ConstantImpactExcitationPEC(ImpactExcitationPEC): + """ + Constant electron impact excitation PEC for test purpose. + """ + + def __init__(self, value): + self.value = value + + def evaluate(self, density, temperature): + + return self.value + + +class ConstantRecombinationPEC(RecombinationPEC): + """ + Constant recombination PEC for test purpose. + """ + + def __init__(self, value): + self.value = value + + def evaluate(self, density, temperature): + + return self.value + + +class ConstantThermalCXPEC(ThermalCXPEC): + """ + Constant recombination PEC for test purpose. + """ + + def __init__(self, value): + self.value = value + + def evaluate(self, electron_density, electron_temperature, donor_temperature): + + return self.value + + +class MockAtomicData(AtomicData): + """Fake atomic data for test purpose.""" + + def impact_excitation_pec(self, ion, charge, transition): + + return ConstantImpactExcitationPEC(1.4e-39) + + def recombination_pec(self, ion, charge, transition): + + return ConstantRecombinationPEC(8.e-41) + + def thermal_cx_pec(self, donor_ion, donor_charge, receiver_ion, receiver_charge, transition): + + return ConstantThermalCXPEC(1.2e-46) + + def wavelength(self, ion, charge, transition): + + return 529.27 + + +class TestExcitationLine(unittest.TestCase): + + def setUp(self): + + self.world = World() + + self.atomic_data = MockAtomicData() + + plasma_species = [(carbon, 5, 2.e18, 800., Vector3D(0, 0, 0))] + self.slab_length = 1.2 + self.plasma = build_constant_slab_plasma(length=self.slab_length, width=1, height=1, + electron_density=1e19, electron_temperature=1000., + plasma_species=plasma_species, b_field=Vector3D(0, 10., 0)) + self.plasma.atomic_data = self.atomic_data + self.plasma.parent = self.world + + def test_default_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [ExcitationLine(line)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + excit_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.impact_excitation_pec(line.element, line.charge, line.transition)(ne, te) + target_species = self.plasma.composition.get(line.element, line.charge) + ni = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * ni * ne * self.slab_length + + gaussian_line = GaussianLine(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = gaussian_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(excit_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='ExcitationLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + def test_custom_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [ExcitationLine(line, lineshape=ZeemanTriplet)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + excit_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.impact_excitation_pec(line.element, line.charge, line.transition)(ne, te) + target_species = self.plasma.composition.get(line.element, line.charge) + ni = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * ni * ne * self.slab_length + + zeeman_line = ZeemanTriplet(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = zeeman_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(excit_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='ExcitationLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + +class TestRecombinationLine(unittest.TestCase): + + def setUp(self): + + self.world = World() + + self.atomic_data = MockAtomicData() + + plasma_species = [(carbon, 6, 1.67e18, 800., Vector3D(0, 0, 0))] + self.slab_length = 1.2 + self.plasma = build_constant_slab_plasma(length=self.slab_length, width=1, height=1, + electron_density=1e19, electron_temperature=1000., + plasma_species=plasma_species, b_field=Vector3D(0, 10., 0)) + self.plasma.atomic_data = self.atomic_data + self.plasma.parent = self.world + + def test_default_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [RecombinationLine(line)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + recomb_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.recombination_pec(line.element, line.charge, line.transition)(ne, te) + target_species = self.plasma.composition.get(line.element, line.charge + 1) + ni = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * ni * ne * self.slab_length + + gaussian_line = GaussianLine(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = gaussian_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(recomb_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='RecombinationLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + def test_custom_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [RecombinationLine(line, lineshape=ZeemanTriplet)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + recomb_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.recombination_pec(line.element, line.charge, line.transition)(ne, te) + target_species = self.plasma.composition.get(line.element, line.charge + 1) + ni = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * ni * ne * self.slab_length + + zeeman_line = ZeemanTriplet(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = zeeman_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(recomb_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='RecombinationLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + +class TestThermalCXLine(unittest.TestCase): + + def setUp(self): + + self.world = World() + + self.atomic_data = MockAtomicData() + + plasma_species = [(carbon, 6, 1.67e18, 800., Vector3D(0, 0, 0)), + (deuterium, 0, 1.e19, 100., Vector3D(0, 0, 0))] + self.slab_length = 1.2 + self.plasma = build_constant_slab_plasma(length=self.slab_length, width=1, height=1, + electron_density=1e19, electron_temperature=1000., + plasma_species=plasma_species, b_field=Vector3D(0, 10., 0)) + self.plasma.atomic_data = self.atomic_data + self.plasma.parent = self.world + + def test_default_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [ThermalCXLine(line)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + thermalcx_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + donor_species = self.plasma.composition.get(deuterium, 0) + donor_density = donor_species.distribution.density(0.5, 0, 0) + donor_temperature = donor_species.distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.thermal_cx_pec(deuterium, 0, line.element, line.charge, line.transition)(ne, te, donor_temperature) + target_species = self.plasma.composition.get(line.element, line.charge + 1) + receiver_density = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * receiver_density * donor_density * self.slab_length + + gaussian_line = GaussianLine(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = gaussian_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(thermalcx_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='ThermalCXLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + def test_custom_lineshape(self): + # setting up the model + line = Line(carbon, 5, (8, 7)) + self.plasma.models = [ThermalCXLine(line, lineshape=ZeemanTriplet)] + wavelength = self.atomic_data.wavelength(line.element, line.charge, line.transition) + + # observing + origin = Point3D(1.5, 0, 0) + direction = Vector3D(-1, 0, 0) + ray = Ray(origin=origin, direction=direction, + min_wavelength=wavelength - 1.5, max_wavelength=wavelength + 1.5, bins=512) + thermalcx_spectrum = ray.trace(self.world) + + # validating + ne = self.plasma.electron_distribution.density(0.5, 0, 0) # constant slab + te = self.plasma.electron_distribution.effective_temperature(0.5, 0, 0) + donor_species = self.plasma.composition.get(deuterium, 0) + donor_density = donor_species.distribution.density(0.5, 0, 0) + donor_temperature = donor_species.distribution.effective_temperature(0.5, 0, 0) + rate = self.atomic_data.thermal_cx_pec(deuterium, 0, line.element, line.charge, line.transition)(ne, te, donor_temperature) + target_species = self.plasma.composition.get(line.element, line.charge + 1) + receiver_density = target_species.distribution.density(0.5, 0, 0) + radiance = 0.25 / np.pi * rate * receiver_density * donor_density * self.slab_length + + zeeman_line = ZeemanTriplet(line, wavelength, target_species, self.plasma) + spectrum = Spectrum(ray.min_wavelength, ray.max_wavelength, ray.bins) + spectrum = zeeman_line.add_line(radiance, Point3D(0.5, 0, 0), direction, spectrum) + + for i in range(ray.bins): + self.assertAlmostEqual(thermalcx_spectrum.samples[i], spectrum.samples[i], delta=1e-8, + msg='ThermalCXLine model gives a wrong value at {} nm.'.format(spectrum.wavelengths[i])) + + +if __name__ == '__main__': + unittest.main() diff --git a/cherab/core/tests/test_lineshapes.py b/cherab/core/tests/test_lineshapes.py index 6f300e59..b6827256 100644 --- a/cherab/core/tests/test_lineshapes.py +++ b/cherab/core/tests/test_lineshapes.py @@ -43,16 +43,21 @@ class TestLineShapes(unittest.TestCase): - plasma_species = [(deuterium, 0, 1.e18, 5., Vector3D(2.e4, 0, 0)), - (nitrogen, 1, 1.e17, 10., Vector3D(1.e4, 5.e4, 0))] - plasma = build_constant_slab_plasma(length=1, width=1, height=1, electron_density=1e19, electron_temperature=20., - plasma_species=plasma_species, b_field=Vector3D(0, 5., 0)) - atomic_data = AtomicData() - beam = Beam() - beam.plasma = plasma - beam.energy = 60000 - beam.temperature = 10 - beam.element = deuterium + def setUp(self): + + plasma_species = [(deuterium, 0, 1.e18, 5., Vector3D(2.e4, 0, 0)), + (nitrogen, 1, 1.e17, 10., Vector3D(1.e4, 5.e4, 0))] + self.plasma = build_constant_slab_plasma(length=1, width=1, height=1, + electron_density=1e19, + electron_temperature=20., + plasma_species=plasma_species, + b_field=Vector3D(0, 5., 0)) + self.atomic_data = AtomicData() + self.beam = Beam() + self.beam.plasma = self.plasma + self.beam.energy = 60000 + self.beam.temperature = 10 + self.beam.element = deuterium def test_gaussian_line(self): # setting up a line shape model diff --git a/cherab/openadas/install.py b/cherab/openadas/install.py index 27e0c5f4..53e68899 100644 --- a/cherab/openadas/install.py +++ b/cherab/openadas/install.py @@ -1,7 +1,7 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2023 Euratom +# Copyright 2016-2023 United Kingdom Atomic Energy Authority +# Copyright 2016-2023 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -21,9 +21,11 @@ import os import urllib.parse import urllib.request +import numpy as np from cherab.openadas import repository from cherab.openadas.parse import * +from cherab.core.atomic import hydrogen from cherab.core.utility import RecursiveDict, PerCm3ToPerM3, Cm3ToM3 @@ -264,19 +266,29 @@ def install_adf15(element, ionisation, file_path, download=False, repository_pat # decode file and write out rates rates, wavelengths = parse_adf15(element, ionisation, path, header_format=header_format) + + if 'thermalcx' in rates: + cx_rates = rates.pop('thermalcx') + # CX rates for Tdon = Trec (2D function of Ne, Te) + # converting to 3D function of Ne, Te, Tdon + repository.update_pec_thermal_cx_rates(_thermalcx_adf15_2dto3d_converter(cx_rates)) + repository.update_pec_rates(rates, repository_path) repository.update_wavelengths(wavelengths, repository_path) def install_adf21(beam_species, target_ion, target_charge, file_path, download=False, repository_path=None, adas_path=None): - # """ - # Adds the rate defined in an ADF21 file to the repository. - # - # :param file_path: Path relative to ADAS root. - # :param download: Attempt to download file if not present (Default=True). - # :param repository_path: Path to the repository in which to install the rates (optional). - # :param adas_path: Path to ADAS files repository (optional). - # """ + """ + Adds the beam stopping rate defined in an ADF21 file to the repository. + + :param beam_species: Beam neutral atom (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param file_path: Path relative to ADAS root. + :param download: Attempt to download file if not present (Default=True). + :param repository_path: Path to the repository in which to install the rates (optional). + :param adas_path: Path to ADAS files repository (optional). + """ print('Installing {}...'.format(file_path)) path = _locate_adas_file(file_path, download, adas_path, repository_path) @@ -289,15 +301,18 @@ def install_adf21(beam_species, target_ion, target_charge, file_path, download=F def install_adf22bmp(beam_species, beam_metastable, target_ion, target_charge, file_path, download=False, repository_path=None, adas_path=None): - pass - # """ - # Adds the rate defined in an ADF21 file to the repository. - # - # :param file_path: Path relative to ADAS root. - # :param download: Attempt to download file if not present (Default=True). - # :param repository_path: Path to the repository in which to install the rates (optional). - # :param adas_path: Path to ADAS files repository (optional). - # """ + """ + Adds the beam population rate defined in an ADF22 BMP file to the repository. + + :param beam_species: Beam neutral atom (Element/Isotope). + :param beam_metastable: Metastable/excitation level of beam neutral atom. + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param file_path: Path relative to ADAS root. + :param download: Attempt to download file if not present (Default=True). + :param repository_path: Path to the repository in which to install the rates (optional). + :param adas_path: Path to ADAS files repository (optional). + """ print('Installing {}...'.format(file_path)) path = _locate_adas_file(file_path, download, adas_path, repository_path) @@ -310,15 +325,18 @@ def install_adf22bmp(beam_species, beam_metastable, target_ion, target_charge, f def install_adf22bme(beam_species, target_ion, target_charge, transition, file_path, download=False, repository_path=None, adas_path=None): - pass - # """ - # Adds the rate defined in an ADF21 file to the repository. - # - # :param file_path: Path relative to ADAS root. - # :param download: Attempt to download file if not present (Default=True). - # :param repository_path: Path to the repository in which to install the rates (optional). - # :param adas_path: Path to ADAS files repository (optional). - # """ + """ + Adds the beam emission rate defined in an ADF22 BME file to the repository. + + :param beam_species: Beam neutral atom (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param transition: Tuple containing (initial level, final level). + :param file_path: Path relative to ADAS root. + :param download: Attempt to download file if not present (Default=True). + :param repository_path: Path to the repository in which to install the rates (optional). + :param adas_path: Path to ADAS files repository (optional). + """ print('Installing {}...'.format(file_path)) path = _locate_adas_file(file_path, download, adas_path, repository_path) @@ -393,3 +411,23 @@ def _notation_adf11_adas2cherab(rate_adas, filetype): return rate_cherab +def _thermalcx_adf15_2dto3d_converter(rates): + """ + Converts thermal CX PEC rates parsed from a standard ADF 15 file + to the format supported by the repository. + + In the standard ADF 15 file, it is assumed that the donor is H0 and Tdon = Trec. + """ + new_rates = RecursiveDict() + for element, charge_states in rates.items(): + for charge, transitions in charge_states.items(): + for transition, rate in transitions.items(): + data = np.empty((len(rate['ne']), len(rate['te']), 2)) + data[:, :, :] = rate['rate'][:, :, None] + new_rate = {'ne': rate['ne'], + 'te': rate['te'], + 'td': np.array([0.01, 10000]), + 'rate': data} + new_rates[hydrogen][0][element][charge + 1][transition] = new_rate + + return new_rates diff --git a/cherab/openadas/openadas.py b/cherab/openadas/openadas.py index 7e850bc4..5d42935a 100644 --- a/cherab/openadas/openadas.py +++ b/cherab/openadas/openadas.py @@ -386,6 +386,48 @@ def recombination_pec(self, ion, charge, transition): return RecombinationPEC(wavelength, data, extrapolate=self._permit_extrapolation) + def thermal_cx_pec(self, donor_element, donor_charge, receiver_element, receiver_charge, transition): + """ + Thermal CX photon emission coefficient for a given species. + + Open-ADAS data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used when permit_extrapolation is True. + + :param donor_element: Element object defining the donor ion type. + :param donor_charge: Charge state of the donor ion. + :param receiver_element: Element object defining the receiver ion type. + :param receiver_charge: Charge state of the receiver ion. + :param transition: Tuple containing (initial level, final level) of the receiver + in charge state receiver_charge - 1. + :return: Thermal charge exchange photon emission coefficient in W.m^3 + as a function of electron density, electron temperature and donor temperature. + """ + + # extract elements from isotopes because there are no isotope rates in ADAS + if isinstance(donor_element, Isotope): + donor_element = donor_element.element + + if isinstance(receiver_element, Isotope): + receiver_element = receiver_element.element + + try: + # read thermal CX rate from json file in the repository + data = repository.get_pec_thermal_cx_rate(donor_element, donor_charge, + receiver_element, receiver_charge, + transition, + repository_path=self._data_path) + + except RuntimeError: + if self._missing_rates_return_null: + return NullThermalCXPEC() + raise + + # obtain isotope's rest wavelength for a given transition + # the wavelength is used ot convert the PEC from photons/s/m3 to W/m3 + wavelength = self.wavelength(receiver_element, receiver_charge - 1, transition) + + return ThermalCXPEC(wavelength, data, extrapolate=self._permit_extrapolation) + def line_radiated_power_rate(self, ion, charge): """ Line radiated power coefficient for a given species. diff --git a/cherab/openadas/parse/adf15.py b/cherab/openadas/parse/adf15.py index 51bbcf6b..269edd0a 100644 --- a/cherab/openadas/parse/adf15.py +++ b/cherab/openadas/parse/adf15.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2023 Euratom +# Copyright 2016-2023 United Kingdom Atomic Energy Authority +# Copyright 2016-2023 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -126,7 +126,7 @@ def _scrape_metadata_hydrogen(file, element, charge): elif rate_type_adas == 'RECOM': rate_type = 'recombination' elif rate_type_adas == 'CHEXC': - rate_type = 'cx_thermal' + rate_type = 'thermalcx' else: raise ValueError("Unrecognised rate type - {}".format(rate_type_adas)) @@ -169,7 +169,7 @@ def _scrape_metadata_hydrogen_like(file, element, charge): elif rate_type_adas == 'RECOM': rate_type = 'recombination' elif rate_type_adas == 'CHEXC': - rate_type = 'cx_thermal' + rate_type = 'thermalcx' else: raise ValueError("Unrecognised rate type - {}".format(rate_type_adas)) @@ -237,7 +237,7 @@ def _scrape_metadata_full(file, element, charge): elif rate_type_adas == 'RECOM': rate_type = 'recombination' elif rate_type_adas == 'CHEXC': - rate_type = 'cx_thermal' + rate_type = 'thermalcx' else: raise ValueError("Unrecognised rate type - {}".format(rate_type_adas)) diff --git a/cherab/openadas/rates/atomic.pyx b/cherab/openadas/rates/atomic.pyx index 04500124..eecb9bf5 100644 --- a/cherab/openadas/rates/atomic.pyx +++ b/cherab/openadas/rates/atomic.pyx @@ -1,7 +1,7 @@ -# Copyright 2016-2021 Euratom -# Copyright 2016-2021 United Kingdom Atomic Energy Authority -# Copyright 2016-2021 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -24,12 +24,26 @@ from raysect.core.math.function.float cimport Interpolator2DArray cdef class IonisationRate(CoreIonisationRate): + """ + Ionisation rate. + + Data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param dict data: Ionisation rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with ionisation rate in m^3.s^-1. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, dict data, extrapolate=False): - """ - :param data: Dictionary containing rate data. - :param extrapolate: Enable extrapolation (default=False). - """ self.raw_data = data @@ -50,11 +64,8 @@ cdef class IonisationRate(CoreIonisationRate): cpdef double evaluate(self, double density, double temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(density), log10(temperature)) @@ -62,7 +73,7 @@ cdef class IonisationRate(CoreIonisationRate): cdef class NullIonisationRate(CoreIonisationRate): """ - A PEC rate that always returns zero. + An ionisation rate that always returns zero. Needed for use cases where the required atomic data is missing. """ @@ -71,12 +82,26 @@ cdef class NullIonisationRate(CoreIonisationRate): cdef class RecombinationRate(CoreRecombinationRate): + """ + Recombination rate. + + Data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param dict data: Recombination rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with recombination rate in m^3.s^-1. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, dict data, extrapolate=False): - """ - :param data: Dictionary containing rate data. - :param extrapolate: Enable extrapolation (default=False). - """ self.raw_data = data @@ -97,11 +122,8 @@ cdef class RecombinationRate(CoreRecombinationRate): cpdef double evaluate(self, double density, double temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(density), log10(temperature)) @@ -109,7 +131,7 @@ cdef class RecombinationRate(CoreRecombinationRate): cdef class NullRecombinationRate(CoreRecombinationRate): """ - A PEC rate that always returns zero. + A recombination rate that always returns zero. Needed for use cases where the required atomic data is missing. """ @@ -118,12 +140,26 @@ cdef class NullRecombinationRate(CoreRecombinationRate): cdef class ThermalCXRate(CoreThermalCXRate): + """ + Thermal charge exchange rate. + + Data is interpolated with cubic spline in log-log space. + Linear extrapolation is used if extrapolate is True. + + :param dict data: CX rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with thermal CX rate in m^3.s^-1. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, dict data, extrapolate=False): - """ - :param data: Dictionary containing rate data. - :param extrapolate: Enable extrapolation (default=False). - """ self.raw_data = data @@ -143,11 +179,8 @@ cdef class ThermalCXRate(CoreThermalCXRate): cpdef double evaluate(self, double density, double temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(density), log10(temperature)) @@ -155,7 +188,7 @@ cdef class ThermalCXRate(CoreThermalCXRate): cdef class NullThermalCXRate(CoreThermalCXRate): """ - A PEC rate that always returns zero. + A thermal CX rate that always returns zero. Needed for use cases where the required atomic data is missing. """ diff --git a/cherab/openadas/rates/beam.pyx b/cherab/openadas/rates/beam.pyx index f40af8dd..30d43664 100644 --- a/cherab/openadas/rates/beam.pyx +++ b/cherab/openadas/rates/beam.pyx @@ -1,6 +1,6 @@ -# Copyright 2016-2021 Euratom -# Copyright 2016-2021 United Kingdom Atomic Energy Authority -# Copyright 2016-2021 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -31,8 +31,26 @@ cdef class BeamStoppingRate(CoreBeamStoppingRate): """ The beam stopping coefficient interpolation class. - :param data: A dictionary holding the beam coefficient data. - :param extrapolate: Set to True to enable extrapolation, False to disable (default). + Data is interpolated with cubic spline in log-log space. + Linear and quadratic extrapolations are used for "sen" and "st" respectively + if extrapolate is True. + + :param dict data: A beam stopping rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n': 1D array of size (M) with target electron density in m^-3, + | 't': 1D array of size (K) with target electron temperature in eV, + | 'sen': 2D array of size (N, M) with beam stopping rate energy component in m^3.s^-1. + | 'st': 1D array of size (K) with beam stopping rate temperature component in m^3.s^-1. + | 'sref': reference beam stopping rate in m^3.s^-1. + | The total beam stopping rate: s = sen * st / sref. + + :param bint extrapolate: Set to True to enable extrapolation, False to disable (default). + + :ivar tuple beam_energy_range: Interaction energy interpolation range. + :ivar tuple density_range: Target electron density interpolation range. + :ivar tuple temperature_range: Target electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. """ @cython.cdivision(True) @@ -78,14 +96,8 @@ cdef class BeamStoppingRate(CoreBeamStoppingRate): """ # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if energy < 1.e-300: - energy = 1.e-300 - - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if energy <= 0 or density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** (self._npl_eb.evaluate(log10(energy), log10(density)) + self._tp.evaluate(log10(temperature))) @@ -93,7 +105,7 @@ cdef class BeamStoppingRate(CoreBeamStoppingRate): cdef class NullBeamStoppingRate(CoreBeamStoppingRate): """ - A beam rate that always returns zero. + A beam stopping rate that always returns zero. Needed for use cases where the required atomic data is missing. """ @@ -105,8 +117,26 @@ cdef class BeamPopulationRate(CoreBeamPopulationRate): """ The beam population coefficient interpolation class. - :param data: A dictionary holding the beam coefficient data. - :param extrapolate: Set to True to enable extrapolation, False to disable (default). + Data is interpolated with cubic spline in log-log space. + Linear and quadratic extrapolations are used for "sen" and "st" respectively + if extrapolate is True. + + :param dict data: Beam population rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n': 1D array of size (M) with target electron density in m^-3, + | 't': 1D array of size (K) with target electron temperature in eV, + | 'sen': 2D array of size (N, M) with dimensionless beam population rate energy component. + | 'st': 1D array of size (K) with dimensionless beam population rate temperature component. + | 'sref': reference dimensionless beam population rate. + | The total beam population rate: s = sen * st / sref. + + :param bint extrapolate: Set to True to enable extrapolation, False to disable (default). + + :ivar tuple beam_energy_range: Interaction energy interpolation range. + :ivar tuple density_range: Target electron density interpolation range. + :ivar tuple temperature_range: Target electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. """ @cython.cdivision(True) @@ -152,14 +182,8 @@ cdef class BeamPopulationRate(CoreBeamPopulationRate): """ # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if energy < 1.e-300: - energy = 1.e-300 - - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if energy <= 0 or density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** (self._npl_eb.evaluate(log10(energy), log10(density)) + self._tp.evaluate(log10(temperature))) @@ -167,7 +191,7 @@ cdef class BeamPopulationRate(CoreBeamPopulationRate): cdef class NullBeamPopulationRate(CoreBeamPopulationRate): """ - A beam rate that always returns zero. + A beam population rate that always returns zero. Needed for use cases where the required atomic data is missing. """ @@ -179,9 +203,26 @@ cdef class BeamEmissionPEC(CoreBeamEmissionPEC): """ The beam emission coefficient interpolation class. - :param data: A dictionary holding the beam coefficient data. - :param wavelength: The natural wavelength of the emission line associated with the rate data in nm. - :param extrapolate: Set to True to enable extrapolation, False to disable (default). + Data is interpolated with cubic spline in log-log space. + Linear and quadratic extrapolations are used for "sen" and "st" respectively + if extrapolate is True. + + :param dict data: Beam emission rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n' 1D array of size (M) with target electron density in m^-3, + | 't' 1D array of size (K) with target electron temperature in eV, + | 'sen' 2D array of size (N, M) with beam emission rate energy component in photon.m^3.s^-1. + | 'st' 1D array of size (K) with beam emission rate temperature component in photon.m^3.s^-1. + | 'sref': reference beam emission rate in photon.m^3.s^-1. + + :param double wavelength: The natural wavelength of the emission line associated with the rate data in nm. + :param bint extrapolate: Set to True to enable extrapolation, False to disable (default). + + :ivar tuple beam_energy_range: Interaction energy interpolation range. + :ivar tuple density_range: Target electron density interpolation range. + :ivar tuple temperature_range: Target electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. """ @cython.cdivision(True) @@ -194,7 +235,7 @@ cdef class BeamEmissionPEC(CoreBeamEmissionPEC): e = data["e"] # eV/amu n = data["n"] # m^-3 t = data["t"] # eV - sen = np.log10(PhotonToJ.to(data["sen"], wavelength)) # W.m^3/s + sen = np.log10(PhotonToJ.to(data["sen"], wavelength)) # W.m^3 st = np.log10(data["st"] / data["sref"]) # dimensionless # store limits of data @@ -228,14 +269,8 @@ cdef class BeamEmissionPEC(CoreBeamEmissionPEC): """ # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if energy < 1.e-300: - energy = 1.e-300 - - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if energy <= 0 or density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** (self._npl_eb.evaluate(log10(energy), log10(density)) + self._tp.evaluate(log10(temperature))) @@ -243,7 +278,7 @@ cdef class BeamEmissionPEC(CoreBeamEmissionPEC): cdef class NullBeamEmissionPEC(CoreBeamEmissionPEC): """ - A beam rate that always returns zero. + A beam emission PEC that always returns zero. Needed for use cases where the required atomic data is missing. """ diff --git a/cherab/openadas/rates/cx.pxd b/cherab/openadas/rates/cx.pxd index 393027d0..4afcadfc 100644 --- a/cherab/openadas/rates/cx.pxd +++ b/cherab/openadas/rates/cx.pxd @@ -25,7 +25,6 @@ cdef class BeamCXPEC(CoreBeamCXPEC): cdef readonly: dict raw_data double wavelength - int donor_metastable Function1D _eb, _ti, _ni, _zeff, _b readonly tuple beam_energy_range readonly tuple density_range diff --git a/cherab/openadas/rates/cx.pyx b/cherab/openadas/rates/cx.pyx index cb827a8b..3e0735d9 100644 --- a/cherab/openadas/rates/cx.pyx +++ b/cherab/openadas/rates/cx.pyx @@ -1,6 +1,6 @@ -# Copyright 2016-2021 Euratom -# Copyright 2016-2021 United Kingdom Atomic Energy Authority -# Copyright 2016-2021 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -26,18 +26,48 @@ from raysect.core.math.function.float cimport Interpolator1DArray, Constant1D cdef class BeamCXPEC(CoreBeamCXPEC): """ - The effective cx rate interpolation class. - - :param donor_metastable: The metastable state of the donor species for which the rate data applies. - :param wavelength: The natural wavelength of the emission line associated with the rate data in nm. - :param data: A dictionary holding the rate data. - :param extrapolate: Set to True to enable extrapolation, False to disable (default). + Effective charge exchange photon emission coefficient. + + The data for "qeb" is interpolated with a cubic spline in log-log space. + The data for "qti", "qni", "qz" and "qb" are interpolated with a cubic spline + in linear space. + + Quadratic extrapolation is used for "qeb" and nearest neighbour extrapolation is used for + "qti", "qni", "qz" and "qb" if extrapolate is True. + + :param int donor_metastable: The metastable state of the donor species for which the rate data applies. + :param double wavelength: The natural wavelength of the emission line associated with the rate data in nm. + :param data: Beam CX PEC dictionary containing the following entries: + + | 'eb': 1D array of size (N) with beam energy in eV/amu, + | 'ti': 1D array of size (M) with receiver ion temperature in eV, + | 'ni': 1D array of size (K) with receiver ion density in m^-3, + | 'z': 1D array of size (L) with receiver Z-effective, + | 'b': 1D array of size (J) with magnetic field strength in Tesla, + | 'qeb': 1D array of size (N) with CX PEC energy component in photon.m^3.s-1, + | 'qti': 1D array of size (M) with CX PEC temperature component in photon.m^3.s-1, + | 'qni': 1D array of size (K) with CX PEC density component in photon.m^3.s-1, + | 'qz': 1D array of size (L) with CX PEC Zeff component in photon.m^3.s-1, + | 'qb': 1D array of size (J) with CX PEC B-field component in photon.m^3.s-1, + | 'qref': reference CX PEC in photon.m^3.s-1. + | The total beam CX PEC: q = qeb * qti * qni * qz * qb / qref^4. + + :param bint extrapolate: Set to True to enable extrapolation, False to disable (default). + + :ivar tuple beam_energy_range: Interaction energy interpolation range. + :ivar tuple density_range: Receiver ion density interpolation range. + :ivar tuple temperature_range: Receiver ion temperature interpolation range. + :ivar tuple zeff_range: Z-effective interpolation range. + :ivar tuple b_field_range: Magnetic field strength interpolation range. + :ivar int donor_metastable: The metastable state of the donor species. + :ivar double wavelength: The natural wavelength of the emission line in nm. + :ivar dict raw_data: Dictionary containing the raw data. """ @cython.cdivision(True) def __init__(self, int donor_metastable, double wavelength, dict data, bint extrapolate=False): - self.donor_metastable = donor_metastable + super().__init__(donor_metastable) self.wavelength = wavelength self.raw_data = data @@ -79,7 +109,7 @@ cdef class BeamCXPEC(CoreBeamCXPEC): :param energy: Interaction energy in eV/amu. :param temperature: Receiver ion temperature in eV. - :param density: Receiver ion density in m^-3 + :param density: Plasma total ion density in m^-3 :param z_effective: Plasma Z-effective. :param b_field: Magnetic field magnitude in Tesla. :return: The effective cx rate in W.m^3 @@ -88,8 +118,8 @@ cdef class BeamCXPEC(CoreBeamCXPEC): cdef double rate # need to handle zeros for log-log interpolation - if energy < 1.e-300: - energy = 1.e-300 + if energy <= 0: + return 0 rate = 10 ** self._eb.evaluate(log10(energy)) diff --git a/cherab/openadas/rates/pec.pxd b/cherab/openadas/rates/pec.pxd index 46a54cbe..57bc7560 100644 --- a/cherab/openadas/rates/pec.pxd +++ b/cherab/openadas/rates/pec.pxd @@ -19,7 +19,7 @@ from cherab.core cimport ImpactExcitationPEC as CoreImpactExcitationPEC from cherab.core cimport RecombinationPEC as CoreRecombinationPEC from cherab.core cimport ThermalCXPEC as CoreThermalCXPEC -from cherab.core.math cimport Function2D +from cherab.core.math cimport Function2D, Function3D cdef class ImpactExcitationPEC(CoreImpactExcitationPEC): @@ -48,8 +48,14 @@ cdef class NullRecombinationPEC(CoreRecombinationPEC): pass -# cdef class CXThermalRate(CoreCXThermalRate): -# pass -# -# cdef class ThermalCXRate(CoreThermalCXRate): -# pass +cdef class ThermalCXPEC(CoreThermalCXPEC): + + cdef: + readonly dict raw_data + readonly double wavelength + readonly tuple density_range, temperature_range, donor_temperature_range + Function3D _rate + + +cdef class NullThermalCXPEC(CoreThermalCXPEC): + pass diff --git a/cherab/openadas/rates/pec.pyx b/cherab/openadas/rates/pec.pyx index eafc6c64..b10b0547 100644 --- a/cherab/openadas/rates/pec.pyx +++ b/cherab/openadas/rates/pec.pyx @@ -1,6 +1,6 @@ -# Copyright 2016-2021 Euratom -# Copyright 2016-2021 United Kingdom Atomic Energy Authority -# Copyright 2016-2021 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -20,18 +20,32 @@ import numpy as np from libc.math cimport INFINITY, log10 -from raysect.core.math.function.float cimport Interpolator2DArray +from raysect.core.math.function.float cimport Interpolator2DArray, Interpolator3DArray from cherab.core.utility.conversion import PhotonToJ cdef class ImpactExcitationPEC(CoreImpactExcitationPEC): + """ + Electron impact excitation photon emission coefficient. + + The data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param double wavelength: Resting wavelength of corresponding emission line in nm. + :param dict data: Excitation PEC dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with excitation PEC in photon.m^3.s^-1. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, double wavelength, dict data, extrapolate=False): - """ - :param wavelength: Resting wavelength of corresponding emission line in nm. - :param data: Dictionary containing rate data. - :param extrapolate: Enable extrapolation (default=False). - """ self.wavelength = wavelength self.raw_data = data @@ -56,11 +70,8 @@ cdef class ImpactExcitationPEC(CoreImpactExcitationPEC): cpdef double evaluate(self, double density, double temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(density), log10(temperature)) @@ -68,7 +79,7 @@ cdef class ImpactExcitationPEC(CoreImpactExcitationPEC): cdef class NullImpactExcitationPEC(CoreImpactExcitationPEC): """ - A PEC rate that always returns zero. + A electron impact excitation PEC rate that always returns zero. Needed for use cases where the required atomic data is missing. """ @@ -77,13 +88,27 @@ cdef class NullImpactExcitationPEC(CoreImpactExcitationPEC): cdef class RecombinationPEC(CoreRecombinationPEC): + """ + Recombination photon emission coefficient. + + The data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param double wavelength: Resting wavelength of corresponding emission line in nm. + :param dict data: Rcombination PEC dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with recombination PEC in photon.m^3.s^-1. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, double wavelength, dict data, extrapolate=False): - """ - :param wavelength: Resting wavelength of corresponding emission line in nm. - :param data: Dictionary containing rate data. - :param extrapolate: Enable extrapolation (default=False). - """ self.wavelength = wavelength self.raw_data = data @@ -108,11 +133,8 @@ cdef class RecombinationPEC(CoreRecombinationPEC): cpdef double evaluate(self, double density, double temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if density < 1.e-300: - density = 1.e-300 - - if temperature < 1.e-300: - temperature = 1.e-300 + if density <= 0 or temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(density), log10(temperature)) @@ -120,9 +142,60 @@ cdef class RecombinationPEC(CoreRecombinationPEC): cdef class NullRecombinationPEC(CoreRecombinationPEC): """ - A PEC rate that always returns zero. + A recombination PEC rate that always returns zero. Needed for use cases where the required atomic data is missing. """ cpdef double evaluate(self, double density, double temperature) except? -1e999: return 0.0 + + +cdef class ThermalCXPEC(CoreThermalCXPEC): + + def __init__(self, double wavelength, dict data, extrapolate=False): + """ + :param wavelength: Resting wavelength of corresponding emission line in nm. + :param data: Dictionary containing rate data. + :param extrapolate: Enable nearest-neighbour extrapolation (default=False). + """ + + self.wavelength = wavelength + self.raw_data = data + + # unpack + ne = data['ne'] + te = data['te'] + td = data['td'] + rate = data['rate'] + + # pre-convert data to W m^3 from Photons s^-1 cm^3 prior to interpolation + rate = np.log10(PhotonToJ.to(rate, wavelength)) + + # store limits of data + self.density_range = ne.min(), ne.max() + self.temperature_range = te.min(), te.max() + self.donor_temperature_range = td.min(), td.max() + + # interpolate rate + # using nearest extrapolation to avoid infinite values at 0 for some rates + extrapolation_type = 'nearest' if extrapolate else 'none' + self._rate = Interpolator3DArray(np.log10(ne), np.log10(te), np.log10(td), rate, 'cubic', extrapolation_type, INFINITY, INFINITY, INFINITY) + + cpdef double evaluate(self, double electron_density, double electron_temperature, double donor_temperature) except? -1e999: + + # need to handle zeros, also density and temperature can become negative due to cubic interpolation + if electron_density <= 0 or electron_temperature <= 0 or donor_temperature <= 0: + return 0 + + # calculate rate and convert from log10 space to linear space + return 10 ** self._rate.evaluate(log10(electron_density), log10(electron_temperature), log10(donor_temperature)) + + +cdef class NullThermalCXPEC(CoreThermalCXPEC): + """ + A PEC rate that always returns zero. + Needed for use cases where the required atomic data is missing. + """ + + cpdef double evaluate(self, double electron_density, double electron_temperature, double donor_temperature) except? -1e999: + return 0.0 diff --git a/cherab/openadas/rates/radiated_power.pyx b/cherab/openadas/rates/radiated_power.pyx index 570ced92..85d4c3d5 100644 --- a/cherab/openadas/rates/radiated_power.pyx +++ b/cherab/openadas/rates/radiated_power.pyx @@ -1,7 +1,7 @@ -# Copyright 2016-2021 Euratom -# Copyright 2016-2021 United Kingdom Atomic Energy Authority -# Copyright 2016-2021 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -24,7 +24,26 @@ from raysect.core.math.function.float cimport Interpolator2DArray cdef class LineRadiationPower(CoreLineRadiationPower): - """Base class for radiated powers.""" + """ + Line radiated power coefficient. + + The data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param Element species: Element object defining the ion type. + :param int ionisation: Charge state of the ion. + :param dict data: Line radiated power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with radiated power rate in W.m^3. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, species, ionisation, dict data, extrapolate=False): @@ -49,11 +68,8 @@ cdef class LineRadiationPower(CoreLineRadiationPower): cdef double evaluate(self, double electron_density, double electron_temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if electron_density < 1.e-300: - electron_density = 1.e-300 - - if electron_temperature < 1.e-300: - electron_temperature = 1.e-300 + if electron_density <= 0 or electron_temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(electron_density), log10(electron_temperature)) @@ -70,7 +86,26 @@ cdef class NullLineRadiationPower(CoreLineRadiationPower): cdef class ContinuumPower(CoreContinuumPower): - """Base class for radiated powers.""" + """ + Recombination continuum radiated power coefficient. + + The data is interpolated with cubic spline in log-log space. + Nearest neighbour extrapolation is used if extrapolate is True. + + :param Element species: Element object defining the ion type. + :param int ionisation: Charge state of the ion. + :param dict data: Recombination continuum radiated power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with radiated power rate in W.m^3. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, species, ionisation, dict data, extrapolate=False): @@ -95,11 +130,8 @@ cdef class ContinuumPower(CoreContinuumPower): cdef double evaluate(self, double electron_density, double electron_temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if electron_density < 1.e-300: - electron_density = 1.e-300 - - if electron_temperature < 1.e-300: - electron_temperature = 1.e-300 + if electron_density <= 0 or electron_temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(electron_density), log10(electron_temperature)) @@ -116,7 +148,26 @@ cdef class NullContinuumPower(CoreContinuumPower): cdef class CXRadiationPower(CoreCXRadiationPower): - """Base class for radiated powers.""" + """ + Charge exchange radiated power coefficient. + + The data is interpolated with cubic spline in log-log space. + Linear extrapolation is used if extrapolate is True. + + :param Element species: Element object defining the ion type. + :param int ionisation: Charge state of the ion. + :param dict data: CX radiated power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with radiated power rate in W.m^3. + + :param bint extrapolate: Enable extrapolation (default=False). + + :ivar tuple density_range: Electron density interpolation range. + :ivar tuple temperature_range: Electron temperature interpolation range. + :ivar dict raw_data: Dictionary containing the raw data. + """ def __init__(self, species, ionisation, dict data, extrapolate=False): @@ -140,11 +191,8 @@ cdef class CXRadiationPower(CoreCXRadiationPower): cdef double evaluate(self, double electron_density, double electron_temperature) except? -1e999: # need to handle zeros, also density and temperature can become negative due to cubic interpolation - if electron_density < 1.e-300: - electron_density = 1.e-300 - - if electron_temperature < 1.e-300: - electron_temperature = 1.e-300 + if electron_density <= 0 or electron_temperature <= 0: + return 0 # calculate rate and convert from log10 space to linear space return 10 ** self._rate.evaluate(log10(electron_density), log10(electron_temperature)) diff --git a/cherab/openadas/repository/atomic.py b/cherab/openadas/repository/atomic.py index c24234f3..cb8add6c 100644 --- a/cherab/openadas/repository/atomic.py +++ b/cherab/openadas/repository/atomic.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -34,8 +34,15 @@ def add_ionisation_rate(species, charge, rate, repository_path=None): function instead. The update function avoids repeatedly opening and closing the rate files. - :param repository_path: - :return: + :param species: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param rate: Ionisation rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with ionisation rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ update_ionisation_rates({ @@ -47,11 +54,21 @@ def add_ionisation_rate(species, charge, rate, repository_path=None): def update_ionisation_rates(rates, repository_path=None): """ - Ionisation rate file structure - - /ionisation/.json + Updates the ionisation rate files `/ionisation/.json` + in atomic data repository. File contains multiple rates, indexed by the ion charge state. + + :param rates: Dictionary in the form {: {: }}, where + + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the ionisation rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with ionisation rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -75,8 +92,15 @@ def add_recombination_rate(species, charge, rate, repository_path=None): function instead. The update function avoids repeatedly opening and closing the rate files. - :param repository_path: - :return: + :param species: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param rate: Recombination rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with recombination rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ update_recombination_rates({ @@ -88,11 +112,21 @@ def add_recombination_rate(species, charge, rate, repository_path=None): def update_recombination_rates(rates, repository_path=None): """ - Ionisation rate file structure - - /recombination/.json + Updates the recombination rate files `/recombination/.json` + in the atomic data repository. File contains multiple rates, indexed by the ion charge state. + + :param rates: Dictionary in the form {: {: }}, where + + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the recombination rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with recombination rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -109,7 +143,6 @@ def update_recombination_rates(rates, repository_path=None): def add_thermal_cx_rate(donor_element, donor_charge, receiver_element, rate, repository_path=None): - """ Adds a single thermal charge exchange rate to the repository. @@ -118,11 +151,16 @@ def add_thermal_cx_rate(donor_element, donor_charge, receiver_element, rate, rep the rate files. :param donor_element: Element donating the electron. - :param donor_charge: Charge of the donating atom/ion - :param receiver_element: Element receiving the electron - :param rate: rates - :param repository_path: - :return: + :param donor_charge: Charge of the donating atom/ion. + :param receiver_element: Element receiving the electron. + :param receiver_charge: Charge of the receiving atom/ion. + :param rate: Thermal CX rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with thermal CX rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ rates2update = RecursiveDict() @@ -133,11 +171,25 @@ def add_thermal_cx_rate(donor_element, donor_charge, receiver_element, rate, rep def update_thermal_cx_rates(rates, repository_path=None): """ - Thermal charge exchange rate file structure - - /thermal_cx///.json + Updates the thermal charge exchange rate files + `/thermal_cx///.json` + in the atomic data repository. File contains multiple rates, indexed by the ion charge state. + + :param rates: Dictionary in the form: + + | { : { : { : { : } } } }, where + | is the element donating the electron. + | is the charge of the donating atom/ion. + | is the element receiving the electron. + | is the charge of the receiving atom/ion. + | is the thermal CX rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with thermal CX rate in m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -203,6 +255,21 @@ def _update_and_write_adf11(species, rate_data, path): def get_ionisation_rate(element, charge, repository_path=None): + """ + Reads the ionisation rate for the given species and charge + from the atomic data repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param repository_path: Path to the atomic data repository. + + :return rate: Ionisation rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with ionisation rate in m^3.s^-1. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -224,6 +291,21 @@ def get_ionisation_rate(element, charge, repository_path=None): def get_recombination_rate(element, charge, repository_path=None): + """ + Reads the recombination rate for the given species and charge + from the atomic data repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param repository_path: Path to the atomic data repository. + + :return rate: Recombination rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with recombination rate in m^3.s^-1. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -245,6 +327,23 @@ def get_recombination_rate(element, charge, repository_path=None): def get_thermal_cx_rate(donor_element, donor_charge, receiver_element, receiver_charge, repository_path=None): + """ + Reads the thermal charge exchange rate for the given species and charge + from the atomic data repository. + + :param donor_element: Element donating the electron. + :param donor_charge: Charge of the donating atom/ion. + :param receiver_element: Element receiving the electron. + :param receiver_charge: Charge of the receiving atom/ion. + :param repository_path: Path to the atomic data repository. + + :return rate: Thermal CX rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with thermal CX rate in m^3.s^-1. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH diff --git a/cherab/openadas/repository/beam/cx.py b/cherab/openadas/repository/beam/cx.py index 65bc3ceb..ef79c102 100644 --- a/cherab/openadas/repository/beam/cx.py +++ b/cherab/openadas/repository/beam/cx.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -29,19 +29,33 @@ def add_beam_cx_rate(donor_ion, donor_metastable, receiver_ion, receiver_charge, transition, rate, repository_path=None): """ - Adds a single beam CX rate to the repository. + Adds a single beam CX PEC to the repository. If adding multiple rate, consider using the update_beam_cx_rates() function instead. The update function avoid repeatedly opening and closing the rate files. - :param donor_ion: - :param donor_metastable: - :param receiver_ion: - :param receiver_charge: - :param rate: - :param repository_path: - :return: + :param donor_ion: Beam neutral atom (Element/Isotope) donating the electron. + :param donor_metastable: Metastable/excited level of beam neutral atom. + :param receiver_ion: Element/Isotope receiving the electron. + :param receiver_charge: Charge of the receiving atom/ion. + :param transition: Tuple containing (initial level, final level). + :param rate: Beam CX PEC dictionary containing the following entries: + + | 'eb': array-like of size (N) with beam energy in eV/amu, + | 'ti': array-like of size (M) with receiver ion temperature in eV, + | 'ni': array-like of size (K) with plasma ion density in m^-3, + | 'z': array-like of size (L) with plasma Z-effective, + | 'b': array-like of size (J) with magnetic field strength in Tesla, + | 'qeb': array-like of size (N) with CX PEC energy component in photon.m^3.s-1, + | 'qti': array-like of size (M) with CX PEC temperature component in photon.m^3.s-1, + | 'qni': array-like of size (K) with CX PEC density component in photon.m^3.s-1, + | 'qz': array-like of size (L) with CX PEC Zeff component in photon.m^3.s-1, + | 'qb': array-like of size (J) with CX PEC B-field component in photon.m^3.s-1, + | 'qref': reference CX PEC in photon.m^3.s-1. + | The total beam CX PEC: q = qeb * qti * qni * qz * qb / qref^4. + + :param repository_path: Path to the atomic data repository. """ update_beam_cx_rates({ @@ -58,10 +72,37 @@ def add_beam_cx_rate(donor_ion, donor_metastable, receiver_ion, receiver_charge, def update_beam_cx_rates(rates, repository_path=None): - # organisation in repository: - # beam/cx/donor_ion/receiver_ion/receiver_charge.json - # inside json file: - # transition: [list of donor_metastables with rates] + """ + Updates the beam CX PEC files + beam/cx///.json + in the atomic data repository. + + File contains multiple metastable-resolved rates, indexed by transition. + + :param rates: Dictionary in the form: + + | { : { : { : { : {: } } } } }, where + | is the beam neutral atom (Element/Isotope) donating the electron. + | is the metastable/excited level of beam neutral atom. + | is the Element/Isotope receiving the electron. + | is the charge of the receiving atom/ion. + | is the tuple containing (initial level, final level). + | is the beam CX PEC dictionary containing the following entries: + | 'eb': array-like of size (N) with beam energy in eV/amu, + | 'ti': array-like of size (M) with receiver ion temperature in eV, + | 'ni': array-like of size (K) with plasma ion density in m^-3, + | 'z': array-like of size (L) with plasma Z-effective, + | 'b': array-like of size (J) with magnetic field strength in Tesla, + | 'qeb': array-like of size (N) with CX PEC energy component in photon.m^3.s-1, + | 'qti': array-like of size (M) with CX PEC temperature component in photon.m^3.s-1, + | 'qni': array-like of size (K) with CX PEC density component in photon.m^3.s-1, + | 'qz': array-like of size (L) with CX PEC Zeff component in photon.m^3.s-1, + | 'qb': array-like of size (J) with CX PEC B-field component in photon.m^3.s-1, + | 'qref': reference CX PEC in photon.m^3.s-1. + | The total beam CX PEC: q = qeb * qti * qni * qz * qb / qref^4. + + :param repository_path: Path to the atomic data repository. + """ def sanitise_and_validate(data, x_key, x_name, y_key, y_name): """ @@ -167,6 +208,32 @@ def sanitise_and_validate(data, x_key, x_name, y_key, y_name): def get_beam_cx_rates(donor_ion, receiver_ion, receiver_charge, transition, repository_path=None): + """ + Reads a single beam CX PEC from the repository. + + :param donor_ion: Beam neutral atom (Element/Isotope) donating the electron. + :param donor_metastable: Metastable/excited level of beam neutral atom. + :param receiver_ion: Element/Isotope receiving the electron. + :param receiver_charge: Charge of the receiving atom/ion. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. + + :return rate: Beam CX PEC dictionary containing the following entries: + + | 'eb': 1D array of size (N) with beam energy in eV/amu, + | 'ti': 1D array of size (M) with receiver ion temperature in eV, + | 'ni': 1D array of size (K) with plasma ion density in m^-3, + | 'z': 1D array of size (L) with plasma Z-effective, + | 'b': 1D array of size (J) with magnetic field strength in Tesla, + | 'qeb': 1D array of size (N) with CX PEC energy component in photon.m^3.s-1, + | 'qti': 1D array of size (M) with CX PEC temperature component in photon.m^3.s-1, + | 'qni': 1D array of size (K) with CX PEC density component in photon.m^3.s-1, + | 'qz': 1D array of size (L) with CX PEC Zeff component in photon.m^3.s-1, + | 'qb': 1D array of size (J) with CX PEC B-field component in photon.m^3.s-1, + | 'qref': reference CX PEC in photon.m^3.s-1. + | The total beam CX PEC: q = qeb * qti * qni * qz * qb / qref^4. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH path = os.path.join(repository_path, 'beam/cx/{}/{}/{}.json'.format(donor_ion.symbol.lower(), receiver_ion.symbol.lower(), receiver_charge)) diff --git a/cherab/openadas/repository/beam/emission.py b/cherab/openadas/repository/beam/emission.py index b6cd14f2..c7c6e2e1 100644 --- a/cherab/openadas/repository/beam/emission.py +++ b/cherab/openadas/repository/beam/emission.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -36,8 +36,24 @@ def add_beam_emission_rate(beam_species, target_ion, target_charge, transition, function instead. The update function avoid repeatedly opening and closing the rate files. - :param repository_path: - :return: + :param beam_species: Beam neutral species (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param transition: Tuple containing (initial level, final level). + :param rate: Beam emission rate dictionary containing the following entries: + + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n' array-like of size (M) with target electron density in m^-3, + | 't' array-like of size (K) with target electron temperature in eV, + | 'sen' array-like of size (N, M) with beam emission rate energy component in photon.m^3.s^-1. + | 'st' array-like of size (K) with beam emission rate temperature component in photon.m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam emission rate in photon.m^3.s^-1. + | The total beam emission rate: s = sen * st / sref. + + :param repository_path: Path to the atomic data repository. """ update_beam_emission_rates({ @@ -53,11 +69,32 @@ def add_beam_emission_rate(beam_species, target_ion, target_charge, transition, def update_beam_emission_rates(rates, repository_path=None): """ - Beam emission rate file structure - + Updates the beam emission rate files: /beam/emission///.json + in the atomic repository. File contains multiple rates, indexed by transition. + + :param rates: Dictionary in the form: + + | { : { : { : {: } } } }, where + | is the beam neutral species (Element/Isotope) + | is the target species (Element/Isotope). + | is the charge of the target species. + | is the tuple containing (initial level, final level). + | Beam emission rate dictionary containing the following entries: + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n' array-like of size (M) with target electron density in m^-3, + | 't' array-like of size (K) with target electron temperature in eV, + | 'sen' array-like of size (N, M) with beam emission rate energy component in photon.m^3.s^-1. + | 'st' array-like of size (K) with beam emission rate temperature component in photon.m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam emission rate in photon.m^3.s^-1. + | The total beam emission rate: s = sen * st / sref. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -136,6 +173,29 @@ def update_beam_emission_rates(rates, repository_path=None): def get_beam_emission_rate(beam_species, target_ion, target_charge, transition, repository_path=None): + """ + Reads a single beam emission rate from the repository. + + :param beam_species: Beam neutral species (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. + + :return rate: Beam emission rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n' 1D array of size (M) with target electron density in m^-3, + | 't' 1D array of size (K) with target electron temperature in eV, + | 'sen' 2D array of size (N, M) with beam emission rate energy component in photon.m^3.s^-1. + | 'st' 1D array of size (K) with beam emission rate temperature component in photon.m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam emission rate in photon.m^3.s^-1. + | The total beam emission rate: s = sen * st / sref. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH path = os.path.join(repository_path, 'beam/emission/{}/{}/{}.json'.format(beam_species.symbol.lower(), target_ion.symbol.lower(), target_charge)) diff --git a/cherab/openadas/repository/beam/population.py b/cherab/openadas/repository/beam/population.py index 54c57df1..9ecfa6eb 100644 --- a/cherab/openadas/repository/beam/population.py +++ b/cherab/openadas/repository/beam/population.py @@ -31,12 +31,24 @@ def add_beam_population_rate(beam_species, beam_metastable, target_ion, target_c """ Adds a single beam population rate to the repository. - :param beam_species: - :param beam_metastable: - :param target_ion: - :param target_charge: - :param rate: - :return: + :param beam_species: Beam neutral species (Element/Isotope). + :param beam_metastable: Metastable level of beam neutral atom. + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param rate: Beam population rate dictionary containing the following entries: + + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n': array-like of size (M) with target electron density in m^-3, + | 't': array-like of size (K) with target electron temperature in eV, + | 'sen': array-like of size (N, M) with dimensionless beam population rate energy component. + | 'st': array-like of size (K) with dimensionless beam population rate temperature component. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference dimensionless beam population rate. + | The total beam population rate: s = sen * st / sref. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -102,11 +114,32 @@ def add_beam_population_rate(beam_species, beam_metastable, target_ion, target_c def update_beam_population_rates(rates, repository_path=None): """ - Beam population rate file structure - + Updates the beam population rate files /beam/population////.json + in the atomic data repository. Each json file contains a single rate, so it can simply be replaced. + + :param rates: Dictionary in the form: + + | { : { : { : {: } } } }, where + | is the beam neutral species (Element/Isotope) + | is the metastable level of beam neutral atom. + | is the target species (Element/Isotope). + | is the charge of the target species. + | is the beam population rate dictionary containing the following fields: + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n': array-like of size (M) with target electron density in m^-3, + | 't': array-like of size (K) with target electron temperature in eV, + | 'sen': array-like of size (N, M) with dimensionless beam population rate energy component. + | 'st': array-like of size (K) with dimensionless beam population rate temperature component. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference dimensionless beam population rate. + | The total beam population rate: s = sen * st / sref. + + :param repository_path: Path to the atomic data repository. """ for beam_species, beam_metastables in rates.items(): @@ -117,6 +150,29 @@ def update_beam_population_rates(rates, repository_path=None): def get_beam_population_rate(beam_species, beam_metastable, target_ion, target_charge, repository_path=None): + """ + Reads a single beam population rate from the repository. + + :param beam_species: Beam neutral species (Element/Isotope). + :param beam_metastable: Metastable level of beam neutral atom. + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param repository_path: Path to the atomic data repository. + + :return rate: Beam population rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n': 1D array of size (M) with target electron density in m^-3, + | 't': 1D array of size (K) with target electron temperature in eV, + | 'sen': 2D array of size (N, M) with dimensionless beam population rate energy component. + | 'st': 1D array of size (K) with dimensionless beam population rate temperature component. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference dimensionless beam population rate. + | The total beam population rate: s = sen * st / sref. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH path = os.path.join(repository_path, 'beam/population/{}/{}/{}/{}.json'.format(beam_species.symbol.lower(), beam_metastable, target_ion.symbol.lower(), target_charge)) diff --git a/cherab/openadas/repository/beam/stopping.py b/cherab/openadas/repository/beam/stopping.py index cf8d1492..46d8caae 100644 --- a/cherab/openadas/repository/beam/stopping.py +++ b/cherab/openadas/repository/beam/stopping.py @@ -31,11 +31,23 @@ def add_beam_stopping_rate(beam_species, target_ion, target_charge, rate, reposi """ Adds a single beam stopping/excitation rate to the repository. - :param beam_species: - :param target_ion: - :param target_charge: - :param rate: - :return: + :param beam_species: Beam neutral atom (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param rate: Beam stopping rate dictionary containing the following entries: + + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n': array-like of size (M) with target electron density in m^-3, + | 't': array-like of size (K) with target electron temperature in eV, + | 'sen': array-like of size (N, M) with beam stopping rate energy component in m^3.s^-1. + | 'st': array-like of size (K) with beam stopping rate temperature component in m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam stopping rate in m^3.s^-1. + | The total beam stopping rate: s = sen * st / sref. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -98,11 +110,30 @@ def add_beam_stopping_rate(beam_species, target_ion, target_charge, rate, reposi def update_beam_stopping_rates(rates, repository_path=None): """ - Beam stopping rate file structure - - /beam/stopping///.json + Updates the beam stopping rate files + /beam/stopping////.json + in the atomic data repository. Each json file contains a single rate, so it can simply be replaced. + + :param rates: Dictionary in the form: + + | { : { : { : {: } } } }, where + | is the beam neutral species (Element/Isotope). + | is the target species (Element/Isotope). + | is the charge of the target species. + | is the beam stopping rate dictionary containing the following entries: + | 'e': array-like of size (N) with interaction energy in eV/amu, + | 'n': array-like of size (M) with target electron density in m^-3, + | 't': array-like of size (K) with target electron temperature in eV, + | 'sen': array-like of size (N, M) with beam stopping rate energy component in m^3.s^-1. + | 'st': array-like of size (K) with beam stopping rate temperature component in m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam stopping rate in m^3.s^-1. + | The total beam stopping rate: s = sen * st / sref. + """ for beam_species, target_ions in rates.items(): @@ -112,6 +143,28 @@ def update_beam_stopping_rates(rates, repository_path=None): def get_beam_stopping_rate(beam_species, target_ion, target_charge, repository_path=None): + """ + Reads a single beam stopping/excitation rate from the repository. + + :param beam_species: Beam neutral atom (Element/Isotope). + :param target_ion: Target species (Element/Isotope). + :param target_charge: Charge of the target species. + :param repository_path: Path to the atomic data repository. + + :return rate: Beam stopping rate dictionary containing the following entries: + + | 'e': 1D array of size (N) with interaction energy in eV/amu, + | 'n': 1D array of size (M) with target electron density in m^-3, + | 't': 1D array of size (K) with target electron temperature in eV, + | 'sen': 2D array of size (N, M) with beam stopping rate energy component in m^3.s^-1. + | 'st': 1D array of size (K) with beam stopping rate temperature component in m^3.s^-1. + | 'eref': reference interaction energy in eV/amu, + | 'nref': reference target electron density in m^-3, + | 'tref': reference target electron temperature in eV, + | 'sref': reference beam stopping rate in m^3.s^-1. + | The total beam stopping rate: s = sen * st / sref. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH path = os.path.join(repository_path, 'beam/stopping/{}/{}/{}.json'.format(beam_species.symbol.lower(), target_ion.symbol.lower(), target_charge)) diff --git a/cherab/openadas/repository/create.py b/cherab/openadas/repository/create.py index 6b9a4f33..b42e5651 100644 --- a/cherab/openadas/repository/create.py +++ b/cherab/openadas/repository/create.py @@ -147,6 +147,7 @@ def populate(download=True, repository_path=None, adas_path=None): (carbon, 0, 'adf15/pec96#c/pec96#c_vsu#c0.dat'), (carbon, 1, 'adf15/pec96#c/pec96#c_vsu#c1.dat'), (carbon, 2, 'adf15/pec96#c/pec96#c_vsu#c2.dat'), + (carbon, 5, 'adf15/pec96#c/pec96#c_pju#c5.dat'), # (neon, 0, 'adf15/pec96#ne/pec96#ne_pju#ne0.dat'), #TODO: OPENADAS DATA CORRUPT # (neon, 1, 'adf15/pec96#ne/pec96#ne_pju#ne1.dat'), #TODO: OPENADAS DATA CORRUPT (nitrogen, 0, 'adf15/pec96#n/pec96#n_vsu#n0.dat'), diff --git a/cherab/openadas/repository/pec.py b/cherab/openadas/repository/pec.py index 8eb867fc..ef933b23 100644 --- a/cherab/openadas/repository/pec.py +++ b/cherab/openadas/repository/pec.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -36,12 +36,16 @@ def add_pec_excitation_rate(element, charge, transition, rate, repository_path=N instead. The update function avoid repeatedly opening and closing the rate files. - :param element: - :param charge: - :param transition: - :param rate: - :param repository_path: - :return: + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param rate: Excitation PEC dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with excitation PEC in photon.m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ update_pec_rates({ @@ -63,12 +67,16 @@ def add_pec_recombination_rate(element, charge, transition, rate, repository_pat instead. The update function avoid repeatedly opening and closing the rate files. - :param element: - :param charge: - :param transition: - :param rate: - :param repository_path: - :return: + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param rate: Recombination PEC dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with recombination PEC in photon.m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ update_pec_rates({ @@ -82,44 +90,59 @@ def add_pec_recombination_rate(element, charge, transition, rate, repository_pat }, repository_path) -def add_pec_thermalcx_rate(element, charge, transition, rate, repository_path=None): +def add_pec_thermal_cx_rate(donor_element, donor_charge, receiver_element, receiver_charge, transition, rate, repository_path=None): """ - Adds a single PEC thermalcx rate to the repository. + Adds a single PEC thermal charge exchange rate to the repository. - If adding multiple rate, consider using the update_pec_rates() function + If adding multiple rate, consider using the update_pec_thermal_rates() function instead. The update function avoid repeatedly opening and closing the rate files. - :param element: - :param charge: - :param transition: - :param rate: - :param repository_path: - :return: + :param donor_element: Electron donor plasma species (Element/Isotope). + :param donor_charge: Electron donor charge. + :param receiver_element: Electron receiver plasma species (Element/Isotope). + :param receiver_charge: Electron receiver charge. + :param transition: Tuple containing (initial level, final level). + :param rate: Thermal CX PEC dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'td': array-like of size (K) with donor temperature in eV, + | 'rate': array-like of size (N, M, K) with thermal CX PEC in photon.m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ + rates2update = RecursiveDict() + rates2update[donor_element][donor_charge][receiver_element][receiver_charge][transition] = rate - update_pec_rates({ - 'thermalcx': { - element: { - charge: { - transition: rate - } - } - } - }, repository_path) + update_pec_thermal_cx_rates(rates2update.freeze(), repository_path) def update_pec_rates(rates, repository_path=None): """ - PEC rate file structure + Updates excitation and recombination PEC files /pec///.json. + in the atomic data repository. - /pec///.json + File contains multiple PECs, indexed by the transition. + + :param rates: Dictionary in the form: + + | { : { : { : { : } } } }, where + | is the one of the following PEC types: 'excitation', 'recombination'. + | is the plasma species (Element/Isotope). + | is the charge of the plasma species. + | is the tuple containing (initial level, final level). + | is the PEC dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with PEC in photon.m^3.s^-1. + + :param repository_path: Path to the atomic data repository. """ valid_classes = [ 'excitation', - 'recombination', - 'thermalcx' + 'recombination' ] repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -183,16 +206,139 @@ def update_pec_rates(rates, repository_path=None): json.dump(content, f, indent=2, sort_keys=True) +def update_pec_thermal_cx_rates(rates, repository_path=None): + """ + Updates thermal CX PEC files /pec/thermal_cx////.json + in the atomic data repository. + + File contains multiple PECs, indexed by the transition. + + :param rates: Dictionary in the form: + + | { : { : { : { : { : } } } } }, where + | is the electron donor species (Element/Isotope). + | is the electron donor charge. + | is the electron receiver species (Element/Isotope). + | is the electron receiver charge. + | is the tuple containing (initial level, final level). + | is the PEC dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'td': array-like of size (K) with donor temperature in eV, + | 'rate': array-like of size (N, M, K) with PEC in photon.m^3.s^-1. + + :param repository_path: Path to the atomic data repository. + """ + repository_path = repository_path or DEFAULT_REPOSITORY_PATH + + for donor_element, donor_charge_states in rates.items(): + for donor_charge, receiver_elements in donor_charge_states.items(): + for receiver_element, receiver_charge_states in receiver_elements.items(): + for receiver_charge, transitions in receiver_charge_states.items(): + + # sanitise and validate + if not isinstance(donor_element, Element): + raise TypeError('The donor element must be an Element object.') + + if not valid_charge(donor_element, donor_charge + 1): + raise ValueError('Donor charge state is equal to or larger than the number of protons in the element.') + + if not isinstance(receiver_element, Element): + raise TypeError('The receiver element must be an Element object.') + + if not valid_charge(receiver_element, receiver_charge): + raise ValueError('Receiver charge state is larger than the number of protons in the element.') + + rate_path = 'pec/thermal_cx/{0}/{1}/{2}/{3}.json'.format(donor_element.symbol.lower(), donor_charge, + receiver_element.symbol.lower(), receiver_charge) + path = os.path.join(repository_path, rate_path) + + # read in any existing rates + try: + with open(path, 'r') as f: + content = RecursiveDict.from_dict(json.load(f)) + except FileNotFoundError: + content = RecursiveDict() + + # add/replace data for a transition + for transition, data in transitions.items(): + key = encode_transition(transition) + + # sanitise/validate data + ne = np.array(data['ne'], np.float64) + te = np.array(data['te'], np.float64) + td = np.array(data['td'], np.float64) + rate = np.array(data['rate'], np.float64) + + if ne.ndim != 1: + raise ValueError('Electron density array must be a 1D array.') + + if te.ndim != 1: + raise ValueError('Electron temperature array must be a 1D array.') + + if td.ndim != 1: + raise ValueError('Donor temperature array must be a 1D array.') + + if (ne.shape[0], te.shape[0], td.shape[0]) != rate.shape: + raise ValueError('Electron density, electron temperature, donor temperature' + ' and rate data arrays have inconsistent sizes.') + + content[key] = { + 'ne': ne.tolist(), + 'te': te.tolist(), + 'td': td.tolist(), + 'rate': rate.tolist() + } + + # create directory structure if missing + directory = os.path.dirname(path) + os.makedirs(directory, exist_ok=True) + + # write new data + with open(path, 'w') as f: + json.dump(content.freeze(), f, indent=2, sort_keys=True) + + def get_pec_excitation_rate(element, charge, transition, repository_path=None): + """ + Reads the excitation PEC from the repository for the given + element, charge and transition. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. + + :return rate: Excitation PEC dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with excitation PEC in photon.m^3.s^-1. + + """ + return _get_pec_rate('excitation', element, charge, transition, repository_path) def get_pec_recombination_rate(element, charge, transition, repository_path=None): - return _get_pec_rate('recombination', element, charge, transition, repository_path) + """ + Reads the recombination PEC from the repository for the given + element, charge and transition. + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. -def get_pec_thermalcx_rate(element, charge, transition, repository_path=None): - return _get_pec_rate('thermalcx', element, charge, transition, repository_path) + :return rate: Recombination PEC dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with recombination PEC in photon.m^3.s^-1. + + """ + + return _get_pec_rate('recombination', element, charge, transition, repository_path) def _get_pec_rate(cls, element, charge, transition, repository_path=None): @@ -214,3 +360,46 @@ def _get_pec_rate(cls, element, charge, transition, repository_path=None): return d + +def get_pec_thermal_cx_rate(donor_element, donor_charge, receiver_element, receiver_charge, transition, repository_path=None): + """ + Reads the thermal charge exchange PEC from the repository for the given + donor element, donor charge, receiver element, receiver charge and transition. + + :param donor_element: Electron donor plasma species (Element/Isotope). + :param donor_charge: Electron donor charge. + :param receiver_element: Electron receiver plasma species (Element/Isotope). + :param receiver_charge: Electron receiver charge. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. + + :return rate: Thermal CX PEC dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'td': 1D array of size (K) with donor temperature in eV, + | 'rate': 2D array of size (N, M, K) with thermal CX PEC in photon.m^3.s^-1. + + """ + + repository_path = repository_path or DEFAULT_REPOSITORY_PATH + + rate_path = 'pec/thermal_cx/{0}/{1}/{2}/{3}.json'.format(donor_element.symbol.lower(), donor_charge, + receiver_element.symbol.lower(), receiver_charge) + path = os.path.join(repository_path, rate_path) + try: + with open(path, 'r') as f: + content = json.load(f) + d = content[encode_transition(transition)] + except (FileNotFoundError, KeyError): + raise RuntimeError('Requested thermal charge-exchange PEC (donor={}, donor charge={}, receiver={}, receiver charge={})' + ' is not available.' + ''.format(donor_element.symbol, donor_charge, receiver_element.symbol, receiver_charge)) + + # convert to numpy arrays + d['ne'] = np.array(d['ne'], np.float64) + d['te'] = np.array(d['te'], np.float64) + d['td'] = np.array(d['td'], np.float64) + d['rate'] = np.array(d['rate'], np.float64) + + return d diff --git a/cherab/openadas/repository/radiated_power.py b/cherab/openadas/repository/radiated_power.py index a6a4b390..7b908cfc 100644 --- a/cherab/openadas/repository/radiated_power.py +++ b/cherab/openadas/repository/radiated_power.py @@ -1,7 +1,7 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -29,13 +29,21 @@ def add_line_power_rate(species, charge, rate, repository_path=None): """ - Adds a single LineRadiationPower rate to the repository. + Adds a single line radiated power rate to the repository. If adding multiple rates, consider using the update_line_power_rates() function instead. The update function avoids repeatedly opening and closing the rate files. - :param repository_path: + :param species: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param rate: Line radiated power rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with line radiated power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ update_line_power_rates({ @@ -47,13 +55,22 @@ def add_line_power_rate(species, charge, rate, repository_path=None): def update_line_power_rates(rates, repository_path=None): """ - Update the repository of LineRadiationPower rates. - - LineRadiationPower rate file structure - + Update the files for the line radiated power rates: /radiated_power/line/.json + in the atomic data repository. File contains multiple rates, indexed by the ion's charge state. + + :param rates: Dictionary in the form {: {: }}, where + + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the line radiated rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with line radiated power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -71,13 +88,21 @@ def update_line_power_rates(rates, repository_path=None): def add_continuum_power_rate(species, charge, rate, repository_path=None): """ - Adds a single ContinuumPower rate to the repository. + Adds a single continuum power rate to the repository. If adding multiple rates, consider using the update_continuum_power_rates() function instead. The update function avoids repeatedly opening and closing the rate files. - :param repository_path: + :param species: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param rate: Continuum power rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with continuum power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ update_line_power_rates({ @@ -89,13 +114,22 @@ def add_continuum_power_rate(species, charge, rate, repository_path=None): def update_continuum_power_rates(rates, repository_path=None): """ - Update the repository of ContinuumPower rates. - - ContinuumPower rate file structure - + Update the files for the continuum power rates: /radiated_power/continuum/.json + in the atomic data repository. File contains multiple rates, indexed by ion's charge state. + + :param rates: Dictionary in the form {: {: }}, where + + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the continuum power rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with continuum power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -113,13 +147,22 @@ def update_continuum_power_rates(rates, repository_path=None): def add_cx_power_rate(species, charge, rate, repository_path=None): """ - Adds a single CXRadiationPower rate to the repository. + Adds a single CX radiation power rate to the repository + (charge exchage with neutral hydrogen). If adding multiple rates, consider using the update_cx_power_rates() function instead. The update function avoids repeatedly opening and closing the rate files. - :param repository_path: + :param species: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param rate: CX power rate dictionary containing the following entries: + + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with CX power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ update_line_power_rates({ @@ -131,13 +174,23 @@ def add_cx_power_rate(species, charge, rate, repository_path=None): def update_cx_power_rates(rates, repository_path=None): """ - Update the repository of CXRadiationPower rates. - - CXRadiationPower rate file structure - + Update the files for the CX radiation power rates + (charge exchage with neutral hydrogen): /radiated_power/cx/.json + in the atomic data repository. File contains multiple rates, indexed by ion's charge state. + + :param rates: Dictionary in the form {: {: }}, where + + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the thermal CX power rate dictionary containing the following entries: + | 'ne': array-like of size (N) with electron density in m^-3, + | 'te': array-like of size (M) with electron temperature in eV, + | 'rate': array-like of size (N, M) with thermal CX power rate in W.m^3. + + :param repository_path: Path to the atomic data repository. """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -199,6 +252,21 @@ def _update_and_write_adf11(species, rate_data, path): def get_line_radiated_power_rate(element, charge, repository_path=None): + """ + Reads the line radiated power rate for the given species and charge + from the atomic data repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param repository_path: Path to the atomic data repository. + + :return rate: Line radiated power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with line radiated power rate in W.m^3. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -220,6 +288,21 @@ def get_line_radiated_power_rate(element, charge, repository_path=None): def get_continuum_radiated_power_rate(element, charge, repository_path=None): + """ + Reads the continuum power rate for the given species and charge + from the atomic data repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param repository_path: Path to the atomic data repository. + + :return rate: Continuum power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with continuum power rate in W.m^3. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -241,6 +324,21 @@ def get_continuum_radiated_power_rate(element, charge, repository_path=None): def get_cx_radiated_power_rate(element, charge, repository_path=None): + """ + Reads the CX radiation power rate for the given species and charge + from the atomic data repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param repository_path: Path to the atomic data repository. + + :return rate: CX radiation power rate dictionary containing the following entries: + + | 'ne': 1D array of size (N) with electron density in m^-3, + | 'te': 1D array of size (M) with electron temperature in eV, + | 'rate': 2D array of size (N, M) with CX radiation power rate in W.m^3. + + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH diff --git a/cherab/openadas/repository/wavelength.py b/cherab/openadas/repository/wavelength.py index 5981e9e6..83758f51 100644 --- a/cherab/openadas/repository/wavelength.py +++ b/cherab/openadas/repository/wavelength.py @@ -1,6 +1,6 @@ -# Copyright 2016-2018 Euratom -# Copyright 2016-2018 United Kingdom Atomic Energy Authority -# Copyright 2016-2018 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas +# Copyright 2016-2024 Euratom +# Copyright 2016-2024 United Kingdom Atomic Energy Authority +# Copyright 2016-2024 Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas # # Licensed under the EUPL, Version 1.1 or – as soon they will be approved by the # European Commission - subsequent versions of the EUPL (the "Licence"); @@ -35,11 +35,11 @@ def add_wavelength(element, charge, transition, wavelength, repository_path=None function instead. The update function avoid repeatedly opening and closing the rate files. - :param element: - :param charge: - :param transition: - :param wavelength: - :param repository_path: + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param wavelength: Transition's wavelength in nm. + :param repository_path: Path to the atomic data repository. """ update_wavelengths({ @@ -52,6 +52,22 @@ def add_wavelength(element, charge, transition, wavelength, repository_path=None def update_wavelengths(wavelengths, repository_path=None): + """ + Updates the wavelength files `/wavelength//.json` + in atomic data repository. + + File contains multiple rates, indexed by the transitions. + + :param wavelengths: Dictionary in the form: + + | { : { : { : } } }, where + | is the plasma species (Element/Isotope), + | is the charge of the plasma species, + | is the tuple containing (initial level, final level), + | is the transition's wavelength in nm. + + :param repository_path: Path to the atomic data repository. + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH @@ -90,6 +106,16 @@ def update_wavelengths(wavelengths, repository_path=None): def get_wavelength(element, charge, transition, repository_path=None): + """ + Reads the wavelength for the given species, charge and transition from the repository. + + :param element: Plasma species (Element/Isotope). + :param charge: Charge of the plasma species. + :param transition: Tuple containing (initial level, final level). + :param repository_path: Path to the atomic data repository. + + :return wavelength: Wavelength in nm. + """ repository_path = repository_path or DEFAULT_REPOSITORY_PATH path = os.path.join(repository_path, 'wavelength/{}/{}.json'.format(element.symbol.lower(), charge)) diff --git a/demos/emission_models/charge_exchange.py b/demos/emission_models/charge_exchange.py index 11d95386..02d61235 100644 --- a/demos/emission_models/charge_exchange.py +++ b/demos/emission_models/charge_exchange.py @@ -122,7 +122,7 @@ integration_step = 0.0025 beam_transform = translate(-0.5, 0.0, 0) * rotate_basis(Vector3D(1, 0, 0), Vector3D(0, 0, 1)) -beam_energy = 50000 # keV +beam_energy = 50000 # eV beam_full = Beam(parent=world, transform=beam_transform) beam_full.plasma = plasma diff --git a/demos/emission_models/thermal_charge_exchange.py b/demos/emission_models/thermal_charge_exchange.py new file mode 100644 index 00000000..42c7b88a --- /dev/null +++ b/demos/emission_models/thermal_charge_exchange.py @@ -0,0 +1,126 @@ + +import matplotlib.pyplot as plt +from scipy.constants import electron_mass, atomic_mass + +from raysect.core import translate, rotate_basis, Point3D, Vector3D +from raysect.primitive import Box +from raysect.optical import World, Ray + +from cherab.core import Plasma, Beam, Maxwellian, Species +from cherab.core.math import ScalarToVectorFunction3D +from cherab.core.atomic import hydrogen, deuterium, carbon, Line +from cherab.core.model import SingleRayAttenuator, BeamCXLine, ThermalCXLine +from cherab.tools.plasmas.slab import NeutralFunction, IonFunction +from cherab.openadas import OpenADAS +from cherab.openadas.install import install_adf15 + + +############### +# Make Plasma # + +width = 1.0 +length = 1.0 +height = 3.0 +peak_density = 5e19 +edge_density = 1e18 +pedestal_top = 1 +neutral_temperature = 0.5 +peak_temperature = 2500 +edge_temperature = 25 +impurities = [(carbon, 6, 0.005)] + +world = World() +adas = OpenADAS(permit_extrapolation=True, missing_rates_return_null=True) + +plasma = Plasma(parent=world) +plasma.atomic_data = adas + +# plasma slab along x direction +plasma.geometry = Box(Point3D(0, -width / 2, -height / 2), Point3D(length, width / 2, height / 2)) + +species = [] + +# make a non-zero velocity profile for the plasma +vy_profile = IonFunction(1E5, 0, pedestal_top=pedestal_top) +velocity_profile = ScalarToVectorFunction3D(0, vy_profile, 0) + +# define neutral species distribution +h0_density = NeutralFunction(peak_density, 0.1, pedestal_top=pedestal_top) +h0_temperature = neutral_temperature +h0_distribution = Maxwellian(h0_density, h0_temperature, velocity_profile, + hydrogen.atomic_weight * atomic_mass) +species.append(Species(hydrogen, 0, h0_distribution)) + +# define hydrogen ion species distribution +h1_density = IonFunction(peak_density, edge_density, pedestal_top=pedestal_top) +h1_temperature = IonFunction(peak_temperature, edge_temperature, pedestal_top=pedestal_top) +h1_distribution = Maxwellian(h1_density, h1_temperature, velocity_profile, + hydrogen.atomic_weight * atomic_mass) +species.append(Species(hydrogen, 1, h1_distribution)) + +# add impurities +if impurities: + for impurity, ionisation, concentration in impurities: + imp_density = IonFunction(peak_density * concentration, edge_density * concentration, pedestal_top=pedestal_top) + imp_temperature = IonFunction(peak_temperature, edge_temperature, pedestal_top=pedestal_top) + imp_distribution = Maxwellian(imp_density, imp_temperature, velocity_profile, + impurity.atomic_weight * atomic_mass) + species.append(Species(impurity, ionisation, imp_distribution)) + +# define the electron distribution +e_density = IonFunction(peak_density, edge_density, pedestal_top=pedestal_top) +e_temperature = IonFunction(peak_temperature, edge_temperature, pedestal_top=pedestal_top) +e_distribution = Maxwellian(e_density, e_temperature, velocity_profile, electron_mass) + +# define species +plasma.electron_distribution = e_distribution +plasma.composition = species + +# add thermal CX emission model +cVI_5_4 = Line(carbon, 5, (5, 4)) +plasma.models = [ThermalCXLine(cVI_5_4)] + +# trace thermal CX spectrum along y direction +ray = Ray(origin=Point3D(0.4, -3.5, 0), direction=Vector3D(0, 1, 0), + min_wavelength=112.3, max_wavelength=112.7, bins=512) +thermal_cx_spectrum = ray.trace(world) + +########################### +# Inject beam into plasma # + +integration_step = 0.0025 +# injected along x direction +beam_transform = translate(-0.5, 0.0, 0) * rotate_basis(Vector3D(1, 0, 0), Vector3D(0, 0, 1)) +beam_energy = 50000 # eV + +beam = Beam(parent=world, transform=beam_transform) +beam.plasma = plasma +beam.atomic_data = adas +beam.energy = beam_energy +beam.power = 3e6 +beam.element = deuterium +beam.sigma = 0.05 +beam.divergence_x = 0.5 +beam.divergence_y = 0.5 +beam.length = 3.0 +beam.attenuator = SingleRayAttenuator(clamp_to_zero=True) +beam.integrator.step = integration_step +beam.integrator.min_samples = 10 + +# remove thermal CX model and add beam CX model +plasma.models = [] +beam.models = [BeamCXLine(cVI_5_4)] + +# trace the spectrum again +beam_cx_spectrum = ray.trace(world) + +# plot the spectra +plt.figure() +plt.plot(thermal_cx_spectrum.wavelengths, thermal_cx_spectrum.samples, label='thermal CX') +plt.plot(beam_cx_spectrum.wavelengths, beam_cx_spectrum.samples, label='beam CX') +plt.legend() +plt.xlabel('Wavelength (nm)') +plt.ylabel('Radiance (W/m^2/str/nm)') +plt.title('Sampled CX spectra') + +plt.show() diff --git a/docs/source/atomic/atomic_data.rst b/docs/source/atomic/atomic_data.rst index 650d89b4..dcff3270 100644 --- a/docs/source/atomic/atomic_data.rst +++ b/docs/source/atomic/atomic_data.rst @@ -7,3 +7,7 @@ Atomic Data emission_lines rate_coefficients gaunt_factors + atomic_data_interface + data_interpolators + repository + openadas diff --git a/docs/source/atomic/atomic_data_interface.rst b/docs/source/atomic/atomic_data_interface.rst new file mode 100644 index 00000000..0d54ff67 --- /dev/null +++ b/docs/source/atomic/atomic_data_interface.rst @@ -0,0 +1,19 @@ + +Atomic Data Interface +===================== + +Abstract (interface) class +-------------------------- + +Abstract atomic data interface. + +.. autoclass:: cherab.core.atomic.interface.AtomicData + :members: + +OpenADAS atomic data source +--------------------------- + +Interface to local atomic data repository. + +.. autoclass:: cherab.openadas.openadas.OpenADAS + :members: diff --git a/docs/source/atomic/data_interpolators.rst b/docs/source/atomic/data_interpolators.rst new file mode 100644 index 00000000..56fa0327 --- /dev/null +++ b/docs/source/atomic/data_interpolators.rst @@ -0,0 +1,89 @@ +Atomic data interpolators +========================= + +The following classes interpolate atomic data defined on a numerical grid. + + +Atomic Processes +^^^^^^^^^^^^^^^^ + +.. autoclass:: cherab.openadas.rates.atomic.IonisationRate + :members: + +.. autoclass:: cherab.openadas.rates.atomic.NullIonisationRate + :members: + +.. autoclass:: cherab.openadas.rates.atomic.RecombinationRate + :members: + +.. autoclass:: cherab.openadas.rates.atomic.NullRecombinationRate + :members: + +.. autoclass:: cherab.openadas.rates.atomic.ThermalCXRate + :members: + +.. autoclass:: cherab.openadas.rates.atomic.NullThermalCXRate + :members: + +Photon Emissivity Coefficients +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +.. autoclass:: cherab.openadas.rates.pec.ImpactExcitationPEC + :members: + +.. autoclass:: cherab.openadas.rates.pec.NullImpactExcitationPEC + :members: + +.. autoclass:: cherab.openadas.rates.pec.RecombinationPEC + :members: + +.. autoclass:: cherab.openadas.rates.pec.NullRecombinationPEC + :members: + +Beam-Plasma Interaction Rates +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +.. autoclass:: cherab.openadas.rates.cx.BeamCXPEC + :members: + +.. autoclass:: cherab.openadas.rates.cx.NullBeamCXPEC + :members: + +.. autoclass:: cherab.openadas.rates.beam.BeamStoppingRate + :members: + +.. autoclass:: cherab.openadas.rates.beam.NullBeamStoppingRate + :members: + +.. autoclass:: cherab.openadas.rates.beam.BeamPopulationRate + :members: + +.. autoclass:: cherab.openadas.rates.beam.NullBeamPopulationRate + :members: + +.. autoclass:: cherab.openadas.rates.beam.BeamEmissionPEC + :members: + +.. autoclass:: cherab.openadas.rates.beam.NullBeamEmissionPEC + :members: + +Radiated Power +^^^^^^^^^^^^^^ + +.. autoclass:: cherab.openadas.rates.radiated_power.LineRadiationPower + :members: + +.. autoclass:: cherab.openadas.rates.radiated_power.NullLineRadiationPower + :members: + +.. autoclass:: cherab.openadas.rates.radiated_power.ContinuumPower + :members: + +.. autoclass:: cherab.openadas.rates.radiated_power.NullContinuumPower + :members: + +.. autoclass:: cherab.openadas.rates.radiated_power.CXRadiationPower + :members: + +.. autoclass:: cherab.openadas.rates.radiated_power.NullCXRadiationPower + :members: diff --git a/docs/source/atomic/openadas.rst b/docs/source/atomic/openadas.rst new file mode 100644 index 00000000..3dde4242 --- /dev/null +++ b/docs/source/atomic/openadas.rst @@ -0,0 +1,29 @@ +Open-ADAS +--------- + +Although a typical Open-ADAS data set is installed to the local atomic data repository +using the `populate()` function, additional atomic data can be installed manually. + +The following functions allow to parse the Open-ADAS files and install the rates of the atomic processes +to the local atomic data repository. + +Parse +^^^^^ + +.. autofunction:: cherab.openadas.parse.adf11.parse_adf11 + +.. autofunction:: cherab.openadas.parse.adf12.parse_adf12 + +.. autofunction:: cherab.openadas.parse.adf15.parse_adf15 + +.. autofunction:: cherab.openadas.parse.adf21.parse_adf21 + +.. autofunction:: cherab.openadas.parse.adf22.parse_adf22bmp + +.. autofunction:: cherab.openadas.parse.adf22.parse_adf22bme + +Install +^^^^^^^ + +.. automodule:: cherab.openadas.install + :members: diff --git a/docs/source/atomic/rate_coefficients.rst b/docs/source/atomic/rate_coefficients.rst index b6547820..82116bf3 100644 --- a/docs/source/atomic/rate_coefficients.rst +++ b/docs/source/atomic/rate_coefficients.rst @@ -15,6 +15,35 @@ provide theoretical equations. Cherab emission models only need to know how to c them after they have been instantiated. +Atomic Processes +^^^^^^^^^^^^^^^^ + +.. autoclass:: cherab.core.atomic.rates.IonisationRate + +.. autoclass:: cherab.core.atomic.rates.RecombinationRate + +.. autoclass:: cherab.core.atomic.rates.ThermalCXRate + +The `IonisationRate`, `RecombinationRate` and `ThermalCXRate` classes all share +the same call signatures. + +.. function:: __call__(density, temperature) + + Returns an effective rate coefficient at the specified plasma conditions. + + This function just wraps the cython evaluate() method. + +.. function:: evaluate(density, temperature) + + an effective recombination rate coefficient at the specified plasma conditions. + + This function needs to be implemented by the atomic data provider. + + :param float density: Electron density in m^-3 + :param float temperature: Electron temperature in eV. + :return: The effective ionisation rate in [m^3.s^-1]. + + Photon Emissivity Coefficients ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ @@ -39,8 +68,8 @@ the same call signatures. This function needs to be implemented by the atomic data provider. - :param float temperature: Receiver ion temperature in eV. :param float density: Receiver ion density in m^-3 + :param float temperature: Receiver ion temperature in eV. :return: The effective PEC rate [Wm^3]. Some example code for requesting PEC objects and sampling them with the __call__() diff --git a/docs/source/atomic/repository.rst b/docs/source/atomic/repository.rst new file mode 100644 index 00000000..3a93085e --- /dev/null +++ b/docs/source/atomic/repository.rst @@ -0,0 +1,83 @@ + +Atomic data repository +---------------------- + +The following functions allow to manipulate the local atomic data repository: +add the rates of the atomic processes, update existing ones or get the data +already present in the repository. + +The default repository is created at `~/.cherab/openadas/repository`. +Cherab supports multiple atomic data repositories. The user can configure different +repositories by setting the `repository_path` parameter. +The data in these repositories can be accessed through the `OpenADAS` atomic data provider +by specifying the `data_path` parameter. + +To create the new atomic data repository at the default location and populate it with a typical +set of rates and wavelengths from Open-ADAS, do: + +.. code-block:: pycon + + >>> from cherab.openadas.repository import populate + >>> populate() + +.. autofunction:: cherab.openadas.repository.create.populate + +Wavelength +^^^^^^^^^^ + +.. automodule:: cherab.openadas.repository.wavelength + :members: + +Ionisation +^^^^^^^^^^ + +.. autofunction:: cherab.openadas.repository.atomic.add_ionisation_rate + +.. autofunction:: cherab.openadas.repository.atomic.get_ionisation_rate + +.. autofunction:: cherab.openadas.repository.atomic.update_ionisation_rates + +Recombination +^^^^^^^^^^^^^ + +.. autofunction:: cherab.openadas.repository.atomic.add_recombination_rate + +.. autofunction:: cherab.openadas.repository.atomic.get_recombination_rate + +.. autofunction:: cherab.openadas.repository.atomic.update_recombination_rates + +Thermal Charge Exchange +^^^^^^^^^^^^^^^^^^^^^^^ + +.. autofunction:: cherab.openadas.repository.atomic.add_thermal_cx_rate + +.. autofunction:: cherab.openadas.repository.atomic.get_thermal_cx_rate + +.. autofunction:: cherab.openadas.repository.atomic.update_thermal_cx_rates + +Photon Emissivity Coefficients +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +.. automodule:: cherab.openadas.repository.pec + :members: + +Radiated Power +^^^^^^^^^^^^^^ + +.. automodule:: cherab.openadas.repository.radiated_power + :members: + +Beam +^^^^ + +.. automodule:: cherab.openadas.repository.beam.cx + :members: + +.. automodule:: cherab.openadas.repository.beam.emission + :members: + +.. automodule:: cherab.openadas.repository.beam.population + :members: + +.. automodule:: cherab.openadas.repository.beam.stopping + :members: diff --git a/docs/source/demonstrations/active_spectroscopy/CXS_multi_sightlines.png b/docs/source/demonstrations/active_spectroscopy/CXS_multi_sightlines.png index b04c94c0..f94e95ab 100644 Binary files a/docs/source/demonstrations/active_spectroscopy/CXS_multi_sightlines.png and b/docs/source/demonstrations/active_spectroscopy/CXS_multi_sightlines.png differ diff --git a/docs/source/demonstrations/active_spectroscopy/CXS_spectrum.png b/docs/source/demonstrations/active_spectroscopy/CXS_spectrum.png index 75b083d7..a258c6ed 100644 Binary files a/docs/source/demonstrations/active_spectroscopy/CXS_spectrum.png and b/docs/source/demonstrations/active_spectroscopy/CXS_spectrum.png differ diff --git a/docs/source/models/basic_line/basic_line_emission.rst b/docs/source/models/basic_line/basic_line_emission.rst index f8c80c60..b71fc03c 100644 --- a/docs/source/models/basic_line/basic_line_emission.rst +++ b/docs/source/models/basic_line/basic_line_emission.rst @@ -5,3 +5,5 @@ Basic Line Emission .. autoclass:: cherab.core.model.plasma.impact_excitation.ExcitationLine .. autoclass:: cherab.core.model.plasma.recombination.RecombinationLine + +.. autoclass:: cherab.core.model.plasma.thermal_cx.ThermalCXLine