# Eigenvalues broken in Version 12.0

Bug introduced in 12.0.

The following code calculates the eigenvalues of a certain complex matrix, which come in pairs of opposite complex numbers. Therefore one can check whether the sum of all eigenvalues is equal to the trace of the matrix, which is zero.

This is indeed the case in Version 10.1 & 11.3 as far as I tested. However, Version 12.0 (Windows, Mac, Linux) gives something seriously wrong.

NN = 374; R = 0.05;
t1 = -1 + Cos[x] - I Sin[x] + I R; t1p = -1 + Cos[x] + I Sin[x] +
I R;
mat[x_] =
DiagonalMatrix[Table[If[EvenQ[n], t1, -1], {n, 0, 2 NN - 1 - 1}],
1] + DiagonalMatrix[
Table[If[EvenQ[n], t1p, -1], {n, 0, 2 NN - 1 - 1}], -1] +
DiagonalMatrix[Table[If[EvenQ[n], -1, 0], {n, 0, 2 NN - 1 - 3}],
3] + DiagonalMatrix[
Table[If[EvenQ[n], -1, 0], {n, 0, 2 NN - 1 - 3}], -3];
mat0 = mat[-0.2 \[Pi]];
Tr@mat0  (* 0. *)
Total@Eigenvalues@mat0  (* 0.394003 - 0.566499 I *)


I would rather switch back to 11.3 for a while. This looks really dangerous...

Original post of a more complex matrix with the same issue:

The code plots the real part of adding each pair. So the correct plot should be just zeros everywhere. This is the case in Version 10.1 & 11.3 as far as I tested (scattered numbers around $$10^{-14}$$ or so). However, Version 12.0 (Windows, Mac, Linux) gives something different as shown below.

NN = 200; R = 0.05;
xlist = Table[x, {x, -0.2 \[Pi], 0.2 \[Pi], 0.01}];
modl[n_] := 2*^-3 (Quotient[n, 2] - NN/2);
t1 = -1 + Cos[x] - I Sin[x] + I R; t1p = -1 + Cos[x] + I Sin[x] + I R;
t2a[n_] := -1 - modl[n]; t2b[n_] := -1 + modl[n];
mat[x_] =
DiagonalMatrix[
Table[If[EvenQ[n], t1, t2a[n]], {n, 0, 2 NN - 1 - 1}], 1] +
DiagonalMatrix[
Table[If[EvenQ[n], t1p, t2a[n]], {n, 0, 2 NN - 1 - 1}], -1] +
DiagonalMatrix[
Table[If[EvenQ[n], t2b[n], 0], {n, 0, 2 NN - 1 - 3}], 3] +
DiagonalMatrix[
Table[If[EvenQ[n], t2b[n], 0], {n, 0, 2 NN - 1 - 3}], -3];
list0 = Sort@Re@Eigenvalues[mat[xlist[]]];
list0p = Table[list0[[i]] + list0[[2 NN - i + 1]], {i, NN}];
ListPlot[Tooltip@list0p, PlotRange -> All] • bugs: "This tag is reserved for questions where the problem has been vetted by this community and the observed behavior is confirmed to be a bug. Please do not use this tag for new questions." – AccidentalFourierTransform May 18 '19 at 0:50
• A simpler test is comparing Tr[matrix] with Total[Eigenvalues[matrix]] which shows a large discrepancy. I think this is worth reporting to support. It happens on OSX as well. – Carl Woll May 18 '19 at 1:14
• @Carl Woll. Interesting: On my machine (Windows 10; version 12.0) Tr[matrix] gives 0. and Total[Eigenvalues[matrix]] gives -2.23821*10^-13 - 2.66454*10^-15 I. – ulvi May 18 '19 at 4:05
• Filing a report (apparently that has not been done to date). – Daniel Lichtblau Jun 10 '19 at 23:14
• @xiaohuamao To my understanding, the problem is with the MKL library. It has been fixed in the development version (which uses a newer version of MKL). – ilian Jul 29 '19 at 14:04

Not a solution but too big for a comment. There seems to be a catastrophic failure in Eigenvalues happening that is not due to the matrix being crazy. As a diagnostic, let's calculate the smallest (by absolute value) eigenvalue of the upper-left $$n\times n$$ part of the matrix

M = mat[xlist[]];


For odd $$n$$ the answer is zero, so let's only do this for even $$n$$. We do this in two ways

1. Calculate all eigenvalues and pick the one with the smallest absolute value:
    e1[n_?EvenQ] := M[[;; n, ;; n]] // Eigenvalues // Abs // Min

1. Calculate only the smallest eigenvalue (by absolute value) with the Arnoldi algorithm:
    e2[n_?EvenQ] := Eigenvalues[M[[;; n, ;; n]], 1,
Method -> {"Arnoldi", "Criteria" -> "Magnitude", "Shift" -> 0}] // First // Abs


Method (2) is very reliable, whereas method (1) breaks down for $$n=358$$ and above: Considering that the Arnoldi algorithm has no problems with this matrix, there seems to be something really strange going on in method 1.

$Version (* 12.0.0 for Mac OS X x86 (64-bit) (April 7, 2019) *)  • This persists for Eigensystem also. – CA Trevillian May 19 '19 at 3:38 I attempted a workaround, to see if Eigensystem had any issues also. It does. This is very unfortunate. (Will we have to wait for 12.1 for the fix (?!)) My code here: e3[n_?EvenQ] := Eigensystem[M[[;; n, ;; n]]][] // Abs // Min  Produces the following, which matches with @Roman shows: (Apologies the colors/styles don't match with the plot from @Roman !!) $Version
(* 12.0.0 for Microsoft Windows (64-bit) (April 6, 2019) *)

• The generated eigenvectors aren't even really eigenvectors of the matrix, so the error is not just in the eigenvalues. – Roman May 19 '19 at 10:32
• If I understand what you mean, this would be some fundamental (internal) issue with the eigensolver being used? Maybe some error in the linear decomposition being used? – CA Trevillian May 19 '19 at 17:48
• Yes it looks like a serious foul in the default "Direct" diagonalizer. This is very surprising for something that's been around for such a long time (coded in LAPACK). – Roman May 19 '19 at 19:35
• No need to fear, there are easy sanity checks, as @CarlWoll has already posted one (trace rule). A sanity check for the eigenvectors is to compute {eval, evec} = Eigensystem[M] and then go evec.Transpose[M].Inverse[evec] - DiagonalMatrix[eval] // Norm which should give zero (or something close to it) if the eigenvectors evec truly diagonalize the matrix M. – Roman May 19 '19 at 20:44
• This sort of breakage of previously working code now happens so frequently on major updates that I've decided to keep one or more old executables of Mathematica on hand. Note: test suites should be catching all of these! – Ralph Dratman May 22 '19 at 2:29

I am not qualified to be at this site because the last time I used eigenvectors was well over half a century ago. The word "stiff matrix" came back to me, so I increased the precision of the author's code by rounding the two real numbers to 50 places. It took forever to compute, but Mathematica solved the problem accurately. That is, R = N[5/100, 50]; and mat0 = mat[N[-2/10 [Pi], 50]]; I ended up with zero to 47 places.

Occasionally--when solving differential equations numerically--I came across stiff systems, so I checked for this condition prior to working with them (I forget the method I used).

Again, sorry for my layman's interjection.

Nick Bagley

• Charles, according to Wikipedia, "[...] the eigenvalue problem for all normal matrices is well-conditioned." So stiffness should not play a role here, unless a bad algorithm is being used. – Roman Jan 26 at 17:09