Determinant and permanent

From Rosetta Code
Task
Determinant and permanent
You are encouraged to solve this task according to the task description, using any language you may know.

For a given matrix, return the determinant and the permanent of the matrix.

The determinant is given by

Failed to parse (Conversion error. Server ("cli") reported: "[INVALID]"): {\displaystyle \det(A) = \sum_\sigma\sgn(\sigma)\prod_{i=1}^n M_{i,\sigma_i}}

while the permanent is given by

In both cases the sum is over the permutations of the permutations of 1, 2, ..., n. (A permutation's sign is 1 if there are an even number of inversions and -1 otherwise; see parity of a permutation.)

More efficient algorithms for the determinant are known: LU decomposition, see for example wp:LU decomposition#Computing the determinant. Efficient methods for calculating the permanent are not known.


Related task



11l[edit]

Translation of: Nim
F s_permutations(seq)
   V items = [[Int]()]
   L(j) seq
      [[Int]] new_items
      L(item) items
         V i = L.index
         I i % 2
            new_items [+]= (0 .. item.len).map(i -> @item[0 .< i] [+] [@j] [+] @item[i ..])
         E
            new_items [+]= (item.len .< -1).step(-1).map(i -> @item[0 .< i] [+] [@j] [+] @item[i ..])
      items = new_items

   R enumerate(items).map((i, item) -> (item, I i % 2 {-1} E 1))

F det(a)
   V result = 0.0
   L(sigma, _sign_) s_permutations(Array(0 .< a.len))
      V x = Float(_sign_)
      L(i) 0 .< a.len
         x *= a[i][sigma[i]]
      result += x
   R result

F perm(a)
   V result = 0.0
   L(sigma, _sign_) s_permutations(Array(0 .< a.len))
      V x = 1.0
      L(i) 0 .< a.len
         x *= a[i][sigma[i]]
      result += x
   R result

V a = [[1.0, 2.0],
       [3.0, 4.0]]

V b = [[Float( 1),  2,  3,  4],
       [Float( 4),  5,  6,  7],
       [Float( 7),  8,  9, 10],
       [Float(10), 11, 12, 13]]

V c = [[Float( 0),  1,  2,  3,  4],
       [Float( 5),  6,  7,  8,  9],
       [Float(10), 11, 12, 13, 14],
       [Float(15), 16, 17, 18, 19],
       [Float(20), 21, 22, 23, 24]]

print(‘perm: ’perm(a)‘ det: ’det(a))
print(‘perm: ’perm(b)‘ det: ’det(b))
print(‘perm: ’perm(c)‘ det: ’det(c))
Output:
perm: 10 det: -2
perm: 29556 det: 0
perm: 6778800 det: 0

360 Assembly[edit]

For maximum compatibility, this program uses only the basic instruction set (S/360) and two ASSIST macros (XDECO,XPRNT) to keep it as short as possible. It works on OS/360 family (MVS,z/OS), on DOS/360 family (z/VSE) use GETVIS,FREEVIS instead of GETMAIN,FREEMAIN.

*        Matrix arithmetic         13/05/2016
MATARI   START
         STM    R14,R12,12(R13)    save caller's registers
         LR     R12,R15            set R12 as base register
         USING  MATARI,R12         notify assembler
         LA     R11,SAVEAREA       get the address of my savearea
         ST     R13,4(R11)         save caller's savearea pointer
         ST     R11,8(R13)         save my savearea pointer
         LR     R13,R11            set R13 to point to my savearea
         LA     R1,TT              @tt
         BAL    R14,DETER          call deter(tt)
         LR     R2,R0              R2=deter(tt)
         LR     R3,R1              R3=perm(tt)
         XDECO  R2,PG1+12          edit determinant
         XPRNT  PG1,80             print determinant
         XDECO  R3,PG2+12          edit permanent
         XPRNT  PG2,80             print permanent
EXITALL  L      R13,SAVEAREA+4     restore caller's savearea address
         LM     R14,R12,12(R13)    restore caller's registers
         XR     R15,R15            set return code to 0
         BR     R14                return to caller
SAVEAREA DS     18F                main savearea
TT       DC     F'3'               matrix size 
         DC     F'2',F'9',F'4',F'7',F'5',F'3',F'6',F'1',F'8' <==input
PG1      DC     CL80'determinant='
PG2      DC     CL80'permanent='
XDEC     DS     CL12
*        recursive function        (R0,R1)=deter(t)   (python style)
DETER    CNOP   0,4                  returns determinant and permanent
         STM    R14,R12,12(R13)    save all registers
         LR     R9,R1              save R1
         L      R2,0(R1)           n
         BCTR   R2,0               n-1
         LR     R11,R2             n-1
         MR     R10,R2             (n-1)*(n-1)
         SLA    R11,2              (n-1)*(n-1)*4
         LA     R11,1(R11)         size of q array
         A      R11,=A(STACKLEN)   R11 storage amount required
         GETMAIN RU,LV=(R11)       allocate storage for stack
         USING  STACK,R10          make storage addressable
         LR     R10,R1             establish stack addressability
         LA     R1,SAVEAREB        get the address of my savearea
         ST     R13,4(R1)          save caller's savearea pointer
         ST     R1,8(R13)          save my savearea pointer
         LR     R13,R1             set R13 to point to my savearea
         LR     R1,R9              restore R1
         LR     R9,R1              @t
         L      R4,0(R9)           t(0)
         ST     R4,N               n=t(0)
IF1      CH     R4,=H'1'           if n=1
         BNE    SIF1               then
         L      R2,4(R9)             t(1)
         ST     R2,R                 r=t(1)
         ST     R2,S                 s=t(1)
         B      EIF1               else
SIF1     L      R2,N                 n
         BCTR   R2,0                 n-1
         ST     R2,Q                 q(0)=n-1
         ST     R2,NM1               nm1=n-1
         LA     R0,1                 1
         ST     R0,SGN               sgn=1
         SR     R0,R0                0
         ST     R0,R                 r=0
         ST     R0,S                 s=0
         LA     R6,1                 k=1
LOOPK    C      R6,N                 do k=1 to n
         BH     ELOOPK               leave k
         SR     R0,R0                  0
         ST     R0,JQ                  jq=0
         ST     R0,KTI                 kti=0
         LA     R7,1                   iq=1
LOOPIQ   C      R7,NM1                 do iq=1 to n-1
         BH     ELOOPIQ                leave iq
         LR     R2,R7                    iq
         LA     R2,1(R2)                 iq+1
         ST     R2,IT                    it=iq+1
         L      R2,KTI                   kti
         A      R2,N                     kti+n
         ST     R2,KTI                   kti=kti+n
         ST     R2,KT                    kt=kti
         LA     R8,1                     jt=1
LOOPJT   C      R8,N                     do jt=1 to n
         BH     ELOOPJT                  leave jt
         L      R2,KT                      kt
         LA     R2,1(R2)                   kt+1
         ST     R2,KT                      kt=kt+1
IF2      CR     R8,R6                      if jt<>k
         BE     EIF2                       then
         L      R2,JQ                        jq
         LA     R2,1(R2)                     jq+1
         ST     R2,JQ                        jq=jq+1
         L      R1,KT                        kt
         SLA    R1,2                         *4
         L      R2,0(R1,R9)                  t(kt)
         L      R1,JQ                        jq
         SLA    R1,2                         *4
         ST     R2,Q(R1)                     q(jq)=t(kt)
EIF2     EQU    *                          end if
         LA     R8,1(R8)                   jt=jt+1
         B      LOOPJT                   next jt
ELOOPJT  LA     R7,1(R7)                 iq=iq+1
         B      LOOPIQ                 next iq
ELOOPIQ  LR     R1,R6                  k
         SLA    R1,2                   *4
         L      R5,0(R1,R9)            t(k)
         LR     R2,R5                  R2,R5=t(k)
         LA     R1,Q                   @q
         BAL    R14,DETER              call deter(q)
         LR     R3,R0                  R3=deter(q)
         ST     R1,P                   p=perm(q)
         MR     R4,R3                  R5=t(k)*deter(q)
         M      R4,SGN                 R5=sgn*t(k)*deter(q)
         A      R5,R                   +r
         ST     R5,R                   r=r+sgn*t(k)*deter(q)
         LR     R5,R2                  t(k)
         M      R4,P                   R5=t(k)*perm(q)
         A      R5,S                   +s
         ST     R5,S                   s=s+t(k)*perm(q)
         L      R2,SGN                 sgn
         LCR    R2,R2                  -sgn
         ST     R2,SGN                 sgn=-sgn
         LA     R6,1(R6)               k=k+1
         B      LOOPK                next k
ELOOPK   EQU    *                    end do
EIF1     EQU    *                  end if 
EXIT     L      R13,SAVEAREB+4     restore caller's savearea address
         L      R2,R               return value (determinant)
         L      R3,S               return value (permanent)
         XR     R15,R15            set return code to 0        
         FREEMAIN A=(R10),LV=(R11) free allocated storage
         LR     R0,R2              first return value
         LR     R1,R3              second return value
         L      R14,12(R13)        restore caller's return address
         LM     R2,R12,28(R13)     restore registers R2 to R12
         BR     R14                return to caller
IT       DS     F                  static area (out of stack)
KT       DS     F                  "
JQ       DS     F                  "
KTI      DS     F                  "
P        DS     F                  "
         DROP   R12                base no longer needed
STACK    DSECT                     dynamic area (stack)
SAVEAREB DS     18F                function savearea
N        DS     F                  n
NM1      DS     F                  n-1
R        DS     F                  determinant accu
S        DS     F                  permanent accu
SGN      DS     F                  sign
STACKLEN EQU    *-STACK
Q        DS     F                  sub matrix q((n-1)*(n-1)+1)
         YREGS  
         END    MATARI
Output:
determinant=        -360
permanent=           900

Arturo[edit]

printMatrix: function [m][
    loop m 'row -> print map row 'val [pad to :string .format:".2f" val 6]
    print "--------------------------------"
]

permutations: function [arr][
    d: 1
    c: array.of: size arr 0
    xs: new arr
    sign: 1

    ret: new @[@[xs, sign]]

    while [true][
        while [d > 1][
            d: d-1
            c\[d]: 0
        ]

        while [c\[d] >= d][
            d: d+1
            if d >= size arr -> return ret
        ]

        i: (1 = and d 1)? -> c\[d] -> 0
        tmp: xs\[i]
        xs\[i]: xs\[d]
        xs\[d]: tmp

        sign: neg sign
        'ret ++ @[new @[xs, sign]]
        c\[d]: c\[d] + 1
    ]

    return ret
]

perm: function [a][
    n: 0..dec size a
    result: new 0.0
    loop permutate n 'sigma [
        x: 1.0
        loop n 'i -> x: x * get a\[i] sigma\[i]
        'result + x
    ]
    return result
]

det: function [a][
    n: 0..dec size a
    result: new 0.0
    loop.with:'i permutations n 'p[
        x: p\1
        loop n 'i -> x: x * get a\[i] p\0\[i]
        'result + x
    ]
    return result
]

A: [[1.0 2.0]
    [3.0 4.0]]

B: [[ 1.0  2  3  4]
    [ 4.0  5  6  7]
    [ 7.0  8  9 10]
    [10.0 11 12 13]]

C: [[ 0.0  1  2  3  4]
    [ 5.0  6  7  8  9]
    [10.0 11 12 13 14]
    [15.0 16 17 18 19]
    [20.0 21 22 23 24]]

print ["A: perm ->" perm A "det ->" det A]
print ["B: perm ->" perm B "det ->" det B]
print ["C: perm ->" perm C "det ->" det C]
Output:
A: perm -> 10.0 det -> -2.0 
B: perm -> 29556.0 det -> 0.0 
C: perm -> 6778800.0 det -> 0.0

C[edit]

C99 code. By no means efficient or reliable. If you need it for serious work, go find a serious library.

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

double det_in(double **in, int n, int perm)
{
	if (n == 1) return in[0][0];

	double sum = 0, *m[--n];
	for (int i = 0; i < n; i++)
		m[i] = in[i + 1] + 1;

	for (int i = 0, sgn = 1; i <= n; i++) {
		sum += sgn * (in[i][0] * det_in(m, n, perm));
		if (i == n) break;

		m[i] = in[i] + 1;
		if (!perm) sgn = -sgn;
	}
	return sum;
}

/* wrapper function */
double det(double *in, int n, int perm)
{
	double *m[n];
	for (int i = 0; i < n; i++)
		m[i] = in + (n * i);

	return det_in(m, n, perm);
}

int main(void)
{
	double x[] = {	0, 1, 2, 3, 4,
			5, 6, 7, 8, 9,
			10, 11, 12, 13, 14,
			15, 16, 17, 18, 19,
			20, 21, 22, 23, 24 };

	printf("det:  %14.12g\n", det(x, 5, 0));
	printf("perm: %14.12g\n", det(x, 5, 1));

	return 0;
}

A method to calculate determinant that might actually be usable:

#include <stdio.h>
#include <stdlib.h>
#include <tgmath.h>

void showmat(const char *s, double **m, int n)
{
	printf("%s:\n", s);
	for (int i = 0; i < n; i++) {
		for (int j = 0; j < n; j++)
			printf("%12.4f", m[i][j]);
			putchar('\n');
	}
}

int trianglize(double **m, int n)
{
	int sign = 1;
	for (int i = 0; i < n; i++) {
		int max = 0;

		for (int row = i; row < n; row++)
			if (fabs(m[row][i]) > fabs(m[max][i]))
				max = row;

		if (max) {
			sign = -sign;
			double *tmp = m[i];
			m[i] = m[max], m[max] = tmp;
		}

		if (!m[i][i]) return 0;

		for (int row = i + 1; row < n; row++) {
			double r = m[row][i] / m[i][i];
			if (!r)	continue;

			for (int col = i; col < n; col ++)
				m[row][col] -= m[i][col] * r;
		}
	}
	return sign;
}

double det(double *in, int n)
{
	double *m[n];
	m[0] = in;

	for (int i = 1; i < n; i++)
		m[i] = m[i - 1] + n;

	showmat("Matrix", m, n);

	int sign = trianglize(m, n);
	if (!sign)
		return 0;

	showmat("Upper triangle", m, n);

	double p = 1;
	for (int i = 0; i < n; i++)
		p *= m[i][i];
	return p * sign;
}

#define N 18
int main(void)
{
	double x[N * N];
	srand(0);
	for (int i = 0; i < N * N; i++)
		x[i] = rand() % N;

	printf("det: %19f\n", det(x, N));
	return 0;
}

C++[edit]

Translation of: Java
#include <iostream>
#include <vector>

template <typename T>
std::ostream &operator<<(std::ostream &os, const std::vector<T> &v) {
    auto it = v.cbegin();
    auto end = v.cend();

    os << '[';
    if (it != end) {
        os << *it;
        it = std::next(it);
    }
    while (it != end) {
        os << ", " << *it;
        it = std::next(it);
    }
    return os << ']';
}

using Matrix = std::vector<std::vector<double>>;

Matrix squareMatrix(size_t n) {
    Matrix m;
    for (size_t i = 0; i < n; i++) {
        std::vector<double> inner;
        for (size_t j = 0; j < n; j++) {
            inner.push_back(nan(""));
        }
        m.push_back(inner);
    }
    return m;
}

Matrix minor(const Matrix &a, int x, int y) {
    auto length = a.size() - 1;
    auto result = squareMatrix(length);
    for (int i = 0; i < length; i++) {
        for (int j = 0; j < length; j++) {
            if (i < x && j < y) {
                result[i][j] = a[i][j];
            } else if (i >= x && j < y) {
                result[i][j] = a[i + 1][j];
            } else if (i < x && j >= y) {
                result[i][j] = a[i][j + 1];
            } else {
                result[i][j] = a[i + 1][j + 1];
            }
        }
    }
    return result;
}

double det(const Matrix &a) {
    if (a.size() == 1) {
        return a[0][0];
    }

    int sign = 1;
    double sum = 0;
    for (size_t i = 0; i < a.size(); i++) {
        sum += sign * a[0][i] * det(minor(a, 0, i));
        sign *= -1;
    }
    return sum;
}

double perm(const Matrix &a) {
    if (a.size() == 1) {
        return a[0][0];
    }

    double sum = 0;
    for (size_t i = 0; i < a.size(); i++) {
        sum += a[0][i] * perm(minor(a, 0, i));
    }
    return sum;
}

void test(const Matrix &m) {
    auto p = perm(m);
    auto d = det(m);

    std::cout << m << '\n';
    std::cout << "Permanent: " << p << ", determinant: " << d << "\n\n";
}

int main() {
    test({ {1, 2}, {3, 4} });
    test({ {1, 2, 3, 4}, {4, 5, 6, 7}, {7, 8, 9, 10}, {10, 11, 12, 13} });
    test({ {0, 1, 2, 3, 4}, {5, 6, 7, 8, 9}, {10, 11, 12, 13, 14}, {15, 16, 17, 18, 19}, {20, 21, 22, 23, 24} });

    return 0;
}
Output:
[[1, 2], [3, 4]]
Permanent: 10, determinant: -2

[[1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10], [10, 11, 12, 13]]
Permanent: 29556, determinant: 0

[[0, 1, 2, 3, 4], [5, 6, 7, 8, 9], [10, 11, 12, 13, 14], [15, 16, 17, 18, 19], [20, 21, 22, 23, 24]]
Permanent: 6.7788e+06, determinant: 0

Common Lisp[edit]

A recursive version, no libraries required, it doesn't use much consing, only for the list of columns to skip

(defun determinant (rows &optional (skip-cols nil))
  (let* ((result 0) (sgn -1))
    (dotimes (col (length (car rows)) result)
      (unless (member col skip-cols)
        (if (null (cdr rows))
          (return-from determinant (elt (car rows) col))
          (incf result (* (setq sgn (- sgn)) (elt (car rows) col) (determinant (cdr rows) (cons col skip-cols)))) )))))

(defun permanent (rows &optional (skip-cols nil))
  (let* ((result 0))
    (dotimes (col (length (car rows)) result)
      (unless (member col skip-cols)
        (if (null (cdr rows))
          (return-from permanent (elt (car rows) col))
          (incf result (* (elt (car rows) col) (permanent (cdr rows) (cons col skip-cols)))) )))))


Test using the first set of definitions (from task description):

(setq m2
  '((1 2)
    (3 4)))

(setq m3
  '((-2 2 -3)
    (-1 1  3)
    ( 2 0 -1)))

(setq m4
  '(( 1  2  3  4)
    ( 4  5  6  7)
    ( 7  8  9 10)
    (10 11 12 13)))

(setq m5
  '(( 0  1  2  3  4)
    ( 5  6  7  8  9)
    (10 11 12 13 14)
    (15 16 17 18 19)
    (20 21 22 23 24)))

(dolist (m (list m2 m3 m4 m5))
  (format t "~a determinant: ~a, permanent: ~a~%" m (determinant m) (permanent m)) )
Output:
((1 2) (3 4)) determinant: -2, permanent: 10
((-2 2 -3) (-1 1 3) (2 0 -1)) determinant: 18, permanent: 10
((1 2 3 4) (4 5 6 7) (7 8 9 10) (10 11 12 13)) determinant: 0, permanent: 29556
((0 1 2 3 4) (5 6 7 8 9) (10 11 12 13 14) (15 16 17 18 19) (20 21 22 23 24)) determinant: 0, permanent: 6778800

D[edit]

This requires the modules from the Permutations and Permutations by swapping tasks.

Translation of: Python
import std.algorithm, std.range, std.traits, permutations2,
       permutations_by_swapping1;

auto prod(Range)(Range r) nothrow @safe @nogc {
    return reduce!q{a * b}(ForeachType!Range(1), r);
}

T permanent(T)(in T[][] a) nothrow @safe
in {
    assert(a.all!(row => row.length == a[0].length));
} body {
    auto r = a.length.iota;
    T tot = 0;
    foreach (const sigma; r.array.permutations)
        tot += r.map!(i => a[i][sigma[i]]).prod;
    return tot;
}

T determinant(T)(in T[][] a) nothrow
in {
    assert(a.all!(row => row.length == a[0].length));
} body {
    immutable n = a.length;
    auto r = n.iota;
    T tot = 0;
    //foreach (sigma, sign; n.spermutations) {
    foreach (const sigma_sign; n.spermutations) {
        const sigma = sigma_sign[0];
        immutable sign = sigma_sign[1];
        tot += sign * r.map!(i => a[i][sigma[i]]).prod;
    }
    return tot;
}

void main() {
    import std.stdio;

    foreach (/*immutable*/ const a; [[[1, 2],
                                      [3, 4]],

                                     [[1, 2, 3, 4],
                                      [4, 5, 6, 7],
                                      [7, 8, 9, 10],
                                      [10, 11, 12, 13]],

                                     [[ 0,  1,  2,  3,  4],
                                      [ 5,  6,  7,  8,  9],
                                      [10, 11, 12, 13, 14],
                                      [15, 16, 17, 18, 19],
                                      [20, 21, 22, 23, 24]]]) {
        writefln("[%([%(%2s, %)],\n %)]]", a);
        writefln("Permanent: %s, determinant: %s\n",
                 a.permanent, a.determinant);
    }
}
Output:
[[ 1,  2],
 [ 3,  4]]
Permanent: 10, determinant: -2

[[ 1,  2,  3,  4],
 [ 4,  5,  6,  7],
 [ 7,  8,  9, 10],
 [10, 11, 12, 13]]
Permanent: 29556, determinant: 0

[[ 0,  1,  2,  3,  4],
 [ 5,  6,  7,  8,  9],
 [10, 11, 12, 13, 14],
 [15, 16, 17, 18, 19],
 [20, 21, 22, 23, 24]]
Permanent: 6778800, determinant: 0

Delphi[edit]

Translation of: Java
program Determinant_and_permanent;

{$APPTYPE CONSOLE}

uses
  System.SysUtils;

type
  TMatrix = TArray<TArray<Double>>;

function Minor(a: TMatrix; x, y: Integer): TMatrix;
begin
  var len := Length(a) - 1;
  SetLength(result, len, len);
  for var i := 0 to len - 1 do
  begin
    for var j := 0 to len - 1 do
    begin
      if ((i < x) and (j < y)) then
      begin
        result[i][j] := a[i][j];
      end
      else if ((i >= x) and (j < y)) then
      begin
        result[i][j] := a[i + 1][j];
      end
      else if ((i < x) and (j >= y)) then
      begin
        result[i][j] := a[i][j + 1];
      end
      else //i>x  and  j>y
        result[i][j] := a[i + 1][j + 1];
    end;
  end;
end;

function det(a: TMatrix): Double;
begin
  if length(a) = 1 then
    exit(a[0][0]);

  var sign := 1;
  result := 0.0;
  for var i := 0 to high(a) do
  begin
    result := result + sign * a[0][i] * det(minor(a, 0, i));
    sign := sign *  - 1;
  end;
end;

function perm(a: TMatrix): Double;
begin
  if Length(a) = 1 then
    exit(a[0][0]);
  Result := 0;
  for var i := 0 to high(a) do
    result := result + a[0][i] * perm(Minor(a, 0, i));
end;

function Readint(Min, Max: Integer; Prompt: string): Integer;
var
  val: string;
  vali: Integer;
begin
  Result := -1;
  repeat
    writeln(Prompt);
    Readln(val);
    if TryStrToInt(val, vali) then
      if (vali < Min) or (vali > Max) then
        writeln(vali, ' is out range [', Min, '...', Max, ']')
      else
        exit(vali)
    else
      writeln(val, ' is not a number valid');
  until false;
end;

function ReadDouble(Min, Max: double; Prompt: string): double;
var
  val: string;
  vali: double;
begin
  Result := -1;
  repeat
    writeln(Prompt);
    Readln(val);
    if TryStrToFloat(val, vali) then
      if (vali < Min) or (vali > Max) then
        writeln(vali, ' is out range [', Min, '...', Max, ']')
      else
        exit(vali)
    else
      writeln(val, ' is not a number valid');
  until false;
end;

procedure ShowMatrix(a: TMatrix);
begin
  var sz := length(a);
  for var i := 0 to sz - 1 do
  begin
    Write('[');
    for var j := 0 to sz - 1 do
      write(a[i][j]: 3: 2, ' ');
    Writeln(']');
  end;
end;

var
  a: TMatrix;
  sz: integer;

begin
  sz := Readint(1, 10, 'Enter with matrix size: ');
  SetLength(a, sz, sz);
  for var i := 0 to sz - 1 do
    for var j := 0 to sz - 1 do
    begin
      a[i][j] := ReadDouble(-1000, 1000, format('Enter a value of position (%d,%d):',
        [i, j]));
    end;

  writeln('Matrix defined: ');
  ShowMatrix(a);
  writeln(#10'Determinant: ', det(a): 3: 2);
  writeln(#10'Permanent: ', perm(a): 3: 2);
  readln;
end.
Output:
Enter with matrix size:
2
Enter a value of position (0,0):
1
Enter a value of position (0,1):
2
Enter a value of position (1,0):
3
Enter a value of position (1,1):
4
Matrix defined:
[1.00 2.00 ]
[3.00 4.00 ]

Determinant: -2.00

Permanent: 10.00

EchoLisp[edit]

This requires the 'list' library for (in-permutations n) and the 'matrix' library for the built-in (determinant M).

(lib 'list)
(lib 'matrix)

;; adapted from Racket
(define (permanent M)
    (let (( n (matrix-row-num M)))
    (for/sum ( (in-permutations n)])
        (for/product ([i n] [σi σ])
            (array-ref M i σi)))))

;; output
(define A (list->array '(1 2 3 4) 2 2))
(array-print A)
  1  2
  3  4
(determinant A)  -2
(permanent A)  10

(define M (list->array (iota 25) 5 5))
(array-print M)
   0   1   2   3   4
   5   6   7   8   9
  10  11  12  13  14
  15  16  17  18  19
  20  21  22  23  24
(determinant M)  0
(permanent M)  6778800

Factor[edit]

USING: fry kernel math.combinatorics math.matrices sequences ;

: permanent ( matrix -- x )
    dup square-matrix? [ "Matrix must be square." throw ] unless
    [ dim first <iota> ] keep
    '[ [ _ nth nth ] map-index product ] map-permutations sum ;

Example output:

IN: scratchpad USE: math.matrices.laplace   ! for determinant
               { { 2 9 4 } { 7 5 3 } { 6 1 8 } }
               [ determinant ] [ permanent ] bi

--- Data stack:
-360
900


Forth[edit]

Works with: gforth version 0.7.9_20170427

Requiring a permute.fs file from the Permutations by swapping task.

S" fsl-util.fs" REQUIRED
S" fsl/dynmem.seq" REQUIRED
[UNDEFINED] defines [IF] SYNONYM defines IS [THEN]
S" fsl/structs.seq" REQUIRED
S" fsl/lufact.seq" REQUIRED
S" fsl/dets.seq" REQUIRED
S" permute.fs" REQUIRED

VARIABLE the-mat
: add-perm ( p0 p1 p2 ... pn n s -- )
  DROP  \ sign
  1E
  1 DO
    the-mat @ SWAP 1- I 1- }} F@ F*
  LOOP
  DROP  \ Dummy element because we're using 1-based indexing
  F+ ;
: permanent ( len mat -- ) ( F: -- perm )
  the-mat !
  0E
  ['] add-perm perms ;

3 SET-PRECISION
2 2 float matrix m2{{
1e 2e  3e 4e  2 2 m2{{ }}fput
lumatrix lmat
3 3 float matrix m3{{
2e 9e 4e  7e 5e 3e  6e 1e 8e  3 3 m3{{ }}fput

lmat 2 lu-malloc
m2{{ lmat lufact
lmat det F. 2 m2{{ permanent F. CR
lmat lu-free

lmat 3 lu-malloc
m3{{ lmat lufact
lmat det F. 3 m3{{ permanent F. CR
lmat lu-free

Fortran[edit]

Please find the compilation and example run at the start of the comments in the f90 source. Thank you.

!-*- mode: compilation; default-directory: "/tmp/" -*-
!Compilation started at Sat May 18 23:25:42
!
!a=./F && make $a && $a < unixdict.txt
!f95 -Wall -ffree-form F.F -o F
! j example, determinant:    7.00000000    
! j example, permanent:      5.00000000    
! maxima, determinant:      -360.000000    
! maxima, permanent:         900.000000    
!
!Compilation finished at Sat May 18 23:25:43



!   NB. example computed by J
!   NB. fixed seed random matrix
!   _2+3 3?.@$5
! 2 _1  1
!_1 _2  1
!_1 _1 _1
!
!   (-/ .*)_2+3 3?.@$5  NB. determinant
!7
!   (+/ .*)_2+3 3?.@$5  NB. permanent
!5

!maxima example
!a: matrix([2, 9, 4], [7, 5, 3], [6, 1, 8])$
!determinant(a);
!-360
! 
!permanent(a);
!900


! compute permanent or determinant
program f
  implicit none
  real, dimension(3,3) :: j, m
  data j/ 2,-1, 1,-1,-2, 1,-1,-1,-1/
  data m/2, 9, 4, 7, 5, 3, 6, 1, 8/
  write(6,*) 'j example, determinant: ',det(j,3,-1)
  write(6,*) 'j example, permanent:   ',det(j,3,1)
  write(6,*) 'maxima, determinant:    ',det(m,3,-1)
  write(6,*) 'maxima, permanent:      ',det(m,3,1)

contains

  recursive function det(a,n,permanent) result(accumulation)
    ! setting permanent to 1 computes the permanent.
    ! setting permanent to -1 computes the determinant.
    real, dimension(n,n), intent(in) :: a
    integer, intent(in) :: n, permanent
    real, dimension(n-1, n-1) :: b
    real :: accumulation
    integer :: i, sgn
    if (n .eq. 1) then
      accumulation = a(1,1)
    else
      accumulation = 0
      sgn = 1
      do i=1, n
        b(:, :(i-1)) = a(2:, :i-1)
        b(:, i:) = a(2:, i+1:)
        accumulation = accumulation + sgn * a(1, i) * det(b, n-1, permanent)
        sgn = sgn * permanent
      enddo
    endif
  end function det

end program f

FreeBASIC[edit]

sub make_S( M() as double, S() as double, i as uinteger, j as uinteger )
    'removes row j, column i from the matrix, stores result in S()
    dim as uinteger ii, jj, size=ubound(M), ix, jx
    for ii = 1 to size-1
        if ii<i then ix = ii else ix = ii + 1
        for jj = 1 to size-1
            if jj<j then jx = jj else jx = jj + 1
            S(ii, jj) = M(ix, jx)
        next jj
    next ii
end sub

function deperminant( M() as double, det as boolean ) as double
    'calculates the determinant or the permanent of a square matrix M
    'det = true  for determinant, false for permanent
    'assumes a square matrix
    dim as uinteger size = ubound(M,1), i
    dim as integer sign
    dim as double S(1 to size-1, 1 to size-1)
    dim as double ret = 0.0, inc
    if size = 1 then return M(1,1)                                       'matrices of size < 3 are easy to calculate
    if size = 2 and det then return M(1,1)*M(2,2) - M(1,2)*M(2,1)
    if size = 2         then return M(1,1)*M(2,2) + M(1,2)*M(2,1)
    for i = 1 to size
            if det then sign = (-1)^(i+1) else sign = 1                 'this bit is what distinguishes a determinant from a permanent
            make_S( M(), S(), i, 1 )
            inc = sign*M(i,1)*deperminant( S(), det )                   'recursively call on submatrices
            ret += inc
    next i
    return ret
end function
    
dim as double A(1 to 2, 1 to 2) = {{1,2},{3,4}}

dim as double B(1 to 4, 1 to 4) = {_
    {1,2,3,4}, {4,5,6,7}, {7,8,9,10}, {10,11,12,13} }

dim as double C(1 to 5, 1 to 5) = {_
     {  0,  1,  2,  3,  4 },_
     {  5,  6,  7,  8,  9 },_
     { 10, 11, 12, 13, 14 },_
     { 15, 16, 17, 18, 19 },_
     { 20, 21, 22, 23, 24 } }

print deperminant( A(), true ), deperminant( A(), false )
print deperminant( B(), true ), deperminant( B(), false )
print deperminant( C(), true ), deperminant( C(), false )
Output:

-2 10 0 29556 0 6778800

FunL[edit]

From the task description:

def sgn( p ) = product( (if s(0) < s(1) xor i(0) < i(1) then -1 else 1) | (s, i) <- p.combinations(2).zip( (0:p.length()).combinations(2) ) )

def perm( m ) = sum( product(m(i, sigma(i)) | i <- 0:m.length()) | sigma <- (0:m.length()).permutations() )

def det( m ) = sum( sgn(sigma)*product(m(i, sigma(i)) | i <- 0:m.length()) | sigma <- (0:m.length()).permutations() )

Laplace expansion (recursive):

def perm( m )
  | m.length() == 1 and m(0).length() == 1 = m(0, 0)
  | otherwise = sum( m(i, 0)*perm(m(0:i, 1:m.length()) + m(i+1:m.length(), 1:m.length())) | i <- 0:m.length() )

def det( m )
  | m.length() == 1 and m(0).length() == 1 = m(0, 0)
  | otherwise = sum( (-1)^i*m(i, 0)*det(m(0:i, 1:m.length()) + m(i+1:m.length(), 1:m.length())) | i <- 0:m.length() )

Test using the first set of definitions (from task description):

matrices = [
  ( (1, 2),
    (3, 4)),
  ( (-2, 2, -3),
    (-1, 1,  3),
    ( 2, 0, -1)),
  ( ( 1,  2,  3,  4),
    ( 4,  5,  6,  7),
    ( 7,  8,  9, 10),
    (10, 11, 12, 13)),
  ( ( 0,  1,  2,  3,  4),
    ( 5,  6,  7,  8,  9),
    (10, 11, 12, 13, 14),
    (15, 16, 17, 18, 19),
    (20, 21, 22, 23, 24)) ]

for m <- matrices
  println( m, 'perm: ' + perm(m), 'det: ' + det(m) )
Output:
((1, 2), (3, 4)), perm: 10, det: -2
((-2, 2, -3), (-1, 1, 3), (2, 0, -1)), perm: 10, det: 18
((1, 2, 3, 4), (4, 5, 6, 7), (7, 8, 9, 10), (10, 11, 12, 13)), perm: 29556, det: 0
((0, 1, 2, 3, 4), (5, 6, 7, 8, 9), (10, 11, 12, 13, 14), (15, 16, 17, 18, 19), (20, 21, 22, 23, 24)), perm: 6778800, det: 0

GLSL[edit]

  mat4 m1 = mat3(1, 2, 3, 4,
                5, 6, 7, 8
                9,10,11,12,
                13,14,15,16);

  float d = det(m1);

Go[edit]

Implementation[edit]

This implements a naive algorithm for each that follows from the definitions. It imports the permute packge from the Permutations by swapping task.

package main

import (
    "fmt"
    "permute"
)

func determinant(m [][]float64) (d float64) {
    p := make([]int, len(m))
    for i := range p {
        p[i] = i
    }
    it := permute.Iter(p)
    for s := it(); s != 0; s = it() {
        pr := 1.
        for i, σ := range p {
            pr *= m[i][σ]
        }
        d += float64(s) * pr
    }
    return
}

func permanent(m [][]float64) (d float64) {
    p := make([]int, len(m))
    for i := range p {
        p[i] = i
    }
    it := permute.Iter(p)
    for s := it(); s != 0; s = it() {
        pr := 1.
        for i, σ := range p {
            pr *= m[i][σ]
        }
        d += pr
    }
    return
}

var m2 = [][]float64{
    {1, 2},
    {3, 4}}

var m3 = [][]float64{
    {2, 9, 4},
    {7, 5, 3},
    {6, 1, 8}}

func main() {
    fmt.Println(determinant(m2), permanent(m2))
    fmt.Println(determinant(m3), permanent(m3))
}
Output:
-2 10
-360 900

Ryser permanent[edit]

package main

import "fmt"

func main() {
    fmt.Println(ryser([][]float64{
        {1, 2},
        {3, 4}}))
    fmt.Println(ryser([][]float64{
        {2, 9, 4},
        {7, 5, 3},
        {6, 1, 8}}))
}

func ryser(m [][]float64) (d float64) {
    gray := 0
    csum := make([]float64, len(m))
    sgn := float64(len(m)&1<<1 - 1)
    n2 := uint32(1) << uint(len(m))
    for i := uint32(1); i < n2; i++ {
        r := [...]byte{
            0, 1, 28, 2, 29, 14, 24, 3, 30, 22, 20, 15, 25, 17, 4, 8,
            31, 27, 13, 23, 21, 19, 16, 7, 26, 12, 18, 6, 11, 5, 10, 9,
        }[i&-i*0x077CB531>>27]
        b := 1 << r
        if gray&b == 0 {
            for c, e := range m[r] {
                csum[c] += e
            }
        } else {
            for c, e := range m[r] {
                csum[c] -= e
            }
        }
        gray ^= b
        p := sgn
        for _, e := range csum {
            p *= e
        }
        d += p
        sgn = -sgn
    }
    return
}
Output:
10
900

Library determinant[edit]

go.matrix:

package main

import (
    "fmt"

    "github.com/skelterjohn/go.matrix"
)

func main() {
    fmt.Println(matrix.MakeDenseMatrixStacked([][]float64{
        {1, 2},
        {3, 4}}).Det())
    fmt.Println(matrix.MakeDenseMatrixStacked([][]float64{
        {2, 9, 4},
        {7, 5, 3},
        {6, 1, 8}}).Det())
}
Output:
-2
-360

gonum/mat:

package main

import (
    "fmt"

    "gonum.org/v1/gonum/mat"
)

func main() {
    fmt.Println(mat.Det(mat.NewDense(2, 2, []float64{
        1, 2,
        3, 4})))
    fmt.Println(mat.Det(mat.NewDense(3, 3, []float64{
        2, 9, 4,
        7, 5, 3,
        6, 1, 8})))
}
Output:
-2
-360.00000000000006

Haskell[edit]

sPermutations :: [a] -> [([a], Int)]
sPermutations = flip zip (cycle [1, -1]) . foldl aux [[]]
  where
    aux items x = do
      (f, item) <- zip (cycle [reverse, id]) items
      f (insertEv x item)
    insertEv x [] = [[x]]
    insertEv x l@(y:ys) = (x : l) : ((y :) <$>) (insertEv x ys)

elemPos :: [[a]] -> Int -> Int -> a
elemPos ms i j = (ms !! i) !! j

prod
  :: Num a
  => ([[a]] -> Int -> Int -> a) -> [[a]] -> [Int] -> a
prod f ms = product . zipWith (f ms) [0 ..]

sDeterminant
  :: Num a
  => ([[a]] -> Int -> Int -> a) -> [[a]] -> [([Int], Int)] -> a
sDeterminant f ms = sum . fmap (\(is, s) -> fromIntegral s * prod f ms is)

determinant
  :: Num a
  => [[a]] -> a
determinant ms =
  sDeterminant elemPos ms . sPermutations $ [0 .. pred . length $ ms]

permanent
  :: Num a
  => [[a]] -> a
permanent ms =
  sum . fmap (prod elemPos ms . fst) . sPermutations $ [0 .. pred . length $ ms]

-- TEST -----------------------------------------------------------------------
result
  :: (Num a, Show a)
  => [[a]] -> String
result ms =
  unlines
    [ "Matrix:"
    , unlines (show <$> ms)
    , "Determinant:"
    , show (determinant ms)
    , "Permanent:"
    , show (permanent ms)
    ]

main :: IO ()
main =
  mapM_
    (putStrLn . result)
    [ [[5]]
    , [[1, 0, 0], [0, 1, 0], [0, 0, 1]]
    , [[0, 0, 1], [0, 1, 0], [1, 0, 0]]
    , [[4, 3], [2, 5]]
    , [[2, 5], [4, 3]]
    , [[4, 4], [2, 2]]
    ]
Output:
Matrix:
[5]

Determinant:
5
Permanent:
5

Matrix:
[1,0,0]
[0,1,0]
[0,0,1]

Determinant:
1
Permanent:
1

Matrix:
[0,0,1]
[0,1,0]
[1,0,0]

Determinant:
-1
Permanent:
1

Matrix:
[4,3]
[2,5]

Determinant:
14
Permanent:
26

Matrix:
[2,5]
[4,3]

Determinant:
-14
Permanent:
26

Matrix:
[4,4]
[2,2]

Determinant:
0
Permanent:
16

Via Cramer's rule[edit]

Here is code for computing the determinant and permanent very inefficiently, via Cramer's rule (for the determinant, as well as its analog for the permanent):

outer :: (a->b->c) -> [a] -> [b] -> [[c]]
outer f [] _       = []
outer f _ []       = []
outer f (h1:t1) x2 = (f h1 <$> x2) : outer f t1 x2

dot [] []           = 0
dot (h1:t1) (h2:t2) = (h1*h2) + (dot t1 t2)

transpose [] = []
transpose ([] : xss) = transpose xss
transpose ((x:xs) : xss)
  = (x : [h | (h:_) <- xss]) : transpose (xs : [ t | (_:t) <- xss])

mul :: Num a => [[a]] -> [[a]] -> [[a]]
mul a b = outer dot a (transpose b)

delRow :: Int -> [a] -> [a]
delRow i v =
  (first ++ rest) where (first, _:rest) = splitAt i v

delCol :: Int -> [[a]] -> [[a]]
delCol j m = (delRow j) <$> m

-- Determinant:
adj :: Num a => [[a]] -> [[a]]
adj [] = []
adj m =
  [
    [(-1)^(i+j) * det (delRow i $ delCol j m)
    | i <- [0.. -1+length m]
    ]
  | j <- [0.. -1+length m]
  ]
det :: Num a => [[a]] -> a
det [] = 1
det m  = (mul m (adj m)) !! 0 !! 0

-- Permanent:
padj :: Num a => [[a]] -> [[a]]
padj [] = []
padj m =
  [
    [perm (delRow i $ delCol j m)
    | i <- [0.. -1+length m]
    ]
  | j <- [0.. -1+length m]
  ]
perm :: Num a => [[a]] -> a
perm [] = 1
perm m  = (mul m (padj m)) !! 0 !! 0

J[edit]

J has a conjunction for defining verbs which can act as determinant (especially -/ .* ). This conjunction is symbolized as a space followed by a dot. And you can get the permanent by replacing - in that definition with +.

For example, given the matrix:

   i. 5 5
 0  1  2  3  4
 5  6  7  8  9
10 11 12 13 14
15 16 17 18 19
20 21 22 23 24

Its determinant is 0. When we use IEEE floating point, we only get an approximation of this result:

   -/ .* i. 5 5
_1.30277e_44

If we use exact (rational) arithmetic, we get a precise result:

   -/ .* i. 5 5x
0

Meanwhile, the permanent does not have this problem in this example (the matrix contains no negative values and permanent does not use subtraction):

   +/ .* i. 5 5
6778800

As an aside, note also that for specific verbs (like -/ .*) J uses an algorithm which is more efficient than the brute force approach implied by the definition of .. (In general, where there are common, useful, concise definitions where special code can improve resource use by more than a factor of 2, the implementors of J try to make sure that that special code gets used for those definitions.)

Java[edit]

import java.util.Scanner;

public class MatrixArithmetic {
	public static double[][] minor(double[][] a, int x, int y){
		int length = a.length-1;
		double[][] result = new double[length][length];
		for(int i=0;i<length;i++) for(int j=0;j<length;j++){
			if(i<x && j<y){
				result[i][j] = a[i][j];
			}else if(i>=x && j<y){
				result[i][j] = a[i+1][j];
			}else if(i<x && j>=y){
				result[i][j] = a[i][j+1];
			}else{ //i>x && j>y
				result[i][j] = a[i+1][j+1];
			}
		}
		return result;
	}
	public static double det(double[][] a){
		if(a.length == 1){
			return a[0][0];
		}else{
			int sign = 1;
			double sum = 0;
			for(int i=0;i<a.length;i++){
				sum += sign * a[0][i] * det(minor(a,0,i));
				sign *= -1;
			}
			return sum;
		}
	}
	public static double perm(double[][] a){
		if(a.length == 1){
			return a[0][0];
		}else{
			double sum = 0;
			for(int i=0;i<a.length;i++){
				sum += a[0][i] * perm(minor(a,0,i));
			}
			return sum;
		}
	}
	public static void main(String args[]){
		Scanner sc = new Scanner(System.in);
		int size = sc.nextInt();
		double[][] a = new double[size][size];
		for(int i=0;i<size;i++) for(int j=0;j<size;j++){
			a[i][j] = sc.nextDouble();
		}
		sc.close();
		System.out.println("Determinant: "+det(a));
		System.out.println("Permanent: "+perm(a));
	}
}

Note that the first input is the size of the matrix.

For example:

2
1 2
3 4
Determinant: -2.0
Permanent: 10.0


5
0 1 2 3 4
5 6 7 8 9
10 11 12 13 14
15 16 17 18 19
20 21 22 23 24
Determinant: 0.0
Permanent: 6778800.0

JavaScript[edit]

const determinant = arr =>
    arr.length === 1 ? (
        arr[0][0]
    ) : arr[0].reduce(
        (sum, v, i) => sum + v * (-1) ** i * determinant(
            arr.slice(1)
            .map(x => x.filter((_, j) => i !== j))
        ), 0
    );

const permanent = arr =>
    arr.length === 1 ? (
        arr[0][0]
    ) : arr[0].reduce(
        (sum, v, i) => sum + v * permanent(
            arr.slice(1)
            .map(x => x.filter((_, j) => i !== j))
        ), 0
    );

const M = [
    [0, 1, 2, 3, 4],
    [5, 6, 7, 8, 9],
    [10, 11, 12, 13, 14],
    [15, 16, 17, 18, 19],
    [20, 21, 22, 23, 24]
];
console.log(determinant(M));
console.log(permanent(M));
Output:
0
6778800

jq[edit]

Works with: jq version 1.4

Recursive definitions[edit]

# Eliminate row i and row j
def except(i;j):
  reduce del(.[i])[] as $row ([]; . + [$row | del(.[j]) ] );
 
def det:
  def parity(i): if i % 2 == 0 then 1 else -1 end;
  if length == 1 and (.[0] | length) == 1 then .[0][0]
  else . as $m
    | reduce range(0; length) as $i
        (0; . + parity($i) * $m[0][$i] * ( $m | except(0;$i) | det) )
  end ;

def perm:
  if length == 1 and (.[0] | length) == 1 then .[0][0]
  else . as $m
    | reduce range(0; length) as $i
        (0; . + $m[0][$i] * ( $m | except(0;$i) | perm) )
  end ;

Examples

def matrices:
  [ [1, 2],
    [3, 4]],

  [ [-2, 2, -3],
    [-1, 1,  3],
    [ 2, 0, -1]],

  [ [ 1,  2,  3,  4],
    [ 4,  5,  6,  7],
    [ 7,  8,  9, 10],
    [10, 11, 12, 13]],

  [ [ 0,  1,  2,  3,  4],
    [ 5,  6,  7,  8,  9],
    [10, 11, 12, 13, 14],
    [15, 16, 17, 18, 19],
    [20, 21, 22, 23, 24]] 
;

"Determinants: ", (matrices | det), 
"Permanents:   ",  (matrices | perm)
Output:
$ jq -n -r -f Matrix_arithmetic.jq
Determinants: 
-2
18
0
0
Permanents:   
10
10
29556
6778800

Determinant via LU Decomposition[edit]

The following uses the jq infrastructure at LU decomposition to achieve an efficient implementation of det/0:

# Requires lup/0
def det:
  def product_diagonal:
    . as $m | reduce range(0;length) as $i (1; . * $m[$i][$i]);
  def tidy: if . == -0 then 0 else . end;
  lup 
  | (.[0]|product_diagonal) as $l
  | if $l == 0 then 0 else $l * (.[1]|product_diagonal) | tidy end ;

Examples

Using matrices/0 as defined above:

matrices | det
Output:
$ /usr/local/bin/jq -M -n -f LU.rc
2
-18
0
0

Julia[edit]

 using LinearAlgebra

The determinant of a matrix A can be computed by the built-in function

det(A)
Translation of: Python

The following function computes the permanent of a matrix A from the definition:

function perm(A)
  m, n = size(A)
  if m != n; throw(ArgumentError("permanent is for square matrices only")); end
  sum(σ -> prod(i -> A[i,σ[i]], 1:n), permutations(1:n))
end

Example output:

julia> A = [2 9 4; 7 5 3; 6 1 8]
julia> det(A), perm(A)
(-360.0,900)

Kotlin[edit]

// version 1.1.2

typealias Matrix = Array<DoubleArray>

fun johnsonTrotter(n: Int): Pair<List<IntArray>, List<Int>> {
    val p = IntArray(n) { it }  // permutation
    val q = IntArray(n) { it }  // inverse permutation
    val d = IntArray(n) { -1 }  // direction = 1 or -1
    var sign = 1
    val perms = mutableListOf<IntArray>()
    val signs = mutableListOf<Int>()

    fun permute(k: Int) {
        if (k >= n) {
            perms.add(p.copyOf())
            signs.add(sign)
            sign *= -1
            return
        } 
        permute(k + 1)
        for (i in 0 until k) {
            val z = p[q[k] + d[k]]
            p[q[k]] = z
            p[q[k] + d[k]] = k
            q[z] = q[k]
            q[k] += d[k]
            permute(k + 1)
        }
        d[k] *= -1
    } 

    permute(0)
    return perms to signs
}

fun determinant(m: Matrix): Double {
    val (sigmas, signs) = johnsonTrotter(m.size)
    var sum = 0.0    
    for ((i, sigma) in sigmas.withIndex()) {
        var prod = 1.0
        for ((j, s) in sigma.withIndex()) prod *= m[j][s]
        sum += signs[i] * prod
    }
    return sum
}

fun permanent(m: Matrix) : Double {
    val (sigmas, _) = johnsonTrotter(m.size)
    var sum = 0.0
    for (sigma in sigmas) {
        var prod = 1.0
        for ((i, s) in sigma.withIndex()) prod *= m[i][s]
        sum += prod
    }
    return sum
}

fun main(args: Array<String>) {
    val m1 = arrayOf(
        doubleArrayOf(1.0)
    )

    val m2 = arrayOf(
        doubleArrayOf(1.0, 2.0),
        doubleArrayOf(3.0, 4.0)
    )

    val m3 = arrayOf(
        doubleArrayOf(2.0, 9.0, 4.0),
        doubleArrayOf(7.0, 5.0, 3.0),
        doubleArrayOf(6.0, 1.0, 8.0)
    )

    val m4 = arrayOf(
        doubleArrayOf( 1.0,  2.0,  3.0,  4.0),
        doubleArrayOf( 4.0,  5.0,  6.0,  7.0),
        doubleArrayOf( 7.0,  8.0,  9.0, 10.0),     
        doubleArrayOf(10.0, 11.0, 12.0, 13.0)
    )      

    val matrices = arrayOf(m1, m2, m3, m4)
    for (m in matrices) {
        println("m${m.size} -> ")
        println("  determinant = ${determinant(m)}")
        println("  permanent   = ${permanent(m)}\n")
    } 
}
Output:
m1 -> 
  determinant = 1.0
  permanent   = 1.0

m2 -> 
  determinant = -2.0
  permanent   = 10.0

m3 -> 
  determinant = -360.0
  permanent   = 900.0

m4 -> 
  determinant = 0.0
  permanent   = 29556.0

Lambdatalk[edit]

{require lib_matrix}

{M.determinant
 {M.new [[1,2,3],
         [4,5,6],
         [7,8,9]]}}
-> 0
{M.permanent
 {M.new [[1,2,3],
         [4,5,6],
         [7,8,9]]}}
-> 450

{M.determinant
 {M.new [[1,2,3],
         [4,5,6],
         [7,8,-9]]}}
-> 54
{M.permanent
 {M.new [[1,2,3],
         [4,5,6],
         [7,8,-9]]}}
-> 216

Lua[edit]

-- Johnson–Trotter permutations generator
_JT={}
function JT(dim)
  local n={ values={}, positions={}, directions={}, sign=1 }
  setmetatable(n,{__index=_JT})
  for i=1,dim do
    n.values[i]=i
    n.positions[i]=i
    n.directions[i]=-1
  end
  return n
end

function _JT:largestMobile()
  for i=#self.values,1,-1 do
    local loc=self.positions[i]+self.directions[i]
    if loc >= 1 and loc <= #self.values and self.values[loc] < i then
      return i
    end
  end
  return 0
end

function _JT:next()
  local r=self:largestMobile()
  if r==0 then return false end
  local rloc=self.positions[r]
  local lloc=rloc+self.directions[r]
  local l=self.values[lloc]
  self.values[lloc],self.values[rloc] = self.values[rloc],self.values[lloc]
  self.positions[l],self.positions[r] = self.positions[r],self.positions[l]
  self.sign=-self.sign
  for i=r+1,#self.directions do self.directions[i]=-self.directions[i] end
  return true
end  

-- matrix class

_MTX={}
function MTX(matrix)
  setmetatable(matrix,{__index=_MTX})
  matrix.rows=#matrix
  matrix.cols=#matrix[1]
  return matrix
end

function _MTX:dump()
  for _,r in ipairs(self) do
    print(unpack(r))
  end
end

function _MTX:perm() return self:det(1) end
function _MTX:det(perm)
  local det=0
  local jt=JT(self.cols)
  repeat
    local pi=perm or jt.sign
    for i,v in ipairs(jt.values) do
      pi=pi*self[i][v]
    end
    det=det+pi
  until not jt:next()
  return det
end

-- test

matrix=MTX
{
  { 7,  2, -2,  4},
  { 4,  4,  1,  7},
  {11, -8,  9, 10},
  {10,  5, 12, 13}
}
matrix:dump();
print("det:",matrix:det(), "permanent:",matrix:perm(),"\n")

matrix2=MTX
{
  {-2, 2,-3},
  {-1, 1, 3},
  { 2, 0,-1}
}
matrix2:dump();
print("det:",matrix2:det(), "permanent:",matrix2:perm())
Output:
7       2       -2      4
4       4       1       7
11      -8      9       10
10      5       12      13
det:    -4319   permanent:      10723

-2      2       -3
-1      1       3
2       0       -1
det:    18      permanent:      10

Maple[edit]

M:=<<2|9|4>,<7|5|3>,<6|1|8>>:

with(LinearAlgebra):

Determinant(M);
Permanent(M);

Output:

                                    -360
                                     900

Mathematica/Wolfram Language[edit]

Determinant is a built in function Det

Permanent is also a built in function, but here is a way it could be implemented:

Permanent[m_List] :=
    With[{v = Array[x, Length[m]]},
      Coefficient[Times @@ (m.v), Times @@ v]
  ]

Maxima[edit]

a: matrix([2, 9, 4], [7, 5, 3], [6, 1, 8])$

determinant(a);
-360

permanent(a);
900

МК-61/52[edit]

П4	ИПE	П2	КИП0	ИП0	П1	С/П	ИП4	/	КП2
L1	06	ИПE	П3	ИП0	П1	Сx	КП2	L1	17
ИП0	ИП2	+	П1	П2	ИП3	-	x#0	34	С/П
ПП	80	БП	21	КИП0	ИП4	С/П	КИП2	-	*
П4	ИП0	П3	x#0	35	Вx	С/П	КИП2	-	<->
/	КП1	L3	45	ИП1	ИП0	+	П3	ИПE	П1
П2	КИП1	/-/	ПП	80	ИП3	+	П3	ИП1	-
x=0	61	ИП0	П1	КИП3	КП2	L1	74	БП	12
ИП0	<->	^	КИП3	*	КИП1	+	КП2	->	L0
82	->	П0	В/О

This program calculates the determinant of the matrix of order <= 5. Prior to startup, РE entered 13, entered the order of the matrix Р0, and the elements are introduced with the launch of the program after one of them, the last on the screen will be determinant. Permanent is calculated in this way.

Nim[edit]

Translation of: Python

Using the permutationsswap module from Permutations by swapping:

import sequtils, permutationsswap

type Matrix[M,N: static[int]] = array[M, array[N, float]]

proc det[M,N](a: Matrix[M,N]): float =
  let n = toSeq 0..a.high
  for sigma, sign in n.permutations:
    var x = sign.float
    for i in n: x *= a[i][sigma[i]]
    result += x

proc perm[M,N](a: Matrix[M,N]): float =
  let n = toSeq 0..a.high
  for sigma, sign in n.permutations:
    var x = 1.0
    for i in n: x *= a[i][sigma[i]]
    result += x

const
  a = [ [1.0, 2.0]
      , [3.0, 4.0]
      ]
  b = [ [ 1.0,  2,  3,  4]
      , [ 4.0,  5,  6,  7]
      , [ 7.0,  8,  9, 10]
      , [10.0, 11, 12, 13]
      ]
  c = [ [ 0.0,  1,  2,  3,  4]
      , [ 5.0,  6,  7,  8,  9]
      , [10.0, 11, 12, 13, 14]
      , [15.0, 16, 17, 18, 19]
      , [20.0, 21, 22, 23, 24]
      ]

echo "perm: ", a.perm, " det: ", a.det
echo "perm: ", b.perm, " det: ", b.det
echo "perm: ", c.perm, " det: ", c.det

Output:

perm: 10.0 det: -2.0
perm: 29556.0 det: 0.0
perm: 6778800.0 det: 0.0

Ol[edit]

; helper function that returns rest of matrix by col/row
(define (rest matrix i j)
   (define (exclude1 l x) (append (take l (- x 1)) (drop l x)))
   (exclude1
      (map exclude1
         matrix (repeat i (length matrix)))
      j))

; superfunction for determinant and permanent
(define (super matrix math)
   (let loop ((n (length matrix)) (matrix matrix))
      (if (eq? n 1)
         (caar matrix)
         (fold (lambda (x a j)
                  (+ x (* a (lref math (mod j 2)) (super (rest matrix j 1) math))))
            0
            (car matrix)
            (iota n 1)))))


; det/per calculators
(define (det matrix) (super matrix '(-1 1)))
(define (per matrix) (super matrix '( 1 1)))

; ---=( testing )=---------------------
(print (det '(
   (1 2)
   (3 4))))
; ==> -2

(print (per '(
   (1 2)
   (3 4))))
; ==> 10


(print (det '(
   ( 1  2  3  1)
   (-1 -1 -1  2)
   ( 1  3  1  1)
   (-2 -2  0 -1))))
; ==> 26

(print (per '(
   ( 1  2  3  1)
   (-1 -1 -1  2)
   ( 1  3  1  1)
   (-2 -2  0 -1))))
; ==> -10


(print (det '(
   ( 0  1  2  3  4)
   ( 5  6  7  8  9)
   (10 11 12 13 14)
   (15 16 17 18 19)
   (20 21 22 23 24))))
; ==> 0

(print (per '(
   ( 0  1  2  3  4)
   ( 5  6  7  8  9)
   (10 11 12 13 14)
   (15 16 17 18 19)
   (20 21 22 23 24))))
; ==> 6778800

PARI/GP[edit]

The determinant is built in:

matdet(M)

and the permanent can be defined as

matperm(M)=my(n=#M,t);sum(i=1,n!,t=numtoperm(n,i);prod(j=1,n,M[j,t[j]]))

For better performance, here's a version using Ryser's formula:

matperm(M)=
{
	my(n=matsize(M)[1],innerSums=vectorv(n));
	if(n==0, return(1));
	sum(x=1,2^n-1,
		my(k=valuation(x,2),s=M[,k+1],gray=bitxor(x, x>>1));
		if(bittest(gray,k),
			innerSums += s;
		,
			innerSums -= s;
		);
		(-1)^hammingweight(gray)*factorback(innerSums)
	)*(-1)^n;
}
Works with: PARI/GP version 2.10.0+

As of version 2.10, the matrix permanent is built in:

matpermanent(M)

Perl[edit]

Translation of: C
#!/usr/bin/perl
use strict;
use warnings;
use PDL;
use PDL::NiceSlice;

sub permanent{
	my $mat = shift;
	my $n = shift // $mat->dim(0);
	return undef if $mat->dim(0) != $mat->dim(1);
	return $mat(0,0) if $n == 1;
	my $sum = 0;
	--$n;
	my $m = $mat(1:,1:)->copy;
	for(my $i = 0; $i <= $n; ++$i){
		$sum += $mat($i,0) * permanent($m, $n);
		last if $i == $n;
		$m($i,:) .= $mat($i,1:);
	}
	return sclr($sum);
}

my $M = pdl([[2,9,4], [7,5,3], [6,1,8]]);
print "M = $M\n";
print "det(M) = " . $M->determinant . ".\n";
print "det(M) = " . $M->det . ".\n";
print "perm(M) = " . permanent($M) . ".\n";

determinant and det are already defined in PDL, see[1]. permanent has to be defined manually.

Output:
M = 
[
 [2 9 4]
 [7 5 3]
 [6 1 8]
]

det(M) = -360.
det(M) = -360.
perm(M) = 900.

Phix[edit]

Translation of: Java
with javascript_semantics
function minor(sequence a, integer x, integer y)
    integer l = length(a)-1
    sequence result = repeat(repeat(0,l),l)
    for i=1 to l do
        for j=1 to l do
            result[i][j] = a[i+(i>=x)][j+(j>=y)]
        end for
    end for
    return result
end function
 
function det(sequence a)
    if length(a)=1 then
        return a[1][1]
    end if
    integer res = 0,
            sgn = 1
    for i=1 to length(a) do
        res += sgn*a[1][i]*det(minor(a,1,i))
        sgn *= -1
    end for
    return res
end function
 
function perm(sequence a)
    if length(a)=1 then
        return a[1][1]
    end if
    integer res = 0
    for i=1 to length(a) do
        res += a[1][i]*perm(minor(a,1,i))
    end for
    return res
end function
 
constant tests = {
{{1,  2},
 {3,  4}},
--Determinant: -2, permanent: 10
{{2, 9, 4},
 {7, 5, 3},
 {6, 1, 8}},
--Determinant: -360, permanent: 900
{{ 1,  2,  3,  4},
 { 4,  5,  6,  7},
 { 7,  8,  9, 10},
 {10, 11, 12, 13}},
--Determinant: 0, permanent: 29556
{{ 0,  1,  2,  3,  4},
 { 5,  6,  7,  8,  9},
 {10, 11, 12, 13, 14},
 {15, 16, 17, 18, 19},
 {20, 21, 22, 23, 24}},
--Determinant: 0, permanent: 6778800
{{5}},
--Determinant: 5, permanent: 5 
{{1,0,0},
 {0,1,0},
 {0,0,1}},
--Determinant: 1, permanent: 1
{{0,0,1},
 {0,1,0},
 {1,0,0}},
--Determinant: -1, Permanent: 1
{{4,3},
 {2,5}},
--Determinant: 14, Permanent: 26
{{2,5},
 {4,3}},
--Determinant: -14, Permanent: 26
{{4,4},
 {2,2}},
--Determinant: 0, Permanent: 16
{{7,    2,      -2,     4},
 {4,    4,      1,      7},
 {11,   -8,     9,      10},
 {10,   5,      12,     13}},
--det:  -4319   permanent:      10723
 
{{-2,   2,      -3},
 {-1,   1,      3},
 {2 ,   0,      -1}}
--det:  18      permanent:      10
}
for i=1 to length(tests) do
    sequence ti = tests[i]
    ?{det(ti),perm(ti)}
end for
Output:
{-2,10}
{-360,900}
{0,29556}
{0,6778800}
{5,5}
{1,1}
{-1,1}
{14,26}
{-14,26}
{0,16}
{-4319,10723}
{18,10}

PowerShell[edit]

function det-perm ($array) {
    if($array) {
        $size = $array.Count
        function prod($A) {
            $prod = 1
            if($A) { $A | foreach{$prod *= $_} }
            $prod
        }
        function generate($sign, $n, $A) {
            if($n -eq 1) {
                $i = 0               
                $prod = prod @($A | foreach{$array[$i++][$_]})
                [pscustomobject]@{det = $sign*$prod; perm = $prod}
            }
            else{
                for($i = 0; $i -lt ($n - 1); $i += 1) {
                    generate $sign ($n - 1) $A
                    if($n % 2 -eq 0){
                        $i1, $i2 = $i, ($n-1)
                        $A[$i1], $A[$i2] = $A[$i2], $A[$i1]
                    }
                    else{
                        $i1, $i2 = 0, ($n-1)
                        $A[$i1], $A[$i2] = $A[$i2], $A[$i1]
                    }
                    $sign *= -1
                }
                generate $sign ($n - 1) $A
            }
        }
        $det = $perm = 0
        generate 1 $size @(0..($size-1)) | foreach{
            $det += $_.det
            $perm += $_.perm
        }
        [pscustomobject]@{det =  "$det"; perm = "$perm"}
    } else {Write-Error "empty array"}
}
det-perm 5
det-perm @(@(1,0,0),@(0,1,0),@(0,0,1))
det-perm @(@(0,0,1),@(0,1,0),@(1,0,0))
det-perm @(@(4,3),@(2,5))
det-perm @(@(2,5),@(4,3))
det-perm @(@(4,4),@(2,2))

Output:

det                                                    perm                                                  
---                                                    ----                                                  
5                                                      5                                                     
1                                                      1                                                     
-1                                                     1                                                     
14                                                     26                                                    
-14                                                    26                                                    
0                                                      16                                                    

Python[edit]

Using the module file spermutations.py from Permutations by swapping. The algorithm for the determinant is a more literal translation of the expression in the task description and the Wikipedia reference.

from itertools import permutations
from operator import mul
from math import fsum
from spermutations import spermutations

def prod(lst):
    return reduce(mul, lst, 1)

def perm(a):
    n = len(a)
    r = range(n)
    s = permutations(r)
    return fsum(prod(a[i][sigma[i]] for i in r) for sigma in s)

def det(a):
    n = len(a)
    r = range(n)
    s = spermutations(n)
    return fsum(sign * prod(a[i][sigma[i]] for i in r)
                for sigma, sign in s)

if __name__ == '__main__':
    from pprint import pprint as pp

    for a in ( 
            [
             [1, 2], 
             [3, 4]], 

            [
             [1, 2, 3, 4],
             [4, 5, 6, 7],
             [7, 8, 9, 10],
             [10, 11, 12, 13]],        

            [
             [ 0,  1,  2,  3,  4],
             [ 5,  6,  7,  8,  9],
             [10, 11, 12, 13, 14],
             [15, 16, 17, 18, 19],
             [20, 21, 22, 23, 24]],
        ):
        print('')
        pp(a)
        print('Perm: %s Det: %s' % (perm(a), det(a)))
Sample output
[[1, 2], [3, 4]]
Perm: 10 Det: -2

[[1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10], [10, 11, 12, 13]]
Perm: 29556 Det: 0

[[0, 1, 2, 3, 4],
 [5, 6, 7, 8, 9],
 [10, 11, 12, 13, 14],
 [15, 16, 17, 18, 19],
 [20, 21, 22, 23, 24]]
Perm: 6778800 Det: 0

The second matrix above is that used in the Tcl example. The third matrix is from the J language example. Note that the determinant seems to be 'exact' using this method of calculation without needing to resort to other than Pythons default numbers.

R[edit]

R has matrix algebra built in, so we do not need to import anything when calculating the determinant. However, we will use a library to generate the permutations of 1:n.

library(combinat)
perm <- function(A)
{
  stopifnot(is.matrix(A))
  n <- nrow(A)
  if(n != ncol(A)) stop("Matrix is not square.")
  if(n < 1) stop("Matrix has a dimension of size 0.")
  sum(sapply(combinat::permn(n), function(sigma) prod(sapply(1:n, function(i) A[i, sigma[i]]))))
}

#We copy our test cases from the Python example.
testData <- list("Test 1" = rbind(c(1, 2), c(3, 4)),
                 "Test 2" = rbind(c(1, 2, 3, 4), c(4, 5, 6, 7), c(7, 8, 9, 10), c(10, 11, 12, 13)),
                 "Test 3" = rbind(c(0, 1, 2, 3, 4), c(5, 6, 7, 8, 9), c(10, 11, 12, 13, 14),
                               c(15, 16, 17, 18, 19), c(20, 21, 22, 23, 24)))
print(sapply(testData, function(x) list(Determinant = det(x), Permanent = perm(x))))
Output:
            Test 1 Test 2       Test 3 
Determinant -2     1.131522e-29 0      
Permanent   10     29556        6778800

Racket[edit]

#lang racket
(require math)
(define determinant matrix-determinant)

(define (permanent M)
  (define n (matrix-num-rows M))
  (for/sum ([σ (in-permutations (range n))])
    (for/product ([i n] [σi σ])
      (matrix-ref M i σi))))

Raku[edit]

(formerly Perl 6)

Works with: Rakudo version 2015.12

Uses the permutations generator from the Permutations by swapping task. This implementation is naive and brute-force (slow) but exact.

sub insert ($x, @xs) { ([flat @xs[0 ..^ $_], $x, @xs[$_ .. *]] for 0 .. @xs) }
sub order ($sg, @xs) { $sg > 0 ?? @xs !! @xs.reverse }

multi σ_permutations ([]) { [] => 1 }

multi σ_permutations ([$x, *@xs]) {
    σ_permutations(@xs).map({ |order($_.value, insert($x, $_.key)) }) Z=> |(1,-1) xx *
}

sub m_arith ( @a, $op ) {
    note "Not a square matrix" and return
      if [||] map { @a.elems cmp @a[$_].elems }, ^@a;
    sum σ_permutations([^@a]).race.map: {
        my $permutation = .key;
        my $term = $op eq 'perm' ?? 1 !! .value;
        for $permutation.kv -> $i, $j { $term *= @a[$i][$j] };
        $term
    }
}

######### helper subs #########
sub hilbert-matrix (\h) {[(1..h).map(-> \n {[(n..^n+h).map: {(1/$_).FatRat}]})]}

sub rat-or-int ($num) {
    return $num unless $num ~~ Rat|FatRat;
    return $num.narrow if $num.narrow.WHAT ~~ Int;
    $num.nude.join: '/';
}

sub say-it ($message, @array) {
    my $max;
    @array.map: {$max max= max $_».&rat-or-int.comb(/\S+/)».chars};
    say "\n$message";
    $_».&rat-or-int.fmt(" %{$max}s").put for @array;
}

########### Testing ###########
my @tests = (
    [
        [ 1, 2 ],
        [ 3, 4 ]
    ],
    [
        [  1,  2,  3,  4 ],
        [  4,  5,  6,  7 ],
        [  7,  8,  9, 10 ],
        [ 10, 11, 12, 13 ]
    ],
    hilbert-matrix 7
);

for @tests -> @matrix {
    say-it 'Matrix:', @matrix;
    say "Determinant:\t", rat-or-int @matrix.&m_arith: <det>;
    say "Permanent:  \t", rat-or-int @matrix.&m_arith: <perm>;
    say '-' x 40;
}

Output

Matrix:
 1  2
 3  4
Determinant:	-2
Permanent:  	10
----------------------------------------

Matrix:
  1   2   3   4
  4   5   6   7
  7   8   9  10
 10  11  12  13
Determinant:	0
Permanent:  	29556
----------------------------------------

Matrix:
    1   1/2   1/3   1/4   1/5   1/6   1/7
  1/2   1/3   1/4   1/5   1/6   1/7   1/8
  1/3   1/4   1/5   1/6   1/7   1/8   1/9
  1/4   1/5   1/6   1/7   1/8   1/9  1/10
  1/5   1/6   1/7   1/8   1/9  1/10  1/11
  1/6   1/7   1/8   1/9  1/10  1/11  1/12
  1/7   1/8   1/9  1/10  1/11  1/12  1/13
Determinant:	1/2067909047925770649600000
Permanent:  	29453515169174062608487/2067909047925770649600000
----------------------------------------

REXX[edit]

/* REXX ***************************************************************
* Test the two functions determinant and permanent
* using the matrix specifications shown for other languages
* 21.05.2013 Walter Pachl
**********************************************************************/
Call test ' 1  2',
          ' 3  4',2

Call test ' 1  2  3  4',
          ' 4  5  6  7',
          ' 7  8  9 10',
          '10 11 12 13',4

Call test ' 0  1  2  3  4',
          ' 5  6  7  8  9',
          '10 11 12 13 14',
          '15 16 17 18 19',
          '20 21 22 23 24',5

Exit

test:
/**********************************************************************
* Show the given matrix and compute and show determinant and permanent
**********************************************************************/
Parse Arg as,n
asc=as
Do i=1 To n
  ol=''
  Do j=1 To n
    Parse Var asc a.i.j asc
    ol=ol right(a.i.j,3)
    End
   Say ol
  End
Say 'determinant='right(determinant(as),7)
Say '  permanent='right(permanent(as),7)
Say copies('-',50)
Return
/* REXX ***************************************************************
* determinant.rex
* compute the determinant of the given square matrix
* Input: as: the representation of the matrix as vector (n**2 elements)
* 21.05.2013 Walter Pachl
**********************************************************************/
  Parse Arg as
  n=sqrt(words(as))
  Do i=1 To n
    Do j=1 To n
      Parse Var as a.i.j as
      End
    End
  Select
    When n=2 Then det=a.1.1*a.2.2-a.1.2*a.2.1
    When n=3 Then det= a.1.1*a.2.2*a.3.3,
                      +a.1.2*a.2.3*a.3.1,
                      +a.1.3*a.2.1*a.3.2,
                      -a.1.3*a.2.2*a.3.1,
                      -a.1.2*a.2.1*a.3.3,
                      -a.1.1*a.2.3*a.3.2
    Otherwise Do
      det=0
      Do k=1 To n
        det=det+((-1)**(k+1))*a.1.k*determinant(subm(k))
        End
      End
    End
  Return det

subm: Procedure Expose a. n
/**********************************************************************
* compute the submatrix resulting when row 1 and column k are removed
* Input: a.*.*, k
* Output: bs the representation of the submatrix as vector
**********************************************************************/
  Parse Arg k
  bs=''
  do i=2 To n
    Do j=1 To n
      If j=k Then Iterate
      bs=bs a.i.j
      End
    End
  Return bs

sqrt: Procedure
/**********************************************************************
* compute and return the (integer) square root of the given argument
* terminate the program if the argument is not a square
**********************************************************************/
  Parse Arg nn
  Do n=1 By 1 while n*n<nn
    End
  If n*n=nn Then
    Return n
  Else Do
    Say 'invalid number of elements:' nn 'is not a square.'
    Exit
    End
/* REXX ***************************************************************
* permanent.rex
* compute the permanent of a matrix
* I found an algorithm here:
* http://www.codeproject.com/Articles/21282/Compute-Permanent-of-a-Matrix-with-Ryser-s-Algorit
* see there for the original author.
* translated it to REXX (hopefully correctly) to REXX
* and believe that I can "publish" it here, on rosettacode
* when I look at the copyright rules shown there:
* http://www.codeproject.com/info/cpol10.aspx
* 20.05.2013 Walter Pachl
**********************************************************************/
Call init arg(1)                 /* initialize the matrix (n and a.* */
sum=0
rowsumprod=0
rowsum=0
chi.=0
c=2**n
Do k=1 To c-1                       /* loop all 2^n submatrices of A */
  rowsumprod = 1
  chis=dec2binarr(k,n)              /* characteristic vector         */
  Do ci=0 By 1 While chis<>''
    Parse Var chis chi.ci chis
    End
  Do m=0 To n-1                     /* loop columns of submatrix #k  */
    rowsum = 0
    Do p=0 To n-1                   /* loop rows and compute rowsum  */
      mnp=m*n+p
      rowsum=rowsum+chi.p*A.mnp
      End
    rowsumprod=rowsumprod*rowsum  /* update product of rowsums     */
                            /* (optional -- use for sparse matrices) */
                            /* if (rowsumprod == 0) break;           */
    End
  sum=sum+((-1)**(n-chi.n))*rowsumprod
  End
Return sum
/**********************************************************************
* Notes
* 1.The submatrices are chosen by use of a characteristic vector chi
* (only the columns are considered, where chi[p] == 1).
* To retrieve the t from Ryser's formula, we need to save the number
* n-t, as is done in chi[n]. Then we get t = n - chi[n].
* 2.The matrix parameter A is expected to be a one-dimensional integer
* array -- should the matrix be encoded row-wise or column-wise?
* -- It doesn't matter. The permanent is invariant under
* row-switching and column-switching, and it is Screenshot
* - per_inv.gif .
* 3.Further enhancements: If any rowsum equals zero,
* the entire rowsumprod becomes zero, and thus the m-loop can be broken.
* Since if-statements are relatively expensive compared to integer
* operations, this might save time only for sparse matrices
* (where most entries are zeros).
* 4.If anyone finds a polynomial algorithm for permanents,
* he will get rich and famous (at least in the computer science world).
**********************************************************************/
/**********************************************************************
* At first, we need to transform a decimal to a binary array
* with an additional element
* (the last one) saving the number of ones in the array:
**********************************************************************/
dec2binarr: Procedure
  Parse Arg n,dim
  ol='n='n 'dim='dim
  res.=0
  pos=dim-1
  Do While n>0
    res.pos=n//2
    res.dim=res.dim+res.pos
    n=n%2
    pos=pos-1
    End
  res_s=''
  Do i=0 To dim
    res_s=res_s res.i
    End
  Return res_s

init: Procedure Expose a. n
/**********************************************************************
* a.* (starting with index 0) contains all array elements
* n is the dimension of the square matrix
**********************************************************************/
Parse Arg as
n=sqrt(words(as))
a.=0
Do ai=0 By 1 While as<>''
  Parse Var as a.ai as
  End
Return

sqrt: Procedure
/**********************************************************************
* compute and return the (integer) square root of the given argument
* terminate the program if the argument is not a square
**********************************************************************/
  Parse Arg nn
  Do n=1 By 1 while n*n<nn
    End
  If n*n=nn Then
    Return n
  Else Do
    Say 'invalid number of elements:' nn 'is not a square.'
    Exit
    End

Output:

   1   2
   3   4
determinant=     -2
  permanent=     10
--------------------------------------------------
   1   2   3   4
   4   5   6   7
   7   8   9  10
  10  11  12  13
determinant=      0
  permanent=  29556
--------------------------------------------------
   0   1   2   3   4
   5   6   7   8   9
  10  11  12  13  14
  15  16  17  18  19
  20  21  22  23  24
determinant=      0
  permanent=6778800
--------------------------------------------------

Ruby[edit]

Matrix in the standard library provides a method for the determinant, but not for the permanent.

require 'matrix'

class Matrix
  # Add "permanent" method to Matrix class
  def permanent
    r = (0...row_count).to_a # [0,1] (first example), [0,1,2,3] (second example)
    r.permutation.inject(0) do |sum, sigma|
       sum += sigma.zip(r).inject(1){|prod, (row, col)| prod *= self[row, col] }
    end
  end
end

m1 = Matrix[[1,2],[3,4]] # testcases from Python version

m2 = Matrix[[1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10], [10, 11, 12, 13]]

m3 = Matrix[[0, 1, 2, 3, 4],
            [5, 6, 7, 8, 9],
            [10, 11, 12, 13, 14],
            [15, 16, 17, 18, 19], 
            [20, 21, 22, 23, 24]]

[m1, m2, m3].each do |m|
  puts "determinant:\t #{m.determinant}", "permanent:\t #{m.permanent}"
  puts 
end
Output:
determinant:	 -2
permanent:	 10

determinant:	 0
permanent:	 29556

determinant:	 0
permanent:	 6778800

Rust[edit]

Translation of: Java
fn main() {
let mut m1: Vec<Vec<f64>> = vec![vec![1.0,2.0],vec![3.0,4.0]];
    let mut r_m1 = &mut m1;
    let rr_m1 = &mut r_m1;
    let mut m2: Vec<Vec<f64>> = vec![vec![1.0, 2.0, 3.0, 4.0], vec![4.0, 5.0, 6.0, 7.0], vec![7.0, 8.0, 9.0, 10.0], vec![10.0, 11.0, 12.0, 13.0]];
    let mut r_m2 = &mut m2;
    let rr_m2 = &mut r_m2;
    let mut m3: Vec<Vec<f64>> = vec![vec![0.0, 1.0, 2.0, 3.0, 4.0],
                                vec![5.0, 6.0, 7.0, 8.0, 9.0],
                                vec![10.0, 11.0, 12.0, 13.0, 14.0],
                                vec![15.0, 16.0, 17.0, 18.0, 19.0], 
                                vec![20.0, 21.0, 22.0, 23.0, 24.0]];
    let mut r_m3 = &mut m3;
    let rr_m3 = &mut r_m3;

    println!("Determinant of m1: {}", determinant(rr_m1));
    println!("Permanent of m1: {}", permanent(rr_m1));

    println!("Determinant of m2: {}", determinant(rr_m2));
    println!("Permanent of m2: {}", permanent(rr_m2));

    println!("Determinant of m3: {}", determinant(rr_m3));
    println!("Permanent of m3: {}", permanent(rr_m3));

}

fn minor( a: &mut Vec<Vec<f64>>, x: usize, y: usize) ->  Vec<Vec<f64>> {
    let mut out_vec: Vec<Vec<f64>> = vec![vec![0.0; a.len() - 1]; a.len() -1];
    for i in 0..a.len()-1 {
        for j in 0..a.len()-1 {
            match () {
                _ if (i < x && j < y) => {
                    out_vec[i][j] = a[i][j];
                },
                _ if (i >= x && j < y) => {
                    out_vec[i][j] = a[i + 1][j];
                },
                _ if (i < x && j >= y) => {
                    out_vec[i][j] = a[i][j + 1];
                },
                _ => { //i > x && j > y 
                    out_vec[i][j] = a[i + 1][j + 1];
                },
            }
        }
    }
    out_vec
}

fn determinant (matrix: &mut Vec<Vec<f64>>) -> f64 {
    match () {
        _ if (matrix.len() == 1) => {
            matrix[0][0]
        },
        _ => {
            let mut sign = 1.0;
            let mut sum = 0.0;
            for i in 0..matrix.len() {
                sum = sum + sign * matrix[0][i] * determinant(&mut minor(matrix, 0, i));
                sign = sign * -1.0;
            }
            sum
        }
    }
}

fn permanent (matrix: &mut Vec<Vec<f64>>) -> f64 {
    match () {
        _ if (matrix.len() == 1) => {
            matrix[0][0]
        },
        _ => {
            let mut sum = 0.0;
            for i in 0..matrix.len() {
                sum = sum + matrix[0][i] * permanent(&mut minor(matrix, 0, i));
            }
            sum
        }
    }
}
Output:
Determinant of m1: -2
Permanent of m1: 10
Determinant of m2: 0
Permanent of m2: 29556
Determinant of m3: 0
Permanent of m3: 6778800

Scala[edit]

def permutationsSgn[T]: List[T] => List[(Int,List[T])] = {
  case Nil => List((1,Nil))
  case xs => {
    for {
      (x, i) <- xs.zipWithIndex
      (sgn,ys) <- permutationsSgn(xs.take(i) ++ xs.drop(1 + i))
    } yield {
      val sgni = sgn * (2 * (i%2) - 1)
      (sgni, (x :: ys))
    }
  }
}

def det(m:List[List[Int]]) = {
  val summands =
    for {
      (sgn,sigma) <- permutationsSgn((0 to m.length - 1).toList).toList
    }
    yield {
      val factors =
        for (i <- 0 to (m.length - 1))
        yield m(i)(sigma(i))
      factors.toList.foldLeft(sgn)({case (x,y) => x * y})
    }
  summands.toList.foldLeft(0)({case (x,y) => x + y})

Sidef[edit]

The `determinant` method is provided by the Array class.

Translation of: Ruby
class Array {
    method permanent {
        var r = @^self.len

        var sum = 0
        r.permutations { |*a|
            var prod = 1
            [a,r].zip {|row,col| prod *= self[row][col] }
            sum += prod
        }

        return sum
    }
}

var m1 = [[1,2],[3,4]]

var m2 = [[1, 2, 3, 4],
          [4, 5, 6, 7],
          [7, 8, 9, 10],
          [10, 11, 12, 13]]

var m3 = [[0, 1, 2, 3, 4],
          [5, 6, 7, 8, 9],
          [10, 11, 12, 13, 14],
          [15, 16, 17, 18, 19],
          [20, 21, 22, 23, 24]]

[m1, m2, m3].each { |m|
  say "determinant:\t #{m.determinant}\npermanent:\t #{m.permanent}\n"
}
Output:
determinant:	 -2
permanent:	 10

determinant:	 0
permanent:	 29556

determinant:	 0
permanent:	 6778800

Simula[edit]

! MATRIX ARITHMETIC ;
BEGIN

    INTEGER PROCEDURE LENGTH(A); ARRAY A;
        LENGTH := UPPERBOUND(A, 1) - LOWERBOUND(A, 1) + 1;

    ! Set MAT to the first minor of A dropping row X and column Y ;
    PROCEDURE MINOR(A, X, Y, MAT); ARRAY A, MAT; INTEGER X, Y;
    BEGIN
        INTEGER I, J, rowA, M; M := LENGTH(A) - 1; ! not a constant;
        FOR I := 1 STEP 1 UNTIL M DO BEGIN
            rowA := IF I < X THEN I ELSE I + 1;
            FOR J := 1 STEP 1 UNTIL M DO
                MAT(I, J) := A(rowA, IF J < Y THEN J else J + 1);
        END
    END MINOR;

    REAL PROCEDURE DET(A); REAL ARRAY A;
    BEGIN
        INTEGER N; N := LENGTH(A);
        IF N = 1 THEN
            DET := A(1, 1)
        ELSE 
        BEGIN
            INTEGER I, SIGN;
            REAL SUM;
            SIGN := 1;
            FOR I := 1 STEP 1 UNTIL N DO
            BEGIN
                REAL ARRAY MAT(1:N-1, 1:N-1);
                MINOR(A, 1, I, MAT);
                SUM := SUM + SIGN * A(1, I) * DET(MAT);
                SIGN := SIGN * -1
            END;
            DET := SUM
        END
    END DET;

    REAL PROCEDURE PERM(A); REAL ARRAY A;
    BEGIN
        INTEGER N; N := LENGTH(A);
        IF N = 1 THEN
            PERM := A(1, 1)
        ELSE 
        BEGIN
            REAL SUM;
            INTEGER I;

            FOR I := 1 STEP 1 UNTIL N DO
            BEGIN
                REAL ARRAY MAT(1:N-1, 1:N-1);
                MINOR(A, 1, I, MAT);
                SUM := SUM  + A(1, I) * PERM(MAT)
            END;
            PERM := SUM
        END
    END PERM;

    INTEGER SIZE;
    SIZE := ININT;
    BEGIN
        REAL ARRAY A(1:SIZE, 1:SIZE);
        INTEGER I, J;

        FOR I := 1 STEP 1 UNTIL SIZE DO BEGIN
            ! may be need here: INIMAGE;
            FOR J := 1 STEP 1 UNTIL SIZE DO
                A(I, J) := INREAL
        END;
        OUTTEXT("DETERMINANT ... : "); OUTREAL(DET (A), 10, 20); OUTIMAGE;
        OUTTEXT("PERMANENT ..... : "); OUTREAL(PERM(A), 10, 20); OUTIMAGE;
    END

    COMMENT  THE FIRST INPUT IS THE SIZE OF THE MATRIX, FOR EXAMPLE:

    ! 2
    ! 1 2
    ! 3 4
    ! DETERMINANT: -2.0
    ! PERMANENT: 10.0 ;
     
    COMMENT
    ! 5
    ! 0 1 2 3 4
    ! 5 6 7 8 9
    ! 10 11 12 13 14
    ! 15 16 17 18 19
    ! 20 21 22 23 24
    ! DETERMINANT: 0.0
    ! PERMANENT: 6778800.0 ;

END

Input:

2
1 2
3 4
Output:
DETERMINANT ... :    -2.000000000&+000
PERMANENT ..... :     1.000000000&+001

Input:

5
0 1 2 3 4
5 6 7 8 9
10 11 12 13 14
15 16 17 18 19
20 21 22 23 24
Output:
DETERMINANT ... :     0.000000000&+000
PERMANENT ..... :     6.778800000&+006

SPAD[edit]

Works with: FriCAS
Works with: OpenAxiom
Works with: Axiom
(1) -> M:=matrix [[2, 9, 4], [7, 5, 3], [6, 1, 8]]

        +2  9  4+
        |       |
   (1)  |7  5  3|
        |       |
        +6  1  8+
                                                        Type: Matrix(Integer)
(2) -> determinant M

   (2)  - 360
                                                                Type: Integer
(3) -> permanent M

   (3)  900
                                                        Type: PositiveInteger

Domain:Matrix(R)

Stata[edit]

Two auxiliary functions: range1(n,i) returns the column vector with numbers 1 to n except i is removed. And submat(a,i,j) returns matrix a with row i and column j removed. For x=-1, the main function sumrec(a,x) computes the determinant of a by developing the determinant along the first column. For x=1, one gets the permanent.

real vector range1(real scalar n, real scalar i) {
	if (i < 1 | i > n) {
		return(1::n)
	} else if (i == 1) {
		return(2::n)
	} else if (i == n) {
		return(1::n-1)
	} else {
		return(1::i-1\i+1::n)
	}
}

real matrix submat(real matrix a, real scalar i, real scalar j) {
	return(a[range1(rows(a), i), range1(cols(a), j)])
}

real scalar sumrec(real matrix a, real scalar x) {
	real scalar n, s, p
	n = rows(a)
	if (n==1) return(a[1,1])
	s = 0
	p = 1
	for (i=1; i<=n; i++) {
		s = s+p*a[i,1]*sumrec(submat(a, i, 1), x)
		p = p*x
	}
	return(s)
}

Example:

: a=1,1,1,0\1,1,0,1\1,0,1,1\0,1,1,1
: a
[symmetric]
       1   2   3   4
    +-----------------+
  1 |  1              |
  2 |  1   1          |
  3 |  1   0   1      |
  4 |  0   1   1   1  |
    +-----------------+

: det(a)
  -3

: sumrec(a,-1)
  -3

: sumrec(a,1)
  9

Tcl[edit]

The determinant is provided by the linear algebra package in Tcllib. The permanent (being somewhat less common) requires definition, but is easily described:

Library: Tcllib (Package: math::linearalgebra)
Library: Tcllib (Package: struct::list)
package require math::linearalgebra
package require struct::list

proc permanent {matrix} {
    for {set plist {};set i 0} {$i<[llength $matrix]} {incr i} {
	lappend plist $i
    }
    foreach p [::struct::list permutations $plist] {
	foreach i $plist j $p {
	    lappend prod [lindex $matrix $i $j]
	}
	lappend sum [::tcl::mathop::* {*}$prod[set prod {}]]
    }
    return [::tcl::mathop::+ {*}$sum]
}

Demonstrating with a sample matrix:

set mat {
    {1 2 3 4}
    {4 5 6 7}
    {7 8 9 10}
    {10 11 12 13}
}
puts [::math::linearalgebra::det $mat]
puts [permanent $mat]
Output:
1.1315223609263888e-29
29556

Visual Basic .NET[edit]

Translation of: Java
Module Module1

    Function Minor(a As Double(,), x As Integer, y As Integer) As Double(,)
        Dim length = a.GetLength(0) - 1
        Dim result(length - 1, length - 1) As Double
        For i = 1 To length
            For j = 1 To length
                If i < x AndAlso j < y Then
                    result(i - 1, j - 1) = a(i - 1, j - 1)
                ElseIf i >= x AndAlso j < y Then
                    result(i - 1, j - 1) = a(i, j - 1)
                ElseIf i < x AndAlso j >= y Then
                    result(i - 1, j - 1) = a(i - 1, j)
                Else
                    result(i - 1, j - 1) = a(i, j)
                End If
            Next
        Next
        Return result
    End Function

    Function Det(a As Double(,)) As Double
        If a.GetLength(0) = 1 Then
            Return a(0, 0)
        Else
            Dim sign = 1
            Dim sum = 0.0
            For i = 1 To a.GetLength(0)
                sum += sign * a(0, i - 1) * Det(Minor(a, 0, i))
                sign *= -1
            Next
            Return sum
        End If
    End Function

    Function Perm(a As Double(,)) As Double
        If a.GetLength(0) = 1 Then
            Return a(0, 0)
        Else
            Dim sum = 0.0
            For i = 1 To a.GetLength(0)
                sum += a(0, i - 1) * Perm(Minor(a, 0, i))
            Next
            Return sum
        End If
    End Function

    Sub WriteLine(a As Double(,))
        For i = 1 To a.GetLength(0)
            Console.Write("[")
            For j = 1 To a.GetLength(1)
                If j > 1 Then
                    Console.Write(", ")
                End If
                Console.Write(a(i - 1, j - 1))
            Next
            Console.WriteLine("]")
        Next
    End Sub

    Sub Test(a As Double(,))
        If a.GetLength(0) <> a.GetLength(1) Then
            Throw New ArgumentException("The dimensions must be equal")
        End If

        WriteLine(a)
        Console.WriteLine("Permanant  : {0}", Perm(a))
        Console.WriteLine("Determinant: {0}", Det(a))
        Console.WriteLine()
    End Sub

    Sub Main()
        Test({{1, 2}, {3, 4}})
        Test({{1, 2, 3, 4}, {4, 5, 6, 7}, {7, 8, 9, 10}, {10, 11, 12, 13}})
        Test({{0, 1, 2, 3, 4}, {5, 6, 7, 8, 9}, {10, 11, 12, 13, 14}, {15, 16, 17, 18, 19}, {20, 21, 22, 23, 24}})
    End Sub

End Module
Output:
[1, 2]
[3, 4]
Permanant  : 10
Determinant: -2

[1, 2, 3, 4]
[4, 5, 6, 7]
[7, 8, 9, 10]
[10, 11, 12, 13]
Permanant  : 29556
Determinant: 0

[0, 1, 2, 3, 4]
[5, 6, 7, 8, 9]
[10, 11, 12, 13, 14]
[15, 16, 17, 18, 19]
[20, 21, 22, 23, 24]
Permanant  : 6778800
Determinant: 0

Press any key to continue . . .

VBA[edit]

Translation of: Phix

As an extra, the results of the built in WorksheetFuction.MDeterm are shown. The latter does not work for scalars.

Option Base 1
Private Function minor(a As Variant, x As Integer, y As Integer) As Variant
    Dim l As Integer: l = UBound(a) - 1
    Dim result() As Double
    If l > 0 Then ReDim result(l, l)
    For i = 1 To l
        For j = 1 To l
            result(i, j) = a(i - (i >= x), j - (j >= y))
        Next j
    Next i
    minor = result
End Function
 
Private Function det(a As Variant)
    If IsArray(a) Then
        If UBound(a) = 1 Then
            On Error GoTo err
            det = a(1, 1)
            Exit Function
        End If
    Else
        det = a
        Exit Function
    End If
    Dim sgn_ As Integer: sgn_ = 1
    Dim res As Integer: res = 0
    Dim i As Integer
    For i = 1 To UBound(a)
        res = res + sgn_ * a(1, i) * det(minor(a, 1, i))
        sgn_ = sgn_ * -1
    Next i
    det = res
    Exit Function
err:
    det = a(1)
End Function
 
Private Function perm(a As Variant) As Double
    If IsArray(a) Then
        If UBound(a) = 1 Then
            On Error GoTo err
            perm = a(1, 1)
            Exit Function
        End If
    Else
        perm = a
        Exit Function
    End If
    Dim res As Double
    Dim i As Integer
    For i = 1 To UBound(a)
        res = res + a(1, i) * perm(minor(a, 1, i))
    Next i
    perm = res
    Exit Function
err:
    perm = a(1)
End Function
 
Public Sub main()
    Dim tests(13) As Variant
    tests(1) = [{1,  2; 3,  4}]
    '--Determinant: -2, permanent: 10
    tests(2) = [{2, 9, 4; 7, 5, 3; 6, 1, 8}]
    '--Determinant: -360, permanent: 900
    tests(3) = [{ 1,  2,  3,  4; 4,  5,  6,  7; 7,  8,  9, 10; 10, 11, 12, 13}]
    '--Determinant: 0, permanent: 29556
    tests(4) = [{ 0,  1,  2,  3,  4; 5,  6,  7,  8,  9; 10, 11, 12, 13, 14; 15, 16, 17, 18, 19; 20, 21, 22, 23, 24}]
    '--Determinant: 0, permanent: 6778800
    tests(5) = [{5}]
    '--Determinant: 5, permanent: 5
    tests(6) = [{1,0,0; 0,1,0; 0,0,1}]
    '--Determinant: 1, permanent: 1
    tests(7) = [{0,0,1; 0,1,0; 1,0,0}]
    '--Determinant: -1, Permanent: 1
    tests(8) = [{4,3; 2,5}]
    '--Determinant: 14, Permanent: 26
    tests(9) = [{2,5; 4,3}]
    '--Determinant: -14, Permanent: 26
    tests(10) = [{4,4; 2,2}]
    '--Determinant: 0, Permanent: 16
    tests(11) = [{7,    2,      -2,     4; 4,    4,      1,      7; 11,   -8,     9,      10; 10,   5,      12,     13}]
    '--det:  -4319   permanent:      10723
    tests(12) = [{-2,   2,      -3; -1,   1,      3; 2 ,   0,      -1}]
    '--det:  18      permanent:      10
    tests(13) = 13
    Debug.Print "Determinant", "Builtin det", "Permanent"
    For i = 1 To 12
        Debug.Print det(tests(i)), WorksheetFunction.MDeterm(tests(i)), perm(tests(i))
    Next i
    Debug.Print det(tests(13)), "error", perm(tests(13))
End Sub
Output:
Determinant   Builtin det   Permanent
-2            -2             10 
-360          -360           900 
 0             0             29556 
 0             0             6778800 
 5             5             5 
 1             1             1 
-1            -1             1 
 14            14            26 
-14           -14            26 
 0             0             16 
-4319         -4319          10723 
 18            18            10 
 13           error          13 

Wren[edit]

Library: Wren-matrix
Library: Wren-fmt
import "/matrix" for Matrix
import "/fmt" for Fmt

var arrays = [
    [ [1, 2],
      [3, 4] ],
 
    [ [-2, 2, -3],
      [-1, 1,  3],
      [ 2, 0, -1] ],
 
    [ [ 1,  2,  3,  4],
      [ 4,  5,  6,  7],
      [ 7,  8,  9, 10],
      [10, 11, 12, 13] ],
 
    [ [ 0,  1,  2,  3,  4],
      [ 5,  6,  7,  8,  9],
      [10, 11, 12, 13, 14],
      [15, 16, 17, 18, 19],
      [20, 21, 22, 23, 24] ]
]

for (array in arrays) {
    var m = Matrix.new(array)
    Fmt.mprint(m, 2, 0)
    System.print("\nDeterminant: %(m.det)")
    System.print("Permanent  : %(m.perm)\n") 
}
Output:
| 1  2|
| 3  4|

Determinant: -2
Permanent  : 10

|-2  2 -3|
|-1  1  3|
| 2  0 -1|

Determinant: 18
Permanent  : 10

| 1  2  3  4|
| 4  5  6  7|
| 7  8  9 10|
|10 11 12 13|

Determinant: 0
Permanent  : 29556

| 0  1  2  3  4|
| 5  6  7  8  9|
|10 11 12 13 14|
|15 16 17 18 19|
|20 21 22 23 24|

Determinant: 0
Permanent  : 6778800

zkl[edit]

var [const] GSL=Import("zklGSL");	// libGSL (GNU Scientific Library)
fcn perm(A){  // should verify A is square
   numRows:=A.rows;
   Utils.Helpers.permute(numRows.toList()).reduce(  // permute(0,1,..numRows)
      'wrap(s,pm){ s + numRows.reduce('wrap(x,i){ x*A[i,pm[i]] },1.0) },
      0.0)
}
test:=fcn(A){
   println(A.format());
   println("Permanent: %.2f, determinant: %.2f".fmt(perm(A),A.det()));
};
A:=GSL.Matrix(2,2).set(1,2, 3,4);
B:=GSL.Matrix(4,4).set(1,2,3,4, 4,5,6,7, 7,8,9,10, 10,11,12,13);
C:=GSL.Matrix(5,5).set( 0, 1, 2, 3, 4,  5, 6, 7, 8, 9, 10,11,12,13,14, 
		       15,16,17,18,19, 20,21,22,23,24);
T(A,B,C).apply2(test);
Output:
      1.00,      2.00
      3.00,      4.00
Permanent: 10.00, determinant: -2.00
      1.00,      2.00,      3.00,      4.00
      4.00,      5.00,      6.00,      7.00
      7.00,      8.00,      9.00,     10.00
     10.00,     11.00,     12.00,     13.00
Permanent: 29556.00, determinant: 0.00
      0.00,      1.00,      2.00,      3.00,      4.00
      5.00,      6.00,      7.00,      8.00,      9.00
     10.00,     11.00,     12.00,     13.00,     14.00
     15.00,     16.00,     17.00,     18.00,     19.00
     20.00,     21.00,     22.00,     23.00,     24.00
Permanent: 6778800.00, determinant: 0.00