# General Solution to Linear Schrodinger equation

I am trying to find a solution to $$\displaystyle \left[-\frac{1}{2}\nabla^2 - \frac{2}{r} + C(r)\right]\phi(r) = E\phi(r)$$ where $C(r)$ is a known function of r. I am just looking for some help on the strategy to solve explicitly for $\phi(r)$. I have read through the solution for hydrogen using separation of variables but

1) I don't really understand it, they make a lot of handwavy arguments about terms disappearing :)

2) How does one handle the extra $C(r)$?

FYI: this is for the Hartree Fock method on Helium so the actual equation is the one posted in the answer here: Hartree Fock equations

P.S. I posted this in math.stackechange first should I close that one?

• Do you have mathematica? – Sponge Bob Nov 28 '15 at 3:02
• Haha no I don't :( but can Mathematica handle it? Also I am trying to implement my own HF algorithm in Python so if I could come up with an analytical solution that would be awesome. – user2879934 Nov 28 '15 at 3:04
• There's no general algorithm that can provide solutions for a generic $C(r)$. If only it were that simple! Wolfram Alpha's 'DSolve' can provide solace for some simple $C(r)$ but will quickly run out of computing steam for more complicated $C(r)$. – Gert Nov 28 '15 at 3:07
• @user2879934: HF-based methods provide numerical solutions, not analytical ones. – Gert Nov 28 '15 at 3:10
• @Gert well kind of this is just one step in HF for helium...we have already turned it into a one body problem that should be solvable by integrating out the effect of the second electron and turned it into an effective field. The step in HF that I am stuck on is I guessed my orbitals as Guassian and I calculated C(r) but how the heck to I generate a new psi using the Hamiltonian? Sorry for typos on phone – user2879934 Nov 28 '15 at 3:14

To do Hartree-Fock to Helium atom, you just need to calculate one orbital, which for helium is spherically symmetric. The Hartree-Fock 'integro-differential' equation for spherically symmetric atom with one eigenstate can be written as

$$\left(-\frac{1}{2}\nabla^2 - \frac{2}{r} + V_{Hx}(r) \right) u(r) = \epsilon u(r),$$

where $u(r) = \psi(r) r$ and

$$V_{Hx}(r) = \frac{1}{2} \int d{\bf r}' \frac{n({\bf r'})}{|{\bf r-r}'|}$$ The factor $\frac{1}{2}$ comes from exchange term cancelling half of the Hartree-potential.

I would recommend shooting method. (Alternatives: gaussian basis set, finite element, finite difference methods).

Shooting method is simple. Guess an eigenvalue and integrate the solution of radial Schrödinger equation from r=0 to r=r_max. Require that boundary condition is fulfilled at r_max (zero). If not, change the eigenvalue accordingly.

https://en.wikipedia.org/wiki/Shooting_method

Solving the radial Hartree-Fock potential equation is also super simple.

$$\int dr' \frac{n(r')}{|r-r'|} = \frac{1}{r} \int_0^r dr' 4\pi n(r) r^2 + \int_r^\infty 4\pi n(r) r$$

The latter part of my answer here explains how this equation works. Dipole in a spherical cavity in an infinite dielectric

Beware, that the Fock-operator is only multiplicative for the Helium 1s orbitals. Rest of the eigenvalue spectrum of the given equation does not match to Fock-operator spectrum.

edit: I just read a comment of yours, which stated that you are using gaussian orbitals. You can solve it easily. You just need to calculate some matrix elements. Given a set of gaussian orbitals $\phi_i$, you calculate two matrices

$$H_{ij} = < \phi_i | H | \phi_j>$$

and

$$S_{ij} = < \phi_i | \phi_j >$$.

Then you solve this generalized eigenvalue equation

$$Hc=\epsilon Sc$$, and you get the new coefficients for the wave function as the eigenvector. Your wave function is now

$$\psi(r) = \sum_n c_n \phi_n(r)$$

The method is called subspace diagonalization.

To solve the generalized eigenvalue equation with octave or Matlab,

octave:1> H = [ -2 0.1 ; 0.1 -3 ];
octave:2> S = [ 1 0.1 ; 0.1 1 ];
octave:5> [psi, lambda] = eig(H,S)
psi =

-0.35088  -0.94180
0.97217  -0.25494

lambda =

Diagonal Matrix

-3.1498        0
0  -1.9209


Here you obtain the eigenvalues (-3.14 is the lowest), and the corresponding eigenvector [-0.35, 0.97]. Now, all that there is left, is to fill matrices H and S with the real matrix elements.

To be more spesific, following integrals must be evaluated:

$$H_{ij} = \int dr \phi_i(r) (-1/2 \nabla^2 + V(r) ) \phi_j(r) \\ S_{ij} = \int dr \phi_i(r) \phi_j(r)$$

• In $$Hc=\epsilon Sc$$ are you solving for $\epsilon$ to update the coefficient? Also how does one solve that? Thanks so much Mikael! – user2879934 Nov 29 '15 at 2:02
• I updated the answer accordingly. One doesn't update epsilon, but obtains directly the correct value. This is called direct diagonalization. For large problems, one uses iterative diagonalization where $\epsilon$ is estimated and converged iteratively. – Mikael Kuisma Nov 29 '15 at 2:09
• oh man! I didn't realize math packages had the ability to solve eigien-value equations. That is awesome thank you! – user2879934 Nov 29 '15 at 2:15
• Actually one last thing...what "does fill matrices H and S with the real matrix elements" mean? Sorry I am just a computer programmer :) – user2879934 Nov 29 '15 at 2:52
• Haven't we already filled them by evaluating the two given integrals? – user2879934 Nov 29 '15 at 14:32

I feel like this question is a little too open ended and I'm sorry for that...but I did find a WONDERFUL HF walk through http://www.phys.sinica.edu.tw/TIGP-NANO/Course/2011_Spring/classnotes/CMS_20110511.pdf that clears everything up.