Compute the determinant
Challenge
A simple challenge: Given a two-dimensional matrix (an array of arrays) of real numbers, compute the determinant.
The determinant of a matrix is a mathematical construct used in many applications, such as solving polynomial equations, identifying the invertibility of the matrix, and finding the scaling factor under a matrix transformation. For more information about it, see this Wikipedia entry.
There are a couple of different ways to compute the determinant, and it is up to you how you implement it.
For instance, you may compute it using the Laplace expansion, a recursive algorithm which goes
- Pick a row or column.
- Start with a sum of zero.
- For each entry of the row/column:
- Create a new matrix with the row and column of the entry removed. This new matrix is a square matrix of size one less that the original.
- Compute the determinant of that smaller matrix.
- Multiply that determinant by the entry.
- If the row index plus the column index is even[1], add it to the sum, otherwise, subtract it.
- The final sum is the determinant.
As an example, here is an ungolfed implementation along the first column.
function laplaceDet(matrix) {
if (matrix.length === 1) return matrix[0][0];
let sum = 0;
for (let rowIndex = 0; rowIndex < matrix.length; ++rowIndex) {
let minorMatrix = matrix.filter((_, index) => index !== rowIndex)
.map(row => row.slice(1));
sum += ((-1) ** rowIndex) * matrix[rowIndex][0] * laplaceDet(minorMatrix);
}
return sum;
}
This is code-golf, so the program with the lowest byte-count wins!
Test cases
$$ \begin{aligned} \det\begin{bmatrix}1&0\\0&1\end{bmatrix}&=1 \\ \det\begin{bmatrix}1.5&2\\-3&4.5\end{bmatrix}&=12.75 \\ \det\begin{bmatrix}3&7\\1&-4\end{bmatrix}&=-19 \\ \det\begin{bmatrix}1&0&0\\0&1&0\\0&0&1\end{bmatrix}&=1 \\ \det\begin{bmatrix}1&2&3\\4&5&6\\7&8&9\end{bmatrix}&=0 \end{aligned} $$
In text form,
[[1,0],[0,1]] -> 1
[[1.5,2],[-3,4.5]] -> 12.75
[[3,7],[1,-4]] -> -19
[[1,0,0],[0,1,0],[0,0,1]] -> 1
[[1,2,3],[4,5,6],[7,8,9]] -> 0
-
Note that it doesn't matter if it is 1-indexed or 0-indexed, as 1+1 and 0+0 are both the same parity. ↩︎
CJam, 45 bytes {:A, …
3y ago
[Python 3.8 (pre-release)], 10 …
3y ago
[Haskell], 75 73 bytes -2 b …
3y ago
C, 147 bytes ```c float d( …
3y ago
Scala, 130 125 bytes Saved …
3y ago
J, 5 bytes ``` -/ . ``` …
3y ago
[Python 3], 29 bytes …
3y ago
[Wolfram Language (Mathematica …
1y ago
Ruby, 35 bytes ->m{requi …
3y ago
9 answers
CJam, 45 bytes
{:A_,{1$_,,.=:+\)/:CAff*A@zf{\f.*::+}..-}/;C}
This implementation is an anonymous block (~function). Online test suite.
Dissection
This implements the Faddeev-LeVerrier algorithm.
The objective is to calculate the coefficients $c_k$ of the characteristic polynomial of the $n\times n$ matrix $A$, $$\begin{eqnarray}p(\lambda)\equiv \det(\lambda I_{n}-A) = \sum_{k=0}^{n}c_{k}\lambda^{k}\end{eqnarray}$$ where, evidently, $c_n = 1$ and $c_0 = (-1)^n \det A$.
The coefficients are determined recursively from the top down, by dint of the auxiliary matrices $M$, $$\begin{aligned}M_{0}&\equiv 0&c_{n}&=1\qquad &(k=0)\\M_{k}&\equiv AM_{k-1}+c_{n-k+1}I\qquad \qquad &c_{n-k}&=-{\frac {1}{k}}\mathrm {tr} (AM_{k})\qquad &k=1,\ldots ,n~.\end{aligned}$$
The code never works directly with $c_{n-k}$ and $M_k$, but always with $(-1)^k c_{n-k}$ and $(-1)^{k+1}AM_k$, so the recurrence is $$\begin{eqnarray*}(-1)^k c_{n-k} &=& \frac1k \mathrm{tr} ((-1)^{k+1} AM_{k}) \\ (-1)^{k+2} AM_{k+1} &=& (-1)^k c_{n-k}A - A((-1)^{k+1}A M_k)\end{eqnarray*}$$
{ e# Define a block
:A e# Store the input matrix in A
_, e# Take the length of a copy
{ e# for i = 0 to n-1
e# Stack: (-1)^{i+2}AM_{i+1} i
1$_,,.=:+ e# Calculate tr((-1)^{i+2}AM_{i+1})
\)/:C e# Divide by (i+1) and store in C
Aff* e# Multiply by A
A@ e# Push a copy of A, bring (-1)^{i+2}AM_{i+1} to the top
zf{\f.*::+} e# Matrix multiplication
..- e# Matrix subtraction
}/
; e# Pop (-1)^{n+2}AM_{n+1} (which incidentally is 0)
C e# Fetch the last stored value of C
}
0 comment threads
Python 3.8 (pre-release), 106 95 bytes
Saved 11 bytes thanks to Peter Taylor!
f=lambda m:sum((-1)**i*x*f([r[:i]+r[i+1:]for r in m[1:]])for i,x in enumerate(m[0]))if m else 1
Here's an attempt at beating Quintec's answer. Only 78 67 more bytes to go!
Haskell, 75 73 bytes
-2 bytes thanks to @user
f[]=1
f(x:y)=foldr(-)0$zipWith(\i->(*f[take i r++drop(i+1)r|r<-y]))[0..]x
C, 147 bytes
float d(float**m,int r){if(r<2)return**m;int i=r,j;float s=0,*n[--r];for(;i--;s+=(i%2?-1:1)**m[i]*d(n,r))for(j=r;j--;)n[j]=m[j+(j>=i)]+1;return s;}
Basically a C version of my example code. Takes an array of pointers and the size.
Prettified and commented version:
float d(float **m, int r) {
if(r < 2) return **m; // Base case: 1x1 matrix
int i = r, j; // Variable initialization, i and j are counters
float s = 0, *n[--r]; // s is the sum, n is the minor matrix
for(; i--; s += (i % 2 ? -1 : 1) * *m[i] * d(n, r)) // Outer loop: loop over the rows of M.
// After each iteration, add the term to the sum
for(j = r; j--; ) // Inner loop to fill in the minor matrix
n[j] = m[j + (j >= i)] + 1; // (j >= i) term skips over the ith row of M.
// Addition of 1 effectively removes the first element
return s;
}
0 comment threads
Scala, 130 125 bytes
Saved 5 bytes by returning 1 like Hakerh400's great answer
def f(m:Seq[Seq[Double]]):Double=if(m.size<1)1 else(m.indices:\.0){(i,a)=>m(0)(i)*f(m.tail.map(r=>r.take(i)++r.drop(i+1)))-a}
A recursive function that uses the first row for the Laplace expansion. Explanation coming soon.
def f(m: Seq[Seq[Double]]): Double =
if (m.size < 1) 1 //If it's empty, return 1
else
//Otherwise, fold right over [0..m.size-1]
//The initial accumulator is 0.0
(m.indices :\ .0) {
(i, a) =>
//Multiply the entry by
m(0)(i) *
//The determinant of the smaller matrix
f(
//Drop the first row
m.tail.map(
//And for each row r, drop column i (0-indexed)
r => r.take(i) ++ r.drop(i + 1)
)
)
//Subtract the accumulator from that to invert
//its sign each time
- a
}
0 comment threads
Python 3, 29 bytes
import numpy
numpy.linalg.det
Yeah, this is really boring, but I'm curious if this can be beat in Python.
1 comment thread