[Numpy-discussion] Stick (line segments) percolation algorithm - graph theory?

[Josè Luis Mietta]

Hi experts!
I wrote an algorithm for study stick percolation (i.e.: networks between line segments that intersect between them). In my algorithm N sticks (line segments) are created inside a rectanglar box of sides 'b' and 'h' and then, one by one, the algorithm explores the intersection between all line segments. This is a Monte Carlo simulation, so the 'experiment' is executed many times (no less than 100 times). Writen like that, very much RAM is consumed:
array_x1=uniform.rvs(loc=-b/2,scale=b,size=N) 
array_y1=uniform.rvs(loc=-h/2,scale=h,size=N)
array_x2=uniform.rvs(loc=-b/2,scale=b,size=N) 
array_y2=uniform.rvs(loc=-h/2,scale=h,size=N) 

M =np.zeros([N,N]) 

foru inxrange(100): ---->This'100'isthe number of experiments.
    forj inxrange(N):
        ifj>0:
            x_A1B1 =array_x2[j]-array_x1[j]
            y_A1B1 =array_y2[j]-array_y1[j]
            x_A1A2 =array_x1[0:j]-array_x1[j]
            y_A1A2 =array_y1[0:j]-array_y1[j]     
            x_A2A1 =-1*x_A1A2
            y_A2A1 =-1*y_A1A2
            x_A2B2 =array_x2[0:j]-array_x1[0:j]
            y_A2B2 =array_y2[0:j]-array_y1[0:j]
            x_A1B2 =array_x2[0:j]-array_x1[j]
            y_A1B2 =array_y2[0:j]-array_y1[j]
            x_A2B1 =array_x2[j]-array_x1[0:j]
            y_A2B1 =array_y2[j]-array_y1[0:j]

            p1 =x_A1B1*y_A1A2 -y_A1B1*x_A1A2
            p2 =x_A1B1*y_A1B2 -y_A1B1*x_A1B2
            p3 =x_A2B2*y_A2B1 -y_A2B2*x_A2B1
            p4 =x_A2B2*y_A2A1 -y_A2B2*x_A2A1

            condition_1=p1*p2
            condition_2=p3*p4                         

            fork inxrange (j):
                ifcondicion_1[k]<=0andcondicion_2[k]<=0:
                    M[j,k]=1

        ifj+1<N+4:
            x_A1B1 =array_x2[j]-array_x1[j]
            y_A1B1 =array_y2[j]-array_y1[j]
            x_A1A2 =array_x1[j+1:]-array_x1[j]
            y_A1A2 =array_y1[j+1:]-array_y1[j]    
            x_A2A1 =-1*x_A1A2
            y_A2A1 =-1*y_A1A2
            x_A2B2 =array_x2[j+1:]-array_x1[j+1:]
            y_A2B2 =array_y2[j+1:]-array_y1[j+1:]
            x_A1B2 =array_x2[j+1:]-array_x1[j]
            y_A1B2 =array_y2[j+1:]-array_y1[j]
            x_A2B1 =array_x2[j]-array_x1[j+1:]
            y_A2B1 =array_y2[j]-array_y1[j+1:]

            p1 =x_A1B1*y_A1A2 -y_A1B1*x_A1A2
            p2 =x_A1B1*y_A1B2 -y_A1B1*x_A1B2
            p3 =x_A2B2*y_A2B1 -y_A2B2*x_A2B1
            p4 =x_A2B2*y_A2A1 -y_A2B2*x_A2A1

            condicion_1=p1*p2
            condicion_2=p3*p4                         

            fork inxrange (N-j-1):
                ifcondicion_1[k]<=0andcondicion_2[k]<=0:
                    M[j,k+j+1]=1

Here, the element Mij=1 if stick i intersect stick j and Mij=0 if not.
How can i optimize my algorithm? Graph theory is usefull in this case?
Waiting for your answers.
Thanks a lot!
Best regards
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20130823/1e8597c1/attachment.html>