Dear sage-support,

I have stumbled into a performance bottleneck in one of my Sage programs. I 
would like to share the relevant problem  here in hope anyone has a 
constructive suggestion for optimizing the given program.

I am given a n x n, (0,1) matrix C  where n < 20. C has up to 30% of 
nonzero elements. 

I define the inner product <u,v> = u (2I - C) v^t and need to compute the 
set S of all (0,1) n-vectors u such that <u,u> = 2 and <u,j> = -1, where j 
is the vector with all ones. Finally I need to record all pairs u,v from S 
such that <u,v> == -1 or 0.

The one optimization that I use is to cache the product u*(2I-C). Other 
than that I don't see any other way to optimize the program hence it looks 
like

====
global vectors
field = RR # this looks like the fastest option

def init_vectors(n):
    global vectors
    vectors = [matrix(field,b) for b in product(*[[0,1]]*n)]
def foo(C,output): 
        
    global vectors

    n = C.nrows()
    # we compute the inverse in QQ just to gain a bit of numerical 
stability 
    D  = (2*identity_matrix(n)-C).inverse()
    D = Matrix(field,D)

    j = matrix(field, tuple([1 for _ in xrange(n)]))
    cur = 1 

    vec2int = {}
    cache = {}

    for bm in vectors:
        val = bm*D
        if abs( (val*bm.transpose())[0,0]-2) <= 0.000001 and 
abs((val*j.transpose())[0,0]+1) <= 0.000001: 
            vec2int[cur] = bm
            output.write(str(el for el in bm[0])+'\n')
            cache[cur] = val 
            cur+=1

    for i in xrange(1, cur):
        for j in xrange(i+1, cur):
            iv = (cache[i]*vec2int[j].transpose())[0,0]
            if abs(iv) <= 0.0000001 or abs(iv+1) <= 0.000001: 
               output.write(str(i) + ' ' + str(j) + '\n')
===

For convenience's sake I have attached the result of %prun on one instance 
of the problem. From the output of prun I suppose it would make sense to 
cache the transpose of the vectors as well, though I am not sure this is 
going to make the crucial difference.

On average it takes 200s to process one matrix and given that I have 
millions of such matrices it would take years to process all the input.  

Other than rewriting the whole thing in C I currently do not see any viable 
solution. Hence I am wondering: Do you guys happen to see any clever 
optimization? Is there a way to compute the named product more efficiently? 

Best,

Jernej

-- 
You received this message because you are subscribed to the Google Groups 
"sage-support" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to [email protected].
To post to this group, send email to [email protected].
Visit this group at http://groups.google.com/group/sage-support.
For more options, visit https://groups.google.com/d/optout.

Attachment: prof.out
Description: Binary data

Reply via email to