Using Computational Essays to Demystify Quantum Mechanics

Mohammad Bahrami, PhD 
Document and Media Systems, Wolfram Reaseach
and
California State University at Los Angeles

Join the Conversation #WolframTechConf

ABSTRACT

Even after 100 years, there is no consensus on the meaning of quantum theory. Many different interpretations and alternative theories have been proposed to resolve this impasse. However, it is quite challenging to understand the implications of different proposals. I will discuss how Wolfram Language provides the technology to address this issue. I will discuss how a computational quantum laboratory can be designed in the form of a computational essay. The core idea of alternative theories (or thought-experiments) can be delivered as computational outputs in an interactive environment. Therefore, one will be able to directly observe the description of alternative models vs. quantum theory, or the consequences of a thought-experiment. Students/researchers will manipulate relevant parameters in the computational essay, to elaborate on differences in more detail.  I will discuss how a computational essay can bridge the gap between the abstract formulation and the qualitative understanding. I will comment how this style of communicating in the form of a computational essay can actively engage someone in the learning process. I will focus on two alternative theories, de Broglie-Bohm pilot wave and the spontaneous collapse models. While in the pilot wave model, particles follow deterministic trajectories, in collapse models particles have a random Brownian motion. I will demonstrate how these characteristic features can communicated effectively within a computational essay.

De Brogile-Bohm pilot-wave model

De Broglie-Bohm in a nutshell

Quantum theory finds its explication in the statistical analysis of Bohmian mechanics following Boltzmann’s ideas. Bohmian mechanics is about the motion of particles. The positions of particles constitute the primitive variables, the primary ontology: whenever you say particle, mean it! The key insight for analyzing Bohmian mechanics lies within the foundations of statistical physics. Typicality is the ticket to get to the statistical import of Bohmian mechanics, which is succinctly captured by ρ=((Abs[ψ])^2)

Deterministic motion of particles (two equations)

Probabilistic reasoning within deterministic physics following Boltzmann’s approach

Probabilities in de Broglie-Bohm model

The objectification of a random measurement result and corresponding statistic (i.e. Born rule) are postulated in standard quantum theory. However, in alternative quantum theory such as de Broglie-Bohm model, postulates are explained (or derived) from underlying dynamics.

Two summarizing questions (once asked by M. von Smoluchowski):

How can it be that the effect of chance is computable, i.e., that random causes can have lawlike effects?

How can randomness arise if all events are reducible to deterministic laws of nature? Or in other words, how can lawlike causes have random effects?

The answer to both is typicality, and this answer has in principle been known since the work of Ludwig Boltzmann: how to reason probabilistically within deterministic physics?

Typicality

The probability theory is the doctrine of typical behavior of physical systems.

The word "typical" is synonymous with "overwhelmingly many". Overwhelmingly many phase points will then mean that the relevant phase points make up a set of very large measure. What we have is a new mechanics, with the old statistical reasoning. Quantum equilibrium refers to the typical behavior of particle positions given the wave function. The particles are in equilibrium given the wave function.

Quantum Equilibrium Hypothesis: for an ensemble of identical systems, each having the wave function ψ, the typical empirical distribution of the configurations of the particles is given approximately by  ρ= |ψ((|)^2). In short, Born' statistical law holds.

Deterministic equations of motion

Bohmian mechanics is the new mechanics for point particles. In the equations for Bohmian mechanics there are parameters (Subscript[m, 1]),...,(Subscript[m, N]) which we shall call masses. A particle of mass (Subscript[m, i]) is described by its position (Subscript[(OverVector[Q]), i]). For an N-particle system, we have:

(∂_t)(Subscript[(OverVector[x]), i])=(ℏ/(Subscript[m, i]))Im[(((Subscript[∇, i])ψ(X,t))/(ψ(X,t)))]=(1/(Subscript[m, i]))(Subscript[∇, i])S(X,t),            i=1,2,...,N

where S is the phase of the wavefunction ψ=R (e^(i S/ℏ)), (Subscript[∇, i])=(∂/(∂(Subscript[(OverVector[x]), i]))) and X={(Subscript[(OverVector[x]), 1]),(Subscript[(OverVector[x]), 2]),...,(Subscript[(OverVector[x]), N])} is the set of configurations.

The wave function ψ follows the Schrodinger equation:

i ℏ (∂_t)ψ(X,t)=(-(∑_(i=1)^N)(((ℏ^2) )/(2 (Subscript[m, i])))(Subscript[∇, i]^2)+V(X))ψ(X,t).

Bohmian trajectories, as described by Eq.(1), are nothing but the flux lines along which the probability gets transported. The velocity field is thus simply the tangent vector field of the flux lines. The vector field is given by

(Subscript[(OverVector[v]), i])=(ℏ/m)Im[(((Subscript[∇, i])ψ)/ψ)]=(Subscript[∇, i])S.

An important implication of Eqs.(1-3) is that the total potential of a quantum system is given by V(X)+Q(X) at X=X(t) which is obtained from Eq.(1) and

Q(X)=-((ℏ^2)/(2m))(((∇^2)R)/R),

where Q(X) is called the quantum potential. Thefore, the local expectation value of Hamilotnian operator along a trajectory is given by:

E=((((∇S))^2)/(2m))+Q+V,

Quantum Tunneling: all mysteries together!

Numerical Evaluation of Schrodinger equation

Time interval (tf) and the spatial dimension of the box (L)

tf=10.;L=100.;

Initial wave function

μ0=-6.;σ0=1.;
ψ0[x_]=Exp[- (x-μ0)^2/((σ0)^2)+2I x];
Row[{Plot[Abs[ψ0[x]]^2,{x,-.1L,.1L},PlotRange->All,Frame->True,FrameLabel->{"x","|ψ((|)^2)"},FrameTicks->None,ImageSize->{400}],Spacer[40],"|ψ(0,x)((|)^2) vs x"}]

Scattering particles off of a potential barrier in one dimension

vL=1;(*barrier spatial size*)
v[x_]=6Exp[-x^2](*HeavisideTheta[x+vL] HeavisideTheta[vL-x]*);
Plot[{v[x],Style[4,Dashed]},{x,-.1L,.1L},PlotRange->All,Exclusions->None,PlotLegends->{"V(x)","Initial KE"},Frame->True,FrameTicks->None]

Schrodinger equation

eqns={I D[u[t,x],t]==-0.5D[u[t,x],x,x]+v[x]u[t,x],u[0,x]==ψ0[x],u[t,L]==0,u[t,-L]==0};

Evolution of the maginary and the real parts of wave funciton: ψ=(Subscript[ψ, R])+(Subscript[iψ, I])

expandedeqns=Flatten[ComplexExpand[ReIm[eqns/.{u[x_,t_]->uR[x,t]+I uI[x,t],Derivative[y__][u][x_,t_]->Derivative[y][uR][x,t]+I Derivative[y][uI][x,t]}
]]];

Numerical solution of Schrodinger equation

sol=Flatten@NDSolve[Evaluate[expandedeqns],{uR,uI},{t,0,tf},{x,-L,L},AccuracyGoal->5]

Probability distribution:  ρ=|ψ((|)^2)

ρ[t_,x_]=(uI[t,x]^2+uR[t,x]^2)/.sol;

Time evolution of the wave packet and its deflection and transition through the potential

Manipulate[
Plot[{ρ[t,x]},{x,-.1L,.1L},PlotRange->{0,ρMax},Ticks-> None,GridLines-> {{-vL,vL},None},GridLinesStyle->Directive[Orange,Dashed],AxesLabel->{"x","|ψ(t,x)((|)^2)"}],
{{t,0,Row[{Spacer[50],"time"}]},0,tf},{{ρMax,1,"Max. ρ"},1,.1,VerticalSlider},ControlPlacement->{Top,Left}]

Bohm trajectory

Obtain the argumnet of the wave function as ψ=((Re)^(iS))

S[t_,x_]=ArcTan[uR[t,x],uI[t,x]]/.sol;

Bohm trajectory starting from the initial point q0 and the final time tf

bohmTrajectory[q0_,tf_]:=NDSolve[{q'[t]==Derivative[0,1][S][t,q[t]],q[0]==q0},q,{t,0,tf}];

Highlight initial positions where the initial ρ is nonzero

highlightρ=Plot[{μ0-2σ0,μ0+2σ0},{t,0,tf},PlotRange->{-.1L,.1L},Filling->{1->{{2},LightGray}},PlotStyle->{{Dashed,White},{Dashed,White}}];

Highlight the potential barrier

highlightV=Plot[{-vL,vL},{t,0,tf},PlotRange->{-.1L,.1L},Filling->{1->{{2},LightOrange}},PlotStyle->{{Dashed,White},{Dashed,White}}];

Obtain Bohm trajectories for initial positions where ρ is nonzero

trajectories=Table[Plot[{Evaluate[q[t]/.bohmTrajectory[x0,tf]]},{t,0,tf},PlotRange->{-.1L,.1L},Frame->True,Ticks->None,PlotStyle->Directive[Red,Opacity[2ρ[0,x0]]]],{x0,μ0-2σ0,μ0+2σ0,.2}];

Show Bohm trajectories

Row[{Show[{highlightρ,highlightV,trajectories},Frame->True,FrameLabel->{"time","position"},FrameTicks-> None,ImageSize->{400}],Spacer[40],"Bohm trajectories"}]

ρBackground=DensityPlot[ρ[t,x],{t,0,tf},{x,-.1L,.1L},PlotPoints->200,ColorFunction->GrayLevel];

Row[{Show[{ρBackground,trajectories},Frame->True,FrameLabel->{"time","position"},FrameTicks-> None,ImageSize->{400}],Spacer[40],"Bohm trajectories in ρ background"}]

Quantum potential

Quantum potential

Q[t_,x_]=-(1/2)((D[(Sqrt[ρ[t,x]]),x,x])/(Sqrt[ρ[t,x]]))/.sol;

Plot quantum potential

qPotential=Plot3D[Q[t,x],{t,0,tf},{x,-10,10},PlotPoints->50,Ticks->None,AxesLabel->{"time","position","Q"},ImageSize->{450}];

Row[{qPotential,"Quantum potential vs time and position"}]

Find quantum potential Q along a trajectory as a function of time for given initial position x0

qP[x0_]:=Block[{posi},
posi[t_]=First@(q[t]/.bohmTrajectory[x0,tf]);
ParametricPlot3D[{t,posi[t],Q[t,posi[t]]},{t,0,tf},PlotStyle->Blue,Exclusions->None]];

Find quantum potential Q along a trajectory as a function of time for a set of given initial positions

seriesOfTr=qP[#]&/@Range[-8,-4,.2];

Find quantum potential Q along a trajectory as a function of time

Show[qPotential,seriesOfTr]

Probabilities: random walk as the ideal Galton board

(Entity["Galton board", Entity["MathWorld","GaltonBoard"], "Entity[\"MathWorld\", \"GaltonBoard\"]", "MathWorld topic"])["TypesetDescription"]

Galton board: relative number of balls in slots vs a normal distrbution for number of balls in slots

Row[{Manipulate[
sd=(data["LastValues"])[[1;;m]];
Show[Histogram[sd,Automatic,AxesOrigin->{0,0},Ticks->{Automatic,None}],
Plot[ 2000 PDF[NormalDistribution[0,(Sqrt[n])],x]//Evaluate,{x,-n,n},PlotRange->All,PlotStyle->Red],PlotRange-> {{-n,n},All}],
{{m,0,"number of balls"},1,1000,1,Appearance->"Labeled"},
{{n,20,"number of slots"},{20,40,80,100}},
Button[Style["new simulation","Label",Small],{data=RandomFunction[RandomWalkProcess[.5],{0,n},1000],m=0}],
Initialization:>{
data=RandomFunction[RandomWalkProcess[.5],{0,20},1000];},TrackedSymbols->{m}
],"Relative number of balls in slots vs a normal distrbution for number of balls in slots"}]

Collapse models

“The abrupt change by measurement ... is the most interesting point of the entire theory.”

E. Schrodinger, in what is often called his “cat paradox” paper

“I think equation guessing might be the best method to proceed to obtain the laws for the part of physics which is presently unknown.”

Richard Feynman, in his 1965 Nobel Prize Lecture

Collapse models in a nutshell

Collapse models are modiffied Schrodinger’s equation, which provides a continuous evolution of the wave function
in all situations. And, in particular, when applied to a measurement situation, the dynamical evolution of the state vector is such that it ends up consistent with the collapse postulate. Collapse models have two equations:

a modified Schrodinger equation

i ℏ ((∂ψ)/(∂t))=((OverHat[H])+i(Subscript[(OverHat[V]), w]))ψ

As a result of the non-Hermiticity of new Hamiltonian, the state vector squared norm |ψ(x,t)((|)^2) changes magnitude as time evolves, depending upon a noise field w(x,t).

probability nature supplies a particular noise field w(x,t), is proportional to |ψ(x,t)((|)^2)

One finds that highly probable w(x,t) (i.e. large |ψ(x,t)((|)^2)), in measurement situations, rapidly, dynamically, bring about wave function collapse. And, not only does the state vector then correspond to one of the possible outcomes we observe, these outcomes occur with the Born probabilities.

Dynamics of state collapse: Gambler’s Ruin Game

Dynamical collapse

Let us consider a quantum state with two possible outcomes for a measurement of observable A: (Subscript[a, 1]) and (Subscript[a, 2]). Normalization implies that the summation of probablaities for each outcome is one:

(Subscript[P, 1])(t)+(Subscript[P, 2])(t)=1

Collapse means that after sufficiently large enough time, we have:

either   (Subscript[P, 1])(t->∞)=1,  (Subscript[P, 2])(t->∞)=0   occures with the probability (Subscript[P, 1])(0)

or   (Subscript[P, 1])(t->∞)=0,  (Subscript[P, 2])(t->∞)=1   occures with the probability (Subscript[P, 2])(0)

Gambler’s Ruin Game

(Entity["gambler’s ruin", Entity["MathWorld","GamblersRuin"], "Entity[\"MathWorld\", \"GamblersRuin\"]", "MathWorld topic"])["TypesetDescription"]

Final result of one run of game depending on the initial number of pennies for the 1st player

results[n1_]:=NestWhile[(#+RandomChoice[{1,-1}])&,n1,!(#==0||#==100)&];

Final result of many games (either winner with $1 or a loser with nothing)

finalResultList[n1_,0]={};
finalResultList[n1_,trials_]:=finalResultList[n1,trials]=Join[results[n1],finalResultList[n1,trials-1]];

Frequency of wins after some games for a given initial bet (Subscript[n, 1])

freqWins[n1_,attempts_]:=Table[{trial,((Count[finalResultList[n1,trial],100])/trial)*100},{trial,attempts}];

ListPlot of wins frequency vs number of times the game played

plot[n1_,trial_]:=ListPlot[freqWins[n1,trial],Joined->True,PlotRange->{All,{0,101}},GridLines->{None, {n1}},GridLinesStyle->Directive[Red, Dashed],Frame->True,FrameLabel->{"# of times that the game is repeated","% win"},PlotRange->{All,{-5,101}}]

ListPlot of wins frequency vs number of times the game played for different initial amount of pennies and different number of times the game is played

Row[{Manipulate[plot[n1,trial],{{n1,30,"initial number of pennies for 1st player"},1,99,1,Appearance->"Labeled"},{{trial,100,"# of times that the game is repeated"},1,1000,1,Appearance->"Labeled"}],Spacer[20],"%wins in a Gambler Ruin Game for different initial bet as a function of number of times that the game is repeated"}]

The simplest possible case: collapse in a 2D Hilbert space

Let us consider the simplest possible case, collapse in a two-state Hilbert space, analogous to the two-gambler game. We are going to set (OverHat[H])=0 for the rest of the present discussion, so as to see the collapse dynamics take place without interference of the usual Hamiltonian dynamics. The solution to the usual Schrodinger equation is (ℏ->1):

(Subscript[ψ, t])=Exp[(∫_0^t)(Subscript[(OverHat[V]), w])dt](Subscript[ψ, 0])=Exp[-(1/(4λt))((((Subscript[w, t])-2λ t(OverHat[A])))^2)](Subscript[ψ, 0])
                                             =(Subscript[c, 1])(e^(-(1/(4λt))((((Subscript[w, t])-2λ t (Subscript[a, 1])))^2)))(Subscript[ψ, (Subscript[a, 1])])+(Subscript[c, 2])(e^(-(1/(4λt))((((Subscript[w, t])-2λ t (Subscript[a, 2])))^2)))(Subscript[ψ, (Subscript[a, 2])])

where

(Subscript[ψ, 0])=(Subscript[c, 1])(Subscript[ψ, (Subscript[a, 1])])+(Subscript[c, 2])(Subscript[ψ, (Subscript[a, 2])])        with ((|(Subscript[c, 1])|)^2)+((|(Subscript[c, 2])|)^2)=1

and (OverHat[A]) is the collapse operator with the eigenvalues (Subscript[a, 1]) and (Subscript[a, 2]), (Subscript[w, t]) is a Brownian motion trajectory (i.e. Wiener process), and λ is the collape rate. The probability that nature will provide a particular Brownian random function (Subscript[w, t]) is:

P((Subscript[w, t]))d(Subscript[w, t])=((d(Subscript[w, t]))/(Sqrt[2π λ t]))((Abs[ψ])^2)=((d(Subscript[w, t]))/(Sqrt[2π λ t]))(((|(Subscript[c, 1])|)^2) (e^(-(1/(2λt))((((Subscript[w, t])-2λ t (Subscript[a, 1])))^2)))+((|(Subscript[c, 2])|)^2) (e^(-(1/(2λt))((((Subscript[w, t])-2λ t (Subscript[a, 2])))^2))))

Since the modified Schrodinger equation is not normalized, one needs the normalized state vector to make predictions:

(Subscript[ψ, t])=(1/N)((Subscript[c, 1])(e^(-(1/(4λt))((((Subscript[w, t])-2λ t (Subscript[a, 1])))^2)))(Subscript[ψ, (Subscript[a, 1])])+(Subscript[c, 2])(e^(-(1/(4λt))((((Subscript[w, t])-2λ t (Subscript[a, 2])))^2)))(Subscript[ψ, (Subscript[a, 2])]))

where N is normalization constant

N=(Sqrt[((|(Subscript[c, 1])|)^2) (e^(-(1/(2λt))((((Subscript[w, t])-2λ t (Subscript[a, 1])))^2)))+((|(Subscript[c, 2])|)^2) (e^(-(1/(2λt))((((Subscript[w, t])-2λ t (Subscript[a, 2])))^2)))])

Numerical simulation of collapse dynamics in 2D Hilbert space

Random motion of probability

Row[{DynamicModule[{data,randomP1},
Manipulate[
ListLinePlot[data,GridLines->{None,{p1}},PlotRange->{All,{0,1}},Frame->True,FrameLabel->{"time","(Subscript[P, 1])(t)"}],{{λ,1,"Collapse rate λ (Hz)"},1,10,Appearance->"Labeled"},{{realization,1},1,10,1,Appearance->"Labeled"},{{p1,.5,"(Subscript[P, 1])(0)"},0,1,Appearance->"Labeled"},Button["New realization of noise field (Subscript[w, t])",data=RandomFunction[randomP1[p1,λ],{0,10,.01},realization]]
],
Initialization:>(randomP1[p1_,λ_]:=TransformedProcess[If[t==0,p1,((p1 (E^(-(1/(2 λ t))(((w[t]-2λ t))^2))))/((p1 (E^(-(1/(2 λ t))(((w[t]-2λ t))^2)))+(1-p1)(E^(-(1/(2λ t))(((w[t]+2λ t))^2))))))],w\[Distributed]WienerProcess[],t];
(*p1 is the initial probability, λ is the collapse rate*)
data=RandomFunction[randomP1[.5,1],{0,10,.01},1];)],"Random change in the probability of observing the first outcome in a 2D Hilbert Space"}]

Histogram of probability as a function of λt

collapse rate

λ=1;

initial probability

p1=.6;

values of A measurements (eigenvalues)

a1=-1;a2=1;

final time

tf=10;

Random process of (Subscript[P, 1])(t)

randomP1=TransformedProcess[If[t==0,p1,((p1 (E^(-(1/(2 λ t))(((w[t]-2a1 λ t))^2))))/(p1 (E^(-(1/(2 λ t))(((w[t]-2a1 λ t))^2)))+(1-p1)(E^(-(1/(2λ t))(((w[t]-2a2 λ t))^2)))))],w\[Distributed]WienerProcess[],t];

Simulation of random process (Subscript[P, 1])(t)

data=RandomFunction[randomP1,{0 ,tf,.01},100];

Histogram of results

Row[{Manipulate[
SmoothHistogram[data["SliceData",λt],.01,PlotRange->{{-.05,1.05},All},Ticks->{Automatic,None},AxesLabel->{"(Subscript[P, 1])","Rel. frequency of (Subscript[P, 1])"}],{λt,0,.4tf,Appearance-> "Labeled"}],"Relative frequency of (Subscript[P, 1]) for differnt values of (Subscript[P, 1]) as a function of λ.t"}]

References

D. Durr, and S. Teufel, Bohmian Mechanics (Springer 2009).

Philip Pearle, Realist Approaches to Quantum Theory (book in preparation)

Acknowledgement

Stephen Wolfram (Wolfram CEO and founder)

Swede White (Wolfram, Public Relations)

Vitaliy Kaurov (Wolfram, Document and Media Systems)

Kyle Keane (MIT and Wolfram Public Relations)

Paul Abbott (University of Western Australia)

Daavid Vaananen (Computational Thinking Initiative)

Andre Kuzniarek (Wolfram, Director of Document and Media Systems)

Danielle Rommel (Wolfram, Director of Outreach and Communications)

Philip Pearle (Physics Department, Hamilton College NYC)