fast sum computation

Question

I would like to compute sum. How is it possible to compute the sum fast? May be with the help of replacing Sum with NSum or NIntegrate?

Here I define the functions

M = 100; 
delt[s_] := 1/2 Sin[s];
U[s_] := 1/2 Cos[s];
p1[s_, k_] := -1 - U[s] - (1 - U[s]) Cos[k];
p2[s_, k_] := (1 - U[s]) Sin[k];
p3[s_] := 1/2 Sin[s];
p[s_, k_] := 
 Sqrt[p1[s, k]^2 + p2[s, k]^2 + p3[s]^2]; 
sp1[s_, k_] := 
 0.5 Sin[s] (1 - Cos[k]);
sp2[s_, k_] := 0.5 Sin[s] Sin[k];
sp3[s_] := 0.5 Cos[s];
spc[s_, k_] := sp1[s, k] + I sp2[s, k];
kp1[s_, k_] := (1 - U[s]) Sin[k];
kp2[s_, k_] := (1 - U[s]) Cos[k];
kpc[s_, k_] := kp1[s, k] + I kp2[s, k];
a1[s_, k_] := 
 Sqrt[1/2 + 
   1/2 p3[s]/
     p[s, k]]; 
b1[s_, k_] := (p1[s, k] + I p2[s, k])/
  Sqrt[2 p[s, k]^2 + 2 p3[s] p[s, k]];
a0[s_, k_] := 
 Sqrt[1/2 - 
   1/2 p3[s]/
     p[s, k]]; 
b0[s_, k_] := -(p1[s, k] + I p2[s, k])/
  Sqrt[2 p[s, k]^2 - 2 p3[s] p[s, k]];
matr[s_, k_] := -1/
    2 ((a1[s, k] sp3[s] + b1[s, k]\[Conjugate] spc[s, k]) a0[s, 
        k] + (a1[s, k] spc[s, k]\[Conjugate] - 
         b1[s, k]\[Conjugate] sp3[s]) b0[s, k])/p[s, k];
berry[s_, k_] := 
  a1[s, k] Derivative[1, 0][a1][s, k] + 
   b1[s, k]\[Conjugate] Derivative[1, 0][b1][s, k] - 
   a0[s, k] Derivative[1, 0][a0][s, k] - 
   b0[s, k]\[Conjugate] Derivative[1, 0][b0][s, k];
ds = 2 Pi/30; ds1 = 2 Pi/50 ; 
dyn[k_] := 
 Sum[(2 p[s, k] ) ds1, {s, 0, 2 Pi, 
   ds1 }];  
dynb[k_] := Sum[berry[s, k] ds1, {s, 0, 2 Pi, ds1 }];
(*a00[s_,k_]:=I Sum[1/2Abs[matr[sk,k]]^2/p[sk,k] ds,{sk,0,s, ds }];*)
a00[s_, k_] := 0;   
a10[s_, k_] := 1/2 I matr[s, k]/p[s, k];
a11[k_] := -1/2 I matr[0, k]/p[0, k];
(*a11[k_]:=0;*)
a[s_, k_] := (1 + G a00[s, k]) a0[s, k] + 
  G (a10[s, k] + Exp[-I dyn[k]/G - dynb[k]] a11[k]) a1[s, 
    k];(*\[Psi]=(a b)^T approximate wave function*)
b[s_, k_] := (1 + G a00[s, k]) b0[s, k] + 
          G (a10[s, k] + Exp[-I dyn[k]/G - dynb[k]] a11[k]) b1[s, k];

And here I evaluate the sum

G = 0.01; 
Sum[-1/G 1/
   M  ds (a[s, k] b[s, k]\[Conjugate] kpc[s, k] + 
    a[s, k]\[Conjugate] b[s, k] kpc[s, k]\[Conjugate])/(Abs[
      a[s, k]]^2 + Abs[b[s, k]]^2)
 , {k, -Pi + 2 Pi/M, Pi, 2 Pi/M}, {s, 0, 2 Pi, ds}] // Timing

The sum over k must be a discreate, but the integration over s I replaced by the sum (the same as in the definition of dyn[k_])

Are you interested in a numeric or symbolic result? What is the value of G? — Thies Heidecke, Dec 10 '18 at 14:52
I would like to compute the plot graph for the sum versus G. That is why i need to compute the sum fast. But for this case G = 0.01; — Chipa-Chipa, Dec 10 '18 at 14:54

Henrik Schumacher · Answer 1 · 2018-12-10T16:03:44.190

5

Here a method with Compile. The actual computation of 1000 of these sums takes about three seconds.

cf = Block[{G, s, k}, 
   With[{code = -1./G 1/
        M ds (a[s, k] b[s, k]\[Conjugate] kpc[s, k] + 
          a[s, k]\[Conjugate] b[s, k] kpc[s, k]\[Conjugate])/(Abs[
            a[s, k]]^2 + Abs[b[s, k]]^2)},
    Compile[{{G, _Real}, {s, _Real}, {k, _Real}},
     code,
     CompilationTarget -> "C",
     RuntimeAttributes -> {Listable},
     Parallelization -> True
     ]
    ]
   ];

{sdata, kdata} = Transpose[Tuples[{Subdivide[0., 2. Pi, 30], Rest@Subdivide[-1. Pi , 1. Pi, M]}]];
F[G_] := Total[cf[G, sdata, kdata]];
Glist = Subdivide[0.1, 10., 1000];
data = F /@ Glist; // AbsoluteTiming // First
ListLinePlot[Transpose[{Glist, data}]]

3.29177

edited Dec 10 '18 at 16:03

answered Dec 10 '18 at 15:07

Henrik Schumacher

106,770
7
179
309

thx very much! Did you create a list of 3-variables function (G,s,k) and simply sum over k and s, right? I obtained the plot, however, I got also 2 errors: 1) CCompilerDriver`CreateLibrary::nocomp: A C compiler cannot be found on your system. Please consult the documentation to learn how to set up suitable compilers. 2) Compile::nogen: A library could not be generated from the compiled function. – Chipa-Chipa Dec 10 '18 at 15:58
Ah. Either you have to install a C-compiler for maximal performance (e.g., cygwin for Windows, gcc for Linux, or XCode with command line tool for macOS) or just set `CompilationTarget -> "WVM" (will be slower though). – Henrik Schumacher Dec 10 '18 at 16:01
The function cf is indeed a 3-argument function. With RuntimeAttributes -> {Listable}, I enfoce that it threads over input lists, so cf[G, sdata, kdata] returns a list of function values that simply have to be summed up. – Henrik Schumacher Dec 10 '18 at 16:04

fast sum computation

1 Answers1