1
# 查看声子谱
2
import yaml
3

4

5
def read_band_yaml(filename):
6
    """读取 band.yaml 文件"""
7
    with open(filename, "r") as file:
8
        data = yaml.safe_load(file)
9

10
    # 提取数据
11
    q_points = []
12
    frequencies = []
13
    for phonon in data["phonon"]:
14
        q_points.append(phonon["q-position"])
15
        freqs = [band["frequency"] for band in phonon["band"]]
16
        frequencies.append(freqs)
17

18
    return np.array(q_points), np.array(frequencies)
19

20
def calculate_path_distances(q_points):
21
    """计算路径距离"""
22
    distances = [0]
23
    for i in range(1, len(q_points)):
24
        delta = np.linalg.norm(q_points[i] - q_points[i-1])
25
        distances.append(distances[-1] + delta)
26
    return np.array(distances)/10
27

28
def calculate_label_positions(config_file):
29
    with open(config_file, 'r') as f:
30
        lines = f.readlines()
31

32
    # 提取 BAND 行
33
    band_line = None
34
    band_labels = None
35
    for line in lines:
36
        if 'BAND_LABELS' in line.split():
37
            band_labels = line.split('=')[1].strip().split()
38
        elif 'BAND' in line.split():
39
            band_line = line.split('=')[1].strip()
40

41
    if not band_line or not band_labels:
42
        return None, None
43
    # 解析坐标
44
    coords = np.array([float(x) for x in band_line.split()]).reshape(-1, 3)
45

46
    # 计算累积距离
47
    distances = [0]
48
    for i in range(1, len(coords)):
49
        delta = np.linalg.norm(coords[i] - coords[i-1])
50
        distances.append(distances[-1] + delta)
51

52
    return distances, band_labels

1
fig,ax= plt.subplots(figsize=(cm(7), cm(5)),dpi=500)
2

3
for i in range(frequencies.shape[1]):
4
    plt.plot(x, frequencies[:, i], color='#1763bb', linewidth=0.5)
5
x_max= np.max(x)
6
distabces=distabces/distabces[-1]*x_max
7
set_ticks(ax,xrange=(distabces[0], distabces[-1], 0.05), yrange=(-1, 13, 2))
8
ax.set_yticks(np.arange(0,12.1, 2))
9
for label, point in zip(band_labels, distabces):
10
    ax.axvline(point, color='dimgray', linewidth=0.5, linestyle='--')
11
    ax.text(point, -1.2, label, ha='center', va='top', fontsize=8)
12
ax.axhline(0, color='dimgray', linewidth=0.5, linestyle='--')
13
# 不显示x轴坐标和刻度
14
ax.set_xticks([])
15
ax.set_xticklabels([])
16
ax.set_ylabel("Frequency (THz)", fontsize=8)
17
half_plot_set(ax)

1
# 计算声子态密度
2
# 读取mesh.yaml文件来计算DOS
3
mesh_data = yaml.safe_load(open(phonopy_dir / "mesh.yaml"))
4

5
# 提取频率数据
6
frequencies_mesh = []
7
weights = []
8
for phonon in mesh_data["phonon"]:
9
    freqs = [band["frequency"] for band in phonon["band"]]
10
    frequencies_mesh.append(freqs)
11
    weights.append(phonon["weight"])
12

13
frequencies_mesh = np.array(frequencies_mesh)
14
weights = np.array(weights)
15

16
# 计算DOS
17
freq_range = np.linspace(-1, 13, 1000)
18
dos = np.zeros_like(freq_range)
19

20
# 使用高斯展宽计算DOS
21
sigma = 0.1  # 展宽参数
22
for i, freq_point in enumerate(freq_range):
23
    dos_set=weights * np.exp(-((freq_point - frequencies_mesh) ** 2) / (2 * sigma**2)).T
24
    dos[i] = np.sum(dos_set)
25

26
# 归一化
27
dos = dos / (sigma * np.sqrt(2 * np.pi) * np.sum(weights))

1
fig, axs = plt.subplots(
2
    1, 2, figsize=(cm(7), cm(5)), dpi=500, gridspec_kw={"width_ratios": [3.5, 1]}
3
)
4
ax = axs[0]
5

6
# 绘制所有声子支
7
for i in range(frequencies.shape[1]):
8
    ax.plot(x, frequencies[:, i], color="#1763bb", linewidth=0.5)
9
x_max = np.max(x)
10
distabces = distabces / distabces[-1] * x_max
11
set_ticks(ax, xrange=(distabces[0], distabces[-1], 0.05), yrange=(-1, 13, 2))
12
ax.set_yticks(np.arange(0, 12.1, 2))
13
i = 0
14
for label, point in zip(band_labels, distabces):
15
    ax.axvline(point, color="dimgray", linewidth=0.5, linestyle="--")
16
    if i == 5:
17
        ax.text(point - 0.005, -1.2, label, ha="center", va="top", fontsize=8)
18
    elif i == 6:
19
        ax.text(point + 0.005, -1.2, label, ha="center", va="top", fontsize=8)
20
    else:
21
        ax.text(point, -1.2, label, ha="center", va="top", fontsize=8)
22
    i += 1
23
ax.axhline(0, color="dimgray", linewidth=0.5, linestyle="--")
24
# 不显示x轴坐标和刻度
25
ax.set_xticks([])
26
ax.set_xticklabels([])
27
ax.set_ylabel("Frequency (THz)", fontsize=10)
28
half_plot_set(ax)
29

30
ax = axs[1]
31
# 绘制DOS
32
ax.plot(dos, freq_range, color=color_map[0], linewidth=0.5)
33
ax.fill_betweenx(freq_range, dos, alpha=0.3, color=color_map[0])
34
set_ticks(ax, xrange=(0, 9, 2), yrange=(-1, 13, 2))
35
ax.axhline(0, color="dimgray", linewidth=0.5, linestyle="--")
36
ax.set_yticklabels([])
37
ax.set_xticks([])
38
ax.set_xticklabels([])
39
half_plot_set(ax)
40
ax.text(4.5, -1.2, "DOS", ha="center", va="top", fontsize=8)
41

42
fig.subplots_adjust(wspace=0)  # 调整子图间距