トップ «前の日記(2007-04-19) 最新 次の日記(2007-05-10)» 編集

U-memo

2006|11|12|
2007|01|02|03|04|05|06|07|08|09|10|11|12|
2008|01|02|03|04|05|06|08|
2009|08|10|
2010|02|03|
2011|11|12|
2012|04|
2016|02|
All= / Today= / Yesterday=

2007-04-20 [長年日記]

_ [Golf][OCaml] あなごる 逆行列

(やはり)kskさんに負けた。

Gauss-Jordan 法あたりで求めるとして、

  • 普通にやると倍精度浮動小数でも誤差が大きく出る
  • 途中で対角要素が 0 のやつがいるので pivot 操作相当が必要

への対処が問題といえる。

387B 版は浮動小数をあきらめて Num (無限精度有理数)を使い、pivot 操作として最下行を足し込むようにした。(400B台後半の頃はまじめにやってたが)

post-mortem (359B) では浮動小数点にして誤差をごまかした。

初期の版数のもさらしておこう。

#load "nums.cma"
open Num
let a=Array.make_matrix 9 18(Int 0);;
for i=0to 8do
 for j=0to 8do
  Scanf.scanf "%d "(fun x->a.(i).(j)<-Int x)
 done;
 a.(i).(i+9)<-(Int 1)
done;;
for i=0to 8do
 if a.(i).(i)=Int 0 then
  for j=i+1to 8do
   for k=0to 17do
    a.(i).(k)<-a.(i).(k)+/a.(j).(k)
   done
  done;
 for j=0to 8do
  let u=a.(i).(i)in
  let c=a.(j).(i)//u in
  for k=0to 17do
   a.(j).(k)<-if i=j then a.(j).(k)//u
              else a.(j).(k)-/a.(i).(k)*/c
  done
 done
done;
for i=0to 8do
 for j=9to 16do
  Printf.printf "%d "(int_of_num a.(i).(j))
 done;
 Printf.printf "%d\n"(int_of_num a.(i).(17))
done
本日のツッコミ(全2件) [ツッコミを入れる]
_ ksk (2007-04-21 00:12)

m.ukaiのさんのコードを基に再投稿させていただきました.

_ うかい (2007-04-27 12:05)

おー素晴らしいですー