Add outer product of self-energies

adga.py

@@ -178,3 +178,20 @@ class Adga(object):
         output_file.close()
         susceptibilities = np.array(susceptibilities)
         return susceptibilities[:, 0, 0, :, :, 0, :]
+
+def outer_product_of_self_energies(a, b):
+    """Return the outer product of ``a`` and ``b`` in the frequency
+    index.
+    
+    :arg a: self-energy as 3D numpy array
+    :arg b: self-energy as 3D numpy array
+    :return: 4D numpy array :math:`c_{(i)(j)kl}=a_{(i)(j)k}b_{(i)(j)l}`
+    """
+    
+    result_shape = list(a.shape)
+    result_shape.append(a.shape[-1])
+    result = np.empty(result_shape)
+    for i in range(len(a)):
+        for j in range(len(a[0])):
+            result[i,j] = np.outer(a[i,j,:], np.conjugate(b[i,j,:]))
+    return result

main.py

@@ -60,7 +60,9 @@ else:
 
 # %%
 # Do the Jackknife estimation
-jackknife = jk.Jackknife(x_generator, n, f, general["output_file_prefix"],
+jackknife = jk.Jackknife(x_generator, n, f,
+                         adga.outer_product_of_self_energies,
+                         general["output_file_prefix"],
                          strtobool(general["store_output_samples"]))
 jackknife.do_estimation()
 jackknife.write_results_to_file()