Stash

notebooks/lib.py

@@ -0,0 +1,69 @@
+import os
+from glob import glob
+from pathlib import Path
+
+import numpy as np
+import pandas as pd
+
+
+def read_xy(path: str):
+    df = pd.read_csv(path, skipinitialspace=True)
+    df['N'] = df.index + 1
+    df['r'] = (df.x ** 2 + df.y ** 2) ** 0.5
+    df['cr'] = df.r.cummax()
+    df['fd'] = np.log(df.N) / np.log(df.cr)
+    df['run'] = os.path.splitext(Path(path).name)[0]
+
+    return df.replace([np.inf, -np.inf], np.nan).dropna()
+
+
+def read_load(load_dir: str, reader=read_xy):
+    paths = glob(f'{load_dir}/*.csv')
+    return pd.concat([reader(path) for path in paths])
+
+
+def read_sp_xy(specific_probability_dir: str):
+    probability = float(Path(specific_probability_dir).name)
+    df = read_load(specific_probability_dir)
+    df['probability'] = probability
+
+    return df
+
+
+def read_sp(sp_dir: str):
+    if not Path(sp_dir).exists():
+        raise Exception("Root does not exist")
+
+    return pd.concat([read_sp_xy(specific_probability_dir) for specific_probability_dir in glob(f'{sp_dir}/*')])
+
+
+def convergent_tail_index(series, tol):
+    diffs = np.abs(np.ediff1d(series))
+    for i in range(0, len(diffs)):
+        if np.max(diffs[i:]) <= tol:
+            return i
+
+    # No convergence found
+    return None
+
+
+def mean_of_tail(series, tol=0.05):
+    tail_index = convergent_tail_index(series, tol)
+    if tail_index is None:
+        raise Exception("No convergence found.")
+
+    return np.mean(series[tail_index:])
+
+
+def std_of_tail(series, tol=0.05):
+    tail_index = convergent_tail_index(series, tol)
+    if tail_index is None:
+        raise Exception("No convergence found.")
+
+    return np.std(series[tail_index:])
+
+
+def fd_stats(dfs):
+    fds = [mean_of_tail(df.fd, 0.1) for df in dfs]
+    fds_clean = [f for f in fds if f < np.inf]
+    return np.mean(fds_clean), np.mean(fds_clean) / np.sqrt(fds_clean.length())

notebooks/nc-fd-convergence.py

@@ -0,0 +1,36 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from notebooks.lib import read_load
+
+alpha = read_load("../data/alpha")
+meaned_by_N = alpha.groupby('N').agg({'fd': ['mean', 'std']}) \
+    .reset_index() \
+    .replace([np.inf, -np.inf], np.nan)
+
+without_prefix = meaned_by_N[50:]
+fig, ax = plt.subplots(figsize=(6, 6))
+
+plt.fill_between(
+    without_prefix.N,
+    # TODO Check error math here
+    (without_prefix['fd']['mean'] - without_prefix['fd']['std'] / np.sqrt(20)),
+    (without_prefix['fd']['mean'] + without_prefix['fd']['std'] / np.sqrt(20)),
+    alpha=0.2, label=f"Standard error band"
+)
+
+plt.plot(
+    without_prefix.N,
+    without_prefix['fd']['mean'],
+    color='tab:blue', label='fd mean, seeds = 20'
+)
+
+plt.plot([50, 10000], [1.71, 1.71], color='red', label='Theory')
+plt.fill_between(without_prefix.N, 1.71 - 0.01, 1.71 + 0.01, alpha=0.2, label='Theory error band')
+
+plt.xlabel("$N_C$")
+plt.ylabel("$fd$ (instantaneous)")
+plt.legend()
+
+plt.savefig('../figures/nc-fd-convergence.svg')
+plt.savefig('../figures/nc-fd-convergence.png')
+plt.show()