#!/usr/bin/env python
import matplotlib.pyplot as plt
import matplotlib as mpl
from numpy import arange, array, append

eV2kJ = 96.4853075

H15 = array(
[[2	,-1.9530201E+02],
[2.25	,-1.9667989E+02],
[2.5	,-1.9723435E+02],
[2.75	,-1.9744930E+02],
[3	,-1.9752748E+02],
[3.25	,-1.9754805E+02],
[3.5	,-1.9754324E+02],
[3.75	,-1.9752771E+02],
[4	,-1.9750826E+02],
[4.25	,-1.9748812E+02],
[4.5	,-1.9746865E+02],
[4.75	,-1.9745049E+02],
[5	,-1.9743393E+02],
[5.25	,-1.9741896E+02],
[5.5	,-1.9740566E+02],
[5.75	,-1.9739415E+02],
[6	,-1.9738442E+02],
[6.25	,-1.9737634E+02],
[6.5	,-1.9736973E+02],
[6.75	,-1.9736445E+02],
[7	,-1.9736048E+02],
[7.25	,-1.9735784E+02],
[7.5	,-1.9735649E+02]]
)
H15[:,1] -= H15[len(H15)-1][1]
H15[:,1] *= eV2kJ

N15 = array(
[[1	,-1.4557077E+02],
[1.25	,-1.7787255E+02],
[1.5	,-1.9103218E+02],
[1.75	,-1.9607795E+02],
[1.9	,-1.9732938E+02],
[2	,-1.9776464E+02],
[2.1	,-1.9800029E+02],
[2.2	,-1.9810721E+02],
[2.25	,-1.9812840E+02],
[2.3	,-1.9813433E+02],
[2.4	,-1.9811412E+02],
[2.5	,-1.9806767E+02],
[2.75	,-1.9791524E+02],
[3	,-1.9777786E+02],
[3.25	,-1.9767577E+02],
[3.5	,-1.9760358E+02],
[3.75	,-1.9755175E+02],
[4	,-1.9751276E+02],
[4.25	,-1.9748220E+02],
[4.5	,-1.9745760E+02],
[4.75	,-1.9743754E+02],
[5	,-1.9742081E+02],
[5.25	,-1.9740648E+02],
[5.5	,-1.9739412E+02],
[5.75	,-1.9738372E+02],
[6	,-1.9737533E+02],
[6.25	,-1.9736877E+02],
[6.5	,-1.9736375E+02],
[6.75	,-1.9736002E+02],
[7	,-1.9735759E+02],
[7.25	,-1.9735642E+02],
[7.5	,-1.9735648E+02]]
)
N15[:,1] -= N15[len(N15)-1][1]
N15[:,1] *= eV2kJ

mpl.rc("text", usetex=True)
fig = plt.figure(figsize=(3.69,3))
plt.rcParams.update({'font.size': 12})
				
plt.plot( H15[:,0], H15[:,1], marker='o', label='H', c='b' )
plt.plot( N15[:,0], N15[:,1], marker='o', label='N', c='r' )
#plt.plot( OH15Ar, OH15AE, marker='o', label='15A OH', c='r' )
#plt.plot( CH315Ar, CH315AE, marker='o', label='15A CH3', c='g' )
#plt.plot( flat30Ar, flat30AE, marker='o', label='30A perpendicular', c='b', markerfacecolor='None', markeredgecolor='b' )
#plt.plot( OH30Ar, OH30AE, marker='o', label='30A OH', c='r', markerfacecolor='None', markeredgecolor='r' )
#plt.plot( CH330Ar, CH330AE, marker='o', label='30A CH3', c='g', markerfacecolor='None', markeredgecolor='g' )
plt.plot( [0.,7.5], [0,0], c='black')
plt.legend(numpoints=1, loc='lower right')
plt.xlim(1.5,7.5)
plt.ylim(-80,10)
plt.xlabel(r'Distance N - surface (\r{A})')
plt.ylabel('Potential energy (kJ/mol)')
plt.tick_params(length=6, width=1)

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