ruby-dev

まつもと ゆきひろです

rubyforgeバグレポートに以下のようなレポートが来てますが、私
には判断できません。どうでしょう。石塚さん？

# コメントは編集してあります。パッチは(3つある)最終版のもの
# です。

---
There is a problem with the Matrix.inverse method.
Given the following matrix

    m = Matrix[[0.417787968829298, 0.542037605627772, 0.729142268138943],
    [0.0, 6.12303176911189e-17, 1.0],
    [-0.542037605627772, 0.417787968829298, -2.55812900589452e-17]]

the result of

    m.inverse * m

is

    [[1.0, 0.0, 0.650423523448542],
    [5.55111512312578e-17, 1.0, 0.84385869292008],
    [0.0, 9.24446373305873e-33, 1.0]]

but it should be

    [[1.0, 0.0, 0.0],
    [0.0, 1.0, 0.0],
    [0.0, 0.0, 1.0]]

or close to this!

----
The patch below should "fix" the inverse method. It adds partial pivoting to the Gauss-Jordan algorithm, making it stable. (I used Ruby 1.8.5 (2006-08-22).)

I don't think there is an error in the algorithm, and the matrix is not ill-conditioned (a condition number of about 2.5 isn't so bad).

Looks to me like BigDecimal uses Gaussian elimination with partial pivoting too, or at least that's what the comments say. It includes some additional scaling, but this may be specific to BigDecimal (precision works differently for BigDecimals and floats).

Peter

--- lib/matrix.rb.old	2007-03-19 01:46:44.000000000 +0100
+++ lib/matrix.rb	2007-03-19 03:26:38.000000000 +0100
@@ -590,15 +590,21 @@
     a = src.to_a
     
     for k in 0..size
-      if (akk = a[k][k]) == 0
-        i = k
-        begin
-          Matrix.Raise ErrNotRegular if (i += 1) > size
-        end while a[i][k] == 0
+      i = k
+      akk = a[k][k].abs
+      for j in (k+1)..size
+        v = a[j][k].abs
+        if v > akk
+          i = j
+          akk = v
+        end
+      end
+      Matrix.Raise ErrNotRegular if akk == 0
+      if i != k
         a[i], a[k] = a[k], a[i]
         @rows[i], @rows[k] = @rows[k], @rows[i]
-        akk = a[k][k]
       end
+      akk = a[k][k]
       
       for i in 0 .. size
         next if i == k

Thread

Prev Next

In This Thread

Prev Next

[#30454] [ANN] CVSサービスの廃止 — Shugo Maeda <shugo@...>

[#30455] a[] += 1 で SEGV — "Hiroki Motoyoshi" <himotoyoshi.ml@...>

[#30465] Ruby 1.8.6 delayed for seven days — "Akinori MUSHA" <knu@...>

[#30476] Re: [ruby-cvs:19211] Ruby:r11975: * lib/fileutils.rb (touch): last #touch change causes error when :mtime option was not given. — Nobuyoshi Nakada <nobu@...>

[#30477] ThreadError on shell.rb — "Akira ODA" <oda.org@...>

[#30482] Ruby 1.8.6 preview3 has been released — "Akinori MUSHA" <knu@...>

[#30486] str.lines[1..-1] — Minero Aoki <aamine@...>

[#30487] Time#to_s format — OHARA Shigeki <os@...>

[#30500] MD5.md5 in RubyGems — Kazuhiro Yoshida <moriq@...>

[#30504] num_tとthread_initという名前について — "Yutaka Kanemoto" <kinpoco@...>

[#30505] [BUG?] SIGTERM の取り扱い — Hidetoshi NAGAI <nagai@...>

[#30507] monitor.rb for ruby-1.9 — keiju@... (Keiju ISHITSUKA)

[#30524] DL::Symbol#call dumps core — sheepman <sheepman@...>

[#30530] dl/import is not thread safe — sheepman <sheepman@...>

[#30531] DL::PtrData#to_a dumps core — sheepman <sheepman@...>

[#30532] DL::PtrData#inspect — sheepman <sheepman@...>

[#30536] Ruby 1.8.6 will be released in 24 hours — "Akinori MUSHA" <knu@...>

[#30537] Microsoft-IIS 上で cgi.rb がレスポンス行を２重に出力する — MIYASAKA Masaru <alkaid@...>

[#30547] dl memory leak — sheepman <sheepman@...>

[#30549] [ANN] Ruby 1.8.6 has been released — "Akinori MUSHA" <knu@...>

[#30551] rb_rescue2 with vm — wanabe <s.wanabe@...>

[#30553] help: lib/shell for ruby 1.9 — keiju@... (Keiju ISHITSUKA)

[#30554] rb_ary2cary in dl.c has a bug — sheepman <sheepman@...>

[#30556] ruby 1.8.6のdateとrails — akira yamada / やまだあきら <akira@...>

[#30563] ext/openssl のメンテナについて — "Akinori MUSHA" <knu@...>

[#30564] release code and revision number — Nobuyoshi Nakada <nobu@...>

[#30575] RSS::Parser.parseでSecurityError — Kazuhiro NISHIYAMA <zn@...>

[#30576] Changelog-1.8.6 — moonwolf@...

[#30581] [ruby-core:10580] Thread後のexecでENOTSUPPエラー対応 — "H.Holon" <holon@...>

[#30586] Etc.getgrgid が getgid の値に固定されています — 植田裕之 <ueda@...>

[#30593] [BUG] waitpid for ruby-trunk — keiju@... (Keiju ISHITSUKA)

[#30612] stdout flush in fork — keiju@... (Keiju ISHITSUKA)

[#30615] Etc モジュールのテスト — 植田 裕之 <ueda@...>

[#30618] possible bug in rb_detach_process() — "KUBO Takehiro" <kubo@...>

[#30620] DL.callback as free function causes errors — sheepman <sheepman@...>

[#30622] dl when $SAFE = 1 — sheepman <sheepman@...>

[#30624] test/rdoc occurs errors. — NAKANO Kouichi <curving@...>

[#30627] [ ruby-Bugs-9360 ] Matrix inverse problem — Yukihiro Matsumoto <matz@...>

[#30629] at_exit — keiju@... (Keiju ISHITSUKA)

[#30637] strlcpy and strlcat — "U.Nakamura" <usa@...>

[#30641] Re: [ruby-cvs:19321] Ruby:r12085: fixed the problem reported in — Urabe Shyouhei <shyouhei@...>

[#30648] YARVでのalias名の取得方法 — Yukihiro Matsumoto <matz@...>

[#30653] [BUG] queue 0x424010 freed with live thread(s) waiting — Masatoshi SEKI <m_seki@...>

[#30655] win32/registry wtime — moonwolf@...

[#30657] ロード済み拡張ライブラリの判別 — arton <artonx@...>

[#30667] $! assignment does not work properly — Yukihiro Matsumoto <matz@...>

[#30671] invalid CFP in rb_protect — wanabe <s.wanabe@...>

[#30672] fail to copy argv on Win32? — wanabe <s.wanabe@...>

[#30674] socket blocking on mingw32 — wanabe <s.wanabe@...>

[#30675] coredump on NetBSD with trunk — "U.Nakamura" <usa@...>

[#30676] excel is not closed in runtime using excel.Quit() — "Md. Ruhul Amin" <amin@...>

[ruby-dev:30627] [ ruby-Bugs-9360 ] Matrix inverse problem

Thread

In This Thread

[#30615] Etc モジュールのテスト — 植田裕之 <ueda@...>