myprojects-hg-reborn

Fork

#! /usr/bin/env python

import scipy as s
import numpy as n
import pylab as p
#from rpy import r
import distance_calc as d_calc
import igraph as ig


box_size=10000.  #the box is now LARGE

threshold=1.06 #threshold to consider to particles as directly connected


test_arr=p.load("cluster_10_monomers")

x_pos=test_arr[:,0]
y_pos=test_arr[:,1]
z_pos=test_arr[:,2]

dist_mat=d_calc.distance_calc(x_pos,y_pos,z_pos, box_size)

#print "dist_mat is, ", dist_mat
p.save("dist_mat.dat", dist_mat)

adjacency=s.where(dist_mat <= threshold, 1, 0)

#print "adjacency is, ", adjacency

p.save("adjacency.csv", adjacency, delimiter=",")
p.save("adjacency.dat", adjacency)

cluster_obj=ig.Graph.Adjacency((dist_mat <= threshold).tolist(),\
					       ig.ADJ_UNDIRECTED)

# Now I start the calculation of the coordination number

coord_list=cluster_obj.degree()  #Now I have a list with the number of 1rst
		         #neughbors of each particle

coord_arr=s.asarray(coord_list) -2

print "coord_arr is, ", coord_arr

print "the coordination number is, ", coord_arr.mean()

p.save("neighbor_list.dat", coord_arr)


cluster_obj.simplify()
clustering=cluster_obj.clusters()

n_clusters=len(clustering)
print "the total number of clusters is,", n_clusters

def eta_0_calc(coord_numbers, n):
    eta_0=1.-(coord_numbers/12.)**(1./n)
    return eta_0


def eta_1_calc(coord_numbers,n):
    eta_1=1.-(coord_numbers/12.)**(1./n)
    return eta_1



my_exp=4.

eta_0=eta_0_calc(coord_arr, my_exp)
eta_1=eta_1_calc(coord_arr, my_exp)

p.save("eta_0", eta_0)
p.save("eta_1", eta_1)




print "So far so good"