1% Atom-field interaction: vacuum Rabi oscillations (Sec. V)
2% Jaynes-Cummings Hamiltonian; the atom starts excited with n photons in the
3% cavity and periodically exchanges a photon with the field.
4% Pe(t) = cos^2(g t sqrt(n+1)), Pg(t) = sin^2(g t sqrt(n+1))
6clear; % Clear memory
7clc; % Clear the command window/screen
8d = 10; % dimension of the cavity field
9hbar = 1;
10W0 = 1; % atomic frequency
11Wf = 1; % cavity field frequency
12g = 0.1; % coupling constant
13A = diag(sqrt(1:d-1), 1); % Annihilation operator
14Ad = A'; % Creation operator
15Sz = [1,0;0,-1]; % sigma_z
16Splus = [0,1;0,0]; % sigma_+
17Sminus = [0,0;1,0]; % sigma_-
18gs = [0;1]; % ground state
19es = [1;0]; % excited state
20I_a = eye(2); % identity operator for the atom
21I_f = eye(d); % identity operator for the field
22Hatom = (1/2)*hbar*W0*kron(Sz,I_f); % atomic Hamiltonian
23Hfield = hbar*Wf*kron(I_a,Ad*A); % field Hamiltonian
24Hint = hbar*g*(kron(Splus,A)+kron(Sminus,Ad)); % interaction Hamiltonian
25H = Hatom + Hfield + Hint; % JC Hamiltonian
26n = 4; % initial number of photons in the cavity
27en = kron(es,I_f(:,n+1)); % atom excited, field has n photons |e,n>
28gn = kron(gs,I_f(:,n+2)); % atom ground, field has n+1 photons |g,n+1>
29Psi = en; % initial state
30dt = 0.1; % time step
31U = expm(-i*H*dt/hbar); % unitary time-evolution operator
32T = 0:dt:30; % total evolution time
33for t = 1:length(T)
34 Pe(t) = norm(en'*Psi)^2; % probability of |e,n>
35 Pg(t) = norm(gn'*Psi)^2; % probability of |g,n+1>
36 Psi = U*Psi; % time-evolved state
37 Psi = Psi/norm(Psi); % renormalize
38end
39plot(T, Pe, 'r', T, Pg, 'k') % Pe in red, Pg in black
40title('Vacuum Rabi oscillations (n = 4, g = 0.1)')
41xlabel('t'); ylabel('Probability'); legend('P_e', 'P_g')