In [1]:

    

import matplotlib.pyplot as plt
%matplotlib inline
from pymatgen.util.plotting_utils import get_publication_quality_plot


    

In [2]:

    

import numpy as np
import pandas as pd


    

In [3]:

    

ls


    

    




Q3.ipynb            Si_50_11_10.26.out  Si_50_13_10.26.out  Si_50_5_10.26.out   Si_50_7_10.26.out   Si_50_9_10.26.out   results.csv
Si_50_10_10.26.out  Si_50_12_10.26.out  Si_50_14_10.26.out  Si_50_6_10.26.out   Si_50_8_10.26.out   analyze.py          results_Q4.csv

In [4]:

    

datafile = pd.read_csv('results_Q4.csv')


    

In [5]:

    

datafile.head(0)


    

    
Out[5]:






  

    

      

      
filename
      
ecut (Ry)
      
nkpts
      
alat
      
energy
      
total_force
      
Force Total (eV/Angstrom)
      
Force Difference (meV/Angstrom)
    
  
  

    

      
0
      
  Si_50_5_10.26.out
      
 50
      
  39
      
 10.26
      
-93.430544
      
 0.150617
      
 3.872520
      
  0.000000
    
    

      
1
      
  Si_50_6_10.26.out
      
 50
      
  68
      
 10.26
      
-93.433690
      
 0.149247
      
 3.837296
      
 35.224129
    
    

      
2
      
  Si_50_7_10.26.out
      
 50
      
 100
      
 10.26
      
-93.434826
      
 0.148690
      
 3.822975
      
 14.321051
    
    

      
3
      
  Si_50_8_10.26.out
      
 50
      
 150
      
 10.26
      
-93.435264
      
 0.148453
      
 3.816881
      
  6.093517
    
    

      
4
      
  Si_50_9_10.26.out
      
 50
      
 205
      
 10.26
      
-93.435442
      
 0.148346
      
 3.814130
      
  2.751082
    
    

      
5
      
 Si_50_10_10.26.out
      
 50
      
 282
      
 10.26
      
-93.435517
      
 0.148296
      
 3.812845
      
  1.285552
    
    

      
6
      
 Si_50_11_10.26.out
      
 50
      
 366
      
 10.26
      
-93.435549
      
 0.148269
      
 3.812151
      
  0.694198
    
    

      
7
      
 Si_50_12_10.26.out
      
 50
      
 476
      
 10.26
      
-93.435564
      
 0.148258
      
 3.811868
      
  0.282821
    
    

      
8
      
 Si_50_13_10.26.out
      
 50
      
 595
      
 10.26
      
-93.435571
      
 0.148250
      
 3.811662
      
  0.205688
    
    

      
9
      
 Si_50_14_10.26.out
      
 50
      
 744
      
 10.26
      
-93.435575
      
 0.148243
      
 3.811482
      
  0.179977
    
  

In [7]:

    

x = datafile['nkpts'].tolist()
x


    

    
Out[7]:





[39, 68, 100, 150, 205, 282, 366, 476, 595, 744]

In [9]:

    

y_energy = datafile['Force Total (eV/Angstrom)'].tolist()
y_conver = datafile['Force Difference (meV/Angstrom)'].tolist()


    

In [33]:

    

fig, (ax0,ax1)=plt.subplots(nrows=2)
fig.set_size_inches(10,10)

ax0.scatter(x, y_energy,marker='o',alpha=0.75,s=60,c='r')
ax0.plot(x,y_energy,c='black')
ax0.set_title('Forces VS K-points',fontsize=30,y=1.05)
ax0.set_xlabel('# of K-points',fontsize=24)
ax0.set_ylabel(r'Force $eV/\AA$',fontsize=20)
ax0.set_ylim(3.81,3970)
ax0.tick_params(labelsize = 18)

ax1.scatter(x,y_conver,marker='o',alpha=0.75, s=60,c='b')
ax1.plot(x,y_conver,c='black')
ax1.set_title('Energy Difference',fontsize=24,y=1.05)
ax1.set_xlabel('Kinetic Energy Cutoff (Ry)',fontsize=24)
ax1.set_ylabel('Force Converge $meV/\AA$',fontsize=18)
# ax1.set_yscale('log')
ax1.set_ylim(-10,100)
ax1.annotate(r'Convergence < 5 $meV/\AA$ Cutoff 30 Ry', xy=(32, 5),  xycoords='data',
                xytext=(0.95, 0.95), textcoords='axes fraction',
                arrowprops=dict(facecolor='black', shrink=0.05),
                horizontalalignment='right', verticalalignment='top',fontsize=16
                )
ax1.tick_params(labelsize = 18)

plt.subplots_adjust(hspace=0.5)
plt.show()


    

In [17]:

    




    

    
Out[17]:





<matplotlib.text.Text at 0x10891e9d0>

In [18]:

    

plt.show()