;; ;; Matrices ;; (define \$M.* (cambda \$ms (foldl M.b.* (car ms) (cdr ms)))) (define \$M.b.* (lambda [%m1 %m2] (with-symbols {j} (. m1~#~j m2_j_#)))) (define \$M.*' (cambda \$ms (foldl M.b.*' (car ms) (cdr ms)))) (define \$M.b.*' (lambda [%m1 %m2] (with-symbols {j} (.' m1~#~j m2_j)))) (define \$M.power (lambda [%m \$n] (repeated-squaring M.* m n))) (define \$M.comm (lambda [%m1 %m2] (with-symbols {i j k} (- (. m1~i~j m2_j_k) (. m2~i~j m1_j_k))))) (define \$M.inverse (lambda [%m] (let {[\$d (M.det m)]} (generate-tensor 2#(match m matrix {[ (if (even? (+ %1 %2)) (/ (M.det (M.join A B C D)) d) (* -1 (/ (M.det (M.join A B C D)) d)))]}) (tensor-size m))))) (define \$trace (lambda [%t] (with-symbols {i} (contract + t~i_i)))) (define \$matrix (matcher {[ [math-expr matrix matrix matrix] {[\$tgt (match (tensor-size tgt) (list integer) {[> {[tgt_1_1 tgt_1_[2 n] tgt_[2 m]_1 tgt_[2 m]_[2 n]]}] [_ {}]})]}] [ [math-expr matrix matrix matrix matrix] {[\$tgt (let* {[\$ns (tensor-size tgt)] [\$m (nth 1 ns)] [\$n (nth 2 ns)]} {[tgt_i_j tgt_[1 (- i 1)]_[1 (- j 1)] tgt_[1 (- i 1)]_[(+ j 1) n] tgt_[(+ i 1) m]_[1 (- j 1)] tgt_[(+ i 1) m]_[(+ j 1) n] ]})]}] [,\$val [] {[\$tgt (if (eq? val tgt) {[]} {})]}] [\$ [something] {[\$tgt {tgt}]}] })) (define \$M.join (lambda [%A %B %C %D] (let* {[\$as (tensor-size A)] [\$a1 (nth 1 as)] [\$a2 (nth 2 as)] [\$bs (tensor-size B)] [\$b1 (nth 1 bs)] [\$b2 (nth 2 bs)] [\$cs (tensor-size C)] [\$c1 (nth 1 cs)] [\$c2 (nth 2 cs)] [\$ds (tensor-size D)] [\$d1 (nth 1 ds)] [\$d2 (nth 2 ds)] [\$m1 (max a1 b1)] [\$m2 (max a2 c2)] [\$n1 (max c1 d1)] [\$n2 (max b2 d2)] } (generate-tensor 2#(match [%1 %2] [integer integer] {[[?(lte? \$ a1) ?(lte? \$ a2)] A_%1_%2] [[?(lte? \$ m1) _] B_%1_(- %2 a2)] [[_ ?(lte? \$ m2)] C_(- %1 a1)_%2] [[_ _] D_(- %1 m1)_(- %2 m2)]}) {(+ m1 n1) (+ m2 n2)})))) ;; ;; Determinant ;; (define \$even-and-odd-permutations (lambda [\$n] (let {[[\$es \$os] (even-and-odd-permutations' n)]} [(map 1#(lambda [\$i] (nth i %1)) es) (map 1#(lambda [\$i] (nth i %1)) os)]))) (define \$even-and-odd-permutations' (lambda [\$n] (match n integer {[,1 [{{1}} {}]] [,2 [{{1 2}} {{2 1}}]] [_ (let* {[[\$es \$os] (even-and-odd-permutations' (- n 1))] [\$es' (map 1#{@%1 n} es)] [\$os' (map 1#{@%1 n} os)]} [{@es' @(concat (map (lambda [\$i] (map (permutate i n \$) os')) (between 1 (- n 1)))) } {@os' @(concat (map (lambda [\$i] (map (permutate i n \$) es')) (between 1 (- n 1)))) } ] )]}))) (define \$permutate (lambda [\$x \$y \$xs] (match xs (list eq) {[>>> {@hs y @ms x @ts}] [>>> {@hs x @ms y @ts}]}))) (define \$M.determinant (lambda [%m] (match (tensor-size m) (list integer) {[>> 1] [>> (let {[[\$es \$os] (even-and-odd-permutations' n)]} (- (sum (map (lambda [\$e] (product (map2 (lambda [\$i \$j] m_i_j) (between 1 n) e))) es)) (sum (map (lambda [\$o] (product (map2 (lambda [\$i \$j] m_i_j) (between 1 n) o))) os))))] [_ undefined]}))) (define \$M.det M.determinant) ;;; ;;; Eigenvalues and eigenvectors ;;; (define \$M.eigenvalues (lambda [%m] (match (tensor-size m) (list integer) {[>> (let {[[\$e1 \$e2] (q-f (M.det (T.- m (scalar-to-tensor x {2 2}))) x)]} {e1 e2})] [_ undefined]}))) (define \$M.eigenvectors (lambda [%m] (match (tensor-size m) (list integer) {[>> (let {[[\$e1 \$e2] (q-f (M.det (T.- m (scalar-to-tensor x {2 2}))) x)]} {[e1 (clear-index (T.- m (scalar-to-tensor e1 {2 2}))_i_1)] [e2 (clear-index (T.- m (scalar-to-tensor e2 {2 2}))_i_1)]}) ] [_ undefined]}))) ;; ;; LU decomposition ;; (define \$M.LU (lambda [%x] (match (tensor-size x) (list integer) {[>> (let* {[\$L (generate-tensor 2#(match (compare %1 %2) ordering {[ 0] [ 1] [ b_%1_%2]}) {2 2})] [\$U (generate-tensor 2#(match (compare %1 %2) ordering {[ 0] [_ c_%1_%2]}) {2 2})] [\$m (M.* L U)] [\$ret (solve {[m_1_1 x_1_1 c_1_1] [m_1_2 x_1_2 c_1_2] [m_2_1 x_2_1 b_2_1] [m_2_2 x_2_2 c_2_2]})]} [(substitute ret L) (substitute ret U)])] [>> (let* {[\$L (generate-tensor 2#(match (compare %1 %2) ordering {[ 0] [ 1] [ b_%1_%2]}) {3 3})] [\$U (generate-tensor 2#(match (compare %1 %2) ordering {[ 0] [_ c_%1_%2]}) {3 3})] [\$m (M.* L U)] [\$ret (solve {[m_1_1 x_1_1 c_1_1] [m_1_2 x_1_2 c_1_2] [m_1_3 x_1_3 c_1_3] [m_2_1 x_2_1 b_2_1] [m_2_2 x_2_2 c_2_2] [m_2_3 x_2_3 c_2_3] [m_3_1 x_3_1 b_3_1] [m_3_2 x_3_2 b_3_2] [m_3_3 x_3_3 c_3_3]})]} [(substitute ret L) (substitute ret U)])] [_ undefined]}))) ;; ;; Utility ;; (define \$generate-matrix-from-quadratic-expr (lambda [\$f \$xs] (generate-tensor 2#(coefficient2 f (nth %1 xs) (nth %2 xs)) {(length xs) (length xs)})))