liaodeyang.py

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
import math


def Hamiltonian3(k):
    # 元胞内的两个C原子坐标
    pos1 = np.array([0, 0])
    pos2 = np.array([0, 1])
    pos = [pos1, pos2]
    # 取晶格基矢
    a1 = np.array([math.sqrt(3), 0])
    a2 = np.array([math.sqrt(3) / 2, 3 / 2])

    # 新建一个名叫 fig的画图窗口
    fig = plt.figure()
    # 准备绘制三维图形
    ax = Axes3D(fig)
    # k点的取值范围
    Kx = np.arange(-math.pi, math.pi, 0.01)
    Ky = np.arange(-math.pi, math.pi, 0.01)
    n = len(Kx)
    E1 = np.mat(np.zeros((n, n)))
    t = -1
    H = np.mat(np.zeros((2, 2)), dtype = complex)
    eps = 0.01
    for i in range(0, 2):
        for j in range(0, 2):
            for m in range(-1, 2):
                for n in range(-1, 2):
                    d = pos[i] + n * a1 + m * a2 - pos[j]
                    if abs(np.linalg.norm(d) - 1) < eps:
                        H[i, j] = H[i, j] + t * math.e ** (1j * np.dot(k, d))
    return H

n=10

for i in range(n):
    for j in range(n):
        # 取k点
        k = np.array([Kx[i], Ky[j]])
        # 计算Hamilton
        H = Hamiltonian3(k)
        # 取本征值
        Es = np.linalg.eig(H)[0]
        # 存储本征值
        E1[i, j] = Es[0]
E2 = -E1
# 画图
Axes3D.ax.plot_surface(Kx, Ky, E1, rstride = 1, cstride = 1, cmap = 'rainbow')
Axes3D.ax.plot_surface(Kx, Ky, E2, rstride = 1, cstride = 1, cmap = 'rainbow')
plt.show()