Is the Wave function a "Smooth" function of the Potential? Consider the Schroedinger equation in a spherically symmetric system. In the unit system under which energy is measured in Hartree and length in Bohr radius $a_0$, the schroedinger equation can be written as
$$
-{d^2 R_{l,\epsilon}(r) \over dr^2} + \left( V_d(r) + {l(l+1)\over 2 r^2} \right)R_{l, \epsilon}(r) = \epsilon R_{l, \epsilon}(r)
$$ 
. In the above equation, $V_d(r)$ is an attractive potential that's dependent on $d$. Specifically, we will consider the following potential. $V_{0}(r) = -Z(r)/r$, where $\lim_{r \rightarrow 0}Z(r) = Z$ and $\lim_{r \rightarrow \infty} Z(r) = 1$. $V_d(r) = V_0(r) + {1 \over r} - {e^{-r/d} \over r}$ ($d>0$). Consider the wave function $R_{l, \epsilon}^{reg}(r)$ that is regular at the origin and which is normalised to unit magnitude at infinity. In other words, 
$$
R_{l, \epsilon}^{reg}(r) \rightarrow \sin(k r + \phi)
$$,
as $r \rightarrow \infty$. $k = \sqrt{2 \epsilon}$. My question is, for at least $\epsilon>0$, is $R_{l, \epsilon}^{reg}(r)$ a "smooth" function of $d$?
 A: The radial Hydrogen equation is
$$
      \frac{1}{r^2}\frac{d}{dr}r^2\frac{dR}{dr} + \left[\frac{2\mu}{\hbar^2}\{E-V(r)\}-\frac{l(l+1)}{r^2}\right]R=0,
$$
where $V$ is the Coulomb potential
$$
          V(r) = -\frac{Ze^2}{4\pi\epsilon_{0}r}.
$$
This equation is considered on $L^2_{r^2}(0,\infty)$ because of the $r^2$ weight arising from the spherical coordinate volume element. Because this is a Math site, I'll get rid of the Physical constants:
$$
               \frac{1}{r^2}\frac{d}{dr}r^2\frac{dR}{dr}+\left[\lambda+\frac{a}{r}-\frac{l(l+1)}{r^2}\right]R=0,
$$
where $\lambda$ is an eigenvalue parameter, $a$ is a fixed constant, and $l$ is an integer coming from the eigenvalues of the (associated) Legendre equation. Of particular importance is the lowest energy solution with no angular momentum corresponding to $l=0$; the solution has the form $R(r)=Ce^{-\mu r}$ where $\mu$ is determined from the equation
$$
     \left[\frac{d^2}{dr^2}+\frac{2}{r}\frac{d}{dr}+\lambda+\frac{a}{r}\right]R=0 \\
         (\mu^2+\lambda)+\frac{1}{r}(-2\mu+a)=0,
$$
which gives a solution if $\mu=a/2$ and $\lambda=-\mu^2=-a^{2}/4$. The smallest eigenvalue is then $\lambda = -a^{2}/4$, with solution $R_1(r)=Ce^{-ar/2}$. The second solution of the equation with the same eigenvalue $\lambda$ is not square integrable with respect to weight function $r^2$. In fact, the second solution has an exponential growth, and is given by variation of parameters to be
$$
 R_2(r)=-\frac{e^{ar/2}}{r}+a\ln(r)e^{ar/2}-a^2e^{-ar/2}\int_{0}^{r}\ln(s)e^{as}ds
$$
Your equation: The functions $u_1(r)=rR_1(r)$, $u_2(r)=rR_2(r)$ are then solutions of
$$
          -\frac{d^2u}{dr^2}+\left(-\lambda-\frac{a}{r}+\frac{l(l+1)}{r^2}\right)u=0,\;\;\lambda = -a^2/4,\; l=0.
$$
This is very close to your equation, but there is no combination of the solutions $u_1$, $u_2$ that has the asymptotic behavior you are expecting, which makes me doubt your expected asymptotic behavior.
However, if $\lambda$ is positive instead, then the perturbation obtained by changing $d$ is integrable on $(0,\infty)$, which would make me think that you will have a nice regularity in $d$. $L^1$ perturbations are very tame.
