Skip to content

Specific Concepts

James Matthews edited this page Sep 26, 2013 · 2 revisions

Stratified Sampling

Originally, Python utilized selected photon frequencies so that each photon bundle had the same energy. This approach is still utilized for calculating detailed spectra. This has the advanatge that one obtains the highest fidelity spectrum in the wavelength range where most of the energy is emitted. However, uniform energy sampling is not ideal for establishin the ionization fractions of ions with ionization potentials that a far from the peak of the energy distribution. These are important, especially in situations, where the energy sources have spectra which diverge significantly from a black body, as is the case for a CV with a boundary layer or an AGN.

To address this problem, Python allows one to force at least certain a fraction of the photons bundles generated to appear in user-definable wavelength/frequency intervals. Within an inteval, photon bundles are generated according to the energy distribution of the various sources of emission. If for example one defines, 3 bands, say 1 micron to 13.6 ev, 13.6 eV to the Helium edge, and the Helium edge to 1000 eV and requires f/3 of the photons to be in each interval, then exaclty 1/3 of the phtons will be generated in the interval. If instead one requires a minimum of 0.1 of the total number of photons to be in each interval, 70% of the photons will be distributed according to the ovrall energy distribution and this will be added to the 10% in each band.

When stratified sampling is employed, the weights (the energy carried by the bundle) is adjusted to reflect the actual energy that should have been released in the interval.

Radiative transfer

Ionization Calculation

Macro Atoms

The scheme in Python currently (Nov 09) is that developed for the YSO application to hydrogen. Although it is designed to be flexible and allow for multiple Macro Atom elements and simultaneously treating both Macro Atom and non-Macro Atom elements, neither of these things have so far been used in a publication and they have not been thoroughly tested:

1) what is a Macro Atom?

Following Lucy 2002, the term "Macro Atom" refers to a numerical scheme designed to deal with the excitation and de-excitation of atomic/ionic states subject to the constraint of statistical equilibrium. In the code, elements can be designated as being Macro Atom or non-Macro Atom - this determines whether or not this scheme is used. The decision is made via the input file and the atomic data set. Elements must currently be either completely Macro Atom or completely non-Macro Atom. If any level of an element is flagged as Macro Atom, any non-Macro Atom for the data should be ignored.

2) where is the action? - ionization cycles

The core of the Macro Atom implementation is in the files matom.c and estimators.c. The macro_gov and it's daughter routines matom and kpkt are the core parts of the implementation for the ionization cycles. Whenever our MC quanta interact with matter (absorbed by bound-bound line, bound-free continuum, free-free) the macro atom machinery kicks in. If the photon energy was removed from the propagating radiation field by excitation or ionization of an atom/ion then the packet is passed to "matom" which deals with the de-exciation / recombination problem to decide how that parcel of excitation/ionization energy should be dealt with. If the photon energy was transfered to the thermal pool (i.e. via electron kinetic energy) then it is sent instead to "kpkt". Conceptually, this means the energy went into heating the gas. Since we want thermal equilibrium, each piece of energy sent to "kpkt" (heating) must be balanced by cooling so kpkt's job is to processes the energy packets sent to it and send their energy back out by selecting a suitable cooling process. In fact, matom and kpkt need to communicate since means of depopulating atomic states (e.g. collisional de-excitation) involve sending energy to the thermal pool and, similarly, some cooling processes (e.g. collisional excitation) involve sending energy to the pool of exciation/ionization energy. This cross talk between macro atoms and kpkts is managed by macro_gov. When macro_gov is done, it will have converted the energy back to a photon packet (r-pkt) which can propagate again. Since we are explicitly forcing thermal/statistical equilibrium this means that all energy packets launched in the simulation eventually leave the computational domain. So there is no loss of packets in the simulation and no need to generate "wind photons" to represent the cooling radiation of the wind (as needed in the non-macro atom version of the code). 3) level populations and estimators The code uses the flight paths of the photons to compute volume-based MC estimators for the radiative rates of bound-bound and bound-free transitions of Macro Atom elements. This is done following Lucy 2003. These rates are needed for both the transition probabilities between states in matom and some of the cooling rates in kpkt. Also, they are used by "macro_pops" which is called during the wind_update part of the code. The function of macro_pops is to determine (more) accurate level populations for the Macro Atom levels, based on an inversion of the rate equations.

4) where is the action? - spectrum cycles

Much of what happens in the spectrum cycles is the same as in the ionization cycles. But there are some differences.

  • the estimators are no longer recomputed but (should be!) held fixed at their values from the last ionization cycle.
  • similarly, the level populations are not re-computed
  • since the spectrum cycles do not cover the entire wavelength range of the calculation, they cannot use the standard version of the macro atom routines - when photons only generated from a narrow range of energies, we cannot simulated the rate at which photons are fluoresced INTO that range from elsewhere in the spectrum. So, instead, a version of macro atom is used which works on a principle similar to the "wind radiation" form the original code version. Based on the last ionization cycle estimators, a calculation is performed to compute the rate at which de-exciting macro atoms and kpkts are contributing to the emission of radiation in the spectral interval we want to know about. These rates and then used to generate "wind photons" in a manner similar to the old way of making photons to account for the thermal radiation of the wind. These photons generated by macro atoms and kpkts and then launched by each grid cell in the wind and simulated along with the photons from the disk (and whatever other radiation sources are included). During the MC simulation, photons which are sent to either macro atom or kpkt are DESTROYED - provided the ionization cycles were converged, the rate of this destruction should balance the rate of macro atom/kpkt wind photon generation so that we have (statistically) the correct equilibrium

5) treatment of non-Macro Atom lines/continua

The code can deal with elements that do not have complete Macro Atom data at the same time as elements that do. As far as possible the treatment of these elements echoes that of the original code but there are some differences - especially for bf continua (which might be an issue right now - see below). For non-MA elements, bound-bound absorption is treated via a two-level atom with an effective escape probability (using the formulae in KSL's Python notes). The two-level atom can either radiatively de-excite (thus a resonance scattering happened) or it can collisional de-excite - in which case the packet is sent to kpkt (i.e. it went into heating the gas). For bound-free the difference from the old method is bigger and more subtle. Bound-free continua are also treated via a two-level approximation: upon absorbing a photon, the packet energy is either sent to kpkt (representing the energy carried off by the photo electron) or to fake_matom_bf which does the two-level approximation. Neglecting collisional recombination, the two-level approximation ONLY allows for a radiative recombination IN THE SAME CONTINUUM. Thus all it does is select a new photon frequency within the continuum and re-emit the photon. For trace elements, I suspect this is not a big deal but for something with important continuum edges (e.g. hydrogen itself) this is clearly a bad thing. So it may not be wise to use MAs ? only for minor species right now - H (at least) should always be considered. Specific issues/problems/concerns to bear in mind for now!

Some facts about the current implementation that one should bear in mind:

  • The point about treatment of non-Macro Atom bound-free continua (paragraph above) is likely important for some cases
  • Macro_pops (the routine which solves for level populations of macro atoms) was designed with hydrogen in mind. It has the problem that it only knows about atomic data that it gets from input. This has the negative consequence that it doesn't have any common sense about missing atomic physics - e.g. if there is no line data for two levels there will be no connection between them - thus, e.g., J-degenerate levels of complex ions will generally only be indirectly coupled via other states whereas we might physically expect them to be close to in Boltzmann equilibrium. There a kluge for this - if you set the lifetimes of macro atom levels to be very long in the input it will try to fudge their populations to be in Boltzmann equilibrium with the ground state.
  • get_matom_f (the routine which computes the macro_atom emissivities for the spectrum cycles) is rather inefficient. One could probably make it quite a bit better (and that may be worth doing since it's pretty slow right now - it's why the code appears to hang for a long time at the end of the ionization cycles).
  • stimulated recombination/stimulated emission - these are not really dealt with fully nor entirely self-consistently. They are generally included as negative absorption but are artificially knocked out in some place to make sure that the rates never go negative. Generally, I think we don't care too much about these processes it should be borne in mind if they are important for any particular application.