# How to find the roots of an equation which is almost singular everywhere

I'm trying to find the roots of the equation as a variable of k:

a = 0.81; b = 1.94;
LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, 1] == 0


The problem is that this function is almost divergent everywhere except only a few points, and FindRoot and Reduce seems to fail on this

FindRoot[LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, 1] == 0, {k, 0}]


FindRoot::nlnum: The function value {ComplexInfinity} is not a list of numbers with dimensions {1} at {k} = {1.49012*10^-8}. >>

Reduce[LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, 1] == 0, k]


Reduce::inex: Reduce was unable to solve the system with inexact coefficients or the system obtained by direct rationalization of inexact numbers present in the system. Since many of the methods used by Reduce require exact input, providing Reduce with an exact version of the system may help. >>

With[{a = 81/100, b = 194/100},
Reduce[LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, 1] == 0, k]]


Reduce::nsmet: This system cannot be solved with the methods available to Reduce. >>

How would I solve this equation? The roots are k = n*a, where n is integer.

Edit

The motivation of this problem is to determine the bound states of a one-dimensional potential well

$$-\frac{\hbar^2}{2m}\frac{d^2\psi}{dx^2}+V(x)\psi=-\frac{\hbar^2k^2}{2m}\psi$$

where

$$V(x)=-\frac{\hbar^2b}{m}sech^2(a x)$$

,$k$ is a parameter.

This equation has solution

ψ[x_]:=A*LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, Tanh[a x]] +
B*LegendreQ[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, Tanh[a x]]

V[x_] := -((\[HBar]^2 b)/m) Sech[a x]^2
-(\[HBar]^2/(2 m)) ψ''[x] + V[x] ψ[x] == -((\[HBar]^2 k^2)/(2 m)) ψ[x] // FullSimplify
(*True*)


Quantum mechanics requires the wave function is finite, which gives a constrain on parameter k and thus the eigen-energies($\sim k^2$) are some discretized values. And I'm trying to find these eigen-energies. Mathematically, this is equivalent to finding the $k's$ that keep associate Legendre polynomials finite in region 0 <= x <= 1

A*LegendreP[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, x]+B*LegendreQ[(-a + Sqrt[8 b + a^2])/(2 a),  k/a, x]


.

Edit2

For the wave function must be normalizable and square-integratable, here is the explanation in the Griffiths's Introduction to Quantum Mechanics:

For some solutions to the Schrodinger equation the integral is infinite; in that case no multiplicative factor is going to make it 1. The same goes for the trivial solution $\Psi = 0$. Such non-normalizable solutions cannot represent particles, and must be rejected. Physically realizable states correspond to the square-integrable solutions to Schrodinger's equation.

For the case of free particle, he also explained why it's not normalizable:

In the case of the free particle, then, the separable solutions do not represent physically realizable states. A free particle cannot exist in a stationary state; or, to put it another way, there is no such thing as a free particle with a definite energy.

-
the Schroedinger equation with this potential has an exact solution. You may find it in the Landau and Lifshitz book, V. 3 Quantum Mechanics, or in my old paper Ferroelectrics V. 98, p. 277 (1989), where it is Eq. (6). Have fun. – Alexei Boulbitch Oct 1 '13 at 15:17
@AlexeiBoulbitch Thanks a lot! I never knew that it was possible to solve the problem exactly. That's really helpful. – xslittlegrass Oct 1 '13 at 16:41
@AlexeiBoulbitch: OP already posted the exact solution. The line of code where OP replaces the wavefunction into Schrodinger's equation proves that it is the exact solution. OP needs help finding the subspace of the space of solutions that also satisfy the boundary conditions for this problem. – Hector Oct 1 '13 at 17:05
@Hector The energies for bound state are $E_n=-\frac{\hbar^2a^2}{8m}\left[\sqrt{1+\frac{8b}{a^2}}-(2n+1)\right]^2$ and n is an non-negative integer and $n<\frac{1}{2}\left(\sqrt{1+\frac{8b}{a^2}}-1\right)$. For this a = 0.81; b = 1.94; there are two bound states at -35.0908eV and -8.62015eV, consistent with the numerical solution at here. It seems that somehow the Legendre P and Q cancel the divergence of each other. – xslittlegrass Oct 1 '13 at 17:29
@xslittlegrass It is always advisable to know, what is written in the Landau&Lifshitz book. – Alexei Boulbitch Oct 2 '13 at 9:09

Quantum mechanics requires the wave function is finite, which gives a constrain on parameter k and thus the eigen-energies(∼k2) are some discretized values.

Not true! The wave function does not need to be finite. It does not even have to be an $L^2$ function. If your teacher insists on it, bring the free-particle wave function as an example.

For the sake of a short argument, let us assume that there is no such a case where the divergence of LegendreP exactly cancels the divergence of LegendreQ. Also, let us assume that $k>0$.

With the replacements {(-a + Sqrt[a^2 + 8 b])/(2 a) -> Δ, k/a -> κ, a x -> 1/ξ} you get a simpler function to analyze:

LegendreP[Δ, κ, Tanh[1/ξ]]


The reasoning behind the following calculation is quite obvious: find the behavior of ψ at infinity.

( Series[LegendreP[Δ, κ, Tanh[1/ξ]], {ξ, 0, 0}] // FullSimplify) /.
(1 - Tanh[1/ξ])^(-κ/2) (1 + Tanh[1/ξ])^(κ/2) -> E^(κ/ξ) /.
ξ -> 1/(a x)


Now, $k>0$ implies $κ>0$ which in turns causes a divergent exponential. That is bad unless the other factors go to zero faster. Thus, either $κ=1$ (from the numerator) or $κ=2,3,4,…$ (from the denominator).

Similarly, the behavior at minus infinity

gives you another constraint (because $\lim_{x\to-∞} sech(x)=0$ and $-κ$ is a negative integer): $Δ=2n$ with $n$ an integer (positive, negative, or zero).

Now, from the physical point of view, this is a rather interesting constraint. If the calculations above are correct, it says that bound states do not happen unless b = 1/2 a^2 n (1 + n).

Lather, rinse, and repeat for LegendreQ. Check that linear combinations of LegendreP and LegendreQ do not cancel each other for special cases.

Update

(1) The procedure above shows that the behavior at +∞ (computed using Series) imply $κ$ must be an integer. Since $κ=k/a$, the original question has been answered.

(2) From the analysis of the behavior at -∞, a relationship between $a$ and $b$ was found. That relationship is necessary for the solutions posted by OP. The following plot shows how using the values provided by the OP, it is impossible to get a linear combination of LegendrePand LegendreQ that does not diverge at -∞. The numerical calculation quoted by OP imply solutions more general than those posted by OP.

On the other hand, modifying $b$ slightly to comply with $Δ=2$, we obtain the following plot where LegendreP is a function that converges at both +∞ and -∞.

(3) Instead of changing the variable from $x$ to $ξ$, you could have used

Series[LegendreP[Δ, κ, Tanh[x]], {x, ∞, 0}]

-
By the way, use Series[LegendreP[Δ, -κ, Tanh[1/ξ]], {ξ, 0, 0}] to get the behavior for negative values of $κ$. Do not replace $κ$ by $-κ$ in the calculations above. Also, the values of $a$ and $b$ that you provided yield $Δ=1.98$. Rounding errors? – Hector Oct 1 '13 at 11:57
Thanks a lot for the answer. But according to Griffiths's book, the wave-function does require finite and square-integrable, see by edit. – xslittlegrass Oct 1 '13 at 15:38
Inextensible strings, massless springs, incompressible water, definite-energy states, and others do not really exist in nature. However, we do not reject them from our mathematical toolbox because they are abstractions that allow us to solve problems easily. You might want to follow the links from here. – Hector Oct 1 '13 at 16:51
Thanks a lot for the analysis, but I believe the eigen energies gave in Landau's book are correct. For example, you can verify that there are two bound states solution satisfy the Schrodinger equation, without requiring relationship between a and b :Table[ \[Psi][x_] := Hypergeometric2F1[Sqrt[1 + (8 b)/a^2] - n, -n, 1/2 + 1/2 Sqrt[1 + (8 b)/a^2] - n, 1/2 (1 - Tanh[a x])] (Sech[a x]^2)^( 1/4 (-1 + Sqrt[1 + (8 b)/a^2] - 2 n)); k = a/2 (Sqrt[1 + (8 b)/a^2] - (2 n + 1)); -(1/2) \[Psi]''[x] + V[x] \[Psi][x] == -(k^2/2) \[Psi][x] // FullSimplify , {n, 0, 1}] – xslittlegrass Oct 2 '13 at 16:03
And the eigen energies and wavefunctions are all consistant with the numerical results I referred. I think maybe we can write Hypergometric functions in terms of Legendre polynomials, so that the divergence of Legendre polynomials would cancel each other. – xslittlegrass Oct 2 '13 at 16:07