Using FindRoot and having problems with $MaxExtraPrecision

Question

I'm trying to find the roots of the function radial, written below. And to make this, I'm using the function \[Omega]ini to define a point closer.

But using any value of prec, the FindRoot returns the error "N::meprec: Internal precision limit $MaxExtraPrecision = 2 prec".

I've tried to N[...,{∞,prec}] as suggested in Mathematica, but it doesn't solve the problem.

An important point is that 0 < Q < 1.

Someone could take a look, please?

radial[{Q_, q_, l_}, ω_, prec_] := 
 Block[{$MinPrecision = prec, $MaxExtraPrecision = 2 prec, λ, 
   rh, rc, ζ, η, ξ, α, β, γ, \
δ, ϵ, z},
  
  rh = 1 + Sqrt[1 - Q^2];
  rc = 1 - Sqrt[1 - Q^2];
  λ = l (l + 1);
  
  ζ = I ω;
  η = I rc (ω rc - q Q)/(rh - rc);
  ξ = -I rh (ω rh - q Q)/(rh - rc);
  
  α = 2 ζ (rh - rc);
  β = 1 + 2 ξ;
  γ = 1 + 2 η;
  δ = 
   2 (rc - rh) ((-2 + rc + 
         rh) ζ^2 - ζ (1 + η + ξ) + (q Q - 
         2 ω) ω);
  ϵ = 
   q^2 Q^2 + (-2 Q^2 - rc^2 + (-2 + rh)^2 + 
       2 rc rh) ζ^2 + η + η^2 - λ + ξ + 
    2 η ξ + ξ^2 - ζ (rc + 2 rc ξ - 
       4 (1 + η + ξ) + rh (3 + 4 η + 2 ξ)) - 
    2 q Q (2 + rc) ω + 4 ω^2 - 2 Q^2 ω^2 + 
    4 rc ω^2;
  
  z = 100 E^(-I (π/2 + Arg[-ω]));
  
  E^((-rc + rh) z ζ) (1 - 
     z)^ξ z^η  HeunC[-δ - ϵ, -δ, \
β, γ, -α, 1 - z]]

ωini[{Q_, q_, l_}, n_] := Block[{s, λ, rh, rc},
  
  s = 0;
  
  λ = (l - s) (l + s + 1);
  rh = 1 + Sqrt[1 - Q^2];
  rc = 1 - Sqrt[1 - Q^2];
  
  (q Q)/rh - 
   I ((1 + 2 n) (rh - rc))/(
    4 rh^2) + ((rh - rc) (2 rh (1 + 2 s + 2 λ) + 
      rc (4 s^2 - 1)))/(16 rh^3 q Q) - 
   I ((1 + 2 n) (3 rc - rh) (4 s^2 - 1) (rh - rc)^2)/(64 rh^4 q^2 Q^2)
    ]

minimum[{Q_, q_, l_},  n_, prec_] := 
 minimum[{Q, q, l}, n, prec] = 
  N[FindRoot[
    radial[{Q, q, l}, ω, prec + 10], {ω, 
     N[ωini[{Q, q, l}, n], prec]}, 
    WorkingPrecision -> prec], {∞, prec}]

Answer 1

The problem is with epsilon in radial[]: It, or it plus delta, equals zero but does not auto-simplify to 0. When it is numericized at an arbitrary precision, the quantity being zero, it is impossible to achieve any (relative) precision.

You can simplify epsilon, e.g. HeunC[-δ - Simplify[ϵ],...] or HeunC[Simplify[-δ - ϵ],...]. FullSimplify[] seems not to be needed.

Or you can numericize Q, q, l, etc. before computing epsilon and replace approximate zeros (underflows) by 0:

Block[{QQ, qq, ll},
 {QQ, qq, ll} = SetPrecision[{Q, q, l}, prec];
 ... (* change `Q` to `QQ` etc *)
 E^((-rc + rh) z ζ) (1 - z)^ξ z^η *
   Hold[HeunC][-δ - ϵ, -δ, β, γ, -α, 1 - z] /.
    x_ /; x == 0 -> 0 // ReleaseHold
]

It takes a long time for minimum[] to run, and I don't know what a good test case is. I get other errors (e.g. FindRoot::lstol) sometimes, and it could be because the test case I try has no solution. So I will leave it at that: You get a N::meprec warning because of numericizing zero to a non-MachinePrecision precision, like this example:

N[Cos[Pi/8] - Root[1 - 8 #1^2 + 8 #1^4 &, 4], 5]