原です。

In [ruby-math:00911]

> 正木です。

> 私が試しに使っている gcd,iroot の code を付記します。

まつもとさんを説得するには、特に数学分野では:-)、コードが大事
みたいです。

ちょっと細かいですが、ツッコミを入れさせてください。

> static int
> int_gcd(i, j)
>     int i, j;
> {
>     if (i < 0) i = -i;
>     if (j < 0) j = -j;
>     if (j == 0) return i;
>     return int_gcd(j, i % j);
> }

組み込みのメソッドにするなら、可能なら再帰関数は使わない方が
いいみたいです。

> static VALUE
> rb_gcd(x, y)
>     VALUE x, y;
> {
>     VALUE z;
> 
>     if (FIXNUM_P(x) && FIXNUM_P(y))
>         return INT2FIX(int_gcd(FIX2INT(x), FIX2INT(y)));
>     if (FIXNUM_P(y) && FIX2INT(y) == 0) return rb_big_abs(x);
>     if (FIXNUM_P(y)) return rb_gcd(y, rb_big_modulo(x, y));
>     y = rb_big_abs(y);
>     while (!FIXNUM_P(y)) {
>         z = rb_big_modulo(x, y);
> 	x = y;
> 	y = z;
>     }
>     return rb_gcd(x, y);
> }

同じ値に対して繰り返し同じ関数(FIXNUM_P)を呼ばないようにすべき。


> static int
> int_pow(i, m)
>     int i, m;
> {
>     if (m == 0) return 1;
>     return i * int_pow(i, m - 1);
> }

power は

     if (m == 0) return 1;
     n = int_pow(i, m/2);
     return (i & 1) ? (n * n * i) : (n * n);

みたいなのが速いです。bignum.c の該当部分は

	yy = NUM2LONG(y);
	if (yy > 0) {
	    VALUE z = x;

	    for (;;) {
		yy -= 1;
		if (yy == 0) break;
		while (yy % 2 == 0) {
		    yy /= 2;
		    x = rb_big_mul(x, x);
		}
		z = rb_big_mul(z, x);
	    }
	    return bignorm(z);
	}

です。

> static int
> int_root(i, m)
>     int i, m;
> {
>     int j, k;
>     if (i < 0) rb_raise(rb_eTypeError,"argument is negative");
>     if (i == 0) return 0;
>     j = int_root(i >> m, m) << 1;
>     k = j | 1;
>     if (int_pow(k,m) > i) return j;
>     return k;
> }

TypeError ではなくて RangeError かな。
このロジックで m 乗根も求まるんですね!


Ruby の一分の隙もないようなコードを見ると、組み込みにするとし
たら、やはりそれに準じるような形にしたくなってしまいます。
コーディングスタイルも気になります。(と、このまえあおきさん
につっこまれたばかりの私。(^^;)



後、アルゴリズムの正当性だけでなく、スピードについても組織的に
ベンチマークを取って最速なものを選びましょう。もちろん保守性を
犠牲にしない範囲で。

この話、今回はなんとかまとめてしまいたいですね。