Roots of a function for eigensystem

Question

I want to find the roots for $\kappa$ for the equation

$$\sqrt{\alpha - 1} \cos{\left (\frac{\sqrt{2} \sqrt{\alpha - 1}}{2 \sqrt{\epsilon}} \right )} \cosh{\left (\frac{\sqrt{2} \sqrt{\alpha + 1}}{2 \sqrt{\epsilon}} \right )} - \sqrt{\alpha - 1} \\ -\frac{1}{\sqrt{\alpha + 1}} \sin{\left (\frac{\sqrt{2} \sqrt{\alpha - 1}}{2 \sqrt{\epsilon}} \right )} \sinh{\left (\frac{\sqrt{2} \sqrt{\alpha + 1}}{2 \sqrt{\epsilon}} \right )} = 0 \enspace ,$$ where $\alpha=\sqrt{1 + 4\epsilon\kappa^2}$. This equation has infinite roots, but I am interested in the first $N$ of them.

One option to solve this problem is to use Newton method, the problem is to define the initial points, since the values of the function can be quite high, as can be seen below

This problem comes from finding the eigenvalues of

$$\left[\frac{d^2 u}{ds^2} - \epsilon \frac{d^4u}{ds^4} + \kappa^2\right]u = 0$$

then I can obtain an approximation using perturbation methods, i.e., the eigenvalues are approximated by

$$\kappa_n^2 = n^2\pi^2 + \epsilon n^4\pi^4$$

for small $\epsilon$. Then, for small values of $\epsilon$, I can use the values $\kappa_n^2$ as initial guesses to the Newton algorithm. But when $\epsilon$ increases these initial guesses fail.

Since I know the original differential equation, I can use FEM or FDM to find the eigenvalues, but I am interested in other methods. Below, you can see the comparison using FDM (1001 points), perturbation solution and Newton method (using the perturbation solution as initial guess). All the curves are for $\epsilon=10^{-3}$, but I couldn't make those guesses to work for greater values of $\epsilon$.

Question: Is there any other method to solve this problem? Maybe some kind of transformation that can be applied to the equation?

Answer 1

I suspect the main problem is the magnitude of the values. If you divide through by $\cosh$ to make all the numbers smaller, then ApproxFun doesn't seem to have a problem finding all the roots.

function lhs(ε, α)
    A, B = sqrt((α-1)/(2ε)), sqrt((α+1)/(2ε))
    return sqrt(α*α-1)*(cos(A) - 1/cosh(B)) - sin(A)*tanh(B)
end

u = Fun(x -> lhs(1e-3, x), [1.0, 10.0], 1024)
r = ApproxFun.roots(u)
r2 = [fzero(x -> lhs(1e-3, x), r[i]*(1+1e-6*rand())) for i in 1:length(r)]

This finds all the roots in $\alpha\in[1,10]$, for $\epsilon=10^{-3}$, $\kappa^2\lesssim 25000$, using an 1024-term Chebyshev series, taking a very small amount of time. This seems to be accurate to about $10^{-8}$ and the roots can be polished off with Newton's method or similar.