#!/usr/bin/env python
import matplotlib.pyplot as plt
import matplotlib as mpl

# ENCUT=400 eV, 3x3 supercell, 4,6,8,10 kpoints
Eb3_400_x = [4,5,6,7,8,9,10]
Eb3_400 = [[44.003,	43.483,	42.888,	42.945],
[46.623,	45.669,	46.252,	46.181],
[47.883,	47.547,	47.391,	47.320],
[45.806,	45.267,	45.350,	45.649],
[45.054,	45.267,	45.190,	44.917],
[44.043,	44.136,	43.942,	43.720],
[44.646,	44.931,	44.584,	44.605]]
				
# ENCUT=400 eV, 4x4 supercell, 4,6,8 kpoints
Eb4_400_x = [4,5,6,7,8,9,10]
Eb4_400 = [[43.993,	42.214,	42.772,	42.767],
[46.007,	45.753,	45.657,	45.314],
[48.232,	46.819,	47.549,	47.419],
[46.065,	45.904,	45.473,	45.648],
[44.668,	45.073,	44.805,	44.513],
[43.916,	44.322,	43.827,	43.604],
[44.825,	44.859,	45.002,	44.754]]

mpl.rc("font", size=15)
fig, ax = plt.subplots(nrows=2, ncols=1, figsize=(8,8))

plt.subplot(2,1,1)
plt.plot( Eb3_400_x, Eb3_400, marker='o' )
#y_avg = (Eb3_400[5][3] + Eb3_400[6][3]) / 2
y_avg = (Eb4_400[5][3] + Eb4_400[6][3]) / 2
plt.plot( [0,12], [y_avg, y_avg], linestyle='--', c='black' )
plt.legend(['4', '6', '8', '10'])
plt.xlim(3.5,10.5)
plt.ylim(42.,49.)
plt.tick_params(labelbottom='off', length=12, width=1)
plt.ylabel('Barrier (kJ/mol)')
plt.annotate('3x3', xy=(8, 47.), size=20)

plt.subplot(2,1,2)
plt.plot( Eb4_400_x, Eb4_400, marker='o' )
#plt.plot( Eb3_400_x, Eb3_400, marker='o' )
#y_avg = (Eb4_400[5][2] + Eb4_400[6][2]) / 2
plt.plot( [0,12], [y_avg, y_avg], linestyle='--', c='black' )
plt.legend(['3', '4', '6', '8'])
plt.xlim(3.5,10.5)
plt.ylim(42.,49.)
plt.annotate('4x4', xy=(8, 47.), size=20)
plt.xlabel('Amount of layers')
plt.ylabel('Barrier (kJ/mol)')
plt.tick_params(length=12, width=1)

plt.tight_layout()
plt.subplots_adjust(wspace=0, hspace=0)
plt.savefig('convergence.pdf')
#plt.show()