# Message Boards

Answer
(Unmark)

GROUPS:

6

# Bootstrapping quantum mechanical systems

Posted 19 days ago

Bootstrapping quantum mechanical systems by D. Berenstein, G. Hulsey Department of Physics, UC Santa Barbara September 2021
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 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 X
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 2 p 2 In a given energy eigenstate of energy E O <[H,O]>0; <HO>E<O>
Recursion between moments These identities allow us to generate relations between expectation values of the operators x,p O m x 0m< m-1 x 2 p 1 4 m-3 x m x dV dx By choosing the operator O < n x 02mE< 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 m1 E 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.
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 ℓ ℓ ℓ 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 n0 n1 ψ(r)∼p(r) -r/a e < 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= 1 8 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 r m-3 r m-2 r We are working on the half line, where the moments < n r n≥-2 m1 E- 1 2 1 r Notably, this already implies that E<0 < 0 r < n r n≥1
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 O K ∑ n c n n r < † O < † 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 KK M nm n+m r c n M M⪰0 K>0 O r Or>0 < † O O T c M where this time, M nm m+n+1 r a n μ a n n r > μ KK M, M K>0 K 7
Bootstrapping in Mathematica
Algorithmic Overview The bootstrap algorithm consists of three steps, which are repeated at each depth K 1 .Generate a moment sequence using the recursion relation derived from the Hamiltonian. 2 .Create Hankel matrices M nm n+m r M nm n+m+1 r 3 .If the energy value satisfies the check at depth K The result is a (possibly disjoint) set 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 The following block of code uses the recursion relation, derived from the Hamiltonian, to generate a sequence of moments from < -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 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*)]; From here the problem is algorithmic. The following procedures search a range of energies for validity under the bootstrap constraints, at fixed depth In[]:= (*thisauxiliaryfunctionjoinsalistofIntervalobjectsintoasingle,disjointIntervalobject*)joinIntervals[intervalList_]:=Module[{result},result=intervalList[[1]];If[Length[intervalList]1,Return[result]];Do[result=IntervalUnion[result,intervalList[[i]]];,{i,2,Length[intervalList]}];Return[result];];(*thisfunctionchecksanenergyintervalforallowedenergiesgiventhebootstrapparametersK,L,withaspecifiednumberofpointstocheck*)CheckEnergies[eInterval_,L_,depth_,npts_,p_]:=Module[{emin,emax,step,pos,evals,bins,intervals,result,joined},emin=N[eInterval[[1]],p];emax=N[eInterval[[2]],p];step=N[Abs[(emin-emax)/npts],p];(*HerewecheckHankelpositivity*)pos=Table[ene*PosCheck[ene,L,depth,p],{ene,N[emin,p],N[emax-step,p],step}];(*tableofenergyvaluesifallowedand0ifdisallowedatgivendepth*)evals=Sort[DeleteCases[pos,0.],Less];If[Length[evals]<1,Return[{}]];(*removezerosandaccountforfindingnopoints*)bins=Split[evals,#2-#1<2*step&];(*splitslistofallowedenergieswhentheyareseparatedbymorethanastep-size*)intervals=Table[Interval[{bins[[i]][[1]]-step,bins[[i]][[-1]]+step}],{i,Length[bins]}];(*createsintervalsofallowedenergies,expandingthembyonestep-sizetoaccountforthefiniteresolution*)intervals=Sort[DeleteCases[intervals,{}],Less];joined=joinIntervals[intervals];(*taketheunionoftheallowedintervals*)Return[joined];];(*theexactenergies,inourconventions*)Henergy[n_]:=-1/2*1/(n+1)^2;energyList[m_]:=N[Table[Henergy[n],{n,0,m}]]; The CheckEnergies function takes an energy interval, quantum number In[]:= (*runabootstrapatfixeddepth*)depthK=14;points=1000;quantumNumberL=0;checked=CheckEnergies[{N[-1,25],0},quantumNumberL,depthK,points,25];NumberLinePlot[checked,PlotRange{-.52,0},PlotLabel"Allowed Energies, K = "<>ToString[depthK],AxesLabel"E",GridLines{energyList[7],None}] Out[]= The above shows the intervals of allowed energies, with the exact energy levels shown in gray. One can see islands split and form around the exact energies. From here, one can proceed to higher and higher depths, utilizing the fact that the set of allowed energies at depth In[]:= multiBandCheck[eInterval_,L_,K_,nptsPerBand_,p_]:=Module[{np=nptsPerBand,ni,resultantIntervals},resultantIntervals={};Do[AppendTo[resultantIntervals,CheckEnergies[eInterval[[i]],L,K,np,p]],{i,Length[eInterval]-1}];(*checkeachintervalexceptthelast*)AppendTo[resultantIntervals,CheckEnergies[eInterval[[-1]],L,K,np*10,p]];(*checkthelastintervalwith10xresolution,becausethisiswherenewdisjointintervalsareliabletoform*)resultantIntervals=DeleteCases[resultantIntervals,{}];Return[joinIntervals[Flatten[resultantIntervals,1]]];];
Bootstrapping Hydrogen The following block of code searches a range of depths by utilizing the “recursive” nature of the algorithm. First, a cursory search is made. This returns a possibly disjoint interval. At successive depths the algorithm searches only this interval, and so on. The result is data of allowed energies at each depth In[]:= KRangeSearch[initialInterval_,kRange_,L_,p_]:=Module[{eint=initialInterval,searched,kmin,kmax,xx,results},kmin=kRange[[1]];kmax=kRange[[2]];(*First,wedoacursorysearchwithalargenumberofpointsatthemindepth*)searched=CheckEnergies[eint,L,kmin,1000,p];(*theresultsare{depthK,Interval[{...}]}*)results={{kmin,searched}};Monitor[Do[(*ateachhigherdepth,searchtheresultsfromthepreviousdepth,re-griddingeachconnectedinterval*)xx=multiBandCheck[searched,L,k,50,p];AppendTo[results,{k,xx}];searched=xx;(*updatethedatatosearch;thisisthe'recursive'step*),{k,kmin+1,kmax}],k];Return[results]];(*returnsalist{{kmin,Interval[{...}]},...{kmax,Interval[{...}]}}*) Run the following to search a range of depths and plot the results. In[]:= KRange={5,15};eps=1/1000;(*inputtingexactly0canleadtonumericalissues;thisisasmallcutoff*)eRange=N[{-1-eps,-1*eps},50];azimuthL=0;exampleSearch=KRangeSearch[eRange,KRange,azimuthL,50];(*thisfunctionprintsthedepthKitiscurrentlychecking*) 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[]:= (*plottheresults*)pkr={5,15};(*K-rangetoplot*)NumberLinePlot[exampleSearch[[pkr[[1]]-KRange[[1]]+1;;pkr[[2]]-KRange[[1]]+1]],PlotRange{eRange[[1]],0},AxesLabel{"E"},PlotLegendsTable[If[Mod[k,5]0,"K = "<>ToString[k],None],{k,pkr[[1]],pkr[[2]]}],PlotLabel"Allowed Energies at Variable Depth",GridLines{energyList[5],None}] Out[]=
In the following we animate the convergence of the bootstrap data just found. In[]:= Manipulate[NumberLinePlot[exampleSearch[[K-KRange[[1]]+1]],PlotLabel"Allowed Energies",AxesLabel{"E"},PlotRangeeRange,GridLines{energyList[5],None}],{K,KRange[[1]],KRange[[2]],1},SaveDefinitionsTrue](*requiresdynamicupdating:Evaluation>DynamicUpdatingEnabled*) Out[]=
Convergence Properties Of course, we would like to see exactly how quickly the bootstrap converges. We can do this by plotting the width of each separated interval of allowed energies versus the depth K In[]:= ConvPlot[searchData_]:=Module[{sd=searchData,kr,ints,results,island,leg},kr={sd[[1]][[1]],sd[[1]][[1]]+Length[sd]-1};(*Krange*)ints=Table[Length[sd[[i]][[2]]],{i,Length[sd]}];(*intervalsfromthesearch*)results={};(*createalistcorrespondingtoeachdistinctislandofitswidthversusK*)Do[island={};Do[If[ints[[q]]≥j,AppendTo[island,{q+kr[[1]]-1,Abs[sd[[q]][[2]][[j]][[2]]-sd[[q]][[2]][[j]][[1]]]}]];,{q,Length[ints]}];AppendTo[results,island],{j,1,Max[ints]}];(*plottheresultsonalogarithmicscale*)Li |