In [1]:

    

import MagicTools


    

In [2]:

    

%pylab inline


    

    




Populating the interactive namespace from numpy and matplotlib

In [3]:

    

MagicTools.Deviation(['0magic0.out','0magic1.out','0magic2.out','0magic3.out'])


    

    




0magic0.out.iter1    S:   5.52561       RDF:   3.37115
0magic0.out.iter2    S:   2.39021       RDF:   0.69364
0magic0.out.iter3    S:   1.48770       RDF:   0.43659
0magic0.out.iter4    S:   1.03896       RDF:   0.45803
0magic0.out.iter5    S:   0.85950       RDF:   0.30883
0magic0.out.iter6    S:   0.74291       RDF:   0.29362
0magic0.out.iter7    S:   0.72854       RDF:   0.27098
0magic0.out.iter8    S:   0.64791       RDF:   0.23766
0magic0.out.iter9    S:   0.61064       RDF:   0.21673
0magic0.out.iter10    S:   0.60868       RDF:   0.18731
0magic1.out.iter1    S:   0.46928       RDF:   0.15734
0magic1.out.iter2    S:   1.16588       RDF:   0.48066
0magic1.out.iter3    S:   1.41678       RDF:   0.57081
0magic1.out.iter4    S:   1.45686       RDF:   0.57529
0magic1.out.iter5    S:   1.31140       RDF:   0.51230
0magic1.out.iter6    S:   0.39708       RDF:   0.14900
0magic1.out.iter7    S:   0.36810       RDF:   0.14673
0magic1.out.iter8    S:   0.38558       RDF:   0.16176
0magic1.out.iter9    S:   0.33607       RDF:   0.15052
0magic1.out.iter10    S:   0.28330       RDF:   0.15601
0magic2.out.iter1    S:   0.28138       RDF:   0.14967
0magic2.out.iter2    S:   0.18970       RDF:   0.08410
0magic2.out.iter3    S:   0.18404       RDF:   0.03791
0magic2.out.iter4    S:   0.18306       RDF:   0.05579
0magic2.out.iter5    S:   0.11599       RDF:   0.04050
0magic2.out.iter6    S:   0.07843       RDF:   0.03290
0magic2.out.iter7    S:   0.25446       RDF:   0.12119
0magic2.out.iter8    S:   0.18376       RDF:   0.08919
0magic2.out.iter9    S:   0.14580       RDF:   0.04161
0magic2.out.iter10    S:   0.11295       RDF:   0.03040
0magic3.out.iter1    S:   0.09846       RDF:   0.04370
0magic3.out.iter2    S:   0.18243       RDF:   0.04362
0magic3.out.iter3    S:   0.17147       RDF:   0.05027
0magic3.out.iter4    S:   0.11582       RDF:   0.03514
0magic3.out.iter5    S:   0.07416       RDF:   0.02440

In [5]:

    

MagicTools.AnalyzeIMCOuput('0magic3.out')


    

    




Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 5
Total number of checks per iteration: 7
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 3
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 4
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 5
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 5
Total number of checks per iteration: 7
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2

In [6]:

    

import MagicTools


    

In [199]:

    

rdf = []
rdf.append(MagicTools.GetDFsFromFile_magic('0magic0.out'));
rdf.append(MagicTools.GetDFsFromFile_magic('0magic1.out'));
rdf.append(MagicTools.GetDFsFromFile_magic('0magic2.out'));
rdf.append(MagicTools.GetDFsFromFile_magic('0magic3.out'));


    

    




Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 10
Total number of checks per iteration: 22
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 3
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 4
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 5
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 6
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 7
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 8
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 9
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 10
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 10
Total number of checks per iteration: 7
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 3
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 4
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 5
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 6
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 7
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 8
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 9
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 10
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 10
Total number of checks per iteration: 7
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 3
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 4
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 5
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 6
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 7
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 8
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 9
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 10
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Automatically defining list of interactions from input file. If this does not work well for your input, edit subroutine GetDFsFromFile_magic in MagicTools.py and manualy add list of your beads to PairNamesList
3  intermolecular pairs detected  ['N-N', 'N-C', 'C-C']
2  pairwise bonds detected  ['N1-C1', 'C1-C2']
1  angles detected  ['N1-C2-C1']
Total number of iterations: 5
Total number of checks per iteration: 7
Iteration 1
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 2
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 3
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 4
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1
Iteration 5
Pair N-N
Pair N-C
Pair C-C
Pair N1-C1
Pair C1-C2
Pair N1-C2-C1






    
Out[199]:





35

In [284]:

    

r2 = np.array(rdf)
rx =[]
its = [10,10,10,5]
for i in xrange(len(its)):
    for j in xrange(its[i]):
        rx.append(r2[i][j][0])
print len(r2)


    

    




4

In [4]:

    

p1=[]
p1.append(MagicTools.GetPotsFromFile_pot('100sC.pot',mcmfile=['M3dmdds.mcm','M3br.mcm'],Ucut= 1000))
p1.append(MagicTools.GetPotsFromFile_pot('100sN.pot',mcmfile=['M3dmdds.mcm','M3br.mcm'],Ucut= 1000))
p1.append(MagicTools.TotalPots(p1[0],70))


    

In [94]:

    

MagicTools.PlotAllDFs([p1[0]])


    

In [8]:

    

import matplotlib.cm as cm
plt.figure(figsize=(15,7), dpi=300)
j=2
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.axis([3, 20, -3,10])
plt.plot(p1[0][0].g[:,0],p1[0][0].g[:,1],'g.-',label='Ni-Ni')
plt.plot(p1[1][0].g[:,0],p1[1][0].g[:,1],'r.-',label='Ni-Ni Incorporated')
plt.plot(p1[2][0].g[:,0],p1[2][0].g[:,1],'b.-',label='Ni-Ni + Electrostatics')
plt.plot(p1[0][2].g[:,0],p1[0][2].g[:,1],'g-',label='Ni-Nn')
plt.plot(p1[1][2].g[:,0],p1[1][2].g[:,1],'r-',label='Ni-Nn Incorporated')
plt.plot(p1[0][7].g[:,0],p1[0][7].g[:,1],'g--',label='Nn-Nn')
plt.plot(p1[1][7].g[:,0],p1[1][7].g[:,1],'r--',label='Nn-Nn Incorporated')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("NNCha.png",bbox_inches='tight', dpi=300)
plt.show()


    

In [7]:

    

import matplotlib.cm as cm
for j in xrange(13):
    plt.figure(figsize=(20,10), dpi=300)
    if (j == 12):
        plt.axis([60, 180, -1,30])
    plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
    plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
    plt.plot(p1[0][j].g[:,0],p1[0][j].g[:,1],'g-',label=p1[0][j].Name)
    plt.plot(p1[1][j].g[:,0],p1[1][j].g[:,1],'r-',label=p1[0][j].Name + ' Incorporated')
    plt.plot(p1[2][j].g[:,0],p1[2][j].g[:,1],'b-',label=p1[0][j].Name + ' + Electrostatics')
    plt.legend( title='Potentials',loc='upper right')
    plt.grid(True)
    savefig("PtsChargedSys{0:}.png".format(j),bbox_inches='tight', dpi=300)
    plt.show()


    

In [3]:

    

p2=[]
p2.append(MagicTools.GetPotsFromFile_pot('45nns.pot',mcmfile=['M3dmdd.mcm','M3osai.mcm','M3na.mcm'],Ucut= 1000))
p2.append(MagicTools.GetPotsFromFile_pot('30n.pot',mcmfile=['M3dmdd.mcm','M3osai.mcm','M3na.mcm'],Ucut= 1000))
p2.append(MagicTools.GetPotsFromFile_pot('100n.pot',mcmfile=['M3dmdd.mcm','M3osai.mcm','M3na.mcm'],Ucut= 1000))


    

In [16]:

    

plt.figure(figsize=(30,5), dpi=300)
plt.subplot2grid((1,3), (0,0))
j=0
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j].g[:,0],p2[0][j].g[:,1],'g-',label=p2[0][j].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,1))
j=1
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j].g[:,0],p2[0][j].g[:,1],'g-',label=p2[0][j].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,2))
j=4
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j-2].g[:,0],p2[0][j-2].g[:,1],'g-',label=p2[0][j-2].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipotsinter.png",bbox_inches='tight', dpi=300)
plt.show()

plt.figure(figsize=(20,5), dpi=300)
plt.subplot2grid((1,2), (0,0))
j=10
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.axis([4, 7, -1,27])
plt.plot(p2[0][j-7].g[:,0],p2[0][j-7].g[:,1],'g-',label='N-C SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label='N-C SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label='N-C LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,2), (0,1))
j=11
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.axis([50, 180, -6,30])
plt.plot(p2[0][j-7].g[:,0],p2[0][j-7].g[:,1],'g-',label='C-C SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label='C-C SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label='C-C LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipotsintra.png",bbox_inches='tight', dpi=300)

plt.figure(figsize=(20,15), dpi=300)
plt.subplot2grid((2,5), (0,0))
j=3
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([4, 7, -1,27])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((2,5), (0,1))
j=4
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((2,5), (0,2))
j=5
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((2,5), (0,3))
j=6
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((2,5), (0,4))
j=7
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((2,5), (1,0))
j=8
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)

plt.subplot2grid((2,5), (1,2))
j=10
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)

plt.subplot2grid((2,5), (1,2))
j=10
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
#plt.axis([50, 180, -6,30])
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipotssi.png",bbox_inches='tight', dpi=300)
plt.show()


    

In [5]:

    

plt.figure(figsize=(30,5), dpi=300)
plt.subplot2grid((1,3), (0,0))
j=0
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j].g[:,0],p2[0][j].g[:,1],'g-',label=p2[0][j].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,1))
j=1
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j].g[:,0],p2[0][j].g[:,1],'g-',label=p2[0][j].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,2))
j=4
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[0][j-2].g[:,0],p2[0][j-2].g[:,1],'g-',label=p2[0][j-2].Name +' SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipots1.png",bbox_inches='tight', dpi=300)
plt.show()

plt.figure(figsize=(30,5), dpi=300)
plt.subplot2grid((1,3), (0,0))
j=2
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,1))
j=3
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,2))
j=5
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipots2.png",bbox_inches='tight', dpi=300)
plt.show()

plt.figure(figsize=(30,5), dpi=300)
plt.subplot2grid((1,3), (0,0))
j=6
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,1))
j=7
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,2))
j=8
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipots3.png",bbox_inches='tight', dpi=300)
plt.show()

plt.figure(figsize=(30,5), dpi=300)
plt.subplot2grid((1,3), (0,0))
j=9
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label=p2[1][j].Name +' SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label=p2[2][j].Name +' LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,1))
j=10
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.axis([4, 7, -1,27])
plt.plot(p2[0][j-7].g[:,0],p2[0][j-7].g[:,1],'g-',label='N-C SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label='N-C SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label='N-C LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
plt.subplot2grid((1,3), (0,2))
j=11
plt.xlabel(r'$r,\  (\AA)$',fontsize=20)
plt.ylabel(r'$\mathcal{U}(r),\ (\frac{kJ}{mol})$',fontsize=20)
plt.axis([50, 180, -6,30])
plt.plot(p2[0][j-7].g[:,0],p2[0][j-7].g[:,1],'g-',label='C-C SR Pure')
plt.plot(p2[1][j].g[:,0],p2[1][j].g[:,1],'r-',label='C-C SR Silica')
plt.plot(p2[2][j].g[:,0],p2[2][j].g[:,1],'b-',label='C-C LR Silica')
plt.legend( title='Potentials',loc='upper right')
plt.grid(True)
savefig("Sipots4.png",bbox_inches='tight', dpi=300)


    

In [ ]: