The Setup: We will take as our space of states H=L2([0,1]) and hamiltonian
H=−d2dx2+V(x).
for some real function V(x) defined on [0,1].
Now fix some large positive integer N. Let ϵ=1/N. We will consider the subspace HN of H spanned by those functions that are constant on the subintervals iϵ,(i+1)ϵ). Such a function is defined (a.e.) by the N values it takes on these intervals, so we may identify HN≅CN as vector spaces. Let us denote elements of HN by ψi for i=0,…,N−1. Thinking of these as the values of a function ψ(x) at x=iϵ, we see that the inner product on HN is given by
(ϕ,ψ)=∫10ˉϕ(x)ψ(x)dx=N−1∑i=0ˉϕiψ(i)
So we see that with these identifications, HN is just CN with the usual hermitian inner product.
Now, we can approximate d/dx with the forward and backward difference operators
(D+ψ)i=ψi+1−ψiϵ(D−ψ)i=ψi−ψi−1ϵ
Note: throughout I will assume periodic boundary conditions to make life easy. In this case we have ψi+N=ψi.
Now consider
(D+ϕ,ψ)=ϵ−1∑i(ˉϕi+1ψi−ˉϕiψi)=ϵ−1∑i(ˉϕiψi−1−ˉϕiψi)=−(ϕ,D−ψ).
Thus we have D∗+=−D−. Similarly, D∗−=−D+. Then the operator D=(D++D−)/2 approximates d/dx and satisfies D∗=−D (as it should!), and the discrete Laplacian D2 is self-adjoint.
Finally, we can form the discrete hamiltonian HN by taking HN=−D2+ˆV, where ˆV is the operator ψi↦V(xi)ψi, where xi=iϵ.
Note: typically one further imposes Dirichlet or Neumann boundary conditions. This corresponds to projecting to a smaller subspace of HN.
I wrote some Sage code to test this. With V(x)=0, and N=500, here is one of the lowest-energy states:
Rather encouraging.