Faster computation of Hybrid PRD coeffs

[Goobley]

This calculation method consumes more memory (temporarily), but reduces the time taken to compute the hybrid PRD coefficients by around a factor of 100 on the cases I am working on. The version in this branch is more RH style than that in my master.

[jaimedelacruz]

You can gain a lot more speed if you only compute the interpolation coeffs for wavelengths that actually correspond to PRD lines. At the moment you are computing it for all lambdas, but the interpolated value of J is only used for PRD/XRD lines. I did that in STiC at least, but it requires some extra book keeping.

[Goobley]

That's a good point which I considered, but I initially chose not to modify the implementation outside of this computation. For me this makes it fast enough for the wavelength grid I'm using (20s -> tens of ms) that it goes from being a significant overhead to insignificant, so I left it at that. Just computing this for the PRD lines would also allow you to save more time and some memory of course, which is never a bad thing! Feel free to close this if it isn't of much interest here -- I just thought I would offer it since it's been helpful for me.

[tiagopereira]

Thank you for this. Sorry for the delay it took me to look at this in more detail. I ported this to the rh15d version as it's easier for me to check there (it's basically the same, just plug the changes on initial_p.c).

I ran a few tests with a 5-level Mg II atom, but unfortunately I am not getting the same results with your fast calculation compared to the normal and "slow" way. In fact, for most cases the iteration either crashes or does not converge, but using snapshots from Bifrost or the FALC model. It is entirely possible I am doing something wrong, but for the time being I'd rather not merge this into the master.

[Goobley]

I completely forgot this was open. Sorry if the changes here didn't work, I was maintaining two forks of RH at the time. As Jaime pointed out we can gain substantially further by only looking at the lines in question. I did it this way in Lightweaver (implementation here), and I believe that this combination of linear searches around the line core to actually be the fastest approach (several orders of magnitude faster than what I provided here -- if this exact version worked or not!). In Lightweaver I also modified the formalism slightly to also recompute the formal solution for wavelengths that (at rest) don't "belong" to a PRD line but are Doppler shifted into the line's grid. I found some cases for which RH's HPRD would never fully converge due to this difference in effective Jbar between the full formal solutions and the PRD sub-iterations. If I am able to come back to this at some point I would probably follow this method. It's been relatively well tested against the full angle-dependent (using Han's RH) and default HPRD (using RH1.5D). I think we can close this for now 😄