def dist(I,J):
    return math.sqrt( (I[0]-J[0])**2 + (I[1]-J[1])**2 )
G=nx.empty_graph(N) #empty graph
I=[(random.random(),random.random()) for i in range(0,N)] #nodes positions
for i in range(0,N-1):
    for j in range(i+1,N):
        xm=0.5*(I[i][0] + I[j][0]) #M is the center of i and j
        ym=0.5*(I[i][1] + I[j][1])
        M=(xm,ym)
        d=min(dist(M,I[i]),dist(M,I[j])) #d is the half-distance between i and j
        possibilites=[k for k in range(0,N) if xm-d<I[k][0]<xm+d and ym-d<I[k][1]<ym+d] #we are only looking around the point M
        inter=[k for k in possibilites if dist(I[k],M)<d]
        if len(inter)==0:
            G.add_edge(i,j)