Update Wilson_loop_in_SSH_model.py

This commit is contained in:
guanjihuan 2024-01-25 23:30:24 +08:00
parent 867d15ae26
commit 9ea3d0a41e

View File

@ -7,7 +7,6 @@ import numpy as np
import cmath import cmath
from math import * from math import *
def hamiltonian(k): # SSH模型哈密顿量 def hamiltonian(k): # SSH模型哈密顿量
gamma = 0.5 gamma = 0.5
lambda0 = 1 lambda0 = 1
@ -18,15 +17,14 @@ def hamiltonian(k): # SSH模型哈密顿量
h[1,0] = gamma+lambda0*cmath.exp(1j*k) h[1,0] = gamma+lambda0*cmath.exp(1j*k)
return h return h
def main(): def main():
Num_k = 100 Num_k = 100
k_array = np.linspace(-pi, pi, Num_k) k_array = np.linspace(-pi, pi, Num_k)
vector_array = [] vector_array = []
for k in k_array: for k in k_array:
vector = get_occupied_bands_vectors(k, hamiltonian) vector = get_occupied_bands_vectors(k, hamiltonian)
vector_array.append(vector) # vector_array.append(vector)
# vector_array.append(vector*cmath.exp(1j*np.random.uniform(0, pi))) # 随机相位测试 vector_array.append(vector*cmath.exp(1j*np.random.uniform(0, pi))) # 随机相位测试
# 波函数固定一个规范 # 波函数固定一个规范
vector_sum = 0 vector_sum = 0
@ -34,11 +32,9 @@ def main():
vector_sum += np.abs(vector_array[i0]) vector_sum += np.abs(vector_array[i0])
index = np.argmax(np.abs(vector_sum)) index = np.argmax(np.abs(vector_sum))
for i0 in range(Num_k): for i0 in range(Num_k):
vector_array[i0] = find_vector_with_fixed_gauge_by_making_one_component_real(vector_array[i0], index=index) angle = cmath.phase(vector_array[i0][index])
vector_array[i0] = vector_array[i0]*cmath.exp(-1j*angle)
# # 波函数固定一个规范 # vector_array[i0] = find_vector_with_fixed_gauge_by_making_one_component_real(vector_array[i0], index=index) # 或者使用这种寻找查找方法,计算效率比较低
# import guan
# vector_array = guan.find_vector_array_with_fixed_gauge_by_making_one_component_real(vector_array, precision=0.005)
# 计算Wilson loop # 计算Wilson loop
W_k = 1 W_k = 1
@ -48,30 +44,27 @@ def main():
nu = np.log(W_k)/2/pi/1j nu = np.log(W_k)/2/pi/1j
# if np.real(nu) < 0: # if np.real(nu) < 0:
# nu += 1 # nu += 1
print('p=', nu, '\n') print('p=', nu)
def get_occupied_bands_vectors(x, matrix): def get_occupied_bands_vectors(x, matrix):
matrix0 = matrix(x) matrix0 = matrix(x)
eigenvalue, eigenvector = np.linalg.eig(matrix0) eigenvalue, eigenvector = np.linalg.eig(matrix0)
vector = eigenvector[:, np.argsort(np.real(eigenvalue))[0]] vector = eigenvector[:, np.argsort(np.real(eigenvalue))[0]]
return vector return vector
# def find_vector_with_fixed_gauge_by_making_one_component_real(vector, precision=0.005, index=None):
def find_vector_with_fixed_gauge_by_making_one_component_real(vector, precision=0.005, index=None): # if index == None:
if index == None: # index = np.argmax(np.abs(vector))
index = np.argmax(np.abs(vector)) # sign_pre = np.sign(np.imag(vector[index]))
sign_pre = np.sign(np.imag(vector[index])) # for phase in np.arange(0, 2*np.pi, precision):
for phase in np.arange(0, 2*np.pi, precision): # sign = np.sign(np.imag(vector[index]*cmath.exp(1j*phase)))
sign = np.sign(np.imag(vector[index]*cmath.exp(1j*phase))) # if np.abs(np.imag(vector[index]*cmath.exp(1j*phase))) < 1e-9 or sign == -sign_pre:
if np.abs(np.imag(vector[index]*cmath.exp(1j*phase))) < 1e-9 or sign == -sign_pre: # break
break # sign_pre = sign
sign_pre = sign # vector = vector*cmath.exp(1j*phase)
vector = vector*cmath.exp(1j*phase) # if np.real(vector[index]) < 0:
if np.real(vector[index]) < 0: # vector = -vector
vector = -vector # return vector
return vector
if __name__ == '__main__': if __name__ == '__main__':
main() main()