当前位置：首页 > news >正文

qutip——玩（3）

news 2026/4/15 1:29:28

十、单粒子激发+多粒子能量传递

前面好几节我们分别讨论了多粒子激发和能量传递，那么如果把他们两个结合在一起结果会怎么样呢。这节我们就分析一下，代码如下，结果随后

N=3   
psi = tensor(basis(2,0),basis(2,0)) # 制备初态a_list, ad_list = [], []
for i in range(N): # 用循环定义每个单格点上的产生湮灭算符op_list = [qeye(2)] * Nop_list[i] = destroy(2)a_list.append(tensor(op_list))op_list[i] = create(2)ad_list.append(tensor(op_list))hal = 0
for i in range(N-1):hal += ad_list[i] + ad_list[i+1] * a_list[i]ttotal = 20
tsteps = 500
tlist = np.linspace(0, ttotal, tsteps)  # 演化时间
res = sesolve(hal,psi,tlist) # 解薛定谔方程
states = [s * s.dag() for s in res.states] # 计算密度矩阵

v2-f17289f7a6dce6c4fa13066b09960aa8_1440w (1)

结果是非常的令人满意的，初态|000>的占据从100%快速衰落，然后中间态各领风骚，最终大趋势不变的是末态|111>，占据趋近于100%。末态持续激发。

十一、多粒子激发+多粒子能量双向传递

上一节我们看到了激发单个粒子，然后能量向右传递达到稳态的结果。如果我们使得能量可以双向传递，结果会怎么样的。代码如下，结果随后

N=3   
psi = tensor(basis(2,0),basis(2,0)) # 制备初态a_list, ad_list = [], []
for i in range(N): # 用循环定义每个单格点上的产生湮灭算符op_list = [qeye(2)] * Nop_list[i] = destroy(2)a_list.append(tensor(op_list))op_list[i] = create(2)ad_list.append(tensor(op_list))hal = 0
for i in range(N-1):hal += ad_list[i] + ad_list[i+1] * a_list[i] + a_list[i+1] * ad_list[i]ttotal = 20
tsteps = 500
tlist = np.linspace(0, ttotal, tsteps)  # 演化时间
res = sesolve(hal,psi,tlist) # 解薛定谔方程
states = [s * s.dag() for s in res.states] # 计算密度矩阵

v2-4bc50a18e34faa14a30db2a4a7f4cc80_1440w (1)

结果也挺不错的，每个态的趋势还是基本一样的。不同的就是，以|111>举例，上一节是平滑上升达到稳态，加入双向能量传递以后可以看到上升过程中会有波动，这也是可以理解的。能量到达右侧粒子后会像水一样反弹回来形成类似于水波纹一样的结果。是振荡激发的。

十二、单粒子激发+多粒子双格点能量双向传递

上节讨论了能量双向传递，如果我们给他加上不同的隧穿系数，能量传递会受到什么影响呢。代码如下，结果随后

N=3
a_list, ad_list = [], []
for i in range(N):op_list = [qeye(2)] * Nop_list[i] = destroy(2)a_list.append(tensor(op_list))op_list[i] = create(2)ad_list.append(tensor(op_list))
sigmap_list, sigmam_list = [], []
for i in range(N):op_list = [qeye(2)] * Nop_list[i] = sigmap()sigmap_list.append(tensor(op_list))op_list[i] = sigmam()sigmam_list.append(tensor(op_list))
ttotal = 40
tsteps = 500
tlist = np.linspace(0, ttotal, tsteps)
psi = tensor(basis(2,0),basis(2,0),basis(2,0))
hal = 0
t_1 = 1
t_2 = 0.1
for i in range(N-1):if i%2 == 0: hal += -t_1 * (ad_list[i+1] * a_list[i] + ad_list[i] * a_list[i+1])else:hal += -t_2 * (ad_list[i+1] * a_list[i] + ad_list[i] * a_list[i+1])hal += ad_list[0]
res = sesolve(hal,psi,tlist)
np.set_printoptions(threshold=np.inf)  # 遇到显示不全的用
states = [s * s.dag() for s in res.states]
# np.set_printoptions(linewidth=1000000)
states = np.real(np.array(states))
rounded_states = np.round(states, decimals=4)
np.set_printoptions(suppress=True)
# print(rounded_states)
psi000 = []
psi001 = []
psi010 = []
psi011 = []
psi100 = []
psi101 = []
psi110 = []
psi111 = []
for i in range(tsteps):psi000.append(rounded_states[i][0][0]) 
for i in range(tsteps):psi001.append(rounded_states[i][1][1])
for i in range(tsteps):psi010.append(rounded_states[i][2][2])
for i in range(tsteps):psi011.append(rounded_states[i][3][3]) 
for i in range(tsteps):psi100.append(rounded_states[i][4][4])
for i in range(tsteps):psi101.append(rounded_states[i][5][5])
for i in range(tsteps):psi110.append(rounded_states[i][6][6]) 
for i in range(tsteps):psi111.append(rounded_states[i][7][7])
plt.plot(tlist,psi000, color='r', label='000')
plt.plot(tlist,psi001, label='001')
plt.plot(tlist,psi010, label='010')
plt.plot(tlist,psi011, label='011')
plt.plot(tlist,psi100, label='100')
plt.plot(tlist,psi101, label='101')
plt.plot(tlist,psi110, label='110')
plt.plot(tlist,psi111, color='b', label='111')
plt.xlabel('Time')
plt.ylabel('Occuption')
text = "Evolution of states"
plt.title(text)
plt.legend()

v2-9c9da1cec5ccdbfdf2032611f49c7a24_1440w