# Find all degenerate eigenvalues of a cubic equation

I have an equation that is cubic in w. The three solutions correspond to bands in a bandstructure, and are a function of wavevector {x, y}. Additionally it is parameterised by 0 <= θ <= π and 0 <= φ <= 2 π. I need a procedure that takes θ and φ as inputs, and can reliably calculate if there are any points at all in wavevector space where either two or all of the bands touch (i.e. are there any values of {x, y} which give rise to a degenerate solution w). Ideally the procedure would then be able to calculate the coordinates in momentum space of these touching points. The code to generate the cubic equation is given below.

M = {{s, ab, ac}, {ab, s, bc}, {ac, bc, s}};
f[e_] := 2 Norm[e]^-3 (1 - 3 Sin[θ]^2 Cos[φ - ArcTan[e[[1]], e[[2]]]]^2)Cos[{x, y}.e];
ab = f[{1, 0}];
ac = f[{0, 1}];
bc = f[{1, 1}] + f[{-1, 1}];
s = f[{2, 0}] + f[{0, 2}];


A slight variation on this question has already been solved here by @bbgodfrey and others and can serve as a more detailed reference. However that question called for a solution when all three bands touch, as opposed to two or more, which is the current question.

One of my ultimate goals is to be able to generate a plot in {θ, φ} space showing what values give rise to a bandstructure with touching points. As an example, in previous work a procedure for a different physical system generated the plot below, where purple is a gapless bandstructure and cream is a gapped bandstructure:

• @bbgodfrey thanks for suggesting I create a new question. Here it is.
– Tom
Jul 14, 2015 at 10:00
• Any reason not to use Eigenvalues[] here? Jul 14, 2015 at 10:08
• It always spits out something in the form of Root[] which I didn't really know how to handle, so I did it this way.
– Tom
Jul 14, 2015 at 10:44
• That's its default behavior; if you want radicals in your output, use ToRadicals[]. Jul 14, 2015 at 10:46
• Thanks. Have modified question to incorporate this.
– Tom
Jul 14, 2015 at 10:55

For a numerical approximation you may try something like:

M = {{s, ab, ac}, {ab, s, bc}, {ac, bc, s}};
disc = Discriminant[CharacteristicPolynomial[M, x], x] // FullSimplify
f[e_] := 2 Norm[e]^-3 (1 - 3 Sin[θ]^2 Cos[φ - ArcTan[e[[1]], e[[2]]]]^2) Cos[{x, y}.e];
ab = f[{1, 0}];
ac = f[{0, 1}];
bc = f[{1, 1}] + f[{-1, 1}];
s = f[{2, 0}] + f[{0, 2}];

Quiet@RegionPlot[With[{θ = t, φ = p},
FindMinValue[disc, {x, y}, AccuracyGoal -> Infinity]] < 10^-6,
{t, 0, Pi}, {p, 0, 2 Pi}, PlotPoints -> 60]


You should take care of the numeric subtleties. Please note:

Quiet@RegionPlot[
FindMinValue[disc /. {θ -> t, φ -> p}, {{x, 0, Pi}, {y, 0, Pi}}] < .001,
{t, 0, Pi}, {p, 0, 2 Pi}, PlotPoints -> 30]


• Oooo interesting thank you very much. I'll have a play with this.
– Tom
Jul 14, 2015 at 20:39