# Fast Plot3D, failing NIntegrate, and reckless surgery

I have a square matrix, m which depends on kx and ky. It isn't Hermitian, but it does have real eigenvalues. I would like to obtain the integral of the sum of the eigenvalues of this matrix over a region in the space of those two variables. After several minutes of evaluation I begin to receive notifications that not all eigenvalues are being found.

Message[Eigenvalues::eival]


Presumably this is for only a few points in the region which are probably of vanishing measure.

In contrast to the trouble I have with NIntegrate, when I plot the sum of eigenvalues, Plot3D returns a highly-resolved and smooth plot in a very modest amount of time. Here's my code:

omega[k_?VectorQ] := Abs[Chop[Eigenvalues[m /. {kx -> k.{1, 0}, ky -> k.{0, 1}}]]]

Plot3D[Total[omega[{kx, ky}]], {kx, ky} \[Element] region,
AspectRatio -> Automatic,
PlotPoints -> 60, Mesh -> All]

NIntegrate[Total[omega[{kx, ky}]], {kx, ky} \[Element] region, PrecisionGoal -> 1]


I thought I might just use the data in the plot to coarsely approximate the integral myself. So I set MaxRecursion->0 and took the average of the points from the plot.

plot = Plot3D[Total[omega[{kx, ky}]], {kx, ky} \[Element] region, AspectRatio -> Automatic, PlotPoints -> 60, Mesh -> All]

dataPoints = Flatten[Cases[plot, x_GraphicsComplex :> First@x, Infinity], 1];
zData = dataPoints[[;; , 3]];
Mean[zData]*Area[region]


Obviously this is the dumbest way to do the integral, but for some matrices, m, it converges to three decimal places with PlotPoints->90 (and with far fewer points in other cases). I think that's evidence that my integrand isn't too pathological.

Shouldn't integration of a smooth function take about as much time as plotting a smooth function? Especially if I can get a sensible integral from such reckless surgery on Plot3D.

So what gives, NIntegrate???

Thanks!

• Please provide also the definition of m. – Henrik Schumacher Jan 28 at 22:28
• Thanks for taking a look, @HenrikSchumacher, the minimal working examples which I could easily construct don't seem to have this problem. I can usually obtain the integral fairly easily (although it is still a factor of 100 slower than a very nice Plot3D. The problematic m is actually a 36x36 matrix. It's relatively sparse and the entries tend to look like (1. I) E^(1/2 I (3 kx - Sqrt[3] ky)). I'd prefer not to post the matrix anyway, for privacy's sake, but I'll continue to try to find a minimum example. – thejacksonjive Jan 28 at 23:22