原です。 ちょっと間が空いてしまいました。すいません。Integer成分の determinantの話の続きです。 In [ruby-dev:27700] >けいじゅ＠いしつかです. >>でも、quoを使った方が現状よりは（クレームの頻度の点で）かなりまし >>じゃないですかね。実際にdetを計算したいのはFloatの行列が多いらしい >>し。 > >ましなのは理解しますが, Rationalをrequireするよりましな理由がわかりま >せんが? 私が誤解しているのか、どうも石塚さんの話がぴんとこないのです が、、、私は、現行ではRationalをrequireするだけでなく、更に 各成分をRational化してdetに放り込まなくてはならないのが問題 だと思うのです。 >det計算中Rationalをrequireする場合, >すべてがIntegerの場合, 結果は必ずInitegerになり, >行列中1要素でもFloatが入っている場合, 結果は必ずFloatになり, >Rationalになることはありません. 今回問題になっているのは、整数行列の行列式はどうあつかうべき か、ですよね？ >ですので, (内部の計算がどうなっているかは別として), それぞれの >クラスに応じた自然な結果が出ることが必ず保証できます. 細かい話ですが、整数行列の成分をRational(x, 1)で全て有理数 化してdetを計算させると、結果はRationalですよね。mathn.rbを 導入するとIntegerですけど。 >あと, パフォーマンスのことを気にするなら, 最初っから要素は >Floatにすべきだと思います. 私がquoを使った方がいいと思うのは、整数行列のdetを計算したと きとんでもない値が出るのを防げるからという理由で、ただその一 点に尽きます。 それに対して石塚さんはrational.rbをrequireして各成分、 Rational()でラップしてくれ、と主張しているわけですよね。 私はそれは面倒すぎるのではないかと、そして今後も忘れて 意図しない結果に悩む人が大勢出るのではないかと思います。 それならちょっとかっこ悪くても「quoでFloat」でいいのでは ないでしょうか。 現行をquo式にしたメリットデメリットまとめると、 -------------------------------------------------------------- 成分 | 現行 | quo式 -------------------+--------------------+--------------------- 全てInteger | NG | 戻値:Float | | 誤差が気になることも -------------------+--------------------+--------------------- 全てInteger | NG | 戻値:Rational rationa.rb導入のみ | | 遅い -------------------+--------------------+--------------------- 全てRational | 戻値:Rational | 戻値:Rational | 遅さは気にならない | 遅さは気にならない -------------------+--------------------+--------------------- 全てFloat | 戻値:Float | 戻値:Float | 誤差は気にならない | 誤差は気にならない -------------------+--------------------+--------------------- となります。「遅さは気にならない」とは、Rational計算である のが明白なので遅いことは覚悟が出来ている、、という意味です。 :-) quo式は、rational.rbをrequireしてるか否かで結果が違うのが確 かに気持ちが悪いのですが、それ以上にとにかく現行の"NG"は良く ないのではないか。 quo式で困るとしたら、rational.rbをrequireしている状態で、整 数行列のdetをFloatで高速に欲しい人（いるのか？）ですが、その 人には各成分Float()してからquo式detを計算してもらうことにし ます。 どうでしょう。 更にもう一つ提案なんですが、整数行列のdetの計算には、互除法 を使ったそこそこ速いのがあります： def determinant_e size = row_size - 1 a = to_a det = 1 0.upto(size) do |k| (k + 1).upto(size) do |i| until a[i][k].zero? q = a[k][k] / a[i][k] k.upto(size) do |j| a[k][j] -= a[i][j] * q end a[k], a[i] = a[i], a[k] det *= -1 end end if a[k][k].zero? return 0 else det *= a[k][k] end end det end これは、整数行列をRatinal()してから、現行detを計算するより かなり速いです。次は15次行列の整数成分の行列式を5回計算させ たのにかかった時間です。 determinant : 4.278sec. determinant_e : 0.329sec. 上段が各成分をRational()して現行のdetによる計算です。 拡張ライブラリのrational.soを使っても1.4秒でした。 このdeterminant_eの面白いのは、IntegerとRationalの混合した 成分を持つ行列でも正しく動作する所です。まとめると、 ------------------------------------------------------------ 成分 | quo式 | determinant_e -------------------+----------------------+----------------- 全てInteger | 戻値:Float | 戻値:Integer | 誤差が気になることも | やや遅い -------------------+----------------------+----------------- 全てInteger | 戻値:Rational | 戻値:Integer rationa.rb導入のみ | 遅い | やや遅い -------------------+----------------------+----------------- 全てRational | 戻値:Rational | 戻値:Rational | 気にならない遅さ | 気にならない遅さ -------------------+----------------------+----------------- 全てFloat | 戻値:Float | NG | 気にならない誤差 | -------------------+----------------------+----------------- そこで提案なんですが、Matrix のdetはquo式にし、と同時に determinant_e（det_iという名前もいいですが）も定義しておいて、 【整数行列を扱う場合の注意】 整数行列の行列式を計算させるときは、Matrix#det_iを使いましょ う。Matrix#detは、rational.rbをrequireしていない状態でFloat、 requireしている状態でRationalを返します。 と、マニュアルに記述するのです。 どうでしょう？ ちなみに、determinant_eは[ruby-math:625]のコードをシンプルに したものです。他にも幾つかバリエーションを試してみたのですが、 結局[ruby-math:625]のコードが速いのでこちらの方がいいかもし れません。 テストスクリプトを載せます。 #!/usr/local/bin/ruby require "matrix" require "rational" class Matrix def determinant_e size = row_size - 1 a = to_a det = 1 0.upto(size) do |k| (k + 1).upto(size) do |i| until a[i][k].zero? q = a[k][k] / a[i][k] k.upto(size) do |j| a[k][j] -= a[i][j] * q end a[k], a[i] = a[i], a[k] det *= -1 end end if a[k][k].zero? return 0 else det *= a[k][k] end end det end alias det_e determinant_e #[ruby-math:625] def determinant_625 return 0 unless square? size = row_size - 1 a = to_a det = 1 k = 0 begin if a[k][k].zero? i = k begin return 0 if (i += 1) > size end while a[i][k].zero? a[i], a[k] = a[k], a[i] det *= -1 end (k + 1).upto(size) do |i| q = a[i][k] / a[k][k] k.upto(size) do |j| a[i][j] -= a[k][j] * q end unless a[i][k].zero? a[i], a[k] = a[k], a[i] det *= -1 redo end end det *= a[k][k] end while (k += 1) <= size det end end class Test attr_reader :value def initialize(name) @name = name @value = nil @time = 0.0 end def stime(v = nil) t = Time.now @value = yield @time += Time.now - t if v && v != @value puts "Test Failed: #{@value} != #{v}" exit end end def to_s sprintf("%-30s: %3.3fsec.", @name, @time) end end def test(dim, num, range = 100) tr = Test.new("determinant(rational)") te = Test.new("determinant_e(int)") t1 = Test.new("determinant_e(int & rational)") tm = Test.new("determinant_625(int)") tf = Test.new("determinant(float)") num.times do a = (0...dim).map{(0...dim).map{rand(2*range+1)-range}} m = Matrix[*a] mr = Matrix[*a.map{|b| b.map{|c| Rational(c, 1)}}] m1 = Matrix[*a.map{|b| b.map{|c| rand(2).zero? ? Rational(c, 1) : c}}] mf = Matrix[*a.map{|b| b.map{|c| Float(c)}}] tr.stime{Matrix[*mr].det} te.stime(tr.value){Matrix[*m].determinant_e} t1.stime(tr.value){Matrix[*m1].determinant_e} tm.stime(tr.value){Matrix[*m].determinant_625} tf.stime(){Matrix[*mf].determinant} end puts "--- #{dim}-dim, #{num}times ---" puts tr, te, t1, tm, tf puts end test(5, 500) test(10, 20) test(15, 5) test(20, 2) 以下、テストの結果： --- 5-dim, 500times --- determinant(rational) : 4.503sec. determinant(rational) : 0.647sec. (with rational.so) determinant_e(int) : 1.158sec. determinant_e(int & rational) : 7.313sec. determinant_625(int) : 0.703sec. determinant(float) : 0.193sec. --- 10-dim, 20times --- determinant(rational) : 3.058sec. determinant(rational) : 0.786sec. (with rational.so) determinant_e(int) : 0.421sec. determinant_e(int & rational) : 4.223sec. determinant_625(int) : 0.251sec. determinant(float) : 0.077sec. --- 15-dim, 5times --- determinant(rational) : 4.278sec. determinant(rational) : 1.390sec. (with rational.so) determinant_e(int) : 0.329sec. determinant_e(int & rational) : 5.347sec. determinant_625(int) : 0.238sec. determinant(float) : 0.027sec. --- 20-dim, 2times --- determinant(rational) : 5.384sec. determinant(rational) : 2.208sec. (with rational.so) determinant_e(int) : 0.884sec. determinant_e(int & rational) : 6.511sec. determinant_625(int) : 0.766sec. determinant(float) : 0.024sec.