Convergence plots.py
From EVOCD
This code generates plots for the convergency study from the SUMMARY
files generated by the Ev curve bash script using data from Code: Quantum Espresso simulations.
This code is part of a group of five codes all used for post-processing data from Quantum Espresso simulations.
- The EvA EvV plot.py code generates plots for the energy vs volume and energy vs lattice parameter curves for any set of EvsA an EvsV files generated by the ev_curve.bash script.
- The Convergence plots.py code generates plots for the convergency study from the
SUMMARY
files. - The EOS plot.py code plots the difference in the energy values found from DFT and from the chosen EOS.
- The EOS comp plot.py code plots the energy vs lattice parameter curve from DFT and EOS fitting for multiple equations of state and compares the difference between each equation of state to values found directly from DFT.
- The Ecut conv.py code plots the values found from the ecut convergence studies.
# This code generates plots for k point convergence studies # It will need to be modified depending on the number and range of # K points tested import numpy as np import matplotlib.pyplot as plt import os # This file will have to be edited to match your file system and # number of convergence runs done. I recommend the following file system # # Number2 (for all data related to convergence stydies # L____________________________________ # \ \ \ \ \ # k7 k9 k11 k13 etc for k# = number of runs # \ # This code expects a SUMMARY file in each of the k# folders # Alter as needed # Getting the path to current working directory. This should be ran from the # Number 2 directory path = os.getcwd() ####################################################### ### User Settings ### ####################################################### # Print to screen/pdf or not write_to_screen = False print_to_pdf = True energy_offset = 2858.8298734 # Set by running simulation with very large lattice parameter #Color options for plots # Best source for colors = https://xkcd.com/color/rgb/ color1 = '#b66325' #for FCC color2 = '#8b2e16' #For BCC color3 = '#00035b' #For HCP mk_size = 18 # Markersize for plots tick_size = 24 ax_label_size = 36 # Font size for axis labels title_font_size = 42 # Font size for table title marker_type = 'o' # Type of marker used for plotting google matplotlib for more options #array with # of k points sized by your number of runs k = [7,9,11,13,15,17,19,21,23,25] # this will need to be altered depending on your convergence study # this also assumes you put each of your runs in a filder called k# where # # is the number of k points you used # the following block of for loops searches each of your folders and compiles # all of the values then writes them to a plot. This will be edited to match your # number of k point values. yes its disghusting. ################################################################# #Initializing arrays EE = [] # Equilibrium energy A0 = [] # array of Lattice Parameter in Angrstrom K0 = [] # array of Energy in eV time = [] # run time lines = [] # for k7 lines = [] with open (path+"/k7/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k9 lines = [] with open (path+"/k9/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k11 lines = [] with open (path+"/k11/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k13 lines = [] with open (path+"/k13/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k15 lines = [] with open (path+"/k15/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k17 lines = [] with open (path+"/k17/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k19 lines = [] with open (path+"/k19/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k21 lines = [] with open (path+"/k21/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k23 lines = [] with open (path+"/k23/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) # for k25 lines = [] with open (path+"/k25/SUMMARY") as in_file: for line in in_file: lines.append(line) EE1 = lines[2].split() EE.append(float(EE1[5])) A01 = lines[3].split() A0.append(float(A01[4])) K01 = lines[4].split() K0.append(float(K01[3])) t1 = lines[5].split() time.append(float(t1[5])) bulkm = np.array(K0) cpu_time = np.array(time) latt_param = np.array(A0) kpts = np.array(k) eq_energy = np.array(EE) eq_energy = eq_energy + energy_offset #shifting the energy values # X axis ticks every 2 steps fom 5 k pts to 27 kpts, modify as needed major_ticks = np.arange(5, 27, 2) eng_spc = np.arange(-10, 0, 2) #creating plots A0_converge, ax = plt.subplots(figsize = (12,10)) ax.plot(k, latt_param, marker_type, markersize = mk_size, color = color2) ax.set_xticks(major_ticks) ax.set_xlabel('Number of k points', fontsize = ax_label_size) ax.set_ylabel('Lattice Parameter ($\AA$)' , fontsize = ax_label_size) ax.set_title('Lattice Parameter vs. k Point Refinement', fontsize = title_font_size) ax.grid(True, linestyle='-.') ax.tick_params(labelcolor='k', labelsize='large', width=3) ax.tick_params(labelsize=tick_size) K0_converge, ax = plt.subplots(figsize = (12,10)) ax.plot(k, bulkm, marker_type,markersize = mk_size, color = color2) ax.set_xticks(major_ticks) ax.set_xlabel('Number of k points', fontsize = ax_label_size) ax.set_ylabel('Bulk Modulus (kbar)' , fontsize = ax_label_size) ax.set_title('Bulk Modulus vs. k Point Refinement', fontsize = title_font_size) ax.grid(True, linestyle='-.') ax.tick_params(labelcolor='k', labelsize='large', width=3) ax.tick_params(labelsize=tick_size) run_time , ax = plt.subplots(figsize = (12,10)) ax.plot(k, cpu_time, marker_type, markersize = mk_size, color = color2) ax.set_xticks(major_ticks) ax.set_xlabel('Number of k points', fontsize = ax_label_size) ax.set_ylabel('CPU time (sec)' , fontsize = ax_label_size) ax.set_title('CPU time vs. k Point Refinement', fontsize = title_font_size) ax.grid(True, linestyle='-.') ax.tick_params(labelcolor='k', labelsize='large', width=3) equilibrium , ax = plt.subplots(figsize = (12,10)) ax.plot(k, eq_energy, marker_type, markersize = mk_size, color = color2) ax.set_xticks(major_ticks) ax.set_xlabel('Number of k points', fontsize = ax_label_size) ax.set_ylabel('Equilibrium Energy ($eV$)' , fontsize = ax_label_size) ax.set_title('Equilibrium Energy vs. k Point Refinement', fontsize = title_font_size) ax.grid(True, linestyle='-.') ax.tick_params(labelcolor='k', labelsize='large', width=3) ax.tick_params(labelsize=tick_size) #consolidating to a single plot god_fig, axarr = plt.subplots(2,2, figsize = (30,25), sharex='col') axarr[0, 0].plot(k, latt_param, marker_type, markersize = mk_size, color = color2) axarr[0, 0].set_xticks(major_ticks) # axarr[0, 0].set_xlabel('Number of k points', fontsize = ax_label_size) axarr[0, 0].set_ylabel('Lattice Parameter ($\AA$)' , fontsize = ax_label_size) axarr[0, 0].set_title('Lattice Parameter vs. k Point Refinement', fontsize = title_font_size) axarr[0, 0].grid(True, linestyle='-.') axarr[0, 0].tick_params(labelcolor='k', labelsize=tick_size, width=3) axarr[0, 1].plot(k, bulkm, marker_type,markersize = mk_size, color = color2) axarr[0, 1].set_xticks(major_ticks) # axarr[0, 1].set_xlabel('Number of k points', fontsize = ax_label_size) axarr[0, 1].set_ylabel('Bulk Modulus (kbar)' , fontsize = ax_label_size) axarr[0, 1].set_title('Bulk Modulus vs. k Point Refinement', fontsize = title_font_size) axarr[0, 1].grid(True, linestyle='-.') axarr[0, 1].tick_params(labelcolor='k', labelsize=tick_size, width=3) axarr[1, 0].plot(k, cpu_time, marker_type, markersize = mk_size, color = color2) axarr[1, 0].set_xticks(major_ticks) axarr[1, 0].set_xlabel('Number of k points', fontsize = ax_label_size) axarr[1, 0].set_ylabel('CPU time (sec)' , fontsize = ax_label_size) axarr[1, 0].set_title('CPU time vs. k Point Refinement', fontsize = title_font_size) axarr[1, 0].grid(True, linestyle='-.') axarr[1, 0].tick_params(labelcolor='k', labelsize=tick_size, width=3) axarr[1, 1].plot(k, eq_energy, marker_type, markersize = mk_size, color = color2) axarr[1, 1].set_xticks(major_ticks) axarr[1, 1].set_xlabel('Number of k points', fontsize = ax_label_size) axarr[1, 1].set_ylabel('Equilibrium Energy ($eV$)' , fontsize = ax_label_size) axarr[1, 1].set_title('Equilibrium Energy vs. k Point Refinement', fontsize = title_font_size) axarr[1, 1].grid(True, linestyle='-.') axarr[1, 1].tick_params(labelcolor='k', labelsize=tick_size, width=3) # Printing plots to a .pdf file if print_to_pdf == True: A0_converge.savefig("Latt_param.pdf", bbox_inches = 'tight') K0_converge.savefig("Bulk_mod.pdf", bbox_inches = 'tight') run_time.savefig("run_time.pdf", bbox_inches = 'tight') equilibrium.savefig("eq_energy.pdf", bbox_inches = 'tight') god_fig.savefig("combined_conv.pdf", bbox_inches = 'tight') # Print to screen if write_to_screen == True: plt.show() |