陈数Chern number的计算（Wilson loop方法，附Python代码）

这是之前关于陈数的计算：

本篇计算方法见参考资料：https://topocondmat.org/w4_haldane/ComputingChern.html。该方法同样不需要对波函数取连续规范。

说明：本篇的代码不支持能带交叉或简并，多条能带同时计算的代码阅读这篇：陈数Chern number的计算（多条能带的Wilson loop方法，附Python代码）。

为了防止原链接失效后无法访问，这里截图如下：

下面仍然以这篇“陈数Chern number的计算（定义法，附Python/Matlab代码）”中的哈密顿量为例子，Python程序如下。

一、初始代码

首先给出一个初始代码，这里计算结果的实部是正确的，但虚部的值比较大，一个可能的原因是因为计算一个small plaquettes的Wilson loop时取点的个数太少，仅取了四个点。

补充说明：Wilson loop方法计算结果之所以存在虚部，根本原因是因为Wilson loop之后没有除以归一化系数。如果考虑了归一化系数，就没有虚部，参考这篇中的公式：陈数Chern number的计算（高效法，附Python/Matlab代码），感谢Song MR同学参与的讨论。另外，微信公众号的这篇文章也给出了原因：https://mp.weixin.qq.com/s/Ci8FCjtj7I93PE5aHmgrzA。如果考虑了归一化系数，这里的表达式和高效法的表达式是完全一致的。个人推荐加上这个归一化系数。

"""
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/18319
"""

"""
Notice that this code is not correct for calculating the Chern number because of the low precision in the calculation of Wilson loop!
"""

import numpy as np
from math import * 
import time
import cmath


def hamiltonian(kx, ky):  # 量子反常霍尔QAH模型（该参数对应的陈数为2）
    t1 = 1.0
    t2 = 1.0
    t3 = 0.5
    m = -1.0
    matrix = np.zeros((2, 2), dtype=complex)
    matrix[0, 1] = 2*t1*cos(kx)-1j*2*t1*cos(ky)
    matrix[1, 0] = 2*t1*cos(kx)+1j*2*t1*cos(ky)
    matrix[0, 0] = m+2*t3*sin(kx)+2*t3*sin(ky)+2*t2*cos(kx+ky)
    matrix[1, 1] = -(m+2*t3*sin(kx)+2*t3*sin(ky)+2*t2*cos(kx+ky))
    return matrix


def main():
    start_time = time.time()
    n = 200  # 积分密度
    delta = 2*pi/n
    chern_number = 0
    for kx in np.arange(-pi, pi, delta):
        for ky in np.arange(-pi, pi, delta):
            H = hamiltonian(kx, ky)
            eigenvalue, eigenvector = np.linalg.eig(H)
            vector = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]  # 价带波函数
            # vector = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]*cmath.exp(1j*np.random.uniform(0, pi))  # 验证规范不依赖性
           
            H_delta_kx = hamiltonian(kx+delta, ky) 
            eigenvalue, eigenvector = np.linalg.eig(H_delta_kx)
            vector_delta_kx = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]   # 略偏离kx的波函数

            H_delta_ky = hamiltonian(kx, ky+delta)  
            eigenvalue, eigenvector = np.linalg.eig(H_delta_ky)
            vector_delta_ky = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]  # 略偏离ky的波函数

            H_delta_kx_ky = hamiltonian(kx+delta, ky+delta)  
            eigenvalue, eigenvector = np.linalg.eig(H_delta_kx_ky)
            vector_delta_kx_ky = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]  # 略偏离kx和ky的波函数

            line_1 = np.dot(vector.transpose().conj(), vector_delta_kx)
            line_2 = np.dot(vector_delta_kx.transpose().conj(), vector_delta_kx_ky)
            line_3 = np.dot(vector_delta_kx_ky.transpose().conj(), vector_delta_ky)
            line_4 = np.dot(vector_delta_ky.transpose().conj(), vector)

            arg = np.log(np.dot(np.dot(np.dot(line_1, line_2), line_3), line_4))/1j
            chern_number = chern_number + arg
    chern_number = chern_number/(2*pi)
    print('Chern number = ', chern_number)
    end_time = time.time()
    print('运行时间（秒）=', end_time-start_time)


if __name__ == '__main__':
    main()

运行结果：

Chern number = (2.0000000000000013+2.6085429982998165j)

二、Wilson loop增加取点的代码

"""
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/18319
"""

import numpy as np
from math import * 
import time
import cmath


def hamiltonian(kx, ky):  # 量子反常霍尔QAH模型（该参数对应的陈数为2）
    t1 = 1.0
    t2 = 1.0
    t3 = 0.5
    m = -1.0
    matrix = np.zeros((2, 2), dtype=complex)
    matrix[0, 1] = 2*t1*cos(kx)-1j*2*t1*cos(ky)
    matrix[1, 0] = 2*t1*cos(kx)+1j*2*t1*cos(ky)
    matrix[0, 0] = m+2*t3*sin(kx)+2*t3*sin(ky)+2*t2*cos(kx+ky)
    matrix[1, 1] = -(m+2*t3*sin(kx)+2*t3*sin(ky)+2*t2*cos(kx+ky))
    return matrix


def main():
    start_time = time.time()
    n1 = 10 # small plaquettes精度
    n2 = 800 # Wilson loop精度
    delta = 2*pi/n1
    chern_number = 0
    for kx in np.arange(-pi, pi, delta):
        for ky in np.arange(-pi, pi, delta):
            vector_array = []
            # line_1
            for i2 in range(n2):
                H_delta = hamiltonian(kx+delta/n2*i2, ky) 
                eigenvalue, eigenvector = np.linalg.eig(H_delta)
                vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]
                # vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]*cmath.exp(1j*np.random.uniform(0, pi))  # 验证规范不依赖性
                vector_array.append(vector_delta)
            # line_2
            for i2 in range(n2):
                H_delta = hamiltonian(kx+delta, ky+delta/n2*i2)  
                eigenvalue, eigenvector = np.linalg.eig(H_delta)
                vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]
                vector_array.append(vector_delta)
            # line_3
            for i2 in range(n2):
                H_delta = hamiltonian(kx+delta-delta/n2*i2, ky+delta)  
                eigenvalue, eigenvector = np.linalg.eig(H_delta)
                vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]
                vector_array.append(vector_delta)
            # line_4
            for i2 in range(n2):
                H_delta = hamiltonian(kx, ky+delta-delta/n2*i2)  
                eigenvalue, eigenvector = np.linalg.eig(H_delta)
                vector_delta = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]
                vector_array.append(vector_delta)
            Wilson_loop = 1
            for i0 in range(len(vector_array)-1):
                Wilson_loop = Wilson_loop*np.dot(vector_array[i0].transpose().conj(), vector_array[i0+1])
            Wilson_loop = Wilson_loop*np.dot(vector_array[len(vector_array)-1].transpose().conj(), vector_array[0])
            arg = np.log(Wilson_loop)/1j
            chern_number = chern_number + arg
    chern_number = chern_number/(2*pi)
    print('Chern number = ', chern_number)
    end_time = time.time()
    print('运行时间（秒）=', end_time-start_time)


if __name__ == '__main__':
    main()

运行结果：

Chern number = (1.9999999999999987+0.0032600066414606387j)

随着Wilson loop计算精度的提高，虚部将趋于0。

其他参考资料：

[1] 数学上陈数（Chern number）或 Berry Phase 有何意义？

[2] 贝里相位、贝里联络、贝里曲率和陈数

[3] B. Andrei Bernevig - Topological Insulators and Topological Superconductors 第30页

7,682 次浏览

【说明：本站主要是个人的一些笔记和代码分享，内容可能会不定期修改。为了使全网显示的始终是最新版本，这里的文章未经同意请勿转载。引用请注明出处：https://www.guanjihuan.com】

25 thoughts on “陈数Chern number的计算（Wilson loop方法，附Python代码）”

您好，如果除以了归一化系数消掉了虚部，那一个wilson loop取多少个点啊？取4个点够吗？同样的哈密顿量同样的delta，是高效法更快还是wilson loop方法更快啊？（因为我要算一个几百×几百的哈密顿量，所以必须要考虑快慢），谢谢

guanjihuan说道：

2024年8月20日 14:06

四个点够的。如果有归一化系数，wilson loop方法取四个点时的计算公式和高效法的计算公式好像是完全一样的，我没看出区别。推荐直接使用高效法。

回复

您好，老师您这个，光子体系能使用吗？

guanjihuan说道：

2024年4月24日 13:59

如果有波函数和布里渊区，那么就可以按公式进行计算。

回复

关老师您好，请问如果不能写出系统的哈密顿量，只能通过数值模拟计算出能带图，那我们能用这个方法去计算某一条能带的陈数吗

guanjihuan说道：

2024年4月16日 20:04

如果有波函数和布里渊区，就可以算，像Vasp计算就可以得到波函数的数值信息。另外，也有一些软件可以计算不变量，例如wanniertool之类的。

回复

关老师，想问一下，在一个系统中，陈数一定是成对存在的吗？例如，一个系统是否存在陈数为0、-1、1和2的情况呢

guanjihuan说道：

2024年3月14日 21:34

不一定是成对存在的，但所有能带的总陈数应该为零，即所有能带的陈数之和为零，这是因为真空能级是拓扑平庸的。

回复

您好，请问一下这种方法中将φ分为n块后，φn中有一个i，exp(iφn)中也有一个i，把φn带入时为什么没有出现负号呢？

guanjihuan说道：

2022年12月14日 15:30

嗯，好像是少了一个负号，应该是笔误。可以考虑为exp(-iφn)，最后公式会多一个负号，但这对是否拓扑没太大影响。可以参考这篇文献中的公式：Band topology in classical waves: Wilson-loop approach to topological numbers and fragile topology。

回复

你好，请问代码里面的量子反常霍尔QAH模型是指具体哪个模型呀

guanjihuan说道：

2022年12月5日 04:48

参考这篇文章以及其参考文献：Route towards Localization for Quantum Anomalous Hall Systems with Chern Number 2。

回复

博主您好。这里的虚部似乎不是个问题，公式中使用的angle（复数的幅角），而代码实现中是通过ln()/i，这样操作得到的实部的确是幅角，但虚部是公式所不要的。

guanjihuan说道：

2022年6月9日 23:27

嗯，最终只要拿到实部就行了。但严格意义上，虚部应该要趋于零，结果才算是正确的，我个人感觉这样用起来会更放心一点。如果有数学上证明说，精度不高时虚部的存在完全不影响实部的值，那么才敢直接把虚部丢掉。

回复

请问计算Wilson loop 的时候，最后为什么要进行这一行操作？
Wilson_loop = ilson_loop*np.dot(vector_array[len(vector_array)-1].transpose().conj(), vector_array[0])。 end的波函数再与begin的波函数dot？

在计算SSH Wilson loop的时候，让begin的波函数与end的波函数相等，并没有这一步操作。

guanjihuan说道：

2022年4月24日 22:34

Wilson loop就是要首尾相连的。在计算SSH的Wilson loop中，也是有首尾相连，对于布里渊区来说，k=-pi和k=pi实际上是同一点，和end的波函数内积，就相当于和begin的波函数内积。为了在数值上避开波函数的规范问题，所以令end的波函数和begin的波函数相等。

回复
1. 匿名说道：
  
  2022年8月25日 10:02
  
  您好像在上面那个不对的代码里没有让line4再与line1做内积，loop没有闭合。原因是在于，首尾两个波函数取为一致，内积为1所以省略掉了吗。
  
  回复
  1. guanjihuan说道：
    
    2022年8月25日 11:30
    
    已经闭合了呀。line1和line4已经是内积的结果了，首尾都是vector。
    
    回复

博主您好，我想问一下，我修改plaquettes精度会影响陈数的值，差几倍的那种，这有可能是什么原因呢？谢谢~

guanjihuan说道：

2022年3月20日 03:11

尽量还是取大一点吧，只要数值能收敛就没有太大问题。但如果增大精度数值一直不收敛，那计算结果可能就不对了。

回复

老师您好，我想请问一下如何计算实陈数呢

guanjihuan说道：

2022年1月1日 13:16

参考前面几篇的方法。这里的方法，我还没明白为什么结果的虚部会这么大。按理来说，应该为无穷小。

回复
1. 匿名说道：
  
  2022年1月1日 14:03
  
  实陈数的化我们是否是在berry curvature那里加一个-i*σy，然后求个trace呢，这样操作后感觉有问题
  
  回复
  1. guanjihuan说道：
    
    2022年1月1日 15:07
    
    没看懂为什么要加。
    
    回复
  2. guanjihuan说道：
    
    2022年1月4日 16:29
    
    找到问题的原因了，是因为计算Wilson loop时，取的点太少了。已更新最新的代码。
    
    回复

陈数Chern number的计算（Wilson loop方法，附Python代码）

一、初始代码

二、Wilson loop增加取点的代码

Published by guanjihuan

25 thoughts on “陈数Chern number的计算（Wilson loop方法，附Python代码）”

发表评论取消回复

一、初始代码

二、Wilson loop增加取点的代码

Published by guanjihuan

25 thoughts on “陈数Chern number的计算（Wilson loop方法，附Python代码）”

发表评论 取消回复

发表评论取消回复