正木です。


|[ruby-math:00672] Re: power
|From: sinara / blade.nagaokaut.ac.jp

|       z = Matrix.identity(self.column_size)
|       while other != 0
|         z *= x if other & 1 != 0
|         x *= x
|         other >>= 1
|       end
|       z

|これが一番きれいなのですが、必ず Identity * x というしなく
|てもいい計算をしてしまいます。そこで先のは n = other - 1 
|みたいなことをしているわけです。
|そこのところ、きれいな解決ってあるでしょうか?

私は次のようにしています。

-----
class Numeric
  def unit
    1
  end
end
class Matrix
  def unit
    raise "not square matrix" if not square?
    Matrix.identity(row_size)
  end
end

Power=Function.new{"|x| 
  Sequence.new([x.unit,x],x){'|n,x| z=self[n>>1]**2;(n&1==0)? z: z*x '}
  "}
----
# Power[x][n]=x**n

普通の関数形式で書けば

def power(x,n)
  return x.unit if n==0
  return x if n==1
  z=power(x,n>>1)**2
  (n&1==0)? z : z*x
end

Ruby では再帰定義を使っても速度的には問題は無いような印象を持っていますが
本当のところはどうなんでしょうか?

因みに上記の花谷さんの code を同じ流儀で書けば

def power(x,n)
  return power(x.inverse,-n) if n<0
  return x.unit if n==0
  return x if n==1
  power(x,n&1)*power(x*x,n>>1)
end

となります。

ところで、冪乗を求める algorithm をそれに使われる乗算の回数だけで評価する
事にすると、上の code は必ずしも最適ではありません。
x**n を計算するのに最低限必要な乗算の回数を f(n) とすると、f(n) 自身を
求める方法は分かりませんが
f(0)=0
f(1)=0
f(a*b) <= f(a)+f(b)
f(a+b) <= f(a)+f(b)+1
から、次の F[n]( >=f(n)) を代わりに使う事にします。
F=Sequence.new([0,0]){" |n|
  ((1..n/2).collect{|i| F[i]+F[n-i]+1}+
  Prime.divisors(n).find_all{|j| 1<j && j**2<=n}.collect{|i| F[i]+F[n/i]}).min
  "}
上の code での乗算の回数は
F1=Sequence.new([0,0]){" |n| self[n>>1]+1+(n&1) "}
で求められるので
for i in 1..n
  p [i,F[i],F1[i]] if F1[i]-F[i] >= m
end
n=50, m=1 としてみると
[15, 5, 6]
[27, 6, 7]
[30, 6, 7]
[31, 7, 8]
[39, 7, 8]
[45, 7, 8]
n=100,m=2 では
[63, 8, 10]
[95, 9, 11]
n=1000,m=4 では
[255, 10, 14]
[495, 11, 15]
[510, 11, 15]
[511, 12, 16]
[735, 12, 16]
[765, 12, 16]
[891, 12, 16]
[959, 13, 17]
[975, 12, 16]
[990, 12, 16]
[991, 13, 17]

となり、これ位なら苦労して最適化する必要は無さそうです。

最初に差が出る 15 で言えば
15=3*(2*2+1)
が最適な分解で
15=(2*(2*(2+1)+1)+1
が我々の algorithm に相当します。