BlackHolePerturbationToolkit / SpinWeightedSpheroidalHarmonics

Tools for computing spin-weighted spheroidal harmonics and their associated eigenvalues
https://bhptoolkit.org/SpinWeightedSpheroidalHarmonics/
MIT License
12 stars 6 forks source link

Branch cuts #5

Open marccasals opened 7 years ago

marccasals commented 7 years ago

The spheroidal eigenvalues possess branch cuts in the complex-frequency plane (see Studies in Applied Mathematics 113, 271, 2004). The package's "SpinWeightedSpheroidalEigenvalue" seems to not only choose the branch cuts going in a direction different from the direction choice in Mathematica’s in-built “SpheroidalEigenvalue” but also it seems to jump to a branch different from Mathematica’s. See, eg, the plot attached for the eigenvalue as a function of arg(\gamma) for spin=0, \ell=1, m=0, |\gamma|=10 (blue curve is the package's "SpinWeightedSpheroidalEigenvalue" and red curve is Mathematica’s “SpheroidalEigenvalue”). It might not necessarily mean that the choices of direction and branch are wrong, but it'd be good to have control over the choices and perhaps convenient to make the same choices as Mathematica’s “SpheroidalEigenvalue”. Plot_lambda_branch.pdf

nielsw2 commented 7 years ago

Thanks for the report, Marc. The package has had less development for complex gamma but we're keen to implement improvements.

We should (i) check that the current version gives a correct answer (but for a different branch cut than Mathematica uses) and (ii) think about how we can implement control over which branch cut is chosen

KwintenF commented 5 years ago

->(i) for the examples I looked at, there was always a nearby choice of l which did match suggesting that the current version does give the right answer (iii) Could this choice affect statements about QNMs? (an initial survey (rather low lying modes) suggests not)

barrywardell commented 5 years ago

The following is an example showing this problem (taken from 0e61cfb):

Module[{l = 2, m = 2, c = 2.28831099800269655020201753359287977218627929687549.74505997240401 + 5.351950938495406262518372386693954467773437550.11405704809964 I}, {N[SpheroidalEigenvalue[l, m, c], 10], N[SpinWeightedSpheroidalEigenvalue[0, l, m, I*c] + 2*m*I*c, 10]}]

barrywardell commented 4 years ago

Marc's originally reported problem was fixed in 250b0531. However, the following plots show that there are still cases where the code is suddenly jumping to a different branch: branches-2.pdf, branches-1.pdf. This does not happen with Mathematica's built-in SpheroidalEigenvalue, so we may be able to learn how to fix it by studying that (via Falloon's code). These plots were generated using the following code:

Module[{l = 2, m = 2,  c = 2.288310998002696550202017533592879772186279296875`49.74505997240401 + 5.3519509384954062625183723866939544677734375`50.11405704809964 I}, Plot[{Re[N[SpinWeightedSpheroidalEigenvalue[0, l, m, I*c Exp[I \[Theta]]] + 2*m*I*c Exp[I \[Theta]], 10]], Re[N[SpinWeightedSpheroidalEigenvalue[0, l + 1, m, I*c Exp[I \[Theta]]] + 2*m*I*c Exp[I \[Theta]], 10]], Re[N[SpheroidalEigenvalue[l, m, c Exp[I \[Theta]]], 10]]}, {\[Theta], -0.2, 0.2}, PlotTheme -> "Detailed", PlotStyle -> {Automatic, Automatic, Directive[Black, Dashed]}]]

and

Module[{l = 2, m = 2,  c = 2.288310998002696550202017533592879772186279296875`49.74505997240401 + 5.3519509384954062625183723866939544677734375`50.11405704809964 I}, Plot[{Im[N[SpinWeightedSpheroidalEigenvalue[0, l, m, I*c Exp[I \[Theta]]] + 2*m*I*c Exp[I \[Theta]], 10]], Im[N[SpinWeightedSpheroidalEigenvalue[0, l + 1, m, I*c Exp[I \[Theta]]] + 2*m*I*c Exp[I \[Theta]], 10]], Im[N[SpheroidalEigenvalue[l, m, c Exp[I \[Theta]]], 10]]}, {\[Theta], -0.2, 0.2}, PlotTheme -> "Detailed", PlotStyle -> {Automatic, Automatic, Directive[Black, Dashed]}]]