Bootstrapping quantum mechanical systems
by D. Berenstein, G. Hulsey
Department of Physics, UC Santa Barbara
September 2021
by D. Berenstein, G. Hulsey
Department of Physics, UC Santa Barbara
September 2021
Introduction
Introduction
Efficient numerical solution of quantum mechanical problems has been a focus of applied computing in physics for decades. Depending on the system in question, certain solvers are better-suited than others. However, newer models, like matrix quantum mechanics, do not admit the same kinds of numerical solution as more traditional, one-particle problems. In light of this, new numerical techniques are needed which apply generally to quantum mechanical systems.
In this interactive notebook, we describe the novel ‘bootstrap’ approach to solving one-dimensional quantum systems. The bootstrap approach has a varied history of other use in physics: conformal field theory, lattice gauge theory, and, more recently, to numerically solve matrix models in the double-scaling limit. We focus on the quantum mechanical bootstrap, and its simplest instance: the Coulomb model of hydrogen. This system is exactly solved, so we can compare the performance of the bootstrap to the known answers.
The result is an algorithm which is exponentially convergent and easily generalizable to new problems. The notebook will allow you to run a bootstrap for the hydrogen model, obtain a numerical solution of the spectrum, and observe the convergence. This is based off our work which appeared in a recent pair of preprints available on the arXiv: https://arxiv.org/abs/2108.08757 and https://arxiv.org/abs/2109.06251.
Our approach follows closely the procedure introduced by Han, Hartnoll and Kruthoff. They applied the method to quantum mechanics and extended it to the matrix quantum mechanics in the large limit. Their work was preceded by a paper of Lin in which the bootstrap method for matrix models in the double scaling limit was developed.
The application to single-particle quantum mechanics has appeared in numerous recent papers. These include the work of Y. Aikawa et al on the harmonic case and the bootstrap on the circle with theta term, J. Bhattacharya et al on hydrogen and the double well, and a recent preprint (contemporary with our second) on the appearance of Bloch band structure for bootstrapping periodic potentials.
To view the in this notebook, make sure to enable Evaluation > Dynamic Updating Enabled.
In this interactive notebook, we describe the novel ‘bootstrap’ approach to solving one-dimensional quantum systems. The bootstrap approach has a varied history of other use in physics: conformal field theory, lattice gauge theory, and, more recently, to numerically solve matrix models in the double-scaling limit
1
The result is an algorithm which is exponentially convergent and easily generalizable to new problems. The notebook will allow you to run a bootstrap for the hydrogen model, obtain a numerical solution of the spectrum, and observe the convergence. This is based off our work which appeared in a recent pair of preprints available on the arXiv: https://arxiv.org/abs/2108.08757 and https://arxiv.org/abs/2109.06251.
Our approach follows closely the procedure introduced by Han, Hartnoll and Kruthoff
2
N
3
The application to single-particle quantum mechanics has appeared in numerous recent papers. These include the work of Y. Aikawa et al
4
5
6
To view the
T
E
XThe method: a theoretical overview
The method: a theoretical overview
Let us consider a very general class of one-dimensional Hamiltonians which describe the time evolution of a quantum particle in a confining potential:
H+V(x)
2
p
2
In a given energy eigenstate of energy , operators obey the following identities of expectation values:
E
O
<[H,O]>0; <HO>E<O>
Recursion between moments
Recursion between moments
These identities allow us to generate relations between expectation values of the operators . For example: we may take and the first relation above to obtain:
x,p
Op
m
x
0m<>+m(m-1)(m-2)<>-<>
m-1
x
2
p
1
4
m-3
x
m
x
dV
dx
By choosing the operator properly, we can generate a few of these such equations and combine them in order to eliminate the momentum operator from inside the expectation values. The result is a recursion relation for the positional moments of energy eigenstates . This recursion depends on the energy and perhaps a few moments:
O
<>
n
x
02mE<>+m(m-1)(m-2)<>-<>-2m<V(x)>
m-1
x
1
4
m-3
x
m
x
dV
dx
m-1
x
If the potential is confining these moments will be well-defined and finite. The case is nothing but the familiar virial theorem:
m1
E<x>+<V(x)>
1
2
dV
dx
To use this we need to initialize the recursion. The number of initialization parameters differs depending on the model. The hydrogen recursion relation is unique, as it requires only the energy to initialize (the only other such example is the harmonic oscillator). All higher radial moments are fully determined by the energy.
The moment recursion can be considered the quantum-mechanical version of the Schwinger-Dyson equations, which in field theory relate different correlation functions. In the matrix model application of the bootstrap, this role is played by the loop equations. Fundamentally it is a set of equations, dictated by averages over dynamics, which relate data that characterize states.
The moment recursion can be considered the quantum-mechanical version of the Schwinger-Dyson equations, which in field theory relate different correlation functions. In the matrix model application of the bootstrap, this role is played by the loop equations. Fundamentally it is a set of equations, dictated by averages over dynamics, which relate data that characterize states.
Recursion for the hydrogen model
Recursion for the hydrogen model
For our treatment of the hydrogen model we take the following radial Hamiltonian:
H+-
2
p
2
ℓ(ℓ+1)
2
2
r
1
r
We work in natural units. This Hamiltonian contains an angular momentum term, and we can choose to vary the azimuthal quantum number as we see fit (Indeed, something you may choose to investigate in this notebook is the spectrum for non-integral values of ). The quantum number is normally quantized when passing from the 3d Schrodinger equation to the radial Hamiltonian we use here.
The system is exactly solved. Eigenstates of fixed energy are indexed by the principal quantum number. The energies, in our conventions, are given by
ℓ
ℓ
ℓ
The system is exactly solved. Eigenstates of fixed energy are indexed by the principal quantum number
n≥0
E
n
1
2
2
(n+1)
where we take the convention that the ground state is , the first excited state is , and so forth. The exact wavefunctions are schematically ; the exact radial moments are thus linear combinations of gamma functions and thus grow approximately factorially: . You can check for yourself:
n0
n1
ψ(r)∼p(r)
-r/a
e
<>∼k!
k
r
In[]:=
(*exactsolution*)hydrogenWF[r_,n_,l_]:=2/n^2*Sqrt[Factorial[n-1-l]/(Factorial[n+l])]*(2r/n)^l*Exp[-r/n]*LaguerreL[n-1-l,2l+1,2r/n];(*youcancheckthattheseareproperlynormalized*)(*computetheradialmoment*)radMom[k_,n_,l_]:=Integrate[r^(2+k)*hydrogenWF[r,n,l]^2,{r,0,Infinity}];(*radialmomentsforsomestate*)nn=2;ll=0;Print["<r^k>:"];radMom[k,nn,ll]//TraditionalForm
<r^k>:
Out[]//TraditionalForm=
This means the moment sequences can be highly sensitive to perturbation of the initial data for the recursion. Such sensitivity may be part of what makes the algorithm converge so quickly.
Using the bootstrapping formulae from above, we obtain a recursion relation for the radial moments :
<>
n
r
08mE<>+(m-1)[m(m-2)-4ℓ(ℓ+1)]<>+4(2m-1)<>
m-1
r
m-3
r
m-2
r
We are working on the half line, where the moments > are well-defined for . The virial theorem ( in the above) reduces to a relation between the energy and the -1-th moment:
<
n
r
n≥-2
m1
E-<>
1
2
1
r
Notably, this already implies that . This relation, in combination with the normalization constraint , initialize the recursion. As a result, we only need guess the value of the energy. That will determine all the radial moments>, .
E<0
<>1
0
r
<
n
r
n≥1
Constraints from positivity
Constraints from positivity
Given some initial data for the recursion, we can compute an arbitrary number of moments. The second step of the bootstrap is determining which of these sequences can actually correspond to the moment sequence of a normalizable, positive measure. This is of course required, since these moment sequences should be derived from a probability distribution associated to a wavefunction. Since we are considering the hydrogen model, we will work on the half-line. Consider the operator . In an energy eigenstate, we must have ; this is the norm of some state in the Hilbert space. This means:
O
K
∑
n
c
n
n
r
<O>≥0
†
O
<O><>≡Mc≥0
†
O
K
∑
n,m
*
c
m
n+m
r
c
n
K
∑
n,m
*
c
m
M
nm
c
n
T
c
where we have defined the matrix <>, 0≤n,m≤K-1 by its elements. This condition, which must be true for all complex numbers , implies that the matrix is positive (semi)-definite: , for all . We can play the same game by defining the operator , where the operator square root is well-defined and unique since . This defines a similar relation:
KK
M
nm
n+m
r
c
n
M
M⪰0
K>0
O
r
Or>0
<>c≥0
†
O
O
T
c
M
where this time, <>. So these two matrices must (both) be positive (semi)-definite, simply for consistency with the positivity of the norm. You can think that this is essentially a unitarity constraint. Positivity of these Hankel matrices is exactly the condition proven in 1894 by Stieltjes in his solution of the classical moment problem on the half line. The classical moment problem of Stieltjes asks: given a sequence of real numbers, does there exist a positive measure on the half line such that the sequence is the moment sequence of the distribution: <? The solution is that each Hankel matrix , defined as above, should be positive (semi)definite for all . We call the “depth” of the bootstrap; it controls the number of constraints. The moment problem was later solved on the real line by Hamburger, and in the trigonometric instance by Carathéodory and Toeplitz. The key is showing that the positivity condition is sufficient to guarantee existence of a measure. With the theoretical description of the bootstrap complete, we move on to actually implementing it.
M
nm
m+n+1
r
a
n
μ
a
n
n
r
>
μ
KK
M,
M
K>0
K
7
Bootstrapping in Mathematica
Bootstrapping in Mathematica
Algorithmic Overview
Algorithmic Overview
The bootstrap algorithm consists of three steps, which are repeated at each depth . Given some set of initial data (in our case, values of the energy in eigenstates), we:
K
1
.Generate a moment sequence using the recursion relation derived from the Hamiltonian.
2
.Create Hankel matrices <>, <> and check whether they are both positive definite.
M
nm
n+m
r
M
nm
n+m+1
r
3
.If the energy value satisfies the check at depth , keep it. If not, reject it.
K
The result is a (possibly disjoint) set of energy values allowed at depth . Conveniently, the energies allowed at depth are necessarily contained within this set: ⊂. This is because the matrix is the principal submatrix of , and so must be positive if the matrix ⪰0. What this means that the same steps can be repeated on only the data returned from the previous depth. From a computational perspective, this allows for a huge speedup. If data is rejected at some depth it must continue to be rejected at all higher depths. This weak notion of convergence allows us to only shrink the search space as we increase the constraint depth. The algorithm which implements this is thoroughly commented inline.
Below, we construct a simple but efficient bootstrap program to numerically solve the hydrogen model. The code blocks are annotated to explain their part to play in the algorithm. Precision can be manually set--for high depths this is quite necessary to avoid numerical artifacts.
E
K
K
K+1
E
K+1
E
K
KK
M
K
M
K+1
M
K+1
Below, we construct a simple but efficient bootstrap program to numerically solve the hydrogen model. The code blocks are annotated to explain their part to play in the algorithm. Precision can be manually set--for high depths
K
Moments & Matrices
Moments & Matrices
The following block of code uses the recursion relation, derived from the Hamiltonian, to generate a sequence of moments from to for some .
<>
-1
r
<>
N
r
N>0
In[]:=
(*thisfunctiongeneratesamomentsequenceoflengthNgivenanenergy,thequantumnumberL,andaprecisionp*)Seq[Ε_,L_,Nn_,p_]:=Module[{A},A=N[{-2*Ε,1},p];(*initialdata,normalization*)Do[AppendTo[A,N[(s(4L(L+1)-(s+1)(s-1))A[[s]]-4(2s+1)A[[s+1]])/(8(s+1)*Ε),p]];(*momentrecursion*),{s,Nn}];Return[A];];(*inputanenergyvaluehere!*)Energy=N[-0.5];Ell=0;length=10;precision=MachinePrecision;ss=Seq[Energy,Ell,length,precision]
Out[]=
{1.,1.,1.5,3.,7.5,22.5,78.75,315.,1417.5,7087.5,38981.3,233888.}
You can see how quickly the sequence tends to grow. Below, the moment sequence is repackaged into the two Hankel matrices, whose positivity will control acceptance/rejection of the initial guess for the energy. Since the sequences grow so quickly, we rescale the matrix elements as /(). This rescaling preserves positivity of the eigenvalues and decreases the absolute size of the matrix elements.
M
nm
M
nm
M
1m
M
n1
In[]:=
(*thisfunctionreturnsthetwoKxKHankelmatricesfromanenergy,l,withprecisionp*)MomentMatrices[Ε_,L_,K_,p_]:=Module[{sequence,s0,s1,M0,M1},sequence=Seq[Ε,L,2*K,p];(*createmomentsequence*)M0=N[Table[sequence[[i+j]]/(sequence[[i+1]]*sequence[[j+1]]),{i,K},{j,K}],p];(*createfirstHankelmatrix,rescalingindividualminorstopreservepositivitybutdecreasenumericalsizeofentries*)M1=N[Table[sequence[[i+j+1]]/(sequence[[i+2]]*sequence[[j+2]]),{i,K},{j,K}],p];(*createsecondHankelmatrix,rescaledasabove*)Return[{M0,M1}]];(*thisfunctioncheckspositivityofbothHankelmatrices*)PosCheck[Ε_,L_,K_,p_]:=Module[{M0,M1,mats,b0,b1},mats=MomentMatrices[Ε,L,K,p];M0=mats[[1]];M1=mats[[2]];(*checkpositivityof(rescaled)Hankels*)b0=Boole[PositiveDefiniteMatrixQ[M0]];b1=Boole[PositiveDefiniteMatrixQ[M1]];Return[b0*b1](*returnthelogical'and'asbinary*)];
Run the following to search a range of depths and plot the results.
To use the rest of the cells, the object exampleSearch must be created by running the cell directly above. Below, you can plot the results. The grid lines are located at the exact energies.
In the following we animate the convergence of the bootstrap data just found.
Convergence Properties
Convergence Properties
A bootstrapped spectrum
A bootstrapped spectrum
At each depth, our data are lists of intervals. To get the most naive guess for the correct energy we can use the midpoints of these intervals; the intervals will naturally furnish error bars. Given the set of results of the bootstrap run above (exampleSearch), we can extract these data at each depth and show another picture of the convergence. This time, it’s characterized as a numerical spectrum of low-lying energies.
The exact energies are marked as gridlines. To give an idea of the precision, print the ground state and excited state values with uncertainties, at the highest depth tested:
Outlook
Outlook
Non-integral L?
Non-integral L?
References
References