hui@yrloc.ipsa.reuter.COM (Roger Hui) (05/19/91)
The determinant has been implemented in J this week. The determinant has
rank 2 (applies to matrices), and is recursively defined by expansion
by minors along column 0:
f . g y is (0{"1 y) f . g f . g minors0 y if 0<#y
identity element of g if 0=#y
The dyad f . g plays a key role in this definition, thus motivating the
choice of the monad f . g to denote the determinant. -/ . * is the
ordinary determinant.
The verb minors0 computes the minors of a matrix along column 0. Its
derivation is an interesting story. It was originally defined using
complementary indexing:
i0 =. <&.>&.> @ { @ (;&0) @ i. @ #
minors0 =. i0 { [
m =. 4 4$'abcdefghijklmnop'
m
abcd
efgh
ijkl
mnop
i0 m
+---------+---------+---------+---------+
|+---+---+|+---+---+|+---+---+|+---+---+|
||+-+|+-+|||+-+|+-+|||+-+|+-+|||+-+|+-+||
|||0|||0|||||1|||0|||||2|||0|||||3|||0|||
||+-+|+-+|||+-+|+-+|||+-+|+-+|||+-+|+-+||
|+---+---+|+---+---+|+---+---+|+---+---+|
+---------+---------+---------+---------+
minors0 m
fgh
jkl
nop
bcd
jkl
nop
bcd
fgh
nop
bcd
fgh
jkl
$minors0 m
4 3 3
Subsequently, Eugene McDonnell proposed the following alternative:
ind =. -."1 0~ @ i. @ #
drop0 =. }."1
minors0a =. ind { drop0
The right tine of minors0a drops the first element of each vector (that is,
drops the first column of the matrix), and the left tine is a matrix where
successive rows successively exclude one index. For example:
ind m
1 2 3
0 2 3
0 1 3
0 1 2
drop0 m
bcd
fgh
jkl
nop
(minors0a m) -: minors0 m
1
Ken Iverson suggested yet another alternative. The phrase x f\.y applies
f to outfixes of y obtained by suppressing successive infixes of size x.
For present purposes, we use the verb "same" ([). Thus:
f =. 1 & ([\.)
minors0b =. drop0 @ f
f m
efgh
ijkl
mnop
abcd
ijkl
mnop
abcd
efgh
mnop
abcd
efgh
ijkl
(minors0b m) -: minors0 m
1
Moreover, this technique can be used to compute the major minors (i.e.
minors of size one less than the size of the matrix), an immediate
precursor to the classical adjoint:
box =. <"2
box minors m
+---+---+---+---+
|fgh|egh|efh|efg|
|jkl|ikl|ijl|ijk|
|nop|mop|mnp|mno|
+---+---+---+---+
|bcd|acd|abd|abc|
|jkl|ikl|ijl|ijk|
|nop|mop|mnp|mno|
+---+---+---+---+
|bcd|acd|abd|abc|
|fgh|egh|efh|efg|
|nop|mop|mnp|mno|
+---+---+---+---+
|bcd|acd|abd|abc|
|fgh|egh|efh|efg|
|jkl|ikl|ijl|ijk|
+---+---+---+---+
f"1"_ applies f to each rank 1 object (i.e. each vector) in the entire
argument; 0 2 1 3&|: transposes axes 1 and 2; and "box" boxes each matrix
to produce a compact display.
The following is a translation of definitions on p. 58 of "Some Uses of
{ and }", APL87 Conference Proceedings. Here, we use the technique of
outfixes in place of complementary indexing, to compute the minors.
(f and minors are repeated from above.)
f =. 1 & ([\.)
minors =. (0 2 1 3&|:) @ (f"1"_) @ f
det =. -/ . *
checker =. */~ @ ($&1 _1) @ #
adjoint =. |: @ (checker * det@minors)
matinv =. adjoint % det
[ n =. _40+?4 4$100
_27 35 5 13
_19 _36 27 27
53 _2 11 43
_37 _35 12 27
det@minors n
_21546 _49272 29140 _21394
_10227 4362 _170024 _55897
20979 _3342 _24952 _44171
33264 _33408 144316 76364
checker n
1 _1 1 _1
_1 1 _1 1
1 _1 1 _1
_1 1 _1 1
adjoint n
_21546 10227 20979 _33264
49272 4362 3342 _33408
29140 170024 _24952 _144316
21394 _55897 44171 76364
det n
2730084
n +/ . * matinv n
1 6.77626e_21 _1.35525e_20 0
_2.71051e_20 1 _4.06576e_20 1.35525e_20
2.71051e_20 2.71051e_20 1 0
_2.71051e_20 0 _5.42101e_20 1
-----------------------------------------------------------------
Roger Hui
Iverson Software Inc., 33 Major Street, Toronto, Ontario M5S 2K9
(416) 925 6096sam@kalessin.jpl.nasa.gov (Sam Sirlin) (05/21/91)
In article <1991May19.053021.7840@yrloc.ipsa.reuter.COM>, hui@yrloc.ipsa.reuter.COM (Roger Hui) writes: |> The determinant has been implemented in J this week. The determinant has |> rank 2 (applies to matrices), and is recursively defined by expansion |> by minors along column 0: |> This sounds interesting, but if I was given a matrix and told to take its determinant, I would use LU or SV decomposition if the matrix was above about 6th order. Isn't the minor calculation ill conditioned? If so then isn't a built in determinant using minors missleading to naive programmers? -- Sam Sirlin Jet Propulsion Laboratory sam@kalessin.jpl.nasa.gov
hui@yrloc.ipsa.reuter.COM (Roger Hui) (05/22/91)
Sam Sirlin writes: >This sounds interesting, but if I was given a matrix and told to take its >determinant, I would use LU or SV decomposition if the matrix was above >about 6th order. Isn't the minor calculation ill conditioned? If so then >isn't a built in determinant using minors missleading to naive programmers? The system will recognize the ordinary determinant -/ . * as a special case, and use an equivalent but O(n^3) algorithm to compute it. Having a determinant operator encourages experimentation with generalized determinants. For example: + / . * permanent <./ . + minimum cost in assignment problem +./ . *. whether a Latin square is covered + / . *. how many Latin squares are covered [ . (,&.>) try it All except the last are from K.E. Iverson, "Determinant-like Functions Produced by the Dot Operator", Sharp APL Technical Notes 42, 1982 4 1. ----------------------------------------------------------------- Roger Hui Iverson Software Inc., 33 Major Street, Toronto, Ontario M5S 2K9 (416) 925 6096