agnwinds / python

This is the repository for Python, the radiative transfer code used to winds in AGN and other syatems
GNU General Public License v3.0
26 stars 25 forks source link

nsh 130718 - Free free issues #29

Closed Higginbottom closed 11 years ago

Higginbottom commented 11 years ago

This issue related to free free heating, and was first observed in python 76a (parallel).

It was noted that the free free heating was rather high in proga models - dominating in large parts of the wind. This seemed odd, since free free cooling did not balance it, rather compton cooling did this at a very high temperature.

Examination of the code has revealed several issues, listed below in order of apparent seriousness.

OK, so I have tracked down the reason for 1 - it turns out to be a rather more serious bug that I initially thought - the sum over ions in pop_kappa_ff_array was being reset in the wrong loop - this means the 'prefactor' went up and up and up. This could potentially affect the benchmark model, so I'm going to fix it, and rerun the benchmark model. Lets hope that the free free was never that important, otherwise we are in trouble.

I'll go ahead and fix bug 2 as well. Since that is simple...

1: The plasma variable kappa_ff_factor seems to be unexpectedly changing by very large amounts after it has been computed . This is the sum of gaunt factors times ionic abundances, used as the 'prefactor' for calculating the free free opacity, and so errors will drastically change the free free opacity, and hence heating. It is computed just before the photons are proagated thruogh the wind at line 1716 in python.c with a call to pop_kappa_ff_array. Examination in the subroutine shows that the value computed is plausible, however, the value contained in the variable seems to change on return to python.c and before trans_phot. This does not happen in earlier, non parallel versions of python, so it seems likely there is an issue of things being overwritten. This bug gives rise to errors of many orders of magnitude.

2: There is (I think) a simple error in the implementation of free free with the sutherland gaunt factors. This affects opacity (heating) but will also affect cooling. In each bit of code, a sum is carried out over ions. This effectively replicates equation 18 in sutherland (1998) in the case of cooling. This summation includes neutral ions as currently implemented - but it shouldn't. It also uses ion[].Z (the atomic number) as Z in the equation. Sutherland uses Z as the ionization stage of the ion being summed over - in other words the charge. So Z=0 for H1, Z=1 for H2 etc. This will give rise to errors of order Z^2 for everythong other than Hydrogen. Given the fact that HII is normally the most abundant ion by two orders of magnitude, the overall error should not be huge, unless you have a very cold, neutral plasma in which case free free should not occur at all. Fortunately, a simple bug to fix.

3: Whilst the use of gffint (the gaunt factor integrated over frequency for a maxwellian distrubution of electrons) is probably fine for the bulk free free emission, I now believe we should be using the gffgu (gaunt factor for a given frequency for a maxwellian distribution of electrons) for opacity and also producing individual photons. This might be a little tricky to implement, and the change should be relatively small, since the gaunt factor tends to be of order 1 - but I'll quantify this better in the future once problems 1 and 2 are fixed!

Higginbottom commented 11 years ago

Fixed the first two issues. Issue 3 remains open, and should be addressed

jhmatthews commented 11 years ago

I just thought I would flag up here how this bug affects what we get out of python. Here is the output spectrum from python 75, with the new fixed version in green and the older version in blue. Barely any difference in the spectra because the lines mainly come from the scattered contribution.

run110e_compsummary

But, interestingly, the thermal emission from the wind is much lower:

run110e_compsources

This is simply because the base of the wind was hotter than it should be due to incorrect ff heating, meaning the thermal emission was higher. This may mean that we are further away from producing BELS than we think.

kslong commented 11 years ago

Nick, I see, you at one point closed this issue and then reopoened and closed it again. Can you explain why you did this, so if we can understand what really is the situation. Is the the 3rd part of the issue fixed, or not necessary, or just deferred. If it is just deferred, then could you open another issue for that and note here that this is why this particular issue is closed. Knox

Higginbottom commented 11 years ago

This bug closed - point 3 not yet addressed and reopened as bug #33