4.1.6. Clus: projected cluster model¶
Warning
This is a complex model, so please read carefully this documentation before using the model
4.1.6.1. Overview¶
This model calculates the spectrum and radial profiles of a cluster of galaxies. The 3-D density, temperature, abundance, turbulence and in/outflow fields are parameterized. The cluster is subdivided into 3D-shells; the number of shells can be adjusted. For each shell, the spectrum is calculated using the collisional equilibrium (CIE) model for the proper parameters.
The shell spectra are projected on the sky, and the user can choose an inner and outer projected radius between which the spectrum is calculated. By choosing different radii, the spectra of different regions of the cluster can be evaluated.
The model calculates the total spectrum, and this spectrum can be fit like any other SPEX model. At this moment, there are two ascii-output options associated with the model.
The first option (“clus”) provides output of radius, hydrogen density, temperature, gas pressure, emission measure, turbulent velocity, outflow velocity and relative abundances. The second option (“clup”) provides the projected radial profile of the number of photons and energy flux, in an energy band specified by the user.
Warning
this model has a lot of free parameters (currently 80). Do not blindly thaw all parameters, but make as much as possible use of other constraints and common sense.
4.1.6.2. Two important radii¶
It is important to distinguish two important radii in the model.
Internally, everything is scaled to , which is the first parameter of the model (input name r500). It should be given in units of meter. The default value is 3 (i.e., m, or approximately (but not exactly!)) 1 Mpc.
The cluster is terminated at a radius (input name: rout). Outside this radius, the density is cut-off to zero, so this means the effective edge of the cluster is at . It should be provided in units of , and its default value is 2 (i.e., or m or approximately 648 kpc).
Warning
If you really want to go to the outermost cluster outskirts, and you increase , then you should take care to adjust possibly the number of shells (default 25) and/or projected annuli (default 50), otherwise the spatial resolution of your model may be insufficient.
The two parameters governing the radial scaling are summarized below.
Acronym |
Variable |
Default value |
Units |
r500 |
3.0 |
m |
|
rout |
2.0 |
4.1.6.3. Radial grid: 3-d shells¶
The 3-d model is evaluated for a logarithmic radius grid, in order to cover the possible large dynamical range between the core of the cluster and the outermost radius. The user specifies (input name nr), the number of shells for which the model is evaluated.
The default value for nr is 25, but it can be changed (obviously, for large values the computation time increases). The logarithmic grid is chosen to be beween and . Here is taken as 1% of , or if that is smaller, 10% of the smallest of the core radii for the density profile (see later). Obviously, for the innermost shell the inner boundary is replaced by 0 instead.
The hydrodynamical parameters like density, temperature etc. are then evaluated for the middle radius of each shell using analytical parameterisation of the profiles (see later), and we assume that these quantities are constant within the shell.
Therefore, the emission measure of each shell is simply evaluated as the product of electron and hydrogen density at the central radius of the shell, multiplied by the volume of the shell; we do not perform finer integration.
In case the user introduces a density discontinuity in the model, the radius of the density discontinuity is also chosen as an additional grid point (if it is close to an existing grid point, that existing grid point is moved to the discontinuity).
4.1.6.4. Radial grid: projected profile¶
While as said we take a logarithmic grid for the 3-d shells of the cluster, to cover the large dynamical range, we take a linear grid for the projected 2-d profile. This is because that projected profile is often closely linked to a real detector, which has a linear scale.
The default number of projected annuli (input name: npro) is 50, but again that number can be changed.
While for the 3-d shells we only evaluate quantities at the shell central radii, in evaluating the projected profiles we do some sub-shell evaluation (typically using 20 points for each shell). This is in order to avid a kind of “saw-tooth” projected profiles, which would result if we simply kept the emission measure constant within the shell. Within the shell, we adopt a power-law behaviour of the emission measure as a function of radius. The power-law slope is simply determined by comparing the emission measures of the neighbouring shells.
Warning
If you adjust the parameters for the model, take care that the radial bin sizes of the outermost shells are compatible in magnitude with the bin sizes of the projected profiles. You can see both grids by using the ascii-output options for the “clus” and “clup” options, respectively.
The parameters governing both the 3-d and projected bins are summarized below.
Acronym |
Variable |
Default value |
Units |
nr |
25 |
||
npro |
npro |
50 |
4.1.6.5. Radial density distribution¶
The hydrogen density distribution is parametrised as the sum of two beta-models, modified by a possible density jump (shock):
with for
where is the central density in units of per cubic meter for both components (input names hd1 and hd2), are the core radii (rc1 and rc2 in the input) in units of and the beta parameters (bet1 and bet2 in the input).
The optional density jump parameterisation is:
where the discontinuity radius (input name rsh) is given again in units of and further is the density jump at (its input name is dfac). Please note that, in the case of no jump (i.e. ), the term is ignored. For the density increases outside relative to the undisturbed model, while for it decreases. Further fine-tuning of the discontinuity can be achieved by playing with the power-law parameter (input name dgad). Positive values means increasing effects at large radii, z zero value (the default) a constant jump, and negative values correspond to a decreasing effect of the jump at large radii.
The parameters describing the density profile are summarized below.
Acronym |
Variable |
Default value |
Units |
hd1 |
5000 |
||
rc1 |
0.1 |
||
bet1 |
0.7 |
||
hd2 |
25000 |
||
rc2 |
0.01 |
||
bet2 |
1.0 |
||
rsh |
2.0 |
||
dfac |
1.0 |
||
dgad |
0.0 |
4.1.6.6. Radial temperature distribution¶
The temperature distribution is based upon Vikhlinin et al. (2006) (Eq. 6) (2006ApJ…640..691V) but with a temperature jump added and described below.
The Vikhlinin parameterisation is the product of two terms. For the first term, describing the inner part of the cluster ( below), we rewrite his expression to an equivalent but slightly different form, as given by Kaastra et al. 2004, eq. 10 (https://www.aanda.org/articles/aa/pdf/2004/02/aah4230.pdf), which in turn is based on Allen et al. 2001.
With this, the temperature profile can be written as
with
The central and outer temperatures and are not the actual temperatures, but the temperatures that would exist without the and terms.
Similar to the density profile, we allow for a possible temperature jump of any magnitude (i.e., increase or decrease) at the radius (which is the same radius where a possible density jump may occur).
Most default parameters are based on Vikhlinin et al. Eq. 8 for the typical temperature profile that they derived.
We assume that the ion temperature equals the electron temperature for each shell. The parameters describing the temperature profile are summarized below.
Acronym |
Variable |
Default value |
Units |
tc |
1.8 |
keV |
|
th |
4.0 |
keV |
|
rtc |
0.045 |
||
mu |
1.9 |
||
rto |
0.6 |
||
at |
0.0 |
||
bt |
2.0 |
||
ct |
2.0 |
||
tfac |
1.0 |
||
gamt |
0.0 |
4.1.6.7. Radial abundance profile¶
The relative metal abundances can be adjusted similar to the CIE model (see the last parts of this chapter). But these abundances can be multiplied by a radial scaling law .
The scaling law, including the default parameters, is taken from Mernier et al. (2017) https://www.aanda.org/articles/aa/pdf/2017/07/aa30075-16.pdf However, we write it in a slightly different, equivalent form as follows:
If the abundances should be kept constant as a function of radius, the user should take care that for all radii. This can be achieved for instance by setting A=0, G=1.
The constant term with G was not included in Mernier et al. (2017), but may be useful for some applications. Its default value is 0.
Note that the radial scaling works the same way for all chemical elements; the abundances themselves (as given by the parameters 01, 02, 03 etc.) can of course be different.
Warning
The radial scaling is only done for elements with nuclear charge 3 or more (i.e., Li and higher). Hydrogen and helium are excluded from the scaling. It is therefore highly recommended to keep the reference atom to its default (Hydrogen) for this model, and not to use another one.
Warning
If you do a spectral fit of a spectrum that is for example dominated by iron lines, take care that the parameter afe and the iron abundance 26 are strongly correlated. In those situations it may be wise to keep one of them frozen.
The parameters describing the abundance profile are summarized below.
Acronym |
Variable |
Default value |
Units |
afe |
A |
1.34 |
|
bfe |
B |
0.021 |
|
cfe |
C |
0.48 |
|
dfe |
D |
0.414 |
|
efe |
E |
0.163 |
|
ffe |
F |
0.0165 |
|
gfe |
G |
0.0. |
4.1.6.8. Turbulence¶
The plasma turbulence can also be parametrised. We use the following description:
with
Here has the same meaning as the vrms parameter of the CIE model of SPEX (see there for more documentation).
At the center, v is given by , at large distances it is given by Note that due to this definition, positive values of means increasing turbulence for larger radii, while negative values means decreasing turbulence for larger radii.
Warning
If you do a spectral fit or evaluate the model, take care that for negative values of with smaller than , the turbulent velocity would become imaginary (). To avoid a crash, SPEX cuts these values off to zero, and continues, but be aware that your results are not realistic anymore in those situations.
Warning
The turbulent velocity here is only the component due to microscopic or macroscopic motions of the gas cells. In addition, the thermal motions of the ions will be added in quadrature using the relevant (ion) temperature.
The parameters describing the turbulence profile are summarized below.
Acronym |
Variable |
Default value |
Units |
av |
100 |
km/s |
|
bv |
0 |
km/s |
|
rv |
1 |
4.1.6.9. Radial velocity fields¶
In addition to turbulence, the model allows for systematic radial motion (inflow or outflow). This is done as follows:
with
Positive values correspond to inflow, negative values to outflow. The parameter corresponds to the flow velocity at the core, is the flow velocity at large distances. Note that the default values for these parameters correspond to no in- or outflow.
Acronym |
Variable |
Default value |
Units |
zc |
0 |
km/s |
|
zh |
0 |
km/s |
|
rz |
1 |
4.1.6.10. Projected spectra and radial profiles¶
The cluster model projects the spectra of all shells onto the sky. But sometimes the user is not interested in the spectrum of the full cluster, but only in the spectrum within a projected annulus.
One can specify the inner and outer radius and for this projected annulus annulus. The spectrum is then only calculated between those boundaries.
Contrary to all other radii, we specify these radii in units of , the outer radius of the cluster. Thus, for and , the full cluster spectrum is obtained.
On the other hand, the “clup” output option allows to display the radial intensity profile of the cluster for a given energy band. These energies (always in keV) can be specified by the parameters and .
It may be that the user wants the spectrum or profile in a more complex region than a simple annulus, for instance the projection of a square detector (pixel). The cluster model provides an option for this by using the parameter with name “azim”. When this parameter has value 0, the full annulus is used. When it has a value of 1, a more complex region can be used. In this latter case, the user must specify the name of an ascii-file (parameter fazi).
This file must have the following contents:
first line: number of points
next lines: for each specified radius, from small to large radii, two numbers, namely the radius and the azimuthal fraction .
The radius must be in units of , the azimuthal fraction should be between 0 (no contribution of this annular radius) and 1 (the full annulus is contributing.) The spacing between the radial points should be such that linear extrapolation between the points gives a fair representation of the shape of the extraction region.
For example, if the extraction region is for example a quarter of the sky counted from the cluster center (e.g, everything between 0 and 90 degrees azimuthal angle), the content of the file could be simply three lines:
2
0 0.25
1 0.25
If the first radius is larger than 0 or the last radius less than 1, we adopt (no emission) for the non-specified regions.
An error message will be given if the file is incorrect, and in that case the model is evaluated as if the flag azim=0 (i.e., everywhere).
We summarise the parameters below.
Acronym |
Variable |
Default value |
Units |
rmin |
0 |
||
rmax |
1 |
||
emin |
0 |
keV |
|
emax |
keV |
||
azim |
0 |
||
fazi |
4.1.6.11. Resonance scattering¶
Resonance scattering of the strong resonance lines can be taken into account, if the user needs this. This is indicated by the parameter with name rsca. Allowed values are 0 (no resonance scattering) or 1 (resonance scattering).
The resonance scattering is calculated using a Monte Carlo approach. It is therefore less suited for spectral fitting, due to the random effect inherently caused by this technique.
The initial 3D emissivity profile for the relevant lines is taken as starting point. From this profile, N x nr random photons are drawn (here nr is the number of shells introduced earlier, and N is an adjustable number that can be set by the user, see the table below).
There are two modes: for mode=2, each shell starts with the same number of photons (the final result is of course weighted with the proper shell emissivity). This has the advantage that the calculation also have a reasonable accuracy in the outer parts, where the emissivity is small. It is therefore most relevant when studying radial line profiles of the resonant lines in the outskirts. For the default mode=1, the number of initial photons is distributed according to emissivity; therefore this case is more suited for the overall spectrum of the full cluster. Experience shows it is slightly more accurate in the part of the cluster that emits the bulk of the photons.
The fate of each individual photon is followed. The calculation for the photon stops when either it is absorbed in the continuum, or when it leaves the cluster. Alternatively, it can be absorbed and then a) re-emitted in a new random direction (the resonance scattering), or b) it decays to a non-ground level, resulting in two or more photons until the atom reaches the ground state again. The fate of these multiple photons is followed until they are destroyed or escape. At the end of the calculation, some statistics is collected on the history of the photon.
When resonance scattering is included, in principle three diagnostic files with information on the resonant lines are
produced. These files are always named cluslin1.asc
, cluslin2.asc
and cluslin3.asc
, they are placed in
the directory from which SPEX is running, and they are overwritten each time the model is being evaluated. When the
user wishes to store these files, a simple spex command like sys exe "cp cluslin1.asc newa1.asc"
is sufficient
to copy it to a file newa1.asc
, for example.
If the user does not want to produce this output, the parameters out1, out2 and out3 should be set to zero (see below). The default value for them is 1, i.e. produce the output. If some of these numbers are put to zero, no output of that kind is produced.
In most other cases of SPEX output, we make use of the asc-output options. However, in particular when the number of shells is large or the number of iterations N is large, the computational time can be very large, and with different random seeds for each run results may not be fully reproducible. For that reason, the resonance scattering line information is given “on the fly”.
The relevant parameters are listed below.
Acronym |
Variable |
Default value |
Units |
rsca |
0 |
||
nit0 |
N |
10000 |
|
rmod |
mode |
1 |
|
out1 |
out1 |
1 |
|
out2 |
out2 |
1 |
|
out3 |
out3 |
1 |
Warning
When running resonance scattering, be aware that the computation times may become very large. For example, a run with 128 radial shells, 256 projected annuli, with cluster parameters adjusted to the Perseus cluster, and with nit0 equal to and thus random drawings per spectral line (and 666 spectral lines included) takes about 6 hours to complete on a single core of a 2.3 GHz processor using the gfortran compiler; it may be faster using other processors or compilers. For less ambitious accuracy a couple of minutes per model evaluation is more realistic.
4.1.6.12. Output files for resonance scattering¶
As noted, in principle three files are produced when resonance scattering is included in the calculations.
These file cluslin1.asc
, cluslin2.asc
and cluslin3.asc
are stored in the directory where SPEX is running,
and are overwritten each time the model is evaluated again. Here we describe the contents of these files.
cluslin1.asc:
This is a summary of all the spectral lines where resonance scattering is taken into account. The file contains a descriptive header (one line) followed by several lines for each spectral line. The columns are as follows:
Line number. This is a unique number for each line, which is used as a label for the
cluslin2.asc
andcluslin3.asc
files to identify the spectral lines. These numbers are always the same for a given transition, but note that the numbering is different for the “old” set of atomic lines and the “new” set of atomic lines.Element and ionisation stage.
Designation of the transition. For the “old” calculation, typically the Mewe et al. (1985) notation, for the “new” calculation, it gives the lower and upper level of the transition. Both notations are the same as for the line-ascii output option.
The energy of the line, in keV.
The oscillator strength of the line. However, for lines not ending at the ground level (the non-resonant lines), we put here a value of zero, even if the true oscillator strength is non-zero. The cluster model assumes that all ions are in their ground state when they absorb a photon.
Branching ratio. This is the probability that an atom that has absorbed a resonant photon falls back to the ground state immediately, rather than to an intermediate energy level.
tau_l: the optical depth of the line at line centre.
tau_c: the continuum optical depth (i.e., the optical depth relevant for absorption of the photon in the continuum).
photons: the number of photons that leave the cluster (either as original, scattered or splitted photon).
lost: the number of photons lost by continuum absorption. This number, added to the number of escaping photons gives the original number of photons for this line for the Monte-Carlo calculation.
ori_flux: the original emitted flux (number of photons per second) in the line, without resonance scattering.
cluslin2.asc:
This file gives the projected properties of all the lines summarized in cluslin1.asc. The file contains a descriptive header (one line) followed by npro lines for each spectral line. The columns are as follows:
line number (corresponding to the first number in
cluslin1.asc
)ir: the projected annulus number
r: the radius of the projected annulus in m. The radius is here the mean of the inner and outer radius of the annulus.
fluxscale: the scaling factor by which the number of photons in the next three columns (5 to 7) needs to be multiplied in order to get fluxes in photons per square meter per second, emitted at the cluster (not to be confused with the flux received at Earth!).
f1: the number of photons that would have been visible in this annulus if no resonance scattering occurs.
f2: the actual number of photons that leave the cluster from this annulus to the observer.
f2-f1: the difference f2-f1.
f2/f1: the ratio between f2 and f1. This is actually the relative flux decrease or increase due to resonance scattering.
# scat: The average number of scatterings that the observed photon has experienced, since creation of the photon.
DelE: the average energy shift of the photons, in eV.
sigma: The r.m.s. line width in eV. This incorporates the original line width and all effects of the resonance scattering. Note that some lines may have Lorentzian wings in reality; those are visible in the final spectrum. Voigt profiles are computed there based on this r.m.s. and the original natural broadening of the line. The Gaussian component of the line therefore ignores any non-gaussianity, but only considers the width of the line, which is accurate enough in the majority of cases.
cluslin3.asc:
This file gives the 3D, radial properties of all the lines summarized in cluslin1.asc
.
The file contains a descriptive header (one line) followed by several lines for each spectral line.
The columns are as follows:
line number (corresponding to the first number in
cluslin1.asc
)ir: the 3D shell number
r: the radius of the shell in m. The radius is here the mean of the inner and outer radius of the shell.
T: the temperature of the shell at mid-radius
fr: the fraction of the total line flux of the cluster emitted from this shell (without resonance scattering)
concentrat: the ion concentration of the emitting ion at this radius.
4.1.6.13. Atomic physics of the resonance lines¶
For the “old” atomic data set, we used a list of 204 selected lines deemed to be relevant in 1996, when the first version of this model was created. For the “new” atomic data set, we selected 666 transitions for H-like to Na-like ions, under the following assumptions:
We use a simple scaling law where is the total column density from the core of the cluster to infinity, in units of per square meter, and the temperature in keV. This matches approximately the Perseus cluster at a temperature of 4 keV. It is based on the simple scaling laws for cluster mass M, radius R and density
, constant, and
For a grid of temperatures T from 0.25 to 16 keV, step size a factor of 2, assuming no turbulence (i.e., the highest possible optical depths) and a column density of 5 times higher than the scaling relation given under 1. above, we calculated the optical depths of all ground-state absorption lines in SPEX using the “hot” absorption model.
From this list, we selected all lines with energy larger than 0.1 keV, optical depth larger than 0.03, and optical depth larger than 1% of the strongest line of that ion at that temperature.
For each resonance line, we only considered the most important decay channels from the excited state (typically stronger than 1% of the total decay rate) in order to limit the number of decay routes to the ground to a manageable number. In practice, the number of decay routes was thus limited to a maximum of 7 (including the main resonance line), and we limited the number of lines to a maximum of 3 per route (ignoring the lowest energy lines in a few cases; these ignored lines were always at much lower energy than the X-ray band).
Transition energies, oscillator strengths and branching ratios were thus obtained from the current atomic data base of SPEX. When this database will be updated, the same transitions will be used for the resonance scattering, but branching ratios, oscillator strengths and line energies will automatically be adjusted; these quantities are initialised each time upon the first call to the cluster model.
4.1.6.14. Non-thermal electron distributions and abundances¶
The effects of non-thermal electron distribution on the spectrum can be included like any SPEX model. See Non-thermal electron distributions for more details.
We refer here to the description of the CIE model for the abundances. We use the same parameters as in that model. See description below.
4.1.6.15. Parameter description¶
The remaining parameters of the model are:
ref
: reference element. Default value 1 (hydrogen). See above for
more details. The value corresponds to the atomic number of the
reference element.01
: Abundance of hydrogen (H, Z=1) in Solar units. Default 1.02
: Abundance of helium (He, Z=2) in Solar units. Default 1.03...30
: Abundance of all elements up to zinc (Zn, Z=30) in Solar units.
Default 1.file
: Filename for the non-thermal electron distribution. If not
present, non-thermal effects are not taken into account (default).Recommended citation: Kaastra et al. (2022) (this manual for the moment).