这是之前的一篇:陈数Chern number的计算(Wilson loop方法,附Python代码),仅仅是计算某条能带的陈数,不支持能带交叉或简并。
本篇给出支持能带交叉或简并的代码,即多带同时计算Wilson loop,并求行列式。公式和之前那篇博文一样,为这个链接:https://topocondmat.org/w4_haldane/ComputingChern.html。
为了验证算法和代码的正确性,这里以这篇文章为例子:方格子紧束缚模型中朗道能级的陈数/霍尔电导(附Python代码)。当Ny为偶数时,存在能带交叉的情况,这时候需要两条带同时计算才能给出正确的陈数结果。之前那篇已经给出计算能带的代码,这里代码略,本篇底部有直接贴上对应的能带图。
代码如下:
"""
This code is supported by the website: https://www.guanjihuan.com
The newest version of this code is on the web page: https://www.guanjihuan.com/archives/23989
"""
import numpy as np
import math
from math import *
import cmath
import functools
def hamiltonian(kx, ky, Ny, B):
h00 = np.zeros((Ny, Ny), dtype=complex)
h01 = np.zeros((Ny, Ny), dtype=complex)
t = 1
for iy in range(Ny-1):
h00[iy, iy+1] = t
h00[iy+1, iy] = t
h00[Ny-1, 0] = t*cmath.exp(1j*ky)
h00[0, Ny-1] = t*cmath.exp(-1j*ky)
for iy in range(Ny):
h01[iy, iy] = t*cmath.exp(-2*np.pi*1j*B*iy)
matrix = h00 + h01*cmath.exp(1j*kx) + h01.transpose().conj()*cmath.exp(-1j*kx)
return matrix
def main():
Ny = 20
H_k = functools.partial(hamiltonian, Ny=Ny, B=1/Ny)
chern_number = calculate_chern_number_for_square_lattice_with_wilson_loop_for_degenerate_case(H_k, index_of_bands=range(int(Ny/2)-1), precision_of_wilson_loop=5)
print('价带:', chern_number)
print()
chern_number = calculate_chern_number_for_square_lattice_with_wilson_loop_for_degenerate_case(H_k, index_of_bands=range(int(Ny/2)+2), precision_of_wilson_loop=5)
print('价带(包含两个交叉能带):', chern_number)
print()
chern_number = calculate_chern_number_for_square_lattice_with_wilson_loop_for_degenerate_case(H_k, index_of_bands=range(Ny), precision_of_wilson_loop=5)
print('所有能带:', chern_number)
# # 函数可通过Guan软件包调用。安装方法:pip install --upgrade guan
# import guan
# chern_number = guan.calculate_chern_number_for_square_lattice_with_wilson_loop_for_degenerate_case(hamiltonian_function, index_of_bands=[0, 1], precision_of_plaquettes=20, precision_of_wilson_loop=5, print_show=0)
def calculate_chern_number_for_square_lattice_with_wilson_loop_for_degenerate_case(hamiltonian_function, index_of_bands=[0, 1], precision_of_plaquettes=20, precision_of_wilson_loop=5, print_show=0):
delta = 2*math.pi/precision_of_plaquettes
chern_number = 0
for kx in np.arange(-math.pi, math.pi, delta):
if print_show == 1:
print(kx)
for ky in np.arange(-math.pi, math.pi, delta):
vector_array = []
# line_1
for i0 in range(precision_of_wilson_loop):
H_delta = hamiltonian_function(kx+delta/precision_of_wilson_loop*i0, ky)
eigenvalue, eigenvector = np.linalg.eig(H_delta)
vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))]
vector_array.append(vector_delta)
# line_2
for i0 in range(precision_of_wilson_loop):
H_delta = hamiltonian_function(kx+delta, ky+delta/precision_of_wilson_loop*i0)
eigenvalue, eigenvector = np.linalg.eig(H_delta)
vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))]
vector_array.append(vector_delta)
# line_3
for i0 in range(precision_of_wilson_loop):
H_delta = hamiltonian_function(kx+delta-delta/precision_of_wilson_loop*i0, ky+delta)
eigenvalue, eigenvector = np.linalg.eig(H_delta)
vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))]
vector_array.append(vector_delta)
# line_4
for i0 in range(precision_of_wilson_loop):
H_delta = hamiltonian_function(kx, ky+delta-delta/precision_of_wilson_loop*i0)
eigenvalue, eigenvector = np.linalg.eig(H_delta)
vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))]
vector_array.append(vector_delta)
wilson_loop = 1
dim = len(index_of_bands)
for i0 in range(len(vector_array)-1):
dot_matrix = np.zeros((dim , dim), dtype=complex)
i01 = 0
for dim1 in index_of_bands:
i02 = 0
for dim2 in index_of_bands:
dot_matrix[i01, i02] = np.dot(vector_array[i0][:, dim1].transpose().conj(), vector_array[i0+1][:, dim2])
i02 += 1
i01 += 1
det_value = np.linalg.det(dot_matrix)
wilson_loop = wilson_loop*det_value
dot_matrix_plus = np.zeros((dim , dim), dtype=complex)
i01 = 0
for dim1 in index_of_bands:
i02 = 0
for dim2 in index_of_bands:
dot_matrix_plus[i01, i02] = np.dot(vector_array[len(vector_array)-1][:, dim1].transpose().conj(), vector_array[0][:, dim2])
i02 += 1
i01 += 1
det_value = np.linalg.det(dot_matrix_plus)
wilson_loop = wilson_loop*det_value
arg = np.log(wilson_loop)/1j
chern_number = chern_number + arg
chern_number = chern_number/(2*math.pi)
return chern_number
if __name__ == '__main__':
main()
运算结果:
价带: (-8.999999999999996+22.620063283010385j)
价带(包含两个交叉能带): (8.000000000000004+22.225601922823458j)
所有能带: (8.988228160493337e-16+4.549167968104345j)
当增加precision_of_wilson_loop值时,虚数会变小。这里似乎对实部结果没影响。
另外说明:之所以存在虚部,很大可能是因为Wilson loop之后没有除以归一化系数。如果考虑了归一化系数,那么虚部就可以特别小,参考高效法中的公式:陈数Chern number的计算(高效法,附Python/Matlab代码)、陈数Chern number的计算(多条能带的高效法,附Python代码)。
以上计算对应的能带图:

【说明:本站主要是个人的一些笔记和代码分享,内容可能会不定期修改。为了使全网显示的始终是最新版本,这里的文章未经同意请勿转载。引用请注明出处:https://www.guanjihuan.com】
您好! 请问您程序里的index_of_bands如果小于能带总数时,是计算的前s条能带的陈数总和吗?似乎这样s无法遍历与所有能带求内积,是不需要吗?
index_of_bands中可以取任意带,数量不限。如果是算体系的陈数,一般是计算费米能以下的价带的陈数总和。不需要算所有能带,所有能带的陈数一定为零。
你好!我想请教您,这代码计算不是交叉能带的话,结果还是正确的吗?
也是正确的,取一条就可以了。可以用原来单条能带的计算方法作为验证。