astro-group-bristol / Gradus.jl

Extensible spacetime agnostic general relativistic ray-tracing (GRRT).
https://astro-group-bristol.github.io/Gradus.jl/dev/
GNU General Public License v3.0
18 stars 2 forks source link

Updated examples.md to use new syntax #71

Closed phajy closed 1 year ago

phajy commented 1 year ago

Updated the examples.md file to use the new syntax. A couple of issues:

Fig2

fjebaker commented 1 year ago

~Hmmm, this looks like an issue in the way the branch split is being performed... let me test a few things.~

Edit: nevermind, I was thinking of the wrong example. I can reproduce this.

fjebaker commented 1 year ago

I seem to get a full error when doing Fig. 1

ERROR: Transfer function integration failed (rₑ=4.0, θ=4.490256166494502).
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] cunningham_transfer_function(m::KerrMetric{Float64}, u::SVector{4, Float64}, d::GeometricThinDisc{Float64}, rₑ::Float64; max_time::Float64, diff_order::Int64, redshift_pf::PointFunction{typeof(Gradus._redshift_guard)}, offset_max::Float64, zero_atol::Float64, N::Int64, tracer_kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ Gradus ~/Developer/devenv/dev/Gradus/src/transfer-functions/cunningham-transfer-functions.jl:103
 [3] top-level scope
   @ ~/Developer/devenv/temp.jl:11

At some point I removed running the examples from the test suite, since they took too long on the CI machine, and subsequently I've tweaked all the defaults so there is a chance I broke something.

fjebaker commented 1 year ago

Found the problem:

https://github.com/astro-group-bristol/Gradus.jl/blob/ec031bb31a9458a97816859b5b52257cb5102e40/src/transfer-functions/precision-solvers.jl#L67

This line is supposed to calculate the radial coordinate on the disc, which the root finder then uses to hone in the on correct impact parameters. I toyed about with these lines in

which was some time after I'd written the example

Changing this line to just gp.u2[2] fixes the integration problem. Although sin(gp.u2[3]) should be 1 in the equitorial plane, there are some numerical ideosyncracies that the root finder is very good at uncovering that ruin this.

The projection was included so the method would work with thick discs too, but there are better ways that this can be implemented so it doesn't break the thin disc examples.

fjebaker commented 1 year ago

Changing this line doesn't seem to break any of the tests (which means I need to add a test case for it), so feel free to change it and add the commit to this PR 👍

phajy commented 1 year ago

I think I've added two changes to the PR. One trivial change to .gitignore and the other to remove sin(gp.u2[3]) as suggested (that's in the "Conflicting files" below).

However, I now get a Roots.ConvergenceFailed("Algorithm failed to converge") error when trying to reproduce the Cunningham transfer functions from Fig. 1 of the Bambi et al. (2017) paper (see code snippet below).

https://github.com/astro-group-bristol/Gradus.jl/blob/8b532bc4ffe35b4ef666dabda73804d564c534af/docs/src/examples/examples.md?plain=1#L350-L365

This might be because I'm "between versions" though. E.g., I don't think the new implementation of cunningham_transfer_function requires the outer radius to be passed, and gstar should be g\varstar.

fjebaker commented 1 year ago

@phajy conflicts resolved, running test suite now.

fjebaker commented 1 year ago

This might be because I'm "between versions" though. E.g., I don't think the new implementation of cunningham_transfer_function requires the outer radius to be passed, and gstar should be g\varstar.

This is most likely the case, as I cannot reproduce the error on main.

codecov-commenter commented 1 year ago

Codecov Report

Merging #71 (5fe05d0) into main (8b532bc) will not change coverage. The diff coverage is n/a.

:mega: This organization is not using Codecov’s GitHub App Integration. We recommend you install it so Codecov can continue to function properly for your repositories. Learn more

@@           Coverage Diff           @@
##             main      #71   +/-   ##
=======================================
  Coverage   71.73%   71.73%           
=======================================
  Files          44       44           
  Lines        1514     1514           
=======================================
  Hits         1086     1086           
  Misses        428      428           
Impacted Files Coverage Δ
src/transfer-functions/precision-solvers.jl 63.82% <ø> (ø)

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

phajy commented 1 year ago

I still seem to get an error using Gradus v0.2.3 (although there was a problem with one of the automated workflows). Here's the error message:

using Gradus
using StaticArrays
m = KerrMetric(M=1.0, a=0.998)
d = GeometricThinDisc(0.0, 100.0, π/2)
angle = 85
u = @SVector [0.0, 1000.0, deg2rad(angle), 0.0]
ctf = cunningham_transfer_function(m, u, d, 4.0)
ERROR: Roots.ConvergenceFailed("Algorithm failed to converge")

I also found that the angle = 74 case does work but has quite a large "step" in the path when g* gets close to 1.

I think I'm running the latest version, but user error is definitely a possibility!

fjebaker commented 1 year ago

tempfig

@phajy is this the step you are describing?

fjebaker commented 1 year ago

The case for angle = 85 as you have above reproduces the error your describe, and the reason for this is quite simple: the initial radius of photons traced by the root finder are missing the disc underneath the equitorial plane. Setting offset_max = 15.0 for example helps the root finder converge fine, but the transfer function has a visible error:

tempfig

This is known behaviour, but it's still a bug and there should be an issue raised for it.

I'll quickly write a 3D plotting script and post it here (pending to add to examples), that will hopefully let you debug the root finder problems a little better.

Thank you for these bug reports though -- they are extremely useful! 👍