burgess_tully_descale#
- fiasco.util.burgess_tully_descale(x, y, energy_ratio, c, scaling_type)[source]#
Convert scaled Burgess-Tully [BT92] parameters to physical quantities.
For a scaled temperature, \(x\) and scaled effective collision strength \(y\), the effective collision strength can be calculated as a function of the scaled energy \(U=k_BT_e/\Delta E_{ij}\) the ratio between the thermal energy and the energy of the transition \(ij\).
There are 6 different scaling types, depending on the type of transition. This scaling is explained in detail in section 5 of Burgess and Tully [BT92]. The scaled temperatures and collision strengths are related to \(U\) and \(\Upsilon\) by,
type 1
\[x = 1 - \frac{\ln C}{\ln{(U + C)}},\quad y = \frac{\Upsilon}{\log(U + e)}\]type 2
\[x = \frac{U}{U + C},\quad y = \Upsilon\]type 3
\[x = \frac{U}{U + C},\quad y = (U + 1)\Upsilon\]type 4
\[x = 1 - \frac{\ln C}{\ln{(U + C)}},\quad y = \frac{\Upsilon}{\log(U + C)}\]type 5
\[x = \frac{U}{U + C},\quad y = \Upsilon U\]type 6
\[x = \frac{U}{U + C},\quad y = \log_{10}\Upsilon\]
where \(C\) is a scaling constant that is different for each transition. Note that Burgess and Tully [BT92] only considered scaling types 1 through 4. Types 5 and 6 correspond to dielectron and proton transitions, respectively.
To “descale” the scaled effective collision strengths that are stored in the database, a spline fit is computed to the new \(x\) as computed from \(U\) and then the relationship between \(\Upsilon\) and \(y\) is inverted to get \(\Upsilon\) as a function of \(U\).
- Parameters:
x (
array-like
) – Scaled temperature. First dimension should have lengthn
, the number of transitions. The second dimension will be the number of spline points, but may be different for each row. If each row hasl
spline points,x
should have shape(n,l)
. If they are not all equal,x
will have shape(n,)
.y (
array-like
) – Scaled collision strength. Must have the same dimensions asx
.energy_ratio (
array-like
) – Ratio between the thermal energy and that of each transition with shape(n,m)
, wherem
is the dimension of the temperature array.c (
array-like
) – Scaling constant for each transition with shape(n,)
scaling_type (
array-like
) – The type of descaling to apply for each transition with shape(n,)
. Must be between 1 and 6
- Returns:
upsilon (
array-like
) – Descaled collision strength or cross-section with the same shape asenergy_ratio
.