Slow Light

Here is the general program for slow light with homogeneous broadening with light entering, slowing down and leaving the medium, but not enough memory to run entire path

Program

***medium between 0 and 7.5**

s[x_]:=UnitStep[x]-UnitStep[x-7.5]

** Eqs. (9.30) coupled with density matrix equations (9.15). f1 is the probe field and g[t] the pump field which is turned off when the pulse is in the medium ***

g[t_]:=If[t<9||t>14,1,Exp[-(t-9)^2]+Exp[-(t-14)^2]]

** pump turned off and back on when light is in medium **

Clear[c1,c2,f1,t]

sol[d_,dp_,c1_,c2_,g2_,g12_,g21_,ti_,tf_,eta_]:=NDSolve[{D[r11[t,x],t]==-I s[x] c1 Conjugate[r12[t,x]] Conjugate[f1[t,x]]+I c1  s[x]f1[t,x] r12[t,x]+g21 r22[t,x],D[r12[t,x],t]==I c1 s[x]Conjugate[f1[t,x] ]r11[t,x]+I d r12[t,x]-g12 r12[t,x]+I c2 g[t] s[x]r13[t,x]-I c1 s[x]Conjugate[f1[t,x] ]r22[t,x],D[r13[t,x],t]==-I c1 s[x] Conjugate[f1[t,x]]Conjugate[r32[t,x]]+I d r13[t,x]-I dp r13[t,x]+I c2 s[x]g[t] r12[t,x],D[r22[t,x],t]==I c1 s[x]Conjugate[r12[t,x]] Conjugate[f1[t,x]]+I c2 s[x]Conjugate[r32[t,x]] g[t]-I c1 s[x]f1[t,x] r12[t,x]-I s[x]c2 g[t] r32[t,x]-g2 r22[t,x],D[r32[t,x],t]==I c2 s[x] g[t] r33[t,x]+I dp r32[t,x]+I c1 s[x] f1[t,x] Conjugate[r13[t,x]]-I c2 s[x] g[t] r22[t,x]-g12 r32[t,x],D[r33[t,x],t]==-I c2 s[x]Conjugate[r32[t,x]] g[t]+I c2  s[x]g[t] r32[t,x]+(g2-g21) r33[t,x], D[f1[t,x],x]+ D[f1[t,x],t]==-I eta UnitStep[x+15-t] Conjugate[r12[t,x]]/c1,r11[ti,x]==s[x],r12[ti,x]==0,r13[ti,x]==0,r22[ti,x]==0,r32[ti,x]==0,r33[ti,x]==0,f1[ti,x]==Exp[-(x+4)^2]},{r11,r12,r13,r22,r32,r33,f1},{t,ti,tf},{x,-9 ,14 },MaxSteps->100000]

sol[0,0,.1,20,1,.5,.5,0,16,400]

Plot3D[Evaluate[f1[t,x] UnitStep[x+15-t]/.First[%]],{t,0,15.9},{x,-9,7},PlotPoints->80,PlotRange->All]

** light enters with almost zero reflection, is converted totally to coherence and then restored **

** *initial pulse ***

Plot[Evaluate[Re[f1[0,x]] /.First[%%]],{x,-6,-2},PlotRange->All]

** pulse in medium

Plot[Evaluate[Re[f1[7,x]] /.First[%%%]],{x,0,6},PlotRange->All]

** *shows initial compression of pulse by a factor of 2 = vg/c

Pulse exiting medium

***look at exit of pulse***

solexit[d_,dp_,c1_,c2_,g2_,g12_,g21_,ti_,tf_,eta_]:=NDSolve[{D[r11[t,x],t]==-I s[x] c1 Conjugate[r12[t,x]] Conjugate[f1[t,x]]+I c1  s[x]f1[t,x] r12[t,x]+g21 r22[t,x],D[r12[t,x],t]==I c1 s[x]Conjugate[f1[t,x] ]r11[t,x]+I d r12[t,x]-g12 r12[t,x]+I c2 g[t] s[x]r13[t,x]-I c1 s[x]Conjugate[f1[t,x] ]r22[t,x],D[r13[t,x],t]==-I c1 s[x] Conjugate[f1[t,x]]Conjugate[r32[t,x]]+I d r13[t,x]-I dp r13[t,x]+I c2 s[x]g[t] r12[t,x],D[r22[t,x],t]==I c1 s[x]Conjugate[r12[t,x]] Conjugate[f1[t,x]]+I c2 s[x]Conjugate[r32[t,x]] g[t]-I c1 s[x]f1[t,x] r12[t,x]-I s[x]c2 g[t] r32[t,x]-g2 r22[t,x],D[r32[t,x],t]==I c2 s[x] g[t] r33[t,x]+I dp r32[t,x]+I c1 s[x] f1[t,x] Conjugate[r13[t,x]]-I c2 s[x] g[t] r22[t,x]-g12 r32[t,x],D[r33[t,x],t]==-I c2 s[x]Conjugate[r32[t,x]] g[t]+I c2  s[x]g[t] r32[t,x]+(g2-g21) r33[t,x], D[f1[t,x],x]+ D[f1[t,x],t]==-I eta UnitStep[x+15-t] Conjugate[r12[t,x]]/c1,r11[ti,x]==s[x],r13[ti,x]==-c1 Exp[-4(x-5)^2] /c2,r12[ti,x]==0,r22[ti,x]==0,r32[ti,x]==0,r33[ti,x]==0,f1[ti,x]==Exp[-4(x-5)^2]},{r11,r12,r13,r22,r32,r33,f1},{t,ti,tf},{x,-9 ,19 },MaxSteps->100000]

solexit[0,0,.1,20,1,.5,.5,16,24,400]

Plot3D[Evaluate[f1[t,x] /.First[%]],{t,16,23.9},{x,0,16},PlotPoints->80,PlotRange->All]

Plot[Evaluate[Re[f1[16,x]] /.First[%%]],{x,0,10},PlotRange->All]

Plot[Evaluate[f1[23,x] /.First[%%%]],{x,2,12},PlotRange->All]

** pulse speeds up and broadens as it leaves medium ***