# Eigenvalues of a Matrix

Given a square matrix, output the matrix's eigenvalues. Each eigenvalue should be repeated a number of times equal to its algebraic multiplicity.

The eigenvalues of a matrix A are scalar values λ such that, for some column vector v, A*v = λ*v. They are also the solutions to the characteristic polynomial of A: det(A - λ*I) = 0 (where I is the identity matrix with the same dimensions as A).

Outputs must be accurate to 3 significant digits. All inputs and outputs will be within the representable range of numerical values for your chosen language.

Builtins are acceptable, but you are encouraged to include solutions that do not use builtins.

## Test Cases

In these test cases, I represents the imaginary unit. Complex numbers are written in the form a + b*I. All outputs have 3 significant digits of precision.

[[42.0]] -> [42.0]
[[1.0, 0.0], [0.0, 1.0]] -> [1.00, 1.00]
[[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]] -> [16.1, -1.12, -1.24e-15]
[[1.2, 3.4, 5.6, 7.8], [6.3, 0.9, -5.4, -2.3], [-12.0, -9.7, 7.3, 5.9], [-2.5, 7.9, 5.3, 4.4]] -> [7.20 + 5.54*I, 7.20 - 5.54*I, -4.35, 3.75]
[[-3.22 - 9.07*I, 0.193 + 9.11*I, 5.59 + 1.33*I, -3.0 - 6.51*I, -3.73 - 6.42*I], [8.49 - 3.46*I, -1.12 + 6.39*I, -8.25 - 0.455*I, 9.37 - 6.43*I, -6.82 + 8.34*I], [-5.26 + 8.07*I, -6.68 + 3.72*I, -3.21 - 5.63*I, 9.31 + 3.86*I, 4.11 - 8.82*I], [-1.24 + 9.04*I, 8.87 - 0.0352*I, 8.35 + 4.5*I, -9.62 - 2.21*I, 1.76 - 5.72*I], [7.0 - 4.79*I, 9.3 - 2.31*I, -2.41 - 7.3*I, -7.77 - 6.85*I, -9.32 + 2.71*I]] -> [5.18 + 16.7*I, -24.9 - 2.01*I, -5.59 - 13.8*I, 0.0438 - 10.6*I, -1.26 + 1.82*I]
[[-30.6 - 73.3*I, 1.03 - 15.6*I, -83.4 + 72.5*I, 24.1 + 69.6*I, 52.3 + 2.68*I, 23.8 + 98.0*I, 96.8 + 49.7*I, -26.2 - 5.87*I, -52.4 + 98.2*I, 78.1 + 6.69*I], [-59.7 - 66.9*I, -26.3 + 65.0*I, 5.71 + 4.75*I, 91.9 + 82.5*I, -94.6 + 51.8*I, 61.7 + 82.3*I, 54.8 - 27.8*I, 45.7 + 59.2*I, -28.3 + 78.1*I, -59.9 - 54.5*I], [-36.0 + 22.9*I, -51.7 + 10.8*I, -46.6 - 88.0*I, -52.8 - 32.0*I, -75.7 - 23.4*I, 96.2 - 71.2*I, -15.3 - 32.7*I, 26.9 + 6.31*I, -59.2 + 25.8*I, -0.836 - 98.3*I], [-65.2 - 90.6*I, 65.6 - 24.1*I, 72.5 + 33.9*I, 1.47 - 93.8*I, -0.143 + 39.0*I, -3.71 - 30.1*I, 60.1 - 42.4*I, 55.6 + 5.65*I, 48.2 - 53.0*I, -3.9 - 33.0*I], [7.04 + 0.0326*I, -12.8 - 50.4*I, 70.1 - 30.3*I, 42.7 - 76.3*I, -3.24 - 64.1*I, 97.3 + 66.8*I, -11.0 + 16.5*I, -40.6 - 90.7*I, 71.5 - 26.2*I, 83.1 - 49.4*I], [-59.5 + 8.08*I, 74.6 + 29.1*I, -65.8 + 26.3*I, -76.7 - 83.2*I, 26.2 + 99.0*I, -54.8 + 33.3*I, 2.79 - 16.6*I, -85.2 - 3.64*I, 98.4 - 12.4*I, -27.6 - 62.3*I], [82.6 - 95.3*I, 55.8 - 73.6*I, -49.9 + 42.1*I, 53.4 + 16.5*I, 80.2 - 43.6*I, -43.3 - 3.9*I, -2.26 - 58.3*I, -19.9 + 98.1*I, 47.2 + 62.4*I, -63.3 - 54.0*I], [-88.7 + 57.7*I, 55.6 + 70.9*I, 84.1 - 52.8*I, 71.3 - 29.8*I, -3.74 - 19.6*I, 29.7 + 1.18*I, -70.6 - 10.5*I, 37.6 + 99.9*I, 87.0 + 19.0*I, -26.1 - 82.0*I], [69.5 - 47.1*I, 11.3 - 59.0*I, -84.3 - 35.1*I, -3.61 - 35.7*I, 88.0 + 88.1*I, -47.5 + 0.956*I, 14.1 + 89.8*I, 51.3 + 0.14*I, -78.5 - 66.5*I, 2.12 - 53.2*I], [0.599 - 71.2*I, 21.7 + 10.8*I, 19.9 - 97.1*I, 20.5 + 37.4*I, 24.7 + 40.6*I, -82.7 - 29.1*I, 77.9 + 12.5*I, 94.1 - 87.4*I, 78.6 - 89.6*I, 82.6 - 69.6*I]] -> [262. - 180.*I, 179. + 117.*I, 10.3 + 214.*I, 102. - 145.*I, -36.5 + 97.7*I, -82.2 + 89.8*I, -241. - 104.*I, -119. - 26.0*I, -140. - 218.*I, -56.0 - 160.*I]

• Sandbox
– user45941
Nov 27, 2017 at 8:59
• Nov 27, 2017 at 9:09
• Related? Probably related? (depends on the approach) Nov 27, 2017 at 12:20

# Haskell, 576 554 532 507 bytes

No built-ins!

import Data.Complex
s=sum
l=length
m=magnitude
i=fromIntegral
(&)=zip
t=zipWith
(x!a)b=x*a+b
a#b=[[s$t(*)x y|y<-foldr(t(:))([]<$b)b]|x<-a]
f a|let c=[1..l a];g(u,d)k|m<-[t(+)a b|(a,b)<-a#u&[[s[d|x==y]|y<-c]|x<-c]]=(m,-s[s[b|(n,b)<-c&a,n==m]|(a,m)<-a#m&c]/i k)=snd<$>scanl g(0<$c<$c,1)c p?x|let f=foldl1(x!);c=l p-1;n=i c;q p=init$t(*)p$i<$>[c,c-1..];o=f(q p)/f p;a|d<-sqrt$(n-1)*(n*(o^2-f(q$q p)/f p)-o^2)=n/last(o-d:[o+d|m(o-d)<m(o+d)])=last$p?(x-a):[x|m a<1e-9] z[a,b]=[-b/a] z p=p?0:z(init$scanl1(p?0!)p)


Try it online!

Many thanks @ØrjanJohansen for a total of -47 bytes!

### Explanation

First this computes the characteristic polynomial with the Faddeev–LeVerrier algorithm which is the function f. Then the function z computes all roots of that polynomial by iterating g which implements Laguerre's Method for finding a root, once a root is found it is removed and g gets called again until the polynomial has degree 1 which is trivially solved by z[a,b]=[-b/a].

### Ungolfed

I re-inlined the functions sum,length,magnitude, fromIntegral, zipWith and (&) as well as the little helper (!). The function faddeevLeVerrier corresponds to f, roots to z and g to laguerre respectively.

-- Transpose a matrix/list
transpose a = foldr (zipWith(:)) (replicate (length a) []) a

-- Straight forward implementation for matrix-matrix multiplication
(#) :: [[Complex Double]] -> [[Complex Double]] -> [[Complex Double]]
a # b = [[sum $zipWith (*) x y | y <- transpose b]|x<-a] -- Faddeev-LeVerrier algorithm faddeevLeVerrier :: [[Complex Double]] -> [Complex Double] faddeevLeVerrier a = snd <$> scanl go (zero,1) [1..n]
where n = length a
zero = replicate n (replicate n 0)
trace m = sum [sum [b|(n,b)<-zip [1..n] a,n==m]|(m,a)<-zip [1..n] m]
diag d = [[sum[d|x==y]|y<-[1..n]]|x<-[1..n]]
add as bs = [[x+y | (x,y) <- zip a b] | (b,a) <- zip as bs]
go (u,d) k = (m, -trace (a#m) / fromIntegral k)
where m = add (diag d) (a#u)

-- Compute roots by succesively removing newly computed roots
roots :: [Complex Double] -> [Complex Double]
roots [a,b] = [-b/a]
roots   p   = root : roots (removeRoot p)
where root = laguerre p 0
removeRoot = init . scanl1 (\a b -> root*a + b)

-- Compute a root of a polynomial p with an initial guess x
laguerre :: [Complex Double] -> Complex Double -> Complex Double
laguerre p x = if magnitude a < 1e-9 then x else laguerre p new_x
where evaluate = foldl1 (\a b -> x*a+b)
order' = length p - 1
order  = fromIntegral $length p - 1 derivative p = init$ zipWith (*) p $map fromIntegral [order',order'-1..] g = evaluate (derivative p) / evaluate p h = (g ** 2 - evaluate (derivative (derivative p)) / evaluate p) d = sqrt$ (order-1) * (order*h - g**2)
ga = g - d
gb = g + d
s = if magnitude ga < magnitude gb then gb else ga
a = order /s
new_x = x - a

• As the only submission that doesn't use builtins, this should be the highest-voted answer. Dec 15, 2017 at 5:46
• +1 for calculating something-related-to-determinant in time less than n!! Dec 15, 2017 at 7:18
• Thanks guys! @user202729: Initially I oversaw the ! and was really confused :D Dec 15, 2017 at 17:34

# Octave, 4 bytes

@eig


Try it online!

Only two bytes more than the MATL golfing language equivalent!

Defines an anonymous function handle to the eig built-in. Interestingly, the MATLAB design philosophy goes against many high-end languages, which like to use DescripteFunctionNamesTakingArguments(), whereas MATLAB and consequently Octave tend to get the shortest unambiguous function name possible. For example, to get a subset of eigenvalues (e.g., the smallest n in absolute magnitude), you use eigs.

As a bonus, here's a function (works in MATLAB, and in theory could work in Octave but their solve is not really up to the task) that doesn't use built-ins, but instead symbolically solves the eigenvalue problem, det(A-λI)=0, and converts it to numerical form using vpa

@(A)vpa(solve(det(A-sym('l')*eye(size(A)))))


# MATL, 2 bytes

Yv


Try it online!

### Explanation

I followed the usual advice in numerical linear algebra: instead of writing your own function, use a built-in that is specifically designed to avoid numerical instabilities.

Incidentally, it's shorter. ¯\_(ツ)_/¯

• This begs the question, how long would it be without Yv? Nov 27, 2017 at 10:33
• @Sanchises I'm not sure. I'd probably go about it by finding the roots (ZQ) of the characteristic polynomial. But explicitly computing the coeffients of the polynomial may be a lot of work Nov 27, 2017 at 11:56

# Mathematica, 11 bytes

Eigenvalues


Try it online!

• Yes, I expected a built-in answer before clicking "1 new answer to this question". Let's wait for some non-builtin answer... / Basically Mathematica solution is often <first word in the title> Nov 27, 2017 at 9:20
• The shortest non-pure built-in I got is First@Eigensystem@#& (20 bytes) Nov 27, 2017 at 13:19
• I actually agree with user202729 here. While it is fun to joke about Mathematica having a builtin for everything, it is very annoying both as a challenge poster and answering to see a builtin answer equivalent to something you tried relatively hard to do. Golfing (IMO) is about attempting to find the shortest algorithm and implementation of said algorithm, but a builtin answer takes that away from the "sport". Nov 27, 2017 at 16:25
• @cairdcoinheringaahing We should really start putting xnor's proposal into practice. Nov 29, 2017 at 14:03

# R, 22 bytes

function(m)eigen(m)\$va


Try it online!

Takes m as a matrix. Frustratingly, the eigen function in R returns an object of class eigen, which has two fields: values, the eigenvalues, and vectors, the eigenvectors.

However, more annoyingly, the optional argument only.values returns a list with two fields, values, containing the eigenvalues, and vectors, set to NULL, but since eigen(m,,T) is also 22 bytes, it's a wash.

# Python + numpy, 33 bytes

from numpy.linalg import*
eigvals


# Pari/GP, 19 bytes

m->mateigen(m,1)[1]


Try it online!

# Pari/GP, 24 bytes

It seems that someone downvoted every answer that uses a built-in. OK, if you don't like built-in:

m->polroots(matdet(x-m))


Try it online!

# Julia 0.6, 12 bytes

n->eig(n)[1]


Try it online!

Unfortunately, eig returns both the eigenvalues and eigenvectors, as a tuple, so we waste another 9 bytes to lambdify it and grab the first item.

• or eigvals for 7 bytes Dec 9, 2021 at 11:18