# NDEigensystem convergence and comparison to DEigensystem

Consider the following eigenvalue differential equation

$$-u_n''(x)+x^2u_n(x)=E_nu_n(x),\qquad x\in(-\infty,\infty).$$

If you know a little bit of physics, you would know that this is the differential equation for a quantum harmonic oscillator, and that the energies are exactly given by

$$E_n=2n+1,\qquad n=0,1,2\cdots$$

We can also confirm this using Mathematica

 {eig, func} = DEigensystem[-u''[x] + x^2*u[x], u[x], {x, -Infinity, Infinity}, 5];
ListPlot[eig]


This is the expected result. But things go wrong quickly if we want to solve this numerically. In that case our code will look like

 {eig, func} = NDEigensystem[-u''[x] + x^2*u[x], u[x], {x, -1000, 1000}, 10];
ListPlot[eig]


We expect the numerical calculation to converge to the exact calculation when we increase the bounds of $$x$$ but that does not happen. As we increase the boundaries, the energy eigenvalues blow up, and they do not converge.

I had the same experience using NDEigensystem for other complicated scenarios (NDEigensystem eigenvalue convergence problem), and therefore I wanted to frame the problem in this simple case. Is there something fundamental, or am I completely butchering NDEigensystem. Any help regarding this question and the linked question would be appreciated.

You're setting so large a domain for approximation of $$\infty$$ that the default grid is too coarse to accurately describe the eigensystem. Try {x, -5, 5} instead:

{eig, func} = NDEigensystem[-u''[x] + x^2*u[x], u[x], {x, -5, 5}, 10];
eig
(* {1.00016, 3.00109, 5.00381, 7.0094, 9.01882, 11.0328, 13.0513, 15.0712,
17.0803, 19.046} *)


Or setting a denser grid:

{eig, func} =
NDEigensystem[-u''[x] + x^2*u[x], u[x], {x, -1000, 1000}, 10,
Method -> {"PDEDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> 0.1}}}];
eig
(* {1., 3., 5.00001, 7.00002, 9.00003, 11.0001, 13.0001, 15.0001, 17.0002,
19.0003} *)