c33ae98Companion to 'Coding closed and open quantum systems in MATLAB' (arXiv:1911.04906)Jeremy Magland 1Nsim = 25; % Number of simulations
2L = 2; % Number of cavities
3wc = 1; % Cavity frequency
4g = 1e-2*wc; % Atom-light coupling
5J = 1e-4*wc; % Coupling between cavities
6Nph = 2; % Number of photons per cavity
7dimFock = Nph+1; % Dimension Fock space for photons
8dimT = 2*dimFock; % Dimension atom+cavity system
9Deltai = 10^(-2)*g; % Initial detuning
10Deltaf = 10^(+2)*g; % Final detuning
11xi = log10(Deltai/g); % Initial detuning in log scale
12xf = log10(Deltaf/g); % Final detuning in log scale
13dx = (xf-xi)/(Nsim-1); % Step dx
14x = xi:dx:xf; % Vector x to plot transition phase
15OP_JC = zeros(size(x)); % Order parameter Jaynes-Cummings model
16OP_R = zeros(size(x)); % Order parameter Rabi model
17A = cell(1,L); % Cell array to storage \hat{a}_i operators
18Sp = cell(1,L); % Cell array to storage \hat{\sigma}_i^+ operators
19N_ex = cell(1,L); % Cell array to storage N_i operators
20Iatom = eye(2); % Identity matrix atom system
21Icav = eye(dimFock ); % Identity matrix cavity system
22Is = eye(2*dimFock); % Identity matrix atom+cavity system
23for i=1:L
24 A{i} = acav(i,L,Nph,Is,Iatom); % Photon operator \hat{a}_i
25 Sp{i} = sigmap(i,L,Is,Icav); % Atom operator \hat{\sigma}_i
26 N_ex{i} = A{i}'*A{i}+Sp{i}*Sp{i}'; % number of excitations per cavity
27end
28Ad = ones(L); % Adyacent matrix
29Ad = triu(Ad)-eye(L); % A(i,j)=1 for j>i
30Hhopp = zeros(dimT^L,dimT^L); % Interaction Hamiltonian
31for i=1:L
32 for j=1:L
33 Hhopp = Hhopp - J*Ad(i,j)*A{i}'*A{j} - J*Ad(i,j)*A{i}*A{j}';
34 end
35end
36Nt = 10000; % Length time vector
37ti = 0.01/J; % Initial time
38tf = 1/J; % Final time
39dt = (tf-ti)/(Nt-1); % Step time dt
40t = ti:dt:tf; % Time vector
41Lambda_R = zeros(Nsim,Nt);
42Lambda_JC = zeros(Nsim,Nt);
43parfor n=1:Nsim
44 D = g*10^(x(n)); % Detuning in each simulation
45 Model = 'Rabi';
46 [OP_R(n), P1m_R] = QuantumSimulationCavityArray(wc,D,g,Nph,L,A,Sp,N_ex,Hhopp,t,Model);
47 Model = 'Jaynes-Cummings';
48 [OP_JC(n),P1m_JC] = QuantumSimulationCavityArray(wc,D,g,Nph,L,A,Sp,N_ex,Hhopp,t,Model);
49 Lambda_R(n,:) = -1/L*log2(P1m_R); % Rate function for the Rabi model
50 Lambda_JC(n,:) = -1/L*log2(P1m_JC); % Rate function for the Jaynes-Cummings model
51end
52figure()
53box on
54hold on
55plot(J*t,Lambda_R(end,:),'b-','Linewidth',3)
56plot(J*t,Lambda_JC(end,:),'r-','Linewidth',3)
57hold off
58xlabel('$Jt$','Interpreter','LaTex','Fontsize', 30)
59ylabel('$\Lambda(t)$','Interpreter','LaTex','Fontsize', 30)
60set(gca,'fontsize',21)
61legend({'$\mbox{RH}$','$\mbox{JCH}$'},'Interpreter','latex','Fontsize', 21,'Location','best')
63figure()
64box on
65hold on
66plot(x,real(OP_R),'.b','Markersize',30)
67plot(x,real(OP_JC),'.r','Markersize',30)
68hold off
69xlabel('$\mbox{Log}_{10}(\Delta/g)$','Interpreter','LaTex','Fontsize', 30)
70ylabel('$\mbox{OP}$','Interpreter','LaTex','Fontsize', 30)
71set(gca,'fontsize',21)
72legend({'$\mbox{RH}$','$\mbox{JCH}$'},'Interpreter','latex','Fontsize', 21,'Location','best')