Add wannier plot

band_comp.py

@@ -9,6 +9,7 @@
 import matplotlib
 matplotlib.use("Agg")
 from matplotlib import pyplot as plt
+import argparse
 
 def get_ef_from_scfout(scfout):
     """
@@ -78,11 +79,11 @@ def get_klabel(wannier_band_gnu):
     return x_list, label_list
 
 
-def main():
+def main(output_dir):
     scfout = "scf.out"
     wannier_band = "pwscf_band.dat"
     wannier_band_gnu = "pwscf_band.gnu"
-    pwscf_band = "./band/bands.out.gnu"
+    pwscf_band = "../band/bands.out.gnu"
 
 
     x, y = get_band_data(wannier_band)
@@ -104,9 +105,14 @@ def main():
     y_max = np.max(y-ef)
     plt.ylim([y_min - 0.05*(y_max-y_min), y_max + 0.05*(y_max-y_min)])
 
-    plt.savefig("./band/band_compare.png", bbox_inches="tight")
-    plt.savefig("./band/band_compare.eps", bbox_inches="tight")
+    plt.savefig(f"{output_dir}/band_compare.png", bbox_inches="tight")
+    plt.savefig(f"{output_dir}/band_compare.eps", bbox_inches="tight")
 
 
 if __name__ == "__main__":
-    main()
+    parser = argparse.ArgumentParser(description="compare band")
+    parser.add_argument("-o", dest="odir")
+
+    args = parser.parse_args()
+    output_dir = args.odir
+    main(output_dir)

cif2qewan.py

@@ -280,6 +280,8 @@ def write_pw2wan(self, pw2wan):
             fp.write(" prefix = 'pwscf'\n")
             fp.write(" seedname = 'pwscf'\n")
             fp.write(" spin_component = 'none'\n")
+            if "wannier_plot_supercell" in info.keys():
+              fp.write(f" wannier_plot_supercell = {info['wannier_plot_supercell']}\n")
             fp.write(" write_mmn = .true.\n")
             fp.write(" write_amn = .true.\n")
             fp.write(" write_unk = {}\n".format(info["write_unk"]))
@@ -348,6 +350,29 @@ def write_proj(self, proj_in):
             fp.write(" Emin = \n")
             fp.write("/\n")
 
+    def write_pp(self, pp_in):
+        so_factor = 1
+        if(self.so or self.mag): so_factor = 2
+        with open(pp_in, "w") as fp:
+            fp.write("! generated by {}\n".format(__file__.split("/")[-1]))
+            fp.write("&inputpp\n")
+            fp.write(" prefix = 'pwscf'\n")
+            fp.write(" outdir = './work'\n")
+            fp.write(" filplot = 'wf_pp'\n")
+            fp.write(" plot_num = 7\n")
+            fp.write(" kpoint = 1\n")
+            fp.write(" kband(1) = {}\n".format(self.nexclude*so_factor + 1))
+            fp.write(" kband(2) = {}\n".format(self.num_wann*so_factor + self.nexclude*so_factor))
+            fp.write(" lsign = .TRUE.\n")
+            fp.write("/\n")
+            fp.write("&plot\n")
+            fp.write(" iflag = 3\n")
+            fp.write(" output_format = 5\n")
+            fp.write(" fileout = '.xsf'\n")
+            fp.write("/\n")
+
+
+
 
 class pseudo_list:
     def __init__(self, pseudo_dir, pp_list_path):
@@ -398,7 +423,6 @@ def read_pp_info(self, pp_list_path):
     qe_wan.calc_bands_seekpath()
     qe_wan.write_pw2wan("pw2wan.in")
     qe_wan.write_wannier("pwscf.win")
-    qe_wan.write_proj("proj.in")
 
     if not os.path.exists("check_wannier"): os.mkdir("check_wannier")
     qe_wan.shift_k_nscf()
@@ -408,3 +432,5 @@ def read_pp_info(self, pp_list_path):
     qe_wan.convert2band()
     qe_wan.write_pwscf_in("band/nscf.in")
     qe_wan.write_band_in("band/band.in")
+    qe_wan.write_proj("band/proj.in")
+    qe_wan.write_pp("band/pp.in")

pp_psl_rrkj.csv

@@ -11,13 +11,13 @@ F,F.pbe-n-rrkjus_psl.0.1,1,p,48.0,316.0
 Ne,Ne.pbe-n-rrkjus_psl.1.0.0,1,p,54.0,264.0
 Na,Na.pbe-spn-rrkjus_psl.0.2,4,s,41.0,163.0
 Mg,Mg.pbe-n-rrkjus_psl.0.3.0,0,s,13.0,87.0
-Al,Al.pbe-n-rrkjus_psl.0.2.2,0,sp,36.0,145.0
+Al,Al.pbe-n-rrkjus_psl.0.1,0,sp,36.0,145.0
 Si,Si.pbe-n-rrkjus_psl.0.1,0,sp,38.0,151.0
 P,P.pbe-n-rrkjus_psl.0.1,1,p,17.0,79.0
 S,S.pbe-n-rrkjus_psl.0.1,1,p,17.0,77.0
 Cl,Cl.pbe-n-rrkjus_psl.0.3.0,1,p,36.0,146.0
 Ar,Ar.pbe-n-rrkjus_psl.0.3.0,1,p,26.0,225.0
-K,K.pbe-spn-rrkjus_psl.0.2.3,4,,34.0,366.0
+K,K.pbe-spn-rrkjus_psl.0.2.3,4,s,34.0,366.0
 Ca,Ca.pbe-spn-rrkjus_psl.0.2.3,4,s,36.0,170.0
 Sc,Sc.pbe-spn-rrkjus_psl.0.2.3,4,sd,46.0,477.0
 Ti,Ti.pbe-spn-rrkjus_psl.0.3.1,4,spd,50.0,738.0

wannier_conv.py

@@ -3,6 +3,7 @@
 import itertools
 import numpy as np
 import os
+import argparse
 
 
 class Hamiltonian(object):
@@ -118,7 +119,7 @@ def __init__(self, nscf_out):
         with open(nscf_out, "r") as fp:
             lines = fp.readlines()
         for i, line in enumerate(lines):
-            if "Fermi energy" in line:
+            if "the Fermi energy" in line:
                 self.ef = float(line[26:35])
             if "number of Kohn-Sham" in line:
                 self.nbnd = int(line[35:])
@@ -138,7 +139,19 @@ def __init__(self, nscf_out):
             kp_str = "k =%7.4f%7.4f%7.4f" % tuple(self.kp_cart[j])
             for i, line in enumerate(lines):
                 if kp_str in line:
-                    self.energy[j, :] = [float(x) for x in ''.join(lines[i+2:i+2+nline]).split()]
+                    tmp_list = ''.join(lines[i+2:i+2+nline]).split()
+                    ene_list = []
+                    for l_tmp in tmp_list:
+                        len_max = 8
+                        if len(l_tmp) > len_max:
+                            #print(len(l_tmp),len(l_tmp)//(len_max+1))
+                            [ene_list.append(l_tmp[s*(len_max+1):(s+1)*(len_max+1)]) for s in range(len(l_tmp)//(len_max+1)) ]
+                        else:
+                            ene_list.append(l_tmp)
+
+                    #print(ene_list)
+                    #self.energy[j, :] = [float(x) for x in ''.join(lines[i+2:i+2+nline]).split()]
+                    self.energy[j, :] = [float(x) for x in ene_list]
 
 
 def get_nexclude(pwscf_win):
@@ -152,13 +165,23 @@ def get_nexclude(pwscf_win):
 
 
 if __name__ == "__main__":
-    nscf_data = Nscfout("./check_wannier/nscf.out")
+    parser = argparse.ArgumentParser(description="check accuracy of Wannier")
+    parser.add_argument("-e", type=float, default=0.0, dest="emax", help="emax")
+    parser.add_argument("-o", dest="odir")
+    parser.add_argument("-i", default="./check_wannier/nscf.out", dest="nscf_file")
+
+    args = parser.parse_args()
+
+    # Energy window for check
+    emin = -200.0
+    emax = args.emax #0.0
+    output_dir = args.odir
+    nscf_file = args.nscf_file
+
+    nscf_data = Nscfout(nscf_file)
     nexclude = get_nexclude("./pwscf.win")
     h = Hamiltonian("./pwscf_hr.dat")
 
-    # Energy window for check
-    emin = -100.0
-    emax = 0.0
 
     # calculate energy difference
     delta_sum = 0
@@ -180,7 +203,8 @@ def get_nexclude(pwscf_win):
         delta_max = max([delta_max, np.max(ediff)])
 
     # output the results
-    with open("check_wannier/CONV", "w") as fp:
+    #with open("check_wannier/CONV", "w") as fp:
+    with open(os.path.join(output_dir, "CONV_{}".format(emax)), "w") as fp:
         fp.write("# energy window [{:>5.2f}:{:>5.2f}]\n".format(emin, emax))
         if nek > 0:
             fp.write("average diff = {:>15.8f}\n".format(np.sqrt(delta_sum/nek)))