How can I speed up this code with 4 variables sum

Question

I have to make a sum over 4 variables. My code is very very slow. I want to know how to speed up this code. This problem is related to but different from one previous problem. Any help or suggestion will be highly appreciated! The code is shown below：

data = Table[   Exp[-((i + j - 20.)/5)^2] Exp[-((i - j)/5)^2], {i, 20}, {j, 20}];
data = Chop[data, 0.00001];
data = data/Sqrt[Sum[(data[[i, j]])^2, {i, 1, 20}, {j, 1, 20}]];
ListDensityPlot[data, InterpolationOrder -> 0, Mesh -> All,  PlotRange -> All, ColorFunction -> (Blend[{Hue[2/3], Hue[0]}, #] &)]
c = 3*10^8; 
Δ = 0.5; 
λ0 = 1500;
CC1[i_, j_, k_, l_, t_] := (data[[i, l]] data[[j, 
  k]] Cos[π*(c/(λ0 - 10 + i*Δ - 
      0.5 Δ) + 
     c/(λ0 - 10 + l*Δ - 
      0.5 Δ)) t] Cos[π*(c/(λ0 - 10 + 
      j*Δ - 0.5 Δ) + 
     c/(λ0 - 10 + k*Δ - 
      0.5 Δ)) t] + 
data[[i, k]] data[[j, 
  l]] Cos[π*(c/(λ0 - 10 + i*Δ - 
      0.5 Δ) + 
     c/(λ0 - 10 + k*Δ - 
      0.5 Δ)) t] Cos[π*(c/(λ0 - 10 + 
      j*Δ - 0.5 Δ) + 
     c/(λ0 - 10 + l*Δ - 
      0.5 Δ)) t] - 
data[[i, j]] data[[k, 
  l]] Sin[π*(c/(λ0 - 10 + i*Δ - 
      0.5 Δ) + 
     c/(λ0 - 10 + j*Δ - 
      0.5 Δ)) t] Sin[π*(c/(λ0 - 10 + 
      k*Δ - 0.5 Δ) + 
     c/(λ0 - 10 + l*Δ - 
      0.5 Δ)) t])^2;
CC2[t_] := \!\(\*UnderoverscriptBox[\(∑\), \(i = 1\), \(20\)]\(\*UnderoverscriptBox[\(∑\), \(j = 1\), \(20\)]\(\*UnderoverscriptBox[\(∑\), \(k = 1\), \(20\)]\(\*UnderoverscriptBox[\(∑\), \(l = 1\), \(20\)]CC1[i, j, k, l,t]\)\)\)\);
ListPlot[Table[{i, CC2[i*0.001]}, {i, -10, 10, 1}], Joined -> True, Axes -> None, PlotRange -> All, Frame -> True, ImageSize -> {400, 250}]
ListPlot[Table[{i, CC2[i*0.001]}, {i, -10, 10, 0.001}],  Joined -> True, Axes -> None, PlotRange -> All, Frame -> True,  ImageSize -> {400, 250}]

Answer 1

The problem is that it reevaluates the sum every single time you call it, recomputing every 20^4 term again and again. You just need to compile the function CC2 so that it performs the summation only once. Using the code you have, it takes my machine about 6 seconds to compute a single data point:

CC2[0.003] // AbsoluteTiming

(*  {6.069311, 1.49893}  *)

But if I compile it first,

CC2comp = Compile[{{t, _Real}},
  Evaluate[
     Sum[CC1[i, j, k, l, t], {i, 20}, {j, 20}, {k, 20}, {l, 20}]
     ]];

I now get the answer in milliseconds

CC2comp[0.003] // AbsoluteTiming

{0.002888, 1.49893}

Now you can make the plot you were going for in about a minute.

ListPlot[Table[{i, CC2comp[i*0.001]}, {i, -10, 10, 0.001}],  
  Joined -> True, Axes -> None, PlotRange -> All, Frame -> True,  
  ImageSize -> {400, 250}] // AbsoluteTiming

Answer 2

Okay, so I tried to write the function like you were saying, taking advantage of the Total built-in function, and there is some improvement.

cc2Comp = Compile[{{t, _Real}},
Evaluate@
Module[{ndim, data, c, Δ, λ0, costable, 
  sintable, table1, table2, table3, table4},
  ndim = 20;
  data = 
  Array[Exp[-((#1 + #2 - 20.)/5)^2] Exp[-((#1 - #2)/5)^2] &
    , {ndim, ndim}];
  data = Chop[data, 0.00001];
  data = data/Norm[data];
  data = SparseArray[data];
  c = 3*10^8;
  Δ = 0.5;
  λ0 = 1500;

  costable = 
  Table[Cos[π t (c/(-10 - 0.5` Δ + 
    j Δ + λ0) + c/(-10 - 0.5` Δ + k Δ + λ0))], {j, ndim}, {k,ndim}];
  sintable = 
  Table[Sin[π t (c/(-10 - 0.5` Δ + 
    i Δ + λ0) + c/(-10 - 0.5` Δ + j Δ + λ0))] , {i, ndim}, {j, ndim}];
  table1 = 
  Transpose[
    Outer[Times, costable data, costable data], {1, 4, 2, 3}];
  table2 = 
  Transpose[
    Outer[Times, costable data, costable data], {1, 3, 2, 4}];
  table3 = Outer[Times, sintable data, sintable data];


  (* This next part can be modified, changing the amount of time needed to 
    compile the code versus the time needed to execute the code  *)
  Total[
  table1^2 + 2 table1 table2 + table2^2 - 2 table1 table3 - 2 table2 table3 + table3^2, 4]
  ]
]; // AbsoluteTiming


cc2Comp[0.003] // AbsoluteTiming
ListPlot[Table[{i, cc2Comp[i*0.001]}, {i, -10, 10, 0.001}], 
  Joined -> True, Axes -> None, PlotRange -> All, Frame -> True, 
  ImageSize -> {400, 250}] // AbsoluteTiming

As you can see, on my machine, this took about 10 seconds to compile (similar to the method in my first answer), about .026 seconds to compute a single data point (ten times longer than above), and 16 seconds to make the plot (one fourth of the time above).

If I change the last line of the Module to read

table4 = SparseArray[table1 + table2 - table3];

Total[(table4)^2, 4]

Then the result is it takes less time to compile but more time to make the plot, and a little more time altogether.

Just straightforwardly applying this method to larger arrays, as in your other question, is problematic. I don't have the RAM to work do this method when the dimension is raised to 100 (since we are talking about 100 million matrix elements in the summation), and it crashes the kernel.