0

I'm trying to plot a function $E_{01}(p)=0.5^{2}(-8+6p)+\frac{1}{N}\displaystyle\sum_{\mathbf{k}}(\sqrt{A_{\mathbf{k}}^{2}-B_{\mathbf{k}}^{2}}-A_{\mathbf{k}})$ where $A_{\mathbf{k}}=0.5(8-6p(1-\frac{1}{3}(\cos{2k_x}+\cos{2k_y}+\cos{2k_z})))$ and $B_{\mathbf{k}}=8*0.5*\cos{k_x}\cos{k_y}\cos{k_z}$ (in code below, I have a=$k_x$, b=$k_y$, c= $k_c$).

Subscript[E, 01][p_] :=0.5^2*(-8 + 6p) +1/10^15Sum[Sqrt[(0.5*(8 -6p (1 - 1/3(Cos[2a]+Cos[2b]+Cos[2c]))))^2 - (80.5*Cos[a]Cos[b]Cos[c])^2] -0.5(8 - 6p (1 - 1/3*(Cos[2a]+Cos[2b]+Cos[2c]))), {a, -1Pi/2, 1Pi/2}, {b, -1Pi/2, 1Pi/2}, {c, -1Pi/2, 1*Pi/2}]
If i plot above code like: Subscript[h, e01] = Plot[Subscript[E, 01][p], {p, 0, 1}]

I get: enter image description here

Now, if I now move from Sum to integral I get problems. Moving from sum to integral is usually done like $\frac{1}{N}\displaystyle\sum_{\mathbf{k}}F(\mathbf{k}) \to \frac{d^{3}}{(2\pi)^3}\int_{-\frac{\pi}{d}}^{\frac{\pi}{d}}\int_{-\frac{\pi}{d}}^{\frac{\pi}{d}}\int_{-\frac{\pi}{d}}^{\frac{\pi}{d}}dk_x dk_y dk_z F(\mathbf{k})$. If we look at a 3D lattice, then N is the number of atoms and d is distance between atoms. I tried numerical integration like:

Subscript[E, I2][p_] := 
 0.5^2*(-8 + 6*p) + 
  10^(-30)*NIntegrate[
    1/(2*\[Pi])^3*
      Sqrt[(0.5*(8 - 
            6*p (1 - 1/3*(Cos[2*a]*Cos[2*b]*Cos[2*c]))))^2 - (8*0.5*
          Cos[a]*Cos[b]*Cos[c])^2] - 
     0.5*(8 - 6*p (1 - 1/3*(Cos[2*a]*Cos[2*b]*Cos[2*c]))), {a, -1*
      Pi/10^(-10), 1*Pi/10^(-10)}, {b, -1*Pi/10^(-10), 1*Pi/10^(-10)}, {c, -1*Pi/10^(-10), 1*Pi/10^(-10)}]

And now if I try to plot this like:

Plot[Subscript[E, I2][p], {p, 0, 1}, PlotRange -> {{0, 1}, {-3, 0}}]

Wolfram mathematica keeps running without stopping.

I tried using different methods like Trapezoidal, NewtonCotes, GaussKronrod and "SymbolicProcesessing" -> 0, but I have not observed any speeding up. I might have used them incorectly and that is way they did not really work (That's why I wrote bare minimum of code here). Is there a way to speed up NIntegrate?

Thanks!

Alen
  • 101
  • 1
  • 1
    Your integrand is periodic (with period 2π) in all integration variables. You can therefore integrate it in the cube $[0,2\pi]^3$ (which is fast with Method -> "LocalAdaptive") and there is no need for integrating over your entire integration volume (unless you care about edge effects). In other words, averaging over a huge volume in $k$-space will give the same result as integrating over a single unit cell. – Roman Sep 09 '22 at 09:54
  • 2
    Syntactically, (1) don't use E because it is Euler's number and (2) don't use subscripts. – Roman Sep 09 '22 at 10:02
  • 1
    If I understood you correctly, instead of of a,b,c that are between $-\frac{\pi}{10^{-10}}$ to $\frac{\pi}{10^{-10}}$ i should use interval [0,2$\pi$] for a,b and c, also I should add Method -> "LocalAdaptive"? – Alen Sep 09 '22 at 10:12
  • Yes. For an integrand $u(a)$ that is $2\pi$-periodic, and $0<d\ll1$, we have $\int_{-\pi/d}^{\pi/d}u(a)da \approx d^{-1}\int_{-\pi}^{\pi}u(a)da$. Do this in 3D and you're almost done. – Roman Sep 09 '22 at 10:29
  • I've added 1NIntegrate[equation, {a, -1Pi/1, 1Pi/1}, {b, -1Pi/1, 1Pi/1}, {c, -1Pi/1, 1*Pi/1}, Method -> "LocalAdaptive"]. Hopefully I did it correctly :D – Alen Sep 09 '22 at 10:35

0 Answers0