How to find the Nth zero of a function

Question

It's almost the same question with "multiple root finding", but the order of roots is important in my case.

Basically, I'm trying to find $n$-th crossing point of the following graph:

Plot[{x BesselJ[1, x], 2 BesselJ[0, x]}, {x, 0, 15}, PlotRange -> All]

For example, I need to find out the 1005th (positive and the order count from the smallest) root of the equation below:

x BesselJ[1, x]-2 BesselJ[0, x]==0

Until now, FindRoot seems to be the fastest and clearest method, but I couldn't find a way to specify the $n$-th root.

FindRoot[x BesselJ[1, x] == A BesselJ[0, x], {x, 8, 12}, Method -> "Brent"]

score 4 · Accepted Answer · edited Aug 22 '17 at 09:09

x BesselJ[1, x] very soon dominates the behavior and you can use BesselJZero[n, k] to constrain FindRoot[]:

f[x_] = x BesselJ[1, x] - 2 BesselJ[0, x]
zero[n_] := BesselJZero[1, n - 1] // N
FindRoot[f[x] == 0, {x, zero[1005] - 1, zero[1005] + 1}]

Of course you need to be a little careful about the interval you search for the root in. A more thorough approach where we aren't guaranteed to have nicely spaced roots would be to step forward one root at a time.

Another approach if you know that roots are at least dx apart would be to observe the sign of your function as you step forward on the x-axis with dx/2 steps and simply count the number of sign changes, then again limit the search with two values - FindRoot[f[x], {x, x1, x2}] - to find the root in the identified interval.

score 3 · Answer 2 · answered Aug 21 '17 at 03:36

Another possibility when trying to find the $k$-th root of a function would be to adapt the NDSolve[]-based event detection strategy by Daniel given in this answer; that is, keep a running count as each zero is crossed, and then bail when you've hit the one you wanted. For instance:

f[x_] := x BesselJ[1, x] - 2 BesselJ[0, x]

k = 0;
j = 1005; (* index of zero *)
sol = Reap[NDSolve[{y'[x] == f'[x], y[0] == f[0], 
                    WhenEvent[y[x] == 0, k++; If[k == j, Sow[x]; "StopIntegration"], 
                              "DetectionMethod" -> "Interpolation",
                              "IntegrateEvent" -> True, "LocationMethod" -> "Brent"]},
                   {}, {x, 0, ∞},
                   AccuracyGoal -> ∞, WorkingPrecision -> 20];][[-1, 1, 1]]

which yields

3154.9449374317543982

You can polish further with FindRoot[] if wanted:

sol = \[FormalX] /. First[FindRoot[f[\[FormalX]], {\[FormalX], sol},
                                   WorkingPrecision -> 20]]
   3154.9449374319666351

If you wish for an entire range (say, the 26th to the 30th zeros), the event inside WhenEvent[] is easily modified:

k = 0;
{kmin, kmax} = {26, 30}
Reap[NDSolve[{y'[x] == D[f[x], x], y[0] == f[0], 
              WhenEvent[y[x] == 0, k++; 
                        If[kmin <= k <= kmax, Sow[x], If[30 < k, StopIntegration"]], 
                        "DetectionMethod" -> "Interpolation",
                        "IntegrateEvent" -> True, "LocationMethod" -> "Brent"]},
             {}, {x, 0, ∞}, AccuracyGoal -> ∞, WorkingPrecision -> 20];][[-1, 1]]
   {79.345692016011214192, 82.486505133453918880, 85.627375391257917317,
    88.768296740115704530, 91.909263963793310359}

and again one can post-process with FindRoot[] if need be.

How to find the Nth zero of a function

2 Answers2

Linked