From: mame@... Date: 2018-10-20T07:57:24+00:00 Subject: [ruby-core:89491] [Ruby trunk Bug#15233][Assigned] Speeding Up Matrix Powers Issue #15233 has been updated by mame (Yusuke Endoh). Status changed from Open to Assigned Assignee set to marcandre (Marc-Andre Lafortune) Interesting. The current implementation accumulates the value in LSB-to-MSB order, and SymPy's implementation does in MSB-to-LSB because. The former is slower than the latter because the latter multiplies the original matrix (which is smaller than an accumulated matrix for each digit). @marcandre, could you review this? ```diff diff --git a/lib/matrix.rb b/lib/matrix.rb index 62852bdad0..fd9115cfad 100644 --- a/lib/matrix.rb +++ b/lib/matrix.rb @@ -1086,12 +1086,12 @@ def **(other) return self.class.identity(self.column_count) if other == 0 other = -other end - z = nil - loop do - z = z ? z * x : x if other[0] == 1 - return z if (other >>= 1).zero? - x *= x + z = self + (other.bit_length - 2).downto(0) do |i| + z *= z + z *= self if other[i] == 1 end + return z when Numeric v, d, v_inv = eigensystem v * self.class.diagonal(*d.each(:diagonal).map{|e| e ** other}) * v_inv ``` ---------------------------------------- Bug #15233: Speeding Up Matrix Powers https://bugs.ruby-lang.org/issues/15233#change-74538 * Author: octonion (Christopher Long) * Status: Assigned * Priority: Normal * Assignee: marcandre (Marc-Andre Lafortune) * Target version: * ruby -v: ruby 2.5.1p57 (2018-03-29 revision 63029) [x86_64-linux-gnu] * Backport: 2.3: UNKNOWN, 2.4: UNKNOWN, 2.5: UNKNOWN ---------------------------------------- I've been looking at SymPy's slow matrix power speed, and noticed that replacing Ruby's matrix power code with the equivalent code from SymPy makes it significantly faster. As this is a recursively defined function, presumably it can be made even faster with a proper iterative version. Base: ~~~ ruby require 'matrix' Q = Matrix[[0,0,1,0,2,0], [0,0,0,1,0,0], [1,0,0,1,0,2], [0,2,1,0,0,0], [1,0,0,0,0,0], [0,0,1,0,0,0]] r = Vector[1,0,0,0,0,0] p = (Q**100000000*r).sum puts p.to_s.size ~~~ Output: ~~~ time ruby main.rb 35950264 real 1m42.250s user 1m41.050s sys 0m1.184s ~~~ Modified: ~~~ ruby require 'matrix' def pow(a,num) if (num==1) return a end if (num%2)==1 return a*pow(a,num-1) end ret = pow(a,(num/2)) ret*ret end Q = Matrix[[0,0,1,0,2,0], [0,0,0,1,0,0], [1,0,0,1,0,2], [0,2,1,0,0,0], [1,0,0,0,0,0], [0,0,1,0,0,0]] r = Vector[1,0,0,0,0,0] p = (pow(Q,100000000)*r).sum puts p.to_s.size ~~~ Output: ~~~ time ruby main.rb 35950264 real 1m9.489s user 1m8.661s sys 0m0.828s ~~~ -- https://bugs.ruby-lang.org/ Unsubscribe: