librarylibs/Matrix/src/Matrix.xtl

Matrix: matrices -- identity, diagonal, trace, product, determinant, inverse, solving linear systems. Import it with an alias of your choice: "mx:" u_se< "Matrix". Put libs/Matrix/src on XETAL_PATH ("just path"); the reference is libs/Matrix/docs. Names with l: are exported; those under h: are private to this file.

A matrix is a rank-2 array, rows first. The numerical functions (d_et, i_nverse, s_olve) work in Floats, by Gaussian elimination with partial pivoting; compare their results with Check's k:n_ear. The transpose is the built-in o_\.

source

Basics

ˡi̲dentity : Num a => Int -> a

function · line 15

i_dentity n: the n by n identity matrix (1 then n zeros, cycled).

ˡi̲dentity ← { n → (n c̲at n) r̲eshape 1 c̲at n r̲eshape 0 }

Used in: ˡi̲nverse

ˡd̲iag : a -> a

function · line 18

d_iag m: the items on the main diagonal.

ˡd̲iag ← { m →
  k ← 'm̲in r̲/ s̲hape m
  (1 + (o̲ffsets k) × 1 + f̲irst r̲ev s̲hape m) s̲elect r̲avel m
}

Used in: ˡt̲race

ˡt̲race : Num a => a -> a

function · line 24

t_race m: the sum of the diagonal.

ˡt̲race ← { m → '+ r̲/ ˡd̲iag m }

ˡm̲ul : Num a => a -> a -> a

function · line 27

a m_ul b: the matrix product (a matrix and a vector work too).

ˡm̲ul ← { a b → a '+ '× i̲nner b }

Used in: b, r

Solving

ʰs̲wapping : Int -> Int -> Int -> Int

function (private) · line 32

The permutation of 1..n that swaps i and j.

ʰs̲wapping ← { n i j →
  r ← r̲ange n
  r + ((r = i) × j − i) + (r = j) × i − j
}

Used in: ʰg̲j, ʰd̲etF

ʰp̲ivot : Num a => Int -> a -> Int

function (private) · line 38

The row (k..n) whose item in column k is largest in size.

ʰp̲ivot ← { k m →
  c ← a̲bs (k − 1) d̲rop k s̲elect₂ m
  k + (f̲irst g̲rade n̲eg c) − 1
}

Used in: ʰg̲j, ʰd̲etF

ʰg̲j : Int -> Int -> Float -> Float

function (private) · line 46

Gauss-Jordan on m (Float) for columns k..n of its first n: each column's pivot row moved up, scaled to 1, and the column cleared from every other row.

ʰg̲j ← { n k m →
  k > n ? m
  m ← ((n ʰs̲wapping k)_ k ʰp̲ivot m) s̲elect m
  r ← r̲ange t̲ally m
  v ← f̲irst k s̲elect k s̲elect₂ m              ⍝ the pivot
  0.000000000001 > a̲bs v ? @ p̲anic< "the matrix is singular: it has no inverse, and its system no unique solution"
  row ← (k s̲elect m) ÷ v
  c ← (k s̲elect₂ m) − f̲loat r = k
  ((n ʰg̲j k + 1)_ m − c '× t̲able row)
}

p̲anic< expands to

(⎕P̲ANIC ("the matrix is singular: it has no inverse, and its system no unique solution"))

Used in: ʰg̲j, ˡs̲olve

ˡs̲olve : Num a => a -> a -> Float

function · line 60

b s_olve a: x with a m_ul x = b, for a square a; b is a vector or a matrix of right-hand sides (one per column). A singular a divides by zero (error[division-by-zero]).

ˡs̲olve ← { b a →
  n ← t̲ally a
  x ← n d̲rop₂ ((n ʰg̲j 1)_ f̲loat a c̲at₂ b)
  (s̲hape b) r̲eshape x
}

Used in: ˡi̲nverse, b, c

ˡi̲nverse : Num a => a -> Float

function · line 67

i_nverse a: the inverse of a square matrix.

ˡi̲nverse ← { a → (ˡi̲dentity t̲ally a) ˡs̲olve a }

Determinant

ʰd̲etF : Float -> Float

function (private) · line 73

The determinant of a Float matrix m, by elimination: the first pivot times the determinant of what is left, negated for a row swap.

ʰd̲etF ← { m →
  n ← t̲ally m
  n = 1 ? f̲irst r̲avel m
  p ← 1 ʰp̲ivot m
  m ← ((n ʰs̲wapping 1)_ p) s̲elect m
  v ← f̲irst r̲avel m
  v = 0 ? 0.0
  s ← f̲loat (p = 1) − p ≠ 1
  r ← (1 d̲rop₂ 1 d̲rop m) − ((1 d̲rop f̲irst₂ m) ÷ v) '× t̲able 1 d̲rop f̲irst m
  s × v × ʰd̲etF r
}

Used in: ʰd̲etF, ˡd̲et

ˡd̲et : Num a => a -> Float

function · line 86

d_et a: the determinant of a square matrix.

ˡd̲et ← { a → ʰd̲etF f̲loat a }