Differences between "simple-ion only " and "standard (weight reduction)" line transfer modes

[jhmatthews]

If one provides Python with standard77 atomic data - i.e. no macro-atom levels and lines - and specifies one of the macro-atom line transfer modes then the model should run in indivisible packet mode. The results of this run should, ideally, be the same as when run in standard line transfer mode, as all one is changing is the philosophy with regards to how we deal with certain continuum opacities.

I am finding there are some differences in the ionization state of a fiducial model when run with these two line transfer modes. Here's an example spectrum comparison for some high inclination viewing angles - the green is simple-ion calculation:

We see a change in the ionization state of the wind, especially at higher velocities which affects the line profiles.

My current suspicion is that this is due to use of Xsections. I believe that what currently happens is that first the code looks for a topbase cross-section, and if it doesn't exist for that ion it defaults to a VFKY one. When calculating heating, you then loop over all Xsections which have threshold frequencies lower than the frequency of the photon. My concern is that this doesn't seem to be done in simple-ion mode, and so I worry that we are not adding the contribution to bf heating from the ions with only VFKY cross-sections, which changes the temperature solution.

If this is the problem, it is clearly linked to #86 - different forms of sigma_phot.

[jhmatthews]

I'm pretty certain this is to do with VFKY cross-sections. If I put a very small factor of multiplication into the sigma_phot routine the result changes, but both modes then agree. The use of cross-sections is fairly wide-ranging so I'm finding it a little tricky to figure out exactly what is happening here - it's not just due to the loop which reduces weights and calculates bf heating.

Whatever the exact difference is, this seems like the appropriate time to implement the suggested change in #86 - basically have one photoionization structure containing a mixture of topbase and tabulated verner data, and one routine sigma_phot which returns the appropriate Xsection. We'll need to be careful with index and ordering here and how the older routines use the Xsections- but to be honest, that's probably something worth reviewing anyway.

[jhmatthews]

Bug!

in tabulate_verner(), we loop over the cross-sections

for (j=0; j < nxphot; j++)
    {   
    xver=&xphot[j];
    xphot_tab[j].z = xphot->z;
    xphot_tab[j].istate = xphot->istate;
    xphot_tab[j].nion = xphot->nion;
    f1=xver->freq_t*(1+very_small); //We need to start our tabulation just a tiny way up from from the threshold, otherwise it is equal to zero.
    f2=xver->freq_max*(1-very_small); //We need to start our tabulation just a tiny way up from from the threshold, otherwise it is equal to zero.
    lf1=log(f1);
    lf2=log(f2);
    dlogf=(lf2-lf1)/(N_VERNER_TAB);

    for (n=0;n<N_VERNER_TAB+1;n++)
        {
        xphot_tab[j].freq[n]=freq=exp(lf1+(n*dlogf));
        xphot_tab[j].x[n]=sigma_phot (xver, freq);
        }

    xphot_tab[j].np=N_VERNER_TAB;
    xphot_tab[j].nlast = -1;

    }

the lines which do something like xphot_tab[j].z = xphot->z; should read xphot_tab[j].z = xver->z; right, @Higginbottom ? This won't affect CV models because they used ML93 ionization and tabulated verner data is only used when we use spectral models. I suspect this doesn't make much difference to the overal ionization state anyway- not sure how much it affects the original issue yet.

[Higginbottom]

Yup - just taken a look at this and I agree its a mistake - and I think the proposed fix is correct!

[kslong]

Hmm — It seems to me as part of this process, we should resort to using either the tabulation or the formulae everywhere. Remind me what we were trying to accomplish with this.

[Higginbottom]

I agree completely. This little bit of code was mainly intended to speed up the execution of python, and in my normal conservative mode I only made changes to the bits of the code I was 'responsible' for. I see no reason why we shouldn't use tabulated data everywhere. This comes back to the question of whether we still do the tabulation in python, or have helper routines to make a generic photoionization data input file. We already tinker heavily with the PI cross sections (extrapolating topbase for example) outside.....

[jhmatthews]

I also agree. If we tabulate it we can have one structure which contains all photoionization data, and one function which you always call. It's simpler, quicker, and we don't really lose anything. I'm writing some code to tabulate outside anyway at the moment, which seems as good as approach as any.

We can discuss further tomorrow, and I think a bit of a review of how we use XSections and their level information would be good.

[jhmatthews]

These two modes now agree pretty well. Slight differences in ionization state but spectra very similar. This was basically originally caused by not looping over the full set of cross-sections in the macro-atom heating and cooling calculations which led to small differences, and there are some small differences between the schemes anyway. Closing.