\ Matrix Determinant
\
\ Copyright (C) 2011 David K端hling <dvdkhlng TA gmx TOD de>
\
\ License: GPL3 or later, NO WARRANTY
\
\ Created: Aug 2011
require ./imatrix.fs
require ./mod.fs
\ todo: factor similar to SNF.
: }}rdet { n M{{ -- result }
\ note: destroys the contents of M{{ during computation
1 { result }
n 0 ?DO
\ search pivot (i.e. swap rows to ensure m(i,i)>0)
n I ?DO
M{{ I J }} @ IF
I J <> IF
M{{ I 0 }} M{{ J 0 }} n iexchange
result rnegate TO result
THEN
LEAVE
THEN
LOOP
\ read pivot, adjust result, precompute inverse
M{{ I I }} @ dup result r* TO result ( s: pivot)
dup 0= IF \ zero pivot
UNLOOP EXIT \ -> abort with determinant zero
THEN
rinv { p^-1 }
\ zero columns I in rows below
n I 1+ ?DO
M{{ I J }} @ ?dup IF
p^-1 r*
n J ?DO
DUP M{{ K I }} @ r* M{{ J I }} r-!
LOOP
drop
THEN
LOOP
LOOP
result ;
\ Customize Emacs
0 [IF]
Local Variables:
compile-command: "gforth rdet-test.fs -e bye"
End:
[THEN]