MATLAB Examples

## Matrix Permanent tests

Here I'll generate a sequence of magic squares, then compute the permanent for each. I'll also show the time required.

Note: Since the result for order 9 and 10 matrices are integers, yet greater than 2^53-1, they will not be exact, since they are floating point numbers too large to fit completely into a double precision flint.

Strong recommendation: Do not go beyond 10x10 matrices to compute a permanent. This computation gets seriously large for bigger problems.

```for n = 2:10 disp(' ') disp(['Permanent of an nxn magic square, with n = ',num2str(n)]) A = magic(n) tic,P = permanent(A),toc end ```
``` Permanent of an nxn magic square, with n = 2 A = 1 3 4 2 P = 14 Elapsed time is 0.000327 seconds. Permanent of an nxn magic square, with n = 3 A = 8 1 6 3 5 7 4 9 2 P = 900 Elapsed time is 0.000253 seconds. Permanent of an nxn magic square, with n = 4 A = 16 2 3 13 5 11 10 8 9 7 6 12 4 14 15 1 P = 153716 Elapsed time is 0.000283 seconds. Permanent of an nxn magic square, with n = 5 A = 17 24 1 8 15 23 5 7 14 16 4 6 13 20 22 10 12 19 21 3 11 18 25 2 9 P = 53131650 Elapsed time is 0.000352 seconds. Permanent of an nxn magic square, with n = 6 A = 35 1 6 26 19 24 3 32 7 21 23 25 31 9 2 22 27 20 8 28 33 17 10 15 30 5 34 12 14 16 4 36 29 13 18 11 P = 34387479996 Elapsed time is 0.000917 seconds. Permanent of an nxn magic square, with n = 7 A = 30 39 48 1 10 19 28 38 47 7 9 18 27 29 46 6 8 17 26 35 37 5 14 16 25 34 36 45 13 15 24 33 42 44 4 21 23 32 41 43 3 12 22 31 40 49 2 11 20 P = 36757635682250 Elapsed time is 0.001102 seconds. Permanent of an nxn magic square, with n = 8 A = 64 2 3 61 60 6 7 57 9 55 54 12 13 51 50 16 17 47 46 20 21 43 42 24 40 26 27 37 36 30 31 33 32 34 35 29 28 38 39 25 41 23 22 44 45 19 18 48 49 15 14 52 53 11 10 56 8 58 59 5 4 62 63 1 P = 6.1755405170642e+16 Elapsed time is 0.006958 seconds. Permanent of an nxn magic square, with n = 9 A = 47 58 69 80 1 12 23 34 45 57 68 79 9 11 22 33 44 46 67 78 8 10 21 32 43 54 56 77 7 18 20 31 42 53 55 66 6 17 19 30 41 52 63 65 76 16 27 29 40 51 62 64 75 5 26 28 39 50 61 72 74 4 15 36 38 49 60 71 73 3 14 25 37 48 59 70 81 2 13 24 35 P = 1.41901767733966e+20 Elapsed time is 0.089605 seconds. Permanent of an nxn magic square, with n = 10 A = 92 99 1 8 15 67 74 51 58 40 98 80 7 14 16 73 55 57 64 41 4 81 88 20 22 54 56 63 70 47 85 87 19 21 3 60 62 69 71 28 86 93 25 2 9 61 68 75 52 34 17 24 76 83 90 42 49 26 33 65 23 5 82 89 91 48 30 32 39 66 79 6 13 95 97 29 31 38 45 72 10 12 94 96 78 35 37 44 46 53 11 18 100 77 84 36 43 50 27 59 P = 4.73325614942977e+23 Elapsed time is 1.233916 seconds. ```

## Permanent applies to non-integer matrices, boolean, or even sparse matrices

Be careful with integer types like uint8, since the permanent can be a very large number, and overflows can easily occurr.

```Philb = permanent(hilb(5)) Pbool = permanent(rand(10) > 0.5) Psparse = permanent(sprand(8,8,.5)) ```
```Philb = 0.068250218962585 Pbool = 2088 Psparse = 0.170886975922394 ```

## Symbolic verification for 2x2, 3x3 and 4x4 matrices

It is fairly simple to test small matices to show the result is correct. I'll guess you'll just need to check that these next results have all the necessary terms, or you can take my word on that claim.

```A2 = sym('A',[2,2]) P2 = permanent(A2) A3 = sym('A',[3,3]) P3 = permanent(A3) A4 = sym('A',[4,4]) P4 = permanent(A4) ```
```A2 = [ A1_1, A1_2] [ A2_1, A2_2] P2 = A1_1*A2_2 + A1_2*A2_1 A3 = [ A1_1, A1_2, A1_3] [ A2_1, A2_2, A2_3] [ A3_1, A3_2, A3_3] P3 = A1_1*A2_2*A3_3 + A1_1*A2_3*A3_2 + A1_2*A2_1*A3_3 + A1_2*A2_3*A3_1 + A1_3*A2_1*A3_2 + A1_3*A2_2*A3_1 A4 = [ A1_1, A1_2, A1_3, A1_4] [ A2_1, A2_2, A2_3, A2_4] [ A3_1, A3_2, A3_3, A3_4] [ A4_1, A4_2, A4_3, A4_4] P4 = A1_1*A2_2*A3_3*A4_4 + A1_1*A2_2*A3_4*A4_3 + A1_1*A2_3*A3_2*A4_4 + A1_1*A2_3*A3_4*A4_2 + A1_1*A2_4*A3_2*A4_3 + A1_1*A2_4*A3_3*A4_2 + A1_2*A2_1*A3_3*A4_4 + A1_2*A2_1*A3_4*A4_3 + A1_2*A2_3*A3_1*A4_4 + A1_2*A2_3*A3_4*A4_1 + A1_2*A2_4*A3_1*A4_3 + A1_2*A2_4*A3_3*A4_1 + A1_3*A2_1*A3_2*A4_4 + A1_3*A2_1*A3_4*A4_2 + A1_3*A2_2*A3_1*A4_4 + A1_3*A2_2*A3_4*A4_1 + A1_3*A2_4*A3_1*A4_2 + A1_3*A2_4*A3_2*A4_1 + A1_4*A2_1*A3_2*A4_3 + A1_4*A2_1*A3_3*A4_2 + A1_4*A2_2*A3_1*A4_3 + A1_4*A2_2*A3_3*A4_1 + A1_4*A2_3*A3_1*A4_2 + A1_4*A2_3*A3_2*A4_1 ```

## Non-Square matrix permanents

Really, the reason I wrote permanent was to answer a question, how to compute the permanent of a non-square matrix. Logically, the extension of a permanent to a NxM matrix, where N~=M is easy enough. (Assume that N < M here, since a matrix transpose would not impact the permanent. And of course, row or column exchanges do not impact the permanent, unlike a determinant.)

So, for N < M, the permanent is simply the sum of all sets of element products, where one element is taken from each of N rows, but no two elements in each product may be taken from the same row.

A simple conclusion is that for an Nx(N+1) matrix, the permanent is the sum of N+1 permanents of NxN matrices, just excluding each column one at a time.

A test of that claim is the difference should be zero here:

```A = sym('A',[3,4]) permanent(A) - (permanent(A(:,[1 2 3])) + permanent(A(:,[1 2 4])) + permanent(A(:,[1 3 4])) + permanent(A(:,[2 3 4]))) ```
```A = [ A1_1, A1_2, A1_3, A1_4] [ A2_1, A2_2, A2_3, A2_4] [ A3_1, A3_2, A3_3, A3_4] ans = 0 ```

## Another way of testing the code is to write the permanent recursively

using minors, like the standard recursive computation for a determinant. Here I'll expand it for a 4x2 symbolic array. The test for correctness is again a zero result:

```A = sym('A',[4,2]) P = permanent(A) - (A(1,1)*permanent(A([2 3 4],2)) + A(2,1)*permanent(A([1 3 4],2)) + A(3,1)*permanent(A([1 2 4],2)) + A(4,1)*permanent(A([1 2 3],2))) simplify(P) ```
```A = [ A1_1, A1_2] [ A2_1, A2_2] [ A3_1, A3_2] [ A4_1, A4_2] P = A1_1*A2_2 - A2_1*(A1_2 + A3_2 + A4_2) - A3_1*(A1_2 + A2_2 + A4_2) - A4_1*(A1_2 + A2_2 + A3_2) - A1_1*(A2_2 + A3_2 + A4_2) + A1_2*A2_1 + A1_1*A3_2 + A1_2*A3_1 + A1_1*A4_2 + A1_2*A4_1 + A2_1*A3_2 + A2_2*A3_1 + A2_1*A4_2 + A2_2*A4_1 + A3_1*A4_2 + A3_2*A4_1 ans = 0 ```

## The permanent of a vector is just the sum of the elements.

```A = sym('A',[1,20]) P = permanent(A) ```
```A = [ A1, A2, A3, A4, A5, A6, A7, A8, A9, A10, A11, A12, A13, A14, A15, A16, A17, A18, A19, A20] P = A1 + A2 + A3 + A4 + A5 + A6 + A7 + A8 + A9 + A10 + A11 + A12 + A13 + A14 + A15 + A16 + A17 + A18 + A19 + A20 ```