add spectrum plot

bfps/test/test_Parseval.py

@@ -8,6 +8,8 @@ import sys
 import bfps
 from bfps import DNS
 
+import matplotlib.pyplot as plt
+
 def main():
     niterations = 10
     c = DNS()
@@ -25,7 +27,20 @@ def main():
              '--wd', './'] +
              sys.argv[1:])
     c.compute_statistics()
-    print((c.statistics['energy(t)'] - c.statistics['renergy(t)']) / c.statistics['renergy(t)'])
+    print((c.statistics['energy(t)']*3 - c.statistics['renergy(t)']) / c.statistics['renergy(t)'])
+    print(list(c.statistics.keys()))
+    energyk = np.mean(c.statistics['energy(t, k)'], axis = 0)
+    nshell = c.get_data_file()['kspace/nshell'].value
+    renergy = np.mean(c.statistics['renergy(t)'])
+    print(renergy, np.trapz(energyk[:-2], c.statistics['kshell'][:-2]))
+
+    f = plt.figure()
+    a = f.add_subplot(111)
+    a.plot(c.statistics['kshell'], energyk)
+    a.plot(c.statistics['kshell'], (energyk / nshell)*(4*np.pi*c.statistics['kshell']**2))
+    a.set_yscale('log')
+    a.set_xscale('log')
+    f.savefig('spectrum.pdf')
     return None
 
 if __name__ == '__main__':