fbartolic / caustics Goto Github PK

Currently caustics supports only the linear limb-darkening law as a model for the intensity profile of the source star disc. This is implemented in _brightness_profile. It would be useful to allow for an arbitrary intensity profile by instead evaluating the intensity of a starry map. This would enable computation of microlensing light curves with spots etc. . The relevant reference is the paper The microlensing signatures of photospheric starspots.

@rodluger would you be interested in this? A basic implementation would involve switching _brightness_profile with a JAX function which takes a vector of spherical harmonic coefficients as an input and return the flux at arbitrary (r, theta) as an output. This seems like a pretty trivial task no?

We can either write the whole thing in JAX or write a custom op which wraps the Boost implementation available here. It should be possible to use the Leibniz integral rule to derive a custom gradient for the definite
integral of f(x) from a(x) to b(x).

It's not clear that this would necessarily be faster than the current approach but we'd get error control.

No idea how to fix this. I've tried to install caustics on an M1 Mac dozens of times and succeeded only once with some specific permutation of compiler, Python and JAX versions. It might just be easier to rewrite the root solver entirely in JAX than to fix this.

The tests from Bozza et. al. 2018 are quite wasteful in the sense that there are tons of false positives. Here's an example:

The red region is where the error is > $10^{-4}$ and the grey region is where the tests says that the full calculation is necessary. I feel like there has to be a way to improve the accuracy of these tests.

Perhaps we should treat this as a machine learning prediction problem and construct physically meaningful features such as
the magnitude of the higher order terms in the multipole expansion and the value of Jacobian.

Related to #19.

This means using $(s,q)$ as binary lens parameters and possibly also shifting the centre of the coordinate system to the centre of mass. Need to look into triple lens parametrizations.

This will require major changes to how the code is organised and the use of bounded while loops. The adaptive Gauss-Kronrod algorithm for estimating the LD integrals (see #13) to a specific precision has to be implemented first.

Being able to specify a target relative precision will substantially improve performance with large sources.

Neither the cusp test nor the "ghost images" test from Bozza et.al. 2018 seem to be working for triple lenses. I must be doing something wrong when computing the derivatives of $f(z)$.

A related issue is that as it is currently implemented, the ghost images tests checks that
$$c_g(|I(z_1)| + |I(z_2)|)<\delta$$
where $I$ is an indicator function that evaluates to true if the condition is satisfied, rather than
$$c_g(|I(z_1)|)<\delta\quad \text{and} \quad c_g(|I(z_2)|)<\delta$$
The prohibited regions in the source plane don't look right with the latter condition.

Only way to fix these issues is to re-derive all the math from scratch and check that it all makes sense.

I currently use jit at the level of smallest functions. Performance could be improved by exploring the effect of JIT at higher level functions.

This is just a double integral in the source plane.

In very rare cases the two lines formed by the endpoints of two segments are almost exactly parallel and due to precision issues the intersection points between the two lines can end up being far away from the endpoints of either segment. I need to think about how to differentiate this situation from a situation when the two lines are almost exactly parallel but they do not belong to the same part of the contour.

I set up the coordinate system such that lens 1 is at $a$ and lens 2 is at $-a$ which is stupid because the coefficients are simpler when one lens is at the origin. This should speed up the computation of the coefficients and hence the magnification calculation.

Installing caustics on Collab with a GPU runtime and running the root solver results in an immediate crash of the kernel. Not sure what's going on here, could be a memory issue perhaps.

The quadrupole and hexadecapole corrections to the point source magnifications are insignificant when approach a fold from the the outside of a caustic. Bozza et.al. 2018 devised a test (Eq. 47 in the paper) to catch when the source is close to a fold and trigger the full calculation at those points. I've implemented this test but the output looks like this:

No idea what these weird curves emerging out of the fold near the top are. Could be a mistake in my implementation.

Compilation times are OK when using caustics on its own but they're completely unacceptable when caustics is incorporated within a numpyro model.

This should live in trajectory.py.

Using JAX for this problem was a mistake. I'm building a Julia version of the caustics code focused exclusively on binary and triple-lens events. This repository will mostly consists of Python bindings to the Julia code.

Hi Frac,

As we discussed, the semi-analytic root solver may substantially faster than the AE solver. Let's see if we can implement it in JAX and see how it works out.

https://github.com/kmzzhang/analytic-lensing
https://arxiv.org/abs/2207.12412

Keming

There's a slight positive bias in the extended source magnification (less than $10^{-4}$ in magnitude):

The plot is error relative to VBB but the hexadecapole approximation confirms that the bias is real.

Equation 24 in https://arxiv.org/abs/astro-ph/9804059

There should be floating point issues with lens equation evaluation when the third lens is far away from the others. Or something like that.