Z2pack 计算chern数

Jan 12, 2022 4 min read z2packpython

z2pack 是一个强大的计算Berryphase相关拓扑数的python工具包，由欧洲老哥 （我们组内对他的趣称，从他的代码里学到不少，现已入职微软)

Z2pack

Z2Pack is a tool for calculating topological invariants. The method is based on tracking the evolution of hybrid Wannier functions, which is equivalent to the computation of the Wilson loop. Originally developed for calculating invariants, it is now also capable of calculating Chern numbers. Moreover, through the use of individual Chern numbers it can be used to identify any kind of topological phase.

A key feature of Z2Pack is its flexibility: It can be used with various kinds of calculations and models. Currently, interfaces to first-principles codes, tight-binding models and models are implemented. However, it can in principle be extended to any kind of system.

关于安装与使用的一些说明

  • z2pack 与第一性原理做接口时需要使用 wannier90, 对VASP来说 wannier90.2.1.0 效果较好。
  • z2pack 计算 TB 模型时，借助另外一个库：tbmodels （也是由欧洲老哥开发，tbmodels可以方便的自己定义TB模型或者从wannier转入，tbmodels可以快速方便地对wannier作对称化，详见站内链接）
  • z2pack 计算 kp 模型，直接在代码里把kp模型写出来即可。

从实践的角度出发

一般而言，只要我们熟悉了 z2pack 计算 TB模型（wannier90_hr.dat）的流程，我们就可以计算无论是模型或者第一性原理计算的结果了。

代码拆解

import module

 1import os
 2import sys
 3import pathlib
 4import symmetry_representation as sr
 5import numpy as np
 6import matplotlib.pyplot as plt
 7import pymatgen as mg
 8import pymatgen.symmetry.analyzer
 9import pymatgen.symmetry.bandstructure
10# import logging
11import itertools
12import z2pack
13import tbmodels

handy function

 1def PrintDict(dict):
 2    # for i in dict.keys():
 3    #     print(i)  
 4    #     print(dict[i])  
 5    # for i in dict.values():
 6    #     print(i)
 7    for key, value in dict.items():
 8        print(key + ": " )
 9        for ivalue in value:
10            print("      "+str(ivalue))

导入 wannier90_hr.dat

 1# load first
 2
 3model = tbmodels.Model.from_wannier_files(
 4    hr_file='wannier90_hr.dat',
 5   # win_file='wannier90.win',
 6   # wsvec_file='wannier90_wsvec.dat',
 7    #xyz_file='wannier90_centres.xyz'
 8)
 9OccupyBand = 3
10SelectBands = OccupyBand
11model.occ = OccupyBand

读取wannier90_hr.dat 时可以根据情况选择性地加入:wannier90.win,wannier90_wsvec.dat, wannier90_centres.xyz等信息

从wannier90_orbital.dat 导入 pos

如果上一步导入了 xyz_file，则无需导入pos信息 如果是自己定义的模型或者用symmetry-adapter wannier 则这里借助 wannier_orbital.dat 导入晶格pos信息 如果不引入pos，则默认所有轨道在原点

wannier90_orbital.dat

10 #这里 0表示spinless 1表示spinful (spinful-会按照wannier90 的轨道顺序，即一个轨道自旋上紧跟着自旋下)
26 #这里表示总轨道数
3s,0.300000012  ,       0.000000000  ,       0.000000000
4s,0.000000000  ,       0.300000012  ,       0.000000000
5s,0.699999988  ,       0.699999988  ,       0.000000000
6s,0.699999988  ,       0.000000000  ,       0.000000000
7s,0.000000000  ,       0.699999988  ,       0.000000000
8s,0.300000012  ,       0.300000012  ,       0.000000000
  • 后面紧跟着 每个轨道位置第一列表示轨道的函数形式，z2pack 用不上，可以全部写s轨道
  • 逗号是必须的
  • 轨道的含义：
    • For the first col, we give a reference here
    • s:
    •    1
    • p:
    •    z

    •    x

    •    y
    • d:
    •    z**2

    •    x * z

    •    y * z

    •    x**2 - y**2

    •    x * y
    • f:
    •    z**3

    •    x * z**2

    •    y * z**2

    •    z * (x**2 - y**2)

    •    x * y * z

    •    x * (x**2 - 3 * y**2)

    •    y * (3 * x**2 - y**2)
    • sp:
    •    1 + x

    •    1 - x
    • sp2:
    •    1 / sqrt(3) - x / sqrt(6) + y / sqrt(2)

    •    1 / sqrt(3) - x / sqrt(6) - y / sqrt(2)

    •    1 / sqrt(3) + 2 * x / sqrt(6)
    • sp3:
    •    1 + x + y + z

    •    1 + x - y  - z

    •    1 - x + y - z

    •    1 - x - y + z
    • sp3d:
    •    1 / sqrt(3) - x / sqrt(6) + y / sqrt(2)

    •    1 / sqrt(3) - x / sqrt(6) - y / sqrt(2)

    •    1 / sqrt(3) + 2 * x / sqrt(6)

    •    z / sqrt(2) + z**2 / sqrt(2)

    •    -z / sqrt(2) + z**2 / sqrt(2)
    • sp3d2:
    •    1 / sqrt(6) - x / sqrt(2) - z**2 / sqrt(12) + (x**2 - y**2) / 2

    •    1 / sqrt(6) + x / sqrt(2) - z**2 / sqrt(12) + (x**2 - y**2) / 2

    •    1 / sqrt(6) - y / sqrt(2) - z**2 / sqrt(12) - (x**2 - y**2) / 2

    •    1 / sqrt(6) + y / sqrt(2) - z**2 / sqrt(12) - (x**2 - y**2) / 2

    •    1 / sqrt(6) - z / sqrt(2) + z**2 / sqrt(3)

    •    1 / sqrt(6) + z / sqrt(2) + z**2 / sqrt(3)

读取

 1# read from wannier90_orbital.dat
 2f = open("wannier90_orbital.dat","r")
 3orbitals = []
 4print('-----------------------------------------------------------')
 5SpinlessOrNot  = int(f.readline().split("#")[0])
 6if SpinlessOrNot == 1:
 7    print('-----------------------    SOC    -----------------------')
 8else:
 9    print('----------------------    w/o SOC    ----------------------')
10Wan_num = int(f.readline())
11print(Wan_num)
12if Wan_num != len(model.pos):
13    print('The second row is the number of wannier orbital')
14Wan_num_str = f.readlines()
15Orb_enforce = [];
16for line in Wan_num_str:
17    temp_str = line.strip('\n')
18    print(temp_str )
19    if temp_str[0] == '#':
20        continue
21    temp_str_list = temp_str.split(',')
22    Orbital_coords = [float(temp_str_list[1]),float(temp_str_list[2]),float(temp_str_list[3])]
23    Orbital_temp_str = temp_str_list[0]
24    #print(Orbital_temp_str)
25    #Orbital_temp_str = temp_str_list[0].strip()
26    try:
27        Orbital_str = sr.WANNIER_ORBITALS[Orbital_temp_str]
28    except (KeyError):
29        Orbital_str = [Orbital_temp_str]
30    except:
31        print('Unknown error')
32    # if sr.WANNIER_ORBITALS[Orbital_temp_str] is None:
33    #     Orbital_str = Orbital_temp_str
34    # else:
35        
36    for iorbstr in Orbital_str:
37        if SpinlessOrNot == 1:
38            for spin in (sr.SPIN_UP, sr.SPIN_DOWN):
39                Orb_enforce.append(Orbital_coords);
40        else:
41            Orb_enforce.append(Orbital_coords);
42if len(orbitals) != len(model.pos):
43    print('The orbitals you construct is different with wannier90_hr given')
44#print(np.array(Orb_enforce))
45model.pos = np.array(Orb_enforce);
46print('--------------------- wan information----------------------')
47print('The lattice uc is : ')
48print(model.uc)
49print('The fractional pos of each wannier_center (Default is zero)')
50print(model.pos)

z2pack set

setting

1settings = {
2    'num_lines': 41,
3    'pos_tol': 1e-2,
4    'gap_tol': 2e-2,
5    'move_tol': 0.1,
6    'iterator': range(8, 47, 2),
7    'min_neighbour_dist': 1e-2,
8}

Construct z2pack obj

1tb_system = z2pack.tb.System(model,dim=3,bands = SelectBands,pos=model.pos)
2#tb_system = z2pack.tb.System(model,bands = SelectBands)
3#tb_system = z2pack.tb.System(model)

z2pack calculation

wcc

运算

1result_wcc = z2pack.surface.run(
2        system=tb_system, surface=lambda s, t: [t, s,0], **settings 
3    )

plot

 1z2pack.plot.wcc(
 2    result_wcc,
 3    axis=ax,
 4    settings=dict(
 5        marker='o',
 6        color='b',
 7        markerfacecolor='b'
 8    )
 9)
10plt.savefig('wcc.pdf', bbox_inches='tight')

后面如有需要可以把wcc的数据弄出来

chern

 1 result_chern = z2pack.surface.run(
 2         system=tb_system, surface=lambda s, t: [t, s,0], **settings
 3     )
 4
 5 fig, ax = plt.subplots()
 6 z2pack.plot.chern(
 7     result_chern,
 8     axis=ax,
 9     settings=dict(
10         marker='o',
11         color='b',
12         markerfacecolor='b'
13     )
14 )
15 plt.savefig('chern.pdf', bbox_inches='tight')
16
17 print(
18     " Chern invariant: {0}".format(
19         z2pack.invariant.chern(result_chern)
20     )
21 )

z2

 1 result_z2 = z2pack.surface.run(
 2         system=tb_system, surface=lambda s, t: [t, s/2,0], **settings
 3     )
 4
 5 fig, ax = plt.subplots()
 6 z2pack.plot.wcc(
 7     result_z2,
 8     axis=ax
 9 )
10 plt.savefig('wcc_z2.pdf', bbox_inches='tight')
11
12 print(
13     " Chern invariant: {0}".format(
14         z2pack.invariant.z2(result_z2)
15     )
16 )

Related Posts

900 Words|This article has been read times