On 2013-05-28, Travis Scrimshaw <tsc...@ucdavis.edu> wrote: > ------=_Part_1805_33045155.1369720622849 > Content-Type: text/plain; charset=ISO-8859-1 > > Hey Stefan, > From just looking at the function, there is only one case in the for > loop that is non-trivial, so I'm thinking the function could be > restructured as: > > cdef bint __exchange(self, long x, long y): > """ > Put element indexed by ``x`` into basis, taking out element ``y``. > Assumptions are that this is a valid basis exchange. > > .. NOTE:: > > Safe for noncommutative rings. > """ > cdef long px, py, r > cdef LeanMatrix A = self._A > px = self._prow[x] > py = self._prow[y] > piv = A.get_unsafe(px, py) > pivi = piv ** (-1) > > # This would be a simple modification of the previous, removing the > for loop > #a = pivi + self._one > #A.row_scale(px, pivi) > #A.set_unsafe(px, py, a) # pivoting without column scaling. > Add extra so column does not need adjusting > # if A and A' are the matrices before and after pivoting, then > # ker[I A] equals ker[I A'] except for the labelling of the columns > #if a: > # A.row_add(px, px, -a) > #A.set_unsafe(px, py, pivi) > > # This is with some other (possible) (micro) optimizations > if pivi + self._one != 0: > A.row_scale(px, -pivi * pivi) > else: > A.row_scale(px, pivi) > A.set_unsafe(px, py, pivi) > > self._prow[y] = px > self._prow[x] = py > BasisExchangeMatroid.__exchange(self, x, y) > > For the second one, I use the fact that > > pivi * A_ij - (pivi + 1) * pivi * A_ij = -pivi^2 * A_ij > > All I've done is unraveled the outcome of the code by following the logic > (interpret that as I haven't tested it), but this should net (again, > provided I didn't make any mistakes) a nice speed increase. > > Also, if you're not really using matrices, such as doing a significant > number of multiplications/using category model/etc., but instead as a nice > data structure,
it looks like something that should be well-known how to implement well, as in the simplex methods for linear programming one does that kind of operations a lot. Also, did you try implementing this using numpy arrays? This would give quite a speedup, IMHO. I imagine there could be a BLAS function (something like daxpy.f) called from numpy, to do the most expensive part (and BLAS can be much faster than C code on good CPUs, due to parallelization). > I'd be more inclined to say to strip out the non > speed-essential methods, making it as lightweight and specific as possible > and converting to proper matrices when needed, and include it into sage > (after perhaps changing the name). How does this sound to everyone? > > Best, > Travis Dima> -- You received this message because you are subscribed to the Google Groups "sage-devel" group. To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscr...@googlegroups.com. To post to this group, send email to sage-devel@googlegroups.com. Visit this group at http://groups.google.com/group/sage-devel?hl=en. For more options, visit https://groups.google.com/groups/opt_out.
