Consider the following one dimensional Schrodinger equation within the complex plane of $z$ $$ −ψ''(z) − (iz)^ N ψ(z) = Eψ(z). $$ where $N$ can be any real number, the boundary condition is $ψ(z) → 0$ as $|z| → ∞$ .
For simplicity, let's replace $z$ to be $x$ and $x$ is real. I set $N=2$, so that I have a harmonic oscillator:
$$
−ψ''(x) + x^ 2 ψ(x) = Eψ(x).
$$
I solve the last equation by using shooting method and implicit runge-kutta integrator. For eigenvalue $E=1$, the corresponding eigenfunction looks like:
where I set $x=4$ as $+∞$ and $x=-4$ as $-∞$ . For eigenvalue $E=5$, the corresponding eigenfunction looks like:
where the eigenfunction roughly still exponentially decays at $x=\pm 4$. However, when eigenvalue $E=17$, the corresponding eigenfunction won't exponentially decays at $x=\pm 4$ as shown on the following Fig.

In fact, the eigenfunction decays at $x=\pm 6$ .

My Question
The numerical infinity should be large (like $x=6$) for large eigenvalue $E$, and small (like $x=4$) for small eigenvalue. If independent from the eigenvalue, I uniformly set numerical infinity as a very large number, it would waste computational resource. So, is there any formula to determine the exact numerical infinity based on different inputted eigenvalue please? Thank you!
It seems OP is essentially asking how to determine which regions the wavefunction $\psi(x)$ is oscillatory (decaying)? In quantum mechanics, these regions are called classically allowed (forbidden) regions, respectively. The turning points are determined by the condition that the total energy $$\tag{1} V(x)~=~E$$ matches the potential energy $$\tag{2}V(x)~:=~x^2.$$ So the turning points can be found as $$\tag{3}|x|~=~\sqrt{E}.$$ If the seed values for $x$ should be in the classically forbidden regions, then one should pick $|x|$ a bit bigger than (3).