savage_core 0.2.0

A primitive computer algebra system (library)
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57

// SPDX-License-Identifier: AGPL-3.0-or-later
// Copyright (C) 2021-2022  Philipp Emanuel Weidmann <pew@worldwidemann.com>

use permutohedron::heap_recursive;
use savage_macros::function;

use crate::{expression::Expression, functions::SquareMatrix, helpers::*};

#[function(
    name = "det",
    description = "determinant of a square matrix",
    examples = r#"[
        ("det([[1, 2], [3, 4]])", "-2"),
        ("det([[a, b], [c, d]])", "a * d - b * c"),
        ("det([])", "1"),
    ]"#,
    categories = r#"[
        "linear algebra",
    ]"#
)]
fn determinant(matrix: SquareMatrix) -> Expression {
    if matrix.is_empty() {
        return int(1);
    }

    let mut indices = (0..matrix.nrows()).collect::<Vec<usize>>();

    let mut products = Vec::new();

    heap_recursive(indices.as_mut_slice(), |permutation| {
        products.push(
            (0..matrix.nrows())
                .map(|i| matrix[(i, permutation[i])].clone())
                .reduce(|a, b| a * b)
                .unwrap(),
        );
    });

    // The first permutation generated by Heap's algorithm is the identity,
    // which has positive sign...
    let mut positive = true;

    products
        .into_iter()
        .reduce(|a, b| {
            // ... and every following permutation differs from its predecessor
            // by exactly one transposition, which flips the sign.
            positive = !positive;

            if positive {
                a + b
            } else {
                a - b
            }
        })
        .unwrap()
}