Landau-Zener effect

Quantum Mechanics
K. Konishi, G. Paffuti

Solve numerically the Schrodinger equation for a Hamiltonian (Landau-Zener effect):

H = (({{(Subscript[E, 1])[t], F}, {(SuperStar[F]), (Subscript[E, 2])[t]}}));     (Subscript[E, 1])[t]=- v (t/2);  (Subscript[E, 2])[t]=+ v (t/2).

Landau-Zener transition

Off[General::"spell",General::"spell1",InterpolatingFunction::"dmval"];
Off[ClearAll::"clloc"];
ClearAll["Global`*"];

Short review of the theory

In the text the Landau Zener approximation for non-adiabatic transition has been exemplified by considering a two state system with Hamiltonian

H = (({{(Subscript[E, 1])[t], F}, {(SuperStar[F]), (Subscript[E, 2])[t]}}));       i ℏ ((∂ψ)/(∂t)) = H ψ

with (Subscript[E, 1])[0]=(Subscript[E, 2])[0]. This almost degenerate system has adiabatic eigenvalues (i.e. eigenvalues of H at fixed t)

(Subscript[λ, ±])[t] = (1/2)((Subscript[E, 1])[t]+(Subscript[E, 2])[t]± (Sqrt[((((Subscript[E, 1])[t]-(Subscript[E, 2])[t]))^2)+ 4|F((|)^2)])).

As a concrete example let' s take

E1[t_]:= -(v/2) t; E2[t_]:=  (v/2) t;
λ1[t_]:=(1/2)(E1[t]+ E2[t] + (Sqrt[(((E1[t]- E2[t]))^2)+ 4  (F^2)]));
λ2[t_]:=(1/2)(E1[t]+ E2[t] - (Sqrt[(((E1[t]- E2[t]))^2)+ 4  (F^2)]));

The behavior of eigenvalues with time is sketched below

v=1;F=0.2;

Plot[{λ1[t],λ2[t],E1[t],E2[t]},{t,-2,2},Axes-> None,
PlotStyle-> {Blue,Blue,{Dashed,Red},{Dashed,Red}},
Epilog-> {Text["(SubPlus[λ])(t)",{-0.9,.7}],Text["(SubMinus[λ])(t)",{-0.9,-0.7}],Text["(Subscript[E, 1])(t)",{-1.4,+0.5}],Text["(Subscript[E, 2])(t)",{-1.4,-0.5}] }]

For large negative t, the state |b⟩= (0,1) is the found state and for an adiabatic evolution it will remain the ground state for large t. At t≃0 a sizable probability exists to switch to higher state and the system will be in the excited state.

Let us stress that at t→∞ the state (0,1) is the excited state (as (Subscript[E, 2])> (Subscript[E, 1])≫ F ), then if U[t] is the evolution operator the amplitude for the switching is

S=⟨b|U[∞]|b⟩;   P = |S((|)^2).

In the text it has been shown that

P = Exp[-(( 2π |F((|)^2))/(ℏ v))].

The Scrödinger equation () can be solved analytically in the simple case at hand but we prefer to give a faster numerical solution.

Byextracting the adiabatic change of phase

ψ = (({{A[t]Exp[-i (∫_(-∞)^t) (Subscript[E, 1])[t]}, {B[t]Exp[-i (∫_(-∞)^t) (Subscript[E, 2])[t]}})),

the Scrödinger equation takes the form (see text)

i ℏ(dA/dt)=B F Exp[-i(∫_(-∞)^t) dτ ω[τ]];     i ℏ(dB/dt)=A (SuperStar[F]) Exp[i(∫_(-∞)^t) dτ ω[τ]];

with

ω[t] = (((Subscript[E, 2])-(Subscript[E, 1]))/ℏ)= ((v t)/ℏ).

The boundary conditions being

A[-∞] = 0;    B[-∞]= 1.

Below we will solve numerically this equation, substituting a large value T for - ∞.

Numerical solution

We will use units ℏ = 1. The numerical solution is straightforward

ClearAll["Global`*"];

eq = {I A'[t]- B[t] f Exp[-I v (((t^2)-(T^2))/2)]== 0,I B'[t]-A[t]f Exp[I v (((t^2)-(T^2))/2)]== 0,
A[-T]== 0, B[-T]== 1};

v=0.2;f=0.2;T = 200.;

sol={A,B}/.First[NDSolve[eq,{A,B},{t,-T,T},MaxSteps->50000]];

expectedP = Exp[-2π (f^2)/v]

Plot[Abs[sol[[2]][t]]^2,{t,-T,T},
AxesLabel-> {"t","P[t]"},
Epilog-> {Red,Dashed,Line[{{0,expectedP},{T,expectedP}}]}]

As expected the evolution is almost adiabatic up to the time of quasi degeneracy of the levels (here t∼0), then a period of oscillation precede a quasi adiabatic evolution in the new state (at t > 0 the state B is the excited state).

A more systematic analysis can be perfomed noticing that the change of variables t = τ/(Sqrt[v]) in eq() bring it in the form

i ℏ(dA/(dτ))=B (F/(Sqrt[v])) Exp[-i(∫_(-∞)^τ) dx (x^2)];     i ℏ(dB/dt)=A ((SuperStar[F])/(Sqrt[v])) Exp[i(∫_(-∞)^τ) dx (x^2)];

so the solution has the form

B = B[τ,((f^2)/v)],

this explain the functional dependence of B[∞,f,v], and suggests that the discriminant between fast and slow transtions depends on the value of (f^2)/v.

In these new variables the equation is practically identical to the previous one, with fv = f/(Sqrt[v]):

ClearAll["Global`*"];

eq = {I A'[τ]- B[τ] fv  Exp[-I  (((τ^2)-(T^2))/2)]== 0,I B'[τ]-A[τ]fv Exp[I  (((τ^2)-(T^2))/2)]== 0,
A[-T]== 0, B[-T]== 1};

We show three different approaches to the asymptotic value of B below

fvValues = {0.2,0.6,1.1};
T = 130;

Table[
fv = fvValues[[k]];
expectedP = Exp[-2π (fv^2)];
sol={A,B}/.First[NDSolve[eq,{A,B},{τ,-T,T},MaxSteps->100000]];
Plot[Abs[sol[[2]][t]]^2,{t,0,T},
AxesLabel-> {"τ","P[τ]"},
Epilog-> {Red,Dashed,Line[{{0,expectedP},{T,expectedP}}]}],{k,1,Length[fvValues]}]

It can be interesting to have an idea of the mean evolution, averaging over fast flutuations. This can be obtained by averaging over bins of length r using the Mathematica function MovingAverage[data,r]. Here is an example

τValues = Range[-T/2,T,0.1];
lengthBins = 20;

Table[
fv = fvValues[[k]];
expectedP = Exp[-2π (fv^2)];
sol={A,B}/.First[NDSolve[eq,{A,B},{τ,-T,T},MaxSteps->100000]];
smoothedData = MovingAverage[Abs[sol[[2]][τValues]]^2,lengthBins];
smoothedTau = MovingAverage[τValues,lengthBins];
ListPlot[Transpose[{smoothedTau,smoothedData}],
AxesLabel-> {"τ","P[τ]"},
Epilog-> {Red,Dashed,Line[{{0,expectedP},{T,expectedP}}]}],
{k,1,Length[fvValues]}]