5
$\begingroup$

Let me explain the problem. I am trying to integrate a one dimensional integral: $$\int {g\left( {{k_x}},parameter1,parameter2,...\right)d{\mkern 1mu} {k_x}} $$

for the sake of clarity, I will give function g in the last.

I choose a specific set of parameter, and do the following integration:

In:=NIntegrate[g[kx, 9, -2, 7, 1, 1, 0.05], {kx, -\[Pi]/3,\[Pi]/3}]//AbsoluteTiming

Out:={422.396160, 0.163126 + 0.103155 I}

This takes me 422s on my computer!!!

I use the following code to show the sampling points used during the integration.

sampp = Reap[
 NIntegrate[#, {kx, -\[Pi]/(3), \[Pi]/(3)}, 
  EvaluationMonitor :> Sow[{kx, #}]]] &[
g[kx, 9, -2, 7, 1, 1, 0.05]] ;
{Length[#], ListPlot[#, AxesOrigin -> {0, 0}, Filling -> Axis]} &[
 relist@sampp[[2, 2, 1]]]

sampling points

it shows that it used 733 sample points. But should that take 400 seconds?

Actually, the evaluation of the integrand at each sample point is acually fast. I use the most simple interpolation method to do the same integration as below, I use 1000 sample points that is evenly distributed and it only takes 1.3 second. Besides the integration result is quite accurate.

In:=(sampp = Table[{kx, 
 g[kx, 9, -2, 7, 1, 1, 0.05]}, {kx, -\[Pi]/3, \[Pi]/3, 
 2 \[Pi]/3/1000}];
Integrate[
Interpolation[sampp, InterpolationOrder -> 2][
x], {x, -\[Pi]/3, \[Pi]/3}]) // AbsoluteTiming
Out:={1.308075, 0.163126 + 0.103155 I}

So what did Mathematica do in the extremely long 400 second? Of course, it needs to do many error estimation, but that should not take too much time.

I scan the documents and related questions in stackexchange. Somebody suggest to set "SymbolicProcessing"->0. OK, I tried

In:=NIntegrate[g[kx, 9, -2, 7, 1, 1, 0.05], {kx, -\[Pi]/3,     \[Pi]/3}, 
Method -> {"GlobalAdaptive", Method -> "GaussKronrodRule", 
"SymbolicProcessing" -> 0}]//AbsoluteTiming
Out:={401.791981, 0.163126 + 0.103155 I}

unfortunately, the result is the same.

Even more peculiar, when I use "Trapezoidal" strategy and limit the MaxRecursion

 In:=NIntegrate[g[kx, 9, -2, 7, 1, 1, 0.05], {kx, -\[Pi]/3, \[Pi]/3}, 
 Method -> "Trapezoidal", MinRecursion -> 1, MaxRecursion -> 4]
Out:={413.894674, 0.163126 + 0.103155 I}

the sample points shows as follows sample points2 it only used 257 sample point, but the time cost is the same??!!It is unbelievable! Anyway, I think the "Trapezoidal" should be the same as evenly sample interpolation method.

I really can't understand what did Mathematica do in the NIntegrate of function g. Can somebody explain it?

Of course, it seems that I can use the simple interpolation integration method manually, because it is fast, and the result seems good in this case. But it is not robust, I insisted using the built in "Adaptive" strategy, because I have to Integrate some other functions which are sharp at specific points.


In the last the form of function g

g[kx_, xx_, e_, width_, ii_, label_, \[Eta]_] := 
Inverse[(e - I \[Eta]) IdentityMatrix[2*width] - 
 armchairibbonmat[kx, width]][[2*(ii - 1) + label, 
2*(ii - 1) + label]] E^(I kx xx)

(*armchairibbonmat is a function used in function g*)

armchairibbonmat[kx_, n_] :=
(
h = Table[0, {i, 1, 2 (n + 2)}, {j, 1, 2 (n + 2)}];
a[m_] := 2 m - 1;
b[m_] := 2 m;
t1 = 1;
aa = 1;
Do[
h[[b[i], a[i]]] = t1 E^(I kx aa);
h[[a[i], b[i]]] = t1 E^(-I kx aa);
h[[b[i + 1], a[i]]] = t1 E^(-I kx aa/2);
h[[b[i - 1], a[i]]] = t1 E^(-I kx aa/2);
h[[a[i - 1], b[i]]] = t1 E^(I kx aa/2);
h[[a[i + 1], b[i]]] = t1 E^(I kx aa/2);
, {i, 2, n + 1}];
h = ArrayPad[h, -2]
)
$\endgroup$
1
  • $\begingroup$ Sorry, everyone, I made a few typoerrors in the definition of function g at first, but now it is fixed. thank you. $\endgroup$
    – matheorem
    Aug 19, 2013 at 8:26

1 Answer 1

10
$\begingroup$

If you replace

g[kx_,

with

g[kx_?NumericQ

Things will work as you expect.

Update: Evaluation of

g[kx, 9, -2, 7, 1, 1, 0.05]

takes a long time and the _?NumericQ absolutely prevents the symbolic evaluation.

$\endgroup$
6
  • 1
    $\begingroup$ I like short answers that get to the point! $\endgroup$
    – Yves Klett
    Aug 19, 2013 at 9:24
  • $\begingroup$ @ruebenko it works!! Thank you!! If I would know the answer is so simple, I would not have written such a long detailed question, HaHa $\endgroup$
    – matheorem
    Aug 19, 2013 at 9:59
  • 1
    $\begingroup$ @ruebenko But why this fixed the problem? Can you explain a little? $\endgroup$
    – matheorem
    Aug 19, 2013 at 10:01
  • $\begingroup$ See update, there many similar examples in this forum, have a look. $\endgroup$
    – user21
    Aug 19, 2013 at 10:36
  • $\begingroup$ @ruebenko but why "SymbolicProcessing" -> 0 didn't work? $\endgroup$
    – matheorem
    Aug 23, 2013 at 9:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.