Comparision steady state and mesolve at infinity

[Hassan Bassereh]

Hello every body

I want to solve master equation in attached-PDF file and calculate P_sink=<n3 | rho(t) | n3> at infinity (with mesolve function) and compare it with steadystate function. Generally, the number in both methods should be equal but sth else is obtained. Can anybody guide me?

All the best,
Hassan


"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""

from qutip import *
from numpy import *
from pylab import *

ee=basis(2,0);g=basis(2,1)

n0=tensor(g,g,g)
n1=tensor(ee,g,g)
n2=tensor(g,ee,g)
n3=tensor(g,g,ee)

H=n1*n1.dag()+n2*n2.dag()+n2*n1.dag()+n1*n2.dag()#Hamiltonian
psi0=n1*n1.dag()#Initial state

gs=0.1#gamma_sink
gd=0.2#gamma_dissipation

sigm1=n0*n1.dag()
sigm2=n0*n2.dag()
sigm32=n3*n2.dag()

c_op_list=[]
c_op_list.append(sqrt(gd)*sigm1)
c_op_list.append(sqrt(gd)*sigm2)
c_op_list.append(sqrt(gd)*sigm32)

L=liouvillian(H,c_op_list)
ss=steadystate(L, method='eigen', use_rcm=False)
ss3=expect(n3*n3.dag(),ss)
print ss3

tlist = linspace(0, 40, 2000)
output=mesolve(H,psi0,tlist,c_op_list,[n3*n3.dag()])
P_sink=output.expect[0]
print P_sink[-1]


""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""

[Hassan Bassereh]

Nobody knows?????