Just for fun I wrote an equivalent program in C and tested the Sage function and an equivalent C program on 1 instance of the problem.
=== $ time ./a.out test real 0m0.366s user 0m0.306s sys 0m0.060s === And Sage === $ time sage test.sage test real 2m6.285s user 2m5.256s sys 0m0.858s === I thought that may be interesting to you as well. Best, Jernej On Wednesday, 3 December 2014 22:33:46 UTC+1, Jernej wrote: > > 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 sage-support+unsubscr...@googlegroups.com. To post to this group, send email to sage-support@googlegroups.com. Visit this group at http://groups.google.com/group/sage-support. For more options, visit https://groups.google.com/d/optout.