Convergence plots.py

From EVOCD
Jump to: navigation, search

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()





Code: Quantum Espresso

Personal tools
Namespaces

Variants
Actions
home
Materials
Material Models
Design
Resources
Projects
Education
Toolbox