rodluger / Limbdark.jl

Analytical transit light curves for limb darkened stars
MIT License
14 stars 4 forks source link

Derivatives of special cases of s_2 #7

Closed ericagol closed 6 years ago

ericagol commented 6 years ago

Since the special cases such as b=0, b=r, and b=1-r only depend on one variable, the computed partial derivatives of s_2 are incorrect.

A couple of options exist: 1). We could add the Taylor series expansion to the evaluation, and then use AutoDiff (but this wouldn't allow arbitrary order to be computed). 2). We could compute the derivatives analytically for these cases (these are already given in Pal 2008).

ericagol commented 6 years ago

I worked out some (all?) of these derivatives for the s_2 special cases. I still need to check that these are correct, and code up. The other cases are either trivial (derivatives are zero), or they can be differentiated with auto-diff (or we could just work out the full, complicated expressions).

image

\frac{d \Lambda}{d b} &= \begin{cases} 0 & b = 0 \\ -\frac{2}{3\pi} & b=r=\tfrac{1}{2} \\ \frac{4r}{3\pi}\mathrm{cel}(k_c,1,-1,k_c^2) & b=r < 1/2\\ -\frac{2}{3\pi} \mathrm{cel}(k_c,1,1,2k_c^2) & b=r > 1/2\\ -\frac{8r}{3\pi}\sqrt{r(1-r)} & b+r = 1 \\ \end{cases}\\ \frac{d \Lambda}{d r} &= \begin{cases} 2 r\sqrt{1-r^2} & b=0\\ \frac{2}{\pi} & b = r = \tfrac{1}{2}\\ \frac{4r}{\pi} E(4r^2) & b=r < 1/2\\ \frac{2}{\pi} \mathrm{cel}(k_c,1,1,0) & b=r > 1/2\\ \frac{8r}{\pi}\sqrt{r(1-r)} & b+r = 1\\ \end{cases}

rodluger commented 6 years ago

Nice!

ericagol commented 6 years ago

I'm puzzled by one thing: in the b=0 case from symmetry considerations, I would expect the derivative of \Lambda with respect to b to be zero, but I'm finding a Mathematica calculation Mathematica Notebook that it is non-zero. I also find from numerical differentiation that the derivative is zero, so I'm fairly sure that the Mathematica result is wrong.

rodluger commented 6 years ago

It has to be zero because you should get the same solution for b = eps as for b = -eps...

ericagol commented 6 years ago

Exactly - that’s why I’m puzzled by the Mathematica result, but it could be wrong due to not treating it as an absolute value.

ericagol commented 6 years ago

For the special cases b=r and b=1-r, it turns out that all I need to do is to add in (b-r) d\Lambda/db and (b+r-1) d\Lambda/db to the computation to give the correct automatic differentiated derivatives. For the b=r=1/2 case, I have to add (b-1/2)d\Lambda/db + (r-1/2)d\Lambda/dr. This gives the correct derivatives, which I show in test_s2_gradient.jl.

ericagol commented 6 years ago

These have now been added to Julia code.