set terminal windows enhanced
set xrange[8:15]
set key box left top
set xlabel "Atomic volume ({\305}^3)"
set ylabel "Energy (Ry)"
plot 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' u 1:2 title "{/Symbol a}-Fe"
rep 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' u 1:4 title "{/Symbol g}-Fe"
rep 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' u 1:3 title "{/Symbol d}-Fe"

ba = 0.0869502
v0a = 11.092
e0a = -2522.82
dbdpa = 4.53622
alpha(x) = ba * x / (dbdpa * (dbdpa - 1)) * (dbdpa * (1 - v0a / x) + (v0a / x) ** dbdpa - 1) + e0a
fit alpha(x) 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' using 1:2 via ba,v0a,e0a,dbdpa
rep alpha(x) lc rgb "red" notitle


bb = 0.1
v0b = 11
e0b = -2522.8
dbdpb = 4
beta(x) = bb * x / (dbdpb * (dbdpb - 1)) * (dbdpb * (1 - v0b / x) + (v0b / x) ** dbdpb - 1) + e0b
fit beta(x) 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' using 1:4 via bb,v0b,e0b,dbdpb
rep beta(x) lc rgb "green" notitle

bc = 0.1
v0c = 11
e0c = -2522.8
dbdpc = 4
delta(x) = bc * x / (dbdpc * (dbdpc - 1)) * (dbdpc * (1 - v0c / x) + (v0c / x) ** dbdpc - 1) + e0c
fit delta(x) 'C:\cygwin\home\Administrator\kkr\cpa2002v009c\project\stableFe\analysis\Energy.dat' using 1:3 via bc,v0c,e0c,dbdpc
rep delta(x) lc rgb "blue" notitle